1 // Copyright (c) 2012, 2017 INRIA Sophia-Antipolis (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/Alpha_shapes_2/include/CGAL/internal/Lazy_alpha_nt_2.h $ 7 // $Id: Lazy_alpha_nt_2.h 254d60f 2019-10-19T15:23:19+02:00 Sébastien Loriot 8 // SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial 9 // 10 // Author(s) : Sébastien Loriot <sebastien.loriot@geometryfactory.com> 11 // Mael Rouxel-Labbé 12 13 #ifndef CGAL_INTERNAL_LAZY_ALPHA_NT_2_H 14 #define CGAL_INTERNAL_LAZY_ALPHA_NT_2_H 15 16 #include <CGAL/license/Alpha_shapes_2.h> 17 18 #include <CGAL/assertions.h> 19 #include <CGAL/internal/Exact_type_selector.h> 20 #include <CGAL/number_type_basic.h> 21 #include <CGAL/Cartesian_converter.h> 22 #include <CGAL/Has_conversion.h> 23 24 #include <boost/mpl/has_xxx.hpp> 25 #include <boost/type_traits.hpp> 26 27 #include <iostream> 28 29 namespace CGAL { 30 31 namespace internal { 32 33 // check whether Cartesian_converter can do the following conversions 34 // -- Input_traits::(Weighted_)Point_2 to K2::(Weighted_)Point_2 35 // -- Input_traits::(Weighted_)Point_2 to K3::(Weighted_)Point_2 36 // 37 template < class Input_traits, class Kernel_approx, class Kernel_exact, 38 class Weighted_tag > 39 class Is_traits_point_convertible_2 40 { 41 typedef typename Kernel_traits<typename Input_traits::Point_2>::Kernel Kernel_input; 42 43 typedef typename Input_traits::Point_2 K1P; 44 typedef typename Kernel_approx::Point_2 K2P; 45 typedef typename Kernel_exact::Point_2 K3P; 46 47 public: 48 static const bool value 49 = (Has_conversion<Kernel_input, Kernel_approx, K1P, K2P>::value && 50 Has_conversion<Kernel_input, Kernel_exact, K1P, K3P>::value); 51 }; 52 53 template < class Input_traits, class Kernel_approx, class Kernel_exact > 54 class Is_traits_point_convertible_2<Input_traits, Kernel_approx, Kernel_exact, 55 ::CGAL::Tag_true /* Weighted_tag */> 56 { 57 typedef typename Kernel_traits<typename Input_traits::Point_2>::Kernel Kernel_input; 58 59 typedef typename Input_traits::Weighted_point_2 K1WP; 60 typedef typename Kernel_approx::Weighted_point_2 K2WP; 61 typedef typename Kernel_exact::Weighted_point_2 K3WP; 62 63 public: 64 static const bool value 65 = (Has_conversion<Kernel_input, Kernel_approx, K1WP, K2WP>::value && 66 Has_conversion<Kernel_input, Kernel_exact, K1WP, K3WP>::value); 67 }; 68 69 template <class T> 70 struct Input_points_for_lazy_alpha_nt_2 71 { 72 int nbpts; 73 const T* p0; 74 const T* p1; 75 const T* p2; 76 }; 77 78 //non-weighted case 79 template <class Weighted_tag, class Input_traits, class Kernel_input, 80 class Kernel_approx, class Kernel_exact> 81 struct Types_for_alpha_nt_2 82 { 83 //Converter types 84 typedef CGAL::Cartesian_converter<Kernel_input, Kernel_approx> To_approx; 85 typedef CGAL::Cartesian_converter<Kernel_input, Kernel_exact> To_exact; 86 87 //Point types 88 typedef typename Kernel_approx::Point_2 Approx_point; 89 typedef typename Kernel_exact::Point_2 Exact_point; 90 typedef typename Input_traits::Point_2 Input_point; 91 92 //Constructions 93 typedef typename Kernel_approx::Compute_squared_radius_2 Approx_squared_radius; 94 typedef typename Kernel_exact::Compute_squared_radius_2 Exact_squared_radius; 95 }; 96 97 //weighted case 98 template <class Input_traits, class Kernel_input, 99 class Kernel_approx, class Kernel_exact> 100 struct Types_for_alpha_nt_2< ::CGAL::Tag_true /* Weighted_tag */, 101 Input_traits, Kernel_input, 102 Kernel_approx, Kernel_exact > 103 { 104 //Converter types 105 typedef CGAL::Cartesian_converter<Kernel_input, Kernel_approx> To_approx; 106 typedef CGAL::Cartesian_converter<Kernel_input, Kernel_exact> To_exact; 107 108 //Point types 109 typedef typename Kernel_approx::Weighted_point_2 Approx_point; 110 typedef typename Kernel_exact::Weighted_point_2 Exact_point; 111 typedef typename Input_traits::Weighted_point_2 Input_point; 112 113 //Constructions 114 typedef typename Kernel_approx::Compute_squared_radius_smallest_orthogonal_circle_2 Approx_squared_radius; 115 typedef typename Kernel_exact::Compute_squared_radius_smallest_orthogonal_circle_2 Exact_squared_radius; 116 }; 117 118 template<class Input_traits, bool mode, class Weighted_tag> 119 class Lazy_alpha_nt_2 120 { 121 //NT & kernels 122 typedef CGAL::Interval_nt<mode> NT_approx; 123 124 //Gmpq or Quotient<MP_float> 125 typedef Exact_field_selector<double>::Type NT_exact; 126 typedef CGAL::Simple_cartesian<NT_approx> Kernel_approx; 127 typedef CGAL::Simple_cartesian<NT_exact> Kernel_exact; 128 typedef typename Kernel_traits<typename Input_traits::Point_2>::Kernel Kernel_input; 129 130 //Helper class for weighted and non-weighted case 131 typedef Types_for_alpha_nt_2<Weighted_tag, Input_traits, Kernel_input, 132 Kernel_approx, Kernel_exact> Types; 133 134 //Converters 135 typedef typename Types::To_approx To_approx; 136 typedef typename Types::To_exact To_exact; 137 138 //Constructions class 139 typedef typename Types::Approx_squared_radius Approx_squared_radius; 140 typedef typename Types::Exact_squared_radius Exact_squared_radius; 141 142 //Point 143 typedef typename Types::Approx_point Approx_point; 144 typedef typename Types::Exact_point Exact_point; 145 typedef typename Types::Input_point Input_point; 146 147 //Convertion functions 148 Approx_point to_approx(const Input_point& wp) const 149 { 150 // The traits class' Point_2 must be convertible using the Cartesian converter 151 CGAL_static_assertion((Is_traits_point_convertible_2< 152 Input_traits, Kernel_approx, Kernel_exact, Weighted_tag>::value)); 153 154 To_approx converter; 155 return converter(wp); 156 } 157 158 Exact_point to_exact(const Input_point& wp) const 159 { 160 // The traits class' Point_2 must be convertible using the Cartesian converter 161 CGAL_static_assertion((Is_traits_point_convertible_2< 162 Input_traits, Kernel_approx, Kernel_exact, Weighted_tag>::value)); 163 164 To_exact converter; 165 return converter(wp); 166 } 167 168 //members 169 //the members can be updated when calling method exact() 170 mutable boost::optional<NT_exact> exact_; 171 mutable NT_approx approx_; 172 173 //private functions 174 typedef Input_points_for_lazy_alpha_nt_2<Input_point> Data_vector; 175 Data_vector input_points; 176 177 const Data_vector& data() const{ return input_points;} 178 Data_vector& data(){ return input_points;} 179 180 static double & relative_precision_of_to_double_internal() 181 { 182 CGAL_STATIC_THREAD_LOCAL_VARIABLE(double, relative_precision_of_to_double, 0.00001); 183 return relative_precision_of_to_double; 184 } 185 186 public: 187 188 static const double & get_relative_precision_of_to_double() 189 { 190 return relative_precision_of_to_double_internal(); 191 } 192 193 static void set_relative_precision_of_to_double(double d) 194 { 195 CGAL_assertion((0 < d) & (d < 1)); 196 relative_precision_of_to_double_internal() = d; 197 } 198 199 typedef NT_exact Exact_nt; 200 typedef NT_approx Approximate_nt; 201 202 void update_exact() const 203 { 204 switch(data().nbpts) { 205 case 1: 206 exact_ = Exact_squared_radius()( to_exact(*data().p0) ); 207 break; 208 case 2: 209 exact_ = Exact_squared_radius()( to_exact(*data().p0),to_exact(*data().p1) ); 210 break; 211 case 3: 212 exact_ = Exact_squared_radius()( to_exact(*data().p0),to_exact(*data().p1),to_exact(*data().p2) ); 213 break; 214 default: 215 CGAL_assertion(false); 216 } 217 } 218 219 void set_approx() 220 { 221 switch(data().nbpts) { 222 case 1: 223 approx_ = Approx_squared_radius()( to_approx(*data().p0) ); 224 break; 225 case 2: 226 approx_ = Approx_squared_radius()( to_approx(*data().p0),to_approx(*data().p1) ); 227 break; 228 case 3: 229 approx_ = Approx_squared_radius()( to_approx(*data().p0),to_approx(*data().p1),to_approx(*data().p2) ); 230 break; 231 default: 232 CGAL_assertion(false); 233 } 234 } 235 236 const NT_exact& exact() const 237 { 238 if (exact_ == boost::none) { 239 update_exact(); 240 approx_=to_interval(*exact_); 241 } 242 return *exact_; 243 } 244 245 const NT_approx& approx() const { 246 return approx_; 247 } 248 249 //Constructors 250 Lazy_alpha_nt_2() 251 : exact_(Exact_nt(0)), approx_(0) 252 { 253 data().nbpts=0; 254 data().p0=nullptr; 255 data().p1=nullptr; 256 data().p2=nullptr; 257 } 258 259 Lazy_alpha_nt_2(double d) 260 : exact_(Exact_nt(d)), approx_(d) 261 { 262 data().nbpts=0; 263 data().p0=nullptr; 264 data().p1=nullptr; 265 data().p2=nullptr; 266 } 267 268 Lazy_alpha_nt_2(const Input_point& wp0) 269 { 270 data().nbpts=1; 271 data().p0=&wp0; 272 data().p1=nullptr; 273 data().p2=nullptr; 274 set_approx(); 275 } 276 277 Lazy_alpha_nt_2(const Input_point& wp0, 278 const Input_point& wp1) 279 { 280 data().nbpts=2; 281 data().p0=&wp0; 282 data().p1=&wp1; 283 data().p2=nullptr; 284 set_approx(); 285 } 286 287 Lazy_alpha_nt_2(const Input_point& wp0, 288 const Input_point& wp1, 289 const Input_point& wp2) 290 { 291 data().nbpts=3; 292 data().p0=&wp0; 293 data().p1=&wp1; 294 data().p2=&wp2; 295 set_approx(); 296 } 297 298 #define CGAL_LANT_COMPARE_FUNCTIONS(CMP) \ 299 bool \ 300 operator CMP (const Lazy_alpha_nt_2<Input_traits,mode,Weighted_tag> &other) const \ 301 { \ 302 Uncertain<bool> res = this->approx() CMP other.approx(); \ 303 if (res.is_certain()) \ 304 return res; \ 305 else \ 306 return this->exact() CMP other.exact(); \ 307 } \ 308 309 CGAL_LANT_COMPARE_FUNCTIONS(<) 310 CGAL_LANT_COMPARE_FUNCTIONS(>) 311 CGAL_LANT_COMPARE_FUNCTIONS(>=) 312 CGAL_LANT_COMPARE_FUNCTIONS(<=) 313 CGAL_LANT_COMPARE_FUNCTIONS(==) 314 CGAL_LANT_COMPARE_FUNCTIONS(!=) 315 316 #undef CGAL_LANT_COMPARE_FUNCTIONS 317 }; 318 319 template<class Input_traits, bool mode, class Weighted_tag> 320 std::ostream& 321 operator<< (std::ostream& os, 322 const Lazy_alpha_nt_2<Input_traits, mode, Weighted_tag>& a){ 323 return os << ::CGAL::to_double(a.approx()); 324 } 325 326 // small classes to select the functors in weighted or unweighted cases 327 template <class GeomTraits, class Weighted_tag> 328 struct iCompute_squared_radius_2; 329 330 template <class GeomTraits> 331 struct iCompute_squared_radius_2<GeomTraits, Tag_false /* Weighted_tag*/> 332 { 333 template <class As> 334 typename GeomTraits::Compute_squared_radius_2 335 operator()(const As& as) const { 336 return static_cast<const typename As::Triangulation&>(as).geom_traits(). 337 compute_squared_radius_2_object(); 338 } 339 }; 340 341 template <class GeomTraits> 342 struct iCompute_squared_radius_2<GeomTraits, Tag_true /* Weighted_tag*/> 343 { 344 template <class As> 345 typename GeomTraits::Compute_squared_radius_smallest_orthogonal_circle_2 346 operator()(const As& as) const { 347 return static_cast<const typename As::Triangulation&>(as).geom_traits(). 348 compute_squared_radius_smallest_orthogonal_circle_2_object(); 349 } 350 }; 351 352 template <class GeomTraits, class Weighted_tag> 353 struct iSide_of_bounded_circle_2; 354 355 template <class GeomTraits> 356 struct iSide_of_bounded_circle_2<GeomTraits, Tag_false /* Weighted_tag*/> 357 { 358 template <class As> 359 typename GeomTraits::Side_of_bounded_circle_2 360 operator()(const As& as) const { 361 return static_cast<const typename As::Triangulation&>(as).geom_traits(). 362 side_of_bounded_circle_2_object(); 363 } 364 }; 365 366 template <class GeomTraits> 367 struct iSide_of_bounded_circle_2<GeomTraits, Tag_true /* Weighted_tag*/> 368 { 369 template <class As> 370 typename GeomTraits::Power_side_of_bounded_power_circle_2 371 operator()(const As& as) const { 372 return static_cast<const typename As::Triangulation&>(as).geom_traits(). 373 power_side_of_bounded_power_circle_2_object(); 374 } 375 }; 376 377 template <class Type_of_alpha, class Point> 378 struct Lazy_compute_squared_radius_2 379 { 380 Type_of_alpha operator()(const Point& p, 381 const Point& q, 382 const Point& r) 383 { return Type_of_alpha(p,q,r); } 384 385 Type_of_alpha operator()(const Point& p, 386 const Point& q) 387 { return Type_of_alpha(p,q); } 388 389 Type_of_alpha operator()(const Point& p) 390 { return Type_of_alpha(p); } 391 }; 392 393 template <class GeomTraits, class ExactAlphaComparisonTag, class Weighted_tag> 394 struct Alpha_nt_selector_impl_2; 395 396 template <class GeomTraits, class Weighted_tag> 397 struct Alpha_nt_selector_impl_2<GeomTraits, 398 Tag_false /* ExactAlphaComparisonTag */, 399 Weighted_tag> 400 { 401 typedef typename GeomTraits::FT Type_of_alpha; 402 typedef iCompute_squared_radius_2<GeomTraits, Weighted_tag> Compute_squared_radius_2; 403 typedef iSide_of_bounded_circle_2<GeomTraits, Weighted_tag> Side_of_bounded_circle_2; 404 }; 405 406 template <class GeomTraits, class Weighted_tag> 407 struct Alpha_nt_selector_impl_2<GeomTraits, 408 Tag_true /* ExactAlphaComparisonTag */, 409 Weighted_tag> 410 { 411 typedef Lazy_alpha_nt_2< 412 GeomTraits, true /* mode */, Tag_false /* Weighted_tag */> Type_of_alpha; 413 typedef Lazy_compute_squared_radius_2< 414 Type_of_alpha, typename GeomTraits::Point_2> Functor; 415 416 struct Compute_squared_radius_2 417 { 418 template<class As> 419 Functor operator()(const As&){return Functor();} 420 }; 421 422 typedef iSide_of_bounded_circle_2<GeomTraits, Weighted_tag> Side_of_bounded_circle_2; 423 }; 424 425 template <class GeomTraits> 426 struct Alpha_nt_selector_impl_2<GeomTraits, 427 Tag_true /* ExactAlphaComparisonTag */, 428 Tag_true /* Weighted_tag */ > 429 { 430 typedef Lazy_alpha_nt_2< 431 GeomTraits, true /* mode */, Tag_true /* Weighted_tag */> Type_of_alpha; 432 433 typedef Lazy_compute_squared_radius_2< 434 Type_of_alpha, typename GeomTraits::Weighted_point_2> Functor; 435 436 struct Compute_squared_radius_2 437 { 438 template<class As> 439 Functor operator()(const As&){return Functor();} 440 }; 441 442 typedef iSide_of_bounded_circle_2<GeomTraits, Tag_true> Side_of_bounded_circle_2; 443 }; 444 445 template <class GeomTraits, class ExactAlphaComparisonTag, class Weighted_tag> 446 struct Alpha_nt_selector_2 447 : public Alpha_nt_selector_impl_2< 448 GeomTraits, 449 // If the base traits is already exact then we don't need to do anything, 450 // and we can simply directly use the traits class 451 Boolean_tag<boost::is_floating_point<typename GeomTraits::FT>::value && 452 ExactAlphaComparisonTag::value >, 453 Weighted_tag> 454 { }; 455 456 } // namespace internal 457 458 template<class Input_traits, bool mode, class Weighted_tag> 459 double to_double(const internal::Lazy_alpha_nt_2<Input_traits, mode, Weighted_tag>& a) 460 { 461 double r; 462 if (fit_in_double(a.approx(), r)) 463 return r; 464 465 // If it isn't precise enough, 466 // we trigger the exact computation first, 467 // which will refine the approximation. 468 if (!has_smaller_relative_precision(a.approx(), a.get_relative_precision_of_to_double())) 469 a.exact(); 470 471 return CGAL_NTS to_double(a.approx()); 472 } 473 474 } // namespace CGAL 475 476 #endif // CGAL_INTERNAL_LAZY_ALPHA_NT_2_H 477 478