packages feed

moonlight-triangulation-0.1.0.0: src-dcel/Moonlight/Triangulation/IntersectionIterator.hs

{-# LANGUAGE BangPatterns #-}

module Moonlight.Triangulation.IntersectionIterator
  ( Intersection (..)
  , lineIntersections
  , lineIntersectionsBetweenVertices
  , foldCorridorBetweenPoints
  , foldCorridorBetweenVertices
  , conflictingEdges
  , segmentIntersectsNonCollinear
  ) where

import Data.List (sortBy)
import Data.Maybe (mapMaybe)
import Data.Ord (comparing)
import Data.Void (Void, absurd)
import Moonlight.Triangulation.Dcel
import Moonlight.Triangulation.Handles.HandleDefs
import Moonlight.Triangulation.Handles.Iterators.FixedIterators (undirectedEdges, vertices)
import Moonlight.Triangulation.Math
import Moonlight.Triangulation.PointLocation
import Moonlight.Triangulation.Types

-- | One crossing: an edge cut, a vertex hit, or a collinear overlap.
data Intersection
  = EdgeIntersection !DirectedEdgeId
  | VertexIntersection !VertexId
  | EdgeOverlap !DirectedEdgeId
  deriving stock (Eq, Ord, Show)

-- | Every crossing between two points, ordered along the segment.
lineIntersections
  :: Triangulation mode vertex directed undirected face
  -> QueryPoint
  -> QueryPoint
  -> [Intersection]
lineIntersections triangulation queryFrom queryTo =
  let !from = queryPointValue queryFrom
      !to = queryPointValue queryTo
   in case firstIntersection triangulation queryFrom queryTo of
        Nothing -> []
        Just first -> walkIntersections triangulation from to first

-- | 'lineIntersections' between two existing vertices.
lineIntersectionsBetweenVertices :: Triangulation mode vertex directed undirected face -> VertexId -> VertexId -> [Intersection]
lineIntersectionsBetweenVertices triangulation fromVertex toVertex =
  let from = vertexPoint triangulation fromVertex
      to = vertexPoint triangulation toVertex
   in walkIntersections triangulation from to (VertexIntersection fromVertex)

-- | The directed edges a crossing list cuts.
conflictingEdges :: [Intersection] -> [DirectedEdgeId]
conflictingEdges = mapMaybe asConflict
 where
  asConflict (EdgeIntersection edge) = Just edge
  asConflict _ = Nothing

-- | Fold the corridor one crossing at a time, stopping the instant the step
-- function answers.
--
-- The list-producing walks cannot stop early: their step budget is only known
-- to have been respected once the walk ends, so the whole corridor is
-- materialized before the first event is visible. A caller whose answer is
-- settled by a prefix — anything asking whether some crossing exists — should
-- not pay for the suffix. 'Nothing' reports a walk that outran its budget and
-- is the caller's signal to fall back to the materialized walk, which
-- substitutes the exact scan.
foldCorridorBetweenVertices
  :: Triangulation mode vertex directed undirected face
  -> VertexId
  -> VertexId
  -> (state -> Intersection -> Either answer state)
  -> state
  -> Maybe (Either answer state)
foldCorridorBetweenVertices triangulation fromVertex toVertex =
  foldCorridor
    triangulation
    (vertexPoint triangulation fromVertex)
    (vertexPoint triangulation toVertex)
    (VertexIntersection fromVertex)

-- | As 'foldCorridorBetweenVertices', for a corridor given by its endpoints.
-- A segment that meets nothing at all folds to the initial state.
foldCorridorBetweenPoints
  :: Triangulation mode vertex directed undirected face
  -> QueryPoint
  -> QueryPoint
  -> (state -> Intersection -> Either answer state)
  -> state
  -> Maybe (Either answer state)
foldCorridorBetweenPoints triangulation queryFrom queryTo step state =
  case firstIntersection triangulation queryFrom queryTo of
    Nothing -> Just (Right state)
    Just first ->
      foldCorridor
        triangulation
        (queryPointValue queryFrom)
        (queryPointValue queryTo)
        first
        step
        state

foldCorridor
  :: Triangulation mode vertex directed undirected face
  -> Point
  -> Point
  -> Intersection
  -> (state -> Intersection -> Either answer state)
  -> state
  -> Maybe (Either answer state)
foldCorridor triangulation from to first step =
  go (2 * numDirectedEdges triangulation + numVertices triangulation + 8) first
 where
  go !remaining !current !state
    | remaining <= 0 = Nothing
    | otherwise = case step state current of
        Left answer -> Just (Left answer)
        Right advanced ->
          case nextIntersection triangulation from to current of
            Nothing -> Just (Right advanced)
            Just following -> go (remaining - 1) following advanced
-- The early-answer 'Either' exists to let a caller stop at the first event it
-- cares about; it should not survive to runtime. Inlining the non-recursive
-- wrapper puts the worker at each call site with 'step' statically known, so
-- the constructor is matched where it is built — and at 'walkIntersections',
-- where the answer is 'Void', the left branch is erased outright.
{-# INLINE foldCorridor #-}

walkIntersections :: Triangulation mode vertex directed undirected face -> Point -> Point -> Intersection -> [Intersection]
walkIntersections triangulation from to first =
  case foldCorridor triangulation from to first collect [] of
    Nothing -> exactIntersectionScan triangulation from to
    Just (Left impossible) -> absurd impossible
    Just (Right events) -> reverse events
 where
  collect :: [Intersection] -> Intersection -> Either Void [Intersection]
  collect accumulated event = Right (event : accumulated)

nextIntersection :: Triangulation mode vertex directed undirected face -> Point -> Point -> Intersection -> Maybe Intersection
nextIntersection triangulation lineFrom lineTo current = case current of
  EdgeIntersection edge -> case traceDirectionOutOfEdge triangulation edge lineFrom lineTo of
    EdgeOutHull -> Nothing
    EdgeOutVertex vertex -> Just (VertexIntersection vertex)
    EdgeOutEdge nextEdge -> Just (EdgeIntersection nextEdge)
    EdgeOutNone -> Nothing
  VertexIntersection vertex
    | vertexPoint triangulation vertex == lineTo -> Nothing
    | otherwise -> case traceDirectionOutOfVertex triangulation vertex lineTo of
        VertexOutHull -> Nothing
        VertexOutOverlap edge -> Just (EdgeOverlap edge)
        VertexOutEdge edge ->
          let from = vertexPoint triangulation (origin triangulation edge)
              to = vertexPoint triangulation (destination triangulation edge)
           in if orient2d from to lineTo == LT then Nothing else Just (EdgeIntersection edge)
  EdgeOverlap edge
    | lineFrom == lineTo -> Nothing
    | onClosedSegment lineFrom lineTo (vertexPoint triangulation (destination triangulation edge)) ->
        Just (VertexIntersection (destination triangulation edge))
    | otherwise -> Nothing

firstIntersection :: Triangulation mode vertex directed undirected face -> QueryPoint -> QueryPoint -> Maybe Intersection
firstIntersection triangulation queryFrom queryTo =
  case locatePoint triangulation queryFrom of
    EmptyTriangulation -> singleVertexHit
    OnVertex vertex -> Just (VertexIntersection vertex)
    OnEdge edge -> Just (classifyStartingEdge edge)
    InFace face -> firstFromFace face
    OutsideConvexHull entry -> firstFromOutside entry
 where
  !lineFrom = queryPointValue queryFrom
  !lineTo = queryPointValue queryTo
  singleVertexHit = case vertices triangulation of
    [vertex]
      | onClosedSegment lineFrom lineTo (vertexPoint triangulation vertex) -> Just (VertexIntersection vertex)
    _ -> Nothing

  classifyStartingEdge edge =
    let a = vertexPoint triangulation (origin triangulation edge)
        b = vertexPoint triangulation (destination triangulation edge)
     in if orient2d lineFrom lineTo a == EQ && orient2d lineFrom lineTo b == EQ
          then EdgeOverlap (orientAlongLine edge)
          else EdgeIntersection (orientTowardTarget edge)

  firstFromFace face = firstEdgeFromRing (faceDirectedEdges triangulation face)

  firstEdgeFromRing [] = Nothing
  firstEdgeFromRing (edge : remaining) =
    let a = vertexPoint triangulation (origin triangulation edge)
        b = vertexPoint triangulation (destination triangulation edge)
     in if segmentIntersectsNonCollinear lineFrom lineTo a b
          then
            if orient2d lineFrom lineTo a == EQ
              then Just (VertexIntersection (origin triangulation edge))
              else
                if orient2d lineFrom lineTo b == EQ
                  then Just (VertexIntersection (destination triangulation edge))
                  else Just (EdgeIntersection (reverseEdge edge))
          else firstEdgeFromRing remaining

  -- Outside the region the segment's first contact with it lies on the ring
  -- the locator's edge sits on, so only that ring can carry the earliest
  -- event. A ring that outran its budget, an absent locator edge, and
  -- endpoints that leave the parameter comparison without a total order all
  -- keep the exact scan.
  firstFromOutside (Just edge)
    | finiteEndpoint lineFrom && finiteEndpoint lineTo =
        case ringEntryEvent triangulation lineFrom lineTo edge of
          Just entry -> eventValue <$> entry
          Nothing -> firstFromScan
  firstFromOutside _ = firstFromScan

  firstFromScan = case exactIntersectionScan triangulation lineFrom lineTo of
    event : _ -> Just event
    [] -> Nothing

  orientAlongLine edge =
    let a = vertexPoint triangulation (origin triangulation edge)
        b = vertexPoint triangulation (destination triangulation edge)
     in if projectionFactor lineFrom lineTo a <= projectionFactor lineFrom lineTo b then edge else reverseEdge edge

  orientTowardTarget edge =
    let a = vertexPoint triangulation (origin triangulation edge)
        b = vertexPoint triangulation (destination triangulation edge)
     in if orient2d a b lineTo == LT then reverseEdge edge else edge

data VertexOut
  = VertexOutHull
  | VertexOutOverlap !DirectedEdgeId
  | VertexOutEdge !DirectedEdgeId

data EdgeOut
  = EdgeOutHull
  | EdgeOutVertex !VertexId
  | EdgeOutEdge !DirectedEdgeId
  | EdgeOutNone

traceDirectionOutOfVertex :: Triangulation mode vertex directed undirected face -> VertexId -> Point -> VertexOut
traceDirectionOutOfVertex triangulation vertex target =
  case vertexOutEdge triangulation vertex of
    Nothing -> VertexOutHull
    Just start ->
      let !startSide = sideOf start
          !rotateCounterClockwise = startSide == GT
       in go rotateCounterClockwise (numDirectedEdges triangulation + 1) start startSide
 where
  go !rotateCounterClockwise !remaining !current !currentSide
    | remaining <= 0 = VertexOutHull
    | currentSide == EQ && projectionFactor currentPoint (edgeTarget current) target >= 0 =
        VertexOutOverlap current
    | otherwise =
        let following = if rotateCounterClockwise then counterClockwise triangulation current else clockwise triangulation current
            followingSide = sideOf following
         in if followingSide == EQ && projectionFactor currentPoint (edgeTarget following) target >= 0
              then VertexOutOverlap following
              else
                let faceBetween = if rotateCounterClockwise then incidentFace triangulation current else incidentFace triangulation following
                 in if faceBetween == outerFace
                      then VertexOutHull
                      else
                        if rotateCounterClockwise == (followingSide == LT)
                          then
                            let segment = if rotateCounterClockwise then next triangulation current else previous triangulation (reverseEdge current)
                             in VertexOutEdge (reverseEdge segment)
                          else go rotateCounterClockwise (remaining - 1) following followingSide

  currentPoint = vertexPoint triangulation vertex
  edgeTarget edge = vertexPoint triangulation (destination triangulation edge)
  sideOf edge = orient2d currentPoint (edgeTarget edge) target

traceDirectionOutOfEdge :: Triangulation mode vertex directed undirected face -> DirectedEdgeId -> Point -> Point -> EdgeOut
traceDirectionOutOfEdge triangulation edge lineFrom lineTo
  | incidentFace triangulation edge == outerFace = EdgeOutHull
  | otherwise =
      case (previousIntersects, nextIntersects) of
        (True, False) -> EdgeOutEdge (reverseEdge edgePrevious)
        (False, True) -> EdgeOutEdge (reverseEdge edgeNext)
        (True, True) -> EdgeOutVertex (origin triangulation edgePrevious)
        (False, False) -> EdgeOutNone
 where
  edgePrevious = previous triangulation edge
  edgeNext = next triangulation edge

  -- The face runs @edge@ A->B, @edgeNext@ B->C, @edgePrevious@ C->A, so the
  -- two candidates share C and each vertex's side of the line is read once
  -- rather than once per candidate. Whether the segment reaches a candidate is
  -- then asked only of one the line already separates, and a walk that entered
  -- across A->B leaves that true for exactly one of the two.
  pointA = vertexPoint triangulation (origin triangulation edge)
  pointB = vertexPoint triangulation (origin triangulation edgeNext)
  pointC = vertexPoint triangulation (origin triangulation edgePrevious)

  sideA = orient2d lineFrom lineTo pointA
  sideB = orient2d lineFrom lineTo pointB
  sideC = orient2d lineFrom lineTo pointC

  previousIntersects = sideC /= sideA && reaches pointC pointA
  nextIntersects = sideB /= sideC && reaches pointB pointC

  reaches from to = orient2d from to lineFrom /= orient2d from to lineTo

-- | Whether two segments properly cross; collinear touching does not.
segmentIntersectsNonCollinear :: Point -> Point -> Point -> Point -> Bool
segmentIntersectsNonCollinear p0 p1 q0 q1 =
  -- Equality admits an endpoint on the opposite segment; four equal sides
  -- reject the collinear case without recomputing either orientation pair.
  p0Side /= p1Side
    && q0Side /= q1Side
 where
  !p0Side = orient2d q0 q1 p0
  !p1Side = orient2d q0 q1 p1
  !q0Side = orient2d p0 p1 q0
  !q1Side = orient2d p0 p1 q1
{-# INLINE segmentIntersectsNonCollinear #-}

exactIntersectionScan :: Triangulation mode vertex directed undirected face -> Point -> Point -> [Intersection]
exactIntersectionScan triangulation lineFrom lineTo =
  map eventValue . sortBy compareEvent $ vertexEvents ++ edgeEvents
 where
  vertexEvents = mapMaybe (vertexEvent triangulation lineFrom lineTo) (vertices triangulation)
  edgeEvents = mapMaybe (edgeEvent triangulation lineFrom lineTo) (undirectedEdges triangulation)

data Event = Event
  { eventParameter :: !Double
  , eventPriority :: {-# UNPACK #-} !Int
  , eventValue :: !Intersection
  }

compareEvent :: Event -> Event -> Ordering
compareEvent = comparing (\event -> (eventParameter event, eventPriority event, eventValue event))

-- Ties in this order are the same 'Intersection', so a traversal that reaches
-- the minimum by a different route still reports the value the sorted scan's
-- head reports.
vertexEvent
  :: Triangulation mode vertex directed undirected face
  -> Point
  -> Point
  -> VertexId
  -> Maybe (Event)
vertexEvent triangulation lineFrom lineTo vertex
  | onClosedSegment lineFrom lineTo point = Just (Event (projectionFactor lineFrom lineTo point) 0 (VertexIntersection vertex))
  | otherwise = Nothing
 where
  point = vertexPoint triangulation vertex

edgeEvent
  :: Triangulation mode vertex directed undirected face
  -> Point
  -> Point
  -> UndirectedEdgeId
  -> Maybe (Event)
edgeEvent triangulation lineFrom lineTo edge
  | oa == EQ && ob == EQ =
      if low < high
        then Just (Event low 1 (EdgeOverlap (if projectedA <= projectedB then directed else reverseEdge directed)))
        else Nothing
  | segmentsIntersect lineFrom lineTo a b && oa /= EQ && ob /= EQ =
      Just (Event (segmentIntersectionParameter lineFrom lineTo a b) 2 (EdgeIntersection oriented))
  | otherwise = Nothing
 where
  directed = normalizedDirected edge
  a = vertexPoint triangulation (origin triangulation directed)
  b = vertexPoint triangulation (destination triangulation directed)
  oa = orient2d lineFrom lineTo a
  ob = orient2d lineFrom lineTo b
  projectedA = projectionFactor lineFrom lineTo a
  projectedB = projectionFactor lineFrom lineTo b
  low = max 0 (min projectedA projectedB)
  high = min 1 (max projectedA projectedB)
  oriented = if orient2d a b lineTo == LT then reverseEdge directed else directed

-- | The earliest event carried by the outer-face ring the given edge sits on,
-- walked from that edge. 'Nothing' reports a ring that outran its budget.
ringEntryEvent
  :: Triangulation mode vertex directed undirected face
  -> Point
  -> Point
  -> DirectedEdgeId
  -> Maybe (Maybe (Event))
ringEntryEvent triangulation lineFrom lineTo entry =
  go (numDirectedEdges triangulation + 1) entry False Nothing
 where
  go !remaining !edge !departed !earliest
    | remaining <= 0 = Nothing
    | departed && edge == entry = Just earliest
    | otherwise =
        let !stepped =
              keepEarliest (edgeEvent triangulation lineFrom lineTo (asUndirected edge)) $
                keepEarliest (vertexEvent triangulation lineFrom lineTo (origin triangulation edge)) earliest
         in go (remaining - 1) (next triangulation edge) True stepped

keepEarliest :: Maybe (Event) -> Maybe (Event) -> Maybe (Event)
keepEarliest Nothing held = held
keepEarliest candidate Nothing = candidate
keepEarliest candidate@(Just proposed) held@(Just incumbent)
  | compareEvent proposed incumbent == LT = candidate
  | otherwise = held

finiteEndpoint :: Point -> Bool
finiteEndpoint point = isFinite (pointX point) && isFinite (pointY point)

segmentIntersectionParameter :: Point -> Point -> Point -> Point -> Double
segmentIntersectionParameter (Point ax ay) (Point bx by) (Point cx cy) (Point dx dy)
  | denominator == 0 = 0
  | otherwise = ((cx - ax) * (dy - cy) - (cy - ay) * (dx - cx)) / denominator
 where
  denominator = (bx - ax) * (dy - cy) - (by - ay) * (dx - cx)