add coordtype (#167)

Closes #128

Currently on a GI.MultiPolygon this takes 4ns. Its not a compile time operation. We could put this in the type parameters of the wrapper geometries at construction so its instant.

A lot of packages will be able to set this absolutely to a single type, e.g. Shapefile.jl, GeoJSON.jl and LibGEOS.jl can return Float64

Co-authored-by: Anshul Singhvi <[email protected]>

Author: rafaqz

src/fallbacks.jl

@@ -8,6 +8,11 @@ const default_coord_names = (:X, :Y, :Z, :M)
 
 coordnames(t::AbstractGeometryTrait, geom) = default_coord_names[1:ncoord(t, geom)]
 
+@inline coordtype(::FeatureCollectionTrait, fc) = coordtype(getfeature(fc, 1))
+@inline coordtype(::FeatureTrait, feature) = coordtype(geometry(feature))
+@inline coordtype(::AbstractGeometryTrait, geom) = coordtype(getgeom(geom, 1))
+@inline coordtype(::PointTrait, geom) = typeof(x(geom))
+
 # Maybe hardcode dimension order? At least for X and Y?
 x(t::AbstractPointTrait, geom) = getcoord(t, geom, findfirst(isequal(:X), coordnames(geom)))
 y(t::AbstractPointTrait, geom) = getcoord(t, geom, findfirst(isequal(:Y), coordnames(geom)))

src/interface.jl

@@ -597,6 +597,15 @@ Returns a geometric object that represents the convex hull of the given `geom`.
 """
 convexhull(geom) = convexhull(geomtrait(geom), geom)
 
+"""
+    coordtype(geom) -> Type
+
+Return the type of the coordinates in the geom. 
+
+Usually `<: AbstractFloat` but can be any `Number`.
+"""
+coordtype(geom) = coordtype(trait(geom), geom)
+
 """
     x(geom) -> Number
 

src/wrappers.jl

@@ -434,6 +434,11 @@ convert(::Type{Point}, ::PointTrait, geom) = Point(geom)
 convert(::Type{Point}, ::PointTrait, geom::Point) = geom
 extent(trait::PointTrait, geom::Point) = extent(trait, parent(geom)) 
 
+coordtype(::Point{<:Any,<:Any,<:AbstractArray{T}}) where T = T 
+coordtype(::Point{<:Any,<:Any,<:NTuple{<:Any,T}}) where T = T 
+coordtype(::Point{<:Any,<:Any,<:NamedTuple{<:Any,<:NTuple{<:Any,T}}}) where T = T 
+coordtype(p::Point) = coordtype(parent(p))
+
 x(trait::PointTrait, geom::Point) = x(trait, parent(geom))
 y(trait::PointTrait, geom::Point) = y(trait, parent(geom))
 z(trait::PointTrait, geom::Point{true}) = z(trait, parent(geom))

test/test_wrappers.jl

@@ -23,6 +23,7 @@ GI.getcoord(point, 1)
 @test !GI.ismeasured(point)
 @test !GI.is3d(point)
 @test GI.ncoord(point) == 2
+@test GI.coordtype(point) == Int
 @test GI.extent(point) == Extent(X=(1, 1), Y=(2, 2))
 @test point == GI.Point(point)
 @test (GI.x(point), GI.y(point)) == (1, 2)
@@ -42,6 +43,7 @@ pointz = GI.Point(1, 2, 3)
 @test !GI.ismeasured(pointz)
 @test GI.is3d(pointz)
 @test GI.ncoord(pointz) == 3
+@test GI.coordtype(pointz) == Int
 @test (GI.x(pointz), GI.y(pointz), GI.z(pointz)) == (1, 2, 3)
 @test GI.testgeometry(pointz)
 @test GI.convert(GI, pointz) === pointz
@@ -55,6 +57,7 @@ pointzm = GI.Point(; X=1, Y=2, Z=3, M=4)
 @test GI.ismeasured(pointzm)
 @test GI.is3d(pointzm)
 @test GI.ncoord(pointzm) == 4
+@test GI.coordtype(pointzm) == Int
 @test pointzm == GI.Point(pointzm)
 @test point != GI.Point(pointzm)
 @test (GI.x(pointzm), GI.y(pointzm), GI.z(pointzm), GI.m(pointzm)) == (1, 2, 3, 4)
@@ -71,6 +74,7 @@ pointm = GI.Point((X=1, Y=2, M=3))
 @test GI.ismeasured(pointm)
 @test !GI.is3d(pointm)
 @test GI.ncoord(pointm) == 3
+@test GI.coordtype(pointm) == Int
 @test pointm == GI.Point(pointm)
 @test point != GI.Point(pointm)
 @test (GI.x(pointm), GI.y(pointm), GI.m(pointm)) == (1, 2, 3)
@@ -88,6 +92,7 @@ pointtm = GI.Point{false,true}(1, 2, 3)
 @test GI.ismeasured(pointtm)
 @test !GI.is3d(pointtm)
 @test GI.ncoord(pointtm) == 3
+@test GI.coordtype(pointtm) == Int
 @test (GI.x(pointtm), GI.y(pointtm), GI.m(pointtm)) == (1, 2, 3)
 @test_throws ArgumentError GI.z(pointtm)
 @test GI.testgeometry(pointtm)
@@ -102,6 +107,7 @@ pointa = GI.Point([1, 2])
 @test !GI.ismeasured(pointa)
 @test !GI.is3d(pointa)
 @test GI.ncoord(pointa) == 2
+@test GI.coordtype(pointa) == Int
 @test (GI.x(pointa), GI.y(pointa)) == (1, 2)
 @test GI.testgeometry(pointa)
 test_display(pointa, "Point{false, false}((1, 2))", "Point((1,2))")
@@ -159,6 +165,7 @@ linestring = GI.LineString([(1, 2), (3, 4)])
 @test GI.testgeometry(linestring)
 @test !GI.is3d(linestring)
 @test GI.ncoord(linestring) == 2
+@test GI.coordtype(linestring) == Int
 test_display(linestring, "LineString{false, false}([(1, 2), (3, 4)])", "LineString([(1,2),(3,4)])")
 @test @inferred(GI.extent(linestring)) == Extent(X=(1, 3), Y=(2, 4))
 @test_throws ArgumentError GI.LineString([(1, 2)])
@@ -242,6 +249,7 @@ collection = GI.GeometryCollection(geoms)
 @test GI.testgeometry(collection)
 @test !GI.is3d(collection)
 @test GI.ncoord(collection) == 2
+@test GI.coordtype(collection) == Int
 @test GI.extent(collection) == reduce(Extents.union, map(GI.extent, geoms))
 test_display(collection, "GeometryCollection{false, false}([Line([(1, 2), (3, 4)]), … (3) … , (1, 2)])",
                     "GeometryCollection([Line([(1,2),(3,4)]),LineString([(1,2),(3,4)]),…(2)…,(1,2)])")
@@ -259,6 +267,7 @@ multicurve = GI.MultiCurve([linestring, linearring])
 @test GI.getgeom(multicurve, 1) === linestring
 @test !GI.is3d(multicurve)
 @test GI.ncoord(multicurve) == 2
+@test GI.coordtype(multicurve) == Int
 @test GI.extent(multicurve) == Extent(X=(1, 5), Y=(2, 6))
 @test_throws ArgumentError GI.MultiCurve([pointz, polygon])
 @test GI.testgeometry(multicurve)
@@ -278,6 +287,7 @@ multipolygon = GI.MultiPolygon([polygon])
 @test GI.getgeom(multipolygon, 1) === polygon
 @test !GI.is3d(multipolygon)
 @test GI.ncoord(multipolygon) == 2
+@test GI.coordtype(multipolygon) == Int
 test_display(multipolygon, "MultiPolygon{false, false}([Polygon([LinearRing([(1, 2), … (2) … , (1, 2)]), LinearRing([(1, 2), … (2) … , (1, 2)])])])",
                             "MultiPolygon([Polygon([LinearRing([(1,2),…(2)…,(1,2)]),LinearRing([(1,2),…(2)…,(1,2)])])])")
 # MultiPolygon extent does not infer, maybe due to nesting
@@ -297,6 +307,7 @@ polyhedralsurface = GI.PolyhedralSurface([polygon, polygon])
 @test GI.PolyhedralSurface(polygon) == GI.PolyhedralSurface(polygon)
 @test !GI.is3d(polyhedralsurface)
 @test GI.ncoord(polyhedralsurface) == 2
+@test GI.coordtype(polyhedralsurface) == Int
 @test @inferred(GI.extent(polyhedralsurface)) == Extent(X=(1, 5), Y=(2, 6))
 @test GI.getgeom(polyhedralsurface, 1) === polygon
 @test collect(GI.getgeom(polyhedralsurface)) == [polygon, polygon]
@@ -315,6 +326,7 @@ test_display(polyhedralsurface_crs, "PolyhedralSurface{false, false}([Polygon([L
 multipolygon_coords = [[[[1, 2], [3, 4], [3, 2], [1, 4]]]]
 multipolygon = GI.MultiPolygon(multipolygon_coords)
 @test GI.coordinates(multipolygon) == multipolygon_coords
+@test GI.coordtype(multipolygon) == Int
 test_display(multipolygon, "MultiPolygon{false, false}([Polygon([LinearRing([[1, 2], [3, 4], [3, 2], [1, 4]])])])",
                         "MultiPolygon([Polygon([LinearRing([[1,2],[3,4],[3,2],[1,4]])])])")
 
@@ -345,6 +357,7 @@ test_display(feature, "Feature(MultiPolygon{false, false}([Polygon([LinearRing([
 @test GI.crs(feature) == EPSG(4326)
 @test GI.extent(feature) == GI.extent(multipolygon) 
 @test GI.testfeature(feature)
+@test GI.coordtype(feature) == Int
 @test_throws ArgumentError GI.Feature(:not_a_feature; properties=(x=1, y=2, z=3))
 @test GI.properties(GI.Feature(multipolygon)) == NamedTuple()
 
@@ -364,6 +377,7 @@ test_display(fc, "FeatureCollection([Feature(MultiPolygon{false, false}([Polygon
 @test_throws ArgumentError GI.FeatureCollection([1])
 vecfc = GI.FeatureCollection([(geometry=(1,2), a=1, b=2)])
 @test GI.getfeature(vecfc, 1) == (geometry=(1,2), a=1, b=2)
+@test GI.coordtype(vecfc) == Int