diff --git a/Kernel_23/include/CGAL/internal/Projection_traits_3.h b/Kernel_23/include/CGAL/internal/Projection_traits_3.h index 4d152acda3e..8aa46bca30b 100644 --- a/Kernel_23/include/CGAL/internal/Projection_traits_3.h +++ b/Kernel_23/include/CGAL/internal/Projection_traits_3.h @@ -225,6 +225,13 @@ public: return Point_2(x(p),y(p)); } + FT alpha(const Point_2& p, const Point_2& source, const Point_2& target) const + { + FT dx = target.x() - source.x(); + FT dy = target.y() - source.y(); + return (CGAL::abs(dx)>CGAL::abs(dy)) ? ( p.x()-source.x() ) / dx : (p.y()-source.y() ) / dy; + } + Object operator()(const Segment_3& s1, const Segment_3& s2) const { Point_2 s1_source = project(s1.source()); @@ -244,15 +251,15 @@ public: const Segment_2* si=CGAL::object_cast(&o); if (si==NULL) return Object(); FT src[3],tgt[3]; - tgt[dim] = src[dim] = FT(0); //the third coordinate is the midpoint between the points on s1 and s2 - - FT z1 = s1.source()[dim] + ( si->source().x()-s1_source.x() ) / (s1_target.x() - s1_source.x()) * ( s1.target()[dim] - s1.source()[dim] ); - FT z2 = s2.source()[dim] + ( si->source().x()-s2_source.x() ) / (s2_target.x() - s2_source.x()) * ( s2.target()[dim] - s2.source()[dim] ); + //the third coordinate is the midpoint between the points on s1 and s2 + FT z1 = s1.source()[dim] + ( alpha(si->source(), s1_source, s1_target) * ( s1.target()[dim] - s1.source()[dim] )); + FT z2 = s2.source()[dim] + ( alpha(si->source(), s2_source, s2_target) * ( s2.target()[dim] - s2.source()[dim] )); src[dim] = (z1+z2) / FT(2); - z1 = s1.source()[dim] + ( si->target().x()-s1_source.x() ) / (s1_target.x() - s1_source.x()) * ( s1.target()[dim] - s1.source()[dim] ); - z2 = s2.source()[dim] + ( si->target().x()-s2_source.x() ) / (s2_target.x() - s2_source.x()) * ( s2.target()[dim] - s2.source()[dim] ); + z1 = s1.source()[dim] + ( alpha(si->target(), s1_source, s1_target) * ( s1.target()[dim] - s1.source()[dim] )); + z2 = s2.source()[dim] + ( alpha(si->target(), s2_source, s2_target) * ( s2.target()[dim] - s2.source()[dim] )); + tgt[dim] = (z1+z2) / FT(2); @@ -264,8 +271,8 @@ public: } FT coords[3]; //compute the third coordinate of the projected intersection point onto 3D segments - FT z1 = s1.source()[dim] + ( pi->x()-s1_source.x() ) / (s1_target.x() - s1_source.x()) * ( s1.target()[dim] - s1.source()[dim] ); - FT z2 = s2.source()[dim] + ( pi->x()-s2_source.x() ) / (s2_target.x() - s2_source.x()) * ( s2.target()[dim] - s2.source()[dim] ); + FT z1 = s1.source()[dim] + ( alpha(*pi, s1_source, s1_target) * ( s1.target()[dim] - s1.source()[dim] )); + FT z2 = s2.source()[dim] + ( alpha(*pi, s2_source, s2_target) * ( s2.target()[dim] - s2.source()[dim] )); coords[dim] = (z1+z2) / FT(2); coords[Projector::x_index] = pi->x(); diff --git a/Kernel_23/test/Kernel_23/test_Projection_traits_xy_3_Intersect_2.cpp b/Kernel_23/test/Kernel_23/test_Projection_traits_xy_3_Intersect_2.cpp new file mode 100644 index 00000000000..cebb0fe760a --- /dev/null +++ b/Kernel_23/test/Kernel_23/test_Projection_traits_xy_3_Intersect_2.cpp @@ -0,0 +1,41 @@ +#include +#include + +typedef CGAL::Exact_predicates_inexact_constructions_kernel Epik; +typedef CGAL::Projection_traits_xy_3 K; + +typedef K::Orientation_2 Orientation_2; +typedef K::Intersect_2 Intersect_2; + +typedef K::Point_2 Point_2; +typedef K::Segment_2 Segment_2; +typedef CGAL::Object Object; + +int main() +{ + Point_2 p(0,0,0), q(1,1,1), r(1,0,0), s(0,1,1), t(0,1,0), u(0,1,-1); + Point_2 v(0.5, 0, 0), w(0.5,1,0); + + Segment_2 pq(p,q), rs(r,s), rt(r,t), ru(r,u), vw(v,w); + + Point_2 pqrs, pqrt, pqru, pqvw; + + Object o = Intersect_2()(pq,rs); + assert(assign(pqrs,o)); + assert(pqrs == Point_2(0.5, 0.5, 0.5)); + + o = Intersect_2()(pq,rt); + assert(assign(pqrt,o)); + assert(pqrt == Point_2(0.5, 0.5, 0.25)); + + o = Intersect_2()(pq,ru); + assert(assign(pqru,o)); + assert(pqru == Point_2(0.5, 0.5, 0)); + + o = Intersect_2()(pq,vw); + assert(assign(pqvw,o)); + assert(pqvw == Point_2(0.5, 0.5, 0.25)); + + std::cerr << "done" << std::endl; + return 0; +}