GeographicLib  1.21
LocalCartesian.hpp
Go to the documentation of this file.
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