1 // Boost.Geometry - gis-projections (based on PROJ4) 2 3 // Copyright (c) 2008-2015 Barend Gehrels, Amsterdam, the Netherlands. 4 5 // This file was modified by Oracle on 2017, 2018, 2019. 6 // Modifications copyright (c) 2017-2019, Oracle and/or its affiliates. 7 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle. 8 9 // Use, modification and distribution is subject to the Boost Software License, 10 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at 11 // http://www.boost.org/LICENSE_1_0.txt) 12 13 // This file is converted from PROJ4, http://trac.osgeo.org/proj 14 // PROJ4 is originally written by Gerald Evenden (then of the USGS) 15 // PROJ4 is maintained by Frank Warmerdam 16 // PROJ4 is converted to Boost.Geometry by Barend Gehrels 17 18 // Last updated version of proj: 5.0.0 19 20 // Original copyright notice: 21 22 // Copyright (c) 2003, 2006 Gerald I. Evenden 23 24 // Permission is hereby granted, free of charge, to any person obtaining a 25 // copy of this software and associated documentation files (the "Software"), 26 // to deal in the Software without restriction, including without limitation 27 // the rights to use, copy, modify, merge, publish, distribute, sublicense, 28 // and/or sell copies of the Software, and to permit persons to whom the 29 // Software is furnished to do so, subject to the following conditions: 30 31 // The above copyright notice and this permission notice shall be included 32 // in all copies or substantial portions of the Software. 33 34 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS 35 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 36 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL 37 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER 38 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING 39 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER 40 // DEALINGS IN THE SOFTWARE. 41 42 #ifndef BOOST_GEOMETRY_PROJECTIONS_ROUSS_HPP 43 #define BOOST_GEOMETRY_PROJECTIONS_ROUSS_HPP 44 45 #include <boost/geometry/srs/projections/impl/base_static.hpp> 46 #include <boost/geometry/srs/projections/impl/base_dynamic.hpp> 47 #include <boost/geometry/srs/projections/impl/projects.hpp> 48 #include <boost/geometry/srs/projections/impl/factory_entry.hpp> 49 #include <boost/geometry/srs/projections/impl/proj_mdist.hpp> 50 51 namespace boost { namespace geometry 52 { 53 54 namespace projections 55 { 56 #ifndef DOXYGEN_NO_DETAIL 57 namespace detail { namespace rouss 58 { 59 template <typename T> 60 struct par_rouss 61 { 62 T s0; 63 T A1, A2, A3, A4, A5, A6; 64 T B1, B2, B3, B4, B5, B6, B7, B8; 65 T C1, C2, C3, C4, C5, C6, C7, C8; 66 T D1, D2, D3, D4, D5, D6, D7, D8, D9, D10, D11; 67 mdist<T> en; 68 }; 69 70 template <typename T, typename Parameters> 71 struct base_rouss_ellipsoid 72 { 73 par_rouss<T> m_proj_parm; 74 75 // FORWARD(e_forward) ellipsoid 76 // Project coordinates from geographic (lon, lat) to cartesian (x, y) fwdboost::geometry::projections::detail::rouss::base_rouss_ellipsoid77 inline void fwd(Parameters const& par, T const& lp_lon, T const& lp_lat, T& xy_x, T& xy_y) const 78 { 79 T s, al, cp, sp, al2, s2; 80 81 cp = cos(lp_lat); 82 sp = sin(lp_lat); 83 s = proj_mdist(lp_lat, sp, cp, this->m_proj_parm.en) - this->m_proj_parm.s0; 84 s2 = s * s; 85 al = lp_lon * cp / sqrt(1. - par.es * sp * sp); 86 al2 = al * al; 87 xy_x = par.k0 * al*(1.+s2*(this->m_proj_parm.A1+s2*this->m_proj_parm.A4)-al2*(this->m_proj_parm.A2+s*this->m_proj_parm.A3+s2*this->m_proj_parm.A5 88 +al2*this->m_proj_parm.A6)); 89 xy_y = par.k0 * (al2*(this->m_proj_parm.B1+al2*this->m_proj_parm.B4)+ 90 s*(1.+al2*(this->m_proj_parm.B3-al2*this->m_proj_parm.B6)+s2*(this->m_proj_parm.B2+s2*this->m_proj_parm.B8)+ 91 s*al2*(this->m_proj_parm.B5+s*this->m_proj_parm.B7))); 92 } 93 94 // INVERSE(e_inverse) ellipsoid 95 // Project coordinates from cartesian (x, y) to geographic (lon, lat) invboost::geometry::projections::detail::rouss::base_rouss_ellipsoid96 inline void inv(Parameters const& par, T const& xy_x, T const& xy_y, T& lp_lon, T& lp_lat) const 97 { 98 T s, al, x = xy_x / par.k0, y = xy_y / par.k0, x2, y2; 99 100 x2 = x * x; 101 y2 = y * y; 102 al = x*(1.-this->m_proj_parm.C1*y2+x2*(this->m_proj_parm.C2+this->m_proj_parm.C3*y-this->m_proj_parm.C4*x2+this->m_proj_parm.C5*y2-this->m_proj_parm.C7*x2*y) 103 +y2*(this->m_proj_parm.C6*y2-this->m_proj_parm.C8*x2*y)); 104 s = this->m_proj_parm.s0 + y*(1.+y2*(-this->m_proj_parm.D2+this->m_proj_parm.D8*y2))+ 105 x2*(-this->m_proj_parm.D1+y*(-this->m_proj_parm.D3+y*(-this->m_proj_parm.D5+y*(-this->m_proj_parm.D7+y*this->m_proj_parm.D11)))+ 106 x2*(this->m_proj_parm.D4+y*(this->m_proj_parm.D6+y*this->m_proj_parm.D10)-x2*this->m_proj_parm.D9)); 107 lp_lat=proj_inv_mdist(s, this->m_proj_parm.en); 108 s = sin(lp_lat); 109 lp_lon=al * sqrt(1. - par.es * s * s)/cos(lp_lat); 110 } 111 get_nameboost::geometry::projections::detail::rouss::base_rouss_ellipsoid112 static inline std::string get_name() 113 { 114 return "rouss_ellipsoid"; 115 } 116 117 }; 118 119 // Roussilhe Stereographic 120 template <typename Parameters, typename T> setup_rouss(Parameters const & par,par_rouss<T> & proj_parm)121 inline void setup_rouss(Parameters const& par, par_rouss<T>& proj_parm) 122 { 123 T N0, es2, t, t2, R_R0_2, R_R0_4; 124 125 if (!proj_mdist_ini(par.es, proj_parm.en)) 126 BOOST_THROW_EXCEPTION( projection_exception(0) ); 127 128 es2 = sin(par.phi0); 129 proj_parm.s0 = proj_mdist(par.phi0, es2, cos(par.phi0), proj_parm.en); 130 t = 1. - (es2 = par.es * es2 * es2); 131 N0 = 1./sqrt(t); 132 R_R0_2 = t * t / par.one_es; 133 R_R0_4 = R_R0_2 * R_R0_2; 134 t = tan(par.phi0); 135 t2 = t * t; 136 proj_parm.C1 = proj_parm.A1 = R_R0_2 / 4.; 137 proj_parm.C2 = proj_parm.A2 = R_R0_2 * (2 * t2 - 1. - 2. * es2) / 12.; 138 proj_parm.A3 = R_R0_2 * t * (1. + 4. * t2)/ ( 12. * N0); 139 proj_parm.A4 = R_R0_4 / 24.; 140 proj_parm.A5 = R_R0_4 * ( -1. + t2 * (11. + 12. * t2))/24.; 141 proj_parm.A6 = R_R0_4 * ( -2. + t2 * (11. - 2. * t2))/240.; 142 proj_parm.B1 = t / (2. * N0); 143 proj_parm.B2 = R_R0_2 / 12.; 144 proj_parm.B3 = R_R0_2 * (1. + 2. * t2 - 2. * es2)/4.; 145 proj_parm.B4 = R_R0_2 * t * (2. - t2)/(24. * N0); 146 proj_parm.B5 = R_R0_2 * t * (5. + 4.* t2)/(8. * N0); 147 proj_parm.B6 = R_R0_4 * (-2. + t2 * (-5. + 6. * t2))/48.; 148 proj_parm.B7 = R_R0_4 * (5. + t2 * (19. + 12. * t2))/24.; 149 proj_parm.B8 = R_R0_4 / 120.; 150 proj_parm.C3 = R_R0_2 * t * (1. + t2)/(3. * N0); 151 proj_parm.C4 = R_R0_4 * (-3. + t2 * (34. + 22. * t2))/240.; 152 proj_parm.C5 = R_R0_4 * (4. + t2 * (13. + 12. * t2))/24.; 153 proj_parm.C6 = R_R0_4 / 16.; 154 proj_parm.C7 = R_R0_4 * t * (11. + t2 * (33. + t2 * 16.))/(48. * N0); 155 proj_parm.C8 = R_R0_4 * t * (1. + t2 * 4.)/(36. * N0); 156 proj_parm.D1 = t / (2. * N0); 157 proj_parm.D2 = R_R0_2 / 12.; 158 proj_parm.D3 = R_R0_2 * (2 * t2 + 1. - 2. * es2) / 4.; 159 proj_parm.D4 = R_R0_2 * t * (1. + t2)/(8. * N0); 160 proj_parm.D5 = R_R0_2 * t * (1. + t2 * 2.)/(4. * N0); 161 proj_parm.D6 = R_R0_4 * (1. + t2 * (6. + t2 * 6.))/16.; 162 proj_parm.D7 = R_R0_4 * t2 * (3. + t2 * 4.)/8.; 163 proj_parm.D8 = R_R0_4 / 80.; 164 proj_parm.D9 = R_R0_4 * t * (-21. + t2 * (178. - t2 * 26.))/720.; 165 proj_parm.D10 = R_R0_4 * t * (29. + t2 * (86. + t2 * 48.))/(96. * N0); 166 proj_parm.D11 = R_R0_4 * t * (37. + t2 * 44.)/(96. * N0); 167 } 168 169 }} // namespace detail::rouss 170 #endif // doxygen 171 172 /*! 173 \brief Roussilhe Stereographic projection 174 \ingroup projections 175 \tparam Geographic latlong point type 176 \tparam Cartesian xy point type 177 \tparam Parameters parameter type 178 \par Projection characteristics 179 - Azimuthal 180 - Ellipsoid 181 \par Example 182 \image html ex_rouss.gif 183 */ 184 template <typename T, typename Parameters> 185 struct rouss_ellipsoid : public detail::rouss::base_rouss_ellipsoid<T, Parameters> 186 { 187 template <typename Params> rouss_ellipsoidboost::geometry::projections::rouss_ellipsoid188 inline rouss_ellipsoid(Params const& , Parameters const& par) 189 { 190 detail::rouss::setup_rouss(par, this->m_proj_parm); 191 } 192 }; 193 194 #ifndef DOXYGEN_NO_DETAIL 195 namespace detail 196 { 197 198 // Static projection BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_rouss,rouss_ellipsoid)199 BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_rouss, rouss_ellipsoid) 200 201 // Factory entry(s) 202 BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(rouss_entry, rouss_ellipsoid) 203 204 BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_BEGIN(rouss_init) 205 { 206 BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(rouss, rouss_entry) 207 } 208 209 } // namespace detail 210 #endif // doxygen 211 212 } // namespace projections 213 214 }} // namespace boost::geometry 215 216 #endif // BOOST_GEOMETRY_PROJECTIONS_ROUSS_HPP 217 218