|
BWAPI
|
00001 // Copyright (c) 1999-2004 Utrecht University (The Netherlands), 00002 // ETH Zurich (Switzerland), Freie Universitaet Berlin (Germany), 00003 // INRIA Sophia-Antipolis (France), Martin-Luther-University Halle-Wittenberg 00004 // (Germany), Max-Planck-Institute Saarbruecken (Germany), RISC Linz (Austria), 00005 // and Tel-Aviv University (Israel). All rights reserved. 00006 // 00007 // This file is part of CGAL (www.cgal.org); you can redistribute it and/or 00008 // modify it under the terms of the GNU Lesser General Public License as 00009 // published by the Free Software Foundation; version 2.1 of the License. 00010 // See the file LICENSE.LGPL distributed with CGAL. 00011 // 00012 // Licensees holding a valid commercial license may use this file in 00013 // accordance with the commercial license agreement provided with the software. 00014 // 00015 // This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE 00016 // WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. 00017 // 00018 // $URL: svn+ssh://scm.gforge.inria.fr/svn/cgal/branches/CGAL-3.5-branch/Homogeneous_kernel/include/CGAL/Homogeneous/function_objects.h $ 00019 // $Id: function_objects.h 50388 2009-07-06 11:36:23Z afabri $ 00020 // 00021 // 00022 // Author(s) : Stefan Schirra, Sylvain Pion, Michael Hoffmann 00023 00024 #ifndef CGAL_HOMOGENEOUS_FUNCTION_OBJECTS_H 00025 #define CGAL_HOMOGENEOUS_FUNCTION_OBJECTS_H 00026 00027 #include <CGAL/Kernel/function_objects.h> 00028 #include <CGAL/Cartesian/function_objects.h> 00029 #include <CGAL/Kernel/Return_base_tag.h> 00030 #include <CGAL/predicates/sign_of_determinant.h> 00031 00032 CGAL_BEGIN_NAMESPACE 00033 00034 namespace HomogeneousKernelFunctors { 00035 00036 using namespace CommonKernelFunctors; 00037 00038 // For lazyness... 00039 using CartesianKernelFunctors::Are_parallel_2; 00040 using CartesianKernelFunctors::Are_parallel_3; 00041 using CartesianKernelFunctors::Compute_squared_area_3; 00042 using CartesianKernelFunctors::Compare_squared_radius_3; 00043 using CartesianKernelFunctors::Collinear_3; 00044 using CartesianKernelFunctors::Construct_line_3; 00045 using CartesianKernelFunctors::Construct_equidistant_line_3; 00046 using CartesianKernelFunctors::Construct_barycenter_2; 00047 using CartesianKernelFunctors::Construct_barycenter_3; 00048 00049 using CartesianKernelFunctors::Compute_approximate_area_3; 00050 using CartesianKernelFunctors::Compute_approximate_squared_length_3; 00051 using CartesianKernelFunctors::Compute_area_divided_by_pi_3; 00052 using CartesianKernelFunctors::Compute_squared_length_divided_by_pi_square_3; 00053 using CartesianKernelFunctors::Construct_radical_plane_3; 00054 00055 template <typename K> 00056 class Angle_2 00057 { 00058 typedef typename K::Point_2 Point_2; 00059 typedef typename K::Construct_vector_2 Construct_vector_2; 00060 Construct_vector_2 c; 00061 public: 00062 typedef typename K::Angle result_type; 00063 00064 Angle_2() {} 00065 Angle_2(const Construct_vector_2& c_) : c(c_) {} 00066 00067 result_type 00068 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 00069 { return enum_cast<Angle>(CGAL_NTS sign(c(q,p) * c(q,r))); } 00070 // FIXME: scalar product 00071 }; 00072 00073 template <typename K> 00074 class Angle_3 00075 { 00076 typedef typename K::Point_3 Point_3; 00077 typedef typename K::Construct_vector_3 Construct_vector_3; 00078 Construct_vector_3 c; 00079 public: 00080 typedef typename K::Angle result_type; 00081 00082 Angle_3() {} 00083 Angle_3(const Construct_vector_3& c_) : c(c_) {} 00084 00085 result_type 00086 operator()(const Point_3& p, const Point_3& q, const Point_3& r) const 00087 { return enum_cast<Angle>(CGAL_NTS sign(c(q,p) * c(q,r))); } 00088 // FIXME: scalar product 00089 }; 00090 00091 00092 template <typename K> 00093 class Bounded_side_2 00094 { 00095 typedef typename K::Point_2 Point_2; 00096 typedef typename K::Circle_2 Circle_2; 00097 typedef typename K::Triangle_2 Triangle_2; 00098 typedef typename K::Iso_rectangle_2 Iso_rectangle_2; 00099 public: 00100 typedef typename K::Bounded_side result_type; 00101 00102 result_type 00103 operator()( const Circle_2& c, const Point_2& p) const 00104 { 00105 typename K::Compute_squared_distance_2 squared_distance; 00106 return enum_cast<Bounded_side>(CGAL_NTS compare(c.squared_radius(), 00107 squared_distance(c.center(),p))); 00108 } 00109 00110 result_type 00111 operator()( const Triangle_2& t, const Point_2& p) const 00112 { 00113 typename K::Collinear_are_ordered_along_line_2 00114 collinear_are_ordered_along_line; 00115 typename K::Orientation_2 orientation; 00116 typename K::Orientation o1 = orientation(t.vertex(0), t.vertex(1), p), 00117 o2 = orientation(t.vertex(1), t.vertex(2), p), 00118 o3 = orientation(t.vertex(2), t.vertex(3), p); 00119 00120 if (o2 == o1 && o3 == o1) 00121 return ON_BOUNDED_SIDE; 00122 return 00123 (o1 == COLLINEAR 00124 && collinear_are_ordered_along_line(t.vertex(0), p, t.vertex(1))) || 00125 (o2 == COLLINEAR 00126 && collinear_are_ordered_along_line(t.vertex(1), p, t.vertex(2))) || 00127 (o3 == COLLINEAR 00128 && collinear_are_ordered_along_line(t.vertex(2), p, t.vertex(3))) 00129 ? ON_BOUNDARY 00130 : ON_UNBOUNDED_SIDE; 00131 } 00132 00133 result_type 00134 operator()( const Iso_rectangle_2& r, const Point_2& p) const 00135 { 00136 return r.rep().bounded_side(p); 00137 } 00138 }; 00139 00140 template <typename K> 00141 class Bounded_side_3 00142 { 00143 typedef typename K::RT RT; 00144 typedef typename K::Point_3 Point_3; 00145 typedef typename K::Vector_3 Vector_3; 00146 typedef typename K::Sphere_3 Sphere_3; 00147 typedef typename K::Tetrahedron_3 Tetrahedron_3; 00148 typedef typename K::Iso_cuboid_3 Iso_cuboid_3; 00149 public: 00150 typedef typename K::Bounded_side result_type; 00151 00152 result_type 00153 operator()( const Sphere_3& s, const Point_3& p) const 00154 { return s.rep().bounded_side(p); } 00155 00156 result_type 00157 operator()( const Tetrahedron_3& t, const Point_3& p) const 00158 { 00159 Vector_3 v1 = t.vertex(1)-t.vertex(0); 00160 Vector_3 v2 = t.vertex(2)-t.vertex(0); 00161 Vector_3 v3 = t.vertex(3)-t.vertex(0); 00162 00163 Vector_3 vp = p - t.vertex(0); 00164 00165 // want to solve alpha*v1 + beta*v2 + gamma*v3 == vp 00166 // let vi' == vi*vi.hw() 00167 // we solve alpha'*v1' + beta'*v2' + gamma'*v3' == vp' / vp.hw() 00168 // muliplied by vp.hw() 00169 // then we have alpha = alpha'*v1.hw() / vp.hw() 00170 // and beta = beta' *v2.hw() / vp.hw() 00171 // and gamma = gamma'*v3.hw() / vp.hw() 00172 00173 const RT & v1x = v1.hx(); 00174 const RT & v1y = v1.hy(); 00175 const RT & v1z = v1.hz(); 00176 const RT & v2x = v2.hx(); 00177 const RT & v2y = v2.hy(); 00178 const RT & v2z = v2.hz(); 00179 const RT & v3x = v3.hx(); 00180 const RT & v3y = v3.hy(); 00181 const RT & v3z = v3.hz(); 00182 const RT & vpx = vp.hx(); 00183 const RT & vpy = vp.hy(); 00184 const RT & vpz = vp.hz(); 00185 00186 RT alpha = determinant( vpx, v2x, v3x, 00187 vpy, v2y, v3y, 00188 vpz, v2z, v3z ); 00189 RT beta = determinant( v1x, vpx, v3x, 00190 v1y, vpy, v3y, 00191 v1z, vpz, v3z ); 00192 RT gamma = determinant( v1x, v2x, vpx, 00193 v1y, v2y, vpy, 00194 v1z, v2z, vpz ); 00195 RT det = determinant( v1x, v2x, v3x, 00196 v1y, v2y, v3y, 00197 v1z, v2z, v3z ); 00198 00199 CGAL_kernel_assertion( det != 0 ); 00200 if (det < 0 ) 00201 { 00202 alpha = - alpha; 00203 beta = - beta ; 00204 gamma = - gamma; 00205 det = - det ; 00206 } 00207 00208 bool t1 = ( alpha < 0 ); 00209 bool t2 = ( beta < 0 ); 00210 bool t3 = ( gamma < 0 ); 00211 // t1 || t2 || t3 == not contained in cone 00212 00213 RT lhs = alpha*v1.hw() + beta*v2.hw() + gamma*v3.hw(); 00214 RT rhs = det * vp.hw(); 00215 00216 bool t4 = ( lhs > rhs ); 00217 // alpha + beta + gamma > 1 ? 00218 bool t5 = ( lhs < rhs ); 00219 // alpha + beta + gamma < 1 ? 00220 bool t6 = ( (alpha > 0) && (beta > 0) && (gamma > 0) ); 00221 00222 if ( t1 || t2 || t3 || t4 ) 00223 { 00224 return ON_UNBOUNDED_SIDE; 00225 } 00226 return (t5 && t6) ? ON_BOUNDED_SIDE : ON_BOUNDARY; 00227 } 00228 00229 result_type 00230 operator()( const Iso_cuboid_3& c, const Point_3& p) const 00231 { return c.rep().bounded_side(p); } 00232 }; 00233 00234 template <typename K> 00235 class Collinear_are_ordered_along_line_2 00236 { 00237 typedef typename K::RT RT; 00238 typedef typename K::Point_2 Point_2; 00239 #ifdef CGAL_kernel_exactness_preconditions 00240 typedef typename K::Collinear_2 Collinear_2; 00241 Collinear_2 c; 00242 #endif // CGAL_kernel_exactness_preconditions 00243 public: 00244 typedef typename K::Boolean result_type; 00245 00246 #ifdef CGAL_kernel_exactness_preconditions 00247 Collinear_are_ordered_along_line_2() {} 00248 Collinear_are_ordered_along_line_2(const Collinear_2& c_) : c(c_) {} 00249 #endif // CGAL_kernel_exactness_preconditions 00250 00251 result_type 00252 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 00253 { 00254 CGAL_kernel_exactness_precondition( c(p, q, r) ); 00255 00256 const RT& phx = p.hx(); 00257 const RT& phy = p.hy(); 00258 const RT& phw = p.hw(); 00259 const RT& qhx = q.hx(); 00260 const RT& qhy = q.hy(); 00261 const RT& qhw = q.hw(); 00262 const RT& rhx = r.hx(); 00263 const RT& rhy = r.hy(); 00264 const RT& rhw = r.hw(); 00265 00266 if ( !(phx * rhw == rhx * phw ) ) // non-vertical ? 00267 { 00268 return !( ( ( phx * qhw < qhx * phw) 00269 &&( rhx * qhw < qhx * rhw)) 00270 ||( ( qhx * phw < phx * qhw) 00271 &&( qhx * rhw < rhx * qhw)) ); 00272 } 00273 else if ( !(phy * rhw == rhy * phw ) ) 00274 { 00275 return !( ( ( phy * qhw < qhy * phw) 00276 &&( rhy * qhw < qhy * rhw)) 00277 ||( ( qhy * phw < phy * qhw) 00278 &&( qhy * rhw < rhy * qhw)) ); 00279 } 00280 else 00281 return (( phx*qhw == qhx*phw) && ( phy*qhw == qhy*phw)); 00282 } 00283 }; 00284 00285 template <typename K> 00286 class Collinear_are_ordered_along_line_3 00287 { 00288 typedef typename K::Point_3 Point_3; 00289 #ifdef CGAL_kernel_exactness_preconditions 00290 typedef typename K::Collinear_3 Collinear_3; 00291 Collinear_3 c; 00292 #endif // CGAL_kernel_exactness_preconditions 00293 public: 00294 typedef typename K::Boolean result_type; 00295 00296 #ifdef CGAL_kernel_exactness_preconditions 00297 Collinear_are_ordered_along_line_3() {} 00298 Collinear_are_ordered_along_line_3(const Collinear_3& c_) : c(c_) {} 00299 #endif // CGAL_kernel_exactness_preconditions 00300 00301 result_type 00302 operator()(const Point_3& p, const Point_3& q, const Point_3& r) const 00303 { 00304 CGAL_kernel_exactness_precondition( c(p, q, r) ); 00305 typedef typename K::RT RT; 00306 const RT & phx = p.hx(); 00307 const RT & phw = p.hw(); 00308 const RT & qhx = q.hx(); 00309 const RT & qhw = q.hw(); 00310 const RT & rhx = r.hx(); 00311 const RT & rhw = r.hw(); 00312 00313 const RT pqx = phx*qhw; 00314 const RT qpx = qhx*phw; 00315 const RT prx = phx*rhw; 00316 const RT qrx = qhx*rhw; 00317 const RT rqx = rhx*qhw; 00318 const RT rpx = rhx*phw; 00319 00320 if ( prx != rpx ) // px != rx 00321 { 00322 // (px <= qx)&&(qx <= rx) || (px >= qx)&&(qx >= rx) 00323 // !(((qx < px)||(rx < qx))&&((px < qx)||(qx < rx))) 00324 return ! ( ((qpx < pqx) || (rqx < qrx)) 00325 && ((pqx < qpx) || (qrx < rqx)) ); 00326 } 00327 00328 const RT & phy = p.hy(); 00329 const RT & qhy = q.hy(); 00330 const RT & rhy = r.hy(); 00331 00332 const RT pqy = phy*qhw; 00333 const RT qpy = qhy*phw; 00334 const RT pry = phy*rhw; 00335 const RT qry = qhy*rhw; 00336 const RT rqy = rhy*qhw; 00337 const RT rpy = rhy*phw; 00338 00339 if ( pry != rpy ) 00340 { 00341 return ! ( ((qpy < pqy) || (rqy < qry)) 00342 && ((pqy < qpy) || (qry < rqy)) ); 00343 } 00344 00345 const RT & phz = p.hz(); 00346 const RT & qhz = q.hz(); 00347 const RT & rhz = r.hz(); 00348 00349 const RT pqz = phz*qhw; 00350 const RT qpz = qhz*phw; 00351 const RT prz = phz*rhw; 00352 const RT qrz = qhz*rhw; 00353 const RT rqz = rhz*qhw; 00354 const RT rpz = rhz*phw; 00355 00356 if ( prz != rpz ) 00357 { 00358 return ! ( ((qpz < pqz) || (rqz < qrz)) 00359 && ((pqz < qpz) || (qrz < rqz)) ); 00360 } 00361 // p == r 00362 return ((rqx == qrx) && (rqy == qry) && (rqz == qrz)); 00363 } 00364 }; 00365 00366 template <typename K> 00367 class Collinear_are_strictly_ordered_along_line_2 00368 { 00369 typedef typename K::Point_2 Point_2; 00370 #ifdef CGAL_kernel_exactness_preconditions 00371 typedef typename K::Collinear_2 Collinear_2; 00372 Collinear_2 c; 00373 #endif // CGAL_kernel_exactness_preconditions 00374 public: 00375 typedef typename K::Boolean result_type; 00376 00377 #ifdef CGAL_kernel_exactness_preconditions 00378 Collinear_are_strictly_ordered_along_line_2() {} 00379 Collinear_are_strictly_ordered_along_line_2(const Collinear_2& c_) : c(c_) 00380 {} 00381 #endif // CGAL_kernel_exactness_preconditions 00382 00383 result_type 00384 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 00385 { 00386 CGAL_kernel_exactness_precondition( c(p, q, r) ); 00387 typedef typename K::RT RT; 00388 00389 const RT& phx = p.hx(); 00390 const RT& phy = p.hy(); 00391 const RT& phw = p.hw(); 00392 const RT& qhx = q.hx(); 00393 const RT& qhy = q.hy(); 00394 const RT& qhw = q.hw(); 00395 const RT& rhx = r.hx(); 00396 const RT& rhy = r.hy(); 00397 const RT& rhw = r.hw(); 00398 00399 if ( !(phx * rhw == rhx * phw ) ) 00400 { 00401 return ( ( phx * qhw < qhx * phw) 00402 &&( qhx * rhw < rhx * qhw)) 00403 ||( ( qhx * phw < phx * qhw) // ( phx * qhw > qhx * phw) 00404 &&( rhx * qhw < qhx * rhw)); // ( qhx * rhw > rhx * qhw) 00405 } 00406 else 00407 { 00408 return ( ( phy * qhw < qhy * phw) 00409 &&( qhy * rhw < rhy * qhw)) 00410 ||( ( qhy * phw < phy * qhw) // ( phy * qhw > qhy * phw) 00411 &&( rhy * qhw < qhy * rhw)); // ( qhy * rhw > rhy * qhw) 00412 } 00413 } 00414 }; 00415 00416 template <typename K> 00417 class Collinear_are_strictly_ordered_along_line_3 00418 { 00419 typedef typename K::Point_3 Point_3; 00420 typedef typename K::Direction_3 Direction_3; 00421 #ifdef CGAL_kernel_exactness_preconditions 00422 typedef typename K::Collinear_3 Collinear_3; 00423 Collinear_3 c; 00424 #endif // CGAL_kernel_exactness_preconditions 00425 public: 00426 typedef typename K::Boolean result_type; 00427 00428 #ifdef CGAL_kernel_exactness_preconditions 00429 Collinear_are_strictly_ordered_along_line_3() {} 00430 Collinear_are_strictly_ordered_along_line_3(const Collinear_3& c_) : c(c_) 00431 {} 00432 #endif // CGAL_kernel_exactness_preconditions 00433 00434 result_type 00435 operator()(const Point_3& p, const Point_3& q, const Point_3& r) const 00436 { 00437 CGAL_kernel_exactness_precondition( c(p, q, r) ); 00438 if ( p == r) return false; 00439 Direction_3 dir_pq = (p - q).direction(); 00440 Direction_3 dir_rq = (r - q).direction(); 00441 return (dir_pq == -dir_rq); 00442 } // FIXME 00443 }; 00444 00445 template <typename K> 00446 class Collinear_has_on_2 00447 { 00448 typedef typename K::Point_2 Point_2; 00449 typedef typename K::Direction_2 Direction_2; 00450 typedef typename K::Ray_2 Ray_2; 00451 typedef typename K::Segment_2 Segment_2; 00452 typedef typename K::Construct_point_on_2 Construct_point_on_2; 00453 typedef typename K::Compare_xy_2 Compare_xy_2; 00454 typedef typename K::Collinear_are_ordered_along_line_2 00455 Collinear_are_ordered_along_line_2; 00456 Collinear_are_ordered_along_line_2 co; 00457 Construct_point_on_2 cp; 00458 Compare_xy_2 cxy; 00459 public: 00460 typedef typename K::Boolean result_type; 00461 00462 Collinear_has_on_2() {} 00463 Collinear_has_on_2(const Construct_point_on_2& cp_, 00464 const Compare_xy_2& cxy_) 00465 : cp(cp_), cxy(cxy_) 00466 {} 00467 00468 result_type 00469 operator()( const Ray_2& r, const Point_2& p) const 00470 { 00471 const Point_2 & source = cp(r,0); 00472 return p == source || Direction_2(p - source) == r.direction(); 00473 } // FIXME 00474 00475 result_type 00476 operator()( const Segment_2& s, const Point_2& p) const 00477 { 00478 return co(cp(s,0), p, cp(s,1)); 00479 } 00480 }; 00481 00482 template <typename K> 00483 class Collinear_2 00484 { 00485 typedef typename K::Point_2 Point_2; 00486 typedef typename K::Orientation_2 Orientation_2; 00487 Orientation_2 o; 00488 public: 00489 typedef typename K::Boolean result_type; 00490 00491 Collinear_2() {} 00492 Collinear_2(const Orientation_2 o_) : o(o_) {} 00493 00494 result_type 00495 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 00496 { 00497 typedef typename K::RT RT; 00498 00499 const RT& phx = p.hx(); 00500 const RT& phy = p.hy(); 00501 const RT& phw = p.hw(); 00502 const RT& qhx = q.hx(); 00503 const RT& qhy = q.hy(); 00504 const RT& qhw = q.hw(); 00505 const RT& rhx = r.hx(); 00506 const RT& rhy = r.hy(); 00507 const RT& rhw = r.hw(); 00508 00509 // | A B | 00510 // | C D | 00511 00512 RT A = phx*rhw - phw*rhx; 00513 RT B = phy*rhw - phw*rhy; 00514 RT C = qhx*rhw - qhw*rhx; 00515 RT D = qhy*rhw - qhw*rhy; 00516 00517 RT det = A*D - B*C; 00518 00519 /* 00520 RT det_old = p.hx() * (q.hy()*r.hw() - q.hw()*r.hy() ) 00521 + p.hy() * (q.hw()*r.hx() - q.hx()*r.hw() ) 00522 + p.hw() * (q.hx()*r.hy() - q.hy()*r.hx() ); 00523 00524 if ( !(CGAL_NTS sign(det) == CGAL_NTS sign(det_old)) ) 00525 { 00526 std::cerr << "det: " << det << " det_old: " << det_old << flush; 00527 } 00528 */ 00529 00530 return CGAL_NTS is_zero(det); 00531 } 00532 }; 00533 00534 template <typename K> 00535 class Compare_angle_with_x_axis_2 00536 { 00537 typedef typename K::Point_2 Point_2; 00538 typedef typename K::Vector_2 Vector_2; 00539 typedef typename K::Direction_2 Direction_2; 00540 public: 00541 typedef typename K::Comparison_result result_type; 00542 00543 result_type 00544 operator()(const Direction_2& d1, const Direction_2& d2) const 00545 { 00546 typedef typename K::RT RT; 00547 CGAL_kernel_precondition( 00548 static_cast<int>(COUNTERCLOCKWISE) == static_cast<int>(LARGER) 00549 && static_cast<int>(COLLINEAR) == static_cast<int>(EQUAL) 00550 && static_cast<int>(CLOCKWISE) == static_cast<int>(SMALLER) ); 00551 00552 const RT RT0(0); 00553 00554 Vector_2 dirvec1(d1.x(), d1.y()); // Added 00555 Point_2 p1 = CGAL::ORIGIN + dirvec1; // Added 00556 Vector_2 dirvec2(d2.x(), d2.y()); // Added 00557 Point_2 p2 = ORIGIN + dirvec2; // Added 00558 // Point_2 p1 = ORIGIN + d1.vector(); // Commented out 00559 // Point_2 p2 = ORIGIN + d2.vector(); // Commented out 00560 00561 CGAL_kernel_precondition( RT0 < p1.hw() ); 00562 CGAL_kernel_precondition( RT0 < p2.hw() ); 00563 00564 int x_sign1 = static_cast<int>(CGAL_NTS sign( p1.hx() )); 00565 int x_sign2 = static_cast<int>(CGAL_NTS sign( p2.hx() )); 00566 int y_sign1 = static_cast<int>(CGAL_NTS sign( p1.hy() )); 00567 int y_sign2 = static_cast<int>(CGAL_NTS sign( p2.hy() )); 00568 00569 if ( y_sign1 * y_sign2 < 0) 00570 { 00571 return (0 < y_sign1 ) ? SMALLER : LARGER; 00572 } 00573 00574 Point_2 origin( RT0 , RT0 ); 00575 00576 if ( 0 < y_sign1 * y_sign2 ) 00577 { 00578 return orientation(origin, p2, p1); 00579 00580 // Precondition on the enums: 00581 // COUNTERCLOCKWISE == LARGER ( == 1 ) 00582 // COLLINEAR == EQUAL ( == 0 ) 00583 // CLOCKWISE == SMALLER ( == -1 ) 00584 } 00585 00586 // ( y_sign1 * y_sign2 == 0 ) 00587 00588 bool b1 = (y_sign1 == 0) && (x_sign1 >= 0); 00589 bool b2 = (y_sign2 == 0) && (x_sign2 >= 0); 00590 00591 if ( b1 ) { return b2 ? EQUAL : SMALLER; } 00592 if ( b2 ) { return b1 ? EQUAL : LARGER; } 00593 if ( y_sign1 == y_sign2 ) // == 0 00594 return EQUAL; 00595 else 00596 return (orientation(origin, p1, p2) == COUNTERCLOCKWISE) ? 00597 SMALLER : LARGER; 00598 } 00599 }; 00600 00601 template <typename K> 00602 class Compare_distance_2 00603 { 00604 typedef typename K::Point_2 Point_2; 00605 public: 00606 typedef typename K::Comparison_result result_type; 00607 00608 result_type 00609 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 00610 { 00611 typedef typename K::RT RT; 00612 00613 const RT & phx = p.hx(); 00614 const RT & phy = p.hy(); 00615 const RT & phw = p.hw(); 00616 const RT & qhx = q.hx(); 00617 const RT & qhy = q.hy(); 00618 const RT & qhw = q.hw(); 00619 const RT & rhx = r.hx(); 00620 const RT & rhy = r.hy(); 00621 const RT & rhw = r.hw(); 00622 00623 RT dosd = // difference of squared distances 00624 00625 // phx * phx * qhw * qhw * rhw * rhw 00626 // -RT(2) * phx * qhx * phw * qhw * rhw * rhw 00627 // + qhx * qhx * phw * phw * rhw * rhw 00628 // 00629 // + phy * phy * qhw * qhw * rhw * rhw 00630 // -RT(2) * phy * qhy * phw * qhw * rhw * rhw 00631 // + qhy * qhy * phw * phw * rhw * rhw 00632 // 00633 // - ( phx * phx * qhw * qhw * rhw * rhw 00634 // -RT(2) * phx * rhx * phw * qhw * qhw * rhw 00635 // + rhx * rhx * phw * phw * qhw * qhw 00636 // 00637 // + phy * phy * qhw * qhw * rhw * rhw 00638 // -RT(2) * phy * rhy * phw * qhw * qhw * rhw 00639 // + rhy * rhy * phw * phw * qhw * qhw 00640 00641 rhw*rhw * ( phw * ( qhx*qhx + qhy*qhy ) 00642 - 2 * qhw * ( phx*qhx + phy*qhy ) 00643 ) 00644 - qhw*qhw * ( phw * ( rhx*rhx + rhy*rhy ) 00645 - 2 * rhw * ( phx*rhx + phy*rhy ) 00646 ); 00647 00648 return CGAL_NTS sign(dosd); 00649 } 00650 }; 00651 00652 template <typename K> 00653 class Compare_squared_distance_2 00654 { 00655 typedef typename K::Point_2 Point_2; 00656 typedef typename K::FT FT; 00657 public: 00658 typedef typename K::Comparison_result result_type; 00659 00660 result_type 00661 operator()(const Point_2& p, const Point_2& q, const FT& d2) const 00662 { 00663 return CGAL_NTS compare(squared_distance(p, q), d2); 00664 } 00665 }; 00666 00667 template <typename K> 00668 class Compare_distance_3 00669 { 00670 typedef typename K::Point_3 Point_3; 00671 public: 00672 typedef typename K::Comparison_result result_type; 00673 00674 result_type 00675 operator()(const Point_3& p, const Point_3& q, const Point_3& r) const 00676 { 00677 typedef typename K::RT RT; 00678 00679 const RT & phx = p.hx(); 00680 const RT & phy = p.hy(); 00681 const RT & phz = p.hz(); 00682 const RT & phw = p.hw(); 00683 const RT & qhx = q.hx(); 00684 const RT & qhy = q.hy(); 00685 const RT & qhz = q.hz(); 00686 const RT & qhw = q.hw(); 00687 const RT & rhx = r.hx(); 00688 const RT & rhy = r.hy(); 00689 const RT & rhz = r.hz(); 00690 const RT & rhw = r.hw(); 00691 00692 RT dosd = // difference of squared distances 00693 00694 rhw*rhw * ( phw * ( qhx*qhx + qhy*qhy + qhz*qhz ) 00695 - 2 * qhw * ( phx*qhx + phy*qhy + phz*qhz ) 00696 ) 00697 - qhw*qhw * ( phw * ( rhx*rhx + rhy*rhy + rhz*rhz ) 00698 - 2 * rhw * ( phx*rhx + phy*rhy + phz*rhz ) 00699 ); 00700 00701 return CGAL_NTS sign(dosd); 00702 } 00703 }; 00704 00705 template <typename K> 00706 class Compare_squared_distance_3 00707 { 00708 typedef typename K::Point_3 Point_3; 00709 typedef typename K::FT FT; 00710 public: 00711 typedef typename K::Comparison_result result_type; 00712 00713 result_type 00714 operator()(const Point_3& p, const Point_3& q, const FT& d2) const 00715 { 00716 return CGAL_NTS compare(squared_distance(p, q), d2); 00717 } 00718 }; 00719 00720 template <typename K> 00721 class Compare_slope_2 00722 { 00723 typedef typename K::Line_2 Line_2; 00724 typedef typename K::Segment_2 Segment_2; 00725 public: 00726 typedef typename K::Comparison_result result_type; 00727 00728 result_type 00729 operator()(const Line_2& l1, const Line_2& l2) const 00730 { 00731 typedef typename K::RT RT; 00732 if (l1.is_horizontal()) 00733 return l2.is_vertical() ? 00734 SMALLER : CGAL_NTS sign(l2.a()) * CGAL_NTS sign(l2.b()); 00735 if (l2.is_horizontal()) 00736 return l1.is_vertical() ? 00737 LARGER : - CGAL_NTS sign(l1.a()) * CGAL_NTS sign(l1.b()); 00738 if (l1.is_vertical()) return l2.is_vertical() ? EQUAL : LARGER; 00739 if (l2.is_vertical()) return SMALLER; 00740 int l1_sign = CGAL_NTS sign(-l1.a() * l1.b()); 00741 int l2_sign = CGAL_NTS sign(-l2.a() * l2.b()); 00742 00743 if (l1_sign < l2_sign) return SMALLER; 00744 if (l1_sign > l2_sign) return LARGER; 00745 00746 if (l1_sign > 0) 00747 return CGAL_NTS compare( CGAL_NTS abs(l1.a() * l2.b()), 00748 CGAL_NTS abs(l2.a() * l1.b()) ); 00749 00750 return CGAL_NTS compare( CGAL_NTS abs(l2.a() * l1.b()), 00751 CGAL_NTS abs(l1.a() * l2.b()) ); 00752 } // FIXME 00753 00754 result_type 00755 operator()(const Segment_2& s1, const Segment_2& s2) const 00756 { 00757 typedef typename K::FT FT; 00758 00759 typename K::Comparison_result cmp_y1 = compare_y(s1.source(), s1.target()); 00760 if (cmp_y1 == EQUAL) // horizontal 00761 { 00762 typename K::Comparison_result cmp_x2 = compare_x(s2.source(), s2.target()); 00763 00764 if (cmp_x2 == EQUAL) return SMALLER; 00765 FT s_hw = s2.source().hw(); 00766 FT t_hw = s2.target().hw(); 00767 return - CGAL_NTS sign(s2.source().hy()*t_hw - s2.target().hy()*s_hw) * 00768 CGAL_NTS sign(s2.source().hx()*t_hw - s2.target().hx()*s_hw); 00769 } 00770 00771 typename K::Comparison_result cmp_y2 = compare_y(s2.source(), s2.target()); 00772 if (cmp_y2 == EQUAL) 00773 { 00774 typename K::Comparison_result cmp_x1 = compare_x(s1.source(), s1.target()); 00775 00776 if (cmp_x1 == EQUAL) return LARGER; 00777 FT s_hw = s1.source().hw(); 00778 FT t_hw = s1.target().hw(); 00779 return CGAL_NTS sign(s1.source().hy()*t_hw - s1.target().hy()*s_hw) * 00780 CGAL_NTS sign(s1.source().hx()*t_hw - s1.target().hx()*s_hw); 00781 } 00782 00783 typename K::Comparison_result cmp_x1 = compare_x(s1.source(), s1.target()); 00784 typename K::Comparison_result cmp_x2 = compare_x(s2.source(), s2.target()); 00785 if (cmp_x1 == EQUAL) 00786 return cmp_x2 == EQUAL ? EQUAL : LARGER; 00787 00788 if (cmp_x2 == EQUAL) return SMALLER; 00789 00790 FT s1_s_hw = s1.source().hw(); 00791 FT s1_t_hw = s1.target().hw(); 00792 FT s2_s_hw = s2.source().hw(); 00793 FT s2_t_hw = s2.target().hw(); 00794 FT s1_xdiff = s1.source().hx()*s1_t_hw - s1.target().hx()*s1_s_hw; 00795 FT s1_ydiff = s1.source().hy()*s1_t_hw - s1.target().hy()*s1_s_hw; 00796 FT s2_xdiff = s2.source().hx()*s2_t_hw - s2.target().hx()*s2_s_hw; 00797 FT s2_ydiff = s2.source().hy()*s2_t_hw - s2.target().hy()*s2_s_hw; 00798 typename K::Sign s1_sign = CGAL_NTS sign(s1_ydiff * s1_xdiff); 00799 typename K::Sign s2_sign = CGAL_NTS sign(s2_ydiff * s2_xdiff); 00800 00801 if (s1_sign < s2_sign) return SMALLER; 00802 if (s1_sign > s2_sign) return LARGER; 00803 00804 if (s1_sign > 0) 00805 return CGAL_NTS sign(CGAL_NTS abs(s1_ydiff * s2_xdiff) - 00806 CGAL_NTS abs(s2_ydiff * s1_xdiff)); 00807 00808 return CGAL_NTS sign(CGAL_NTS abs(s2_ydiff * s1_xdiff) - 00809 CGAL_NTS abs(s1_ydiff * s2_xdiff)); 00810 } 00811 }; 00812 00813 template <typename K> 00814 class Compare_x_at_y_2 00815 { 00816 typedef typename K::Point_2 Point_2; 00817 typedef typename K::Line_2 Line_2; 00818 public: 00819 typedef typename K::Comparison_result result_type; 00820 00821 result_type 00822 operator()( const Point_2& p, const Line_2& h) const 00823 { 00824 typedef typename K::RT RT; 00825 CGAL_kernel_precondition( ! h.is_horizontal() ); 00826 typename K::Oriented_side ors = h.oriented_side( p ); 00827 if ( h.a() < RT(0) ) 00828 ors = -ors; 00829 if ( ors == ON_POSITIVE_SIDE ) 00830 return LARGER; 00831 return ( ors == ON_NEGATIVE_SIDE ) ? SMALLER : EQUAL; 00832 } // FIXME 00833 00834 result_type 00835 operator()( const Point_2& p, const Line_2& h1, const Line_2& h2) const 00836 { return CGAL_NTS compare(h1.x_at_y( p.y() ), h2.x_at_y( p.y() )); } 00837 // FIXME 00838 00839 result_type 00840 operator()( const Line_2& l1, const Line_2& l2, const Line_2& h) const 00841 { return compare_x_at_y( gp_linear_intersection( l1, l2 ), h); } 00842 // FIXME 00843 00844 result_type 00845 operator()( const Line_2& l1, const Line_2& l2, 00846 const Line_2& h1, const Line_2& h2) const 00847 { return compare_x_at_y( gp_linear_intersection( l1, l2 ), h1, h2 ); } 00848 // FIXME 00849 }; 00850 00851 template <typename K> 00852 class Compare_xyz_3 00853 { 00854 typedef typename K::Point_3 Point_3; 00855 public: 00856 typedef typename K::Comparison_result result_type; 00857 00858 result_type 00859 operator()( const Point_3& p, const Point_3& q) const 00860 { 00861 typedef typename K::RT RT; 00862 RT pV = p.hx()*q.hw(); 00863 RT qV = q.hx()*p.hw(); 00864 if ( pV < qV ) 00865 { 00866 return SMALLER; 00867 } 00868 if ( qV < pV ) // ( pV > qV ) 00869 { 00870 return LARGER; 00871 } 00872 // same x 00873 pV = p.hy()*q.hw(); 00874 qV = q.hy()*p.hw(); 00875 if ( pV < qV ) 00876 { 00877 return SMALLER; 00878 } 00879 if ( qV < pV ) // ( pV > qV ) 00880 { 00881 return LARGER; 00882 } 00883 // same x and y 00884 pV = p.hz()*q.hw(); 00885 qV = q.hz()*p.hw(); 00886 return CGAL_NTS compare(pV, qV); 00887 } 00888 }; 00889 00890 template <typename K> 00891 class Compare_xy_2 00892 { 00893 typedef typename K::Point_2 Point_2; 00894 public: 00895 typedef typename K::Comparison_result result_type; 00896 00897 result_type 00898 operator()( const Point_2& p, const Point_2& q) const 00899 { 00900 typedef typename K::RT RT; 00901 00902 const RT& phx = p.hx(); 00903 const RT& phy = p.hy(); 00904 const RT& phw = p.hw(); 00905 const RT& qhx = q.hx(); 00906 const RT& qhy = q.hy(); 00907 const RT& qhw = q.hw(); 00908 00909 RT pV = phx*qhw; 00910 RT qV = qhx*phw; 00911 if ( pV == qV ) 00912 { 00913 pV = phy*qhw; 00914 qV = qhy*phw; 00915 } 00916 return CGAL_NTS compare(pV, qV); 00917 } 00918 }; 00919 00920 template <typename K> 00921 class Compare_yx_2 00922 { 00923 typedef typename K::Point_2 Point_2; 00924 public: 00925 typedef typename K::Comparison_result result_type; 00926 00927 result_type 00928 operator()( const Point_2& p, const Point_2& q) const 00929 { 00930 typedef typename K::RT RT; 00931 00932 const RT& phx = p.hx(); 00933 const RT& phy = p.hy(); 00934 const RT& phw = p.hw(); 00935 const RT& qhx = q.hx(); 00936 const RT& qhy = q.hy(); 00937 const RT& qhw = q.hw(); 00938 00939 RT pV = phy*qhw; 00940 RT qV = qhy*phw; 00941 if ( pV == qV ) 00942 { 00943 pV = phx*qhw; 00944 qV = qhx*phw; 00945 } 00946 return CGAL_NTS compare(pV, qV); 00947 } 00948 }; 00949 00950 template <typename K> 00951 class Compare_xy_3 00952 { 00953 typedef typename K::Point_3 Point_3; 00954 public: 00955 typedef typename K::Comparison_result result_type; 00956 00957 result_type 00958 operator()( const Point_3& p, const Point_3& q) const 00959 { 00960 typedef typename K::RT RT; 00961 RT pV = p.hx()*q.hw(); 00962 RT qV = q.hx()*p.hw(); 00963 if ( pV < qV ) 00964 { 00965 return SMALLER; 00966 } 00967 if ( qV < pV ) // ( pV > qV ) 00968 { 00969 return LARGER; 00970 } 00971 // same x 00972 pV = p.hy()*q.hw(); 00973 qV = q.hy()*p.hw(); 00974 return CGAL_NTS compare(pV, qV); 00975 } 00976 }; 00977 00978 template <typename K> 00979 class Compare_x_2 00980 { 00981 typedef typename K::Point_2 Point_2; 00982 typedef typename K::Line_2 Line_2; 00983 public: 00984 typedef typename K::Comparison_result result_type; 00985 00986 result_type 00987 operator()( const Point_2& p, const Point_2& q) const 00988 { 00989 return CGAL_NTS compare(p.hx()*q.hw(), q.hx()*p.hw()); 00990 } 00991 00992 result_type 00993 operator()( const Point_2& p, const Line_2& l1, const Line_2& l2) const 00994 { 00995 Point_2 ip = gp_linear_intersection( l1, l2 ); 00996 return this->operator()(p, ip); 00997 } // FIXME 00998 00999 result_type 01000 operator()( const Line_2& l, const Line_2& h1, const Line_2& h2) const 01001 { 01002 return this->operator()(l, h1, l, h2); 01003 } // FIXME 01004 01005 result_type 01006 operator()( const Line_2& l1, const Line_2& l2, 01007 const Line_2& h1, const Line_2& h2) const 01008 { 01009 Point_2 lip = gp_linear_intersection( l1, l2 ); 01010 Point_2 hip = gp_linear_intersection( h1, h2 ); 01011 return this->operator()(lip, hip); 01012 } // FIXME 01013 }; 01014 01015 template <typename K> 01016 class Compare_x_3 01017 { 01018 typedef typename K::Point_3 Point_3; 01019 public: 01020 typedef typename K::Comparison_result result_type; 01021 01022 result_type 01023 operator()( const Point_3& p, const Point_3& q) const 01024 { return CGAL_NTS compare(p.hx() * q.hw(), q.hx() * p.hw() ); } 01025 }; 01026 01027 template <typename K> 01028 class Compare_y_at_x_2 01029 { 01030 typedef typename K::Point_2 Point_2; 01031 typedef typename K::Line_2 Line_2; 01032 typedef typename K::Segment_2 Segment_2; 01033 public: 01034 typedef typename K::Comparison_result result_type; 01035 01036 result_type 01037 operator()( const Point_2& p, const Line_2& h) const 01038 { 01039 typedef typename K::RT RT; 01040 CGAL_kernel_precondition( ! h.is_vertical() ); 01041 typename K::Oriented_side ors = h.oriented_side( p ); 01042 if ( h.b() < 0 ) 01043 ors = -ors; 01044 return ors; 01045 } // FIXME 01046 01047 result_type 01048 operator()( const Point_2& p, const Line_2& h1, const Line_2& h2) const 01049 { return CGAL_NTS compare(h1.y_at_x( p.x() ), h2.y_at_x( p.x() )); } 01050 // FIXME 01051 01052 result_type 01053 operator()( const Line_2& l1, const Line_2& l2, const Line_2& h) const 01054 { return compare_y_at_x( gp_linear_intersection( l1, l2 ), h); } 01055 // FIXME 01056 01057 result_type 01058 operator()( const Line_2& l1, const Line_2& l2, 01059 const Line_2& h1, const Line_2& h2) const 01060 { return compare_y_at_x( gp_linear_intersection( l1, l2 ), h1, h2 ); } 01061 // FIXME 01062 01063 result_type 01064 operator()( const Point_2& p, const Segment_2& s) const 01065 { 01066 // compares the y-coordinates of p and the vertical projection of p on s. 01067 // Precondition : p is in the x-range of s. 01068 01069 if (compare_x(s.source(), s.target()) == SMALLER) { 01070 CGAL_kernel_precondition(compare_x(s.source(), p) != LARGER 01071 && compare_x(p, s.target()) != LARGER); 01072 return (Comparison_result) orientation(p, s.source(), s.target()); 01073 } 01074 else if (compare_x(s.source(), s.target()) == LARGER) { 01075 CGAL_kernel_precondition(compare_x(s.target(), p) != LARGER 01076 && compare_x(p, s.source()) != LARGER); 01077 return (Comparison_result) orientation(p, s.target(), s.source()); 01078 } 01079 else { 01080 CGAL_kernel_precondition(compare_x(s.target(), p) == EQUAL); 01081 if (compare_y(p, s.source()) == SMALLER && 01082 compare_y(p, s.target()) == SMALLER) 01083 return SMALLER; 01084 if (compare_y(p, s.source()) == LARGER && 01085 compare_y(p, s.target()) == LARGER) 01086 return LARGER; 01087 return EQUAL; 01088 } 01089 } // FIXME 01090 01091 result_type 01092 operator()( const Point_2& p, 01093 const Segment_2& s1, const Segment_2& s2) const 01094 { 01095 // compares the y-coordinates of the vertical projections 01096 // of p on s1 and s2 01097 // Precondition : p is in the x-range of s1 and s2. 01098 // - if one or two segments are vertical : 01099 // - if the segments intersect, return EQUAL 01100 // - if not, return the obvious SMALLER/LARGER. 01101 01102 typedef typename K::FT FT; 01103 FT px = p.x(); 01104 FT s1sx = s1.source().x(); 01105 FT s1sy = s1.source().y(); 01106 FT s1tx = s1.target().x(); 01107 FT s1ty = s1.target().y(); 01108 FT s2sx = s2.source().x(); 01109 FT s2sy = s2.source().y(); 01110 FT s2tx = s2.target().x(); 01111 FT s2ty = s2.target().y(); 01112 01113 CGAL_kernel_precondition(px >= (CGAL::min)(s1sx, s1tx) && 01114 px <= (CGAL::max)(s1sx, s1tx)); 01115 CGAL_kernel_precondition(px >= (CGAL::min)(s2sx, s2tx) && 01116 px <= (CGAL::max)(s2sx, s2tx)); 01117 01118 if (s1sx != s1tx && s2sx != s2tx) { 01119 FT s1stx = s1sx-s1tx; 01120 FT s2stx = s2sx-s2tx; 01121 01122 return CGAL_NTS compare(s1sx, s1tx) * 01123 CGAL_NTS compare(s2sx, s2tx) * 01124 CGAL_NTS compare(-(s1sx-px)*(s1sy-s1ty)*s2stx, 01125 (s2sy-s1sy)*s2stx*s1stx 01126 -(s2sx-px)*(s2sy-s2ty)*s1stx); 01127 } 01128 else { 01129 if (s1sx == s1tx) { // s1 is vertical 01130 typename K::Comparison_result c1, c2; 01131 c1 = compare_y_at_x(s1.source(), s2); 01132 c2 = compare_y_at_x(s1.target(), s2); 01133 if (c1 == c2) 01134 return c1; 01135 return EQUAL; 01136 } 01137 // s2 is vertical 01138 typename K::Comparison_result c3, c4; 01139 c3 = compare_y_at_x(s2.source(), s1); 01140 c4 = compare_y_at_x(s2.target(), s1); 01141 if (c3 == c4) 01142 return -c3; 01143 return EQUAL; 01144 } 01145 } // FIXME 01146 }; 01147 01148 template <typename K> 01149 class Compare_y_2 01150 { 01151 typedef typename K::Point_2 Point_2; 01152 typedef typename K::Line_2 Line_2; 01153 public: 01154 typedef typename K::Comparison_result result_type; 01155 01156 result_type 01157 operator()( const Point_2& p, const Point_2& q) const 01158 { 01159 typedef typename K::RT RT; 01160 01161 const RT& phy = p.hy(); 01162 const RT& phw = p.hw(); 01163 const RT& qhy = q.hy(); 01164 const RT& qhw = q.hw(); 01165 return CGAL_NTS compare(phy * qhw, qhy * phw); 01166 } 01167 01168 result_type 01169 operator()( const Point_2& p, const Line_2& l1, const Line_2& l2) const 01170 { 01171 Point_2 ip = gp_linear_intersection( l1, l2 ); 01172 return compare_y( p, ip ); 01173 } // FIXME 01174 01175 result_type 01176 operator()( const Line_2& l, const Line_2& h1, const Line_2& h2) const 01177 { 01178 return this->operator()(l, h1, l, h2); 01179 } 01180 01181 result_type 01182 operator()( const Line_2& l1, const Line_2& l2, 01183 const Line_2& h1, const Line_2& h2) const 01184 { 01185 Point_2 lip = gp_linear_intersection( l1, l2 ); 01186 Point_2 hip = gp_linear_intersection( h1, h2 ); 01187 return this->operator()( lip, hip ); 01188 } // FIXME 01189 }; 01190 01191 template <typename K> 01192 class Compare_y_3 01193 { 01194 typedef typename K::Point_3 Point_3; 01195 public: 01196 typedef typename K::Comparison_result result_type; 01197 01198 result_type 01199 operator()( const Point_3& p, const Point_3& q) const 01200 { return CGAL_NTS compare(p.hy() * q.hw(), q.hy() * p.hw() ); } 01201 }; 01202 01203 template <typename K> 01204 class Compare_z_3 01205 { 01206 typedef typename K::Point_3 Point_3; 01207 public: 01208 typedef typename K::Comparison_result result_type; 01209 01210 result_type 01211 operator()( const Point_3& p, const Point_3& q) const 01212 { return CGAL_NTS compare(p.hz() * q.hw(), q.hz() * p.hw() ); } 01213 }; 01214 01215 template <typename K> 01216 class Compute_area_2 01217 { 01218 typedef typename K::RT RT; 01219 typedef typename K::FT FT; 01220 typedef typename K::Iso_rectangle_2 Iso_rectangle_2; 01221 typedef typename K::Triangle_2 Triangle_2; 01222 typedef typename K::Point_2 Point_2; 01223 typedef typename K::Vector_2 Vector_2; 01224 typedef typename K::Construct_vector_2 Construct_vector_2; 01225 Construct_vector_2 co; 01226 public: 01227 typedef FT result_type; 01228 01229 FT 01230 operator()( const Point_2& p, const Point_2& q, const Point_2& r ) const 01231 { 01232 Vector_2 v1 = co(p, q); 01233 Vector_2 v2 = co(p, r); 01234 01235 RT num = v1.hx()*v2.hy() - v2.hx()*v1.hy(); 01236 RT den = RT(2) * v1.hw() * v2.hw(); 01237 return FT(num)/FT(den); 01238 } 01239 01240 FT 01241 operator()( const Iso_rectangle_2& r ) const 01242 { return (r.xmax()-r.xmin()) * (r.ymax()-r.ymin()); } 01243 01244 FT 01245 operator()( const Triangle_2& t ) const 01246 { return t.area(); } 01247 }; 01248 01249 template <typename K> 01250 class Compute_determinant_2 01251 { 01252 typedef typename K::FT FT; 01253 typedef typename K::Vector_2 Vector_2; 01254 public: 01255 typedef FT result_type; 01256 01257 result_type 01258 operator()(const Vector_2& v, const Vector_2& w) const 01259 { 01260 return determinant(v.hx(), v.hy(), 01261 w.hx(), w.hy()) / FT(v.hw() * w.hw()); 01262 } 01263 }; 01264 01265 template <typename K> 01266 class Compute_determinant_3 01267 { 01268 typedef typename K::FT FT; 01269 typedef typename K::Vector_3 Vector_3; 01270 public: 01271 typedef FT result_type; 01272 01273 result_type 01274 operator()(const Vector_3& v, const Vector_3& w, const Vector_3& t) const 01275 { 01276 return determinant(v.hx(), v.hy(), v.hz(), 01277 w.hx(), w.hy(), w.hz(), 01278 t.hx(), t.hy(), t.hz()) 01279 / FT(v.hw() * w.hw() * t.hw()); 01280 } 01281 }; 01282 01283 template <typename K> 01284 class Compute_scalar_product_2 01285 { 01286 typedef typename K::RT RT; 01287 typedef typename K::FT FT; 01288 typedef typename K::Vector_2 Vector_2; 01289 public: 01290 typedef FT result_type; 01291 01292 FT 01293 operator()(const Vector_2& v, const Vector_2& w) const 01294 { 01295 return FT( RT(v.hx()*w.hx() + v.hy()*w.hy()) ) / 01296 FT( RT(v.hw()*w.hw() ) ); 01297 } 01298 }; 01299 01300 template <typename K> 01301 class Compute_scalar_product_3 01302 { 01303 typedef typename K::RT RT; 01304 typedef typename K::FT FT; 01305 typedef typename K::Vector_3 Vector_3; 01306 public: 01307 typedef FT result_type; 01308 01309 FT 01310 operator()(const Vector_3& v, const Vector_3& w) const 01311 { 01312 return FT( RT(v.hx()*w.hx() + v.hy()*w.hy()) + v.hz()*w.hz() ) / 01313 FT( RT(v.hw()*w.hw() ) ); 01314 } 01315 }; 01316 01317 // TODO ... 01318 template <typename K> 01319 class Compute_squared_radius_2 01320 { 01321 typedef typename K::FT FT; 01322 typedef typename K::Point_2 Point_2; 01323 typedef typename K::Circle_2 Circle_2; 01324 public: 01325 typedef FT result_type; 01326 01327 FT 01328 operator()( const Circle_2& c) const 01329 { return c.rep().squared_radius(); } 01330 01331 FT 01332 operator()( const Point_2& p) const 01333 { return FT(0); } 01334 01335 FT 01336 operator()( const Point_2& p, const Point_2& q) const 01337 { 01338 typedef typename K::FT FT; 01339 return squared_distance(p, q)/FT(4); 01340 } // FIXME 01341 01342 FT 01343 operator()( const Point_2& p, const Point_2& q, const Point_2& r) const 01344 { return squared_distance(p, circumcenter(p, q, r)); } 01345 // FIXME 01346 }; 01347 01348 template <typename K> 01349 class Compute_squared_radius_3 01350 { 01351 typedef typename K::FT FT; 01352 typedef typename K::Point_3 Point_3; 01353 typedef typename K::Sphere_3 Sphere_3; 01354 public: 01355 typedef FT result_type; 01356 01357 FT 01358 operator()( const Sphere_3& s) const 01359 { return s.rep().squared_radius(); } 01360 01361 FT 01362 operator()( const Point_3& p) const 01363 { return FT(0); } 01364 01365 FT 01366 operator()( const Point_3& p, const Point_3& q) const 01367 { 01368 typedef typename K::FT FT; 01369 return squared_distance(p, q) / FT(4); 01370 } // FIXME 01371 01372 FT 01373 operator()( const Point_3& p, const Point_3& q, const Point_3& r) const 01374 { 01375 return squared_distance(p, circumcenter(p, q, r)); 01376 } // FIXME 01377 01378 FT 01379 operator()( const Point_3& p, const Point_3& q, 01380 const Point_3& r, const Point_3& s) const 01381 { 01382 return squared_distance(p, circumcenter(p, q, r, s)); 01383 } // FIXME 01384 }; 01385 01386 template <typename K> 01387 class Compute_volume_3 01388 { 01389 typedef typename K::FT FT; 01390 typedef typename K::Point_3 Point_3; 01391 typedef typename K::Vector_3 Vector_3; 01392 typedef typename K::Tetrahedron_3 Tetrahedron_3; 01393 typedef typename K::Iso_cuboid_3 Iso_cuboid_3; 01394 public: 01395 typedef FT result_type; 01396 01397 FT 01398 operator()(const Point_3& p0, const Point_3& p1, 01399 const Point_3& p2, const Point_3& p3) const 01400 { 01401 Vector_3 vec1 = p1 - p0; 01402 Vector_3 vec2 = p2 - p0; 01403 Vector_3 vec3 = p3 - p0; 01404 01405 // first compute (vec1.hw * vec2.hw * vec3.hw * det(vec1, vec2, vec3)) 01406 // then divide by (6 * vec1.hw * vec2.hw * vec3.hw) 01407 const FT w123 (vec1.hw() * vec2.hw() * vec3.hw()); 01408 const FT& hx1 = vec1.hx(); 01409 const FT& hy1 = vec1.hy(); 01410 const FT& hz1 = vec1.hz(); 01411 const FT& hx2 = vec2.hx(); 01412 const FT& hy2 = vec2.hy(); 01413 const FT& hz2 = vec2.hz(); 01414 const FT& hx3 = vec3.hx(); 01415 const FT& hy3 = vec3.hy(); 01416 const FT& hz3 = vec3.hz(); 01417 01418 return ( (hx1 * (hy2 * hz3 - hy3 * hz2)) 01419 - (hy1 * (hx2 * hz3 - hx3 * hz2)) 01420 + (hz1 * (hx2 * hy3 - hx3 * hy2)))/ (6 * w123); 01421 } 01422 01423 FT 01424 operator()( const Tetrahedron_3& t ) const 01425 { 01426 return this->operator()(t.vertex(0), t.vertex(1), 01427 t.vertex(2), t.vertex(3)); 01428 } 01429 01430 FT 01431 operator()( const Iso_cuboid_3& c ) const 01432 { return c.rep().volume(); } 01433 }; 01434 01435 01436 template <typename K> 01437 class Compute_x_2 01438 { 01439 typedef typename K::FT FT; 01440 typedef typename K::Point_2 Point_2; 01441 typedef typename K::Vector_2 Vector_2; 01442 01443 public: 01444 typedef FT result_type; 01445 01446 FT 01447 operator()(const Point_2& p) const 01448 { 01449 return p.rep().x(); 01450 } 01451 01452 FT 01453 operator()(const Vector_2& v) const 01454 { 01455 return v.rep().x(); 01456 } 01457 }; 01458 01459 template <typename K> 01460 class Compute_x_3 01461 { 01462 typedef typename K::FT FT; 01463 typedef typename K::Point_3 Point_3; 01464 typedef typename K::Vector_3 Vector_3; 01465 01466 public: 01467 typedef FT result_type; 01468 01469 FT 01470 operator()(const Point_3& p) const 01471 { 01472 return p.rep().x(); 01473 } 01474 01475 FT 01476 operator()(const Vector_3& v) const 01477 { 01478 return v.rep().x(); 01479 } 01480 }; 01481 01482 template <typename K> 01483 class Compute_y_2 01484 { 01485 typedef typename K::FT FT; 01486 typedef typename K::Point_2 Point_2; 01487 typedef typename K::Vector_2 Vector_2; 01488 01489 public: 01490 typedef FT result_type; 01491 01492 FT 01493 operator()(const Point_2& p) const 01494 { 01495 return p.rep().y(); 01496 } 01497 01498 FT 01499 operator()(const Vector_2& v) const 01500 { 01501 return v.rep().y(); 01502 } 01503 }; 01504 01505 template <typename K> 01506 class Compute_y_3 01507 { 01508 typedef typename K::FT FT; 01509 typedef typename K::Point_3 Point_3; 01510 typedef typename K::Vector_3 Vector_3; 01511 01512 public: 01513 typedef FT result_type; 01514 01515 FT 01516 operator()(const Point_3& p) const 01517 { 01518 return p.rep().y(); 01519 } 01520 01521 FT 01522 operator()(const Vector_3& v) const 01523 { 01524 return v.rep().y(); 01525 } 01526 }; 01527 01528 template <typename K> 01529 class Compute_z_3 01530 { 01531 typedef typename K::FT FT; 01532 typedef typename K::Point_3 Point_3; 01533 typedef typename K::Vector_3 Vector_3; 01534 01535 public: 01536 typedef FT result_type; 01537 01538 FT 01539 operator()(const Point_3& p) const 01540 { 01541 return p.rep().z(); 01542 } 01543 01544 FT 01545 operator()(const Vector_3& v) const 01546 { 01547 return v.rep().z(); 01548 } 01549 }; 01550 01551 template <typename K> 01552 class Compute_dx_2 : public Has_qrt 01553 { 01554 typedef typename K::RT RT; 01555 typedef typename K::Direction_2 Direction_2; 01556 01557 public: 01558 typedef RT result_type; 01559 01560 const result_type & 01561 operator()(const Direction_2& d) const 01562 { 01563 return d.rep().dx(); 01564 } 01565 }; 01566 01567 template <typename K> 01568 class Compute_dx_3 : public Has_qrt 01569 { 01570 typedef typename K::RT RT; 01571 typedef typename K::Direction_3 Direction_3; 01572 01573 public: 01574 typedef RT result_type; 01575 01576 const result_type & 01577 operator()(const Direction_3& d) const 01578 { 01579 return d.rep().dx(); 01580 } 01581 }; 01582 01583 template <typename K> 01584 class Compute_dy_2 : public Has_qrt 01585 { 01586 typedef typename K::RT RT; 01587 typedef typename K::Direction_2 Direction_2; 01588 01589 public: 01590 typedef RT result_type; 01591 01592 const result_type & 01593 operator()(const Direction_2& d) const 01594 { 01595 return d.rep().dy(); 01596 } 01597 }; 01598 01599 template <typename K> 01600 class Compute_dy_3 : public Has_qrt 01601 { 01602 typedef typename K::RT RT; 01603 typedef typename K::Direction_3 Direction_3; 01604 01605 public: 01606 typedef RT result_type; 01607 01608 const result_type & 01609 operator()(const Direction_3& d) const 01610 { 01611 return d.rep().dy(); 01612 } 01613 }; 01614 01615 template <typename K> 01616 class Compute_dz_3 : public Has_qrt 01617 { 01618 typedef typename K::RT RT; 01619 typedef typename K::Direction_3 Direction_3; 01620 01621 public: 01622 typedef RT result_type; 01623 01624 const result_type & 01625 operator()(const Direction_3& d) const 01626 { 01627 return d.rep().dz(); 01628 } 01629 }; 01630 01631 template <typename K> 01632 class Compute_hx_2 : public Has_qrt 01633 { 01634 typedef typename K::FT FT; 01635 typedef typename K::RT RT; 01636 typedef typename K::Point_2 Point_2; 01637 typedef typename K::Vector_2 Vector_2; 01638 01639 public: 01640 typedef RT result_type; 01641 01642 const result_type & 01643 operator()(const Point_2& p) const 01644 { 01645 return p.rep().hx(); 01646 } 01647 01648 const result_type & 01649 operator()(const Vector_2& v) const 01650 { 01651 return v.rep().hx(); 01652 } 01653 }; 01654 01655 template <typename K> 01656 class Compute_hx_3 : public Has_qrt 01657 { 01658 typedef typename K::FT FT; 01659 typedef typename K::RT RT; 01660 typedef typename K::Point_3 Point_3; 01661 typedef typename K::Vector_3 Vector_3; 01662 01663 public: 01664 typedef RT result_type; 01665 01666 const result_type & 01667 operator()(const Point_3& p) const 01668 { 01669 return p.rep().hx(); 01670 } 01671 01672 const result_type & 01673 operator()(const Vector_3& v) const 01674 { 01675 return v.rep().hx(); 01676 } 01677 }; 01678 01679 template <typename K> 01680 class Compute_hy_2 : public Has_qrt 01681 { 01682 typedef typename K::FT FT; 01683 typedef typename K::RT RT; 01684 typedef typename K::Point_2 Point_2; 01685 typedef typename K::Vector_2 Vector_2; 01686 01687 public: 01688 typedef RT result_type; 01689 01690 const result_type & 01691 operator()(const Point_2& p) const 01692 { 01693 return p.rep().hy(); 01694 } 01695 01696 const result_type & 01697 operator()(const Vector_2& v) const 01698 { 01699 return v.rep().hy(); 01700 } 01701 }; 01702 01703 template <typename K> 01704 class Compute_hy_3 : public Has_qrt 01705 { 01706 typedef typename K::FT FT; 01707 typedef typename K::RT RT; 01708 typedef typename K::Point_3 Point_3; 01709 typedef typename K::Vector_3 Vector_3; 01710 01711 public: 01712 typedef RT result_type; 01713 01714 const result_type & 01715 operator()(const Point_3& p) const 01716 { 01717 return p.rep().hy(); 01718 } 01719 01720 const result_type & 01721 operator()(const Vector_3& v) const 01722 { 01723 return v.rep().hy(); 01724 } 01725 }; 01726 01727 template <typename K> 01728 class Compute_hz_3 : public Has_qrt 01729 { 01730 typedef typename K::FT FT; 01731 typedef typename K::RT RT; 01732 typedef typename K::Point_3 Point_3; 01733 typedef typename K::Vector_3 Vector_3; 01734 01735 public: 01736 typedef RT result_type; 01737 01738 const result_type & 01739 operator()(const Point_3& p) const 01740 { 01741 return p.rep().hz(); 01742 } 01743 01744 const result_type & 01745 operator()(const Vector_3& v) const 01746 { 01747 return v.rep().hz(); 01748 } 01749 }; 01750 01751 template <typename K> 01752 class Compute_hw_2 : public Has_qrt 01753 { 01754 typedef typename K::FT FT; 01755 typedef typename K::RT RT; 01756 typedef typename K::Point_2 Point_2; 01757 typedef typename K::Vector_2 Vector_2; 01758 01759 public: 01760 typedef RT result_type; 01761 01762 const result_type & 01763 operator()(const Point_2& p) const 01764 { 01765 return p.rep().hw(); 01766 } 01767 01768 const result_type & 01769 operator()(const Vector_2& v) const 01770 { 01771 return v.rep().hw(); 01772 } 01773 }; 01774 01775 template <typename K> 01776 class Compute_hw_3 : public Has_qrt 01777 { 01778 typedef typename K::FT FT; 01779 typedef typename K::RT RT; 01780 typedef typename K::Point_3 Point_3; 01781 typedef typename K::Vector_3 Vector_3; 01782 01783 public: 01784 typedef RT result_type; 01785 01786 const result_type & 01787 operator()(const Point_3& p) const 01788 { 01789 return p.rep().hw(); 01790 } 01791 01792 const result_type & 01793 operator()(const Vector_3& v) const 01794 { 01795 return v.rep().hw(); 01796 } 01797 }; 01798 01799 template <typename K> 01800 class Construct_base_vector_3 01801 { 01802 typedef typename K::Vector_3 Vector_3; 01803 typedef typename K::Plane_3 Plane_3; 01804 typedef typename K::RT RT; 01805 typedef typename K::Construct_orthogonal_vector_3 01806 Construct_orthogonal_vector_3; 01807 Construct_orthogonal_vector_3 co; 01808 public: 01809 typedef Vector_3 result_type; 01810 01811 Construct_base_vector_3() {} 01812 Construct_base_vector_3(const Construct_orthogonal_vector_3& co_) 01813 : co(co_) 01814 {} 01815 01816 Vector_3 01817 operator()( const Plane_3& h, int index ) const 01818 { 01819 if (index == 1) { 01820 // point(): 01821 // a() != RT0 : Point_3( -d(), RT0, RT0, a() ); 01822 // b() != RT0 : Point_3( RT0, -d(), RT0, b() ); 01823 // : Point_3( RT0, RT0, -d(), c() ); 01824 // point1(): 01825 // a() != RT0 : Point_3( -b()-d(), a(), RT0, a() ); 01826 // b() != RT0 : Point_3( RT0, -c()-d(), b(), b() ); 01827 // : Point_3( c(), RT0, -a()-d(), c() ); 01828 01829 const RT RT0(0); 01830 if ( h.a() != RT0 ) 01831 { 01832 return Vector_3( -h.b(), h.a(), RT0, h.a() ); 01833 } 01834 if ( h.b() != RT0 ) 01835 { 01836 return Vector_3( RT0, -h.c(), h.b(), h.b() ); 01837 } 01838 CGAL_kernel_assertion ( h.c() != RT(0) ); 01839 return Vector_3( h.c(), RT0, -h.a(), h.c() ); 01840 } else { 01841 Vector_3 a = co(h); 01842 Vector_3 b = this->operator()(h, 1); 01843 return Vector_3(a.hy()*b.hz() - a.hz()*b.hy(), 01844 a.hz()*b.hx() - a.hx()*b.hz(), 01845 a.hx()*b.hy() - a.hy()*b.hx(), 01846 a.hw()*b.hw() ); 01847 } 01848 } 01849 }; 01850 01851 template <typename K> 01852 class Construct_bbox_2 01853 { 01854 typedef typename K::Point_2 Point_2; 01855 typedef typename K::Segment_2 Segment_2; 01856 typedef typename K::Iso_rectangle_2 Iso_rectangle_2; 01857 typedef typename K::Triangle_2 Triangle_2; 01858 typedef typename K::Circle_2 Circle_2; 01859 public: 01860 typedef Bbox_2 result_type; 01861 01862 Bbox_2 01863 operator()( const Point_2& p) const 01864 { 01865 Interval_nt<> ihx = CGAL_NTS to_interval(p.hx()); 01866 Interval_nt<> ihy = CGAL_NTS to_interval(p.hy()); 01867 Interval_nt<> ihw = CGAL_NTS to_interval(p.hw()); 01868 01869 Interval_nt<> ix = ihx/ihw; 01870 Interval_nt<> iy = ihy/ihw; 01871 01872 return Bbox_2(ix.inf(), iy.inf(), ix.sup(), iy.sup()); 01873 } 01874 01875 Bbox_2 01876 operator()( const Segment_2& s) const 01877 { return s.source().bbox() + s.target().bbox(); } 01878 01879 Bbox_2 01880 operator()( const Triangle_2& t) const 01881 { 01882 typename K::Construct_bbox_2 construct_bbox_2; 01883 return construct_bbox_2(t.vertex(0)) 01884 + construct_bbox_2(t.vertex(1)) 01885 + construct_bbox_2(t.vertex(2)); 01886 } 01887 01888 Bbox_2 01889 operator()( const Iso_rectangle_2& r) const 01890 { 01891 typename K::Construct_bbox_2 construct_bbox_2; 01892 return construct_bbox_2((r.min)()) + construct_bbox_2((r.max)()); 01893 } 01894 01895 Bbox_2 01896 operator()( const Circle_2& c) const 01897 { 01898 typename K::Construct_bbox_2 construct_bbox_2; 01899 Bbox_2 b = construct_bbox_2(c.center()); 01900 01901 Interval_nt<> x (b.xmin(), b.xmax()); 01902 Interval_nt<> y (b.ymin(), b.ymax()); 01903 01904 Interval_nt<> sqr = CGAL_NTS to_interval(c.squared_radius()); 01905 Interval_nt<> r = CGAL::sqrt(sqr); 01906 Interval_nt<> minx = x-r; 01907 Interval_nt<> maxx = x+r; 01908 Interval_nt<> miny = y-r; 01909 Interval_nt<> maxy = y+r; 01910 01911 return Bbox_2(minx.inf(), miny.inf(), maxx.sup(), maxy.sup()); } 01912 }; 01913 01914 01915 template <typename K> 01916 class Construct_bbox_3 01917 { 01918 typedef typename K::Point_3 Point_3; 01919 typedef typename K::Segment_3 Segment_3; 01920 typedef typename K::Triangle_3 Triangle_3; 01921 typedef typename K::Tetrahedron_3 Tetrahedron_3; 01922 typedef typename K::Iso_cuboid_3 Iso_cuboid_3; 01923 typedef typename K::Sphere_3 Sphere_3; 01924 public: 01925 typedef Bbox_3 result_type; 01926 01927 Bbox_3 01928 operator()(const Point_3& p) const 01929 { 01930 Interval_nt<> ihx = CGAL_NTS to_interval(p.hx()); 01931 Interval_nt<> ihy = CGAL_NTS to_interval(p.hy()); 01932 Interval_nt<> ihz = CGAL_NTS to_interval(p.hz()); 01933 Interval_nt<> ihw = CGAL_NTS to_interval(p.hw()); 01934 01935 Interval_nt<> ix = ihx/ihw; 01936 Interval_nt<> iy = ihy/ihw; 01937 Interval_nt<> iz = ihz/ihw; 01938 01939 return Bbox_3(ix.inf(), iy.inf(), iz.inf(), 01940 ix.sup(), iy.sup(), iz.sup()); 01941 } 01942 01943 Bbox_3 01944 operator()(const Segment_3& s) const 01945 { return s.source().bbox() + s.target().bbox(); } 01946 01947 Bbox_3 01948 operator()(const Triangle_3& t) const 01949 { 01950 typename K::Construct_bbox_3 construct_bbox; 01951 return construct_bbox(t.vertex(0)) 01952 + construct_bbox(t.vertex(1)) 01953 + construct_bbox(t.vertex(2)); 01954 } 01955 01956 Bbox_3 01957 operator()(const Iso_cuboid_3& r) const 01958 { 01959 typename K::Construct_bbox_3 construct_bbox; 01960 return construct_bbox((r.min)()) + construct_bbox((r.max)()); 01961 } 01962 01963 Bbox_3 01964 operator()(const Tetrahedron_3& t) const 01965 { 01966 typename K::Construct_bbox_3 construct_bbox_3; 01967 return construct_bbox_3(t.vertex(0)) + construct_bbox_3(t.vertex(1)) 01968 + construct_bbox_3(t.vertex(2)) + construct_bbox_3(t.vertex(3)); 01969 } 01970 01971 Bbox_3 01972 operator()(const Sphere_3& s) const 01973 { 01974 Bbox_3 b = s.center().bbox(); 01975 01976 Interval_nt<> x (b.xmin(), b.xmax()); 01977 Interval_nt<> y (b.ymin(), b.ymax()); 01978 Interval_nt<> z (b.zmin(), b.zmax()); 01979 01980 Interval_nt<> sqr = CGAL_NTS to_interval(s.squared_radius()); 01981 Interval_nt<> r = CGAL::sqrt(sqr); 01982 Interval_nt<> minx = x-r; 01983 Interval_nt<> maxx = x+r; 01984 Interval_nt<> miny = y-r; 01985 Interval_nt<> maxy = y+r; 01986 Interval_nt<> minz = z-r; 01987 Interval_nt<> maxz = z+r; 01988 01989 return Bbox_3(minx.inf(), miny.inf(), minz.inf(), 01990 maxx.sup(), maxy.sup(), maxz.sup()); 01991 } 01992 }; 01993 01994 01995 template <typename K> 01996 class Construct_bisector_2 01997 { 01998 typedef typename K::RT RT; 01999 typedef typename K::FT FT; 02000 typedef typename K::Point_2 Point_2; 02001 typedef typename K::Line_2 Line_2; 02002 public: 02003 typedef Line_2 result_type; 02004 02005 Line_2 02006 operator()(const Point_2& p, const Point_2& q) const 02007 { 02008 // Bisector equation is based on equation 02009 // ( X - p.x())^2 + (Y - p.y())^2 == ( X - q.x())^2 + (Y - q.y()) 02010 // and x() = hx()/hw() ... 02011 02012 const RT &phx = p.hx(); 02013 const RT &phy = p.hy(); 02014 const RT &phw = p.hw(); 02015 const RT &qhx = q.hx(); 02016 const RT &qhy = q.hy(); 02017 const RT &qhw = q.hw(); 02018 02019 RT a = RT(2) * ( phx*phw*qhw*qhw - qhx*qhw*phw*phw ); 02020 RT b = RT(2) * ( phy*phw*qhw*qhw - qhy*qhw*phw*phw ); 02021 RT c = qhx*qhx*phw*phw + qhy*qhy*phw*phw 02022 - phx*phx*qhw*qhw - phy*phy*qhw*qhw; 02023 02024 return Line_2( a, b, c ); 02025 } 02026 02027 Line_2 02028 operator()(const Line_2& p, const Line_2& q) const 02029 { 02030 RT a, b, c; 02031 bisector_of_linesC2(p.a(), p.b(), p.c(), 02032 q.a(), q.b(), q.c(), 02033 a, b, c); 02034 return Line_2(a, b, c); 02035 } 02036 }; 02037 02038 template <typename K> 02039 class Construct_bisector_3 02040 { 02041 typedef typename K::RT RT; 02042 typedef typename K::FT FT; 02043 typedef typename K::Point_3 Point_3; 02044 typedef typename K::Plane_3 Plane_3; 02045 public: 02046 typedef Plane_3 result_type; 02047 02048 Plane_3 02049 operator()(const Point_3& p, const Point_3& q) const 02050 { 02051 // Bisector equation is based on equation 02052 // ( X - p.x())^2 + (Y - p.y())^2 == ( X - q.x())^2 + (Y - q.y()) 02053 // and x() = hx()/hw() ... 02054 02055 const RT& phx = p.hx(); 02056 const RT& phy = p.hy(); 02057 const RT& phz = p.hz(); 02058 const RT& phw = p.hw(); 02059 const RT& qhx = q.hx(); 02060 const RT& qhy = q.hy(); 02061 const RT& qhz = q.hz(); 02062 const RT& qhw = q.hw(); 02063 02064 RT a = RT(2) * ( phx*phw*qhw*qhw - qhx*qhw*phw*phw ); 02065 RT b = RT(2) * ( phy*phw*qhw*qhw - qhy*qhw*phw*phw ); 02066 RT c = RT(2) * ( phz*phw*qhw*qhw - qhz*qhw*phw*phw ); 02067 RT d = qhx*qhx*phw*phw + qhy*qhy*phw*phw + qhz*qhz*phw*phw 02068 - phx*phx*qhw*qhw - phy*phy*qhw*qhw - phz*phz*qhw*qhw; 02069 02070 return Plane_3( a, b, c, d ); 02071 } 02072 02073 Plane_3 02074 operator()(const Plane_3& p, const Plane_3& q) const 02075 { 02076 RT a, b, c, d; 02077 bisector_of_planesC3(p.a(), p.b(), p.c(), p.d(), 02078 q.a(), q.b(), q.c(), q.d(), 02079 a, b, c, d); 02080 return Plane_3(a, b, c, d); 02081 } 02082 }; 02083 02084 02085 template <typename K> 02086 class Construct_centroid_2 02087 { 02088 typedef typename K::FT FT; 02089 typedef typename K::Point_2 Point_2; 02090 typedef typename K::Triangle_2 Triangle_2; 02091 public: 02092 typedef Point_2 result_type; 02093 02094 Point_2 02095 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 02096 { 02097 typedef typename K::RT RT; 02098 const RT phw(p.hw()); 02099 const RT qhw(q.hw()); 02100 const RT rhw(r.hw()); 02101 RT hx(p.hx()*qhw*rhw + q.hx()*phw*rhw + r.hx()*phw*qhw); 02102 RT hy(p.hy()*qhw*rhw + q.hy()*phw*rhw + r.hy()*phw*qhw); 02103 RT hw( phw*qhw*rhw * 3); 02104 return Point_2(hx, hy, hw); 02105 } 02106 02107 Point_2 02108 operator()(const Triangle_2& t) const 02109 { 02110 return this->operator()(t.vertex(0), t.vertex(1), t.vertex(2)); 02111 } 02112 02113 Point_2 02114 operator()(const Point_2& p, const Point_2& q, 02115 const Point_2& r, const Point_2& s) const 02116 { 02117 typedef typename K::RT RT; 02118 const RT phw(p.hw()); 02119 const RT qhw(q.hw()); 02120 const RT rhw(r.hw()); 02121 const RT shw(s.hw()); 02122 RT hx(p.hx()*qhw*rhw*shw + q.hx()*phw*rhw*shw + r.hx()*phw*qhw*shw 02123 + s.hx()*phw*qhw*rhw); 02124 RT hy(p.hy()*qhw*rhw*shw + q.hy()*phw*rhw*shw + r.hy()*phw*qhw*shw 02125 + s.hy()*phw*qhw*rhw); 02126 RT hw( phw*qhw*rhw*shw * 4); 02127 return Point_2(hx, hy, hw); 02128 } 02129 }; 02130 02131 template <typename K> 02132 class Construct_centroid_3 02133 { 02134 typedef typename K::RT RT; 02135 typedef typename K::Point_3 Point_3; 02136 typedef typename K::Triangle_3 Triangle_3; 02137 typedef typename K::Tetrahedron_3 Tetrahedron_3; 02138 public: 02139 typedef Point_3 result_type; 02140 02141 Point_3 02142 operator()(const Point_3& p, const Point_3& q, const Point_3& r) const 02143 { 02144 const RT& phw = p.hw(); 02145 const RT& qhw = q.hw(); 02146 const RT& rhw = r.hw(); 02147 RT hx(p.hx()*qhw*rhw + q.hx()*phw*rhw + r.hx()*phw*qhw); 02148 RT hy(p.hy()*qhw*rhw + q.hy()*phw*rhw + r.hy()*phw*qhw); 02149 RT hz(p.hz()*qhw*rhw + q.hz()*phw*rhw + r.hz()*phw*qhw); 02150 RT hw( phw*qhw*rhw * RT(3)); 02151 return Point_3(hx, hy, hz, hw); 02152 } 02153 02154 Point_3 02155 operator()(const Point_3& p, const Point_3& q, 02156 const Point_3& r, const Point_3& s) const 02157 { 02158 const RT& phw = p.hw(); 02159 const RT& qhw = q.hw(); 02160 const RT& rhw = r.hw(); 02161 const RT& shw = s.hw(); 02162 RT hx(p.hx()*qhw*rhw*shw + q.hx()*phw*rhw*shw + r.hx()*phw*qhw*shw 02163 + s.hx()*phw*qhw*rhw); 02164 RT hy(p.hy()*qhw*rhw*shw + q.hy()*phw*rhw*shw + r.hy()*phw*qhw*shw 02165 + s.hy()*phw*qhw*rhw); 02166 RT hz(p.hz()*qhw*rhw*shw + q.hz()*phw*rhw*shw + r.hz()*phw*qhw*shw 02167 + s.hz()*phw*qhw*rhw); 02168 RT hw( phw*qhw*rhw*shw * RT(4)); 02169 return Point_3(hx, hy, hz, hw); 02170 } 02171 02172 Point_3 02173 operator()(const Triangle_3& t) const 02174 { 02175 return this->operator()(t.vertex(0), t.vertex(1), t.vertex(2)); 02176 } 02177 02178 Point_3 02179 operator()(const Tetrahedron_3& t) const 02180 { 02181 return this->operator()(t.vertex(0), t.vertex(1), 02182 t.vertex(2), t.vertex(3)); 02183 } 02184 }; 02185 02186 template <typename K> 02187 class Construct_circumcenter_2 02188 { 02189 typedef typename K::FT FT; 02190 typedef typename K::Point_2 Point_2; 02191 typedef typename K::Triangle_2 Triangle_2; 02192 public: 02193 typedef Point_2 result_type; 02194 02195 Point_2 02196 operator()(const Point_2& p, const Point_2& q) const 02197 { 02198 typename K::Construct_midpoint_2 construct_midpoint_2; 02199 return construct_midpoint_2(p, q); 02200 } 02201 02202 Point_2 02203 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 02204 { 02205 typedef typename K::RT RT; 02206 const RT & phx = p.hx(); 02207 const RT & phy = p.hy(); 02208 const RT & phw = p.hw(); 02209 const RT & qhx = q.hx(); 02210 const RT & qhy = q.hy(); 02211 const RT & qhw = q.hw(); 02212 const RT & rhx = r.hx(); 02213 const RT & rhy = r.hy(); 02214 const RT & rhw = r.hw(); 02215 02216 #ifdef CGAL_EXPANDED_CIRCUMCENTER_COMPUTATION 02217 RT vvx = 02218 ( qhy*qhw*phw*phw - phy*phw*qhw*qhw ) 02219 *( phx*phx*rhw*rhw + phy*phy*rhw*rhw - 02220 rhx*rhx*phw*phw - rhy*rhy*phw*phw ) 02221 - ( rhy*rhw*phw*phw - phy*phw*rhw*rhw ) 02222 *( phx*phx*qhw*qhw + phy*phy*qhw*qhw - 02223 qhx*qhx*phw*phw - qhy*qhy*phw*phw ); 02224 02225 RT vvy = 02226 - ( qhx*qhw*phw*phw - phx*phw*qhw*qhw ) 02227 *( phx*phx*rhw*rhw + phy*phy*rhw*rhw - 02228 rhx*rhx*phw*phw - rhy*rhy*phw*phw ) 02229 + ( rhx*rhw*phw*phw - phx*phw*rhw*rhw ) 02230 *( phx*phx*qhw*qhw + phy*phy*qhw*qhw - 02231 qhx*qhx*phw*phw - qhy*qhy*phw*phw ); 02232 02233 RT vvw = RT(2) * 02234 ( ( qhx*qhw*phw*phw - phx*phw*qhw*qhw ) 02235 *( rhy*rhw*phw*phw - phy*phw*rhw*rhw ) 02236 - ( rhx*rhw*phw*phw - phx*phw*rhw*rhw ) 02237 *( qhy*qhw*phw*phw - phy*phw*qhw*qhw ) ); 02238 #endif // CGAL_EXPANDED_CIRCUMCENTER_COMPUTATION 02239 02240 RT qy_py = ( qhy*qhw*phw*phw - phy*phw*qhw*qhw ); 02241 RT qx_px = ( qhx*qhw*phw*phw - phx*phw*qhw*qhw ); 02242 RT rx_px = ( rhx*rhw*phw*phw - phx*phw*rhw*rhw ); 02243 RT ry_py = ( rhy*rhw*phw*phw - phy*phw*rhw*rhw ); 02244 02245 RT px2_py2_rx2_ry_2 = 02246 phx*phx*rhw*rhw + phy*phy*rhw*rhw - 02247 rhx*rhx*phw*phw - rhy*rhy*phw*phw ; 02248 RT px2_py2_qx2_qy_2 = 02249 phx*phx*qhw*qhw + phy*phy*qhw*qhw - 02250 qhx*qhx*phw*phw - qhy*qhy*phw*phw ; 02251 02252 RT vvx = qy_py * px2_py2_rx2_ry_2 - ry_py * px2_py2_qx2_qy_2; 02253 RT vvy = rx_px * px2_py2_qx2_qy_2 - qx_px * px2_py2_rx2_ry_2; 02254 RT vvw = RT(2) * ( qx_px * ry_py - rx_px * qy_py ); 02255 02256 return Point_2( vvx, vvy, vvw ); 02257 } 02258 02259 Point_2 02260 operator()(const Triangle_2& t) const 02261 { 02262 return this->operator()(t.vertex(0), t.vertex(1), t.vertex(2)); 02263 } 02264 }; 02265 02266 template <typename K> 02267 class Construct_circumcenter_3 02268 { 02269 typedef typename K::FT FT; 02270 typedef typename K::Point_3 Point_3; 02271 typedef typename K::Triangle_3 Triangle_3; 02272 typedef typename K::Tetrahedron_3 Tetrahedron_3; 02273 typedef typename K::Plane_3 Plane_3; 02274 public: 02275 typedef Point_3 result_type; 02276 02277 Point_3 02278 operator()(const Point_3& p, const Point_3& q) const 02279 { 02280 typename K::Construct_midpoint_3 construct_midpoint_3; 02281 return construct_midpoint_3(p, q); 02282 } 02283 02284 Point_3 02285 operator()(const Point_3& p, const Point_3& q, const Point_3& r) const 02286 { 02287 return gp_linear_intersection( Plane_3(p,q,r), 02288 bisector(p,q), 02289 bisector(p,r)); 02290 } // FIXME 02291 02292 Point_3 02293 operator()(const Triangle_3& t) const 02294 { 02295 return this->operator()(t.vertex(0), t.vertex(1), t.vertex(2)); 02296 } 02297 02298 Point_3 02299 operator()(const Point_3& p, const Point_3& q, 02300 const Point_3& r, const Point_3& s) const 02301 { 02302 typedef typename K::RT RT; 02303 02304 RT phw( p.hw() ); 02305 RT qhw( q.hw() ); 02306 RT rhw( r.hw() ); 02307 RT shw( s.hw() ); 02308 02309 RT phx( p.hx() ); 02310 RT phy( p.hy() ); 02311 RT phz( p.hz() ); 02312 RT qhx( q.hx() ); 02313 RT qhy( q.hy() ); 02314 RT qhz( q.hz() ); 02315 RT rhx( r.hx() ); 02316 RT rhy( r.hy() ); 02317 RT rhz( r.hz() ); 02318 RT shx( s.hx() ); 02319 RT shy( s.hy() ); 02320 RT shz( s.hz() ); 02321 02322 RT pssq( phx*phx + phy*phy + phz*phz ); 02323 RT qssq( qhx*qhx + qhy*qhy + qhz*qhz ); 02324 RT rssq( rhx*rhx + rhy*rhy + rhz*rhz ); 02325 RT sssq( shx*shx + shy*shy + shz*shz ); 02326 02327 phx *= phw; 02328 phy *= phw; 02329 phz *= phw; 02330 phw *= phw; 02331 qhx *= qhw; 02332 qhy *= qhw; 02333 qhz *= qhw; 02334 qhw *= qhw; 02335 rhx *= rhw; 02336 rhy *= rhw; 02337 rhz *= rhw; 02338 rhw *= rhw; 02339 shx *= shw; 02340 shy *= shw; 02341 shz *= shw; 02342 shw *= shw; 02343 02344 RT chx = determinant(phy, phz, pssq, phw, 02345 qhy, qhz, qssq, qhw, 02346 rhy, rhz, rssq, rhw, 02347 shy, shz, sssq, shw ); 02348 RT chy = determinant(phx, phz, pssq, phw, 02349 qhx, qhz, qssq, qhw, 02350 rhx, rhz, rssq, rhw, 02351 shx, shz, sssq, shw ); 02352 RT chz = determinant(phx, phy, pssq, phw, 02353 qhx, qhy, qssq, qhw, 02354 rhx, rhy, rssq, rhw, 02355 shx, shy, sssq, shw ); 02356 RT chw = determinant(phx, phy, phz, phw, 02357 qhx, qhy, qhz, qhw, 02358 rhx, rhy, rhz, rhw, 02359 shx, shy, shz, shw ); 02360 02361 return Point_3( chx, -chy, chz, RT(2)*chw); 02362 } 02363 02364 Point_3 02365 operator()(const Tetrahedron_3& t) const 02366 { 02367 return this->operator()(t.vertex(0), t.vertex(1), 02368 t.vertex(2), t.vertex(3)); 02369 } 02370 }; 02371 02372 template <typename K> 02373 class Construct_cross_product_vector_3 02374 { 02375 typedef typename K::Vector_3 Vector_3; 02376 public: 02377 typedef Vector_3 result_type; 02378 02379 Vector_3 02380 operator()(const Vector_3& a, const Vector_3& b) const 02381 { 02382 return Vector_3(a.hy()*b.hz() - a.hz()*b.hy(), 02383 a.hz()*b.hx() - a.hx()*b.hz(), 02384 a.hx()*b.hy() - a.hy()*b.hx(), 02385 a.hw()*b.hw() ); 02386 } 02387 }; 02388 02389 template <typename K> 02390 class Construct_difference_of_vectors_2 02391 { 02392 typedef typename K::Vector_2 Vector_2; 02393 public: 02394 typedef Vector_2 result_type; 02395 02396 Vector_2 02397 operator()(const Vector_2& v, const Vector_2& w) const 02398 { 02399 return Vector_2( v.hx()*w.hw() - w.hx()*v.hw(), 02400 v.hy()*w.hw() - w.hy()*v.hw(), 02401 v.hw()*w.hw() ); 02402 } 02403 }; 02404 02405 template <typename K> 02406 class Construct_difference_of_vectors_3 02407 { 02408 typedef typename K::Vector_3 Vector_3; 02409 public: 02410 typedef Vector_3 result_type; 02411 02412 Vector_3 02413 operator()(const Vector_3& v, const Vector_3& w) const 02414 { 02415 return Vector_3( v.hx()*w.hw() - w.hx()*v.hw(), 02416 v.hy()*w.hw() - w.hy()*v.hw(), 02417 v.hz()*w.hw() - w.hz()*v.hw(), 02418 v.hw()*w.hw() ); 02419 } 02420 }; 02421 02422 02423 template <typename K> 02424 class Construct_direction_2 02425 { 02426 typedef typename K::Direction_2 Direction_2; 02427 typedef typename Direction_2::Rep Rep; 02428 typedef typename K::Point_2 Point_2; 02429 typedef typename K::Vector_2 Vector_2; 02430 typedef typename K::Line_2 Line_2; 02431 typedef typename K::Ray_2 Ray_2; 02432 typedef typename K::Segment_2 Segment_2; 02433 typedef typename K::RT RT; 02434 02435 public: 02436 typedef Direction_2 result_type; 02437 02438 Rep // Direction_2 02439 operator()(Return_base_tag, const RT& x, const RT& y) const 02440 { return Rep(x, y); } 02441 02442 Rep // Direction_2 02443 operator()(Return_base_tag, const Vector_2& v) const 02444 { return Rep(v.hx(),v.hy()); } 02445 02446 Rep // Direction_2 02447 operator()(Return_base_tag, const Line_2& l) const 02448 { return Rep(l.b(), -l.a()); } 02449 02450 Rep // Direction_2 02451 operator()(Return_base_tag, const Ray_2& r) const 02452 { return this->operator()(Return_base_tag(), r.source(), r.second_point()); } 02453 02454 Rep // Direction_2 02455 operator()(Return_base_tag, const Segment_2& s) const 02456 { return this->operator()(Return_base_tag(), s.source(), s.target()); } 02457 02458 Rep // Direction_2 02459 operator()(Return_base_tag, const Point_2& q, const Point_2& p) const 02460 { 02461 return Rep( p.hx()*q.hw() - q.hx()*p.hw(), 02462 p.hy()*q.hw() - q.hy()*p.hw() ); 02463 } 02464 02465 02466 Direction_2 02467 operator()(const RT& x, const RT& y) const 02468 { return this->operator()(Return_base_tag(), x, y); } 02469 02470 Direction_2 02471 operator()(const Vector_2& v) const 02472 { return this->operator()(Return_base_tag(), v); } 02473 02474 Direction_2 02475 operator()(const Line_2& l) const 02476 { return this->operator()(Return_base_tag(), l); } 02477 02478 Direction_2 02479 operator()(const Ray_2& r) const 02480 { return this->operator()(Return_base_tag(), r); } 02481 02482 Direction_2 02483 operator()(const Segment_2& s) const 02484 { return this->operator()(Return_base_tag(), s); } 02485 02486 Direction_2 02487 operator()(const Point_2& q, const Point_2& p) const 02488 { 02489 return this->operator()(Return_base_tag(), p, q); 02490 } 02491 }; 02492 02493 template <typename K> 02494 class Construct_direction_3 02495 { 02496 typedef typename K::Direction_3 Direction_3; 02497 typedef typename K::Vector_3 Vector_3; 02498 typedef typename K::Line_3 Line_3; 02499 typedef typename K::Ray_3 Ray_3; 02500 typedef typename K::Segment_3 Segment_3; 02501 typedef typename K::RT RT; 02502 typedef typename Direction_3::Rep Rep; 02503 public: 02504 typedef Direction_3 result_type; 02505 02506 Rep // Direction_3 02507 operator()(Return_base_tag, const RT& x, const RT& y, const RT& z) const 02508 { return Rep(x, y, z); } 02509 02510 Rep // Direction_3 02511 operator()(Return_base_tag, const Vector_3& v) const 02512 { return Rep(v.hx(), v.hy(), v.hz()); } 02513 02514 Rep // Direction_3 02515 operator()(Return_base_tag, const Line_3& l) const 02516 { return Rep(l); } 02517 02518 Rep // Direction_3 02519 operator()(Return_base_tag, const Ray_3& r) const 02520 { return Rep(r); } 02521 02522 Rep // Direction_3 02523 operator()(Return_base_tag, const Segment_3& s) const 02524 { return Rep(s); } 02525 02526 02527 Direction_3 02528 operator()(const RT& x, const RT& y, const RT& z) const 02529 { return this->operator()(Return_base_tag(), x, y, z); } 02530 02531 Direction_3 02532 operator()(const Vector_3& v) const 02533 { return this->operator()(Return_base_tag(), v); } 02534 02535 Direction_3 02536 operator()(const Line_3& l) const 02537 { return this->operator()(Return_base_tag(), l); } 02538 02539 Direction_3 02540 operator()(const Ray_3& r) const 02541 { return this->operator()(Return_base_tag(), r); } 02542 02543 Direction_3 02544 operator()(const Segment_3& s) const 02545 { return this->operator()(Return_base_tag(), s); } 02546 }; 02547 02548 02549 template <typename K> 02550 class Construct_sum_of_vectors_2 02551 { 02552 typedef typename K::Vector_2 Vector_2; 02553 public: 02554 typedef Vector_2 result_type; 02555 02556 Vector_2 02557 operator()(const Vector_2& v, const Vector_2& w) const 02558 { 02559 return Vector_2(v.hx()*w.hw() + w.hx()*v.hw(), 02560 v.hy()*w.hw() + w.hy()*v.hw(), 02561 v.hw()*w.hw()); 02562 } 02563 }; 02564 02565 template <typename K> 02566 class Construct_sum_of_vectors_3 02567 { 02568 typedef typename K::Vector_3 Vector_3; 02569 public: 02570 typedef Vector_3 result_type; 02571 02572 Vector_3 02573 operator()(const Vector_3& v, const Vector_3& w) const 02574 { 02575 return Vector_3(v.hx()*w.hw() + w.hx()*v.hw(), 02576 v.hy()*w.hw() + w.hy()*v.hw(), 02577 v.hz()*w.hw() + w.hz()*v.hw(), 02578 v.hw()*w.hw()); 02579 } 02580 }; 02581 02582 template <typename K> 02583 class Construct_divided_vector_2 02584 { 02585 typedef typename K::FT FT; 02586 typedef typename K::RT RT; 02587 typedef typename K::Vector_2 Vector_2; 02588 public: 02589 typedef Vector_2 result_type; 02590 02591 Vector_2 02592 operator()(const Vector_2& v, const FT& f ) const 02593 { 02594 return Vector_2( v.hx()*f.denominator(), v.hy()*f.denominator(), 02595 v.hw()*f.numerator() ); 02596 } 02597 02598 Vector_2 02599 operator()(const Vector_2& v, const RT& f ) const 02600 { 02601 return Vector_2( v.hx(), v.hy(), v.hw()*f ); 02602 } 02603 }; 02604 02605 template <typename K> 02606 class Construct_divided_vector_3 02607 { 02608 typedef typename K::FT FT; 02609 typedef typename K::RT RT; 02610 typedef typename K::Vector_3 Vector_3; 02611 public: 02612 typedef Vector_3 result_type; 02613 02614 Vector_3 02615 operator()(const Vector_3& v, const FT& f ) const 02616 { 02617 return Vector_3( v.hx()*f.denominator(), v.hy()*f.denominator(), 02618 v.hz()*f.denominator(), v.hw()*f.numerator() ); 02619 } 02620 02621 Vector_3 02622 operator()(const Vector_3& v, const RT& f ) const 02623 { 02624 return Vector_3( v.hx(), v.hy(), v.hz(), v.hw()*f ); 02625 } 02626 }; 02627 02628 02629 template <typename K> 02630 class Construct_iso_rectangle_2 02631 { 02632 typedef typename K::RT RT; 02633 typedef typename K::FT FT; 02634 typedef typename K::Point_2 Point_2; 02635 typedef typename K::Iso_rectangle_2 Iso_rectangle_2; 02636 typedef typename Iso_rectangle_2::Rep Rep; 02637 02638 public: 02639 typedef Iso_rectangle_2 result_type; 02640 02641 Rep // Iso_rectangle_2 02642 operator()(Return_base_tag, const Point_2& p, const Point_2& q, int) const 02643 { 02644 // I have to remove the assertions, because of Cartesian_converter. 02645 // CGAL_kernel_assertion(p.x()<=q.x()); 02646 // CGAL_kernel_assertion(p.y()<=q.y()); 02647 return Rep(p, q, 0); 02648 } 02649 02650 Rep // Iso_rectangle_2 02651 operator()(Return_base_tag, const Point_2& p, const Point_2& q) const 02652 { 02653 bool px_g_qx = ( p.hx()*q.hw() > q.hx()*p.hw() ); 02654 bool py_g_qy = ( p.hy()*q.hw() > q.hy()*p.hw() ); 02655 02656 if ( px_g_qx || py_g_qy) 02657 { 02658 if ( px_g_qx && py_g_qy ) 02659 { 02660 return Rep(q, p, 0); 02661 } 02662 else 02663 { 02664 if ( px_g_qx ) 02665 { 02666 return Rep( 02667 Point_2(q.hx()*p.hw(), p.hy()*q.hw(), q.hw()*p.hw() ), 02668 Point_2(p.hx()*q.hw(), q.hy()*p.hw(), q.hw()*p.hw() ), 0); 02669 } 02670 if ( py_g_qy ) 02671 { 02672 return Rep( 02673 Point_2(p.hx()*q.hw(), q.hy()*p.hw(), q.hw()*p.hw() ), 02674 Point_2(q.hx()*p.hw(), p.hy()*q.hw(), q.hw()*p.hw() ), 0); 02675 } 02676 } 02677 } 02678 return Rep(p, q, 0); 02679 } 02680 02681 Rep // Iso_rectangle_2 02682 operator()(Return_base_tag, const Point_2 &left, const Point_2 &right, 02683 const Point_2 &bottom, const Point_2 &top) const 02684 { 02685 CGAL_kernel_assertion_code(typename K::Less_x_2 less_x;) 02686 CGAL_kernel_assertion_code(typename K::Less_y_2 less_y;) 02687 CGAL_kernel_assertion(!less_x(right, left)); 02688 CGAL_kernel_assertion(!less_y(top, bottom)); 02689 return Rep(Point_2(left.hx() * bottom.hw(), 02690 bottom.hy() * left.hw(), 02691 left.hw() * bottom.hw()), 02692 Point_2(right.hx() * top.hw(), 02693 top.hy() * right.hw(), 02694 right.hw() * top.hw()), 0); 02695 } 02696 02697 Rep // Iso_rectangle_2 02698 operator()(Return_base_tag, const RT& min_hx, const RT& min_hy, 02699 const RT& max_hx, const RT& max_hy) const 02700 { 02701 CGAL_kernel_precondition(min_hx <= max_hx); 02702 CGAL_kernel_precondition(min_hy <= max_hy); 02703 return Rep(Point_2(min_hx, min_hy), 02704 Point_2(max_hx, max_hy), 0); 02705 } 02706 02707 Rep // Iso_rectangle_2 02708 operator()(Return_base_tag, const RT& min_hx, const RT& min_hy, 02709 const RT& max_hx, const RT& max_hy, const RT& hw) const 02710 { 02711 return Rep(Point_2(min_hx, min_hy, hw), 02712 Point_2(max_hx, max_hy, hw), 0); 02713 } 02714 02715 02716 Iso_rectangle_2 02717 operator()(const Point_2& p, const Point_2& q, int i) const 02718 { 02719 return this->operator()(Return_base_tag(), p, q, i); 02720 } 02721 02722 Iso_rectangle_2 02723 operator()(const Point_2& p, const Point_2& q) const 02724 { 02725 return this->operator()(Return_base_tag(), p, q); 02726 } 02727 02728 Iso_rectangle_2 02729 operator()(const Point_2 &left, const Point_2 &right, 02730 const Point_2 &bottom, const Point_2 &top) const 02731 { 02732 return this->operator()(Return_base_tag(), left, right, bottom, top); 02733 } 02734 02735 Iso_rectangle_2 02736 operator()(const RT& min_hx, const RT& min_hy, 02737 const RT& max_hx, const RT& max_hy) const 02738 { 02739 return this->operator()(Return_base_tag(), min_hx, min_hy, max_hx, max_hy); 02740 } 02741 02742 Iso_rectangle_2 02743 operator()(const RT& min_hx, const RT& min_hy, 02744 const RT& max_hx, const RT& max_hy, const RT& hw) const 02745 { 02746 return this->operator()(Return_base_tag(), min_hx, min_hy, max_hx, max_hy, hw); 02747 } 02748 }; 02749 02750 02751 template <typename K> 02752 class Construct_lifted_point_3 02753 { 02754 typedef typename K::RT RT; 02755 typedef typename K::Point_2 Point_2; 02756 typedef typename K::Point_3 Point_3; 02757 typedef typename K::Plane_3 Plane_3; 02758 public: 02759 typedef Point_3 result_type; 02760 02761 Point_3 02762 operator()(const Plane_3& h, const Point_2& p) const 02763 { 02764 Point_3 hp( p.hx(), p.hy(), RT(0.0), p.hw()); 02765 return hp.transform( h.transform_to_2d().inverse() ); 02766 } 02767 }; 02768 02769 template <typename K> 02770 class Construct_line_2 02771 { 02772 typedef typename K::RT RT; 02773 typedef typename K::FT FT; 02774 typedef typename K::Point_2 Point_2; 02775 typedef typename K::Vector_2 Vector_2; 02776 typedef typename K::Direction_2 Direction_2; 02777 typedef typename K::Segment_2 Segment_2; 02778 typedef typename K::Ray_2 Ray_2; 02779 typedef typename K::Line_2 Line_2; 02780 typedef typename Line_2::Rep Rep; 02781 typedef typename K::Construct_point_on_2 Construct_point_on_2; 02782 Construct_point_on_2 cp; 02783 public: 02784 typedef Line_2 result_type; 02785 02786 Construct_line_2() {} 02787 Construct_line_2(const Construct_point_on_2& cp_) : cp(cp_) {} 02788 02789 Rep // Line_2 02790 operator()(Return_base_tag, const RT& a, const RT& b, const RT& c) const 02791 { return Rep(a, b, c); } 02792 02793 Rep // Line_2 02794 operator()(Return_base_tag, const Point_2& p, const Point_2& q) const 02795 { 02796 return Rep( 02797 // a() * X + b() * Y + c() * W() == 0 02798 // | X Y W | 02799 // | p.hx() p.hy() p.hw() | 02800 // | q.hx() q.hy() q.hw() | 02801 02802 p.hy()*q.hw() - p.hw()*q.hy(), 02803 p.hw()*q.hx() - p.hx()*q.hw(), 02804 p.hx()*q.hy() - p.hy()*q.hx() ); 02805 } 02806 02807 Rep // Line_2 02808 operator()(Return_base_tag, const Point_2& p, const Vector_2& v) const 02809 { 02810 Point_2 q = p + v; 02811 return Rep( p.hy()*q.hw() - p.hw()*q.hy(), 02812 p.hw()*q.hx() - p.hx()*q.hw(), 02813 p.hx()*q.hy() - p.hy()*q.hx() ); 02814 } 02815 02816 Rep // Line_2 02817 operator()(Return_base_tag, const Point_2& p, const Direction_2& d) const 02818 { 02819 Point_2 q = p + d.to_vector(); 02820 return Rep( p.hy()*q.hw() - p.hw()*q.hy(), 02821 p.hw()*q.hx() - p.hx()*q.hw(), 02822 p.hx()*q.hy() - p.hy()*q.hx() ); 02823 } 02824 02825 Rep // Line_2 02826 operator()(Return_base_tag, const Segment_2& s) const 02827 { return this->operator()(Return_base_tag(), cp(s, 0), cp(s, 1)); } 02828 02829 Rep // Line_2 02830 operator()(Return_base_tag, const Ray_2& r) const 02831 { return this->operator()(Return_base_tag(), cp(r, 0), cp(r, 1)); } 02832 02833 02834 Line_2 02835 operator()(const RT& a, const RT& b, const RT& c) const 02836 { return this->operator()(Return_base_tag(), a, b, c); } 02837 02838 Line_2 02839 operator()(const Point_2& p, const Point_2& q) const 02840 { return this->operator()(Return_base_tag(), p, q); } 02841 02842 Line_2 02843 operator()(const Point_2& p, const Vector_2& v) const 02844 { return this->operator()(Return_base_tag(), p, v); } 02845 02846 Line_2 02847 operator()(const Point_2& p, const Direction_2& d) const 02848 { return this->operator()(Return_base_tag(), p, d); } 02849 02850 Line_2 02851 operator()(const Segment_2& s) const 02852 { return this->operator()(Return_base_tag(), s); } 02853 02854 Line_2 02855 operator()(const Ray_2& r) const 02856 { return this->operator()(Return_base_tag(), r); } 02857 }; 02858 02859 template <typename K> 02860 class Construct_midpoint_2 02861 { 02862 typedef typename K::FT FT; 02863 typedef typename K::Point_2 Point_2; 02864 public: 02865 typedef Point_2 result_type; 02866 02867 Point_2 02868 operator()(const Point_2& p, const Point_2& q) const 02869 { 02870 typedef typename K::RT RT; 02871 const RT& phw = p.hw(); 02872 const RT& qhw = q.hw(); 02873 return Point_2( p.hx()*qhw + q.hx()*phw, 02874 p.hy()*qhw + q.hy()*phw, 02875 phw * qhw * RT( 2)); 02876 } 02877 }; 02878 02879 template <typename K> 02880 class Construct_midpoint_3 02881 { 02882 typedef typename K::FT FT; 02883 typedef typename K::Point_3 Point_3; 02884 public: 02885 typedef Point_3 result_type; 02886 02887 Point_3 02888 operator()(const Point_3& p, const Point_3& q) const 02889 { 02890 typedef typename K::RT RT; 02891 RT phw = p.hw(); 02892 RT qhw = q.hw(); 02893 return Point_3( p.hx()*qhw + q.hx()*phw, 02894 p.hy()*qhw + q.hy()*phw, 02895 p.hz()*qhw + q.hz()*phw, 02896 RT(2) * phw * qhw ); 02897 } 02898 }; 02899 02900 // TODO ... 02901 template <typename K> 02902 class Construct_opposite_vector_2 02903 { 02904 typedef typename K::Vector_2 Vector_2; 02905 public: 02906 typedef Vector_2 result_type; 02907 02908 Vector_2 02909 operator()( const Vector_2& v) const 02910 { return Vector_2(-v.hx(), -v.hy(), v.hw()); } 02911 }; 02912 02913 template <typename K> 02914 class Construct_opposite_vector_3 02915 { 02916 typedef typename K::Vector_3 Vector_3; 02917 public: 02918 typedef Vector_3 result_type; 02919 02920 Vector_3 02921 operator()( const Vector_3& v) const 02922 { return Vector_3(-v.hx(), -v.hy(), -v.hz(), v.hw()); } 02923 }; 02924 02925 template <typename K> 02926 class Construct_orthogonal_vector_3 02927 { 02928 typedef typename K::Point_3 Point_3; 02929 typedef typename K::Vector_3 Vector_3; 02930 typedef typename K::Plane_3 Plane_3; 02931 public: 02932 typedef Vector_3 result_type; 02933 02934 Vector_3 02935 operator()( const Plane_3& p ) const 02936 { return p.rep().orthogonal_vector(); } 02937 02938 Vector_3 02939 operator()( const Point_3& p, const Point_3& q, const Point_3& r ) const 02940 { 02941 return operator()(Plane_3(p, q, r)); 02942 } 02943 }; 02944 02945 template <typename K> 02946 class Construct_perpendicular_vector_2 02947 { 02948 typedef typename K::Vector_2 Vector_2; 02949 public: 02950 typedef Vector_2 result_type; 02951 02952 Vector_2 02953 operator()( const Vector_2& v, Orientation o) const 02954 { 02955 CGAL_kernel_precondition( o != COLLINEAR ); 02956 if (o == COUNTERCLOCKWISE) 02957 return K().construct_vector_2_object()(-v.hy(), v.hx(), v.hw()); 02958 else 02959 return K().construct_vector_2_object()(v.hy(), -v.hx(), v.hw()); 02960 } 02961 }; 02962 02963 template <typename K> 02964 class Construct_perpendicular_direction_2 02965 { 02966 typedef typename K::Direction_2 Direction_2; 02967 public: 02968 typedef Direction_2 result_type; 02969 02970 Direction_2 02971 operator()( const Direction_2& d, Orientation o) const 02972 { 02973 CGAL_kernel_precondition(o != COLLINEAR); 02974 if (o == COUNTERCLOCKWISE) { 02975 return Direction_2(-d.dy(), d.dx()); 02976 } else { 02977 return Direction_2(d.dy(), -d.dx()); 02978 } 02979 } 02980 }; 02981 02982 02983 template <typename K> 02984 class Construct_perpendicular_line_2 02985 { 02986 typedef typename K::Line_2 Line_2; 02987 typedef typename K::Point_2 Point_2; 02988 public: 02989 typedef Line_2 result_type; 02990 02991 Line_2 02992 operator()( const Line_2& l, const Point_2& p) const 02993 { return typename K::Line_2( -l.b()*p.hw(), l.a()*p.hw(), l.b()*p.hx() - l.a()*p.hy()); } 02994 }; 02995 02996 02997 template <typename K> 02998 class Construct_point_2 02999 { 03000 typedef typename K::RT RT; 03001 typedef typename K::FT FT; 03002 typedef typename K::Point_2 Point_2; 03003 typedef typename K::Vector_2 Vector_2; 03004 typedef typename K::Line_2 Line_2; 03005 typedef typename Point_2::Rep Rep; 03006 public: 03007 typedef Point_2 result_type; 03008 03009 Rep // Point_2 03010 operator()(Return_base_tag, Origin o) const 03011 { return Rep(o); } 03012 03013 template < typename Tx, typename Ty > 03014 Rep // Point_2 03015 operator()(Return_base_tag, const Tx & x, const Ty & y) const 03016 { return Rep(x, y); } 03017 03018 Rep // Point_2 03019 operator()(Return_base_tag, const RT& x, const RT& y, const RT& w) const 03020 { return Rep(x, y, w); } 03021 03022 Point_2 03023 operator()(const Line_2& l) const 03024 { 03025 CGAL_kernel_precondition( ! l.is_degenerate() ); 03026 if (l.is_vertical() ) 03027 { 03028 return Rep(-l.c(), RT(0) , l.a() ); 03029 } else { 03030 return Rep(RT(0) , -l.c(), l.b() ); 03031 } 03032 } 03033 03034 Point_2 03035 operator()(const Line_2& l, int i) const 03036 { 03037 Point_2 p = K().construct_point_2_object()(l); 03038 Vector_2 v = K().construct_vector_2_object()(l); 03039 return K().construct_translated_point_2_object() 03040 (p, K().construct_scaled_vector_2_object()(v, RT(i))); 03041 } 03042 03043 03044 Point_2 03045 operator()(Origin o) const 03046 { return this->operator()(Return_base_tag(), o); } 03047 03048 template < typename Tx, typename Ty > 03049 Point_2 03050 operator()(const Tx & x, const Ty & y) const 03051 { return this->operator()(Return_base_tag(), x, y); } 03052 03053 Point_2 03054 operator()(const RT& x, const RT& y, const RT& w) const 03055 { return this->operator()(Return_base_tag(), x, y, w); } 03056 }; 03057 03058 template <typename K> 03059 class Construct_point_3 03060 { 03061 typedef typename K::RT RT; 03062 typedef typename K::FT FT; 03063 typedef typename K::Point_3 Point_3; 03064 typedef typename Point_3::Rep Rep; 03065 public: 03066 typedef Point_3 result_type; 03067 03068 Rep // Point_3 03069 operator()(Return_base_tag, Origin o) const 03070 { return Rep(o); } 03071 03072 template < typename Tx, typename Ty, typename Tz > 03073 Rep // Point_3 03074 operator()(Return_base_tag, const Tx& x, const Ty& y, const Tz& z) const 03075 { return Rep(x, y, z); } 03076 03077 Rep // Point_3 03078 operator()(Return_base_tag, const FT& x, const FT& y, const FT& z) const 03079 { return Rep(x, y, z); } 03080 03081 Rep // Point_3 03082 operator()(Return_base_tag, const RT& x, const RT& y, const RT& z, const RT& w) const 03083 { return Rep(x, y, z, w); } 03084 03085 03086 Point_3 03087 operator()(Origin o) const 03088 { return this->operator()(Return_base_tag(), o); } 03089 03090 template < typename Tx, typename Ty, typename Tz > 03091 Point_3 03092 operator()(const Tx& x, const Ty& y, const Tz& z) const 03093 { return this->operator()(Return_base_tag(), x, y, z); } 03094 03095 Point_3 03096 operator()(const FT& x, const FT& y, const FT& z) const 03097 { return this->operator()(Return_base_tag(), x, y, z); } 03098 03099 Point_3 03100 operator()(const RT& x, const RT& y, const RT& z, const RT& w) const 03101 { return this->operator()(Return_base_tag(), x, y, z, w); } 03102 }; 03103 03104 03105 template <typename K> 03106 class Construct_projected_point_2 03107 { 03108 typedef typename K::Point_2 Point_2; 03109 typedef typename K::Direction_2 Direction_2; 03110 typedef typename K::Line_2 Line_2; 03111 public: 03112 typedef Point_2 result_type; 03113 03114 Point_2 03115 operator()( const Line_2& l, const Point_2& p ) const 03116 { 03117 CGAL_kernel_precondition( ! l.is_degenerate() ); 03118 Line_2 l2( p, Direction_2( l.a(), l.b() )); 03119 return Point_2( l.b()*l2.c() - l2.b()*l.c(), 03120 l2.a()*l.c() - l.a()*l2.c(), 03121 l.a()*l2.b() - l2.a()*l.b() ); 03122 } 03123 }; 03124 03125 template <typename K> 03126 class Construct_projected_point_3 03127 { 03128 typedef typename K::RT RT; 03129 typedef typename K::Point_3 Point_3; 03130 typedef typename K::Plane_3 Plane_3; 03131 typedef typename K::Line_3 Line_3; 03132 typedef typename K::Vector_3 Vector_3; 03133 public: 03134 typedef Point_3 result_type; 03135 03136 Point_3 03137 operator()( const Line_3& l, const Point_3& p ) const 03138 { 03139 if ( l.has_on(p) ) 03140 return p; 03141 Vector_3 v = p - l.point(); 03142 const RT& vx = v.hx(); 03143 const RT& vy = v.hy(); 03144 const RT& vz = v.hz(); 03145 const RT& vw = v.hw(); 03146 Vector_3 dir = l.to_vector(); 03147 const RT& dx = dir.hx(); 03148 const RT& dy = dir.hy(); 03149 const RT& dz = dir.hz(); 03150 const RT& dw = dir.hw(); 03151 03152 RT lambda_num = (vx*dx + vy*dy + vz*dz)*dw; // *dw 03153 RT lambda_den = (dx*dx + dy*dy + dz*dz)*vw; // *dw 03154 03155 return l.point() + ( (lambda_num * dir)/lambda_den ); 03156 } 03157 03158 Point_3 03159 operator()( const Plane_3& h, const Point_3& p ) const 03160 { return h.rep().projection(p); } 03161 }; 03162 03163 template <class K> 03164 class Construct_radical_line_2 03165 { 03166 typedef typename K::Line_2 Line_2; 03167 typedef typename K::Circle_2 Circle_2; 03168 typedef typename K::RT RT; 03169 typedef typename K::FT FT; 03170 03171 public: 03172 03173 typedef Line_2 result_type; 03174 03175 result_type 03176 operator() (const Circle_2 & c1, const Circle_2 & c2) const 03177 { 03178 // Concentric Circles don't have radical line 03179 CGAL_kernel_precondition (c1.center() != c2.center()); 03180 const FT a = 2*(c2.center().x() - c1.center().x()); 03181 const FT b = 2*(c2.center().y() - c1.center().y()); 03182 const FT c = CGAL::square(c1.center().x()) + 03183 CGAL::square(c1.center().y()) - c1.squared_radius() - 03184 CGAL::square(c2.center().x()) - 03185 CGAL::square(c2.center().y()) + c2.squared_radius(); 03186 const RT aa = a.numerator() * b.denominator() * c.denominator(); 03187 const RT bb = a.denominator() * b.numerator() * c.denominator(); 03188 const RT cc = a.denominator() * b.denominator() * c.numerator(); 03189 return Line_2(aa, bb, cc); 03190 } 03191 }; 03192 03193 03194 template <typename K> 03195 class Construct_scaled_vector_2 03196 { 03197 typedef typename K::RT RT; 03198 typedef typename K::FT FT; 03199 typedef typename K::Vector_2 Vector_2; 03200 public: 03201 typedef Vector_2 result_type; 03202 03203 Vector_2 03204 operator()( const Vector_2& v, const RT& c) const 03205 { 03206 return Vector_2(c * v.hx(), c * v.hy(), v.hw()); 03207 } 03208 03209 Vector_2 03210 operator()( const Vector_2& v, const FT& c) const 03211 { 03212 return Vector_2( v.hx()*c.numerator(), v.hy()*c.numerator(), 03213 v.hw()*c.denominator() ); } 03214 }; 03215 03216 template <typename K> 03217 class Construct_scaled_vector_3 03218 { 03219 typedef typename K::RT RT; 03220 typedef typename K::FT FT; 03221 typedef typename K::Vector_3 Vector_3; 03222 public: 03223 typedef Vector_3 result_type; 03224 03225 Vector_3 03226 operator()( const Vector_3& v, const RT& c) const 03227 { 03228 return Vector_3(c * v.hx(), c * v.hy(), c * v.hz(), v.hw()); 03229 } 03230 03231 Vector_3 03232 operator()( const Vector_3& v, const FT& c) const 03233 { 03234 return Vector_3( v.hx()*c.numerator(), v.hy()*c.numerator(), 03235 v.hz()*c.numerator(), v.hw()*c.denominator() ); } 03236 }; 03237 03238 template <typename K> 03239 class Construct_translated_point_2 03240 { 03241 typedef typename K::Point_2 Point_2; 03242 typedef typename K::Vector_2 Vector_2; 03243 public: 03244 typedef Point_2 result_type; 03245 03246 Point_2 03247 operator()( const Point_2& p, const Vector_2& v) const 03248 { 03249 return Point_2( p.hx()*v.hw() + v.hx()*p.hw(), 03250 p.hy()*v.hw() + v.hy()*p.hw(), 03251 p.hw()*v.hw() ); 03252 } 03253 03254 Point_2 03255 operator()( const Origin&, const Vector_2& v) const 03256 { 03257 return Point_2( v.hx(), v.hy(), v.hw() ); 03258 } 03259 }; 03260 03261 template <typename K> 03262 class Construct_translated_point_3 03263 { 03264 typedef typename K::Point_3 Point_3; 03265 typedef typename K::Vector_3 Vector_3; 03266 public: 03267 typedef Point_3 result_type; 03268 03269 Point_3 03270 operator()( const Point_3& p, const Vector_3& v) const 03271 { 03272 return Point_3(p.hx()*v.hw() + v.hx()*p.hw(), 03273 p.hy()*v.hw() + v.hy()*p.hw(), 03274 p.hz()*v.hw() + v.hz()*p.hw(), 03275 p.hw()*v.hw() ); 03276 } 03277 03278 Point_3 03279 operator()( const Origin&, const Vector_3& v) const 03280 { 03281 return Point_3( v.hx(), v.hy(), v.hz(), v.hw() ); 03282 } 03283 }; 03284 03285 template <typename K> 03286 class Construct_vector_2 03287 { 03288 typedef typename K::RT RT; 03289 typedef typename K::FT FT; 03290 typedef typename K::Segment_2 Segment_2; 03291 typedef typename K::Ray_2 Ray_2; 03292 typedef typename K::Line_2 Line_2; 03293 typedef typename K::Vector_2 Vector_2; 03294 typedef typename K::Point_2 Point_2; 03295 typedef typename K::Direction_2 Direction_2; 03296 typedef typename Vector_2::Rep Rep; 03297 public: 03298 typedef Vector_2 result_type; 03299 03300 Rep // Vector_2 03301 operator()(Return_base_tag, const Point_2& p, const Point_2& q) const 03302 { 03303 return Rep( q.hx()*p.hw() - p.hx()*q.hw(), 03304 q.hy()*p.hw() - p.hy()*q.hw(), 03305 p.hw()*q.hw() ); 03306 } 03307 03308 Rep // Vector_2 03309 operator()(Return_base_tag, const Origin& , const Point_2& q) const 03310 { 03311 return Rep( q.hx(), q.hy(), q.hw() ); 03312 } 03313 03314 Rep // Vector_2 03315 operator()(Return_base_tag, const Point_2& p, const Origin& ) const 03316 { 03317 return Rep( - p.hx(), - p.hy(), p.hw() ); 03318 } 03319 03320 Rep // Vector_2 03321 operator()(Return_base_tag, const Direction_2& d ) const 03322 { return Rep(d.dx(), d.dy()); } 03323 03324 Rep // Vector_2 03325 operator()(Return_base_tag, const Segment_2& s) const 03326 { return K().construct_vector_2_object()(s.source(), s.target()); } 03327 03328 Rep // Vector_2 03329 operator()(Return_base_tag, const Ray_2& r) const 03330 { return K().construct_vector_2_object()(r.source(), r.point(1)); } 03331 03332 Rep // Vector_2 03333 operator()(Return_base_tag, const Line_2& l) const 03334 { return K().construct_vector_2_object()( l.b(), -l.a()); } 03335 03336 Rep // Vector_2 03337 operator()(Return_base_tag, Null_vector) const 03338 { return Rep(RT(0), RT(0), RT(1)); } 03339 03340 template < typename Tx, typename Ty > 03341 Rep // Vector_2 03342 operator()(Return_base_tag, const Tx & x, const Ty & y) const 03343 { return Rep(x, y); } 03344 03345 Rep // Vector_2 03346 operator()(Return_base_tag, const RT& x, const RT& y, const RT& w) const 03347 { return Rep(x, y, w); } 03348 03349 03350 Vector_2 03351 operator()( const Point_2& p, const Point_2& q) const 03352 { return this->operator()(Return_base_tag(), p, q); } 03353 03354 Vector_2 03355 operator()( const Origin& o, const Point_2& q) const 03356 { return this->operator()(Return_base_tag(), o, q); } 03357 03358 Vector_2 03359 operator()( const Point_2& p, const Origin& o) const 03360 { return this->operator()(Return_base_tag(), p, o); } 03361 03362 Vector_2 03363 operator()( const Direction_2& d ) const 03364 { return this->operator()(Return_base_tag(), d); } 03365 03366 Vector_2 03367 operator()( const Segment_2& s) const 03368 { return this->operator()(Return_base_tag(), s); } 03369 03370 Vector_2 03371 operator()( const Ray_2& r) const 03372 { return this->operator()(Return_base_tag(), r); } 03373 03374 Vector_2 03375 operator()( const Line_2& l) const 03376 { return this->operator()(Return_base_tag(), l); } 03377 03378 Vector_2 03379 operator()( Null_vector n) const 03380 { return this->operator()(Return_base_tag(), n); } 03381 03382 template < typename Tx, typename Ty > 03383 Vector_2 03384 operator()(const Tx & x, const Ty & y) const 03385 { return this->operator()(Return_base_tag(), x, y); } 03386 03387 Vector_2 03388 operator()( const RT& x, const RT& y, const RT& w) const 03389 { return this->operator()(Return_base_tag(), x, y, w); } 03390 }; 03391 03392 template <typename K> 03393 class Construct_vector_3 03394 { 03395 typedef typename K::RT RT; 03396 typedef typename K::FT FT; 03397 typedef typename K::Direction_3 Direction_3; 03398 typedef typename K::Segment_3 Segment_3; 03399 typedef typename K::Ray_3 Ray_3; 03400 typedef typename K::Line_3 Line_3; 03401 typedef typename K::Vector_3 Vector_3; 03402 typedef typename K::Point_3 Point_3; 03403 typedef typename Vector_3::Rep Rep; 03404 public: 03405 typedef Vector_3 result_type; 03406 03407 Rep // Vector_3 03408 operator()(Return_base_tag, const Point_3& p, const Point_3& q) const 03409 { 03410 return Rep(q.hx()*p.hw() - p.hx()*q.hw(), 03411 q.hy()*p.hw() - p.hy()*q.hw(), 03412 q.hz()*p.hw() - p.hz()*q.hw(), 03413 q.hw()*p.hw() ); 03414 } 03415 03416 Rep // Vector_3 03417 operator()(Return_base_tag, const Origin&, const Point_3& q) const 03418 { 03419 return Rep( q.hx(), q.hy(), q.hz(), q.hw()); 03420 } 03421 03422 Rep // Vector_3 03423 operator()(Return_base_tag, const Point_3& p, const Origin& ) const 03424 { 03425 return Rep( - p.hx(), - p.hy(), - p.hz(), p.hw() ); 03426 } 03427 03428 Rep // Vector_3 03429 operator()(Return_base_tag, const Direction_3& d) const 03430 { return d.rep().to_vector(); } 03431 03432 Rep // Vector_3 03433 operator()(Return_base_tag, const Segment_3& s) const 03434 { return s.rep().to_vector(); } 03435 // { return this->operator()(s.start(), s.end()); } 03436 03437 Rep // Vector_3 03438 operator()(Return_base_tag, const Ray_3& r) const 03439 { return r.rep().to_vector(); } 03440 03441 Rep // Vector_3 03442 operator()(Return_base_tag, const Line_3& l) const 03443 { return l.rep().to_vector(); } 03444 03445 Rep // Vector_3 03446 operator()(Return_base_tag, const Null_vector&) const 03447 { return Rep(RT(0), RT(0), RT(0), RT(1)); } 03448 03449 template < typename Tx, typename Ty, typename Tz > 03450 Rep // Vector_3 03451 operator()(Return_base_tag, const Tx & x, const Ty & y, const Tz & z) const 03452 { return Rep(x, y, z); } 03453 03454 Rep // Vector_3 03455 operator()(Return_base_tag, const RT& x, const RT& y, const RT& z, const RT& w) const 03456 { return Rep(x, y, z, w); } 03457 03458 03459 Vector_3 03460 operator()( const Point_3& p, const Point_3& q) const 03461 { return this->operator()(Return_base_tag(), p, q); } 03462 03463 Vector_3 03464 operator()( const Origin& o, const Point_3& q) const 03465 { return this->operator()(Return_base_tag(), o, q); } 03466 03467 Vector_3 03468 operator()( const Point_3& p, const Origin& q) const 03469 { return this->operator()(Return_base_tag(), p, q); } 03470 03471 Vector_3 03472 operator()( const Direction_3& d) const 03473 { return this->operator()(Return_base_tag(), d); } 03474 03475 Vector_3 03476 operator()( const Segment_3& s) const 03477 { return this->operator()(Return_base_tag(), s); } 03478 03479 Vector_3 03480 operator()( const Ray_3& r) const 03481 { return this->operator()(Return_base_tag(), r); } 03482 03483 Vector_3 03484 operator()( const Line_3& l) const 03485 { return this->operator()(Return_base_tag(), l); } 03486 03487 Vector_3 03488 operator()( const Null_vector& n) const 03489 { return this->operator()(Return_base_tag(), n); } 03490 03491 template < typename Tx, typename Ty, typename Tz > 03492 Vector_3 03493 operator()(const Tx & x, const Ty & y, const Tz & z) const 03494 { return this->operator()(Return_base_tag(), x, y, z); } 03495 03496 Vector_3 03497 operator()( const RT& x, const RT& y, const RT& z, const RT& w) const 03498 { return this->operator()(Return_base_tag(), x, y, z, w); } 03499 }; 03500 03501 template <typename K> 03502 class Construct_vertex_2 03503 { 03504 typedef typename K::Point_2 Point_2; 03505 typedef typename K::Segment_2 Segment_2; 03506 typedef typename K::Iso_rectangle_2 Iso_rectangle_2; 03507 typedef typename K::Triangle_2 Triangle_2; 03508 public: 03509 typedef Point_2 result_type; 03510 03511 const Point_2 & 03512 operator()( const Segment_2& s, int i) const 03513 { return s.vertex(i); } 03514 03515 const Point_2 & 03516 operator()( const Triangle_2& t, int i) const 03517 { return t.rep().vertex(i); } 03518 03519 const Point_2 03520 operator()( const Iso_rectangle_2& r, int i) const 03521 { 03522 switch (i%4) { 03523 case 0: return (r.min)(); 03524 case 1: 03525 return Point_2( (r.max)().hx()*(r.min)().hw(), 03526 (r.min)().hy()*(r.max)().hw(), 03527 (r.min)().hw()*(r.max)().hw() ); 03528 case 2: return (r.max)(); 03529 default: return Point_2( (r.min)().hx()*(r.max)().hw(), 03530 (r.max)().hy()*(r.min)().hw(), 03531 (r.min)().hw()*(r.max)().hw() ); 03532 } 03533 } 03534 }; 03535 03536 } //namespace HomogeneousKernelFunctors 03537 03538 03539 #ifndef CGAL_CFG_DONT_OVERLOAD_TOO_MUCH 03540 template < typename K> 03541 struct Qualified_result_of<HomogeneousKernelFunctors::Construct_vertex_2<K>, typename K::Segment_2, int > 03542 { 03543 typedef typename K::Point_2 const & type; 03544 }; 03545 03546 template < typename K> 03547 struct Qualified_result_of<HomogeneousKernelFunctors::Construct_vertex_2<K>, typename K::Triangle_2, int > 03548 { 03549 typedef typename K::Point_2 const & type; 03550 }; 03551 #endif 03552 03553 // For Iso_rectangle the non specialized template will do the right thing, namely return a copy of a point 03554 03555 namespace HomogeneousKernelFunctors { 03556 03557 template <typename K> 03558 class Coplanar_orientation_3 03559 { 03560 typedef typename K::Point_3 Point_3; 03561 #ifdef CGAL_kernel_exactness_preconditions 03562 typedef typename K::Coplanar_3 Coplanar_3; 03563 typedef typename K::Collinear_3 Collinear_3; 03564 Coplanar_3 cp; 03565 Collinear_3 cl; 03566 #endif // CGAL_kernel_exactness_preconditions 03567 public: 03568 typedef typename K::Orientation result_type; 03569 03570 #ifdef CGAL_kernel_exactness_preconditions 03571 Coplanar_orientation_3() {} 03572 Coplanar_orientation_3(const Coplanar_3& cp_, const Collinear_3& cl_) 03573 : cp(cp_), cl(cl_) {} 03574 #endif // CGAL_kernel_exactness_preconditions 03575 03576 result_type 03577 operator()(const Point_3& p, const Point_3& q, const Point_3& r) const 03578 { 03579 Orientation oxy_pqr = sign_of_determinant(p.hx(), p.hy(), p.hw(), 03580 q.hx(), q.hy(), q.hw(), 03581 r.hx(), r.hy(), r.hw()); 03582 if (oxy_pqr != COLLINEAR) 03583 return oxy_pqr; 03584 03585 Orientation oyz_pqr = sign_of_determinant(p.hy(), p.hz(), p.hw(), 03586 q.hy(), q.hz(), q.hw(), 03587 r.hy(), r.hz(), r.hw()); 03588 if (oyz_pqr != COLLINEAR) 03589 return oyz_pqr; 03590 03591 return sign_of_determinant(p.hx(), p.hz(), p.hw(), 03592 q.hx(), q.hz(), q.hw(), 03593 r.hx(), r.hz(), r.hw()); 03594 } 03595 03596 result_type 03597 operator()( const Point_3& p, const Point_3& q, 03598 const Point_3& r, const Point_3& s) const 03599 { 03600 // p,q,r,s supposed to be coplanar 03601 // p,q,r supposed to be non collinear 03602 // tests whether s is on the same side of p,q as r 03603 // returns : 03604 // COLLINEAR if pqr collinear 03605 // POSITIVE if qrp and qrs have the same orientation 03606 // NEGATIVE if qrp and qrs have opposite orientations 03607 CGAL_kernel_exactness_precondition( ! cl(p, q, r) ); 03608 CGAL_kernel_exactness_precondition( cp(p, q, r, s) ); 03609 03610 // compute orientation of p,q,s in this plane P: 03611 Orientation save; 03612 if ( (save = sign_of_determinant(p.hy(), p.hz(), p.hw(), 03613 q.hy(), q.hz(), q.hw(), 03614 r.hy(), r.hz(), r.hw())) != COLLINEAR) 03615 { return save * sign_of_determinant(p.hy(), p.hz(), p.hw(), 03616 q.hy(), q.hz(), q.hw(), 03617 s.hy(), s.hz(), s.hw()); 03618 } 03619 if ( (save = sign_of_determinant(p.hx(), p.hz(), p.hw(), 03620 q.hx(), q.hz(), q.hw(), 03621 r.hx(), r.hz(), r.hw())) != COLLINEAR) 03622 { return save * sign_of_determinant(p.hx(), p.hz(), p.hw(), 03623 q.hx(), q.hz(), q.hw(), 03624 s.hx(), s.hz(), s.hw()); 03625 } 03626 if ( (save = sign_of_determinant(p.hx(), p.hy(), p.hw(), 03627 q.hx(), q.hy(), q.hw(), 03628 r.hx(), r.hy(), r.hw())) != COLLINEAR) 03629 { return save * sign_of_determinant( p.hx(), p.hy(), p.hw(), 03630 q.hx(), q.hy(), q.hw(), 03631 s.hx(), s.hy(), s.hw()); 03632 } 03633 CGAL_kernel_assertion( false); 03634 return COLLINEAR; 03635 } 03636 }; 03637 03638 template <typename K> 03639 class Coplanar_side_of_bounded_circle_3 03640 { 03641 typedef typename K::Point_3 Point_3; 03642 #ifdef CGAL_kernel_exactness_preconditions 03643 typedef typename K::Coplanar_3 Coplanar_3; 03644 typedef typename K::Collinear_3 Collinear_3; 03645 Coplanar_3 cp; 03646 Collinear_3 cl; 03647 #endif // CGAL_kernel_exactness_preconditions 03648 public: 03649 typedef typename K::Bounded_side result_type; 03650 03651 #ifdef CGAL_kernel_exactness_preconditions 03652 Coplanar_side_of_bounded_circle_3() {} 03653 Coplanar_side_of_bounded_circle_3(const Coplanar_3& cp_, 03654 const Collinear_3& cl_) 03655 : cp(cp_), cl(cl_) {} 03656 #endif // CGAL_kernel_exactness_preconditions 03657 03658 result_type 03659 operator()( const Point_3& p, const Point_3& q, 03660 const Point_3& r, const Point_3& t) const 03661 { 03662 // p,q,r,t are supposed to be coplanar. 03663 // p,q,r determine an orientation of this plane (not collinear). 03664 // returns the equivalent of side_of_bounded_circle(p,q,r,t) 03665 // in this plane 03666 CGAL_kernel_exactness_precondition( cp(p,q,r,t) ); 03667 CGAL_kernel_exactness_precondition( !cl(p,q,r) ); 03668 return enum_cast<Bounded_side>( 03669 side_of_oriented_sphere(p, q, r, t+cross_product(q-p, r-p), t)); 03670 } // FIXME 03671 }; 03672 03673 template <typename K> 03674 class Equal_xy_3 03675 { 03676 typedef typename K::Point_3 Point_3; 03677 public: 03678 typedef typename K::Boolean result_type; 03679 03680 result_type 03681 operator()( const Point_3& p, const Point_3& q) const 03682 { 03683 return (p.hx() * q.hw() == q.hx() * p.hw() ) 03684 && (p.hy() * q.hw() == q.hy() * p.hw() ); 03685 } 03686 }; 03687 03688 template <typename K> 03689 class Equal_x_2 03690 { 03691 typedef typename K::Point_2 Point_2; 03692 public: 03693 typedef typename K::Boolean result_type; 03694 03695 result_type 03696 operator()( const Point_2& p, const Point_2& q) const 03697 { return p.hx()*q.hw() == q.hx()*p.hw(); } 03698 }; 03699 03700 template <typename K> 03701 class Equal_x_3 03702 { 03703 typedef typename K::Point_3 Point_3; 03704 public: 03705 typedef typename K::Boolean result_type; 03706 03707 result_type 03708 operator()( const Point_3& p, const Point_3& q) const 03709 { return p.hx()*q.hw() == q.hx()*p.hw(); } 03710 }; 03711 03712 template <typename K> 03713 class Equal_y_2 03714 { 03715 typedef typename K::Point_2 Point_2; 03716 public: 03717 typedef typename K::Boolean result_type; 03718 03719 result_type 03720 operator()( const Point_2& p, const Point_2& q) const 03721 { return p.hy()*q.hw() == q.hy()*p.hw(); } 03722 }; 03723 03724 template <typename K> 03725 class Equal_y_3 03726 { 03727 typedef typename K::Point_3 Point_3; 03728 public: 03729 typedef typename K::Boolean result_type; 03730 03731 result_type 03732 operator()( const Point_3& p, const Point_3& q) const 03733 { return p.hy()*q.hw() == q.hy()*p.hw(); } 03734 }; 03735 03736 template <typename K> 03737 class Equal_z_3 03738 { 03739 typedef typename K::Point_3 Point_3; 03740 public: 03741 typedef typename K::Boolean result_type; 03742 03743 result_type 03744 operator()( const Point_3& p, const Point_3& q) const 03745 { return p.hz()*q.hw() == q.hz()*p.hw(); } 03746 }; 03747 03748 template <typename K> 03749 class Has_on_3 03750 { 03751 typedef typename K::Point_3 Point_3; 03752 typedef typename K::Line_3 Line_3; 03753 typedef typename K::Ray_3 Ray_3; 03754 typedef typename K::Segment_3 Segment_3; 03755 typedef typename K::Plane_3 Plane_3; 03756 typedef typename K::Triangle_3 Triangle_3; 03757 typedef typename K::Tetrahedron_3 Tetrahedron_3; 03758 public: 03759 typedef typename K::Boolean result_type; 03760 03761 result_type 03762 operator()( const Line_3& l, const Point_3& p) const 03763 { return l.rep().has_on(p); } 03764 03765 result_type 03766 operator()( const Ray_3& r, const Point_3& p) const 03767 { return r.rep().has_on(p); } 03768 03769 result_type 03770 operator()( const Segment_3& s, const Point_3& p) const 03771 { return s.has_on(p); } 03772 03773 result_type 03774 operator()( const Plane_3& pl, const Point_3& p) const 03775 { return pl.rep().has_on(p); } 03776 03777 result_type 03778 operator()( const Plane_3& pl, const Line_3& l) const 03779 { return pl.rep().has_on(l); } 03780 03781 result_type 03782 operator()( const Triangle_3& t, const Point_3& p) const 03783 { 03784 if (!t.is_degenerate() ) 03785 { 03786 Plane_3 sup_pl = t.supporting_plane(); 03787 if ( !sup_pl.has_on(p) ) 03788 { 03789 return false; 03790 } 03791 Tetrahedron_3 tetrapak( t.vertex(0), 03792 t.vertex(1), 03793 t.vertex(2), 03794 t.vertex(0) + sup_pl.orthogonal_vector()); 03795 return tetrapak.has_on_boundary(p); 03796 } 03797 Point_3 minp( t.vertex(0) ); 03798 Point_3 maxp( t.vertex(1) ); 03799 if (lexicographically_xyz_smaller(t.vertex(1),t.vertex(0)) ) 03800 { 03801 minp = t.vertex(1); 03802 maxp = t.vertex(0); 03803 } 03804 if (lexicographically_xyz_smaller(t.vertex(2),minp ) ) 03805 { 03806 minp = t.vertex(2); 03807 } 03808 if (lexicographically_xyz_smaller(maxp, t.vertex(2)) ) 03809 { 03810 maxp = t.vertex(2); 03811 } 03812 if (minp == maxp) 03813 { 03814 return (p == maxp); 03815 } 03816 Segment_3 s(minp,maxp); 03817 return s.has_on(p); 03818 } 03819 }; 03820 03821 template <typename K> 03822 class Less_distance_to_point_2 03823 { 03824 typedef typename K::Point_2 Point_2; 03825 public: 03826 typedef typename K::Boolean result_type; 03827 03828 result_type 03829 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 03830 { 03831 typedef typename K::RT RT; 03832 03833 const RT & phx = p.hx(); 03834 const RT & phy = p.hy(); 03835 const RT & phw = p.hw(); 03836 const RT & qhx = q.hx(); 03837 const RT & qhy = q.hy(); 03838 const RT & qhw = q.hw(); 03839 const RT & rhx = r.hx(); 03840 const RT & rhy = r.hy(); 03841 const RT & rhw = r.hw(); 03842 03843 RT dosd = // difference of squared distances 03844 03845 // phx * phx * qhw * qhw * rhw * rhw 03846 // -RT(2) * phx * qhx * phw * qhw * rhw * rhw 03847 // + qhx * qhx * phw * phw * rhw * rhw 03848 // 03849 // + phy * phy * qhw * qhw * rhw * rhw 03850 // -RT(2) * phy * qhy * phw * qhw * rhw * rhw 03851 // + qhy * qhy * phw * phw * rhw * rhw 03852 // 03853 // - ( phx * phx * qhw * qhw * rhw * rhw 03854 // -RT(2) * phx * rhx * phw * qhw * qhw * rhw 03855 // + rhx * rhx * phw * phw * qhw * qhw 03856 // 03857 // + phy * phy * qhw * qhw * rhw * rhw 03858 // -RT(2) * phy * rhy * phw * qhw * qhw * rhw 03859 // + rhy * rhy * phw * phw * qhw * qhw 03860 03861 rhw*rhw * ( phw * ( qhx*qhx + qhy*qhy ) 03862 - 2 * qhw * ( phx*qhx + phy*qhy ) 03863 ) 03864 - qhw*qhw * ( phw * ( rhx*rhx + rhy*rhy ) 03865 - 2 * rhw * ( phx*rhx + phy*rhy ) 03866 ); 03867 03868 03869 return dosd < 0; 03870 } 03871 }; 03872 03873 template <typename K> 03874 class Less_distance_to_point_3 03875 { 03876 typedef typename K::Point_3 Point_3; 03877 public: 03878 typedef typename K::Boolean result_type; 03879 03880 result_type 03881 operator()(const Point_3& p, const Point_3& q, const Point_3& r) const 03882 { 03883 typedef typename K::RT RT; 03884 03885 const RT & phx = p.hx(); 03886 const RT & phy = p.hy(); 03887 const RT & phz = p.hz(); 03888 const RT & phw = p.hw(); 03889 const RT & qhx = q.hx(); 03890 const RT & qhy = q.hy(); 03891 const RT & qhz = q.hz(); 03892 const RT & qhw = q.hw(); 03893 const RT & rhx = r.hx(); 03894 const RT & rhy = r.hy(); 03895 const RT & rhz = r.hz(); 03896 const RT & rhw = r.hw(); 03897 03898 RT dosd = // difference of squared distances 03899 rhw*rhw * ( phw * ( qhx*qhx + qhy*qhy + qhz*qhz ) 03900 - 2 * qhw * ( phx*qhx + phy*qhy + phz*qhz ) 03901 ) 03902 - qhw*qhw * ( phw * ( rhx*rhx + rhy*rhy + rhz*rhz ) 03903 - 2 * rhw * ( phx*rhx + phy*rhy + phz*rhz ) 03904 ); 03905 03906 return dosd < 0; 03907 } 03908 }; 03909 03910 template <typename K> 03911 class Less_signed_distance_to_line_2 03912 { 03913 typedef typename K::Point_2 Point_2; 03914 typedef typename K::Line_2 Line_2; 03915 public: 03916 typedef typename K::Boolean result_type; 03917 03918 result_type 03919 operator()(const Point_2& p, const Point_2& q, 03920 const Point_2& r, const Point_2& s) const 03921 { 03922 typedef typename K::RT RT; 03923 03924 const RT & phx= p.hx(); 03925 const RT & phy= p.hy(); 03926 const RT & phw= p.hw(); 03927 const RT & qhx= q.hx(); 03928 const RT & qhy= q.hy(); 03929 const RT & qhw= q.hw(); 03930 const RT & rhx= r.hx(); 03931 const RT & rhy= r.hy(); 03932 const RT & rhw= r.hw(); 03933 const RT & shx= s.hx(); 03934 const RT & shy= s.hy(); 03935 const RT & shw= s.hw(); 03936 03937 RT scaled_dist_r_minus_scaled_dist_s = 03938 ( rhx*shw - shx*rhw ) * (phy*qhw - qhy*phw) 03939 - ( rhy*shw - shy*rhw ) * (phx*qhw - qhx*phw); 03940 03941 return scaled_dist_r_minus_scaled_dist_s < 0; 03942 } 03943 03944 result_type 03945 operator()(const Line_2& l, const Point_2& p, 03946 const Point_2& q) const 03947 { 03948 typedef typename K::RT RT; 03949 03950 const RT & la = l.a(); 03951 const RT & lb = l.b(); 03952 const RT & phx= p.hx(); 03953 const RT & phy= p.hy(); 03954 const RT & phw= p.hw(); 03955 const RT & qhx= q.hx(); 03956 const RT & qhy= q.hy(); 03957 const RT & qhw= q.hw(); 03958 03959 RT scaled_dist_p_minus_scaled_dist_q = 03960 la*( phx*qhw - qhx*phw ) 03961 + lb*( phy*qhw - qhy*phw ); 03962 03963 return scaled_dist_p_minus_scaled_dist_q < 0; 03964 } 03965 }; 03966 03967 template <typename K> 03968 class Less_signed_distance_to_plane_3 03969 { 03970 typedef typename K::RT RT; 03971 typedef typename K::Point_3 Point_3; 03972 typedef typename K::Plane_3 Plane_3; 03973 typedef typename K::Construct_plane_3 Construct_plane_3; 03974 public: 03975 typedef typename K::Boolean result_type; 03976 03977 result_type 03978 operator()( const Plane_3& pl, const Point_3& p, const Point_3& q) const 03979 { 03980 const RT & pla = pl.a(); 03981 const RT & plb = pl.b(); 03982 const RT & plc = pl.c(); 03983 const RT & phx = p.hx(); 03984 const RT & phy = p.hy(); 03985 const RT & phz = p.hz(); 03986 const RT & phw = p.hw(); 03987 const RT & qhx = q.hx(); 03988 const RT & qhy = q.hy(); 03989 const RT & qhz = q.hz(); 03990 const RT & qhw = q.hw(); 03991 03992 RT scaled_dist_p_minus_scaled_dist_q = 03993 pla*( phx*qhw - qhx*phw ) 03994 + plb*( phy*qhw - qhy*phw ) 03995 + plc*( phz*qhw - qhz*phw ); 03996 03997 return scaled_dist_p_minus_scaled_dist_q < 0; 03998 } 03999 04000 result_type 04001 operator()(const Point_3& plp, const Point_3& plq, const Point_3& plr, 04002 const Point_3& p, const Point_3& q) const 04003 { 04004 Construct_plane_3 construct_plane_3; 04005 return operator()(construct_plane_3(plp, plq, plr), p, q); 04006 } 04007 }; 04008 04009 template <typename K> 04010 class Less_xyz_3 04011 { 04012 typedef typename K::Point_3 Point_3; 04013 typedef typename K::Compare_xyz_3 Compare_xyz_3; 04014 Compare_xyz_3 c; 04015 public: 04016 typedef typename K::Boolean result_type; 04017 04018 Less_xyz_3() {} 04019 Less_xyz_3(const Compare_xyz_3& c_) : c(c_) {} 04020 04021 result_type 04022 operator()( const Point_3& p, const Point_3& q) const 04023 { return c(p, q) == SMALLER; } 04024 }; 04025 04026 template <typename K> 04027 class Less_xy_2 04028 { 04029 typedef typename K::Point_2 Point_2; 04030 typedef typename K::Compare_xy_2 Compare_xy_2; 04031 Compare_xy_2 c; 04032 public: 04033 typedef typename K::Boolean result_type; 04034 04035 Less_xy_2() {} 04036 Less_xy_2(const Compare_xy_2& c_) : c(c_) {} 04037 04038 result_type 04039 operator()( const Point_2& p, const Point_2& q) const 04040 { return c(p, q) == SMALLER; } 04041 }; 04042 04043 template <typename K> 04044 class Less_xy_3 04045 { 04046 typedef typename K::Point_3 Point_3; 04047 typedef typename K::Compare_xy_3 Compare_xy_3; 04048 Compare_xy_3 c; 04049 public: 04050 typedef typename K::Boolean result_type; 04051 04052 Less_xy_3() {} 04053 Less_xy_3(const Compare_xy_3& c_) : c(c_) {} 04054 04055 result_type 04056 operator()( const Point_3& p, const Point_3& q) const 04057 { return c(p, q) == SMALLER; } 04058 }; 04059 04060 template <typename K> 04061 class Less_x_2 04062 { 04063 typedef typename K::Point_2 Point_2; 04064 public: 04065 typedef typename K::Boolean result_type; 04066 04067 result_type 04068 operator()( const Point_2& p, const Point_2& q) const 04069 { return ( p.hx()*q.hw() < q.hx()*p.hw() ); } 04070 }; 04071 04072 template <typename K> 04073 class Less_x_3 04074 { 04075 typedef typename K::Point_3 Point_3; 04076 public: 04077 typedef typename K::Boolean result_type; 04078 04079 result_type 04080 operator()( const Point_3& p, const Point_3& q) const 04081 { return ( p.hx()*q.hw() < q.hx()*p.hw() ); } 04082 }; 04083 04084 template <typename K> 04085 class Less_yx_2 04086 { 04087 typedef typename K::Point_2 Point_2; 04088 public: 04089 typedef typename K::Boolean result_type; 04090 04091 result_type 04092 operator()( const Point_2& p, const Point_2& q) const 04093 { 04094 typedef typename K::RT RT; 04095 04096 const RT& phx = p.hx(); 04097 const RT& phy = p.hy(); 04098 const RT& phw = p.hw(); 04099 const RT& qhx = q.hx(); 04100 const RT& qhy = q.hy(); 04101 const RT& qhw = q.hw(); 04102 04103 RT pV = phy * qhw; 04104 RT qV = qhy * phw; 04105 if ( qV < pV ) 04106 { 04107 return false; 04108 } 04109 else if ( pV < qV ) 04110 { 04111 return true; 04112 } 04113 pV = phx * qhw; 04114 qV = qhx * phw; 04115 return ( pV < qV ); 04116 } 04117 }; 04118 04119 template <typename K> 04120 class Less_y_2 04121 { 04122 typedef typename K::Point_2 Point_2; 04123 public: 04124 typedef typename K::Boolean result_type; 04125 04126 result_type 04127 operator()( const Point_2& p, const Point_2& q) const 04128 { return ( p.hy()*q.hw() < q.hy()*p.hw() ); } 04129 }; 04130 04131 template <typename K> 04132 class Less_y_3 04133 { 04134 typedef typename K::Point_3 Point_3; 04135 public: 04136 typedef typename K::Boolean result_type; 04137 04138 result_type 04139 operator()( const Point_3& p, const Point_3& q) const 04140 { return ( p.hy()*q.hw() < q.hy()*p.hw() ); } 04141 }; 04142 04143 template <typename K> 04144 class Less_z_3 04145 { 04146 typedef typename K::Point_3 Point_3; 04147 public: 04148 typedef typename K::Boolean result_type; 04149 04150 result_type 04151 operator()( const Point_3& p, const Point_3& q) const 04152 { return (p.hz() * q.hw() < q.hz() * p.hw() ); } 04153 }; 04154 04155 template <typename K> 04156 class Orientation_2 04157 { 04158 typedef typename K::Point_2 Point_2; 04159 typedef typename K::Vector_2 Vector_2; 04160 typedef typename K::Circle_2 Circle_2; 04161 public: 04162 typedef typename K::Orientation result_type; 04163 04164 result_type 04165 operator()(const Point_2& p, const Point_2& q, const Point_2& r) const 04166 { 04167 typedef typename K::RT RT; 04168 04169 const RT& phx = p.hx(); 04170 const RT& phy = p.hy(); 04171 const RT& phw = p.hw(); 04172 const RT& qhx = q.hx(); 04173 const RT& qhy = q.hy(); 04174 const RT& qhw = q.hw(); 04175 const RT& rhx = r.hx(); 04176 const RT& rhy = r.hy(); 04177 const RT& rhw = r.hw(); 04178 04179 // | A B | 04180 // | C D | 04181 04182 RT A = phx*rhw - phw*rhx; 04183 RT B = phy*rhw - phw*rhy; 04184 RT C = qhx*rhw - qhw*rhx; 04185 RT D = qhy*rhw - qhw*rhy; 04186 04187 return CGAL_NTS compare(A*D, B*C); 04188 } 04189 04190 result_type 04191 operator()(const Vector_2& u, const Vector_2& v) const 04192 { 04193 return sign_of_determinant(u.hx(), u.hy(), 04194 v.hx(), v.hy()); 04195 } 04196 04197 result_type 04198 operator()(const Circle_2& c) const 04199 { 04200 return c.rep().orientation(); 04201 } 04202 }; 04203 04204 template <typename K> 04205 class Orientation_3 04206 { 04207 typedef typename K::Point_3 Point_3; 04208 typedef typename K::Vector_3 Vector_3; 04209 typedef typename K::Tetrahedron_3 Tetrahedron_3; 04210 typedef typename K::Sphere_3 Sphere_3; 04211 public: 04212 typedef typename K::Orientation result_type; 04213 04214 result_type 04215 operator()( const Point_3& p, const Point_3& q, 04216 const Point_3& r, const Point_3& s) const 04217 { 04218 // Two rows are switched, because of the homogeneous column. 04219 return sign_of_determinant( p.hx(), p.hy(), p.hz(), p.hw(), 04220 r.hx(), r.hy(), r.hz(), r.hw(), 04221 q.hx(), q.hy(), q.hz(), q.hw(), 04222 s.hx(), s.hy(), s.hz(), s.hw()); 04223 } 04224 04225 result_type 04226 operator()( const Vector_3& u, const Vector_3& v, const Vector_3& w) const 04227 { 04228 return sign_of_determinant( u.hx(), u.hy(), u.hz(), 04229 v.hx(), v.hy(), v.hz(), 04230 w.hx(), w.hy(), w.hz()); 04231 } 04232 04233 result_type 04234 operator()( const Tetrahedron_3& t) const 04235 { 04236 return t.rep().orientation(); 04237 } 04238 04239 result_type 04240 operator()(const Sphere_3& s) const 04241 { 04242 return s.rep().orientation(); 04243 } 04244 }; 04245 04246 04247 template <typename K> 04248 class Oriented_side_2 04249 { 04250 typedef typename K::RT RT; 04251 typedef typename K::Point_2 Point_2; 04252 typedef typename K::Circle_2 Circle_2; 04253 typedef typename K::Line_2 Line_2; 04254 typedef typename K::Triangle_2 Triangle_2; 04255 public: 04256 typedef typename K::Oriented_side result_type; 04257 04258 result_type 04259 operator()( const Circle_2& c, const Point_2& p) const 04260 { return Oriented_side(c.bounded_side(p) * c.orientation()); } 04261 04262 result_type 04263 operator()( const Line_2& l, const Point_2& p) const 04264 { 04265 CGAL_kernel_precondition( ! l.is_degenerate() ); 04266 RT v = l.a()*p.hx() + l.b()*p.hy() + l.c()*p.hw(); 04267 return CGAL_NTS sign(v); 04268 } 04269 04270 result_type 04271 operator()( const Triangle_2& t, const Point_2& p) const 04272 { 04273 typename K::Collinear_are_ordered_along_line_2 04274 collinear_are_ordered_along_line; 04275 typename K::Orientation_2 orientation; 04276 // depends on the orientation of the vertices 04277 Orientation o1 = orientation(t.vertex(0), t.vertex(1), p), 04278 o2 = orientation(t.vertex(1), t.vertex(2), p), 04279 o3 = orientation(t.vertex(2), t.vertex(3), p), 04280 ot = orientation(t.vertex(0), t.vertex(1), t.vertex(2)); 04281 04282 if (o1 == ot && o2 == ot && o3 == ot) // ot cannot be COLLINEAR 04283 return ot; 04284 return 04285 (o1 == COLLINEAR 04286 && collinear_are_ordered_along_line(t.vertex(0), p, t.vertex(1))) || 04287 (o2 == COLLINEAR 04288 && collinear_are_ordered_along_line(t.vertex(1), p, t.vertex(2))) || 04289 (o3 == COLLINEAR 04290 && collinear_are_ordered_along_line(t.vertex(2), p, t.vertex(3))) 04291 ? ON_ORIENTED_BOUNDARY 04292 : -ot; 04293 } 04294 }; 04295 04296 04297 template <typename K> 04298 class Side_of_bounded_circle_2 04299 { 04300 typedef typename K::Point_2 Point_2; 04301 public: 04302 typedef typename K::Bounded_side result_type; 04303 04304 result_type 04305 operator()( const Point_2& p, const Point_2& q, const Point_2& t) const 04306 { 04307 typedef typename K::RT RT; 04308 04309 const RT& phx = p.hx(); 04310 const RT& phy = p.hy(); 04311 const RT& phw = p.hw(); 04312 const RT& qhx = q.hx(); 04313 const RT& qhy = q.hy(); 04314 const RT& qhw = q.hw(); 04315 const RT& thx = t.hx(); 04316 const RT& thy = t.hy(); 04317 const RT& thw = t.hw(); 04318 04319 return enum_cast<Bounded_side>( 04320 CGAL_NTS compare((thx*phw-phx*thw)*(qhx*thw-thx*qhw), 04321 (thy*phw-phy*thw)*(thy*qhw-qhy*thw)) ); 04322 } 04323 04324 result_type 04325 operator()( const Point_2& q, const Point_2& r, 04326 const Point_2& s, const Point_2& t) const 04327 { 04328 typedef typename K::RT RT; 04329 04330 const RT& qhx = q.hx(); 04331 const RT& qhy = q.hy(); 04332 const RT& qhw = q.hw(); 04333 const RT& rhx = r.hx(); 04334 const RT& rhy = r.hy(); 04335 const RT& rhw = r.hw(); 04336 const RT& shx = s.hx(); 04337 const RT& shy = s.hy(); 04338 const RT& shw = s.hw(); 04339 const RT& thx = t.hx(); 04340 const RT& thy = t.hy(); 04341 const RT& thw = t.hw(); 04342 04343 CGAL_kernel_precondition( ! collinear(q,r,s) ); 04344 04345 // compute sign of |qx qy qx^2+qy^2 1 | | a b c d | 04346 // | -- r -- | = | e f g h | 04347 // determinant | -- s -- | = | i j k l | 04348 // | -- t -- | | m n o p | 04349 // where 04350 04351 RT a = qhx*qhw; 04352 RT b = qhy*qhw; 04353 RT c = qhx*qhx + qhy*qhy; 04354 RT d = qhw*qhw; 04355 04356 RT e = rhx*rhw; 04357 RT f = rhy*rhw; 04358 RT g = rhx*rhx + rhy*rhy; 04359 RT h = rhw*rhw; 04360 04361 RT i = shx*shw; 04362 RT j = shy*shw; 04363 RT k = shx*shx + shy*shy; 04364 RT l = shw*shw; 04365 04366 RT m = thx*thw; 04367 RT n = thy*thw; 04368 RT o = thx*thx + thy*thy; 04369 RT p = thw*thw; 04370 04371 RT det = a * ( f*(k*p - l*o) + j*(h*o - g*p) + n*(g*l - h*k) ) 04372 - e * ( b*(k*p - l*o) + j*(d*o - c*p) + n*(c*l - d*k) ) 04373 + i * ( b*(g*p - h*o) + f*(d*o - c*p) + n*(c*h - d*g) ) 04374 - m * ( b*(g*l - h*k) + f*(d*k - c*l) + j*(c*h - d*g) ); 04375 04376 if ( det == 0 ) 04377 return ON_BOUNDARY; 04378 else 04379 { 04380 if (orientation(q,r,s) == CLOCKWISE) 04381 det = -det; 04382 return (0 < det ) ? ON_BOUNDED_SIDE : ON_UNBOUNDED_SIDE; 04383 } 04384 } 04385 }; 04386 04387 template <typename K> 04388 class Side_of_bounded_sphere_3 04389 { 04390 typedef typename K::Point_3 Point_3; 04391 public: 04392 typedef typename K::Bounded_side result_type; 04393 04394 result_type 04395 operator()( const Point_3& p, const Point_3& q, const Point_3& t) const 04396 { 04397 typedef typename K::RT RT; 04398 04399 const RT& phx = p.hx(); 04400 const RT& phy = p.hy(); 04401 const RT& phz = p.hz(); 04402 const RT& phw = p.hw(); 04403 const RT& qhx = q.hx(); 04404 const RT& qhy = q.hy(); 04405 const RT& qhz = q.hz(); 04406 const RT& qhw = q.hw(); 04407 const RT& thx = t.hx(); 04408 const RT& thy = t.hy(); 04409 const RT& thz = t.hz(); 04410 const RT& thw = t.hw(); 04411 04412 return enum_cast<Bounded_side>( 04413 CGAL_NTS sign((thx*phw-phx*thw)*(qhx*thw-thx*qhw) 04414 + (thy*phw-phy*thw)*(qhy*thw-thy*qhw) 04415 + (thz*phw-phz*thw)*(qhz*thw-thz*qhw))); 04416 } 04417 04418 result_type 04419 operator()( const Point_3& p, const Point_3& q, 04420 const Point_3& r, const Point_3& t) const 04421 { 04422 Point_3 center = circumcenter(p, q, r); 04423 return enum_cast<Bounded_side>( compare_distance_to_point(center, p, t) ); 04424 } // FIXME 04425 04426 result_type 04427 operator()( const Point_3& p, const Point_3& q, const Point_3& r, 04428 const Point_3& s, const Point_3& test) const 04429 { 04430 Oriented_side oside = side_of_oriented_sphere(p,q,r,s,test); 04431 if ( are_positive_oriented( p,q,r,s) ) 04432 { 04433 switch (oside) 04434 { 04435 case ON_POSITIVE_SIDE : return ON_BOUNDED_SIDE; 04436 case ON_ORIENTED_BOUNDARY: return ON_BOUNDARY; 04437 case ON_NEGATIVE_SIDE : return ON_UNBOUNDED_SIDE; 04438 } 04439 } 04440 else 04441 { 04442 switch (oside) 04443 { 04444 case ON_POSITIVE_SIDE : return ON_UNBOUNDED_SIDE; 04445 case ON_ORIENTED_BOUNDARY: return ON_BOUNDARY; 04446 case ON_NEGATIVE_SIDE : return ON_BOUNDED_SIDE; 04447 } 04448 } 04449 return ON_BOUNDARY; // Pls, no warnings anylonger 04450 } // FIXME 04451 }; 04452 04453 template <typename K> 04454 class Side_of_oriented_circle_2 04455 { 04456 typedef typename K::Point_2 Point_2; 04457 public: 04458 typedef typename K::Oriented_side result_type; 04459 04460 result_type 04461 operator()( const Point_2& q, const Point_2& r, 04462 const Point_2& s, const Point_2& t) const 04463 { 04464 typedef typename K::RT RT; 04465 04466 const RT& qhx = q.hx(); 04467 const RT& qhy = q.hy(); 04468 const RT& qhw = q.hw(); 04469 const RT& rhx = r.hx(); 04470 const RT& rhy = r.hy(); 04471 const RT& rhw = r.hw(); 04472 const RT& shx = s.hx(); 04473 const RT& shy = s.hy(); 04474 const RT& shw = s.hw(); 04475 const RT& thx = t.hx(); 04476 const RT& thy = t.hy(); 04477 const RT& thw = t.hw(); 04478 04479 CGAL_kernel_precondition( ! collinear(q,r,s) ); 04480 04481 // compute sign of |qx qy qx^2+qy^2 1 | | a b c d | 04482 // | -- r -- | = | e f g h | 04483 // determinant | -- s -- | = | i j k l | 04484 // | -- t -- | | m n o p | 04485 // where 04486 04487 RT a = qhx*qhw; 04488 RT b = qhy*qhw; 04489 RT c = qhx*qhx + qhy*qhy; 04490 RT d = qhw*qhw; 04491 04492 RT e = rhx*rhw; 04493 RT f = rhy*rhw; 04494 RT g = rhx*rhx + rhy*rhy; 04495 RT h = rhw*rhw; 04496 04497 RT i = shx*shw; 04498 RT j = shy*shw; 04499 RT k = shx*shx + shy*shy; 04500 RT l = shw*shw; 04501 04502 RT m = thx*thw; 04503 RT n = thy*thw; 04504 RT o = thx*thx + thy*thy; 04505 RT p = thw*thw; 04506 04507 RT det = a * ( f*(k*p - l*o) + j*(h*o - g*p) + n*(g*l - h*k) ) 04508 - e * ( b*(k*p - l*o) + j*(d*o - c*p) + n*(c*l - d*k) ) 04509 + i * ( b*(g*p - h*o) + f*(d*o - c*p) + n*(c*h - d*g) ) 04510 - m * ( b*(g*l - h*k) + f*(d*k - c*l) + j*(c*h - d*g) ); 04511 04512 return CGAL_NTS sign(det); 04513 } 04514 }; 04515 04516 template <typename K> 04517 class Side_of_oriented_sphere_3 04518 { 04519 typedef typename K::Point_3 Point_3; 04520 public: 04521 typedef typename K::Oriented_side result_type; 04522 04523 result_type 04524 operator()( const Point_3& p, const Point_3& q, const Point_3& r, 04525 const Point_3& s, const Point_3& t) const 04526 { 04527 typedef typename K::RT RT; 04528 CGAL_kernel_precondition( !coplanar(p,q,r,s) ); 04529 const RT & phx = p.hx(); 04530 const RT & phy = p.hy(); 04531 const RT & phz = p.hz(); 04532 const RT & phw = p.hw(); 04533 const RT phw2 = phw*phw; 04534 04535 const RT & qhx = q.hx(); 04536 const RT & qhy = q.hy(); 04537 const RT & qhz = q.hz(); 04538 const RT & qhw = q.hw(); 04539 const RT qhw2 = qhw*qhw; 04540 04541 const RT & rhx = r.hx(); 04542 const RT & rhy = r.hy(); 04543 const RT & rhz = r.hz(); 04544 const RT & rhw = r.hw(); 04545 const RT rhw2 = rhw*rhw; 04546 04547 const RT & shx = s.hx(); 04548 const RT & shy = s.hy(); 04549 const RT & shz = s.hz(); 04550 const RT & shw = s.hw(); 04551 const RT shw2 = shw*shw; 04552 04553 const RT & thx = t.hx(); 04554 const RT & thy = t.hy(); 04555 const RT & thz = t.hz(); 04556 const RT & thw = t.hw(); 04557 const RT thw2 = thw*thw; 04558 04559 const RT det = determinant<RT>( 04560 phx*phw, phy*phw, phz*phw, phx*phx + phy*phy + phz*phz, phw2, 04561 qhx*qhw, qhy*qhw, qhz*qhw, qhx*qhx + qhy*qhy + qhz*qhz, qhw2, 04562 rhx*rhw, rhy*rhw, rhz*rhw, rhx*rhx + rhy*rhy + rhz*rhz, rhw2, 04563 shx*shw, shy*shw, shz*shw, shx*shx + shy*shy + shz*shz, shw2, 04564 thx*thw, thy*thw, thz*thw, thx*thx + thy*thy + thz*thz, thw2); 04565 return - CGAL_NTS sign(det); 04566 } 04567 }; 04568 04569 } // namespace HomogeneousKernelFunctors 04570 04571 CGAL_END_NAMESPACE 04572 04573 #endif // CGAL_HOMOGENEOUS_FUNCTION_OBJECTS_H
1.7.6.1
