• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
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 // Permission is hereby granted, free of charge, to any person obtaining a
23 // copy of this software and associated documentation files (the "Software"),
24 // to deal in the Software without restriction, including without limitation
25 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
26 // and/or sell copies of the Software, and to permit persons to whom the
27 // Software is furnished to do so, subject to the following conditions:
28 
29 // The above copyright notice and this permission notice shall be included
30 // in all copies or substantial portions of the Software.
31 
32 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
33 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
34 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
35 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
36 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
37 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
38 // DEALINGS IN THE SOFTWARE.
39 
40 #ifndef BOOST_GEOMETRY_PROJECTIONS_GN_SINU_HPP
41 #define BOOST_GEOMETRY_PROJECTIONS_GN_SINU_HPP
42 
43 #include <boost/geometry/srs/projections/impl/aasincos.hpp>
44 #include <boost/geometry/srs/projections/impl/base_static.hpp>
45 #include <boost/geometry/srs/projections/impl/base_dynamic.hpp>
46 #include <boost/geometry/srs/projections/impl/factory_entry.hpp>
47 #include <boost/geometry/srs/projections/impl/pj_mlfn.hpp>
48 #include <boost/geometry/srs/projections/impl/pj_param.hpp>
49 #include <boost/geometry/srs/projections/impl/projects.hpp>
50 
51 #include <boost/geometry/util/math.hpp>
52 
53 namespace boost { namespace geometry
54 {
55 
56 namespace projections
57 {
58     #ifndef DOXYGEN_NO_DETAIL
59     namespace detail { namespace gn_sinu
60     {
61 
62             static const double epsilon10 = 1e-10;
63             static const int max_iter = 8;
64             static const double loop_tol = 1e-7;
65 
66             template <typename T>
67             struct par_gn_sinu_e
68             {
69                 detail::en<T> en;
70             };
71 
72             template <typename T>
73             struct par_gn_sinu_s
74             {
75                 T m, n, C_x, C_y;
76             };
77 
78             /* Ellipsoidal Sinusoidal only */
79 
80             template <typename T, typename Parameters>
81             struct base_gn_sinu_ellipsoid
82             {
83                 par_gn_sinu_e<T> m_proj_parm;
84 
85                 // FORWARD(e_forward)  ellipsoid
86                 // Project coordinates from geographic (lon, lat) to cartesian (x, y)
fwdboost::geometry::projections::detail::gn_sinu::base_gn_sinu_ellipsoid87                 inline void fwd(Parameters const& par, T const& lp_lon, T const& lp_lat, T& xy_x, T& xy_y) const
88                 {
89                     T s, c;
90 
91                     xy_y = pj_mlfn(lp_lat, s = sin(lp_lat), c = cos(lp_lat), this->m_proj_parm.en);
92                     xy_x = lp_lon * c / sqrt(1. - par.es * s * s);
93                 }
94 
95                 // INVERSE(e_inverse)  ellipsoid
96                 // Project coordinates from cartesian (x, y) to geographic (lon, lat)
invboost::geometry::projections::detail::gn_sinu::base_gn_sinu_ellipsoid97                 inline void inv(Parameters const& par, T const& xy_x, T const& xy_y, T& lp_lon, T& lp_lat) const
98                 {
99                     static const T half_pi = detail::half_pi<T>();
100 
101                     T s;
102 
103                     if ((s = fabs(lp_lat = pj_inv_mlfn(xy_y, par.es, this->m_proj_parm.en))) < half_pi) {
104                         s = sin(lp_lat);
105                         lp_lon = xy_x * sqrt(1. - par.es * s * s) / cos(lp_lat);
106                     } else if ((s - epsilon10) < half_pi)
107                         lp_lon = 0.;
108                     else
109                         BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
110                 }
111                 /* General spherical sinusoidals */
112 
get_nameboost::geometry::projections::detail::gn_sinu::base_gn_sinu_ellipsoid113                 static inline std::string get_name()
114                 {
115                     return "gn_sinu_ellipsoid";
116                 }
117 
118             };
119 
120             template <typename T, typename Parameters>
121             struct base_gn_sinu_spheroid
122             {
123                 par_gn_sinu_s<T> m_proj_parm;
124 
125                 // FORWARD(s_forward)  sphere
126                 // Project coordinates from geographic (lon, lat) to cartesian (x, y)
fwdboost::geometry::projections::detail::gn_sinu::base_gn_sinu_spheroid127                 inline void fwd(Parameters const& , T const& lp_lon, T lp_lat, T& xy_x, T& xy_y) const
128                 {
129                     if (this->m_proj_parm.m == 0.0)
130                         lp_lat = this->m_proj_parm.n != 1. ? aasin(this->m_proj_parm.n * sin(lp_lat)): lp_lat;
131                     else {
132                         T k, V;
133                         int i;
134 
135                         k = this->m_proj_parm.n * sin(lp_lat);
136                         for (i = max_iter; i ; --i) {
137                             lp_lat -= V = (this->m_proj_parm.m * lp_lat + sin(lp_lat) - k) /
138                                 (this->m_proj_parm.m + cos(lp_lat));
139                             if (fabs(V) < loop_tol)
140                                 break;
141                         }
142                         if (!i) {
143                             BOOST_THROW_EXCEPTION( projection_exception(error_tolerance_condition) );
144                         }
145                     }
146                     xy_x = this->m_proj_parm.C_x * lp_lon * (this->m_proj_parm.m + cos(lp_lat));
147                     xy_y = this->m_proj_parm.C_y * lp_lat;
148                 }
149 
150                 // INVERSE(s_inverse)  sphere
151                 // Project coordinates from cartesian (x, y) to geographic (lon, lat)
invboost::geometry::projections::detail::gn_sinu::base_gn_sinu_spheroid152                 inline void inv(Parameters const& , T const& xy_x, T xy_y, T& lp_lon, T& lp_lat) const
153                 {
154                     xy_y /= this->m_proj_parm.C_y;
155                     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) :
156                         ( this->m_proj_parm.n != 1. ? aasin(sin(xy_y) / this->m_proj_parm.n) : xy_y );
157                     lp_lon = xy_x / (this->m_proj_parm.C_x * (this->m_proj_parm.m + cos(xy_y)));
158                 }
159 
get_nameboost::geometry::projections::detail::gn_sinu::base_gn_sinu_spheroid160                 static inline std::string get_name()
161                 {
162                     return "gn_sinu_spheroid";
163                 }
164 
165             };
166 
167             template <typename Parameters, typename T>
setup(Parameters & par,par_gn_sinu_s<T> & proj_parm)168             inline void setup(Parameters& par, par_gn_sinu_s<T>& proj_parm)
169             {
170                 par.es = 0;
171 
172                 proj_parm.C_x = (proj_parm.C_y = sqrt((proj_parm.m + 1.) / proj_parm.n))/(proj_parm.m + 1.);
173             }
174 
175 
176             // General Sinusoidal Series
177             template <typename Params, typename Parameters, typename T>
setup_gn_sinu(Params const & params,Parameters & par,par_gn_sinu_s<T> & proj_parm)178             inline void setup_gn_sinu(Params const& params, Parameters& par, par_gn_sinu_s<T>& proj_parm)
179             {
180                 if (pj_param_f<srs::spar::n>(params, "n", srs::dpar::n, proj_parm.n)
181                  && pj_param_f<srs::spar::m>(params, "m", srs::dpar::m, proj_parm.m)) {
182                     if (proj_parm.n <= 0 || proj_parm.m < 0)
183                         BOOST_THROW_EXCEPTION( projection_exception(error_invalid_m_or_n) );
184                 } else
185                     BOOST_THROW_EXCEPTION( projection_exception(error_invalid_m_or_n) );
186 
187                 setup(par, proj_parm);
188             }
189 
190             // Sinusoidal (Sanson-Flamsteed)
191             template <typename Parameters, typename T>
setup_sinu(Parameters const & par,par_gn_sinu_e<T> & proj_parm)192             inline void setup_sinu(Parameters const& par, par_gn_sinu_e<T>& proj_parm)
193             {
194                 proj_parm.en = pj_enfn<T>(par.es);
195             }
196 
197             // Sinusoidal (Sanson-Flamsteed)
198             template <typename Parameters, typename T>
setup_sinu(Parameters & par,par_gn_sinu_s<T> & proj_parm)199             inline void setup_sinu(Parameters& par, par_gn_sinu_s<T>& proj_parm)
200             {
201                 proj_parm.n = 1.;
202                 proj_parm.m = 0.;
203                 setup(par, proj_parm);
204             }
205 
206             // Eckert VI
207             template <typename Parameters, typename T>
setup_eck6(Parameters & par,par_gn_sinu_s<T> & proj_parm)208             inline void setup_eck6(Parameters& par, par_gn_sinu_s<T>& proj_parm)
209             {
210                 proj_parm.m = 1.;
211                 proj_parm.n = 2.570796326794896619231321691;
212                 setup(par, proj_parm);
213             }
214 
215             // McBryde-Thomas Flat-Polar Sinusoidal
216             template <typename Parameters, typename T>
setup_mbtfps(Parameters & par,par_gn_sinu_s<T> & proj_parm)217             inline void setup_mbtfps(Parameters& par, par_gn_sinu_s<T>& proj_parm)
218             {
219                 proj_parm.m = 0.5;
220                 proj_parm.n = 1.785398163397448309615660845;
221                 setup(par, proj_parm);
222             }
223 
224     }} // namespace detail::gn_sinu
225     #endif // doxygen
226 
227     /*!
228         \brief General Sinusoidal Series projection
229         \ingroup projections
230         \tparam Geographic latlong point type
231         \tparam Cartesian xy point type
232         \tparam Parameters parameter type
233         \par Projection characteristics
234          - Pseudocylindrical
235          - Spheroid
236         \par Projection parameters
237          - m (real)
238          - n (real)
239         \par Example
240         \image html ex_gn_sinu.gif
241     */
242     template <typename T, typename Parameters>
243     struct gn_sinu_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<T, Parameters>
244     {
245         template <typename Params>
gn_sinu_spheroidboost::geometry::projections::gn_sinu_spheroid246         inline gn_sinu_spheroid(Params const& params, Parameters & par)
247         {
248             detail::gn_sinu::setup_gn_sinu(params, par, this->m_proj_parm);
249         }
250     };
251 
252     /*!
253         \brief Sinusoidal (Sanson-Flamsteed) projection
254         \ingroup projections
255         \tparam Geographic latlong point type
256         \tparam Cartesian xy point type
257         \tparam Parameters parameter type
258         \par Projection characteristics
259          - Pseudocylindrical
260          - Spheroid
261          - Ellipsoid
262         \par Example
263         \image html ex_sinu.gif
264     */
265     template <typename T, typename Parameters>
266     struct sinu_ellipsoid : public detail::gn_sinu::base_gn_sinu_ellipsoid<T, Parameters>
267     {
268         template <typename Params>
sinu_ellipsoidboost::geometry::projections::sinu_ellipsoid269         inline sinu_ellipsoid(Params const& , Parameters & par)
270         {
271             detail::gn_sinu::setup_sinu(par, this->m_proj_parm);
272         }
273     };
274 
275     /*!
276         \brief Sinusoidal (Sanson-Flamsteed) projection
277         \ingroup projections
278         \tparam Geographic latlong point type
279         \tparam Cartesian xy point type
280         \tparam Parameters parameter type
281         \par Projection characteristics
282          - Pseudocylindrical
283          - Spheroid
284          - Ellipsoid
285         \par Example
286         \image html ex_sinu.gif
287     */
288     template <typename T, typename Parameters>
289     struct sinu_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<T, Parameters>
290     {
291         template <typename Params>
sinu_spheroidboost::geometry::projections::sinu_spheroid292         inline sinu_spheroid(Params const& , Parameters & par)
293         {
294             detail::gn_sinu::setup_sinu(par, this->m_proj_parm);
295         }
296     };
297 
298     /*!
299         \brief Eckert VI projection
300         \ingroup projections
301         \tparam Geographic latlong point type
302         \tparam Cartesian xy point type
303         \tparam Parameters parameter type
304         \par Projection characteristics
305          - Pseudocylindrical
306          - Spheroid
307         \par Example
308         \image html ex_eck6.gif
309     */
310     template <typename T, typename Parameters>
311     struct eck6_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<T, Parameters>
312     {
313         template <typename Params>
eck6_spheroidboost::geometry::projections::eck6_spheroid314         inline eck6_spheroid(Params const& , Parameters & par)
315         {
316             detail::gn_sinu::setup_eck6(par, this->m_proj_parm);
317         }
318     };
319 
320     /*!
321         \brief McBryde-Thomas Flat-Polar Sinusoidal projection
322         \ingroup projections
323         \tparam Geographic latlong point type
324         \tparam Cartesian xy point type
325         \tparam Parameters parameter type
326         \par Projection characteristics
327          - Pseudocylindrical
328          - Spheroid
329         \par Example
330         \image html ex_mbtfps.gif
331     */
332     template <typename T, typename Parameters>
333     struct mbtfps_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<T, Parameters>
334     {
335         template <typename Params>
mbtfps_spheroidboost::geometry::projections::mbtfps_spheroid336         inline mbtfps_spheroid(Params const& , Parameters & par)
337         {
338             detail::gn_sinu::setup_mbtfps(par, this->m_proj_parm);
339         }
340     };
341 
342     #ifndef DOXYGEN_NO_DETAIL
343     namespace detail
344     {
345 
346         // Static projection
BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_gn_sinu,gn_sinu_spheroid)347         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_gn_sinu, gn_sinu_spheroid)
348         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI2(srs::spar::proj_sinu, sinu_spheroid, sinu_ellipsoid)
349         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_eck6, eck6_spheroid)
350         BOOST_GEOMETRY_PROJECTIONS_DETAIL_STATIC_PROJECTION_FI(srs::spar::proj_mbtfps, mbtfps_spheroid)
351 
352         // Factory entry(s)
353         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(gn_sinu_entry, gn_sinu_spheroid)
354         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI2(sinu_entry, sinu_spheroid, sinu_ellipsoid)
355         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(eck6_entry, eck6_spheroid)
356         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_ENTRY_FI(mbtfps_entry, mbtfps_spheroid)
357 
358         BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_BEGIN(gn_sinu_init)
359         {
360             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(gn_sinu, gn_sinu_entry);
361             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(sinu, sinu_entry);
362             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(eck6, eck6_entry);
363             BOOST_GEOMETRY_PROJECTIONS_DETAIL_FACTORY_INIT_ENTRY(mbtfps, mbtfps_entry);
364         }
365 
366     } // namespace detail
367     #endif // doxygen
368 
369 } // namespace projections
370 
371 }} // namespace boost::geometry
372 
373 #endif // BOOST_GEOMETRY_PROJECTIONS_GN_SINU_HPP
374 
375