1 // Copyright (c) 2016-2018 INRIA Sophia Antipolis, INRIA Nancy (France). 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/Periodic_4_hyperbolic_triangulation_2/include/CGAL/Periodic_4_hyperbolic_Delaunay_triangulation_traits_2.h $ 7 // $Id: Periodic_4_hyperbolic_Delaunay_triangulation_traits_2.h 5c8df66 2020-09-25T14:25:14+02:00 Jane Tournois 8 // SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial 9 // 10 // Author(s) : Iordan Iordanov 11 // Monique Teillaud 12 // 13 14 #ifndef CGAL_PERIODIC_4_HYPERBOLIC_DELAUNAY_TRIANGULATION_TRAITS_2_H 15 #define CGAL_PERIODIC_4_HYPERBOLIC_DELAUNAY_TRIANGULATION_TRAITS_2_H 16 17 #include <CGAL/license/Periodic_4_hyperbolic_triangulation_2.h> 18 19 #include <CGAL/Exact_predicates_exact_constructions_kernel_with_sqrt.h> 20 #include <CGAL/Hyperbolic_octagon_translation.h> 21 #include <CGAL/Hyperbolic_Delaunay_triangulation_traits_2.h> 22 23 #include <CGAL/Bbox_2.h> 24 #include <CGAL/Cartesian.h> 25 26 #include <boost/tuple/tuple.hpp> 27 #include <boost/variant.hpp> 28 29 #include <utility> 30 31 namespace CGAL { 32 33 namespace internal { 34 35 36 template <typename K, typename Predicate_> 37 class Hyperbolic_traits_with_translations_2_adaptor 38 : Predicate_ 39 { 40 typedef K Kernel; 41 typedef Predicate_ Predicate; 42 43 typedef typename Kernel::FT FT; 44 typedef typename Kernel::Point_2 Point; 45 typedef typename Kernel::Hyperbolic_translation Hyperbolic_translation; 46 47 // Use the construct_point_2 predicate from the kernel to convert the periodic points to Euclidean points 48 typedef typename Kernel::Construct_hyperbolic_point_2 Construct_hyperbolic_point_2; 49 50 public: 51 typedef typename Predicate::result_type result_type; 52 53 using Predicate::operator(); 54 Predicate_(pred)55 Hyperbolic_traits_with_translations_2_adaptor(const Predicate_ pred = Predicate_()) : Predicate_(pred) {} 56 operator()57 result_type operator()(const Point& p0, const Point& p1, 58 const Hyperbolic_translation& o0, const Hyperbolic_translation& o1) const 59 { 60 return operator()(pp(p0, o0), pp(p1, o1)); 61 } 62 operator()63 result_type operator()(const Point& p0, const Point& p1, const Point& p2, 64 const Hyperbolic_translation& o0, const Hyperbolic_translation& o1, const Hyperbolic_translation& o2) const 65 { 66 return operator()(pp(p0, o0), pp(p1, o1), pp(p2, o2)); 67 } 68 operator()69 result_type operator()(const Point& p0, const Point& p1, 70 const Point& p2, const Point& p3, 71 const Hyperbolic_translation& o0, const Hyperbolic_translation& o1, 72 const Hyperbolic_translation& o2, const Hyperbolic_translation& o3) const 73 { 74 return operator()(pp(p0, o0), pp(p1, o1), pp(p2, o2), pp(p3, o3)); 75 } 76 77 // Added for Side_of_hyperbolic_triangle_2 78 template <typename T> operator()79 result_type operator()(const Point& p0, const Point& p1, 80 const Point& p2, const Point& p3, T& t) const 81 { 82 return operator()(p0, p1, p2, p3, t); 83 } 84 85 template <typename T> operator()86 result_type operator()(const Point& p0, const Point& p1, 87 const Point& p2, const Point& p3, 88 const Hyperbolic_translation& o0, const Hyperbolic_translation& o1, 89 const Hyperbolic_translation& o2, const Hyperbolic_translation& o3, T& t) const 90 { 91 return operator()(pp(p0, o0), pp(p1, o1), pp(p2, o2), pp(p3, o3), t); 92 } 93 private: pp(const Point & p,const Hyperbolic_translation & o)94 Point pp(const Point &p, const Hyperbolic_translation &o) const 95 { 96 return Construct_hyperbolic_point_2()(p, o); 97 } 98 }; 99 100 template < typename K, typename Construct_point_base> 101 class Periodic_4_construct_hyperbolic_point_2 102 : public Construct_point_base 103 { 104 105 private: 106 typedef K Kernel; 107 typedef typename Kernel::FT NT; 108 typedef typename Kernel::Point_2 Point; 109 typedef typename Kernel::Hyperbolic_translation Hyperbolic_translation; 110 111 public: 112 typedef Point result_type; 113 114 using Construct_point_base::operator(); 115 Periodic_4_construct_hyperbolic_point_2()116 Periodic_4_construct_hyperbolic_point_2() { } 117 operator()118 Point operator()(const Point& pt, const Hyperbolic_translation& tr) const 119 { 120 if(tr.is_identity()) 121 return operator()(pt); 122 123 Point p = operator()(pt); 124 std::pair<NT, NT> a, b; 125 a = tr.alpha(); 126 b = tr.beta(); 127 128 // Prepare variables 129 NT ax(a.first); 130 NT bx(b.first); 131 NT zx(p.x()); 132 NT ay(a.second); 133 NT by(b.second); 134 NT zy(p.y()); 135 136 // Compute parts of fraction 137 NT rn(ax*zx - ay*zy + bx); // real part of numerator 138 NT in(ay*zx + ax*zy + by); // imaginary part of numerator 139 NT rd(bx*zx + by*zy + ax); // real part of denominator 140 NT id(bx*zy - by*zx - ay); // imaginary part of denominator 141 142 // The denominator cannot be zero 143 NT den(rd*rd + id*id); 144 CGAL_assertion(den != NT(0)); 145 146 // Compute real and imaginary part of result 147 NT rx((rn*rd + in*id)); 148 NT ry((in*rd - rn*id)); 149 150 return operator()(rx/den, ry/den); 151 } 152 operator()153 Point operator() (const Point& p) const 154 { 155 return p; 156 } 157 158 }; 159 160 161 template <typename Traits> 162 class Compute_approximate_hyperbolic_diameter 163 { 164 165 typedef typename Traits:: Euclidean_line_2 Euclidean_line_2; 166 typedef typename Traits::Circle_2 Circle_2; 167 typedef typename Traits::Hyperbolic_point_2 Hyperbolic_point_2; 168 public: 169 typedef double result_type; 170 _gt(gt)171 Compute_approximate_hyperbolic_diameter(const Traits& gt = Traits()) : _gt(gt) {} 172 operator()173 result_type operator()(const Hyperbolic_point_2& p1, 174 const Hyperbolic_point_2& p2, 175 const Hyperbolic_point_2& p3) 176 { 177 178 Circle_2 c = _gt.construct_circle_2_object()(p1, p2, p3); 179 180 Hyperbolic_point_2 p0 = _gt.construct_point_2_object()(0, 0); 181 Circle_2 c0 = _gt.construct_circle_2_object()(p0, 1); 182 Euclidean_line_2 ell = _gt.construct_line_2_object()(p0, c.center()); 183 184 if(ell.is_degenerate()) 185 return 5.; // for the Bolza surface, this just needs to be larger than half the systole (which is ~1.5) 186 187 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> res1 = _gt.construct_inexact_intersection_2_object()(c0, ell); 188 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> res2 = _gt.construct_inexact_intersection_2_object()(c , ell); 189 190 Hyperbolic_point_2 a = res1.first; 191 Hyperbolic_point_2 b = res1.second; 192 193 Hyperbolic_point_2 p = res2.first; 194 Hyperbolic_point_2 q = res2.second; 195 196 double aq = CGAL::sqrt(CGAL::to_double(CGAL::squared_distance(a, q))); 197 double pb = CGAL::sqrt(CGAL::to_double(CGAL::squared_distance(p, b))); 198 double ap = CGAL::sqrt(CGAL::to_double(CGAL::squared_distance(a, p))); 199 double qb = CGAL::sqrt(CGAL::to_double(CGAL::squared_distance(q, b))); 200 201 double hyperdist = std::fabs(std::log(CGAL::to_double((aq*pb)/(ap*qb)))); 202 203 return hyperdist; 204 } 205 206 private: 207 const Traits& _gt; 208 }; 209 210 211 template <typename Traits> 212 class Construct_hyperbolic_line_2 213 { 214 typedef typename Traits::Weighted_point_2 Weighted_point_2; 215 typedef typename Traits::Hyperbolic_point_2 Hyperbolic_point_2; 216 typedef typename Traits::Hyperbolic_segment_2 Hyperbolic_segment_2; 217 typedef typename Traits::Point_2 Bare_point; 218 typedef typename Traits::Circle_2 Circle_2; 219 typedef typename Traits::Circular_arc_2 Circular_arc_2; 220 typedef typename Traits::FT FT; 221 222 public: _gt(gt)223 Construct_hyperbolic_line_2(const Traits gt = Traits()) : _gt(gt) { } 224 operator()225 Hyperbolic_segment_2 operator()(const Hyperbolic_point_2& p, const Hyperbolic_point_2& q) 226 { 227 Origin o; 228 if(_gt.collinear_2_object()(p, q, _gt.construct_hyperbolic_point_2_object(o))) 229 { 230 std::pair<Hyperbolic_point_2,Hyperbolic_point_2> inters = _gt.construct_intersection_2_object()(_gt.construct_line_2_object()(p,q), _gt.construct_circle_2_object()(_gt.construct_hyperbolic_point_2_object()(0,0),1)); 231 return Euclidean_segment_2(inters.first, inters.second); 232 } 233 234 Weighted_point_2 wp(p); 235 Weighted_point_2 wq(q); 236 Weighted_point_2 wo(_gt.construct_hyperbolic_point_2_object()(o), FT(1)); // Poincaré circle 237 238 Bare_point center = _gt.construct_weighted_circumcenter_2_object()(wp, wo, wq); 239 FT sq_radius = _gt.compute_squared_Euclidean_distance_2_object()(p, center); 240 241 Circle_2 circle = _gt.construct_circle_2_object()(center, sq_radius); 242 std::pair<Hyperbolic_point_2,Hyperbolic_point_2> inters = _gt.construct_intersection_2_object()(circle, _gt.construct_circle_2_object()(_gt.construct_hyperbolic_point_2_object()(0,0),1)); 243 244 if(_gt.orientation_2_object()(circle.center(), inters.first, inters.second) == POSITIVE) { 245 return Circular_arc_2(circle, inters.first, inters.second); 246 } else { 247 return Circular_arc_2(circle, inters.second, inters.first); 248 } 249 } 250 251 private: 252 const Traits& _gt; 253 }; 254 255 256 257 template <typename Traits> 258 class Construct_inexact_intersection_2 259 { 260 typedef typename Traits::FT FT; 261 typedef typename Traits::Circle_2 Circle_2; 262 typedef typename Traits::Euclidean_line_2 Euclidean_line_2; 263 typedef typename Traits::Hyperbolic_point_2 Hyperbolic_point_2; 264 typedef typename Traits::Hyperbolic_segment_2 Hyperbolic_segment_2; 265 typedef typename Traits::Euclidean_segment_2 Euclidean_segment_2; 266 typedef typename Traits::Circular_arc_2 Circular_arc_2; 267 268 public: _gt(gt)269 Construct_inexact_intersection_2(const Traits& gt = Traits()) : _gt(gt) {} 270 operator()271 Hyperbolic_point_2 operator()(const Euclidean_line_2& ell1, const Euclidean_line_2& ell2) 272 { 273 if(std::fabs(CGAL::to_double(ell1.b())) < 1e-16) 274 std::swap(ell1, ell2); 275 276 double a1 = CGAL::to_double(ell1.a()), b1 = CGAL::to_double(ell1.b()), c1 = CGAL::to_double(ell1.c()); 277 double a2 = CGAL::to_double(ell2.a()), b2 = CGAL::to_double(ell2.b()), c2 = CGAL::to_double(ell2.c()); 278 279 CGAL_assertion(std::fabs(b1) > 1e-16); 280 if(std::fabs(b2) > 1e-16) 281 CGAL_assertion(std::fabs(a1/b1 - a2/b2) > 1e-16); 282 283 double lambda1 = -a1/b1; 284 double mu1 = -c1/b1; 285 double x = (-c2 - mu1*b2)/(a2 + lambda1*b2); 286 double y = lambda1*x + mu1; 287 return Hyperbolic_point_2(x, y); 288 } 289 operator()290 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> operator()(const Euclidean_line_2& ell, 291 const Circle_2& cc) 292 { 293 double a = CGAL::to_double(ell.a()), 294 b = CGAL::to_double(ell.b()), 295 c = CGAL::to_double(ell.c()); 296 double p = CGAL::to_double(cc.center().x()), 297 q = CGAL::to_double(cc.center().y()), 298 r2 = CGAL::to_double(cc.squared_radius()); 299 300 double A, B, C, D; 301 double x1, y1, x2, y2; 302 if(std::fabs(a) < 1e-16) 303 { 304 y1 = -c/b; y2 = -c/b; 305 A = b*p; 306 D = -b*b*q*q + b*b*r2 - 2.*b*c*q - c*c; 307 x1 = (A + CGAL::sqrt(D))/b; 308 x2 = (A - CGAL::sqrt(D))/b; 309 } 310 else if(std::fabs(b) < 1e-16) 311 { 312 x1 = -c/a; x2 = -c/a; 313 A = q*a; 314 D = -a*a*p*p + r2*a*a - 2.*a*c*p - c*c; 315 y1 = (A + CGAL::sqrt(D))/a; 316 y2 = (A - CGAL::sqrt(D))/a; 317 } 318 else 319 { 320 A = a*a*q - a*b*p-b*c; 321 C = (-b*q - c)*a*a + b*b*p*a; 322 D = -a*a*(b*b*q*q + 2.*q*(p*a + c)*b - b*b*r2 + (p*p - r2)*a*a + 2.*a*c*p + c*c); 323 B = a*a + b*b; 324 325 y1 = (A + CGAL::sqrt(D))/B; 326 y2 = (A - CGAL::sqrt(D))/B; 327 x1 = (C - b*CGAL::sqrt(D))/(a*(a*a + b*b)); 328 x2 = (C + b*CGAL::sqrt(D))/(a*(a*a + b*b)); 329 } 330 331 Hyperbolic_point_2 p1(x1, y1); 332 Hyperbolic_point_2 p2(x2, y2); 333 334 return std::make_pair(p1, p2); 335 } 336 operator()337 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> operator()(const Circle_2& c, 338 const Euclidean_line_2& ell) 339 { 340 return operator()(ell, c); 341 } 342 operator()343 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> operator()(const Circle_2& c1, 344 const Circle_2& c2) 345 { 346 double xa = CGAL::to_double(c1.center().x()), 347 ya = CGAL::to_double(c1.center().y()); 348 double xb = CGAL::to_double(c2.center().x()), 349 yb = CGAL::to_double(c2.center().y()); 350 double d2 = (xa-xb)*(xa-xb) + (ya-yb)*(ya-yb); 351 352 double ra = CGAL::sqrt(CGAL::to_double(c1.squared_radius())); 353 double rb = CGAL::sqrt(CGAL::to_double(c2.squared_radius())); 354 double K = CGAL::sqrt(((ra+rb)*(ra+rb)-d2)*(d2-(ra-rb)*(ra-rb)))/4.; 355 356 double xbase = (xb + xa)/2. + (xb - xa)*(ra*ra - rb*rb)/d2/2.; 357 double xdiff = 2.*(yb - ya)*K/d2; 358 double x1 = xbase + xdiff; 359 double x2 = xbase - xdiff; 360 361 double ybase = (yb + ya)/2. + (yb - ya)*(ra*ra - rb*rb)/d2/2.; 362 double ydiff = -2.*(xb - xa)*K/d2; 363 double y1 = ybase + ydiff; 364 double y2 = ybase - ydiff; 365 366 Hyperbolic_point_2 res1(x1, y1); 367 Hyperbolic_point_2 res2(x2, y2); 368 return std::make_pair(res1, res2); 369 } 370 operator()371 Hyperbolic_point_2 operator()(const Hyperbolic_segment_2& s1, 372 const Hyperbolic_segment_2& s2) 373 { 374 if(Circular_arc_2* c1 = boost::get<Circular_arc_2>(&s1)) 375 { 376 if(Circular_arc_2* c2 = boost::get<Circular_arc_2>(&s2)) 377 { 378 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> res = operator()(c1->circle(), 379 c2->circle()); 380 Hyperbolic_point_2 p1 = res.first; 381 if(p1.x()*p1.x() + p1.y()*p1.y() < FT(1)) 382 return p1; 383 384 Hyperbolic_point_2 p2 = res.second; 385 CGAL_assertion(p2.x()*p2.x() + p2.y()*p2.y() < FT(1)); 386 return p2; 387 } 388 else 389 { 390 Euclidean_segment_2* ell2 = boost::get<Euclidean_segment_2>(&s2); 391 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> res = operator()(c1->circle(), 392 ell2->supporting_line()); 393 Hyperbolic_point_2 p1 = res.first; 394 if(p1.x()*p1.x() + p1.y()*p1.y() < FT(1)) 395 return p1; 396 397 Hyperbolic_point_2 p2 = res.second; 398 CGAL_assertion(p2.x()*p2.x() + p2.y()*p2.y() < FT(1)); 399 return p2; 400 } 401 } 402 else 403 { 404 Euclidean_segment_2* ell1 = boost::get<Euclidean_segment_2>(&s1); 405 if(Circular_arc_2* c2 = boost::get<Circular_arc_2>(&s2)) 406 { 407 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> res = operator()(ell1->supporting_line(), 408 c2->circle()); 409 Hyperbolic_point_2 p1 = res.first; 410 if(p1.x()*p1.x() + p1.y()*p1.y() < FT(1)) 411 return p1; 412 413 Hyperbolic_point_2 p2 = res.second; 414 CGAL_assertion(p2.x()*p2.x() + p2.y()*p2.y()) < FT(1); 415 return p2; 416 } 417 else 418 { 419 Euclidean_segment_2* ell2 = boost::get<Euclidean_segment_2>(&s2); 420 Hyperbolic_point_2 p1 = operator()(ell1->supporting_line(), ell2->supporting_line()); 421 CGAL_assertion(p1.x()*p1.x() + p1.y()*p1.y()) < FT(1); 422 return p1; 423 } 424 } 425 } 426 427 private: 428 const Traits& _gt; 429 }; 430 431 432 433 template <typename Traits> 434 class Construct_inexact_hyperbolic_circumcenter_2 435 { 436 typedef typename Traits::FT FT; 437 typedef typename Traits::Hyperbolic_Voronoi_point_2 Hyperbolic_Voronoi_point_2; 438 typedef typename Traits::Euclidean_circle_or_line_2 Euclidean_circle_or_line_2; 439 typedef typename Traits::Hyperbolic_point_2 Hyperbolic_point_2; 440 typedef typename Traits::Circle_2 Circle_2; 441 typedef typename Traits::Euclidean_line_2 Euclidean_line_2; 442 443 public: _gt(gt)444 Construct_inexact_hyperbolic_circumcenter_2(const Traits& gt = Traits()) : _gt(gt) {} 445 446 typedef Hyperbolic_Voronoi_point_2 result_type; 447 operator()448 Hyperbolic_Voronoi_point_2 operator()(const Hyperbolic_point_2& p, 449 const Hyperbolic_point_2& q, 450 const Hyperbolic_point_2& r) 451 { 452 Origin o; 453 Hyperbolic_point_2 po = _gt.construct_hyperbolic_point_2_object()(o); 454 Circle_2 l_inf = _gt.construct_circle_2_object()(po, FT(1)); 455 456 // Check if |p,O| = |q,O| = |r,O| -- then the circumcenter is the origin O 457 if(_gt.compare_distance_2_object()(po,p,q) == EQUAL && _gt.compare_distance_2_object()(po,p,r) == EQUAL) 458 return po; 459 460 Euclidean_circle_or_line_2 bis_pq = _gt.construct_circle_or_line_supporting_bisector_2_object()(p,q); 461 Euclidean_circle_or_line_2 bis_qr = _gt.construct_circle_or_line_supporting_bisector_2_object()(q,r); 462 463 // now supporting objects cannot both be Euclidean lines 464 465 Euclidean_line_2* l; 466 Circle_2* c; 467 468 if(Circle_2* c_pq = boost::get<Circle_2>(&bis_pq)) 469 { 470 if(Circle_2* c_qr = boost::get<Circle_2>(&bis_qr)) 471 { 472 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> inters = _gt.construct_inexact_intersection_2_object()(*c_pq, *c_qr); 473 474 if(_gt.has_on_bounded_side_2_object()(l_inf, inters.first)) 475 return inters.first; 476 477 return inters.second; 478 } 479 // here bis_qr is a line 480 l = boost::get<Euclidean_line_2>(&bis_qr); 481 c = c_pq; 482 } 483 else 484 { 485 // here bis_pq is a line 486 l = boost::get<Euclidean_line_2>(&bis_pq); 487 c = boost::get<Circle_2>(&bis_qr); 488 } 489 490 std::pair<Hyperbolic_point_2, Hyperbolic_point_2> inters = _gt.construct_inexact_intersection_2_object()(*c, *l); 491 492 if(_gt.has_on_bounded_side_2_object()(l_inf, inters.first)) 493 return inters.first; 494 495 return inters.second; 496 } 497 498 template <typename Face_handle> operator()499 Hyperbolic_Voronoi_point_2 operator()(const Face_handle fh) 500 { 501 return operator()(_gt.construct_hyperbolic_point_2_object()(fh.vertex(0)->point(), fh.translation(0)), 502 _gt.construct_hyperbolic_point_2_object()(fh.vertex(1)->point(), fh.translation(1)), 503 _gt.construct_hyperbolic_point_2_object()(fh.vertex(2)->point(), fh.translation(2)) ); 504 } 505 506 private: 507 const Traits& _gt; 508 }; // end Construct_inexact_hyperbolic_circumcenter_2 509 510 511 512 template <typename Traits> 513 class Side_of_original_octagon 514 { 515 typedef typename Traits::FT FT; 516 typedef typename Traits::Circle_2 Circle_2; 517 typedef typename Traits::Hyperbolic_point_2 Hyperbolic_point_2; 518 519 public: _gt(gt)520 Side_of_original_octagon(const Traits& gt = Traits()) : _gt(gt) {} 521 522 template <class Point_2_template> operator()523 CGAL::Bounded_side operator()(const Point_2_template& p) 524 { 525 FT n2(2); 526 FT n4(4); 527 FT sq2(CGAL::sqrt(n2)); 528 FT p14(CGAL::sqrt(CGAL::sqrt(n2*n2*n2))); // 2^{3/4} 529 FT p34(CGAL::sqrt(CGAL::sqrt(n2*n2*n2))); // 2^{3/4} 530 FT s1 (CGAL::sqrt(n2 + CGAL::sqrt(n2))); // CGAL::sqrt(2 + CGAL::sqrt(2)) 531 FT s2 (CGAL::sqrt(n2 - CGAL::sqrt(n2))); // CGAL::sqrt(2 - CGAL::sqrt(2)) 532 FT s3 (CGAL::sqrt(sq2 - FT(1))); // CGAL::sqrt(CGAL::sqrt(2) - 1) 533 _gt.construct_point_2_object()(s1*p34/n4, -s2*p34/n4); 534 Hyperbolic_point_2 V0 = _gt.construct_hyperbolic_point_2_object()(_gt.construct_point_2_object()(s1*p34/n4, -s2*p34/n4)); 535 Hyperbolic_point_2 V1 = _gt.construct_hyperbolic_point_2_object()(_gt.construct_point_2_object()(p14*(s1+s2)/n4, p14*(s1-s2)/n4)); 536 Hyperbolic_point_2 V2 = _gt.construct_hyperbolic_point_2_object()(_gt.construct_point_2_object()(s2*p34/n4, s1*p34/n4)); 537 538 Hyperbolic_point_2 M0 = _gt.construct_hyperbolic_point_2_object()(_gt.construct_point_2_object()(s3, FT(0))); 539 Hyperbolic_point_2 M1 = _gt.construct_hyperbolic_point_2_object()(_gt.construct_point_2_object()(sq2*s3/n2, sq2*s3/n2)); 540 541 // Poincare disk (i.e., unit Euclidean disk) 542 Circle_2 Poincare = _gt.construct_circle_2_object()(_gt.construct_hyperbolic_point_2_object()(_gt.construct_point_2_object()(0, 0)), 1*1); 543 544 // Euclidean circle supporting the side [V0, V1] 545 Circle_2 SuppCircle1 = _gt.construct_circle_2_object()(M0, V0, V1); 546 547 // Euclidean circle supporting the side [V1, V2] 548 Circle_2 SuppCircle2 = _gt.construct_circle_2_object()(M1, V1, V2); 549 550 // This transformation brings the point in the first quadrant (positive x, positive y) 551 FT x(FT(p.x()) > FT(0) ? p.x() : -p.x()); 552 FT y(FT(p.y()) > FT(0) ? p.y() : -p.y()); 553 554 // This brings the point in the first octant (positive x and y < x) 555 if(y > x) 556 { 557 FT tmp = x; 558 x = y; 559 y = tmp; 560 } 561 562 // This tells us whether the point is on the side of the open boundary 563 bool on_open_side = ((p.y() + std::tan(CGAL_PI / 8.) * p.x()) < 0.0); 564 565 Hyperbolic_point_2 t = _gt.construct_hyperbolic_point_2_object()(x, y); 566 567 CGAL::Bounded_side PoincareSide = _gt.bounded_side_2_object()(Poincare, t); 568 CGAL::Bounded_side Circ1Side = _gt.bounded_side_2_object()(SuppCircle1, t); 569 CGAL::Bounded_side Circ2Side = _gt.bounded_side_2_object()(SuppCircle2, t); 570 571 // First off, the point needs to be inside the Poincare disk. if not, there's no hope. 572 if(PoincareSide == CGAL::ON_BOUNDED_SIDE) 573 { 574 // Inside the Poincare disk, but still outside the original domain 575 if(Circ1Side == CGAL::ON_BOUNDED_SIDE || Circ2Side == CGAL::ON_BOUNDED_SIDE) 576 return CGAL::ON_UNBOUNDED_SIDE; 577 578 // Inside the Poincare disk and inside the original domain 579 if(Circ1Side == CGAL::ON_UNBOUNDED_SIDE && Circ2Side == CGAL::ON_UNBOUNDED_SIDE) 580 return CGAL::ON_BOUNDED_SIDE; 581 582 // This is boundary, but we only consider the upper half. The lower half means outside. 583 if(on_open_side) 584 return CGAL::ON_UNBOUNDED_SIDE; 585 else 586 return CGAL::ON_BOUNDED_SIDE; 587 } else 588 { 589 return CGAL::ON_UNBOUNDED_SIDE; 590 } 591 } 592 593 private: 594 const Traits& _gt; 595 }; 596 597 598 } // end namespace internal 599 600 601 template <typename Kernel = Exact_predicates_exact_constructions_kernel_with_sqrt, 602 template<typename> class Translation_type = Hyperbolic_octagon_translation> 603 class Periodic_4_hyperbolic_Delaunay_triangulation_traits_2 604 : public Hyperbolic_Delaunay_triangulation_traits_2<Kernel> 605 { 606 typedef Hyperbolic_Delaunay_triangulation_traits_2<Kernel> Base; 607 typedef Periodic_4_hyperbolic_Delaunay_triangulation_traits_2<Kernel, Translation_type> Self; 608 609 public: 610 typedef typename Base::FT FT; 611 typedef Translation_type<FT> Hyperbolic_translation; 612 613 typedef typename Base::Hyperbolic_point_2 Hyperbolic_point_2; 614 typedef Hyperbolic_point_2 Hyperbolic_Voronoi_point_2; 615 typedef typename Base::Line_2 Euclidean_line_2; 616 typedef typename Base::Euclidean_circle_or_line_2 Euclidean_circle_or_line_2; 617 typedef typename Base::Circular_arc_2 Circular_arc_2; 618 typedef typename Base::Euclidean_segment_2 Euclidean_segment_2; //only used internally here 619 typedef typename Base::Hyperbolic_segment_2 Hyperbolic_segment_2; 620 typedef Euclidean_segment_2 Line_segment_2; // kept for demo 621 typedef Hyperbolic_segment_2 Segment_2; // kept for demo 622 623 // the following types are only used internally in this traits class, 624 // so they need not be documented, and they don't need _object() 625 typedef typename Base::Construct_Euclidean_bisector_2 Construct_Euclidean_bisector_2; 626 typedef typename Base::Construct_intersection_2 Construct_intersection_2; 627 typedef typename Base::Construct_hyperbolic_bisector_2 Construct_hyperbolic_bisector_2; 628 typedef typename Base::Construct_circle_or_line_supporting_bisector Construct_circle_or_line_supporting_bisector; 629 typedef typename Base::Euclidean_collinear_2 Euclidean_collinear_2; 630 typedef typename Base::Compute_squared_Euclidean_distance_2 Compute_squared_Euclidean_distance_2; 631 typedef typename Base::Has_on_bounded_side_2 Has_on_bounded_side_2; 632 633 // Wrappers for the translation adapter 634 typedef internal::Hyperbolic_traits_with_translations_2_adaptor< 635 Self, typename Base::Orientation_2> Orientation_2; 636 637 typedef internal::Hyperbolic_traits_with_translations_2_adaptor< 638 Self, typename Base::Side_of_oriented_circle_2> Side_of_oriented_circle_2; 639 640 typedef internal::Hyperbolic_traits_with_translations_2_adaptor< 641 Self, typename Base::Construct_hyperbolic_circumcenter_2> Construct_hyperbolic_circumcenter_2; 642 643 typedef internal::Hyperbolic_traits_with_translations_2_adaptor< 644 Self, typename Base::Construct_hyperbolic_segment_2> Construct_hyperbolic_segment_2; 645 646 typedef internal::Hyperbolic_traits_with_translations_2_adaptor< 647 Self, typename Base::Construct_segment_2> Construct_segment_2; 648 649 typedef internal::Hyperbolic_traits_with_translations_2_adaptor< 650 Self, typename Base::Construct_triangle_2> Construct_triangle_2; 651 652 typedef internal::Hyperbolic_traits_with_translations_2_adaptor< 653 Self, typename Base::Compare_distance_2> Compare_distance_2; 654 655 typedef internal::Periodic_4_construct_hyperbolic_point_2< 656 Self, typename Kernel::Construct_point_2> Construct_hyperbolic_point_2; 657 658 typedef internal::Hyperbolic_traits_with_translations_2_adaptor< 659 Self, typename Base::Side_of_oriented_hyperbolic_segment_2> Side_of_oriented_hyperbolic_segment_2; 660 661 662 typedef internal::Compute_approximate_hyperbolic_diameter<Self> Compute_approximate_hyperbolic_diameter; 663 typedef internal::Construct_hyperbolic_line_2<Self> Construct_hyperbolic_line_2; 664 typedef internal::Construct_inexact_intersection_2<Self> Construct_inexact_intersection_2; 665 typedef internal::Construct_inexact_hyperbolic_circumcenter_2<Self> Construct_inexact_hyperbolic_circumcenter_2; 666 typedef internal::Side_of_original_octagon<Self> Side_of_original_octagon; 667 668 669 public: 670 671 Compute_approximate_hyperbolic_diameter compute_approximate_hyperbolic_diameter_object()672 compute_approximate_hyperbolic_diameter_object() const 673 { return Compute_approximate_hyperbolic_diameter(*this); } 674 675 676 Construct_hyperbolic_line_2 construct_hyperbolic_line_2_object()677 construct_hyperbolic_line_2_object() const 678 { return Construct_hyperbolic_line_2(*this); } 679 680 681 Construct_inexact_intersection_2 construct_inexact_intersection_2_object()682 construct_inexact_intersection_2_object() const 683 { return Construct_inexact_intersection_2(*this); } 684 685 686 Construct_inexact_hyperbolic_circumcenter_2 construct_inexact_hyperbolic_circumcenter_2_object()687 construct_inexact_hyperbolic_circumcenter_2_object() const 688 { return Construct_inexact_hyperbolic_circumcenter_2(*this); } 689 690 691 Side_of_original_octagon side_of_original_octagon_object()692 side_of_original_octagon_object() const 693 { return Side_of_original_octagon(*this); } 694 695 696 Side_of_oriented_hyperbolic_segment_2 side_of_oriented_hyperbolic_segment_2_object()697 side_of_oriented_hyperbolic_segment_2_object() const 698 { return Side_of_oriented_hyperbolic_segment_2(*this); } 699 700 701 Construct_triangle_2 construct_triangle_2_object()702 construct_triangle_2_object() const 703 { return Construct_triangle_2(); } 704 705 Construct_hyperbolic_point_2 construct_hyperbolic_point_2_object()706 construct_hyperbolic_point_2_object() const 707 { return Construct_hyperbolic_point_2(); } 708 709 Construct_hyperbolic_segment_2 construct_hyperbolic_segment_2_object()710 construct_hyperbolic_segment_2_object() const 711 { return Construct_hyperbolic_segment_2(); } 712 713 Construct_segment_2 construct_segment_2_object()714 construct_segment_2_object() const 715 { return Construct_segment_2(); } 716 717 Orientation_2 orientation_2_object()718 orientation_2_object() const 719 { return Orientation_2(); } 720 721 Side_of_oriented_circle_2 side_of_oriented_circle_2_object()722 side_of_oriented_circle_2_object() const 723 { return Side_of_oriented_circle_2(); } 724 725 726 public: Base(base)727 Periodic_4_hyperbolic_Delaunay_triangulation_traits_2(const Base& base = Base()) : Base(base) {} 728 729 730 }; // class Periodic_4_hyperbolic_Delaunay_triangulation_traits_2 731 732 } // namespace CGAL 733 734 #endif // CGAL_PERIODIC_4_HYPERBOLIC_DELAUNAY_TRIANGULATION_TRAITS_2_H 735