[wip] antimeridian cutting with trees for regridding

[asinghvi17]

implement an antimeridian cutting fix ported from the python package

use traittarget in application_level for fix

Merge branch 'as/antimeridian' into as/trees-and-antimeridian

[asinghvi17]

Warning

This pull request is not mergeable via GitHub because a downstack PR is open. Once all requirements are satisfied, merge this PR as a stack on Graphite.
Learn more

• [wip] antimeridian cutting with trees for regridding #283 👈 (View in Graphite)
• Tree based acceleration for polygon clipping / boolean ops #274 : 1 other dependent PR (#278 )
• Algorithm + manifold and accelerator passthrough in clipping functions #273
• Spatial tree interface #297
• Convert the orient predicate to AdaptivePredicates + minor predicate changes #275
• Add a TGGeometry extension for fast inexact geometric predicates #271 : 1 other dependent PR (#291 )
• [Core] Add an algorithm interface + StableTasks #269 : 1 other dependent PR (#282 )
• main

This stack of pull requests is managed by Graphite. Learn more about stacking.

[asinghvi17]

This is what the method will ultimately look like...it takes me ~3.5s on my machine, single threaded, to cross reference two ~80,000 element grids.

It's possible to multithread this but not the easiest thing...might need some kind of master-worker setup. We need to multithread within the dual tree query, so maybe we let the query run on one thread and spawn nthreads()*2-2 worker threads that can consume the output from a channel. There would likely also need to be a writer task that is single-threadedly writing to A.

[milankl]

I'm getting the following error when I try to install this branch

julia> using GeometryOps
Precompiling GeometryOps...
Info Given GeometryOps was explicitly requested, output will be shown live
WARNING: could not import GeometryOpsCore.AutoManifold into GeometryOps
WARNING: could not import GeometryOpsCore.WrongManifoldException into GeometryOps
WARNING: could not import GeometryOpsCore.Algorithm into GeometryOps
WARNING: could not import GeometryOpsCore.AutoAlgorithm into GeometryOps
WARNING: could not import GeometryOpsCore.ManifoldIndependentAlgorithm into GeometryOps
WARNING: could not import GeometryOpsCore.SingleManifoldAlgorithm into GeometryOps
WARNING: could not import GeometryOpsCore.NoAlgorithm into GeometryOps
ERROR: LoadError: UndefVarError: `SingleManifoldAlgorithm` not defined in `GeometryOpsCore`yOps
Stacktrace:
 [1] getproperty(x::Module, f::Symbol)
   @ Base ./Base.jl:42
 [2] top-level scope
   @ ~/.julia/packages/GeometryOps/4octE/src/types.jl:34
 [3] include(mod::Module, _path::String)
   @ Base ./Base.jl:557
 [4] include(x::String)
   @ GeometryOps ~/.julia/packages/GeometryOps/4octE/src/GeometryOps.jl:3
 [5] top-level scope
   @ ~/.julia/packages/GeometryOps/4octE/src/GeometryOps.jl:37
 [6] include
   @ ./Base.jl:557 [inlined]
 [7] include_package_for_output(pkg::Base.PkgId, input::String, depot_path::Vector{String}, dl_load_path::Vector{String}, load_path::Vector{String}, concrete_deps::Vector{Pair{Base.PkgId, UInt128}}, source::Nothing)
   @ Base ./loading.jl:2881
 [8] top-level scope
   @ stdin:6
in expression starting at /Users/milan/.julia/packages/GeometryOps/4octE/src/types.jl:34
in expression starting at /Users/milan/.julia/packages/GeometryOps/4octE/src/GeometryOps.jl:3
in expression starting at stdin:6
  ✗ GeometryOps
  0 dependencies successfully precompiled in 3 seconds. 56 already precompiled.

ERROR: The following 1 direct dependency failed to precompile:

GeometryOps

[asinghvi17]

Ah, you need to add GeometryOpsCore on this branch as well

If you have it dev'ed, dev ~/.julia/dev/GeometryOps/GeometryOpsCore should do it

[milankl]

This seems to work 🥳 @simone-silvestri

and quick check whether it's conserving

looks pretty good I would say!!!

[milankl]

With some combinations of grids and resolutions, e.g.

grid1 = randn(OctaHEALPixGrid, 24)
grid2 = randn(OctaminimalGaussianGrid, 12)

I'm hitting an error though

┌ Error: Intersection failed!
│   i1 = 1152
│   i2 = 289
└ @ Main In[135]:25

BoundsError: attempt to access 4-element Vector{GeometryOps.PolyNode{Float64}} at index [6]

