diff --git a/FastEnvelope/examples/FastEnvelope/fast.cpp b/FastEnvelope/examples/FastEnvelope/fast.cpp index 9d3bfdf8a75..18781ae18a5 100644 --- a/FastEnvelope/examples/FastEnvelope/fast.cpp +++ b/FastEnvelope/examples/FastEnvelope/fast.cpp @@ -366,16 +366,12 @@ struct Envelope { for (int i = 0; i < cutp.size(); i++){ const Plane& plane_i = prism[cutp[i]]; - CGAL::cpp11::result_of::type - result = CGAL::intersection(line, plane_i.eplane); - if(! result){ + boost::optional op = CGAL::intersection_point(line, plane_i.eplane); + if(! op){ std::cout << "there must be an intersection 2" << std::endl; } - const ePoint_3* ipp = boost::get(&*result); - CGAL_assertion(ipp != nullptr); - - const ePoint_3& ip = *ipp; + const ePoint_3& ip = *op; for(int j = 0; j < cutp.size(); j++) { if (i == j){ @@ -412,18 +408,14 @@ struct Envelope { const std::vector &prismindex, const int &jump, int &id) const { CGAL::Oriented_side ori; - CGAL::cpp11::result_of::type - result = CGAL::intersection(eline, plane.eplane); + boost::optional op = CGAL::intersection_point(eline, plane.eplane); #ifdef TRACE - if(! result){ + if(! op){ std::cout << "there must be an intersection 3" << std::endl; } #endif - const ePoint_3* ipp = boost::get(&*result); - CGAL_assertion(ipp != nullptr); - - const ePoint_3& ip = *ipp; + const ePoint_3& ip = *op; for (int i = 0; i < prismindex.size(); i++){ if (prismindex[i] == jump){ continue; @@ -551,10 +543,16 @@ struct Envelope { */ ) const { + + Point_3 itp(CGAL::approx(tp).x().inf(), CGAL::approx(tp).y().inf(), CGAL::approx(tp).z().inf()); + Point_3 itq(CGAL::approx(tq).x().inf(), CGAL::approx(tq).y().inf(), CGAL::approx(tq).z().inf()); + Point_3 itr(CGAL::approx(tr).x().inf(), CGAL::approx(tr).y().inf(), CGAL::approx(tr).z().inf()); + Point_3 in = itp + CGAL::cross_product((itp - itq), (itp - itr)); + ePoint_3 n(in.x(),in.y(),in.z()); // return 2; // todo: what do we test here with n ? // todo : do this not with Epeck - ePoint_3 n = tp + CGAL::cross_product((tp - tq), (tp - tr)); + // ePoint_3 n = tp + CGAL::cross_product((tp - tq), (tp - tr)); if (CGAL::orientation(n, tp, tq, tr) == 0){ std::cout << "todo degeneration handling" << std::endl; @@ -583,7 +581,13 @@ struct Envelope { const ePoint_3& tr, const ePoint_3& ip) const { - ePoint_3 n = tp + CGAL::cross_product((tp - tq), (tp - tr)); + Point_3 itp(CGAL::approx(tp).x().inf(), CGAL::approx(tp).y().inf(), CGAL::approx(tp).z().inf()); + Point_3 itq(CGAL::approx(tq).x().inf(), CGAL::approx(tq).y().inf(), CGAL::approx(tq).z().inf()); + Point_3 itr(CGAL::approx(tr).x().inf(), CGAL::approx(tr).y().inf(), CGAL::approx(tr).z().inf()); + Point_3 in = itp + CGAL::cross_product((itp - itq), (itp - itr)); + ePoint_3 n(in.x(),in.y(),in.z()); + + //ePoint_3 n = tp + CGAL::cross_product((tp - tq), (tp - tr)); #if 0 if (Predicates::orient_3d(n, tp, q, tr) == 0) { @@ -725,17 +729,14 @@ struct Envelope { // std::cout << *(seg0[k]) << " " << *(seg1[k]) << std::endl; eLine_3 eline(*(seg0[k]), *(seg1[k])); - CGAL::cpp11::result_of::type - result = CGAL::intersection(eline, plane_i.eplane); - if(! result){ + boost::optional op = CGAL::intersection_point(eline, plane_i.eplane); + if(! op){ #ifdef TRACE std::cout << "there must be an intersection 6" << std::endl; #endif } - const ePoint_3* ipp = boost::get(&*result); - CGAL_assertion(ipp != nullptr); - const ePoint_3& ip = *ipp; + const ePoint_3& ip = *op; for (int j = 0; j < cutp.size(); j++){ if (i == j){ @@ -900,18 +901,14 @@ struct Envelope { { eLine_3 eline(segpoint0,segpoint1); // todo replace parameter of function - CGAL::cpp11::result_of::type - result = CGAL::intersection(eline, eplane); - if(! result){ + boost::optional op = CGAL::intersection_point(eline, eplane); + if(! op){ #ifdef TRACE std::cout << "there must be an intersection 9" << std::endl; #endif } - const ePoint_3* ipp = boost::get(&*result); - CGAL_assertion(ipp != nullptr); - - const ePoint_3& ip = *ipp; + const ePoint_3& ip = *op; int tot, fid, ori; for (int i = 0; i < prismindex.size(); i++){ if (prismindex[i] == jump){ @@ -985,18 +982,14 @@ struct Envelope { { eLine_3 eline(segpoint0,segpoint1); // todo replace parameter of function - CGAL::cpp11::result_of::type - result = CGAL::intersection(eline, eplane); - if(! result){ + boost::optional op = CGAL::intersection_point(eline, eplane); + if(! op){ #ifdef TRACE std::cout << "there must be an intersection 10" << std::endl; #endif } - const ePoint_3* ipp = boost::get(&*result); - CGAL_assertion(ipp != nullptr); - - const ePoint_3& ip = *ipp; + const ePoint_3& ip = *op; int tot, ori, fid; @@ -1472,20 +1465,15 @@ struct Envelope { // AF: We moved the intersection here // In case there is no intersection point we continue - CGAL::cpp11::result_of::type - result = CGAL::intersection(etriangle_eplane, - halfspace[jump1][intersect_face[queue[i]][k]].eplane, - halfspace[jump2][intersect_face[queue[j]][h]].eplane); - if(! result){ + boost::optional + op = CGAL::intersection_point(etriangle_eplane, + halfspace[jump1][intersect_face[queue[i]][k]].eplane, + halfspace[jump2][intersect_face[queue[j]][h]].eplane); + if(! op){ continue; } - const ePoint_3* ipp = boost::get(&*result); - if(ipp == nullptr){ - continue; - } - - const ePoint_3& ip = *ipp; + const ePoint_3& ip = *op; cut = is_3_triangle_cut(etriangle[0], etriangle[1], etriangle[2], ip); @@ -1938,17 +1926,17 @@ int main(int argc, char* argv[]) std::ofstream inside("inside.txt"); std::ofstream outside("outside.txt"); - for(int i = 0; i < env_vertices.size() ; i+=10){ + for(int i = 0; i < env_vertices.size() ; i+=10){ for(int j = i+1; j < env_vertices.size(); j+= 10){ for(int k = j+1; k < env_vertices.size(); k+=10){ if( ( i != j) && (i != k) && (j != k)){ if(! CGAL::collinear(env_vertices[i], env_vertices[j],env_vertices[k])){ if(envelope(env_vertices[i], env_vertices[j], env_vertices[k])){ inside_count++; - inside << i << " " << j << " "<< k < P2(CGAL_FE_TONEAREST); - boost::optional oep = ec(CGAL::exact(l1),CGAL::exact(l2),CGAL::exact(l3)); + boost::optional oep = ec(CGAL::exact(l1), CGAL::exact(l2), CGAL::exact(l3)); if(oep == boost::none){ return boost::none; } - typedef Lazy_rep_0 LazyPointRep; + typedef Lazy_rep_0 LazyPointRep; + const typename EK::Point_3 ep = *oep; + LazyPointRep *rep = new LazyPointRep(ep); + typename LK::Point_3 lp(rep); + return boost::make_optional(lp); + } + // AF can we get here?? + return boost::none; + } + + // for Intersect_point_3 with Plane_3 Line_3 + template + result_type operator()(const L1& l1, const L2& l2) const + { + Protect_FPU_rounding P; + + try { + boost::optional oap = ac(CGAL::approx(l1),CGAL::approx(l2)); + if(oap == boost::none){ + return boost::none; + } + // Now we have to construct a rep for a lazy point with the three lazy planes. + typedef Lazy_rep_optional_n LazyPointRep; + const typename AK::Point_3 ap = *oap; + LazyPointRep *rep = new LazyPointRep(ap, ec, l1, l2); + typename LK::Point_3 lp(rep); + return boost::make_optional(lp); + + + } catch (Uncertain_conversion_exception&) { + Protect_FPU_rounding P2(CGAL_FE_TONEAREST); + boost::optional oep = ec(CGAL::exact(l1), CGAL::exact(l2)); + if(oep == boost::none){ + return boost::none; + } + typedef Lazy_rep_0 LazyPointRep; const typename EK::Point_3 ep = *oep; LazyPointRep *rep = new LazyPointRep(ep); typename LK::Point_3 lp(rep); diff --git a/Intersections_3/include/CGAL/Intersections_3/Line_3_Plane_3.h b/Intersections_3/include/CGAL/Intersections_3/Line_3_Plane_3.h index 0eebb3eb50f..5dc93e65edd 100644 --- a/Intersections_3/include/CGAL/Intersections_3/Line_3_Plane_3.h +++ b/Intersections_3/include/CGAL/Intersections_3/Line_3_Plane_3.h @@ -22,6 +22,22 @@ namespace CGAL { CGAL_INTERSECTION_FUNCTION(Plane_3, Line_3, 3) CGAL_DO_INTERSECT_FUNCTION(Plane_3, Line_3, 3) + +template < class K > +inline +boost::optional +intersection_point(const Plane_3& plane, const Line_3& line) +{ + return K().intersect_point_3_object()(plane, line); +} + + template < class K > +inline +boost::optional + intersection_point(const Line_3& line, const Plane_3& plane) +{ + return K().intersect_point_3_object()(plane, line); +} } #endif // CGAL_INTERSECTIONS_3_LINE_PLANE_3_H diff --git a/Intersections_3/include/CGAL/Intersections_3/internal/intersection_3_1_impl.h b/Intersections_3/include/CGAL/Intersections_3/internal/intersection_3_1_impl.h index c0822ae73a0..8206961896d 100644 --- a/Intersections_3/include/CGAL/Intersections_3/internal/intersection_3_1_impl.h +++ b/Intersections_3/include/CGAL/Intersections_3/internal/intersection_3_1_impl.h @@ -79,6 +79,33 @@ do_intersect(const Plane_3 &plane1, const Plane_3 &plane2, namespace Intersections { namespace internal { +template +boost::optional +intersection_point(const typename K::Plane_3 &plane, + const typename K::Line_3 &line, + const K& /*k*/) +{ + typedef typename K::Point_3 Point_3; + typedef typename K::Direction_3 Direction_3; + typedef typename K::RT RT; + + const Point_3 &line_pt = line.point(); + const Direction_3 &line_dir = line.direction(); + + RT num = plane.a()*line_pt.hx() + plane.b()*line_pt.hy() + + plane.c()*line_pt.hz() + wmult_hw((K*)0, plane.d(), line_pt); + RT den = plane.a()*line_dir.dx() + plane.b()*line_dir.dy() + + plane.c()*line_dir.dz(); + if (den == 0) { + return boost::none; + } + return boost::make_optional(Point_3(den*line_pt.hx()-num*line_dir.dx(), + den*line_pt.hy()-num*line_dir.dy(), + den*line_pt.hz()-num*line_dir.dz(), + wmult_hw((K*)0, den, line_pt))); +} + + template typename Intersection_traits::result_type intersection(const typename K::Plane_3 &plane, diff --git a/Kernel_23/include/CGAL/Kernel/function_objects.h b/Kernel_23/include/CGAL/Kernel/function_objects.h index 110be9f6eb8..a788522379e 100644 --- a/Kernel_23/include/CGAL/Kernel/function_objects.h +++ b/Kernel_23/include/CGAL/Kernel/function_objects.h @@ -3573,6 +3573,7 @@ namespace CommonKernelFunctors { { public: typedef typename K::Point_3 Point_3; + typedef typename K::Line_3 Line_3; typedef typename K::Plane_3 Plane_3; typedef typename boost::optional result_type; @@ -3581,6 +3582,12 @@ namespace CommonKernelFunctors { { return Intersections::internal::intersection_point(pl1, pl2, pl3, K() ); } + + result_type + operator()(const Plane_3& plane, const Line_3& line) const + { + return Intersections::internal::intersection_point(plane, line, K() ); + } };