File indexing completed on 2025-01-18 09:36:41
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010 #ifndef BOOST_GEOMETRY_SRS_TRANSFORMATION_HPP
0011 #define BOOST_GEOMETRY_SRS_TRANSFORMATION_HPP
0012
0013
0014 #include <string>
0015 #include <type_traits>
0016
0017 #include <boost/range/size.hpp>
0018 #include <boost/throw_exception.hpp>
0019
0020 #include <boost/geometry/algorithms/convert.hpp>
0021
0022 #include <boost/geometry/core/coordinate_dimension.hpp>
0023 #include <boost/geometry/core/static_assert.hpp>
0024
0025 #include <boost/geometry/geometries/point.hpp>
0026 #include <boost/geometry/geometries/multi_point.hpp>
0027 #include <boost/geometry/geometries/linestring.hpp>
0028 #include <boost/geometry/geometries/ring.hpp>
0029 #include <boost/geometry/geometries/segment.hpp>
0030
0031 #include <boost/geometry/srs/projection.hpp>
0032 #include <boost/geometry/srs/projections/grids.hpp>
0033 #include <boost/geometry/srs/projections/impl/pj_transform.hpp>
0034
0035 #include <boost/geometry/views/detail/indexed_point_view.hpp>
0036
0037
0038 namespace boost { namespace geometry
0039 {
0040
0041 namespace projections { namespace detail
0042 {
0043
0044 template <typename T1, typename T2>
0045 inline bool same_object(T1 const& , T2 const& )
0046 {
0047 return false;
0048 }
0049
0050 template <typename T>
0051 inline bool same_object(T const& o1, T const& o2)
0052 {
0053 return boost::addressof(o1) == boost::addressof(o2);
0054 }
0055
0056 template
0057 <
0058 typename PtIn,
0059 typename PtOut,
0060 bool SameUnits = std::is_same
0061 <
0062 typename geometry::detail::cs_angular_units<PtIn>::type,
0063 typename geometry::detail::cs_angular_units<PtOut>::type
0064 >::value
0065 >
0066 struct transform_geometry_point_coordinates
0067 {
0068 static inline void apply(PtIn const& in, PtOut & out, bool )
0069 {
0070 geometry::set<0>(out, geometry::get<0>(in));
0071 geometry::set<1>(out, geometry::get<1>(in));
0072 }
0073 };
0074
0075 template <typename PtIn, typename PtOut>
0076 struct transform_geometry_point_coordinates<PtIn, PtOut, false>
0077 {
0078 static inline void apply(PtIn const& in, PtOut & out, bool enable_angles)
0079 {
0080 if (enable_angles)
0081 {
0082 geometry::set_from_radian<0>(out, geometry::get_as_radian<0>(in));
0083 geometry::set_from_radian<1>(out, geometry::get_as_radian<1>(in));
0084 }
0085 else
0086 {
0087 geometry::set<0>(out, geometry::get<0>(in));
0088 geometry::set<1>(out, geometry::get<1>(in));
0089 }
0090 }
0091 };
0092
0093 template <typename Geometry, typename CT>
0094 struct transform_geometry_point
0095 {
0096 typedef typename geometry::point_type<Geometry>::type point_type;
0097
0098 typedef geometry::model::point
0099 <
0100 typename select_most_precise
0101 <
0102 typename geometry::coordinate_type<point_type>::type,
0103 CT
0104 >::type,
0105 geometry::dimension<point_type>::type::value,
0106 typename geometry::coordinate_system<point_type>::type
0107 > type;
0108
0109 template <typename PtIn, typename PtOut>
0110 static inline void apply(PtIn const& in, PtOut & out, bool enable_angles)
0111 {
0112 transform_geometry_point_coordinates<PtIn, PtOut>::apply(in, out, enable_angles);
0113 projections::detail::copy_higher_dimensions<2>(in, out);
0114 }
0115 };
0116
0117 template <typename Geometry, typename CT>
0118 struct transform_geometry_range_base
0119 {
0120 struct convert_strategy
0121 {
0122 convert_strategy(bool enable_angles)
0123 : m_enable_angles(enable_angles)
0124 {}
0125
0126 template <typename PtIn, typename PtOut>
0127 void apply(PtIn const& in, PtOut & out)
0128 {
0129 transform_geometry_point<Geometry, CT>::apply(in, out, m_enable_angles);
0130 }
0131
0132 bool m_enable_angles;
0133 };
0134
0135 template <typename In, typename Out>
0136 static inline void apply(In const& in, Out & out, bool enable_angles)
0137 {
0138
0139
0140
0141
0142
0143 geometry::detail::conversion::range_to_range
0144 <
0145 In,
0146 Out,
0147 geometry::point_order<In>::value != geometry::point_order<Out>::value
0148 >::apply(in, out, convert_strategy(enable_angles));
0149 }
0150 };
0151
0152 template
0153 <
0154 typename Geometry,
0155 typename CT,
0156 typename Tag = typename geometry::tag<Geometry>::type
0157 >
0158 struct transform_geometry
0159 {};
0160
0161 template <typename Point, typename CT>
0162 struct transform_geometry<Point, CT, point_tag>
0163 : transform_geometry_point<Point, CT>
0164 {};
0165
0166 template <typename Segment, typename CT>
0167 struct transform_geometry<Segment, CT, segment_tag>
0168 {
0169 typedef geometry::model::segment
0170 <
0171 typename transform_geometry_point<Segment, CT>::type
0172 > type;
0173
0174 template <typename In, typename Out>
0175 static inline void apply(In const& in, Out & out, bool enable_angles)
0176 {
0177 apply<0>(in, out, enable_angles);
0178 apply<1>(in, out, enable_angles);
0179 }
0180
0181 private:
0182 template <std::size_t Index, typename In, typename Out>
0183 static inline void apply(In const& in, Out & out, bool enable_angles)
0184 {
0185 geometry::detail::indexed_point_view<In const, Index> in_pt(in);
0186 geometry::detail::indexed_point_view<Out, Index> out_pt(out);
0187 transform_geometry_point<Segment, CT>::apply(in_pt, out_pt, enable_angles);
0188 }
0189 };
0190
0191 template <typename MultiPoint, typename CT>
0192 struct transform_geometry<MultiPoint, CT, multi_point_tag>
0193 : transform_geometry_range_base<MultiPoint, CT>
0194 {
0195 typedef model::multi_point
0196 <
0197 typename transform_geometry_point<MultiPoint, CT>::type
0198 > type;
0199 };
0200
0201 template <typename LineString, typename CT>
0202 struct transform_geometry<LineString, CT, linestring_tag>
0203 : transform_geometry_range_base<LineString, CT>
0204 {
0205 typedef model::linestring
0206 <
0207 typename transform_geometry_point<LineString, CT>::type
0208 > type;
0209 };
0210
0211 template <typename Ring, typename CT>
0212 struct transform_geometry<Ring, CT, ring_tag>
0213 : transform_geometry_range_base<Ring, CT>
0214 {
0215 typedef model::ring
0216 <
0217 typename transform_geometry_point<Ring, CT>::type,
0218 geometry::point_order<Ring>::value == clockwise,
0219 geometry::closure<Ring>::value == closed
0220 > type;
0221 };
0222
0223
0224 template
0225 <
0226 typename OutGeometry,
0227 typename CT,
0228 bool EnableTemporary = ! std::is_same
0229 <
0230 typename select_most_precise
0231 <
0232 typename geometry::coordinate_type<OutGeometry>::type,
0233 CT
0234 >::type,
0235 typename geometry::coordinate_type<OutGeometry>::type
0236 >::value
0237 >
0238 struct transform_geometry_wrapper
0239 {
0240 typedef transform_geometry<OutGeometry, CT> transform;
0241 typedef typename transform::type type;
0242
0243 template <typename InGeometry>
0244 transform_geometry_wrapper(InGeometry const& in, OutGeometry & out, bool input_angles)
0245 : m_out(out)
0246 {
0247 transform::apply(in, m_temp, input_angles);
0248 }
0249
0250 type & get() { return m_temp; }
0251 void finish() { geometry::convert(m_temp, m_out); }
0252
0253 private:
0254 type m_temp;
0255 OutGeometry & m_out;
0256 };
0257
0258 template
0259 <
0260 typename OutGeometry,
0261 typename CT
0262 >
0263 struct transform_geometry_wrapper<OutGeometry, CT, false>
0264 {
0265 typedef transform_geometry<OutGeometry, CT> transform;
0266 typedef OutGeometry type;
0267
0268 transform_geometry_wrapper(OutGeometry const& in, OutGeometry & out, bool input_angles)
0269 : m_out(out)
0270 {
0271 if (! same_object(in, out))
0272 transform::apply(in, out, input_angles);
0273 }
0274
0275 template <typename InGeometry>
0276 transform_geometry_wrapper(InGeometry const& in, OutGeometry & out, bool input_angles)
0277 : m_out(out)
0278 {
0279 transform::apply(in, out, input_angles);
0280 }
0281
0282 OutGeometry & get() { return m_out; }
0283 void finish() {}
0284
0285 private:
0286 OutGeometry & m_out;
0287 };
0288
0289 template <typename CT>
0290 struct transform_range
0291 {
0292 template
0293 <
0294 typename Proj1, typename Par1,
0295 typename Proj2, typename Par2,
0296 typename RangeIn, typename RangeOut,
0297 typename Grids
0298 >
0299 static inline bool apply(Proj1 const& proj1, Par1 const& par1,
0300 Proj2 const& proj2, Par2 const& par2,
0301 RangeIn const& in, RangeOut & out,
0302 Grids const& grids1, Grids const& grids2)
0303 {
0304
0305 bool const input_angles = !par1.is_geocent && par1.is_latlong;
0306
0307 transform_geometry_wrapper<RangeOut, CT> wrapper(in, out, input_angles);
0308
0309 bool res = true;
0310 try
0311 {
0312 res = pj_transform(proj1, par1, proj2, par2, wrapper.get(), grids1, grids2);
0313 }
0314 catch (projection_exception const&)
0315 {
0316 res = false;
0317 }
0318 catch(...)
0319 {
0320 BOOST_RETHROW
0321 }
0322
0323 wrapper.finish();
0324
0325 return res;
0326 }
0327 };
0328
0329 template <typename Policy>
0330 struct transform_multi
0331 {
0332 template
0333 <
0334 typename Proj1, typename Par1,
0335 typename Proj2, typename Par2,
0336 typename MultiIn, typename MultiOut,
0337 typename Grids
0338 >
0339 static inline bool apply(Proj1 const& proj1, Par1 const& par1,
0340 Proj2 const& proj2, Par2 const& par2,
0341 MultiIn const& in, MultiOut & out,
0342 Grids const& grids1, Grids const& grids2)
0343 {
0344 if (! same_object(in, out))
0345 range::resize(out, boost::size(in));
0346
0347 return apply(proj1, par1, proj2, par2,
0348 boost::begin(in), boost::end(in),
0349 boost::begin(out),
0350 grids1, grids2);
0351 }
0352
0353 private:
0354 template
0355 <
0356 typename Proj1, typename Par1,
0357 typename Proj2, typename Par2,
0358 typename InIt, typename OutIt,
0359 typename Grids
0360 >
0361 static inline bool apply(Proj1 const& proj1, Par1 const& par1,
0362 Proj2 const& proj2, Par2 const& par2,
0363 InIt in_first, InIt in_last, OutIt out_first,
0364 Grids const& grids1, Grids const& grids2)
0365 {
0366 bool res = true;
0367 for ( ; in_first != in_last ; ++in_first, ++out_first )
0368 {
0369 if ( ! Policy::apply(proj1, par1, proj2, par2, *in_first, *out_first, grids1, grids2) )
0370 {
0371 res = false;
0372 }
0373 }
0374 return res;
0375 }
0376 };
0377
0378 template
0379 <
0380 typename Geometry,
0381 typename CT,
0382 typename Tag = typename geometry::tag<Geometry>::type
0383 >
0384 struct transform
0385 : not_implemented<Tag>
0386 {};
0387
0388 template <typename Point, typename CT>
0389 struct transform<Point, CT, point_tag>
0390 {
0391 template
0392 <
0393 typename Proj1, typename Par1,
0394 typename Proj2, typename Par2,
0395 typename PointIn, typename PointOut,
0396 typename Grids
0397 >
0398 static inline bool apply(Proj1 const& proj1, Par1 const& par1,
0399 Proj2 const& proj2, Par2 const& par2,
0400 PointIn const& in, PointOut & out,
0401 Grids const& grids1, Grids const& grids2)
0402 {
0403
0404 bool const input_angles = !par1.is_geocent && par1.is_latlong;
0405
0406 transform_geometry_wrapper<PointOut, CT> wrapper(in, out, input_angles);
0407
0408 typedef typename transform_geometry_wrapper<PointOut, CT>::type point_type;
0409 point_type * ptr = boost::addressof(wrapper.get());
0410
0411 std::pair<point_type *, point_type *> range = std::make_pair(ptr, ptr + 1);
0412
0413 bool res = true;
0414 try
0415 {
0416 res = pj_transform(proj1, par1, proj2, par2, range, grids1, grids2);
0417 }
0418 catch (projection_exception const&)
0419 {
0420 res = false;
0421 }
0422 catch(...)
0423 {
0424 BOOST_RETHROW
0425 }
0426
0427 wrapper.finish();
0428
0429 return res;
0430 }
0431 };
0432
0433 template <typename MultiPoint, typename CT>
0434 struct transform<MultiPoint, CT, multi_point_tag>
0435 : transform_range<CT>
0436 {};
0437
0438 template <typename Segment, typename CT>
0439 struct transform<Segment, CT, segment_tag>
0440 {
0441 template
0442 <
0443 typename Proj1, typename Par1,
0444 typename Proj2, typename Par2,
0445 typename SegmentIn, typename SegmentOut,
0446 typename Grids
0447 >
0448 static inline bool apply(Proj1 const& proj1, Par1 const& par1,
0449 Proj2 const& proj2, Par2 const& par2,
0450 SegmentIn const& in, SegmentOut & out,
0451 Grids const& grids1, Grids const& grids2)
0452 {
0453
0454 bool const input_angles = !par1.is_geocent && par1.is_latlong;
0455
0456 transform_geometry_wrapper<SegmentOut, CT> wrapper(in, out, input_angles);
0457
0458 typedef typename geometry::point_type
0459 <
0460 typename transform_geometry_wrapper<SegmentOut, CT>::type
0461 >::type point_type;
0462
0463 point_type points[2];
0464
0465 geometry::detail::assign_point_from_index<0>(wrapper.get(), points[0]);
0466 geometry::detail::assign_point_from_index<1>(wrapper.get(), points[1]);
0467
0468 std::pair<point_type*, point_type*> range = std::make_pair(points, points + 2);
0469
0470 bool res = true;
0471 try
0472 {
0473 res = pj_transform(proj1, par1, proj2, par2, range, grids1, grids2);
0474 }
0475 catch (projection_exception const&)
0476 {
0477 res = false;
0478 }
0479 catch(...)
0480 {
0481 BOOST_RETHROW
0482 }
0483
0484 geometry::detail::assign_point_to_index<0>(points[0], wrapper.get());
0485 geometry::detail::assign_point_to_index<1>(points[1], wrapper.get());
0486
0487 wrapper.finish();
0488
0489 return res;
0490 }
0491 };
0492
0493 template <typename Linestring, typename CT>
0494 struct transform<Linestring, CT, linestring_tag>
0495 : transform_range<CT>
0496 {};
0497
0498 template <typename MultiLinestring, typename CT>
0499 struct transform<MultiLinestring, CT, multi_linestring_tag>
0500 : transform_multi<transform_range<CT> >
0501 {};
0502
0503 template <typename Ring, typename CT>
0504 struct transform<Ring, CT, ring_tag>
0505 : transform_range<CT>
0506 {};
0507
0508 template <typename Polygon, typename CT>
0509 struct transform<Polygon, CT, polygon_tag>
0510 {
0511 template
0512 <
0513 typename Proj1, typename Par1,
0514 typename Proj2, typename Par2,
0515 typename PolygonIn, typename PolygonOut,
0516 typename Grids
0517 >
0518 static inline bool apply(Proj1 const& proj1, Par1 const& par1,
0519 Proj2 const& proj2, Par2 const& par2,
0520 PolygonIn const& in, PolygonOut & out,
0521 Grids const& grids1, Grids const& grids2)
0522 {
0523 bool r1 = transform_range
0524 <
0525 CT
0526 >::apply(proj1, par1, proj2, par2,
0527 geometry::exterior_ring(in),
0528 geometry::exterior_ring(out),
0529 grids1, grids2);
0530 bool r2 = transform_multi
0531 <
0532 transform_range<CT>
0533 >::apply(proj1, par1, proj2, par2,
0534 geometry::interior_rings(in),
0535 geometry::interior_rings(out),
0536 grids1, grids2);
0537 return r1 && r2;
0538 }
0539 };
0540
0541 template <typename MultiPolygon, typename CT>
0542 struct transform<MultiPolygon, CT, multi_polygon_tag>
0543 : transform_multi
0544 <
0545 transform
0546 <
0547 typename boost::range_value<MultiPolygon>::type,
0548 CT,
0549 polygon_tag
0550 >
0551 >
0552 {};
0553
0554
0555 }}
0556
0557 namespace srs
0558 {
0559
0560
0561
0562
0563
0564
0565
0566
0567
0568
0569 template
0570 <
0571 typename Proj1 = srs::dynamic,
0572 typename Proj2 = srs::dynamic,
0573 typename CT = double
0574 >
0575 class transformation
0576 {
0577 typedef typename projections::detail::promote_to_double<CT>::type calc_t;
0578
0579 public:
0580
0581 transformation()
0582 {}
0583
0584
0585 template
0586 <
0587 typename Parameters1,
0588 std::enable_if_t
0589 <
0590 std::is_same<Proj1, srs::dynamic>::value
0591 && projections::dynamic_parameters<Parameters1>::is_specialized,
0592 int
0593 > = 0
0594 >
0595 explicit transformation(Parameters1 const& parameters1)
0596 : m_proj1(parameters1)
0597 {}
0598
0599
0600 explicit transformation(Proj1 const& parameters1)
0601 : m_proj1(parameters1)
0602 {}
0603
0604
0605 template
0606 <
0607 typename Parameters1, typename Parameters2,
0608 std::enable_if_t
0609 <
0610 std::is_same<Proj1, srs::dynamic>::value
0611 && std::is_same<Proj2, srs::dynamic>::value
0612 && projections::dynamic_parameters<Parameters1>::is_specialized
0613 && projections::dynamic_parameters<Parameters2>::is_specialized,
0614 int
0615 > = 0
0616 >
0617 transformation(Parameters1 const& parameters1,
0618 Parameters2 const& parameters2)
0619 : m_proj1(parameters1)
0620 , m_proj2(parameters2)
0621 {}
0622
0623
0624 template
0625 <
0626 typename Parameters1,
0627 std::enable_if_t
0628 <
0629 std::is_same<Proj1, srs::dynamic>::value
0630 && projections::dynamic_parameters<Parameters1>::is_specialized,
0631 int
0632 > = 0
0633 >
0634 transformation(Parameters1 const& parameters1,
0635 Proj2 const& parameters2)
0636 : m_proj1(parameters1)
0637 , m_proj2(parameters2)
0638 {}
0639
0640
0641 template
0642 <
0643 typename Parameters2,
0644 std::enable_if_t
0645 <
0646 std::is_same<Proj2, srs::dynamic>::value
0647 && projections::dynamic_parameters<Parameters2>::is_specialized,
0648 int
0649 > = 0
0650 >
0651 transformation(Proj1 const& parameters1,
0652 Parameters2 const& parameters2)
0653 : m_proj1(parameters1)
0654 , m_proj2(parameters2)
0655 {}
0656
0657
0658 transformation(Proj1 const& parameters1,
0659 Proj2 const& parameters2)
0660 : m_proj1(parameters1)
0661 , m_proj2(parameters2)
0662 {}
0663
0664 template <typename GeometryIn, typename GeometryOut>
0665 bool forward(GeometryIn const& in, GeometryOut & out) const
0666 {
0667 return forward(in, out, transformation_grids<detail::empty_grids_storage>());
0668 }
0669
0670 template <typename GeometryIn, typename GeometryOut>
0671 bool inverse(GeometryIn const& in, GeometryOut & out) const
0672 {
0673 return inverse(in, out, transformation_grids<detail::empty_grids_storage>());
0674 }
0675
0676 template <typename GeometryIn, typename GeometryOut, typename GridsStorage>
0677 bool forward(GeometryIn const& in, GeometryOut & out,
0678 transformation_grids<GridsStorage> const& grids) const
0679 {
0680 BOOST_GEOMETRY_STATIC_ASSERT(
0681 (projections::detail::same_tags<GeometryIn, GeometryOut>::value),
0682 "Not supported combination of Geometries.",
0683 GeometryIn, GeometryOut);
0684
0685 return projections::detail::transform
0686 <
0687 GeometryOut,
0688 calc_t
0689 >::apply(m_proj1.proj(), m_proj1.proj().params(),
0690 m_proj2.proj(), m_proj2.proj().params(),
0691 in, out,
0692 grids.src_grids,
0693 grids.dst_grids);
0694 }
0695
0696 template <typename GeometryIn, typename GeometryOut, typename GridsStorage>
0697 bool inverse(GeometryIn const& in, GeometryOut & out,
0698 transformation_grids<GridsStorage> const& grids) const
0699 {
0700 BOOST_GEOMETRY_STATIC_ASSERT(
0701 (projections::detail::same_tags<GeometryIn, GeometryOut>::value),
0702 "Not supported combination of Geometries.",
0703 GeometryIn, GeometryOut);
0704
0705 return projections::detail::transform
0706 <
0707 GeometryOut,
0708 calc_t
0709 >::apply(m_proj2.proj(), m_proj2.proj().params(),
0710 m_proj1.proj(), m_proj1.proj().params(),
0711 in, out,
0712 grids.dst_grids,
0713 grids.src_grids);
0714 }
0715
0716 template <typename GridsStorage>
0717 inline transformation_grids<GridsStorage> initialize_grids(GridsStorage & grids_storage) const
0718 {
0719 transformation_grids<GridsStorage> result(grids_storage);
0720
0721 using namespace projections::detail;
0722 pj_gridlist_from_nadgrids(m_proj1.proj().params(),
0723 result.src_grids);
0724 pj_gridlist_from_nadgrids(m_proj2.proj().params(),
0725 result.dst_grids);
0726
0727 return result;
0728 }
0729
0730 private:
0731 projections::proj_wrapper<Proj1, CT> m_proj1;
0732 projections::proj_wrapper<Proj2, CT> m_proj2;
0733 };
0734
0735
0736 }
0737
0738
0739 }}
0740
0741
0742 #endif