[Proj] Spheroidal gnomonic projection
Noel Zinn
ndzinn at comcast.net
Mon Jun 14 08:25:37 EST 2010
Charles,
Here's an algorithm for an ellipsoidal gnomonic that is a direct perspective
from the geocenter through the ellipsoid onto the tangential plane via ECEF,
topocentric coordinates and a little matrix algebra. No intermediate
mappings. No truncations. No approximations.
Choose an oblique point of tangency (lat, lon, hgt=0). Compute the ECEF XYZ
as you have described. Define that as (0,0,0) topocentric UVW. Chose a
point to the east (100000,0,0) and a point to the north (0,100000,0).
Convert these points from topocentric UVW to ECEF XYZ. EPSG Guidance Note
7-2 Section 2.2.2 will do (www.epsg.org) or, if Bugayevskiy and Snyder is
your slug of vodka, use the equations on page 3 (less elegant). You now
have three points in ECEF XYZ that define the tangential plane. Any point
on the ellipsoid can be converted to ECEF XYZ. The geocenter is (0,0,0) in
ECEF XYZ. Those two points define the perspective ray. Now, we need the
point of intersection of the ray with the plane. Wikipedia
(http://en.wikipedia.org/wiki/Line-plane_intersection) provides the answer,
which requires the inversion of a 3x3 matrix, but we've already computed all
the required xyz points. Multiply the parameter 't' (throw away 'u' and
'v') by the ECEF XYZ coordinates of the point on the ellipsoid to get the
ECEF XYZ coordinates on the tangential plane. Convert those XYZ coordinates
to topocentric UVW. The W will be zero (because we're on the plane). The U
and V coordinates are our ellipsoidal gnomonic Easting and Northing.
On the 'downloads' page of my website (link below) I've posted a
presentation containing a mapping of some simple North American coastline
culture in the ellipsoidal gnomonic perspective projection and the essential
Matlab code that accomplishes this. You may prefer to use Octave.
Assessing the linearity of a geodesic on this implementation is left as an
exercise for the interested.
Regards,
Noel Zinn, Principal
Hydrometronics LLC, Consultancy and Technical Software
+1-832-539-1472, noel.zinn at hydrometronics.com
http://www.hydrometronics.com
-----Original Message-----
From: proj-bounces at lists.maptools.org
[mailto:proj-bounces at lists.maptools.org] On Behalf Of Karney, Charles
Sent: Sunday, June 13, 2010 9:49 PM
To: proj at lists.maptools.org
Cc: deakin at rmit.edu.au; knud.poder at mail.dk; Craig.M.Rollins at nga.mil
Subject: Re: [Proj] Spheroidal gnomonic projection
I've tracked down another ellipsoidal gnomonic projection
I. G. Letoval'tsev, Generalization of the Gnomonic Projection for a
Spheroid and the Principal Geodetic Problems Involved in the Alignment
of Surface Routes, Geodesy and Aerophotography (5), 271-274 (1963).
As described by Bugayevskiy and Snyder, this consists of a conformal
projection to a sphere followed by a central projection onto a tangent
plane. However, on reading the paper I see that there's an important
wrinkle: the conformal projection is offset in latitude in order to make
the geographic and conformal latitudes match at the center point
(instead of the equator), i.e.,
asinh(tan(beta)) = asinh(tan(phi))
- e * ( atanh( e * sin(phi) ) - atanh( e * sin(phi0) ) )
With this projection, normal sections are mapped to straight lines. The
maximum error in the direction of a geodesic within a radius r of the
center point scales as f*(r/a)^2.
In summary, we have 3 ellipsoidal gnomonic projections
* center projection, Williams (1998)
great ellipses are straight
error in azimuth for geodesics = f*(r/a) = 108" for r = 1000km
* conformal mapping, Letoval'tsev (1963)
normal sections through center point are straight
error in azimuth for geodesics = f*(r/a)^2 = 17.4" for r = 1000km
* geodesic mapping, Karney (2010)
geodesics through center point are straight
error in azimuth for geodesics = 0.4*f*(r/a)^3 = 1.04" for r = 1000km
Postscript: The affiliation of Letoval'tsev is listed as "Moscow
Institute of Railroad Transportation Engineers". I wonder how often you
see "Railroad" appearing in the addresses of articles on geodesy
nowadays.
--
Charles Karney <ckarney at sarnoff.com>
Sarnoff Corporation, Princeton, NJ 08543-5300
Tel: +1 609 734 2312
Fax: +1 609 734 2662
_______________________________________________
Proj mailing list
Proj at lists.maptools.org
http://lists.maptools.org/mailman/listinfo/proj
More information about the Proj
mailing list