package org.opentrafficsim.base.geometry;
   
   import java.awt.geom.Line2D;
   
   import org.djunits.value.vdouble.scalar.Direction;
   import org.djutils.draw.line.PolyLine2d;
   import org.djutils.draw.point.Point2d;
   import org.djutils.exceptions.Throw;
   import org.djutils.math.AngleUtil;
  
  /**
   * Fractional projection helper for {@link OtsLine2d}.
   * <p>
   * Copyright (c) 2026-2026 Delft University of Technology, PO Box 5, 2600 AA, Delft, the Netherlands. All rights reserved. <br>
   * BSD-style license. See <a href="https://opentrafficsim.org/docs/license.html">OpenTrafficSim License</a>.
   * </p>
   * @author Wouter Schakel
   * @see OtsLine2d#projectFractionalAt
   */
  public final class FractionalProjectionHelper
  {
  
      /** Numerical precision for fractional projection comparisons. */
      private static final double FRAC_PROJ_PRECISION = 2e-5;
  
      /** Epsilon to consider two lines near-parallel in intersection tests. */
      private static final double INTERSECTION_EPS = 1e-7;
  
      /** Owning line (used to avoid duplicating points and cumulative lengths). */
      private final OtsLine2d line;
  
      /** Number of segments. */
      private final int n;
  
      /** Whether fixed helper 1..n-2 is a center, other wise a direction. */
      private final boolean[] isCenter; // length n; interior entries set, edges ignored
  
      /** Helper center x coordinate. */
      private final double[] centerX; // length n; only valid when isCenter[i] == true
  
      /** Helper center y coordinate. */
      private final double[] centerY; // length n
  
      /** Helper direction x component. */
      private final double[] dirX; // length n; only valid when isCenter[i] == false
  
      /** Helper direction y component. */
      private final double[] dirY; // length n
  
      /** Intersection of first two unit-offset segments (start direction-independent). */
      private Point2d firstOffsetIntersection;
  
      /** Intersection of last two unit-offset segments (end direction-independent). */
      private Point2d lastOffsetIntersection;
  
      /**
       * Cached helpers for {@code null} start and end direction. Often used by {@link OtsLine2d#radiusAtFraction} and
       * {@link OtsLine2d#radiusAtVertex(int)}.
       */
      private EdgeHelpers edgeNullNull; // cached helpers for (null,null)
  
      /** Last none-both-{@code null} directions. */
      private DirectionKey lastKey; // last used (may include null on one side)
  
      /** Helpers of last none-both-{@code null} directions. */
      private EdgeHelpers edgeLast; // cached helpers for lastKey
  
      /**
       * Constructor.
       * @param line line
       */
      FractionalProjectionHelper(final OtsLine2d line)
      {
          this.line = line;
          this.n = this.line.size() - 1;
  
          this.isCenter = new boolean[this.n];
          this.centerX = new double[this.n];
          this.centerY = new double[this.n];
          this.dirX = new double[this.n];
          this.dirY = new double[this.n];
  
          // Pre-compute interior (none-edge) helpers and the two direction-independent offset intersections
          precomputeFixedHelpers();
      }
  
      /**
       * Pre-computes all fixed helpers.
       */
      private void precomputeFixedHelpers()
      {
          if (this.n < 2)
          {
              // No interior segments; nothing to pre-compute
              return;
          }
  
          // First two unit-offset segments (direction independent)
          PolyLine2d prevOfs = unitOffsetSegment(0);
         PolyLine2d nextOfs = unitOffsetSegment(1);
 
         Point2d parStart = intersectionOrFallbackMid(prevOfs.get(0), prevOfs.get(1), nextOfs.get(0), nextOfs.get(1));
         this.firstOffsetIntersection = parStart;
 
         // Special case: exactly two segments -> lastOffsetIntersection equals parStart
         if (this.n == 2)
         {
             this.lastOffsetIntersection = parStart;
             return; // no interior helpers to compute
         }
 
         // Build for interior segments i = 1 .. n-2
         for (int i = 1; i <= this.n - 2; i++)
         {
             prevOfs = nextOfs;
             if (i + 1 <= this.n - 1)
             {
                 nextOfs = unitOffsetSegment(i + 1);
             }
             final Point2d parEnd = intersectionOrFallbackMid(prevOfs.get(0), prevOfs.get(1), nextOfs.get(0), nextOfs.get(1));
 
             // Intersection of helper lines: (vertex i -> parStart) and (vertex i+1 -> parEnd)
             final Point2d c = intersectionOrNull(this.line.get(i), parStart, this.line.get(i + 1), parEnd);
             if (c != null)
             {
                 this.isCenter[i] = true;
                 this.centerX[i] = c.x;
                 this.centerY[i] = c.y;
             }
             else
             {
                 this.isCenter[i] = false;
                 this.dirX[i] = parStart.x - this.line.get(i).x;
                 this.dirY[i] = parStart.y - this.line.get(i).y;
             }
 
             parStart = parEnd;
             if (i == this.n - 2)
             {
                 this.lastOffsetIntersection = parStart;
             }
         }
     }
 
     /**
      * Unit-offset line to the LEFT of segment i (distance 1.0).
      * @param i segment
      * @return offset line
      */
     private PolyLine2d unitOffsetSegment(final int i)
     {
         return new PolyLine2d(this.line.get(i), this.line.get(i + 1)).offsetLine(1.0);
     }
 
     /**
      * Fractionally project a point on the polyline using the fractional helper logic. Falls back via the given strategy when
      * fractional projection is not applicable.
      * @param start direction in first point
      * @param end direction in last point
      * @param x x-coordinate of point to project
      * @param y y-coordinate of point to project
      * @param fallback fallback method for when fractional projection fails
      * @return fractional position along this line of the fractional projection on that line of a point
      * @see OtsLine2d#projectFractionalAt
      */
     public synchronized double projectFractionalAt(final Direction start, final Direction end, final double x, final double y,
             final FractionalFallback fallback)
     {
         Throw.whenNull(fallback, "fallback");
 
         // Determine edge helpers for these directions (2-slot cache: (null,null) + lastKey).
         final EdgeHelpers edges = getEdgeHelpers(start, end);
 
         // Compute distances to segments; consider only those near the minimum distance.
         double minD = Double.POSITIVE_INFINITY;
         final double[] segDist = new double[this.n];
         for (int i = 0; i < this.n; i++)
         {
             Point2d a = this.line.get(i);
             Point2d b = this.line.get(i + 1);
             segDist[i] = Line2D.ptSegDist(a.x, a.y, b.x, b.y, x, y);
             if (segDist[i] < minD)
             {
                 minD = segDist[i];
             }
         }
 
         double bestDistance = Double.POSITIVE_INFINITY;
         int bestSeg = -1;
         double bestSegFrac = 0.0;
 
         final Point2d ext = new Point2d(x, y);
 
         for (int i = 0; i < this.n; i++)
         {
             if (segDist[i] > minD + FRAC_PROJ_PRECISION)
             {
                 continue;
             }
             final Helper h = helperForSegment(i, edges);
             final Point2d p = intersectProjection(i, h, ext);
             if (p == null)
             {
                 continue;
             }
 
             // Ensure intersection lies on segment (within tolerance)
             Point2d a = this.line.get(i);
             Point2d b = this.line.get(i + 1);
             final double segLen = a.distance(b) + FRAC_PROJ_PRECISION;
             if (p.distance(a) > segLen || p.distance(b) > segLen)
             {
                 continue;
             }
 
             // Prefer the nearest intersection to the external point
             final double dist = p.distance(ext);
             if (dist < bestDistance)
             {
                 bestDistance = dist;
                 // Compute fraction within segment
                 final double segActual = a.distance(b);
                 final double along = a.distance(p);
                 bestSegFrac = Math.min(1.0, Math.max(0.0, segActual > 0.0 ? (along / segActual) : 0.0));
                 bestSeg = i;
             }
         }
 
         if (bestSeg < 0)
         {
             // Fractional projection not applicable; fallback
             return fallback.getFraction(this, x, y);
         }
 
         // Convert segment-local fraction to global fraction using line.lengthAtIndex
         final double segStartLen = this.line.lengthAtIndex(bestSeg);
         final double segEndLen = this.line.lengthAtIndex(bestSeg + 1);
         final double abs = segStartLen + bestSegFrac * (segEndLen - segStartLen);
         final double total = this.line.lengthAtIndex(this.line.size() - 1);
         return abs / total;
     }
 
     /**
      * Get the helper (center or direction) at a vertex index using (start,end) dependent edges. Useful for curvature / radius
      * logic. Supply ({@code null}, {@code null}) if that is the desired configuration.
      * @param vertexIndex 0 .. size-1
      * @param start direction at start (can be {@code null})
      * @param end direction at end (can be {@code null})
      * @return helper
      */
     public synchronized Helper helperAtVertex(final int vertexIndex, final Direction start, final Direction end)
     {
         Throw.when(vertexIndex < 0 || vertexIndex > this.line.size() - 1, IndexOutOfBoundsException.class,
                 "vertexIndex %s out of bounds [0..%s]", vertexIndex, this.line.size() - 1);
         if (vertexIndex < 0 || vertexIndex > this.line.size() - 1)
         {
             throw new IndexOutOfBoundsException("vertexIndex out of bounds");
         }
         final EdgeHelpers edges = getEdgeHelpers(start, end);
         if (vertexIndex == 0)
         {
             return edges.first;
         }
         if (vertexIndex >= this.n - 1)
         {
             return edges.last;
         }
         // Interior vertex i corresponds to segment i (between i and i+1)
         return helperForInteriorSegment(vertexIndex);
     }
 
     /**
      * Fallback strategies for when fractional projection is not applicable.
      */
     public enum FractionalFallback
     {
 
         /** Orthogonal projection clamped to [0,1]. */
         ORTHOGONAL
         {
             @Override
             double getFraction(final FractionalProjectionHelper helper, final double x, final double y)
             {
                 return helper.line.projectOrthogonalSnapAt(x, y);
             }
         },
 
         /** Orthogonal projection allowing extension beyond end-points. */
         ORTHOGONAL_EXTENDED
         {
             @Override
             double getFraction(final FractionalProjectionHelper helper, final double x, final double y)
             {
                 return helper.line.projectOrthogonalSnapAt(x, y, false);
             }
         },
 
         /** Nearest end-point as fraction &lt;0 before start, &gt;1 after end. */
         ENDPOINT
         {
             @Override
             double getFraction(final FractionalProjectionHelper helper, final double x, final double y)
             {
                 final Point2d p = new Point2d(x, y);
                 final Point2d a = helper.line.get(0);
                 final Point2d b = helper.line.get(helper.line.size() - 1);
                 final double dStart = p.distance(a);
                 final double dEnd = p.distance(b);
                 final double total = helper.line.lengthAtIndex(helper.line.size() - 1);
                 return (dStart < dEnd) ? (-dStart / total) : ((total + dEnd) / total);
             }
         },
 
         /** Return NaN. */
         NaN
         {
             @Override
             double getFraction(final FractionalProjectionHelper h, final double x, final double y)
             {
                 return Double.NaN;
             }
         };
 
         /**
          * Returns fraction for when fractional projection fails as the point is beyond the line or from numerical limitations.
          * @param helper helper
          * @param x x coordinate of point
          * @param y y coordinate of point
          * @return fraction for when fractional projection fails
          */
         abstract double getFraction(FractionalProjectionHelper helper, double x, double y);
     }
 
     /**
      * Value object: per-segment helper, either center or direction.
      * @param hasCenter whether this is a center helper
      * @param cx center x coordinate
      * @param cy center y coordinate
      * @param dx direction x component
      * @param dy direction y component
      */
     public record Helper(boolean hasCenter, double cx, double cy, double dx, double dy)
     {
         /**
          * Factory for a center-based helper.
          * @param cx center x coordinate
          * @param cy center y coordinate
          * @return center-based helper
          */
         public static Helper center(final double cx, final double cy)
         {
             return new Helper(true, cx, cy, Double.NaN, Double.NaN);
         }
 
         /**
          * Factory for a direction-based helper.
          * @param dx direction x component
          * @param dy direction y component
          * @return direction-based helper
          */
         public static Helper direction(final double dx, final double dy)
         {
             return new Helper(false, Double.NaN, Double.NaN, dx, dy);
         }
     }
 
     /**
      * Value object for the two edge helpers.
      * @param first helper
      * @param last helper
      */
     private record EdgeHelpers(Helper first, Helper last)
     {
     }
 
     /**
      * Cache key for directions (normalized angles; nulls allowed).
      * @param start start direction
      * @param end end direction
      */
     public record DirectionKey(Double start, Double end)
     {
         /**
          * Create key with normalized directions, which may be {@code null}.
          * @param start start direction
          * @param end end direction
          * @return key
          */
         public static DirectionKey of(final Direction start, final Direction end)
         {
             Double sa = (start == null) ? null : AngleUtil.normalizeAroundZero(start.si);
             Double ea = (end == null) ? null : AngleUtil.normalizeAroundZero(end.si);
             return new DirectionKey(sa, ea);
         }
     }
 
     /**
      * Returns helper for interior (non-edge) segment.
      * @param segIndex index
      * @return helper for interior (non-edge) segment
      */
     private Helper helperForInteriorSegment(final int segIndex)
     {
         if (this.isCenter[segIndex])
         {
             return Helper.center(this.centerX[segIndex], this.centerY[segIndex]);
         }
         return Helper.direction(this.dirX[segIndex], this.dirY[segIndex]);
     }
 
     /**
      * Returns helper for segment.
      * @param segIndex index
      * @param edges edge helpers
      * @return helper for segment
      */
     private Helper helperForSegment(final int segIndex, final EdgeHelpers edges)
     {
         if (segIndex == 0)
         {
             return edges.first;
         }
         if (segIndex == this.n - 1)
         {
             return edges.last;
         }
         return helperForInteriorSegment(segIndex);
     }
 
     /**
      * Returns edge helpers based on directions.
      * @param start start direction
      * @param end end direction
      * @return edge helpers based on directions
      */
     private EdgeHelpers getEdgeHelpers(final Direction start, final Direction end)
     {
         if (start == null && end == null)
         {
             if (this.edgeNullNull == null)
             {
                 this.edgeNullNull = computeEdgeHelpers(null, null);
             }
             return this.edgeNullNull;
         }
         final DirectionKey key = DirectionKey.of(start, end);
         if (key.equals(this.lastKey) && this.edgeLast != null)
         {
             return this.edgeLast;
         }
         this.edgeLast = computeEdgeHelpers(start, end);
         this.lastKey = key;
         return this.edgeLast;
     }
 
     /**
      * Compute edge helpers.
      * @param start start direction
      * @param end end direction
      * @return edge helpers
      */
     private EdgeHelpers computeEdgeHelpers(final Direction start, final Direction end)
     {
         // Angles (default to segment direction if null)
         final double startAng = (start == null)
                 ? Math.atan2(this.line.get(1).y - this.line.get(0).y, this.line.get(1).x - this.line.get(0).x) : start.si;
         final double endAng = (end == null) ? Math.atan2(this.line.get(this.n).y - this.line.get(this.n - 1).y,
                 this.line.get(this.n).x - this.line.get(this.n - 1).x) : end.si;
 
         // Unit offset points at start and end of line
         final Point2d p1 = new Point2d(this.line.get(0).x + Math.cos(startAng + Math.PI / 2.0),
                 this.line.get(0).y + Math.sin(startAng + Math.PI / 2.0));
         final Point2d p2 = new Point2d(this.line.get(this.n).x + Math.cos(endAng + Math.PI / 2.0),
                 this.line.get(this.n).y + Math.sin(endAng + Math.PI / 2.0));
 
         Helper hFirst, hLast;
         if (this.n == 1)
         {
             // Single segment: edges collapse; center is intersection of offset rays
             final Point2d c = intersectionOrNull(this.line.get(0), p1, this.line.get(1), p2);
             if (c != null)
             {
                 hFirst = Helper.center(c.x, c.y);
             }
             else
             {
                 hFirst = Helper.direction(p1.x - this.line.get(0).x, p1.y - this.line.get(0).y);
             }
             // Same for last (same segment)
             hLast = hFirst;
         }
         else
         {
             // Use direction-independent offset intersections
             Point2d cFirst = intersectionOrNull(this.line.get(0), p1, this.line.get(1), this.firstOffsetIntersection);
             if (cFirst != null)
             {
                 hFirst = Helper.center(cFirst.x, cFirst.y);
             }
             else
             {
                 hFirst = Helper.direction(p1.x - this.line.get(0).x, p1.y - this.line.get(0).y);
             }
 
             Point2d cLast =
                     intersectionOrNull(this.line.get(this.n - 1), this.lastOffsetIntersection, this.line.get(this.n), p2);
             if (cLast != null)
             {
                 hLast = Helper.center(cLast.x, cLast.y);
             }
             else
             {
                 hLast = Helper.direction(p2.x - this.line.get(this.n).x, p2.y - this.line.get(this.n).y);
             }
         }
         return new EdgeHelpers(hFirst, hLast);
     }
 
     /**
      * Projects external point to segment.
      * @param segIndex index
      * @param helper helper of the segment
      * @param ext external point
      * @return projection to segment, or {@code null} if no valid projection
      */
     private Point2d intersectProjection(final int segIndex, final Helper helper, final Point2d ext)
     {
         final Point2d a = this.line.get(segIndex);
         final Point2d b = this.line.get(segIndex + 1);
 
         if (helper.hasCenter)
         {
             // Intersection of (center -> ext) ray with segment line
             final Point2d c = new Point2d(helper.cx, helper.cy);
             final Point2d p = intersectionOrNull(c, ext, a, b);
             if (p == null)
             {
                 return null;
             }
             // Ensure center is not between ext and intersection:
             final double v1x = p.x - c.x, v1y = p.y - c.y;
             final double v2x = ext.x - c.x, v2y = ext.y - c.y;
             final double dot = v1x * v2x + v1y * v2y;
             if (dot <= FRAC_PROJ_PRECISION)
             {
                 return null;
             }
             return p;
         }
         // Parallel helper lines: project along stored direction from ext
         final Point2d off = new Point2d(ext.x + helper.dx, ext.y + helper.dy);
         return intersectionOrNull(ext, off, a, b);
     }
 
     /**
      * Intersection of infinite lines (p1->p2) and (p3->p4). Returns null if near-parallel.
      * @param p1 first point of first line
      * @param p2 second point of first line
      * @param p3 first point of second line
      * @param p4 second point of second line
      * @return intersection
      */
     // TODO: can be replaced with djutils version with eps once it is published in djutils
     private static Point2d intersectionOrNull(final Point2d p1, final Point2d p2, final Point2d p3, final Point2d p4)
     {
         final double x1 = p1.x, y1 = p1.y;
         final double x2 = p2.x, y2 = p2.y;
         final double x3 = p3.x, y3 = p3.y;
         final double x4 = p4.x, y4 = p4.y;
 
         final double dx1 = x2 - x1, dy1 = y2 - y1;
         final double dx2 = x4 - x3, dy2 = y4 - y3;
 
         final double denom = dx1 * dy2 - dy1 * dx2;
         if (Math.abs(denom) < INTERSECTION_EPS)
         {
             return null; // near-parallel
         }
         final double t = ((x3 - x1) * dy2 - (y3 - y1) * dx2) / denom;
         return new Point2d(x1 + t * dx1, y1 + t * dy1);
     }
 
     /**
      * Returns intersection for parallel segments, or midpoint fallback when intersection is null or unstable.
      * @param a1 first point of first parallel segment
      * @param a2 second point of first parallel segment
      * @param b1 first point of second parallel segment
      * @param b2 second point of second parallel segment
      * @return intersection
      */
     private static Point2d intersectionOrFallbackMid(final Point2d a1, final Point2d a2, final Point2d b1, final Point2d b2)
     {
         final Point2d inter = intersectionOrNull(a1, a2, b1, b2);
         if (inter == null)
         {
             return midpoint(a2, b1);
         }
         final double dStraight = a2.distance(b1);
         if (dStraight < Math.min(a2.distance(inter), b1.distance(inter)))
         {
             return midpoint(a2, b1);
         }
         return inter;
     }
 
     /**
      * Computes the mid-point.
      * @param p first point
      * @param q second point
      * @return mid-point
      */
     private static Point2d midpoint(final Point2d p, final Point2d q)
     {
         return new Point2d(0.5 * (p.x + q.x), 0.5 * (p.y + q.y));
     }
 
 }