diff --git a/benchmarks/benchmark_plots.jl b/benchmarks/benchmark_plots.jl index 1819116f1..e9255120c 100644 --- a/benchmarks/benchmark_plots.jl +++ b/benchmarks/benchmark_plots.jl @@ -89,7 +89,6 @@ function plot_trials( end tag_attrs = capture_tag_attrs(results.tags) - tag_theme = merge(theme, tag_attrs) lp = if legend_position isa Makie.Automatic gp.layout[gp.span.rows, gp.span.cols, TopRight()] @@ -101,12 +100,15 @@ function plot_trials( error() end - return Makie.with_theme(tag_theme) do + ax = Makie.with_theme(theme) do ax = Axis( gp; tag_attrs.Axis..., xlabel = "Number of points", ylabel = "Time to calculate", xscale = log10, yscale = log10, ytickformat = _prettytime, + xticksvisible = true, xticklabelsvisible = true, + yticks = Makie.LogTicks(Makie.WilkinsonTicks(7; k_min = 4)), + ygridwidth = 0.75, ) plots = [scatterlines!(ax, x, y; label = label) for (x, y, label) in zip(xs, ys, labels)] setproperty!.(getindex.(getproperty.(plots, :plots), 1), :alpha, 0.1) @@ -118,13 +120,11 @@ function plot_trials( valign = legend_valign, orientation = legend_orientation ) - ax.xticksvisible[] = true ax.xtickcolor[] = ax.xgridcolor[] - ax.xticklabelsvisible[] = true - ax.yticks[] = Makie.LogTicks(Makie.WilkinsonTicks(7; k_min = 4)) - ax.ygridwidth[] = 0.75 - return ax + ax end + + return ax end const _tag_includelist = ["title", "subtitle"] diff --git a/benchmarks/polygon_scaling.jl b/benchmarks/polygon_scaling.jl new file mode 100644 index 000000000..f38cfcc23 --- /dev/null +++ b/benchmarks/polygon_scaling.jl @@ -0,0 +1,19 @@ +using CairoMakie, Chairmarks, BenchmarkTools +import GeometryOps as GO, GeoInterface as GI, LibGEOS as LG +include("benchmark_plots.jl") + +p25 = GI.Polygon([[(0.0, 0.0), (5.0, 0.0), (5.0, 8.0), (0.0, 8.0), (0.0, 0.0)], [(4.0, 0.5), (4.5, 0.5), (4.5, 3.5), (4.0, 3.5), (4.0, 0.5)], [(2.0, 4.0), (4.0, 4.0), (4.0, 6.0), (2.0, 6.0), (2.0, 4.0)]]) +p26 = GI.Polygon([[(3.0, 1.0), (8.0, 1.0), (8.0, 7.0), (3.0, 7.0), (3.0, 5.0), (6.0, 5.0), (6.0, 3.0), (3.0, 3.0), (3.0, 1.0)], [(3.5, 5.5), (6.0, 5.5), (6.0, 6.5), (3.5, 6.5), (3.5, 5.5)], [(5.5, 1.5), (5.5, 2.5), (3.5, 2.5), (3.5, 1.5), (5.5, 1.5)]]) + +suite = BenchmarkGroup(["title:Polygon intersection timing","subtitle:Single polygon, densified"]) + +for max_distance in exp10.(LinRange(-1, 1.5, 10)) + p25s = GO.segmentize(p25; max_distance) + p26s = GO.segmentize(p26; max_distance) + n_verts = GI.npoint(p25s) + suite["GeometryOps"][n_verts] = @be GO.intersection($p25s, $p26s; target = $(GI.PolygonTrait()), fix_multipoly = $nothing) + suite["LibGEOS"][n_verts] = @be LG.intersection($(GI.convert(LG, p25s)), $(GI.convert(LG, p26s))) +end + + +plot_trials(suite) \ No newline at end of file diff --git a/src/GeometryOps.jl b/src/GeometryOps.jl index 80f50b5f0..a4a774e38 100644 --- a/src/GeometryOps.jl +++ b/src/GeometryOps.jl @@ -5,8 +5,9 @@ module GeometryOps using GeoInterface using GeometryBasics import Tables -using ExactPredicates, LinearAlgebra, Statistics +using LinearAlgebra, Statistics import GeometryBasics.StaticArrays +import ExactPredicates import Base.@kwdef using GeoInterface.Extents: Extents diff --git a/src/methods/clipping/intersection.jl b/src/methods/clipping/intersection.jl index 0fe302574..e5ffd18ce 100644 --- a/src/methods/clipping/intersection.jl +++ b/src/methods/clipping/intersection.jl @@ -263,8 +263,8 @@ function _intersection_point(::Type{T}, (a1, a2)::Edge, (b1, b2)::Edge) where T r_cross_s = rx * sy - ry * sx if Predicates.isparallel((rx, ry), (sx, sy)) != 0 # non-parallel lines # Calculate α ratio if lines cross or touch - a1_orient = ExactPredicates.orient(b1, b2, a1) - a2_orient = ExactPredicates.orient(b1, b2, a2) + a1_orient = Predicates.orient(b1, b2, a1) + a2_orient = Predicates.orient(b1, b2, a2) α = if a1_orient == 0 # α = 0 zero(T) elseif a2_orient == 0 # α = 1 @@ -276,8 +276,8 @@ function _intersection_point(::Type{T}, (a1, a2)::Edge, (b1, b2)::Edge) where T return line_orient, intr1, intr2 end # Calculate β ratio if lines touch or cross - b1_orient = ExactPredicates.orient(a1, a2, b1) - b2_orient = ExactPredicates.orient(a1, a2, b2) + b1_orient = Predicates.orient(a1, a2, b1) + b2_orient = Predicates.orient(a1, a2, b2) β = if b1_orient == 0 # β = 0 zero(T) elseif b2_orient == 0 # β = 1 @@ -295,10 +295,10 @@ function _intersection_point(::Type{T}, (a1, a2)::Edge, (b1, b2)::Edge) where T line_orient = (α == 0 || α == 1 || β == 0 || β == 1) ? line_hinge : line_cross elseif Predicates.iscollinear((Δqp_x, Δqp_y), (sx, sy)) == 0 # collinear parallel lines # Determine if lines touch or overlap and with what α and β values - a1_side = ExactPredicates.sameside(a1, b1, b2) - a2_side = ExactPredicates.sameside(a2, b1, b2) - b1_side = ExactPredicates.sameside(b1, a1, a2) - b2_side = ExactPredicates.sameside(b2, a1, a2) + a1_side = Predicates.sameside(a1, b1, b2) + a2_side = Predicates.sameside(a2, b1, b2) + b1_side = Predicates.sameside(b1, a1, a2) + b2_side = Predicates.sameside(b2, a1, a2) # Lines touch or overlap if endpoints of line a are on/in line b and visa versa r_dot_s = rx * sx + ry * sy # Determine which endpoints start and end the overlapping region diff --git a/src/methods/clipping/predicates.jl b/src/methods/clipping/predicates.jl index f038e247d..41dcd4c3a 100644 --- a/src/methods/clipping/predicates.jl +++ b/src/methods/clipping/predicates.jl @@ -6,6 +6,9 @@ module Predicates using ExactPredicates: ext, inp + orient(args...) = ExactPredicates.orient(args...) + sameside(args...) = ExactPredicates.sameside(args...) + # This is the implementation of r_cross_s from 2 points in the `_intersection_points` file. # 0 == parallel function isparallel(a1, a2, b1, b2)