7

There are few ways of articulating my question. It starts with a known circle on a sphere, i.e. a point-radius, and one that does not cover a pole. I would like to calculate the point furthest west on this circle and the point furthest east. These points are the tangent points of lines of longitude which encircle the poles and are vertical. These points are opposite each other along the vertical axis so there is one hard problem to solve from which the other is trivially determined. This is highly related to calculating the minimum bounding lat-lon box of a circle. I know how to determine the longitude delta from the center to the edges, which is:

asin( sin( dist_rad / cos(center_lat_rad) ) )

whereas dist_rad is the circle's radial distance is radians, and center_lat_rad is the latitude point of the circle's center in radians. This equation was elusive to me until I discovered it with a wonderful explanation here. It contains a wonderful illustration:

wonderful illustration.

The illustration shows T1 and T2 which are the points I am trying to calculate. What is left to figure out is at what latitude the eastern and western boundary meridians intersect this circle. If the circle is at the equator, then it's zero, and if it's north of the equator, then this latitude is going to be higher than the circle's center latitude.

4
  • 1
    Would a numerical solution suffice? In other words could something like Newton's method be used to make successive guesses at the latitude be made. For each guess you could find the geodesic distance to M, the guess with the minimum distance wins. Commented Jan 24, 2012 at 21:07
  • In some cases that's a good approach, @Kirk, but here it won't be very precise. If you were to graph the distance to M as a function of position along the extreme meridian through T2, you would find that it hardly changes close to T2: the function is quadratic there. This is precisely the kind of situation where methods like Newton's can break down or lose a great deal of precision. (It's related to the reason we use a Haversine law rather than a law of Cosines for finding spherical distances.)
    – whuber
    Commented Jan 24, 2012 at 21:18
  • Would this problem be simpler if you converted your Lat/Long to X/Y/Z space? Then you're just solving for an (x,y,z) that lies on the two spheres expressed by XX+YY+ZZ=RR? Commented Jan 24, 2012 at 23:37
  • That's always a reasonable option @Bicycle. Given that the spherical trig solution is so simple (and that all spherical trig identities can be derived from your suggestion), one would expect that doing that math would yield the same answer. (My guess is that it would be the same answer in a more complicated form due to the conversions from spherical to 3D cartesian coordinates and back again.)
    – whuber
    Commented Jan 24, 2012 at 23:45

2 Answers 2

10

Questions like this can often be answered by solving triangles using laws of spherical trigonometry. You can look these up: they are referenced below. Normally, when you know three of the six parts of a spherical triangle, the laws let you find the other three in terms of the sines and cosines of the parts you know. The trick is to draw a useful triangle. When you're dealing with meridians, think of using a pole (which is common to all meridians) as one of the triangle's vertices. This suggests focusing on the triangles N-M-T1 or, equivalently, N-M-T2. Here are the details, a general solution, and a worked example.

Let the point in question be at latitude phi and the circle's radius be r (expressed as an angle, not as a distance. The angle in radians is the radial distance divided by the sphere's radius). We start by finding the difference between the local longitude and the extreme longitudes; let this difference be lambda.

Because the extreme meridians are tangent to the circle, angles M-T1-N and M-T2-N are right angles. This is nice because the sine and cosine of a right angle are easy to find and simple (they equal 1 and 0, respectively). The spherical Law of Sines says

sin(lambda) / sin(r) = sin(90 degrees) / sin(90 - phi) 
= 1 / cos(phi).

Solving this gives

lambda = ArcSin( sin(r) / cos(phi) ).

For example, let r = 45 degrees and phi = -30 degrees. The solution is lambda = 54.7356 degrees.

This much you already found out: it's the formula quoted in the question. Now we can apply the spherical Law of Cosines to relate the latitude phi' of the extremal points T1 and T2 to the radius r: this is what we are looking for.

sin(phi) = sin(phi') cos(r).

The solution for the common latitude of the tangent points is

phi' = ArcSin( sin(phi) / cos(r) )

(for a circle centered at latitude phi of radius r radians.)


The remainder of the reply discusses and illustrates this result.

Notice in particular that when we start at the equator, phi = 0, easily giving phi' = 0 as expected. Otherwise, the sine of phi' is larger (in size) than the sine of phi, implying the solution is closer to the nearest pole, again as expected. In the example, phi' works out to -45 degrees.

Figure

The south pole is near the bottom of this pseudo-3D image of a sphere. The meridians start at the circle's center and go in 15 degree increments to either side, showing +-15, +-30, and +-45 degree longitudinal displacements, then stop at the limiting values of +-54.7356 degrees. The red dots are situated along these limiting meridians at latitudes of -45 degrees.

These formulas work for any circle that does not include either pole.

Figure 2

In this image, the circle's center (black dot) is just north of 60 degrees north latitude. Approximate circles of latitude are shown as dashed curves in 10 degree increments. The radius, equal almost to 30 degrees of arc, takes this circle almost to the north pole (where the meridians are converging at the top). The extremal meridians are therefore spread far apart (lambda is about 70 degrees) and, accordingly, the points of tangency are also close to the north pole, near 80 degrees latitude. This shows why the points of tangency (red dots) are usually closer to the nearest pole than the circle's center is.

11
  • If, as you say, "angles M-T1-N and M-T2-N are right angles". Then the line M-T1 must follow the Latitude line (all lines of latitude cross each meridian at a right angle (in a sphere)). Therefore phi is solved - it is the latitude of M, just figure out lambda now.
    – brenth
    Commented Jan 24, 2012 at 21:22
  • I'm afraid that argument is incorrect, @Brenth. It assumes lines of latitude are straight, but (except at the equator) they are not. (This also helps us see why your argument does give the correct answer when the center lies on the equator.) For a detailed explanation of the phenomenon, please refer to the answer I posted at gis.stackexchange.com/questions/6822/….
    – whuber
    Commented Jan 24, 2012 at 21:26
  • 1
    Sorry brenth, but whuber is correct. Spherical geometry is a bitch, and you are unfortunately making some incorrect assumptions. As whuber points out, lines of latitude are not the spherical equivalent of "straight" lines: only great circles are.
    – nmpeterson
    Commented Jan 24, 2012 at 22:27
  • 2
    This answer is stupendous; it really shows the value of stack exchange. And I got the answer fast. I took the formula at the end of the long explanation and put it into my system which has an extensive test suite to effectively prove the answer is correct. There were some NaN result edge-cases when the circle touches the pole in which I mapped the answer to +90 or -90 depending on wether the circle's center is above or below the equator. Commented Jan 26, 2012 at 16:28
  • 1
    Thanks @whuber , I was trying to have a bounding box around a point, now it is clear, and also test is matching the result. For anyone who is interested Python code in gist.github.com/alexcpn/f95ae83a7ee0293a5225 Commented Feb 20, 2016 at 10:52
0

Struggled a bit to understand first the answer. Here is the python port of whubers code; Putting it in an answer, for better visibility for others as earlier it was buried in comments

Python code for bounding box with test for the above

Here is more amateur friendly drawing to explain the same

enter image description here

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.