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

buffer: add end round and join miter

Browse Source
This commit is contained in:
Barend Gehrels 2022-09-09 08:59:13 +02:00
parent b18db001b8
commit 475684732c
12 changed files with 612 additions and 107 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
include/boost/geometry/strategies/cartesian/buffer_end_flat.hpp
View File

@ -63,20 +63,13 @@ public :
DistanceStrategy const& distance,
RangeOut& range_out) const
{
typedef typename coordinate_type<Point>::type coordinate_type;
auto const dist_left = distance.apply(penultimate_point, ultimate_point, buffer_side_left);
auto const dist_right = distance.apply(penultimate_point, ultimate_point, buffer_side_right);
typedef typename geometry::select_most_precise
<
coordinate_type,
double
>::type promoted_type;
bool const reversed =
(side == buffer_side_left && dist_right < 0 && -dist_right > dist_left)
|| (side == buffer_side_right && dist_left < 0 && -dist_left > dist_right);
promoted_type const dist_left = distance.apply(penultimate_point, ultimate_point, buffer_side_left);
promoted_type const dist_right = distance.apply(penultimate_point, ultimate_point, buffer_side_right);
bool reversed = (side == buffer_side_left && dist_right < 0 && -dist_right > dist_left)
|| (side == buffer_side_right && dist_left < 0 && -dist_left > dist_right)
;
if (reversed)
{
range_out.push_back(perp_right_point);
@ -89,7 +82,7 @@ public :
}
// Don't add the ultimate_point (endpoint of the linestring).
// The buffer might be generated completely at one side.
// In other cases it does no harm but is further useless
// In other cases it does no harm but it is useless
}
template <typename NumericType>

2
include/boost/geometry/strategies/cartesian/buffer_end_round.hpp
View File

@ -106,7 +106,7 @@ public :
#ifndef DOXYGEN_SHOULD_SKIP_THIS
//! Fills output_range with a flat end
//! Fills output_range with a round end
template <typename Point, typename RangeOut, typename DistanceStrategy>
inline void apply(Point const& penultimate_point,
Point const& perp_left_point,

163
include/boost/geometry/strategies/geographic/buffer_end_round.hpp Normal file
View File

@ -0,0 +1,163 @@
// Boost.Geometry
// Copyright (c) 2022 Barend Gehrels, Amsterdam, the Netherlands.
// Use, modification and distribution is subject to the Boost Software License,
// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
// http://www.boost.org/LICENSE_1_0.txt)
#ifndef BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_BUFFER_END_ROUND_HPP
#define BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_BUFFER_END_ROUND_HPP
#include <boost/range/value_type.hpp>
#include <boost/geometry/core/radian_access.hpp>
#include <boost/geometry/srs/spheroid.hpp>
#include <boost/geometry/strategies/buffer.hpp>
#include <boost/geometry/strategies/geographic/parameters.hpp>
#include <boost/geometry/util/math.hpp>
#include <boost/geometry/util/select_calculation_type.hpp>
namespace boost { namespace geometry
{
namespace strategy { namespace buffer
{
template
<
typename FormulaPolicy = strategy::andoyer,
typename Spheroid = srs::spheroid<double>,
typename CalculationType = void
>
class geographic_end_round
{
static bool const enable_azimuth = true;
static bool const enable_coordinates = true;
template <typename T>
using inverse = typename FormulaPolicy::template inverse
<
T, false, enable_azimuth, false, false, false
>;
template <typename T>
using direct = typename FormulaPolicy::template direct
<
T, enable_coordinates, false, false, false
>;
public :
//! \brief Constructs the strategy
//! \param points_per_circle points which would be used for a full circle
//! (if points_per_circle is smaller than 4, it is internally set to 4)
explicit inline geographic_end_round(std::size_t points_per_circle = 90)
: m_points_per_circle((points_per_circle < 4u) ? 4u : points_per_circle)
{}
#ifndef DOXYGEN_SHOULD_SKIP_THIS
template <typename T, typename RangeOut>
inline void generate(T lon_rad, T lat_rad, T distance, T azimuth, RangeOut& range_out) const
{
using point_t = typename boost::range_value<RangeOut const>::type;
std::size_t const n = m_points_per_circle / 2;
T const angle_diff = geometry::math::pi<T>() / n;
T azi = math::wrap_azimuth_in_radian(azimuth + angle_diff);
// Generate points between 0 and n, not including them
// because left and right are inserted before and after this range.
for (std::size_t i = 1; i < n; i++)
{
auto const d = direct<T>::apply(lon_rad, lat_rad, distance, azi, m_spheroid);
point_t point;
set_from_radian<0>(point, d.lon2);
set_from_radian<1>(point, d.lat2);
range_out.emplace_back(point);
azi = math::wrap_azimuth_in_radian(azi + angle_diff);
}
}
//! Fills output_range with a round end
template <typename Point, typename RangeOut, typename DistanceStrategy>
inline void apply(Point const& penultimate_point, Point const& perp_left_point,
Point const& ultimate_point, Point const& perp_right_point,
buffer_side_selector side, DistanceStrategy const& distance,
RangeOut& range_out) const
{
using calc_t = typename select_calculation_type
<
Point,
typename boost::range_value<RangeOut>::type,
CalculationType
>::type;
calc_t const dist_left = distance.apply(penultimate_point, ultimate_point, buffer_side_left);
calc_t const dist_right = distance.apply(penultimate_point, ultimate_point, buffer_side_right);
bool const reversed = (side == buffer_side_left && dist_right < 0 && -dist_right > dist_left)
|| (side == buffer_side_right && dist_left < 0 && -dist_left > dist_right)
;
calc_t const lon_rad = get_as_radian<0>(ultimate_point);
calc_t const lat_rad = get_as_radian<1>(ultimate_point);
calc_t const lon1_rad = get_as_radian<0>(perp_left_point);
calc_t const lat1_rad = get_as_radian<1>(perp_left_point);
auto const azimuth = inverse<calc_t>::apply(lon_rad, lat_rad, lon1_rad, lat1_rad, m_spheroid).azimuth;
if (reversed)
{
range_out.push_back(perp_right_point);
// generate
range_out.push_back(perp_left_point);
}
else
{
range_out.push_back(perp_left_point);
if (geometry::math::equals(dist_left, dist_right))
{
generate(lon_rad, lat_rad, dist_left, azimuth, range_out);
}
else
{
static calc_t const two = 2.0;
calc_t const dist_average = (dist_left + dist_right) / two;
calc_t const dist_half
= (side == buffer_side_right
? (dist_right - dist_left)
: (dist_left - dist_right)) / two;
auto const shifted = direct<calc_t>::apply(lon_rad, lat_rad, dist_half, azimuth, m_spheroid);
generate(shifted.lon2, shifted.lat2, dist_average, azimuth, range_out);
}
range_out.push_back(perp_right_point);
}
}
template <typename NumericType>
static inline NumericType max_distance(NumericType const& distance)
{
return distance;
}
//! Returns the piece_type (flat end)
static inline piece_type get_piece_type()
{
return buffered_round_end;
}
#endif // DOXYGEN_SHOULD_SKIP_THIS
private :
std::size_t m_points_per_circle;
Spheroid m_spheroid;
};
}} // namespace strategy::buffer
}} // namespace boost::geometry
#endif // BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_BUFFER_END_ROUND_HPP

155
include/boost/geometry/strategies/geographic/buffer_join_miter.hpp Normal file
View File

@ -0,0 +1,155 @@
// Boost.Geometry
// Copyright (c) 2022 Barend Gehrels, Amsterdam, the Netherlands.
// Use, modification and distribution is subject to the Boost Software License,
// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
// http://www.boost.org/LICENSE_1_0.txt)
#ifndef BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_BUFFER_JOIN_MITER_HPP
#define BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_BUFFER_JOIN_MITER_HPP
#include <boost/range/value_type.hpp>
#include <boost/geometry/core/radian_access.hpp>
#include <boost/geometry/srs/spheroid.hpp>
#include <boost/geometry/strategies/buffer.hpp>
#include <boost/geometry/strategies/geographic/parameters.hpp>
#include <boost/geometry/util/math.hpp>
#include <boost/geometry/util/select_calculation_type.hpp>
namespace boost { namespace geometry
{
namespace strategy { namespace buffer
{
template
<
typename FormulaPolicy = strategy::andoyer,
typename Spheroid = srs::spheroid<double>,
typename CalculationType = void
>
class geographic_join_miter
{
static bool const enable_azimuth = true;
static bool const enable_coordinates = true;
template <typename T>
using inverse = typename FormulaPolicy::template inverse
<
T, false, enable_azimuth, false, false, false
>;
template <typename T>
using direct = typename FormulaPolicy::template direct
<
T, enable_coordinates, false, false, false
>;
public :
//! \brief Constructs the strategy
//! \param miter_limit The miter limit, to avoid excessively long miters around sharp corners
explicit inline geographic_join_miter(double miter_limit = 5.0)
: m_miter_limit(valid_limit(miter_limit))
{}
#ifndef DOXYGEN_SHOULD_SKIP_THIS
//! Fills output_range with a sharp shape around a vertex
template <typename Point, typename DistanceType, typename RangeOut>
inline bool apply(Point const& , Point const& vertex,
Point const& perp1, Point const& perp2,
DistanceType const& buffer_distance,
RangeOut& range_out) const
{
using calc_t = typename select_calculation_type
<
Point,
typename boost::range_value<RangeOut>::type,
CalculationType
>::type;
calc_t const lon_rad = get_as_radian<0>(vertex);
calc_t const lat_rad = get_as_radian<1>(vertex);
calc_t const lon1_rad = get_as_radian<0>(perp1);
calc_t const lat1_rad = get_as_radian<1>(perp1);
calc_t const lon2_rad = get_as_radian<0>(perp2);
calc_t const lat2_rad = get_as_radian<1>(perp2);
// Calculate angles from vertex to perp1/perp2
auto const inv1 = inverse<calc_t>::apply(lon_rad, lat_rad, lon1_rad, lat1_rad, m_spheroid);
auto const inv2 = inverse<calc_t>::apply(lon_rad, lat_rad, lon2_rad, lat2_rad, m_spheroid);
// For a sharp corner, perpendicular points are nearly opposite and the
// angle between the two azimuths can be nearly 180, but not more.
calc_t const two_pi = geometry::math::two_pi<calc_t>();
bool const wrapped = inv2.azimuth < inv1.azimuth;
calc_t const angle_diff = wrapped
? ((two_pi + inv2.azimuth) - inv1.azimuth)
: inv2.azimuth - inv1.azimuth;
if (angle_diff < 0 || angle_diff > geometry::math::pi<calc_t>())
{
// Defensive check with asserts
BOOST_GEOMETRY_ASSERT(angle_diff >= 0);
BOOST_GEOMETRY_ASSERT(angle_diff <= geometry::math::pi<calc_t>());
return false;
}
calc_t const half = 0.5;
calc_t const half_angle_diff = half * angle_diff;
calc_t const cos_angle = std::cos(half_angle_diff);
calc_t const max_distance = m_miter_limit * geometry::math::abs(buffer_distance);
if (cos_angle == 0)
{
// It is opposite, perp1==perp2, do not generate a miter cap
return false;
}
// If it is sharp (angle close to 0), the distance will become too high and will be capped.
calc_t const distance = (std::min)(max_distance, buffer_distance / cos_angle);
calc_t const azi = math::wrap_azimuth_in_radian(inv1.azimuth + half_angle_diff);
Point point;
auto const d = direct<calc_t>::apply(lon_rad, lat_rad, distance, azi, m_spheroid);
set_from_radian<0>(point, d.lon2);
set_from_radian<1>(point, d.lat2);
range_out.push_back(perp1);
range_out.push_back(point);
range_out.push_back(perp2);
return true;
}
template <typename NumericType>
inline NumericType max_distance(NumericType const& distance) const
{
return distance * m_miter_limit;
}
#endif // DOXYGEN_SHOULD_SKIP_THIS
private :
double valid_limit(double miter_limit) const
{
if (miter_limit < 1.0)
{
// It should always exceed the buffer distance
miter_limit = 1.0;
}
return miter_limit;
}
double m_miter_limit;
Spheroid m_spheroid;
};
}} // namespace strategy::buffer
}} // namespace boost::geometry
#endif // BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_BUFFER_JOIN_MITER_HPP

72
include/boost/geometry/strategies/geographic/buffer_join_round.hpp
View File

@ -9,8 +9,6 @@
#ifndef BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_BUFFER_JOIN_ROUND_HPP
#define BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_BUFFER_JOIN_ROUND_HPP
#include <cstddef>
#include <boost/range/value_type.hpp>
#include <boost/geometry/core/radian_access.hpp>
@ -36,6 +34,20 @@ template
>
class geographic_join_round
{
static bool const enable_azimuth = true;
static bool const enable_coordinates = true;
template <typename T>
using inverse = typename FormulaPolicy::template inverse
<
T, false, enable_azimuth, false, false, false
>;
template <typename T>
using direct = typename FormulaPolicy::template direct
<
T, enable_coordinates, false, false, false
>;
public :
//! \brief Constructs the strategy
@ -44,17 +56,13 @@ public :
: m_points_per_circle(points_per_circle)
{}
#ifndef DOXYGEN_SHOULD_SKIP_THIS
//! Fills output_range with a rounded shape around a vertex
template
<
typename Point,
typename DistanceType,
typename RangeOut
>
inline bool apply(Point const &ip, Point const &vertex,
Point const &perp1, Point const &perp2,
DistanceType const &buffer_distance,
RangeOut &range_out) const
template <typename Point, typename DistanceType, typename RangeOut>
inline bool apply(Point const& ip, Point const& vertex,
Point const& perp1, Point const& perp2,
DistanceType const& buffer_distance,
RangeOut& range_out) const
{
using calc_t = typename select_calculation_type
<
@ -70,38 +78,25 @@ public :
calc_t const lon2_rad = get_as_radian<0>(perp2);
calc_t const lat2_rad = get_as_radian<1>(perp2);
// Define the types for the Formulas to calculate a point
// at a certain distance (using <direct> with coordinates)
// and to calculate the angle between two specified points
// (using <inverse> with azimuth)
// See also boost/geometry/strategies/geographic/parameters.hpp
constexpr bool enable_azimuth = true;
constexpr bool enable_coordinates = true;
using direct_t = typename FormulaPolicy::template direct
<
calc_t, enable_coordinates, false, false, false
>;
using inverse_t = typename FormulaPolicy::template inverse
<
calc_t, false, enable_azimuth, false, false, false
>;
calc_t const two_pi = geometry::math::two_pi<calc_t>();
// Calculate angles from vertex to perp1/perp2
auto const inv1 = inverse_t::apply(lon_rad, lat_rad, lon1_rad, lat1_rad, m_spheroid);
auto const inv2 = inverse_t::apply(lon_rad, lat_rad, lon2_rad, lat2_rad, m_spheroid);
auto const inv1 = inverse<calc_t>::apply(lon_rad, lat_rad, lon1_rad, lat1_rad, m_spheroid);
auto const inv2 = inverse<calc_t>::apply(lon_rad, lat_rad, lon2_rad, lat2_rad, m_spheroid);
// For a sharp corner, perpendicular points are nearly opposite and the
// angle between the two azimuths can be nearly 180, but not more.
calc_t const two_pi = geometry::math::two_pi<calc_t>();
bool const wrapped = inv2.azimuth < inv1.azimuth;
calc_t const angle_diff = wrapped
? ((two_pi + inv2.azimuth) - inv1.azimuth)
: inv2.azimuth - inv1.azimuth;
BOOST_GEOMETRY_ASSERT(angle_diff >= 0);
BOOST_GEOMETRY_ASSERT(angle_diff <= geometry::math::pi<calc_t>());
if (angle_diff < 0 || angle_diff > geometry::math::pi<calc_t>())
{
// Defensive check with asserts
BOOST_GEOMETRY_ASSERT(angle_diff >= 0);
BOOST_GEOMETRY_ASSERT(angle_diff <= geometry::math::pi<calc_t>());
return false;
}
calc_t const circle_fraction = angle_diff / two_pi;
std::size_t const n = (std::max)(static_cast<std::size_t>(
@ -116,7 +111,7 @@ public :
// because perp1 and perp2 are inserted before and after this range.
for (std::size_t i = 1; i < n; i++)
{
auto const d = direct_t::apply(lon_rad, lat_rad, buffer_distance, azi, m_spheroid);
auto const d = direct<calc_t>::apply(lon_rad, lat_rad, buffer_distance, azi, m_spheroid);
Point p;
set_from_radian<0>(p, d.lon2);
set_from_radian<1>(p, d.lat2);
@ -135,12 +130,13 @@ public :
return distance;
}
#endif // DOXYGEN_SHOULD_SKIP_THIS
private :
std::size_t m_points_per_circle;
Spheroid m_spheroid;
};
}} // namespace strategy::buffer
}} // namespace boost::geometry

