GeographicLib
1.21
|
00001 /** 00002 * \file LocalCartesian.hpp 00003 * \brief Header for GeographicLib::LocalCartesian class 00004 * 00005 * Copyright (c) Charles Karney (2008-2011) <charles@karney.com> and licensed 00006 * under the MIT/X11 License. For more information, see 00007 * http://geographiclib.sourceforge.net/ 00008 **********************************************************************/ 00009 00010 #if !defined(GEOGRAPHICLIB_LOCALCARTESIAN_HPP) 00011 #define GEOGRAPHICLIB_LOCALCARTESIAN_HPP \ 00012 "$Id: 133014a30695adf3bdea9f4b52613fab458c505a $" 00013 00014 #include <GeographicLib/Geocentric.hpp> 00015 #include <GeographicLib/Constants.hpp> 00016 00017 namespace GeographicLib { 00018 00019 /** 00020 * \brief Local Cartesian coordinates 00021 * 00022 * Convert between geodetic coordinates latitude = \e lat, longitude = \e 00023 * lon, height = \e h (measured vertically from the surface of the ellipsoid) 00024 * to local cartesian coordinates (\e x, \e y, \e z). The origin of local 00025 * cartesian coordinate system is at \e lat = \e lat0, \e lon = \e lon0, \e h 00026 * = \e h0. The \e z axis is normal to the ellipsoid; the \e y axis points 00027 * due north. The plane \e z = - \e h0 is tangent to the ellipsoid. 00028 * 00029 * The conversions all take place via geocentric coordinates using a 00030 * Geocentric object (by default Geocentric::WGS84). 00031 * 00032 * Example of use: 00033 * \include example-LocalCartesian.cpp 00034 * 00035 * <a href="CartConvert.1.html">CartConvert</a> is a command-line utility 00036 * providing access to the functionality of Geocentric and LocalCartesian. 00037 **********************************************************************/ 00038 00039 class GEOGRAPHIC_EXPORT LocalCartesian { 00040 private: 00041 typedef Math::real real; 00042 static const size_t dim_ = 3; 00043 static const size_t dim2_ = dim_ * dim_; 00044 Geocentric _earth; 00045 real _lat0, _lon0, _h0; 00046 real _x0, _y0, _z0, _r[dim2_]; 00047 void IntForward(real lat, real lon, real h, real& x, real& y, real& z, 00048 real M[dim2_]) const throw(); 00049 void IntReverse(real x, real y, real z, real& lat, real& lon, real& h, 00050 real M[dim2_]) const throw(); 00051 void MatrixMultiply(real M[dim2_]) const throw(); 00052 public: 00053 00054 /** 00055 * Constructor setting the origin. 00056 * 00057 * @param[in] lat0 latitude at origin (degrees). 00058 * @param[in] lon0 longitude at origin (degrees). 00059 * @param[in] h0 height above ellipsoid at origin (meters); default 0. 00060 * @param[in] earth Geocentric object for the transformation; default 00061 * Geocentric::WGS84. 00062 **********************************************************************/ 00063 LocalCartesian(real lat0, real lon0, real h0 = 0, 00064 const Geocentric& earth = Geocentric::WGS84) throw() 00065 : _earth(earth) 00066 { Reset(lat0, lon0, h0); } 00067 00068 /** 00069 * Default constructor. 00070 * 00071 * @param[in] earth Geocentric object for the transformation; default 00072 * Geocentric::WGS84. 00073 * 00074 * Sets \e lat0 = 0, \e lon0 = 0, \e h0 = 0. 00075 **********************************************************************/ 00076 explicit LocalCartesian(const Geocentric& earth = Geocentric::WGS84) 00077 throw() 00078 : _earth(earth) 00079 { Reset(real(0), real(0), real(0)); } 00080 00081 /** 00082 * Reset the origin. 00083 * 00084 * @param[in] lat0 latitude at origin (degrees). 00085 * @param[in] lon0 longitude at origin (degrees). 00086 * @param[in] h0 height above ellipsoid at origin (meters); default 0. 00087 **********************************************************************/ 00088 void Reset(real lat0, real lon0, real h0 = 0) 00089 throw(); 00090 00091 /** 00092 * Convert from geodetic to local cartesian coordinates. 00093 * 00094 * @param[in] lat latitude of point (degrees). 00095 * @param[in] lon longitude of point (degrees). 00096 * @param[in] h height of point above the ellipsoid (meters). 00097 * @param[out] x local cartesian coordinate (meters). 00098 * @param[out] y local cartesian coordinate (meters). 00099 * @param[out] z local cartesian coordinate (meters). 00100 * 00101 * \e lat should be in the range [-90, 90]; \e lon and \e lon0 should be in 00102 * the range [-180, 360]. 00103 **********************************************************************/ 00104 void Forward(real lat, real lon, real h, real& x, real& y, real& z) 00105 const throw() { 00106 IntForward(lat, lon, h, x, y, z, NULL); 00107 } 00108 00109 /** 00110 * Convert from geodetic to local cartesian coordinates and return rotation 00111 * matrix. 00112 * 00113 * @param[in] lat latitude of point (degrees). 00114 * @param[in] lon longitude of point (degrees). 00115 * @param[in] h height of point above the ellipsoid (meters). 00116 * @param[out] x local cartesian coordinate (meters). 00117 * @param[out] y local cartesian coordinate (meters). 00118 * @param[out] z local cartesian coordinate (meters). 00119 * @param[out] M if the length of the vector is 9, fill with the rotation 00120 * matrix in row-major order. 00121 * 00122 * Let \e v be a unit vector located at (\e lat, \e lon, \e h). We can 00123 * express \e v as \e column vectors in one of two ways 00124 * - in east, north, up coordinates (where the components are relative to a 00125 * local coordinate system at (\e lat, \e lon, \e h)); call this 00126 * representation \e v1. 00127 * - in \e x, \e y, \e z coordinates (where the components are relative to 00128 * the local coordinate system at (\e lat0, \e lon0, \e h0)); call this 00129 * representation \e v0. 00130 * . 00131 * Then we have \e v0 = \e M . \e v1. 00132 **********************************************************************/ 00133 void Forward(real lat, real lon, real h, real& x, real& y, real& z, 00134 std::vector<real>& M) 00135 const throw() { 00136 if (M.end() == M.begin() + dim2_) { 00137 real t[dim2_]; 00138 IntForward(lat, lon, h, x, y, z, t); 00139 copy(t, t + dim2_, M.begin()); 00140 } else 00141 IntForward(lat, lon, h, x, y, z, NULL); 00142 } 00143 00144 /** 00145 * Convert from local cartesian to geodetic coordinates. 00146 * 00147 * @param[in] x local cartesian coordinate (meters). 00148 * @param[in] y local cartesian coordinate (meters). 00149 * @param[in] z local cartesian coordinate (meters). 00150 * @param[out] lat latitude of point (degrees). 00151 * @param[out] lon longitude of point (degrees). 00152 * @param[out] h height of point above the ellipsoid (meters). 00153 * 00154 * The value of \e lon returned is in the range [-180, 180). 00155 **********************************************************************/ 00156 void Reverse(real x, real y, real z, real& lat, real& lon, real& h) 00157 const throw() { 00158 IntReverse(x, y, z, lat, lon, h, NULL); 00159 } 00160 00161 /** 00162 * Convert from local cartesian to geodetic coordinates and return rotation 00163 * matrix. 00164 * 00165 * @param[in] x local cartesian coordinate (meters). 00166 * @param[in] y local cartesian coordinate (meters). 00167 * @param[in] z local cartesian coordinate (meters). 00168 * @param[out] lat latitude of point (degrees). 00169 * @param[out] lon longitude of point (degrees). 00170 * @param[out] h height of point above the ellipsoid (meters). 00171 * @param[out] M if the length of the vector is 9, fill with the rotation 00172 * matrix in row-major order. 00173 * 00174 * Let \e v be a unit vector located at (\e lat, \e lon, \e h). We can 00175 * express \e v as \e column vectors in one of two ways 00176 * - in east, north, up coordinates (where the components are relative to a 00177 * local coordinate system at (\e lat, \e lon, \e h)); call this 00178 * representation \e v1. 00179 * - in \e x, \e y, \e z coordinates (where the components are relative to 00180 * the local coordinate system at (\e lat0, \e lon0, \e h0)); call this 00181 * representation \e v0. 00182 * . 00183 * Then we have \e v1 = \e M^T . \e v0, where \e M^T is the transpose of \e 00184 * M. 00185 **********************************************************************/ 00186 void Reverse(real x, real y, real z, real& lat, real& lon, real& h, 00187 std::vector<real>& M) 00188 const throw() { 00189 if (M.end() == M.begin() + dim2_) { 00190 real t[dim2_]; 00191 IntReverse(x, y, z, lat, lon, h, t); 00192 copy(t, t + dim2_, M.begin()); 00193 } else 00194 IntReverse(x, y, z, lat, lon, h, NULL); 00195 } 00196 00197 /** \name Inspector functions 00198 **********************************************************************/ 00199 ///@{ 00200 /** 00201 * @return latitude of the origin (degrees). 00202 **********************************************************************/ 00203 Math::real LatitudeOrigin() const throw() { return _lat0; } 00204 00205 /** 00206 * @return longitude of the origin (degrees). 00207 **********************************************************************/ 00208 Math::real LongitudeOrigin() const throw() { return _lon0; } 00209 00210 /** 00211 * @return height of the origin (meters). 00212 **********************************************************************/ 00213 Math::real HeightOrigin() const throw() { return _h0; } 00214 00215 /** 00216 * @return \e a the equatorial radius of the ellipsoid (meters). This is 00217 * the value of \e a inherited from the Geocentric object used in the 00218 * constructor. 00219 **********************************************************************/ 00220 Math::real MajorRadius() const throw() { return _earth.MajorRadius(); } 00221 00222 /** 00223 * @return \e f the flattening of the ellipsoid. This is the value 00224 * inherited from the Geocentric object used in the constructor. 00225 **********************************************************************/ 00226 Math::real Flattening() const throw() { return _earth.Flattening(); } 00227 ///@} 00228 00229 /// \cond SKIP 00230 /** 00231 * <b>DEPRECATED</b> 00232 * @return \e r the inverse flattening of the ellipsoid. 00233 **********************************************************************/ 00234 Math::real InverseFlattening() const throw() 00235 { return _earth.InverseFlattening(); } 00236 /// \endcond 00237 }; 00238 00239 } // namespace GeographicLib 00240 00241 #endif // GEOGRAPHICLIB_LOCALCARTESIAN_HPP