module Algorithms.Geometry.LineSegmentIntersection.BentleyOttmann where
import Algorithms.Geometry.LineSegmentIntersection.Types
import Control.Lens hiding (contains)
import qualified Data.BalBST as SS
import Data.Ext
import Data.Function (on)
import Data.Geometry.Interval
import Data.Geometry.Line
import Data.Geometry.LineSegment
import Data.Geometry.Point
import Data.Geometry.Properties
import qualified Data.List as L
import Data.List.NonEmpty (NonEmpty(..))
import qualified Data.List.NonEmpty as NonEmpty
import qualified Data.Map as M
import Data.Maybe
import Data.Ord (Down(..), comparing)
import Data.Semigroup
import qualified Data.Set as EQ
import Data.Vinyl
import Frames.CoRec
import Debug.Trace
intersections :: (Ord r, Fractional r)
=> [LineSegment 2 p r] -> Intersections p r
intersections ss = merge $ sweep pts (SS.empty $ ordAtNav undefined)
where
pts = EQ.fromAscList . groupStarts . L.sort . concatMap f $ ss
f s = let [p,q] = L.sortBy ordPoints [s^.start.core,s^.end.core]
in [Event p (Start $ s :| []), Event q (End s)]
merge :: Ord r => [IntersectionPoint p r] -> Intersections p r
merge = foldr (\(IntersectionPoint p a) -> M.insertWith (<>) p a) M.empty
groupStarts :: Eq r => [Event p r] -> [Event p r]
groupStarts [] = []
groupStarts (Event p (Start s) : es) = Event p (Start ss) : groupStarts rest
where
(ss',rest) = L.span sameStart es
ss = let (x:|xs) = s in x :| (xs ++ concatMap startSegs ss')
sameStart (Event q (Start _)) = p == q
sameStart _ = False
groupStarts (e : es) = e : groupStarts es
data EventType s = Start !(NonEmpty s)| Intersection | End !s deriving (Show)
instance Eq (EventType s) where
a == b = a `compare` b == EQ
instance Ord (EventType s) where
(Start _) `compare` (Start _) = EQ
(Start _) `compare` _ = LT
Intersection `compare` (Start _) = GT
Intersection `compare` Intersection = EQ
Intersection `compare` (End _) = LT
(End _) `compare` (End _) = EQ
(End _) `compare` _ = GT
data Event p r = Event { eventPoint :: !(Point 2 r)
, eventType :: !(EventType (LineSegment 2 p r))
} deriving (Show,Eq)
instance Ord r => Ord (Event p r) where
(Event p s) `compare` (Event q t) = case ordPoints p q of
EQ -> s `compare` t
x -> x
ordPoints :: Ord r => Point 2 r -> Point 2 r -> Ordering
ordPoints a b = let f p = (Down $ p^.yCoord, p^.xCoord) in comparing f a b
startSegs :: Event p r -> [LineSegment 2 p r]
startSegs e = case eventType e of
Start ss -> NonEmpty.toList ss
_ -> []
ordAtNav :: (Ord r, Fractional r) => r -> SS.TreeNavigator r (LineSegment 2 p r)
ordAtNav y = SS.Nav (\s x -> h s <= x) (min `on` h)
where
h s = match (s `intersect` horizontalLine y) $
(H $ \NoIntersection -> error "ordAtNav: No intersection")
:& (H $ \p -> p^.xCoord)
:& (H $ \_ -> rightEndpoint s)
:& RNil
type EventQueue p r = EQ.Set (Event p r)
type StatusStructure p r = SS.BalBST r (LineSegment 2 p r)
sweep :: (Ord r, Fractional r)
=> EventQueue p r -> StatusStructure p r -> [IntersectionPoint p r]
sweep eq ss = case EQ.minView eq of
Nothing -> []
Just (e,eq') -> handle e eq' ss
handle :: (Ord r, Fractional r)
=> Event p r -> EventQueue p r -> StatusStructure p r
-> [IntersectionPoint p r]
handle e@(eventPoint -> p) eq ss = toReport <> sweep eq' ss'
where
starts = startSegs e
(before,contains',after) = extractContains p ss
(ends,contains) = L.partition (endsAt p) contains'
toReport = case starts ++ contains' of
(_:_:_) -> [IntersectionPoint p $ associated (starts <> ends) contains]
_ -> []
ss' = before `SS.join` newSegs `SS.join` after
newSegs = toStatusStruct p $ starts ++ contains
eq' = foldr EQ.insert eq es
es | SS.null newSegs = maybeToList $ app (findNewEvent p) sl sr
| otherwise = let s' = fst <$> SS.minView newSegs
s'' = fst <$> SS.maxView newSegs
in catMaybes [ app (findNewEvent p) sl s'
, app (findNewEvent p) s'' sr
]
sl = fst <$> SS.maxView before
sr = fst <$> SS.minView after
app f x y = do { x' <- x ; y' <- y ; f x' y'}
extractContains :: (Fractional r, Ord r)
=> Point 2 r -> StatusStructure p r
-> (StatusStructure p r, [LineSegment 2 p r], StatusStructure p r)
extractContains p ss = (before, mid1 ++ mid2, after)
where
n = ordAtNav (p^.yCoord)
SS.Split before (mid1,mid2) after = SS.splitExtract pred' sel $ ss { SS.nav = n}
pred' s = not $ SS.goLeft n s (p^.xCoord)
sel s = p `onSegment` s
toStatusStruct :: (Fractional r, Ord r)
=> Point 2 r -> [LineSegment 2 p r] -> StatusStructure p r
toStatusStruct p xs = ss { SS.nav = ordAtNav $ p^.yCoord } `SS.join` hors
where
(hors',rest) = L.partition isHorizontal xs
ss = SS.fromList (ordAtNav $ maxY xs) rest
hors = SS.fromList (SS.ordNavBy rightEndpoint) hors'
isHorizontal s = s^.start.core.yCoord == s^.end.core.yCoord
maxY = maximum . filter (< p^.yCoord)
. concatMap (\s -> [s^.start.core.yCoord,s^.end.core.yCoord])
rightEndpoint :: Ord r => LineSegment 2 p r -> r
rightEndpoint s = (s^.start.core.xCoord) `max` (s^.end.core.xCoord)
endsAt :: Ord r => Point 2 r -> LineSegment 2 p r -> Bool
endsAt p (LineSegment' a b) = all (\q -> ordPoints (q^.core) p /= GT) [a,b]
findNewEvent :: (Ord r, Fractional r)
=> Point 2 r -> LineSegment 2 p r -> LineSegment 2 p r
-> Maybe (Event p r)
findNewEvent p l r = match (l `intersect` r) $
(H $ \NoIntersection -> Nothing)
:& (H $ \q -> if ordPoints q p == GT then Just (Event q Intersection)
else Nothing)
:& (H $ \_ -> Nothing)
:& RNil