detail/relate/areal_areal.hpp

0001 // Boost.Geometry (aka GGL, Generic Geometry Library)
0002
0003 // Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
0004
0005 // This file was modified by Oracle on 2013-2022.
0006 // Modifications copyright (c) 2013-2022 Oracle and/or its affiliates.
0007
0008 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
0009
0010 // Use, modification and distribution is subject to the Boost Software License,
0011 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
0012 // http://www.boost.org/LICENSE_1_0.txt)
0013
0014 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_AREAL_AREAL_HPP
0015 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_AREAL_AREAL_HPP
0016
0017 #include <boost/geometry/core/topological_dimension.hpp>
0018
0019 #include <boost/geometry/util/condition.hpp>
0020 #include <boost/geometry/util/range.hpp>
0021
0022 #include <boost/geometry/algorithms/num_interior_rings.hpp>
0023 #include <boost/geometry/algorithms/detail/point_on_border.hpp>
0024 #include <boost/geometry/algorithms/detail/sub_range.hpp>
0025 #include <boost/geometry/algorithms/detail/single_geometry.hpp>
0026
0027 #include <boost/geometry/algorithms/detail/relate/point_geometry.hpp>
0028 #include <boost/geometry/algorithms/detail/relate/turns.hpp>
0029 #include <boost/geometry/algorithms/detail/relate/boundary_checker.hpp>
0030 #include <boost/geometry/algorithms/detail/relate/follow_helpers.hpp>
0031
0032 #include <boost/geometry/geometries/helper_geometry.hpp>
0033
0034 #include <boost/geometry/util/numeric_cast.hpp>
0035
0036 namespace boost { namespace geometry
0037 {
0038
0039 #ifndef DOXYGEN_NO_DETAIL
0040 namespace detail { namespace relate {
0041
0042 // WARNING!
0043 // TODO: In the worst case calling this Pred in a loop for MultiPolygon/MultiPolygon may take O(NM)
0044 // Use the rtree in this case!
0045
0046 // may be used to set EI and EB for an Areal geometry for which no turns were generated
0047 template
0048 <
0049     typename OtherAreal,
0050     typename Result,
0051     typename PointInArealStrategy,
0052     bool TransposeResult
0053 >
0054 class no_turns_aa_pred
0055 {
0056 public:
0057     no_turns_aa_pred(OtherAreal const& other_areal,
0058                      Result & res,
0059                      PointInArealStrategy const& point_in_areal_strategy)
0060         : m_result(res)
0061         , m_point_in_areal_strategy(point_in_areal_strategy)
0062         , m_other_areal(other_areal)
0063         , m_flags(0)
0064     {
0065         // check which relations must be analysed
0066
0067         if ( ! may_update<interior, interior, '2', TransposeResult>(m_result)
0068           && ! may_update<boundary, interior, '1', TransposeResult>(m_result)
0069           && ! may_update<exterior, interior, '2', TransposeResult>(m_result) )
0070         {
0071             m_flags |= 1;
0072         }
0073
0074         if ( ! may_update<interior, exterior, '2', TransposeResult>(m_result)
0075           && ! may_update<boundary, exterior, '1', TransposeResult>(m_result) )
0076         {
0077             m_flags |= 2;
0078         }
0079     }
0080
0081     template <typename Areal>
0082     bool operator()(Areal const& areal)
0083     {
0084         using detail::within::point_in_geometry;
0085
0086         // if those flags are set nothing will change
0087         if ( m_flags == 3 )
0088         {
0089             return false;
0090         }
0091
0092         using point_type = typename geometry::point_type<Areal>::type;
0093         typename helper_geometry<point_type>::type pt;
0094         bool const ok = geometry::point_on_border(pt, areal);
0095
0096         // TODO: for now ignore, later throw an exception?
0097         if ( !ok )
0098         {
0099             return true;
0100         }
0101
0102         // check if the areal is inside the other_areal
0103         // TODO: This is O(N)
0104         // Run in a loop O(NM) - optimize!
0105         int const pig = point_in_geometry(pt,
0106                                           m_other_areal,
0107                                           m_point_in_areal_strategy);
0108         //BOOST_GEOMETRY_ASSERT( pig != 0 );
0109
0110         // inside
0111         if ( pig > 0 )
0112         {
0113             update<interior, interior, '2', TransposeResult>(m_result);
0114             update<boundary, interior, '1', TransposeResult>(m_result);
0115             update<exterior, interior, '2', TransposeResult>(m_result);
0116             m_flags |= 1;
0117
0118             // TODO: OPTIMIZE!
0119             // Only the interior rings of other ONE single geometry must be checked
0120             // NOT all geometries
0121
0122             // Check if any interior ring is outside
0123             ring_identifier ring_id(0, -1, 0);
0124             std::size_t const irings_count = geometry::num_interior_rings(areal);
0125             for ( ; static_cast<std::size_t>(ring_id.ring_index) < irings_count ;
0126                     ++ring_id.ring_index )
0127             {
0128                 typename detail::sub_range_return_type<Areal const>::type
0129                     range_ref = detail::sub_range(areal, ring_id);
0130
0131                 if ( boost::empty(range_ref) )
0132                 {
0133                     // TODO: throw exception?
0134                     continue; // ignore
0135                 }
0136
0137                 // TODO: O(N)
0138                 // Optimize!
0139                 int const hpig = point_in_geometry(range::front(range_ref),
0140                                                    m_other_areal,
0141                                                    m_point_in_areal_strategy);
0142
0143                 // hole outside
0144                 if ( hpig < 0 )
0145                 {
0146                     update<interior, exterior, '2', TransposeResult>(m_result);
0147                     update<boundary, exterior, '1', TransposeResult>(m_result);
0148                     m_flags |= 2;
0149                     break;
0150                 }
0151             }
0152         }
0153         // outside
0154         else
0155         {
0156             update<interior, exterior, '2', TransposeResult>(m_result);
0157             update<boundary, exterior, '1', TransposeResult>(m_result);
0158             m_flags |= 2;
0159
0160             // Check if any interior ring is inside
0161             ring_identifier ring_id(0, -1, 0);
0162             std::size_t const irings_count = geometry::num_interior_rings(areal);
0163             for ( ; static_cast<std::size_t>(ring_id.ring_index) < irings_count ;
0164                     ++ring_id.ring_index )
0165             {
0166                 typename detail::sub_range_return_type<Areal const>::type
0167                     range_ref = detail::sub_range(areal, ring_id);
0168
0169                 if ( boost::empty(range_ref) )
0170                 {
0171                     // TODO: throw exception?
0172                     continue; // ignore
0173                 }
0174
0175                 // TODO: O(N)
0176                 // Optimize!
0177                 int const hpig = point_in_geometry(range::front(range_ref),
0178                                                    m_other_areal,
0179                                                    m_point_in_areal_strategy);
0180
0181                 // hole inside
0182                 if ( hpig > 0 )
0183                 {
0184                     update<interior, interior, '2', TransposeResult>(m_result);
0185                     update<boundary, interior, '1', TransposeResult>(m_result);
0186                     update<exterior, interior, '2', TransposeResult>(m_result);
0187                     m_flags |= 1;
0188                     break;
0189                 }
0190             }
0191         }
0192
0193         return m_flags != 3 && !m_result.interrupt;
0194     }
0195
0196 private:
0197     Result & m_result;
0198     PointInArealStrategy const& m_point_in_areal_strategy;
0199     OtherAreal const& m_other_areal;
0200     int m_flags;
0201 };
0202
0203 // The implementation of an algorithm calculating relate() for A/A
0204 template <typename Geometry1, typename Geometry2>
0205 struct areal_areal
0206 {
0207     // check Linear / Areal
0208     BOOST_STATIC_ASSERT(topological_dimension<Geometry1>::value == 2
0209                      && topological_dimension<Geometry2>::value == 2);
0210
0211     static const bool interruption_enabled = true;
0212
0213     template <typename Result, typename Strategy>
0214     static inline void apply(Geometry1 const& geometry1, Geometry2 const& geometry2,
0215                              Result & result,
0216                              Strategy const& strategy)
0217     {
0218 // TODO: If Areal geometry may have infinite size, change the following line:
0219
0220         update<exterior, exterior, result_dimension<Geometry2>::value>(result);// FFFFFFFFd, d in [1,9] or T
0221
0222         if ( BOOST_GEOMETRY_CONDITION(result.interrupt) )
0223             return;
0224
0225         // get and analyse turns
0226         using turn_type = typename turns::get_turns
0227             <
0228                 Geometry1, Geometry2
0229             >::template turn_info_type<Strategy>::type;
0230         std::vector<turn_type> turns;
0231
0232         interrupt_policy_areal_areal<Result> interrupt_policy(geometry1, geometry2, result);
0233
0234         turns::get_turns<Geometry1, Geometry2>::apply(turns, geometry1, geometry2, interrupt_policy, strategy);
0235         if ( BOOST_GEOMETRY_CONDITION(result.interrupt) )
0236             return;
0237
0238         no_turns_aa_pred<Geometry2, Result, Strategy, false>
0239             pred1(geometry2, result, strategy);
0240         for_each_disjoint_geometry_if<0, Geometry1>::apply(turns.begin(), turns.end(), geometry1, pred1);
0241         if ( BOOST_GEOMETRY_CONDITION(result.interrupt) )
0242             return;
0243
0244         no_turns_aa_pred<Geometry1, Result, Strategy, true>
0245             pred2(geometry1, result, strategy);
0246         for_each_disjoint_geometry_if<1, Geometry2>::apply(turns.begin(), turns.end(), geometry2, pred2);
0247         if ( BOOST_GEOMETRY_CONDITION(result.interrupt) )
0248             return;
0249
0250         if ( turns.empty() )
0251             return;
0252
0253         if ( may_update<interior, interior, '2'>(result)
0254           || may_update<interior, exterior, '2'>(result)
0255           || may_update<boundary, interior, '1'>(result)
0256           || may_update<boundary, exterior, '1'>(result)
0257           || may_update<exterior, interior, '2'>(result) )
0258         {
0259             // sort turns
0260             using less_t = turns::less<0, turns::less_op_areal_areal<0>, Strategy>;
0261             std::sort(turns.begin(), turns.end(), less_t());
0262
0263             /*if ( may_update<interior, exterior, '2'>(result)
0264               || may_update<boundary, exterior, '1'>(result)
0265               || may_update<boundary, interior, '1'>(result)
0266               || may_update<exterior, interior, '2'>(result) )*/
0267             {
0268                 // analyse sorted turns
0269                 turns_analyser<turn_type, 0> analyser;
0270                 analyse_each_turn(result, analyser, turns.begin(), turns.end(), strategy);
0271
0272                 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) )
0273                     return;
0274             }
0275
0276             if ( may_update<interior, interior, '2'>(result)
0277               || may_update<interior, exterior, '2'>(result)
0278               || may_update<boundary, interior, '1'>(result)
0279               || may_update<boundary, exterior, '1'>(result)
0280               || may_update<exterior, interior, '2'>(result) )
0281             {
0282                 // analyse rings for which turns were not generated
0283                 // or only i/i or u/u was generated
0284                 uncertain_rings_analyser<0, Result, Geometry1, Geometry2, Strategy>
0285                     rings_analyser(result, geometry1, geometry2, strategy);
0286                 analyse_uncertain_rings<0>::apply(rings_analyser, turns.begin(), turns.end());
0287
0288                 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) )
0289                     return;
0290             }
0291         }
0292
0293         if ( may_update<interior, interior, '2', true>(result)
0294           || may_update<interior, exterior, '2', true>(result)
0295           || may_update<boundary, interior, '1', true>(result)
0296           || may_update<boundary, exterior, '1', true>(result)
0297           || may_update<exterior, interior, '2', true>(result) )
0298         {
0299             // sort turns
0300             using less_t = turns::less<1, turns::less_op_areal_areal<1>, Strategy>;
0301             std::sort(turns.begin(), turns.end(), less_t());
0302
0303             /*if ( may_update<interior, exterior, '2', true>(result)
0304               || may_update<boundary, exterior, '1', true>(result)
0305               || may_update<boundary, interior, '1', true>(result)
0306               || may_update<exterior, interior, '2', true>(result) )*/
0307             {
0308                 // analyse sorted turns
0309                 turns_analyser<turn_type, 1> analyser;
0310                 analyse_each_turn(result, analyser, turns.begin(), turns.end(), strategy);
0311
0312                 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) )
0313                     return;
0314             }
0315
0316             if ( may_update<interior, interior, '2', true>(result)
0317               || may_update<interior, exterior, '2', true>(result)
0318               || may_update<boundary, interior, '1', true>(result)
0319               || may_update<boundary, exterior, '1', true>(result)
0320               || may_update<exterior, interior, '2', true>(result) )
0321             {
0322                 // analyse rings for which turns were not generated
0323                 // or only i/i or u/u was generated
0324                 uncertain_rings_analyser<1, Result, Geometry2, Geometry1, Strategy>
0325                     rings_analyser(result, geometry2, geometry1, strategy);
0326                 analyse_uncertain_rings<1>::apply(rings_analyser, turns.begin(), turns.end());
0327
0328                 //if ( result.interrupt )
0329                 //    return;
0330             }
0331         }
0332     }
0333
0334     // interrupt policy which may be passed to get_turns to interrupt the analysis
0335     // based on the info in the passed result/mask
0336     template <typename Result>
0337     class interrupt_policy_areal_areal
0338     {
0339     public:
0340         static bool const enabled = true;
0341
0342         interrupt_policy_areal_areal(Geometry1 const& geometry1,
0343                                      Geometry2 const& geometry2,
0344                                      Result & result)
0345             : m_result(result)
0346             , m_geometry1(geometry1)
0347             , m_geometry2(geometry2)
0348         {}
0349
0350         template <typename Range>
0351         inline bool apply(Range const& turns)
0352         {
0353             for (auto it = boost::begin(turns) ; it != boost::end(turns) ; ++it)
0354             {
0355                 per_turn<0>(*it);
0356                 per_turn<1>(*it);
0357             }
0358
0359             return m_result.interrupt;
0360         }
0361
0362     private:
0363         template <std::size_t OpId, typename Turn>
0364         inline void per_turn(Turn const& turn)
0365         {
0366             //static const std::size_t other_op_id = (OpId + 1) % 2;
0367             static const bool transpose_result = OpId != 0;
0368
0369             overlay::operation_type const op = turn.operations[OpId].operation;
0370
0371             if ( op == overlay::operation_union )
0372             {
0373                 // ignore u/u
0374                 /*if ( turn.operations[other_op_id].operation != overlay::operation_union )
0375                 {
0376                     update<interior, exterior, '2', transpose_result>(m_result);
0377                     update<boundary, exterior, '1', transpose_result>(m_result);
0378                 }*/
0379
0380                 update<boundary, boundary, '0', transpose_result>(m_result);
0381             }
0382             else if ( op == overlay::operation_intersection )
0383             {
0384                 // ignore i/i
0385                 /*if ( turn.operations[other_op_id].operation != overlay::operation_intersection )
0386                 {
0387                     // not correct e.g. for G1 touching G2 in a point where a hole is touching the exterior ring
0388                     // in this case 2 turns i/... and u/u will be generated for this IP
0389                     //update<interior, interior, '2', transpose_result>(m_result);
0390
0391                     //update<boundary, interior, '1', transpose_result>(m_result);
0392                 }*/
0393
0394                 update<boundary, boundary, '0', transpose_result>(m_result);
0395             }
0396             else if ( op == overlay::operation_continue )
0397             {
0398                 update<boundary, boundary, '1', transpose_result>(m_result);
0399                 update<interior, interior, '2', transpose_result>(m_result);
0400             }
0401             else if ( op == overlay::operation_blocked )
0402             {
0403                 update<boundary, boundary, '1', transpose_result>(m_result);
0404                 update<interior, exterior, '2', transpose_result>(m_result);
0405             }
0406         }
0407
0408         Result & m_result;
0409         Geometry1 const& m_geometry1;
0410         Geometry2 const& m_geometry2;
0411     };
0412
0413     // This analyser should be used like Input or SinglePass Iterator
0414     // IMPORTANT! It should be called also for the end iterator - last
0415     template <typename TurnInfo, std::size_t OpId>
0416     class turns_analyser
0417     {
0418         typedef typename TurnInfo::point_type turn_point_type;
0419
0420         static const std::size_t op_id = OpId;
0421         static const std::size_t other_op_id = (OpId + 1) % 2;
0422         static const bool transpose_result = OpId != 0;
0423
0424     public:
0425         turns_analyser()
0426             : m_previous_turn_ptr(0)
0427             , m_previous_operation(overlay::operation_none)
0428             , m_enter_detected(false)
0429             , m_exit_detected(false)
0430         {}
0431
0432         template <typename Result, typename TurnIt, typename EqPPStrategy>
0433         void apply(Result & result, TurnIt it, EqPPStrategy const& strategy)
0434         {
0435             //BOOST_GEOMETRY_ASSERT( it != last );
0436
0437             overlay::operation_type const op = it->operations[op_id].operation;
0438
0439             if ( op != overlay::operation_union
0440               && op != overlay::operation_intersection
0441               && op != overlay::operation_blocked
0442               && op != overlay::operation_continue )
0443             {
0444                 return;
0445             }
0446
0447             segment_identifier const& seg_id = it->operations[op_id].seg_id;
0448             //segment_identifier const& other_id = it->operations[other_op_id].seg_id;
0449
0450             const bool first_in_range = m_seg_watcher.update(seg_id);
0451
0452             if ( m_previous_turn_ptr )
0453             {
0454                 if ( m_exit_detected /*m_previous_operation == overlay::operation_union*/ )
0455                 {
0456                     // real exit point - may be multiple
0457                     if ( first_in_range
0458                       || ! turn_on_the_same_ip<op_id>(*m_previous_turn_ptr, *it, strategy) )
0459                     {
0460                         update_exit(result);
0461                         m_exit_detected = false;
0462                     }
0463                     // fake exit point, reset state
0464                     else if ( op != overlay::operation_union )
0465                     {
0466                         m_exit_detected = false;
0467                     }
0468                 }
0469                 /*else*/
0470                 if ( m_enter_detected /*m_previous_operation == overlay::operation_intersection*/ )
0471                 {
0472                     // real entry point
0473                     if ( first_in_range
0474                       || ! turn_on_the_same_ip<op_id>(*m_previous_turn_ptr, *it, strategy) )
0475                     {
0476                         update_enter(result);
0477                         m_enter_detected = false;
0478                     }
0479                     // fake entry point, reset state
0480                     else if ( op != overlay::operation_intersection )
0481                     {
0482                         m_enter_detected = false;
0483                     }
0484                 }
0485             }
0486
0487             if ( op == overlay::operation_union )
0488             {
0489                 // already set in interrupt policy
0490                 //update<boundary, boundary, '0', transpose_result>(m_result);
0491
0492                 // ignore u/u
0493                 //if ( it->operations[other_op_id].operation != overlay::operation_union )
0494                 {
0495                     m_exit_detected = true;
0496                 }
0497             }
0498             else if ( op == overlay::operation_intersection )
0499             {
0500                 // ignore i/i
0501                 if ( it->operations[other_op_id].operation != overlay::operation_intersection )
0502                 {
0503                     // this was set in the interrupt policy but it was wrong
0504                     // also here it's wrong since it may be a fake entry point
0505                     //update<interior, interior, '2', transpose_result>(result);
0506
0507                     // already set in interrupt policy
0508                     //update<boundary, boundary, '0', transpose_result>(result);
0509                     m_enter_detected = true;
0510                 }
0511             }
0512             else if ( op == overlay::operation_blocked )
0513             {
0514                 // already set in interrupt policy
0515             }
0516             else // if ( op == overlay::operation_continue )
0517             {
0518                 // already set in interrupt policy
0519             }
0520
0521             // store ref to previously analysed (valid) turn
0522             m_previous_turn_ptr = boost::addressof(*it);
0523             // and previously analysed (valid) operation
0524             m_previous_operation = op;
0525         }
0526
0527         // it == last
0528         template <typename Result>
0529         void apply(Result & result)
0530         {
0531             //BOOST_GEOMETRY_ASSERT( first != last );
0532
0533             if ( m_exit_detected /*m_previous_operation == overlay::operation_union*/ )
0534             {
0535                 update_exit(result);
0536                 m_exit_detected = false;
0537             }
0538
0539             if ( m_enter_detected /*m_previous_operation == overlay::operation_intersection*/ )
0540             {
0541                 update_enter(result);
0542                 m_enter_detected = false;
0543             }
0544         }
0545
0546         template <typename Result>
0547         static inline void update_exit(Result & result)
0548         {
0549             update<interior, exterior, '2', transpose_result>(result);
0550             update<boundary, exterior, '1', transpose_result>(result);
0551         }
0552
0553         template <typename Result>
0554         static inline void update_enter(Result & result)
0555         {
0556             update<interior, interior, '2', transpose_result>(result);
0557             update<boundary, interior, '1', transpose_result>(result);
0558             update<exterior, interior, '2', transpose_result>(result);
0559         }
0560
0561     private:
0562         segment_watcher<same_ring> m_seg_watcher;
0563         TurnInfo * m_previous_turn_ptr;
0564         overlay::operation_type m_previous_operation;
0565         bool m_enter_detected;
0566         bool m_exit_detected;
0567     };
0568
0569     // call analyser.apply() for each turn in range
0570     // IMPORTANT! The analyser is also called for the end iterator - last
0571     template
0572     <
0573         typename Result,
0574         typename Analyser,
0575         typename TurnIt,
0576         typename EqPPStrategy
0577     >
0578     static inline void analyse_each_turn(Result & res,
0579                                          Analyser & analyser,
0580                                          TurnIt first, TurnIt last,
0581                                          EqPPStrategy const& strategy)
0582     {
0583         if ( first == last )
0584             return;
0585
0586         for ( TurnIt it = first ; it != last ; ++it )
0587         {
0588             analyser.apply(res, it, strategy);
0589
0590             if ( BOOST_GEOMETRY_CONDITION(res.interrupt) )
0591                 return;
0592         }
0593
0594         analyser.apply(res);
0595     }
0596
0597     template
0598     <
0599         std::size_t OpId,
0600         typename Result,
0601         typename Geometry,
0602         typename OtherGeometry,
0603         typename PointInArealStrategy
0604     >
0605     class uncertain_rings_analyser
0606     {
0607         static const bool transpose_result = OpId != 0;
0608         static const int other_id = (OpId + 1) % 2;
0609
0610     public:
0611         inline uncertain_rings_analyser(Result & result,
0612                                         Geometry const& geom,
0613                                         OtherGeometry const& other_geom,
0614                                         PointInArealStrategy const& point_in_areal_strategy)
0615             : geometry(geom)
0616             , other_geometry(other_geom)
0617             , interrupt(result.interrupt) // just in case, could be false as well
0618             , m_result(result)
0619             , m_point_in_areal_strategy(point_in_areal_strategy)
0620             , m_flags(0)
0621         {
0622             // check which relations must be analysed
0623             // NOTE: 1 and 4 could probably be connected
0624
0625             if ( ! may_update<interior, interior, '2', transpose_result>(m_result) )
0626             {
0627                 m_flags |= 1;
0628             }
0629
0630             if ( ! may_update<interior, exterior, '2', transpose_result>(m_result)
0631               && ! may_update<boundary, exterior, '1', transpose_result>(m_result) )
0632             {
0633                 m_flags |= 2;
0634             }
0635
0636             if ( ! may_update<boundary, interior, '1', transpose_result>(m_result)
0637               && ! may_update<exterior, interior, '2', transpose_result>(m_result) )
0638             {
0639                 m_flags |= 4;
0640             }
0641         }
0642
0643         inline void no_turns(segment_identifier const& seg_id)
0644         {
0645             // if those flags are set nothing will change
0646             if ( m_flags == 7 )
0647             {
0648                 return;
0649             }
0650
0651             auto const& sub_range = detail::sub_range(geometry, seg_id);
0652
0653             if ( boost::empty(sub_range) )
0654             {
0655                 // TODO: throw an exception?
0656                 return; // ignore
0657             }
0658
0659             // TODO: possible optimization
0660             // if the range is an interior ring we may use other IPs generated for this single geometry
0661             // to know which other single geometries should be checked
0662
0663             // TODO: optimize! e.g. use spatial index
0664             // O(N) - running it in a loop gives O(NM)
0665             using detail::within::point_in_geometry;
0666             int const pig = point_in_geometry(range::front(sub_range),
0667                                               other_geometry,
0668                                               m_point_in_areal_strategy);
0669
0670             //BOOST_GEOMETRY_ASSERT(pig != 0);
0671             if ( pig > 0 )
0672             {
0673                 update<interior, interior, '2', transpose_result>(m_result);
0674                 m_flags |= 1;
0675
0676                 update<boundary, interior, '1', transpose_result>(m_result);
0677                 update<exterior, interior, '2', transpose_result>(m_result);
0678                 m_flags |= 4;
0679             }
0680             else
0681             {
0682                 update<boundary, exterior, '1', transpose_result>(m_result);
0683                 update<interior, exterior, '2', transpose_result>(m_result);
0684                 m_flags |= 2;
0685             }
0686
0687 // TODO: break if all things are set
0688 // also some of them could be checked outside, before the analysis
0689 // In this case we shouldn't relay just on dummy flags
0690 // Flags should be initialized with proper values
0691 // or the result should be checked directly
0692 // THIS IS ALSO TRUE FOR OTHER ANALYSERS! in L/L and L/A
0693
0694             interrupt = m_flags == 7 || m_result.interrupt;
0695         }
0696
0697         template <typename TurnIt>
0698         inline void turns(TurnIt first, TurnIt last)
0699         {
0700             // if those flags are set nothing will change
0701             if ( (m_flags & 6) == 6 )
0702             {
0703                 return;
0704             }
0705
0706             bool found_ii = false;
0707             bool found_uu = false;
0708
0709             for ( TurnIt it = first ; it != last ; ++it )
0710             {
0711                 if ( it->operations[0].operation == overlay::operation_intersection
0712                   && it->operations[1].operation == overlay::operation_intersection )
0713                 {
0714                     found_ii = true;
0715                 }
0716                 else if ( it->operations[0].operation == overlay::operation_union
0717                        && it->operations[1].operation == overlay::operation_union )
0718                 {
0719                     found_uu = true;
0720                 }
0721                 else // ignore
0722                 {
0723                     return; // don't interrupt
0724                 }
0725             }
0726
0727             // only i/i was generated for this ring
0728             if ( found_ii )
0729             {
0730                 update<interior, interior, '2', transpose_result>(m_result);
0731                 m_flags |= 1;
0732
0733                 //update<boundary, boundary, '0', transpose_result>(m_result);
0734
0735                 update<boundary, interior, '1', transpose_result>(m_result);
0736                 update<exterior, interior, '2', transpose_result>(m_result);
0737                 m_flags |= 4;
0738             }
0739
0740             // only u/u was generated for this ring
0741             if ( found_uu )
0742             {
0743                 update<boundary, exterior, '1', transpose_result>(m_result);
0744                 update<interior, exterior, '2', transpose_result>(m_result);
0745                 m_flags |= 2;
0746             }
0747
0748             interrupt = m_flags == 7 || m_result.interrupt; // interrupt if the result won't be changed in the future
0749         }
0750
0751         Geometry const& geometry;
0752         OtherGeometry const& other_geometry;
0753         bool interrupt;
0754
0755     private:
0756         Result & m_result;
0757         PointInArealStrategy const& m_point_in_areal_strategy;
0758         int m_flags;
0759     };
0760
0761     template <std::size_t OpId>
0762     class analyse_uncertain_rings
0763     {
0764     public:
0765         template <typename Analyser, typename TurnIt>
0766         static inline void apply(Analyser & analyser, TurnIt first, TurnIt last)
0767         {
0768             if ( first == last )
0769                 return;
0770
0771             for_preceding_rings(analyser, *first);
0772             //analyser.per_turn(*first);
0773
0774             TurnIt prev = first;
0775             for ( ++first ; first != last ; ++first, ++prev )
0776             {
0777                 // same multi
0778                 if ( prev->operations[OpId].seg_id.multi_index
0779                   == first->operations[OpId].seg_id.multi_index )
0780                 {
0781                     // same ring
0782                     if ( prev->operations[OpId].seg_id.ring_index
0783                       == first->operations[OpId].seg_id.ring_index )
0784                     {
0785                         //analyser.per_turn(*first);
0786                     }
0787                     // same multi, next ring
0788                     else
0789                     {
0790                         //analyser.end_ring(*prev);
0791                         analyser.turns(prev, first);
0792
0793                         //if ( prev->operations[OpId].seg_id.ring_index + 1
0794                         //   < first->operations[OpId].seg_id.ring_index)
0795                         {
0796                             for_no_turns_rings(analyser,
0797                                                *first,
0798                                                prev->operations[OpId].seg_id.ring_index + 1,
0799                                                first->operations[OpId].seg_id.ring_index);
0800                         }
0801
0802                         //analyser.per_turn(*first);
0803                     }
0804                 }
0805                 // next multi
0806                 else
0807                 {
0808                     //analyser.end_ring(*prev);
0809                     analyser.turns(prev, first);
0810                     for_following_rings(analyser, *prev);
0811                     for_preceding_rings(analyser, *first);
0812                     //analyser.per_turn(*first);
0813                 }
0814
0815                 if ( analyser.interrupt )
0816                 {
0817                     return;
0818                 }
0819             }
0820
0821             //analyser.end_ring(*prev);
0822             analyser.turns(prev, first); // first == last
0823             for_following_rings(analyser, *prev);
0824         }
0825
0826     private:
0827         template <typename Analyser, typename Turn>
0828         static inline void for_preceding_rings(Analyser & analyser, Turn const& turn)
0829         {
0830             segment_identifier const& seg_id = turn.operations[OpId].seg_id;
0831
0832             for_no_turns_rings(analyser, turn, -1, seg_id.ring_index);
0833         }
0834
0835         template <typename Analyser, typename Turn>
0836         static inline void for_following_rings(Analyser & analyser, Turn const& turn)
0837         {
0838             segment_identifier const& seg_id = turn.operations[OpId].seg_id;
0839
0840             signed_size_type
0841                 count = util::numeric_cast<signed_size_type>(
0842                             geometry::num_interior_rings(
0843                                 detail::single_geometry(analyser.geometry, seg_id)));
0844
0845             for_no_turns_rings(analyser, turn, seg_id.ring_index + 1, count);
0846         }
0847
0848         template <typename Analyser, typename Turn>
0849         static inline void for_no_turns_rings(Analyser & analyser,
0850                                               Turn const& turn,
0851                                               signed_size_type first,
0852                                               signed_size_type last)
0853         {
0854             segment_identifier seg_id = turn.operations[OpId].seg_id;
0855
0856             for ( seg_id.ring_index = first ; seg_id.ring_index < last ; ++seg_id.ring_index )
0857             {
0858                 analyser.no_turns(seg_id);
0859             }
0860         }
0861     };
0862 };
0863
0864 }} // namespace detail::relate
0865 #endif // DOXYGEN_NO_DETAIL
0866
0867 }} // namespace boost::geometry
0868
0869 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_AREAL_AREAL_HPP