ceph/src/boost/boost/geometry/algorithms/detail/disjoint/segment_box.hpp

   1 // Boost.Geometry (aka GGL, Generic Geometry Library)
   2
   3 // Copyright (c) 2007-2014 Barend Gehrels, Amsterdam, the Netherlands.
   4 // Copyright (c) 2008-2014 Bruno Lalande, Paris, France.
   5 // Copyright (c) 2009-2014 Mateusz Loskot, London, UK.
   6 // Copyright (c) 2013-2014 Adam Wulkiewicz, Lodz, Poland.
   7
   8 // This file was modified by Oracle on 2013-2021.
   9 // Modifications copyright (c) 2013-2021, Oracle and/or its affiliates.
  10
  11 // Contributed and/or modified by Vissarion Fysikopoulos, on behalf of Oracle
  12 // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
  13 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
  14
  15 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
  16 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
  17
  18 // Use, modification and distribution is subject to the Boost Software License,
  19 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  20 // http://www.boost.org/LICENSE_1_0.txt)
  21
  22 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP
  23 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP
  24
  25 #include <cstddef>
  26
  27 #include <boost/geometry/core/tags.hpp>
  28 #include <boost/geometry/core/radian_access.hpp>
  29
  30 #include <boost/geometry/algorithms/detail/assign_indexed_point.hpp>
  31 #include <boost/geometry/algorithms/detail/disjoint/point_box.hpp>
  32 #include <boost/geometry/algorithms/detail/disjoint/box_box.hpp>
  33 #include <boost/geometry/algorithms/detail/envelope/segment.hpp>
  34 #include <boost/geometry/algorithms/detail/normalize.hpp>
  35 #include <boost/geometry/algorithms/dispatch/disjoint.hpp>
  36
  37 #include <boost/geometry/formulas/vertex_longitude.hpp>
  38
  39 #include <boost/geometry/geometries/box.hpp>
  40
  41 // Temporary, for envelope_segment_impl
  42 #include <boost/geometry/strategy/spherical/envelope_segment.hpp>
  43
  44 namespace boost { namespace geometry
  45 {
  46
  47
  48 #ifndef DOXYGEN_NO_DETAIL
  49 namespace detail { namespace disjoint
  50 {
  51
  52 template <typename CS_Tag>
  53 struct disjoint_segment_box_sphere_or_spheroid
  54 {
  55     struct disjoint_info
  56     {
  57         enum type
  58         {
  59             intersect,
  60             disjoint_no_vertex,
  61             disjoint_vertex
  62         };
  63         disjoint_info(type t) : m_(t){}
  64         operator type () const {return m_;}
  65         type m_;
  66     private :
  67         //prevent automatic conversion for any other built-in types
  68         template <typename T>
  69         operator T () const;
  70     };
  71
  72     template
  73     <
  74         typename Segment, typename Box,
  75         typename AzimuthStrategy,
  76         typename NormalizeStrategy,
  77         typename DisjointPointBoxStrategy,
  78         typename DisjointBoxBoxStrategy
  79     >
  80     static inline bool apply(Segment const& segment,
  81                              Box const& box,
  82                              AzimuthStrategy const& azimuth_strategy,
  83                              NormalizeStrategy const& normalize_strategy,
  84                              DisjointPointBoxStrategy const& disjoint_point_box_strategy,
  85                              DisjointBoxBoxStrategy const& disjoint_box_box_strategy)
  86     {
  87         typedef typename point_type<Segment>::type segment_point;
  88         segment_point vertex;
  89         return apply(segment, box, vertex,
  90                      azimuth_strategy,
  91                      normalize_strategy,
  92                      disjoint_point_box_strategy,
  93                      disjoint_box_box_strategy) != disjoint_info::intersect;
  94     }
  95
  96     template
  97     <
  98         typename Segment, typename Box,
  99         typename P,
 100         typename AzimuthStrategy,
 101         typename NormalizeStrategy,
 102         typename DisjointPointBoxStrategy,
 103         typename DisjointBoxBoxStrategy
 104     >
 105     static inline disjoint_info apply(Segment const& segment,
 106                                       Box const& box,
 107                                       P& vertex,
 108                                       AzimuthStrategy const& azimuth_strategy,
 109                                       NormalizeStrategy const& ,
 110                                       DisjointPointBoxStrategy const& disjoint_point_box_strategy,
 111                                       DisjointBoxBoxStrategy const& disjoint_box_box_strategy)
 112     {
 113         assert_dimension_equal<Segment, Box>();
 114
 115         typedef typename point_type<Segment>::type segment_point_type;
 116
 117         segment_point_type p0, p1;
 118         geometry::detail::assign_point_from_index<0>(segment, p0);
 119         geometry::detail::assign_point_from_index<1>(segment, p1);
 120
 121         //vertex not computed here
 122         disjoint_info disjoint_return_value = disjoint_info::disjoint_no_vertex;
 123
 124         // Simplest cases first
 125
 126         // Case 1: if box contains one of segment's endpoints then they are not disjoint
 127         if ( ! disjoint_point_box(p0, box, disjoint_point_box_strategy)
 128           || ! disjoint_point_box(p1, box, disjoint_point_box_strategy) )
 129         {
 130             return disjoint_info::intersect;
 131         }
 132
 133         // Case 2: disjoint if bounding boxes are disjoint
 134
 135         typedef typename coordinate_type<segment_point_type>::type CT;
 136
 137         segment_point_type p0_normalized;
 138         NormalizeStrategy::apply(p0, p0_normalized);
 139         segment_point_type p1_normalized;
 140         NormalizeStrategy::apply(p1, p1_normalized);
 141
 142         CT lon1 = geometry::get_as_radian<0>(p0_normalized);
 143         CT lat1 = geometry::get_as_radian<1>(p0_normalized);
 144         CT lon2 = geometry::get_as_radian<0>(p1_normalized);
 145         CT lat2 = geometry::get_as_radian<1>(p1_normalized);
 146
 147         if (lon1 > lon2)
 148         {
 149             std::swap(lon1, lon2);
 150             std::swap(lat1, lat2);
 151         }
 152
 153         geometry::model::box<segment_point_type> box_seg;
 154
 155         strategy::envelope::detail::envelope_segment_impl
 156             <
 157                 CS_Tag
 158             >::template apply<geometry::radian>(lon1, lat1,
 159                                                 lon2, lat2,
 160                                                 box_seg,
 161                                                 azimuth_strategy);
 162
 163         if (disjoint_box_box(box, box_seg, disjoint_box_box_strategy))
 164         {
 165             return disjoint_return_value;
 166         }
 167
 168         // Case 3: test intersection by comparing angles
 169
 170         CT alp1, a_b0, a_b1, a_b2, a_b3;
 171
 172         CT b_lon_min = geometry::get_as_radian<geometry::min_corner, 0>(box);
 173         CT b_lat_min = geometry::get_as_radian<geometry::min_corner, 1>(box);
 174         CT b_lon_max = geometry::get_as_radian<geometry::max_corner, 0>(box);
 175         CT b_lat_max = geometry::get_as_radian<geometry::max_corner, 1>(box);
 176
 177         azimuth_strategy.apply(lon1, lat1, lon2, lat2, alp1);
 178         azimuth_strategy.apply(lon1, lat1, b_lon_min, b_lat_min, a_b0);
 179         azimuth_strategy.apply(lon1, lat1, b_lon_max, b_lat_min, a_b1);
 180         azimuth_strategy.apply(lon1, lat1, b_lon_min, b_lat_max, a_b2);
 181         azimuth_strategy.apply(lon1, lat1, b_lon_max, b_lat_max, a_b3);
 182
 183         int s0 = formula::azimuth_side_value(alp1, a_b0);
 184         int s1 = formula::azimuth_side_value(alp1, a_b1);
 185         int s2 = formula::azimuth_side_value(alp1, a_b2);
 186         int s3 = formula::azimuth_side_value(alp1, a_b3);
 187
 188         if (s0 == 0 || s1 == 0 || s2 == 0 || s3 == 0)
 189         {
 190             return disjoint_info::intersect;
 191         }
 192
 193         bool s0_positive = s0 > 0;
 194         bool s1_positive = s1 > 0;
 195         bool s2_positive = s2 > 0;
 196         bool s3_positive = s3 > 0;
 197
 198         bool all_positive = s0_positive && s1_positive && s2_positive && s3_positive;
 199         bool all_non_positive = !(s0_positive || s1_positive || s2_positive || s3_positive);
 200         bool vertex_north = lat1 + lat2 > 0;
 201
 202         if ((all_positive && vertex_north) || (all_non_positive && !vertex_north))
 203         {
 204             return disjoint_info::disjoint_no_vertex;
 205         }
 206
 207         if (!all_positive && !all_non_positive)
 208         {
 209             return disjoint_info::intersect;
 210         }
 211
 212         // Case 4: The only intersection case not covered above is when all four
 213         // points of the box are above (below) the segment in northern (southern)
 214         // hemisphere. Then we have to compute the vertex of the segment
 215
 216         CT vertex_lat;
 217
 218         if ((lat1 < b_lat_min && vertex_north)
 219                 || (lat1 > b_lat_max && !vertex_north))
 220         {
 221             CT b_lat_below; //latitude of box closest to equator
 222
 223             if (vertex_north)
 224             {
 225                 vertex_lat = geometry::get_as_radian<geometry::max_corner, 1>(box_seg);
 226                 b_lat_below = b_lat_min;
 227             } else {
 228                 vertex_lat = geometry::get_as_radian<geometry::min_corner, 1>(box_seg);
 229                 b_lat_below = b_lat_max;
 230             }
 231
 232             //optimization TODO: computing the spherical longitude should suffice for
 233             // the majority of cases
 234             CT vertex_lon = geometry::formula::vertex_longitude<CT, CS_Tag>
 235                                     ::apply(lon1, lat1,
 236                                             lon2, lat2,
 237                                             vertex_lat,
 238                                             alp1,
 239                                             azimuth_strategy);
 240
 241             geometry::set_from_radian<0>(vertex, vertex_lon);
 242             geometry::set_from_radian<1>(vertex, vertex_lat);
 243             disjoint_return_value = disjoint_info::disjoint_vertex; //vertex_computed
 244
 245             // Check if the vertex point is within the band defined by the
 246             // minimum and maximum longitude of the box; if yes, then return
 247             // false if the point is above the min latitude of the box; return
 248             // true in all other cases
 249             if (vertex_lon >= b_lon_min && vertex_lon <= b_lon_max
 250                     && std::abs(vertex_lat) > std::abs(b_lat_below))
 251             {
 252                 return disjoint_info::intersect;
 253             }
 254         }
 255
 256         return disjoint_return_value;
 257     }
 258 };
 259
 260 struct disjoint_segment_box
 261 {
 262     template <typename Segment, typename Box, typename Strategy>
 263     static inline bool apply(Segment const& segment,
 264                              Box const& box,
 265                              Strategy const& strategy)
 266     {
 267         return strategy.disjoint(segment, box).apply(segment, box);
 268     }
 269 };
 270
 271 }} // namespace detail::disjoint
 272 #endif // DOXYGEN_NO_DETAIL
 273
 274
 275 #ifndef DOXYGEN_NO_DISPATCH
 276 namespace dispatch
 277 {
 278
 279
 280 template <typename Segment, typename Box, std::size_t DimensionCount>
 281 struct disjoint<Segment, Box, DimensionCount, segment_tag, box_tag, false>
 282         : detail::disjoint::disjoint_segment_box
 283 {};
 284
 285
 286 } // namespace dispatch
 287 #endif // DOXYGEN_NO_DISPATCH
 288
 289
 290 }} // namespace boost::geometry
 291
 292
 293 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP