boost_mirror/geometry
1
0
Fork 0
mirror of https://github.com/boostorg/geometry.git synced 2025-05-09 23:24:02 +00:00
Code Issues Packages Projects Releases Wiki Activity

[closest_points] Several fixes: indentation, formating, pass correct types, activate tests

Browse Source
This commit is contained in:
Vissarion Fisikopoulos 2022-01-05 14:21:10 +02:00
parent 588fe9e751
commit 1b0a1dfe61
11 changed files with 148 additions and 289 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

18
doc/readme.txt
View File

@ -32,14 +32,14 @@ src: examples used in documentation and tools (doxygen_xml2qbk)
Per new algorithm (e.g. foo), one should add: Per new algorithm (e.g. foo), one should add:
1) in file boost/geometry/algorithms/foo.hpp, include a "\ingroup foo" in the doxygen comments 1) in file boost/geometry/algorithms/foo.hpp, include a "\ingroup foo" in the doxygen comments
(if a directory is created e.g. boost/geometry/algorithms/detail/foo the path should be added in doc/doxy/Doxyfile) 2) if a directory is created e.g. boost/geometry/algorithms/detail/foo the path should be added in doc/doxy/Doxyfile
2) in file doc/doxy/doxygen_input/groups/groups.hpp, define the group "foo" 3) in file doc/doxy/doxygen_input/groups/groups.hpp, define the group "foo"
3) in file doc/make_qbk.py, include the algorithm "foo" 4) in file doc/make_qbk.py, include the algorithm "foo"
4) in file doc/reference.qbk, include the foo.qbk ([include generated/foo.qbk]) 5) in file doc/reference.qbk, include the foo.qbk ([include generated/foo.qbk])
5) in file doc/quickref.xml, include a section on foo conform other sections 6) in file doc/quickref.xml, include a section on foo conform other sections
6) in file doc/src/docutils/tools/support_status/support_status.cpp include the algorithm (3 places) (optionally) 7) in file doc/src/docutils/tools/support_status/support_status.cpp include the algorithm (3 places) (optionally)
7) in file doc/reference/foo.qbk (to be created), include the support status and write other text, and include examples (optionally) 8) in file doc/reference/foo.qbk (to be created), include the support status and write other text, and include examples (optionally)
8) in file doc/imports.qbk, include the example foo.cpp (if any) 9) in file doc/imports.qbk, include the example foo.cpp (if any)
9) create file doc/src/examples/algorithm/foo.cpp (optional) 10) create file doc/src/examples/algorithm/foo.cpp (optional)

2
include/boost/geometry/algorithms/detail/buffer/buffered_piece_collection.hpp
View File

@ -672,7 +672,7 @@ struct buffered_piece_collection
{ {
// Make sure the closing point is identical (they are calculated // Make sure the closing point is identical (they are calculated
// separately by different pieces) // separately by different pieces)
range::back(added) = range::front(added); range::back(added) = range::added);
} }
for (std::size_t i = first_piece_index; i < boost::size(m_pieces); i++) for (std::size_t i = first_piece_index; i < boost::size(m_pieces); i++)

40
include/boost/geometry/algorithms/detail/closest_points/linear_or_areal_to_areal.hpp
View File

@ -34,23 +34,20 @@ struct linear_to_areal
Segment& shortest_seg, Segment& shortest_seg,
Strategies const& strategies) Strategies const& strategies)
{ {
using most_precise_type = typename select_coordinate_type using most_precise_type = typename select_coordinate_type<Linear, Areal>::type;
<
Linear,
Areal
>::type;
using point_type = typename std::conditional using point_type = typename std::conditional
< <
std::is_same<coordinate_type<Linear>, most_precise_type>::value, std::is_same<typename coordinate_type<Linear>::type, most_precise_type>::value,
typename point_type<Linear>::type, typename point_type<Linear>::type,
typename point_type<Areal>::type typename point_type<Areal>::type
>::type; >::type;
using linestring_type = geometry::model::linestring<point_type>; using linestring_type = geometry::model::linestring<point_type>;
/* TODO: currently intersection does not support tupled input /* TODO: currently intersection does not support some cases of tupled input
* this should be implemented directly with dynamic geometries * such as linestring - multipolygon
* this could be implemented directly with dynamic geometries
using polygon_type = geometry::model::polygon<point_type>; using polygon_type = geometry::model::polygon<point_type>;
std::tuple std::tuple
< <
@ -111,8 +108,7 @@ struct segment_to_areal
Strategies const& strategies, Strategies const& strategies,
bool = false) bool = false)
{ {
using linestring_type = geometry::model::linestring using linestring_type = geometry::model::linestring<typename point_type<Segment>::type>;
<typename point_type<Segment>::type>;
linestring_type linestring; linestring_type linestring;
convert(segment, linestring); convert(segment, linestring);
linear_to_areal::apply(linestring, areal, shortest_seg, strategies); linear_to_areal::apply(linestring, areal, shortest_seg, strategies);
@ -141,18 +137,14 @@ struct areal_to_areal
Segment& shortest_seg, Segment& shortest_seg,
Strategies const& strategies) Strategies const& strategies)
{ {
using most_precise_type = typename select_coordinate_type using most_precise_type = typename select_coordinate_type<Areal1, Areal2>::type;
<
Areal1,
Areal2
>::type;
using point_type = typename std::conditional using point_type = typename std::conditional
< <
std::is_same<coordinate_type<Areal1>, most_precise_type>::value, std::is_same<typename coordinate_type<Areal1>::type, most_precise_type>::value,
typename point_type<Areal1>::type, typename point_type<Areal1>::type,
typename point_type<Areal2>::type typename point_type<Areal2>::type
>::type; >::type;
using linestring_type = geometry::model::linestring<point_type>; using linestring_type = geometry::model::linestring<point_type>;
using polygon_type = geometry::model::polygon<point_type>; using polygon_type = geometry::model::polygon<point_type>;

27
include/boost/geometry/algorithms/detail/closest_points/segment_to_segment.hpp
View File

@ -87,26 +87,15 @@ public:
closest_points::creturn_t<Segment1, Segment2, Strategies> d[4]; closest_points::creturn_t<Segment1, Segment2, Strategies> d[4];
//TODO: this could be simplified by getting distance strategy in a unique way auto const cds = strategies::distance::detail::make_comparable(strategies)
d[0] = strategy::distance::services::get_comparable .distance(detail::dummy_point(), detail::dummy_point());
<
decltype(strategies.distance(cp0, q[0])) d[0] = cds.apply(cp0, q[0]);
>::apply(strategies.distance(cp0, q[0])).apply(cp0, q[0]); d[1] = cds.apply(cp1, q[1]);
d[1] = strategy::distance::services::get_comparable d[2] = cds.apply(p[0], cp2);
< d[3] = cds.apply(p[1], cp3);
decltype(strategies.distance(cp1, q[1]))
>::apply(strategies.distance(cp1, q[1])).apply(cp1, q[1]);
d[2] = strategy::distance::services::get_comparable
<
decltype(strategies.distance(p[0], cp2))
>::apply(strategies.distance(p[0], cp2)).apply(p[0], cp2);
d[3] = strategy::distance::services::get_comparable
<
decltype(strategies.distance(p[1], cp3))
>::apply(strategies.distance(p[1], cp3)).apply(p[1], cp3);
std::size_t imin = std::distance(boost::addressof(d[0]), std::size_t imin = std::distance(boost::addressof(d[0]), std::min_element(d, d + 4));
std::min_element(d, d + 4));
switch (imin) switch (imin)
{ {

51
include/boost/geometry/strategies/closest_points/geographic.hpp
View File

@ -27,10 +27,6 @@ namespace boost { namespace geometry
namespace strategies { namespace closest_points namespace strategies { namespace closest_points
{ {
#ifndef DOXYGEN_NO_DETAIL
namespace detail
{
template template
< <
typename FormulaPolicy = geometry::strategy::andoyer, typename FormulaPolicy = geometry::strategy::andoyer,
@ -43,6 +39,7 @@ class geographic
using base_t = strategies::distance::geographic<FormulaPolicy, Spheroid, CalculationType>; using base_t = strategies::distance::geographic<FormulaPolicy, Spheroid, CalculationType>;
public: public:
geographic() = default; geographic() = default;
explicit geographic(Spheroid const& spheroid) explicit geographic(Spheroid const& spheroid)
@ -50,8 +47,8 @@ public:
{} {}
template <typename Geometry1, typename Geometry2> template <typename Geometry1, typename Geometry2>
auto closest_points(Geometry1 const&, Geometry2 const&, auto closest_points(Geometry1 const&, Geometry2 const&,
distance::detail::enable_if_ps_t<Geometry1, Geometry2> * = nullptr) const distance::detail::enable_if_ps_t<Geometry1, Geometry2> * = nullptr) const
{ {
return strategy::closest_points::geographic_cross_track return strategy::closest_points::geographic_cross_track
< <
@ -60,40 +57,6 @@ public:
CalculationType CalculationType
>(base_t::m_spheroid); >(base_t::m_spheroid);
} }
};
} // namespace detail
#endif // DOXYGEN_NO_DETAIL
template
<
typename FormulaPolicy = strategy::andoyer,
typename Spheroid = srs::spheroid<double>,
typename CalculationType = void
>
class geographic
: public strategies::closest_points::detail::geographic
<
FormulaPolicy,
Spheroid,
CalculationType
>
{
using base_t = strategies::closest_points::detail::geographic
<
FormulaPolicy,
Spheroid,
CalculationType
>;
public:
geographic() = default;
explicit geographic(Spheroid const& spheroid)
: base_t(spheroid)
{}
}; };
@ -101,13 +64,7 @@ namespace services
{ {
template <typename Geometry1, typename Geometry2> template <typename Geometry1, typename Geometry2>
struct default_strategy struct default_strategy<Geometry1, Geometry2, geographic_tag, geographic_tag>
<
Geometry1,
Geometry2,
geographic_tag,
geographic_tag
>
{ {
using type = strategies::closest_points::geographic<>; using type = strategies::closest_points::geographic<>;
}; };

52
include/boost/geometry/strategies/closest_points/spherical.hpp
View File

@ -10,15 +10,6 @@
#ifndef BOOST_GEOMETRY_STRATEGIES_CLOSEST_POINTS_SPHERICAL_HPP #ifndef BOOST_GEOMETRY_STRATEGIES_CLOSEST_POINTS_SPHERICAL_HPP
#define BOOST_GEOMETRY_STRATEGIES_CLOSEST_POINTS_SPHERICAL_HPP #define BOOST_GEOMETRY_STRATEGIES_CLOSEST_POINTS_SPHERICAL_HPP
//#include <boost/geometry/strategies/cartesian/azimuth.hpp>
//#include <boost/geometry/strategies/cartesian/distance_projected_point.hpp>
//#include <boost/geometry/strategies/cartesian/distance_pythagoras.hpp>
//#include <boost/geometry/strategies/spherical/distance_pythagoras_box_box.hpp>
//#include <boost/geometry/strategies/spherical/distance_pythagoras_point_box.hpp>
//#include <boost/geometry/strategies/spherical/distance_segment_box.hpp>
#include <boost/geometry/strategies/spherical/closest_points_pt_seg.hpp> #include <boost/geometry/strategies/spherical/closest_points_pt_seg.hpp>
#include <boost/geometry/strategies/detail.hpp> #include <boost/geometry/strategies/detail.hpp>
@ -36,15 +27,15 @@ namespace boost { namespace geometry
namespace strategies { namespace closest_points namespace strategies { namespace closest_points
{ {
#ifndef DOXYGEN_NO_DETAIL template
namespace detail <
{ typename RadiusTypeOrSphere = double,
typename CalculationType = void
template <typename RadiusTypeOrSphere, typename CalculationType> >
class spherical class spherical
: public strategies::distance::detail::spherical<RadiusTypeOrSphere, CalculationType> : public strategies::distance::spherical<RadiusTypeOrSphere, CalculationType>
{ {
using base_t = strategies::distance::detail::spherical<RadiusTypeOrSphere, CalculationType>; using base_t = strategies::distance::spherical<RadiusTypeOrSphere, CalculationType>;
public: public:
spherical() = default; spherical() = default;
@ -56,35 +47,10 @@ public:
template <typename Geometry1, typename Geometry2> template <typename Geometry1, typename Geometry2>
auto closest_points(Geometry1 const&, Geometry2 const&, auto closest_points(Geometry1 const&, Geometry2 const&,
distance::detail::enable_if_ps_t<Geometry1, Geometry2> * = nullptr) const distance::detail::enable_if_ps_t<Geometry1, Geometry2> * = nullptr) const
{ {
return strategy::closest_points::cross_track<CalculationType>( return strategy::closest_points::cross_track<CalculationType>(base_t::radius());
base_t::radius());
} }
};
} // namespace detail
#endif // DOXYGEN_NO_DETAIL
template
<
typename RadiusTypeOrSphere = double,
typename CalculationType = void
>
class spherical
: public strategies::closest_points::detail::spherical<RadiusTypeOrSphere, CalculationType>
{
using base_t = strategies::closest_points::detail::spherical<RadiusTypeOrSphere, CalculationType>;
public:
spherical() = default;
template <typename RadiusOrSphere>
explicit spherical(RadiusOrSphere const& radius_or_sphere)
: base_t(radius_or_sphere)
{}
}; };

15
include/boost/geometry/strategies/geographic/closest_points_pt_seg.hpp
View File

@ -85,6 +85,19 @@ class geographic_cross_track
false, false,
true true
>; >;
template <typename Point, typename PointOfSegment>
struct calculation_type
: promote_floating_point
<
typename select_calculation_type
<
Point,
PointOfSegment,
CalculationType
>::type
>
{};
public : public :
explicit geographic_cross_track(Spheroid const& spheroid = Spheroid()) explicit geographic_cross_track(Spheroid const& spheroid = Spheroid())
@ -103,7 +116,7 @@ public :
model::point model::point
< <
typename base_t::template return_type<Point, PointOfSegment>::type, typename calculation_type<Point, PointOfSegment>::type,
dimension<PointOfSegment>::value, dimension<PointOfSegment>::value,
typename coordinate_system<PointOfSegment>::type typename coordinate_system<PointOfSegment>::type
> cp; > cp;

175
include/boost/geometry/strategies/geographic/distance_cross_track.hpp
View File

@ -65,7 +65,7 @@ namespace strategy { namespace distance
namespace detail namespace detail
{ {
template <bool> template <bool EnableClosestPoint>
struct set_result struct set_result
{ {
template <typename CT, typename ResultType> template <typename CT, typename ResultType>
@ -165,8 +165,7 @@ public:
strategy::default_order<FormulaPolicy>::value strategy::default_order<FormulaPolicy>::value
>; >;
return meridian_inverse::meridian_not_crossing_pole_dist(lat1, lat2, return meridian_inverse::meridian_not_crossing_pole_dist(lat1, lat2, m_spheroid);
m_spheroid);
} }
private: private:
@ -192,16 +191,14 @@ private:
if (g4 < -1.25*pi)//close to -270 if (g4 < -1.25*pi)//close to -270
{ {
#ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK #ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK
std::cout << "g4=" << g4 * math::r2d<CT>() << ", close to -270" std::cout << "g4=" << g4 * math::r2d<CT>() << ", close to -270" << std::endl;
<< std::endl;
#endif #endif
return g4 + 1.5 * pi; return g4 + 1.5 * pi;
} }
else if (g4 > 1.25*pi)//close to 270 else if (g4 > 1.25*pi)//close to 270
{ {
#ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK #ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK
std::cout << "g4=" << g4 * math::r2d<CT>() << ", close to 270" std::cout << "g4=" << g4 * math::r2d<CT>() << ", close to 270" << std::endl;
<< std::endl;
#endif #endif
der = -der; der = -der;
return - g4 + 1.5 * pi; return - g4 + 1.5 * pi;
@ -209,8 +206,7 @@ private:
else if (g4 < 0 && g4 > -0.75*pi)//close to -90 else if (g4 < 0 && g4 > -0.75*pi)//close to -90
{ {
#ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK #ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK
std::cout << "g4=" << g4 * math::r2d<CT>() << ", close to -90" std::cout << "g4=" << g4 * math::r2d<CT>() << ", close to -90" << std::endl;
<< std::endl;
#endif #endif
der = -der; der = -der;
return -g4 - pi/2; return -g4 - pi/2;
@ -227,23 +223,10 @@ private:
CT const& s14_start, CT const& a12, CT const& s14_start, CT const& a12,
result_type<CT>& result) result_type<CT>& result)
{ {
using direct_distance_type = typename FormulaPolicy::template direct using direct_distance_type =
< typename FormulaPolicy::template direct<CT, true, false, false, false>;
CT, using inverse_distance_type =
true, typename FormulaPolicy::template inverse<CT, true, false, false, false, false>;
false,
false,
false
>;
using inverse_distance_type = typename FormulaPolicy::template inverse
<
CT,
true,
false,
false,
false,
false
>;
geometry::formula::result_direct<CT> res14; geometry::formula::result_direct<CT> res14;
@ -263,40 +246,34 @@ private:
#ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK #ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK
std::cout << "dist(pl,p3)=" std::cout << "dist(pl,p3)="
<< inverse_distance_type::apply(lon3, lat3, << inverse_distance_type::apply(lon3, lat3, pr_lon, pr_lat, spheroid).distance
pr_lon, pr_lat,
spheroid).distance
<< std::endl; << std::endl;
std::cout << "dist(pr,p3)=" std::cout << "dist(pr,p3)="
<< inverse_distance_type::apply(lon3, lat3, << inverse_distance_type::apply(lon3, lat3, pr_lon, pr_lat, spheroid).distance
pr_lon, pr_lat,
spheroid).distance
<< std::endl; << std::endl;
#endif #endif
CT const dist_l = inverse_distance_type::apply CT const dist_l =
(lon3, lat3, pl_lon, pl_lat, spheroid).distance; inverse_distance_type::apply(lon3, lat3, pl_lon, pl_lat, spheroid).distance;
CT const dist_r = inverse_distance_type::apply CT const dist_r =
(lon3, lat3, pr_lon, pr_lat, spheroid).distance; inverse_distance_type::apply(lon3, lat3, pr_lon, pr_lat, spheroid).distance;
if (dist_l < dist_r) if (dist_l < dist_r)
{ {
s14 -= inverse_distance_type::apply(res14.lon2, res14.lat2, s14 -= inverse_distance_type::apply(res14.lon2, res14.lat2, pl_lon,
pl_lon, pl_lat, pl_lat, spheroid).distance/2;
spheroid).distance/2;
pr_lon = res14.lon2; pr_lon = res14.lon2;
pr_lat = res14.lat2; pr_lat = res14.lat2;
} }
else else
{ {
s14 += inverse_distance_type::apply(res14.lon2, res14.lat2, s14 += inverse_distance_type::apply(res14.lon2, res14.lat2, pr_lon,
pr_lon, pr_lat, pr_lat, spheroid).distance/2;
spheroid).distance/2;
pl_lon = res14.lon2; pl_lon = res14.lon2;
pl_lat = res14.lat2; pl_lat = res14.lat2;
} }
CT new_distance = inverse_distance_type::apply(lon3, lat3, CT new_distance = inverse_distance_type::apply(lon3, lat3, res14.lon2, res14.lat2,
res14.lon2, res14.lat2, spheroid).distance;
spheroid).distance;
dist_improve = new_distance != result.distance; dist_improve = new_distance != result.distance;
@ -316,16 +293,14 @@ private:
<< std::endl<< std::endl; << std::endl<< std::endl;
if (!dist_improve) if (!dist_improve)
{ {
std::cout << "Stop msg: res34.distance >= prev_distance" std::cout << "Stop msg: res34.distance >= prev_distance" << std::endl;
<< std::endl;
} }
if (counter == BOOST_GEOMETRY_DETAIL_POINT_SEGMENT_DISTANCE_MAX_STEPS) if (counter == BOOST_GEOMETRY_DETAIL_POINT_SEGMENT_DISTANCE_MAX_STEPS)
{ {
std::cout << "Stop msg: counter" << std::endl; std::cout << "Stop msg: counter" << std::endl;
} }
#endif #endif
set_result<EnableClosestPoint>::apply(new_distance, res14.lon2, set_result<EnableClosestPoint>::apply(new_distance, res14.lon2, res14.lat2, result);
res14.lat2, result);
} while (dist_improve && counter++ } while (dist_improve && counter++
< BOOST_GEOMETRY_DETAIL_POINT_SEGMENT_DISTANCE_MAX_STEPS); < BOOST_GEOMETRY_DETAIL_POINT_SEGMENT_DISTANCE_MAX_STEPS);
@ -339,14 +314,14 @@ private:
CT const& s14_start, CT const& a12, CT const& s14_start, CT const& a12,
result_type<CT>& result) result_type<CT>& result)
{ {
using inverse_distance_azimuth_quantities_type = typename using inverse_distance_azimuth_quantities_type =
FormulaPolicy::template inverse<CT, true, true, false, true, true>; typename FormulaPolicy::template inverse<CT, true, true, false, true, true>;
using inverse_dist_azimuth_type = typename FormulaPolicy::template using inverse_dist_azimuth_type =
inverse<CT, false, true, false, false, false>; typename FormulaPolicy::template inverse<CT, false, true, false, false, false>;
using direct_distance_type = typename FormulaPolicy::template using direct_distance_type =
direct<CT, true, false, false, false>; typename FormulaPolicy::template direct<CT, true, false, false, false>;
CT const half_pi = math::pi<CT>() / CT(2); CT const half_pi = math::pi<CT>() / CT(2);
geometry::formula::result_direct<CT> res14; geometry::formula::result_direct<CT> res14;
@ -370,10 +345,10 @@ private:
// g4 is the angle between segment (p1,p2) and segment (p3,p4) // g4 is the angle between segment (p1,p2) and segment (p3,p4)
// that meet on p4 (GEO) // that meet on p4 (GEO)
CT a4 = inverse_dist_azimuth_type CT a4 = inverse_dist_azimuth_type::apply(res14.lon2, res14.lat2,
::apply(res14.lon2, res14.lat2, lon2, lat2, spheroid).azimuth; lon2, lat2, spheroid).azimuth;
res34 = inverse_distance_azimuth_quantities_type res34 = inverse_distance_azimuth_quantities_type::apply(res14.lon2, res14.lat2,
::apply(res14.lon2, res14.lat2, lon3, lat3, spheroid); lon3, lat3, spheroid);
g4 = res34.azimuth - a4; g4 = res34.azimuth - a4;
// cos(s14/earth_radius) is the spherical limit // cos(s14/earth_radius) is the spherical limit
@ -391,13 +366,13 @@ private:
if (dist_improve) if (dist_improve)
{ {
set_result<EnableClosestPoint>::apply(res34.distance, set_result<EnableClosestPoint>::apply(res34.distance, res14.lon2, res14.lat2,
res14.lon2, res14.lat2, result); result);
} }
else else
{ {
set_result<EnableClosestPoint>::apply(prev_distance, set_result<EnableClosestPoint>::apply(prev_distance, prev_res.lon2, prev_res.lat2,
prev_res.lon2, prev_res.lat2, result); result);
} }
#ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK #ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK
@ -420,8 +395,7 @@ private:
} }
if (!dist_improve) if (!dist_improve)
{ {
std::cout << "Stop msg: res34.distance >= prev_distance" std::cout << "Stop msg: res34.distance >= prev_distance" << std::endl;
<< std::endl;
} }
if (delta_g4 == 0) if (delta_g4 == 0)
{ {
@ -436,8 +410,7 @@ private:
} while (g4 != half_pi } while (g4 != half_pi
&& dist_improve && dist_improve
&& delta_g4 != 0 && delta_g4 != 0
&& counter++ && counter++ < BOOST_GEOMETRY_DETAIL_POINT_SEGMENT_DISTANCE_MAX_STEPS);
< BOOST_GEOMETRY_DETAIL_POINT_SEGMENT_DISTANCE_MAX_STEPS);
#ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK #ifdef BOOST_GEOMETRY_DEBUG_GEOGRAPHIC_CROSS_TRACK
std::cout << "distance=" << res34.distance << std::endl; std::cout << "distance=" << res34.distance << std::endl;
@ -449,22 +422,25 @@ private:
<< res14.lat2 * math::r2d<CT>() << ")" << res14.lat2 * math::r2d<CT>() << ")"
<< std::endl; << std::endl;
CT s31 = inverse_distance_azimuth_quantities_type CT s31 = inverse_distance_azimuth_quantities_type::apply(lon3, lat3, lon1, lat1, spheroid)
::apply(lon3, lat3, lon1, lat1, spheroid).distance; .distance;
CT s32 = inverse_distance_azimuth_quantities_type CT s32 = inverse_distance_azimuth_quantities_type::apply(lon3, lat3, lon2, lat2, spheroid)
::apply(lon3, lat3, lon2, lat2, spheroid).distance; .distance;
CT a4 = inverse_dist_azimuth_type CT a4 = inverse_dist_azimuth_type::apply(res14.lon2, res14.lat2, lon2, lat2, spheroid)
::apply(res14.lon2, res14.lat2, lon2, lat2, spheroid).azimuth; .azimuth;
geometry::formula::result_direct<CT> res4 = geometry::formula::result_direct<CT> res4 =
direct_distance_type::apply(res14.lon2, res14.lat2, .04, a4, spheroid); direct_distance_type::apply(res14.lon2, res14.lat2, .04, a4, spheroid);
CT p4_plus = inverse_distance_azimuth_quantities_type
::apply(res4.lon2, res4.lat2, lon3, lat3, spheroid).distance; CT p4_plus = inverse_distance_azimuth_quantities_type::apply(res4.lon2, res4.lat2, lon3,
lat3, spheroid).distance;
geometry::formula::result_direct<CT> res1 geometry::formula::result_direct<CT> res1 =
= direct_distance_type::apply(lon1, lat1, s14-.04, a12, spheroid); direct_distance_type::apply(lon1, lat1, s14-.04, a12, spheroid);
CT p4_minus = inverse_distance_azimuth_quantities_type
::apply(res1.lon2, res1.lat2, lon3, lat3, spheroid).distance; CT p4_minus = inverse_distance_azimuth_quantities_type::apply(res1.lon2, res1.lat2, lon3,
lat3, spheroid).distance;
std::cout << "s31=" << s31 << "\ns32=" << s32 std::cout << "s31=" << s31 << "\ns32=" << s32
<< "\np4_plus=" << p4_plus << ", p4=(" << "\np4_plus=" << p4_plus << ", p4=("
@ -522,11 +498,11 @@ protected:
CT const& lo3, CT const& la3, //query point p3 CT const& lo3, CT const& la3, //query point p3
Spheroid const& spheroid) Spheroid const& spheroid)
{ {
using inverse_dist_azimuth_type = typename FormulaPolicy::template using inverse_dist_azimuth_type =
inverse<CT, true, true, false, false, false>; typename FormulaPolicy::template inverse<CT, true, true, false, false, false>;
using inverse_dist_azimuth_reverse_type = typename FormulaPolicy:: using inverse_dist_azimuth_reverse_type =
template inverse<CT, true, true, true, false, false>; typename FormulaPolicy::template inverse<CT, true, true, true, false, false>;
CT const earth_radius = geometry::formula::mean_radius<CT>(spheroid); CT const earth_radius = geometry::formula::mean_radius<CT>(spheroid);
@ -572,8 +548,7 @@ protected:
//segment on equator //segment on equator
//Note: antipodal points on equator does not define segment on equator //Note: antipodal points on equator does not define segment on equator
//but pass by the pole //but pass by the pole
CT diff = geometry::math::longitude_distance_signed<geometry CT diff = geometry::math::longitude_distance_signed<geometry::radian>(lon1, lon2);
::radian>(lon1, lon2);
using meridian_inverse = typename formula::meridian_inverse<CT>; using meridian_inverse = typename formula::meridian_inverse<CT>;
@ -622,20 +597,16 @@ protected:
#endif #endif
CT sign_non_zero = lat3 >= c0 ? 1 : -1; CT sign_non_zero = lat3 >= c0 ? 1 : -1;
auto res13 = apply(lon1, lat1, lon1, half_pi * sign_non_zero, lon3, auto res13 = apply(lon1, lat1, lon1, half_pi * sign_non_zero, lon3, lat3, spheroid);
lat3, spheroid);
auto res23 = apply(lon2, lat2, lon2, half_pi * sign_non_zero, lon3, auto res23 = apply(lon2, lat2, lon2, half_pi * sign_non_zero, lon3, lat3, spheroid);
lat3, spheroid);
return (res13.distance) < (res23.distance) ? res13 : res23; return (res13.distance) < (res23.distance) ? res13 : res23;
} }
auto res12 = inverse_dist_azimuth_reverse_type::apply(lon1, lat1, lon2, auto res12 = inverse_dist_azimuth_reverse_type::apply(lon1, lat1, lon2, lat2, spheroid);
lat2, spheroid);
auto res13 = inverse_dist_azimuth_type::apply(lon1, lat1, lon3, lat3, auto res13 = inverse_dist_azimuth_type::apply(lon1, lat1, lon3, lat3, spheroid);
spheroid);
if (geometry::math::equals(res12.distance, c0)) if (geometry::math::equals(res12.distance, c0))
{ {
@ -644,8 +615,7 @@ protected:
std::cout << "distance between points=" std::cout << "distance between points="
<< res13.distance << std::endl; << res13.distance << std::endl;
#endif #endif
typename meridian_inverse::result res = meridian_inverse::apply( auto res = meridian_inverse::apply(lon1, lat1, lon3, lat3, spheroid);
lon1, lat1, lon3, lat3, spheroid);
return non_iterative_case(lon3, lat3, lon1, lat2, return non_iterative_case(lon3, lat3, lon1, lat2,
res.meridian ? res.distance : res13.distance); res.meridian ? res.distance : res13.distance);
@ -657,7 +627,7 @@ protected:
// TODO: meridian case optimization // TODO: meridian case optimization
if (geometry::math::equals(a312, c0) && meridian_not_crossing_pole) if (geometry::math::equals(a312, c0) && meridian_not_crossing_pole)
{ {
std::tuple<CT,CT> minmax_elem = std::minmax(lat1, lat2); auto const minmax_elem = std::minmax(lat1, lat2);
if (lat3 >= std::get<0>(minmax_elem) && if (lat3 >= std::get<0>(minmax_elem) &&
lat3 <= std::get<1>(minmax_elem)) lat3 <= std::get<1>(minmax_elem))
@ -690,8 +660,7 @@ protected:
return non_iterative_case(lon3, lat3, lon1, lat1, spheroid); return non_iterative_case(lon3, lat3, lon1, lat1, spheroid);
} }
auto res23 = inverse_dist_azimuth_type::apply(lon2, lat2, lon3, lat3, auto res23 = inverse_dist_azimuth_type::apply(lon2, lat2, lon3, lat3, spheroid);
spheroid);
CT a321 = res23.azimuth - res12.reverse_azimuth + pi; CT a321 = res23.azimuth - res12.reverse_azimuth + pi;
CT projection2 = cos( a321 ) * res23.distance / res12.distance; CT projection2 = cos( a321 ) * res23.distance / res12.distance;
@ -733,14 +702,13 @@ protected:
geometry::strategy::distance::haversine<CT> str(earth_radius); geometry::strategy::distance::haversine<CT> str(earth_radius);
CT s13_sph = str.apply(p1, p3); CT s13_sph = str.apply(p1, p3);
//CT s14 = acos( cos(s13/earth_radius) / cos(s34/earth_radius) ) //CT s14 = acos( cos(s13/earth_radius) / cos(s34/earth_radius) ) * earth_radius;
// * earth_radius;
CT cos_frac = cos(s13_sph / earth_radius) / cos(s34_sph / earth_radius); CT cos_frac = cos(s13_sph / earth_radius) / cos(s34_sph / earth_radius);
CT s14_sph = cos_frac >= 1 ? CT(0) CT s14_sph = cos_frac >= 1 ? CT(0)
: cos_frac <= -1 ? pi * earth_radius : cos_frac <= -1 ? pi * earth_radius
: acos(cos_frac) * earth_radius; : acos(cos_frac) * earth_radius;
CT a12_sph = geometry::formula::spherical_azimuth<>(lon1, lat1, lon2, lat2); CT const a12_sph = geometry::formula::spherical_azimuth<>(lon1, lat1, lon2, lat2);
auto res = geometry::formula::spherical_direct<true, false>(lon1, lat1, auto res = geometry::formula::spherical_direct<true, false>(lon1, lat1,
s14_sph, a12_sph, srs::sphere<CT>(earth_radius)); s14_sph, a12_sph, srs::sphere<CT>(earth_radius));
@ -772,8 +740,7 @@ protected:
CT CT
>::apply(lon1, lat1, res.lon2, res.lat2, spheroid); >::apply(lon1, lat1, res.lon2, res.lat2, spheroid);
newton(lon1, lat1, lon2, lat2, lon3, lat3, spheroid, s14_start, newton(lon1, lat1, lon2, lat2, lon3, lat3, spheroid, s14_start, res12.azimuth, result);
res12.azimuth, result);
} }
return result; return result;

25
include/boost/geometry/strategies/spherical/distance_cross_track.hpp
View File

@ -68,7 +68,7 @@ namespace detail
CalculationType lon = geometry::get_as_radian<0>(p); CalculationType lon = geometry::get_as_radian<0>(p);
CalculationType lat = geometry::get_as_radian<1>(p); CalculationType lat = geometry::get_as_radian<1>(p);
CalculationType crs_AD = geometry::formula::spherical_azimuth CalculationType const crs_AD = geometry::formula::spherical_azimuth
< <
CalculationType, CalculationType,
false false
@ -549,28 +549,7 @@ template
class cross_track class cross_track
{ {
public : public :
using equals_point_point_strategy_type = within::spherical_point_point;
using relate_segment_segment_strategy_type = intersection::spherical_segments
<
CalculationType
>;
static inline relate_segment_segment_strategy_type get_relate_segment_segment_strategy()
{
return relate_segment_segment_strategy_type();
}
using point_in_geometry_strategy_type = within::spherical_winding
<
void, void, CalculationType
>;
static inline point_in_geometry_strategy_type get_point_in_geometry_strategy()
{
return point_in_geometry_strategy_type();
}
template <typename Point, typename PointOfSegment> template <typename Point, typename PointOfSegment>
struct return_type struct return_type
: promote_floating_point : promote_floating_point

12
test/algorithms/closest_points/empty_geometry.hpp
View File

@ -27,8 +27,7 @@ void test_empty_input(Geometry1 const& geometry1,
{ {
return; return;
} }
BOOST_CHECK_MESSAGE(false, BOOST_CHECK_MESSAGE(false, "A empty_input_exception should have been thrown");
"A empty_input_exception should have been thrown");
try try
{ {
@ -38,8 +37,7 @@ void test_empty_input(Geometry1 const& geometry1,
{ {
return; return;
} }
BOOST_CHECK_MESSAGE(false, BOOST_CHECK_MESSAGE(false, "A empty_input_exception should have been thrown");
"A empty_input_exception should have been thrown");
try try
{ {
@ -49,8 +47,7 @@ void test_empty_input(Geometry1 const& geometry1,
{ {
return; return;
} }
BOOST_CHECK_MESSAGE(false, BOOST_CHECK_MESSAGE(false, "A empty_input_exception should have been thrown");
"A empty_input_exception should have been thrown");
try try
{ {
@ -60,8 +57,7 @@ void test_empty_input(Geometry1 const& geometry1,
{ {
return; return;
} }
BOOST_CHECK_MESSAGE(false, BOOST_CHECK_MESSAGE(false, "A empty_input_exception should have been thrown");
"A empty_input_exception should have been thrown");
} }

20
test/algorithms/closest_points/l_l.cpp
View File

@ -85,11 +85,11 @@ void test_closest_points_segment_linestring(Strategies const& strategies)
strategies, true, true); strategies, true, true);
//geometries intersect //geometries intersect
// tester::apply("SEGMENT(0 0,1 1)", tester::apply("SEGMENT(0 0,1 1)",
// "LINESTRING(0 2,0 1,1 0,2 0)", "LINESTRING(0 2,0 1,1 0,2 0)",
// "SEGMENT(0.5 0.5,0.5 0.5)", "SEGMENT(0.5 0.5,0.5 0.5)",
// "SEGMENT(0.5 0.500057,0.5 0.500057)", "SEGMENT(0.5 0.500057,0.5 0.500057)",
// strategies); strategies);
} }
//=========================================================================== //===========================================================================
@ -167,11 +167,11 @@ void test_closest_points_segment_multi_linestring(Strategies const& strategies)
strategies, true, true); strategies, true, true);
//geometries intersect //geometries intersect
// tester::apply("SEGMENT(0 0,1 1)", tester::apply("SEGMENT(0 0,1 1)",
// "MULTILINESTRING((0 2,0 1,1 0,2 0)(1 0,3 0))", "MULTILINESTRING((0 2,0 1,1 0,2 0)(1 0,3 0))",
// "SEGMENT(0.5 0.5,0.5 0.5)", "SEGMENT(0.5 0.5,0.5 0.5)",
// "SEGMENT(0.5 0.500057,0.5 0.500057)", "SEGMENT(0.5 0.500057,0.5 0.500057)",
// strategies); strategies);
} }
//=========================================================================== //===========================================================================