From bd2dc54e8c212f00a7c9fbede5563810ea04ef8a Mon Sep 17 00:00:00 2001 From: Adam Wulkiewicz Date: Tue, 12 Feb 2019 21:48:44 +0100 Subject: [PATCH] [extensions][io] Check rings point-order and assign inner rings. --- .../extensions/gis/io/shapefile/read.hpp | 227 ++++++++++++++---- 1 file changed, 174 insertions(+), 53 deletions(-) diff --git a/include/boost/geometry/extensions/gis/io/shapefile/read.hpp b/include/boost/geometry/extensions/gis/io/shapefile/read.hpp index cf9cfcb01..aa5c1d716 100644 --- a/include/boost/geometry/extensions/gis/io/shapefile/read.hpp +++ b/include/boost/geometry/extensions/gis/io/shapefile/read.hpp @@ -9,6 +9,7 @@ #ifndef BOOST_GEOMETRY_EXTENSIONS_GIS_IO_SHAPEFILE_READ_HPP #define BOOST_GEOMETRY_EXTENSIONS_GIS_IO_SHAPEFILE_READ_HPP + #include #include @@ -16,6 +17,8 @@ #include #include +#include +#include #include #include #include @@ -27,9 +30,16 @@ #include #include +// TEMP - only here for convenience, for now +#include +#include +#include + + namespace boost { namespace geometry { + class read_shapefile_exception : public geometry::exception { public: @@ -145,14 +155,6 @@ struct shape_type }; }; -template -inline typename boost::range_reference::type -push_back(Range & rng, typename boost::range_value::type const& v) -{ - range::push_back(rng, v); - return range::back(rng); -} - template inline void read_and_push_back_point(IStream & is, Range & rng) { @@ -197,8 +199,8 @@ inline void read_parts(IStream & is, struct read_point_policy { - template - static inline void apply(IStream & is, Points & points) + template + static inline void apply(IStream & is, Points & points, Strategy const&) { boost::int32_t type; read_little(is, type); @@ -218,8 +220,8 @@ struct read_point_policy struct read_multipoint_policy { - template - static inline void apply(IStream & is, Points & points) + template + static inline void apply(IStream & is, Points & points, Strategy const&) { typedef typename boost::range_value::type pt_type; @@ -256,8 +258,8 @@ struct read_multipoint_policy struct read_polyline_policy { - template - static inline void apply(IStream &is, Linestrings & linestrings) + template + static inline void apply(IStream &is, Linestrings & linestrings, Strategy const&) { typedef typename boost::range_value::type ls_type; typedef typename boost::range_value::type pt_type; @@ -300,7 +302,8 @@ struct read_polyline_policy BOOST_THROW_EXCEPTION(read_shapefile_exception("Invalid part number")); } - ls_type & ls = push_back(linestrings, ls_type()); + range::push_back(linestrings, ls_type()); + ls_type & ls = range::back(linestrings); std::size_t ls_size = l - f; @@ -321,8 +324,8 @@ struct read_polyline_policy struct read_polygon_policy { - template - static inline void apply(IStream &is, Polygons & polygons) + template + static inline void apply(IStream &is, Polygons & polygons, Strategy const& strategy) { typedef typename boost::range_value::type poly_type; typedef typename geometry::point_type::type pt_type; @@ -331,6 +334,11 @@ struct read_polygon_policy static const bool is_ccw = geometry::point_order::value == geometry::counterclockwise; static const bool is_open = geometry::closure::value == geometry::open; + typename Strategy::point_order_strategy_type + order_strategy = strategy.get_point_order_strategy(); + typename Strategy::template point_in_geometry_strategy::type + within_strategy = strategy.template get_point_in_geometry_strategy(); + boost::int32_t type; //double min_x, min_y, max_x, max_y; boost::int32_t num_parts; @@ -359,7 +367,8 @@ struct read_polygon_policy read_parts(is, parts, num_parts); - poly_type & poly = push_back(polygons, poly_type()); + std::vector outer_rings; + std::vector inner_rings; for (boost::int32_t i = 0; i < num_parts; ++i) { @@ -371,9 +380,7 @@ struct read_polygon_policy BOOST_THROW_EXCEPTION(read_shapefile_exception("Invalid part number")); } - ring_type & ring = (i == 0) - ? geometry::exterior_ring(poly) - : push_back(geometry::interior_rings(poly), ring_type()); + ring_type ring; std::size_t ring_size = l - f - (is_open ? 1 : 0); @@ -398,26 +405,119 @@ struct read_polygon_policy e = boost::end(ring); std::reverse(++b, is_open ? e : (--e)); } + + // TODO: support rval references in range::push_back() + // and/or implement range::emplace_back() + // implement and call move_or_copy(ring) + + // assume outer ring + if (num_parts == 1) + range::push_back(outer_rings, ring); // order could be checked here too + // check order + else if ( is_outer_ring(ring, order_strategy) ) + range::push_back(outer_rings, ring); + else + range::push_back(inner_rings, ring); } + if (inner_rings.empty()) // no inner rings + { + for (size_t i = 0; i < outer_rings.size(); ++i) + { + poly_type poly; + geometry::exterior_ring(poly) = outer_rings[i]; // TODO: move + range::push_back(polygons, poly); // TODO: move + } + } + else if (! outer_rings.empty()) // outer and inner rings + { + std::vector assigned(inner_rings.size(), false); + std::size_t assigned_count = 0; + for (size_t i = 0; i < outer_rings.size(); ++i) + { + poly_type poly; + geometry::exterior_ring(poly) = outer_rings[i]; // TODO: move + for (size_t j = 0; j < inner_rings.size(); ++j) + { + if (! assigned[j]) + { + if (is_inner_ring(inner_rings[j], + geometry::exterior_ring(poly), + within_strategy)) + { + range::push_back(geometry::interior_rings(poly), inner_rings[j]); // TODO: move + ++assigned_count; + assigned[j] = true; + } + // just in case, ignore empty + else if (boost::empty(inner_rings[j])) + { + ++assigned_count; + assigned[j] = true; + } + } + } + range::push_back(polygons, poly); // TODO: move + } + + // check if all interior rings were assigned + if (assigned_count != inner_rings.size()) + { + BOOST_THROW_EXCEPTION(read_shapefile_exception("Not all interior rings were assigned to polygons.")); + } + } + else // inner rings but no outer rings + { + // unexpected, file corrupted, bug or numerical error + BOOST_THROW_EXCEPTION(read_shapefile_exception("Exterior ring expected")); + } + if (!is.good()) { BOOST_THROW_EXCEPTION(read_shapefile_exception("Read error")); } } + + template + static inline bool is_outer_ring(Ring const& ring, Strategy const& order_strategy) + { + geometry::order_selector result = detail::calculate_point_order(ring, order_strategy); + return result == geometry::point_order::value + || result == geometry::order_undetermined; + } + + template + static inline bool is_inner_ring(InnerRing const& inner_ring, + OuterRing const& outer_ring, + Strategy const& within_strategy) + { + // NOTE: The worst case complexity is O(N^2) and best O(N) + // R-tree or partition could be used (~O(NlogN)) + // but in most cases this version should be faster. + typedef typename boost::range_iterator::type iter_type; + for (iter_type it = boost::begin(inner_ring); it != boost::end(inner_ring); ++it) + { + if (detail::within::within_point_geometry(*it, outer_ring, within_strategy)) + { + // at least one point of potentially interior ring found in the interior of exterior ring + return true; + } + } + return false; + } }; -template -inline void add_records(IStream & is, Range & rng) +template +inline void add_records(IStream & is, Range & rng, Strategy const& strategy) { while (read_record_header(is)) { - Policy::apply(is, rng); + Policy::apply(is, rng, strategy); } } -template -inline void add_records_as_new_element(IStream & is, Range & rng) +template +inline void add_records_as_new_element(IStream & is, Range & rng, Strategy const& strategy) { typedef typename boost::range_value::type val_type; @@ -426,25 +526,27 @@ inline void add_records_as_new_element(IStream & is, Range & rng) return; } - val_type & elem = shapefile::push_back(rng, val_type()); + range::push_back(rng, val_type()); + val_type & elem = range::back(rng); do { - Policy::apply(is, elem); + Policy::apply(is, elem, strategy); } while (read_record_header(is)); } -template -inline void add_records_as_new_elements(IStream & is, Range & rng) +template +inline void add_records_as_new_elements(IStream & is, Range & rng, Strategy const& strategy) { typedef typename boost::range_value::type val_type; while (read_record_header(is)) { - val_type & elem = shapefile::push_back(rng, val_type()); + range::push_back(rng, val_type()); + val_type & elem = range::back(rng); - Policy::apply(is, elem); + Policy::apply(is, elem, strategy); } } @@ -465,8 +567,8 @@ struct read_shapefile template struct read_shapefile { - template - static inline void apply(IStream &is, Points & points) + template + static inline void apply(IStream &is, Points & points, Strategy const& strategy) { namespace shp = detail::shapefile; @@ -474,11 +576,11 @@ struct read_shapefile if (type == shp::shape_type::point) { - shp::add_records(is, points); + shp::add_records(is, points, strategy); } else if (type == detail::shapefile::shape_type::multipoint) { - shp::add_records(is, points); + shp::add_records(is, points, strategy); } } }; @@ -486,8 +588,8 @@ struct read_shapefile template struct read_shapefile { - template - static inline void apply(IStream &is, MultiPoints & multi_points) + template + static inline void apply(IStream &is, MultiPoints & multi_points, Strategy const& strategy) { namespace shp = detail::shapefile; @@ -495,11 +597,11 @@ struct read_shapefile if (type == shp::shape_type::point) { - shp::add_records_as_new_element(is, multi_points); + shp::add_records_as_new_element(is, multi_points, strategy); } else if (type == detail::shapefile::shape_type::multipoint) { - shp::add_records_as_new_elements(is, multi_points); + shp::add_records_as_new_elements(is, multi_points, strategy); } } }; @@ -507,8 +609,8 @@ struct read_shapefile template struct read_shapefile { - template - static inline void apply(IStream &is, Linestrings & linestrings) + template + static inline void apply(IStream &is, Linestrings & linestrings, Strategy const& strategy) { namespace shp = detail::shapefile; @@ -516,7 +618,7 @@ struct read_shapefile if (type == shp::shape_type::polyline) { - shp::add_records(is, linestrings); + shp::add_records(is, linestrings, strategy); } } }; @@ -524,8 +626,8 @@ struct read_shapefile template struct read_shapefile { - template - static inline void apply(IStream &is, MultiLinestrings & multi_linestrings) + template + static inline void apply(IStream &is, MultiLinestrings & multi_linestrings, Strategy const& strategy) { namespace shp = detail::shapefile; @@ -533,7 +635,7 @@ struct read_shapefile if (type == shp::shape_type::polyline) { - shp::add_records_as_new_elements(is, multi_linestrings); + shp::add_records_as_new_elements(is, multi_linestrings, strategy); } } }; @@ -541,8 +643,8 @@ struct read_shapefile template struct read_shapefile { - template - static inline void apply(IStream &is, Polygons & polygons) + template + static inline void apply(IStream &is, Polygons & polygons, Strategy const& strategy) { namespace shp = detail::shapefile; @@ -550,7 +652,7 @@ struct read_shapefile if (type == shp::shape_type::polygon) { - shp::add_records(is, polygons); + shp::add_records(is, polygons, strategy); } } }; @@ -558,8 +660,8 @@ struct read_shapefile template struct read_shapefile { - template - static inline void apply(IStream &is, MultiPolygons & multi_polygons) + template + static inline void apply(IStream &is, MultiPolygons & multi_polygons, Strategy const& strategy) { namespace shp = detail::shapefile; @@ -567,7 +669,7 @@ struct read_shapefile if (type == shp::shape_type::polygon) { - shp::add_records_as_new_elements(is, multi_polygons); + shp::add_records_as_new_elements(is, multi_polygons, strategy); } } }; @@ -577,12 +679,31 @@ struct read_shapefile // Note: if an exception is thrown the output range may contain partial data +template +inline void read_shapefile(IStream &is, RangeOfGeometries & range_of_geometries, + Strategy const& strategy) +{ + typedef typename boost::range_value::type geometry_type; + + geometry::concepts::check(); + + dispatch::read_shapefile::apply(is, range_of_geometries, + strategy); +} + template inline void read_shapefile(IStream &is, RangeOfGeometries & range_of_geometries) { typedef typename boost::range_value::type geometry_type; + typedef typename strategy::io::services::default_strategy + < + typename cs_tag::type + >::type strategy_type; + geometry::concepts::check(); - dispatch::read_shapefile::apply(is, range_of_geometries); + + dispatch::read_shapefile::apply(is, range_of_geometries, + strategy_type()); }