detail/disjoint/segment_box.hpp

0001 // Boost.Geometry (aka GGL, Generic Geometry Library)
0002
0003 // Copyright (c) 2007-2014 Barend Gehrels, Amsterdam, the Netherlands.
0004 // Copyright (c) 2008-2014 Bruno Lalande, Paris, France.
0005 // Copyright (c) 2009-2014 Mateusz Loskot, London, UK.
0006 // Copyright (c) 2013-2014 Adam Wulkiewicz, Lodz, Poland.
0007
0008 // This file was modified by Oracle on 2013-2021.
0009 // Modifications copyright (c) 2013-2021, Oracle and/or its affiliates.
0010
0011 // Contributed and/or modified by Vissarion Fysikopoulos, on behalf of Oracle
0012 // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
0013 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
0014
0015 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
0016 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
0017
0018 // Use, modification and distribution is subject to the Boost Software License,
0019 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
0020 // http://www.boost.org/LICENSE_1_0.txt)
0021
0022 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP
0023 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP
0024
0025 #include <cstddef>
0026
0027 #include <boost/geometry/core/tags.hpp>
0028 #include <boost/geometry/core/radian_access.hpp>
0029
0030 #include <boost/geometry/algorithms/detail/assign_indexed_point.hpp>
0031 #include <boost/geometry/algorithms/detail/disjoint/point_box.hpp>
0032 #include <boost/geometry/algorithms/detail/disjoint/box_box.hpp>
0033 #include <boost/geometry/algorithms/detail/envelope/segment.hpp>
0034 #include <boost/geometry/algorithms/detail/normalize.hpp>
0035 #include <boost/geometry/algorithms/dispatch/disjoint.hpp>
0036
0037 #include <boost/geometry/formulas/vertex_longitude.hpp>
0038
0039 #include <boost/geometry/geometries/box.hpp>
0040
0041 // Temporary, for envelope_segment_impl
0042 #include <boost/geometry/strategy/spherical/envelope_segment.hpp>
0043
0044 namespace boost { namespace geometry
0045 {
0046
0047
0048 #ifndef DOXYGEN_NO_DETAIL
0049 namespace detail { namespace disjoint
0050 {
0051
0052 template <typename CS_Tag>
0053 struct disjoint_segment_box_sphere_or_spheroid
0054 {
0055     struct disjoint_info
0056     {
0057         enum type
0058         {
0059             intersect,
0060             disjoint_no_vertex,
0061             disjoint_vertex
0062         };
0063         disjoint_info(type t) : m_(t){}
0064         operator type () const {return m_;}
0065         type m_;
0066     private :
0067         //prevent automatic conversion for any other built-in types
0068         template <typename T>
0069         operator T () const;
0070     };
0071
0072     template
0073     <
0074         typename Segment, typename Box,
0075         typename AzimuthStrategy,
0076         typename NormalizeStrategy,
0077         typename DisjointPointBoxStrategy,
0078         typename DisjointBoxBoxStrategy
0079     >
0080     static inline bool apply(Segment const& segment,
0081                              Box const& box,
0082                              AzimuthStrategy const& azimuth_strategy,
0083                              NormalizeStrategy const& normalize_strategy,
0084                              DisjointPointBoxStrategy const& disjoint_point_box_strategy,
0085                              DisjointBoxBoxStrategy const& disjoint_box_box_strategy)
0086     {
0087         point_type_t<Segment> vertex;
0088         return apply(segment, box, vertex,
0089                      azimuth_strategy,
0090                      normalize_strategy,
0091                      disjoint_point_box_strategy,
0092                      disjoint_box_box_strategy) != disjoint_info::intersect;
0093     }
0094
0095     template
0096     <
0097         typename Segment, typename Box,
0098         typename P,
0099         typename AzimuthStrategy,
0100         typename NormalizeStrategy,
0101         typename DisjointPointBoxStrategy,
0102         typename DisjointBoxBoxStrategy
0103     >
0104     static inline disjoint_info apply(Segment const& segment,
0105                                       Box const& box,
0106                                       P& vertex,
0107                                       AzimuthStrategy const& azimuth_strategy,
0108                                       NormalizeStrategy const& ,
0109                                       DisjointPointBoxStrategy const& disjoint_point_box_strategy,
0110                                       DisjointBoxBoxStrategy const& disjoint_box_box_strategy)
0111     {
0112         assert_dimension_equal<Segment, Box>();
0113
0114         using segment_point_type = point_type_t<Segment>;
0115
0116         segment_point_type p0, p1;
0117         geometry::detail::assign_point_from_index<0>(segment, p0);
0118         geometry::detail::assign_point_from_index<1>(segment, p1);
0119
0120         // Vertex is not computed here
0121         disjoint_info disjoint_return_value = disjoint_info::disjoint_no_vertex;
0122
0123         // Simplest cases first
0124
0125         // Case 1: if box contains one of segment's endpoints then they are not disjoint
0126         if ( ! disjoint_point_box(p0, box, disjoint_point_box_strategy)
0127           || ! disjoint_point_box(p1, box, disjoint_point_box_strategy) )
0128         {
0129             return disjoint_info::intersect;
0130         }
0131
0132         // Case 2: disjoint if bounding boxes are disjoint
0133
0134         using coor_t = coordinate_type_t<segment_point_type>;
0135
0136         segment_point_type p0_normalized;
0137         NormalizeStrategy::apply(p0, p0_normalized);
0138         segment_point_type p1_normalized;
0139         NormalizeStrategy::apply(p1, p1_normalized);
0140
0141         coor_t lon1 = geometry::get_as_radian<0>(p0_normalized);
0142         coor_t lat1 = geometry::get_as_radian<1>(p0_normalized);
0143         coor_t lon2 = geometry::get_as_radian<0>(p1_normalized);
0144         coor_t lat2 = geometry::get_as_radian<1>(p1_normalized);
0145
0146         if (lon1 > lon2)
0147         {
0148             std::swap(lon1, lon2);
0149             std::swap(lat1, lat2);
0150         }
0151
0152         geometry::model::box<segment_point_type> box_seg;
0153
0154         strategy::envelope::detail::envelope_segment_impl
0155             <
0156                 CS_Tag
0157             >::template apply<geometry::radian>(lon1, lat1,
0158                                                 lon2, lat2,
0159                                                 box_seg,
0160                                                 azimuth_strategy);
0161
0162         if (disjoint_box_box(box, box_seg, disjoint_box_box_strategy))
0163         {
0164             return disjoint_return_value;
0165         }
0166
0167         // Case 3: test intersection by comparing angles
0168
0169         coor_t alp1, a_b0, a_b1, a_b2, a_b3;
0170
0171         coor_t b_lon_min = geometry::get_as_radian<geometry::min_corner, 0>(box);
0172         coor_t b_lat_min = geometry::get_as_radian<geometry::min_corner, 1>(box);
0173         coor_t b_lon_max = geometry::get_as_radian<geometry::max_corner, 0>(box);
0174         coor_t b_lat_max = geometry::get_as_radian<geometry::max_corner, 1>(box);
0175
0176         azimuth_strategy.apply(lon1, lat1, lon2, lat2, alp1);
0177         azimuth_strategy.apply(lon1, lat1, b_lon_min, b_lat_min, a_b0);
0178         azimuth_strategy.apply(lon1, lat1, b_lon_max, b_lat_min, a_b1);
0179         azimuth_strategy.apply(lon1, lat1, b_lon_min, b_lat_max, a_b2);
0180         azimuth_strategy.apply(lon1, lat1, b_lon_max, b_lat_max, a_b3);
0181
0182         int s0 = formula::azimuth_side_value(alp1, a_b0);
0183         int s1 = formula::azimuth_side_value(alp1, a_b1);
0184         int s2 = formula::azimuth_side_value(alp1, a_b2);
0185         int s3 = formula::azimuth_side_value(alp1, a_b3);
0186
0187         if (s0 == 0 || s1 == 0 || s2 == 0 || s3 == 0)
0188         {
0189             return disjoint_info::intersect;
0190         }
0191
0192         bool s0_positive = s0 > 0;
0193         bool s1_positive = s1 > 0;
0194         bool s2_positive = s2 > 0;
0195         bool s3_positive = s3 > 0;
0196
0197         bool all_positive = s0_positive && s1_positive && s2_positive && s3_positive;
0198         bool all_non_positive = !(s0_positive || s1_positive || s2_positive || s3_positive);
0199         bool vertex_north = lat1 + lat2 > 0;
0200
0201         if ((all_positive && vertex_north) || (all_non_positive && !vertex_north))
0202         {
0203             return disjoint_info::disjoint_no_vertex;
0204         }
0205
0206         if (!all_positive && !all_non_positive)
0207         {
0208             return disjoint_info::intersect;
0209         }
0210
0211         // Case 4: The only intersection case not covered above is when all four
0212         // points of the box are above (below) the segment in northern (southern)
0213         // hemisphere. Then we have to compute the vertex of the segment
0214
0215         coor_t vertex_lat;
0216
0217         if ((lat1 < b_lat_min && vertex_north)
0218                 || (lat1 > b_lat_max && !vertex_north))
0219         {
0220             coor_t b_lat_below; //latitude of box closest to equator
0221
0222             if (vertex_north)
0223             {
0224                 vertex_lat = geometry::get_as_radian<geometry::max_corner, 1>(box_seg);
0225                 b_lat_below = b_lat_min;
0226             } else {
0227                 vertex_lat = geometry::get_as_radian<geometry::min_corner, 1>(box_seg);
0228                 b_lat_below = b_lat_max;
0229             }
0230
0231             //optimization TODO: computing the spherical longitude should suffice for
0232             // the majority of cases
0233             coor_t vertex_lon = geometry::formula::vertex_longitude<coor_t, CS_Tag>
0234                                     ::apply(lon1, lat1,
0235                                             lon2, lat2,
0236                                             vertex_lat,
0237                                             alp1,
0238                                             azimuth_strategy);
0239
0240             geometry::set_from_radian<0>(vertex, vertex_lon);
0241             geometry::set_from_radian<1>(vertex, vertex_lat);
0242             disjoint_return_value = disjoint_info::disjoint_vertex; //vertex_computed
0243
0244             // Check if the vertex point is within the band defined by the
0245             // minimum and maximum longitude of the box; if yes, then return
0246             // false if the point is above the min latitude of the box; return
0247             // true in all other cases
0248             if (vertex_lon >= b_lon_min && vertex_lon <= b_lon_max
0249                     && std::abs(vertex_lat) > std::abs(b_lat_below))
0250             {
0251                 return disjoint_info::intersect;
0252             }
0253         }
0254
0255         return disjoint_return_value;
0256     }
0257 };
0258
0259 struct disjoint_segment_box
0260 {
0261     template <typename Segment, typename Box, typename Strategy>
0262     static inline bool apply(Segment const& segment,
0263                              Box const& box,
0264                              Strategy const& strategy)
0265     {
0266         return strategy.disjoint(segment, box).apply(segment, box);
0267     }
0268 };
0269
0270 }} // namespace detail::disjoint
0271 #endif // DOXYGEN_NO_DETAIL
0272
0273
0274 #ifndef DOXYGEN_NO_DISPATCH
0275 namespace dispatch
0276 {
0277
0278
0279 template <typename Segment, typename Box, std::size_t DimensionCount>
0280 struct disjoint<Segment, Box, DimensionCount, segment_tag, box_tag, false>
0281         : detail::disjoint::disjoint_segment_box
0282 {};
0283
0284
0285 } // namespace dispatch
0286 #endif // DOXYGEN_NO_DISPATCH
0287
0288
0289 }} // namespace boost::geometry
0290
0291
0292 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP