Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Integrate better with GeoInterface.jl and Proj.jl to plot any geoms in any projections #161

Open
rafaqz opened this issue Apr 18, 2023 · 16 comments
Labels
enhancement New feature or request

Comments

@rafaqz
Copy link

rafaqz commented Apr 18, 2023

GeoMakie.jl could accept any geointerface geometries in any projection.

Some geometries have the projection attached already, so we could do the conversion automaticaly. Otherwise plotting onto a GeoAxis could accept a source_crs keyword and do the conversions.

If the geometries are in different projections to the map and/or each other we can reproject them with JuliaGeo/Proj.jl#79

@asinghvi17
Copy link
Member

hmm, that sounds hard to implement if it's the object which has the projection, since we can't really dispatch on that...

@rafaqz
Copy link
Author

rafaqz commented Apr 18, 2023

We can dispatch on GeoInterface.crs(object) - it will be a GeoFormatTypes.jl wrapper. Users could also pass in crs as a keyword like poly!(geoaxis, geometries; source_crs=EPSG(4326)) when its not attached to the object.

So it would need specialisation fot those methods on GeoAxis.

Then we can use Proj.jl to reproject the geometries to the axis projection.

Its quite easily doable ;)

@asinghvi17
Copy link
Member

We can do that...but I don't think Makie can. Maybe in the poly recipe, using a convert to a geointerface wrapper or custom wrapper type? But not generally.

@felixcremer
Copy link

Wouldn't that also be a use case for the CRS trait discussed in JuliaGeo/GeoInterface.jl#119 ?

@asinghvi17
Copy link
Member

This is a bit orthogonal - the problem lies in extracting the CRS, not in dispatching on it. That's already implemented in the code anyway.

@rafaqz
Copy link
Author

rafaqz commented Jun 7, 2024

What is the problem extracting it exactly?

@asinghvi17
Copy link
Member

The question, I guess, was where to perform the extraction - but I guess it's actually possible to do this in the same place we pop the source kwarg! So it's probably doable then, at least for single-argument plots.

@asinghvi17
Copy link
Member

GeoMakie.jl/src/geoaxis.jl

Lines 817 to 830 in 8941fe9

function Makie.plot!(axis::GeoAxis, plot::Makie.AbstractPlot)
source = pop!(plot.kw, :source, axis.source)
transformfunc = lift(create_transform, axis.dest, source)
trans = Transformation(transformfunc; get(plot.kw, :transformation, Attributes())...)
plot.kw[:transformation] = trans
Makie.plot!(axis.scene, plot)
# some area-like plots basically always look better if they cover the whole plot area.
# adjust the limit margins in those cases automatically.
Makie.needs_tight_limits(plot) && Makie.tightlimits!(axis)
if Makie.is_open_or_any_parent(axis.scene)
Makie.reset_limits!(axis)
end
return plot
end

source would also have to check plot.args to see if they have a CRS. This won't work in recipes, but for user facing plots it's better than nothing.

@rafaqz
Copy link
Author

rafaqz commented Jun 7, 2024

Ahh no it will already be gone right? Converted to GeometryBasics with no crs?

Or is the original object still around?

@asinghvi17
Copy link
Member

The original object should still be around in plot.args...

@rafaqz
Copy link
Author

rafaqz commented Jun 7, 2024

Ok cool, I've no idea how it works internally

@rafaqz
Copy link
Author

rafaqz commented Jun 7, 2024

You should be able to call crs on it. Or if it's a Vector, crs on the first object.

Mixed crs can be undefined 😂

@asinghvi17
Copy link
Member

Yeah that's what I was thinking as well :D

@asinghvi17
Copy link
Member

asinghvi17 commented Jun 7, 2024

Here's what I'm thinking:

get_crs(obj) = get_crs(GI.trait(obj), obj)
get_crs(::GI.AbstractGeometryTrait, obj) = GI.crs(obj)
function get_crs(::Nothing, obj)
	if obj isa AbstractVector
		if GI.trait(first(obj)) isa GI.AbstractGeometryTrait
			return get_crs_vector_checked(obj)
		end
	end
	return GI.crs(obj)
end
function get_crs_vector_checked(obj)
	crs_es = GI.crs.(obj)
	if length(unique(crs_es)) == 1
		return first(crs_es)
	else
		@warn("You have tried to plot objects with conflicting CRS in the same plot call!  Assuming `crs=nothing` for now.")
		return nothing
	end
end

@asinghvi17 asinghvi17 added the enhancement New feature or request label Jun 7, 2024
@asinghvi17
Copy link
Member

This is priority after #261 and #245

@asinghvi17
Copy link
Member

It turned out that this caused some form of infinite loop, I'm still not sure why. Will keep looking.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement New feature or request
Projects
None yet
Development

No branches or pull requests

3 participants