1 // Copyright (c) 2015  Università della Svizzera italiana.
2 // All rights reserved.
3 //
4 // This file is part of CGAL (www.cgal.org).
5 //
6 // $URL: https://github.com/CGAL/cgal/blob/v5.3/Segment_Delaunay_graph_Linf_2/include/CGAL/Segment_Delaunay_graph_Linf_2/Voronoi_vertex_ring_C2.h $
7 // $Id: Voronoi_vertex_ring_C2.h 0779373 2020-03-26T13:31:46+01:00 Sébastien Loriot
8 // SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
9 //
10 //
11 // Author(s)     : Panagiotis Cheilaris, Sandeep Kumar Dey, Evanthia Papadopoulou
12 //philaris@gmail.com, sandeep.kr.dey@gmail.com, evanthia.papadopoulou@usi.ch
13 
14 #ifndef CGAL_SEGMENT_DELAUNAY_GRAPH_LINF_2_VORONOI_VERTEX_RING_C2_H
15 #define CGAL_SEGMENT_DELAUNAY_GRAPH_LINF_2_VORONOI_VERTEX_RING_C2_H
16 
17 #include <CGAL/license/Segment_Delaunay_graph_Linf_2.h>
18 
19 
20 
21 #include <CGAL/Segment_Delaunay_graph_Linf_2/Basic_predicates_C2.h>
22 #include <CGAL/Segment_Delaunay_graph_2/Are_same_points_C2.h>
23 #include <CGAL/Segment_Delaunay_graph_2/Are_same_segments_C2.h>
24 #include <CGAL/Segment_Delaunay_graph_2/Compare_x_2.h>
25 #include <CGAL/Segment_Delaunay_graph_2/Compare_y_2.h>
26 #include <CGAL/Side_of_bounded_square_2.h>
27 #include <CGAL/Segment_Delaunay_graph_Linf_2/Bisector_Linf.h>
28 #include <CGAL/tuple.h>
29 
30 
31 namespace CGAL {
32 
33 namespace SegmentDelaunayGraphLinf_2 {
34 
35 #define sdg_tuple_maker std::forward_as_tuple
36 
37 template<class K>
38 class Voronoi_vertex_ring_C2
39   : public Basic_predicates_C2<K>
40 {
41 public:
42   typedef Basic_predicates_C2<K> Base;
43 
44   typedef enum {PPP = 0, PPS, PSS, SSS} vertex_t;
45   struct PPP_Type {};
46   struct PPS_Type {};
47   struct PSS_Type {};
48   struct SSS_Type {};
49 
50   typedef typename Base::Point_2             Point_2;
51   typedef typename Base::Segment_2           Segment_2;
52   typedef typename Base::Line_2              Line_2;
53   typedef typename Base::Direction_2         Direction_2;
54   typedef typename Base::Site_2              Site_2;
55   typedef typename Base::FT                  FT;
56   typedef typename Base::RT                  RT;
57 
58   typedef typename Base::Homogeneous_point_2 Homogeneous_point_2;
59 
60   typedef typename Base::Orientation         Orientation;
61   typedef typename Base::Comparison_result   Comparison_result;
62   typedef typename Base::Oriented_side       Oriented_side;
63   typedef typename Base::Bounded_side        Bounded_side;
64   typedef typename Base::Sign                Sign;
65 
66   typedef typename Base::Polychainline_2     Polychainline_2;
67 
68   typedef typename Base::Bearing Bearing;
69 
70   using Base::compute_supporting_line;
71   using Base::compute_linf_projection_hom;
72   using Base::compute_linf_projection_nonhom;
73   using Base::oriented_side_of_line;
74   using Base::opposite_line;
75   using Base::compute_line_from_to;
76   using Base::compute_horizontal_projection;
77   using Base::compute_vertical_projection;
78   using Base::compute_linf_perpendicular;
79   using Base::is_site_horizontal;
80   using Base::is_site_vertical;
81   using Base::is_site_h_or_v;
82   using Base::is_line_h_or_v;
83   using Base::test_star;
84   using Base::compute_neg_45_line_at;
85   using Base::compute_pos_45_line_at;
86   using Base::compute_hor_line_at;
87   using Base::compute_ver_line_at;
88   using Base::has_positive_slope;
89   using Base::are_in_same_open_halfspace_of;
90   using Base::horseg_y_coord;
91   using Base::verseg_x_coord;
92   using Base::hvseg_coord;
93   using Base::compute_intersection_of_lines;
94   using Base::is_orth_dist_smaller_than_pt_dist;
95   using Base::touch_same_side;
96   using Base::coord_at;
97   using Base::orient_lines_linf;
98   using Base::compute_line_dir;
99   using Base::bisector_linf_line;
100   using Base::is_endpoint_of;
101   using Base::orient_line_endp;
102   using Base::orient_line_nonendp;
103   using Base::bearing;
104   using Base::bearing_diff;
105   using Base::center_from_corner_and_pt;
106   using Base::points_inside_touching_sides_v;
107   using Base::center_from_opposite_corners;
108   using Base::center_from_same_side_corners;
109   using Base::is_on_hv_seg_line;
110 
111 private:
112   typedef SegmentDelaunayGraph_2::Are_same_points_C2<K>
113           Are_same_points_2;
114   typedef SegmentDelaunayGraph_2::Are_same_segments_C2<K>
115           Are_same_segments_2;
116   typedef Side_of_bounded_square_2<K>    Side_of_bounded_square_2_Type;
117   typedef Bisector_Linf<K>               Bisector_Linf_Type;
118 
119   typedef SegmentDelaunayGraph_2::Compare_x_2<K> Compare_x_2_Sites_Type;
120   typedef SegmentDelaunayGraph_2::Compare_y_2<K> Compare_y_2_Sites_Type;
121 
122   Are_same_points_2                same_points;
123   Are_same_segments_2              same_segments;
124   Side_of_bounded_square_2_Type    side_of_bounded_square;
125   Bisector_Linf_Type               bisector_linf;
126   Compare_x_2_Sites_Type           scmpx;
127   Compare_y_2_Sites_Type           scmpy;
128 
129 private:
130   //--------------------------------------------------------------------------
131 
132   inline void
compute_v_if_not_computed()133   compute_v_if_not_computed() const {
134     if (! is_v_computed) {
135       compute_vertex(p_, q_, r_);
136       is_v_computed = true;
137     }
138   }
139 
140   void
compute_ppp(const Site_2 & sp,const Site_2 & sq,const Site_2 & sr)141   compute_ppp(const Site_2& sp, const Site_2& sq, const Site_2& sr)
142   const
143   {
144     CGAL_precondition( sp.is_point() && sq.is_point() &&
145                        sr.is_point() );
146 
147     CGAL_SDG_DEBUG(std::cout << "debug vring entering compute_ppp"
148         << std::endl;);
149 
150     Point_2 p = sp.point(), q = sq.point(), r = sr.point();
151 
152     compute_ppp(p, q, r);
153   }
154 
155   void
compute_ppp(const Point_2 & p,const Point_2 & q,const Point_2 & r)156   compute_ppp(const Point_2& p, const Point_2& q, const Point_2& r)
157   const
158   {
159     RT x_min, x_max, y_min, y_max;
160     RT two_x_center, two_y_center;
161     RT two(2);
162 
163     bool is_set_x_center(false);
164     bool is_set_y_center(false);
165     bool is_set_x_max(false);
166     bool is_set_y_max(false);
167     bool is_set_x_min(false);
168     bool is_set_y_min(false);
169 
170     Comparison_result cmpxqp = CGAL::compare(q.x(), p.x());
171 
172     if (cmpxqp == SMALLER) { // q.x() < p.x()
173       x_min = q.x();
174       x_max = p.x();
175     } else if (cmpxqp == LARGER) { // q.x() > p.x()
176       x_min = p.x();
177       x_max = q.x();
178     } else { // q.x() = p.x()
179       x_min = p.x();
180       x_max = p.x();
181       two_y_center = p.y() + q.y();
182       is_set_y_center = true;
183 
184       CGAL_SDG_DEBUG(std::cout << "debug set " <<
185         " py, qy =" << p.y() << ' ' << q.y() <<
186         " two_y_center=" << two_y_center << std::endl;);
187 
188       Comparison_result cmpxrothers = CGAL::compare(r.x(), p.x());
189       if (cmpxrothers == SMALLER) {
190         CGAL_SDG_DEBUG(std::cout << "debug r is left of p, q" << std::endl;);
191         Comparison_result cmpyrp = CGAL::compare(r.y(), p.y());
192         Comparison_result cmpyrq = CGAL::compare(r.y(), q.y());
193         if (((cmpyrp == LARGER)  && (cmpyrq == LARGER)) ||
194             ((cmpyrp == SMALLER) && (cmpyrq == SMALLER))
195            ) {
196           // do fix
197           if (cmpyrp == LARGER) {
198             y_min = two_y_center - r.y();
199             is_set_y_min = true;
200             CGAL_SDG_DEBUG(std::cout << "debug set y_min=" <<
201               y_min << std::endl;);
202           } else {
203             y_max = two_y_center - r.y();
204             is_set_y_max = true;
205             CGAL_SDG_DEBUG(std::cout << "debug set y_max=" <<
206               y_max << std::endl;);
207           }
208         }
209       } else if (cmpxrothers == LARGER) {
210         CGAL_SDG_DEBUG(std::cout << "debug r is right of p, q" << std::endl;);
211         Comparison_result cmpyrp = CGAL::compare(r.y(), p.y());
212         Comparison_result cmpyrq = CGAL::compare(r.y(), q.y());
213         if (((cmpyrp == LARGER)  && (cmpyrq == LARGER)) ||
214             ((cmpyrp == SMALLER) && (cmpyrq == SMALLER))
215            ) {
216           // do fix
217           if (cmpyrp == LARGER) {
218             y_min = two_y_center - r.y();
219             is_set_y_min = true;
220             CGAL_SDG_DEBUG(std::cout << "debug set y_min=" <<
221               y_min << std::endl;);
222           } else {
223             y_max = two_y_center - r.y();
224             is_set_y_max = true;
225             CGAL_SDG_DEBUG(std::cout << "debug set y_max=" <<
226               y_max << std::endl;);
227           }
228         }
229       } else {
230         // not possible
231       }
232     }
233 
234     Comparison_result cmpyqp = CGAL::compare(q.y(), p.y());
235 
236     if (cmpyqp == SMALLER) { // q.y() < p.y()
237       if (! is_set_y_min) {
238         y_min = q.y();
239       }
240       if (! is_set_y_max) {
241         y_max = p.y();
242       }
243     } else if (cmpyqp == LARGER) { // q.y() > p.y()
244       if (! is_set_y_min) {
245         y_min = p.y();
246       }
247       if (! is_set_y_max) {
248         y_max = q.y();
249       }
250     } else { //  q.y() = p.y()
251       if (! is_set_y_min) {
252         y_min = p.y();
253       }
254       if (! is_set_y_max) {
255         y_max = p.y();
256       }
257       two_x_center = p.x() + q.x();
258       is_set_x_center = true;
259 
260       Comparison_result cmpyrothers = CGAL::compare(r.y(), p.y());
261       if (cmpyrothers == SMALLER) {
262         Comparison_result cmpxrp = CGAL::compare(r.x(), p.x());
263         Comparison_result cmpxrq = CGAL::compare(r.x(), q.x());
264         if (((cmpxrp == LARGER)  && (cmpxrq == LARGER)) ||
265             ((cmpxrp == SMALLER) && (cmpxrq == SMALLER))
266            ) {
267           // do fix
268           if (cmpxrp == LARGER) {
269             x_min = two_x_center - r.x();
270             is_set_x_min = true;
271           } else {
272             x_max = two_x_center - r.x();
273             is_set_x_max = true;
274           }
275         }
276       } else if (cmpyrothers == LARGER) {
277         Comparison_result cmpxrp = CGAL::compare(r.x(), p.x());
278         Comparison_result cmpxrq = CGAL::compare(r.x(), q.x());
279         if (((cmpxrp == LARGER)  && (cmpxrq == LARGER)) ||
280             ((cmpxrp == SMALLER) && (cmpxrq == SMALLER))
281            ) {
282           // do fix
283           if (cmpxrp == LARGER) {
284             x_min = two_x_center - r.x();
285             is_set_x_min = true;
286           } else {
287             x_max = two_x_center - r.x();
288             is_set_x_max = true;
289           }
290         }
291       } else {
292         // not possible
293       }
294 
295     }
296 
297     Comparison_result cmpxrmin = CGAL::compare(r.x(), x_min);
298     Comparison_result cmpxrmax = CGAL::compare(r.x(), x_max);
299     if (cmpxrmin == SMALLER) {
300         // here r.x() < x_min <= x_max
301         if (! is_set_x_min) {
302           x_min = r.x();
303         }
304     } else if (cmpxrmin == LARGER) {
305       // here r.x() > x_min
306       if (cmpxrmax == LARGER) {
307         // here x_min <= x_max < r.x()
308         if (! is_set_x_max) {
309           x_max = r.x();
310         }
311       } else if (cmpxrmax == SMALLER) {
312         // x_min < r.x() < x_max
313         // do nothing
314       } else { // r.x() = x_max
315         // r.x() = p.x() || r.x() = q.x()
316         if (CGAL::compare(r.x(), p.x()) == EQUAL) {
317           two_y_center = p.y() + r.y();
318           //Comparison_result cmpyqp = CGAL::compare(q.y(),p.y());
319           Comparison_result cmpyqr = CGAL::compare(q.y(),r.y());
320           if ((cmpyqp == LARGER) && (cmpyqr == LARGER)) {
321             y_min = two_y_center - q.y();
322             is_set_y_min = true;
323           }
324           if ((cmpyqp == SMALLER) && (cmpyqr == SMALLER)) {
325             y_max = two_y_center - q.y();
326             is_set_y_max = true;
327           }
328         } else {
329           two_y_center = q.y() + r.y();
330           Comparison_result cmpypq = CGAL::compare(p.y(),q.y());
331           Comparison_result cmpypr = CGAL::compare(p.y(),r.y());
332           if ((cmpypq == LARGER) && (cmpypr == LARGER)) {
333             y_min = two_y_center - p.y();
334             is_set_y_min = true;
335           }
336           if ((cmpypq == SMALLER) && (cmpypr == SMALLER)) {
337             y_max = two_y_center - p.y();
338             is_set_y_max = true;
339           }
340         }
341         is_set_y_center = true;
342       }
343     } else {
344       // here r.x() = x_min
345       // r.x() = p.x() || r.x() = q.x()
346       if (CGAL::compare(r.x(), p.x()) == EQUAL) {
347         CGAL_SDG_DEBUG(std::cout << "debug r.x = p.x" << std::endl;);
348         // r.x() = p.x()
349         two_y_center = p.y() + r.y();
350         //Comparison_result cmpyqp = CGAL::compare(q.y(),p.y());
351         Comparison_result cmpyqr = CGAL::compare(q.y(),r.y());
352         if ((cmpyqp == LARGER) && (cmpyqr == LARGER)) {
353           CGAL_SDG_DEBUG(std::cout << "debug q is above p, r" << std::endl;);
354           y_min = two_y_center - q.y();
355           is_set_y_min = true;
356         }
357         if ((cmpyqp == SMALLER) && (cmpyqr == SMALLER)) {
358           CGAL_SDG_DEBUG(std::cout << "debug q is below p, r" << std::endl;);
359           y_max = two_y_center - q.y();
360           is_set_y_max = true;
361         }
362       } else {
363         // r.x() = q.x()
364         CGAL_SDG_DEBUG(std::cout << "debug r.x = q.x" << std::endl;);
365         two_y_center = q.y() + r.y();
366         Comparison_result cmpypq = CGAL::compare(p.y(),q.y());
367         Comparison_result cmpypr = CGAL::compare(p.y(),r.y());
368         if ((cmpypq == LARGER) && (cmpypr == LARGER)) {
369           CGAL_SDG_DEBUG(std::cout << "debug p is above q, r" << std::endl;);
370           y_min = two_y_center - p.y();
371           is_set_y_min = true;
372         }
373         if ((cmpypq == SMALLER) && (cmpypr == SMALLER)) {
374           CGAL_SDG_DEBUG(std::cout << "debug p is below q, r" << std::endl;);
375           y_max = two_y_center - p.y();
376           is_set_y_max = true;
377         }
378       }
379       is_set_y_center = true;
380     }
381 
382     Comparison_result cmpyrmin = CGAL::compare(r.y(), y_min);
383     Comparison_result cmpyrmax = CGAL::compare(r.y(), y_max);
384     if (cmpyrmin == SMALLER) {
385       // here r.y() < y_min <= y_max
386       if (! is_set_y_min) {
387         y_min = r.y();
388       }
389     } else if (cmpyrmin == LARGER) {
390       // here r.y() > y_min
391       if (cmpyrmax == LARGER) {
392         // here y_min <= y_max < r.y()
393         if (! is_set_y_max) {
394           y_max = r.y();
395         }
396       } else if (cmpyrmax == SMALLER) {
397         // y_min < r.y() < y_max
398         // do nothing
399       } else { // r.y() = y_max
400         // r.y() = p.y() || r.y() = q.y()
401         if (CGAL::compare(r.y(), p.y()) == EQUAL) {
402           two_x_center = p.x() + r.x();
403           //Comparison_result cmpxqp = CGAL::compare(q.x(),p.x());
404           Comparison_result cmpxqr = CGAL::compare(q.x(),r.x());
405           if ((cmpxqp == LARGER) && (cmpxqr == LARGER)) {
406             x_min = two_x_center - q.x();
407             is_set_x_min = true;
408           }
409           if ((cmpxqp == SMALLER) && (cmpxqr == SMALLER)) {
410             x_max = two_x_center - q.x();
411             is_set_x_max = true;
412           }
413         } else {
414           two_x_center = q.x() + r.x();
415           Comparison_result cmpxpq = CGAL::compare(p.x(),q.x());
416           Comparison_result cmpxpr = CGAL::compare(p.x(),r.x());
417           if ((cmpxpq == LARGER) && (cmpxpr == LARGER)) {
418             x_min = two_x_center - p.x();
419             is_set_x_min = true;
420           }
421           if ((cmpxpq == SMALLER) && (cmpxpr == SMALLER)) {
422             x_max = two_x_center - p.x();
423             is_set_x_max = true;
424           }
425         }
426         is_set_x_center = true;
427       }
428     } else {
429       // here r.y() = y_min
430       // r.y() = p.y() || r.y() = q.y()
431       if (CGAL::compare(r.y(), p.y()) == EQUAL) {
432         two_x_center = p.x() + r.x();
433         //Comparison_result cmpxqp = CGAL::compare(q.x(),p.x());
434         Comparison_result cmpxqr = CGAL::compare(q.x(),r.x());
435         if ((cmpxqp == LARGER) && (cmpxqr == LARGER)) {
436           x_min = two_x_center - q.x();
437           is_set_x_min = true;
438         }
439         if ((cmpxqp == SMALLER) && (cmpxqr == SMALLER)) {
440           x_max = two_x_center - q.x();
441           is_set_x_max = true;
442         }
443       } else {
444         two_x_center = q.x() + r.x();
445         Comparison_result cmpxpq = CGAL::compare(p.x(),q.x());
446         Comparison_result cmpxpr = CGAL::compare(p.x(),r.x());
447         if ((cmpxpq == LARGER) && (cmpxpr == LARGER)) {
448           x_min = two_x_center - p.x();
449           is_set_x_min = true;
450         }
451         if ((cmpxpq == SMALLER) && (cmpxpr == SMALLER)) {
452           x_max = two_x_center - p.x();
453           is_set_x_max = true;
454         }
455       }
456       is_set_x_center = true;
457     }
458 
459     Comparison_result cmpsides =
460         CGAL::compare(x_max - x_min, y_max - y_min);
461 
462     // if bounding box is non-square and points are not
463     // on corners of it, then grow it to become square
464     switch(cmpsides) {
465       case SMALLER:
466         CGAL_SDG_DEBUG(std::cout
467             << "debug vring rectangle has to be made fatter" << std::endl;);
468         // make rectangle fatter
469         if (is_set_x_center) {
470           CGAL_SDG_DEBUG(std::cout
471               << "debug vring x_center already set" << std::endl;);
472           // grow in both sides
473           break;
474         }
475         // grow only if any point is inside vertical sides
476         if (((CGAL::compare(p.x(), x_min) == EQUAL)   &&
477              (CGAL::compare(p.y(), y_max) == SMALLER) &&
478              (CGAL::compare(p.y(), y_min) == LARGER)     ) ||
479             ((CGAL::compare(q.x(), x_min) == EQUAL)   &&
480              (CGAL::compare(q.y(), y_max) == SMALLER) &&
481              (CGAL::compare(q.y(), y_min) == LARGER)     ) ||
482             ((CGAL::compare(r.x(), x_min) == EQUAL)   &&
483              (CGAL::compare(r.y(), y_max) == SMALLER) &&
484              (CGAL::compare(r.y(), y_min) == LARGER)     )   )
485         { // grow rectangle to the right
486           CGAL_SDG_DEBUG(std::cout << "debug vring grow right" << std::endl;);
487           x_max = x_min + y_max - y_min;
488         } else
489         { // grow rectangle to the left
490           CGAL_SDG_DEBUG(std::cout << "debug vring grow left" << std::endl;);
491           x_min = x_max - y_max + y_min;
492         }
493         break;
494       case LARGER:
495         CGAL_SDG_DEBUG(std::cout
496             << "debug vring rectangle has to be made taller" << std::endl;);
497         // make rectangle taller
498         if (is_set_y_center) {
499           // grow in both sides
500           CGAL_SDG_DEBUG(std::cout << "y_center already set" << std::endl;);
501           break;
502         }
503         // grow only if any point is inside horizontal sides
504         if (((CGAL::compare(p.y(), y_min) == EQUAL)   &&
505              (CGAL::compare(p.x(), x_max) == SMALLER) &&
506              (CGAL::compare(p.x(), x_min) == LARGER)     ) ||
507             ((CGAL::compare(q.y(), y_min) == EQUAL)   &&
508              (CGAL::compare(q.x(), x_max) == SMALLER) &&
509              (CGAL::compare(q.x(), x_min) == LARGER)     ) ||
510             ((CGAL::compare(r.y(), y_min) == EQUAL)   &&
511              (CGAL::compare(r.x(), x_max) == SMALLER) &&
512              (CGAL::compare(r.x(), x_min) == LARGER)     )   )
513         { // grow rectangle upwards
514           CGAL_SDG_DEBUG(std::cout
515               << "debug vring grow upwards" << std::endl;);
516           y_max = y_min + x_max - x_min;
517         } else
518         { // grow rectangle downwards
519           CGAL_SDG_DEBUG(std::cout
520               << "debug vring grow downwards" << std::endl;);
521           y_min = y_max - x_max + x_min;
522         }
523         break;
524       case EQUAL:
525         // do nothing
526         break;
527     }
528 
529     ux_ = x_min + x_max;
530     uy_ = y_min + y_max;
531     uz_ = RT(2) ;
532   }
533 
534   //--------------------------------------------------------------------------
535 
536   void
compute_pss(const Site_2 & p,const Site_2 & q,const Site_2 & r)537   compute_pss(const Site_2& p, const Site_2& q, const Site_2& r)
538   const
539   {
540     CGAL_precondition( p.is_point() && q.is_segment() &&
541                        r.is_segment() );
542 
543     CGAL_SDG_DEBUG(std::cout << "debug: compute_pss entering p=" << p
544        << " q=" << q << " r=" << r << std::endl;);
545 
546     const bool pq =
547       same_points(p, q.source_site()) || same_points(p, q.target_site());
548     const bool pr =
549       same_points(p, r.source_site()) || same_points(p, r.target_site());
550 
551     if ( pq && pr ) {
552       // philaris: result should be point p
553       const Point_2 pp = p.point();
554       ux_ = pp.hx();
555       uy_ = pp.hy();
556       uz_ = pp.hw();
557       return;
558     }
559     const bool is_q_hor = is_site_horizontal(q);
560     const bool is_q_ver = is_site_vertical(q);
561     const bool is_r_hor = is_site_horizontal(r);
562     const bool is_r_ver = is_site_vertical(r);
563     const bool is_q_hv = is_q_hor || is_q_ver;
564     const bool is_r_hv = is_r_hor || is_r_ver;
565     if (is_q_hv && is_r_hv) {
566       compute_pss_both_hv(p, q, r, is_q_hor, is_r_hor, pq, pr);
567     } else {
568       if (pq || pr) {
569         compute_pss_endp(p, q, r,
570             is_q_hv, is_q_hor, pq, is_r_hv, is_r_hor, pr);
571       } else {
572         compute_pss_nonendp(p, q, r,
573             is_q_hv, is_q_hor, is_r_hv, is_r_hor);
574       }
575     }
576   }
577 
578   // PSS case when not both segments are axis-parallel and p is
579   // not an endpoint of any of the segments q and r
580   inline void
compute_pss_nonendp(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool is_q_hv,const bool is_q_hor,const bool is_r_hv,const bool is_r_hor)581   compute_pss_nonendp(const Site_2& p, const Site_2& q, const Site_2& r,
582       const bool is_q_hv, const bool is_q_hor,
583       const bool is_r_hv, const bool is_r_hor)
584   const
585   {
586     CGAL_USE(is_q_hv);
587     CGAL_USE(is_q_hor);
588     CGAL_USE(is_r_hv);
589     CGAL_USE(is_r_hor);
590     const Line_2 lq = orient_line_nonendp(p, q);
591     const Line_2 lr = orient_line_nonendp(p, r);
592     const unsigned int bq = bearing(lq);
593     const unsigned int br = bearing(lr);
594     const unsigned int bdiff = bearing_diff(bq, br);
595     CGAL_assertion( bdiff != 0 );
596     CGAL_assertion( bdiff != 7 );
597 
598     if (bdiff == 1) {
599       compute_pss_corner_and_pt(p, q, r, lq, lr, bq, br);
600     } else if (bdiff == 2) {
601       compute_pss_nonhv_consecutive(p, q, r, lq, lr, bq, br);
602     } else if ((bdiff == 3) || (bdiff == 4)) {
603       compute_pss_ortho_wedge(p, q, r, lq, lr, bq, br);
604     } else if (bdiff == 5) {
605       compute_pss_side_p_known(p, q, r, lq, lr, bq, br);
606     } else if (bdiff == 6) {
607       compute_pss_lines_side(p, lq, lr, (br+1)%8);
608     } else {
609       CGAL_assertion( false );
610     }
611     CGAL_assertion_code( is_v_computed = true );
612     CGAL_assertion( oriented_side_of_line(lq, this->point()) != ZERO );
613     CGAL_assertion( oriented_side_of_line(lr, this->point()) != ZERO );
614   }
615 
616   inline void
compute_pss_lines_side(const Site_2 & p,const Line_2 & lq,const Line_2 & lr,const Bearing bside)617   compute_pss_lines_side(const Site_2& p,
618       const Line_2& lq, const Line_2 & lr,
619       const Bearing bside)
620   const
621   {
622     CGAL_precondition(bside % 2 == 1);
623     const bool side_ver = (bside % 4 == 1);
624     const FT pcoord = (side_ver) ? p.point().x() : p.point().y();
625     const FT qcoord = coord_at(lq, pcoord, side_ver);
626     const FT rcoord = coord_at(lr, pcoord, side_ver);
627     const FT sidelen = CGAL::abs(qcoord-rcoord);
628     const int sgn = (bside < 4) ? -1 : +1;
629     const RT two(2);
630     if(side_ver) {
631       ux_ = two*pcoord + sgn*sidelen;
632       uy_ = qcoord+rcoord;
633     } else {
634       ux_ = qcoord+rcoord;
635       uy_ = two*pcoord + sgn*sidelen;
636     }
637     uz_ = two;
638   }
639 
640   inline void
compute_pss_side_p_known(const Site_2 & p,const Site_2 & q,const Site_2 & r,const Line_2 & lq,const Line_2 & lr,const Bearing bq,const Bearing br)641   compute_pss_side_p_known(
642       const Site_2& p, const Site_2& q, const Site_2& r,
643       const Line_2& lq, const Line_2 & lr,
644       const Bearing bq, const Bearing br)
645   const
646   {
647     CGAL_USE(q);
648     CGAL_USE(r);
649     CGAL_USE(bq);
650     const Bearing bside = (br + ((br % 2 == 0) ? 1 : 2)) % 8;
651     const bool l_compute_y = (bside % 4 == 1) ? true : false;
652     const FT pcoord = l_compute_y ? p.point().x() : p.point().y();
653     const FT qcoord = coord_at(lq, pcoord, l_compute_y);
654     const FT rcoord = coord_at(lr, pcoord, l_compute_y);
655     const Point_2 qcorner =
656       l_compute_y ? Point_2(pcoord, qcoord) : Point_2(qcoord, pcoord);
657     const Point_2 rcorner =
658       l_compute_y ? Point_2(pcoord, rcoord) : Point_2(rcoord, pcoord);
659     const Point_2 v =
660       center_from_same_side_corners(rcorner, qcorner, bside);
661     ux_ = v.hx();
662     uy_ = v.hy();
663     uz_ = v.hw();
664   }
665 
666   inline void
compute_pss_ortho_wedge(const Site_2 & p,const Site_2 & q,const Site_2 & r,const Line_2 & lq,const Line_2 & lr,const Bearing bq,const Bearing br)667   compute_pss_ortho_wedge(
668       const Site_2& p, const Site_2& q, const Site_2& r,
669       const Line_2& lq, const Line_2 & lr,
670       const Bearing bq, const Bearing br)
671   const
672   {
673     CGAL_USE(q);
674     CGAL_USE(r);
675     const FT xp = p.point().x();
676     const FT yp = p.point().y();
677     const bool lq_compute_y = ((bq / 2) % 2 == 0) ? false : true;
678     const FT & lq_from_p = lq_compute_y ? xp : yp;
679     const FT & lr_from_p = lq_compute_y ? yp : xp;
680     const FT qcoord = coord_at(lq, lq_from_p, lq_compute_y);
681     const FT rcoord = coord_at(lr, lr_from_p, ! lq_compute_y);
682     const FT qdist = (bq < 4) ? qcoord - lr_from_p :
683                                 lr_from_p - qcoord;
684     CGAL_assertion(CGAL::sign(qdist) == POSITIVE);
685     const FT rdist = (bq <= 1) || (bq >= 6) ? rcoord - lq_from_p :
686                                               lq_from_p - rcoord;
687     CGAL_assertion(CGAL::sign(rdist) == POSITIVE);
688     const Comparison_result cmpqr = CGAL::compare(qdist, rdist);
689     const bool q_closer = (cmpqr == SMALLER);
690     const Point_2 corner =
691       q_closer ?
692       (lq_compute_y ? Point_2(xp, qcoord) : Point_2(qcoord, yp)) :
693       (lq_compute_y ? Point_2(rcoord, yp) : Point_2(xp, rcoord)) ;
694     const Bearing bnonhv = (bq % 2 == 1) ? br : bq;
695     CGAL_assertion(bnonhv % 2 == 0);
696     const Line_2 lcorner = (bnonhv % 4 == 0)?
697         compute_neg_45_line_at(corner) :
698         compute_pos_45_line_at(corner) ;
699     const Line_2 & lother = q_closer ? lr : lq;
700     RT hx, hy, hw;
701     compute_intersection_of_lines(lother, lcorner, hx, hy, hw);
702     const Point_2 v =
703       center_from_opposite_corners(Point_2(hx, hy, hw), corner);
704     ux_ = v.hx();
705     uy_ = v.hy();
706     uz_ = v.hw();
707   }
708 
709   inline void
compute_pss_nonhv_consecutive(const Site_2 & p,const Site_2 & q,const Site_2 & r,const Line_2 & lq,const Line_2 & lr,const Bearing bq,const Bearing br)710   compute_pss_nonhv_consecutive(
711       const Site_2& p, const Site_2& q, const Site_2& r,
712       const Line_2& lq, const Line_2 & lr,
713       const Bearing bq, const Bearing br)
714   const
715   {
716     const Bearing bqr = (bq+1)%8;
717     return (bqr % 4) == 1 ?
718       compute_pss_x_consecutive(p, q, r, lq, lr, bq, br, bqr) :
719       compute_pss_y_consecutive(p, q, r, lq, lr, bq, br, bqr) ;
720   }
721 
722   inline void
compute_pss_x_consecutive(const Site_2 & p,const Site_2 & q,const Site_2 & r,const Line_2 & lq,const Line_2 & lr,const Bearing bq,const Bearing br,const Bearing bqr)723   compute_pss_x_consecutive(
724       const Site_2& p, const Site_2& q, const Site_2& r,
725       const Line_2& lq, const Line_2 & lr,
726       const Bearing bq, const Bearing br,
727       const Bearing bqr)
728   const
729   {
730     CGAL_precondition((bqr == 1) || (bqr == 5));
731     CGAL_USE(q);
732     CGAL_USE(r);
733     CGAL_USE(bq);
734     CGAL_USE(br);
735     const FT xp = p.point().x();
736     const FT x =
737       (lr.b()*(lq.b()*xp + lq.c()) - lq.b()*lr.c()) /
738       (lr.b()*(lq.b() -lq.a()) + lq.b()*lr.a()) ;
739     const FT yq = (lq.a()*x + lq.c())/(-lq.b());
740     const FT yr = (lr.a()*x + lr.c())/(-lr.b());
741 
742     const FT yp = p.point().y();
743     if (CGAL::compare(yp, yq) == ((bqr == 1) ? SMALLER : LARGER)) {
744       // p close to q
745       const FT xs = coord_at(lq, yp, false);
746       const FT ys = coord_at(lr, xs, true);
747       ux_ = RT(2)*xs + (yp - ys);
748       uy_ = yp + ys;
749     } else if (CGAL::compare(yp, yr) == ((bqr == 1) ? LARGER : SMALLER)) {
750       // p close to r
751       const FT xs = coord_at(lr, yp, false);
752       const FT ys = coord_at(lq, xs, true);
753       ux_ = RT(2)*xs + (ys - yp);
754       uy_ = yp + ys;
755     } else {
756       // p on opposite side of two lines (or on its corners)
757       ux_ = xp + x;
758       uy_ = yq + yr;
759     }
760     uz_ = RT(2);
761   }
762 
763   inline void
compute_pss_y_consecutive(const Site_2 & p,const Site_2 & q,const Site_2 & r,const Line_2 & lq,const Line_2 & lr,const Bearing bq,const Bearing br,const Bearing bqr)764   compute_pss_y_consecutive(
765       const Site_2& p, const Site_2& q, const Site_2& r,
766       const Line_2& lq, const Line_2 & lr,
767       const Bearing bq, const Bearing br,
768       const Bearing bqr)
769   const
770   {
771     CGAL_precondition((bqr == 3) || (bqr == 7));
772     CGAL_USE(q);
773     CGAL_USE(r);
774     CGAL_USE(bq);
775     CGAL_USE(br);
776     const FT yp = p.point().y();
777     const FT y =
778       (lr.a()*(lq.a()*yp - lq.c()) + lq.a()*lr.c()) /
779       (lr.a()*(lq.a() + lq.b()) - lq.a()*lr.b()) ;
780     const FT xq = (lq.b()*y + lq.c())/(-lq.a());
781     const FT xr = (lr.b()*y + lr.c())/(-lr.a());
782 
783     const FT xp = p.point().x();
784     if (CGAL::compare(xp, xq) == ((bqr == 3) ? LARGER : SMALLER)) {
785       // p close to q
786       const FT ys = coord_at(lq, xp, true);
787       const FT xs = coord_at(lr, ys, false);
788       ux_ = xp + xs;
789       uy_ = RT(2)*ys + (xs - xp);
790     } else if (CGAL::compare(xp, xr) == ((bqr == 3) ? SMALLER : LARGER)) {
791       // p close to r
792       const FT ys = coord_at(lr, xp, true);
793       const FT xs = coord_at(lq, ys, false);
794       ux_ = xp + xs;
795       uy_ = RT(2)*ys + (xp - xs);
796     } else {
797       // p on opposite side of two lines (or on its corners)
798       ux_ = xq + xr;
799       uy_ = yp + y;
800     }
801     uz_ = RT(2);
802   }
803 
804   inline void
compute_pss_corner_and_pt(const Site_2 & p,const Site_2 & q,const Site_2 & r,const Line_2 & lq,const Line_2 & lr,const Bearing bq,const Bearing br)805   compute_pss_corner_and_pt(const Site_2& p, const Site_2& q, const Site_2& r,
806       const Line_2& lq, const Line_2 & lr,
807       const Bearing bq, const Bearing br)
808   const
809   {
810     const Bearing cb = (bq % 2 == 0) ? bq : br;
811     Point_2 v;
812     const bool is_qsrc_r = is_endpoint_of(q.source_site(), r);
813     if (is_qsrc_r) {
814       v = center_from_corner_and_pt(q.source(), cb, p.point());
815     } else {
816       const bool is_qtrg_r = is_endpoint_of(q.target_site(), r);
817       if (is_qtrg_r) {
818         v = center_from_corner_and_pt(q.target(), cb, p.point());
819       } else {
820         RT cx, cy, cw;
821         compute_intersection_of_lines(lq, lr, cx, cy, cw);
822         v = center_from_corner_and_pt(Point_2(cx, cy, cw), cb, p.point());
823       }
824     }
825     ux_ = v.hx();
826     uy_ = v.hy();
827     uz_ = v.hw();
828   }
829 
830 
831   // PSS case when not both segments are axis-parallel and p is
832   // an endpoint of one of the segments
833   inline void
compute_pss_endp(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool is_q_hv,const bool is_q_hor,const bool pq,const bool is_r_hv,const bool is_r_hor,const bool pr)834   compute_pss_endp(const Site_2& p, const Site_2& q, const Site_2& r,
835       const bool is_q_hv, const bool is_q_hor, const bool pq,
836       const bool is_r_hv, const bool is_r_hor, const bool pr)
837   const
838   {
839     CGAL_precondition(pq || pr);
840     CGAL_USE(pr);
841     const Line_2 lendp = orient_line_endp(p, (pq ? q : r), pq);
842     const Line_2 lnon = orient_line_nonendp(p, (pq ? r : q));
843     CGAL_SDG_DEBUG(std::cout << "debug compute_pss_endp lendp="
844         << lendp.a() << ' ' << lendp.b() << ' ' << lendp.c() << ' '
845         << std::endl ; );
846     CGAL_SDG_DEBUG(std::cout << "debug compute_pss_endp lnon="
847         << lnon.a() << ' ' << lnon.b() << ' ' << lnon.c() << ' '
848         << std::endl ; );
849     const Line_2 llbis = bisector_linf_line(
850         (pq ? q : r), (pq ? r : q), lendp, lnon);
851     Line_2 lperp;
852     const bool is_hv = pq ? is_q_hv : is_r_hv;
853     if (is_hv) {
854       const bool is_hor = pq ? is_q_hor : is_r_hor;
855       lperp = is_hor ? compute_ver_line_at(p.point()) :
856                        compute_hor_line_at(p.point()) ;
857     } else {
858       lperp = has_positive_slope(pq ? q : r) ?
859         compute_neg_45_line_at(p.point()) :
860         compute_pos_45_line_at(p.point()) ;
861     }
862     CGAL_SDG_DEBUG(std::cout << "debug compute_pss_endp llbis="
863         << llbis.a() << ' ' << llbis.b() << ' ' << llbis.c() << ' '
864         << std::endl ; );
865     CGAL_SDG_DEBUG(std::cout << "debug compute_pss_endp lperp="
866         << lperp.a() << ' ' << lperp.b() << ' ' << lperp.c() << ' '
867         << std::endl ; );
868     compute_intersection_of_lines(llbis, lperp, ux_, uy_, uz_);
869     CGAL_SDG_DEBUG(std::cout << "debug compute_pss_endp vertex="
870         << ux_ << ' ' << uy_ << ' ' << uz_ << ' '
871         << Point_2(ux_, uy_, uz_) << std::endl ; );
872     CGAL_assertion_code( is_v_computed = true );
873     CGAL_assertion( oriented_side_of_line(lendp, this->point()) != ZERO );
874     CGAL_assertion( oriented_side_of_line(lnon, this->point()) != ZERO );
875   }
876 
877   // both segments are axis-parallel
878   inline void
compute_pss_both_hv(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool is_q_hor,const bool is_r_hor,const bool pq,const bool pr)879   compute_pss_both_hv(const Site_2& p, const Site_2& q, const Site_2& r,
880       const bool is_q_hor, const bool is_r_hor,
881       const bool pq, const bool pr)
882   const
883   {
884     CGAL_precondition(! (pq && pr));
885     if (is_q_hor == is_r_hor) {
886       // parallel segments
887       const RT q_coord = hvseg_coord(q, is_q_hor);
888       const RT r_coord = hvseg_coord(r, is_r_hor);
889       RT & upar = is_q_hor ? ux_ : uy_;
890       RT & uort = is_q_hor ? uy_ : ux_;
891       upar = RT(2)*(is_q_hor ? p.point().x() : p.point().y())
892         + (( pq || pr ) ? RT(0) :
893                           RT(is_q_hor ? +1 : -1)*(r_coord - q_coord));
894       uort = q_coord + r_coord;
895       uz_ = RT(2);
896     } else {
897       return compute_pss_both_hv_nonpar(
898           p, q, r, is_q_hor, is_r_hor, pq, pr);
899     }
900   }
901 
902   // one segment is horizontal and the other is vertical
903   inline void
compute_pss_both_hv_nonpar(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool is_q_hor,const bool is_r_hor,const bool pq,const bool pr)904   compute_pss_both_hv_nonpar(
905       const Site_2& p, const Site_2& q, const Site_2& r,
906       const bool is_q_hor, const bool is_r_hor,
907       const bool pq, const bool pr)
908   const
909   {
910     CGAL_precondition(is_q_hor != is_r_hor);
911     if (pq || pr) {
912       const RT q_coord = hvseg_coord(q, is_q_hor);
913       const RT r_coord = hvseg_coord(r, is_r_hor);
914       const bool is_touched_hor = pq ? is_q_hor : is_r_hor;
915       const RT coord_c = is_touched_hor ? p.point().x() : p.point().y();
916       const RT radius = CGAL::abs(coord_c - (pq ? r_coord : q_coord));
917       RT & upar = is_touched_hor ? ux_ : uy_;
918       RT & uort = is_touched_hor ? uy_ : ux_;
919       const Site_2 & sother =
920         pq ? (same_points(p, q.source_site()) ?
921               q.target_site() : q.source_site()) :
922              (same_points(p, r.source_site()) ?
923               r.target_site() : r.source_site());
924       const bool test = is_touched_hor ?
925         (scmpx(p, sother) == LARGER) : (scmpy(p, sother) == SMALLER);
926       const RT sgn = RT( (pq ? +1: -1)* (test ? -1 : +1) );
927       upar = coord_c;
928       uort = (pq ? q_coord : r_coord) + sgn*radius;
929       uz_ = RT(1);
930       CGAL_SDG_DEBUG(std::cout << "debug: vring compute_pss vv="
931           << Point_2(ux_, uy_, uz_) << " radius=" << radius << std::endl;);
932       CGAL_assertion_code( const Point_2 pother = sother.point() );
933       CGAL_assertion(pq ?
934           CGAL::left_turn(Point_2(ux_,uy_,uz_), p.point(), pother) :
935           CGAL::left_turn(pother, p.point(), Point_2(ux_,uy_,uz_)) );
936       return;
937     } else {
938       return compute_pss_both_hv_nonpar_nonendp(
939                p, q, r, is_q_hor, is_r_hor, pq, pr);
940     }
941   }
942 
943   // one segment is horizontal and the other is vertical and
944   // the point p is not an endpoint of the segments
945   inline void
compute_pss_both_hv_nonpar_nonendp(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool is_q_hor,const bool is_r_hor,const bool pq,const bool pr)946   compute_pss_both_hv_nonpar_nonendp(
947       const Site_2& p, const Site_2& q, const Site_2& r,
948       const bool is_q_hor, const bool is_r_hor,
949       const bool pq, const bool pr)
950   const
951   {
952     CGAL_precondition(! (pq || pr));
953     CGAL_USE(pq);
954     CGAL_USE(pr);
955     const RT q_coord = hvseg_coord(q, is_q_hor);
956     const RT r_coord = hvseg_coord(r, is_r_hor);
957     const RT p_coord_q = is_q_hor ? p.point().y() : p.point().x();
958     const RT p_coord_r = is_r_hor ? p.point().y() : p.point().x();
959     const RT sdistq = p_coord_q - q_coord;
960     const RT sdistr = p_coord_r - r_coord;
961     const RT distq = CGAL::abs(sdistq);
962     const RT distr = CGAL::abs(sdistr);
963     const RT & dx = is_r_hor ? distq : distr;
964     const RT & dy = is_r_hor ? distr : distq;
965     const Comparison_result cmp = CGAL::compare(dx, dy);
966     if (cmp == LARGER) {
967       ux_ = is_q_hor ? (r_coord + p_coord_r) : (q_coord + p_coord_q);
968       uy_ = is_q_hor ? (RT(2)*q_coord + (int)CGAL::sign(sdistq)*dx) :
969                        (RT(2)*r_coord + (int)CGAL::sign(sdistr)*dx) ;
970     } else if (cmp == SMALLER) {
971       uy_ = is_r_hor ? (r_coord + p_coord_r) : (q_coord + p_coord_q);
972       ux_ = is_r_hor ? (RT(2)*q_coord + (int)CGAL::sign(sdistq)*dy) :
973                        (RT(2)*r_coord + (int)CGAL::sign(sdistr)*dy) ;
974     } else {
975       ux_ = is_q_hor ? (r_coord + p_coord_r) : (q_coord + p_coord_q);
976       uy_ = is_q_hor ? (q_coord + p_coord_q) : (r_coord + p_coord_r);
977     }
978     uz_ = RT(2);
979   }
980 
981   inline void
compute_pss_bisectors(const Site_2 & p,const Site_2 & q,const Site_2 & r)982   compute_pss_bisectors(const Site_2& p, const Site_2& q, const Site_2& r)
983   const
984   {
985     const bool pq =
986       same_points(p, q.source_site()) || same_points(p, q.target_site());
987     Polychainline_2 goodbisector;
988     if (pq) {
989       goodbisector = bisector_linf(p, q);
990       CGAL_SDG_DEBUG(std::cout << "debug: vring compute_pss bpq p=" << p
991           << " q=" << q << std::endl;);
992       CGAL_SDG_DEBUG(std::cout << "debug: vring compute_pss bpq ="
993           << goodbisector << std::endl;);
994     } else {
995       goodbisector = bisector_linf(r, p);
996       CGAL_SDG_DEBUG(std::cout << "debug: vring compute_pss brp r=" << r
997           << " p=" << p << std::endl;);
998       CGAL_SDG_DEBUG(std::cout << "debug: vring compute_pss brp ="
999           << goodbisector << std::endl;);
1000     }
1001 
1002     Polychainline_2 bqr = bisector_linf(q, r);
1003     CGAL_SDG_DEBUG(std::cout << "debug: vring compute_pss bqr q=" << q
1004         << " r=" << r << std::endl;);
1005     CGAL_SDG_DEBUG(std::cout << "debug: vring compute_pss bqr ="
1006         << bqr << std::endl;);
1007 
1008     Point_2 vv = goodbisector.first_intersection_point_with(bqr);
1009     CGAL_SDG_DEBUG(std::cout
1010         << "debug: vring compute_pss vv=" << vv << std::endl;);
1011 
1012     ux_ = vv.hx();
1013     uy_ = vv.hy();
1014     uz_ = vv.hw();
1015 
1016   }
1017 
1018 
1019 
1020   //--------------------------------------------------------------------------
1021 
1022   inline void
compute_pps_endp_hv(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool p_endp_r,const bool is_r_horizontal)1023   compute_pps_endp_hv(const Site_2& p, const Site_2& q, const Site_2& r,
1024                    const bool p_endp_r, const bool is_r_horizontal)
1025   const
1026   {
1027     CGAL_USE(r);
1028     const Site_2 & A = p_endp_r ? p : q;
1029     const Site_2 & B = p_endp_r ? q : p;
1030     const RT Apar = is_r_horizontal ? A.point().x() : A.point().y();
1031     const RT Aort = is_r_horizontal ? A.point().y() : A.point().x();
1032     const RT Bpar = is_r_horizontal ? B.point().x() : B.point().y();
1033     const RT Bort = is_r_horizontal ? B.point().y() : B.point().x();
1034     const RT dpar = Apar - Bpar;
1035     const RT dort = Aort - Bort;
1036     const RT absdpar = CGAL::abs(dpar);
1037 
1038     RT & upar = is_r_horizontal ? ux_ : uy_;
1039     RT & uort = is_r_horizontal ? uy_ : ux_;
1040 
1041     if (2*absdpar < CGAL::abs(dort)) {
1042       upar = RT(2)*Apar;
1043       uort = RT(2)*Aort - dort;
1044       uz_ = RT(2);
1045     } else {
1046       upar = Apar;
1047       uort = Aort - (int)CGAL::sign(dort)*absdpar;
1048       uz_ = RT(1);
1049     }
1050   }
1051 
1052   inline void
compute_pps_nonendp_hv_samecoord(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool is_r_horizontal)1053   compute_pps_nonendp_hv_samecoord(
1054       const Site_2& p, const Site_2& q, const Site_2& r,
1055       const bool is_r_horizontal)
1056   const
1057   {
1058     const RT ppar = is_r_horizontal ? p.point().x() : p.point().y();
1059     const RT port = is_r_horizontal ? p.point().y() : p.point().x();
1060     const RT qort = is_r_horizontal ? q.point().y() : q.point().x();
1061     RT & upar = is_r_horizontal ? ux_ : uy_;
1062     RT & uort = is_r_horizontal ? uy_ : ux_;
1063     const RT segort = (is_r_horizontal)?
1064       horseg_y_coord(r) : verseg_x_coord(r);
1065     const RT sumort = port + qort;
1066     uort = sumort;
1067     const int vhsign = is_r_horizontal ? +1 : -1;
1068     const int distsign = CGAL::abs(segort-qort) < CGAL::abs(segort-port) ?
1069       +1: -1;
1070     upar = RT(2)*ppar - RT(vhsign*distsign)*(RT(2)*segort-sumort);
1071     uz_ = RT(2);
1072   }
1073 
1074   /* compute pps vertex when the points p, q are not endpoints of
1075    * segment r, and when segment r is axis-parallel
1076    */
1077   inline void
compute_pps_nonendp_hv(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool is_r_horizontal)1078   compute_pps_nonendp_hv(const Site_2& p, const Site_2& q, const Site_2& r,
1079                          const bool is_r_horizontal)
1080   const
1081   {
1082     if ((is_r_horizontal     && (scmpx(p, q) == EQUAL)) ||
1083         ((! is_r_horizontal) && (scmpy(p, q) == EQUAL))   ) {
1084       return compute_pps_nonendp_hv_samecoord(p, q, r, is_r_horizontal);
1085     }
1086 
1087     // here, the segment is axis-parallel and the two points:
1088     // either: are not sharing a coordinate
1089     // or: they share a coordinate AND
1090     //     the line connecting the two points is parallel to the segment
1091     const Point_2 pp = p.point();
1092     const Point_2 qq = q.point();
1093     Line_2 l = compute_supporting_line(r);
1094     if (oriented_side_of_line(l, pp) == NEGATIVE) {
1095       l = opposite_line(l);
1096     }
1097     CGAL_assertion(oriented_side_of_line(l, pp) == POSITIVE);
1098     CGAL_assertion(oriented_side_of_line(l, qq) == POSITIVE);
1099 
1100     const Comparison_result perpcomp =
1101       is_r_horizontal ? scmpy(p, q) : scmpx(p, q);
1102 
1103     const RT coordr = hvseg_coord(r, is_r_horizontal);
1104 
1105     RT & upar = is_r_horizontal ? ux_ : uy_;
1106     RT & uort = is_r_horizontal ? uy_ : ux_;
1107 
1108     if (perpcomp == EQUAL) {
1109       const RT pqdist = is_r_horizontal ?
1110         CGAL::abs(pp.x()-qq.x()) : CGAL::abs(pp.y()-qq.y());
1111       const RT signrdist = (is_r_horizontal ? pp.y() : pp.x()) - coordr;
1112       Comparison_result comp = CGAL::compare(pqdist, CGAL::abs(signrdist));
1113       upar = is_r_horizontal ? pp.x() + qq.x() : pp.y() + qq.y();
1114       if (comp == LARGER) {
1115         uort = RT(2)*coordr + (int)CGAL::sign(signrdist)*pqdist;
1116       } else {
1117         uort = coordr + (is_r_horizontal ? pp.y() : pp.x());
1118       }
1119       uz_ = RT(2);
1120       return;
1121     }
1122     // here, perpcomp is not EQUAL
1123     const Sign signla = CGAL::sign(l.a());
1124     const Sign signlb = CGAL::sign(l.b());
1125     const Sign & testsign = is_r_horizontal ? signlb : signla;
1126     CGAL_assertion(testsign != ZERO);
1127     const Point_2 & farp =
1128       testsign == POSITIVE ?
1129         (perpcomp == SMALLER ? qq : pp) :
1130         (perpcomp == SMALLER ? pp : qq) ;
1131     CGAL_assertion(Base::compare_linf_distances_to_line(
1132           l, farp == pp ? qq : pp , farp) == SMALLER);
1133     const RT pqdist = (CGAL::max)(
1134       CGAL::abs(pp.x()-qq.x()), CGAL::abs(pp.y()-qq.y()));
1135 
1136     const RT sdistf = (is_r_horizontal ? farp.y() : farp.x()) - coordr;
1137     CGAL_assertion(CGAL::sign(sdistf) == testsign);
1138 
1139     if (CGAL::compare(CGAL::abs(sdistf), pqdist) == LARGER) {
1140       const bool is_p_farthest = farp == pp;
1141       const Point_2 & closep = (is_p_farthest)? qq : pp;
1142       upar = RT(2)* (is_r_horizontal ? closep.x() : closep.y()) +
1143         (is_r_horizontal? -1 : +1) * (is_p_farthest? -1 : +1) * sdistf;
1144       uort = RT(2)*coordr + sdistf;
1145     } else {
1146       upar = is_r_horizontal ? pp.x() + qq.x() : pp.y() + qq.y();
1147       uort = RT(2)*coordr + (int)CGAL::sign(sdistf)*pqdist;
1148     }
1149     uz_ = RT(2);
1150     return;
1151   }
1152 
1153   inline void
compute_pps_endp_slope(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool p_endp_r,const bool pos_slope)1154   compute_pps_endp_slope(const Site_2& p, const Site_2& q, const Site_2& r,
1155                    const bool p_endp_r, const bool pos_slope)
1156   const
1157   {
1158     CGAL_USE(r);
1159     const Site_2 & A = p_endp_r ? p : q;
1160     const Site_2 & B = p_endp_r ? q : p;
1161     const RT Ax = A.point().x();
1162     const RT Ay = A.point().y();
1163     const RT Bx = B.point().x();
1164     const RT By = B.point().y();
1165     const RT dx = Ax - Bx;
1166     const RT dy = Ay - By;
1167     const RT absdx = CGAL::abs(dx);
1168     const RT absdy = CGAL::abs(dy);
1169 
1170     if (absdx > absdy) {
1171       ux_ = RT(2)*Ax - dx;
1172       uy_ = RT(2)*Ay - RT(pos_slope? -1 : +1)*dx;
1173     } else {
1174       ux_ = RT(2)*Ax - RT(pos_slope? -1 : +1)*dy;
1175       uy_ = RT(2)*Ay - dy;
1176     }
1177     uz_ = RT(2);
1178   }
1179 
1180   inline void
compute_pps_endp(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool p_endp_r)1181   compute_pps_endp(const Site_2& p, const Site_2& q, const Site_2& r,
1182                    const bool p_endp_r)
1183   const
1184   {
1185     const bool is_r_horizontal = is_site_horizontal(r);
1186     if (is_r_horizontal || is_site_vertical(r)) {
1187       return compute_pps_endp_hv(p, q, r, p_endp_r, is_r_horizontal);
1188     } else {
1189       const bool pos_slope = has_positive_slope(r);
1190       return compute_pps_endp_slope(p, q, r, p_endp_r, pos_slope);
1191     }
1192   }
1193 
1194   /* compute pps vertex when the points p, q are not endpoints of
1195    * segment r, when segment r is not axis-parallel, and when the
1196    * two points p and q do not share any coordinate
1197    */
1198   inline void
compute_pps_nonendp_nonhv_nonsamec(const Site_2 & p,const Site_2 & q,const Site_2 & r)1199   compute_pps_nonendp_nonhv_nonsamec
1200   (const Site_2& p, const Site_2& q, const Site_2& r)
1201   const
1202   {
1203     Line_2 l = compute_supporting_line(r);
1204     if (oriented_side_of_line(l, p.point()) == NEGATIVE) {
1205       l = opposite_line(l);
1206     }
1207     CGAL_assertion(oriented_side_of_line(l, p.point()) == POSITIVE);
1208     CGAL_assertion(oriented_side_of_line(l, q.point()) == POSITIVE);
1209 
1210     const bool pos_slope = has_positive_slope(r);
1211     const Comparison_result first_comp =
1212       (pos_slope) ? scmpy(p, q) : scmpx(p, q);
1213     const Comparison_result second_comp =
1214       (pos_slope) ? scmpx(p, q) : scmpy(p, q);
1215 
1216     const Sign signla = CGAL::sign(l.a());
1217     const Sign signlb = CGAL::sign(l.b());
1218     const Comparison_result first_value =
1219       (signlb == POSITIVE)? SMALLER : LARGER;
1220 
1221     const Comparison_result second_value =
1222       (signla == NEGATIVE)? SMALLER : LARGER;
1223 
1224     if (first_comp == first_value) {
1225       const RT pcoord = pos_slope ? p.point().x() : p.point().y();
1226       const RT lineval = coord_at(l, pcoord, pos_slope);
1227       const Point_2 corner = pos_slope?
1228         Point_2(pcoord, lineval) : Point_2(lineval, pcoord);
1229       const RT sidelen = (CGAL::max)(CGAL::abs(corner.x() - q.point().x()),
1230                                      CGAL::abs(corner.y() - q.point().y()));
1231       ux_ = RT(2)*corner.x() + (int)signla*sidelen;
1232       uy_ = RT(2)*corner.y() + (int)signlb*sidelen;
1233       uz_ = RT(2);
1234       return;
1235     }
1236     if (second_comp == second_value) {
1237       const RT qcoord = pos_slope ? q.point().y() : q.point().x();
1238       const RT lineval = coord_at(l, qcoord, ! pos_slope);
1239       const Point_2 corner = pos_slope?
1240         Point_2(lineval, qcoord) : Point_2(qcoord, lineval);
1241       const RT sidelen = (CGAL::max)(CGAL::abs(corner.x() - p.point().x()),
1242                                      CGAL::abs(corner.y() - p.point().y()));
1243       ux_ = RT(2)*corner.x() + (int)signla*sidelen;
1244       uy_ = RT(2)*corner.y() + (int)signlb*sidelen;
1245       uz_ = RT(2);
1246       return;
1247     }
1248 
1249     CGAL_assertion((first_comp  == -first_value ) &&
1250                    (second_comp == -second_value)    );
1251 
1252     const RT px = p.point().x();
1253     const RT py = p.point().y();
1254     const RT qx = q.point().x();
1255     const RT qy = q.point().y();
1256     const RT pqdist = (CGAL::max)(CGAL::abs(px - qx), CGAL::abs(py - qy));
1257 
1258     CGAL_SDG_DEBUG(std::cout
1259         << "debug: vring pqdist=" << pqdist << std::endl;);
1260 
1261     const RT & pcoord = pos_slope ? px : py;
1262     const RT plineval = coord_at(l, pcoord, pos_slope);
1263     const RT & pothercoord = pos_slope ? py : px;
1264     const RT plen = CGAL::abs(plineval -  pothercoord);
1265     CGAL_SDG_DEBUG(std::cout
1266         << "debug: vring plen=" << plen << std::endl;);
1267     if (CGAL::compare(pqdist, plen) != SMALLER) {
1268       // here, appropriate projection of p on supporting line of segment r
1269       // is shorter than Linf p, q distance
1270       const Point_2 corner = pos_slope?
1271         Point_2(pcoord, plineval) : Point_2(plineval, pcoord);
1272       ux_ = RT(2)*corner.x() + (int)signla*pqdist;
1273       uy_ = RT(2)*corner.y() + (int)signlb*pqdist;
1274       uz_ = RT(2);
1275       return;
1276     }
1277 
1278     const RT & qcoord = pos_slope ? qy : qx;
1279     const RT qlineval = coord_at(l, qcoord, ! pos_slope);
1280     const RT & qothercoord = pos_slope ? qx : qy;
1281     const RT qlen = CGAL::abs(qlineval -  qothercoord);
1282     CGAL_SDG_DEBUG(std::cout
1283         << "debug: vring qlen=" << qlen << std::endl;);
1284     if (CGAL::compare(pqdist, qlen) != SMALLER) {
1285       // here, appropriate projection of q on supporting line of segment r
1286       // is shorter than Linf p, q distance
1287       const Point_2 corner = pos_slope?
1288         Point_2(qlineval, qcoord) : Point_2(qcoord, qlineval);
1289       ux_ = RT(2)*corner.x() + (int)signla*pqdist;
1290       uy_ = RT(2)*corner.y() + (int)signlb*pqdist;
1291       uz_ = RT(2);
1292       return;
1293     }
1294 
1295     CGAL_assertion((pqdist < plen) && (pqdist < qlen));
1296 
1297     // here, compute corner opposite of corner on line of segment r
1298     const Point_2 opposite_corner = pos_slope ?
1299       Point_2(qx, py) : Point_2(px, qy);
1300     CGAL_SDG_DEBUG(std::cout << "debug: vring opposite_corner="
1301         << opposite_corner << std::endl;);
1302 
1303     const Point_2 corner =
1304       compute_linf_projection_nonhom(l, opposite_corner);
1305 
1306     ux_ = corner.x() + opposite_corner.x();
1307     uy_ = corner.y() + opposite_corner.y();
1308     uz_ = RT(2);
1309   }
1310 
1311   inline void
compute_pps_nonendp_nonhv(const Site_2 & p,const Site_2 & q,const Site_2 & r)1312   compute_pps_nonendp_nonhv(const Site_2& p, const Site_2& q, const Site_2& r)
1313   const
1314   {
1315     const bool samexpq = scmpx(p, q) == EQUAL;
1316     const bool sameypq = (samexpq)? false : make_certain(scmpy(p, q) == EQUAL);
1317     if (! (samexpq || sameypq)) {
1318       return compute_pps_nonendp_nonhv_nonsamec(p, q, r);
1319     } else {
1320       // samexpq || sameypq
1321       CGAL_assertion(samexpq != sameypq);
1322       Line_2 l = compute_supporting_line(r);
1323       const FT common_coord = (samexpq) ? p.point().x() : p.point().y();
1324       const FT sumdiffpq = (samexpq) ?
1325         p.point().y() + q.point().y() :
1326         p.point().x() + q.point().x();
1327       const bool pos_slope = has_positive_slope(r);
1328       FT vsamecoord;
1329       if (touch_same_side(p, q, l, samexpq, pos_slope)) {
1330         vsamecoord = common_coord +
1331           (pos_slope? +1: -1)*
1332             (coord_at(l, common_coord, samexpq) - (sumdiffpq/FT(2)));
1333       } else {
1334         const FT closest_coord =
1335           (samexpq)? ((pos_slope)? q.point().y() : p.point().y()):
1336                      ((pos_slope)? p.point().x() : q.point().x());
1337         if (is_orth_dist_smaller_than_pt_dist(
1338               closest_coord, l, p, q, samexpq)) {
1339           vsamecoord =
1340             coord_at(l, closest_coord, sameypq) +
1341             (((samexpq) ? (q.point().y() - p.point().y()) :
1342                           (p.point().x() - q.point().x())  ) / FT(2)) ;
1343         } else {
1344           const Line_2 lc (
1345               (samexpq) ? 1 : (2 * ((pos_slope) ? +1 : -1)),
1346               (samexpq) ? (2 * ((pos_slope) ? +1 : -1)) : 1 ,
1347               ((pos_slope)? -1 : +1 ) * sumdiffpq - common_coord);
1348           RT hx, hy, hz;
1349           compute_intersection_of_lines(l, lc, hx, hy, hz);
1350           vsamecoord = ((samexpq ? hx/hz : hy/hz) + common_coord)/ FT(2);
1351         }
1352       }
1353       const FT vdiffcoord = sumdiffpq/FT(2);
1354       ux_ = (samexpq) ? vsamecoord : vdiffcoord;
1355       uy_ = (samexpq) ? vdiffcoord : vsamecoord;
1356       uz_ = RT(1);
1357     }
1358   }
1359 
1360   inline void
compute_pps_nonendp(const Site_2 & p,const Site_2 & q,const Site_2 & r)1361   compute_pps_nonendp(const Site_2& p, const Site_2& q, const Site_2& r)
1362   const
1363   {
1364     const bool is_r_horizontal = is_site_horizontal(r);
1365     if (is_r_horizontal || is_site_vertical(r)) {
1366       return compute_pps_nonendp_hv(p, q, r, is_r_horizontal);
1367     } else {
1368       return compute_pps_nonendp_nonhv(p, q, r);
1369     }
1370   }
1371 
1372   void
compute_pps(const Site_2 & p,const Site_2 & q,const Site_2 & r)1373   compute_pps(const Site_2& p, const Site_2& q, const Site_2& r)
1374   const
1375   {
1376     CGAL_precondition( p.is_point() && q.is_point() &&
1377                        r.is_segment() );
1378 
1379     CGAL_SDG_DEBUG(std::cout
1380         << "debug: vring compute_pps entering p=" << p
1381         << " q=" << q << " r=" << r << std::endl;);
1382 
1383     bool p_endp_r = is_endpoint_of(p, r);
1384     bool q_endp_r = is_endpoint_of(q, r);
1385 
1386     if (p_endp_r || q_endp_r) {
1387       return compute_pps_endp(p, q, r, p_endp_r);
1388     }
1389     CGAL_assertion(are_in_same_open_halfspace_of(p, q, r));
1390     return compute_pps_nonendp(p, q, r);
1391   }
1392 
1393 
1394   void
compute_pps_bisectors(const Site_2 & p,const Site_2 & q,const Site_2 & r)1395   compute_pps_bisectors(const Site_2& p, const Site_2& q, const Site_2& r)
1396   const
1397   {
1398     CGAL_precondition( p.is_point() && q.is_point() &&
1399                        r.is_segment() );
1400 
1401     CGAL_SDG_DEBUG(std::cout
1402         << "debug: vring compute_pps_bisectors entering p=" << p
1403         << " q=" << q << " r=" << r << std::endl;);
1404 
1405     bool p_endp_r = is_endpoint_of(p, r);
1406     bool q_endp_r = is_endpoint_of(q, r);
1407     Polychainline_2 bpq = bisector_linf(p, q);
1408     CGAL_SDG_DEBUG(std::cout << "debug: bpq p="
1409         << p << " q=" << q << std::endl;);
1410     CGAL_SDG_DEBUG(std::cout << "debug: bpq =" << bpq << std::endl;);
1411 
1412     CGAL_assertion(! (p_endp_r && q_endp_r));
1413 
1414     bool samexpq = (scmpx(p, q) == EQUAL);
1415     bool sameypq = (scmpy(p, q) == EQUAL);
1416 
1417     bool samecoordpq = samexpq || sameypq ;
1418 
1419     Polychainline_2 goodbisector;
1420     if (p_endp_r) {
1421       goodbisector = bisector_linf(r, p);
1422       CGAL_SDG_DEBUG(std::cout << "debug: vring brp r="
1423           << r << " p=" << p << std::endl;);
1424       CGAL_SDG_DEBUG(std::cout << "debug: vring brp res="
1425           << goodbisector << std::endl;);
1426     } else if (q_endp_r) {
1427       goodbisector = bisector_linf(q, r);
1428       CGAL_SDG_DEBUG(std::cout << "debug: vring bqr q="
1429           << q << " r=" << r << std::endl;);
1430       CGAL_SDG_DEBUG(std::cout << "debug: vring bqr res="
1431           << goodbisector << std::endl;);
1432     } else if (samecoordpq) {
1433       CGAL_SDG_DEBUG(std::cout << "debug vring PPS samecoordpq"
1434           << std::endl;);
1435 
1436       // check which of points p, q is closer to segment r
1437 
1438       bool use_bqr;
1439 
1440       Line_2 l = compute_supporting_line(r.supporting_site());
1441 
1442       if (((CGAL::sign(l.a()) == ZERO) && sameypq) ||
1443           ((CGAL::sign(l.b()) == ZERO) && samexpq)   )  {
1444         // here l is horizontal or vertical and parallel to pq;
1445         // bqr or brp are equally good
1446         use_bqr = true;
1447       } else {
1448         // here l and segment are neither hor. nor ver.
1449         Point_2 proj;
1450 
1451         // philaris: tofix with RT
1452         FT projft, pft, qft;
1453         if (samexpq) {
1454           // compute vertical projection
1455           proj = compute_vertical_projection(l, p.point());
1456           projft = proj.y();
1457           pft = p.point().y();
1458           qft = q.point().y();
1459 
1460         } else {
1461           CGAL_assertion(sameypq);
1462           // compute horizontal projection
1463           proj = compute_horizontal_projection(l, p.point());
1464           projft = proj.x();
1465           pft = p.point().x();
1466           qft = q.point().x();
1467         }
1468         Comparison_result cpq, cqproj;
1469         cpq    = CGAL::compare(pft, qft);
1470         cqproj = CGAL::compare(qft, projft);
1471         if (cpq == cqproj) {
1472           use_bqr = true;
1473         } else {
1474           use_bqr = false;
1475         }
1476       } // end of case of neither hor nor ver segment
1477 
1478       if (use_bqr) {
1479         goodbisector = bisector_linf(q, r);
1480         CGAL_SDG_DEBUG(std::cout << "debug: vring bqr q="
1481             << q << " r=" << r << std::endl;);
1482         CGAL_SDG_DEBUG(std::cout << "debug: vring bqr res="
1483             << goodbisector << std::endl;);
1484       } else {
1485         goodbisector = bisector_linf(r, p);
1486         CGAL_SDG_DEBUG(std::cout << "debug: vring brp r="
1487             << r << " p=" << p << std::endl;);
1488         CGAL_SDG_DEBUG(std::cout << "debug: vring brp res="
1489             << goodbisector << std::endl;);
1490       }
1491     } else {
1492       goodbisector = bisector_linf(q, r);
1493       CGAL_SDG_DEBUG(std::cout << "debug: vring bqr q="
1494           << q << " r=" << r << std::endl;);
1495       CGAL_SDG_DEBUG(std::cout << "debug: vring bqr res="
1496           << goodbisector << std::endl;);
1497     }
1498 
1499 
1500     Point_2 vv = bpq.first_intersection_point_with(goodbisector);
1501 
1502     CGAL_SDG_DEBUG(std::cout << "debug: PPS returns with vv="
1503         << vv << std::endl;);
1504 
1505     ux_ = vv.hx();
1506     uy_ = vv.hy();
1507     uz_ = vv.hw();
1508   }
1509 
1510   //--------------------------------------------------------------------------
1511 
1512   void
compute_sss(const Site_2 & p,const Site_2 & q,const Site_2 & r)1513   compute_sss(const Site_2& p, const Site_2& q, const Site_2& r)
1514   const
1515   {
1516     CGAL_precondition( p.is_segment() && q.is_segment() &&
1517                        r.is_segment() );
1518     const bool is_psrc_q = is_endpoint_of(p.source_site(), q);
1519     const bool is_psrc_r = is_endpoint_of(p.source_site(), r);
1520     const bool is_ptrg_q = is_endpoint_of(p.target_site(), q);
1521     const bool is_ptrg_r = is_endpoint_of(p.target_site(), r);
1522 
1523     if (is_psrc_q && is_psrc_r) {
1524       ux_ = p.source().hx();
1525       uy_ = p.source().hy();
1526       uz_ = p.source().hw();
1527     } else if (is_ptrg_q && is_ptrg_r) {
1528       ux_ = p.target().hx();
1529       uy_ = p.target().hy();
1530       uz_ = p.target().hw();
1531     } else {
1532       // here, not all segments have a common point
1533 
1534       const bool is_p_hor = is_site_horizontal(p);
1535       const bool is_q_hor = is_site_horizontal(q);
1536       const bool is_r_hor = is_site_horizontal(r);
1537 
1538       const bool is_p_hv = is_p_hor || is_site_vertical(p);
1539       const bool is_q_hv = is_q_hor || is_site_vertical(q);
1540       const bool is_r_hv = is_r_hor || is_site_vertical(r);
1541 
1542       if (is_p_hv && is_q_hv && is_r_hv) {
1543         return compute_sss_hv(p, q, r, is_p_hor, is_q_hor, is_r_hor);
1544       }
1545 
1546       Line_2 lines[3];
1547       orient_lines_linf(p, q, r, lines);
1548 
1549       compute_sss_bisectors(p, q, r, lines);
1550       CGAL_assertion_code( is_v_computed = true );
1551       CGAL_assertion( oriented_side_of_line(lines[0], this->point()) != ZERO );
1552       CGAL_assertion( oriented_side_of_line(lines[1], this->point()) != ZERO );
1553       CGAL_assertion( oriented_side_of_line(lines[2], this->point()) != ZERO );
1554     }
1555   }
1556 
1557   // SSS: all sites are axis-parallel
1558   inline void
compute_sss_hv(const Site_2 & p,const Site_2 & q,const Site_2 & r,const bool is_p_hor,const bool is_q_hor,const bool is_r_hor)1559   compute_sss_hv(const Site_2 & p, const Site_2 & q, const Site_2 & r,
1560       const bool is_p_hor, const bool is_q_hor, const bool is_r_hor)
1561   const
1562   {
1563     CGAL_precondition(! (is_p_hor && is_q_hor && is_r_hor));
1564     CGAL_precondition(is_p_hor || is_q_hor || is_r_hor);
1565     const unsigned int num_hor =
1566       (is_p_hor ? 1 : 0) + (is_q_hor ? 1 : 0) + (is_r_hor ? 1 : 0);
1567     CGAL_assertion((num_hor == 1) || (num_hor == 2));
1568     const bool are_common_hor = num_hor == 2;
1569     const bool is_odd_hor = ! are_common_hor;
1570 
1571     const Site_2 & odd = (is_odd_hor) ?
1572       (is_p_hor ? p : (is_q_hor ? q : r)) :
1573       (is_p_hor ? (is_q_hor ? r : q) : p);
1574     CGAL_assertion( (! (num_hor == 1)) || is_site_horizontal(odd) );
1575     CGAL_assertion( (! (num_hor == 2)) || is_site_vertical(odd) );
1576     const Site_2 & prev = (is_odd_hor) ?
1577       (is_p_hor ? r : (is_q_hor ? p : q)) :
1578       (is_p_hor ? (is_q_hor ? q : p) : r);
1579     CGAL_assertion( (! (num_hor == 1)) || is_site_vertical(prev) );
1580     CGAL_assertion( (! (num_hor == 2)) || is_site_horizontal(prev) );
1581     const Site_2 & next = (is_odd_hor) ?
1582       (is_p_hor ? q : (is_q_hor ? r : p)) :
1583       (is_p_hor ? (is_q_hor ? p : r) : q);
1584     CGAL_assertion( (! (num_hor == 1)) || is_site_vertical(next) );
1585     CGAL_assertion( (! (num_hor == 2)) || is_site_horizontal(next) );
1586 
1587     const RT prevc = hvseg_coord(prev, are_common_hor);
1588     const RT nextc = hvseg_coord(next, are_common_hor);
1589     RT & umid = is_odd_hor ? ux_ : uy_;
1590     RT & udis = is_odd_hor ? uy_ : ux_;
1591     umid = prevc + nextc;
1592     udis = RT(2)*hvseg_coord(odd, is_odd_hor) +
1593            RT(are_common_hor ? +1 : -1) * (prevc - nextc);
1594     uz_ = RT(2);
1595   }
1596 
1597   inline void
compute_sss_bisectors(const Site_2 & p,const Site_2 & q,const Site_2 & r,const Line_2 lines[])1598   compute_sss_bisectors(const Site_2& p, const Site_2& q, const Site_2& r,
1599       const Line_2 lines[])
1600   const
1601   {
1602     CGAL_SDG_DEBUG(std::cout << "debug vring compute_sss_bisectors"
1603         << " p=" << p << " q=" << q  << " r=" << r << std::endl;);
1604     Line_2 bpq = bisector_linf_line(p, q, lines[0], lines[1]);
1605     Line_2 bqr = bisector_linf_line(q, r, lines[1], lines[2]);
1606     compute_intersection_of_lines(bpq, bqr, ux_, uy_, uz_);
1607   }
1608 
1609   inline void
compute_sss_bisectors_old(const Site_2 & p,const Site_2 & q,const Site_2 & r)1610   compute_sss_bisectors_old(
1611       const Site_2& p, const Site_2& q, const Site_2& r)
1612   const
1613   {
1614     CGAL_SDG_DEBUG(std::cout << "debug vring compute_sss_bisectors_old"
1615         << " p=" << p << " q=" << q  << " r=" << r << std::endl;);
1616     Polychainline_2 bpq = bisector_linf(p, q);
1617     Polychainline_2 bqr = bisector_linf(q, r);
1618     Point_2 vv = bpq.first_intersection_point_with(bqr);
1619     ux_ = vv.hx();
1620     uy_ = vv.hy();
1621     uz_ = vv.hw();
1622   }
1623 
1624   //--------------------------------------------------------------------------
1625 
1626   void
analyze_vertex(const Site_2 & s1,const Site_2 & s2,const Site_2 & s3)1627   analyze_vertex(const Site_2& s1, const Site_2& s2, const Site_2& s3)
1628   {
1629     if ( s1.is_point() && s2.is_point() && s3.is_point() ) {
1630       v_type = PPP;
1631     } else if ( s1.is_segment() && s2.is_point() && s3.is_point() ) {
1632       pps_idx = 1;
1633       v_type = PPS;
1634     } else if ( s1.is_point() && s2.is_segment() && s3.is_point() ) {
1635       pps_idx = 2;
1636       v_type = PPS;
1637     } else if ( s1.is_point() && s2.is_point() && s3.is_segment() ) {
1638       pps_idx = 0;
1639       v_type = PPS;
1640     } else if ( s1.is_point() && s2.is_segment() && s3.is_segment() ) {
1641       v_type = PSS;
1642     } else if ( s1.is_segment() && s2.is_point() && s3.is_segment() ) {
1643       v_type = PSS;
1644     } else if ( s1.is_segment() && s2.is_segment() && s3.is_point() ) {
1645       v_type = PSS;
1646     } else {
1647       v_type = SSS;
1648     }
1649   }
1650 
1651   void
compute_vertex(const Site_2 & s1,const Site_2 & s2,const Site_2 & s3)1652   compute_vertex(const Site_2& s1, const Site_2& s2, const Site_2& s3)
1653   const
1654   {
1655     CGAL_assertion( ! is_v_computed );
1656     CGAL_SDG_DEBUG(std::cout << "debug vring compute_vertex "
1657         << s1 << ' ' << s2 << ' ' << s3 << std::endl;);
1658     if ( v_type == PPP ) {
1659       compute_ppp(s1, s2, s3);
1660 
1661     } else if ( s1.is_segment() && s2.is_point() && s3.is_point() ) {
1662       compute_pps(s2, s3, s1);
1663     } else if ( s1.is_point() && s2.is_segment() && s3.is_point() ) {
1664       compute_pps(s3, s1, s2);
1665     } else if ( s1.is_point() && s2.is_point() && s3.is_segment() ) {
1666       compute_pps(s1, s2, s3);
1667 
1668     } else if ( s1.is_point() && s2.is_segment() && s3.is_segment() ) {
1669       compute_pss(s1, s2, s3);
1670     } else if ( s1.is_segment() && s2.is_point() && s3.is_segment() ) {
1671       compute_pss(s2, s3, s1);
1672     } else if ( s1.is_segment() && s2.is_segment() && s3.is_point() ) {
1673       compute_pss(s3, s1, s2);
1674     } else {
1675       compute_sss(s1, s2, s3);
1676     }
1677     is_v_computed = true;
1678   }
1679 
1680 
1681   //--------------------------------------------------------------------------
1682   //--------------------------------------------------------------------------
1683   //                           the orientation test
1684   //--------------------------------------------------------------------------
1685   //--------------------------------------------------------------------------
1686 
1687   Orientation
orientation(const Line_2 & l,PPP_Type)1688   orientation(const Line_2& l, PPP_Type) const
1689   {
1690     Sign s_uz = CGAL::sign(uz_);
1691     Sign s_l =
1692       CGAL::sign(l.a() * ux_ + l.b() * uy_ + l.c() * uz_);
1693 
1694     return s_uz * s_l;
1695   }
1696 
1697 
1698   //--------------------------------------------------------------------------
1699 
1700   Orientation
orientation(const Line_2 & l,PPS_Type)1701   orientation(const Line_2& l, PPS_Type) const
1702   {
1703     Sign s_uz = CGAL::sign(uz_);
1704     Sign s_l = CGAL::sign(l.a() * ux_ + l.b() * uy_ + l.c() * uz_);
1705 
1706     return s_uz * s_l;
1707   }
1708 
1709 
1710   //--------------------------------------------------------------------------
1711 
1712   // the cases PSS and SSS are identical
1713   template<class Type>
1714   Orientation
orientation(const Line_2 & l,Type)1715   orientation(const Line_2& l, Type) const
1716   {
1717     Sign s_uz = CGAL::sign(uz_);
1718     Sign s_l = CGAL::sign(l.a() * ux_ + l.b() * uy_ + l.c() * uz_);
1719 
1720     return s_uz * s_l;
1721   }
1722 
1723 
1724   //--------------------------------------------------------------------------
1725   //--------------------------------------------------------------------------
1726   //                              the incircle test
1727   //--------------------------------------------------------------------------
1728   //--------------------------------------------------------------------------
1729 
1730   //--------------------------------------------------------------------------
1731   //  the incircle test when the fourth site is a point
1732   //--------------------------------------------------------------------------
1733 
1734   //--------------------------------------------------------------------------
1735 
check_easy_degeneracies(const Site_2 & t,PPS_Type,bool & use_result)1736   Sign check_easy_degeneracies(const Site_2& t, PPS_Type,
1737                                bool& use_result) const
1738   {
1739     CGAL_precondition( t.is_point() );
1740 
1741     use_result = false;
1742     if (  ( p_.is_point() && same_points(p_, t) ) ||
1743           ( q_.is_point() && same_points(q_, t) ) ||
1744           ( r_.is_point() && same_points(r_, t) )  ) {
1745       use_result = true;
1746       return ZERO;
1747     }
1748 
1749     if (
1750       ( p_.is_segment() && is_endpoint_of(t, p_) ) ||
1751       ( q_.is_segment() && is_endpoint_of(t, q_) ) ||
1752       ( r_.is_segment() && is_endpoint_of(t, r_) )  )
1753     {
1754       use_result = true;
1755       return POSITIVE;
1756     }
1757 
1758     if (
1759       ( p_.is_segment() && is_on_hv_seg_line(t, p_) ) ||
1760       ( q_.is_segment() && is_on_hv_seg_line(t, q_) ) ||
1761       ( r_.is_segment() && is_on_hv_seg_line(t, r_) )  )
1762     {
1763       use_result = true;
1764       return POSITIVE;
1765     }
1766 
1767     return ZERO;
1768   }
1769 
1770   inline
check_easy_degeneracies(const Site_2 & t,PSS_Type,bool & use_result)1771   Sign check_easy_degeneracies(const Site_2& t, PSS_Type,
1772                                bool& use_result) const
1773   {
1774     CGAL_precondition( t.is_point() );
1775 
1776     return check_easy_degeneracies(t, PPS_Type(), use_result);
1777   }
1778 
1779   inline
check_easy_degeneracies(const Site_2 & t,SSS_Type,bool & use_result)1780   Sign check_easy_degeneracies(const Site_2& t, SSS_Type,
1781                                bool& use_result) const
1782   {
1783     CGAL_precondition( t.is_point() );
1784     use_result = false;
1785     if (is_endpoint_of(t, p_) || is_endpoint_of(t, q_) ||
1786         is_endpoint_of(t, r_) ) {
1787       use_result = true;
1788       return POSITIVE;
1789     }
1790     return ZERO;
1791   }
1792 
1793   //--------------------------------------------------------------------------
1794 
1795   template<class Type>
1796   inline
incircle_p(const Site_2 & st,Type type)1797   Sign incircle_p(const Site_2& st, Type type) const
1798   {
1799     CGAL_precondition( st.is_point() );
1800 
1801     bool use_result(false);
1802     Sign s = check_easy_degeneracies(st, type, use_result);
1803     if ( use_result ) { return s; }
1804 
1805     return incircle_p_no_easy(st, type);
1806   }
1807 
1808   //--------------------------------------------------------------------------
1809 
incircle_p(const Site_2 & st,PPP_Type)1810   Sign incircle_p(const Site_2& st, PPP_Type) const
1811   {
1812     CGAL_precondition( st.is_point() );
1813 
1814     Point_2 t = st.point();
1815 
1816     Bounded_side bs =
1817       side_of_bounded_square(p_.point(), q_.point(), r_.point(), t);
1818 
1819     switch(bs) {
1820       case ON_UNBOUNDED_SIDE:
1821         CGAL_SDG_DEBUG(std::cout
1822             << "debug incircle_p returns POSITIVE" << std::endl;);
1823         return POSITIVE;
1824       case ON_BOUNDED_SIDE:
1825         CGAL_SDG_DEBUG(std::cout
1826             << "debug incircle_p returns NEGATIVE" << std::endl;);
1827         return NEGATIVE;
1828       default:
1829         CGAL_SDG_DEBUG(std::cout
1830             << "debug incircle_p returns ZERO" << std::endl;);
1831         return ZERO;
1832     }
1833   }
1834 
1835   //--------------------------------------------------------------------------
1836 
incircle_p_no_easy(const Site_2 & st,PPS_Type)1837   Sign incircle_p_no_easy(const Site_2& st, PPS_Type ) const
1838   {
1839     CGAL_precondition( st.is_point() );
1840     compute_v_if_not_computed();
1841 
1842     CGAL_SDG_DEBUG(std::cout << "debug vring incircle_p_no_easy PPS p="
1843       << p_ << " q=" << q_  << " r=" << r_ << " t=" << st
1844       << std::endl;);
1845 
1846     Point_2 t = st.point();
1847 
1848     Point_2 pointref = p_ref().point();
1849 
1850     CGAL_SDG_DEBUG(std::cout
1851         << "debug vring incircle_p_no_easy PPS pointref="
1852         << pointref << std::endl;);
1853 
1854     RT vx = ux_ - pointref.x() * uz_;
1855     RT vy = uy_ - pointref.y() * uz_;
1856 
1857     RT Rs =
1858       (CGAL::max) ( CGAL::abs(vx), CGAL::abs(vy) );
1859 
1860     RT scalediffdvtx = ux_ - t.x() * uz_;
1861     RT scalediffdvty = uy_ - t.y() * uz_;
1862 
1863     RT Rs1 =
1864       (CGAL::max)(
1865         CGAL::abs(scalediffdvtx),
1866         CGAL::abs(scalediffdvty) );
1867 
1868 
1869     Sign crude = CGAL::sign(Rs1 - Rs);
1870 
1871     if (crude != ZERO) {
1872       return crude;
1873     } else {
1874       CGAL_SDG_DEBUG(std::cout
1875           << "debug vring refining in incircle_p_no_easy PPS pqr=("
1876           << p_ << ", " << q_ << ", " << r_ << "), "
1877           << "t=" << t
1878           << std::endl;);
1879       // here crude == ZERO, so
1880       // we might have to refine
1881 
1882       // tocheck
1883 
1884       const std::tuple<
1885         const Site_2 &, const Site_2 &, const Site_2 &> sites =
1886          r_.is_segment() ? sdg_tuple_maker(p_, q_, r_) :
1887         (p_.is_segment() ? sdg_tuple_maker(q_, r_, p_) :
1888                            sdg_tuple_maker(r_, p_, q_) );
1889 
1890       const Site_2 & p1 = std::get<0>(sites);
1891       const Site_2 & p2 = std::get<1>(sites);
1892       const Site_2 & s  = std::get<2>(sites);
1893 
1894       const RT d_fine = (CGAL::min)(CGAL::abs(scalediffdvtx),
1895                                     CGAL::abs(scalediffdvty));
1896       for (size_t i = 0; i < 2; ++i) {
1897         const Site_2 & cur = (i == 0) ? p1 : p2;
1898         const Point_2 pref = cur.point();
1899         const RT scalediffdvpx = ux_ - pref.x() * uz_;
1900         const RT scalediffdvpy = uy_ - pref.y() * uz_;
1901 
1902         Comparison_result sidecmp = EQUAL;
1903         const bool p_t_samex =
1904           CGAL::compare(scalediffdvpx, scalediffdvtx) == EQUAL;
1905         const bool p_t_on_same_ver_side =
1906           (p_t_samex) &&
1907           (CGAL::compare(CGAL::abs(scalediffdvpx), Rs1) == EQUAL) ;
1908         if (p_t_on_same_ver_side) {
1909           sidecmp = CGAL::compare(d_fine, CGAL::abs(scalediffdvpy));
1910         }
1911         const bool p_t_samey =
1912           CGAL::compare(scalediffdvpy, scalediffdvty) == EQUAL;
1913         const bool p_t_on_same_hor_side =
1914           (p_t_samey) &&
1915           (CGAL::compare(CGAL::abs(scalediffdvpy), Rs1) == EQUAL) ;
1916         if (p_t_on_same_hor_side) {
1917           sidecmp = CGAL::compare(d_fine, CGAL::abs(scalediffdvpx));
1918         }
1919         CGAL_SDG_DEBUG(std::cout << "vring test with p=" << cur
1920             << ", sidecmp=" << sidecmp << std::endl; );
1921         if (sidecmp == SMALLER) {
1922           return NEGATIVE;
1923         } else if (sidecmp == LARGER) {
1924           return POSITIVE;
1925         }
1926       }
1927 
1928       const bool is_s_hor = is_site_horizontal(s);
1929       const bool is_s_ver = is_site_vertical(s);
1930 
1931       const bool is_p1_endp_of_s = is_endpoint_of(p1, s);
1932       const bool is_p2_endp_of_s = is_endpoint_of(p2, s);
1933 
1934       // check for p, q with same coordinate and r non-hv segment
1935       if ((! (is_s_hor || is_s_ver)) &&
1936           (! (is_p1_endp_of_s || is_p2_endp_of_s))
1937          ) {
1938         CGAL_SDG_DEBUG(std::cout << "debug vring seg=" << s
1939             << " is non-axis parallel"
1940             << " and no points are its endpoints" << std::endl;);
1941         const bool pqsamex = scmpx(p1, p2) == EQUAL;
1942         bool pqsamey (false);
1943         if (pqsamex) {
1944           CGAL_SDG_DEBUG(std::cout << "debug vring points have same x, "
1945               << " might be on same Linf vertical side"
1946               << std::endl;);
1947         } else {
1948           pqsamey = scmpy(p1, p2) == EQUAL;
1949           if (pqsamey) {
1950             CGAL_SDG_DEBUG(std::cout << "debug vring points have same y, "
1951                 << " might be on same Linf horizontal side" << std::endl;);
1952           }
1953         }
1954         if (pqsamex || pqsamey) {
1955           Line_2 lr = compute_supporting_line(s.supporting_site());
1956           Homogeneous_point_2 rref = compute_linf_projection_hom(lr, point());
1957           if ( pqsamex && (CGAL::compare(rref.x(), pointref.x()) == EQUAL) ) {
1958             RT scalediffdvry = uy_ - rref.y() * uz_;
1959             if (CGAL::sign(scalediffdvry) == CGAL::sign(scalediffdvty)) {
1960               if (CGAL::compare(CGAL::abs(scalediffdvtx),
1961                                 CGAL::abs(scalediffdvty)) == SMALLER) {
1962                 return NEGATIVE;
1963               }
1964             }
1965           } // end of pqsamex and rref same x case
1966           if ( pqsamey && (CGAL::compare(rref.y(), pointref.y()) == EQUAL) ) {
1967             RT scalediffdvrx = ux_ - rref.x() * uz_;
1968             if (CGAL::sign(scalediffdvrx) == CGAL::sign(scalediffdvtx)) {
1969               if (CGAL::compare(CGAL::abs(scalediffdvty),
1970                                 CGAL::abs(scalediffdvtx)) == SMALLER) {
1971                 return NEGATIVE;
1972               }
1973             }
1974           } // end of pqsamey and rref same y case
1975         } // end of case: pqsamex || pqsamey
1976       } // end of non-hv segment r case with p, q non-endpoints of r
1977 
1978       return ZERO;
1979 
1980       //FT radius_fine = linf_fine_radius(vv, p, q, r, type);
1981 
1982       //return CGAL::compare(d_fine, radius_fine);
1983 
1984     }
1985 
1986   }
1987 
1988   //--------------------------------------------------------------------------
1989 
incircle_p_no_easy(const Site_2 & st,PSS_Type)1990   Sign incircle_p_no_easy(const Site_2& st, PSS_Type ) const
1991   {
1992     CGAL_precondition( st.is_point() );
1993     compute_v_if_not_computed();
1994     const Point_2 t = st.point();
1995 
1996     const Point_2 pref = p_ref().point();
1997 
1998     const RT xref = pref.x();
1999     const RT yref = pref.y();
2000 
2001     const RT vx = ux_ - xref * uz_;
2002     const RT vy = uy_ - yref * uz_;
2003 
2004     const RT Rs =
2005       (CGAL::max) ( CGAL::abs(vx), CGAL::abs(vy) );
2006 
2007     const RT tx = t.x() ;
2008     const RT ty = t.y() ;
2009 
2010     const RT scalediffdvtx = ux_ - tx * uz_;
2011     const RT scalediffdvty = uy_ - ty * uz_;
2012 
2013     const RT Rs1 =
2014       (CGAL::max)(
2015         CGAL::abs(scalediffdvtx),
2016         CGAL::abs(scalediffdvty) );
2017 
2018     const Sign s_Q = CGAL::sign(Rs1 - Rs);
2019 
2020     if (s_Q != ZERO) {
2021       return s_Q;
2022     } else {
2023       CGAL_SDG_DEBUG(std::cout
2024           << "debug vring refining in incircle_p_no_easy PSS pqr=("
2025           << p_ << ", " << q_ << ", " << r_ << "), "
2026           << "t=" << t
2027           << std::endl;);
2028       // here crude == ZERO, so
2029       // we might have to refine
2030 
2031       // tocheck
2032 
2033       const std::tuple<
2034         const Site_2 &, const Site_2 &, const Site_2 &> sites =
2035          p_.is_point() ? sdg_tuple_maker(p_, q_, r_) :
2036         (q_.is_point() ? sdg_tuple_maker(q_, r_, p_) :
2037                          sdg_tuple_maker(r_, p_, q_) );
2038 
2039       const Site_2 & pt_site = std::get<0>(sites);
2040       const Site_2 & s1 = std::get<1>(sites);
2041       const Site_2 & s2 = std::get<2>(sites);
2042 
2043       const bool is_s1src_s2 = is_endpoint_of(s1.source_site(), s2);
2044       const bool is_s1trg_s2 = is_endpoint_of(s1.target_site(), s2);
2045 
2046       if (is_s1src_s2 || is_s1trg_s2) {
2047         if ((is_site_h_or_v(s1) && (! is_site_h_or_v(s2))) ||
2048             (is_site_h_or_v(s2) && (! is_site_h_or_v(s1)))   )
2049         {
2050           CGAL_SDG_DEBUG(std::cout << "debug vring "
2051               << "s1, s2 candidates" << std::endl; );
2052           if (is_site_horizontal(s1) || is_site_horizontal(s2)) {
2053             Site_2 s1test = is_s1src_s2?
2054                        (s1.source_site()):
2055                        (s1.target_site());
2056             if (scmpx(s1test, st)
2057                 == EQUAL)
2058             {
2059               // return NEGATIVE or ZERO
2060               Point_2 s1ref =
2061                       (is_s1src_s2?
2062                        s1.source_site(): s1.target_site())
2063                       .point();
2064               RT scalediffdvs1y = uy_ - s1ref.y() * uz_;
2065               Comparison_result test =
2066                 CGAL::compare(
2067                     CGAL::abs(scalediffdvty),
2068                     CGAL::abs(scalediffdvs1y));
2069               return (test == SMALLER) ? NEGATIVE : ZERO;
2070             }
2071           } else { // one of q, r is vertical
2072             if (scmpy(is_s1src_s2?
2073                       s1.source_site(): s1.target_site(), st)
2074                 == EQUAL)
2075             {
2076               // return NEGATIVE or ZERO
2077               CGAL_SDG_DEBUG(std::cout << "debug vring "
2078                   << "vertical case" << std::endl; );
2079               Point_2 s1ref =
2080                       (is_s1src_s2?
2081                        s1.source_site(): s1.target_site())
2082                       .point();
2083               RT scalediffdvs1x = ux_ - s1ref.x() * uz_;
2084               CGAL_SDG_DEBUG(std::cout << "debug vring "
2085                   << "scalediffdvs1x=" << scalediffdvs1x
2086                   << " scalediffdvtx=" << scalediffdvtx << std::endl; );
2087               Comparison_result test =
2088                 CGAL::compare(
2089                     CGAL::abs(scalediffdvtx),
2090                     CGAL::abs(scalediffdvs1x));
2091               return (test == SMALLER) ? NEGATIVE : ZERO;
2092             }
2093           }
2094         }
2095       }
2096 
2097       const RT d_fine = (CGAL::min)(CGAL::abs(scalediffdvtx),
2098                                     CGAL::abs(scalediffdvty));
2099       const Point_2 pref = pt_site.point();
2100       const RT scalediffdvpx = ux_ - pref.x() * uz_;
2101       const RT scalediffdvpy = uy_ - pref.y() * uz_;
2102       Comparison_result sidecmp = EQUAL;
2103       const bool p_t_samex =
2104         CGAL::compare(scalediffdvpx, scalediffdvtx) == EQUAL;
2105       const bool p_t_on_same_ver_side =
2106         (p_t_samex) &&
2107         (CGAL::compare(CGAL::abs(scalediffdvpx), Rs1) == EQUAL) ;
2108       if (p_t_on_same_ver_side) {
2109         sidecmp = CGAL::compare(d_fine, CGAL::abs(scalediffdvpy));
2110       }
2111       const bool p_t_samey =
2112         CGAL::compare(scalediffdvpy, scalediffdvty) == EQUAL;
2113       const bool p_t_on_same_hor_side =
2114         (p_t_samey) &&
2115         (CGAL::compare(CGAL::abs(scalediffdvpy), Rs1) == EQUAL) ;
2116       if (p_t_on_same_hor_side) {
2117         sidecmp = CGAL::compare(d_fine, CGAL::abs(scalediffdvpx));
2118       }
2119       CGAL_SDG_DEBUG(std::cout << "debug: PSS temporary sidecmp = "
2120           << sidecmp << std::endl;);
2121       if (sidecmp == SMALLER) {
2122         return NEGATIVE;
2123       } else if (sidecmp == LARGER) {
2124         return POSITIVE;
2125       }
2126 
2127       if (p_t_samex || p_t_samey) {
2128         return ZERO;
2129       }
2130 
2131       if (! is_site_h_or_v(s1)) {
2132         // here s1 is non-axis parallel
2133         // therefore, it touches the square at a corner
2134         if (points_inside_touching_sides_v(
2135               s1, pt_site, s2, st, this->point())) {
2136           return NEGATIVE;
2137         }
2138       }
2139 
2140       if (! is_site_h_or_v(s2)) {
2141         // here s2 is non-axis parallel
2142         // therefore, it touches the square at a corner
2143         if (points_inside_touching_sides_v(
2144               s2, pt_site, s1, st, this->point())) {
2145           return NEGATIVE;
2146         }
2147       }
2148 
2149       CGAL_SDG_DEBUG(std::cout
2150           << "debug vring PSS P return final ZERO"
2151           << std::endl;);
2152       return ZERO;
2153     }
2154 
2155   }
2156 
2157   //--------------------------------------------------------------------------
2158 
incircle_p_no_easy(const Site_2 & st,SSS_Type)2159   Sign incircle_p_no_easy(const Site_2& st, SSS_Type ) const
2160   {
2161     CGAL_precondition( st.is_point() );
2162     compute_v_if_not_computed();
2163 
2164     Point_2 t = st.point();
2165 
2166     Line_2 l = compute_supporting_line(p_.supporting_site());
2167     Homogeneous_point_2 hp = compute_linf_projection_hom(l, point());
2168 
2169     RT dup =
2170       (CGAL::max)(CGAL::abs(ux_ - hp.x() * uz_),
2171                   CGAL::abs(uy_ - hp.y() * uz_));
2172 
2173     RT scalediffdvtx = ux_ - t.x() * uz_;
2174     RT scalediffdvty = uy_ - t.y() * uz_;
2175 
2176     RT dut =
2177       (CGAL::max)(
2178         CGAL::abs(scalediffdvtx),
2179         CGAL::abs(scalediffdvty) );
2180 
2181     Sign crude_sign = CGAL::sign(dut - dup);
2182     if (crude_sign != ZERO) {
2183       return crude_sign;
2184     } else {
2185       CGAL_SDG_DEBUG(std::cout
2186           << "debug vring refining in incircle_p_no_easy SSS pqr=("
2187           << p_ << ", " << q_ << ", " << r_ << "), "
2188           << "t=" << t
2189           << std::endl;);
2190 
2191       const Site_2 * s1_ptr;
2192       const Site_2 * s2_ptr;
2193 
2194       const Site_2 * s1ptr_arr[] = {&p_, &q_, &r_};
2195       const Site_2 * s2ptr_arr[] = {&q_, &r_, &p_};
2196 
2197       for (int i = 0; i < 3; i++) {
2198         s1_ptr = s1ptr_arr[i];
2199         s2_ptr = s2ptr_arr[i];
2200 
2201         CGAL_SDG_DEBUG(std::cout << "debug vring check for candidates"
2202             << "(s1, s2) = " << *s1_ptr << ", " << *s2_ptr << std::endl; );
2203 
2204         bool is_s1src_s2 = is_endpoint_of((*s1_ptr).source_site(), *s2_ptr);
2205         bool is_s1trg_s2 = is_endpoint_of((*s1_ptr).target_site(), *s2_ptr);
2206 
2207         if (is_s1src_s2 || is_s1trg_s2) {
2208           if ((is_site_h_or_v(*s1_ptr) && (! is_site_h_or_v(*s2_ptr))) ||
2209               (is_site_h_or_v(*s2_ptr) && (! is_site_h_or_v(*s1_ptr)))   )
2210           {
2211             CGAL_SDG_DEBUG(std::cout << "debug vring "
2212                 << "s1, s2 candidates" << std::endl; );
2213             if (is_site_horizontal(*s1_ptr) || is_site_horizontal(*s2_ptr)) {
2214               Site_2 s1test = is_s1src_s2?
2215                 ((*s1_ptr).source_site()):
2216                 ((*s1_ptr).target_site());
2217               if (scmpx(s1test, st)
2218                   == EQUAL)
2219               {
2220                 // return NEGATIVE or ZERO
2221                 Point_2 s1ref =
2222                   (is_s1src_s2?
2223                    (*s1_ptr).source_site(): (*s1_ptr).target_site())
2224                   .point();
2225                 RT scalediffdvs1y = uy_ - s1ref.y() * uz_;
2226                 Comparison_result test =
2227                   CGAL::compare(
2228                       CGAL::abs(scalediffdvty),
2229                       CGAL::abs(scalediffdvs1y));
2230                 return (test == SMALLER) ? NEGATIVE : ZERO;
2231               }
2232             } else { // one of q, r is vertical
2233               if (scmpy(is_s1src_s2?
2234                     (*s1_ptr).source_site(): (*s1_ptr).target_site(), st)
2235                   == EQUAL)
2236               {
2237                 // return NEGATIVE or ZERO
2238                 CGAL_SDG_DEBUG(std::cout << "debug vring "
2239                     << "vertical case" << std::endl; );
2240                 Point_2 s1ref =
2241                   (is_s1src_s2?
2242                    (*s1_ptr).source_site(): (*s1_ptr).target_site())
2243                   .point();
2244                 RT scalediffdvs1x = ux_ - s1ref.x() * uz_;
2245                 CGAL_SDG_DEBUG(std::cout << "debug vring "
2246                     << "scalediffdvs1x=" << scalediffdvs1x
2247                     << " scalediffdvtx=" << scalediffdvtx << std::endl; );
2248                 Comparison_result test =
2249                   CGAL::compare(
2250                       CGAL::abs(scalediffdvtx),
2251                       CGAL::abs(scalediffdvs1x));
2252                 return (test == SMALLER) ? NEGATIVE : ZERO;
2253               }
2254             }
2255           }
2256         }
2257       } // end for
2258 
2259 
2260 
2261 
2262 
2263       return ZERO;
2264     }
2265   }
2266 
2267 
2268 
2269   //--------------------------------------------------------------------------
2270   //--------------------------------------------------------------------------
2271 
incircle_p(const Site_2 & t)2272   Sign incircle_p(const Site_2& t) const
2273   {
2274     CGAL_SDG_DEBUG(std::cout << "debug: entering vring incircle_p with "
2275       << "v_type=" << v_type << " p="
2276       << p_ << " q=" << q_ << " r=" << r_ << " t=" << t
2277       << std::endl;);
2278 
2279     if ( is_degenerate_Voronoi_circle() ) {
2280       return POSITIVE;
2281     }
2282 
2283     Sign s(ZERO);
2284     switch ( v_type ) {
2285     case PPP:
2286       s = incircle_p(t, PPP_Type());
2287       break;
2288     case PPS:
2289       s = incircle_p(t, PPS_Type());
2290       break;
2291     case PSS:
2292       s = incircle_p(t, PSS_Type());
2293       break;
2294     case SSS:
2295       s = incircle_p(t, SSS_Type());
2296       break;
2297     }
2298 
2299     return s;
2300   }
2301 
incircle_p_no_easy(const Site_2 & t)2302   Sign incircle_p_no_easy(const Site_2& t) const
2303   {
2304     Sign s(ZERO);
2305     switch ( v_type ) {
2306     case PPP:
2307       s = incircle_p(t, PPP_Type());
2308       break;
2309     case PPS:
2310       s = incircle_p_no_easy(t, PPS_Type());
2311       break;
2312     case PSS:
2313       s = incircle_p_no_easy(t, PSS_Type());
2314       break;
2315     case SSS:
2316       s = incircle_p_no_easy(t, SSS_Type());
2317       break;
2318     }
2319 
2320     return s;
2321   }
2322 
2323   //--------------------------------------------------------------------------
2324 
2325   //--------------------------------------------------------------------------
2326   //  the incircle test when the fourth site is a segment
2327   //--------------------------------------------------------------------------
2328 
2329   //--------------------------------------------------------------------------
2330 
2331   Oriented_side
oriented_side_l2(const Line_2 & l,const Point_2 & p,PPP_Type)2332   oriented_side_l2(const Line_2& l, const Point_2& p, PPP_Type) const
2333   {
2334     Sign s_uz = CGAL::sign(uz_);
2335 
2336     RT px = uz_ * p.x() - ux_;
2337     RT py = uz_ * p.y() - uy_;
2338 
2339     Sign s1 = CGAL::sign(l.b() * px - l.a() * py);
2340 
2341     return s_uz * s1;
2342   }
2343 
2344   Oriented_side
oriented_side_l2(const Line_2 & l,const Point_2 & p,PPS_Type)2345   oriented_side_l2(const Line_2& l, const Point_2& p, PPS_Type) const
2346   {
2347     RT dx = ux_ - uz_ * p.x();
2348     RT dy = uy_ - uz_ * p.y();
2349 
2350     return CGAL::sign(uz_) * CGAL::sign(dy * l.a() - dx * l.b());
2351   }
2352 
2353   // the cases PSS and SSS are identical
2354   template<class Type>
2355   Oriented_side
oriented_side_l2(const Line_2 & l,const Point_2 & p,Type)2356   oriented_side_l2(const Line_2& l, const Point_2& p, Type) const
2357   {
2358     RT px = p.x();
2359     RT py = p.y();
2360 
2361     RT dx = ux_ - px * uz_;
2362     RT dy = uy_ - py * uz_;
2363 
2364     RT a = l.a();
2365     RT b = l.b();
2366 
2367     return CGAL::sign(uz_) * CGAL::sign(a * dy - b * dx);
2368   }
2369 
2370 
2371   // philaris:
2372   template<class Type>
2373   Oriented_side
oriented_side_linf(const Line_2 & l,const Point_2 & p,Type)2374   oriented_side_linf(const Line_2& l, const Point_2& p, Type) const
2375   {
2376     CGAL_SDG_DEBUG(std::cout << "debug oriented_side_linf " << std::endl;);
2377 
2378     Point_2 vv (ux_, uy_, uz_);
2379 
2380     Line_2 l1 = compute_linf_perpendicular(l, vv);
2381 
2382     return oriented_side_of_line(l1, p);
2383   }
2384 
2385 
2386   //--------------------------------------------------------------------------
2387   //--------------------------------------------------------------------------
2388 
2389 
incircle(const Line_2 & l,PPP_Type)2390   Sign incircle(const Line_2& l, PPP_Type) const
2391   {
2392 
2393     Point_2 pref = p_ref().point();
2394     Homogeneous_point_2 hp = compute_linf_projection_hom(l, point());
2395 
2396     CGAL_SDG_DEBUG(std::cout << "debug incircle l PPP: pref="
2397       << pref << std::endl;);
2398 
2399     RT dul = (CGAL::max)(
2400         CGAL::abs(ux_ - hp.x() * uz_),
2401         CGAL::abs(uy_ - hp.y() * uz_));
2402 
2403     RT dupref = (CGAL::max)(
2404         CGAL::abs(ux_ - pref.x() * uz_),
2405         CGAL::abs(uy_ - pref.y() * uz_));
2406 
2407     Comparison_result cr = CGAL::compare(dul, dupref);
2408 
2409     if ( cr == LARGER ) { return POSITIVE; }
2410     if ( cr == SMALLER ) { return NEGATIVE; }
2411 
2412     // here cr == EQUAL == ZERO, so
2413     // we might have to refine
2414     CGAL_SDG_DEBUG(std::cout
2415       << "debug vring refining in incircle l PPP pqr=("
2416       << p_ << ", " << q_ << ", " << r_ << "), "
2417       << "hp(x,y)=" << hp.x() << ' ' << hp.y()
2418       << ", l: " << l.a() << ' ' << l.b() << ' ' <<  l.c()
2419       << ", u(x,y,z)= " << ux_ << ' ' << uy_ << ' ' << uz_
2420       << std::endl;);
2421 
2422     Comparison_result other = linf_refine(l, hp);
2423 
2424     if (cr != other) {
2425       CGAL_SDG_DEBUG(std::cout
2426           << "incircle l PPP instead of 0 returning " << other
2427           << std::endl;);
2428     }
2429 
2430     return other;
2431 
2432   }
2433 
2434   //--------------------------------------------------------------------------
2435 
incircle(const Line_2 & l,PPS_Type)2436   Sign incircle(const Line_2& l, PPS_Type) const
2437   {
2438     Point_2 pref = p_ref().point();
2439 
2440     CGAL_SDG_DEBUG(std::cout << "debug incircle l PPS: pref="
2441       << pref << std::endl;);
2442 
2443     RT vx = ux_ - pref.x() * uz_;
2444     RT vy = uy_ - pref.y() * uz_;
2445 
2446     RT dupref = (CGAL::max)(CGAL::abs(vx), CGAL::abs(vy));
2447 
2448     Homogeneous_point_2 hp = compute_linf_projection_hom(l, point());
2449 
2450     RT dul = (CGAL::max)(
2451         CGAL::abs(ux_ - hp.x() * uz_),
2452         CGAL::abs(uy_ - hp.y() * uz_));
2453 
2454     Sign cr = CGAL::sign(dul - dupref);
2455 
2456     if (cr != ZERO) {
2457       return cr;
2458     }
2459 
2460     // here cr == EQUAL == ZERO, so
2461     // we might have to refine
2462     CGAL_SDG_DEBUG(std::cout
2463       << "debug vring refining in incircle l PPS pqr=("
2464       << p_ << ", " << q_ << ", " << r_ << "), "
2465       << "hp(x,y)=" << hp.x() << ' ' << hp.y()
2466       << ", l: " << l.a() << ' ' << l.b() << ' ' <<  l.c()
2467       << ", u(x,y,z)= " << ux_ << ' ' << uy_ << ' ' << uz_
2468       << std::endl;);
2469 
2470     Comparison_result other = linf_refine(l, hp);
2471 
2472     if (cr != other) {
2473       CGAL_SDG_DEBUG(std::cout
2474           << "incircle l PPS instead of 0 returning " << other
2475           << std::endl;);
2476     }
2477 
2478     return other;
2479   }
2480 
2481 
2482   //--------------------------------------------------------------------------
2483 
incircle(const Line_2 & l,PSS_Type)2484   Sign incircle(const Line_2& l, PSS_Type) const
2485   {
2486     Point_2 pref = p_ref().point();
2487 
2488     CGAL_SDG_DEBUG(std::cout << "debug incircle l PSS: pref="
2489       << pref << std::endl;);
2490 
2491     RT vx = ux_ - (pref.x() ) * uz_;
2492     RT vy = uy_ - (pref.y() ) * uz_;
2493 
2494     RT dupref = (CGAL::max)(CGAL::abs(vx), CGAL::abs(vy));
2495 
2496     Homogeneous_point_2 lhp = compute_linf_projection_hom(l, point());
2497 
2498     RT dul = (CGAL::max)(
2499         CGAL::abs(ux_ - lhp.x() * uz_),
2500         CGAL::abs(uy_ - lhp.y() * uz_));
2501 
2502     Sign cr = CGAL::sign(dul - dupref);
2503 
2504     if (cr != ZERO) {
2505       return cr;
2506     }
2507 
2508     // here cr == EQUAL == ZERO, so
2509     // we might have to refine
2510     CGAL_SDG_DEBUG(std::cout
2511       << "debug vring refining in incircle l PSS pqr=("
2512       << p_ << ", " << q_ << ", " << r_ << "), "
2513       << "hp(x,y)=" << lhp.x() << ' ' << lhp.y()
2514       << ", l: " << l.a() << ' ' << l.b() << ' ' <<  l.c()
2515       << ", u(x,y,z)= " << ux_ << ' ' << uy_ << ' ' << uz_
2516       << std::endl;);
2517 
2518     Comparison_result other = linf_refine(l, lhp);
2519 
2520     if (cr != other) {
2521       CGAL_SDG_DEBUG(std::cout
2522           << "incircle l PSS instead of 0 returning " << other
2523           << std::endl;);
2524     }
2525 
2526     return other;
2527 
2528   }
2529 
2530 
2531   //--------------------------------------------------------------------------
2532 
incircle(const Line_2 & l,SSS_Type)2533   Sign incircle(const Line_2& l, SSS_Type) const
2534   {
2535     Line_2 lref = compute_supporting_line(p_.supporting_site());
2536     Homogeneous_point_2 lrefhp =
2537       compute_linf_projection_hom(lref, point());
2538     RT dulref = (CGAL::max)(
2539         CGAL::abs(ux_ - lrefhp.x() * uz_),
2540         CGAL::abs(uy_ - lrefhp.y() * uz_));
2541 
2542     Homogeneous_point_2 lhp = compute_linf_projection_hom(l, point());
2543     RT dul = (CGAL::max)(
2544         CGAL::abs(ux_ - lhp.x() * uz_),
2545         CGAL::abs(uy_ - lhp.y() * uz_));
2546 
2547     Sign cr = CGAL::sign(dul - dulref);
2548 
2549     if (cr != ZERO) {
2550       return cr;
2551     }
2552 
2553     // here cr == EQUAL == ZERO, so
2554     // we might have to refine
2555     CGAL_SDG_DEBUG(std::cout
2556       << "debug vring refining in incircle l PSS pqr=("
2557       << p_ << ", " << q_ << ", " << r_ << "), "
2558       << "lhp(x,y)=" << lhp.x() << ' ' << lhp.y()
2559       << ", l: " << l.a() << ' ' << l.b() << ' ' <<  l.c()
2560       << ", u(x,y,z)= " << ux_ << ' ' << uy_ << ' ' << uz_
2561       << std::endl;);
2562 
2563     Comparison_result other = linf_refine(l, lhp);
2564 
2565     if (cr != other) {
2566       CGAL_SDG_DEBUG(std::cout
2567           << "incircle l SSS instead of 0 returning " << other
2568           << std::endl;);
2569     }
2570 
2571     return other;
2572   }
2573 
2574   //--------------------------------------------------------------------------
2575   //--------------------------------------------------------------------------
2576 
2577   template<class Type>
incircle_s(const Site_2 & t,Type type)2578   Sign incircle_s(const Site_2& t, Type type) const
2579   {
2580     CGAL_precondition( t.is_segment() );
2581 
2582     if ( v_type == PPP || v_type == PPS ) {
2583       if (  p_.is_point() && q_.is_point() &&
2584             is_endpoint_of(p_, t) && is_endpoint_of(q_, t)  ) {
2585         return NEGATIVE;
2586       }
2587 
2588       if (  p_.is_point() && r_.is_point() &&
2589             is_endpoint_of(p_, t) && is_endpoint_of(r_, t)  ){
2590         return NEGATIVE;
2591       }
2592 
2593       if (  q_.is_point() && r_.is_point() &&
2594             is_endpoint_of(q_, t) && is_endpoint_of(r_, t)  ){
2595         return NEGATIVE;
2596       }
2597     }
2598 
2599     if ( v_type == PSS ) {
2600       if ( p_.is_segment() &&
2601            same_segments(p_.supporting_site(),
2602                          t.supporting_site()) ) {
2603         return POSITIVE;
2604       }
2605       if ( q_.is_segment() &&
2606            same_segments(q_.supporting_site(),
2607                          t.supporting_site()) ) {
2608         return POSITIVE;
2609       }
2610       if ( r_.is_segment() &&
2611            same_segments(r_.supporting_site(),
2612                          t.supporting_site()) ) {
2613         return POSITIVE;
2614       }
2615     }
2616 
2617     Sign retval = incircle_s_no_easy(t, type);
2618 
2619     CGAL_SDG_DEBUG(std::cout
2620         << "debug vring incircle_s: about to return retval of"
2621         << " incircle_s_no_easy = " << retval << std::endl;);
2622 
2623     return retval;
2624   }
2625 
2626   template<class Type>
incircle_s_no_easy(const Site_2 & t,Type type)2627   Sign incircle_s_no_easy(const Site_2& t, Type type) const
2628   {
2629     CGAL_SDG_DEBUG(std::cout << "debug vring fn incircle_s_no_easy pqrt= ("
2630         << p_ << ") (" << q_ << ") (" << r_ << ") (" << t << ")"
2631         << std::endl;);
2632 
2633     bool is_p_point = p_.is_point();
2634     bool is_q_point = q_.is_point();
2635     bool is_r_point = r_.is_point();
2636 
2637     unsigned int numpts_in_pqr =
2638       ((is_p_point)? 1 : 0) +
2639       ((is_q_point)? 1 : 0) +
2640       ((is_r_point)? 1 : 0)  ;
2641 
2642     bool is_p_tsrc(false);
2643     CGAL_assertion_code( bool has_p_endp_tsrc(false); )
2644     if ( is_p_point ) {
2645       if ( same_points(p_, t.source_site()) ) {
2646         is_p_tsrc = true;
2647       }
2648     }
2649 #ifndef CGAL_NO_ASSERTIONS
2650     else { // p is segment
2651       if (same_points(p_.source_site(), t.source_site()) ||
2652           same_points(p_.target_site(), t.source_site())   ) {
2653         has_p_endp_tsrc = true;
2654       }
2655     }
2656 #endif
2657 
2658     bool is_q_tsrc(false);
2659     CGAL_assertion_code( bool has_q_endp_tsrc(false); )
2660     if ( is_q_point ) {
2661       if ( same_points(q_, t.source_site()) ) {
2662         is_q_tsrc = true;
2663       }
2664     }
2665 #ifndef CGAL_NO_ASSERTIONS
2666     else { // q is segment
2667       if (same_points(q_.source_site(), t.source_site()) ||
2668           same_points(q_.target_site(), t.source_site())   ) {
2669         has_q_endp_tsrc = true;
2670       }
2671     }
2672 #endif
2673 
2674     bool is_r_tsrc(false);
2675     CGAL_assertion_code( bool has_r_endp_tsrc(false); )
2676     if ( is_r_point ) {
2677       if ( same_points(r_, t.source_site()) ) {
2678         is_r_tsrc = true;
2679       }
2680     }
2681 #ifndef CGAL_NO_ASSERTIONS
2682     else { // r is segment
2683       if (same_points(r_.source_site(), t.source_site()) ||
2684           same_points(r_.target_site(), t.source_site())   ) {
2685         has_r_endp_tsrc = true;
2686       }
2687     }
2688 #endif
2689 
2690 #ifndef CGAL_NO_ASSERTIONS
2691     unsigned int num_common_endp_tsrc =
2692       ((has_p_endp_tsrc)? 1 : 0) +
2693       ((has_q_endp_tsrc)? 1 : 0) +
2694       ((has_r_endp_tsrc)? 1 : 0)  ;
2695     CGAL_USE(num_common_endp_tsrc);
2696     CGAL_SDG_DEBUG(
2697     std::cout << "debug num_common_endp_tsrc="
2698       << num_common_endp_tsrc << std::endl;
2699     );
2700 #endif
2701 
2702     unsigned int numendpts_of_t = 0;
2703 
2704     Sign d1, d2;
2705     if ( is_p_tsrc || is_q_tsrc || is_r_tsrc ) {
2706       d1 = ZERO;
2707       ++numendpts_of_t;
2708     } else {
2709       d1 = incircle_p(t.source_site());
2710     }
2711 
2712     CGAL_SDG_DEBUG(std::cout << "debug incircle_s_no_easy d1="
2713         << d1 << " with tsrc=" << t.source_site() << std::endl;);
2714 
2715     if ( d1 == NEGATIVE ) { return NEGATIVE; }
2716 
2717 
2718     bool is_p_ttrg(false);
2719     CGAL_assertion_code( bool has_p_endp_ttrg(false); )
2720     if ( is_p_point ) {
2721       if ( same_points(p_, t.target_site()) ) {
2722         is_p_ttrg = true;
2723       }
2724     }
2725 #ifndef CGAL_NO_ASSERTIONS
2726     else { // p is segment
2727       if (same_points(p_.source_site(), t.target_site()) ||
2728           same_points(p_.target_site(), t.target_site())   ) {
2729         has_p_endp_ttrg = true;
2730       }
2731     }
2732 #endif
2733 
2734     bool is_q_ttrg(false);
2735     CGAL_assertion_code( bool has_q_endp_ttrg(false); )
2736     if ( is_q_point ) {
2737       if ( same_points(q_, t.target_site()) ) {
2738         is_q_ttrg = true;
2739       }
2740     }
2741 #ifndef CGAL_NO_ASSERTIONS
2742     else { // q is segment
2743       if (same_points(q_.source_site(), t.target_site()) ||
2744           same_points(q_.target_site(), t.target_site())   ) {
2745         has_q_endp_ttrg = true;
2746       }
2747     }
2748 #endif
2749 
2750     bool is_r_ttrg(false);
2751     CGAL_assertion_code( bool has_r_endp_ttrg(false); )
2752     if ( is_r_point ) {
2753       if ( same_points(r_, t.target_site()) ) {
2754         is_r_ttrg = true;
2755       }
2756     }
2757 #ifndef CGAL_NO_ASSERTIONS
2758     else { // r is segment
2759       if (same_points(r_.source_site(), t.target_site()) ||
2760           same_points(r_.target_site(), t.target_site())   ) {
2761         has_r_endp_ttrg = true;
2762       }
2763     }
2764 #endif
2765 
2766 #ifndef CGAL_NO_ASSERTIONS
2767     unsigned int num_common_endp_ttrg =
2768       ((has_p_endp_ttrg)? 1 : 0) +
2769       ((has_q_endp_ttrg)? 1 : 0) +
2770       ((has_r_endp_ttrg)? 1 : 0)  ;
2771     CGAL_USE(num_common_endp_ttrg);
2772     CGAL_SDG_DEBUG(
2773     std::cout << "debug num_common_endp_ttrg="
2774       << num_common_endp_ttrg << std::endl;
2775     );
2776 #endif
2777 
2778     if ( is_p_ttrg || is_q_ttrg || is_r_ttrg ) {
2779       d2 = ZERO;
2780       ++numendpts_of_t;
2781     } else {
2782       d2 = incircle_p(t.target_site());
2783     }
2784     if ( d2 == NEGATIVE ) { return NEGATIVE; }
2785 
2786     CGAL_SDG_DEBUG(std::cout << "debug incircle_s_no_easy d2="
2787         << d2 << std::endl;);
2788 
2789     CGAL_assertion(numendpts_of_t < 2);
2790 
2791     CGAL_SDG_DEBUG(std::cout << "debug incircle_s_no_easy numendpts_of_t= "
2792       << numendpts_of_t << std::endl;);
2793 
2794     compute_v_if_not_computed();
2795 
2796     if (numendpts_of_t > 0) {
2797       bool is_t_horizontal = is_site_horizontal(t);
2798       bool is_t_vertical   = is_site_vertical(t);
2799 
2800       if (is_t_horizontal || is_t_vertical) {
2801         CGAL_assertion(numendpts_of_t == 1);
2802 
2803         // set endp to endpoint in {p,q,r}
2804         Site_2 endp;
2805         if ( is_p_tsrc || is_q_tsrc || is_r_tsrc ) {
2806           endp = t.source_site();
2807         } else {
2808           endp = t.target_site();
2809         }
2810 
2811         // numothers will be the number of segments
2812         // in {p,q,r} that have endp as an endpoint
2813         unsigned int numothers = 0;
2814 
2815         // a possible segment in {p,q,r} which has endpoint endp
2816         Site_2 other;
2817 
2818         // if there is a segment in {p,q,r}, try its endpoints
2819         if (numpts_in_pqr < 3) {
2820 
2821           if ((! is_p_point) && is_endpoint_of(endp, p_)) {
2822             numothers++;
2823             other = p_;
2824           }
2825 
2826           if ((! is_q_point) && is_endpoint_of(endp, q_)) {
2827             numothers++;
2828             other = q_;
2829           }
2830 
2831           if ((! is_r_point) && is_endpoint_of(endp, r_)) {
2832             numothers++;
2833             other = r_;
2834           }
2835 
2836         } // end of case: numpts_in_pqr < 3
2837 
2838         CGAL_assertion(numothers < 2);
2839 
2840         if (numothers == 1) {
2841           bool is_other_horizontal = is_site_horizontal(other);
2842           bool is_other_vertical = is_site_vertical(other);
2843 
2844           if ((is_t_horizontal && is_other_horizontal) ||
2845               (is_t_vertical && is_other_vertical)       ) {
2846             return POSITIVE;
2847           }
2848         } else {
2849           CGAL_assertion(numothers == 0);
2850           Point_2 vv(ux_, uy_,uz_);
2851 
2852           Comparison_result ptcmpxve =
2853             CGAL::compare(vv.x(), endp.point().x());
2854           Comparison_result ptcmpyve =
2855             CGAL::compare(vv.y(), endp.point().y());
2856 
2857           CGAL_SDG_DEBUG(std::cout << "debug vv = " << vv << std::endl;);
2858 
2859           if ( ( (ptcmpxve == EQUAL) && is_t_horizontal ) ||
2860                ( (ptcmpyve == EQUAL) && is_t_vertical   )    ) {
2861             return ZERO;
2862           }
2863 
2864         } // end of case numothers == 0
2865       }  // endif (is_t_horizontal || is_t_vertical)
2866     } // endif ((numendpts_of_t > 0) && (numpts_in_pqr < 3))
2867 
2868     bool is_tsrc_endp_of_p (false);
2869     bool is_tsrc_endp_of_q (false);
2870     bool is_tsrc_endp_of_r (false);
2871     bool is_ttrg_endp_of_p (false);
2872     bool is_ttrg_endp_of_q (false);
2873     bool is_ttrg_endp_of_r (false);
2874 
2875     if (! is_p_point) {
2876       is_tsrc_endp_of_p = same_points(t.source_site(), p_.source_site())
2877                        || same_points(t.source_site(), p_.target_site());
2878       is_ttrg_endp_of_p = same_points(t.target_site(), p_.source_site())
2879                        || same_points(t.target_site(), p_.target_site());
2880     }
2881     if (! is_q_point) {
2882       is_tsrc_endp_of_q = same_points(t.source_site(), q_.source_site())
2883                        || same_points(t.source_site(), q_.target_site());
2884       is_ttrg_endp_of_q = same_points(t.target_site(), q_.source_site())
2885                        || same_points(t.target_site(), q_.target_site());
2886     }
2887     if (! is_r_point) {
2888       is_tsrc_endp_of_r = same_points(t.source_site(), r_.source_site())
2889                        || same_points(t.source_site(), r_.target_site());
2890       is_ttrg_endp_of_r = same_points(t.target_site(), r_.source_site())
2891                        || same_points(t.target_site(), r_.target_site());
2892     }
2893 
2894     if (is_tsrc_endp_of_p && is_tsrc_endp_of_q) {
2895       if (test_star(t.source_site(), p_, q_, t)) {
2896         return NEGATIVE;
2897       }
2898     }
2899     if (is_ttrg_endp_of_p && is_ttrg_endp_of_q) {
2900       if (test_star(t.target_site(), p_, q_, t)) {
2901         return NEGATIVE;
2902       }
2903     }
2904     if (is_tsrc_endp_of_q && is_tsrc_endp_of_r) {
2905       if (test_star(t.source_site(), q_, r_, t)) {
2906         return NEGATIVE;
2907       }
2908     }
2909     if (is_ttrg_endp_of_q && is_ttrg_endp_of_r) {
2910       if (test_star(t.target_site(), q_, r_, t)) {
2911         return NEGATIVE;
2912       }
2913     }
2914     if (is_tsrc_endp_of_r && is_tsrc_endp_of_p) {
2915       if (test_star(t.source_site(), r_, p_, t)) {
2916         return NEGATIVE;
2917       }
2918     }
2919     if (is_ttrg_endp_of_r && is_ttrg_endp_of_p) {
2920       if (test_star(t.target_site(), r_, p_, t)) {
2921         return NEGATIVE;
2922       }
2923     }
2924 
2925     Line_2 l = compute_supporting_line(t.supporting_site());
2926     Sign sl = incircle(l, type);
2927 
2928     CGAL_SDG_DEBUG(std::cout
2929         << "debug vring incircle_s_no_easy: incircle l returned "
2930         << sl << std::endl;);
2931 
2932     if ( sl == POSITIVE ) { return sl; }
2933 
2934     CGAL_SDG_DEBUG(std::cout << "debug incircle_s_no_easy sl=" << sl <<
2935       " d1=" << d1 << " d2=" << d2 << std::endl;);
2936 
2937     CGAL_SDG_DEBUG(std::cout << "debug vring numpts_in_pqr="
2938         << numpts_in_pqr << std::endl;);
2939 
2940     // philaris: here we have a serious change related to L2
2941     if ( sl == ZERO && (d1 == ZERO || d2 == ZERO) ) {
2942 
2943       // if some site in {p,q,r} is a point and it is also:
2944       // an endpoint of t and an endpoint of another site in {p,q,r}
2945 
2946       Site_2 sqpnt, other_t, other_seg;
2947 
2948       if (compute_helper(p_, q_, r_, t, sqpnt, other_t, other_seg)) {
2949 
2950         CGAL_assertion(sqpnt.is_point());
2951         CGAL_assertion(other_t.is_point());
2952         CGAL_assertion(other_seg.is_point());
2953 
2954         Point_2 vv (ux_, uy_, uz_);
2955 
2956         CGAL_SDG_DEBUG(std::cout << "debug vring incircle_s_no_easy "
2957             << "compute_helper true, "
2958             << "  vv=" << vv << "  sqpnt= " << sqpnt
2959             << "  other_t=" << other_t
2960             << "  other_seg=" << other_seg
2961             << std::endl;);
2962 
2963         Line_2 lvs =
2964           compute_line_from_to(vv, sqpnt.point());
2965 
2966         Oriented_side os_t =
2967           oriented_side_of_line(lvs, other_t.point());
2968         Oriented_side os_s =
2969           oriented_side_of_line(lvs, other_seg.point());
2970 
2971         CGAL_assertion(os_s != ON_ORIENTED_BOUNDARY);
2972 
2973         if (os_t == os_s) {
2974           Line_2 lseg =
2975             compute_line_from_to(sqpnt.point(), other_seg.point());
2976 
2977           Oriented_side os_seg_vv =
2978            oriented_side_of_line(lseg, vv);
2979           Oriented_side os_seg_t =
2980            oriented_side_of_line(lseg, other_t.point());
2981 
2982           if (os_seg_t == os_seg_vv) {
2983             return NEGATIVE;
2984           } else {
2985             if (os_seg_t == ON_ORIENTED_BOUNDARY) {
2986               return ZERO;
2987             } else {
2988               return POSITIVE;
2989             }
2990           }
2991         } // end of case: os_t == os_s
2992       } // end of case where
2993 
2994       return ZERO;
2995     }
2996 
2997     Oriented_side os1 = oriented_side_linf(l, t.source(), type);
2998     Oriented_side os2 = oriented_side_linf(l, t.target(), type);
2999 
3000     CGAL_SDG_DEBUG(std::cout << "debug vring incircle_s_no_easy: os1="
3001         << os1 << " os2="
3002         << os2 << std::endl;);
3003 
3004     if ( sl == ZERO ) {
3005       if (os1 == ON_ORIENTED_BOUNDARY || os2 == ON_ORIENTED_BOUNDARY) {
3006         return ZERO;
3007       }
3008       return ( os1 == os2 ) ? POSITIVE : ZERO;
3009     }
3010 
3011     CGAL_SDG_DEBUG(std::cout << "debug incircle_s_no_easy non-zero sl: os1="
3012       << os1 << " os2=" << os2 << std::endl;);
3013 
3014     return (os1 == os2) ? POSITIVE : NEGATIVE;
3015   }
3016 
3017   inline
3018   bool
compute_helper(const Site_2 & p,const Site_2 & q,const Site_2 & r,const Site_2 & t,Site_2 & sqpnt,Site_2 & other_of_t,Site_2 & other_of_seg)3019   compute_helper(const Site_2& p, const Site_2& q, const Site_2& r,
3020       const Site_2& t,
3021       Site_2& sqpnt, Site_2& other_of_t, Site_2& other_of_seg)
3022   const
3023   {
3024     CGAL_assertion(t.is_segment());
3025 
3026     const bool is_p_point = p.is_point();
3027     const bool is_q_point = q.is_point();
3028     const bool is_r_point = r.is_point();
3029 
3030     const unsigned int numpts =
3031       ((is_p_point)? 1 : 0) +
3032       ((is_q_point)? 1 : 0) +
3033       ((is_r_point)? 1 : 0)  ;
3034 
3035     CGAL_SDG_DEBUG(std::cout << "debug vring compute_helper #pts="
3036         << numpts << std::endl;);
3037 
3038     if (numpts == 3) {
3039       return false;
3040     }
3041 
3042     // here and on, there are 1 or 2 points in {p,q,r}
3043 
3044 
3045     bool is_p_tsrc(false);
3046     bool is_p_endp_of_t(false);
3047     if (is_p_point) {
3048       is_p_tsrc = same_points(p, t.source_site());
3049       const bool is_p_ttrg = same_points(p, t.target_site());
3050       is_p_endp_of_t = is_p_tsrc || is_p_ttrg;
3051 
3052       if (is_p_endp_of_t) {
3053         sqpnt = p;
3054       }
3055     }
3056 
3057     bool is_q_tsrc(false);
3058     bool is_q_endp_of_t(false);
3059     if (is_q_point) {
3060       is_q_tsrc = same_points(q, t.source_site());
3061       const bool is_q_ttrg = same_points(q, t.target_site());
3062       is_q_endp_of_t = is_q_tsrc || is_q_ttrg;
3063       if (is_q_endp_of_t) {
3064         sqpnt = q;
3065       }
3066     }
3067 
3068     bool is_r_tsrc(false);
3069     bool is_r_endp_of_t(false);
3070 
3071     if (is_r_point) {
3072       is_r_tsrc = same_points(r, t.source_site());
3073       const bool is_r_ttrg = same_points(r, t.target_site());
3074       is_r_endp_of_t = is_r_tsrc || is_r_ttrg;
3075       if (is_r_endp_of_t) {
3076         sqpnt = r;
3077       }
3078     }
3079 
3080     const unsigned int numendpts_of_t =
3081       ((is_p_endp_of_t)? 1 : 0) +
3082       ((is_q_endp_of_t)? 1 : 0) +
3083       ((is_r_endp_of_t)? 1 : 0)  ;
3084 
3085     CGAL_SDG_DEBUG(std::cout << "debug compute_helper #endpts_of_t=" <<
3086       numendpts_of_t << std::endl;);
3087 
3088     if (numendpts_of_t == 0) {
3089 
3090       bool have_common_p_tsrc(false),
3091            have_common_p_ttrg(false),
3092            have_common_p_t(false);
3093 
3094       if (! is_p_point) {
3095         CGAL_assertion( ! same_segments(p, t) );
3096         const bool is_psrc_tsrc =
3097           same_points(p.source_site(), t.source_site());
3098         const bool is_ptrg_tsrc =
3099           same_points(p.target_site(), t.source_site());
3100         const bool is_psrc_ttrg =
3101           same_points(p.source_site(), t.target_site());
3102         const bool is_ptrg_ttrg =
3103           same_points(p.target_site(), t.target_site());
3104         have_common_p_tsrc = is_psrc_tsrc || is_ptrg_tsrc;
3105         have_common_p_ttrg = is_psrc_ttrg || is_ptrg_ttrg;
3106         have_common_p_t = have_common_p_tsrc || have_common_p_ttrg;
3107       }
3108 
3109       bool have_common_q_tsrc(false),
3110            have_common_q_ttrg(false),
3111            have_common_q_t(false);
3112 
3113       if (! is_q_point) {
3114         CGAL_assertion( ! same_segments(q, t) );
3115         const bool is_qsrc_tsrc = same_points(q.source_site(), t.source_site());
3116         const bool is_qtrg_tsrc = same_points(q.target_site(), t.source_site());
3117         const bool is_qsrc_ttrg = same_points(q.source_site(), t.target_site());
3118         const bool is_qtrg_ttrg = same_points(q.target_site(), t.target_site());
3119         have_common_q_tsrc = is_qsrc_tsrc || is_qtrg_tsrc;
3120         have_common_q_ttrg = is_qsrc_ttrg || is_qtrg_ttrg;
3121         have_common_q_t = have_common_q_tsrc || have_common_q_ttrg;
3122       }
3123 
3124       bool have_common_r_tsrc(false),
3125            have_common_r_ttrg(false),
3126            have_common_r_t(false);
3127 
3128       if (! is_r_point) {
3129         CGAL_assertion( ! same_segments(r, t) );
3130         const bool is_rsrc_tsrc = same_points(r.source_site(), t.source_site());
3131         const bool is_rtrg_tsrc = same_points(r.target_site(), t.source_site());
3132         const bool is_rsrc_ttrg = same_points(r.source_site(), t.target_site());
3133         const bool is_rtrg_ttrg = same_points(r.target_site(), t.target_site());
3134         have_common_r_tsrc = is_rsrc_tsrc || is_rtrg_tsrc;
3135         have_common_r_ttrg = is_rsrc_ttrg || is_rtrg_ttrg;
3136         have_common_r_t = have_common_r_tsrc || have_common_r_ttrg;
3137       }
3138 
3139       const unsigned int numcommon =
3140       ((have_common_p_t)? 1 : 0) +
3141       ((have_common_q_t)? 1 : 0) +
3142       ((have_common_r_t)? 1 : 0)  ;
3143 
3144       CGAL_SDG_DEBUG(std::cout << "debug compute_helper #numcommon=" <<
3145         numcommon << std::endl;);
3146 
3147       CGAL_assertion(numcommon < 3);
3148 
3149       if (numcommon < 2) {
3150         return false;
3151       }
3152 
3153       // here, numcommon == 2
3154 
3155       const unsigned int numcommon_tsrc =
3156       ((have_common_p_tsrc)? 1 : 0) +
3157       ((have_common_q_tsrc)? 1 : 0) +
3158       ((have_common_r_tsrc)? 1 : 0)  ;
3159 
3160       const unsigned int numcommon_ttrg =
3161       ((have_common_p_ttrg)? 1 : 0) +
3162       ((have_common_q_ttrg)? 1 : 0) +
3163       ((have_common_r_ttrg)? 1 : 0)  ;
3164 
3165       CGAL_assertion( numcommon_tsrc + numcommon_ttrg == 2 );
3166 
3167       if ( numcommon_tsrc == numcommon_ttrg )  { // both equal 1
3168         return false;
3169       }
3170 
3171       // here either numcommon_tsrc==2 or numcommon_ttrg==2
3172 
3173       if (numcommon_tsrc > 0) {
3174         // here, numcommon_tsrc == 2
3175         sqpnt = t.source_site();
3176         other_of_t = t.target_site();
3177       } else {
3178         // here, numcommon_ttrg == 2
3179         sqpnt = t.target_site();
3180         other_of_t = t.source_site();
3181       }
3182 
3183       if (have_common_p_t && have_common_q_t) {
3184         compute_helper_two_seg(p, q, sqpnt, other_of_seg);
3185       } else if (have_common_q_t && have_common_r_t) {
3186         compute_helper_two_seg(q, r, sqpnt, other_of_seg);
3187       } else if (have_common_r_t && have_common_p_t) {
3188         compute_helper_two_seg(r, p, sqpnt, other_of_seg);
3189       } else {
3190         CGAL_assertion(false);
3191       }
3192 
3193       return true;
3194     }
3195 
3196     // philaris: tocheck
3197     CGAL_assertion( numendpts_of_t == 1 );
3198 
3199     if (is_p_tsrc || is_q_tsrc || is_r_tsrc) {
3200       other_of_t = t.target_site();
3201     } else {
3202       other_of_t = t.source_site();
3203     }
3204 
3205     if (is_p_endp_of_t) {
3206       if (q.is_segment()) {
3207         bool is_p_qsrc = same_points(p, q.source_site());
3208         bool is_p_qtrg = same_points(p, q.target_site());
3209         if (is_p_qsrc || is_p_qtrg) {
3210           other_of_seg = is_p_qsrc ? q.target_site() : q.source_site();
3211           return true;
3212         }
3213       }
3214       if (r.is_segment()) {
3215         bool is_p_rsrc = same_points(p, r.source_site());
3216         bool is_p_rtrg = same_points(p, r.target_site());
3217         if (is_p_rsrc || is_p_rtrg) {
3218           other_of_seg = is_p_rsrc ? r.target_site() : r.source_site();
3219           return true;
3220         }
3221       }
3222 
3223     } // end of case: is_p_endp_of_t
3224 
3225     if (is_q_endp_of_t) {
3226       if (r.is_segment()) {
3227         bool is_q_rsrc = same_points(q, r.source_site());
3228         bool is_q_rtrg = same_points(q, r.target_site());
3229         if (is_q_rsrc || is_q_rtrg) {
3230           other_of_seg = is_q_rsrc ? r.target_site() : r.source_site();
3231           return true;
3232         }
3233       }
3234 
3235       if (p.is_segment()) {
3236         bool is_q_psrc = same_points(q, p.source_site());
3237         bool is_q_ptrg = same_points(q, p.target_site());
3238         if (is_q_psrc || is_q_ptrg) {
3239           other_of_seg = is_q_psrc ? p.target_site() : p.source_site();
3240           return true;
3241         }
3242       }
3243 
3244     } // end of case: is_q_endp_of_t
3245 
3246     if (is_r_endp_of_t) {
3247       if (p.is_segment()) {
3248         bool is_r_psrc = same_points(r, p.source_site());
3249         bool is_r_ptrg = same_points(r, p.target_site());
3250         if (is_r_psrc || is_r_ptrg) {
3251           other_of_seg = is_r_psrc ? p.target_site() : p.source_site();
3252           return true;
3253         }
3254       }
3255 
3256       if (q.is_segment()) {
3257         bool is_r_qsrc = same_points(r, q.source_site());
3258         bool is_r_qtrg = same_points(r, q.target_site());
3259         if (is_r_qsrc || is_r_qtrg) {
3260           other_of_seg = is_r_qsrc ? q.target_site() : q.source_site();
3261           return true;
3262         }
3263       }
3264 
3265     } // end of case: is_r_endp_of_t
3266 
3267     CGAL_SDG_DEBUG(std::cout << "debug compute_helper about to return false"
3268         << std::endl;);
3269     return false;
3270 
3271   }
3272 
3273   inline
3274   void
compute_helper_two_seg(const Site_2 & a,const Site_2 & b,const Site_2 & common_site,Site_2 & other_of_seg)3275   compute_helper_two_seg(
3276       const Site_2& a, const Site_2& b,
3277       const Site_2& common_site, Site_2& other_of_seg)
3278   const
3279   {
3280     CGAL_assertion(a.is_segment());
3281     CGAL_assertion(b.is_segment());
3282 
3283     CGAL_SDG_DEBUG(std::cout
3284         << "debug compute_helper_two_seg entering with "
3285         << a << " and " << b << " having common "
3286         << common_site << std::endl;);
3287 
3288     if (is_site_h_or_v(a)) {
3289       if ( same_points(common_site, b.source_site()) ) {
3290         other_of_seg = b.target_site();
3291       } else {
3292         other_of_seg = b.source_site();
3293       }
3294     } else {
3295       CGAL_assertion(is_site_h_or_v(b));
3296 
3297       if ( same_points(common_site, a.source_site()) ) {
3298         other_of_seg = a.target_site();
3299       } else {
3300         other_of_seg = a.source_site();
3301       }
3302 
3303     }
3304   } // end of compute_helper_two_seg
3305 
3306 
3307   //--------------------------------------------------------------------------
3308 
3309 
3310 
incircle_s(const Site_2 & t)3311   Sign incircle_s(const Site_2& t) const
3312   {
3313     CGAL_precondition( t.is_segment() );
3314 
3315     CGAL_SDG_DEBUG(std::cout << "debug incircle_s (pqrt) = "
3316       << "(" << p_ << ") (" << q_ << ") (" << r_ << ") "
3317       << "(" << t << ")" << std::endl;);
3318 
3319     if ( is_degenerate_Voronoi_circle() ) {
3320       // case 1: the new segment is not adjacent to the center of the
3321       //         degenerate Voronoi circle
3322       if (  !same_points( p_ref(), t.source_site() ) &&
3323             !same_points( p_ref(), t.target_site() )  ) {
3324         return POSITIVE;
3325       }
3326 
3327       CGAL_assertion( v_type == PSS );
3328 
3329       if ( p_.is_segment() &&
3330            same_segments(p_.supporting_site(),
3331                          t.supporting_site()) ) {
3332         return ZERO;
3333       }
3334 
3335       if ( q_.is_segment() &&
3336            same_segments(q_.supporting_site(),
3337                          t.supporting_site()) ) {
3338         return ZERO;
3339       }
3340 
3341       if ( r_.is_segment() &&
3342            same_segments(r_.supporting_site(),
3343                          t.supporting_site()) ) {
3344         return ZERO;
3345       }
3346 
3347       Site_2 pr;
3348       Site_2 sp, sq;
3349       if ( p_.is_point() ) {
3350         CGAL_assertion( q_.is_segment() && r_.is_segment() );
3351         pr = p_;
3352         sp = q_;
3353         sq = r_;
3354       } else if ( q_.is_point() ) {
3355         CGAL_assertion( r_.is_segment() && p_.is_segment() );
3356         pr = q_;
3357         sp = r_;
3358         sq = p_;
3359       } else {
3360         CGAL_assertion( p_.is_segment() && q_.is_segment() );
3361         pr = r_;
3362         sp = p_;
3363         sq = q_;
3364       }
3365 
3366       Point_2 pq = sq.source(), pp = sp.source(), pt = t.source();
3367 
3368       if ( same_points(sp.source_site(), pr) ) { pp = sp.target(); }
3369       if ( same_points(sq.source_site(), pr) ) { pq = sq.target(); }
3370       if ( same_points( t.source_site(), pr) ) { pt =  t.target(); }
3371 
3372       Point_2 pr_ = pr.point();
3373 
3374       if ( CGAL::orientation(pr_, pp, pt) == LEFT_TURN &&
3375            CGAL::orientation(pr_, pq, pt) == RIGHT_TURN ) {
3376         return NEGATIVE;
3377       }
3378       return ZERO;
3379     } // if ( is_degenerate_Voronoi_circle() )
3380 
3381     Sign s(ZERO);
3382     switch ( v_type ) {
3383     case PPP:
3384       s = incircle_s(t, PPP_Type());
3385       break;
3386     case PPS:
3387       s = incircle_s(t, PPS_Type());
3388       break;
3389     case PSS:
3390       s = incircle_s(t, PSS_Type());
3391       break;
3392     case SSS:
3393       s = incircle_s(t, SSS_Type());
3394       break;
3395     }
3396 
3397     return s;
3398   }
3399 
incircle_s_no_easy(const Site_2 & t)3400   Sign incircle_s_no_easy(const Site_2& t) const
3401   {
3402     Sign s(ZERO);
3403     switch ( v_type ) {
3404     case PPP:
3405       s = incircle_s_no_easy(t, PPP_Type());
3406       break;
3407     case PPS:
3408       s = incircle_s_no_easy(t, PPS_Type());
3409       break;
3410     case PSS:
3411       s = incircle_s_no_easy(t, PSS_Type());
3412       break;
3413     case SSS:
3414       s = incircle_s_no_easy(t, SSS_Type());
3415       break;
3416     }
3417 
3418     return s;
3419   }
3420 
3421   //--------------------------------------------------------------------------
3422   //  subpredicates for the incircle test
3423   //--------------------------------------------------------------------------
3424 
3425 
3426 public:
is_degenerate_Voronoi_circle()3427   bool is_degenerate_Voronoi_circle() const
3428   {
3429     if ( v_type != PSS ) { return false; }
3430 
3431     if ( p_.is_point() ) {
3432       return ( is_endpoint_of(p_, q_) && is_endpoint_of(p_, r_) );
3433     } else if ( q_.is_point() ) {
3434       return ( is_endpoint_of(q_, p_) && is_endpoint_of(q_, r_) );
3435     } else {
3436       CGAL_assertion( r_.is_point() );
3437       return ( is_endpoint_of(r_, p_) && is_endpoint_of(r_, q_) );
3438     }
3439   }
3440 
3441 
3442   //--------------------------------------------------------------------------
3443 
3444 private:
3445 
3446   //--------------------------------------------------------------------------
3447   //  the reference point (valid if v_type != SSS)
3448   //--------------------------------------------------------------------------
3449 
p_ref()3450   Site_2 p_ref() const
3451   {
3452     CGAL_precondition ( v_type != SSS );
3453 
3454 
3455     if ( v_type == PPS ) {
3456       //CGAL_SDG_DEBUG(std::cout << "debug p_ref pps_idx="
3457       //    << pps_idx << std::endl;);
3458 
3459       if ( pps_idx == 0 ) {
3460         CGAL_assertion( p_.is_point());
3461         return p_;
3462       }
3463 
3464       if ( pps_idx == 1 ) {
3465         CGAL_assertion( q_.is_point());
3466         return q_;
3467       }
3468 
3469       //CGAL_SDG_DEBUG(std::cout << "debug p_ref about to return r="
3470       //    << r_ << std::endl;);
3471 
3472       CGAL_assertion( r_.is_point());
3473       return r_;
3474     }
3475 
3476     if ( p_.is_point() ) {
3477       return p_;
3478     } else if ( q_.is_point() ) {
3479       return q_;
3480     } else {
3481       CGAL_assertion( r_.is_point() );
3482       return r_;
3483     }
3484   }
3485 
3486 
3487 public:
3488   //--------------------------------------------------------------------------
3489   //--------------------------------------------------------------------------
3490   //                           access methods
3491   //--------------------------------------------------------------------------
3492   //--------------------------------------------------------------------------
3493 
x(Integral_domain_without_division_tag)3494   inline FT x(Integral_domain_without_division_tag) const {
3495     return FT(CGAL::to_double(hx()) / CGAL::to_double(hw()));
3496   }
y(Integral_domain_without_division_tag)3497   inline FT y(Integral_domain_without_division_tag) const {
3498     return FT(CGAL::to_double(hy()) / CGAL::to_double(hw()));
3499   }
3500 
x(Field_tag)3501   inline FT x(Field_tag) const { return hx() / hw(); }
y(Field_tag)3502   inline FT y(Field_tag) const { return hy() / hw(); }
3503 
x()3504   inline FT x() const {
3505       typedef Algebraic_structure_traits<FT> AST;
3506       return x(typename AST::Algebraic_category());
3507   }
3508 
y()3509   inline FT y() const {
3510       typedef Algebraic_structure_traits<FT> AST;
3511       return y(typename AST::Algebraic_category());
3512   }
3513 
hx()3514   FT hx() const {
3515     CGAL_assertion( is_v_computed );
3516     return ux_;
3517   }
3518 
hy()3519   FT hy() const {
3520     CGAL_assertion( is_v_computed );
3521     return uy_;
3522   }
3523 
hw()3524   FT hw() const {
3525     CGAL_assertion( is_v_computed );
3526     return uz_;
3527   }
3528 
radius()3529   FT radius() const {
3530     switch (v_type) {
3531     case PPP:    case PPS:    case PSS:
3532       {
3533         Point_2 pref = p_ref().point();
3534         //FT absdx = CGAL::abs(x() - pref.x());
3535         //FT absdy = CGAL::abs(y() - pref.y());
3536         return (CGAL::max)( CGAL::abs(x() - pref.x()),
3537                             CGAL::abs(y() - pref.y()) );
3538       }
3539       break;
3540     case SSS:
3541       {
3542         Line_2 l = compute_supporting_line(p_.supporting_site());
3543         Homogeneous_point_2 q = compute_linf_projection_hom(l, point());
3544 
3545         FT dx = CGAL::abs(x() - q.x());
3546         FT dy = CGAL::abs(y() - q.y());
3547         return (CGAL::max)(dx, dy);
3548       }
3549       break;
3550     default:
3551       return FT(0);
3552     }
3553   }
3554 
point()3555   Point_2 point() const {
3556     if ( is_degenerate_Voronoi_circle() ) {
3557       return degenerate_point();
3558     }
3559     compute_v_if_not_computed();
3560     return Point_2(x(), y());
3561   }
3562 
3563 
degenerate_point()3564   Point_2 degenerate_point() const
3565   {
3566     CGAL_precondition( is_degenerate_Voronoi_circle() );
3567     return p_ref().point();
3568   }
3569 
3570   // philaris: the circle is in fact an Iso_rectangle_2
circle()3571   typename K::Iso_rectangle_2 circle() const
3572   {
3573     typedef typename K::Iso_rectangle_2  Iso_rectangle_2;
3574     Point_2 pleftbot (point().x()-radius(), point().y()-radius());
3575     Point_2 prghttop (point().x()+radius(), point().y()+radius());
3576     return Iso_rectangle_2(pleftbot, prghttop);
3577   }
3578 
type()3579   vertex_t type() const { return v_type; }
3580 
3581 public:
Voronoi_vertex_ring_C2(const Site_2 & p,const Site_2 & q,const Site_2 & r)3582   Voronoi_vertex_ring_C2(const Site_2& p,
3583                          const Site_2& q,
3584                          const Site_2& r)
3585     : p_(p), q_(q), r_(r), is_v_computed(false)
3586   {
3587     CGAL_SDG_DEBUG(std::cout << "Voronoi_vertex_ring_C2()" << std::endl;);
3588     analyze_vertex(p, q, r);
3589   }
3590 
3591   //--------------------------------------------------------------------------
3592 
incircle(const Site_2 & t)3593   Sign incircle(const Site_2& t) const
3594   {
3595     Sign s;
3596 
3597     CGAL_SDG_DEBUG(std::cout << "debug ring incircle t=" << t << std::endl;);
3598 
3599     if ( t.is_point() ) {
3600       s = incircle_p(t);
3601     } else {
3602       CGAL_SDG_DEBUG(std::cout << "debug about to run incircle_s with t="
3603         << t << std::endl;);
3604       s = incircle_s(t);
3605     }
3606 
3607     return s;
3608   }
3609 
incircle_no_easy(const Site_2 & t)3610   Sign incircle_no_easy(const Site_2& t) const
3611   {
3612     Sign s;
3613 
3614     if ( t.is_point() ) {
3615       s = incircle_p_no_easy(t);
3616     } else {
3617       s = incircle_s_no_easy(t);
3618     }
3619 
3620     return s;
3621   }
3622 
3623   //--------------------------------------------------------------------------
3624 
3625 
orientation(const Line_2 & l)3626   Orientation orientation(const Line_2& l) const
3627   {
3628     CGAL_assertion( is_v_computed);
3629     Orientation o(COLLINEAR);
3630     switch ( v_type ) {
3631     case PPP:
3632       o = orientation(l, PPP_Type());
3633       break;
3634     case PPS:
3635       o = orientation(l, PPS_Type());
3636       break;
3637     case PSS:
3638       o = orientation(l, PSS_Type());
3639       break;
3640     case SSS:
3641       o = orientation(l, SSS_Type());
3642       break;
3643     }
3644 
3645     return o;
3646   }
3647 
oriented_side(const Line_2 & l)3648   Oriented_side oriented_side(const Line_2& l) const
3649   {
3650     compute_v_if_not_computed();
3651     Orientation o = orientation(l);
3652 
3653     if ( o == COLLINEAR ) { return ON_ORIENTED_BOUNDARY; }
3654     return ( o == LEFT_TURN ) ? ON_POSITIVE_SIDE : ON_NEGATIVE_SIDE;
3655   }
3656 
3657 
3658   // L_inf refinement
3659   inline
3660   Comparison_result
linf_refine(const Line_2 & l,Homogeneous_point_2 & lrefhp)3661   linf_refine( const Line_2& l, Homogeneous_point_2& lrefhp ) const
3662   {
3663     CGAL_assertion( is_v_computed );
3664     Point_2 vv ( ux_, uy_, uz_ );
3665 
3666     bool is_l_h_or_v = is_line_h_or_v(l);
3667 
3668     Comparison_result compare_p(EQUAL);
3669     Comparison_result compare_q(EQUAL);
3670     Comparison_result compare_r(EQUAL);
3671 
3672     FT difxvl = vv.x() - lrefhp.x();
3673     FT difyvl = vv.y() - lrefhp.y();
3674     FT absdifxvl = CGAL::abs(difxvl);
3675     FT absdifyvl = CGAL::abs(difyvl);
3676     Comparison_result cmplabsxy = CGAL::compare(absdifxvl, absdifyvl);
3677 
3678     // lref is corner if and only if the line is not axis-parallel
3679     CGAL_assertion( (cmplabsxy == EQUAL) == (! is_l_h_or_v) );
3680 
3681     Oriented_side oslvv (ON_ORIENTED_BOUNDARY);
3682     if ((p_.is_segment() || q_.is_segment() || r_.is_segment()) &&
3683         is_l_h_or_v) {
3684       oslvv = oriented_side_of_line(l, vv);
3685       CGAL_assertion(oslvv != ON_ORIENTED_BOUNDARY);
3686     }
3687 
3688     // The following boolean variables are used to:
3689     // compute whether lref agrees with the x and y coordinates
3690     // of some of the points in p, q, r
3691     bool corner_agree_pt_x(false);
3692     bool corner_agree_pt_y(false);
3693 
3694     if (p_.is_point()) {
3695       Point_2 pp = p_.point();
3696       FT difxvp = vv.x() - pp.x();
3697       FT difyvp = vv.y() - pp.y();
3698       FT absdifxvp = CGAL::abs(difxvp);
3699       FT absdifyvp = CGAL::abs(difyvp);
3700       Comparison_result cmppabsxy = CGAL::compare(absdifxvp, absdifyvp);
3701       if (! ( (cmplabsxy == SMALLER) && (cmppabsxy == SMALLER) ))
3702       {
3703         if (CGAL::compare(difxvl, difxvp) == EQUAL) {
3704           compare_p = CGAL::compare(absdifyvl, absdifyvp);
3705           corner_agree_pt_x = true;
3706         }
3707       }
3708       if (! ( (cmplabsxy == LARGER ) && (cmppabsxy == LARGER ) ))
3709       {
3710         if (CGAL::compare(difyvl, difyvp) == EQUAL) {
3711           CGAL_assertion(compare_p == EQUAL);
3712           compare_p = CGAL::compare(absdifxvl, absdifxvp);
3713           corner_agree_pt_y = true;
3714         }
3715       }
3716     } else {
3717       if (is_l_h_or_v) {
3718         Oriented_side oslpsrc =
3719           oriented_side_of_line(l, p_.source_site().point());
3720         Oriented_side oslptrg =
3721           oriented_side_of_line(l, p_.target_site().point());
3722         if (((oslpsrc != oslvv) && (oslptrg != oslvv)) &&
3723             ((oslpsrc != ON_ORIENTED_BOUNDARY) ||
3724              (oslptrg != ON_ORIENTED_BOUNDARY)   )         ) {
3725           compare_p = SMALLER;
3726         }
3727       }
3728     }
3729 
3730     if (q_.is_point()) {
3731       Point_2 qq = q_.point();
3732       FT difxvq = vv.x() - qq.x();
3733       FT difyvq = vv.y() - qq.y();
3734       FT absdifxvq = CGAL::abs(difxvq);
3735       FT absdifyvq = CGAL::abs(difyvq);
3736       Comparison_result cmpqabsxy = CGAL::compare(absdifxvq, absdifyvq);
3737       if (! ( (cmplabsxy == SMALLER) && (cmpqabsxy == SMALLER) ))
3738       {
3739         if (CGAL::compare(difxvl, difxvq) == EQUAL) {
3740           compare_q = CGAL::compare(absdifyvl, absdifyvq);
3741           corner_agree_pt_x = true;
3742         }
3743       }
3744       if (! ( (cmplabsxy == LARGER ) && (cmpqabsxy == LARGER ) ))
3745       {
3746         if (CGAL::compare(difyvl, difyvq) == EQUAL) {
3747           CGAL_assertion(compare_q == EQUAL);
3748           compare_q = CGAL::compare(absdifxvl, absdifxvq);
3749           corner_agree_pt_y = true;
3750         }
3751       }
3752     } else {
3753       if (is_l_h_or_v) {
3754         Oriented_side oslqsrc =
3755           oriented_side_of_line(l, q_.source_site().point());
3756         Oriented_side oslqtrg =
3757           oriented_side_of_line(l, q_.target_site().point());
3758         if (((oslqsrc != oslvv) && (oslqtrg != oslvv)) &&
3759             ((oslqsrc != ON_ORIENTED_BOUNDARY) ||
3760              (oslqtrg != ON_ORIENTED_BOUNDARY)   )         ) {
3761           compare_q = SMALLER;
3762         }
3763       }
3764     }
3765 
3766     if (r_.is_point()) {
3767       Point_2 rr = r_.point();
3768       FT difxvr = vv.x() - rr.x();
3769       FT difyvr = vv.y() - rr.y();
3770       FT absdifxvr = CGAL::abs(difxvr);
3771       FT absdifyvr = CGAL::abs(difyvr);
3772       Comparison_result cmprabsxy = CGAL::compare(absdifxvr, absdifyvr);
3773       if (! ( (cmplabsxy == SMALLER) && (cmprabsxy == SMALLER) ))
3774       {
3775         if (CGAL::compare(difxvl, difxvr) == EQUAL) {
3776           compare_r = CGAL::compare(absdifyvl, absdifyvr);
3777           corner_agree_pt_x = true;
3778         }
3779       }
3780       if (! ( (cmplabsxy == LARGER ) && (cmprabsxy == LARGER ) ))
3781       {
3782         if (CGAL::compare(difyvl, difyvr) == EQUAL) {
3783           CGAL_assertion(compare_r == EQUAL);
3784           compare_r = CGAL::compare(absdifxvl, absdifxvr);
3785           corner_agree_pt_y = true;
3786         }
3787       }
3788     } else {
3789       if (is_l_h_or_v) {
3790         Oriented_side oslrsrc =
3791           oriented_side_of_line(l, r_.source_site().point());
3792         Oriented_side oslrtrg =
3793           oriented_side_of_line(l, r_.target_site().point());
3794         if (((oslrsrc != oslvv) && (oslrtrg != oslvv)) &&
3795             ((oslrsrc != ON_ORIENTED_BOUNDARY) ||
3796              (oslrtrg != ON_ORIENTED_BOUNDARY)   )         ) {
3797           compare_r = SMALLER;
3798         }
3799       }
3800     }
3801 
3802     CGAL_SDG_DEBUG(std::cout << "debug linf_refine compare p q r = "
3803       << compare_p << " " << compare_q << " " << compare_r << std::endl;);
3804 
3805     if ((compare_p == SMALLER) ||
3806         (compare_q == SMALLER) ||
3807         (compare_r == SMALLER)   ) {
3808       return SMALLER;
3809     }
3810     if (corner_agree_pt_x && corner_agree_pt_y) {
3811       CGAL_assertion(! is_l_h_or_v);
3812       const unsigned int count_larger =
3813         (compare_p ? 1 : 0) +
3814         (compare_q ? 1 : 0) +
3815         (compare_r ? 1 : 0) ;
3816       if (count_larger >= 2) {
3817         // two points (among p, q, r) hide the line l from the vertex vv
3818         return LARGER;
3819       }
3820     }
3821     return EQUAL;
3822 
3823   }
3824 
3825   inline
3826   Comparison_result
linf_refinement(Homogeneous_point_2 & lrefhp)3827   linf_refinement( Homogeneous_point_2& lrefhp ) const
3828   {
3829     CGAL_assertion( is_v_computed );
3830     Point_2 vv ( ux_, uy_, uz_ );
3831 
3832     Comparison_result compare_p(EQUAL);
3833     Comparison_result compare_q(EQUAL);
3834     Comparison_result compare_r(EQUAL);
3835 
3836     FT difxvl = vv.x() - lrefhp.x();
3837     FT difyvl = vv.y() - lrefhp.y();
3838     FT absdifxvl = CGAL::abs(difxvl);
3839     FT absdifyvl = CGAL::abs(difyvl);
3840     Comparison_result cmplabsxy = CGAL::compare(absdifxvl, absdifyvl);
3841 
3842     // philaris: (cmplabsxy == EQUAL) means that lref is
3843     // one of the corners of the square with center vv
3844 
3845     if (p_.is_point()) {
3846       Point_2 pp = p_.point();
3847       FT difxvp = vv.x() - pp.x();
3848       FT difyvp = vv.y() - pp.y();
3849       FT absdifxvp = CGAL::abs(difxvp);
3850       FT absdifyvp = CGAL::abs(difyvp);
3851       Comparison_result cmppabsxy = CGAL::compare(absdifxvp, absdifyvp);
3852       if (! ( (cmplabsxy == SMALLER) && (cmppabsxy == SMALLER) ))
3853       {
3854         if (CGAL::compare(difxvl, difxvp) == EQUAL) {
3855           compare_p = CGAL::compare(absdifyvl, absdifyvp);
3856         }
3857       }
3858       if (! ( (cmplabsxy == LARGER ) && (cmppabsxy == LARGER ) ))
3859       {
3860         if (CGAL::compare(difyvl, difyvp) == EQUAL) {
3861           CGAL_assertion(compare_p == EQUAL);
3862           compare_p = CGAL::compare(absdifxvl, absdifxvp);
3863         }
3864       }
3865     }
3866 
3867     if (q_.is_point()) {
3868       Point_2 qq = q_.point();
3869       FT difxvq = vv.x() - qq.x();
3870       FT difyvq = vv.y() - qq.y();
3871       FT absdifxvq = CGAL::abs(difxvq);
3872       FT absdifyvq = CGAL::abs(difyvq);
3873       Comparison_result cmpqabsxy = CGAL::compare(absdifxvq, absdifyvq);
3874       if (! ( (cmplabsxy == SMALLER) && (cmpqabsxy == SMALLER) ))
3875       {
3876         if (CGAL::compare(difxvl, difxvq) == EQUAL) {
3877           compare_q = CGAL::compare(absdifyvl, absdifyvq);
3878         }
3879       }
3880       if (! ( (cmplabsxy == LARGER ) && (cmpqabsxy == LARGER ) ))
3881       {
3882         if (CGAL::compare(difyvl, difyvq) == EQUAL) {
3883           CGAL_assertion(compare_q == EQUAL);
3884           compare_q = CGAL::compare(absdifxvl, absdifxvq);
3885         }
3886       }
3887     }
3888 
3889     if (r_.is_point()) {
3890       Point_2 rr = r_.point();
3891       FT difxvr = vv.x() - rr.x();
3892       FT difyvr = vv.y() - rr.y();
3893       FT absdifxvr = CGAL::abs(difxvr);
3894       FT absdifyvr = CGAL::abs(difyvr);
3895       Comparison_result cmprabsxy = CGAL::compare(absdifxvr, absdifyvr);
3896       if (! ( (cmplabsxy == SMALLER) && (cmprabsxy == SMALLER) ))
3897       {
3898         if (CGAL::compare(difxvl, difxvr) == EQUAL) {
3899           compare_r = CGAL::compare(absdifyvl, absdifyvr);
3900         }
3901       }
3902       if (! ( (cmplabsxy == LARGER ) && (cmprabsxy == LARGER ) ))
3903       {
3904         if (CGAL::compare(difyvl, difyvr) == EQUAL) {
3905           CGAL_assertion(compare_r == EQUAL);
3906           compare_r = CGAL::compare(absdifxvl, absdifxvr);
3907         }
3908       }
3909     }
3910 
3911     CGAL_SDG_DEBUG(std::cout << "debug compare p q r = "
3912       << compare_p << " " << compare_q << " " << compare_r << std::endl;);
3913 
3914     if ((compare_p == SMALLER) ||
3915         (compare_q == SMALLER) ||
3916         (compare_r == SMALLER)   ) {
3917       return SMALLER;
3918     }
3919     /*
3920     if ((compare_p == LARGER) ||
3921         (compare_q == LARGER) ||
3922         (compare_r == LARGER)   ) {
3923       // tocheck
3924       return LARGER;
3925     }
3926     */
3927     return EQUAL;
3928 
3929   }
3930 
3931 
3932 
3933   //--------------------------------------------------------------------------
3934 
3935 private:
3936   const Site_2& p_, q_, r_;
3937 
3938   vertex_t v_type;
3939 
3940   // index that indicates the refence point for the case PPS
3941   short pps_idx;
3942 
3943   // philaris: different types are not needed any more
3944   // the case ppp
3945   //RT ux_ppp, uy_ppp, uz_ppp;
3946 
3947   // the case pps
3948   //Sqrt_1 ux_pps, uy_pps, uz_pps;
3949 
3950   // the case pss and sss
3951   //Sqrt_3 ux, uy, uz;
3952 
3953   mutable bool is_v_computed;
3954 
3955   mutable RT ux_, uy_, uz_;
3956 };
3957 
3958 
3959 } //namespace SegmentDelaunayGraphLinf_2
3960 
3961 } //namespace CGAL
3962 
3963 
3964 #endif // CGAL_SEGMENT_DELAUNAY_GRAPH_LINF_2_VORONOI_VERTEX_RING_C2_H
3965