Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2025-01-18 09:35:46 

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_SCONICS_HPP
0041 #define BOOST_GEOMETRY_PROJECTIONS_SCONICS_HPP
0042 
0043 
0044 #include <boost/geometry/util/math.hpp>
0045 #include <boost/math/special_functions/hypot.hpp>
0046 
0047 #include <boost/geometry/srs/projections/impl/base_static.hpp>
0048 #include <boost/geometry/srs/projections/impl/base_dynamic.hpp>
0049 #include <boost/geometry/srs/projections/impl/factory_entry.hpp>
0050 #include <boost/geometry/srs/projections/impl/pj_param.hpp>
0051 #include <boost/geometry/srs/projections/impl/projects.hpp>
0052 
0053 namespace boost { namespace geometry
0054 {
0055 
0056 namespace projections
0057 {
0058     #ifndef DOXYGEN_NO_DETAIL
0059     namespace detail { namespace sconics
0060     {
0061 
0062             enum proj_type {
0063                 proj_euler  = 0,
0064                 proj_murd1  = 1,
0065                 proj_murd2  = 2,
0066                 proj_murd3  = 3,
0067                 proj_pconic = 4,
0068                 proj_tissot = 5,
0069                 proj_vitk1  = 6
0070             };
0071             static const double epsilon10 = 1.e-10;
0072             static const double epsilon = 1e-10;
0073 
0074             template <typename T>
0075             struct par_sconics
0076             {
0077                 T   n;
0078                 T   rho_c;
0079                 T   rho_0;
0080                 T   sig;
0081                 T   c1, c2;
0082                 proj_type type;
0083             };
0084 
0085             /* get common factors for simple conics */
0086             template <typename Params, typename T>
0087             inline int phi12(Params const& params, par_sconics<T>& proj_parm, T *del)
0088             {
0089                 T p1, p2;
0090                 int err = 0;
0091 
0092                 if (!pj_param_r<srs::spar::lat_1>(params, "lat_1", srs::dpar::lat_1, p1) ||
0093                     !pj_param_r<srs::spar::lat_2>(params, "lat_2", srs::dpar::lat_2, p2)) {
0094                     err = -41;
0095                 } else {
0096                     //p1 = pj_get_param_r(par.params, "lat_1"); // set above
0097                     //p2 = pj_get_param_r(par.params, "lat_2"); // set above
0098                     *del = 0.5 * (p2 - p1);
0099                     proj_parm.sig = 0.5 * (p2 + p1);
0100                     err = (fabs(*del) < epsilon || fabs(proj_parm.sig) < epsilon) ? -42 : 0;
0101                 }
0102                 return err;
0103             }
0104 
0105             template <typename T, typename Parameters>
0106             struct base_sconics_spheroid
0107             {
0108                 par_sconics<T> m_proj_parm;
0109 
0110                 // FORWARD(s_forward)  spheroid
0111                 // Project coordinates from geographic (lon, lat) to cartesian (x, y)
0112                 inline void fwd(Parameters const& , T lp_lon, T const& lp_lat, T& xy_x, T& xy_y) const
0113                 {
0114                     T rho;
0115 
0116                     switch (this->m_proj_parm.type) {
0117                     case proj_murd2:
0118                         rho = this->m_proj_parm.rho_c + tan(this->m_proj_parm.sig - lp_lat);
0119                         break;
0120                     case proj_pconic:
0121                         rho = this->m_proj_parm.c2 * (this->m_proj_parm.c1 - tan(lp_lat - this->m_proj_parm.sig));
0122                         break;
0123                     default:
0124                         rho = this->m_proj_parm.rho_c - lp_lat;
0125                         break;
0126                     }
0127                     xy_x = rho * sin( lp_lon *= this->m_proj_parm.n );
0128                     xy_y = this->m_proj_parm.rho_0 - rho * cos(lp_lon);
0129                 }
0130 
0131                 // INVERSE(s_inverse)  ellipsoid & spheroid
0132                 // Project coordinates from cartesian (x, y) to geographic (lon, lat)
0133                 inline void inv(Parameters const& , T xy_x, T xy_y, T& lp_lon, T& lp_lat) const
0134                 {
0135                     T rho;
0136 
0137                     rho = boost::math::hypot(xy_x, xy_y = this->m_proj_parm.rho_0 - xy_y);
0138                     if (this->m_proj_parm.n < 0.) {
0139                         rho = - rho;
0140                         xy_x = - xy_x;
0141                         xy_y = - xy_y;
0142                     }
0143 
0144                     lp_lon = atan2(xy_x, xy_y) / this->m_proj_parm.n;
0145 
0146                     switch (this->m_proj_parm.type) {
0147                     case proj_pconic:
0148                         lp_lat = atan(this->m_proj_parm.c1 - rho / this->m_proj_parm.c2) + this->m_proj_parm.sig;
0149                         break;
0150                     case proj_murd2:
0151                         lp_lat = this->m_proj_parm.sig - atan(rho - this->m_proj_parm.rho_c);
0152                         break;
0153                     default:
0154                         lp_lat = this->m_proj_parm.rho_c - rho;
0155                     }
0156                 }
0157 
0158                 static inline std::string get_name()
0159                 {
0160                     return "sconics_spheroid";
0161                 }
0162 
0163             };
0164 
0165             template <typename Params, typename Parameters, typename T>
0166             inline void setup(Params const& params, Parameters& par, par_sconics<T>& proj_parm, proj_type type)
0167             {
0168                 static const T half_pi = detail::half_pi<T>();
0169 
0170                 T del, cs;
0171                 int err;
0172 
0173                 proj_parm.type = type;
0174 
0175                 err = phi12(params, proj_parm, &del);
0176                 if(err)
0177                     BOOST_THROW_EXCEPTION( projection_exception(err) );
0178 
0179                 switch (proj_parm.type) {
0180                 case proj_tissot:
0181                     proj_parm.n = sin(proj_parm.sig);
0182                     cs = cos(del);
0183                     proj_parm.rho_c = proj_parm.n / cs + cs / proj_parm.n;
0184                     proj_parm.rho_0 = sqrt((proj_parm.rho_c - 2 * sin(par.phi0))/proj_parm.n);
0185                     break;
0186                 case proj_murd1:
0187                     proj_parm.rho_c = sin(del)/(del * tan(proj_parm.sig)) + proj_parm.sig;
0188                     proj_parm.rho_0 = proj_parm.rho_c - par.phi0;
0189                     proj_parm.n = sin(proj_parm.sig);
0190                     break;
0191                 case proj_murd2:
0192                     proj_parm.rho_c = (cs = sqrt(cos(del))) / tan(proj_parm.sig);
0193                     proj_parm.rho_0 = proj_parm.rho_c + tan(proj_parm.sig - par.phi0);
0194                     proj_parm.n = sin(proj_parm.sig) * cs;
0195                     break;
0196                 case proj_murd3:
0197                     proj_parm.rho_c = del / (tan(proj_parm.sig) * tan(del)) + proj_parm.sig;
0198                     proj_parm.rho_0 = proj_parm.rho_c - par.phi0;
0199                     proj_parm.n = sin(proj_parm.sig) * sin(del) * tan(del) / (del * del);
0200                     break;
0201                 case proj_euler:
0202                     proj_parm.n = sin(proj_parm.sig) * sin(del) / del;
0203                     del *= 0.5;
0204                     proj_parm.rho_c = del / (tan(del) * tan(proj_parm.sig)) + proj_parm.sig;
0205                     proj_parm.rho_0 = proj_parm.rho_c - par.phi0;
0206                     break;
0207                 case proj_pconic:
0208                     proj_parm.n = sin(proj_parm.sig);
0209                     proj_parm.c2 = cos(del);
0210                     proj_parm.c1 = 1./tan(proj_parm.sig);
0211                     if (fabs(del = par.phi0 - proj_parm.sig) - epsilon10 >= half_pi)
0212                         BOOST_THROW_EXCEPTION( projection_exception(error_lat_0_half_pi_from_mean) );
0213                     proj_parm.rho_0 = proj_parm.c2 * (proj_parm.c1 - tan(del));
0214                     break;
0215                 case proj_vitk1:
0216                     proj_parm.n = (cs = tan(del)) * sin(proj_parm.sig) / del;
0217                     proj_parm.rho_c = del / (cs * tan(proj_parm.sig)) + proj_parm.sig;
0218                     proj_parm.rho_0 = proj_parm.rho_c - par.phi0;
0219                     break;
0220                 }
0221 
0222                 par.es = 0;
0223             }
0224 
0225 
0226             // Euler
0227             template <typename Params, typename Parameters, typename T>
0228             inline void setup_euler(Params const& params, Parameters& par, par_sconics<T>& proj_parm)
0229             {
0230                 setup(params, par, proj_parm, proj_euler);
0231             }
0232 
0233             // Tissot
0234             template <typename Params, typename Parameters, typename T>
0235             inline void setup_tissot(Params const& params, Parameters& par, par_sconics<T>& proj_parm)
0236             {
0237                 setup(params, par, proj_parm, proj_tissot);
0238             }
0239 
0240             // Murdoch I
0241             template <typename Params, typename Parameters, typename T>
0242             inline void setup_murd1(Params const& params, Parameters& par, par_sconics<T>& proj_parm)
0243             {
0244                 setup(params, par, proj_parm, proj_murd1);
0245             }
0246 
0247             // Murdoch II
0248             template <typename Params, typename Parameters, typename T>
0249             inline void setup_murd2(Params const& params, Parameters& par, par_sconics<T>& proj_parm)
0250             {
0251                 setup(params, par, proj_parm, proj_murd2);
0252             }
0253 
0254             // Murdoch III
0255             template <typename Params, typename Parameters, typename T>
0256             inline void setup_murd3(Params const& params, Parameters& par, par_sconics<T>& proj_parm)
0257             {
0258                 setup(params, par, proj_parm, proj_murd3);
0259             }
0260 
0261             // Perspective Conic
0262             template <typename Params, typename Parameters, typename T>
0263             inline void setup_pconic(Params const& params, Parameters& par, par_sconics<T>& proj_parm)
0264             {
0265                 setup(params, par, proj_parm, proj_pconic);
0266             }
0267 
0268             // Vitkovsky I
0269             template <typename Params, typename Parameters, typename T>
0270             inline void setup_vitk1(Params const& params, Parameters& par, par_sconics<T>& proj_parm)
0271             {
0272                 setup(params, par, proj_parm, proj_vitk1);
0273             }
0274 
0275     }} // namespace detail::sconics
0276     #endif // doxygen
0277 
0278     /*!
0279         \brief Tissot projection
0280         \ingroup projections
0281         \tparam Geographic latlong point type
0282         \tparam Cartesian xy point type
0283         \tparam Parameters parameter type
0284         \par Projection characteristics
0285          - Conic
0286          - Spheroid
0287         \par Projection parameters
0288          - lat_1: Latitude of first standard parallel
0289          - lat_2: Latitude of second standard parallel
0290         \par Example
0291         \image html ex_tissot.gif
0292     */
0293     template <typename T, typename Parameters>
0294     struct tissot_spheroid : public detail::sconics::base_sconics_spheroid<T, Parameters>
0295     {
0296         template <typename Params>
0297         inline tissot_spheroid(Params const& params, Parameters& par)
0298         {
0299             detail::sconics::setup_tissot(params, par, this->m_proj_parm);
0300         }
0301     };
0302 
0303     /*!
0304         \brief Murdoch I projection
0305         \ingroup projections
0306         \tparam Geographic latlong point type
0307         \tparam Cartesian xy point type
0308         \tparam Parameters parameter type
0309         \par Projection characteristics
0310          - Conic
0311          - Spheroid
0312         \par Projection parameters
0313          - lat_1: Latitude of first standard parallel
0314          - lat_2: Latitude of second standard parallel
0315         \par Example
0316         \image html ex_murd1.gif
0317     */
0318     template <typename T, typename Parameters>
0319     struct murd1_spheroid : public detail::sconics::base_sconics_spheroid<T, Parameters>
0320     {
0321         template <typename Params>
0322         inline murd1_spheroid(Params const& params, Parameters& par)
0323         {
0324             detail::sconics::setup_murd1(params, par, this->m_proj_parm);
0325         }
0326     };
0327 
0328     /*!
0329         \brief Murdoch II projection
0330         \ingroup projections
0331         \tparam Geographic latlong point type
0332         \tparam Cartesian xy point type
0333         \tparam Parameters parameter type
0334         \par Projection characteristics
0335          - Conic
0336          - Spheroid
0337         \par Projection parameters
0338          - lat_1: Latitude of first standard parallel
0339          - lat_2: Latitude of second standard parallel
0340         \par Example
0341         \image html ex_murd2.gif
0342     */
0343     template <typename T, typename Parameters>
0344     struct murd2_spheroid : public detail::sconics::base_sconics_spheroid<T, Parameters>
0345     {
0346         template <typename Params>
0347         inline murd2_spheroid(Params const& params, Parameters& par)
0348         {
0349             detail::sconics::setup_murd2(params, par, this->m_proj_parm);
0350         }
0351     };
0352 
0353     /*!
0354         \brief Murdoch III projection
0355         \ingroup projections
0356         \tparam Geographic latlong point type
0357         \tparam Cartesian xy point type
0358         \tparam Parameters parameter type
0359         \par Projection characteristics
0360          - Conic
0361          - Spheroid
0362         \par Projection parameters
0363          - lat_1: Latitude of first standard parallel
0364          - lat_2: Latitude of second standard parallel
0365         \par Example
0366         \image html ex_murd3.gif
0367     */
0368     template <typename T, typename Parameters>
0369     struct murd3_spheroid : public detail::sconics::base_sconics_spheroid<T, Parameters>
0370     {
0371         template <typename Params>
0372         inline murd3_spheroid(Params const& params, Parameters& par)
0373         {
0374             detail::sconics::setup_murd3(params, par, this->m_proj_parm);
0375         }
0376     };
0377 
0378     /*!
0379         \brief Euler projection
0380         \ingroup projections
0381         \tparam Geographic latlong point type
0382         \tparam Cartesian xy point type
0383         \tparam Parameters parameter type
0384         \par Projection characteristics
0385          - Conic
0386          - Spheroid
0387         \par Projection parameters
0388          - lat_1: Latitude of first standard parallel
0389          - lat_2: Latitude of second standard parallel
0390         \par Example
0391         \image html ex_euler.gif
0392     */
0393     template <typename T, typename Parameters>
0394     struct euler_spheroid : public detail::sconics::base_sconics_spheroid<T, Parameters>
0395     {
0396         template <typename Params>
0397         inline euler_spheroid(Params const& params, Parameters& par)
0398         {
0399             detail::sconics::setup_euler(params, par, this->m_proj_parm);
0400         }
0401     };
0402 
0403     /*!
0404         \brief Perspective Conic projection
0405         \ingroup projections
0406         \tparam Geographic latlong point type
0407         \tparam Cartesian xy point type
0408         \tparam Parameters parameter type
0409         \par Projection characteristics
0410          - Conic
0411          - Spheroid
0412         \par Projection parameters
0413          - lat_1: Latitude of first standard parallel
0414          - lat_2: Latitude of second standard parallel
0415         \par Example
0416         \image html ex_pconic.gif
0417     */
0418     template <typename T, typename Parameters>
0419     struct pconic_spheroid : public detail::sconics::base_sconics_spheroid<T, Parameters>
0420     {
0421         template <typename Params>
0422         inline pconic_spheroid(Params const& params, Parameters& par)
0423         {
0424             detail::sconics::setup_pconic(params, par, this->m_proj_parm);
0425         }
0426     };
0427 
0428     /*!
0429         \brief Vitkovsky I projection
0430         \ingroup projections
0431         \tparam Geographic latlong point type
0432         \tparam Cartesian xy point type
0433         \tparam Parameters parameter type
0434         \par Projection characteristics
0435          - Conic
0436          - Spheroid
0437         \par Projection parameters
0438          - lat_1: Latitude of first standard parallel
0439          - lat_2: Latitude of second standard parallel
0440         \par Example
0441         \image html ex_vitk1.gif
0442     */
0443     template <typename T, typename Parameters>
0444     struct vitk1_spheroid : public detail::sconics::base_sconics_spheroid<T, Parameters>
0445     {
0446         template <typename Params>
0447         inline vitk1_spheroid(Params const& params, Parameters& par)
0448         {
0449             detail::sconics::setup_vitk1(params, par, this->m_proj_parm);
0450         }
0451     };
0452 
0453     #ifndef DOXYGEN_NO_DETAIL
0454     namespace detail
0455     {
0456 
0457         // Static projection
0458         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_euler, euler_spheroid)
0459         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_murd1, murd1_spheroid)
0460         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_murd2, murd2_spheroid)
0461         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_murd3, murd3_spheroid)
0462         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_pconic, pconic_spheroid)
0463         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_tissot, tissot_spheroid)
0464         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_vitk1, vitk1_spheroid)
0465 
0466         // Factory entry(s)
0467         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(euler_entry, euler_spheroid)
0468         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(murd1_entry, murd1_spheroid)
0469         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(murd2_entry, murd2_spheroid)
0470         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(murd3_entry, murd3_spheroid)
0471         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(pconic_entry, pconic_spheroid)
0472         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(tissot_entry, tissot_spheroid)
0473         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(vitk1_entry, vitk1_spheroid)
0474 
0475         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_BEGIN(sconics_init)
0476         {
0477             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(euler, euler_entry)
0478             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(murd1, murd1_entry)
0479             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(murd2, murd2_entry)
0480             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(murd3, murd3_entry)
0481             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(pconic, pconic_entry)
0482             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(tissot, tissot_entry)
0483             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(vitk1, vitk1_entry)
0484         }
0485 
0486     } // namespace detail
0487     #endif // doxygen
0488 
0489 } // namespace projections
0490 
0491 }} // namespace boost::geometry
0492 
0493 #endif // BOOST_GEOMETRY_PROJECTIONS_SCONICS_HPP
0494
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!