land() function raises error

[paulmelis]

Not sure what's going on, but seems to be a Nothing ends up getting used somewhere. However, directly calling NaturalEarth.naturalearth returns a valid set of geometry.

melis@blackbox 12:26:~/projects/harmonie$ j --project=.
               _
   _       _ _(_)_     |  Documentation: https://docs.julialang.org
  (_)     | (_) (_)    |
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 1.10.4 (2024-06-04)
 _/ |\__'_|_|_|\__'_|  |  Official https://julialang.org/ release
|__/                   |

julia> using GeoMakie
[ Info: Precompiling GeoMakie [db073c08-6b98-4ee5-b6a4-5efafb3259c6]
[ Info: Skipping precompilation since __precompile__(false). Importing GeoMakie [db073c08-6b98-4ee5-b6a4-5efafb3259c6].

julia> land = GeoMakie.land(50)
ERROR: MethodError: no method matching to_multipoly(::Nothing, ::Vector{GeometryBasics.MultiPolygon{2, Float32, GeometryBasics.Polygon{2, Float32, GeometryBasics.Point{…}, GeometryBasics.LineString{…}, Vector{…}}, Vector{GeometryBasics.Polygon{…}}}})

Closest candidates are:
  to_multipoly(::GeoInterface.MultiPolygonTrait, ::Any)
   @ GeoMakie ~/.julia/packages/GeoMakie/y0Tfj/src/geojson.jl:103
  to_multipoly(::GeoInterface.PolygonTrait, ::Any)
   @ GeoMakie ~/.julia/packages/GeoMakie/y0Tfj/src/geojson.jl:102
  to_multipoly(::GeoInterface.GeometryCollectionTrait, ::Any)
   @ GeoMakie ~/.julia/packages/GeoMakie/y0Tfj/src/geojson.jl:105
  ...

Stacktrace:
 [1] _broadcast_getindex_evalf
   @ ./broadcast.jl:709 [inlined]
 [2] _broadcast_getindex
   @ ./broadcast.jl:682 [inlined]
 [3] getindex
   @ ./broadcast.jl:636 [inlined]
 [4] copyto_nonleaf!(dest::Vector{…}, bc::Base.Broadcast.Broadcasted{…}, iter::Base.OneTo{…}, state::Int64, count::Int64)
   @ Base.Broadcast ./broadcast.jl:1098
 [5] copy
   @ ./broadcast.jl:950 [inlined]
 [6] materialize
   @ ./broadcast.jl:903 [inlined]
 [7] to_multipoly(geom::Vector{Any})
   @ GeoMakie ~/.julia/packages/GeoMakie/y0Tfj/src/geojson.jl:101
 [8] land(scale::Int64)
   @ GeoMakie ~/.julia/packages/GeoMakie/y0Tfj/src/data.jl:63
 [9] top-level scope
   @ REPL[2]:1
Some type information was truncated. Use `show(err)` to see complete types.

julia> using NaturalEarth

julia> NaturalEarth.naturalearth("land", 50)
FeatureCollection with 1420 Features

(harmonie) pkg> st
Status `~/projects/harmonie/Project.toml`
  [13f3f980] CairoMakie v0.12.5
  [b16dfd50] GRIB v0.4.0
  [db073c08] GeoMakie v0.7.2
  [436b0209] NaturalEarth v0.1.0
  [91a5bcdd] Plots v1.40.5

[paulmelis]

Ah,

GeoMakie.jl/src/geojson.jl

Lines 35 to 36 in 4e39e7e

	elseif GeoInterface.isfeaturecollection(input)
	elseif input isa AbstractArray && GeoInterface.isgeometry(eltype(input))

returns Nothing, because it has no body:

    if GeoInterface.isgeometry(input)
        return geo2basic(GeoInterface.geomtrait(input), input)
    elseif GeoInterface.isfeature(input)
        return geo2basic(GeoInterface.getgeom(input))
    elseif GeoInterface.isfeaturecollection(input)
        # No return value
    elseif input isa AbstractArray && GeoInterface.isgeometry(eltype(input))
        return geo2basic.(input)
    else
        @error("Input of type $(typeof(input)) does not support GeoInterface!")
    end

[asinghvi17]

Hmm, you should be able to get away without to_multipoly these days

[paulmelis]

Hmm, you should be able to get away without to_multipoly these days

But the problem is in geo2basic()?

Edit: I tried passing the output of naturalearth() directly, e.g.:

    land = NaturalEarth.naturalearth("land", 50)

    fig = Figure()
    ax = GeoAxis(fig[1,1],
        title=pname,
        source="+proj=longlat +datum=WGS84",
        dest="+proj=lcc +lon_0=5 +lat_1=50 +lat_2=55")
    xlims!(ax, 1, 9)
    ylims!(ax, 50, 54.5)
    surf = surface!(ax, lons, lats, values; shading=NoShading, colormap=:GnBu)
    println("lines!()")
    lines!(ax, land)

But it just seems to get stuck on the lines!() call, drawing CPU

[asinghvi17]

Will take a crack at this and potentially replace the to_multipoly implementation with GeoInterfaceMakie's version once I merge #261 and #245.

[asinghvi17]

I can't actually replicate the land error locally @paulmelis - though I do see the same "freezing up" from lines(naturalearth("land", 50)).

[asinghvi17]

The freezing up is caused by GeoInterfaceMakie not handling feature collections, we just need to provide a better error message there. lines(naturalearth("land", 50).geometry) causes a different error because GeometryBasics multipolygons are not supported in lines, which I will fix in GeoMakie.