Solving the direct geodesic problem with Proj.jl

Hi!
I want to calculate the end point after moving some distance in some direction on an ellipsoid.
I used to use the Proj4.jl function geod_direct. As this higher level function has disappeared in Proj.jl, I’m trying to deal with the lower level version. Taking inspiration from the C example in the proj docs, I end up with this code:

g = Ref{Proj.geod_geodesic}()
a = 6378137.; f = 1/298.257223563
Proj.geod_init(g, a, f)

lon, lat, distance, azimuth = 4., 50., 111321., 90.
outlon, outlat = Ref{Cdouble}(), Ref{Cdouble}()
Proj.geod_direct(g, lat, lon, distance, azimuth, outlat, outlon, Ref(0.))

And I get: outlon[], outlat[] == 4.0024797121980935, 50.000253128588675, which is wrong, I should get something like 5., 50..

I’m not used to interfacing C in Julia, so maybe I’m doing something stupid here. I also tried to get inspiration from the implementation in ´Proj4.jl´, but the proj API has changed since then.

Any hint where this wrong result could come from?

Thank you for your time!

My guess is your problem is that distance and azimuth are in the wrong order. At least if I calculate the endpoint from travelling 90 m along an azimuth of 111,321° == 81°, I get approximately the same result you obtain above.

It may be useful to you to note that GeographicLib.jl provides a Julia implementation of Charles Karney’s algorithm (used by Proj), although it is an unregistered package.

The contents of GeographicLib.jl will be folded into Geodesics.jl fairly soon and then the latter will be registered.

2 Likes

See the geod in GMT.jl
And a cool application in the Longest sail

2 Likes

@anowacki
Wow, sorry, I thought I had really thought through the problem before posting, but apparently I didn’t sufficiently at all, I feel quite stupid :face_in_clouds:.

Yes I would prefer to use GeographicLib.jl but the main stumbling block to me is that it isn’t registered. I’m really happy to know that it will be taken care of soon, thank you!

@joa-quim
Thank you for the suggestion. GMT.jl seems pretty cool, but I think that’s a quite heavy package, and I just want to use this function at the moment. But I’ll definitely consider using it if I need more of such features!