diff --git a/LICENSE b/LICENSE new file mode 100644 --- /dev/null +++ b/LICENSE @@ -0,0 +1,20 @@ +Copyright (c) 2015 Peter France + +Permission is hereby granted, free of charge, to any person obtaining +a copy of this software and associated documentation files (the +"Software"), to deal in the Software without restriction, including +without limitation the rights to use, copy, modify, merge, publish, +distribute, sublicense, and/or sell copies of the Software, and to +permit persons to whom the Software is furnished to do so, subject to +the following conditions: + +The above copyright notice and this permission notice shall be included +in all copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. +IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY +CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, +TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE +SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. diff --git a/cbits/noticehandlers.c b/cbits/noticehandlers.c new file mode 100644 --- /dev/null +++ b/cbits/noticehandlers.c @@ -0,0 +1,28 @@ +#include +#include +#include +#include + +void +geos_notice_handler(const char *fmt, ...) +{ + va_list ap; + va_start(ap, fmt); + vprintf(fmt, ap); + va_end(ap); + +} + +GEOSContextHandle_t init_GEOS() { + +#if GEOS_VERSION_MAJOR > 3 && GEOS_VERSION_MINOR > 4 + GEOSContextHandle_t handle = GEOS_init_r(); + GEOSContext_setNoticeHandler_r(handle, geos_notice_handler); + GEOSContext_setErrorHandler_r(handle, geos_notice_handler); + return handle; +#else + return initGEOS_r(geos_notice_handler, geos_notice_handler); +#endif +} + + diff --git a/cbits/noticehandlers.h b/cbits/noticehandlers.h new file mode 100644 --- /dev/null +++ b/cbits/noticehandlers.h @@ -0,0 +1,6 @@ +#include + +void geos_notice_handler(const char *fmt, ...); +void geos_error_handler(const char *fmt, ...); +GEOSContextHandle_t init_GEOS(); + diff --git a/geos.cabal b/geos.cabal new file mode 100644 --- /dev/null +++ b/geos.cabal @@ -0,0 +1,115 @@ +name: geos + +version: 0.1.0.0 + +synopsis: Bindings for GEOS. + +description: This is a Haskell binding to Geos, the open-source geometry library, which includes geometry types, predicate functions and other operations, spatially indexed geometries, and parsers for WKB and WKT formats. + +license: MIT + +license-file: LICENSE + +author: Peter France + +maintainer: pfrance@gmail.com + +category: Geometry + +build-type: Simple +Extra-source-files: cbits/noticehandlers.h +cabal-version: >=1.18.0 + +Library + Includes: geos_c.h + Extra-libraries: geos_c + cc-options: -fPIC + c-sources: cbits/noticehandlers.c + include-dirs: cbits + GHC-options: -Wall + build-depends: + base <= 4.10.1 + , bytestring + , vector + , transformers + , mtl + , cassava + build-tools: hsc2hs + exposed-modules: + Data.Geometry.Geos.Topology + , Data.Geometry.Geos.Geometry + , Data.Geometry.Geos.Types + , Data.Geometry.Geos.Serialize + , Data.Geometry.Geos.Buffer + , Data.Geometry.Geos.Prepared + + other-modules: + Data.Geometry.Geos.Raw.Base + , Data.Geometry.Geos.Raw.CoordSeq + , Data.Geometry.Geos.Raw.Geometry + , Data.Geometry.Geos.Raw.Serialize + , Data.Geometry.Geos.Raw.Topology + , Data.Geometry.Geos.Raw.Buffer + , Data.Geometry.Geos.Raw.Prepared + , Data.Geometry.Geos.Raw.Internal + + hs-source-dirs: src + frameworks: GEOS + default-extensions: + LambdaCase + , GADTs + , ExistentialQuantification + , DeriveDataTypeable + , RankNTypes + , DataKinds + , KindSignatures + --, ScopedTypeVariables + + default-language: Haskell2010 + +test-suite test + Includes: geos_c.h + Extra-libraries: geos_c + cc-options: -fPIC + c-sources: cbits/noticehandlers.c + include-dirs: cbits + type: + exitcode-stdio-1.0 + default-language: Haskell2010 + hs-source-dirs: tests src + frameworks: GEOS + build-tools: hsc2hs + main-is: Main.hs + build-depends: + base + , bytestring + , vector + , mtl + , geos + , hspec + , cassava + default-extensions: + LambdaCase + , GADTs + , ExistentialQuantification + , DeriveDataTypeable + , RankNTypes + , DataKinds + , KindSignatures + other-modules: + Data.Geometry.Geos.Geometry + , Data.Geometry.Geos.Raw.Internal + , Data.Geometry.Geos.Raw.Base + , Data.Geometry.Geos.Raw.CoordSeq + , Data.Geometry.Geos.Raw.Geometry + , Data.Geometry.Geos.Raw.Serialize + , Data.Geometry.Geos.Serialize + , Data.Geometry.Geos.Types + , ParsingSpec + , RawGeometrySpec + , SpatialOperationsSpec + , SpecSampleData + +source-repository head + type: git + location: https://github.com/ewestern/geos diff --git a/src/Data/Geometry/Geos/Buffer.hs b/src/Data/Geometry/Geos/Buffer.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Buffer.hs @@ -0,0 +1,76 @@ +{-# LANGUAGE ScopedTypeVariables #-} + +module Data.Geometry.Geos.Buffer ( + buffer + , defaultBufferParams + , bufferWithParams + , BufferParams (..) + , BufferJoinStyle (..) + , BufferCapStyle (..) +) where +import qualified Data.Geometry.Geos.Raw.Buffer as R +import Data.Geometry.Geos.Raw.Base +import Data.Geometry.Geos.Types +import Data.Geometry.Geos.Geometry +import qualified Data.Geometry.Geos.Raw.Geometry as RG + +data BufferCapStyle = + RoundCap + | SquareCap + | FlatCap + +data BufferJoinStyle = + RoundJoin + | MitreJoin + | BevelJoin + +convertCapStyle :: BufferCapStyle -> R.BufferCapStyle +convertCapStyle = \case + RoundCap -> R.capRound + SquareCap -> R.capSquare + FlatCap -> R.capFlat + +convertJoinStyle :: BufferJoinStyle -> R.BufferJoinStyle +convertJoinStyle = \case + RoundJoin -> R.joinRound + MitreJoin -> R.joinMitre + BevelJoin -> R.joinBevel + +data BufferParams = BufferParams { + joinStyle :: BufferJoinStyle + , capStyle :: BufferCapStyle + , mitreLimit :: Double + -- | The default number of facets into which to divide a fillet of 90 degrees. + , quadrantSegments :: Int + -- | True for single-sided, False otherwise + , singleSided :: Bool +} + +defaultBufferParams :: BufferParams +defaultBufferParams = BufferParams { + joinStyle = RoundJoin + , capStyle = RoundCap + , mitreLimit = 2.0 + , quadrantSegments = 8 + , singleSided = False +} + +-- | Returns a Geometry that represents all points whose distance from this geometry is less than or equal to the given width. The quadsegs argument sets the number of segments used to approximate a quarter circle. + +buffer :: Geometry a -> Double -> Int -> Some Geometry +buffer g width quadsegs = runGeos $ do + rg :: RG.Geom <- convertGeometryToRaw g + rg' <- R.buffer rg width quadsegs + convertGeometryFromRaw rg' + +bufferWithParams :: Geometry a -> Double -> BufferParams -> Some Geometry +bufferWithParams g width bp = runGeos $ do + rg :: RG.Geom <- convertGeometryToRaw g + rbp <- R.createBufferParams + R.setEndCapStyle rbp (convertCapStyle $ capStyle bp) + R.setJoinStyle rbp (convertJoinStyle $ joinStyle bp) + R.setMitreLimit rbp (mitreLimit bp) + R.setQuadrantSegments rbp (quadrantSegments bp) + R.setSingleSided rbp (singleSided bp) + rg' <- R.bufferWithParams rg rbp width + convertGeometryFromRaw rg' diff --git a/src/Data/Geometry/Geos/Geometry.hs b/src/Data/Geometry/Geos/Geometry.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Geometry.hs @@ -0,0 +1,329 @@ +{-# LANGUAGE LambdaCase, ScopedTypeVariables #-} +{-# LANGUAGE ExistentialQuantification #-} +{-# LANGUAGE FlexibleContexts #-} +{-# OPTIONS_GHC -fno-warn-orphans #-} + +module Data.Geometry.Geos.Geometry ( + convertGeometryFromRaw + , convertGeometryToRaw + , convertMultiPolygonFromRaw + , ensurePoint + , ensureLineString + , ensureLinearRing + , ensurePolygon + , ensureMultiPoint + , ensureMultiPolygon + , ensureMultiLineString + , interpolate + , interpolateNormalized + , project + , projectNormalized + , covers + , coveredBy + , equalsExact + , equals + , overlaps + , contains + , within + , crosses + , touches + , disjoint + , area + , geometryLength + , distance + , hausdorffDistance + , nearestPoints + +) where + +import Data.Geometry.Geos.Types +import qualified Data.Vector as V +import qualified Data.Geometry.Geos.Raw.Geometry as R +import qualified Data.Geometry.Geos.Raw.CoordSeq as RC +import Data.Geometry.Geos.Raw.Base +import Control.Monad + +-- | Returns the distance from the origin of LineString to the point projected on the geometry (that is to a point of the line the closest to the given point). +project :: Geometry LineString -> Geometry Point -> Double +project g1 g2 = runGeos $ do + g1' :: R.Geom <- convertGeometryToRaw g1 + g2' :: R.Geom <- convertGeometryToRaw g2 + R.project g1' g2' + +-- | Like @project@, but returns the distance as a Double between 0 and 1. +projectNormalized :: Geometry LineString -> Geometry Point -> Double +projectNormalized g1 g2 = runGeos $ do + g1' :: R.Geom <- convertGeometryToRaw g1 + g2' :: R.Geom <- convertGeometryToRaw g2 + R.projectNormalized g1' g2' + +-- | Given a distance, returns the point (or closest point) within the geometry LineString that distance. +interpolate :: Geometry LineString -> Double -> Geometry Point +interpolate g d = runGeos $ do + g' :: R.Geom <- convertGeometryToRaw g + s <- R.getSRID g' + p <- convertPointFromRaw =<< (R.interpolate g' $ realToFrac d) + return $ PointGeometry p s + + +-- | Like @interpolate@, but takes the distance as a double between 0 and 1. +interpolateNormalized :: Geometry LineString -> Double -> Geometry Point +interpolateNormalized g d = runGeos $ do + g' :: R.Geom <- convertGeometryToRaw g + s <- R.getSRID g' + p <- convertPointFromRaw =<< (R.interpolateNormalized g' $ realToFrac d) + return $ PointGeometry p s + +binaryPredicate_ :: (R.GeomConst -> R.GeomConst -> Geos Bool) + -> Geometry a + -> Geometry b + -> Bool +binaryPredicate_ f g1 g2 = runGeos . join $ (f <$> convertGeometryToRaw g1 <*> convertGeometryToRaw g2) + +instance Relatable (Geometry a) where + disjoint = binaryPredicate_ R.disjoint + touches = binaryPredicate_ R.touches + intersects = binaryPredicate_ R.intersects + contains = binaryPredicate_ R.contains + within = binaryPredicate_ R.within + crosses = binaryPredicate_ R.crosses + overlaps = binaryPredicate_ R.overlaps + covers = binaryPredicate_ R.covers + coveredBy = binaryPredicate_ R.coveredBy + +-- | Returns True if the DE-9IM intersection matrix for the two geometries is T*F**FFF*. +equals :: Geometry a -> Geometry a -> Bool +equals = binaryPredicate_ R.equals + +-- | Returns True if the two geometries are exactly equal, up to a specified tolerance. The tolerance value should be a floating point number representing the error tolerance in the comparison, e.g., @equalsExact g1 g2 0.001 @ will compare equality to within one thousandth of a unit. +equalsExact :: Geometry a -> Geometry a -> Double -> Bool +equalsExact g1 g2 d = binaryPredicate_ (\g1' g2' -> R.equalsExact g1' g2' d) g1 g2 + +convertGeometryToRaw :: (R.Geometry a, R.CoordSeqInput a ~ cb, RC.CoordinateSequence cb) => Geometry b -> Geos a +convertGeometryToRaw = \case + PointGeometry pg s -> convertPointToRaw pg s + LineStringGeometry lsg s -> convertLineStringToRaw lsg s + LinearRingGeometry lg s -> convertLinearRingToRaw lg s + PolygonGeometry pg s -> convertPolygonToRaw pg s + MultiPointGeometry mp s -> convertMultiPointToRaw mp s + MultiLineStringGeometry ml s -> convertMultiLineStringToRaw ml s + MultiPolygonGeometry mp s -> convertMultiPolygonToRaw mp s + +convertPointToRaw :: R.Geometry b => Point -> Maybe Int -> Geos b +convertPointToRaw (Point c) s = do + cs <- RC.createEmptyCoordinateSequence 1 (dimensionsCoordinate c) + setCoordinateSequence cs 0 c + R.createPoint cs >>= R.setSRID s + + +convertLinearRingToRaw :: R.Geometry b => LinearRing -> SRID -> Geos b +convertLinearRingToRaw (LinearRing cs) s = do + csr <- RC.createEmptyCoordinateSequence len (dimensionsCoordinateSequence cs) + V.zipWithM_ (setCoordinateSequence csr) (V.enumFromN 0 len) cs + R.createLinearRing csr >>= R.setSRID s + where + len = V.length cs + +convertLineStringToRaw :: R.Geometry b => LineString -> SRID -> Geos b +convertLineStringToRaw (LineString cs) s = do + csr <- RC.createEmptyCoordinateSequence len (dimensionsCoordinateSequence cs) + V.zipWithM_ (setCoordinateSequence csr) ( V.enumFromN 0 len) cs + R.createLineString csr >>= R.setSRID s + where + len = V.length cs + +convertPolygonToRaw :: R.Geometry a => Polygon -> SRID -> Geos a +convertPolygonToRaw (Polygon lrs) s = do + ext <- convertLinearRingToRaw (V.head lrs) s + inn <- (flip convertLinearRingToRaw $ s) `V.mapM` V.tail lrs + R.createPolygon ext (V.toList inn) >>= R.setSRID s + +convertMultiPointToRaw :: R.Geometry a => MultiPoint -> SRID -> Geos a +convertMultiPointToRaw (MultiPoint vp) s = do + vr <- (flip convertPointToRaw $ s) `V.mapM` vp + R.createMultiPoint (V.toList vr) >>= R.setSRID s + +convertMultiLineStringToRaw :: R.Geometry a => MultiLineString -> SRID -> Geos a +convertMultiLineStringToRaw (MultiLineString vl) s = do + vr <- (flip convertLineStringToRaw $ s) `V.mapM` vl + R.createMultiLineString (V.toList vr) >>= R.setSRID s + +convertMultiPolygonToRaw :: R.Geometry a => MultiPolygon -> SRID -> Geos a +convertMultiPolygonToRaw (MultiPolygon vp) s = do + vr <- (flip convertPolygonToRaw $ s) `V.mapM` vp + R.createMultiPolygon (V.toList vr) >>= R.setSRID s + +setCoordinateSequence :: RC.CoordinateSequence a => a -> Int -> Coordinate -> Geos () +setCoordinateSequence cs i (Coordinate2 x y) = + RC.setCoordinateSequenceX cs i x >> RC.setCoordinateSequenceY cs i y + +setCoordinateSequence cs i (Coordinate3 x y z) = + RC.setCoordinateSequenceX cs i x >> RC.setCoordinateSequenceY cs i y >> RC.setCoordinateSequenceZ cs i z + +--- Conversions +-- +convertGeometryFromRaw :: (R.Geometry a, R.CoordSeqInput a ~ cb, RC.CoordinateSequence cb) => a -> Geos (Some Geometry) +convertGeometryFromRaw rg = do + s <- R.getSRID rg + tid <- R.getTypeId rg + case tid of + R.PointTypeId -> do + p <- convertPointFromRaw rg + return $ Some $ (PointGeometry p s) + R.LineStringTypeId -> do + l <- convertLineStringFromRaw rg + return $ Some (LineStringGeometry l s) + R.LinearRingTypeId -> do + l <- convertLinearRingFromRaw rg + return $ Some (LinearRingGeometry l s) + R.PolygonTypeId -> do + p <- convertPolygonFromRaw rg + return $ Some (PolygonGeometry p s) + R.MultiPointTypeId -> do + mp <- convertMultiPointFromRaw rg + return $ Some (MultiPointGeometry mp s) + R.MultiLineStringTypeId -> do + ml <- convertMultiLineStringFromRaw rg + return $ Some (MultiLineStringGeometry ml s) + R.MultiPolygonTypeId -> do + mp <- convertMultiPolygonFromRaw rg + return $ Some (MultiPolygonGeometry mp s) + R.GeometryCollectionTypeId -> error "GeometryCollection currently unsupported" + + +-- The following methods are useful when the type of a (Some Geometry) is known a priori +-- (i.e. the result of calling centroid is always a point) + +ensurePoint :: Some Geometry -> Geometry Point +ensurePoint g = withSomeGeometry g $ \g' -> case g' of + PointGeometry _ _ -> g' + _ -> error "This geometry was expected to be a Point" + +ensureLineString :: Some Geometry -> Geometry LineString +ensureLineString g = withSomeGeometry g $ \g' -> case g' of + LineStringGeometry _ _ -> g' + _ -> error "This geometry was expected to be a LineString" + +ensureLinearRing :: Some Geometry -> Geometry LinearRing +ensureLinearRing g = withSomeGeometry g $ \g' -> case g' of + LinearRingGeometry _ _ -> g' + _ -> error "This geometry was expected to be a LinearRing" + +ensurePolygon :: Some Geometry -> Geometry Polygon +ensurePolygon g = withSomeGeometry g $ \g' -> case g' of + PolygonGeometry _ _ -> g' + _ -> error "This geometry was expected to be a Polygon" + +ensureMultiPoint :: Some Geometry -> Geometry MultiPoint +ensureMultiPoint g = withSomeGeometry g $ \p' -> case p' of + MultiPointGeometry _ _ -> p' + _ -> error "This geometry was expected to be a MultiPoint" + +ensureMultiLineString :: Some Geometry -> Geometry MultiLineString +ensureMultiLineString g = withSomeGeometry g $ \p' -> case p' of + MultiLineStringGeometry _ _ -> p' + _ -> error "This geometry was expected to be a MultiLineString" + +ensureMultiPolygon :: Some Geometry -> Geometry MultiPolygon +ensureMultiPolygon g = withSomeGeometry g $ \p' -> case p' of + MultiPolygonGeometry _ _ -> p' + _ -> error "This geometry was expected to be a MultiPolygon" + + +getPosition :: RC.CoordinateSequence a => a -> Int -> Geos Coordinate +getPosition cs i = do + dim <- RC.getCoordinateSequenceDimensions cs + x <- RC.getCoordinateSequenceX cs i + y <- RC.getCoordinateSequenceY cs i + z <- if dim == 3 + then Just <$> RC.getCoordinateSequenceZ cs i + else return Nothing + case z of + Nothing -> return $ Coordinate2 x y + Just z' -> return $ Coordinate3 x y z' + +convertPointFromRaw :: (R.Geometry a, R.CoordSeqInput a ~ ca, RC.CoordinateSequence ca) => a -> Geos Point +convertPointFromRaw g = do + cs <- R.getCoordinateSequence g + Point <$> getPosition cs 0 + + +convertSequenceFromRaw :: (R.Geometry a, R.CoordSeqInput a ~ ca, RC.CoordinateSequence ca) + => a -> Geos CoordinateSequence +convertSequenceFromRaw g = do + cs <- R.getCoordinateSequence g + size <- R.getNumCoordinates g + V.generateM size (getPosition cs) + +convertLineStringFromRaw :: (R.Geometry a, R.CoordSeqInput a ~ ca, RC.CoordinateSequence ca) + => a -> Geos LineString +convertLineStringFromRaw g = LineString <$> convertSequenceFromRaw g + +convertLinearRingFromRaw :: (R.Geometry a, R.CoordSeqInput a ~ ca, RC.CoordinateSequence ca) + => a -> Geos LinearRing +convertLinearRingFromRaw g = LinearRing <$> convertSequenceFromRaw g + +convertPolygonFromRaw :: R.Geometry a => a -> Geos Polygon +convertPolygonFromRaw g = do + is <- R.getNumInteriorRings g + ext :: R.GeomConst <- R.getExteriorRing g + ins :: V.Vector R.GeomConst <- V.generateM is (R.getInteriorRingN g) + Polygon <$> convertLinearRingFromRaw `V.mapM` (ext `V.cons` ins) + +{- +Enforces using GeomConst for following functions + +-} +getGeometryN :: R.Geometry a => a -> Int -> Geos R.GeomConst +getGeometryN = R.getGeometryN + +convertMultiPointFromRaw :: (R.Geometry a, R.CoordSeqInput a ~ ca, RC.CoordinateSequence ca) => a -> Geos MultiPoint +convertMultiPointFromRaw g = do + ng <- R.getNumGeometries g + MultiPoint <$> V.generateM ng (\i -> convertPointFromRaw =<< getGeometryN g i ) + +convertMultiLineStringFromRaw :: (R.Geometry a, R.CoordSeqInput a ~ ca, RC.CoordinateSequence ca) => a -> Geos MultiLineString +convertMultiLineStringFromRaw g = do + ng <- R.getNumGeometries g + MultiLineString <$> V.generateM ng (\i -> convertLineStringFromRaw =<< getGeometryN g i) + +convertMultiPolygonFromRaw :: (R.Geometry a, R.CoordSeqInput a ~ ca, RC.CoordinateSequence ca) => a -> Geos MultiPolygon +convertMultiPolygonFromRaw g = do + ng <- R.getNumGeometries g + MultiPolygon <$> V.generateM ng (\i -> convertPolygonFromRaw =<< getGeometryN g i) + +area :: Geometry a -> Double +area g = runGeos $ do + r :: R.Geom <- convertGeometryToRaw g + d <- R.area r + return d + +-- | Returns the length of this geometry (e.g., 0 for a Point, the length of a LineString, or the circumference of a Polygon). +geometryLength :: Geometry a -> Double +geometryLength g = runGeos $ do + r :: R.Geom <- convertGeometryToRaw g + l <- R.geometryLength r + return l + +-- | NOTE: Data.Geometry.Geos distance calculations are linear – in other words, Data.Geometry.Geos does not perform a spherical calculation even if the SRID specifies a geographic coordinate system. +distance :: Geometry a -> Geometry a -> Double +distance p g = runGeos $ do + p' :: R.Geom <- convertGeometryToRaw p + g' :: R.Geom <- convertGeometryToRaw g + R.distance p' g' + +hausdorffDistance :: Geometry a -> Geometry a -> Double +hausdorffDistance p g = runGeos $ do + p' :: R.Geom <- convertGeometryToRaw p + g' :: R.Geom <- convertGeometryToRaw g + R.hausdorffDistance p' g' + +-- | Returns the closest points of the two geometries. The first point comes from g1 geometry and the second point comes from g2. +nearestPoints :: Geometry a -> Geometry a -> (Coordinate, Coordinate) +nearestPoints g1 g2 = runGeos $ do + g1' :: R.Geom <- convertGeometryToRaw g1 + g2' :: R.Geom <- convertGeometryToRaw g2 + cs :: RC.CoordSeq <- R.nearestPoints g1' g2' + p1 <- getPosition cs 0 + p2 <- getPosition cs 1 + return (p1, p2) diff --git a/src/Data/Geometry/Geos/Prepared.hs b/src/Data/Geometry/Geos/Prepared.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Prepared.hs @@ -0,0 +1,67 @@ +{-# OPTIONS_GHC -fno-warn-orphans #-} +{-# LANGUAGE ScopedTypeVariables #-} +{-| +Module : Data.Geometry.Geos.Prepared + +An interface for classes which prepare Geometrys in order to optimize the performance of repeated calls to specific geometric operations. + +A given implementation may provide optimized implementations for only some of the specified methods, and delegate the remaining methods to the original Geometry operations. An implementation may also only optimize certain situations, and delegate others. See the implementing classes for documentation about which methods and situations they optimize. + +-} + +module Data.Geometry.Geos.Prepared ( + prepare + , contains + , containsProperly + , coveredBy + , covers + , crosses + , disjoint + , intersects + , overlaps + , touches + , within +) where + +import qualified Data.Geometry.Geos.Raw.Prepared as RP +import qualified Data.Geometry.Geos.Raw.Geometry as RG +import Data.Geometry.Geos.Raw.Base +import qualified Data.Geometry.Geos.Geometry as G +import Data.Geometry.Geos.Types + +prepare :: Geometry a -> RP.PreparedGeometry +prepare g = runGeos $ do + r :: RG.Geom <- G.convertGeometryToRaw g + RP.prepare r + + +queryPrepared :: (RP.PreparedGeometry -> RG.GeomConst -> Geos Bool) + -> RP.PreparedGeometry + -> Geometry b + -> Bool +queryPrepared f pg g = runGeos $ G.convertGeometryToRaw g >>= (f pg) + +instance Relatable (RP.PreparedGeometry) where + contains = queryPrepared RP.contains + coveredBy = queryPrepared RP.coveredBy + covers = queryPrepared RP.covers + crosses = queryPrepared RP.crosses + disjoint = queryPrepared RP.disjoint + intersects = queryPrepared RP.intersects + overlaps = queryPrepared RP.overlaps + touches = queryPrepared RP.touches + within = queryPrepared RP.within + +-- | The containsProperly predicate has the following equivalent definitions: + +-- | Every point of the other geometry is a point of this geometry's interior. In other words, if the test geometry has any interaction with the boundary of the target geometry the result of containsProperly is false. This is different semantics to the @contains@ predicate, in which test geometries can intersect the target's boundary and still be contained. + +-- | The advantage of using this predicate is that it can be computed efficiently, since it avoids the need to compute the full topological relationship of the input boundaries in cases where they intersect. + +-- | An example use case is computing the intersections of a set of geometries with a large polygonal geometry. Since intersection is a fairly slow operation, it can be more efficient to use containsProperly to filter out test geometries which lie wholly inside the area. In these cases the intersection is known a priori to be exactly the original test geometry. + + +containsProperly :: RP.PreparedGeometry + -> Geometry a + -> Bool +containsProperly = queryPrepared RP.containsProperly diff --git a/src/Data/Geometry/Geos/Raw/Base.hs b/src/Data/Geometry/Geos/Raw/Base.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/Base.hs @@ -0,0 +1,43 @@ +{-# LANGUAGE GeneralizedNewtypeDeriving, RankNTypes #-} + +module Data.Geometry.Geos.Raw.Base ( + Geos + , runGeos + , throwIfZero + , withGeos + , mkErrorMessage +) where + +import qualified Data.Geometry.Geos.Raw.Internal as I +import Foreign +import Data.Monoid ((<>)) +import System.IO.Unsafe +import qualified Control.Concurrent.MVar as MV +import Control.Monad.Reader + +newtype GEOSHandle = GEOSHandle (MV.MVar (ForeignPtr I.GEOSContextHandle_HS)) + +-- idea: stores list of (Handle -> IO) actions +newtype Geos a = Geos { unGeos :: ReaderT GEOSHandle IO a } + deriving (MonadReader GEOSHandle, Monad, Functor, Applicative) +-- don't derive MonadIO to restrict user from performing arbitrary IO + +withGeos :: (I.GEOSContextHandle_t -> IO a) -> Geos a +withGeos f = Geos . ReaderT $ \(GEOSHandle mv) -> MV.withMVar mv $ \fp -> withForeignPtr fp f + +runGeos :: Geos a -> a +runGeos g = unsafePerformIO $ do + ptrC <- I.geos_init + fptr <- newForeignPtr I.geos_finish ptrC + mv <- MV.newMVar fptr + v <- runReaderT (unGeos g) $ GEOSHandle mv + return v + + +throwIfZero :: (Eq a, Num a) => (a -> String) -> IO a -> IO a +throwIfZero = throwIf ((==) 0) + +mkErrorMessage :: Show a => String -> (a -> String) +mkErrorMessage s = \n -> s <> " has thrown an error: " <> show n + + diff --git a/src/Data/Geometry/Geos/Raw/Buffer.hs b/src/Data/Geometry/Geos/Raw/Buffer.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/Buffer.hs @@ -0,0 +1,103 @@ +module Data.Geometry.Geos.Raw.Buffer ( + BufferParams + , withBufferParams + , BufferCapStyle + , BufferJoinStyle + , buffer + , createBufferParams + , setEndCapStyle + , setJoinStyle + , setMitreLimit + , setQuadrantSegments + , setSingleSided + , bufferWithStyle + , bufferWithParams + , offsetCurve + , capRound + , capFlat + , capSquare + , joinRound + , joinMitre + , joinBevel +) where +import Data.Geometry.Geos.Raw.Internal +import qualified Data.Geometry.Geos.Raw.Geometry as RG +import Data.Geometry.Geos.Raw.Base +import Foreign + +newtype BufferParams = BufferParams (ForeignPtr GEOSBufferParams) deriving (Show, Eq) + + +withBufferParams :: BufferParams -> (Ptr GEOSBufferParams -> IO a ) -> IO a +withBufferParams (BufferParams g) f = withForeignPtr g f + + +-- | Create a buffer around a geometry, where quadsegs is the number of line segments to use to approximate a quarter of a circle. +buffer :: RG.Geometry a => a -> Double -> Int -> Geos a +buffer geo width quadsegs = do + withGeos $ \h -> + RG.withGeometry geo $ \gp -> do + g' <- geos_Buffer h gp (realToFrac width) $ fromIntegral quadsegs + RG.constructGeometry h g' + + +createBufferParams :: Geos BufferParams +createBufferParams = withGeos $ \h -> do + bp <- geos_BufferParamsCreate h + fp <- newForeignPtrEnv geos_BufferParamsDestroy h bp + return $ BufferParams fp + + +setEndCapStyle :: BufferParams -> BufferCapStyle -> Geos () +setEndCapStyle b s = withGeos $ \h -> do + _ <- throwIfZero (mkErrorMessage "setEndCapStyle") $ withBufferParams b $ \bp -> + geos_BufferParamsSetEndCapStyle h bp $ unBufferCapStyle s + return () + +setJoinStyle :: BufferParams -> BufferJoinStyle -> Geos () +setJoinStyle b s = withGeos $ \h -> do + _ <- throwIfZero (mkErrorMessage "setJoinStyle") $ withBufferParams b $ \bp -> + geos_BufferParamsSetJoinStyle h bp $ unBufferJoinStyle s + return () + + +setMitreLimit :: BufferParams -> Double -> Geos () +setMitreLimit b d = withGeos $ \h -> do + _ <- throwIfZero (mkErrorMessage "setJoinStyle") $ withBufferParams b $ \bp -> + geos_BufferParamsSetMitreLimit h bp $ realToFrac d + return () + +setQuadrantSegments :: BufferParams -> Int -> Geos () +setQuadrantSegments b i = withGeos $ \h -> do + _ <- throwIfZero (mkErrorMessage "setJoinStyle") $ withBufferParams b $ \bp -> + geos_BufferParamsSetQuadrantSegments h bp $ fromIntegral i + return () + +setSingleSided :: BufferParams -> Bool -> Geos () +setSingleSided bp b = withGeos $ \h -> do + _ <- throwIfZero (mkErrorMessage "setSingleSided") $ withBufferParams bp $ \bpp -> + geos_BufferParamsSetSingleSided h bpp $ fromBool b + return () + +bufferWithParams :: RG.Geometry a => a -> BufferParams -> Double -> Geos a +bufferWithParams g b width = + withGeos $ \h -> do + RG.withGeometry g $ \gp -> + withBufferParams b $ \bp -> do + g' <- throwIfNull "bufferWithParams" $ geos_BufferWithParams h gp bp $ realToFrac width + RG.constructGeometry h g' + +bufferWithStyle :: RG.Geometry a => a -> Double -> Int -> BufferCapStyle -> BufferJoinStyle -> Double -> Geos a +bufferWithStyle g width quadsegs capstyle joinstyle mitrelimit = + withGeos $ \h -> + RG.withGeometry g $ \gp -> do + ptr <- throwIfNull "bufferWithStyle" $ geos_BufferWithStyle h gp (realToFrac width) (fromIntegral quadsegs) (unBufferCapStyle capstyle) (unBufferJoinStyle joinstyle) $ realToFrac mitrelimit + RG.constructGeometry h ptr + +-- | Will only accept LineString geometries. For the 'width' parameter, negative doubles represent a right-side offset, and positive doubles represent a left-side offset. +offsetCurve :: RG.Geometry a => a -> Double -> Int -> BufferJoinStyle -> Double -> Geos a +offsetCurve g width quadsegs joinstyle mitrelimit = + withGeos $ \h -> + RG.withGeometry g $ \gp -> do + ptr <- throwIfNull "offsetCurve" $ geos_OffsetCurve h gp (realToFrac width) (fromIntegral quadsegs) (unBufferJoinStyle joinstyle) $ realToFrac mitrelimit + RG.constructGeometry h ptr diff --git a/src/Data/Geometry/Geos/Raw/CoordSeq.hs b/src/Data/Geometry/Geos/Raw/CoordSeq.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/CoordSeq.hs @@ -0,0 +1,173 @@ +{-# LANGUAGE FlexibleInstances #-} + +module Data.Geometry.Geos.Raw.CoordSeq ( + CoordinateSequence (..) + , CoordSeqConst (CoordSeqConst) + , CoordSeq (CoordSeq) + , getCoordinateSequenceX + , getCoordinateSequenceY + , getCoordinateSequenceZ + , getCoordinateSequenceSize + , getCoordinateSequenceDimensions + , setCoordinateSequenceX + , setCoordinateSequenceY + , setCoordinateSequenceZ + , setCoordinateSequenceOrd +) where + +import qualified Data.Geometry.Geos.Raw.Internal as I +import Data.Geometry.Geos.Raw.Base +import Foreign +import Foreign.C.Types +import Control.Monad + +class CoordinateSequence a where + withCoordinateSequence :: a -> (Ptr I.GEOSCoordSequence -> IO b) -> IO b + createEmptyCoordinateSequence :: Int -> Int -> Geos a + createCoordinateSequence :: Ptr I.GEOSCoordSequence -> Geos a + + + +newtype CoordSeq = CoordSeq { + _unCoordSeq :: ForeignPtr I.GEOSCoordSequence +} +newtype CoordSeqConst = CoordSeqConst { + _unCoordSeqConst :: Ptr I.GEOSCoordSequence +} + +instance Eq CoordSeq where + (==) = coordSeqEq + +instance Eq CoordSeqConst where + (==) = coordSeqEq + +coordSeqEq :: CoordinateSequence a => a -> a -> Bool +coordSeqEq a b = runGeos $ do + sa <- getCoordinateSequenceSize a + sb <- getCoordinateSequenceSize b + da <- getCoordinateSequenceDimensions a + db <- getCoordinateSequenceDimensions b + if (sa == sb) && (da == db) + then foldM (comp (da == 3)) True [0..(sa-1)] + else return False + where + comp zdim acc i = do + xa <- getCoordinateSequenceX a i + ya <- getCoordinateSequenceY a i + xb <- getCoordinateSequenceX b i + yb <- getCoordinateSequenceY b i + zd <- if zdim + then do + za <- getCoordinateSequenceZ a i + zb <- getCoordinateSequenceZ b i + return $ za == zb + else return True + return $ (xa == xb) && (ya == yb) && acc && zd + +instance Show CoordSeq where + show = coordSeqShow + +instance Show CoordSeqConst where + show = coordSeqShow + +coordSeqShow :: CoordinateSequence a => a -> String +coordSeqShow a = runGeos $ do + sa <- getCoordinateSequenceSize a + unlines `fmap` mapM show' [0..(sa-1)] + where + show' i = do + xa <- getCoordinateSequenceX a i + ya <- getCoordinateSequenceY a i + return . show $ (xa, ya) + + +instance CoordinateSequence CoordSeq where + withCoordinateSequence (CoordSeq fp) f = withForeignPtr fp f + createEmptyCoordinateSequence size dim = do + ptr <- withGeos $ \h -> + throwIfNull "createEmptyCoordinateSequence" $ I.geos_CoordSeqCreate h (fromIntegral size) (fromIntegral dim) + createCoordinateSequence ptr + + createCoordinateSequence ptr = withGeos $ \h -> do + fptr <- newForeignPtrEnv I.geos_CoordSeqDestroy h ptr + return $ CoordSeq fptr + +instance CoordinateSequence CoordSeqConst where + withCoordinateSequence (CoordSeqConst p) f = f p + createEmptyCoordinateSequence size dim = do + ptr <- withGeos $ \h -> + throwIfNull "createEmptyCoordinateSequence" $ I.geos_CoordSeqCreate h (fromIntegral size) (fromIntegral dim) + createCoordinateSequence ptr + createCoordinateSequence ptr = return $ CoordSeqConst ptr + + + + +getCoordinateSequenceD_ :: CoordinateSequence a + => (I.GEOSContextHandle_t -> Ptr I.GEOSCoordSequence -> CUInt -> Ptr CDouble -> IO CInt) + -> a + -> Int + -> Geos Double +getCoordinateSequenceD_ f cs idx = withGeos $ \h -> + alloca $ \dptr -> do + _ <- throwIfZero (mkErrorMessage "getCoordiniateSequenceN") $ + withCoordinateSequence cs $ \pcs -> f h pcs (fromIntegral idx) dptr + d <- peek dptr + return $ realToFrac d + +getCoordinateSequenceX :: CoordinateSequence a => a -> Int -> Geos Double +getCoordinateSequenceX = getCoordinateSequenceD_ I.geos_CoordSeqGetX + +getCoordinateSequenceY :: CoordinateSequence a => a -> Int -> Geos Double +getCoordinateSequenceY = getCoordinateSequenceD_ I.geos_CoordSeqGetY + +getCoordinateSequenceZ :: CoordinateSequence a => a -> Int -> Geos Double +getCoordinateSequenceZ = getCoordinateSequenceD_ I.geos_CoordSeqGetZ + +getCoordinateSequenceSize :: CoordinateSequence a => a -> Geos Int +getCoordinateSequenceSize c = withGeos $ \h -> + alloca $ \ptr -> do + _ <- throwIfZero (mkErrorMessage "getCoordinateSequenceSize") $ + withCoordinateSequence c $ \pc -> + I.geos_CoordSeqGetSize h pc ptr + s <- peek ptr + return $ fromIntegral s + +getCoordinateSequenceDimensions :: CoordinateSequence a => a -> Geos Int +getCoordinateSequenceDimensions c = withGeos $ \h -> + alloca $ \ptr -> do + _ <- throwIfZero (mkErrorMessage "getCoordinateSeqenceDimensions") $ + withCoordinateSequence c $ \pc -> + I.geos_CoordSeqGetDimensions h pc ptr + s <- peek ptr + return $ fromIntegral s + +--- +setCoordinateSequence_ :: CoordinateSequence a + => (I.GEOSContextHandle_t -> Ptr I.GEOSCoordSequence -> CUInt -> CDouble -> IO CInt) + -> a + -> Int + -> Double + -> Geos () +setCoordinateSequence_ f cs idx val = withGeos $ \h -> do + _ <- throwIfZero (mkErrorMessage "setCoordinateSEquenceN") $ + withCoordinateSequence cs $ \pcs -> + f h pcs (fromIntegral idx) (realToFrac val) + return () + + +setCoordinateSequenceX :: CoordinateSequence a => a -> Int -> Double -> Geos () +setCoordinateSequenceX = setCoordinateSequence_ I.geos_CoordSeqSetX + +setCoordinateSequenceY :: CoordinateSequence a => a -> Int -> Double -> Geos () +setCoordinateSequenceY = setCoordinateSequence_ I.geos_CoordSeqSetY + +setCoordinateSequenceZ :: CoordinateSequence a => a -> Int -> Double -> Geos () +setCoordinateSequenceZ = setCoordinateSequence_ I.geos_CoordSeqSetZ + +setCoordinateSequenceOrd :: CoordinateSequence a => a -> Int -> Int -> Double -> Geos Int +setCoordinateSequenceOrd cs idx dim v = withGeos $ \h -> do + i <- throwIfZero (mkErrorMessage "setCoordinateSequenceN") $ + withCoordinateSequence cs $ \pcs -> + I.geos_CoordSeqSetOrdinate h pcs (fromIntegral idx) (fromIntegral dim) (realToFrac v) + return $ fromIntegral i diff --git a/src/Data/Geometry/Geos/Raw/Geometry.hs b/src/Data/Geometry/Geos/Raw/Geometry.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/Geometry.hs @@ -0,0 +1,416 @@ +{-# OPTIONS_GHC -fno-warn-orphans #-} +{-# LANGUAGE TypeFamilies #-} +{-# LANGUAGE FlexibleInstances #-} + + +{-| +Module : Data.Geometry.Geos.Raw.Geometery + +Light wrappers around Geos functions. Must be run within the Geos monad. + +-} +module Data.Geometry.Geos.Raw.Geometry ( + Geom (..) + , GeomConst (..) + , Geometry (..) + , GeomTypeId (..) + , withMaybeGeometry + , getSRID + , setSRID + , getTypeName + , getTypeId + , getCoordinateSequence + , getNumCoordinates + , getNumInteriorRings + , getNumGeometries + , getInteriorRingN + , getGeometryN + , getExteriorRing + , createPoint + , createLinearRing + , createLineString + , createPolygon + , createMultiPoint + , createMultiLineString + , createMultiPolygon + , createCollection + + , project + , projectNormalized + , interpolate + , interpolateNormalized + , disjoint + , touches + , intersects + , crosses + , within + , contains + , overlaps + , equals + , equalsExact + , covers + , coveredBy + -- Misc functions + , area + , geometryLength + , distance + , hausdorffDistance + , nearestPoints + , normalize +) where +import qualified Data.Geometry.Geos.Raw.Internal as I +import Data.Geometry.Geos.Raw.Base +import Data.Geometry.Geos.Raw.CoordSeq +import Foreign +import Foreign.Ptr (nullPtr) +import Foreign.C.Types +import Foreign.C.String + +newtype Geom = Geom (ForeignPtr I.GEOSGeometry) + +newtype GeomConst = GeomConst ( Ptr I.GEOSGeometry) + + +class Geometry a where + type CoordSeqInput a + + withGeometry :: a + -> (Ptr I.GEOSGeometry -> IO b ) + -> IO b + constructGeometry :: I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> IO a + +instance Geometry Geom where + type CoordSeqInput Geom = CoordSeqConst + withGeometry (Geom g) f = withForeignPtr g f + constructGeometry h geo = do + fptr <- newForeignPtrEnv I.geos_GeomDestroy h geo + return $ Geom fptr + +instance Geometry GeomConst where + type CoordSeqInput GeomConst = CoordSeq + withGeometry (GeomConst p) f = f p + constructGeometry _ geo = return $ GeomConst geo + +withMaybeGeometry :: Geometry a => Maybe a -> (Ptr I.GEOSGeometry -> IO b) -> IO b +withMaybeGeometry mg f = case mg of + Just g -> withGeometry g f + Nothing -> f nullPtr + +createGeometryFromCoords :: Geometry b + => (I.GEOSContextHandle_t -> Ptr I.GEOSCoordSequence -> IO (Ptr I.GEOSGeometry)) + -> CoordSeqConst + -> Geos b +createGeometryFromCoords f c = + withGeos $ \h -> + withCoordinateSequence c $ \pcs -> do + ptr <- f h pcs + constructGeometry h ptr + +geomEq :: (Eq ca, Geometry a, CoordSeqInput a ~ ca, CoordinateSequence ca) => a -> a -> Bool +geomEq a b = runGeos $ do + sa <- getSRID a + sb <- getSRID b + ta <- getTypeId a + tb <- getTypeId b + if (sa == sb) && (ta == tb) + then do + csa <- getCoordinateSequence a + csb <- getCoordinateSequence b + return $ csa == csb + else return False + + +instance Eq Geom where + a == b = geomEq a b + +instance Eq GeomConst where + a == b = geomEq a b + +getSRID :: Geometry a => a -> Geos (Maybe Int) +getSRID g = withGeos $ \h -> do + s <- withGeometry g $ I.geos_GetSRID h + case fromIntegral s of + 0 -> return Nothing + i -> return (Just i) + +setSRID :: Geometry a => (Maybe Int) -> a -> Geos a +setSRID Nothing g = return g +setSRID (Just i) g = withGeos $ \h -> do + withGeometry g $ \gp -> I.geos_SetSRID h gp $ fromIntegral i + return g + +data GeomTypeId = PointTypeId | LineStringTypeId | LinearRingTypeId | PolygonTypeId + | MultiPointTypeId | MultiLineStringTypeId | MultiPolygonTypeId | GeometryCollectionTypeId deriving (Eq,Show) + +geomTypeId :: Integer -> GeomTypeId +geomTypeId 0 = PointTypeId +geomTypeId 1 = LineStringTypeId +geomTypeId 2 = LinearRingTypeId +geomTypeId 3 = PolygonTypeId +geomTypeId 4 = MultiPointTypeId +geomTypeId 5 = MultiLineStringTypeId +geomTypeId 6 = MultiPolygonTypeId +geomTypeId 7 = GeometryCollectionTypeId +geomTypeId i = error $ "Not a valid geometry type " ++ (show i) + +getTypeName :: Geometry a => a -> Geos String +getTypeName g = withGeos $ \h -> do + s <- throwIfNull "getType" $ + withGeometry g $ I.geos_GeomType h + return =<< peekCString s + +getTypeId ::Geometry a => a -> Geos GeomTypeId +getTypeId g = withGeos $ \h -> do + i <- throwIfNeg (mkErrorMessage "getTypeId") $ withGeometry g $ I.geos_GeomTypeId h + return $ geomTypeId (fromIntegral i) + + +getCoordinateSequence :: Geometry a => a -> Geos CoordSeq +getCoordinateSequence g = withGeos $ \h -> + withGeometry g $ \gptr -> do + cptr <- throwIfNull "getCoordinateSequence" $ I.geos_GetCoordSeq h gptr + cptr' <- I.geos_CoordSeqClone h cptr + fptr <- newForeignPtrEnv I.geos_CoordSeqDestroy h cptr' + return $ CoordSeq fptr + + +getNum_ :: Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> IO CInt) + -> a + -> Geos Int +getNum_ f g = withGeos $ \h -> do + i <- throwIfNeg (mkErrorMessage "getNumCoordinates") $ + withGeometry g $ f h + return $ fromIntegral i + +getNumCoordinates :: Geometry a => a -> Geos Int +getNumCoordinates = getNum_ I.geos_GetNumCoordinates + +---- Polygons +getNumInteriorRings :: Geometry a => a -> Geos Int +getNumInteriorRings = getNum_ I.geos_GetNumInteriorRings + +--- multi geometries +getNumGeometries :: Geometry a => a -> Geos Int +getNumGeometries = getNum_ I.geos_GetNumGeometries + +getN_ :: (Geometry a , Geometry b) + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> CInt -> IO (Ptr I.GEOSGeometry)) + -> a + -> Int + -> Geos b +getN_ f g i = + withGeos $ \h -> withGeometry g $ \gp -> do + gp' <- throwIfNull "getN" $ f h gp $ fromIntegral i + constructGeometry h gp' + + +getGeometryN :: (Geometry a, Geometry b) => a -> Int -> Geos b +getGeometryN = getN_ I.geos_GetGeometryN + +getExteriorRing :: (Geometry a, Geometry b) => a -> Geos b +getExteriorRing g = do + withGeos $ \h -> do + withGeometry g $ \gp -> do + gp' <- throwIfNull "getExteriorRing" $ I.geos_GetExteriorRing h gp + constructGeometry h gp' + + +getInteriorRingN :: (Geometry a, Geometry b) => a -> Int -> Geos b +getInteriorRingN = getN_ I.geos_GetInteriorRingN + +normalize :: Geometry a => a -> Geos a +normalize g = do + cloned <- cloneGeometry g + withGeos $ \h -> do + _ <- throwIfNeg (mkErrorMessage "normalize") $ withGeometry cloned $ I.geos_Normalize h + return () + return cloned +-- + +cloneGeometry :: Geometry a => a -> Geos a +cloneGeometry g = do + withGeos $ \h -> withGeometry g $ \gp -> + I.geos_GeomClone h gp >>= constructGeometry h + + {-Geometry Constructors.-} + {-GEOSCoordSequence* arguments will become ownership of the returned object.-} + {-All functions return NULL on exception.-} + + +-- Geometry Constructors +{-| +The following require CoordSeqConst as arguments since coordinate sequences become owned by the Geometry object. + +-} +createPoint ::Geometry b => CoordSeqConst -> Geos b +createPoint = createGeometryFromCoords I.geos_GeomCreatePoint + +createLinearRing :: Geometry a => CoordSeqConst -> Geos a +createLinearRing = createGeometryFromCoords I.geos_GeomCreateLinearRing + +createLineString ::Geometry b => CoordSeqConst -> Geos b +createLineString = createGeometryFromCoords I.geos_GeomCreateLineString + +-- TODO: Make this take a vector argument + +-- | The second argument is a list of geometries, +-- | NOTE. geometries become owned by caller. +createPolygon :: Geometry a => GeomConst -> [GeomConst] -> Geos a +createPolygon o hs = do + withGeos $ \h -> do + ptrs <- mapM (\v -> withGeometry v $ return) hs + withGeometry o $ \op -> + withArray ptrs $ \ph -> do + g' <- I.geos_GeomCreatePolygon h op ph $ fromIntegral $ length hs + constructGeometry h g' + + +createMulti_ :: Geometry a => I.GEOSGeomType -> [GeomConst] -> Geos a +createMulti_ t gs = do + withGeos $ \h -> do + ptrs <- mapM (\v -> withGeometry v $ return) gs + withArray ptrs $ \ph -> do + g' <- I.geos_GeomCreateCollection h (I.unGEOSGeomType t) ph $ fromIntegral $ length gs + constructGeometry h g' + +createMultiPoint :: Geometry a => [GeomConst] -> Geos a +createMultiPoint = createMulti_ I.multiPointId + +createMultiLineString :: Geometry a => [GeomConst] -> Geos a +createMultiLineString = createMulti_ I.multiLineStringId + +createMultiPolygon :: Geometry a => [GeomConst] -> Geos a +createMultiPolygon = createMulti_ I.multiPolygonId + +createCollection :: Geometry a => [GeomConst] -> Geos a +createCollection = createMulti_ I.geometryCollectionId + + +--- Linear Referencing +---------------------- +geo_2_ :: Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> Ptr I.GEOSGeometry -> IO CDouble) + -> a + -> a + -> Geos Double +geo_2_ f g p = withGeos $ \h -> do + d <- withGeometry g $ \gp -> + withGeometry p $ f h gp + return . realToFrac $ d + +-- | @project p g@ returns the distance of point @p@ projected on @g@ from origin of @g@. Geometry @g@ must be a lineal geometry +-- +project :: Geometry a => a -> a -> Geos Double +project = geo_2_ I.geos_Project + +projectNormalized :: Geometry a => a -> a -> Geos Double +projectNormalized = geo_2_ I.geos_ProjectNormalized + + +geo_1_d :: Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> CDouble -> IO (Ptr I.GEOSGeometry)) + -> a + -> Double + -> Geos Geom +geo_1_d f g d = do + withGeos $ \h -> withGeometry g $ \gp -> do + gp' <- f h gp $ realToFrac d + constructGeometry h gp' + +-- | Return the closest point to given distance within geometry. Geometry must be a LineString +-- +interpolate :: Geometry a => a -> Double -> Geos Geom +interpolate = geo_1_d I.geos_Interpolate + +interpolateNormalized :: Geometry a => a -> Double -> Geos Geom +interpolateNormalized = geo_1_d I.geos_InterpolateNormalized + +--Binary Predicates +-------------------- +binaryPredicate_ :: Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> Ptr I.GEOSGeometry -> IO CChar) + -> String + -> a + -> a + -> Geos Bool +binaryPredicate_ f s g1 g2 = withGeos $ \h -> do + b <- throwIf ((==) 2) (mkErrorMessage s) $ + withGeometry g1 $ \gp1 -> + withGeometry g2 $ f h gp1 + return . toBool $ b + +disjoint :: Geometry a => a-> a -> Geos Bool +disjoint = binaryPredicate_ I.geos_Disjoint "disjoint" + +touches :: Geometry a => a -> a -> Geos Bool +touches = binaryPredicate_ I.geos_Touches "touches" + +intersects :: Geometry a => a -> a -> Geos Bool +intersects = binaryPredicate_ I.geos_Intersects "intersects" + +crosses :: Geometry a => a -> a -> Geos Bool +crosses = binaryPredicate_ I.geos_Crosses "crosses" + +within :: Geometry a => a -> a -> Geos Bool +within = binaryPredicate_ I.geos_Within "within" + +contains :: Geometry a => a -> a -> Geos Bool +contains = binaryPredicate_ I.geos_Contains "contains" + +overlaps :: Geometry a => a -> a -> Geos Bool +overlaps = binaryPredicate_ I.geos_Overlaps "overlaps" + +equals :: Geometry a => a -> a -> Geos Bool +equals = binaryPredicate_ I.geos_Equals "equals" + +equalsExact :: Geometry a => a -> a -> Double -> Geos Bool +equalsExact g1' g2' d = binaryPredicate_ (\h g1 g2 -> I.geos_EqualsExact h g1 g2 (realToFrac d)) "equalsExact" g1' g2' + +covers :: Geometry a => a -> a -> Geos Bool +covers = binaryPredicate_ I.geos_Covers "covers" + +coveredBy :: Geometry a => a -> a -> Geos Bool +coveredBy = binaryPredicate_ I.geos_CoveredBy "coveredBy" + +-- Misc functions + +geo_1 :: Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> Ptr CDouble -> IO CInt) + -> a + -> Geos Double +geo_1 f g = withGeos $ \h -> alloca $ \dptr -> do + _ <- throwIfZero (mkErrorMessage "geo_1" ) $ withGeometry g $ \gp -> + f h gp dptr + s <- peek dptr + pure $ realToFrac s + +area :: Geometry a => a -> Geos Double +area = geo_1 I.geos_Area + +geometryLength :: Geometry a => a -> Geos Double +geometryLength = geo_1 I.geos_Length + +geo_2_d :: Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> Ptr I.GEOSGeometry -> Ptr CDouble -> IO CInt) + -> a + -> a + -> Geos Double +geo_2_d f g p = withGeos $ \h -> alloca $ \dptr -> do + _ <- throwIfZero (mkErrorMessage "geo_2") $ withGeometry g $ \gp -> + withGeometry p $ \pp -> + f h gp pp dptr + d <- peek dptr + pure . realToFrac $ d + +distance :: Geometry a => a -> a -> Geos Double +distance = geo_2_d I.geos_Distance + +hausdorffDistance :: Geometry a => a -> a -> Geos Double +hausdorffDistance = geo_2_d I.geos_HausdorffDistance + +nearestPoints :: (Geometry a, CoordinateSequence b) => a -> a -> Geos b +nearestPoints g p = do + ptr <- withGeos $ \h -> withGeometry g $ \gp -> + withGeometry p $ I.geos_NearestPoints h gp + createCoordinateSequence ptr diff --git a/src/Data/Geometry/Geos/Raw/Internal.hs b/src/Data/Geometry/Geos/Raw/Internal.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/Internal.hs @@ -0,0 +1,613 @@ +{-# LINE 1 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} +{-# LANGUAGE CPP, ForeignFunctionInterface, EmptyDataDecls #-} +{-# LINE 2 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} +module Data.Geometry.Geos.Raw.Internal where + +import Foreign +import Foreign.C + + +{-# LINE 8 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} + +{-# LINE 9 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} + +newtype GEOSGeomType = GEOSGeomType { unGEOSGeomType :: CInt } + deriving (Eq,Show) + +pointId :: GEOSGeomType +pointId = GEOSGeomType 0 +lineStringId :: GEOSGeomType +lineStringId = GEOSGeomType 1 +polygonId :: GEOSGeomType +polygonId = GEOSGeomType 3 +multiPointId :: GEOSGeomType +multiPointId = GEOSGeomType 4 +multiLineStringId :: GEOSGeomType +multiLineStringId = GEOSGeomType 4 +multiPolygonId :: GEOSGeomType +multiPolygonId = GEOSGeomType 6 +geometryCollectionId :: GEOSGeomType +geometryCollectionId = GEOSGeomType 7 + +{-# LINE 22 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} + +newtype GEOSByteOrder = GEOSByteOrder { unGEOSByteOrder :: CInt} + deriving (Eq, Show) + +bigEndian :: GEOSByteOrder +bigEndian = GEOSByteOrder 0 +littleEndian :: GEOSByteOrder +littleEndian = GEOSByteOrder 1 + +{-# LINE 30 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} + +data GEOSContextHandle_HS +type GEOSContextHandle_t = Ptr GEOSContextHandle_HS +data GEOSGeometry +data GEOSCoordSequence +data GEOSSTRTree +data GEOSBufferParams +data GEOSPreparedGeometry +type GEOSMessageHandler = FunPtr (CString -> (Ptr ()) -> IO ()) +type GEOSMessageHandler_r = FunPtr (CString -> (Ptr ()) -> IO ()) +type GEOSQueryCallback = FunPtr (Ptr () -> Ptr () -> IO ()) + +-- read/write +data GEOSWKBWriter +data GEOSWKBReader +data GEOSWKTWriter +data GEOSWKTReader + +foreign import ccall + "stdlib.h &free" + p_free :: FunPtr (Ptr a -> IO ()) + +foreign import ccall + "notice_handlers.h init_GEOS" + geos_init :: IO GEOSContextHandle_t + +foreign import ccall + "notice_handlers.h &finalise_GEOS" + geos_finish :: FunPtr (GEOSContextHandle_t -> IO ()) + +foreign import ccall + "geos_c.h GEOSContext_setNoticeHandler_r" + geos_setNoticeHandler :: GEOSContextHandle_t -> GEOSMessageHandler -> IO GEOSMessageHandler + +foreign import ccall + "geos_c.h GEOSContext_setErrorHandler_r" + geos_setErrorHandler :: GEOSContextHandle_t -> GEOSMessageHandler -> IO GEOSMessageHandler + +-- Info + +foreign import ccall + "geos_c.h GEOSGetSRID_r" + geos_GetSRID :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + +foreign import ccall + + "geos_c.h GEOSSetSRID_r" + geos_SetSRID :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CInt -> IO () + +foreign import ccall + "geos_c.h GEOSGeom_getCoordSeq_r" + geos_GetCoordSeq :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSCoordSequence) + +foreign import ccall + "geos_c.h GEOSGeomTypeId_r" + geos_GeomTypeId :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + +foreign import ccall + "geos_c.h GEOSGeomType_r" + geos_GeomType :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr CChar) + +foreign import ccall + "geos_c.h GEOSNormalize_r" + geos_Normalize :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + + -- polygons +foreign import ccall + "geos_c.h GEOSGetInteriorRingN_r" + geos_GetInteriorRingN :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CInt -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGetExteriorRing_r" + geos_GetExteriorRing :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGetNumInteriorRings_r" + geos_GetNumInteriorRings :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt +-- multigeometries +foreign import ccall + "geos_c.h GEOSGetNumGeometries_r" + geos_GetNumGeometries :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + +foreign import ccall + "geos_c.h GEOSGetGeometryN_r" + geos_GetGeometryN :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CInt -> IO (Ptr GEOSGeometry) + + +-- Linear Referencing +foreign import ccall + "geos_c.h GEOSProject_r" + geos_Project :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CDouble +--- +foreign import ccall + "geos_c.h GEOSInterpolate_r" + geos_Interpolate :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSProjectNormalized_r" + geos_ProjectNormalized :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CDouble + +foreign import ccall + "geos_c.h GEOSInterpolateNormalized_r" + geos_InterpolateNormalized :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> IO (Ptr GEOSGeometry) + + -- Coord Sequence -- return 0 on exception + -- +foreign import ccall + "geos_c.h GEOSCoordSeq_create_r" + geos_CoordSeqCreate :: GEOSContextHandle_t -> CUInt -> CUInt -> IO (Ptr GEOSCoordSequence) + + +foreign import ccall + "geos_c.h GEOSCoordSeq_clone_r" + geos_CoordSeqClone :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO (Ptr GEOSCoordSequence ) + +foreign import ccall + "geos_c.h &GEOSCoordSeq_destroy_r" + geos_CoordSeqDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO ()) + +foreign import ccall + "geos_c.h GEOSCoordSeq_setX_r" + geos_CoordSeqSetX :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_setY_r" + geos_CoordSeqSetY :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_setZ_r" + geos_CoordSeqSetZ :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_setOrdinate_r" + geos_CoordSeqSetOrdinate :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CUInt -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getX_r" + geos_CoordSeqGetX :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getY_r" + geos_CoordSeqGetY :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getZ_r" + geos_CoordSeqGetZ :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getOrdinate_r" + geos_CoordSeqGetOrdinate :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CUInt -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getSize_r" + geos_CoordSeqGetSize :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> Ptr CUInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSGetNumCoordinates_r" + geos_GetNumCoordinates :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getDimensions_r" + geos_CoordSeqGetDimensions :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> Ptr CUInt -> IO CInt + +foreign import ccall + "geos_c.h &GEOSGeom_destroy_r" + geos_GeomDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSGeometry -> IO ()) + +foreign import ccall + "geos_c.h GEOSGeom_clone_r" + geos_GeomClone :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +-- Geometry Constructors. +-- * GEOSCoordSequence* arguments will become ownership of the returned object. +-- * All functions return NULL on exception. + +foreign import ccall + "geos_c.h GEOSGeom_createPoint_r" + geos_GeomCreatePoint :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createEmptyPoint_r" + geos_GeomCreateEmptyPoint :: GEOSContextHandle_t -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createLinearRing_r" + geos_GeomCreateLinearRing :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createLineString_r" + geos_GeomCreateLineString :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createEmptyLineString_r" + geos_GeomCreateEmptyLineString :: GEOSContextHandle_t -> IO (Ptr GEOSGeometry) + +--Second argument is an array of GEOSGeometry* objects. +--The caller remains owner of the array, but pointed-to +--objects become ownership of the returned GEOSGeometry. + +foreign import ccall + "geos_c.h GEOSGeom_createEmptyPolygon_r" + geos_GeomCreateEmptyPolygon :: GEOSContextHandle_t -> IO (Ptr GEOSGeometry) + + -- todo: this might not work as a plain pointer +foreign import ccall + "geos_c.h GEOSGeom_createPolygon_r" + geos_GeomCreatePolygon :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr (Ptr GEOSGeometry) -> CUInt -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createCollection_r" + geos_GeomCreateCollection :: GEOSContextHandle_t -> CInt -> Ptr (Ptr GEOSGeometry) -> CUInt -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createEmptyCollection_r" + geos_GeomCreateEmptyCollection :: GEOSContextHandle_t -> CInt -> IO (Ptr GEOSGeometry) + +---------- +--- Buffer +--------- + +newtype BufferCapStyle = BufferCapStyle { unBufferCapStyle :: CInt } + deriving (Eq,Show) + +capRound :: BufferCapStyle +capRound = BufferCapStyle 1 +capFlat :: BufferCapStyle +capFlat = BufferCapStyle 2 +capSquare :: BufferCapStyle +capSquare = BufferCapStyle 3 + +{-# LINE 258 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} + +newtype BufferJoinStyle = BufferJoinStyle { unBufferJoinStyle :: CInt } + deriving (Eq, Show) + +joinRound :: BufferJoinStyle +joinRound = BufferJoinStyle 1 +joinMitre :: BufferJoinStyle +joinMitre = BufferJoinStyle 2 +joinBevel :: BufferJoinStyle +joinBevel = BufferJoinStyle 3 + +{-# LINE 267 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} +foreign import ccall + "geos_c.h GEOSBuffer_r" + geos_Buffer :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> CInt -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSBufferParams_create_r" + geos_BufferParamsCreate :: GEOSContextHandle_t -> IO (Ptr GEOSBufferParams) + +foreign import ccall + "geos_c.h &GEOSBufferParams_destroy_r" + geos_BufferParamsDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSBufferParams -> IO ()) + +foreign import ccall + "geos_c.h GEOSBufferParams_setEndCapStyle_r" + geos_BufferParamsSetEndCapStyle :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferParams_setJoinStyle_r" + geos_BufferParamsSetJoinStyle :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferParams_setMitreLimit_r" + geos_BufferParamsSetMitreLimit :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferParams_setQuadrantSegments_r" + geos_BufferParamsSetQuadrantSegments :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferParams_setSingleSided_r" + geos_BufferParamsSetSingleSided :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferWithParams_r" + geos_BufferWithParams :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSBufferParams -> CDouble -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSBufferWithStyle_r" + geos_BufferWithStyle :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> CInt -> CInt -> CInt -> CDouble -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSOffsetCurve_r" + geos_OffsetCurve :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> CInt -> CInt -> CDouble -> IO (Ptr GEOSGeometry) + +----------- +--- Topology +---------- +foreign import ccall + "geos_c.h GEOSEnvelope_r" + geos_Envelope :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSIntersection_r" + geos_Intersection :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSConvexHull_r" + geos_ConvexHull :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSDifference_r" + geos_Difference :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSSymDifference_r" + geos_SymDifference :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSBoundary_r" + geos_Boundary :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSUnion_r" + geos_Union :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSUnaryUnion_r" + geos_UnaryUnion :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGetCentroid_r" + geos_GetCentroid :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSPointOnSurface_r" + geos_PointsOnSurface :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSNode_r" + geos_Node :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + + +foreign import ccall + "geos_c.h GEOSDelaunayTriangulation_r" + geos_DelaunayTriangulation :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> CInt -> IO (Ptr GEOSGeometry) + + + +{-# LINE 365 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} +foreign import ccall + "geos_c.h GEOSVoronoiDiagram_r" + geos_VoronoiDiagram :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> CDouble -> CInt -> IO (Ptr GEOSGeometry) + +{-# LINE 369 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} + +----- +--Binary Predicates. +----- + +foreign import ccall + "geos_c.h GEOSDisjoint_r" + geos_Disjoint :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSTouches_r" + geos_Touches :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSCrosses_r" + geos_Crosses :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSIntersects_r" + geos_Intersects :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSWithin_r" + geos_Within :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSContains_r" + geos_Contains :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSOverlaps_r" + geos_Overlaps :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSEquals_r" + geos_Equals :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSEqualsExact_r" + geos_EqualsExact :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> CDouble -> IO CChar + +foreign import ccall + "geos_c.h GEOSCovers_r" + geos_Covers :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSCoveredBy_r" + geos_CoveredBy :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +--- prepared Geometries + +--- Readers / Writers +foreign import ccall + "geos_c.h GEOSWKBReader_create_r" + geos_WKBReaderCreate :: GEOSContextHandle_t -> IO (Ptr GEOSWKBReader) + +foreign import ccall + "geos_c.h &GEOSWKBReader_destroy_r" + geos_WKBReaderDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSWKBReader -> IO ()) + +foreign import ccall + "geos_c.h GEOSWKBReader_read_r" + geos_WKBReaderRead :: GEOSContextHandle_t -> Ptr GEOSWKBReader -> CString -> CSize -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSWKBReader_readHEX_r" + geos_WKBReaderReadHex :: GEOSContextHandle_t -> Ptr GEOSWKBReader -> CString -> CSize -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSWKBWriter_create_r" + geos_WKBWriterCreate :: GEOSContextHandle_t -> IO (Ptr GEOSWKBWriter) + +foreign import ccall + "geos_c.h &GEOSWKBWriter_destroy_r" + geos_WKBWriterDestroy :: FunPtr ( GEOSContextHandle_t -> Ptr GEOSWKBWriter -> IO ()) + + -- caller owns results +foreign import ccall + "geos_c.h GEOSWKBWriter_write_r" + geos_WKBWriterWrite :: GEOSContextHandle_t -> Ptr GEOSWKBWriter -> Ptr GEOSGeometry -> Ptr CSize -> IO CString + + -- caller owns results +foreign import ccall + "geos_c.h GEOSWKBWriter_writeHEX_r" + geos_WKBWriterWriteHex :: GEOSContextHandle_t -> Ptr GEOSWKBWriter -> Ptr GEOSGeometry -> Ptr CSize -> IO CString + +foreign import ccall + "geos_c.h GEOSWKBWriter_setIncludeSRID_r" + geos_WKBWriterSetIncludeSRID :: GEOSContextHandle_t -> Ptr GEOSWKBWriter -> CChar -> IO () + +foreign import ccall + "geos_c.h GEOSWKTReader_create_r" + geos_WKTReaderCreate :: GEOSContextHandle_t -> IO (Ptr GEOSWKTReader) + +foreign import ccall + "geos_c.h &GEOSWKTReader_destroy_r" + geos_WKTReaderDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSWKTReader -> IO ()) + +foreign import ccall + "geos_c.h GEOSWKTReader_read_r" + geos_WKTReaderRead :: GEOSContextHandle_t -> Ptr GEOSWKTReader -> CString -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSWKTWriter_create_r" + geos_WKTWriterCreate :: GEOSContextHandle_t -> IO (Ptr GEOSWKTWriter) + +foreign import ccall + "geos_c.h &GEOSWKTWriter_destroy_r" + geos_WKTWriterDestroy :: FunPtr ( GEOSContextHandle_t -> Ptr GEOSWKTWriter -> IO ()) + + -- caller owns results +foreign import ccall + "geos_c.h GEOSWKTWriter_write_r" + geos_WKTWriterWrite :: GEOSContextHandle_t -> Ptr GEOSWKTWriter -> Ptr GEOSGeometry -> IO CString + + +--TODO: 1085 finish writer methods + +-- following return 0 on exception, 1 otherwise +foreign import ccall + "geos_c.h GEOSArea_r" + geos_Area :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSLength_r" + geos_Length :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSDistance_r" + geos_Distance :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSHausdorffDistance_r" + geos_HausdorffDistance :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSHausdorffDistanceDensify_r" + geos_HausdorffDistanceDensify :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> CDouble -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSGeomGetLength_r" + geos_GeomGetLength :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSNearestPoints_r" + geos_NearestPoints :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSCoordSequence) + +-- | Prepared Geometries +foreign import ccall + "geos_c.h GEOSPrepare_r" + geos_Prepare :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSPreparedGeometry) + +foreign import ccall + "geos_c.h &GEOSPreparedGeom_destroy_r" + geos_PreparedGeomDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> IO ()) + +foreign import ccall + "geos_c.h GEOSPreparedContains_r" + geos_PreparedContains :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedContainsProperly_r" + geos_PreparedContainsProperly :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedCoveredBy_r" + geos_PreparedCoveredBy :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedCovers_r" + geos_PreparedCovers :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedCrosses_r" + geos_PreparedCrosses :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedDisjoint_r" + geos_PreparedDisjoint :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedIntersects_r" + geos_PreparedIntersects :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedOverlaps_r" + geos_PreparedOverlaps :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedTouches_r" + geos_PreparedTouches :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedWithin_r" + geos_PreparedWithin :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +-- STRTREE + + +{-# LINE 568 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} +foreign import ccall + "geos_c.h GEOSSTRtree_create_r" + geos_STRTreeCreate :: GEOSContextHandle_t -> CSize -> IO (Ptr GEOSSTRTree) + + +foreign import ccall + "geos_c.h GEOSSTRtree_insert_r" + geos_STRTreeInsert :: GEOSContextHandle_t -> Ptr GEOSSTRTree -> Ptr GEOSGeometry -> Ptr () -> IO () + +foreign import ccall + "geos_c.h GEOSSTRtree_query_r" + geos_STRTreeQuery :: GEOSContextHandle_t -> Ptr GEOSSTRTree -> Ptr GEOSGeometry -> GEOSQueryCallback -> Ptr () -> IO () + +foreign import ccall + "geos_c.h GEOSSTRtree_iterate_r" + geos_STRTreeIterate :: GEOSContextHandle_t -> Ptr GEOSSTRTree -> GEOSQueryCallback -> Ptr () -> IO () + +foreign import ccall + "geos_c.h &GEOSSTRtree_destroy_r" + geos_STRTreeDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSSTRTree -> IO ()) + + +{-# LINE 590 "src/Data/Geometry/Geos/Raw/Internal.hsc" #-} + diff --git a/src/Data/Geometry/Geos/Raw/Internal.hsc b/src/Data/Geometry/Geos/Raw/Internal.hsc new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/Internal.hsc @@ -0,0 +1,591 @@ +{-# LANGUAGE CPP, ForeignFunctionInterface, EmptyDataDecls #-} + +module Data.Geometry.Geos.Raw.Internal where + +import Foreign +import Foreign.C + +#include +#include "noticehandlers.h" + +newtype GEOSGeomType = GEOSGeomType { unGEOSGeomType :: CInt } + deriving (Eq,Show) + +#{ enum GEOSGeomType, GEOSGeomType, + pointId = GEOS_POINT + , lineStringId = GEOS_LINESTRING + , polygonId = GEOS_POLYGON + , multiPointId = GEOS_MULTIPOINT + , multiLineStringId = GEOS_MULTILINESTRING + , multiPolygonId = GEOS_MULTIPOLYGON + , geometryCollectionId = GEOS_GEOMETRYCOLLECTION +} + +newtype GEOSByteOrder = GEOSByteOrder { unGEOSByteOrder :: CInt} + deriving (Eq, Show) + +#{ enum GEOSByteOrder, GEOSByteOrder, + bigEndian = 0 + , littleEndian = 1 +} + +data GEOSContextHandle_HS +type GEOSContextHandle_t = Ptr GEOSContextHandle_HS +data GEOSGeometry +data GEOSCoordSequence +data GEOSSTRTree +data GEOSBufferParams +data GEOSPreparedGeometry +type GEOSMessageHandler = FunPtr (CString -> (Ptr ()) -> IO ()) +type GEOSMessageHandler_r = FunPtr (CString -> (Ptr ()) -> IO ()) +type GEOSQueryCallback = FunPtr (Ptr () -> Ptr () -> IO ()) + +-- read/write +data GEOSWKBWriter +data GEOSWKBReader +data GEOSWKTWriter +data GEOSWKTReader + + +foreign import ccall + "notice_handlers.h init_GEOS" + geos_init :: IO GEOSContextHandle_t + +foreign import ccall + "geos_c.h &finishGEOS_r" + geos_finish :: FunPtr (GEOSContextHandle_t -> IO ()) + +foreign import ccall + "geos_c.h GEOSContext_setNoticeHandler_r" + geos_setNoticeHandler :: GEOSContextHandle_t -> GEOSMessageHandler -> IO GEOSMessageHandler + +foreign import ccall + "geos_c.h GEOSContext_setErrorHandler_r" + geos_setErrorHandler :: GEOSContextHandle_t -> GEOSMessageHandler -> IO GEOSMessageHandler + +-- Info + +foreign import ccall + "geos_c.h GEOSGetSRID_r" + geos_GetSRID :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + +foreign import ccall + + "geos_c.h GEOSSetSRID_r" + geos_SetSRID :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CInt -> IO () + +foreign import ccall + "geos_c.h GEOSGeom_getCoordSeq_r" + geos_GetCoordSeq :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSCoordSequence) + +foreign import ccall + "geos_c.h GEOSGeomTypeId_r" + geos_GeomTypeId :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + +foreign import ccall + "geos_c.h GEOSGeomType_r" + geos_GeomType :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr CChar) + +foreign import ccall + "geos_c.h GEOSNormalize_r" + geos_Normalize :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + + -- polygons +foreign import ccall + "geos_c.h GEOSGetInteriorRingN_r" + geos_GetInteriorRingN :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CInt -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGetExteriorRing_r" + geos_GetExteriorRing :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGetNumInteriorRings_r" + geos_GetNumInteriorRings :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt +-- multigeometries +foreign import ccall + "geos_c.h GEOSGetNumGeometries_r" + geos_GetNumGeometries :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + +foreign import ccall + "geos_c.h GEOSGetGeometryN_r" + geos_GetGeometryN :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CInt -> IO (Ptr GEOSGeometry) + + +-- Linear Referencing +foreign import ccall + "geos_c.h GEOSProject_r" + geos_Project :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CDouble +--- +foreign import ccall + "geos_c.h GEOSInterpolate_r" + geos_Interpolate :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSProjectNormalized_r" + geos_ProjectNormalized :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CDouble + +foreign import ccall + "geos_c.h GEOSInterpolateNormalized_r" + geos_InterpolateNormalized :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> IO (Ptr GEOSGeometry) + + -- Coord Sequence -- return 0 on exception + -- +foreign import ccall + "geos_c.h GEOSCoordSeq_create_r" + geos_CoordSeqCreate :: GEOSContextHandle_t -> CUInt -> CUInt -> IO (Ptr GEOSCoordSequence) + + +foreign import ccall + "geos_c.h GEOSCoordSeq_clone_r" + geos_CoordSeqClone :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO (Ptr GEOSCoordSequence ) + +foreign import ccall + "geos_c.h &GEOSCoordSeq_destroy_r" + geos_CoordSeqDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO ()) + +foreign import ccall + "geos_c.h GEOSCoordSeq_setX_r" + geos_CoordSeqSetX :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_setY_r" + geos_CoordSeqSetY :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_setZ_r" + geos_CoordSeqSetZ :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_setOrdinate_r" + geos_CoordSeqSetOrdinate :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CUInt -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getX_r" + geos_CoordSeqGetX :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getY_r" + geos_CoordSeqGetY :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getZ_r" + geos_CoordSeqGetZ :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getOrdinate_r" + geos_CoordSeqGetOrdinate :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> CUInt -> CUInt -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getSize_r" + geos_CoordSeqGetSize :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> Ptr CUInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSGetNumCoordinates_r" + geos_GetNumCoordinates :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CInt + +foreign import ccall + "geos_c.h GEOSCoordSeq_getDimensions_r" + geos_CoordSeqGetDimensions :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> Ptr CUInt -> IO CInt + +foreign import ccall + "geos_c.h &GEOSGeom_destroy_r" + geos_GeomDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSGeometry -> IO ()) + +foreign import ccall + "geos_c.h GEOSGeom_clone_r" + geos_GeomClone :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +-- Geometry Constructors. +-- * GEOSCoordSequence* arguments will become ownership of the returned object. +-- * All functions return NULL on exception. + +foreign import ccall + "geos_c.h GEOSGeom_createPoint_r" + geos_GeomCreatePoint :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createEmptyPoint_r" + geos_GeomCreateEmptyPoint :: GEOSContextHandle_t -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createLinearRing_r" + geos_GeomCreateLinearRing :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createLineString_r" + geos_GeomCreateLineString :: GEOSContextHandle_t -> Ptr GEOSCoordSequence -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createEmptyLineString_r" + geos_GeomCreateEmptyLineString :: GEOSContextHandle_t -> IO (Ptr GEOSGeometry) + +--Second argument is an array of GEOSGeometry* objects. +--The caller remains owner of the array, but pointed-to +--objects become ownership of the returned GEOSGeometry. + +foreign import ccall + "geos_c.h GEOSGeom_createEmptyPolygon_r" + geos_GeomCreateEmptyPolygon :: GEOSContextHandle_t -> IO (Ptr GEOSGeometry) + + -- todo: this might not work as a plain pointer +foreign import ccall + "geos_c.h GEOSGeom_createPolygon_r" + geos_GeomCreatePolygon :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr (Ptr GEOSGeometry) -> CUInt -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createCollection_r" + geos_GeomCreateCollection :: GEOSContextHandle_t -> CInt -> Ptr (Ptr GEOSGeometry) -> CUInt -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGeom_createEmptyCollection_r" + geos_GeomCreateEmptyCollection :: GEOSContextHandle_t -> CInt -> IO (Ptr GEOSGeometry) + +---------- +--- Buffer +--------- + +newtype BufferCapStyle = BufferCapStyle { unBufferCapStyle :: CInt } + deriving (Eq,Show) + +#{ enum BufferCapStyle, BufferCapStyle, + capRound = 1 + , capFlat = 2 + , capSquare = 3 +} + +newtype BufferJoinStyle = BufferJoinStyle { unBufferJoinStyle :: CInt } + deriving (Eq, Show) + +#{ enum BufferJoinStyle, BufferJoinStyle, + joinRound = 1 + , joinMitre = 2 + , joinBevel = 3 +} +foreign import ccall + "geos_c.h GEOSBuffer_r" + geos_Buffer :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> CInt -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSBufferParams_create_r" + geos_BufferParamsCreate :: GEOSContextHandle_t -> IO (Ptr GEOSBufferParams) + +foreign import ccall + "geos_c.h &GEOSBufferParams_destroy_r" + geos_BufferParamsDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSBufferParams -> IO ()) + +foreign import ccall + "geos_c.h GEOSBufferParams_setEndCapStyle_r" + geos_BufferParamsSetEndCapStyle :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferParams_setJoinStyle_r" + geos_BufferParamsSetJoinStyle :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferParams_setMitreLimit_r" + geos_BufferParamsSetMitreLimit :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferParams_setQuadrantSegments_r" + geos_BufferParamsSetQuadrantSegments :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferParams_setSingleSided_r" + geos_BufferParamsSetSingleSided :: GEOSContextHandle_t -> Ptr GEOSBufferParams -> CInt -> IO CInt + +foreign import ccall + "geos_c.h GEOSBufferWithParams_r" + geos_BufferWithParams :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSBufferParams -> CDouble -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSBufferWithStyle_r" + geos_BufferWithStyle :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> CInt -> CInt -> CInt -> CDouble -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSOffsetCurve_r" + geos_OffsetCurve :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> CInt -> CInt -> CDouble -> IO (Ptr GEOSGeometry) + +----------- +--- Topology +---------- +foreign import ccall + "geos_c.h GEOSEnvelope_r" + geos_Envelope :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSIntersection_r" + geos_Intersection :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSConvexHull_r" + geos_ConvexHull :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSDifference_r" + geos_Difference :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSSymDifference_r" + geos_SymDifference :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSBoundary_r" + geos_Boundary :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSUnion_r" + geos_Union :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSUnaryUnion_r" + geos_UnaryUnion :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSGetCentroid_r" + geos_GetCentroid :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSPointOnSurface_r" + geos_PointsOnSurface :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSNode_r" + geos_Node :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSGeometry) + + +foreign import ccall + "geos_c.h GEOSDelaunayTriangulation_r" + geos_DelaunayTriangulation :: GEOSContextHandle_t -> Ptr GEOSGeometry -> CDouble -> CInt -> IO (Ptr GEOSGeometry) + + +#if GEOS_VERSION_MAJOR >= 3 && GEOS_VERSION_MINOR > 4 +foreign import ccall + "geos_c.h GEOSVoronoiDiagram_r" + geos_VoronoiDiagram :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> CDouble -> CInt -> IO (Ptr GEOSGeometry) +#endif + +----- +--Binary Predicates. +----- + +foreign import ccall + "geos_c.h GEOSDisjoint_r" + geos_Disjoint :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSTouches_r" + geos_Touches :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSCrosses_r" + geos_Crosses :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSIntersects_r" + geos_Intersects :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSWithin_r" + geos_Within :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSContains_r" + geos_Contains :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSOverlaps_r" + geos_Overlaps :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSEquals_r" + geos_Equals :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSEqualsExact_r" + geos_EqualsExact :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> CDouble -> IO CChar + +foreign import ccall + "geos_c.h GEOSCovers_r" + geos_Covers :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSCoveredBy_r" + geos_CoveredBy :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO CChar + +--- prepared Geometries + +--- Readers / Writers +foreign import ccall + "geos_c.h GEOSWKBReader_create_r" + geos_WKBReaderCreate :: GEOSContextHandle_t -> IO (Ptr GEOSWKBReader) + +foreign import ccall + "geos_c.h &GEOSWKBReader_destroy_r" + geos_WKBReaderDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSWKBReader -> IO ()) + +foreign import ccall + "geos_c.h GEOSWKBReader_read_r" + geos_WKBReaderRead :: GEOSContextHandle_t -> Ptr GEOSWKBReader -> CString -> CSize -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSWKBReader_readHEX_r" + geos_WKBReaderReadHex :: GEOSContextHandle_t -> Ptr GEOSWKBReader -> CString -> CSize -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSWKBWriter_create_r" + geos_WKBWriterCreate :: GEOSContextHandle_t -> IO (Ptr GEOSWKBWriter) + +foreign import ccall + "geos_c.h &GEOSWKBWriter_destroy_r" + geos_WKBWriterDestroy :: FunPtr ( GEOSContextHandle_t -> Ptr GEOSWKBWriter -> IO ()) + + -- caller owns results +foreign import ccall + "geos_c.h GEOSWKBWriter_write_r" + geos_WKBWriterWrite :: GEOSContextHandle_t -> Ptr GEOSWKBWriter -> Ptr GEOSGeometry -> Ptr CSize -> IO CString + + -- caller owns results +foreign import ccall + "geos_c.h GEOSWKBWriter_writeHEX_r" + geos_WKBWriterWriteHex :: GEOSContextHandle_t -> Ptr GEOSWKBWriter -> Ptr GEOSGeometry -> Ptr CSize -> IO CString + +foreign import ccall + "geos_c.h GEOSWKBWriter_setIncludeSRID_r" + geos_WKBWriterSetIncludeSRID :: GEOSContextHandle_t -> Ptr GEOSWKBWriter -> CChar -> IO () + +foreign import ccall + "geos_c.h GEOSWKTReader_create_r" + geos_WKTReaderCreate :: GEOSContextHandle_t -> IO (Ptr GEOSWKTReader) + +foreign import ccall + "geos_c.h &GEOSWKTReader_destroy_r" + geos_WKTReaderDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSWKTReader -> IO ()) + +foreign import ccall + "geos_c.h GEOSWKTReader_read_r" + geos_WKTReaderRead :: GEOSContextHandle_t -> Ptr GEOSWKTReader -> CString -> IO (Ptr GEOSGeometry) + +foreign import ccall + "geos_c.h GEOSWKTWriter_create_r" + geos_WKTWriterCreate :: GEOSContextHandle_t -> IO (Ptr GEOSWKTWriter) + +foreign import ccall + "geos_c.h &GEOSWKTWriter_destroy_r" + geos_WKTWriterDestroy :: FunPtr ( GEOSContextHandle_t -> Ptr GEOSWKTWriter -> IO ()) + + -- caller owns results +foreign import ccall + "geos_c.h GEOSWKTWriter_write_r" + geos_WKTWriterWrite :: GEOSContextHandle_t -> Ptr GEOSWKTWriter -> Ptr GEOSGeometry -> IO CString + + +--TODO: 1085 finish writer methods + +-- following return 0 on exception, 1 otherwise +foreign import ccall + "geos_c.h GEOSArea_r" + geos_Area :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSLength_r" + geos_Length :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSDistance_r" + geos_Distance :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSHausdorffDistance_r" + geos_HausdorffDistance :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSHausdorffDistanceDensify_r" + geos_HausdorffDistanceDensify :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> CDouble -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSGeomGetLength_r" + geos_GeomGetLength :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr CDouble -> IO CInt + +foreign import ccall + "geos_c.h GEOSNearestPoints_r" + geos_NearestPoints :: GEOSContextHandle_t -> Ptr GEOSGeometry -> Ptr GEOSGeometry -> IO (Ptr GEOSCoordSequence) + +-- | Prepared Geometries +foreign import ccall + "geos_c.h GEOSPrepare_r" + geos_Prepare :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO (Ptr GEOSPreparedGeometry) + +foreign import ccall + "geos_c.h &GEOSPreparedGeom_destroy_r" + geos_PreparedGeomDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> IO ()) + +foreign import ccall + "geos_c.h GEOSPreparedContains_r" + geos_PreparedContains :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedContainsProperly_r" + geos_PreparedContainsProperly :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedCoveredBy_r" + geos_PreparedCoveredBy :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedCovers_r" + geos_PreparedCovers :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedCrosses_r" + geos_PreparedCrosses :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedDisjoint_r" + geos_PreparedDisjoint :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedIntersects_r" + geos_PreparedIntersects :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedOverlaps_r" + geos_PreparedOverlaps :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedTouches_r" + geos_PreparedTouches :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +foreign import ccall + "geos_c.h GEOSPreparedWithin_r" + geos_PreparedWithin :: GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar + +-- STRTREE + +#if GEOS_VERSION_MAJOR >= 3 && GEOS_VERSION_MINOR > 4 +foreign import ccall + "geos_c.h GEOSSTRtree_create_r" + geos_STRTreeCreate :: GEOSContextHandle_t -> CSize -> IO (Ptr GEOSSTRTree) + + +foreign import ccall + "geos_c.h GEOSSTRtree_insert_r" + geos_STRTreeInsert :: GEOSContextHandle_t -> Ptr GEOSSTRTree -> Ptr GEOSGeometry -> Ptr () -> IO () + +foreign import ccall + "geos_c.h GEOSSTRtree_query_r" + geos_STRTreeQuery :: GEOSContextHandle_t -> Ptr GEOSSTRTree -> Ptr GEOSGeometry -> GEOSQueryCallback -> Ptr () -> IO () + +foreign import ccall + "geos_c.h GEOSSTRtree_iterate_r" + geos_STRTreeIterate :: GEOSContextHandle_t -> Ptr GEOSSTRTree -> GEOSQueryCallback -> Ptr () -> IO () + +foreign import ccall + "geos_c.h &GEOSSTRtree_destroy_r" + geos_STRTreeDestroy :: FunPtr (GEOSContextHandle_t -> Ptr GEOSSTRTree -> IO ()) + +#endif + +foreign import ccall + "geos_c.h GEOSisValid_r" + geos_isValid :: GEOSContextHandle_t -> Ptr GEOSGeometry -> IO CChar diff --git a/src/Data/Geometry/Geos/Raw/Prepared.hs b/src/Data/Geometry/Geos/Raw/Prepared.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/Prepared.hs @@ -0,0 +1,79 @@ +module Data.Geometry.Geos.Raw.Prepared ( + PreparedGeometry (..) + , withPreparedGeometry + , prepare + , contains + , containsProperly + , coveredBy + , covers + , crosses + , disjoint + , intersects + , overlaps + , touches + , within +) where +import Data.Geometry.Geos.Raw.Internal +import Data.Geometry.Geos.Raw.Base +import qualified Data.Geometry.Geos.Raw.Geometry as RG +import Foreign.Marshal.Utils +import Foreign +import Foreign.C.Types + + +newtype PreparedGeometry = PreparedGeometry { + unPreparedGeometry :: ForeignPtr GEOSPreparedGeometry +} deriving (Show, Eq) + +withPreparedGeometry :: PreparedGeometry -> (Ptr GEOSPreparedGeometry -> IO a ) -> IO a +withPreparedGeometry (PreparedGeometry g) f = withForeignPtr g f + +prepare :: RG.Geometry a => a -> Geos PreparedGeometry +prepare g = withGeos $ \h -> do + g' <- RG.withGeometry g $ \gp -> + geos_Prepare h gp + fp <- newForeignPtrEnv geos_PreparedGeomDestroy h g' + return $ PreparedGeometry fp + +queryPrepared :: RG.Geometry a + => (GEOSContextHandle_t -> Ptr GEOSPreparedGeometry -> Ptr GEOSGeometry -> IO CChar) + -> String + -> PreparedGeometry + -> a + -> Geos Bool +queryPrepared f s pg g = withGeos $ \h -> do + b <- throwIf (\v -> v == 2) (mkErrorMessage s) $ + RG.withGeometry g $ \gp -> + withPreparedGeometry pg $ \pp -> + f h pp gp + return . toBool $ b + +contains :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +contains = queryPrepared geos_PreparedContains "contains" + +containsProperly :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +containsProperly = queryPrepared geos_PreparedContainsProperly "containsProperly" + +coveredBy :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +coveredBy = queryPrepared geos_PreparedCoveredBy "coveredBy" + +covers :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +covers = queryPrepared geos_PreparedCovers "covers" + +crosses :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +crosses = queryPrepared geos_PreparedCrosses "crosses" + +disjoint :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +disjoint = queryPrepared geos_PreparedDisjoint "disjoint" + +intersects :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +intersects = queryPrepared geos_PreparedIntersects "intersects" + +overlaps :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +overlaps = queryPrepared geos_PreparedOverlaps "overlaps" + +touches :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +touches = queryPrepared geos_PreparedTouches "touches" + +within :: RG.Geometry a => PreparedGeometry -> a-> Geos Bool +within = queryPrepared geos_PreparedWithin "within" diff --git a/src/Data/Geometry/Geos/Raw/Serialize.hs b/src/Data/Geometry/Geos/Raw/Serialize.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/Serialize.hs @@ -0,0 +1,119 @@ +{-# LANGUAGE ScopedTypeVariables #-} +module Data.Geometry.Geos.Raw.Serialize ( + createReader + , createWriter + , createWktReader + , createWktWriter + , read + , readHex + , readWkt + , write + , writeHex + , writeWkt +) where +import Prelude hiding (read) +import qualified Data.Geometry.Geos.Raw.Internal as I +import Data.Geometry.Geos.Raw.Base +import Data.Geometry.Geos.Raw.Geometry +import Foreign +import Foreign.C.String +import Foreign.C.Types +import qualified Data.ByteString.Char8 as BC +import Control.Exception (onException) + + +newtype Reader = Reader { _unReader :: ForeignPtr I.GEOSWKBReader } +newtype Writer = Writer { _unWriter :: ForeignPtr I.GEOSWKBWriter } + +newtype WktReader = WktReader { _unWktReader :: ForeignPtr I.GEOSWKTReader } +newtype WktWriter = WktWriter { _unWktWriter :: ForeignPtr I.GEOSWKTWriter } + +createReader :: Geos Reader +createReader = withGeos $ \h -> do + ptr <- throwIfNull "Create Reader" $ I.geos_WKBReaderCreate h + Reader <$> newForeignPtrEnv I.geos_WKBReaderDestroy h ptr + +createWktReader :: Geos WktReader +createWktReader = withGeos $ \h -> do + ptr <- throwIfNull "Create WKT Reader" $ I.geos_WKTReaderCreate h + WktReader <$> newForeignPtrEnv I.geos_WKTReaderDestroy h ptr + +read_ :: (I.GEOSContextHandle_t -> Ptr I.GEOSWKBReader -> CString -> CSize -> IO (Ptr I.GEOSGeometry)) + -> Reader + -> BC.ByteString + -> Geos (Maybe Geom) +read_ f (Reader r) bs = withGeos $ \h -> + onException (readBlock h) (pure Nothing) + where + readBlock h = do + ptr <- withForeignPtr r $ \rp -> + BC.useAsCStringLen bs $ + \(cs, l) -> f h rp cs $ fromIntegral l + g <- wrapUpGeom h ptr + pure g + +wrapUpGeom :: I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> IO (Maybe Geom) +wrapUpGeom h ptr + | ptr == nullPtr = pure Nothing + | otherwise = Just . Geom <$> newForeignPtrEnv I.geos_GeomDestroy h ptr + +read :: Reader -> BC.ByteString -> Geos (Maybe Geom) +read = readWkb + +readWkb :: Reader -> BC.ByteString -> Geos (Maybe Geom) +readWkb = read_ I.geos_WKBReaderRead + +readHex :: Reader -> BC.ByteString -> Geos (Maybe Geom) +readHex = read_ I.geos_WKBReaderReadHex + +readWkt :: WktReader -> BC.ByteString -> Geos (Maybe Geom) +readWkt (WktReader r) bs = do + withGeos $ \h -> do + ptr <- withForeignPtr r $ \rp -> + BC.useAsCString bs $ \cs -> + I.geos_WKTReaderRead h rp cs + g <- wrapUpGeom h ptr + pure g + +createWriter :: Geos Writer +createWriter = withGeos $ \h -> do + ptr <- throwIfNull "CreateWriter" $ I.geos_WKBWriterCreate h + I.geos_WKBWriterSetIncludeSRID h ptr $ fromBool True + fp <- newForeignPtrEnv I.geos_WKBWriterDestroy h ptr + return $ Writer fp + +createWktWriter :: Geos WktWriter +createWktWriter = withGeos $ \h -> do + ptr <- throwIfNull "CreateWktWriter" $ I.geos_WKTWriterCreate h + fp <- newForeignPtrEnv I.geos_WKTWriterDestroy h ptr + return $ WktWriter fp + +write_ :: Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSWKBWriter -> Ptr I.GEOSGeometry -> Ptr CSize -> IO CString) + -> Writer + -> a + -> Geos BC.ByteString +write_ f (Writer w) g = withGeos $ \h -> do + clen <- withForeignPtr w $ \wp -> + withGeometry g $ \gp -> + alloca $ \lp -> do + cs <- f h wp gp lp + vl <- fromIntegral `fmap` (peek lp) + return (cs, vl) + bs <- BC.packCStringLen clen + return bs + +write :: Geometry a => Writer -> a -> Geos BC.ByteString +write = write_ I.geos_WKBWriterWrite + +writeHex :: Geometry a => Writer -> a -> Geos BC.ByteString +writeHex = write_ I.geos_WKBWriterWriteHex + +writeWkt :: Geometry a => WktWriter -> a -> Geos BC.ByteString +writeWkt (WktWriter w) g = withGeos $ \h -> do + wkt <- withForeignPtr w $ \wp -> + withGeometry g $ \gp -> do + cs <- I.geos_WKTWriterWrite h wp gp + return cs + bs <- BC.packCString wkt + return bs diff --git a/src/Data/Geometry/Geos/Raw/Topology.hs b/src/Data/Geometry/Geos/Raw/Topology.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Raw/Topology.hs @@ -0,0 +1,97 @@ +{-# LANGUAGE CPP #-} +module Data.Geometry.Geos.Raw.Topology ( + envelope + , intersection + , convexHull + , difference + , symmetricDifference + , boundary + , union + , unaryUnion + , pointOnSurface + , centroid + , node + , delaunayTriangulation +) where +import qualified Data.Geometry.Geos.Raw.Internal as I +import Data.Geometry.Geos.Raw.Base +import qualified Data.Geometry.Geos.Raw.Geometry as R +import Foreign + +geo_1 :: R.Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> IO (Ptr I.GEOSGeometry)) + -> String + -> a + -> Geos a +geo_1 f s g = + withGeos $ \h -> + R.withGeometry g $ \gp -> do + ptr <- throwIfNull s $ f h gp + R.constructGeometry h ptr + +geo_2 :: R.Geometry a + => (I.GEOSContextHandle_t -> Ptr I.GEOSGeometry -> Ptr I.GEOSGeometry -> IO (Ptr I.GEOSGeometry)) + -> String + -> a + -> a + -> Geos a +geo_2 f s g1 g2 = do + withGeos $ \h -> + R.withGeometry g1 $ \gp -> + R.withGeometry g2 $ \gp2 -> do + ptr <- throwIfNull s $ f h gp gp2 + R.constructGeometry h ptr + + +envelope :: R.Geometry a => a -> Geos a +envelope = geo_1 I.geos_Envelope "envelope" + +intersection :: R.Geometry a => a -> a -> Geos a +intersection = geo_2 I.geos_Intersection "intersection" + +convexHull :: R.Geometry a => a -> Geos a +convexHull = geo_1 I.geos_ConvexHull "convexHull" + +difference :: R.Geometry a => a -> a -> Geos a +difference = geo_2 I.geos_Difference "difference" + +symmetricDifference :: R.Geometry a => a -> a-> Geos a +symmetricDifference = geo_2 I.geos_SymDifference "symmetricDifference" + +boundary :: R.Geometry a => a -> Geos a +boundary = geo_1 I.geos_Boundary "boundary" + +union :: R.Geometry a => a -> a -> Geos a +union = geo_2 I.geos_Union "union" + +unaryUnion :: R.Geometry a => a -> Geos a +unaryUnion = geo_1 I.geos_UnaryUnion "unaryUnion" + +pointOnSurface :: R.Geometry a => a -> Geos a +pointOnSurface = geo_1 I.geos_PointsOnSurface "pointOnSurface" + +centroid :: R.Geometry a => a -> Geos a +centroid = geo_1 I.geos_GetCentroid "getCentroid" + +node :: R.Geometry a => a -> Geos a +node = geo_1 I.geos_Node "node" + +-- | Return a Delaunay triangulation of the vertex of the given geometry @g@, where @tol@ is the snapping tolerance to use. +delaunayTriangulation :: R.Geometry a => a -> Double -> Geos a +delaunayTriangulation g tol = do + withGeos $ \h -> do + R.withGeometry g $ \gp -> do + ptr <- throwIfNull "delaunayTriangulation" $ I.geos_DelaunayTriangulation h gp (realToFrac tol) $ fromBool True + R.constructGeometry h ptr + +#if GEOS_VERSION_MAJOR > 3 && GEOS_VERSION_MINOR > 4 +-- | +-- TODO: make env Maybe Geometry +voronoiDiagram :: R.Geometry a => a -> Maybe a -> Double -> Bool -> Geos a +voronoiDiagram g menv tol oe = do + withGeos $ \h -> + R.withGeometry g $ \gp -> + R.withMaybeGeometry env $ \ep -> + ptr <- throwIfNull "voronoiDiagram" $ I.geos_VoronoiDiagram hp gp ep (realToFrac tol) $ fromBool oe + R.constructGeometry h ptr +#endif diff --git a/src/Data/Geometry/Geos/Serialize.hs b/src/Data/Geometry/Geos/Serialize.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Serialize.hs @@ -0,0 +1,54 @@ +{-# LANGUAGE OverloadedStrings #-} +{-# LANGUAGE ScopedTypeVariables #-} + +module Data.Geometry.Geos.Serialize ( + readHex + , readLotsOfHex + , writeHex + , readWkt + , writeWkt +) where + +import Data.Geometry.Geos.Raw.Base +import Data.Geometry.Geos.Geometry +import qualified Data.Geometry.Geos.Raw.Serialize as S +import qualified Data.Geometry.Geos.Raw.Geometry as R +import Data.Geometry.Geos.Types +import qualified Data.ByteString.Char8 as BC +import Data.Maybe (catMaybes) + +readHex :: BC.ByteString -> Maybe (Some Geometry) +readHex bs = runGeos $ do + r <- S.createReader + g <- S.readHex r bs + case g of + Just g' -> Just <$> convertGeometryFromRaw g' + Nothing -> return Nothing + +readLotsOfHex :: [BC.ByteString] -> [Some Geometry] +readLotsOfHex bs = runGeos $ do + r <- S.createReader + x <- traverse (S.readHex r) bs + traverse convertGeometryFromRaw $ catMaybes x + +writeHex :: Geometry a -> BC.ByteString +writeHex g = runGeos $ do + w <- S.createWriter + r :: R.Geom <- convertGeometryToRaw g + S.writeHex w r + +readWkt :: Maybe Int -> BC.ByteString -> Maybe (Some Geometry) +readWkt srid bs = runGeos $ do + r <- S.createWktReader + g <- S.readWkt r bs + case g of + Just g' -> do + g'' <- R.setSRID srid g' + Just <$> convertGeometryFromRaw g'' + Nothing -> return Nothing + +writeWkt :: Geometry a -> BC.ByteString +writeWkt g = runGeos $ do + w <- S.createWktWriter + r :: R.Geom <- convertGeometryToRaw g + S.writeWkt w r diff --git a/src/Data/Geometry/Geos/Topology.hs b/src/Data/Geometry/Geos/Topology.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Topology.hs @@ -0,0 +1,102 @@ +{-# LANGUAGE CPP #-} +{-# LANGUAGE ScopedTypeVariables #-} +module Data.Geometry.Geos.Topology ( + envelope + , intersection + , convexHull + , difference + , symmetricDifference + , boundary + , union + , unaryUnion + , pointOnSurface + , centroid + , node + , delaunayTriangulation +) where +import Data.Geometry.Geos.Types +import Data.Geometry.Geos.Raw.Base +import Data.Geometry.Geos.Geometry +import qualified Data.Geometry.Geos.Raw.Topology as R +import qualified Data.Geometry.Geos.Raw.Geometry as RG + + +geo_1 ::( RG.GeomConst -> Geos RG.GeomConst ) + -> Geometry a + -> Some Geometry +geo_1 f g = runGeos $ do + geo <- convertGeometryToRaw g + convertGeometryFromRaw =<< f geo + + +geo_2 :: (RG.GeomConst -> RG.GeomConst -> Geos RG.GeomConst ) + -> Geometry a + -> Geometry b + -> Some Geometry +geo_2 f g1 g2 = runGeos $ do + g1' <- convertGeometryToRaw g1 + g2' <- convertGeometryToRaw g2 + convertGeometryFromRaw =<< f g1' g2' + +-- | Returns a Polygon that represents the bounding envelope of this geometry. Note that it can also return a Point if the input geometry is a point. +envelope :: Geometry a -> Some Geometry +envelope = geo_1 R.envelope + +-- | Returns a Geometry representing the points shared by both geometries. +intersection :: Geometry a -> Geometry b -> Some Geometry +intersection = geo_2 R.intersection + +-- | Returns the smallest Polygon that contains all the points in the geometry. +convexHull :: Geometry a -> Geometry Polygon +convexHull g = ensurePolygon $ geo_1 R.convexHull g + +-- | Returns a Geometry representing the points making up this geometry that do not make up other. +difference :: Geometry a -> Geometry b -> Some Geometry +difference = geo_2 R.difference + +-- | Returns a Geometry combining the points in this geometry not in other, and the points in other not in this geometry. +symmetricDifference :: Geometry a -> Geometry b -> Some Geometry +symmetricDifference = geo_2 R.symmetricDifference + +boundary :: Geometry a -> Some Geometry +boundary = geo_1 R.boundary + +-- | Returns a Geometry representing all the points in both geometries. + +{-| Computes the union of all the elements of this geometry. Heterogeneous GeometryCollections are fully supported. + +The result obeys the following contract: + +Unioning a set of LineStrings has the effect of fully noding and dissolving the linework. +Unioning a set of Polygons will always return a Polygonal geometry (unlike {link #union(Geometry)}, which may return geometrys of lower dimension if a topology collapse occurred. +-} +union :: Geometry a -> Geometry b -> Some Geometry +union = geo_2 R.union + + +unaryUnion :: Geometry a -> Some Geometry +unaryUnion = geo_1 R.unaryUnion + +-- | Computes and returns a Point guaranteed to be on the interior of this geometry. +pointOnSurface :: Geometry a -> Geometry Point +pointOnSurface g = ensurePoint $ geo_1 R.pointOnSurface g + +-- | Returns a Point object representing the geometric center of the geometry. The point is not guaranteed to be on the interior of the geometry. +centroid :: Geometry a -> Geometry Point +centroid g = ensurePoint $ geo_1 R.centroid g + +node :: Geometry a -> Some Geometry +node = geo_1 R.node + +-- | Return a Delaunay triangulation of the vertex of the given geometry @g@, where @tol@ is the snapping tolerance to use. +delaunayTriangulation :: Geometry a -> Double -> Geometry MultiLineString +delaunayTriangulation g d = ensureMultiLineString (geo_1 (flip R.delaunayTriangulation $ d) g) + +#if GEOS_VERSION_MAJOR > 3 && GEOS_VERSION_MINOR > 4 +voronoiDiagram :: Geometry a -> Geometry b -> Double -> Bool -> Some Geometry +voronoiDiagram g env tol onlyEdges = runGeos $ do + g' <- convertGeometryToRaw g + env' <- convertGeometryToRaw env + convertGeometryFromRaw =<< R.voronoiDiagram g' env' tol onlyEndges + +#endif diff --git a/src/Data/Geometry/Geos/Types.hs b/src/Data/Geometry/Geos/Types.hs new file mode 100644 --- /dev/null +++ b/src/Data/Geometry/Geos/Types.hs @@ -0,0 +1,96 @@ +{-# LANGUAGE RankNTypes, KindSignatures, GADTs, DeriveDataTypeable, StandaloneDeriving, FlexibleInstances #-} + +module Data.Geometry.Geos.Types where +import qualified Data.Vector as V +import Data.Monoid +import Data.Data + +class Relatable a where + contains :: a -> Geometry b -> Bool + coveredBy :: a -> Geometry b -> Bool + covers :: a -> Geometry b -> Bool + -- | Returns @True@ if the DE-9IM intersection matrix for the two Geometries is T*T****** (for a point and a curve,a point and an area or a line and an area) 0******** (for two curves). + crosses :: a -> Geometry b -> Bool + -- | Returns @True@ if the DE-9IM intersection matrix for the two geometries is FF*FF****. + disjoint :: a -> Geometry b -> Bool + -- | Returns @True@ if @disjoint@ is False. + intersects :: a -> Geometry b -> Bool + -- | Returns true if the DE-9IM intersection matrix for the two geometries is T*T***T** (for two points or two surfaces) 1*T***T** (for two curves). + overlaps :: a -> Geometry b -> Bool + -- | Returns True if the DE-9IM intersection matrix for the two geometries is FT*******, F**T***** or F***T****. + touches :: a -> Geometry b -> Bool + -- | Returns True if the DE-9IM intersection matrix for the two geometries is T*F**F***. + within :: a -> Geometry b -> Bool + + +type SRID = Maybe Int + +data Some :: (* -> *) -> * where + Some :: f a -> Some f + +withSomeGeometry :: Some Geometry -> (forall a . Geometry a -> b) -> b +withSomeGeometry (Some p) f = f p + +instance Show (Some Geometry) where + show (Some a) = "Some (" <> show a <> ")" + +data Geometry a where + PointGeometry :: Point -> SRID -> Geometry Point + LineStringGeometry :: LineString -> SRID -> Geometry LineString + LinearRingGeometry :: LinearRing -> SRID -> Geometry LinearRing + PolygonGeometry :: Polygon -> SRID -> Geometry Polygon + MultiPointGeometry :: MultiPoint -> SRID -> Geometry MultiPoint + MultiLineStringGeometry :: MultiLineString -> SRID -> Geometry MultiLineString + MultiPolygonGeometry :: MultiPolygon -> SRID -> Geometry MultiPolygon + {-CollectionGeometry :: GeometryCollection -> Geometry GeometryCollection-} + +deriving instance Eq (Geometry a) +deriving instance Show (Geometry a) + +data Coordinate = + Coordinate2 {-# UNPACK #-} !Double {-# UNPACK #-} !Double + | Coordinate3 {-# UNPACK #-} !Double {-# UNPACK #-} !Double {-# UNPACK #-} !Double deriving (Read, Ord, Show, Eq, Data, Typeable) + + +dimensionsCoordinate :: Coordinate -> Int +dimensionsCoordinate = length . gmapQ (const ()) + +type CoordinateSequence = V.Vector Coordinate + +dimensionsCoordinateSequence :: CoordinateSequence -> Int +dimensionsCoordinateSequence = dimensionsCoordinate . V.head + +newtype Point = Point Coordinate + deriving (Read, Ord, Show, Eq, Data, Typeable) + +-- A LinearRing is a LineString that is closed +newtype LinearRing = LinearRing CoordinateSequence + deriving (Read, Ord, Show, Eq, Data, Typeable) + +instance Monoid LinearRing where + mempty = LinearRing V.empty + mappend (LinearRing a) (LinearRing b) = LinearRing (a <> b) + +newtype LineString = LineString CoordinateSequence + deriving (Read, Ord, Show, Eq, Data, Typeable) + +instance Monoid LineString where + mempty = LineString V.empty + mappend (LineString a) (LineString b) = LineString (a <> b) + +-- | In a polygon, the fist LinearRing is the shell, and any following are holes. +newtype Polygon = Polygon (V.Vector LinearRing) + deriving (Read, Ord, Show, Eq, Data, Typeable) + +newtype MultiPoint = MultiPoint (V.Vector Point) + deriving (Read, Ord, Show, Eq, Data, Typeable) + +instance Monoid MultiPoint where + mempty = MultiPoint V.empty + mappend (MultiPoint a) (MultiPoint b) = MultiPoint (a <> b) + +newtype MultiLineString = MultiLineString (V.Vector LineString) + deriving (Read, Ord, Show, Eq, Data, Typeable) + +newtype MultiPolygon = MultiPolygon (V.Vector Polygon) + deriving (Read, Ord, Show, Eq, Data, Typeable) diff --git a/tests/Main.hs b/tests/Main.hs new file mode 100644 --- /dev/null +++ b/tests/Main.hs @@ -0,0 +1,25 @@ +{-# LANGUAGE OverloadedStrings, LambdaCase, ScopedTypeVariables #-} + +module Main where +import qualified Data.ByteString as BS +import Test.Hspec +import Control.Exception +import Data.Geometry.Geos.Types +import Data.Geometry.Geos.Serialize +import Data.Geometry.Geos.Raw.Base +import Data.Geometry.Geos.Geometry +import qualified Data.Geometry.Geos.Raw.CoordSeq as RC +import qualified Data.Geometry.Geos.Raw.Serialize as RS +import qualified Data.Geometry.Geos.Raw.Geometry as R +import qualified Data.Vector as V +import Debug.Trace + +import ParsingSpec +import RawGeometrySpec +import SpatialOperationsSpec + +main :: IO () +main = hspec $ do + rawGeometrySpecs + parsingSpecs + spatialOpsSpecs diff --git a/tests/ParsingSpec.hs b/tests/ParsingSpec.hs new file mode 100644 --- /dev/null +++ b/tests/ParsingSpec.hs @@ -0,0 +1,74 @@ +{-# LANGUAGE OverloadedStrings, LambdaCase, ScopedTypeVariables #-} + +module ParsingSpec where + +import Prelude hiding (read) +import qualified Data.ByteString as BS +import Test.Hspec +import Control.Exception +import Data.Geometry.Geos.Types +import Data.Geometry.Geos.Serialize +import Data.Geometry.Geos.Raw.Base +import Data.Geometry.Geos.Geometry +import qualified Data.Geometry.Geos.Raw.CoordSeq as RC +import qualified Data.Geometry.Geos.Raw.Serialize as RS +import qualified Data.Geometry.Geos.Raw.Geometry as R +import qualified Data.Vector as V +import Data.Maybe (fromJust) + +import SpecSampleData + +parsingSpecs = describe "Tests Serialization" $ do + it "Parses a bytestring to a linestring" $ do + ensureLineString <$> readHex linestringBS `shouldBe` Just linestring + it "Parse a multipolygon" $ do + ensureMultiPolygon <$> readHex multiPolygonStringBS `shouldBe` Just multiPolygon + it "Serializes a LineString into a bytestring" $ do + linestringBS `shouldBe` writeHex linestring + it "Serializes a LineString to WKT" $ do + linestringWkt `shouldBe` (writeWkt linestring) + it "Can parse WKT" $ do + ensureLineString <$> readWkt (Just 4326) linestringWkt `shouldBe` Just linestring + it "can parse lots of things" $ do + polygons <- (fmap ensurePolygon) <$> loadThingsFromFile "tests/sampledata/polygons.csv" + (length polygons) `shouldBe` 98 + points <- (fmap ensurePoint) <$> loadThingsFromFile "tests/sampledata/points.csv" + (length points) `shouldBe` 34582 + it "can read a long list" $ do + points <- (fmap ensurePoint) <$> loadThingsFromFile "tests/sampledata/points.csv" + (length points) `shouldBe` 34582 + -- Because we're in Australia (below the equator) all the points in this set should have a latitude < 0 + let filteredPoints = filter (\(PointGeometry (Point (Coordinate2 lat long)) _) -> long > 0) points + -- running this test multiple times will generate a random length list + length filteredPoints `shouldBe` 0 + it "can read a long list in one op" $ do + points <- (fmap ensurePoint) <$> loadThingsFromFile' "tests/sampledata/points.csv" + (length points) `shouldBe` 34582 + -- Because we're in Australia (below the equator) all the points in this set should have a latitude < 0 + let filteredPoints = filter (\(PointGeometry (Point (Coordinate2 lat long)) _) -> long > 0) points + -- running this test multiple times will generate a random length list + length filteredPoints `shouldBe` 0 + it "handles nonesense properly" $ do + let bs = "2D50491A5DC024275C1ED5DE4040" :: BS.ByteString + Nothing `shouldBe` (fmap ensurePoint $ readHex bs) + + it "can write/read with zero loss of information" $ do + let bs = "0105000020C46E000001000000010200000004000000EE2DE25E859A20410829F224364B5A41BD757FB6679A204115C6D8E9304B5A4177111862469A2041A03F04E32C4B5A41BC1D3ADD459A2041793C8CD92C4B5A41" :: BS.ByteString + srid = Just 28356 + points = [[(544066.685319362, 6892760.57728029), + (544051.856441192, 6892739.65385582), + (544035.191589876, 6892723.54713431), + (544034.932084016, 6892723.39918434)]] + expectedGeo = MultiLineStringGeometry (makeMultiLineString points) srid + parsedGeo = fromJust (fmap ensureMultiLineString $ readHex bs) + parsedGeo `shouldBe` expectedGeo + let encodedHex = writeHex expectedGeo + encodedHex `shouldBe` bs + + + + + +searchPoints :: [Geometry Point] -> Geometry Polygon -> [Geometry Point] +searchPoints points polygon = filter f points + where f point = contains polygon point diff --git a/tests/RawGeometrySpec.hs b/tests/RawGeometrySpec.hs new file mode 100644 --- /dev/null +++ b/tests/RawGeometrySpec.hs @@ -0,0 +1,98 @@ +{-# LANGUAGE OverloadedStrings, LambdaCase, ScopedTypeVariables #-} +module RawGeometrySpec where + +import qualified Data.ByteString as BS +import Test.Hspec +import Control.Exception +import Data.Geometry.Geos.Types +import Data.Geometry.Geos.Serialize +import Data.Geometry.Geos.Raw.Base +import Data.Geometry.Geos.Geometry +import qualified Data.Geometry.Geos.Raw.CoordSeq as RC +import qualified Data.Geometry.Geos.Raw.Serialize as RS +import qualified Data.Geometry.Geos.Raw.Geometry as R +import qualified Data.Vector as V +import Data.Maybe (fromJust) + +import SpecSampleData + +rawGeometrySpecs = describe "raw geometry" $ do + + it "Creates a Coordinate Sequence" $ do + let (size, dim) = runGeos $ do + cs :: RC.CoordSeq <- RC.createEmptyCoordinateSequence 2 2 + size <- RC.getCoordinateSequenceSize cs + dim <- RC.getCoordinateSequenceDimensions cs + return (size, dim) + size `shouldBe` (2 :: Int) + dim `shouldBe` (2 :: Int) + it "Sets a Coordinate Sequence" $ do + let (d1, d2) = runGeos $ do + c :: RC.CoordSeq <- RC.createEmptyCoordinateSequence 2 2 + RC.setCoordinateSequenceX c 0 5.0 + RC.setCoordinateSequenceY c 0 10.0 + d1 <- RC.getCoordinateSequenceX c 0 + d2 <- RC.getCoordinateSequenceY c 0 + return (d1, d2) + d1 `shouldBe` (5.0 :: Double) + d2 `shouldBe` (10.0 :: Double) + + it "Creates a LineString " $ do + let tid = runGeos $ do + cs :: RC.CoordSeqConst <- RC.createEmptyCoordinateSequence 2 2 + ls :: R.Geom <- R.createLineString cs + R.getTypeId ls + tid `shouldBe` R.LineStringTypeId + it "Creates a Geometry" $ do + pending + + it "Converts a LineString" $ do + let (srid, tid, tn) = runGeos $ do + l :: R.Geom <- convertGeometryToRaw linestring + t <- R.getTypeId l + s <- R.getSRID l + tn <- R.getTypeName l + return (s, t, tn) + tid `shouldBe` R.LineStringTypeId + srid `shouldBe` (Just 4326) + tn `shouldBe` "LineString" + it "Converts a Polygon" $ do + let t = runGeos $ do + rp :: R.Geom <- convertGeometryToRaw $ PolygonGeometry polygon1 Nothing + R.getTypeId rp + t `shouldBe` R.PolygonTypeId + it "Converts a MultiPolygon from Raw" $ do + let (x,y) = runGeos $ do + r <- RS.createReader + mpr <- RS.readHex r multiPolygonStringBS + MultiPolygon vps <- convertMultiPolygonFromRaw $ fromJust mpr + let (Polygon vlr) = vps V.! 0 + (LinearRing vc) = vlr V.! 0 + (Coordinate2 x y) = vc V.! 0 + return (x, y) + + x `shouldBe` 153.160525 + y `shouldBe` -27.377412 + it "Tests disjoint" $ do + (disjoint (PolygonGeometry polygon1 Nothing) (PolygonGeometry polygon2 Nothing)) `shouldBe` False + it "Tests intersects" $ do + intersects (PolygonGeometry polygon1 Nothing) (PolygonGeometry polygon2 Nothing) `shouldBe` True + it "Tests raw parsing" $ do + let (x,y) = runGeos $ do + r <- RS.createReader + g <- RS.readHex r multiPolygonStringBS + gi :: R.GeomConst <- R.getGeometryN (fromJust g) 0 + ir :: R.GeomConst <- R.getExteriorRing gi + cs <- R.getCoordinateSequence ir + x <- RC.getCoordinateSequenceX cs 0 + y <- RC.getCoordinateSequenceY cs 0 + return (x,y) + x `shouldBe` 153.160525 + y `shouldBe` -27.377412 + it "Tests Distance" $ do + let a = makePointGeo (0,0) + b = makePointGeo (0,2) + let d = distance a b + hd = hausdorffDistance a b + d `shouldBe` 2.0 + hd `shouldBe` 2.0 diff --git a/tests/SpatialOperationsSpec.hs b/tests/SpatialOperationsSpec.hs new file mode 100644 --- /dev/null +++ b/tests/SpatialOperationsSpec.hs @@ -0,0 +1,68 @@ +{-# LANGUAGE OverloadedStrings #-} +{-# LANGUAGE ScopedTypeVariables #-} + +module SpatialOperationsSpec where + +import Test.Hspec +import qualified Data.ByteString as BS +import qualified Data.Vector as V +import Data.Geometry.Geos.Types +import Data.Geometry.Geos.Geometry +import Data.Geometry.Geos.Topology + +import SpecSampleData + +spatialOpsSpecs = describe "Tests Contains" $ do + it "Does polygon to point comparison" $ do + let polygon = makePolygonGeo [[(0,0), (0,1), (1,1), (1,0), (0,0)]] + polygonWithHole = makePolygonGeo [ [(0,0), (0,1), (1,1), (1,0), (0,0)], [(0.1,0.1),(0.1,0.9),(0.9,0.1),(0.1,0.1)] ] + pointIn = makePointGeo (0.6, 0.6) + pointOut = makePointGeo (1.5, 0.5) + pointInHole = makePointGeo (0.4, 0.4) + (contains polygon pointIn) `shouldBe` True + (contains polygon pointOut) `shouldBe` False + (contains polygonWithHole pointIn) `shouldBe` True + (contains polygonWithHole pointOut) `shouldBe` False + (contains polygonWithHole pointInHole) `shouldBe` False + + it "Does simple polygon to polygon comparison" $ do + let polygonBig = makePolygonGeo [ [(0,0), (0,2), (2,2), (2,0), (0,0)] ] + polygonSmall = makePolygonGeo [ [(0,0), (0,1), (1,1), (1,0), (0,0)] ] + polygonIntersect = makePolygonGeo [ [(0,0), (0,2.5), (1,1), (1,0), (0,0)] ] + (contains polygonBig polygonSmall) `shouldBe` True + (contains polygonBig polygonBig) `shouldBe` True + (contains polygonBig polygonIntersect) `shouldBe` False + + it "Does multi polygon to point comparison" $ do + let polygonBig = [ [(0,0), (0,2), (2,2), (2,0), (0,0)]] + polygonSmall = [ [(3,0), (3,1), (4,1), (4,0), (3,0)]] + multiPoly = makeMultiPolygonGeo [polygonBig, polygonSmall] + pointIn = makePointGeo (0.5, 0.5) + pointOut = makePointGeo (2.5, 0.5) + (contains multiPoly pointIn) `shouldBe` True + (contains multiPoly pointOut) `shouldBe` False + it "Projects a point against a linestring" $ do + let lr = makeLineStringGeo [(0,0), (0, 1), (1, 1)] + p = makePointGeo (0.5, 1.0) +-- The point on this line string nearest (0.5, 1.0) is 1.5 units from the origin. i.e., halfway between the second and third point. + project lr p `shouldBe` 1.5 + interpolate lr 1.5 `shouldBe` p + projectNormalized lr p `shouldBe` 0.75 + interpolateNormalized lr 0.75 `shouldBe` p + it "Tests disjoint geometries" $ do + let poly = makePolygonGeo [[(0,0), (0,1), (1,1), (1,0), (0,0)]] + p1 = makePointGeo (2,2) + p2 = makePointGeo (0.5, 0.5) + disjoint poly p1 `shouldBe` True + disjoint poly p2 `shouldBe` False + it "creates an envelope / boundary of a geometry" $ do + let poly1 = makePolygonGeo [[(0,0), (0,1), (1,1), (1.5, 1.5), (1,0), (0,0)]] + env1 = makePolygonGeo [[(0.0, 0.0), (1.5, 0.0), (1.5, 1.5), (0.0, 1.5), (0.0, 0.0)]] + poly2 = makePolygonGeo [ [(0,0), (0,1), (1,1), (1,0), (0,0)], [(0.1,0.1),(0.1,0.9),(0.9,0.1),(0.1,0.1)] ] + env2 = makeMultiLineStringGeo [[(0.0, 0.0), (0.0, 1.0), (1.0, 1.0), (1.0, 0.0), (0.0, 0.0)], [(0.1, 0.1), (0.1, 0.9), (0.9, 0.1), (0.1, 0.1)]] + env3 = makePolygonGeo [ [(0,0), (0,1), (1,1), (1,0), (0,0)] ] + point = makePointGeo (2.5, 0.5) + (ensurePolygon $ envelope poly1) `shouldBe` env1 + (ensurePoint $ envelope point) `shouldBe` point + (ensureMultiLineString $ boundary poly2) `shouldBe` env2 + convexHull poly2 `shouldBe` env3 diff --git a/tests/SpecSampleData.hs b/tests/SpecSampleData.hs new file mode 100644 --- /dev/null +++ b/tests/SpecSampleData.hs @@ -0,0 +1,54 @@ +{-# LANGUAGE OverloadedStrings #-} + +module SpecSampleData where + +import qualified Data.ByteString as BS +import qualified Data.ByteString.Char8 as BS8 +import qualified Data.Vector as V +import Data.Geometry.Geos.Types +import Data.Geometry.Geos.Serialize +import Data.Geometry.Geos.CSV +import Data.Maybe (catMaybes) + +point = PointGeometry (Point $ Coordinate2 36.1 (-119.1)) (Just 4326) +makePoint (c1, c2) = Point $ Coordinate2 c1 c2 +makePointGeo c = PointGeometry (makePoint c) Nothing +makeLinearRing cs = LinearRing . V.map (uncurry Coordinate2) . V.fromList $ cs +makeLinearRingGeo cs = LinearRingGeometry (makeLinearRing cs) Nothing +makeLineString cs = LineString . V.map (uncurry Coordinate2) . V.fromList $ cs +makeLineStringGeo cs = LineStringGeometry (makeLineString cs ) Nothing +makePolygon = Polygon . V.fromList . (fmap makeLinearRing) +makePolygonGeo cs = PolygonGeometry (makePolygon cs) Nothing +makeMultiLineString = MultiLineString . V.fromList . (fmap makeLineString) +makeMultiLineStringGeo lss = MultiLineStringGeometry (makeMultiLineString lss) Nothing + +-- Expects a file containing one geom per line WKB(Hex) encoded +loadThingsFromFile :: FilePath -> IO [Some Geometry] +loadThingsFromFile fp = do + rows <- BS8.readFile fp + return $ catMaybes $ readHex <$> (BS8.lines rows) + +loadThingsFromFile' :: FilePath -> IO [Some Geometry] +loadThingsFromFile' fp = do + rows <- BS8.readFile fp + return $ readLotsOfHex (BS8.lines rows) + +makeMultiPolygon = MultiPolygon . V.fromList . (fmap makePolygon) +makeMultiPolygonGeo polygons = MultiPolygonGeometry (makeMultiPolygon polygons) Nothing + +polygon1 = makePolygon [[(0, 0), (0, 1), (1, 1), (1, 0), (0, 0)]] +polygon2 = makePolygon [[(0, 0), (0, 1), (1, 1), (1, 0), (0, 0)]] + +polygonBS = "0103000020E6100000010000000C00000073E92D50491A5DC024275C1ED5DE404076E933474C1A5DC02C279CD7DBDE40406EE9A178431A5DC034271059E6DE40406DE9851C431A5DC034271C7AE6DE40406CE9F7AA421A5DC03427A07DE6DE40406CE92B3B421A5DC03427E05BE6DE40406CE955E4411A5DC03427A023E6DE404066E9E9FB3B1A5DC0282718AAD8DE404062E905D4381A5DC02027E877D1DE40406CE95FF8411A5DC014274C4CC6DE40406CE9B5EC421A5DC01427C8A2C6DE404073E92D50491A5DC024275C1ED5DE4040" :: BS.ByteString + +linestringBS = "0102000020E610000037000000603A02BC406A5DC020650408151C4140643A94BD436A5DC020659862171C4140693A684D496A5DC01865E0D70E1C4140703AD69C4F6A5DC0146520320B1C4140753A7462556A5DC0146570780B1C41407B3A6CFB5A6A5DC0186534680C1C4140823A5295626A5DC0186500D80C1C4140863AC848666A5DC018659C390F1C4140883A5E29686A5DC01C65985A101C41408A3A6C466A6A5DC0186538520F1C41408E3A0AF06D6A5DC01865A0B00C1C4140923A4619726A5DC01865B4CA0B1C41409A3A30547A6A5DC01065D8E7061C41409D3AEEDE7C6A5DC01065A855051C41409F3AA2FB7E6A5DC010650095041C4140A03A6258806A5DC010656CDD041C4140A23A2EA9816A5DC01065484D041C4140A43AB895836A5DC00C658C83021C4140A63A78A6856A5DC00C6594E4021C4140AC3A24218C6A5DC00C65040C021C4140AE3AAE598D6A5DC00C6594E4021C4140AF3AF6798E6A5DC0106530CD051C4140B03A5CD68F6A5DC01465887F0A1C4140B23A103F916A5DC01865A8980D1C4140B43ACAFB926A5DC0186568A90F1C4140B53A90AC946A5DC01C656432121C4140B63A9060956A5DC02065087C151C4140B63A8A0C956A5DC02465E0DC181C4140B33A9A77926A5DC02865FC161D1C4140B23A76E7916A5DC02C654C00201C4140B23AF856916A5DC0306560D9231C4140B03A5636906A5DC03065383A271C4140AF3A026E8E6A5DC0386514852D1C4140AE3AC6F58D6A5DC03C6504662F1C4140AE3ADEDD8D6A5DC0406538F0341C4140AE3AAE0D8E6A5DC04065F800371C4140B23A76E7916A5DC0546534FF471C4140B83AB0C5976A5DC05C653063511C4140BA3A345E996A5DC06465B84D571C4140BC3A12AB9B6A5DC06C6510D05E1C4140C23A3A41A16A5DC07C6524F76E1C4140C43A5452A36A5DC08065B49F721C4140C63A329FA56A5DC080651C71751C4140CC3A12C9AA6A5DC0846584AA791C4140CF3A7A4EAE6A5DC08865EC7B7C1C4140D43A5424B36A5DC08C658C4E821C4140D53A8AFCB36A5DC09065883F861C4140D63A3811B56A5DC0946538C9881C4140D73AA40DB66A5DC0946580E9891C4140D83A28A6B76A5DC094655C59891C4140DA3A8E02B96A5DC09465E468881C4140DC3AAC00BB6A5DC0946538C9881C4140DE3A027CBC6A5DC09465D828891C4140DF3AB224BE6A5DC09465FC5E891C4140E33ABCBCC16A5DC094659C0A891C4140" :: BS.ByteString + +multiPolygonStringBS = "0106000020e6100000020000000103000000010000002e0000001826530523256340321ea5129e603bc0e7e26f7b8223634005a568e55e703bc0cade52ce97216340534145d5af703bc0d38558fd9121634051f4c0c760793bc04eb857e62d21634004711e4e60763bc09a42e735f620634071ab2006ba7a3bc09e616a4b1d206340906ad8ef89793bc0f3e84658d41f63402159c0046e7d3bc0d3d9c9e0a8206340709692e524803bc0226e4e25032063404f57772cb6813bc048c2be9d4420634023a46e675f853bc0994869360f1f6340c6f7c5a52a813bc0e605d847271f6340d28f8653e68a3bc08ac745b5081e634025b1a4dc7d863bc01ec539eae81c6340406ce9d1548b3bc0b84082e2c71c6340ed832c0b269a3bc022aaf067f81b63405c72dc291d943bc07a7077d66e1b63400d1afa27b8983bc035b22b2d231b6340882cd2c43b883bc0c828cfbc9c1a63401616dc0f78883bc03a5cab3dec1a6340ddcd531d728b3bc0378b170b431a634024d5777e518e3bc0cbf6216fb9176340de3d40f7e5903bc050508a56ee1763401cb5c2f4bd863bc0f4311f1068186340772b4b7496853bc09b560a815c17634064ac36ffaf7e3bc0e15e99b7ea17634077f35487dc703bc0fce4284094176340a9177c9a935f3bc0a530ef71a6166340bc9179e40f5a3bc03eca880bc0156340087250c24c5f3bc0f8e0b54b1b166340c4b0c398f4573bc05e2a36e6f5176340af08feb792513bc07a3881e9b41863400a2fc1a90f583bc0dc813ae55118634008e753c72a5d3bc04208c897d0186340ff3d78edd2663bc030849cf7ff196340ed647094bc663bc0d653abaf2e1b6340a69718cbf46b3bc0bb809719b61c63405f46b1dcd26a3bc03eeaaf57581d634000ac8e1ce9643bc0fa5e43705c1e634003eb387ea8683bc0e9edcf45431f6340b073d3669c663bc05f27f565e91f6340211ff46c564d3bc09f909db731226340e3a9471adc463bc0b212f3aca4226340d0b7054b75553bc0f14a92e7fa2463406571ff91e95c3bc01826530523256340321ea5129e603bc00103000000010000000b00000003b16ce610236340779fe3a3c5413bc01b2c9ca479226340dae731ca33433bc0221add416c22634083fa96395d3a3bc0ee5d83be74216340acaa97df693a3bc04b5af10d052163403140a20914353bc062dba2cc86216340dc2a8881ae313bc06953758fec22634078b988efc4343bc0569e40d8a9236340e2afc91af5303bc0478e7406c6236340807d74eaca3b3bc054e4107173236340d7868a71fe423bc003b16ce610236340779fe3a3c5413bc0" :: BS.ByteString + +linestringWkt = "LINESTRING (-117.6602010747678833 34.2193918248615319, -117.6603845546283651 34.2194636577717120, -117.6607240217068835 34.2192029805830202, -117.6611091701577152 34.2190916689088738, -117.6614614615480860 34.2191000508120453, -117.6618031079214433 34.2191286331018887, -117.6622670462620874 34.2191419603279314, -117.6624929385526173 34.2192146314284287, -117.6626075191690006 34.2192490810504921, -117.6627365166588390 34.2192175650945387, -117.6629600620164808 34.2191372664621554, -117.6632140336826353 34.2191098576387844, -117.6637163611398194 34.2189607435813059, -117.6638715101675672 34.2189127990951647, -117.6640004238383739 34.2188898326804747, -117.6640835723178498 34.2188984660407414, -117.6641638709502615 34.2188812831392397, -117.6642813014137232 34.2188267169495646, -117.6644073652374516 34.2188382839759413, -117.6648028234291701 34.2188124677141730, -117.6648773385483935 34.2188382839759413, -117.6649460701544143 34.2189270483305563, -117.6650291348148585 34.2190703788748181, -117.6651151331414269 34.2191649267426214, -117.6652211642165753 34.2192279586544714, -117.6653243454446311 34.2193053236207732, -117.6653672607888836 34.2194056550017649, -117.6653472280403037 34.2195087524108033, -117.6651896482606361 34.2196377499006417, -117.6651552824576186 34.2197265980742884, -117.6651208328355835 34.2198440285377501, -117.6650520174105168 34.2199471259467600, -117.6649432203073360 34.2201391553484768, -117.6649145541984751 34.2201964875661986, -117.6649088545043185 34.2203655505531970, -117.6649202538926318 34.2204285824650469, -117.6651552824576186 34.2209471708144122, -117.6655134411802237 34.2212337480839039, -117.6656108388951054 34.2214142942782757, -117.6657512357732571 34.2216434555110425, -117.6660922115943606 34.2221363952366744, -117.6662183592371207 34.2222480421869477, -117.6663587561152724 34.2223340405134877, -117.6666739156746075 34.2224629541842944, -117.6668888276719684 34.2225489525108628, -117.6671839544827094 34.2227266488581279, -117.6672355031872286 34.2228469291686679, -117.6673014687652028 34.2229243779540013, -117.6673616508299887 34.2229587437570046, -117.6674590485448562 34.2229415608555030, -117.6675421132053145 34.2229128947466563, -117.6676637346203620 34.2229243779540013, -117.6677541753556113 34.2229357773423146, -117.6678554287459377 34.2229422314077567, -117.6680747831519938 34.2229321731239509)" :: BS.ByteString + +multiPolygonWkt = "MULTIPOLYGON(((153.160525 -27.377412,153.109678 -27.438948,153.049781 -27.440183,153.049071 -27.474133,153.036853 -27.462407,153.030055 -27.479401,153.003576 -27.474761,152.994671 -27.48996,153.020615 -27.500563,153.000384 -27.506686,153.008376 -27.520987,152.970607 -27.504557,152.973545 -27.542577,152.938563 -27.525358,152.903432 -27.544263,152.8994 -27.602143,152.874073 -27.57857,152.85728 -27.59656,152.848044 -27.532162,152.831633 -27.533082,152.841338 -27.54471,152.820684 -27.555931,152.741386 -27.566009,152.747844 -27.526336,152.762703 -27.521827,152.730042 -27.494873,152.747402 -27.440865,152.736847 -27.373346,152.707818 -27.351805,152.679693 -27.372265,152.690832 -27.343576,152.748767 -27.318645,152.772084 -27.343989,152.759997 -27.363934,152.775463 -27.401656,152.812496 -27.401315,152.849449 -27.421704,152.897229 -27.41728,152.917034 -27.394182,152.948784 -27.408821,152.976962 -27.400824,152.997241 -27.3021,153.068569 -27.276796,153.082602 -27.333821,153.155628 -27.362939,153.160525 -27.377412)),((153.095813 -27.256922,153.077349 -27.262509,153.075715 -27.227985,153.045501 -27.228178,153.031867 -27.207337,153.047705 -27.194069,153.091377 -27.20613,153.114483 -27.19124,153.117923 -27.233565,153.107842 -27.261695,153.095813 -27.256922)))" :: BS.ByteString + +linestring = LineStringGeometry (LineString $ V.fromList [Coordinate2 (-117.66020107476788) 34.21939182486153,Coordinate2 (-117.66038455462837) 34.21946365777171,Coordinate2 (-117.66072402170688) 34.21920298058302,Coordinate2 (-117.66110917015772) 34.219091668908874,Coordinate2 (-117.66146146154809) 34.219100050812045,Coordinate2 (-117.66180310792144) 34.21912863310189,Coordinate2 (-117.66226704626209) 34.21914196032793,Coordinate2 (-117.66249293855262) 34.21921463142843,Coordinate2 (-117.662607519169) 34.21924908105049,Coordinate2 (-117.66273651665884) 34.21921756509454,Coordinate2 (-117.66296006201648) 34.219137266462155,Coordinate2 (-117.66321403368264) 34.219109857638784,Coordinate2 (-117.66371636113982) 34.218960743581306,Coordinate2 (-117.66387151016757) 34.218912799095165,Coordinate2 (-117.66400042383837) 34.218889832680475,Coordinate2 (-117.66408357231785) 34.21889846604074,Coordinate2 (-117.66416387095026) 34.21888128313924,Coordinate2 (-117.66428130141372) 34.218826716949565,Coordinate2 (-117.66440736523745) 34.21883828397594,Coordinate2 (-117.66480282342917) 34.21881246771417,Coordinate2 (-117.6648773385484) 34.21883828397594,Coordinate2 (-117.66494607015441) 34.218927048330556,Coordinate2 (-117.66502913481486) 34.21907037887482,Coordinate2 (-117.66511513314143) 34.21916492674262,Coordinate2 (-117.66522116421658) 34.21922795865447,Coordinate2 (-117.66532434544463) 34.21930532362077,Coordinate2 (-117.66536726078888) 34.219405655001765,Coordinate2 (-117.6653472280403) 34.2195087524108,Coordinate2 (-117.66518964826064) 34.21963774990064,Coordinate2 (-117.66515528245762) 34.21972659807429,Coordinate2 (-117.66512083283558) 34.21984402853775,Coordinate2 (-117.66505201741052) 34.21994712594676,Coordinate2 (-117.66494322030734) 34.22013915534848,Coordinate2 (-117.66491455419848) 34.2201964875662,Coordinate2 (-117.66490885450432) 34.2203655505532,Coordinate2 (-117.66492025389263) 34.22042858246505,Coordinate2 (-117.66515528245762) 34.22094717081441,Coordinate2 (-117.66551344118022) 34.221233748083904,Coordinate2 (-117.6656108388951) 34.221414294278276,Coordinate2 (-117.66575123577326) 34.22164345551104,Coordinate2 (-117.66609221159436) 34.222136395236674,Coordinate2 (-117.66621835923712) 34.22224804218695,Coordinate2 (-117.66635875611527) 34.22233404051349,Coordinate2 (-117.66667391567461) 34.222462954184294,Coordinate2 (-117.66688882767197) 34.22254895251086,Coordinate2 (-117.66718395448271) 34.22272664885813,Coordinate2 (-117.66723550318723) 34.22284692916867,Coordinate2 (-117.6673014687652) 34.222924377954,Coordinate2 (-117.66736165082999) 34.222958743757005,Coordinate2 (-117.66745904854486) 34.2229415608555,Coordinate2 (-117.66754211320531) 34.222912894746656,Coordinate2 (-117.66766373462036) 34.222924377954,Coordinate2 (-117.66775417535561) 34.222935777342315,Coordinate2 (-117.66785542874594) 34.22294223140776,Coordinate2 (-117.668074783152) 34.22293217312395]) (Just 4326) + +multiPolygon = MultiPolygonGeometry (MultiPolygon $ V.fromList [Polygon $ V.fromList [LinearRing $ V.fromList [Coordinate2 153.160525 (-27.377412),Coordinate2 153.109678 (-27.438948),Coordinate2 153.049781 (-27.440183),Coordinate2 153.049071 (-27.474133),Coordinate2 153.036853 (-27.462407),Coordinate2 153.030055 (-27.479401),Coordinate2 153.003576 (-27.474761),Coordinate2 152.994671 (-27.48996),Coordinate2 153.020615 (-27.500563),Coordinate2 153.000384 (-27.506686),Coordinate2 153.008376 (-27.520987),Coordinate2 152.970607 (-27.504557),Coordinate2 152.973545 (-27.542577),Coordinate2 152.938563 (-27.525358),Coordinate2 152.903432 (-27.544263),Coordinate2 152.8994 (-27.602143),Coordinate2 152.874073 (-27.57857),Coordinate2 152.85728 (-27.59656),Coordinate2 152.848044 (-27.532162),Coordinate2 152.831633 (-27.533082),Coordinate2 152.841338 (-27.54471),Coordinate2 152.820684 (-27.555931),Coordinate2 152.741386 (-27.566009),Coordinate2 152.747844 (-27.526336),Coordinate2 152.762703 (-27.521827),Coordinate2 152.730042 (-27.494873),Coordinate2 152.747402 (-27.440865),Coordinate2 152.736847 (-27.373346),Coordinate2 152.707818 (-27.351805),Coordinate2 152.679693 (-27.372265),Coordinate2 152.690832 (-27.343576),Coordinate2 152.748767 (-27.318645),Coordinate2 152.772084 (-27.343989),Coordinate2 152.759997 (-27.363934),Coordinate2 152.775463 (-27.401656),Coordinate2 152.812496 (-27.401315),Coordinate2 152.849449 (-27.421704),Coordinate2 152.897229 (-27.41728),Coordinate2 152.917034 (-27.394182),Coordinate2 152.948784 (-27.408821),Coordinate2 152.976962 (-27.400824),Coordinate2 152.997241 (-27.3021),Coordinate2 153.068569 (-27.276796),Coordinate2 153.082602 (-27.333821),Coordinate2 153.155628 (-27.362939),Coordinate2 153.160525 (-27.377412)]],Polygon $ V.fromList [LinearRing $ V.fromList [Coordinate2 153.095813 (-27.256922),Coordinate2 153.077349 (-27.262509),Coordinate2 153.075715 (-27.227985),Coordinate2 153.045501 (-27.228178),Coordinate2 153.031867 (-27.207337),Coordinate2 153.047705 (-27.194069),Coordinate2 153.091377 (-27.20613),Coordinate2 153.114483 (-27.19124),Coordinate2 153.117923 (-27.233565),Coordinate2 153.107842 (-27.261695),Coordinate2 153.095813 (-27.256922)]]]) (Just 4326)