algorithmgpsspherical-geometry
I need to make an algorithm:
:
I have a point O center of a circle of radius r, on this circle there is two point P1 and P2.
I want to know the latitude and longitude of point P2 on the same circle, I know the angle a in degrees between the line OP1 and OP2.
How to solve this problem?
The gps points O,P1 and P2 are situated in the south of France.
For this purpose simple approximations are more than good enough. North or south, one degree is about 69 miles but east or west, it is only 69*cos(latitude) miles. Because latitudes do not change much over a ten mile span, you can safely use the cosine of the central latitude of the "square." Therefore the desired coordinates for square vertices at distance r miles from a central location (f,l), given as lat-lon, are computed as
df = r/69 // North-south distance in degrees
dl = df / cos(f) // East-west distance in degrees
{(f-df,l-dl), (f+df,l-dl), (f+df,l+dl), (f-df,l+dl)} // List of vertices
For example, suppose r = 10 miles and the central location is at latitude 50 degrees north, longitude 1 degree west, so that (f,l) = (50,-1) degrees. Then
df = 10/69 = 0.145
dl = 0.145 / cos(50 degrees) = 0.145 / 0.6428 = 0.225
f - df = 50 - 0.145 = 49.855 (southernmost latitude)
f + df = 50 + 0.145 = 50.145 (northernmost latitude)
l - dl = -1 - 0.225 = -1.225 (western longitude)
l + dl = -1 + 0.225 = -0.775 (eastern longitude)
and the coordinates are (49.855,-1.225), (50.145,-1.225), (50.145, -0.775), and (49.855, -0.775) as you march clockwise around the square starting at its southwestern corner.
Don't use this approximation near the poles or for squares larger than a few degrees on a side. Also, depending on limitations of the GIS, some care might be needed around the global cut in longitude, usually taken at +-180 degrees.
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.
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.
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.
Best Answer
Point P1 is irrelevant, if it's a circle r will always be the same, this together with the bearing theta and the origin point O are what matters. To find theta, simply trig your way from your provided angle a.
Supposing your surface is a sphere:
In the above formula, lat, lon and theta are to be provided in radians. Theta is, as in your example, the bearing clockwise from the North. Use absolute values if you coordinates are negative.
This is an adaptation of the Haversine equation, if you want to look it up. It'll be as precise as the Earth radius you choose.
Supposing your surface is an elipsoid:
Ok, this is way more complicated, so I'll just link you to Vincenty's formula. You want to use the Direct Formula to find your P2. The same precautions apply here, concerning negative coordinates.