Geo Library

0.5

Author:: Dean Nelson

Version:

0.1 functional beta (01/28/02)

0.2 Functional beta (02/02/02) with test routines

0.3 Functional beta (11/1/03) with geomagnetism routines

0.4 Functional beta (11/20/04) tuned for DLL and VB support

0.5 Functional Beta (01/08/05) Updated World Magnetic Model to 2005

Home

Geodetic Library routines augmented by astrometry routines

Table of Contents

This geodetic library has functions for dealing with many geodesy-based problems found in positioning, pointing, and surveying situations. It is useful to determine absolute position on the earth, pointing vectors, coordinate transformations, and deg/min/sec conversions. It has the following features:

ANSI C code
Sun position
Accurate Azimuth, Elevation, and Range calculation
Cartesian to Polar conversions
It has 23 ellipsoid definitions (can be used worldwide)
Links with the Naval Observatory's Novas library for astronomic calculations
Calculates the earth's magnetic declination of any location and time(2005-2010)

Compiling Instructions

Geo Library was designed to be as flexible as possible. It does not come as a linkable library. You compile the code on your command line or place it in your own library. Geo Library contains these files:

geo.h - Main Geo Library include file
geoEllips.c - Ellipsoid definitions and location descriptor routines
geoPoint.c - Coordinate transformations, polar and Cartesian conversions
geoAstro.c - Wrapper routines for the Novas library used in Sun Position
geoTest.c - test routines for the Geo Library
NOVAS Library - This is the Naval Observatory Vector Astrometry Subroutine library.
geoMag.c - Geo Magnetism routines for the World Magnetism Model of 2005

Novas is located in a subdirectory called novas in the main code tree. When using Novas, you will need to compile in the following files:

There are also header files in that directory that you may want to make your compiler aware of. If, for some reason, Novas is not supplied with your version of Geo Library, you may download it from http://aa.usno.navy.mil/software/novas/novas_info.html

Example compilation with GNU C++

with novas already in a library:

gcc geoTest.c geoAstro.c geoEllips.c geoPoint.c geoMag.c novas\novas.a

with Novas not in a library:

gcc geoTest.c geoAstro.c geoEllips.c geoPoint.c geoMag.c \
      novas\novas.c novas\novascon.c novas\readeph0.c novas\solsys3.c

Example compilation with Borland C++

with novas already in a library:

bcc32 geoTest.c geoAstro.c geoEllips.c geoPoint.c geoMag.c novas\novas.lib

with Novas not in a library:

bcc32 geoTest.c geoAstro.c geoEllips.c geoPoint.c geoMag.c \
      novas\novas.c novas\novascon.c novas\readeph0.c novas\solsys3.c

Example compilation with Microsoft C++

with novas already in a library:

cl geoTest.c geoAstro.c geoEllips.c geoPoint.c geoMag.c novas\novas.lib

with Novas not in a library:

cl geoTest.c geoAstro.c geoEllips.c geoPoint.c geoMag.c \
      novas\novas.c novas\novascon.c novas\readeph0.c novas\solsys3.c

Calling the geo.dll from Visual Basic

Example:

Private Declare Function MagD Lib "geo.dll" Alias "_geoMagGetDecRet@36" _
   (ByVal Latitude As Double, ByVal Longitude As Double, ByVal Altitude As Double, _
    ByVal Month As Integer, ByVal Day As Integer, ByVal Year As Integer) As Double

' Although the path to "geo.dll" could be included, if it is in the system path, the local
' directory, c:\windows or \windows\system32, the path is not required, as it will find it in going 
' through the defaults.

Dim dlat As Double          'Latitude as decimal.  Plus for Northern Latitudes
Dim dlong As Double         'Longitude as decimal.  Minus for West Longitude's
Dim alt As Double           'Altitude in Meters  An average value for a state can be used
Dim mo As Integer           'Month
Dim da As Integer           'Day
Dim yr As Integer           'Year


'_______________________________________________________________
Private Sub cmdDeclination_Click()          'They pushed the calculate button'
'Note, the following values would normally be generated by some program steps
'and checked to be sure they were numeric etc.
dlat = 37.1668611
dlong = -93.35138889
alt = 350
mo = 11
da = 2
yr = 2004

txtResult.Text = MagD(dlat, dlong, alt, mo, da, yr)   'The result is displayed in a text box txtResult

End Sub
'_______________________________________________________________

More examples and definintions are included in geo.bas and geoTest.bas.

Home

Formulas

Basic Geodesy

Ellipsoids

Contrary to popular belief, the Earth is not round. It, like most people, is a bit bigger around the girth. That girth is called the equator. The equator is one of the two axis of the Earth. The other axis is the polar axis.

