removing duplicate code

hgeometry/src/HGeometry/LowerEnvelope/Connected/FromVertexForm.hs

@@ -235,8 +235,7 @@ sortAlongBoundary face = case mv0 of
         u' = projectPoint $ ivLoc u
         v' = projectPoint $ ivLoc v
 
--- | given an edge (u,v) that has h to its left, and all remaining vertices of the face,
--- sorted in CCW order around the face starting with *v*, compute all its edgesf
+-- | given the vertices of the face in CCW order, compute all its edges of the face.
 faceToEdges       :: forall plane r. (Plane_ plane r, Ord r, Fractional r, Ord plane
                                      , Show plane, Show r
                                      )
@@ -297,14 +296,14 @@ twoVertices (Vtx h ui up uDefs) (Vtx _ vi vp vDefs) =
   case extractEdgeDefs h up uDefs vp vDefs >>= withUVOrder of
       Nothing  -> error "twoVertices. absurd, h and h' don't define an edge!?"
       Just (EdgeDefs hr hu hv, uBeforeV) -- -> traceShowWith ("twoVertices",u,v,uBeforeV,) $
-        | uBeforeV -> NonEmpty.fromList [ (vi, Edge ui                h hr)
-                                        , (ui, Edge unboundedVertexId h hu)
-                                        ]
+        | uBeforeV  -> NonEmpty.fromList [ (vi, Edge ui                h hr)
+                                         , (ui, Edge unboundedVertexId h hu)
+                                         ]
                       -- the edge must be oriented from v to u so that h is on the left
-        | otherwise  -> NonEmpty.fromList [ (ui, Edge vi                h hr)
-                                          , (vi, Edge unboundedVertexId h hv)
-                                          ]
-                        -- the edge must be oriented from u to v
+        | otherwise -> NonEmpty.fromList [ (ui, Edge vi                h hr)
+                                         , (vi, Edge unboundedVertexId h hv)
+                                         ]
+                       -- the edge must be oriented from u to v
   where
     -- determine if u lies before v in the order of u and v along the intersection line of
     -- h and hr.
@@ -315,64 +314,6 @@ twoVertices (Vtx h ui up uDefs) (Vtx _ vi vp vDefs) =
               -- v if v does not lie in the left haflpalne
       ) <$> intersectionLine' h hr
 
