projections/proj/gn_sinu.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.
0006 // Modifications copyright (c) 2017-2019, Oracle and/or its affiliates.
0007 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle.
0008
0009 // Use, modification and distribution is subject to the Boost Software License,
0010 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
0011 // http://www.boost.org/LICENSE_1_0.txt)
0012
0013 // This file is converted from PROJ4, http://trac.osgeo.org/proj
0014 // PROJ4 is originally written by Gerald Evenden (then of the USGS)
0015 // PROJ4 is maintained by Frank Warmerdam
0016 // PROJ4 is converted to Boost.Geometry by Barend Gehrels
0017
0018 // Last updated version of proj: 5.0.0
0019
0020 // Original copyright notice:
0021
0022 // Permission is hereby granted, free of charge, to any person obtaining a
0023 // copy of this software and associated documentation files (the "Software"),
0024 // to deal in the Software without restriction, including without limitation
0025 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
0026 // and/or sell copies of the Software, and to permit persons to whom the
0027 // Software is furnished to do so, subject to the following conditions:
0028
0029 // The above copyright notice and this permission notice shall be included
0030 // in all copies or substantial portions of the Software.
0031
0032 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
0033 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
0034 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
0035 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
0036 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
0037 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
0038 // DEALINGS IN THE SOFTWARE.
0039
0040 #ifndef BOOST_GEOMETRY_PROJECTIONS_GN_SINU_HPP
0041 #define BOOST_GEOMETRY_PROJECTIONS_GN_SINU_HPP
0042
0043 #include <boost/geometry/srs/projections/impl/aasincos.hpp>
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_mlfn.hpp>
0048 #include <boost/geometry/srs/projections/impl/pj_param.hpp>
0049 #include <boost/geometry/srs/projections/impl/projects.hpp>
0050
0051 #include <boost/geometry/util/math.hpp>
0052
0053 namespace boost { namespace geometry
0054 {
0055
0056 namespace projections
0057 {
0058     #ifndef DOXYGEN_NO_DETAIL
0059     namespace detail { namespace gn_sinu
0060     {
0061
0062             static const double epsilon10 = 1e-10;
0063             static const int max_iter = 8;
0064             static const double loop_tol = 1e-7;
0065
0066             template <typename T>
0067             struct par_gn_sinu_e
0068             {
0069                 detail::en<T> en;
0070             };
0071
0072             template <typename T>
0073             struct par_gn_sinu_s
0074             {
0075                 T m, n, C_x, C_y;
0076             };
0077
0078             /* Ellipsoidal Sinusoidal only */
0079
0080             template <typename T, typename Parameters>
0081             struct base_gn_sinu_ellipsoid
0082             {
0083                 par_gn_sinu_e<T> m_proj_parm;
0084
0085                 // FORWARD(e_forward)  ellipsoid
0086                 // Project coordinates from geographic (lon, lat) to cartesian (x, y)
0087                 inline void fwd(Parameters const& par, T const& lp_lon, T const& lp_lat, T& xy_x, T& xy_y) const
0088                 {
0089                     T s, c;
0090
0091                     xy_y = pj_mlfn(lp_lat, s = sin(lp_lat), c = cos(lp_lat), this->m_proj_parm.en);
0092                     xy_x = lp_lon * c / sqrt(1. - par.es * s * s);
0093                 }
0094
0095                 // INVERSE(e_inverse)  ellipsoid
0096                 // Project coordinates from cartesian (x, y) to geographic (lon, lat)
0097                 inline void inv(Parameters const& par, T const& xy_x, T const& xy_y, T& lp_lon, T& lp_lat) const
0098                 {
0099                     static const T half_pi = detail::half_pi<T>();
0100
0101                     T s;
0102
0103                     if ((s = fabs(lp_lat = pj_inv_mlfn(xy_y, par.es, this->m_proj_parm.en))) < half_pi) {
0104                         s = sin(lp_lat);
0105                         lp_lon = xy_x * sqrt(1. - par.es * s * s) / cos(lp_lat);
0106                     } else if ((s - epsilon10) < half_pi)
0107                         lp_lon = 0.;
0108                     else
0109                         BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
0110                 }
0111                 /* General spherical sinusoidals */
0112
0113                 static inline std::string get_name()
0114                 {
0115                     return "gn_sinu_ellipsoid";
0116                 }
0117
0118             };
0119
0120             template <typename T, typename Parameters>
0121             struct base_gn_sinu_spheroid
0122             {
0123                 par_gn_sinu_s<T> m_proj_parm;
0124
0125                 // FORWARD(s_forward)  sphere
0126                 // Project coordinates from geographic (lon, lat) to cartesian (x, y)
0127                 inline void fwd(Parameters const& , T const& lp_lon, T lp_lat, T& xy_x, T& xy_y) const
0128                 {
0129                     if (this->m_proj_parm.m == 0.0)
0130                         lp_lat = this->m_proj_parm.n != 1. ? aasin(this->m_proj_parm.n * sin(lp_lat)): lp_lat;
0131                     else {
0132                         T k, V;
0133                         int i;
0134
0135                         k = this->m_proj_parm.n * sin(lp_lat);
0136                         for (i = max_iter; i ; --i) {
0137                             lp_lat -= V = (this->m_proj_parm.m * lp_lat + sin(lp_lat) - k) /
0138                                 (this->m_proj_parm.m + cos(lp_lat));
0139                             if (fabs(V) < loop_tol)
0140                                 break;
0141                         }
0142                         if (!i) {
0143                             BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
0144                         }
0145                     }
0146                     xy_x = this->m_proj_parm.C_x * lp_lon * (this->m_proj_parm.m + cos(lp_lat));
0147                     xy_y = this->m_proj_parm.C_y * lp_lat;
0148                 }
0149
0150                 // INVERSE(s_inverse)  sphere
0151                 // Project coordinates from cartesian (x, y) to geographic (lon, lat)
0152                 inline void inv(Parameters const& , T const& xy_x, T xy_y, T& lp_lon, T& lp_lat) const
0153                 {
0154                     xy_y /= this->m_proj_parm.C_y;
0155                     lp_lat = (this->m_proj_parm.m != 0.0) ? aasin((this->m_proj_parm.m * xy_y + sin(xy_y)) / this->m_proj_parm.n) :
0156                         ( this->m_proj_parm.n != 1. ? aasin(sin(xy_y) / this->m_proj_parm.n) : xy_y );
0157                     lp_lon = xy_x / (this->m_proj_parm.C_x * (this->m_proj_parm.m + cos(xy_y)));
0158                 }
0159
0160                 static inline std::string get_name()
0161                 {
0162                     return "gn_sinu_spheroid";
0163                 }
0164
0165             };
0166
0167             template <typename Parameters, typename T>
0168             inline void setup(Parameters& par, par_gn_sinu_s<T>& proj_parm)
0169             {
0170                 par.es = 0;
0171
0172                 proj_parm.C_x = (proj_parm.C_y = sqrt((proj_parm.m + 1.) / proj_parm.n))/(proj_parm.m + 1.);
0173             }
0174
0175
0176             // General Sinusoidal Series
0177             template <typename Params, typename Parameters, typename T>
0178             inline void setup_gn_sinu(Params const& params, Parameters& par, par_gn_sinu_s<T>& proj_parm)
0179             {
0180                 if (pj_param_f<srs::spar::n>(params, "n", srs::dpar::n, proj_parm.n)
0181                  && pj_param_f<srs::spar::m>(params, "m", srs::dpar::m, proj_parm.m)) {
0182                     if (proj_parm.n <= 0 || proj_parm.m < 0)
0183                         BOOST_THROW_EXCEPTION( projection_exception(error_invalid_m_or_n) );
0184                 } else
0185                     BOOST_THROW_EXCEPTION( projection_exception(error_invalid_m_or_n) );
0186
0187                 setup(par, proj_parm);
0188             }
0189
0190             // Sinusoidal (Sanson-Flamsteed)
0191             template <typename Parameters, typename T>
0192             inline void setup_sinu(Parameters const& par, par_gn_sinu_e<T>& proj_parm)
0193             {
0194                 proj_parm.en = pj_enfn<T>(par.es);
0195             }
0196
0197             // Sinusoidal (Sanson-Flamsteed)
0198             template <typename Parameters, typename T>
0199             inline void setup_sinu(Parameters& par, par_gn_sinu_s<T>& proj_parm)
0200             {
0201                 proj_parm.n = 1.;
0202                 proj_parm.m = 0.;
0203                 setup(par, proj_parm);
0204             }
0205
0206             // Eckert VI
0207             template <typename Parameters, typename T>
0208             inline void setup_eck6(Parameters& par, par_gn_sinu_s<T>& proj_parm)
0209             {
0210                 proj_parm.m = 1.;
0211                 proj_parm.n = 2.570796326794896619231321691;
0212                 setup(par, proj_parm);
0213             }
0214
0215             // McBryde-Thomas Flat-Polar Sinusoidal
0216             template <typename Parameters, typename T>
0217             inline void setup_mbtfps(Parameters& par, par_gn_sinu_s<T>& proj_parm)
0218             {
0219                 proj_parm.m = 0.5;
0220                 proj_parm.n = 1.785398163397448309615660845;
0221                 setup(par, proj_parm);
0222             }
0223
0224     }} // namespace detail::gn_sinu
0225     #endif // doxygen
0226
0227     /*!
0228         \brief General Sinusoidal Series projection
0229         \ingroup projections
0230         \tparam Geographic latlong point type
0231         \tparam Cartesian xy point type
0232         \tparam Parameters parameter type
0233         \par Projection characteristics
0234          - Pseudocylindrical
0235          - Spheroid
0236         \par Projection parameters
0237          - m (real)
0238          - n (real)
0239         \par Example
0240         \image html ex_gn_sinu.gif
0241     */
0242     template <typename T, typename Parameters>
0243     struct gn_sinu_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<T, Parameters>
0244     {
0245         template <typename Params>
0246         inline gn_sinu_spheroid(Params const& params, Parameters & par)
0247         {
0248             detail::gn_sinu::setup_gn_sinu(params, par, this->m_proj_parm);
0249         }
0250     };
0251
0252     /*!
0253         \brief Sinusoidal (Sanson-Flamsteed) projection
0254         \ingroup projections
0255         \tparam Geographic latlong point type
0256         \tparam Cartesian xy point type
0257         \tparam Parameters parameter type
0258         \par Projection characteristics
0259          - Pseudocylindrical
0260          - Spheroid
0261          - Ellipsoid
0262         \par Example
0263         \image html ex_sinu.gif
0264     */
0265     template <typename T, typename Parameters>
0266     struct sinu_ellipsoid : public detail::gn_sinu::base_gn_sinu_ellipsoid<T, Parameters>
0267     {
0268         template <typename Params>
0269         inline sinu_ellipsoid(Params const& , Parameters & par)
0270         {
0271             detail::gn_sinu::setup_sinu(par, this->m_proj_parm);
0272         }
0273     };
0274
0275     /*!
0276         \brief Sinusoidal (Sanson-Flamsteed) projection
0277         \ingroup projections
0278         \tparam Geographic latlong point type
0279         \tparam Cartesian xy point type
0280         \tparam Parameters parameter type
0281         \par Projection characteristics
0282          - Pseudocylindrical
0283          - Spheroid
0284          - Ellipsoid
0285         \par Example
0286         \image html ex_sinu.gif
0287     */
0288     template <typename T, typename Parameters>
0289     struct sinu_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<T, Parameters>
0290     {
0291         template <typename Params>
0292         inline sinu_spheroid(Params const& , Parameters & par)
0293         {
0294             detail::gn_sinu::setup_sinu(par, this->m_proj_parm);
0295         }
0296     };
0297
0298     /*!
0299         \brief Eckert VI projection
0300         \ingroup projections
0301         \tparam Geographic latlong point type
0302         \tparam Cartesian xy point type
0303         \tparam Parameters parameter type
0304         \par Projection characteristics
0305          - Pseudocylindrical
0306          - Spheroid
0307         \par Example
0308         \image html ex_eck6.gif
0309     */
0310     template <typename T, typename Parameters>
0311     struct eck6_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<T, Parameters>
0312     {
0313         template <typename Params>
0314         inline eck6_spheroid(Params const& , Parameters & par)
0315         {
0316             detail::gn_sinu::setup_eck6(par, this->m_proj_parm);
0317         }
0318     };
0319
0320     /*!
0321         \brief McBryde-Thomas Flat-Polar Sinusoidal projection
0322         \ingroup projections
0323         \tparam Geographic latlong point type
0324         \tparam Cartesian xy point type
0325         \tparam Parameters parameter type
0326         \par Projection characteristics
0327          - Pseudocylindrical
0328          - Spheroid
0329         \par Example
0330         \image html ex_mbtfps.gif
0331     */
0332     template <typename T, typename Parameters>
0333     struct mbtfps_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<T, Parameters>
0334     {
0335         template <typename Params>
0336         inline mbtfps_spheroid(Params const& , Parameters & par)
0337         {
0338             detail::gn_sinu::setup_mbtfps(par, this->m_proj_parm);
0339         }
0340     };
0341
0342     #ifndef DOXYGEN_NO_DETAIL
0343     namespace detail
0344     {
0345
0346         // Static projection
0347         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_gn_sinu, gn_sinu_spheroid)
0348         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI2(srs::spar::proj_sinu, sinu_spheroid, sinu_ellipsoid)
0349         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_eck6, eck6_spheroid)
0350         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_mbtfps, mbtfps_spheroid)
0351
0352         // Factory entry(s)
0353         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(gn_sinu_entry, gn_sinu_spheroid)
0354         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI2(sinu_entry, sinu_spheroid, sinu_ellipsoid)
0355         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(eck6_entry, eck6_spheroid)
0356         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(mbtfps_entry, mbtfps_spheroid)
0357
0358         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_BEGIN(gn_sinu_init)
0359         {
0360             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(gn_sinu, gn_sinu_entry);
0361             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(sinu, sinu_entry);
0362             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(eck6, eck6_entry);
0363             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(mbtfps, mbtfps_entry);
0364         }
0365
0366     } // namespace detail
0367     #endif // doxygen
0368
0369 } // namespace projections
0370
0371 }} // namespace boost::geometry
0372
0373 #endif // BOOST_GEOMETRY_PROJECTIONS_GN_SINU_HPP
0374