projections/impl/pj_datum_set.hpp

0001 // Boost.Geometry (aka GGL, Generic Geometry Library)
0002 // This file is manually converted from PROJ4
0003
0004 // Copyright (c) 2008-2012 Barend Gehrels, Amsterdam, the Netherlands.
0005 // Copyright (c) 2023 Adam Wulkiewicz, Lodz, Poland.
0006
0007 // This file was modified by Oracle on 2017-2020.
0008 // Modifications copyright (c) 2017-2020, Oracle and/or its affiliates.
0009 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
0010
0011 // Use, modification and distribution is subject to the Boost Software License,
0012 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
0013 // http://www.boost.org/LICENSE_1_0.txt)
0014
0015 // This file is converted from PROJ4, http://trac.osgeo.org/proj
0016 // PROJ4 is originally written by Gerald Evenden (then of the USGS)
0017 // PROJ4 is maintained by Frank Warmerdam
0018 // PROJ4 is converted to Geometry Library by Barend Gehrels (Geodan, Amsterdam)
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_IMPL_PJ_DATUM_SET_HPP
0041 #define BOOST_GEOMETRY_PROJECTIONS_IMPL_PJ_DATUM_SET_HPP
0042
0043
0044 #include <string>
0045 #include <vector>
0046
0047 #include <boost/algorithm/string/classification.hpp>
0048 #include <boost/algorithm/string/split.hpp>
0049
0050 #include <boost/geometry/srs/projections/dpar.hpp>
0051 #include <boost/geometry/srs/projections/exception.hpp>
0052 #include <boost/geometry/srs/projections/impl/projects.hpp>
0053 #include <boost/geometry/srs/projections/impl/pj_datums.hpp>
0054 #include <boost/geometry/srs/projections/impl/pj_param.hpp>
0055 #include <boost/geometry/srs/projections/proj4.hpp>
0056 #include <boost/geometry/srs/projections/spar.hpp>
0057
0058
0059 namespace boost { namespace geometry { namespace projections {
0060
0061 namespace detail {
0062
0063
0064 /* SEC_TO_RAD = Pi/180/3600 */
0065 template <typename T>
0066 inline T sec_to_rad() { return 4.84813681109535993589914102357e-6; }
0067
0068 /************************************************************************/
0069 /*                        pj_datum_find_datum()                         */
0070 /************************************************************************/
0071
0072 template <typename T>
0073 inline const pj_datums_type<T>* pj_datum_find_datum(srs::detail::proj4_parameters const& params)
0074 {
0075     std::string name = pj_get_param_s(params, "datum");
0076     if(! name.empty())
0077     {
0078         /* find the datum definition */
0079         const pj_datums_type<T>* pj_datums = pj_get_datums<T>().first;
0080         const int n = pj_get_datums<T>().second;
0081         int index = -1;
0082         for (int i = 0; i < n && index == -1; i++)
0083         {
0084             if(pj_datums[i].id == name)
0085             {
0086                 index = i;
0087             }
0088         }
0089
0090         if (index != -1)
0091         {
0092             return pj_datums + index;
0093         }
0094         else
0095         {
0096             BOOST_THROW_EXCEPTION( projection_exception(error_unknown_ellp_param) );
0097         }
0098     }
0099
0100     return NULL;
0101 }
0102
0103 template <typename T>
0104 inline const pj_datums_type<T>* pj_datum_find_datum(srs::dpar::parameters<T> const& params)
0105 {
0106     typename srs::dpar::parameters<T>::const_iterator
0107         it = pj_param_find(params, srs::dpar::datum);
0108
0109     if (it != params.end())
0110     {
0111         const pj_datums_type<T>* pj_datums = pj_get_datums<T>().first;
0112         const int n = pj_get_datums<T>().second;
0113         int i = it->template get_value<int>();
0114         if (i >= 0 && i < n)
0115         {
0116             return pj_datums + i;
0117         }
0118         else
0119         {
0120             BOOST_THROW_EXCEPTION( projection_exception(error_unknown_ellp_param) );
0121         }
0122     }
0123
0124     return NULL;
0125 }
0126
0127 template
0128 <
0129     typename Params,
0130     typename Param = typename geometry::tuples::find_if
0131         <
0132             Params,
0133             srs::spar::detail::is_param_tr<srs::spar::detail::datum_traits>::pred
0134         >::type,
0135     bool IsFound = geometry::tuples::is_found<Param>::value
0136 >
0137 struct pj_datum_find_datum_static
0138 {
0139     template <typename T>
0140     static const pj_datums_type<T>* apply(Params const& )
0141     {
0142         const pj_datums_type<T>* pj_datums = pj_get_datums<T>().first;
0143         const int n = pj_get_datums<T>().second;
0144         const int i = srs::spar::detail::datum_traits<Param>::id;
0145         if (i >= 0 && i < n)
0146         {
0147             return pj_datums + i;
0148         }
0149         else
0150         {
0151             // TODO: Implemnt as BOOST_GEOMETRY_STATIC_ASSERT instead
0152             BOOST_THROW_EXCEPTION( projection_exception(error_unknown_ellp_param) );
0153         }
0154     }
0155 };
0156 template <typename Params, typename Param>
0157 struct pj_datum_find_datum_static<Params, Param, false>
0158 {
0159     template <typename T>
0160     static const pj_datums_type<T>* apply(Params const& )
0161     {
0162         return NULL;
0163     }
0164 };
0165
0166 template <typename T, typename ...Ps>
0167 inline const pj_datums_type<T>* pj_datum_find_datum(srs::spar::parameters<Ps...> const& params)
0168 {
0169     return pj_datum_find_datum_static
0170         <
0171             srs::spar::parameters<Ps...>
0172         >::template apply<T>(params);
0173 }
0174
0175 /************************************************************************/
0176 /*                        pj_datum_find_nadgrids()                      */
0177 /************************************************************************/
0178
0179 inline bool pj_datum_find_nadgrids(srs::detail::proj4_parameters const& params,
0180                                    srs::detail::nadgrids & out)
0181 {
0182     std::string snadgrids = pj_get_param_s(params, "nadgrids");
0183     if (! snadgrids.empty())
0184     {
0185         for (std::string::size_type i = 0 ; i < snadgrids.size() ; )
0186         {
0187             std::string::size_type end = snadgrids.find(',', i);
0188             std::string name = snadgrids.substr(i, end - i);
0189
0190             i = end;
0191             if (end != std::string::npos)
0192                 ++i;
0193
0194             if (! name.empty())
0195                 out.push_back(name);
0196         }
0197     }
0198
0199     return ! out.empty();
0200 }
0201
0202 template <typename T>
0203 inline bool pj_datum_find_nadgrids(srs::dpar::parameters<T> const& params,
0204                                    srs::detail::nadgrids & out)
0205 {
0206     typename srs::dpar::parameters<T>::const_iterator
0207         it = pj_param_find(params, srs::dpar::nadgrids);
0208     if (it != params.end())
0209     {
0210         out = it->template get_value<srs::detail::nadgrids>();
0211     }
0212
0213     return ! out.empty();
0214 }
0215
0216 template
0217 <
0218     typename Params,
0219     int I = geometry::tuples::find_index_if
0220         <
0221             Params,
0222             srs::spar::detail::is_param<srs::spar::nadgrids>::pred
0223         >::value,
0224     int N = geometry::tuples::size<Params>::value
0225 >
0226 struct pj_datum_find_nadgrids_static
0227 {
0228     static void apply(Params const& params, srs::detail::nadgrids & out)
0229     {
0230         out = geometry::tuples::get<I>(params);
0231     }
0232 };
0233 template <typename Params, int N>
0234 struct pj_datum_find_nadgrids_static<Params, N, N>
0235 {
0236     static void apply(Params const& , srs::detail::nadgrids & )
0237     {}
0238 };
0239
0240 template <typename ...Ps>
0241 inline bool pj_datum_find_nadgrids(srs::spar::parameters<Ps...> const& params,
0242                                    srs::detail::nadgrids & out)
0243 {
0244     pj_datum_find_nadgrids_static
0245         <
0246             srs::spar::parameters<Ps...>
0247         >::apply(params, out);
0248
0249     return ! out.empty();
0250 }
0251
0252 /************************************************************************/
0253 /*                        pj_datum_find_towgs84()                       */
0254 /************************************************************************/
0255
0256 template <typename T>
0257 inline bool pj_datum_find_towgs84(srs::detail::proj4_parameters const& params,
0258                                   srs::detail::towgs84<T> & out)
0259 {
0260     std::string towgs84 = pj_get_param_s(params, "towgs84");
0261     if(! towgs84.empty())
0262     {
0263         std::vector<std::string> parm;
0264         boost::split(parm, towgs84, boost::is_any_of(" ,"));
0265
0266         std::size_t n = (std::min<std::size_t>)(parm.size(), 7);
0267         std::size_t z = n <= 3 ? 3 : 7;
0268
0269         /* parse out the pvalues */
0270         for (std::size_t i = 0 ; i < n; ++i)
0271         {
0272             out.push_back(geometry::str_cast<T>(parm[i]));
0273         }
0274         for (std::size_t i = out.size() ; i < z; ++i)
0275         {
0276             out.push_back(T(0));
0277         }
0278     }
0279
0280     return ! out.empty();
0281 }
0282
0283 template <typename T>
0284 inline bool pj_datum_find_towgs84(srs::dpar::parameters<T> const& params,
0285                                   srs::detail::towgs84<T> & out)
0286 {
0287     typename srs::dpar::parameters<T>::const_iterator
0288         it = pj_param_find(params, srs::dpar::towgs84);
0289
0290     if (it != params.end())
0291     {
0292         srs::detail::towgs84<T> const&
0293             towgs84 = it->template get_value<srs::detail::towgs84<T> >();
0294
0295         std::size_t n = (std::min<std::size_t>)(towgs84.size(), 7u);
0296         std::size_t z = n <= 3 ? 3 : 7;
0297
0298         for (std::size_t i = 0 ; i < n; ++i)
0299         {
0300             out.push_back(towgs84[i]);
0301         }
0302         for (std::size_t i = out.size() ; i < z; ++i)
0303         {
0304             out.push_back(T(0));
0305         }
0306     }
0307
0308     return ! out.empty();
0309 }
0310
0311 template
0312 <
0313     typename Params,
0314     int I = geometry::tuples::find_index_if
0315         <
0316             Params,
0317             srs::spar::detail::is_param_t<srs::spar::towgs84>::pred
0318         >::value,
0319     int N = geometry::tuples::size<Params>::value
0320 >
0321 struct pj_datum_find_towgs84_static
0322 {
0323     template <typename T>
0324     static void apply(Params const& params, srs::detail::towgs84<T> & out)
0325     {
0326         typename geometry::tuples::element<I, Params>::type const&
0327             towgs84 = geometry::tuples::get<I>(params);
0328
0329         std::size_t n = (std::min<std::size_t>)(towgs84.size(), 7u);
0330         std::size_t z = n <= 3 ? 3 : 7;
0331
0332         for (std::size_t i = 0 ; i < n; ++i)
0333         {
0334             out.push_back(towgs84[i]);
0335         }
0336         for (std::size_t i = out.size() ; i < z; ++i)
0337         {
0338             out.push_back(T(0));
0339         }
0340     }
0341 };
0342 template <typename Params, int N>
0343 struct pj_datum_find_towgs84_static<Params, N, N>
0344 {
0345     template <typename T>
0346     static void apply(Params const& , srs::detail::towgs84<T> & )
0347     {}
0348 };
0349
0350 template <typename T, typename ...Ps>
0351 inline bool pj_datum_find_towgs84(srs::spar::parameters<Ps...> const& params,
0352                                   srs::detail::towgs84<T> & out)
0353 {
0354     pj_datum_find_towgs84_static
0355         <
0356             srs::spar::parameters<Ps...>
0357         >::apply(params, out);
0358
0359     return ! out.empty();
0360 }
0361
0362 /************************************************************************/
0363 /*                        pj_datum_prepare_towgs84()                    */
0364 /************************************************************************/
0365
0366 template <typename T>
0367 inline bool pj_datum_prepare_towgs84(srs::detail::towgs84<T> & towgs84)
0368 {
0369     if( towgs84.size() == 7
0370      && (towgs84[3] != 0.0
0371       || towgs84[4] != 0.0
0372       || towgs84[5] != 0.0
0373       || towgs84[6] != 0.0) )
0374     {
0375         static const T sec_to_rad = detail::sec_to_rad<T>();
0376
0377         /* transform from arc seconds to radians */
0378         towgs84[3] *= sec_to_rad;
0379         towgs84[4] *= sec_to_rad;
0380         towgs84[5] *= sec_to_rad;
0381         /* transform from parts per million to scaling factor */
0382         towgs84[6] = (towgs84[6]/1000000.0) + 1;
0383         return true;
0384     }
0385     else
0386     {
0387         return false;
0388     }
0389 }
0390
0391 /************************************************************************/
0392 /*                            pj_datum_init()                           */
0393 /************************************************************************/
0394
0395 // This function works differently than the original pj_datum_set().
0396 // It doesn't push parameters defined in datum into params list.
0397 // Instead it tries to use nadgrids and towgs84 and only then
0398 // falls back to nadgrid or towgs84 defiend in datum parameter.
0399 template <typename Params, typename T>
0400 inline void pj_datum_init(Params const& params, parameters<T> & projdef)
0401 {
0402     projdef.datum_type = datum_unknown;
0403
0404     // Check for nadgrids parameter.
0405     if(pj_datum_find_nadgrids(params, projdef.nadgrids))
0406     {
0407         // NOTE: It's different than in the original proj4.
0408         // Nadgrids names are stored in projection definition.
0409
0410         projdef.datum_type = datum_gridshift;
0411     }
0412     // Check for towgs84 parameter.
0413     else if(pj_datum_find_towgs84(params, projdef.datum_params))
0414     {
0415         if (pj_datum_prepare_towgs84(projdef.datum_params))
0416         {
0417             projdef.datum_type = datum_7param;
0418         }
0419         else
0420         {
0421             projdef.datum_type = datum_3param;
0422         }
0423
0424         /* Note that pj_init() will later switch datum_type to
0425            PJD_WGS84 if shifts are all zero, and ellipsoid is WGS84 or GRS80 */
0426     }
0427     // Check for datum parameter.
0428     else
0429     {
0430         const pj_datums_type<T>* datum = pj_datum_find_datum<T>(params);
0431         if (datum != NULL)
0432         {
0433             if (! datum->nadgrids.empty())
0434             {
0435                 projdef.nadgrids = datum->nadgrids;
0436                 projdef.datum_type = datum_gridshift;
0437             }
0438             else if ( ! datum->towgs84.empty() )
0439             {
0440                 projdef.datum_params = datum->towgs84;
0441                 if (pj_datum_prepare_towgs84(projdef.datum_params))
0442                 {
0443                     projdef.datum_type = datum_7param;
0444                 }
0445                 else
0446                 {
0447                     projdef.datum_type = datum_3param;
0448                 }
0449             }
0450         }
0451     }
0452 }
0453
0454 } // namespace detail
0455 }}} // namespace boost::geometry::projections
0456
0457 #endif // BOOST_GEOMETRY_PROJECTIONS_IMPL_PJ_DATUM_SET_HPP