X-Git-Url: https://git.mxchange.org/?a=blobdiff_plain;f=simgear%2Fmath%2Fsg_geodesy.hxx;h=2ee46e0826cc87deb24c4c5ae0aff6dd68893401;hb=aa859c488f21d0b1f8d54c657e23a738dcadacca;hp=b9f7d6db17ff9ced97838bfc5e78025ac41b756d;hpb=7d7e41dacc78c3e997b77a8b0be6b4fed8c1adf2;p=simgear.git diff --git a/simgear/math/sg_geodesy.hxx b/simgear/math/sg_geodesy.hxx index b9f7d6db..2ee46e08 100644 --- a/simgear/math/sg_geodesy.hxx +++ b/simgear/math/sg_geodesy.hxx @@ -1,175 +1,156 @@ -// sg_geodesy.hxx -- routines to convert between geodetic and geocentric -// coordinate systems. -// -// Copied and adapted directly from LaRCsim/ls_geodesy.c -// -// See below for the complete original LaRCsim comments. -// -// $Id$ - - #ifndef _SG_GEODESY_HXX #define _SG_GEODESY_HXX +#include "SGMath.hxx" -#ifndef __cplusplus -# error This library requires C++ -#endif - - -#include -#include - - -// sgGeocToGeod(lat_geoc, radius, *lat_geod, *alt, *sea_level_r) -// INPUTS: -// lat_geoc Geocentric latitude, radians, + = North -// radius C.G. radius to earth center (meters) -// -// OUTPUTS: -// lat_geod Geodetic latitude, radians, + = North -// alt C.G. altitude above mean sea level (meters) -// sea_level_r radius from earth center to sea level at -// local vertical (surface normal) of C.G. (meters) - -void sgGeocToGeod( double lat_geoc, double radius, double - *lat_geod, double *alt, double *sea_level_r ); - - -// sgGeodToGeoc( lat_geod, alt, *sl_radius, *lat_geoc ) -// INPUTS: -// lat_geod Geodetic latitude, radians, + = North -// alt C.G. altitude above mean sea level (meters) -// -// OUTPUTS: -// sl_radius SEA LEVEL radius to earth center (meters) -// (add Altitude to get true distance from earth center. -// lat_geoc Geocentric latitude, radians, + = North -// - -void sgGeodToGeoc( double lat_geod, double alt, double *sl_radius, - double *lat_geoc ); - - -// convert a geodetic point lon(radians), lat(radians), elev(meter) to -// a cartesian point - -inline Point3D sgGeodToCart(const Point3D& geod) { - double gc_lon, gc_lat, sl_radius; - - // printf("A geodetic point is (%.2f, %.2f, %.2f)\n", - // geod[0], geod[1], geod[2]); +// Compatibility header. +// Please use the SGGeodesy and SGMath functions directly. - gc_lon = geod.lon(); - sgGeodToGeoc(geod.lat(), geod.radius(), &sl_radius, &gc_lat); - - // printf("A geocentric point is (%.2f, %.2f, %.2f)\n", gc_lon, - // gc_lat, sl_radius+geod[2]); - - Point3D pp = Point3D( gc_lon, gc_lat, sl_radius + geod.radius()); - return sgPolarToCart3d(pp); +/** + * Convert from geodetic coordinates to geocentric coordinates. + * WARNING: this function is non-reversible. Due to the fact that + * "up" is a different direction for geocentric and geodetic frames, + * you can not simply add your "alt" parameter to the "sl_radius" + * result and get back (via sgGeodToGeoc()) to the coordinates you + * started with. The error under normal conditions will be of + * centimeter order; whether that is important or not is application + * dependent. Consider using sgGeodToCart() instead. + * + * @param lat_geod (in) Geodetic latitude, radians, + = North + * @param alt (in) C.G. altitude above mean sea level (meters) + * @param sl_radius (out) SEA LEVEL radius to earth center (meters) + * @param lat_geoc (out) Geocentric latitude, radians, + = North + */ +inline void sgGeodToGeoc(double lat_geod, double alt, + double *sl_radius, double *lat_geoc) +{ + SGVec3 cart; + SGGeod geod = SGGeod::fromRadM(0, lat_geod, alt); + SGGeodesy::SGGeodToCart(geod, cart); + SGGeoc geoc; + SGGeodesy::SGCartToGeoc(cart, geoc); + *lat_geoc = geoc.getLatitudeRad(); + *sl_radius = SGGeodesy::SGGeodToSeaLevelRadius(geod); } -// given, alt, lat1, lon1, az1 and distance (s), calculate lat2, lon2 -// and az2. Lat, lon, and azimuth are in degrees. distance in meters -int geo_direct_wgs_84 ( double alt, double lat1, double lon1, double az1, - double s, double *lat2, double *lon2, double *az2 ); - - -// given alt, lat1, lon1, lat2, lon2, calculate starting and ending -// az1, az2 and distance (s). Lat, lon, and azimuth are in degrees. -// distance in meters -int geo_inverse_wgs_84( double alt, double lat1, double lon1, double lat2, - double lon2, double *az1, double *az2, double *s ); - - -/*************************************************************************** - - TITLE: ls_geodesy - ----------------------------------------------------------------------------- - - FUNCTION: Converts geocentric coordinates to geodetic positions - ----------------------------------------------------------------------------- - - MODULE STATUS: developmental - ----------------------------------------------------------------------------- - - GENEALOGY: Written as part of LaRCSim project by E. B. Jackson - ----------------------------------------------------------------------------- - - DESIGNED BY: E. B. Jackson - - CODED BY: E. B. Jackson - - MAINTAINED BY: E. B. Jackson - ----------------------------------------------------------------------------- - - MODIFICATION HISTORY: - - DATE PURPOSE BY - - 930208 Modified to avoid singularity near polar region. EBJ - 930602 Moved backwards calcs here from ls_step. EBJ - 931214 Changed erroneous Latitude and Altitude variables to - *lat_geod and *alt in routine ls_geoc_to_geod. EBJ - 940111 Changed header files from old ls_eom.h style to ls_types, - and ls_constants. Also replaced old DATA type with new - SCALAR type. EBJ - - CURRENT RCS HEADER: - -$Header$ - - * Revision 1.5 1994/01/11 18:47:05 bjax - * Changed include files to use types and constants, not ls_eom.h - * Also changed DATA type to SCALAR type. +/** + * Convert a cartesian point to a geodetic lat/lon/altitude. * - * Revision 1.4 1993/12/14 21:06:47 bjax - * Removed global variable references Altitude and Latitude. EBJ - * - * Revision 1.3 1993/06/02 15:03:40 bjax - * Made new subroutine for calculating geodetic to geocentric; changed name - * of forward conversion routine from ls_geodesy to ls_geoc_to_geod. - * - ----------------------------------------------------------------------------- - - REFERENCES: - - [ 1] Stevens, Brian L.; and Lewis, Frank L.: "Aircraft - Control and Simulation", Wiley and Sons, 1992. - ISBN 0-471-61397-5 - - ----------------------------------------------------------------------------- - - CALLED BY: ls_aux - ----------------------------------------------------------------------------- - - CALLS TO: - ----------------------------------------------------------------------------- - - INPUTS: - lat_geoc Geocentric latitude, radians, + = North - radius C.G. radius to earth center, ft - ----------------------------------------------------------------------------- + * @param xyz (in) Pointer to cartesian point. + * @param lat (out) Latitude, in radians + * @param lon (out) Longitude, in radians + * @param alt (out) Altitude, in meters above the WGS84 ellipsoid + */ +inline void sgCartToGeod(const double* xyz, double* lat, double* lon, double* alt) +{ + SGGeod geod; + SGGeodesy::SGCartToGeod(SGVec3(xyz), geod); + *lat = geod.getLatitudeRad(); + *lon = geod.getLongitudeRad(); + *alt = geod.getElevationM(); +} - OUTPUTS: - lat_geod Geodetic latitude, radians, + = North - alt C.G. altitude above mean sea level, ft - sea_level_r radius from earth center to sea level at - local vertical (surface normal) of C.G. +/** + * Convert a geodetic lat/lon/altitude to a cartesian point. + * + * @param lat (in) Latitude, in radians + * @param lon (in) Longitude, in radians + * @param alt (in) Altitude, in meters above the WGS84 ellipsoid + * @param xyz (out) Pointer to cartesian point. + */ +inline void sgGeodToCart(double lat, double lon, double alt, double* xyz) +{ + SGVec3 cart; + SGGeodesy::SGGeodToCart(SGGeod::fromRadM(lon, lat, alt), cart); + xyz[0] = cart(0); + xyz[1] = cart(1); + xyz[2] = cart(2); +} ---------------------------------------------------------------------------*/ +/** + * Given a starting position and an offset radial and distance, + * calculate an ending positon on a wgs84 ellipsoid. + * @param alt (in) meters (unused) + * @param lat1 (in) degrees + * @param lon1 (in) degrees + * @param az1 (in) degrees + * @param s (in) distance in meters + * @param lat2 (out) degrees + * @param lon2 (out) degrees + * @param az2 (out) return course in degrees + */ +inline int geo_direct_wgs_84 ( double lat1, double lon1, double az1, + double s, double *lat2, double *lon2, + double *az2 ) +{ + SGGeod p2; + if (!SGGeodesy::direct(SGGeod::fromDeg(lon1, lat1), az1, s, p2, *az2)) + return 1; + *lat2 = p2.getLatitudeDeg(); + *lon2 = p2.getLongitudeDeg(); + return 0; +} +inline int geo_direct_wgs_84 ( double alt, double lat1, + double lon1, double az1, + double s, double *lat2, double *lon2, + double *az2 ) +{ return geo_direct_wgs_84(lat1, lon1, az1, s, lat2, lon2, az2); } + +/** + * Given a starting position and an offset radial and distance, + * calculate an ending positon on a wgs84 ellipsoid. + * @param p1 (in) geodetic position + * @param az1 (in) degrees + * @param s (in) distance in meters + * @param p2 (out) geodetic position + * @param az2 (out) return course in degrees + */ +inline int geo_direct_wgs_84(const SGGeod& p1, double az1, + double s, SGGeod& p2, double *az2 ) +{ + return !SGGeodesy::direct(p1, az1, s, p2, *az2); +} +/** + * Given an altitude and two sets of (lat, lon) calculate great circle + * distance between them as well as the starting and ending azimuths. + * @param alt (in) meters (unused) + * @param lat1 (in) degrees + * @param lon1 (in) degrees + * @param lat2 (in) degrees + * @param lon2 (in) degrees + * @param az1 (out) start heading degrees + * @param az2 (out) end heading degrees + * @param s (out) distance meters + */ +inline int geo_inverse_wgs_84( double lat1, double lon1, double lat2, + double lon2, double *az1, double *az2, + double *s ) +{ + return !SGGeodesy::inverse(SGGeod::fromDeg(lon1, lat1), + SGGeod::fromDeg(lon2, lat2), *az1, *az2, *s); +} +inline int geo_inverse_wgs_84( double alt, double lat1, + double lon1, double lat2, + double lon2, double *az1, double *az2, + double *s ) +{ return geo_inverse_wgs_84(lat1, lon1, lat2, lon2, az1, az2, s); } + + +/** + * Given an altitude and two sets of (lat, lon) calculate great circle + * distance between them as well as the starting and ending azimuths. + * @param p1 (in) first position + * @param p2 (in) fsecond position + * @param az1 (out) start heading degrees + * @param az2 (out) end heading degrees + * @param s (out) distance meters + */ +inline int geo_inverse_wgs_84(const SGGeod& p1, const SGGeod& p2, + double *az1, double *az2, double *s ) +{ + return !SGGeodesy::inverse(p1, p2, *az1, *az2, *s); +} #endif // _SG_GEODESY_HXX