projection.hpp 17 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592
  1. // Boost.Geometry (aka GGL, Generic Geometry Library)
  2. // Copyright (c) 2008-2012 Barend Gehrels, Amsterdam, the Netherlands.
  3. // This file was modified by Oracle on 2017.
  4. // Modifications copyright (c) 2017, Oracle and/or its affiliates.
  5. // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
  6. // Use, modification and distribution is subject to the Boost Software License,
  7. // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  8. // http://www.boost.org/LICENSE_1_0.txt)
  9. #ifndef BOOST_GEOMETRY_SRS_PROJECTION_HPP
  10. #define BOOST_GEOMETRY_SRS_PROJECTION_HPP
  11. #include <string>
  12. #include <boost/geometry/algorithms/convert.hpp>
  13. #include <boost/geometry/algorithms/detail/convert_point_to_point.hpp>
  14. #include <boost/geometry/core/coordinate_dimension.hpp>
  15. #include <boost/geometry/srs/projections/exception.hpp>
  16. #include <boost/geometry/srs/projections/factory.hpp>
  17. #include <boost/geometry/srs/projections/impl/base_dynamic.hpp>
  18. #include <boost/geometry/srs/projections/impl/base_static.hpp>
  19. #include <boost/geometry/srs/projections/impl/pj_init.hpp>
  20. #include <boost/geometry/srs/projections/invalid_point.hpp>
  21. #include <boost/geometry/srs/projections/par4.hpp>
  22. #include <boost/geometry/srs/projections/proj4.hpp>
  23. #include <boost/geometry/views/detail/indexed_point_view.hpp>
  24. #include <boost/mpl/assert.hpp>
  25. #include <boost/mpl/if.hpp>
  26. #include <boost/smart_ptr/shared_ptr.hpp>
  27. #include <boost/throw_exception.hpp>
  28. #include <boost/type_traits/is_integral.hpp>
  29. #include <boost/type_traits/is_same.hpp>
  30. namespace boost { namespace geometry
  31. {
  32. namespace projections
  33. {
  34. #ifndef DOXYGEN_NO_DETAIL
  35. namespace detail
  36. {
  37. template <typename G1, typename G2>
  38. struct same_tags
  39. {
  40. static const bool value = boost::is_same
  41. <
  42. typename geometry::tag<G1>::type,
  43. typename geometry::tag<G2>::type
  44. >::value;
  45. };
  46. template <typename CT>
  47. struct promote_to_double
  48. {
  49. typedef typename boost::mpl::if_c
  50. <
  51. boost::is_integral<CT>::value || boost::is_same<CT, float>::value,
  52. double, CT
  53. >::type type;
  54. };
  55. // Copy coordinates of dimensions >= MinDim
  56. template <std::size_t MinDim, typename Point1, typename Point2>
  57. inline void copy_higher_dimensions(Point1 const& point1, Point2 & point2)
  58. {
  59. static const std::size_t dim1 = geometry::dimension<Point1>::value;
  60. static const std::size_t dim2 = geometry::dimension<Point2>::value;
  61. static const std::size_t lesser_dim = dim1 < dim2 ? dim1 : dim2;
  62. BOOST_MPL_ASSERT_MSG((lesser_dim >= MinDim),
  63. THE_DIMENSION_OF_POINTS_IS_TOO_SMALL,
  64. (Point1, Point2));
  65. geometry::detail::conversion::point_to_point
  66. <
  67. Point1, Point2, MinDim, lesser_dim
  68. > ::apply(point1, point2);
  69. // TODO: fill point2 with zeros if dim1 < dim2 ?
  70. // currently no need because equal dimensions are checked
  71. }
  72. struct forward_point_projection_policy
  73. {
  74. template <typename LL, typename XY, typename Proj>
  75. static inline bool apply(LL const& ll, XY & xy, Proj const& proj)
  76. {
  77. return proj.forward(ll, xy);
  78. }
  79. };
  80. struct inverse_point_projection_policy
  81. {
  82. template <typename XY, typename LL, typename Proj>
  83. static inline bool apply(XY const& xy, LL & ll, Proj const& proj)
  84. {
  85. return proj.inverse(xy, ll);
  86. }
  87. };
  88. template <typename PointPolicy>
  89. struct project_point
  90. {
  91. template <typename P1, typename P2, typename Proj>
  92. static inline bool apply(P1 const& p1, P2 & p2, Proj const& proj)
  93. {
  94. // (Geographic -> Cartesian) will be projected, rest will be copied.
  95. // So first copy third or higher dimensions
  96. projections::detail::copy_higher_dimensions<2>(p1, p2);
  97. if (! PointPolicy::apply(p1, p2, proj))
  98. {
  99. // For consistency with transformation
  100. set_invalid_point(p2);
  101. return false;
  102. }
  103. return true;
  104. }
  105. };
  106. template <typename PointPolicy>
  107. struct project_range
  108. {
  109. template <typename Proj>
  110. struct convert_policy
  111. {
  112. explicit convert_policy(Proj const& proj)
  113. : m_proj(proj)
  114. , m_result(true)
  115. {}
  116. template <typename Point1, typename Point2>
  117. inline void apply(Point1 const& point1, Point2 & point2)
  118. {
  119. if (! project_point<PointPolicy>::apply(point1, point2, m_proj) )
  120. m_result = false;
  121. }
  122. bool result() const
  123. {
  124. return m_result;
  125. }
  126. private:
  127. Proj const& m_proj;
  128. bool m_result;
  129. };
  130. template <typename R1, typename R2, typename Proj>
  131. static inline bool apply(R1 const& r1, R2 & r2, Proj const& proj)
  132. {
  133. return geometry::detail::conversion::range_to_range
  134. <
  135. R1, R2,
  136. geometry::point_order<R1>::value != geometry::point_order<R2>::value
  137. >::apply(r1, r2, convert_policy<Proj>(proj)).result();
  138. }
  139. };
  140. template <typename Policy>
  141. struct project_multi
  142. {
  143. template <typename G1, typename G2, typename Proj>
  144. static inline bool apply(G1 const& g1, G2 & g2, Proj const& proj)
  145. {
  146. range::resize(g2, boost::size(g1));
  147. return apply(boost::begin(g1), boost::end(g1),
  148. boost::begin(g2),
  149. proj);
  150. }
  151. private:
  152. template <typename It1, typename It2, typename Proj>
  153. static inline bool apply(It1 g1_first, It1 g1_last, It2 g2_first, Proj const& proj)
  154. {
  155. bool result = true;
  156. for ( ; g1_first != g1_last ; ++g1_first, ++g2_first )
  157. {
  158. if (! Policy::apply(*g1_first, *g2_first, proj))
  159. {
  160. result = false;
  161. }
  162. }
  163. return result;
  164. }
  165. };
  166. template
  167. <
  168. typename Geometry,
  169. typename PointPolicy,
  170. typename Tag = typename geometry::tag<Geometry>::type
  171. >
  172. struct project_geometry
  173. {};
  174. template <typename Geometry, typename PointPolicy>
  175. struct project_geometry<Geometry, PointPolicy, point_tag>
  176. : project_point<PointPolicy>
  177. {};
  178. template <typename Geometry, typename PointPolicy>
  179. struct project_geometry<Geometry, PointPolicy, multi_point_tag>
  180. : project_range<PointPolicy>
  181. {};
  182. template <typename Geometry, typename PointPolicy>
  183. struct project_geometry<Geometry, PointPolicy, segment_tag>
  184. {
  185. template <typename G1, typename G2, typename Proj>
  186. static inline bool apply(G1 const& g1, G2 & g2, Proj const& proj)
  187. {
  188. bool r1 = apply<0>(g1, g2, proj);
  189. bool r2 = apply<1>(g1, g2, proj);
  190. return r1 && r2;
  191. }
  192. private:
  193. template <std::size_t Index, typename G1, typename G2, typename Proj>
  194. static inline bool apply(G1 const& g1, G2 & g2, Proj const& proj)
  195. {
  196. geometry::detail::indexed_point_view<G1 const, Index> pt1(g1);
  197. geometry::detail::indexed_point_view<G2, Index> pt2(g2);
  198. return project_point<PointPolicy>::apply(pt1, pt2, proj);
  199. }
  200. };
  201. template <typename Geometry, typename PointPolicy>
  202. struct project_geometry<Geometry, PointPolicy, linestring_tag>
  203. : project_range<PointPolicy>
  204. {};
  205. template <typename Geometry, typename PointPolicy>
  206. struct project_geometry<Geometry, PointPolicy, multi_linestring_tag>
  207. : project_multi< project_range<PointPolicy> >
  208. {};
  209. template <typename Geometry, typename PointPolicy>
  210. struct project_geometry<Geometry, PointPolicy, ring_tag>
  211. : project_range<PointPolicy>
  212. {};
  213. template <typename Geometry, typename PointPolicy>
  214. struct project_geometry<Geometry, PointPolicy, polygon_tag>
  215. {
  216. template <typename G1, typename G2, typename Proj>
  217. static inline bool apply(G1 const& g1, G2 & g2, Proj const& proj)
  218. {
  219. bool r1 = project_range
  220. <
  221. PointPolicy
  222. >::apply(geometry::exterior_ring(g1),
  223. geometry::exterior_ring(g2),
  224. proj);
  225. bool r2 = project_multi
  226. <
  227. project_range<PointPolicy>
  228. >::apply(geometry::interior_rings(g1),
  229. geometry::interior_rings(g2),
  230. proj);
  231. return r1 && r2;
  232. }
  233. };
  234. template <typename MultiPolygon, typename PointPolicy>
  235. struct project_geometry<MultiPolygon, PointPolicy, multi_polygon_tag>
  236. : project_multi
  237. <
  238. project_geometry
  239. <
  240. typename boost::range_value<MultiPolygon>::type,
  241. PointPolicy,
  242. polygon_tag
  243. >
  244. >
  245. {};
  246. } // namespace detail
  247. #endif // DOXYGEN_NO_DETAIL
  248. template <typename Params, typename CT>
  249. struct dynamic_parameters
  250. {
  251. BOOST_MPL_ASSERT_MSG((false),
  252. NOT_IMPLEMENTED_FOR_THESE_PARAMETERS,
  253. (Params));
  254. };
  255. template <typename CT>
  256. struct dynamic_parameters<srs::proj4, CT>
  257. {
  258. static inline projections::parameters<CT> apply(srs::proj4 const& params)
  259. {
  260. return projections::detail::pj_init_plus<CT>(srs::dynamic(), params.str);
  261. }
  262. };
  263. // proj_wrapper class and its specializations wrapps the internal projection
  264. // representation and implements transparent creation of projection object
  265. template <typename Proj, typename CT>
  266. class proj_wrapper
  267. {
  268. BOOST_MPL_ASSERT_MSG((false),
  269. NOT_IMPLEMENTED_FOR_THIS_PROJECTION,
  270. (Proj));
  271. };
  272. template <typename CT>
  273. class proj_wrapper<srs::dynamic, CT>
  274. {
  275. // Some projections do not work with float -> wrong results
  276. // select <double> from int/float/double and else selects T
  277. typedef typename projections::detail::promote_to_double<CT>::type calc_t;
  278. typedef projections::parameters<calc_t> parameters_type;
  279. typedef projections::detail::base_v<calc_t, parameters_type> vprj_t;
  280. public:
  281. template <typename Params>
  282. proj_wrapper(Params const& params)
  283. : m_ptr(create(projections::dynamic_parameters<Params, calc_t>::apply(params)))
  284. {}
  285. vprj_t const& proj() const { return *m_ptr; }
  286. vprj_t & mutable_proj() { return *m_ptr; }
  287. private:
  288. static vprj_t* create(parameters_type const& pj_params)
  289. {
  290. vprj_t* result = projections::detail::create_new(pj_params);
  291. if (result == NULL)
  292. {
  293. if (pj_params.name.empty())
  294. {
  295. BOOST_THROW_EXCEPTION(projection_not_named_exception());
  296. }
  297. else
  298. {
  299. BOOST_THROW_EXCEPTION(projection_unknown_id_exception(pj_params.name));
  300. }
  301. }
  302. return result;
  303. }
  304. boost::shared_ptr<vprj_t> m_ptr;
  305. };
  306. template <typename StaticParameters, typename CT>
  307. class static_proj_wrapper_base
  308. {
  309. typedef typename projections::detail::promote_to_double<CT>::type calc_t;
  310. typedef projections::parameters<calc_t> parameters_type;
  311. typedef typename srs::par4::detail::pick_proj_tag
  312. <
  313. StaticParameters
  314. >::type proj_tag;
  315. typedef typename srs::par4::detail::pick_ellps
  316. <
  317. StaticParameters
  318. >::type ellps_type;
  319. typedef typename projections::detail::static_projection_type
  320. <
  321. proj_tag,
  322. typename geometry::tag
  323. <
  324. typename srs::par4::detail::ellps_traits
  325. <
  326. ellps_type
  327. >::model_type
  328. >::type,
  329. StaticParameters,
  330. calc_t,
  331. parameters_type
  332. >::type projection_type;
  333. public:
  334. projection_type const& proj() const { return m_proj; }
  335. projection_type & mutable_proj() { return m_proj; }
  336. protected:
  337. explicit static_proj_wrapper_base(StaticParameters const& s_params,
  338. bool use_defaults = true)
  339. : m_proj(get_parameters(s_params, "", use_defaults))
  340. {}
  341. static_proj_wrapper_base(StaticParameters const& s_params,
  342. srs::proj4 const& params,
  343. bool use_defaults = true)
  344. : m_proj(get_parameters(s_params, params.str, use_defaults))
  345. {}
  346. private:
  347. static parameters_type get_parameters(StaticParameters const& s_params,
  348. std::string const& params_str,
  349. bool use_defaults)
  350. {
  351. return projections::detail::pj_init_plus<calc_t>(s_params, params_str, use_defaults);
  352. }
  353. projection_type m_proj;
  354. };
  355. template <BOOST_GEOMETRY_PROJECTIONS_DETAIL_TYPENAME_PX, typename CT>
  356. class proj_wrapper<srs::static_proj4<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>, CT>
  357. : public static_proj_wrapper_base<srs::static_proj4<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>, CT>
  358. {
  359. typedef srs::static_proj4<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>
  360. static_parameters_type;
  361. typedef static_proj_wrapper_base
  362. <
  363. static_parameters_type,
  364. CT
  365. > base_t;
  366. public:
  367. proj_wrapper()
  368. : base_t(static_parameters_type())
  369. {}
  370. proj_wrapper(static_parameters_type const& s_params)
  371. : base_t(s_params)
  372. {}
  373. proj_wrapper(srs::proj4 const& params)
  374. : base_t(static_parameters_type(), params)
  375. {}
  376. proj_wrapper(static_parameters_type const& s_params,
  377. srs::proj4 const& params)
  378. : base_t(s_params, params)
  379. {}
  380. };
  381. // projection class implements transparent forward/inverse projection interface
  382. template <typename Proj, typename CT>
  383. class projection
  384. : private proj_wrapper<Proj, CT>
  385. {
  386. typedef proj_wrapper<Proj, CT> base_t;
  387. public:
  388. projection()
  389. {}
  390. template <typename Params>
  391. explicit projection(Params const& params)
  392. : base_t(params)
  393. {}
  394. template <typename SParams, typename Params>
  395. projection(SParams const& s_params, Params const& params)
  396. : base_t(s_params, params)
  397. {}
  398. /// Forward projection, from Latitude-Longitude to Cartesian
  399. template <typename LL, typename XY>
  400. inline bool forward(LL const& ll, XY& xy) const
  401. {
  402. BOOST_MPL_ASSERT_MSG((projections::detail::same_tags<LL, XY>::value),
  403. NOT_SUPPORTED_COMBINATION_OF_GEOMETRIES,
  404. (LL, XY));
  405. concepts::check_concepts_and_equal_dimensions<LL const, XY>();
  406. return projections::detail::project_geometry
  407. <
  408. LL,
  409. projections::detail::forward_point_projection_policy
  410. >::apply(ll, xy, base_t::proj());
  411. }
  412. /// Inverse projection, from Cartesian to Latitude-Longitude
  413. template <typename XY, typename LL>
  414. inline bool inverse(XY const& xy, LL& ll) const
  415. {
  416. BOOST_MPL_ASSERT_MSG((projections::detail::same_tags<XY, LL>::value),
  417. NOT_SUPPORTED_COMBINATION_OF_GEOMETRIES,
  418. (XY, LL));
  419. concepts::check_concepts_and_equal_dimensions<XY const, LL>();
  420. return projections::detail::project_geometry
  421. <
  422. XY,
  423. projections::detail::inverse_point_projection_policy
  424. >::apply(xy, ll, base_t::proj());
  425. }
  426. };
  427. } // namespace projections
  428. namespace srs
  429. {
  430. /*!
  431. \brief Representation of projection
  432. \details Either dynamic or static projection representation
  433. \ingroup projection
  434. \tparam Proj default_dynamic or static projection parameters
  435. \tparam CT calculation type used internally
  436. */
  437. template
  438. <
  439. typename Proj = srs::dynamic,
  440. typename CT = double
  441. >
  442. class projection
  443. {
  444. BOOST_MPL_ASSERT_MSG((false),
  445. NOT_IMPLEMENTED_FOR_THIS_PROJECTION,
  446. (Proj));
  447. };
  448. template <typename CT>
  449. class projection<srs::dynamic, CT>
  450. : public projections::projection<srs::dynamic, CT>
  451. {
  452. typedef projections::projection<srs::dynamic, CT> base_t;
  453. public:
  454. /*!
  455. \ingroup projection
  456. \brief Initializes a projection as a string, using the format with + and =
  457. \details The projection can be initialized with a string (with the same format as the PROJ4 package) for
  458. convenient initialization from, for example, the command line
  459. \par Example
  460. <tt>+proj=labrd +ellps=intl +lon_0=46d26'13.95E +lat_0=18d54S +azi=18d54 +k_0=.9995 +x_0=400000 +y_0=800000</tt>
  461. for the Madagascar projection.
  462. \note Parameters are described in the group
  463. */
  464. template <typename Params>
  465. projection(Params const& params)
  466. : base_t(params)
  467. {}
  468. };
  469. template <BOOST_GEOMETRY_PROJECTIONS_DETAIL_TYPENAME_PX, typename CT>
  470. class projection<srs::static_proj4<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>, CT>
  471. : public projections::projection<srs::static_proj4<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>, CT>
  472. {
  473. typedef projections::projection<srs::static_proj4<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX>, CT> base_t;
  474. public:
  475. projection()
  476. {}
  477. projection(srs::static_proj4<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX> const& params)
  478. : base_t(params)
  479. {}
  480. #ifdef BOOST_GEOMETRY_SRS_ENABLE_STATIC_PROJECTION_HYBRID_INTERFACE
  481. projection(srs::proj4 const& params)
  482. : base_t(params)
  483. {}
  484. projection(srs::static_proj4<BOOST_GEOMETRY_PROJECTIONS_DETAIL_PX> const& s_params,
  485. srs::proj4 const& params)
  486. : base_t(s_params, params)
  487. {}
  488. #endif
  489. };
  490. } // namespace srs
  491. }} // namespace boost::geometry
  492. #endif // BOOST_GEOMETRY_SRS_PROJECTION_HPP