add macro to test all implementations

src/methods/clipping/clipping_processor.jl

@@ -643,6 +643,7 @@ function _add_holes_to_polys!(::Type{T}, return_polys, hole_iterator, remove_pol
     for i in 1:n_polys
         n_new_per_poly = 0
         for curr_hole in Iterators.map(tuples, hole_iterator) # loop through all holes
+            curr_hole = _linearring(curr_hole)
             # loop through all pieces of original polygon (new pieces added to end of list)
             for j in Iterators.flatten((i:i, (n_polys + 1):(n_polys + n_new_per_poly)))
                 curr_poly = return_polys[j]

src/methods/clipping/union.jl

@@ -62,9 +62,9 @@ function _union(
     if n_pieces == 0 # no crossing points, determine if either poly is inside the other
         a_in_b, b_in_a = _find_non_cross_orientation(a_list, b_list, ext_a, ext_b; exact)
         if a_in_b
-            push!(polys, GI.Polygon([tuples(ext_b)]))
+            push!(polys, GI.Polygon([_linearring(tuples(ext_b))]))
         elseif b_in_a
-            push!(polys,  GI.Polygon([tuples(ext_a)]))
+            push!(polys,  GI.Polygon([_linearring(tuples(ext_a))]))
         else
             push!(polys, tuples(poly_a))
             push!(polys, tuples(poly_b))
@@ -124,6 +124,7 @@ function _add_union_holes!(polys, a_in_b, b_in_a, poly_a, poly_b; exact)
         current_poly = n_a_holes > 0 ? ext_poly_b : poly_a
         # Loop over all holes in both original polygons
         for (i, ih) in enumerate(Iterators.flatten((GI.gethole(poly_a), GI.gethole(poly_b))))
+            ih = _linearring(ih)
             in_ext, _, _ = _line_polygon_interactions(ih, curr_exterior_poly; exact, closed_line = true)
             if !in_ext
                 #= if the hole isn't in the overlapping region between the two polygons, add

src/primitives.jl

@@ -260,6 +260,19 @@ end
     geoms = _maptasks(apply_to_geom, 1:GI.ngeom(geom), threaded)
     return _apply_inner(geom, geoms, crs, calc_extent)
 end
+@inline function _apply(f::F, target::TraitTarget{<:PointTrait}, trait::GI.PolygonTrait, geom;
+    crs=GI.crs(geom), calc_extent=_False(), threaded
+)::(GI.geointerface_geomtype(trait)) where F
+    # We need to force rebuilding a LinearRing not a LineString
+    geoms = _maptasks(1:GI.ngeom(geom), threaded) do i
+        lr = GI.getgeom(geom, i)
+        points = map(GI.getgeom(lr)) do p
+            _apply(f, target, p; crs, calc_extent, threaded=_False())
+        end
+        _linearring(_apply_inner(lr, points, crs, calc_extent))
+    end
+    return _apply_inner(geom, geoms, crs, calc_extent)
+end
 function _apply_inner(geom, geoms, crs, calc_extent::_True)
     # Calculate the extent of the sub geometries
     extent = mapreduce(GI.extent, Extents.union, geoms)

src/transformations/correction/intersecting_polygons.jl

@@ -141,4 +141,4 @@ function (::DiffIntersectingPolygons)(::GI.MultiPolygonTrait, multipoly)
         keepat!(diff_multipoly.geom, keep_idx)
     end
     return diff_multipoly
-end
+end

src/utils.jl

@@ -126,3 +126,6 @@ function _point_in_extent(p, extent::Extents.Extent)
     (x1, x2), (y1, y2) = extent.X, extent.Y
     return x1 ≤ GI.x(p) ≤ x2 && y1 ≤ GI.y(p) ≤ y2
 end
+
+_linearring(geom::GI.LineString) = GI.LinearRing(parent(geom); extent=geom.extent, crs=geom.crs)
+_linearring(geom::GI.LinearRing) = geom

test/extensions/flexijoins.jl

@@ -1,22 +1,23 @@
-import GeometryOps as GO, GeoInterface as GI
-using FlexiJoins, DataFrames
+using Test
+using FlexiJoins
+using DataFrames
+import GeometryOps as GO
+import GeoInterface as GI
+using ..TestHelpers
+
+points = GI.MultiPoint(tuple.(rand(100), rand(100)))
 
 pl = GI.Polygon([GI.LinearRing([(0, 0), (1, 0), (1, 1), (0, 0)])])
 pu = GI.Polygon([GI.LinearRing([(0, 0), (0, 1), (1, 1), (0, 0)])])
-poly_df = DataFrame(geometry = [pl, pu], color = [:red, :blue])
-
-points = tuple.(rand(100), rand(100))
-points_df = DataFrame(geometry = points)
 
-@testset "Basic integration" begin
-
-    @test_nowarn joined_df = FlexiJoins.innerjoin((poly_df, points_df), by_pred(:geometry, GO.contains, :geometry))
+@testset_implementations "Polygon DataDrame" begin
+    points_df = DataFrame(geometry=collect(GI.getpoint($points)))
+    poly_df = DataFrame(geometry=[$pl, $pu], color=[:red, :blue])
     # Test that the join happened correctly
     joined_df = FlexiJoins.innerjoin((poly_df, points_df), by_pred(:geometry, GO.contains, :geometry))
-    @test all(GO.contains.((pl,), joined_df.geometry_1[joined_df.color .== :red]))
-    @test all(GO.contains.((pu,), joined_df.geometry_1[joined_df.color .== :blue]))
+    @test all(GO.contains.(($pl,), joined_df.geometry_1[joined_df.color .== :red]))
+    @test all(GO.contains.(($pu,), joined_df.geometry_1[joined_df.color .== :blue]))
     # Test that within also works
     @test_nowarn joined_df = FlexiJoins.innerjoin((points_df, poly_df), by_pred(:geometry, GO.within, :geometry))
-
 end
 

test/helpers.jl

@@ -0,0 +1,119 @@
+module TestHelpers
+
+using Test, GeoInterface, ArchGDAL, GeometryBasics, LibGEOS
+
+export @test_implementations, @testset_implementations
+
+const TEST_MODULES = [GeoInterface, ArchGDAL, GeometryBasics, LibGEOS]
+
+# Monkey-patch GeometryBasics to have correct methods.
+# TODO: push this up to GB!
+
+@eval GeometryBasics begin
+    # MultiGeometry ncoord implementations
+    GeoInterface.ncoord(::GeoInterface.MultiPolygonTrait, ::GeometryBasics.MultiPolygon{N}) where N = N
+    GeoInterface.ncoord(::GeoInterface.MultiLineStringTrait, ::GeometryBasics.MultiLineString{N}) where N = N
+    GeoInterface.ncoord(::GeoInterface.MultiPointTrait, ::GeometryBasics.MultiPoint{N}) where N = N
+    # LinearRing and LineString confusion
+    GeometryBasics.geointerface_geomtype(::GeoInterface.LinearRingTrait) = LineString
+    function GeoInterface.convert(::Type{LineString}, ::GeoInterface.LinearRingTrait, geom)
+        return GeoInterface.convert(LineString, GeoInterface.LineStringTrait(), geom) # forward to the linestring conversion method
+    end
+    # Line interface
+    GeometryBasics.geointerface_geomtype(::GeoInterface.LineTrait) = Line
+    function GeoInterface.convert(::Type{Line}, ::GeoInterface.LineTrait, geom)
+        p1, p2 = GeoInterface.getpoint(geom)
+        return Line(GeoInterface.convert(Point, GeoInterface.PointTrait(), p1), GeoInterface.convert(Point, GeoInterface.PointTrait(), p2))
+    end
+    # GeometryCollection interface - currently just a large Union
+    const _ALL_GB_GEOM_TYPES = Union{Point, Line, LineString, Polygon, MultiPolygon, MultiLineString, MultiPoint}
+    GeometryBasics.geointerface_geomtype(::GeoInterface.GeometryCollectionTrait) = Vector{_ALL_GB_GEOM_TYPES}
+    function GeoInterface.convert(::Type{Vector{_ALL_GB_GEOM_TYPES}}, ::GeoInterface.GeometryCollectionTrait, geoms)
+        return _ALL_GB_GEOM_TYPES[GeoInterface.convert(GeometryBasics, g) for g in GeoInterface.getgeom(geoms)]
+    end
+end
+
+
+# Macro to run a block of `code` for multiple modules,
+# using GeoInterface.convert for each var in `args`
+macro test_implementations(code::Expr)
+    _test_implementations_inner(TEST_MODULES, code)
+end
+macro test_implementations(modules::Union{Expr,Vector}, code::Expr)
+    _test_implementations_inner(modules, code)
+end
+
+function _test_implementations_inner(modules::Union{Expr,Vector}, code::Expr)
+    vars = Dict{Symbol,Symbol}()
+    code1 = _quasiquote!(code, vars)
+    modules1 = modules isa Expr ? modules.args : modules
+    tests = Expr(:block)
+
+    for mod in modules1
+        expr = Expr(:block)
+        for (var, genkey) in pairs(vars)
+            push!(expr.args, :($genkey = $GeoInterface.convert($mod, $var)))
+        end
+        push!(expr.args, :(@test $code1))
+        push!(tests.args, expr)
+    end
+
+    return esc(tests)
+end
+
+# Macro to run a block of `code` for multiple modules,
+# using GeoInterface.convert for each var in `args`
+macro testset_implementations(code::Expr)
+    _testset_implementations_inner("", TEST_MODULES, code)
+end
+macro testset_implementations(arg, code::Expr)
+    if arg isa String || arg isa Expr && arg.head == :string
+        _testset_implementations_inner(arg, TEST_MODULES, code)
+    else
+        _testset_implementations_inner("", arg, code)
+    end
+end
+macro testset_implementations(title, modules::Union{Expr,Vector}, code::Expr)
+    _testset_implementations_inner(title, modules, code)
+end
+
+function _testset_implementations_inner(title, modules::Union{Expr,Vector}, code::Expr)
+    vars = Dict{Symbol,Symbol}()
+    code1 = _quasiquote!(code, vars)
+    modules1 = modules isa Expr ? modules.args : modules
+    testsets = Expr(:block)
+
+    for mod in modules1
+        expr = Expr(:block)
+        for (var, genkey) in pairs(vars)
+            push!(expr.args, :($genkey = $GeoInterface.convert($mod, $var)))
+        end
+        label = title == "" ? "$mod" : "$title: $mod"
+        push!(expr.args, :(@testset $label $code1))
+        push!(testsets.args, expr)
+    end
+
+    return esc(testsets)
+end
+
+# Taken from BenchmarkTools.jl
+_quasiquote!(ex, vars) = ex
+function _quasiquote!(ex::Expr, vars::Dict)
+    if ex.head === :($)
+        v = ex.args[1]
+        gen = if v isa Symbol 
+            haskey(vars, v) ? vars[v] : gensym(v)
+        else
+            gensym()
+        end
+        vars[v] = gen
+        return v 
+    elseif ex.head !== :quote
+        for i in 1:length(ex.args)
+            ex.args[i] = _quasiquote!(ex.args[i], vars)
+        end
+    end
+    return ex
+end
+
+end

test/methods/angles.jl

@@ -1,5 +1,9 @@
 using Test
-import GeoInterface as GI, GeometryOps as GO, LibGEOS as LG
+import GeoInterface as GI
+import GeometryOps as GO 
+import GeometryBasics as GB 
+import LibGEOS as LG
+using ..TestHelpers
 
 pt1 = GI.Point((0.0, 0.0))
 mpt1 = GI.MultiPoint([pt1, pt1])
@@ -9,39 +13,43 @@ concave_coords = [(0.0, 0.0), (0.0, 1.0), (-1.0, 1.0), (-1.0, 2.0), (2.0, 2.0),
 l2 = GI.LineString(concave_coords)
 l3 = GI.LineString(concave_coords[1:(end - 1)])
 r1 = GI.LinearRing(concave_coords)
-r2 = GI.LinearRing(concave_coords[1:(end - 1)])
-r3 = GI.LinearRing([(1.0, 1.0), (1.0, 1.5), (1.5, 1.5), (1.5, 1.0), (1.0, 1.0)])
+r2 = GI.LinearRing([(1.0, 1.0), (1.0, 1.5), (1.5, 1.5), (1.5, 1.0), (1.0, 1.0)])
 concave_angles = [90.0, 270.0, 90.0, 90.0, 90.0, 90.0]
 
-p1 = GI.Polygon([r3])
+p1 = GI.Polygon([r2])
 p2 = GI.Polygon([[(0.0, 0.0), (0.0, 4.0), (3.0, 0.0), (0.0, 0.0)]])
 p3 = GI.Polygon([[(-3.0, -2.0), (0.0,0.0), (5.0, 0.0), (-3.0, -2.0)]])
 p4 = GI.Polygon([r1])
-p5 = GI.Polygon([r1, r3])
+p5 = GI.Polygon([r1, r2])
 
 mp1 = GI.MultiPolygon([p2, p3])
 c1 = GI.GeometryCollection([pt1, l2, p2])
 
-# Points and lines
-@test isempty(GO.angles(pt1))
-@test isempty(GO.angles(mpt1))
-@test isempty(GO.angles(l1))
-
-# LineStrings and Linear Rings
-@test all(isapprox.(GO.angles(l2), concave_angles, atol = 1e-3))
-@test all(isapprox.(GO.angles(l3), concave_angles[2:(end - 1)], atol = 1e-3))
-@test all(isapprox.(GO.angles(r1), concave_angles, atol = 1e-3))
-@test all(isapprox.(GO.angles(r2), concave_angles, atol = 1e-3))
-
-# Polygons
-p2_angles = [90.0, 36.8699, 53.1301]
-p3_angles = [19.6538, 146.3099, 14.0362]
-@test all(isapprox.(GO.angles(p1), [90.0 for _ in 1:4], atol = 1e-3))
-@test all(isapprox.(GO.angles(p2), p2_angles, atol = 1e-3))
-@test all(isapprox.(GO.angles(p3), p3_angles, atol = 1e-3))
-@test all(isapprox.(GO.angles(p4), concave_angles, atol = 1e-3))
-@test all(isapprox.(GO.angles(p5), vcat(concave_angles, [270.0 for _ in 1:4]), atol = 1e-3))
-
-# Multi-geometries
-@test all(isapprox.(GO.angles(mp1), [p2_angles; p3_angles], atol = 1e-3))
-@test all(isapprox.(GO.angles(c1), [concave_angles; p2_angles], atol = 1e-3))
+# Line is not a widely available geometry type
+@testset_implementations "line angles" [GB, GI] begin
+    @test isempty(GO.angles($l1))
+end
+
+@testset_implementations "angles" begin
+    # Points and lines
+    @test isempty(GO.angles($pt1))
+    @test isempty(GO.angles($mpt1))
+
+    # LineStrings and Linear Rings
+    @test all(isapprox.(GO.angles($l2), concave_angles, atol = 1e-3))
+    @test all(isapprox.(GO.angles($l3), concave_angles[2:(end - 1)], atol = 1e-3))
+    @test all(isapprox.(GO.angles($r1), concave_angles, atol = 1e-3))
+
+    # Polygons
+    p2_angles = [90.0, 36.8699, 53.1301]
+    p3_angles = [19.6538, 146.3099, 14.0362]
+    @test all(isapprox.(GO.angles($p1), [90.0 for _ in 1:4], atol = 1e-3))
+    @test all(isapprox.(GO.angles($p2), p2_angles, atol = 1e-3))
+    @test all(isapprox.(GO.angles($p3), p3_angles, atol = 1e-3))
+    @test all(isapprox.(GO.angles($p4), concave_angles, atol = 1e-3))
+    @test all(isapprox.(GO.angles($p5), vcat(concave_angles, [270.0 for _ in 1:4]), atol = 1e-3))
+
+    # Multi-geometries
+    @test all(isapprox.(GO.angles($mp1), [p2_angles; p3_angles], atol = 1e-3))
+    @test all(isapprox.(GO.angles(c1), [concave_angles; p2_angles], atol = 1e-3))
+end