Stacktrace:
  [1] throw_boundserror(A::Vector{GeometryOps.PolyNode{Float64}}, I::Tuple{Int64})
    @ Base ./essentials.jl:14
  [2] getindex
    @ ./essentials.jl:916 [inlined]
  [3] _classify_crossing!(alg::GeometryOps.FosterHormannClipping{GeometryOpsCore.Planar, GeometryOps.AutoAccelerator}, ::Type{Float64}, a_list::Vector{GeometryOps.PolyNode{Float64}}, b_list::Vector{GeometryOps.PolyNode{Float64}}; exact::GeometryOpsCore.True)
    @ GeometryOps ~/.julia/packages/GeometryOps/4octE/src/methods/clipping/clipping_processor.jl:623
  [4] _classify_crossing!
    @ ~/.julia/packages/GeometryOps/4octE/src/methods/clipping/clipping_processor.jl:606 [inlined]
  [5] _build_ab_list(alg::GeometryOps.FosterHormannClipping{GeometryOpsCore.Planar, GeometryOps.AutoAccelerator}, ::Type{Float64}, poly_a::GeoInterface.Wrappers.LinearRing{false, false, Vector{Tuple{Float64, Float64}}, Nothing, Nothing}, poly_b::GeoInterface.Wrappers.LinearRing{false, false, Vector{Tuple{Float64, Float64}}, Nothing, Nothing}, delay_cross_f::typeof(GeometryOps._inter_delay_cross_f), delay_bounce_f::typeof(GeometryOps._inter_delay_bounce_f); exact::GeometryOpsCore.True)
    @ GeometryOps ~/.julia/packages/GeometryOps/4octE/src/methods/clipping/clipping_processor.jl:149
  [6] _build_ab_list
    @ ~/.julia/packages/GeometryOps/4octE/src/methods/clipping/clipping_processor.jl:143 [inlined]
  [7] _intersection(alg::GeometryOps.FosterHormannClipping{GeometryOpsCore.Planar, GeometryOps.AutoAccelerator}, ::GeometryOpsCore.TraitTarget{GeoInterface.PolygonTrait}, ::Type{Float64}, ::GeoInterface.PolygonTrait, poly_a::GeoInterface.Wrappers.Polygon{false, false, Vector{GeoInterface.Wrappers.LinearRing{false, false, Vector{Tuple{Float64, Float64}}, Nothing, Nothing}}, Nothing, Nothing}, ::GeoInterface.PolygonTrait, poly_b::GeoInterface.Wrappers.Polygon{false, false, Vector{GeoInterface.Wrappers.LinearRing{false, false, Vector{Tuple{Float64, Float64}}, Nothing, Nothing}}, Nothing, Nothing}; exact::GeometryOpsCore.True, kwargs::@Kwargs{})
    @ GeometryOps ~/.julia/packages/GeometryOps/4octE/src/methods/clipping/intersection.jl:86
  [8] _intersection
    @ ~/.julia/packages/GeometryOps/4octE/src/methods/clipping/intersection.jl:76 [inlined]
  [9] #intersection#102
    @ ~/.julia/packages/GeometryOps/4octE/src/methods/clipping/intersection.jl:45 [inlined]
 [10] intersection

[asinghvi17]

hmm, those look like this

but computing the intersection seems to work fine for me manually.

[asinghvi17]

hmm - looks like this is an antimeridian cutting issue, this is the polygon from grid 2:

edit: fixed

[simone-silvestri]

Hi @asinghvi17, when I try to use this branch in combination with ConservativeRegridding.jl I get the following error:

ERROR: The following 2 direct dependencies failed to precompile:

GeometryOps

Failed to precompile GeometryOps [3251bfac-6a57-4b6d-aa61-ac1fef2975ab] to "/Users/simonesilvestri/.julia/compiled/v1.10/GeometryOps/jl_BD8CaS".
ERROR: LoadError: UndefVarError: `CutAtAntimeridianAndPoles` not defined
Stacktrace:
 [1] top-level scope
   @ ~/.julia/packages/GeometryOps/LoXLj/src/transformations/correction/cut_at_antimeridian.jl:80
 [2] include(mod::Module, _path::String)
   @ Base ./Base.jl:495
 [3] include(x::String)
   @ GeometryOps ~/.julia/packages/GeometryOps/LoXLj/src/GeometryOps.jl:3
 [4] top-level scope
   @ ~/.julia/packages/GeometryOps/LoXLj/src/GeometryOps.jl:89
 [5] include
   @ ./Base.jl:495 [inlined]
 [6] include_package_for_output(pkg::Base.PkgId, input::String, depot_path::Vector{String}, dl_load_path::Vector{String}, load_path::Vector{String}, concrete_deps::Vector{Pair{Base.PkgId, UInt128}}, source::Nothing)
   @ Base ./loading.jl:2292
 [7] top-level scope
   @ stdin:4
in expression starting at /Users/simonesilvestri/.julia/packages/GeometryOps/LoXLj/src/transformations/correction/cut_at_antimeridian.jl:53
in expression starting at /Users/simonesilvestri/.julia/packages/GeometryOps/LoXLj/src/GeometryOps.jl:3
in expression starting at stdin:4
ConservativeRegridding

Failed to precompile ConservativeRegridding [8e50ac2c-eb48-49bc-a402-07c87b949343] to "/Users/simonesilvestri/.julia/compiled/v1.10/ConservativeRegridding/jl_UH89QO".
ERROR: LoadError: UndefVarError: `CutAtAntimeridianAndPoles` not defined
Stacktrace:
 [1] top-level scope
   @ ~/.julia/packages/GeometryOps/LoXLj/src/transformations/correction/cut_at_antimeridian.jl:80
 [2] include(mod::Module, _path::String)
   @ Base ./Base.jl:495
 [3] include(x::String)
   @ GeometryOps ~/.julia/packages/GeometryOps/LoXLj/src/GeometryOps.jl:3
 [4] top-level scope
   @ ~/.julia/packages/GeometryOps/LoXLj/src/GeometryOps.jl:89
 [5] include
   @ ./Base.jl:495 [inlined]

[asinghvi17]

Ah my bad, looks like I pushed an incomplete change set. If you can go back one commit it should work, and I'll push a fix in half an hour or so.

[simone-silvestri]

No problem, no rush, I imagined you guys were still working on it 😄

[asinghvi17]

OK this should work now!