diff --git a/README.md b/README.md index a3294b983..d9c020d9e 100644 --- a/README.md +++ b/README.md @@ -1,7 +1,7 @@ # HGeometry ![GitHub Workflow Status]( -https://img.shields.io/github/actions/workflow/status/noinia/hgeometry/haskell-ci.yml?branch=hgeom1_again) +https://img.shields.io/github/actions/workflow/status/noinia/hgeometry/haskell-ci.yml?branch=master) [![Hackage](https://img.shields.io/hackage/v/hgeometry.svg?color=success)](https://hackage.haskell.org/package/hgeometry) [![API docs coverage](https://img.shields.io/endpoint?url=https%3A%2F%2Fnoinia.github.io%2Fhgeometry%2Fhaddock_badge.json)](https://noinia.github.io/hgeometry/haddocks) diff --git a/hgeometry-combinatorial/src-quickcheck/HGeometry/Combinatorial/Instances.hs b/hgeometry-combinatorial/src-quickcheck/HGeometry/Combinatorial/Instances.hs index f36e00b78..21f0b02d6 100644 --- a/hgeometry-combinatorial/src-quickcheck/HGeometry/Combinatorial/Instances.hs +++ b/hgeometry-combinatorial/src-quickcheck/HGeometry/Combinatorial/Instances.hs @@ -27,7 +27,7 @@ import HGeometry.Tree.Binary.Static instance (Arbitrary c, Arbitrary e) => Arbitrary (c :+ e) where arbitrary = (:+) <$> arbitrary <*> arbitrary - + shrink = genericShrink -------------------------------------------------------------------------------- @@ -38,9 +38,15 @@ instance (Arbitrary c, Arbitrary e) => Arbitrary (c :+ e) where instance (Arbitrary a, Num a, Eq a) => Arbitrary (GRatio a) where arbitrary = (/) <$> arbitrary <*> (arbitrary `suchThat` (/= 0)) + shrink r = 0 : 1 : [ a' % b' + | a' <- shrink $ numerator r + , b' <- fromInteger 1 : shrink (denominator r) + , b' /= 0 + ] instance KnownNat p => Arbitrary (RealNumber p) where arbitrary = fromFixed <$> arbitrary + shrink = genericShrink instance Arbitrary Sign.Sign where arbitrary = (\b -> if b then Sign.Positive else Sign.Negative) <$> arbitrary @@ -69,6 +75,7 @@ instance (Arbitrary a, Arbitrary v) => Arbitrary (BinLeafTree v a) where | otherwise = do l <- choose (0,n-1) Node <$> f l <*> arbitrary <*> f (n-l-1) + -- shrink = genericShrink instance Arbitrary a => Arbitrary (BinaryTree a) where arbitrary = sized f @@ -76,3 +83,4 @@ instance Arbitrary a => Arbitrary (BinaryTree a) where | otherwise = do l <- choose (0,n-1) Internal <$> f l <*> arbitrary <*> f (n-l-1) + -- shrink = genericShrink diff --git a/hgeometry/data/test-with-ipe/golden/PlaneGraph/smallPlaneGraph.ipe b/hgeometry/data/test-with-ipe/golden/PlaneGraph/smallPlaneGraph.ipe index 8988561db..5ae5285d5 100644 --- a/hgeometry/data/test-with-ipe/golden/PlaneGraph/smallPlaneGraph.ipe +++ b/hgeometry/data/test-with-ipe/golden/PlaneGraph/smallPlaneGraph.ipe @@ -133,7 +133,7 @@ -Point2 0.0 0.0 :+ 0Point2 200.0 200.0 :+ 1Point2 200.0 0.0 :+ 2Point2 (-100.0) 400.0 :+ 30.0 1.0 m +Point2 0.0 0.0 :+ 0Point2 200.0 200.0 :+ 1Point2 200.0 0.0 :+ 2Point2 (-100.0) 400.0 :+ 30.0 1.0 m 200.0 1.0 l "0->2"200.0 -1.0 m 0.0 -1.0 l @@ -153,7 +153,17 @@ 201.0 0.0 l "1->2"199.0 0.0 m 199.0 200.0 l -"2->1"0.0 0.0 m +"2->1"0.0 0.0 m +200.0 0.0 l +Dart (Arc 0) +10.0 0.0 m +200.0 200.0 l +Dart (Arc 1) +10.0 0.0 m +-100.0 400.0 l +Dart (Arc 2) +1200.0 200.0 m +-100.0 400.0 l +Dart (Arc 3) +1200.0 200.0 m +200.0 0.0 l +Dart (Arc 4) +10.0 0.0 m 200.0 0.0 l 200.0 200.0 l h diff --git a/hgeometry/hgeometry.cabal b/hgeometry/hgeometry.cabal index 40865b006..5bdf68b01 100644 --- a/hgeometry/hgeometry.cabal +++ b/hgeometry/hgeometry.cabal @@ -70,7 +70,8 @@ common all-setup , file-io >= 0.1 && < 1 , indexed-traversable >= 0.1.3 && < 1 , vector-builder >= 0.3.8 && < 1 - , HsYAML >= 0.2 && < 1 + , HsYAML >= 0.2 && < 1 + , semialign >= 1.3 && < 1.4 , ghc-typelits-natnormalise >= 0.7.7 && < 1 , ghc-typelits-knownnat >= 0.7.6 && < 1 @@ -100,14 +101,12 @@ common point-setup common vector-setup build-depends: - semialign >= 1.2 && < 2 - , these >= 1.1 && < 2 + these >= 1.1 && < 2 common kernel-setup build-depends: hgeometry:vector , hgeometry:point - , semialign common hgeometry-setup build-depends: @@ -381,6 +380,7 @@ library -- HGeometry.HalfPlane.CommonIntersection + HGeometry.Polygon.Triangulation HGeometry.Polygon.Triangulation.TriangulateMonotone HGeometry.Polygon.Triangulation.MakeMonotone @@ -664,6 +664,7 @@ test-suite with-ipe-hspec LineSegmentSpec HalfLineSpec Polygon.Triangulation.TriangulateMonotoneSpec + Polygon.Triangulation.TriangulateSpec LineSegment.Intersection.BentleyOttmannSpec PlaneGraph.RenderSpec diff --git a/hgeometry/kernel/src-quickcheck/HGeometry/Kernel/Instances.hs b/hgeometry/kernel/src-quickcheck/HGeometry/Kernel/Instances.hs index 11c08eb89..1ec4f8463 100644 --- a/hgeometry/kernel/src-quickcheck/HGeometry/Kernel/Instances.hs +++ b/hgeometry/kernel/src-quickcheck/HGeometry/Kernel/Instances.hs @@ -13,6 +13,7 @@ module HGeometry.Kernel.Instances where import Control.Lens hiding (cons) +import Data.Semialign import GHC.TypeLits import HGeometry.Ball import HGeometry.Box @@ -21,7 +22,7 @@ import HGeometry.HalfSpace import HGeometry.HyperPlane (HyperPlane(..)) import HGeometry.HyperPlane.NonVertical (NonVerticalHyperPlane(..)) import HGeometry.Interval -import HGeometry.Interval.EndPoint() +import HGeometry.Interval.EndPoint () import HGeometry.Line.LineEQ import HGeometry.Line.PointAndVector import HGeometry.LineSegment @@ -41,25 +42,39 @@ import Test.QuickCheck instance Arbitrary r => Arbitrary (EndPoint ep r) where arbitrary = EndPoint <$> arbitrary + shrink (EndPoint p) = EndPoint <$> shrink p instance Arbitrary EndPointType where arbitrary = (\b -> if b then Open else Closed) <$> arbitrary + shrink = \case + Open -> [Closed] + Closed -> [] instance Arbitrary r => Arbitrary (AnEndPoint r) where arbitrary = AnEndPoint <$> arbitrary <*> arbitrary - + shrink = genericShrink instance ( Arbitrary (endPoint r) , Eq (endPoint r), Ord r, IxValue (endPoint r) ~ r, EndPoint_ (endPoint r) ) => Arbitrary (Interval endPoint r) where arbitrary = do p <- arbitrary - q <- arbitrary `suchThat` (isValid p) + q <- arbitrary `suchThat` (isValidInterval p) pure $ buildInterval p q - where - isValid p q = p /= q && ((p^._endPoint == q^._endPoint) `implies` bothClosed p q) - bothClosed p q = endPointType p == Closed && endPointType q == Closed - implies p q = not p || q + shrink i = [ buildInterval p q + | p <- shrink $ i^.startPoint + , q <- shrink $ i^.endPoint + , isValidInterval p q + ] + +isValidInterval :: (Eq (endPoint r), Ord r, IxValue (endPoint r) ~ r, EndPoint_ (endPoint r)) + => endPoint r -> endPoint r -> Bool +isValidInterval p q = p /= q && ((p^._endPoint == q^._endPoint) `implies` bothClosed p q) +bothClosed :: EndPoint_ (endPoint r) => endPoint r -> endPoint r -> Bool +bothClosed p q = endPointType p == Closed && endPointType q == Closed + +implies :: Bool -> Bool -> Bool +implies p q = not p || q instance ( Arbitrary (endPoint point) , IsEndPoint (endPoint point) (endPoint point) @@ -69,6 +84,11 @@ instance ( Arbitrary (endPoint point) arbitrary = do p <- arbitrary q <- arbitrary `suchThat` (\q' -> q'^._endPoint /= p^._endPoint) pure $ LineSegment p q + shrink s = [ LineSegment p q + | p <- shrink $ s^.startPoint + , q <- shrink $ s^.endPoint + , q^._endPoint /= p^._endPoint + ] instance ( Arbitrary point , Arbitrary (NumType point) @@ -77,6 +97,11 @@ instance ( Arbitrary point ) => Arbitrary (Ball point) where arbitrary = Ball <$> arbitrary <*> (arbitrary `suchThat` (> 0)) + shrink (Ball c r) = [ Ball c' r' + | c' <- shrink c + , r' <- 1 : shrink r + , r' > 0 + ] instance ( Arbitrary point , Point_ point 2 r, Num r, Ord r @@ -86,6 +111,8 @@ instance ( Arbitrary point b <- arbitrary `suchThat` (/= a) c <- arbitrary `suchThat` (\c' -> c' /= a && c' /= b && ccw a b c' /= CoLinear) pure $ Triangle a b c + shrink = genericShrink + instance Arbitrary r => Arbitrary (LineEQ r) where arbitrary = LineEQ <$> arbitrary <*> arbitrary @@ -103,10 +130,16 @@ instance ( Arbitrary point , Arbitrary r , Point_ point d r , Num r - , Ord (Vector d r) + , Ord r + , Zip (Vector d) ) => Arbitrary (Box point) where arbitrary = (\p v -> Box p (p .+^ v)) <$> arbitrary - <*> arbitrary `suchThat` (> zero) + <*> arbitrary `suchThat` (allOf components (> 0)) + shrink b = [ Box p (p .+^ v) + | p <- shrink $ b^.minPoint + , v <- shrink $ size b + , allOf components (> 0) v + ] instance ( Has_ Additive_ m r , Has_ Vector_ n (Vector m r) diff --git a/hgeometry/kernel/src/HGeometry/Interval/EndPoint.hs b/hgeometry/kernel/src/HGeometry/Interval/EndPoint.hs index 466c06cd5..06bb118d5 100644 --- a/hgeometry/kernel/src/HGeometry/Interval/EndPoint.hs +++ b/hgeometry/kernel/src/HGeometry/Interval/EndPoint.hs @@ -22,6 +22,7 @@ module HGeometry.Interval.EndPoint import Control.Lens import Data.Foldable1 +import GHC.Generics (Generic) import HGeometry.Properties import Text.Read @@ -40,7 +41,7 @@ class IsEndPoint endPoint endPoint => EndPoint_ endPoint where mkEndPoint :: IxValue endPoint -> endPoint -- | Possible endpoint types; open or closed -data EndPointType = Open | Closed deriving (Show,Read,Eq,Ord,Enum,Bounded) +data EndPointType = Open | Closed deriving (Show,Read,Eq,Ord,Enum,Bounded,Generic) -- testV :: Vector 2 (Point 2 Double) @@ -53,7 +54,7 @@ data EndPointType = Open | Closed deriving (Show,Read,Eq,Ord,Enum,Bounded) -- | EndPoint with a type safe tag newtype EndPoint (et :: EndPointType) r = EndPoint r - deriving stock (Eq,Ord,Functor,Foldable,Traversable) + deriving stock (Eq,Ord,Functor,Foldable,Traversable,Generic) instance Show r => Show (EndPoint Closed r) where showsPrec = showsPrecImpl "ClosedE" @@ -117,7 +118,7 @@ pattern OpenE x = EndPoint x -- | Data type modelling an endpoint that can both be open and closed. data AnEndPoint r = AnEndPoint {-# UNPACK #-} !EndPointType !r - deriving (Show,Read,Eq,Ord,Functor,Foldable,Traversable) + deriving (Show,Read,Eq,Ord,Functor,Foldable,Traversable,Generic) type instance NumType (AnEndPoint r) = r type instance IxValue (AnEndPoint r) = r diff --git a/hgeometry/kernel/src/HGeometry/LineSegment/Internal.hs b/hgeometry/kernel/src/HGeometry/LineSegment/Internal.hs index f62092be6..c1c98a90f 100644 --- a/hgeometry/kernel/src/HGeometry/LineSegment/Internal.hs +++ b/hgeometry/kernel/src/HGeometry/LineSegment/Internal.hs @@ -27,7 +27,6 @@ import Data.Kind (Type) import HGeometry.Box.Boxable import HGeometry.Intersection import HGeometry.Interval -import HGeometry.Interval.Class import HGeometry.Line.Class import HGeometry.Line.PointAndVector import HGeometry.LineSegment.Class diff --git a/hgeometry/kernel/src/HGeometry/Triangle.hs b/hgeometry/kernel/src/HGeometry/Triangle.hs index 3f2b0ea04..4c541ab79 100644 --- a/hgeometry/kernel/src/HGeometry/Triangle.hs +++ b/hgeometry/kernel/src/HGeometry/Triangle.hs @@ -15,22 +15,22 @@ module HGeometry.Triangle ) where import Control.Lens +import GHC.Generics (Generic) import HGeometry.Box.Boxable --- import HGeometry.HalfSpace --- import HGeometry.HyperPlane -import HGeometry.Point import HGeometry.Intersection +import HGeometry.Point import HGeometry.Properties --- import HGeometry.Transformation +import HGeometry.Transformation import HGeometry.Triangle.Class import HGeometry.Vector -import Text.Read import Hiraffe.Graph +import Text.Read -------------------------------------------------------------------------------- -- | Triangles in d-dimensional space newtype Triangle point = MkTriangle (Vector 3 point) + deriving (Generic) -- | Construct a triangle from its three points pattern Triangle :: point -> point -> point -> Triangle point @@ -64,11 +64,9 @@ instance HasVertices (Triangle point) (Triangle point') where instance HasPoints (Triangle point) (Triangle point') point point' where allPoints = _TriangleVector.components --- instance ( DefaultTransformByConstraints (Triangle point) d r --- , Point_ point d r --- , d > 0 --- ) => IsTransformable (Triangle point) - +instance ( DefaultTransformByConstraints (Triangle point) d r + , Point_ point d r + ) => IsTransformable (Triangle point) instance (Show point) => Show (Triangle point) where showsPrec k (Triangle a b c ) = showParen (k > appPrec) $ diff --git a/hgeometry/point/src-quickcheck/HGeometry/Point/Instances.hs b/hgeometry/point/src-quickcheck/HGeometry/Point/Instances.hs index c958b63f0..fbdf43367 100644 --- a/hgeometry/point/src-quickcheck/HGeometry/Point/Instances.hs +++ b/hgeometry/point/src-quickcheck/HGeometry/Point/Instances.hs @@ -19,3 +19,4 @@ import Test.QuickCheck instance Arbitrary v => Arbitrary (PointF v) where arbitrary = Point <$> arbitrary + shrink (Point v) = Point <$> shrink v diff --git a/hgeometry/point/src/HGeometry/Point/Class.hs b/hgeometry/point/src/HGeometry/Point/Class.hs index c00873c01..1d34eaa56 100644 --- a/hgeometry/point/src/HGeometry/Point/Class.hs +++ b/hgeometry/point/src/HGeometry/Point/Class.hs @@ -24,6 +24,7 @@ module HGeometry.Point.Class -- , projectPoint -- , PointFor , HasPoints(..), HasPoints' + , NoDefault(..) ) where import Control.Lens @@ -31,15 +32,12 @@ import Data.Default.Class import Data.Function (on) import qualified Data.List.NonEmpty as NonEmpty import Data.Proxy (Proxy(..)) +import GHC.Generics (Generic) import GHC.TypeNats import HGeometry.Ext import HGeometry.Properties import HGeometry.Vector import qualified Linear.Affine as Linear --- import Linear.V2 (V2(..)) --- import Linear.V3 (V3(..)) --- import Linear.V4 (V4(..)) - -- $setup -- >>> import HGeometry.Point @@ -361,3 +359,12 @@ instance Affine_ point d r => Affine_ (point :+ extra) d r where instance (Point_ point d r, Default extra) => Point_ (point :+ extra) d r where {-# SPECIALIZE instance Point_ point d r => Point_ (point :+ ()) d r #-} fromVector v = fromVector v :+ def + + +-- | A newtype that can discharge the Default constraint in an unsafe way, if you really +-- sure that you'll never actually need the default +newtype NoDefault extra = NoDefault extra + deriving newtype (Show,Read,Eq,Ord,Enum,Num,Bounded,Real,Fractional,RealFrac,Generic) + +instance Default (NoDefault extra) where + def = error "NoDefault does not have an actual default. So something went wrong" diff --git a/hgeometry/src-quickcheck/HGeometry/Polygon/Instances.hs b/hgeometry/src-quickcheck/HGeometry/Polygon/Instances.hs index 3ad2b2ce6..ae0b2fe06 100644 --- a/hgeometry/src-quickcheck/HGeometry/Polygon/Instances.hs +++ b/hgeometry/src-quickcheck/HGeometry/Polygon/Instances.hs @@ -118,6 +118,11 @@ instance Arbitrary (SimplePolygon (Point 2 (RealNumber (p::Nat)))) where -- instance Arbitrary (MultiPolygon () Rational) where -- arbitrary = elements allMultiPolygons' + + + + + simplifyP :: SimplePolygon (Point 2 Rational) -> [SimplePolygon (Point 2 Rational)] simplifyP pg -- Scale up polygon such that each coordinate is a whole number. @@ -127,13 +132,12 @@ simplifyP pg -- Scale down polygon maintaining each coordinate as a whole number | gcdP /= 1 = [ pg&vertices %~ divP gcdP ] - -- unsafeFromCircularVector $ CV.map (over core (divP gcdP)) vs] | minX /= 0 || minY /= 0 = [ pg&vertices %~ align ] -- = [unsafeFromCircularVector $ CV.map (over core align) vs] | otherwise = let pg' = pg&vertices %~ _div2 -- unsafeFromCircularVector $ CV.map (over core _div2) vs - in [ pg' | not (hasSelfIntersections pg') ] + in [ pg' | hasNoSelfIntersections $ pg'^..vertices ] -- otherwise = [] where minX = first1Of (minimumVertexBy (comparing (^.xCoord)).xCoord) pg @@ -154,8 +158,6 @@ simplifyP pg divP v (Point2 c d) = Point2 (c/v) (d/v) _div2 (Point2 a b) = Point2 (numerator a `div` 2 % 1) (numerator b `div` 2 % 1) -hasSelfIntersections = const True --- FIXME! lcmPoint :: SimplePolygon (Point 2 Rational) -> Rational lcmPoint p = realToFrac t diff --git a/hgeometry/src/HGeometry/LineSegment/Intersection/BentleyOttmann.hs b/hgeometry/src/HGeometry/LineSegment/Intersection/BentleyOttmann.hs index 911e23873..d4af6f33c 100644 --- a/hgeometry/src/HGeometry/LineSegment/Intersection/BentleyOttmann.hs +++ b/hgeometry/src/HGeometry/LineSegment/Intersection/BentleyOttmann.hs @@ -33,7 +33,6 @@ module HGeometry.LineSegment.Intersection.BentleyOttmann import Control.Lens hiding (contains) import Data.Coerce -import Data.Default.Class import qualified Data.Foldable as F import Data.Function (on) import qualified Data.List as List @@ -49,9 +48,8 @@ import qualified Data.Vector as Vector import HGeometry.Foldable.Sort import HGeometry.Intersection import HGeometry.Interval.Class -import HGeometry.Interval.EndPoint +-- import HGeometry.Interval.EndPoint import HGeometry.LineSegment -import HGeometry.LineSegment.Class import HGeometry.LineSegment.Intersection.Types import HGeometry.Point import HGeometry.Properties (NumType, Dimension) diff --git a/hgeometry/src/HGeometry/PlaneGraph/Class.hs b/hgeometry/src/HGeometry/PlaneGraph/Class.hs index 02ac8edba..313008613 100644 --- a/hgeometry/src/HGeometry/PlaneGraph/Class.hs +++ b/hgeometry/src/HGeometry/PlaneGraph/Class.hs @@ -22,15 +22,19 @@ module HGeometry.PlaneGraph.Class , dartSegments , edgeSegments + , interiorFacePolygonAt , interiorFacePolygons ) where import Control.Lens +import Data.Coerce +import Data.Default.Class +import Data.Foldable1 import Data.Functor.Apply -import Data.Kind (Constraint) import Data.Maybe (fromMaybe) import Data.Monoid (Endo(..)) import Data.Ord (comparing) +import HGeometry.Ext import HGeometry.LineSegment import HGeometry.Point import HGeometry.Polygon.Simple @@ -49,7 +53,21 @@ class ( PlanarGraph_ planeGraph , NumType vertex ~ NumType planeGraph ) => PlaneGraph_ planeGraph vertex | planeGraph -> vertex where - {-# MINIMAL #-} + {-# MINIMAL fromEmbedding #-} + + -- | Build a graph from its embedding; i.e. for each vertex we expect its adjacencies in + -- CCW order. + -- + -- If the, in the list of neighbours of vertex u we see a vertex v + -- that itself does not appear in the adjacencylist, we may drop + -- it. In other words if u has a neighbour v, then v better have a + -- specification of its neighbours somewhere. + fromEmbedding :: ( Foldable1 f, Functor f, Foldable h, Functor h + , vi ~ VertexIx planeGraph + , v ~ Vertex planeGraph + , e ~ Edge planeGraph + , GraphFromAdjListExtraConstraints planeGraph h + ) => f (vi, v, h (vi, e)) -> planeGraph -- | Getter to access the outer face outerFace :: Eq (FaceIx planeGraph) @@ -191,23 +209,67 @@ interiorFacePolygons :: forall planeGraph vertex r. , Ord r, Num r , Eq (FaceIx planeGraph) ) - => IndexedFold (FaceIx planeGraph) planeGraph (SimplePolygon vertex) + => IndexedFold (FaceIx planeGraph) + planeGraph + (SimplePolygon (vertex :+ VertexIx planeGraph)) interiorFacePolygons = theFold where theFold :: forall p f. - ( Indexable (FaceIx planeGraph) p, Applicative f, Contravariant f) - => p (SimplePolygon vertex) (f (SimplePolygon vertex)) + ( PlaneGraph_ planeGraph vertex + , Indexable (FaceIx planeGraph) p, Applicative f, Contravariant f) + => p (SimplePolygon (vertex :+ VertexIx planeGraph)) + (f (SimplePolygon (vertex :+ VertexIx planeGraph))) -> planeGraph -> f planeGraph theFold pPolyFPoly g = interiorFaces (Indexed draw) g where draw :: FaceIx planeGraph -> Face planeGraph -> f (Face planeGraph) - draw fi _ = let poly = uncheckedFromCCWPoints . fmap (\vi -> g^?!vertexAt vi) - $ boundaryVertices fi g + draw fi _ = let poly = polygonFromFace g fi in poly >$ indexed pPolyFPoly fi poly + +polygonFromFace :: forall planeGraph vertex r.( PlaneGraph_ planeGraph vertex + , Point_ vertex 2 r + ) + => planeGraph -> FaceIx planeGraph + -> SimplePolygon (vertex :+ VertexIx planeGraph) +polygonFromFace gr fi = poly'&vertices.extra %~ coerce + where + poly' :: SimplePolygon (vertex :+ NoDefault (VertexIx planeGraph)) + poly' = uncheckedFromCCWPoints + . fmap (\vi -> gr^?!vertexAt vi :+ NoDefault vi) + $ boundaryVertices fi gr -- note that this is safe, since boundaryVerticesOf guarantees that for -- interior faces, the vertices are returned in CCW order. + -- TODO: why can't I just coerce poyl' to the right type? + + +-- | Renders a single interior face as a simple polygon. +-- +-- Note that this is a fold rather than a getter for the same reason faceAt is a traversal +-- rather than a lens: i.e. if you pass some nonsensical FaceIx the face may not exist. +interiorFacePolygonAt :: forall planeGraph vertex. + ( PlaneGraph_ planeGraph vertex + , Point_ vertex 2 (NumType vertex) + ) + => FaceIx planeGraph + -> IndexedFold (FaceIx planeGraph) + planeGraph + (SimplePolygon (vertex :+ VertexIx planeGraph)) +interiorFacePolygonAt fi = theFold + where + theFold pPolyFPoly gr = faceAt fi draw gr + where + -- draw :: Face planeGraph -> f (Face planeGraph) + draw _ = let poly = polygonFromFace gr fi + in poly >$ indexed pPolyFPoly fi poly + + +newtype NoDefault e = NoDefault e + +instance Default (NoDefault e) where + def = undefined + -------------------------------------------------------------------------------- diff --git a/hgeometry/src/HGeometry/PlaneGraph/Type.hs b/hgeometry/src/HGeometry/PlaneGraph/Type.hs index 52d8942e7..33e861c69 100644 --- a/hgeometry/src/HGeometry/PlaneGraph/Type.hs +++ b/hgeometry/src/HGeometry/PlaneGraph/Type.hs @@ -19,14 +19,19 @@ module HGeometry.PlaneGraph.Type import Control.Lens hiding (holes, holesOf, (.=)) import Data.Coerce import Data.Foldable1 +import qualified Data.Map as Map +import qualified Data.Vector.NonEmpty as Vector import Data.YAML import GHC.Generics (Generic) import HGeometry.Box +import HGeometry.Foldable.Sort (sortBy ) import HGeometry.PlaneGraph.Class import HGeometry.Point import HGeometry.Properties +import HGeometry.Transformation import Hiraffe.PlanarGraph import qualified Hiraffe.PlanarGraph as PG + -------------------------------------------------------------------------------- -- * The PlaneGraph type @@ -89,7 +94,11 @@ instance HasFaces (PlaneGraph s v e f) (PlaneGraph s v e f') where ---------------------------------------- instance DiGraph_ (PlaneGraph s v e f) where type DiGraphFromAdjListExtraConstraints (PlaneGraph s v e f) h = (f ~ (), Foldable1 h) + + -- | The vertices are expected to have their adjacencies in CCW order. diGraphFromAdjacencyLists = PlaneGraph . diGraphFromAdjacencyLists + -- TODO: we should probably use some toEmbedding here as well I think + endPoints (PlaneGraph g) = endPoints g twinDartOf d = twinOf d . to Just outgoingDartsOf v = _PlanarGraph.outgoingDartsOf v @@ -99,14 +108,37 @@ instance BidirGraph_ (PlaneGraph s v e f) where getPositiveDart (PlaneGraph g) e = getPositiveDart g e -instance Graph_ (PlaneGraph s v e f) where +-- | Computes the cyclic order of adjacencies around each vertex. +-- +-- \(O(n \log n)\) +toEmbedding :: ( Foldable1 g, Functor g, Foldable h, Functor h + , vi ~ VertexIx (PlaneGraph s v e f) + , v ~ Vertex (PlaneGraph s v e f) + , e ~ Edge (PlaneGraph s v e f) + , GraphFromAdjListExtraConstraints (PlaneGraph s v e f) h + , Point_ v 2 r, Ord r, Num r + ) => g (vi, v, h (vi, e)) -> g (vi, v, Vector.NonEmptyVector (vi, e)) +toEmbedding vs = fmap sortAround' vs + where + vertexLocs = foldMap (\(vi,v,_) -> Map.singleton vi v) vs + sortAround' (vi,v,adjs) = (vi,v, Vector.unsafeFromVector $ sortBy (ccwCmpAround' v) adjs) + ccwCmpAround' v (ui,_) (wi,_) = ccwCmpAround v (vertexLocs Map.! ui) (vertexLocs Map.! wi) + + + +instance ( Point_ v 2 (NumType v) + , Ord (NumType v), Num (NumType v) + ) => Graph_ (PlaneGraph s v e f) where type GraphFromAdjListExtraConstraints (PlaneGraph s v e f) h = (f ~ (), Foldable1 h) - fromAdjacencyLists = PlaneGraph . fromAdjacencyLists + + fromAdjacencyLists = fromEmbedding . toEmbedding neighboursOf u = _PlanarGraph.neighboursOf u incidentEdgesOf u = _PlanarGraph.incidentEdgesOf u -instance PlanarGraph_ (PlaneGraph s v e f) where +instance ( Point_ v 2 (NumType v) + , Ord (NumType v), Num (NumType v) + ) => PlanarGraph_ (PlaneGraph s v e f) where type DualGraphOf (PlaneGraph s v e f) = PlanarGraph s Dual f e v dualGraph = dualGraph . coerce @_ @(PlanarGraph s Primal v e f) @@ -123,7 +155,8 @@ instance PlanarGraph_ (PlaneGraph s v e f) where instance ( Point_ v 2 (NumType v) , Ord (NumType v), Num (NumType v) - ) => PlaneGraph_ (PlaneGraph s v e f) v + ) => PlaneGraph_ (PlaneGraph s v e f) v where + fromEmbedding = PlaneGraph . fromAdjacencyLists instance ( Point_ v 2 r, Point_ v' 2 r' ) => HasPoints (PlaneGraph s v e f) @@ -134,5 +167,10 @@ instance ( Point_ v 2 r , Ord r, Num r ) => IsBoxable (PlaneGraph s v e f) +instance ( Point_ v 2 r + , DefaultTransformByConstraints (PlaneGraph s v e f) 2 r + ) => IsTransformable (PlaneGraph s v e f) + + -- boundingBox = boundingBoxList' . F.toList . fmap (^._2.location) . vertices diff --git a/hgeometry/src/HGeometry/Polygon/Class.hs b/hgeometry/src/HGeometry/Polygon/Class.hs index e73158769..e5aa6897b 100644 --- a/hgeometry/src/HGeometry/Polygon/Class.hs +++ b/hgeometry/src/HGeometry/Polygon/Class.hs @@ -74,13 +74,15 @@ class HasVertices polygon polygon => HasOuterBoundary polygon where -- -- -- running time :: \(O(n)\) - outerBoundaryEdges :: IndexedFold1 (VertexIx polygon) polygon (Vertex polygon, Vertex polygon) + outerBoundaryEdges :: IndexedFold1 (VertexIx polygon,VertexIx polygon) polygon + (Vertex polygon, Vertex polygon) default outerBoundaryEdges :: Enum (VertexIx polygon) - => IndexedFold1 (VertexIx polygon) polygon (Vertex polygon, Vertex polygon) + => IndexedFold1 (VertexIx polygon,VertexIx polygon) polygon (Vertex polygon, Vertex polygon) outerBoundaryEdges = ifolding1 $ - \pg -> ( \(i,u) -> (i,(u, pg ^.outerBoundaryVertexAt (succ i))) ) - <$> itoNonEmptyOf outerBoundary pg + \pg -> ( \(i,u) -> let (j,v) = pg ^.outerBoundaryVertexAt (succ i).withIndex + in ((i,j) , (u,v)) + ) <$> itoNonEmptyOf outerBoundary pg -- \pg -> fmap ( \(i,u) -> (i,(u, pg ^.outerBoundaryVertexAt (succ i))) ) -- . NonEmpty.fromList -- $ itoListOf outerBoundary pg @@ -91,16 +93,17 @@ class HasVertices polygon polygon => HasOuterBoundary polygon where -- -- running time: \(O(1)\) outerBoundaryEdgeAt :: VertexIx polygon - -> IndexedGetter (VertexIx polygon) polygon + -> IndexedGetter (VertexIx polygon, VertexIx polygon) polygon (Vertex polygon, Vertex polygon) -- default implementation of outerBoundaryEdge. It achieves the -- desired running time when indexing is indeed constant. default outerBoundaryEdgeAt :: Enum (VertexIx polygon) => VertexIx polygon - -> IndexedGetter (VertexIx polygon) polygon + -> IndexedGetter (VertexIx polygon, VertexIx polygon) polygon (Vertex polygon, Vertex polygon) outerBoundaryEdgeAt i = ito $ - \pg -> (i, (pg^.outerBoundaryVertexAt i, pg^.outerBoundaryVertexAt (succ i))) + \pg -> let (j,v) = pg^.outerBoundaryVertexAt (succ i).withIndex + in ( (i,j), (pg^.outerBoundaryVertexAt i, v)) -------------------------------------------------------------------------------- @@ -170,8 +173,9 @@ outerBoundaryEdgeSegmentAt :: ( HasOuterBoundary polygon , Point_ point 2 r ) => VertexIx polygon - -> IndexedGetter (VertexIx polygon) polygon - (ClosedLineSegment point) + -> IndexedGetter (VertexIx polygon, VertexIx polygon) + polygon + (ClosedLineSegment point) outerBoundaryEdgeSegmentAt i = outerBoundaryEdgeAt i. to (uncurry ClosedLineSegment) -- | Get the line segments representing the outer boundary of the polygon. @@ -180,25 +184,27 @@ outerBoundaryEdgeSegments :: forall polygon point r. , Vertex polygon ~ point , Point_ point 2 r ) - => IndexedFold1 (VertexIx polygon) polygon (ClosedLineSegment point) + => IndexedFold1 (VertexIx polygon,VertexIx polygon) + polygon + (ClosedLineSegment point) outerBoundaryEdgeSegments = outerBoundaryEdges . to (uncurry ClosedLineSegment) -- | A fold that associates each vertex on the boundary with its -- predecessor and successor (in that order) along the boundary. outerBoundaryWithNeighbours :: ( HasOuterBoundary polygon - , Enum (VertexIx polygon) + , VertexIx polygon ~ Int ) - => IndexedFold1 (VertexIx polygon) + => IndexedFold1 (VertexIx polygon, + (VertexIx polygon, VertexIx polygon)) polygon (Vertex polygon, (Vertex polygon, Vertex polygon)) outerBoundaryWithNeighbours = ifolding1 $ - \pg -> (\(i,u) -> (i, f pg i u)) <$> itoNonEmptyOf outerBoundary pg + \pg -> f pg (numVertices pg) <$> itoNonEmptyOf outerBoundary pg where - f pg i u = let v = pg^.outerBoundaryVertexAt (pred i) - w = pg^.outerBoundaryVertexAt (succ i) - in (u, (v, w)) - + f pg n (i,u) = let (j,v) = pg^.outerBoundaryVertexAt ((i-1) `mod` n).withIndex + (k,w) = pg^.outerBoundaryVertexAt ((i+1) `mod` n).withIndex + in ( (i, (j,k)) , (u, (v, w)) ) -------------------------------------------------------------------------------- diff --git a/hgeometry/src/HGeometry/Polygon/Simple.hs b/hgeometry/src/HGeometry/Polygon/Simple.hs index 5d1c943db..e1e66c69a 100644 --- a/hgeometry/src/HGeometry/Polygon/Simple.hs +++ b/hgeometry/src/HGeometry/Polygon/Simple.hs @@ -14,11 +14,12 @@ module HGeometry.Polygon.Simple ( SimplePolygon_(..) , SimplePolygon - , SimplePolygonF + , SimplePolygonF(..) , toCyclic , VertexContainer , HasInPolygon(..) , inSimplePolygon + , hasNoSelfIntersections ) where import Control.DeepSeq (NFData) @@ -26,6 +27,7 @@ import Control.Lens import qualified Data.Foldable as F import Data.Functor.Classes import qualified Data.List.NonEmpty as NonEmpty +import qualified Data.Map as Map import Data.Semigroup.Foldable import Data.Vector.NonEmpty.Internal (NonEmptyVector(..)) import GHC.Generics @@ -34,6 +36,7 @@ import HGeometry.Box import HGeometry.Cyclic import HGeometry.Foldable.Util import HGeometry.Intersection +import HGeometry.LineSegment.Intersection.BentleyOttmann import HGeometry.Point import HGeometry.Polygon.Class import HGeometry.Polygon.Simple.Class @@ -123,8 +126,15 @@ instance ( Point_ point 2 r , VertexContainer f point ) => Polygon_ (SimplePolygonF f point) point r where area = areaSimplePolygon - ccwPredecessorOf u = singular $ vertexAt (pred u) - ccwSuccessorOf u = singular $ vertexAt (succ u) + ccwPredecessorOf u = \pvFv pg -> let n = numVertices pg + p = (pred u) `mod` n + l = singular $ vertexAt p + in l pvFv pg + -- make sure to wrap the index to make sure we report the right index. + ccwSuccessorOf u = \pvFv pg -> let n = numVertices pg + s = (succ u) `mod` n + l = singular $ vertexAt s + in l pvFv pg instance ( Point_ point 2 r , VertexContainer f point @@ -197,3 +207,16 @@ instance ( SimplePolygon_ (SimplePolygonF f point) point r q `intersect` pg | q `intersects` pg = Just q | otherwise = Nothing -- this implementation is a bit silly but ok + + +-------------------------------------------------------------------------------- + +-- | verify that some sequence of points has no self intersecting edges. +hasNoSelfIntersections :: forall f point r. + (Foldable f, Functor f, Point_ point 2 r, Ord r, Real r) + => f point -> Bool +hasNoSelfIntersections vs = let vs' = (\p -> (p^.asPoint)&coordinates %~ toRational) <$> vs + pg :: SimplePolygon (Point 2 Rational) + pg = uncheckedFromCCWPoints vs' + in Map.null $ interiorIntersections $ pg^..outerBoundaryEdgeSegments + -- outerBoundaryEdgeSegments interiorIntersections pg diff --git a/hgeometry/src/HGeometry/Polygon/Simple/InPolygon.hs b/hgeometry/src/HGeometry/Polygon/Simple/InPolygon.hs index ca73aebfd..fbf300610 100644 --- a/hgeometry/src/HGeometry/Polygon/Simple/InPolygon.hs +++ b/hgeometry/src/HGeometry/Polygon/Simple/InPolygon.hs @@ -116,11 +116,11 @@ q `inSimplePolygon` pg = case ifoldMapOf outerBoundaryEdges countAbove pg of -- (e.g. the <=> does not hold in that case). If we discover q lies -- on an edge e, we actually report that. where - countAbove i (u,v) = case (u^.xCoord) `compare` (v^.xCoord) of - EQ | onVerticalEdge u v -> OnEdge i - | otherwise -> mempty - LT -> countAbove' i (AnClosedE u) (AnOpenE v) - GT -> countAbove' i (AnOpenE v) (AnClosedE u) + countAbove (i,_) (u,v) = case (u^.xCoord) `compare` (v^.xCoord) of + EQ | onVerticalEdge u v -> OnEdge i + | otherwise -> mempty + LT -> countAbove' i (AnClosedE u) (AnOpenE v) + GT -> countAbove' i (AnOpenE v) (AnClosedE u) -- count the edge if q countAbove' i l r | (Point1 $ q^.xCoord) `intersects` (view xCoord <$> Interval l r) diff --git a/hgeometry/src/HGeometry/Polygon/Triangulation.hs b/hgeometry/src/HGeometry/Polygon/Triangulation.hs new file mode 100644 index 000000000..c606ad733 --- /dev/null +++ b/hgeometry/src/HGeometry/Polygon/Triangulation.hs @@ -0,0 +1,77 @@ +{-# LANGUAGE PartialTypeSignatures #-} +-------------------------------------------------------------------------------- +-- | +-- Module : HGeometry.Polygon.Triangulation +-- Copyright : (C) Frank Staals +-- License : see the LICENSE file +-- Maintainer : Frank Staals +-------------------------------------------------------------------------------- +module HGeometry.Polygon.Triangulation + ( triangulate + , computeDiagonals + + , PolygonEdgeType(..) + , PolygonFaceData(..) + , Diagonal + ) where + +import Control.Lens +import Data.Coerce +import HGeometry.Ext +import HGeometry.PlaneGraph +import HGeometry.Polygon.Simple +import qualified HGeometry.Polygon.Triangulation.MakeMonotone as MM +import qualified HGeometry.Polygon.Triangulation.TriangulateMonotone as TM +import HGeometry.Polygon.Triangulation.Types +import HGeometry.Vector +import Hiraffe.PlanarGraph + +-------------------------------------------------------------------------------- + +-- | Triangulates a polygon of \(n\) vertices +-- +-- running time: \(O(n \log n)\) +triangulate :: forall s polygon point r. + (SimplePolygon_ polygon point r, Ord r, Num r) + => polygon + -> PlaneGraph s point PolygonEdgeType PolygonFaceData +triangulate pg = constructGraph pg (computeDiagonals pg) + +-- | Computes a set of diagaonals that together triangulate the input polygon +-- of \(n\) vertices. +-- +-- running time: \(O(n \log n)\) +computeDiagonals :: forall polygon point r. + (SimplePolygon_ polygon point r, Ord r, Num r) + => polygon -> [Diagonal polygon] +computeDiagonals pg = monotoneDiags <> extraDiags + where + monotoneSubdiv :: PlaneGraph () point PolygonEdgeType PolygonFaceData + monotoneSubdiv = MM.makeMonotone @() pg + -- use some arbitrary proxy type + + -- get the existing diagonals + monotoneDiags :: [Diagonal polygon] + monotoneDiags = map (\(d,_) -> let (u,v) = monotoneSubdiv ^. endPointsOf d.asIndex + in coerce $ Vector2 u v + ) + . filter ((== Diagonal) . snd) + $ monotoneSubdiv ^.. edges.withIndex + + -- and compute the diagonals in each interior y-monotone polygon + extraDiags :: [Diagonal polygon] + extraDiags = foldMap (collectDiags . fst) + $ filter ((== Inside) . snd) + $ monotoneSubdiv^..interiorFaces.withIndex + + -- collectDiags :: FaceIx planeGraph -> [Diagonal polygon] + collectDiags i = let yMonotonePoly = monotoneSubdiv ^?! interiorFacePolygonAt i + in map (coerce . withOriginalId yMonotonePoly) $ + TM.computeDiagonals yMonotonePoly + +withOriginalId :: ( TM.YMonotonePolygon_ yMonotonePolygon (point :+ i) r + ) => yMonotonePolygon -> Diagonal yMonotonePolygon -> Vector 2 i +withOriginalId yMonotonePoly = fmap (\j -> yMonotonePoly^?!vertexAt j.extra) + + + -- getOriginalID :: VertexId monotonePolygon -> VertexId polygon diff --git a/hgeometry/src/HGeometry/Polygon/Triangulation/MakeMonotone.hs b/hgeometry/src/HGeometry/Polygon/Triangulation/MakeMonotone.hs index 820cf7eb4..ed0133bc9 100644 --- a/hgeometry/src/HGeometry/Polygon/Triangulation/MakeMonotone.hs +++ b/hgeometry/src/HGeometry/Polygon/Triangulation/MakeMonotone.hs @@ -6,7 +6,8 @@ -- Maintainer : Frank Staals -------------------------------------------------------------------------------- module HGeometry.Polygon.Triangulation.MakeMonotone - ( computeDiagonals + ( makeMonotone + , computeDiagonals , classifyVertices , VertexType(..) @@ -22,8 +23,10 @@ import qualified Data.Vector as Vector import HGeometry.Ext import HGeometry.Foldable.Sort import HGeometry.LineSegment +import HGeometry.PlaneGraph import HGeometry.Point import HGeometry.Polygon.Class +import HGeometry.Polygon.Simple.Class import HGeometry.Polygon.Triangulation.Types import qualified HGeometry.Set.Util as SS import HGeometry.Vector @@ -32,6 +35,17 @@ import qualified VectorBuilder.Vector as Builder ---------------------------------------------------------------------------------- +-- | Given a polygon, computes a plane subdivision representing a split of the polygon +-- into y-monotone subpolygons. +-- +-- running time: \(O(n\log n)\) +makeMonotone :: forall s polygon point r. + (SimplePolygon_ polygon point r, Ord r, Num r) + => polygon + -> PlaneGraph s point PolygonEdgeType PolygonFaceData +makeMonotone pg = constructGraph pg (computeDiagonals pg) + + -- | Given a polygon, find a set of non-intersecting diagonals that partition -- the polygon into y-monotone pieces. -- @@ -77,7 +91,8 @@ p `cmpSweep` q = comparing (^.yCoord) p q <> comparing (Down . (^.xCoord)) p q -- | Handle an event handle :: forall polygon point r. ( Polygon_ polygon point r , Point_ point 2 r, Num r, Ord r - , Default (VertexIx polygon), Ord (VertexIx polygon) + , Default (VertexIx polygon) + , Ord (VertexIx polygon) ) => polygon -> (StatusSruct polygon, [Diagonal polygon]) -> Event polygon diff --git a/hgeometry/src/HGeometry/Polygon/Triangulation/Triangulate.hs b/hgeometry/src/HGeometry/Polygon/Triangulation/Triangulate.hs deleted file mode 100644 index 553734f88..000000000 --- a/hgeometry/src/HGeometry/Polygon/Triangulation/Triangulate.hs +++ /dev/null @@ -1,75 +0,0 @@ --------------------------------------------------------------------------------- --- | --- Module : HGeometry.Polygon.Triangulation.MakeMonotone --- Copyright : (C) Frank Staals --- License : see the LICENSE file --- Maintainer : Frank Staals --------------------------------------------------------------------------------- -module HGeometry.Polygon.Triangulation.Triangulate - ( - ) where - -import Control.Lens -import Data.Either (lefts) -import qualified Data.Foldable as F -import HGeometry.Ext -import HGeometry.LineSegment -import HGeometry.PlanarSubdivision.Basic -import HGeometry.Polygon -import qualified HGeometry.Polygon.Triangulation.MakeMonotone as MM -import qualified HGeometry.Polygon.Triangulation.TriangulateMonotone.TriangulateMonotone as TM -import HGeometry.Polygon.Triangulation.Types - --------------------------------------------------------------------------------- - --- | Triangulates a polygon of \(n\) vertices --- --- running time: \(O(n \log n)\) -triangulate :: forall s t p r. (Ord r, Fractional r) - => Polygon t p r -> PlanarSubdivision s p PolygonEdgeType PolygonFaceData r -triangulate pg' = constructSubdivision e es diags - where - (pg, diags) = computeDiagonals' pg' - (e:es) = listEdges pg - - --- | Triangulates a polygon of \(n\) vertices --- --- running time: \(O(n \log n)\) -triangulate' :: forall s t p r. (Ord r, Fractional r) - => Polygon t p r -> PlaneGraph s p PolygonEdgeType PolygonFaceData r -triangulate' pg' = constructGraph e es diags - where - (pg, diags) = computeDiagonals' pg' - (e:es) = listEdges pg - - --- | Computes a set of diagaonals that together triangulate the input polygon --- of \(n\) vertices. --- --- running time: \(O(n \log n)\) -computeDiagonals :: (Ord r, Fractional r) => Polygon t p r -> [LineSegment 2 p r] -computeDiagonals = snd . computeDiagonals' - --- | Computes a set of diagaonals that together triangulate the input polygon --- of \(n\) vertices. Returns a copy of the input polygon, whose boundaries are --- oriented in counter clockwise order, as well. --- --- running time: \(O(n \log n)\) -computeDiagonals' :: (Ord r, Fractional r) - => Polygon t p r -> (Polygon t p r, [LineSegment 2 p r]) -computeDiagonals' pg' = (pg, monotoneDiags <> extraDiags) - where - pg = toCounterClockWiseOrder pg' - monotoneP = MM.makeMonotone @() pg -- use some arbitrary proxy type - -- outerFaceId' = outerFaceId monotoneP - - monotoneDiags = map (^._2.core) . filter (\e' -> e'^._2.extra == Diagonal) - . F.toList . edgeSegments $ monotoneP - extraDiags = concatMap (TM.computeDiagonals . toCounterClockWiseOrder') - . lefts . map (^._2.core) - . filter (\mp -> mp^._2.extra == Inside) -- triangulate only the insides - -- . filter (\f -> f^._1 /= outerFaceId') - . F.toList . internalFacePolygons $ monotoneP - -- FIXME: we should already get all polygons in CCW order, so no - -- need for the toClockwiseOrder' call diff --git a/hgeometry/src/HGeometry/Polygon/Triangulation/TriangulateMonotone.hs b/hgeometry/src/HGeometry/Polygon/Triangulation/TriangulateMonotone.hs index 23c2cae0d..424472e3d 100644 --- a/hgeometry/src/HGeometry/Polygon/Triangulation/TriangulateMonotone.hs +++ b/hgeometry/src/HGeometry/Polygon/Triangulation/TriangulateMonotone.hs @@ -7,8 +7,7 @@ -------------------------------------------------------------------------------- module HGeometry.Polygon.Triangulation.TriangulateMonotone ( YMonotonePolygon_ - -- , triangulate - -- , triangulate' + , triangulate , computeDiagonals -- , LR(..) -- , P @@ -32,7 +31,7 @@ import HGeometry.Combinatorial.Util import HGeometry.Ext import HGeometry.Vector (Vector(Vector2)) -- import HGeometry.PlanarSubdivision.Basic (PlanarSubdivision, PolygonFaceData) --- import HGeometry.PlaneGraph (PlaneGraph) +import HGeometry.PlaneGraph (PlaneGraph) import HGeometry.Point import HGeometry.Polygon.Class import HGeometry.Polygon.Simple.Class @@ -48,35 +47,18 @@ import HGeometry.Polygon.Triangulation.Types type YMonotonePolygon_ = SimplePolygon_ data LR = L | R deriving (Show,Eq) -{- - -type PlanarSubdivision s point e f r = () -- | Triangulates a polygon of \(n\) vertices -- -- running time: \(O(n \log n)\) -triangulate :: ( YMonotonePolygon_ yMonotonePolygon point r - , Ord r, Fractional r) +triangulate :: forall s yMonotonePolygon point r. + (YMonotonePolygon_ yMonotonePolygon point r, Ord r, Num r + ) => yMonotonePolygon - -> PlanarSubdivision s point PolygonEdgeType PolygonFaceData r -triangulate pg = constructSubdivision e es (computeDiagonals pg) - where - (e:es) = listEdges pg + -> PlaneGraph s point PolygonEdgeType PolygonFaceData +triangulate pg = constructGraph pg (computeDiagonals pg) -- TODO: Find a way to construct the graph in O(n) time. --- | Triangulates a polygon of \(n\) vertices --- --- running time: \(O(n \log n)\) -triangulate' :: forall s p r. (Ord r, Fractional r) - => yMonotonePolygon-> PlaneGraph s p PolygonEdgeType PolygonFaceData r -triangulate' pg' = constructGraph e es (computeDiagonals pg) - where - pg = toCounterClockWiseOrder pg' - (e:es) = listEdges pg - -- TODO: Find a way to construct the graph in O(n) time. - --} - -- | Given a y-monotone polygon in counter clockwise order computes the diagonals -- to add to triangulate the polygon -- diff --git a/hgeometry/src/HGeometry/Polygon/Triangulation/Types.hs b/hgeometry/src/HGeometry/Polygon/Triangulation/Types.hs index 41c726c58..47f5ecffa 100644 --- a/hgeometry/src/HGeometry/Polygon/Triangulation/Types.hs +++ b/hgeometry/src/HGeometry/Polygon/Triangulation/Types.hs @@ -8,23 +8,21 @@ module HGeometry.Polygon.Triangulation.Types where import Control.Lens -import Control.Monad (forM_) -import qualified Data.Foldable as F +import Data.Coerce import Data.List.NonEmpty (NonEmpty(..)) -import qualified Data.List.NonEmpty as NonEmpty -import qualified Data.Vector as V --- import qualified Data.Vector.Mutable as MV -import HGeometry.Ext -import HGeometry.LineSegment --- import HGeometry.PlanarSubdivision.Basic --- import qualified HGeometry.PlaneGraph as PG -import Hiraffe.Graph.Class +import qualified Data.Map as Map +import Data.Maybe (fromMaybe) +import HGeometry.Lens.Util import HGeometry.PlaneGraph -import HGeometry.Polygon.Class import HGeometry.Point +import HGeometry.Polygon.Class +import HGeometry.Polygon.Simple.Class import HGeometry.Vector +import Hiraffe.Graph.Class import Hiraffe.PlanarGraph as PlanarGraph +-- import qualified Data.Foldable as F +-- import Debug.Trace -------------------------------------------------------------------------------- -- | After triangulation, edges are either from the original polygon or a new diagonal. @@ -90,48 +88,84 @@ constructSubdivision e origs diags = fromPlaneGraph $ constructGraph e origs dia -} +-- fromConnectedSegments :: ( Foldable1 f +-- , LineSegment_ lineSegment vertex +-- , Point_ vertex 2 r, Ord r, Num r +-- ) +-- => f lineSegment -> PlaneGraph s (NonEmpty vertex) lineSegment () +-- fromConnectedSegments segs = fromAdjacencyLists adjLists +-- where +-- adjLists = Map.fromListWith (<>) . concatMap f . zipWith g [0..] . F.toList $ ss + +-- mkVertex + + + +-- fromadjacencylists +-- where +-- planarGraph dts & PG.vertexData .~ vxData +-- where +-- pts = M.fromListWith (<>) . concatMap f . zipWith g [0..] . F.toList $ ss +-- f (s :+ e) = [ ( s^.start.core +-- , SP (sing $ s^.start.extra) [(s^.end.core) :+ h Positive e]) +-- , ( s^.end.core +-- , SP (sing $ s^.end.extra) [(s^.start.core) :+ h Negative e]) +-- ] +-- g i (s :+ e) = s :+ (Arc i :+ e) +-- h d (a :+ e) = (Dart a d, e) + +-- sing x = x NonEmpty.:| [] + +-- vts = map (\(p,sp) -> (p,map (^.extra) . sortAround' (ext p) <$> sp)) +-- . M.assocs $ pts +-- -- vertex Data +-- vxData = V.fromList . map (\(p,sp) -> VertexData p (sp^._1)) $ vts +-- -- The darts +-- dts = map (^._2._2) vts + + -- | Given a list of original edges and a list of diagonals, creates a --- planar-subdivision +-- planar-subdivision. The vertexId's will remain unchanged. -- -- -- running time: \(O(n\log n)\) -constructGraph :: (Polygon_ polygon point r, Point_ point 2 r +constructGraph :: forall s polygon point r f. + ( SimplePolygon_ polygon point r, Point_ point 2 r + , Foldable f, Ord r, Num r ) => polygon -> f (Diagonal polygon) -> PlaneGraph s point PolygonEdgeType PolygonFaceData -constructGraph pg diags = undefined -{- -e origs diags = - subdiv & PG.vertexData.traverse %~ NonEmpty.head - & PG.faceData .~ faceData' - & PG.rawDartData.traverse %~ snd +constructGraph pg diags = gr&faces %@~ computeFaceLabel +-- constructGraph pg diags = gr&faces %@~ computeFaceLabel where - subdiv :: PG.PlaneGraph s (NonEmpty p) (Bool,PolygonEdgeType) () r - subdiv = PG.fromConnectedSegments $ e' : origs' <> diags' - - diags' = (:+ (True, Diagonal)) <$> diags - origs' = (:+ (False,Original)) <$> origs - e' = e :+ (True, Original) - - -- the darts incident to internal faces - queryDarts = concatMap shouldQuery . F.toList . PG.edges' $ subdiv - shouldQuery d = case subdiv^.dataOf d of - (True, Original) -> [d] - (True, Diagonal) -> [d, twin d] - _ -> [] - - -- the interior faces - intFaces = flip PG.leftFace subdiv <$> queryDarts - - faceData' :: V.Vector PolygonFaceData - faceData' = V.create $ do - v' <- MV.replicate (PG.numFaces subdiv) Outside - forM_ intFaces $ \(PG.FaceId (PG.VertexId f)) -> - MV.write v' f Inside - pure v' --} + -- | Note that we use fromAdjacencyLists + gr = fromAdjacencyLists adjLists :: PlaneGraph s point PolygonEdgeType () + + adjLists = uncurry collectDiags <$> itoNonEmptyOf outerBoundaryWithNeighbours pg + + collectDiags :: (VertexIx polygon, (VertexIx polygon, VertexIx polygon)) + -> (point, (point, point)) + -> ( VertexIdIn Primal s + , point, NonEmpty (VertexIdIn Primal s,PolygonEdgeType) + ) + collectDiags (u, (v,w)) (p,_) = coerce (u,p, (v, Original) :| (w, Original) : diagonalsOf u) + + -- get the diagonals incident to vertex u + diagonalsOf u = fromMaybe [] $ Map.lookup u diags' + -- associate every diagonal with its endpoints + diags' :: Map.Map (VertexIx polygon) [(VertexIx polygon, PolygonEdgeType)] + diags' = foldr (\(Vector2 u v) -> Map.insertWith (<>) u [(v, Diagonal)] + . Map.insertWith (<>) v [(u, Diagonal)] + ) Map.empty diags + + theOuterFaceId = outerFaceId gr + computeFaceLabel fi _ + | fi == theOuterFaceId = Outside + | otherwise = Inside + + -- -constructSubdivision px e origs diags = -- - subdiv & planeGraph.PG.vertexData.traverse %~ NonEmpty.head @@ -160,3 +194,6 @@ e origs diags = -- - forM_ intFaces $ \(PG.FaceId (PG.VertexId f)) -> -- - MV.write v' f (FaceData [] Inside) -- - pure v' + + +-------------------------------------------------------------------------------- diff --git a/hgeometry/test-with-ipe/test/PlaneGraph/RenderSpec.hs b/hgeometry/test-with-ipe/test/PlaneGraph/RenderSpec.hs index 96268f6df..caa710ff2 100644 --- a/hgeometry/test-with-ipe/test/PlaneGraph/RenderSpec.hs +++ b/hgeometry/test-with-ipe/test/PlaneGraph/RenderSpec.hs @@ -2,6 +2,11 @@ {-# LANGUAGE OverloadedStrings #-} module PlaneGraph.RenderSpec ( spec + , drawGraph + , drawVertex + , drawDart + , drawFace + , drawEdge ) where import Control.Lens @@ -14,7 +19,6 @@ import HGeometry.Boundary import HGeometry.Ext import HGeometry.Line import HGeometry.LineSegment -import HGeometry.Number.Radical import HGeometry.PlaneGraph import HGeometry.PlaneGraph.Class import HGeometry.Point @@ -94,13 +98,16 @@ spec = describe "render planegraph tests" $ do drawGraph :: ( PlaneGraph_ planeGraph vertex , IsTransformable vertex - , Point_ vertex 2 r, Ord r, Radical r, Fractional r, Show r, Eq (FaceIx planeGraph) + , Point_ vertex 2 r, Ord r, Real r + , Fractional r, Show r, Eq (FaceIx planeGraph) , Show (Vertex planeGraph), Show (Dart planeGraph), Show (Face planeGraph) + , Show (EdgeIx planeGraph) ) => planeGraph -> [IpeObject r] drawGraph gr = theVertices <> theEdges <> theFaces where theVertices = ifoldMapOf vertices drawVertex gr theEdges = ifoldMapOf dartSegments (drawDart gr) gr + <> ifoldMapOf edgeSegments (drawEdge gr) gr theFaces = ifoldMapOf interiorFacePolygons (drawFace gr) gr drawVertex :: ( Point_ vertex 2 r, Show vertex) @@ -110,8 +117,18 @@ drawVertex _ v = [ iO $ ipeDiskMark (v^.asPoint) ! attr SLayer "vertex" -- ! attr SStroke Ipe.red ] +drawEdge :: ( PlaneGraph_ planeGraph vertex, Point_ vertex 2 r, IsTransformable vertex + , Show (EdgeIx planeGraph), Fractional r, Real r) + => planeGraph -> EdgeIx planeGraph -> ClosedLineSegment vertex -> [IpeObject r] +drawEdge gr d s = [ iO $ ipeLineSegment s ! attr SLayer "edges" + , iO $ ipeLabel (tshow d :+ c) ! attr SLayer "edgeLabel" + ] + where + c = interpolate 0.5 s ^. asPoint + + drawDart :: ( PlaneGraph_ planeGraph vertex, Point_ vertex 2 r, IsTransformable vertex - , Show (Dart planeGraph), Radical r, Fractional r) + , Show (Dart planeGraph), Fractional r, Real r) => planeGraph -> DartIx planeGraph -> ClosedLineSegment vertex -> [IpeObject r] drawDart gr d s = [ iO $ ipeLineSegment (offset s) ! attr SArrow normalArrow @@ -123,23 +140,27 @@ drawDart gr d s = [ iO $ ipeLineSegment (offset s) c = interpolate 0.5 s ^. asPoint -- computes the midpoint of the segment. +offset :: forall lineSegment point r. + (LineSegment_ lineSegment point, IsTransformable lineSegment + , HasSupportingLine lineSegment + , Point_ point 2 r, Real r, Fractional r) + => lineSegment -> lineSegment offset s = translateBy theOffset s where - theOffset = negated $ signorm v + theOffset :: Vector 2 r + theOffset = fmap realToFrac . negated $ signorm (realToFrac @_ @Double <$> v) + v :: Vector 2 r v = perpendicularTo (supportingLine s) ^. direction - - - - drawFace :: ( PlaneGraph_ planeGraph vertex, Point_ vertex 2 r , Show (Face planeGraph), Ord r, Fractional r) - => planeGraph -> FaceIx planeGraph -> SimplePolygon vertex -> [IpeObject r] -drawFace gr f pg = [ iO $ ipePolygon (pg&vertices %~ (^.asPoint)) ! attr SLayer "face" + => planeGraph -> FaceIx planeGraph -> SimplePolygon (vertex :+ VertexIx planeGraph) -> [IpeObject r] +drawFace gr f pg = [ iO $ ipePolygon pg' ! attr SLayer "face" , iO $ ipeLabel (tshow (gr^?!faceAt f) :+ c) ! attr SLayer "faceLabel" ] where - c = centroid pg ^.asPoint + pg' = pg&vertices %~ (^.core.asPoint) + c = centroid pg' tshow :: Show a => a -> Text.Text diff --git a/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateMonotoneSpec.hs b/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateMonotoneSpec.hs index 0102146d8..4604cd321 100644 --- a/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateMonotoneSpec.hs +++ b/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateMonotoneSpec.hs @@ -2,53 +2,90 @@ {-# LANGUAGE QuasiQuotes #-} module Polygon.Triangulation.TriangulateMonotoneSpec (spec) where -import Control.Lens -import Data.Maybe -import Golden -import HGeometry -import HGeometry.Ext -import HGeometry.Number.Real.Rational -import HGeometry.Polygon.Class -import HGeometry.Polygon.Simple -import HGeometry.Polygon.Triangulation.TriangulateMonotone -import Ipe -import System.OsPath -import Test.Hspec -import Test.Util - +import Control.Lens +import Data.Maybe +import Golden +import HGeometry +import HGeometry.Ext +import HGeometry.Number.Real.Rational +import HGeometry.Polygon.Class +import HGeometry.Polygon.Simple +import HGeometry.PlaneGraph +import HGeometry.Polygon.Triangulation +import Ipe +import System.OsPath +import Test.Hspec +import Test.Util + +import qualified HGeometry.Polygon.Triangulation.TriangulateMonotone as TM + +import Debug.Trace -------------------------------------------------------------------------------- type R = RealNumber 5 spec :: Spec -spec = do testCases [osp|test-with-ipe/Polygon/Triangulation/monotone.ipe|] - testCases [osp|test-with-ipe/Polygon/Triangulation/simplepolygon6.ipe|] - -testCases :: OsPath -> Spec -testCases fp = (runIO $ readInput =<< getDataFileName fp) >>= \case +spec = do it "triangulate triangle" $ do + let tri :: SimplePolygon (Point 2 R) + tri = uncheckedFromCCWPoints $ [origin, Point2 10 0, Point2 10 10] + gr = TM.triangulate @() tri + extractDiagonals gr `shouldBe` [] + testCases toMonotoneSpec [osp|test-with-ipe/Polygon/Triangulation/monotone.ipe|] + testCases toSpec [osp|test-with-ipe/Polygon/Triangulation/monotone.ipe|] + testCases toSpec [osp|test-with-ipe/Polygon/Triangulation/simplepolygon6.ipe|] + + +testCases :: (TestCase R -> Spec) -> OsPath -> Spec +testCases toSpec' fp = (runIO $ readInput =<< getDataFileName fp) >>= \case Left e -> it "reading TriangulateMonotone file" $ expectationFailure $ "Failed to read ipe file " ++ show e - Right tcs -> mapM_ toSpec tcs - + Right tcs -> mapM_ toSpec' tcs data TestCase r = TestCase { _polygon :: SimplePolygon (Point 2 r) :+ IpeColor r , _solution :: [ClosedLineSegment (Point 2 r)] } deriving (Show,Eq) +toMonotoneSpec :: (Num r, Ord r, Show r) => TestCase r -> Spec +toMonotoneSpec (TestCase (poly :+ c) sol) = do + describe ("testing monotone polygons of color " ++ show c) $ do + it "comparing monotone diagonals with manual solution" $ do + let algSol = TM.computeDiagonals poly + (naiveSet . map toSeg $ algSol) `shouldBe` naiveSet sol + it "comparing monotone diagonals from the triangualtion with manual solution" $ do + let gr = TM.triangulate @() poly + algSol = traceShow ("edges", gr^..edges.withIndex) + $ extractDiagonals gr + (naiveSet algSol) `shouldBe` naiveSet sol + where + toSeg (Vector2 i j) = ClosedLineSegment (poly^?!vertexAt i) (poly^?!vertexAt j) toSpec :: (Num r, Ord r, Show r) => TestCase r -> Spec -toSpec (TestCase (poly :+ c) sol) = - describe ("testing polygions of color " ++ show c) $ do - it "comparing with manual solution" $ do +toSpec (TestCase (poly :+ c) sol) = do + describe ("testing polygons of color " ++ show c) $ do + it "comparing diagonals with manual solution" $ do let algSol = computeDiagonals poly (naiveSet . map toSeg $ algSol) `shouldBe` naiveSet sol + it "comparing diagonals from the triangualtion with manual solution" $ do + let gr = triangulate @() poly + algSol = traceShow ("edges", gr^..edges.withIndex) + $ extractDiagonals gr + (naiveSet algSol) `shouldBe` naiveSet sol where - naiveSet = NaiveSet . map S toSeg (Vector2 i j) = ClosedLineSegment (poly^?!vertexAt i) (poly^?!vertexAt j) -newtype S r = S (ClosedLineSegment (Point 2 r)) deriving (Show) +naiveSet = NaiveSet . map S + + + +extractDiagonals gr = mapMaybe (\(d,edgeType) -> case edgeType of + Diagonal -> gr^?edgeSegmentAt d + _ -> Nothing + ) $ gr^..edges.withIndex + + +newtype S r = S (ClosedLineSegment (Point 2 r)) deriving newtype (Show) instance (Eq r) => Eq (S r) where (S s) == (S y) = (s^.start == y^.start && s^.end == y^.end) diff --git a/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateSpec.hs b/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateSpec.hs new file mode 100644 index 000000000..b6b803d96 --- /dev/null +++ b/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateSpec.hs @@ -0,0 +1,111 @@ +{-# LANGUAGE OverloadedStrings #-} +{-# LANGUAGE QuasiQuotes #-} +module Polygon.Triangulation.TriangulateSpec (spec) where + +import Control.Lens +-- import Debug.Trace +import HGeometry +import HGeometry.Ext +import HGeometry.Number.Real.Rational +import HGeometry.PlaneGraph +import HGeometry.Polygon.Class +import HGeometry.Polygon.Simple +import HGeometry.Polygon.Triangulation +import qualified HGeometry.Polygon.Triangulation.MakeMonotone as MM +import qualified HGeometry.Polygon.Triangulation.TriangulateMonotone as TM +import HGeometry.Transformation +import Ipe +import PlaneGraph.RenderSpec (drawVertex, drawDart, drawEdge) +-- import System.OsPath +import Test.Hspec +-- import Test.QuickCheck +-- import Test.QuickCheck.Instances () + +-------------------------------------------------------------------------------- + +type R = RealNumber 5 + +spec :: Spec +spec = describe "triangulateSpec" $ do + it "buggy polygon monotone area" $ do + let g = TM.triangulate @() buggyPolygon + trigs = graphPolygons g + sum (map area trigs) `shouldBe` area buggyPolygon + -- it "all monotone polygons" + it "monotone diags of buggypolygon2'" $ + MM.computeDiagonals buggyPolygon2' `shouldBe` [Vector2 0 3] + it "buggyPolygon2' diagonals" $ + computeDiagonals buggyPolygon2' `shouldBe` [Vector2 0 3,Vector2 1 3] + -- let g = triangulate @() buggyPolygon2' + -- trigs = graphPolygons g + -- runIO $ do + -- putStrLn "==============g buggy2 ====================" + -- print g + -- putStrLn "==============dual g buggy 2 ====================" + -- print $ dualGraph g + -- putStrLn "==================================" + -- writeIpeFile [osp|/tmp/out.ipe|] . singlePageFromContent $ drawGraph $ scaleUniformlyBy 10 g + it "buggyPolygon2' monotone area" $ do + let g = triangulate @() buggyPolygon2' + trigs = graphPolygons g + sum (map area trigs) `shouldBe` area buggyPolygon2' + it "buggyPolygon2 monotone area" $ do + let g = triangulate @() buggyPolygon2 + trigs = graphPolygons g + sum (map area trigs) `shouldBe` area buggyPolygon2 + +-------------------------------------------------------------------------------- + +_drawGraph :: ( PlaneGraph_ planeGraph vertex + , IsTransformable vertex + , Point_ vertex 2 r, Ord r, Real r, Fractional r, Show r, Eq (FaceIx planeGraph) + , Show (Vertex planeGraph), Show (Dart planeGraph), Show (Face planeGraph) + , Show (EdgeIx planeGraph) + ) => planeGraph -> [IpeObject r] +_drawGraph gr = theVertices <> theEdges <> theFaces + where + theVertices = ifoldMapOf vertices drawVertex gr + theEdges = ifoldMapOf dartSegments (drawDart gr) gr + <> ifoldMapOf edgeSegments (drawEdge gr) gr + theFaces = [] -- ifoldMapOf interiorFacePolygons (drawFace gr) gr + +graphPolygons :: (Ord r, Num r, Point_ point 2 r) + => PlaneGraph s point PolygonEdgeType PolygonFaceData + -> [SimplePolygon (Point 2 r)] +graphPolygons gr = map (&vertices %~ view (core.asPoint)) $ gr^..interiorFacePolygons + +-------------------------------------------------------------------------------- +buggyPolygon :: SimplePolygon (Point 2 R) +buggyPolygon = read "SimplePolygon [Point2 9 9,Point2 3 6,Point2 0 3,Point2 2 3,Point2 0 0]" + + +buggyPolygon2 :: SimplePolygon (Point 2 R) +buggyPolygon2 = read "SimplePolygon [Point2 0 11,Point2 52 0,Point2 68 13,Point2 38 12,Point2 15 16]" + +buggyPolygon2' :: SimplePolygon (Point 2 R) +buggyPolygon2' = read "SimplePolygon [Point2 0 11,Point2 52 0,Point2 68 15,Point2 38 12,Point2 15 16]" + + +-------------------------------------------------------------------------------- + +-- spec :: Spec +-- spec = do testCases [osp|test-with-ipe/Polygon/Triangulation/monotone.ipe|] +-- -- testCases [osp|test-with-ipe/Polygon/Triangulation/simplepolygon6.ipe|] + +-- testCases :: OsPath -> Spec +-- testCases fp = (runIO $ readInput =<< getDataFileName fp) >>= \case +-- Left e -> it "reading TriangulateMonotone file" $ +-- expectationFailure $ "Failed to read ipe file " ++ show e +-- Right tcs -> mapM_ toSpec tcs + +-- data TestCase r = TestCase { _polygon :: SimplePolygon (Point 2 r) :+ IpeColor r +-- , _solution :: [ClosedLineSegment (Point 2 r)] +-- } +-- deriving (Show,Eq) + +-- toSpec :: (Num r, Ord r, Show r) => TestCase r -> Spec +-- toSpec (TestCase (poly :+ c) sol) = +-- describe ("testing polygions of color " ++ show c) $ do +-- it "comparing with manual solution" $ do +-- let algSol = triangulate @() poly +-- (naiveSet . map toSeg $ algSol) `shouldBe` naiveSet sol diff --git a/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateSpec_.hs b/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateSpec_.hs deleted file mode 100644 index 56079ebdc..000000000 --- a/hgeometry/test-with-ipe/test/Polygon/Triangulation/TriangulateSpec_.hs +++ /dev/null @@ -1,32 +0,0 @@ -{-# LANGUAGE OverloadedStrings #-} -module HGeometry.Polygon.Triangulation.TriangulateSpec (spec) where - -import Control.Lens -import qualified Data.Vector as V -import HGeometry -import HGeometry.Ext -import HGeometry.PlanarSubdivision (PolygonFaceData) -import HGeometry.PlaneGraph -import HGeometry.PolygonSpec () -import HGeometry.PolygonTriangulation.Triangulate -import HGeometry.PolygonTriangulation.Types -import Test.Hspec -import Test.QuickCheck -import Test.QuickCheck.Instances () - -spec :: Spec -spec = do - it "sum (map area (triangulate polygon)) == area polygon" $ do - property $ \(poly :: SimplePolygon () Rational) -> - let g = triangulate' @() poly - trigs = graphPolygons g - in sum (map area trigs) === area poly - it "all isTriangle . triangulate" $ do - property $ \(poly :: SimplePolygon () Rational) -> - let g = triangulate' @() poly - trigs = graphPolygons g - in all isTriangle trigs - -graphPolygons :: (Ord r, Fractional r) - => PlaneGraph s p PolygonEdgeType PolygonFaceData r -> [SimplePolygon p r] -graphPolygons g = map (^._2.core) . V.toList . snd $ facePolygons (outerFaceId g) g diff --git a/hgeometry/test/Polygon/TriangulateSpec.hs b/hgeometry/test/Polygon/TriangulateSpec.hs index a9413cc63..cda5c7c14 100644 --- a/hgeometry/test/Polygon/TriangulateSpec.hs +++ b/hgeometry/test/Polygon/TriangulateSpec.hs @@ -1,36 +1,71 @@ {-# LANGUAGE OverloadedStrings #-} module Polygon.TriangulateSpec (spec) where --- import Control.Lens --- import qualified Data.Vector as V --- import HGeometry --- import HGeometry.Ext --- import HGeometry.PlanarSubdivision (PolygonFaceData) --- import HGeometry.PlaneGraph --- import HGeometry.PolygonSpec () --- import HGeometry.PolygonTriangulation.Triangulate --- import HGeometry.PolygonTriangulation.Types +import Control.Lens +import HGeometry +import HGeometry.Ext +import HGeometry.Number.Real.Rational +import HGeometry.PlaneGraph +import HGeometry.Polygon.Class +import HGeometry.Polygon.Instances () +import HGeometry.Polygon.Simple +import HGeometry.Polygon.Triangulation +import qualified HGeometry.Polygon.Triangulation.TriangulateMonotone as TM import Test.Hspec --- import Test.Hspec.QuickCheck +import Test.Hspec.QuickCheck +import Test.QuickCheck import Test.QuickCheck.Instances () + + +-- import Debug.Trace -------------------------------------------------------------------------------- +type R = RealNumber 5 + + + spec :: Spec spec = do - pure () - --- prop "sum (map area (triangulate polygon)) == area polygon" $ --- \(poly :: SimplePolygon () Rational) -> --- let g = triangulate' @() poly --- trigs = graphPolygons g --- in sum (map area trigs) === area poly --- prop "all isTriangle . triangulate" $ --- \(poly :: SimplePolygon () Rational) -> --- let g = triangulate' @() poly --- trigs = graphPolygons g --- in all isTriangle trigs - --- graphPolygons :: (Ord r, Fractional r) --- => PlaneGraph s p PolygonEdgeType PolygonFaceData r -> [SimplePolygon p r] --- graphPolygons g = map (^._2.core) . V.toList . snd $ facePolygons (outerFaceId g) g + it "buggy polygon diagionals" $ + computeDiagonals buggyPolygon `shouldBe` [ Vector2 3 1 + , Vector2 3 0 + ] + it "buggy polygon monotone area" $ + let g = TM.triangulate @() buggyPolygon + trigs = graphPolygons g + in sum (map area trigs) `shouldBe` area buggyPolygon + + it "buggy polygon area" $ + let g = triangulate @() buggyPolygon + trigs = graphPolygons g + in sum (map area trigs) `shouldBe` area buggyPolygon + prop "sum (map area (triangulate polygon)) == area polygon" $ + \(poly :: SimplePolygon (Point 2 R)) -> + let g = triangulate @() poly + trigs = graphPolygons g + in counterexample (show g) $ sum (map area trigs) === area poly + prop "all isTriangle . triangulate" $ + \(poly :: SimplePolygon (Point 2 R)) -> + let g = triangulate @() poly + trigs = graphPolygons g + in all isTriangle trigs + + -- prop "creating the graph does not create additional diagionals" $ + -- \poly -> + -- let gr = triangulate @() poly + -- algSol = extractDiagonals gr + -- in naiveSet algoSol `shouldBe` naiveSet (computeDiagionals poly) + + + +buggyPolygon :: SimplePolygon (Point 2 R) +buggyPolygon = read "SimplePolygon [Point2 9 13,Point2 3 6,Point2 0 3,Point2 2 3,Point2 0 0]" + +isTriangle :: SimplePolygon (Point 2 R) -> Bool +isTriangle = (== 3) . numVertices + +graphPolygons :: (Ord r, Num r, Point_ point 2 r) + => PlaneGraph s point PolygonEdgeType PolygonFaceData + -> [SimplePolygon (Point 2 r)] +graphPolygons gr = map (&vertices %~ view (core.asPoint)) $ gr^..interiorFacePolygons diff --git a/hgeometry/vector/src-quickcheck/HGeometry/Vector/Instances.hs b/hgeometry/vector/src-quickcheck/HGeometry/Vector/Instances.hs index 8c55ec0f1..99c374630 100644 --- a/hgeometry/vector/src-quickcheck/HGeometry/Vector/Instances.hs +++ b/hgeometry/vector/src-quickcheck/HGeometry/Vector/Instances.hs @@ -13,16 +13,18 @@ -------------------------------------------------------------------------------- module HGeometry.Vector.Instances where --- import Control.Lens --- import D --- import HGeometry.Vector.Class -import Test.QuickCheck +import Control.Lens import HGeometry.Vector +import Test.QuickCheck -------------------------------------------------------------------------------- instance (Arbitrary r, Vector_ (Vector d r) d r) => Arbitrary (Vector d r) where arbitrary = generateA (const arbitrary) + shrink v = [ v&component' i .~ x' + | (i,x) <- v^..components.withIndex + , x' <- shrink x + ] -- instance ( forall r. VectorLike_ (Vector d r) -- ) => Arbitrary1 (Vector d) where