projections/proj/merc.hpp

0001 // Boost.Geometry - gis-projections (based on PROJ4)
0002
0003 // Copyright (c) 2008-2015 Barend Gehrels, Amsterdam, the Netherlands.
0004
0005 // This file was modified by Oracle on 2017, 2018, 2019, 2022.
0006 // Modifications copyright (c) 2017-2022, Oracle and/or its affiliates.
0007 // Contributed and/or modified by Vissarion Fysikopoulos, on behalf of Oracle.
0008 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle.
0009
0010 // Use, modification and distribution is subject to the Boost Software License,
0011 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
0012 // http://www.boost.org/LICENSE_1_0.txt)
0013
0014 // This file is converted from PROJ4, http://trac.osgeo.org/proj
0015 // PROJ4 is originally written by Gerald Evenden (then of the USGS)
0016 // PROJ4 is maintained by Frank Warmerdam
0017 // PROJ4 is converted to Boost.Geometry by Barend Gehrels
0018
0019 // Last updated version of proj: 5.0.0
0020
0021 // Original copyright notice:
0022
0023 // Permission is hereby granted, free of charge, to any person obtaining a
0024 // copy of this software and associated documentation files (the "Software"),
0025 // to deal in the Software without restriction, including without limitation
0026 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
0027 // and/or sell copies of the Software, and to permit persons to whom the
0028 // Software is furnished to do so, subject to the following conditions:
0029
0030 // The above copyright notice and this permission notice shall be included
0031 // in all copies or substantial portions of the Software.
0032
0033 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
0034 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
0035 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
0036 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
0037 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
0038 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
0039 // DEALINGS IN THE SOFTWARE.
0040
0041 #ifndef BOOST_GEOMETRY_PROJECTIONS_MERC_HPP
0042 #define BOOST_GEOMETRY_PROJECTIONS_MERC_HPP
0043
0044 #include <boost/geometry/srs/projections/impl/base_static.hpp>
0045 #include <boost/geometry/srs/projections/impl/base_dynamic.hpp>
0046 #include <boost/geometry/srs/projections/impl/factory_entry.hpp>
0047 #include <boost/geometry/srs/projections/impl/pj_msfn.hpp>
0048 #include <boost/geometry/srs/projections/impl/pj_param.hpp>
0049 #include <boost/geometry/srs/projections/impl/pj_phi2.hpp>
0050 #include <boost/geometry/srs/projections/impl/pj_tsfn.hpp>
0051 #include <boost/geometry/srs/projections/impl/projects.hpp>
0052
0053 #include <boost/geometry/util/math.hpp>
0054
0055 namespace boost { namespace geometry
0056 {
0057
0058 namespace projections
0059 {
0060     #ifndef DOXYGEN_NO_DETAIL
0061     namespace detail { namespace merc
0062     {
0063
0064             static const double epsilon10 = 1.e-10;
0065
0066             template <typename T, typename Parameters>
0067             struct base_merc_ellipsoid
0068             {
0069                 // FORWARD(e_forward)  ellipsoid
0070                 // Project coordinates from geographic (lon, lat) to cartesian (x, y)
0071                 inline void fwd(Parameters const& par, T const& lp_lon, T const& lp_lat, T& xy_x, T& xy_y) const
0072                 {
0073                     static const T half_pi = detail::half_pi<T>();
0074
0075                     if (fabs(fabs(lp_lat) - half_pi) <= epsilon10) {
0076                         BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
0077                     }
0078                     xy_x = par.k0 * lp_lon;
0079                     xy_y = - par.k0 * log(pj_tsfn(lp_lat, sin(lp_lat), par.e));
0080                 }
0081
0082                 // INVERSE(e_inverse)  ellipsoid
0083                 // Project coordinates from cartesian (x, y) to geographic (lon, lat)
0084                 inline void inv(Parameters const& par, T const& xy_x, T const& xy_y, T& lp_lon, T& lp_lat) const
0085                 {
0086                     if ((lp_lat = pj_phi2(exp(- xy_y / par.k0), par.e)) == HUGE_VAL) {
0087                         BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
0088                     }
0089                     lp_lon = xy_x / par.k0;
0090                 }
0091
0092                 static inline std::string get_name()
0093                 {
0094                     return "merc_ellipsoid";
0095                 }
0096
0097             };
0098
0099             template <typename T, typename Parameters>
0100             struct base_merc_spheroid
0101             {
0102                 // FORWARD(s_forward)  spheroid
0103                 // Project coordinates from geographic (lon, lat) to cartesian (x, y)
0104                 inline void fwd(Parameters const& par, T const& lp_lon, T const& lp_lat, T& xy_x, T& xy_y) const
0105                 {
0106                     static const T half_pi = detail::half_pi<T>();
0107                     static const T fourth_pi = detail::fourth_pi<T>();
0108
0109                     if (fabs(fabs(lp_lat) - half_pi) <= epsilon10) {
0110                         BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
0111                     }
0112                     xy_x = par.k0 * lp_lon;
0113                     xy_y = par.k0 * log(tan(fourth_pi + .5 * lp_lat));
0114                 }
0115
0116                 // INVERSE(s_inverse)  spheroid
0117                 // Project coordinates from cartesian (x, y) to geographic (lon, lat)
0118                 inline void inv(Parameters const& par, T const& xy_x, T const& xy_y, T& lp_lon, T& lp_lat) const
0119                 {
0120                     static const T half_pi = detail::half_pi<T>();
0121
0122                     lp_lat = half_pi - 2. * atan(exp(-xy_y / par.k0));
0123                     lp_lon = xy_x / par.k0;
0124                 }
0125
0126                 static inline std::string get_name()
0127                 {
0128                     return "merc_spheroid";
0129                 }
0130
0131             };
0132
0133             // Mercator
0134             template <typename Params, typename Parameters>
0135             inline void setup_merc(Params const& params, Parameters& par)
0136             {
0137                 typedef typename Parameters::type calc_t;
0138                 static const calc_t half_pi = detail::half_pi<calc_t>();
0139
0140                 calc_t phits=0.0;
0141                 int is_phits;
0142
0143                 if( (is_phits = pj_param_r<srs::spar::lat_ts>(params, "lat_ts", srs::dpar::lat_ts, phits)) ) {
0144                     phits = fabs(phits);
0145                     if (phits >= half_pi)
0146                         BOOST_THROW_EXCEPTION( projection_exception(error_lat_ts_larger_than_90) );
0147                 }
0148                 if (par.es != 0.0) { /* ellipsoid */
0149                     if (is_phits)
0150                         par.k0 = pj_msfn(sin(phits), cos(phits), par.es);
0151                 } else { /* sphere */
0152                     if (is_phits)
0153                         par.k0 = cos(phits);
0154                 }
0155             }
0156
0157     }} // namespace detail::merc
0158     #endif // doxygen
0159
0160     /*!
0161         \brief Mercator projection
0162         \ingroup projections
0163         \tparam Geographic latlong point type
0164         \tparam Cartesian xy point type
0165         \tparam Parameters parameter type
0166         \par Projection characteristics
0167          - Cylindrical
0168          - Spheroid
0169          - Ellipsoid
0170         \par Projection parameters
0171          - lat_ts: Latitude of true scale (degrees)
0172         \par Example
0173         \image html ex_merc.gif
0174     */
0175     template <typename T, typename Parameters>
0176     struct merc_ellipsoid : public detail::merc::base_merc_ellipsoid<T, Parameters>
0177     {
0178         template <typename Params>
0179         inline merc_ellipsoid(Params const& params, Parameters & par)
0180         {
0181             detail::merc::setup_merc(params, par);
0182         }
0183     };
0184
0185     /*!
0186         \brief Mercator projection
0187         \ingroup projections
0188         \tparam Geographic latlong point type
0189         \tparam Cartesian xy point type
0190         \tparam Parameters parameter type
0191         \par Projection characteristics
0192          - Cylindrical
0193          - Spheroid
0194          - Ellipsoid
0195         \par Projection parameters
0196          - lat_ts: Latitude of true scale (degrees)
0197         \par Example
0198         \image html ex_merc.gif
0199     */
0200     template <typename T, typename Parameters>
0201     struct merc_spheroid : public detail::merc::base_merc_spheroid<T, Parameters>
0202     {
0203         template <typename Params>
0204         inline merc_spheroid(Params const& params, Parameters & par)
0205         {
0206             detail::merc::setup_merc(params, par);
0207         }
0208     };
0209
0210     /*!
0211         \brief Web Mercator projection
0212         \ingroup projections
0213         \tparam Geographic latlong point type
0214         \tparam Cartesian xy point type
0215         \tparam Parameters parameter type
0216         \par Projection characteristics
0217          - Cylindrical
0218          - Spheroid
0219          - Ellipsoid
0220     */
0221     template <typename T, typename Parameters>
0222     struct webmerc_spheroid : public detail::merc::base_merc_spheroid<T, Parameters>
0223     {
0224         template <typename Params>
0225         inline webmerc_spheroid(Params const&, Parameters & par)
0226         {
0227             par.k0 = 1;
0228         }
0229     };
0230
0231     #ifndef DOXYGEN_NO_DETAIL
0232     namespace detail
0233     {
0234
0235         // Static projection
0236         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI2(srs::spar::proj_merc, merc_spheroid,
0237             merc_ellipsoid)
0238         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_webmerc,
0239             webmerc_spheroid)
0240
0241         // Factory entry(s)
0242         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI2(merc_entry, merc_spheroid,
0243             merc_ellipsoid)
0244         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(webmerc_entry, webmerc_spheroid)
0245
0246         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_BEGIN(merc_init)
0247         {
0248             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(merc, merc_entry)
0249             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(webmerc, webmerc_entry)
0250         }
0251
0252     } // namespace detail
0253     #endif // doxygen
0254
0255 } // namespace projections
0256
0257 }} // namespace boost::geometry
0258
0259 #endif // BOOST_GEOMETRY_PROJECTIONS_MERC_HPP