3

I am trying to find a way to calculate the direct problem on the WGS-84 ellipsoid but instead of using the initial Azimuth (e.g. lat1, lon1, distance, fwd azi), I want to compute it with the reverse Azimuth (e.g. lat1, lon1, distance, point2-reverse-azimuth).

Are there any tools/libraries out there to do that?

My purpose is to find the tangential latitudes and longitudes to a circle of radius X from any given point -- to calculate an accurate lat/long bounding box that contains all points within X distance.

3
  • 1
    Welcome to GIS SE. I understand the first part of your question (and don't know the solution). The second part seems relatively trivial: covert a difference in latitude distance into a difference in latitude angle. I don't understand the connection between parts one and two, however. Can you edit to clarify?
    – Martin F
    Commented Jan 24, 2015 at 18:39
  • The questioner's problem is to find not only the latitude extent of the bounding box (relatively simple, as you note), but also the longitude extent (a more interesting problem); see below for my answer.
    – cffk
    Commented Jan 25, 2015 at 15:43
  • OK, thanks to cffk, i understood, and slightly edited, your Q.
    – Martin F
    Commented Jan 25, 2015 at 19:46

1 Answer 1

4

You are trying to solve the ellipsoidal triangle where

φ1, α2, s12

are given. (See Figure 1 of the Wikipedia article, Geodesics on an ellipsoid, for the notation.) In your case, we have α2 = ±½π because of the requirement of tangency. This is Problem 7 in §10 of my paper

Geodesics on an ellipsoid of revolution (Feb. 2011)

The solution entails assuming a value of α1 (perhaps by solving an equivalent spherical problem), solving the ellipsoidal triangle given

φ1, α1, α2

(Problem 2), determining the resulting s12, and correcting α1 using Newton's method. The derivative needed for Newton's method, ds12/dα1 is given by Eq. (79). Problem 2 is trivial to solve: convert to the auxiliary sphere (β1, α1, α2), determine σ12 using some subset of Eqs. (9) thru (18), and find s12 by using the routine Geodesic::ArcDirect in my library GeographicLib. (Here β is the parametric latitude and σ is the arc distance on the auxiliary sphere.)

Assuming that your initial guess is sufficiently good, this will converge to the desired solution quadratically (the number of correct digits in α1 will double on each iteration).

2
  • You are a computational geodesist of the highest order! I wish you'd participate even more often.
    – Martin F
    Commented Jan 25, 2015 at 19:40
  • Not only that, but very generous: I simultaneously emailed my question and received a response within hours. Thank you again!
    – JustinHK
    Commented Jan 27, 2015 at 14:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.