merc.hpp 9.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247
  1. // Boost.Geometry - gis-projections (based on PROJ4)
  2. // Copyright (c) 2008-2015 Barend Gehrels, Amsterdam, the Netherlands.
  3. // This file was modified by Oracle on 2017, 2018.
  4. // Modifications copyright (c) 2017-2018, Oracle and/or its affiliates.
  5. // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle.
  6. // Use, modification and distribution is subject to the Boost Software License,
  7. // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  8. // http://www.boost.org/LICENSE_1_0.txt)
  9. // This file is converted from PROJ4, http://trac.osgeo.org/proj
  10. // PROJ4 is originally written by Gerald Evenden (then of the USGS)
  11. // PROJ4 is maintained by Frank Warmerdam
  12. // PROJ4 is converted to Boost.Geometry by Barend Gehrels
  13. // Last updated version of proj: 5.0.0
  14. // Original copyright notice:
  15. // Permission is hereby granted, free of charge, to any person obtaining a
  16. // copy of this software and associated documentation files (the "Software"),
  17. // to deal in the Software without restriction, including without limitation
  18. // the rights to use, copy, modify, merge, publish, distribute, sublicense,
  19. // and/or sell copies of the Software, and to permit persons to whom the
  20. // Software is furnished to do so, subject to the following conditions:
  21. // The above copyright notice and this permission notice shall be included
  22. // in all copies or substantial portions of the Software.
  23. // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
  24. // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
  25. // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
  26. // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
  27. // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
  28. // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
  29. // DEALINGS IN THE SOFTWARE.
  30. #ifndef BOOST_GEOMETRY_PROJECTIONS_MERC_HPP
  31. #define BOOST_GEOMETRY_PROJECTIONS_MERC_HPP
  32. #include <boost/geometry/srs/projections/impl/base_static.hpp>
  33. #include <boost/geometry/srs/projections/impl/base_dynamic.hpp>
  34. #include <boost/geometry/srs/projections/impl/factory_entry.hpp>
  35. #include <boost/geometry/srs/projections/impl/pj_msfn.hpp>
  36. #include <boost/geometry/srs/projections/impl/pj_param.hpp>
  37. #include <boost/geometry/srs/projections/impl/pj_phi2.hpp>
  38. #include <boost/geometry/srs/projections/impl/pj_tsfn.hpp>
  39. #include <boost/geometry/srs/projections/impl/projects.hpp>
  40. #include <boost/geometry/util/math.hpp>
  41. namespace boost { namespace geometry
  42. {
  43. namespace projections
  44. {
  45. #ifndef DOXYGEN_NO_DETAIL
  46. namespace detail { namespace merc
  47. {
  48. static const double epsilon10 = 1.e-10;
  49. // template class, using CRTP to implement forward/inverse
  50. template <typename T, typename Parameters>
  51. struct base_merc_ellipsoid
  52. : public base_t_fi<base_merc_ellipsoid<T, Parameters>, T, Parameters>
  53. {
  54. inline base_merc_ellipsoid(const Parameters& par)
  55. : base_t_fi<base_merc_ellipsoid<T, Parameters>, T, Parameters>(*this, par)
  56. {}
  57. // FORWARD(e_forward) ellipsoid
  58. // Project coordinates from geographic (lon, lat) to cartesian (x, y)
  59. inline void fwd(T const& lp_lon, T const& lp_lat, T& xy_x, T& xy_y) const
  60. {
  61. static const T half_pi = detail::half_pi<T>();
  62. if (fabs(fabs(lp_lat) - half_pi) <= epsilon10) {
  63. BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
  64. }
  65. xy_x = this->m_par.k0 * lp_lon;
  66. xy_y = - this->m_par.k0 * log(pj_tsfn(lp_lat, sin(lp_lat), this->m_par.e));
  67. }
  68. // INVERSE(e_inverse) ellipsoid
  69. // Project coordinates from cartesian (x, y) to geographic (lon, lat)
  70. inline void inv(T const& xy_x, T const& xy_y, T& lp_lon, T& lp_lat) const
  71. {
  72. if ((lp_lat = pj_phi2(exp(- xy_y / this->m_par.k0), this->m_par.e)) == HUGE_VAL) {
  73. BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
  74. }
  75. lp_lon = xy_x / this->m_par.k0;
  76. }
  77. static inline std::string get_name()
  78. {
  79. return "merc_ellipsoid";
  80. }
  81. };
  82. // template class, using CRTP to implement forward/inverse
  83. template <typename T, typename Parameters>
  84. struct base_merc_spheroid
  85. : public base_t_fi<base_merc_spheroid<T, Parameters>, T, Parameters>
  86. {
  87. inline base_merc_spheroid(const Parameters& par)
  88. : base_t_fi<base_merc_spheroid<T, Parameters>, T, Parameters>(*this, par)
  89. {}
  90. // FORWARD(s_forward) spheroid
  91. // Project coordinates from geographic (lon, lat) to cartesian (x, y)
  92. inline void fwd(T const& lp_lon, T const& lp_lat, T& xy_x, T& xy_y) const
  93. {
  94. static const T half_pi = detail::half_pi<T>();
  95. static const T fourth_pi = detail::fourth_pi<T>();
  96. if (fabs(fabs(lp_lat) - half_pi) <= epsilon10) {
  97. BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
  98. }
  99. xy_x = this->m_par.k0 * lp_lon;
  100. xy_y = this->m_par.k0 * log(tan(fourth_pi + .5 * lp_lat));
  101. }
  102. // INVERSE(s_inverse) spheroid
  103. // Project coordinates from cartesian (x, y) to geographic (lon, lat)
  104. inline void inv(T const& xy_x, T const& xy_y, T& lp_lon, T& lp_lat) const
  105. {
  106. static const T half_pi = detail::half_pi<T>();
  107. lp_lat = half_pi - 2. * atan(exp(-xy_y / this->m_par.k0));
  108. lp_lon = xy_x / this->m_par.k0;
  109. }
  110. static inline std::string get_name()
  111. {
  112. return "merc_spheroid";
  113. }
  114. };
  115. // Mercator
  116. template <typename Params, typename Parameters>
  117. inline void setup_merc(Params const& params, Parameters& par)
  118. {
  119. typedef typename Parameters::type calc_t;
  120. static const calc_t half_pi = detail::half_pi<calc_t>();
  121. calc_t phits=0.0;
  122. int is_phits;
  123. if( (is_phits = pj_param_r<srs::spar::lat_ts>(params, "lat_ts", srs::dpar::lat_ts, phits)) ) {
  124. phits = fabs(phits);
  125. if (phits >= half_pi)
  126. BOOST_THROW_EXCEPTION( projection_exception(error_lat_ts_larger_than_90) );
  127. }
  128. if (par.es != 0.0) { /* ellipsoid */
  129. if (is_phits)
  130. par.k0 = pj_msfn(sin(phits), cos(phits), par.es);
  131. } else { /* sphere */
  132. if (is_phits)
  133. par.k0 = cos(phits);
  134. }
  135. }
  136. }} // namespace detail::merc
  137. #endif // doxygen
  138. /*!
  139. \brief Mercator projection
  140. \ingroup projections
  141. \tparam Geographic latlong point type
  142. \tparam Cartesian xy point type
  143. \tparam Parameters parameter type
  144. \par Projection characteristics
  145. - Cylindrical
  146. - Spheroid
  147. - Ellipsoid
  148. \par Projection parameters
  149. - lat_ts: Latitude of true scale (degrees)
  150. \par Example
  151. \image html ex_merc.gif
  152. */
  153. template <typename T, typename Parameters>
  154. struct merc_ellipsoid : public detail::merc::base_merc_ellipsoid<T, Parameters>
  155. {
  156. template <typename Params>
  157. inline merc_ellipsoid(Params const& params, Parameters const& par)
  158. : detail::merc::base_merc_ellipsoid<T, Parameters>(par)
  159. {
  160. detail::merc::setup_merc(params, this->m_par);
  161. }
  162. };
  163. /*!
  164. \brief Mercator projection
  165. \ingroup projections
  166. \tparam Geographic latlong point type
  167. \tparam Cartesian xy point type
  168. \tparam Parameters parameter type
  169. \par Projection characteristics
  170. - Cylindrical
  171. - Spheroid
  172. - Ellipsoid
  173. \par Projection parameters
  174. - lat_ts: Latitude of true scale (degrees)
  175. \par Example
  176. \image html ex_merc.gif
  177. */
  178. template <typename T, typename Parameters>
  179. struct merc_spheroid : public detail::merc::base_merc_spheroid<T, Parameters>
  180. {
  181. template <typename Params>
  182. inline merc_spheroid(Params const& params, Parameters const& par)
  183. : detail::merc::base_merc_spheroid<T, Parameters>(par)
  184. {
  185. detail::merc::setup_merc(params, this->m_par);
  186. }
  187. };
  188. #ifndef DOXYGEN_NO_DETAIL
  189. namespace detail
  190. {
  191. // Static projection
  192. BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION(srs::spar::proj_merc, merc_spheroid, merc_ellipsoid)
  193. // Factory entry(s)
  194. BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI2(merc_entry, merc_spheroid, merc_ellipsoid)
  195. BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_BEGIN(merc_init)
  196. {
  197. BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(merc, merc_entry)
  198. }
  199. } // namespace detail
  200. #endif // doxygen
  201. } // namespace projections
  202. }} // namespace boost::geometry
  203. #endif // BOOST_GEOMETRY_PROJECTIONS_MERC_HPP