I have a little computational geometry project I'm struggling with for a non-commercial "art" installation. It is driving me crazy and I'd happily pay for an implementable algorithm/solution (apologies if Stack Exchange doesn't allow this - I looked but didn't find any restriction on offering a cash bounty).

Given a starting location P1 on the surface of the earth (latitude phi1, longitude L1), compute the exit point P2 (latitude phi2, longitude L2) of a ray that departs from the starting point at local azimuthal angle alpha1 (from true north) with downward altitude angle beta1.

In less formal terms, imagine sitting on a rotating stool. Turn the stool from true north by some angle alpha1. Point your arm straight out. Now drop your arm by angle beta1. Imagine a ray that extends from your pointed finger into the ground. Imagine this ray continuing straight through the earth. What point P2 does this ray emerge from the earth?

Notes and Suggestions:

1) As a first approach, suggest using a spherical earth model. Here are some useful resources, http://www.movable-type.co.uk/scripts/latlong.html

2) The "direct Vincenty method", http://en.wikipedia.org/wiki/Vincenty%27s_formulae, is potentially useful. It provides an iterative approach to the related problem of determining the azimuthal angles when point P2 is known. The method uses an oblate spheroid model of the earth (same as WGS84)

3) Determine how much accuracy is improved by 2 versus 1.

4) Determine if incorporating height above sea level of P1 or P2 has significant impact on answer.

  • $\begingroup$ Are you considering the Earth as a homogeneous sphere, or are you really interested in considering, say, the direction earthquake waves take through a realistically structured Earth with its many layers? $\endgroup$ Feb 10 '15 at 4:44
  • $\begingroup$ @WolfgangBangerth, for the purposes of my project, the Earth can be considered homogeneous. I'm only concerned with a notionally straight ray and not the propagation of physical waves. $\endgroup$ Feb 10 '15 at 9:04

In spherical coordinates, this is not much harder than solving a quadratic equation. Let the coordinates be $(r,\theta,\varphi)$, (radial, longitude, latitude). If the starting point is $\mathbf{R}=R\mathbf{e}_r$, then after rotating by $\alpha$, the observer is pointing in the direction of the unit vector $$ \mathbf{e}_\alpha = \mathbf{e}_\varphi \cos\alpha + \mathbf{e}_\theta\sin\alpha, $$ and after rotating by $\beta$, the observer is pointing in the direction of $$ \mathbf{e}_\beta = \mathbf{e}_\alpha\cos\beta - \mathbf{e}_r\sin\beta. $$

So on the line passing through $\mathbf{R}$ with direction $\mathbf{e}_\beta$, if $t$ is a parameter, the point where the line crosses the Earth's sphere in a second point is given by a quadratic equation in $t$: $$ \|R\mathbf{e}_r + t\mathbf{e}_\beta\|^2 = R^2. $$ Expanding the vector norm, solving the equation to get the other root $t$, finding the point corresponding to that $t$, and converting that point to spherical coordinates gives you the coordinates of the other intersection.

If you want to use an ellipsoid instead, the basis vectors and the equation of intersection will change. If you want to use terrain height (call it $h$), the change in angles will be on the order of $h/R$ where $R$ is Earth's radius, and $h/R$ will be quite tiny.

  • $\begingroup$ I believe that the ray is given in terms of "local" az/el so there also needs to be a change of coordinates from local az/el to global coordinates. Not difficult, but just one minor addition. $\endgroup$ Feb 9 '15 at 2:29
  • $\begingroup$ @DougLipinski Well, $e_\varphi$ points north and $e_\theta$ points east, so $e_\alpha$ points an angle $\alpha$ to north, so I think that's right, same for $e_\beta$, what change do you have in mind? $\endgroup$
    – Kirill
    Feb 9 '15 at 2:41
  • $\begingroup$ Looks like maybe I just didn't read your answer carefully enough. I'll have to think a bit about the angles involve, but this sounds good. $\endgroup$ Feb 9 '15 at 15:07
  • $\begingroup$ @Kirill, don't we need a pre and post transform to account for the fact that, in general, our starting point isn't going to be the north pole (R*e_r) but some other point (r,theta,phi) on the round Earth's surface? $\endgroup$ Feb 10 '15 at 9:01
  • $\begingroup$ @DarranEdmundson If the observer is at $(r,\theta,\varphi)$ (or $(x,y,z)$), the spherical coordinates basis vectors are, by definition, such that $e_r$ points from the origin to the observer, $e_\varphi$ points north (or south, with a different convention) from the observer at the observer's location on the sphere around the origin with radius $R$, etc. Use this: mathworld.wolfram.com/SphericalCoordinates.html $\endgroup$
    – Kirill
    Feb 10 '15 at 9:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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