2
include/boost/geometry/strategies/strategies.hpp
View File

@ -84,6 +84,8 @@
#include <boost/geometry/strategies/spherical/ssf.hpp>
#include <boost/geometry/strategies/geographic/azimuth.hpp>
#include <boost/geometry/strategies/geographic/buffer_end_round.hpp>
#include <boost/geometry/strategies/geographic/buffer_join_miter.hpp>
#include <boost/geometry/strategies/geographic/buffer_join_round.hpp>
#include <boost/geometry/strategies/geographic/buffer_point_circle.hpp>
#include <boost/geometry/strategies/geographic/buffer_side_straight.hpp>

2
test/algorithms/buffer/Jamfile
View File

@ -21,10 +21,12 @@ test-suite boost-geometry-algorithms-buffer
[ run buffer_point.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_point ]
[ run buffer_point_geo.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_point_geo ]
[ run buffer_linestring.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_linestring ]
[ run buffer_linestring_geo.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_linestring_geo ]
[ run buffer_ring.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_ring ]
[ run buffer_polygon.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_polygon ]
[ run buffer_multi_point.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_multi_point ]
[ run buffer_multi_linestring.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_multi_linestring ]
[ run buffer_multi_linestring_geo.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_multi_linestring_geo ]
[ run buffer_multi_polygon.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_multi_polygon ]
[ run buffer_linestring_aimes.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE : algorithms_buffer_linestring_aimes ]
[ run buffer_linestring.cpp : : : <define>BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE <define>BOOST_GEOMETRY_ROBUSTNESS_ALTERNATIVE : algorithms_buffer_linestring_alternative ]