-
---------------------------------------------------------------------------------
-
--- | Vertices in of the lower envelope in adjacencylist form.
-type BoundedVertex = BoundedVertexF Seq.Seq
-
-
--- | Given a location of a vertex v, a pair of planes h1,h2 and the
--- remaining defining planes of v, computes the outgoing half-line
--- from v on which h1,h2 are the lowest (if such an halfline exists).
-outgoingUnboundedEdge                   :: ( Plane_ plane r, Ord r, Fractional r
-                                           , Foldable1 f
-                                           )
-                                        => Point 3 r -- ^ the location of the vertex v
-                                        -> Two plane -- ^ the pair of planes for which to compute
-                                        -- the halfine
-                                        -> f plane -- ^ the other planes intersecting at v
-                                        -> Maybe (HalfLine (Point 2 r) :+ EdgeDefiners plane)
-outgoingUnboundedEdge v (Two h1 h2) h3s =
-  intersectionLineWithDefiners h1 h2 >>= toHalfLineFrom (projectPoint v) h3s
-  -- todo, if there are more planes, I guess we should check if the hl is not dominated by the other
-  -- planes either.
-
--- | Given :
---
--- v : the projected location of the vertex
--- hs : the remaining planes defining v (typically just one plane h3)
--- l   : the (projection of the) line l in which planes h1 and h2 intersect (containing v)
---
--- we compute the half-line eminating from v in which h1 and h2 define
--- an edge incident to v.
-toHalfLineFrom                  :: (Plane_ plane r, Foldable1 f, Fractional r, Ord r)
-                                => Point 2 r -- ^ vertex v
-                                -> f plane     -- ^ the remaining plane(s) hs
-                                -> LinePV 2 r :+ EdgeDefiners plane -- ^ the line l
-                                -> Maybe (HalfLine (Point 2 r) :+ EdgeDefiners plane)
-toHalfLineFrom v hs ((LinePV _ w) :+ defs@(EdgeDefiners h1 h2)) =
-    validate w defs <|> validate ((-1) *^ w) (EdgeDefiners h2 h1)
-    -- We try both directions. Note that if we reverse the direction
-    -- of the line, what the left/right plane is changes.
-    --
-    -- If neither direction works, then h1,h2 do not define a good
-    -- direction. This should happen only when there are more than 3
-    -- planes intersecting in v, i.e. in degernate sitatuions
-  where
-    -- | test if direction d is a good direction, i.e. if towards direction d
-    -- h1 (and thus also h2) is actually lower than all remaining defining planes.
-    validate d defs' = let zVal = evalAt (v .+^ d)
-                       in if all (\h3 -> zVal h1 < zVal h3) hs
-                          then Just (HalfLine v d :+ defs') else Nothing
-
-
-
-
-
-
-
-
 --------------------------------------------------------------------------------
 -- * Convenience functions
 
@@ -394,7 +335,7 @@ maxBy cmp a b = case cmp a b of
                   _  -> a
 
 
--- | shift the vector by d
+-- | shift the vector by d items to the right
 shiftR     :: Int -> NonEmptyVector v -> NonEmptyVector v
 shiftR d v = let n = length v
              in NonEmptyV.generate1 n $ \i -> v NonEmptyV.! ((i+n+d) `mod` n)

hgeometry/src/HGeometry/LowerEnvelope/Connected/Type.hs

@@ -17,6 +17,7 @@ module HGeometry.LowerEnvelope.Connected.Type
 
   , singleton
 
+  , BoundedVertex
   , BoundedVertexF(Vertex)
   , location, definers, location2
 
@@ -28,7 +29,7 @@ module HGeometry.LowerEnvelope.Connected.Type
   , projectedEdgeGeometries, projectedEdgeGeometry
 
 
-
+  , outgoingUnboundedEdge
   , edgeIntersectionLine
   , intersectionLine'
   , intersectionLineWithDefiners
@@ -84,16 +85,16 @@ instance (Ord (NumType plane), Num (NumType plane)) => IsBoxable (LowerEnvelope'
   boundingBox env = boundingBox . NonEmpty.fromList $ env^..boundedVertices.traverse.location
     -- the fromList is safe since there is alwasy at least one vertex
 
--- instance Functor LowerEnvelope' where
--- instance Foldable LowerEnvelope' where
--- instance Traversable LowerEnvelope' where
-
--- | Traversal of the planes in the lower envelope
+-- | Traversal of the planes in the lower envelope. Since we need an Ord constraint we
+-- can't make LowerEnvelope' an instance of Traversable.
+--
+-- Be aware that this may destroy some of the invariants. So use this function with care.
 traverseLowerEnvelope                         :: ( Applicative f, NumType plane ~ NumType plane'
                                                  , Ord plane'
                                                  )
                                               => (plane -> f plane')
-                                              -> LowerEnvelope' plane -> f (LowerEnvelope' plane')
+                                              -> LowerEnvelope' plane
+                                              -> f (LowerEnvelope' plane')
 traverseLowerEnvelope f (LowerEnvelope v0 vs) =
     LowerEnvelope <$> traverse f v0 <*> traverse (traverseBoundedV f) vs
 
@@ -121,18 +122,6 @@ singleton v = LowerEnvelope v0 (Seq.singleton v')
 
 --------------------------------------------------------------------------------
 
-
---------------------------------------------------------------------------------
-
-
--- | Helper type for edgedefs
-data EdgeDefs plane = EdgeDefs { common :: plane
-                               , uNeigh :: plane
-                               , vNeigh :: plane
-                               } deriving (Show,Eq)
-
---------------------------------------------------------------------------------
-
 -- | The unbounded vertex.
 newtype UnboundedVertex plane = UnboundedVertex { _incidentEdgesU :: Seq.Seq (LEEdge plane) }
                               deriving (Show,Eq,Functor,Foldable,Traversable)
@@ -270,8 +259,9 @@ edgeIntersectionLine e = case intersectionLine' (e^.leftPlane) (e^.rightPlane) o
 
 --------------------------------------------------------------------------------
 
--- |
+-- | Types that have an IncidentEdges' field.
 class HasIncidentEdges t  where
+  -- | Lens to access the incident edges field.
   incidentEdges' :: Lens' (t plane) (Seq.Seq (LEEdge plane))
 
 instance HasIncidentEdges UnboundedVertex where
@@ -450,17 +440,4 @@ projectedEdgeGeometry env (ui,vi) e = case env^?!vertexAt ui of
     -- the unbounded vertex into the bounded vertex that means that to construct the
     -- halfline we actually wish to flip the orientation of the line
 
-
-
---       case
--- --      traceShowWith ("unbEdge", e, "h=",e^.rightPlane, "h'=",e^.leftPlane, "line'",) $
---         intersectionLine' (e^.rightPlane) (e^.leftPlane) of
---       Just l  -> Left $ HalfLine (v^.location2) (l^.direction)
---                  -- edge e is oriented from the lowerId towards the higherId,
---                  -- so in particular *from* the unbounded vertex into the bounded vertex
---                  -- that means that to construct the halfline we actually wish to
---                  -- flip the left and right plane, so that the halfline is directed outwards.
---       Nothing -> error "projectedEdgeGeometry: absurd, no intersection between the planes"
-
-
 --------------------------------------------------------------------------------

hgeometry/src/HGeometry/LowerEnvelope/Naive.hs

@@ -10,22 +10,14 @@ module HGeometry.LowerEnvelope.Naive
 
 --------------------------------------------------------------------------------
 
-import           Control.Lens
 import           Control.Monad (guard)
-import qualified Data.Sequence as Seq
 import qualified Data.Set as Set
 import           HGeometry.Combinatorial.Util
-import           HGeometry.Foldable.Sort
 import           HGeometry.HyperPlane.Class
 import           HGeometry.HyperPlane.NonVertical
-import           HGeometry.Intersection
-import           HGeometry.Line
-import           HGeometry.Line.LineEQ
 import           HGeometry.LowerEnvelope.AdjListForm (LowerEnvelope, fromVertexForm)
 import           HGeometry.LowerEnvelope.VertexForm
 import           HGeometry.Point
-import           HGeometry.Properties
-import           Hiraffe.Graph
 
 --------------------------------------------------------------------------------
 
@@ -61,10 +53,7 @@ lowerEnvelope = fromVertexForm . lowerEnvelopeVertexForm
 lowerEnvelopeVertexForm    :: ( Plane_ plane r
                               , Ord r, Fractional r, Foldable f, Ord plane
                               ) => f plane -> VertexForm plane
-lowerEnvelopeVertexForm hs = foldMap (\t -> case asVertex hs t of
-                                              Nothing -> mempty
-                                              Just v  -> singleton v
-                                     ) $ uniqueTriplets hs
+lowerEnvelopeVertexForm hs = foldMap (maybe mempty singleton . asVertex hs) $ uniqueTriplets hs
 
 -- | Given all planes, and a triple, computes if the triple defines a
 -- vertex of the lower envelope, and if so returns it.
@@ -75,8 +64,8 @@ asVertex hs t@(Three h1 h2 h3) = do v <- intersectionPoint t
                                     pure $ LEVertex v (Set.fromList [h1,h2,h3])
 
 -- | test if v lies below (or on) all the planes in hs
-belowAll      :: (Plane_ plane r, Ord r, Num r, Foldable f) => Point 3 r -> f plane -> Bool
-belowAll v hs = all (\h -> onSideTest v h /= GT) hs
+belowAll   :: (Plane_ plane r, Ord r, Num r, Foldable f) => Point 3 r -> f plane -> Bool
+belowAll v = all (\h -> onSideTest v h /= GT)
 {-# INLINE belowAll #-}
 
 

hgeometry/src/HGeometry/VoronoiDiagram/ViaLowerEnvelope.hs

@@ -20,7 +20,6 @@ module HGeometry.VoronoiDiagram.ViaLowerEnvelope
 
 import           Control.Lens
 import           Data.Default.Class
-import qualified Data.Map as Map
 import           HGeometry.Box
 import           HGeometry.Duality
 import           HGeometry.Ext
@@ -31,7 +30,6 @@ import           HGeometry.LowerEnvelope.Naive (lowerEnvelopeVertexForm)
 import           HGeometry.LowerEnvelope.VertexForm (VertexForm)
 import           HGeometry.Point
 import           HGeometry.Properties
-import           HGeometry.Vector
 import           Hiraffe.Graph
 
 --------------------------------------------------------------------------------

todo.org

@@ -105,11 +105,20 @@ intersecting halflines with boxes seems to go wrong somehow.
 *** TODO plane graph
 
 * TODO 3d-lower-envelope/convex hull
-** DONE naive
+** TODO naive
+*** TODO triangulated envelope
+*** TODO handle degeneracies
+*** TODO handle all colinear
+*** TODO cyclic sorting of the edges
+
 ** TODO define tests
-*** TODO correctly render a lower envelope/vd with 1 vertex and 3 unbounded edges
-*** TODO correctly render bounded edges of some larger set of points
-*** TODO correctly render unbounded edges of some larger set of points
+*** DONE correctly render a lower envelope/vd with 1 vertex and 3 unbounded edges
+*** DONE correctly render bounded edges of some larger set of points
+*** DONE correctly render unbounded edges of some larger set of points
+*** TODO properyt test that every (bounded) face is convex
+
+** TODO some sort of benchmarking for the naive algorithm
+
 
 ** TODO Set3 type to clean up and/or speed up the fromVertexForm code ?