hgeometry-0.4.0.0: src/Algorithms/Geometry/ConvexHull/GrahamScan.hs
module Algorithms.Geometry.ConvexHull.GrahamScan( ConvexHull(..)
, DegenerateCH
, convexHull
, upperHull
, lowerHull
) where
import Control.Lens((^.))
import Data.Ext
import Data.Geometry.Point
import Data.Geometry.Polygon
import qualified Data.List as L
import Data.Monoid
-- | Two dimensional convex hulls
newtype ConvexHull p r = ConvexHull { _hull :: (SimplePolygon p r) }
-- deriving (Show,Eq)
type DegenerateCH p r = Maybe (Point 2 r :+ p)
-- | O(n log n) time ConvexHull using Graham-Scan
convexHull :: (Ord r, Num r)
=> [Point 2 r :+ p] -> Either (DegenerateCH p r) (ConvexHull p r)
convexHull [] = Left Nothing
convexHull [p] = Left (Just p)
convexHull ps = let ps' = L.sortBy incXdecY ps
uh = tail . hull' $ ps'
lh = tail . hull' $ reverse ps'
in Right . ConvexHull . fromPoints $ lh ++ uh
upperHull :: (Ord r, Num r)
=> [Point 2 r :+ p] -> Either (DegenerateCH p r) [Point 2 r :+ p]
upperHull = hull id
lowerHull :: (Ord r, Num r)
=> [Point 2 r :+ p] -> Either (DegenerateCH p r) [Point 2 r :+ p]
lowerHull = hull reverse
-- | Helper function so that that can compute both the upper or the lower hull, depending
-- on the function f
hull :: (Ord r, Num r)
=> ([Point 2 r :+ p] -> [Point 2 r :+ p])
-> [Point 2 r :+ p]
-> Either (DegenerateCH p r) [Point 2 r :+ p]
hull _ [] = Left Nothing
hull _ [p] = Left (Just p)
hull f ps = Right . hull' . f . L.sortBy incXdecY $ ps
incXdecY :: Ord a => Ext t (Point 2 a) -> Ext t1 (Point 2 a) -> Ordering
incXdecY (Point2 px py :+ _) (Point2 qx qy :+ _) =
compare px qx <> compare qy py
-- | Precondition: The list of input points is sorted
hull' :: (Ord r, Num r) => [Point 2 r :+ p] -> [Point 2 r :+ p]
hull' (a:b:ps) = hull'' [b,a] ps
where
hull'' h [] = h
hull'' h (p:ps) = hull'' (cleanMiddle (p:h)) ps
cleanMiddle [b,a] = [b,a]
cleanMiddle h@(c:b:a:rest)
| rightTurn (a^.core) (b^.core) (c^.core) = h
| otherwise = cleanMiddle (c:a:rest)
rightTurn :: (Ord r, Num r) => Point 2 r -> Point 2 r -> Point 2 r -> Bool
rightTurn a b c = ccw a b c == CW