From 6d49003140f43fe7dcb47bb94c2520ce94b2ba0c Mon Sep 17 00:00:00 2001 From: rafaqz Date: Mon, 1 Jan 2024 20:19:55 +0100 Subject: [PATCH 1/3] return tuple points from linestring in geointerface --- src/geo_interface.jl | 32 +++++++++++++++++++++++++++----- 1 file changed, 27 insertions(+), 5 deletions(-) diff --git a/src/geo_interface.jl b/src/geo_interface.jl index 244d834..3b954fc 100644 --- a/src/geo_interface.jl +++ b/src/geo_interface.jl @@ -47,11 +47,33 @@ 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) + refs = Ref{Float64}(), Ref{Float64}(), Ref{Float64}() # 3 Refs is faster than a Vector + seq = getCoordSeq(geom::Union{LineString, LinearRing}) + return _get_tuple_point(geom, seq, refs, i) +end +function GeoInterface.getgeom(::AbstractGeometryTrait, geom::Union{LineString,LinearRing}) + context = get_context(geom) + seq = getCoordSeq(geom, context) + n = getSize(seq, context) # Faster thatn `GI.ngeom(geom)` when we already have `seq` + # Preallocate refse + refs = Ref{Float64}(), Ref{Float64}(), Ref{Float64}() + # Pretetermin if there is a z coordinate + hasz = hasZ(geom, context) + # We specify Unit32 in the generator to avoid unnecesary conversion later + return (_get_tuple_point(geom, seq, refs, i, context, hasz) for i in UInt32(1):UInt32(n)) +end + +function _get_tuple_point(geom, seq, (x, y, z), i, context=get_context(geom), hasz=hasZ(geom, context)) + if hasz + GEOSCoordSeq_getXYZ_r(context, seq, i - UInt32(1), x, y, z) + return x[], y[], Z[] + else + GEOSCoordSeq_getXY_r(context, seq, i - UInt32(1), x, y) + return x[], y[] + end +end + GeoInterface.getgeom(t::AbstractPointTrait, geom::PreparedGeometry) = nothing function GeoInterface.getgeom(::PolygonTrait, geom::Polygon, i::Int) if i == 1 From dbec28cfc34e22e56d22acbfa22087a3324d1d00 Mon Sep 17 00:00:00 2001 From: rafaqz Date: Sun, 7 Jan 2024 15:10:59 +0100 Subject: [PATCH 2/3] fix tests --- test/test_geo_interface.jl | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/test/test_geo_interface.jl b/test/test_geo_interface.jl index 41c3be5..412117f 100644 --- a/test/test_geo_interface.jl +++ b/test/test_geo_interface.jl @@ -82,9 +82,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 GeoInterface.testgeometry(ls) + @test collect(GeoInterface.getpoint(ls)) == [(8, 1), (9, 1), (9, 2), (8, 2)] @test GeoInterface.is3d(ls) == false Plots.plot(ls) Makie.plot(ls) @@ -121,10 +122,11 @@ 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.is3d(lr) == false @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 # Plots.plot(lr) # Makie.plot(lr) From 1e1fff381945a356015aacd54b3e4e32441b4efd Mon Sep 17 00:00:00 2001 From: rafaqz Date: Sun, 7 Jan 2024 15:13:16 +0100 Subject: [PATCH 3/3] fix ambiguity --- src/geo_interface.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/geo_interface.jl b/src/geo_interface.jl index 3b954fc..9c0a3f4 100644 --- a/src/geo_interface.jl +++ b/src/geo_interface.jl @@ -47,7 +47,7 @@ GeoInterface.getgeom(::MultiLineStringTrait, geom::MultiLineString, i) = getGeometry(geom, i)::LineString GeoInterface.getgeom(::MultiPolygonTrait, geom::MultiPolygon, i) = getGeometry(geom, i)::Polygon -function GeoInterface.getgeom(::AbstractGeometryTrait, geom::Union{LineString,LinearRing}, i) +function GeoInterface.getgeom(::Union{GI.LineStringTrait,GI.LinearRingTrait}, geom::Union{LineString,LinearRing}, i) refs = Ref{Float64}(), Ref{Float64}(), Ref{Float64}() # 3 Refs is faster than a Vector seq = getCoordSeq(geom::Union{LineString, LinearRing}) return _get_tuple_point(geom, seq, refs, i)