1 // Boost.Geometry (aka GGL, Generic Geometry Library)
2 // This file is manually converted from PROJ4
3 
4 // Copyright (c) 2008-2012 Barend Gehrels, Amsterdam, the Netherlands.
5 
6 // This file was modified by Oracle on 2017, 2018, 2019.
7 // Modifications copyright (c) 2017-2019, Oracle and/or its affiliates.
8 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
9 
10 // Use, modification and distribution is subject to the Boost Software License,
11 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
12 // http://www.boost.org/LICENSE_1_0.txt)
13 
14 // This file is converted from PROJ4, http://trac.osgeo.org/proj
15 // PROJ4 is originally written by Gerald Evenden (then of the USGS)
16 // PROJ4 is maintained by Frank Warmerdam
17 // PROJ4 is converted to Geometry Library by Barend Gehrels (Geodan, Amsterdam)
18 
19 // Original copyright notice:
20 
21 // Permission is hereby granted, free of charge, to any person obtaining a
22 // copy of this software and associated documentation files (the "Software"),
23 // to deal in the Software without restriction, including without limitation
24 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
25 // and/or sell copies of the Software, and to permit persons to whom the
26 // Software is furnished to do so, subject to the following conditions:
27 
28 // The above copyright notice and this permission notice shall be included
29 // in all copies or substantial portions of the Software.
30 
31 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
32 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
33 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
34 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
35 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
36 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
37 // DEALINGS IN THE SOFTWARE.
38 
39 #ifndef BOOST_GEOMETRY_PROJECTIONS_IMPL_PJ_DATUM_SET_HPP
40 #define BOOST_GEOMETRY_PROJECTIONS_IMPL_PJ_DATUM_SET_HPP
41 
42 
43 #include <string>
44 #include <vector>
45 
46 #include <boost/algorithm/string.hpp>
47 
48 #include <boost/geometry/srs/projections/dpar.hpp>
49 #include <boost/geometry/srs/projections/exception.hpp>
50 #include <boost/geometry/srs/projections/impl/projects.hpp>
51 #include <boost/geometry/srs/projections/impl/pj_datums.hpp>
52 #include <boost/geometry/srs/projections/impl/pj_param.hpp>
53 #include <boost/geometry/srs/projections/proj4.hpp>
54 #include <boost/geometry/srs/projections/spar.hpp>
55 
56 
57 namespace boost { namespace geometry { namespace projections {
58 
59 namespace detail {
60 
61 
62 /* SEC_TO_RAD = Pi/180/3600 */
63 template <typename T>
sec_to_rad()64 inline T sec_to_rad() { return 4.84813681109535993589914102357e-6; }
65 
66 /************************************************************************/
67 /*                        pj_datum_find_datum()                         */
68 /************************************************************************/
69 
70 template <typename T>
pj_datum_find_datum(srs::detail::proj4_parameters const & params)71 inline const pj_datums_type<T>* pj_datum_find_datum(srs::detail::proj4_parameters const& params)
72 {
73     std::string name = pj_get_param_s(params, "datum");
74     if(! name.empty())
75     {
76         /* find the datum definition */
77         const pj_datums_type<T>* pj_datums = pj_get_datums<T>().first;
78         const int n = pj_get_datums<T>().second;
79         int index = -1;
80         for (int i = 0; i < n && index == -1; i++)
81         {
82             if(pj_datums[i].id == name)
83             {
84                 index = i;
85             }
86         }
87 
88         if (index != -1)
89         {
90             return pj_datums + index;
91         }
92         else
93         {
94             BOOST_THROW_EXCEPTION( projection_exception(error_unknown_ellp_param) );
95         }
96     }
97 
98     return NULL;
99 }
100 
101 template <typename T>
pj_datum_find_datum(srs::dpar::parameters<T> const & params)102 inline const pj_datums_type<T>* pj_datum_find_datum(srs::dpar::parameters<T> const& params)
103 {
104     typename srs::dpar::parameters<T>::const_iterator
105         it = pj_param_find(params, srs::dpar::datum);
106 
107     if (it != params.end())
108     {
109         const pj_datums_type<T>* pj_datums = pj_get_datums<T>().first;
110         const int n = pj_get_datums<T>().second;
111         int i = it->template get_value<int>();
112         if (i >= 0 && i < n)
113         {
114             return pj_datums + i;
115         }
116         else
117         {
118             BOOST_THROW_EXCEPTION( projection_exception(error_unknown_ellp_param) );
119         }
120     }
121 
122     return NULL;
123 }
124 
125 template
126 <
127     typename Params,
128     typename Param = typename srs::spar::detail::tuples_find_if
129         <
130             Params,
131             srs::spar::detail::is_param_tr<srs::spar::detail::datum_traits>::pred
132         >::type,
133     bool IsFound = srs::spar::detail::tuples_is_found<Param>::value
134 >
135 struct pj_datum_find_datum_static
136 {
137     template <typename T>
applyboost::geometry::projections::detail::pj_datum_find_datum_static138     static const pj_datums_type<T>* apply(Params const& )
139     {
140         const pj_datums_type<T>* pj_datums = pj_get_datums<T>().first;
141         const int n = pj_get_datums<T>().second;
142         const int i = srs::spar::detail::datum_traits<Param>::id;
143         if (i >= 0 && i < n)
144         {
145             return pj_datums + i;
146         }
147         else
148         {
149             // TODO: Implemnt as MPL_ASSERT instead
150             BOOST_THROW_EXCEPTION( projection_exception(error_unknown_ellp_param) );
151         }
152     }
153 };
154 template <typename Params, typename Param>
155 struct pj_datum_find_datum_static<Params, Param, false>
156 {
157     template <typename T>
applyboost::geometry::projections::detail::pj_datum_find_datum_static158     static const pj_datums_type<T>* apply(Params const& )
159     {
160         return NULL;
161     }
162 };
163 
164 template <typename T, BOOST_GEOMETRY_PROJECTIONS_DETAIL_TYPENAME_PX>
pj_datum_find_datum(srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX> const & params)165 inline const pj_datums_type<T>* pj_datum_find_datum(srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX> const& params)
166 {
167     return pj_datum_find_datum_static
168         <
169             srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>
170         >::template apply<T>(params);
171 }
172 
173 /************************************************************************/
174 /*                        pj_datum_find_nadgrids()                      */
175 /************************************************************************/
176 
pj_datum_find_nadgrids(srs::detail::proj4_parameters const & params,srs::detail::nadgrids & out)177 inline bool pj_datum_find_nadgrids(srs::detail::proj4_parameters const& params,
178                                    srs::detail::nadgrids & out)
179 {
180     std::string snadgrids = pj_get_param_s(params, "nadgrids");
181     if (! snadgrids.empty())
182     {
183         for (std::string::size_type i = 0 ; i < snadgrids.size() ; )
184         {
185             std::string::size_type end = snadgrids.find(',', i);
186             std::string name = snadgrids.substr(i, end - i);
187 
188             i = end;
189             if (end != std::string::npos)
190                 ++i;
191 
192             if (! name.empty())
193                 out.push_back(name);
194         }
195     }
196 
197     return ! out.empty();
198 }
199 
200 template <typename T>
pj_datum_find_nadgrids(srs::dpar::parameters<T> const & params,srs::detail::nadgrids & out)201 inline bool pj_datum_find_nadgrids(srs::dpar::parameters<T> const& params,
202                                    srs::detail::nadgrids & out)
203 {
204     typename srs::dpar::parameters<T>::const_iterator
205         it = pj_param_find(params, srs::dpar::nadgrids);
206     if (it != params.end())
207     {
208         out = it->template get_value<srs::detail::nadgrids>();
209     }
210 
211     return ! out.empty();
212 }
213 
214 template
215 <
216     typename Params,
217     int I = srs::spar::detail::tuples_find_index_if
218         <
219             Params,
220             srs::spar::detail::is_param<srs::spar::nadgrids>::pred
221         >::value,
222     int N = boost::tuples::length<Params>::value
223 >
224 struct pj_datum_find_nadgrids_static
225 {
applyboost::geometry::projections::detail::pj_datum_find_nadgrids_static226     static void apply(Params const& params, srs::detail::nadgrids & out)
227     {
228         out = boost::tuples::get<I>(params);
229     }
230 };
231 template <typename Params, int N>
232 struct pj_datum_find_nadgrids_static<Params, N, N>
233 {
applyboost::geometry::projections::detail::pj_datum_find_nadgrids_static234     static void apply(Params const& , srs::detail::nadgrids & )
235     {}
236 };
237 
238 template <BOOST_GEOMETRY_PROJECTIONS_DETAIL_TYPENAME_PX>
pj_datum_find_nadgrids(srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX> const & params,srs::detail::nadgrids & out)239 inline bool pj_datum_find_nadgrids(srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX> const& params,
240                                    srs::detail::nadgrids & out)
241 {
242     pj_datum_find_nadgrids_static
243         <
244             srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>
245         >::apply(params, out);
246 
247     return ! out.empty();
248 }
249 
250 /************************************************************************/
251 /*                        pj_datum_find_towgs84()                       */
252 /************************************************************************/
253 
254 template <typename T>
pj_datum_find_towgs84(srs::detail::proj4_parameters const & params,srs::detail::towgs84<T> & out)255 inline bool pj_datum_find_towgs84(srs::detail::proj4_parameters const& params,
256                                   srs::detail::towgs84<T> & out)
257 {
258     std::string towgs84 = pj_get_param_s(params, "towgs84");
259     if(! towgs84.empty())
260     {
261         std::vector<std::string> parm;
262         boost::split(parm, towgs84, boost::is_any_of(" ,"));
263 
264         std::size_t n = (std::min<std::size_t>)(parm.size(), 7);
265         std::size_t z = n <= 3 ? 3 : 7;
266 
267         /* parse out the pvalues */
268         for (std::size_t i = 0 ; i < n; ++i)
269         {
270             out.push_back(geometry::str_cast<T>(parm[i]));
271         }
272         for (std::size_t i = out.size() ; i < z; ++i)
273         {
274             out.push_back(T(0));
275         }
276     }
277 
278     return ! out.empty();
279 }
280 
281 template <typename T>
pj_datum_find_towgs84(srs::dpar::parameters<T> const & params,srs::detail::towgs84<T> & out)282 inline bool pj_datum_find_towgs84(srs::dpar::parameters<T> const& params,
283                                   srs::detail::towgs84<T> & out)
284 {
285     typename srs::dpar::parameters<T>::const_iterator
286         it = pj_param_find(params, srs::dpar::towgs84);
287 
288     if (it != params.end())
289     {
290         srs::detail::towgs84<T> const&
291             towgs84 = it->template get_value<srs::detail::towgs84<T> >();
292 
293         std::size_t n = (std::min<std::size_t>)(towgs84.size(), 7u);
294         std::size_t z = n <= 3 ? 3 : 7;
295 
296         for (std::size_t i = 0 ; i < n; ++i)
297         {
298             out.push_back(towgs84[i]);
299         }
300         for (std::size_t i = out.size() ; i < z; ++i)
301         {
302             out.push_back(T(0));
303         }
304     }
305 
306     return ! out.empty();
307 }
308 
309 template
310 <
311     typename Params,
312     int I = srs::spar::detail::tuples_find_index_if
313         <
314             Params,
315             srs::spar::detail::is_param_t<srs::spar::towgs84>::pred
316         >::value,
317     int N = boost::tuples::length<Params>::value
318 >
319 struct pj_datum_find_towgs84_static
320 {
321     template <typename T>
applyboost::geometry::projections::detail::pj_datum_find_towgs84_static322     static void apply(Params const& params, srs::detail::towgs84<T> & out)
323     {
324         typename boost::tuples::element<I, Params>::type const&
325             towgs84 = boost::tuples::get<I>(params);
326 
327         std::size_t n = (std::min<std::size_t>)(towgs84.size(), 7u);
328         std::size_t z = n <= 3 ? 3 : 7;
329 
330         for (std::size_t i = 0 ; i < n; ++i)
331         {
332             out.push_back(towgs84[i]);
333         }
334         for (std::size_t i = out.size() ; i < z; ++i)
335         {
336             out.push_back(T(0));
337         }
338     }
339 };
340 template <typename Params, int N>
341 struct pj_datum_find_towgs84_static<Params, N, N>
342 {
343     template <typename T>
applyboost::geometry::projections::detail::pj_datum_find_towgs84_static344     static void apply(Params const& , srs::detail::towgs84<T> & )
345     {}
346 };
347 
348 template <typename T, BOOST_GEOMETRY_PROJECTIONS_DETAIL_TYPENAME_PX>
pj_datum_find_towgs84(srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX> const & params,srs::detail::towgs84<T> & out)349 inline bool pj_datum_find_towgs84(srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX> const& params,
350                                   srs::detail::towgs84<T> & out)
351 {
352     pj_datum_find_towgs84_static
353         <
354             srs::spar::parameters<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>
355         >::apply(params, out);
356 
357     return ! out.empty();
358 }
359 
360 /************************************************************************/
361 /*                        pj_datum_prepare_towgs84()                    */
362 /************************************************************************/
363 
364 template <typename T>
pj_datum_prepare_towgs84(srs::detail::towgs84<T> & towgs84)365 inline bool pj_datum_prepare_towgs84(srs::detail::towgs84<T> & towgs84)
366 {
367     if( towgs84.size() == 7
368      && (towgs84[3] != 0.0
369       || towgs84[4] != 0.0
370       || towgs84[5] != 0.0
371       || towgs84[6] != 0.0) )
372     {
373         static const T sec_to_rad = detail::sec_to_rad<T>();
374 
375         /* transform from arc seconds to radians */
376         towgs84[3] *= sec_to_rad;
377         towgs84[4] *= sec_to_rad;
378         towgs84[5] *= sec_to_rad;
379         /* transform from parts per million to scaling factor */
380         towgs84[6] = (towgs84[6]/1000000.0) + 1;
381         return true;
382     }
383     else
384     {
385         return false;
386     }
387 }
388 
389 /************************************************************************/
390 /*                            pj_datum_init()                           */
391 /************************************************************************/
392 
393 // This function works differently than the original pj_datum_set().
394 // It doesn't push parameters defined in datum into params list.
395 // Instead it tries to use nadgrids and towgs84 and only then
396 // falls back to nadgrid or towgs84 defiend in datum parameter.
397 template <typename Params, typename T>
pj_datum_init(Params const & params,parameters<T> & projdef)398 inline void pj_datum_init(Params const& params, parameters<T> & projdef)
399 {
400     projdef.datum_type = datum_unknown;
401 
402     // Check for nadgrids parameter.
403     if(pj_datum_find_nadgrids(params, projdef.nadgrids))
404     {
405         // NOTE: It's different than in the original proj4.
406         // Nadgrids names are stored in projection definition.
407 
408         projdef.datum_type = datum_gridshift;
409     }
410     // Check for towgs84 parameter.
411     else if(pj_datum_find_towgs84(params, projdef.datum_params))
412     {
413         if (pj_datum_prepare_towgs84(projdef.datum_params))
414         {
415             projdef.datum_type = datum_7param;
416         }
417         else
418         {
419             projdef.datum_type = datum_3param;
420         }
421 
422         /* Note that pj_init() will later switch datum_type to
423            PJD_WGS84 if shifts are all zero, and ellipsoid is WGS84 or GRS80 */
424     }
425     // Check for datum parameter.
426     else
427     {
428         const pj_datums_type<T>* datum = pj_datum_find_datum<T>(params);
429         if (datum != NULL)
430         {
431             if (! datum->nadgrids.empty())
432             {
433                 projdef.nadgrids = datum->nadgrids;
434                 projdef.datum_type = datum_gridshift;
435             }
436             else if ( ! datum->towgs84.empty() )
437             {
438                 projdef.datum_params = datum->towgs84;
439                 if (pj_datum_prepare_towgs84(projdef.datum_params))
440                 {
441                     projdef.datum_type = datum_7param;
442                 }
443                 else
444                 {
445                     projdef.datum_type = datum_3param;
446                 }
447             }
448         }
449     }
450 }
451 
452 } // namespace detail
453 }}} // namespace boost::geometry::projections
454 
455 #endif // BOOST_GEOMETRY_PROJECTIONS_IMPL_PJ_DATUM_SET_HPP
456