Ellipsoids contain several defining values:

Major Axis - This is the equatorial axis or radius
Minor Axis - This is the polar axis or radius
Flattening - This is a ratio describing the amount of flattening that is present on the ellipsoid. Flattening is represented by
$f = \frac{a - b}{a}$
Flattening is usually a small number (i.e. 0.0033900753) and is often expressed as a fraction (i.e. $\frac{1}{297}$ ). Inverse Flattening $f^{-1}$ is, of course, the denominator of the fraction.
Eccentricity - also called Eccentricity Squared is represented by
$e^2 = \frac{a^2 - b^2}{a^2}$
Eccentricity Prime - also called Eccentricity Squared Prime is represented by
$e^2_p = \frac{a^2 - b^2}{b^2}$

Ellipsoid Parameters

The Geo Library contains definitions of the following ellipsoids:

       Ellipsoid Name,                       ID,  Major Axis, Inverse Flattening
  ---------------------------------------------------------------------------
 Airy 1830,                                  AA, 6377563.396, 299.3249646
 Australian National,                        AN, 6378160    , 298.25
 Bessel 1841,                                BR, 6377397.155, 299.1528128
 Bessel 1841 (Namibia),                      BN, 6377483.865, 299.1528128
 Clarke 1866,                                CC, 6378206.4  , 294.9786982
 Clarke 1880,                                CD, 6378249.145, 293.465
 Everest (Brunei, E. Malaysia ),             EB, 6377298.556, 300.8017
 Everest 1830,                               EA, 6377276.345, 300.8017
 Everest 1956 (India and Nepal),             EC, 6377301.243, 300.8017
 Everest (Pakistan),                         EF, 6377309.613, 300.8017
 Everest 1948 (W. Malaysia and Singapore),   EE, 6377304.063, 300.8017
 Everest 1969 (W. Malaysia),                 ED, 6377295.664, 300.8017
 Geodetic Reference System 1980,             RF, 6378137    , 298.257222101
 Helmert 1906,                               HE, 6378200    , 298.3
 Hough 1960,                                 HO, 6378270    , 297.0
 Indonesian 1974,                            ID, 6378160    , 298.247
 International 1924,                         IN, 6378388    , 297.0
 Krassovsky 1940,                            KA, 6378245    , 298.3
 Modified Airy,	                             AM, 6377340.189, 299.3249646
 Modified Fischer 1960,                      FA, 6378155    , 298.3
 South American 1969,                        SA, 6378160    , 298.25
 WGS 1972,                                   WD, 6378135    , 298.26
 WGS 1984,                                   WE, 6378137    , 298.257223563

Home

Coordinate Systems

There are several coordinate systems used in the Geo Library.

Latitude ( $\phi$ ), Longitude ( $\lambda$ ), Height () This coordinate system references a point on the Earth's ellipsoid. Latitude and longitude are measured in degrees. Height is measured in meters.

Latitude, Longitude, and height

Earth Fixed Geodetic Coordinates This Cartesian coordinate system uses the center of the earth as a reference point. The coordinates represent an X,Y,Z coordinate system. However, in Geo Library, these coordinates are not called X,Y,Z. They are called E,F,G. This is to help distinguish between the tangential plane X,Y,Z coordinate systems.

Earth Fixed Geodetic Coordinates

Range (), Azimuth ( $\alpha$ ), Elevation ( $\epsilon$ ) This polar coordinate system originates from a point on the Earth's ellipsoid. Range () is measured in meters and represents the slant range from the origin to the point. Azimuth ( $\alpha$ ) is measured clockwise from north (zero degrees) in degrees. $\epsilon$ is measured from the tangential plane upward. Also in degrees.

Range, Azimuth, and Elevation Coordinates

X, Y, Z (also known as ENU) This coordinate system also originates at point on the Earth's ellipsoid. X (also known as East) is measured in meters from the point eastward. Y (also known as North) is measured in meters from the point northward. Z (also known as Up) is measured in meters, in the tangential plane, in an upward direction.

X,Y,Z Cartesian Coordinates

Home

Conversion Formulas

Range , Azimuth $\alpha$ , and Elevation $\epsilon$ to X,Y,Z If the slant range , azimuth $\alpha$ , and elevation $\epsilon$ from the origin to a point are known, topocentric Cartesian coordinates can be computed in the following manner:

$X = r \; cos \, \epsilon \; sin \, \alpha$

$Y = r \; cos \, \epsilon \; cos \, \alpha$

$Z = r \; sin \, \epsilon$

Conversion of Range, Azimuth, and Elevation to E,F,G Geocentric Coordinates

