return tuple points from GeoInterface methods

JuliaGeo · Jan 1, 2024 · 519051a · 519051a
1 parent feb0055
commit 519051a
Show file tree

Hide file tree

Showing 2 changed files with 38 additions and 7 deletions.
diff --git a/src/geo_interface.jl b/src/geo_interface.jl
@@ -43,11 +43,40 @@ GeoInterface.getgeom(::MultiLineStringTrait, geom::MultiLineString, i) =
     getGeometry(geom, i)::LineString
 GeoInterface.getgeom(::MultiPolygonTrait, geom::MultiPolygon, i) =
     getGeometry(geom, i)::Polygon
-GeoInterface.getgeom(
-    ::Union{LineStringTrait,LinearRingTrait},
-    geom::Union{LineString,LinearRing},
-    i,
-) = getPoint(geom, i)
+function GeoInterface.getgeom(::AbstractGeometryTrait, geom::Union{LineString,LinearRing}, i)
+    ref = Ref{Float64}()
+    seq = getCoordSeq(geom::Union{LineString, LinearRing})
+    _get_tuple_point(geom, seq, ref, i)
+end
+function GeoInterface.getgeom(::AbstractGeometryTrait, geom::Union{LineString,LinearRing})
+    n = GeoInterface.ngeom(geom)
+    seq = getCoordSeq(geom::Union{LineString, LinearRing})
+    ref = Ref{Float64}()
+    return (_get_tuple_point(geom, seq, ref, i) for i in 1:n)
+ end
+
+function _get_tuple_point(geom, seq, ref, i)
+    context = get_global_context()
+    if GeoInterface.is3d(geom)
+        GC.@preserve ref begin
+            GEOSCoordSeq_getX_r(context, seq, i - 1, ref)
+            x = ref[]
+            GEOSCoordSeq_getY_r(context, seq, i - 1, ref)
+            y = ref[]
+            GEOSCoordSeq_getZ_r(context, seq, i - 1, ref)
+            z = ref[]
+        end
+        return (x, y, z)
+    else
+        GC.@preserve ref begin
+            GEOSCoordSeq_getX_r(context, seq, i - 1, ref)
+            x = ref[]
+            GEOSCoordSeq_getY_r(context, seq, i - 1, ref)
+            y = ref[]
+        end
+        return (x, y)
+    end
+end
 GeoInterface.getgeom(t::AbstractPointTrait, geom::PreparedGeometry) = nothing
 function GeoInterface.getgeom(::PolygonTrait, geom::Polygon, i::Int)
     if i == 1

diff --git a/test/test_geo_interface.jl b/test/test_geo_interface.jl
@@ -71,9 +71,10 @@ const LG = LibGEOS
     @test GeoInterface.geomtrait(ls) == LineStringTrait()
     p = GeoInterface.ngeom(ls) == 4
     p = GeoInterface.getgeom(ls, 3)
-    @test p isa LibGEOS.Point
+    @test p isa Tuple{Float64,Float64}
     @test GeoInterface.coordinates(p) == [9, 2]
     @test GeoInterface.testgeometry(ls)
+    @test collect(GeoInterface.getpoint(ls)) == [(8, 1), (9, 1), (9, 2), (8, 2)]
     plot(ls)
 
     @inferred GeoInterface.ncoord(ls)
@@ -104,9 +105,10 @@ const LG = LibGEOS
     @test GeoInterface.geomtrait(lr) == LinearRingTrait()
     @test GeoInterface.ngeom(lr) == 5
     p = GeoInterface.getgeom(lr, 3)
-    @test p isa LibGEOS.Point
+    @test p isa Tuple{Float64,Float64}
     @test GeoInterface.coordinates(p) == [9, 2]
     @test GeoInterface.testgeometry(lr)
+    @test collect(GeoInterface.getpoint(lr)) == [(8, 1), (9, 1), (9, 2), (8, 2), (8, 1)]
     # Cannot convert LinearRingTrait to series data for plotting
     # plot(lr)