88
test/algorithms/buffer/buffer_linestring_geo.cpp
View File

@ -10,20 +10,65 @@
#include "test_buffer_geo.hpp"
#include "aimes_cases.hpp"
// Road with one segment from west to east
static std::string const simplex = "LINESTRING(10.3965628 63.4276786,10.3953134 63.4299634)";
// Simple road from west to east and from south to north
// Road from west to east and from south to north
static std::string const road = "LINESTRING(10.3966569 63.4276957,10.3998059 63.4279182,10.4003964 63.4288424)";
template <bool Clockwise, typename PointType>
// Road with first a sharp angle to the right, then a rectangular angle to the left, then two obtuse angles
static std::string const sharp = "LINESTRING(10.3939684 63.4255808, 10.3948159 63.4263661, 10.3947354 63.4255918, 10.3950209 63.4255588, 10.3960570 63.4255423, 10.3973017 63.4255570)";
// Linestring forwards and backwards (close to 180 corner)
static std::string const opposite = "LINESTRING(10.4000988 63.4283885,10.4029694 63.4280726,10.4006988 63.4283397)";
template <typename Formula, bool Clockwise, typename PointType>
void test_linestring()
{
using linestring = bg::model::linestring<PointType>;
using polygon = bg::model::polygon<PointType, Clockwise>;
ut_settings settings(0.1, false, 360);
// Because areas can change significantly when another formula is used,
// use a high tolerance.
int const points_per_circle = 360;
ut_settings settings(0.5);
test_one_geo<linestring, polygon>("simplex_5_8", simplex, 2622.1129, 5.0, settings);
test_one_geo<linestring, polygon>("road_5_8", road, 1972.6796, 5.0, settings, 2.5);
#if defined(BOOST_GEOMETRY_TEST_FAILURES)
bool const thomas_skip = false;
#else
// TODO: some cases are missing one or more turns
bool const thomas_skip = std::is_same<Formula, bg::strategy::thomas>::value;
#endif
bg::strategies::buffer::geographic<Formula> strategy;
bg::strategy::buffer::geographic_side_straight<Formula> side;
bg::strategy::buffer::geographic_join_miter<Formula> join_miter;
bg::strategy::buffer::geographic_join_miter<Formula> join_miter25(2.5);
bg::strategy::buffer::geographic_join_round<Formula> join_round(points_per_circle);
bg::strategy::buffer::geographic_end_round<Formula> end_round(points_per_circle);
bg::strategy::buffer::geographic_point_circle<Formula> circle(points_per_circle);
bg::strategy::buffer::end_flat end_flat;
test_one_geo<linestring, polygon>("simplex_5_8", simplex, strategy, side, circle, join_round, end_flat, 2622.0, 5.0, settings);
test_one_geo<linestring, polygon>("road_5_flat", road, strategy, side, circle, join_round, end_flat, 2644.0, 5.0, settings);
test_one_geo<linestring, polygon>("road_5_25_round", road, strategy, side, circle, join_round, end_round, 2016.0, 5.0, settings, 2.5);
test_one_geo<linestring, polygon>("sharp_5_round", sharp, strategy, side, circle, join_round, end_round, 3090.0, 5.0, settings);
test_one_geo<linestring, polygon>("sharp_5_miter", sharp, strategy, side, circle, join_miter, end_round, 3181.0, 5.0, settings);
test_one_geo<linestring, polygon>("sharp_5_miter25", sharp, strategy, side, circle, join_miter25, end_round, 3121.0, 5.0, settings);
if (! BOOST_GEOMETRY_CONDITION(thomas_skip))
{
// Misses an intersection point when using thomas
test_one_geo<linestring, polygon>("opposite", opposite, strategy, side, circle, join_round, end_round, 1658.0, 5.0, settings);
}
{
auto specific = settings;
specific.fraction_buffered_points_too_close = 0.3;
test_one_geo<linestring, polygon>("opposite_miter", opposite, strategy, side, circle, join_miter, end_flat, 1705.0, 5.0, specific);
test_one_geo<linestring, polygon>("opposite_miter25", opposite, strategy, side, circle, join_miter25, end_flat, 1642.0, 5.0, specific);
}
settings.test_area = false;
auto const n = sizeof(testcases_aimes) / sizeof(testcases_aimes[0]);
@ -40,6 +85,8 @@ void test_linestring()
std::set<int> const curved_cases_max_area{5, 95, 119, 142, 192};
// Cases which are rounded such that it results in a large area
std::set<int> const round_cases_max_area{196};
// Cases which are not yet valid or false negatives
std::set<int> const round_cases_invalid{143};
for (auto i = 0; i < n; i++)
{
@ -53,10 +100,32 @@ void test_linestring()
: ut_settings().multiplier_max_area;
settings.fraction_buffered_points_too_close
= cases_with_artefacts.count(i) > 0 ? 0.15
= cases_with_artefacts.count(i) > 0 ? 0.20
: ut_settings().fraction_buffered_points_too_close;
test_one_geo<linestring, polygon>("aimes_" + std::to_string(i), testcases_aimes[i], -1, 25.0, settings);
// With miter/flat, the artefacts are more pronounced and there can
// be more points close than expected.
auto settings_mf = settings;
settings_mf.fraction_buffered_points_too_close *= 4.0;
// With rounded cases, both rounded ends overap, adapt the expected area
auto settings_rr = settings;
settings_rr.multiplier_min_area
= round_cases_max_area.count(i) > 0 ? 0.75
: settings.multiplier_min_area;
settings_rr.set_test_validity(round_cases_invalid.count(i) == 0);
if (i != 181)
{
// 181 fails, it should generate a hole, but instead that is the outer ring now.
test_one_geo<linestring, polygon>("aimes_" + std::to_string(i) + "_rr", testcases_aimes[i],
strategy, side, circle, join_round, end_round, -1, 25.0, settings_rr);
}
test_one_geo<linestring, polygon>("aimes_" + std::to_string(i) + "_rf", testcases_aimes[i],
strategy, side, circle, join_round, end_flat, -1, 25.0, settings);
test_one_geo<linestring, polygon>("aimes_" + std::to_string(i) + "_mf", testcases_aimes[i],
strategy, side, circle, join_miter, end_flat, -1, 25.0, settings_mf);
}
}
@ -64,10 +133,11 @@ int test_main(int, char* [])
{
BoostGeometryWriteTestConfiguration();
test_linestring<true, bg::model::point<default_test_type, 2, bg::cs::geographic<bg::degree> > >();
test_linestring<bg::strategy::andoyer, true, bg::model::point<default_test_type, 2, bg::cs::geographic<bg::degree> > >();
test_linestring<bg::strategy::thomas, true, bg::model::point<default_test_type, 2, bg::cs::geographic<bg::degree> > >();
#if ! defined(BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE)
test_linestring<true, bg::model::point<long double, 2, bg::cs::geographic<bg::degree> > >();
test_linestring<bg::strategy::andoyer, true, bg::model::point<long double, 2, bg::cs::geographic<bg::degree> > >();
#endif
return 0;

86
test/algorithms/buffer/buffer_multi_linestring_geo.cpp Normal file
View File

@ -0,0 +1,86 @@
// Boost.Geometry
// Unit Test
// Copyright (c) 2022 Barend Gehrels, Amsterdam, the Netherlands.
// Use, modification and distribution is subject to the Boost Software License,
// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
// http://www.boost.org/LICENSE_1_0.txt)
#include "test_buffer_geo.hpp"
// Case constructed such that an interior ring will be formed, depending on the buffer distance (11..18)
static std::string const trondheim = "MULTILINESTRING((10.3988246 63.4332944, 10.3991235 63.4341424, 10.3991956 63.4346470), (10.3985108 63.4336039, 10.3985342 63.4339071, 10.3984854 63.4342548, 10.3984154 63.4345156), (10.3973490 63.4332584, 10.4001666 63.4332774), (10.3985596 63.4346152, 10.3991129 63.4345792))";
static std::string const wrangel = "MULTILINESTRING((178.6083 71.0669, 179.9387 71.5239, -178.3429 71.5283, -177.4560 71.2289))";
template <typename Formula, bool Clockwise, typename PointType>
void test_geometry()
{
using ml = bg::model::multi_linestring<bg::model::linestring<PointType>>;
using polygon = bg::model::polygon<PointType, Clockwise>;
int const points_per_circle = 360;
ut_settings settings(0.5);
settings.fraction_buffered_points_too_close = 0.35;
bg::strategies::buffer::geographic<Formula> strategy;
bg::strategy::buffer::geographic_side_straight<Formula> side;
bg::strategy::buffer::geographic_join_miter<Formula> join_miter;
bg::strategy::buffer::geographic_join_round<Formula> join_round(points_per_circle);
bg::strategy::buffer::geographic_end_round<Formula> end_round(points_per_circle);
bg::strategy::buffer::geographic_point_circle<Formula> circle(points_per_circle);
bg::strategy::buffer::end_flat end_flat;
#if defined(BOOST_GEOMETRY_TEST_FAILURES)
bool const andoyer_skip = false;
bool const thomas_skip = false;
#else
// TODO: some cases are missing one or more turns.
bool const andoyer_skip = std::is_same<Formula, bg::strategy::andoyer>::value;
bool const thomas_skip = std::is_same<Formula, bg::strategy::thomas>::value;
#endif
test_one_geo<ml, polygon>("trondheim05_rr", trondheim, strategy, side, circle, join_round, end_round, 4398.0, 5.0, settings);
if (! BOOST_GEOMETRY_CONDITION(thomas_skip))
{
test_one_geo<ml, polygon>("trondheim10_rr", trondheim, strategy, side, circle, join_round, end_round, 8994.0, 10.0, settings);
}
test_one_geo<ml, polygon>("trondheim12_rr", trondheim, strategy, side, circle, join_round, end_round, 10790.0, 12.0, settings);
if (! BOOST_GEOMETRY_CONDITION(thomas_skip) && ! BOOST_GEOMETRY_CONDITION(andoyer_skip))
{
test_one_geo<ml, polygon>("trondheim15_rr", trondheim, strategy, side, circle, join_round, end_round, 13358.0, 15.0, settings);
}
if (! BOOST_GEOMETRY_CONDITION(thomas_skip))
{
test_one_geo<ml, polygon>("trondheim17_rr", trondheim, strategy, side, circle, join_round, end_round, 14824.0, 17.0, settings);
}
test_one_geo<ml, polygon>("trondheim20_rr", trondheim, strategy, side, circle, join_round, end_round, 17055.0, 20.0, settings);
test_one_geo<ml, polygon>("trondheim25_rr", trondheim, strategy, side, circle, join_round, end_round, 20657.0, 25.0, settings);
test_one_geo<ml, polygon>("trondheim05_mf", trondheim, strategy, side, circle, join_miter, end_flat, 4190.0, 5.0, settings);
test_one_geo<ml, polygon>("trondheim10_mf", trondheim, strategy, side, circle, join_miter, end_flat, 8196.0, 10.0, settings);
test_one_geo<ml, polygon>("trondheim12_mf", trondheim, strategy, side, circle, join_miter, end_flat, 9706.0, 12.0, settings);
test_one_geo<ml, polygon>("trondheim15_mf", trondheim, strategy, side, circle, join_miter, end_flat, 11708.0, 15.0, settings);
test_one_geo<ml, polygon>("trondheim17_mf", trondheim, strategy, side, circle, join_miter, end_flat, 12896.0, 17.0, settings);
test_one_geo<ml, polygon>("trondheim20_mf", trondheim, strategy, side, circle, join_miter, end_flat, 14486.0, 20.0, settings);
test_one_geo<ml, polygon>("trondheim25_mf", trondheim, strategy, side, circle, join_miter, end_flat, 17025.0, 25.0, settings);
}
int test_main(int, char* [])
{
BoostGeometryWriteTestConfiguration();
test_geometry<bg::strategy::andoyer, true, bg::model::point<default_test_type, 2, bg::cs::geographic<bg::degree> > >();
test_geometry<bg::strategy::thomas, true, bg::model::point<default_test_type, 2, bg::cs::geographic<bg::degree> > >();
#if ! defined(BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE)
test_geometry<true, bg::model::point<long double, 2, bg::cs::geographic<bg::degree> > >();
#endif
return 0;
}

46
test/algorithms/buffer/buffer_point_geo.cpp
View File

@ -11,33 +11,53 @@
static std::string const simplex = "POINT(4.9 52.0)";
template <bool Clockwise, typename PointType>
template <typename Formula, bool Clockwise, typename PointType>
void test_point()
{
using polygon = bg::model::polygon<PointType, Clockwise>;
// NOTE: for now do not test with a radius less than 2 meter, because is not precise yet (in double)
test_one_geo<PointType, polygon>("simplex_5_8", simplex, 70.7107, 5.0, ut_settings(0.1, false, 8));
test_one_geo<PointType, polygon>("simplex_5_16", simplex, 76.5437, 5.0, ut_settings(0.1, false, 16));
// * Result is different for clang/VCC. Specified expectation is in between, and tolerance higher
test_one_geo<PointType, polygon>("simplex_5_32", simplex, 77.9640, 5.0, ut_settings(0.2, false, 32));
// Test with a high tolerance to account for possible differences in andoyer/thomas.
auto make_settings = [](int points_per_circle)
{
ut_settings result;
result.tolerance = 0.5;
result.points_per_circle = points_per_circle;
return result;
};
// Test with a high tolerance to account for possible differences in andoyer/thomas.
auto make_circle = [](int points_per_circle)
{
return bg::strategy::buffer::geographic_point_circle<Formula>(points_per_circle);
};
bg::strategies::buffer::geographic<Formula> strategy;
bg::strategy::buffer::geographic_join_miter<Formula> join;
bg::strategy::buffer::geographic_side_straight<Formula> side;
bg::strategy::buffer::end_flat end;
test_one_geo<PointType, polygon>("simplex_1_16", simplex, strategy, side, make_circle(360), join, end, 3.1415, 1.0, make_settings(360));
test_one_geo<PointType, polygon>("simplex_5_8", simplex, strategy, side, make_circle(8), join, end, 70.7107, 5.0, make_settings(8));
test_one_geo<PointType, polygon>("simplex_5_16", simplex, strategy, side, make_circle(16), join, end, 76.5437, 5.0, make_settings(16));
test_one_geo<PointType, polygon>("simplex_5_32", simplex, strategy, side, make_circle(32), join, end, 77.9640, 5.0, make_settings(32));
// The more points used for the buffer, the more the area approaches 10*PI square meters
test_one_geo<PointType, polygon>("simplex_10_8", simplex, 282.8430, 10.0, ut_settings(0.1, false, 8));
test_one_geo<PointType, polygon>("simplex_10_16", simplex, 306.1471, 10.0, ut_settings(0.1, false, 16));
test_one_geo<PointType, polygon>("simplex_10_32", simplex, 312.1450, 10.0, ut_settings(0.1, false, 32));
// * Same here
test_one_geo<PointType, polygon>("simplex_10_180", simplex, 313.9051, 10.0, ut_settings(0.2, false, 180));
test_one_geo<PointType, polygon>("simplex_10_8", simplex, strategy, side, make_circle(8), join, end, 282.8430, 10.0, make_settings(8));
test_one_geo<PointType, polygon>("simplex_10_16", simplex, strategy, side, make_circle(16), join, end, 306.1471, 10.0, make_settings(16));
test_one_geo<PointType, polygon>("simplex_10_32", simplex, strategy, side, make_circle(32), join, end, 312.1450, 10.0, make_settings(32));
test_one_geo<PointType, polygon>("simplex_10_180", simplex, strategy, side, make_circle(180), join, end, 313.9051, 10.0, make_settings(180));
test_one_geo<PointType, polygon>("simplex_10_180", simplex, strategy, side, make_circle(360), join, end, 314.15, 10.0, make_settings(360));
}
int test_main(int, char* [])
{
BoostGeometryWriteTestConfiguration();
test_point<true, bg::model::point<default_test_type, 2, bg::cs::geographic<bg::degree> > >();
test_point<bg::strategy::andoyer, true, bg::model::point<default_test_type, 2, bg::cs::geographic<bg::degree> > >();
test_point<bg::strategy::thomas, true, bg::model::point<default_test_type, 2, bg::cs::geographic<bg::degree> > >();
#if ! defined(BOOST_GEOMETRY_TEST_ONLY_ONE_TYPE)
test_point<true, bg::model::point<long double, 2, bg::cs::geographic<bg::degree> > >();
test_point<bg::strategy::andoyer, true, bg::model::point<long double, 2, bg::cs::geographic<bg::degree> > >();
#endif
return 0;

23
test/algorithms/buffer/test_buffer.hpp
View File

@ -67,27 +67,38 @@ template <typename JoinStrategy>
struct JoinTestProperties
{
static std::string name() { return "join_unknown"; }
static bool is_miter() { return false; }
};
template<> struct JoinTestProperties<boost::geometry::strategy::buffer::join_round>
{
static std::string name() { return "round"; }
static bool is_miter() { return false; }
};
template<typename F, typename S, typename CT>
struct JoinTestProperties<boost::geometry::strategy::buffer::geographic_join_round<F, S, CT> >
{
static std::string name() { return "geo_round"; }
static bool is_miter() { return false; }
};
template<> struct JoinTestProperties<boost::geometry::strategy::buffer::join_miter>
{
static std::string name() { return "miter"; }
static bool is_miter() { return true; }
};
template<typename F, typename S, typename CT>
struct JoinTestProperties<boost::geometry::strategy::buffer::geographic_join_miter<F, S, CT> >
{
static std::string name() { return "geo_miter"; }
static bool is_miter() { return true; }
};
template<> struct JoinTestProperties<boost::geometry::strategy::buffer::join_round_by_divide>
{
static std::string name() { return "divide"; }
static bool is_miter() { return false; }
};
@ -98,11 +109,20 @@ struct EndTestProperties { };
template<> struct EndTestProperties<boost::geometry::strategy::buffer::end_round>
{
static std::string name() { return "round"; }
static bool is_round() { return true; }
};
template<typename F, typename S, typename CT>
struct EndTestProperties<boost::geometry::strategy::buffer::geographic_end_round<F, S, CT>>
{
static std::string name() { return "geo_round"; }
static bool is_round() { return true; }
};
template<> struct EndTestProperties<boost::geometry::strategy::buffer::end_flat>
{
static std::string name() { return "flat"; }
static bool is_round() { return false; }
};
struct ut_settings : public ut_base_settings
@ -133,7 +153,10 @@ struct ut_settings : public ut_base_settings
double tolerance;
bool test_area = true;
bool use_ln_area = false;
// Number of points in a circle. Not used for geo tests.
int points_per_circle;
double multiplier_min_area = 0.95;
double multiplier_max_area = 1.05;
double fraction_buffered_points_too_close = 0.10;

61
test/algorithms/buffer/test_buffer_geo.hpp
View File

@ -17,10 +17,21 @@
#include "test_buffer.hpp"
template<typename Geometry, typename GeometryOut>
template<typename Geometry,
typename GeometryOut,
typename UmbrellaStrategy,
typename SideStrategy,
typename PointStrategy,
typename JoinStrategy,
typename EndStrategy>
void test_one_geo(std::string const& caseid,
std::string const& wkt,
int expected_count, int expected_holes_count, double expected_area,
UmbrellaStrategy const& strategy,
SideStrategy const& side_strategy,
PointStrategy const& circle_strategy,
JoinStrategy const& join_strategy,
EndStrategy const& end_strategy,
double expected_area,
double distance_left, ut_settings settings = ut_settings(),
double distance_right = same_distance)
{
@ -28,31 +39,18 @@ void test_one_geo(std::string const& caseid,
bg::read_wkt(wkt, input_geometry);
bg::correct(input_geometry);
const bool symmetric = bg::math::equals(distance_right, same_distance);
bool const symmetric = bg::math::equals(distance_right, same_distance);
if (symmetric)
{
distance_right = distance_left;
}
auto const mean_distance = (distance_left + distance_right) / 2.0;
bg::strategy::buffer::distance_asymmetric
<
typename bg::coordinate_type<Geometry>::type
> distance_strategy(distance_left, distance_right);
// Use the appropriate strategies for geographic points
bg::strategy::buffer::geographic_side_straight<> side_strategy;
bg::strategy::buffer::geographic_point_circle<> circle_strategy(settings.points_per_circle);
bg::strategy::buffer::geographic_join_round<> join_strategy(settings.points_per_circle);
// TODO
bg::strategy::buffer::end_flat end_strategy;
bg::strategies::buffer::geographic
<
bg::strategy::thomas, bg::srs::spheroid<long double>, long double
> strategy;
bg::model::multi_polygon<GeometryOut> buffer;
test_buffer<GeometryOut>
@ -60,7 +58,7 @@ void test_one_geo(std::string const& caseid,
join_strategy, end_strategy,
distance_strategy, side_strategy, circle_strategy,
strategy,
expected_count, expected_holes_count, expected_area,
-1, -1, expected_area,
settings);
if (symmetric && distance_left > 0.0)
@ -84,18 +82,26 @@ void test_one_geo(std::string const& caseid,
const double f = too_close / static_cast<double>(total);
BOOST_CHECK_MESSAGE(f < settings.fraction_buffered_points_too_close,
caseid << " has too many points too close " << too_close << " " << f);
BOOST_CHECK_MESSAGE(too_far == 0,
caseid << " has too far " << too_far);
if (!JoinTestProperties<JoinStrategy>::is_miter())
{
BOOST_CHECK_MESSAGE(too_far == 0,
caseid << " has too far " << too_far);
}
}
if (expected_area < 0 && bg::util::is_linear<Geometry>::value)
{
// Calculate the area of a linear feature using its length and the buffer distance.
// For a straight line, with flat ends, this expectation is perfect.
// For round ends, add the area of a circle (two halves at both ends).
// For a straight line this expectation is perfect.
// For a curved line it might be too large.
// Therefore the default is 95% of it, and it can be modified with a setting.
const auto addition = EndTestProperties<EndStrategy>::is_round()
? mean_distance * mean_distance * bg::math::pi<double>() : 0.0;
const auto area = bg::area(buffer);
const auto estimated_area = bg::length(input_geometry) * (distance_left + distance_right);
const auto estimated_area = addition + bg::length(input_geometry) * (distance_left + distance_right);
const auto min_area = settings.multiplier_min_area * estimated_area;
const auto max_area = settings.multiplier_max_area * estimated_area;
BOOST_CHECK_MESSAGE(area > min_area,
@ -107,15 +113,4 @@ void test_one_geo(std::string const& caseid,
}
}
template<typename Geometry, typename GeometryOut>
void test_one_geo(std::string const& caseid, std::string const& wkt,
double expected_area,
double distance_left, ut_settings const& settings = ut_settings(),
double distance_right = same_distance)
{
test_one_geo<Geometry, GeometryOut>(caseid, wkt,
-1 ,-1, expected_area,
distance_left, settings, distance_right);
}
#endif