[Proj] Geo Transfoms and calculations

Gerald I. Evenden geraldi.evenden at gmail.com
Mon Dec 29 19:55:18 EST 2008

Your timing could not be much better.

I have almost completed a library set 'geod_*' that covers geodesic 
computations for both the spherical and elliptical earth.  I expect to have a 
first release out by this weekend as I am only doing finishing touches on the 
ellps= section (a user program version works right now with the libproj4 
library).  The ellipse method is based upon the NOAA routines rather than the 
USNOO process that I used in 'geod'---accuracy is also a little better.

The release comes with Doxygen html documentation (~2mb), of all the code and 
a self documenting user program 'geodesic'---that is: you type "geodesic" and 
at the prompt "geod:" type "help."  Maybe a little more help on the Doxygen 
front page. ;-)

Also, see note imbedded below and geodesic example at end.

On Monday 29 December 2008 7:07:28 pm Till wrote:
> Hi all,
> We're quite new to geographic fundamentals, but it sounds like PROJ.4 is a
> pretty good base where to find proper 'geographic related' libraries.
> We're currently working on graphic (2D/3D) tools for displaying accurate
> GPS information and our goal is to rely on standard systems, avoiding to
> 'reinvent the wheel', feeling it should be much of a benefit if we could
> share with the community.
> By now, we were able to properly use the 'proj' API in order to match
> different projection schemes into 'Cartesian' (3D space) coordinates, and
> also mix 2D representations.
> However, distance calculation seems to be another issue.
> While we could find ways to implement a reasonable calculation
> (for instance :
> gpsDist(PointF P1, PointF P2)
> {
This is a common equation for spherical geodesic but it is not very good for 
small values of the argument because acos looses precision.  All geodesic 
code I have used employs Snyder's recommended variation as in the geod_* 

> return Acos( sin(P1.y)*sin(P2.y) + cos(P1.y)*cos(P2.y)*cos(P1.x-P2.x) ) *
> earthRadius;
> }),
> our approach considers it could be concentrated on the same system, in
> order to provide consistency (at least sharing 'earthRadius' values).
> In that manner, we failed to use 'geod' or affiliate as API.
> We have to say that we're developing on Windows platform (VS 2005), and
> we're lacking sample/tutorial/help for such a use.
> Any help in that manner would be welcome.
> PS: Is there a proper (recommended) way to increase 'ellipsoid' definitions
> (increasing or modifying  'pj_ellps' table, for instance changing
> "sphere"'s radius)?

Note that both traditional proj DMS and GRASS colon format geographic 
coordinates are read.

gie at charon:~/CARTO/geodesic/src$ ./geodesic
geod: help
command: earth -- set earth shape parameters
command: show -- display [arc|earth] structure contents
command: p1 -- set point 1 [lon lat]
command: p2 -- set point 2 [lon lat]
command: open -- open output file
command: vector -- set vector [distance azimuth
command: line -- compute geodesic line
command: circle -- compute circle or arc [2nd azi. del azi.]
command: prec -- set display precision
command: quit -- exit program
command: help -- show this list

For more detailed assistance type "<command> help" for
more detailed information on that command.
geod: earth ellps=WGS84
geod: p1 100d14'33.2W 30d15
geod: p2 -90:17:55 44d7'30n
    Point 1 Lon:  100d14'33.20000"W  Lat:   30d15'00.00000"N  height: 0.000
    Point 2 Lon:   90d17'55.00000"W  Lat:   44d07'30.00000"N  height: 0.000
    Azimuth p1->p2: 26d55'6.60"
    Distance: 1771916.625
    Azimuth p2->p1: 212d58'58.59"
geod: quit
geodesic at completion.   G'day mate!
gie at charon:~/CARTO/geodesic/src$

The whole religious complexion of the modern world is due
to the absence from Jerusalem of a lunatic asylum.
-- Havelock Ellis (1859-1939) British psychologist

More information about the Proj mailing list