X,Y,Z to Range , Azimuth $\alpha$ , and Elevation $\epsilon$ If Cartesian coordinates X,Y,Z are known, then slant range , azimuth $\alpha$ , and elevation $\epsilon$ can be determined this way:

$r = \sqrt{X^2 + Y^2 + Z^2}$

$\alpha = tan^{-1}\left\{\frac{X}{Y}\right\} (if \; \alpha < 0, \; then \; \alpha=\alpha +2\pi)$

$\epsilon = tan^{-1}\left\{\frac{Z}{ \sqrt{X^2 + Y^2}}\right\}$

Latitude $\phi$ , Longitude $\lambda$ , Height to E,F,G Geodetic coordinates, Latitude, longitude, and height can be converted to geocentric coordinates through this method:

$Radius \; of \; curvature: \; N=\frac{a}{\sqrt{1-e^2sin^2\phi}}$

$E=(N+h) \; cos\phi \; cos\lambda$

$F=(N+h) \; cos\phi \; sin\lambda$

$G=[N(1-e^2)+h] \; sin\phi$

E,F,G to Latitude $\phi$ , Longitude $\lambda$ , Height Geocentric coordinates E,F,G can be converted to latitude $\phi$ , longitude $\lambda$ , and height in the following manner:

$p = \sqrt{E^2+F^2}$

$u = tan^{-1}\left\{ \frac{G}{p} \; \frac{a}{b} \right\}$

$\phi = tan^{-1} \left\{ \frac{G+e^2_p \; b \; sin^3u}{p - e^2 \; a \; cos^3u} \right\}$

$u_p = tan^{-1} (( 1-f)tan \; \phi)$

$h = \sqrt{(p-a \; cos \; u_p)^2 + (G-b \; sin \; u_p)^2}$

$\lambda = tan^{-1} \left\{ \frac{E}{F} \right\}$

Note:: This routine will be exact only for WGS84 coordinates. All other datums will be slightly off.

Geomagnetic Routines

The modeling of the earth's magnetic field is called geomagnetism. The most current application of this model is magnetic declination - deflection of the compass needle off true north. The earth's magnetic field is in constant motion. Thus, a stable model needed to be developed so that one might be able to correct for the changes of the field in time.

By using the geomagnetic routines in this library, true north can be determined quite easily.

How do I ....

Here are some common questions:

How do I convert Range, Azimuth, and Elevation to XYZ coordinates?
- use geoRae2Xyz()
How do I convert XYZ coordinates to Range, Azimuth, and Elevation?
- use geoXyz2Rae()
How do I get a Range, Azimuth, and Elevation from one point to another and all I have are their geodetic coordinates?
- use geoInitLocation() to establish a base location and then again to establish the point location
- use geoEfg2XyzDiff() to get the XYZ coordinates of the point (with the base as the origin).
- use geoXyz2Rae() to get your polar coordinates
How do I get X,Y,Z from one point to another and all I have are their geodetic coordinates?
- use geoInitLocation() to establish a base location and then again to establish the point location
- use geoEfg2XyzDiff() to get the XYZ coordinates of the point (with the base as the origin).
How do I get earth centered coordinates (EFG) from Range, Azimuth, and Elevation?
- use geoInitLocation() to establish a base location
- use geoRae2Efg()
How do I get the Latitude, Longitude and height of the point indicated by Range, Azimuth, and Elevation?
- use geoInitLocation() to establish a base location
- use geoRae2Efg() to get earth centered coordinates
- use geoEfg2Llh() to get Latitude, Longitude, and height.
How do I convert radians to degrees, minutes, and seconds?
- use geoRads2Dms()
How do I convert degrees, minutes, and seconds to radians?
- use geoDms2Rads()
How do I find out the magnetic declination given a Latitude and Longitude?
- use geoMag() to get every blessed geomagnetic variable under the sun, or...
- use geoMagGetDec() to get just the magnetic declination, or ...
- use geoMagGetDecRet() to return a double precision value declination (also used in VB)
How do I get the position of the sun right now?
- use geoInitLocation() to establish a base location
- use geoSunAzElNow() to get the Azimuth and Elevation
How do I get the sun position for any day?
- include time.h
- use geoInitLocation() to establish a base location
- fill in a tm structure with the date and time you want
- use geoSunAzEltm()
How do I get the sun position for any day when I only have a Julian date?
- use geoInitLocation() to establish a base location
- compute or obtain the Julian date/time
- use geoSunAzElJD() to get the Azimuth and Elevation
How do I get a date with this Julian?
- get a life!

Home

Generated on Wed Mar 2 20:30:43 2005 for GeoStarslib - Geodetic Library by

1.4.1