From 20b90ceff09fc758f11fd512dedcc185df367da7 Mon Sep 17 00:00:00 2001 From: rafaqz Date: Tue, 17 Jan 2023 01:05:00 +0000 Subject: [PATCH 1/3] fix base performance and correctness, and test --- src/base.jl | 79 +++++++++++++++++++++++++-------- test/test_primitives.jl | 96 +++++++++++++++++++++++++++++++++++++---- 2 files changed, 150 insertions(+), 25 deletions(-) diff --git a/src/base.jl b/src/base.jl index 788b704c..16584990 100644 --- a/src/base.jl +++ b/src/base.jl @@ -1,47 +1,89 @@ # Implementation of GeoInterface for Base Types -const PointTuple2 = Tuple{<:Real,<:Real} -const PointTuple3 = Tuple{<:Real,<:Real,<:Real} -const PointTuple4 = Tuple{<:Real,<:Real,<:Real,<:Real} -const PointTuple = Union{PointTuple2,PointTuple3,PointTuple4} +# AbstractVector{<:Real} length 2 - 4 (where length checking is possible and efficient) GeoInterface.isgeometry(::Type{<:AbstractVector{<:Real}}) = true GeoInterface.geomtrait(::AbstractVector{<:Real}) = PointTrait() -GeoInterface.ncoord(::PointTrait, geom::AbstractVector{<:Real}) = Base.length(geom) +GeoInterface.ncoord(::PointTrait, geom::AbstractVector{<:Real}) = Base.length(geom) # TODO should this error for length > 4 ? GeoInterface.getcoord(::PointTrait, geom::AbstractVector{<:Real}, i) = getindex(geom, i) +GeoInterface.is3d(::PointTrait, geom::AbstractVector{<:Real}) = Base.length(geom) in (3, 4) +GeoInterface.ismeasured(::PointTrait, geom::AbstractVector{<:Real}) = Base.length(geom) == 4 # Measured must have Z + +# x/y/z/m methods were 100x slower without these custom definitions +# Also allow @inbounds to make them slightly faster when you know the sizes +Base.@propagate_inbounds function GeoInterface.x(::PointTrait, geom::AbstractVector{<:Real}) + @boundscheck Base.length(geom) in (2, 3, 4) || _xy_error(Base.length(geom)) + @inbounds geom[begin] +end +Base.@propagate_inbounds function GeoInterface.y(::PointTrait, geom::AbstractVector{<:Real}) + @boundscheck Base.length(geom) in (2, 3, 4) || _xy_error(Base.length(geom)) + @inbounds geom[begin + 1] +end +Base.@propagate_inbounds function GeoInterface.z(::PointTrait, geom::AbstractVector{<:Real}) + @boundscheck Base.length(geom) in (3, 4) || _z_error(Base.length(geom)) + @inbounds geom[begin + 2] +end +Base.@propagate_inbounds function GeoInterface.m(::PointTrait, geom::AbstractVector{<:Real}) + @boundscheck Base.length(geom) == 4 || _m_error(Base.length(geom)) + @inbounds geom[end] +end + +@noinline _xy_error(l) = throw(ArgumentError("Length of point must be 2, 3 or 4 to use `GeoInterface.x(point)` or `GeoInterface.y(point)`, got $l")) +@noinline _z_error(l) = throw(ArgumentError("Length of point must be 3 or 4 to use `GeoInterface.z(point)`, got $l")) +@noinline _m_error(l) = throw(ArgumentError("Length of point must be 4 to use `GeoInterface.m(point)`, got $l")) + + +# Tuple length 2 - 4 + +const PointTuple2 = Tuple{<:Real,<:Real} +const PointTuple3 = Tuple{<:Real,<:Real,<:Real} +const PointTuple4 = Tuple{<:Real,<:Real,<:Real,<:Real} +const PointTuple = Union{PointTuple2,PointTuple3,PointTuple4} GeoInterface.isgeometry(::Type{<:PointTuple}) = true GeoInterface.geomtrait(::PointTuple) = PointTrait() GeoInterface.ncoord(::PointTrait, geom::PointTuple) = Base.length(geom) GeoInterface.getcoord(::PointTrait, geom::PointTuple, i) = getindex(geom, i) +GeoInterface.is3d(::PointTrait, geom::PointTuple2) = false +GeoInterface.is3d(::PointTrait, geom::Union{PointTuple3,PointTuple4}) = true +GeoInterface.ismeasured(::PointTrait, geom::Union{PointTuple2,PointTuple3}) = false +GeoInterface.ismeasured(::PointTrait, geom::PointTuple4) = true + +GeoInterface.x(::PointTrait, geom::PointTuple) = geom[1] +GeoInterface.y(::PointTrait, geom::PointTuple) = geom[2] +GeoInterface.z(::PointTrait, geom::PointTuple2) = _z_error(Base.length(geom)) +GeoInterface.z(::PointTrait, geom::Union{PointTuple3,PointTuple4}) = geom[3] +GeoInterface.m(::PointTrait, geom::PointTuple4) = geom[4] +GeoInterface.m(::PointTrait, geom::Union{PointTuple2,PointTuple3}) = _m_error(Base.length(geom)) -for (i, pointtype) in enumerate((PointTuple2, PointTuple3, PointTuple4)) - keys = default_coord_names[1:i+1] - sig = NamedTuple{keys,<:pointtype} - @eval GeoInterface.isgeometry(::Type{<:$sig}) = true - @eval GeoInterface.geomtrait(::$sig) = PointTrait() - @eval GeoInterface.ncoord(::PointTrait, geom::$sig) = $i + 1 - @eval GeoInterface.getcoord(::PointTrait, geom::$sig, i) = getindex(geom, i) -end +# NamedTuple +# with X/Y/Z/M names in any order # Define all possible NamedTuple points -const NamedTuplePoint = Union{ +const NamedTuplePointXY = Union{ NamedTuple{(:X, :Y),<:PointTuple2}, NamedTuple{(:Y, :X),<:PointTuple2}, +} + +const NamedTuplePointZ = Union{ NamedTuple{(:X, :Y, :Z),<:PointTuple3}, NamedTuple{(:X, :Z, :Y),<:PointTuple3}, NamedTuple{(:Z, :Y, :X),<:PointTuple3}, NamedTuple{(:Z, :X, :Y),<:PointTuple3}, NamedTuple{(:Y, :X, :Z),<:PointTuple3}, NamedTuple{(:Y, :Z, :X),<:PointTuple3}, +} +const NamedTuplePointM = Union{ NamedTuple{(:X, :Y, :M),<:PointTuple3}, NamedTuple{(:X, :M, :Y),<:PointTuple3}, NamedTuple{(:M, :Y, :X),<:PointTuple3}, NamedTuple{(:M, :X, :Y),<:PointTuple3}, NamedTuple{(:Y, :X, :Z),<:PointTuple3}, NamedTuple{(:Y, :Z, :X),<:PointTuple3}, +} +const NamedTuplePointZM = Union{ NamedTuple{(:X, :Y, :Z, :M),<:PointTuple4}, NamedTuple{(:X, :Y, :M, :Z),<:PointTuple4}, NamedTuple{(:X, :Z, :Y, :M),<:PointTuple4}, @@ -68,16 +110,19 @@ const NamedTuplePoint = Union{ NamedTuple{(:M, :X, :Y, :Z),<:PointTuple4}, } +const NamedTuplePoint = Union{NamedTuplePointXY,NamedTuplePointZ,NamedTuplePointM,NamedTuplePointZM} + _keys(::Type{<:NamedTuple{K}}) where K = K -GeoInterface.isgeometry(::Type{T}) where {T<:NamedTuplePoint} = all(in(default_coord_names), _keys(T)) +GeoInterface.isgeometry(::Type{T}) where {T<:NamedTuplePoint} = true GeoInterface.geomtrait(::NamedTuplePoint) = PointTrait() GeoInterface.ncoord(::PointTrait, geom::NamedTuplePoint) = Base.length(geom) GeoInterface.getcoord(::PointTrait, geom::NamedTuplePoint, i) = getindex(geom, i) GeoInterface.coordnames(::PointTrait, geom::NamedTuplePoint) = _keys(typeof(geom)) GeoInterface.x(::PointTrait, geom::NamedTuplePoint) = geom.X GeoInterface.y(::PointTrait, geom::NamedTuplePoint) = geom.Y -GeoInterface.z(::PointTrait, geom::NamedTuplePoint) = geom.Z -GeoInterface.m(::PointTrait, geom::NamedTuplePoint) = geom.M +GeoInterface.z(::PointTrait, geom::Union{NamedTuplePointZ,NamedTuplePointZM}) = geom.Z +GeoInterface.z(::PointTrait, geom::NamedTuplePointXY) = throw(ArgumentError("NamedTuple point has no Z field")) +GeoInterface.m(::PointTrait, geom::Union{NamedTuplePointXY,NamedTuplePointZ}) = throw(ArgumentError("NamedTuple point has no M field")) # Default features using NamedTuple and AbstractArray diff --git a/test/test_primitives.jl b/test/test_primitives.jl index 8227ff99..ec007dd8 100644 --- a/test/test_primitives.jl +++ b/test/test_primitives.jl @@ -308,10 +308,41 @@ end @testset "Vector" begin geom = [1, 2] + @test !GeoInterface.is3d(geom) + @test !GeoInterface.ismeasured(geom) @test testgeometry(geom) - @test GeoInterface.x(geom) == 1 + @test collect(GeoInterface.getcoord(geom)) == geom @test GeoInterface.ncoord(geom) == 2 + @test GeoInterface.x(geom) == 1 + @test GeoInterface.y(geom) == 2 + @test_throws ArgumentError GeoInterface.z(geom) + @test_throws ArgumentError GeoInterface.m(geom) + geom = [1, 2, 3] + @test testgeometry(geom) @test collect(GeoInterface.getcoord(geom)) == geom + @test GeoInterface.is3d(geom) + @test !GeoInterface.ismeasured(geom) + @test GeoInterface.x(geom) == 1 + @test GeoInterface.y(geom) == 2 + @test GeoInterface.z(geom) == 3 + @test_throws ArgumentError GeoInterface.m(geom) + geom = [1, 2, 3, 4] + @test testgeometry(geom) + @test collect(GeoInterface.getcoord(geom)) == geom + @test GeoInterface.is3d(geom) + @test GeoInterface.ismeasured(geom) + @test GeoInterface.ncoord(geom) == 4 + @test GeoInterface.x(geom) == 1 + @test GeoInterface.y(geom) == 2 + @test GeoInterface.z(geom) == 3 + @test GeoInterface.m(geom) == 4 + geom = [1, 2, 3, 4, 5] + @test !GeoInterface.is3d(geom) + @test !GeoInterface.ismeasured(geom) + @test_throws ArgumentError GeoInterface.x(geom) + @test_throws ArgumentError GeoInterface.y(geom) + @test_throws ArgumentError GeoInterface.z(geom) + @test_throws ArgumentError GeoInterface.m(geom) end @testset "Tuple" begin @@ -319,10 +350,42 @@ end @test GeoInterface.trait(geom) isa PointTrait @test GeoInterface.geomtrait(geom) isa PointTrait @test testgeometry(geom) - @test GeoInterface.x(geom) == 1 - @test GeoInterface.y(geom) == 2.0f0 + @test GeoInterface.x(geom) === 1 + @test GeoInterface.y(geom) === 2.0f0 + @test_throws ArgumentError GeoInterface.z(geom) + @test_throws ArgumentError GeoInterface.m(geom) @test GeoInterface.ncoord(geom) == 2 @test collect(GeoInterface.getcoord(geom)) == [1, 2] + geom = (1, 2, 3.0) + @test GeoInterface.trait(geom) isa PointTrait + @test GeoInterface.geomtrait(geom) isa PointTrait + @test testgeometry(geom) + @test GeoInterface.x(geom) === 1 + @test GeoInterface.y(geom) === 2 + @test GeoInterface.z(geom) === 3.0 + @test_throws ArgumentError GeoInterface.m(geom) + @test GeoInterface.ncoord(geom) == 3 + @test collect(GeoInterface.getcoord(geom)) == [1, 2, 3] + geom = (1, 2, 3, 4.0) + @test GeoInterface.trait(geom) isa PointTrait + @test GeoInterface.geomtrait(geom) isa PointTrait + @test testgeometry(geom) + @test GeoInterface.x(geom) === 1 + @test GeoInterface.y(geom) === 2 + @test GeoInterface.z(geom) === 3 + @test GeoInterface.m(geom) === 4.0 + @test GeoInterface.ncoord(geom) == 4 + @test collect(GeoInterface.getcoord(geom)) == [1, 2, 3, 4] + geom = (1, 2, 3, 4.0, 5) + @test GeoInterface.isgeometry(geom) == false + @test GeoInterface.trait(geom) isa Nothing + @test GeoInterface.geomtrait(geom) isa Nothing + @test_throws MethodError GeoInterface.x(geom) + @test_throws MethodError GeoInterface.y(geom) + @test_throws MethodError GeoInterface.z(geom) + @test_throws MethodError GeoInterface.m(geom) + @test_throws MethodError GeoInterface.ncoord(geom) == 4 + @test_throws MethodError GeoInterface.getcoord(geom) end @testset "NamedTuple" begin @@ -334,8 +397,20 @@ end @test testgeometry(geom) @test GeoInterface.x(geom) == 1 @test GeoInterface.y(geom) == 2.0 + @test_throws ArgumentError GeoInterface.z(geom) + @test_throws ArgumentError GeoInterface.m(geom) @test collect(GeoInterface.getcoord(geom)) == [1, 2] + geom = (; X=1.0, Y=2.0, Z=0x03) + @test testgeometry(geom) + @test GeoInterface.is3d(geom) + @test GeoInterface.ismeasured(geom) == false + @test GeoInterface.coordnames(geom) == (:X, :Y, :Z) + @test GeoInterface.x(geom) == 1 + @test GeoInterface.y(geom) == 2.0 + @test GeoInterface.z(geom) == 0x03 + @test_throws ArgumentError GeoInterface.m(geom) + geom = (; X=1, Y=2, Z=3.0f0, M=4.0) @test GeoInterface.trait(geom) isa PointTrait @test GeoInterface.geomtrait(geom) isa PointTrait @@ -355,15 +430,20 @@ end @test GeoInterface.getcoord(geom, 4) === 4.0 @test collect(GeoInterface.getcoord(geom)) == [1, 2, 3, 4] - geom = (; X=1.0, Y=2.0, Z=0x03) - @test testgeometry(geom) - @test GeoInterface.is3d(geom) - @test GeoInterface.ismeasured(geom) == false - @test GeoInterface.coordnames(geom) == (:X, :Y, :Z) geom = (; Z=3, X=1, Y=2, M=4) @test testgeometry(geom) @test collect(GeoInterface.getcoord(geom)) == [3, 1, 2, 4] @test GeoInterface.coordnames(geom) == (:Z, :X, :Y, :M) + + geom = (; A=0, X=1, Y=2, Z=3.0f0, M=4.0) + @test GeoInterface.trait(geom) isa Nothing + @test GeoInterface.geomtrait(geom) isa Nothing + @test_throws MethodError GeoInterface.x(geom) + @test_throws MethodError GeoInterface.y(geom) + @test_throws MethodError GeoInterface.z(geom) + @test_throws MethodError GeoInterface.m(geom) + @test_throws MethodError GeoInterface.ncoord(geom) == 4 + @test_throws MethodError GeoInterface.getcoord(geom) end @testset "NamedTupleFeature" begin From 6abf141c28f017d7beb0d6e3283f6be91747e1d7 Mon Sep 17 00:00:00 2001 From: rafaqz Date: Tue, 17 Jan 2023 20:58:56 +0000 Subject: [PATCH 2/3] more clearly document x/y/z/m for points --- src/interface.jl | 14 +++++++++++++- 1 file changed, 13 insertions(+), 1 deletion(-) diff --git a/src/interface.jl b/src/interface.jl index cd495e2c..633cff4c 100644 --- a/src/interface.jl +++ b/src/interface.jl @@ -561,6 +561,8 @@ convexhull(geom) = convexhull(geomtrait(geom), geom) Return the :X coordinate of the given `geom`. Note that this is only valid for [`AbstractPointTrait`](@ref)s. + +For `Tuple` and `Vector` points, the first value is returned. """ x(geom) = x(geomtrait(geom), geom) @@ -569,6 +571,8 @@ x(geom) = x(geomtrait(geom), geom) Return the :Y coordinate of the given `geom`. Note that this is only valid for [`AbstractPointTrait`](@ref)s. + +For `Tuple` and `Vector` points, the second value is returned. """ y(geom) = y(geomtrait(geom), geom) @@ -577,14 +581,22 @@ y(geom) = y(geomtrait(geom), geom) Return the :Z coordinate of the given `geom`. Note that this is only valid for [`AbstractPointTrait`](@ref)s. + +For length 3 `Tuple` and `Vector` points, the third value is returned. """ z(geom) = z(geomtrait(geom), geom) """ m(geom) -> Number -Return the :M coordinate of the given `geom`. +Return the :M (measured) coordinate of the given `geom`. Note that this is only valid for [`AbstractPointTrait`](@ref)s. + +For length 4 `Tuple` and `Vector` points, the fouth value +is returned. + +Length 3 `Tuple` and `Vector` points can *not* represent measured points, +and will throw an `ArgumentError`. """ m(geom) = m(geomtrait(geom), geom) From a40dccda9f1ec10500eeed7ce2e840df257d0756 Mon Sep 17 00:00:00 2001 From: rafaqz Date: Tue, 17 Jan 2023 21:04:00 +0000 Subject: [PATCH 3/3] test ismeasured and is3d for Tuple --- test/test_primitives.jl | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/test/test_primitives.jl b/test/test_primitives.jl index ec007dd8..c328f00f 100644 --- a/test/test_primitives.jl +++ b/test/test_primitives.jl @@ -350,6 +350,8 @@ end @test GeoInterface.trait(geom) isa PointTrait @test GeoInterface.geomtrait(geom) isa PointTrait @test testgeometry(geom) + @test !GeoInterface.is3d(geom) + @test !GeoInterface.ismeasured(geom) @test GeoInterface.x(geom) === 1 @test GeoInterface.y(geom) === 2.0f0 @test_throws ArgumentError GeoInterface.z(geom) @@ -360,6 +362,8 @@ end @test GeoInterface.trait(geom) isa PointTrait @test GeoInterface.geomtrait(geom) isa PointTrait @test testgeometry(geom) + @test GeoInterface.is3d(geom) + @test !GeoInterface.ismeasured(geom) @test GeoInterface.x(geom) === 1 @test GeoInterface.y(geom) === 2 @test GeoInterface.z(geom) === 3.0 @@ -370,6 +374,8 @@ end @test GeoInterface.trait(geom) isa PointTrait @test GeoInterface.geomtrait(geom) isa PointTrait @test testgeometry(geom) + @test GeoInterface.is3d(geom) + @test GeoInterface.ismeasured(geom) @test GeoInterface.x(geom) === 1 @test GeoInterface.y(geom) === 2 @test GeoInterface.z(geom) === 3