Merge pull request #994 from xarthurx/dev

fix segments_intersection bug.

include/igl/segment_segment_intersect.cpp

@@ -11,55 +11,55 @@
 
 template<typename DerivedSource, typename DerivedDir>
 IGL_INLINE bool igl::segments_intersect(
-        const Eigen::PlainObjectBase <DerivedSource> &p,
-        const Eigen::PlainObjectBase <DerivedDir> &r,
-        const Eigen::PlainObjectBase <DerivedSource> &q,
-        const Eigen::PlainObjectBase <DerivedDir> &s,
-        double &a_t,
-        double &a_u,
-        double eps
+  const Eigen::PlainObjectBase <DerivedSource> &p,
+  const Eigen::PlainObjectBase <DerivedDir> &r,
+  const Eigen::PlainObjectBase <DerivedSource> &q,
+  const Eigen::PlainObjectBase <DerivedDir> &s,
+  double &a_t,
+  double &a_u,
+  double eps
 )
 {
-    // http://stackoverflow.com/questions/563198/how-do-you-detect-where-two-line-segments-intersect
-    // Search intersection between two segments
-    // p + t*r :  t \in [0,1]
-    // q + u*s :  u \in [0,1]
+  // http://stackoverflow.com/questions/563198/how-do-you-detect-where-two-line-segments-intersect
+  // Search intersection between two segments
+  // p + t*r :  t \in [0,1]
+  // q + u*s :  u \in [0,1]
 
-    // p + t * r = q + u * s  // x s
-    // t(r x s) = (q - p) x s
-    // t = (q - p) x s / (r x s)
+  // p + t * r = q + u * s  // x s
+  // t(r x s) = (q - p) x s
+  // t = (q - p) x s / (r x s)
 
-    // (r x s) ~ 0 --> directions are parallel, they will never cross
-    Eigen::RowVector3d rxs = r.cross(s);
-    if (rxs.norm() <= eps)
-        return false;
+  // (r x s) ~ 0 --> directions are parallel, they will never cross
+  Eigen::RowVector3d rxs = r.cross(s);
+  if (rxs.norm() <= eps)
+    return false;
 
-    int sign;
+  int sign;
 
-    double u;
-    // u = (q − p) × r / (r × s)
-    Eigen::RowVector3d u1 = (q - p).cross(r);
-    sign = ((u1.dot(rxs)) > 0) ? 1 : -1;
-    u = u1.norm() / rxs.norm();
-    u = u * sign;
+  double u;
+  // u = (q − p) × r / (r × s)
+  Eigen::RowVector3d u1 = (q - p).cross(r);
+  sign = ((u1.dot(rxs)) > 0) ? 1 : -1;
+  u = u1.norm() / rxs.norm();
+  u = u * sign;
 
-    if ((u - 1.) > eps || u < -eps)
-        return false;
+  double t;
+  // t = (q - p) x s / (r x s)
+  Eigen::RowVector3d t1 = (q - p).cross(s);
+  sign = ((t1.dot(rxs)) > 0) ? 1 : -1;
+  t = t1.norm() / rxs.norm();
+  t = t * sign;
 
-    double t;
-    // t = (q - p) x s / (r x s)
-    Eigen::RowVector3d t1 = (q - p).cross(s);
-    sign = ((t1.dot(rxs)) > 0) ? 1 : -1;
-    t = t1.norm() / rxs.norm();
-    t = t * sign;
+  a_t = t;
+  a_u = u;
 
-    if (t < -eps || fabs(t) < eps)
-        return false;
+  if ((u - 1.) > eps || u < -eps)
+    return false;
 
-    a_t = t;
-    a_u = u;
+  if ((t - 1.) > eps || t < -eps)
+    return false;
 
-    return true;
+  return true;
 };
 
 #ifdef IGL_STATIC_LIBRARY

tests/include/igl/segment_segment_intersect.cpp

@@ -0,0 +1,39 @@
+#include <test_common.h>
+#include <igl/segment_segment_intersect.h>
+
+TEST_CASE("segment_segment_intersect: examples", "[igl]")
+{
+  // example 1
+  // https://github.com/libigl/libigl/issues/965
+
+  auto s1 = Eigen::RowVector3d(581, 388, -54);
+  auto dir1 = Eigen::RowVector3d(0.75, -4.36, 0.64);
+
+  auto s2 = Eigen::RowVector3d(636, 77, -10);
+  auto dir2 = Eigen::RowVector3d(-2.79, 0.96, 0.39);
+
+  double a1, a2;
+  bool sect1 = igl::segments_intersect(s1, dir1, s2, dir2, a1, a2);
+
+  bool intersectCondition1 = (a1 >= 0 && a1 <= 1 && a2 >= 0 && a2 <= 1);
+
+  REQUIRE(sect1 == false);
+  REQUIRE(intersectCondition1 == false);
+
+  // example 2
+  // https://github.com/libigl/libigl/issues/957
+
+  auto s3 = Eigen::RowVector3d(-56.6, 0, -201.7);
+  auto dir3 = Eigen::RowVector3d(0, 0, 1);
+
+  auto s4 = Eigen::RowVector3d(-65.6, 0, 258.8);
+  auto dir4 = Eigen::RowVector3d(15.4, 0, 1.9);
+
+  double a3, a4;
+  bool sect2 = igl::segments_intersect(s3, dir3, s4, dir4, a3, a4);
+
+  bool intersectCondition2 = (a3 >= 0 && a3 <= 1 && a4 >= 0 && a4 <= 1);
+
+  REQUIRE(sect2 == false);
+  REQUIRE(intersectCondition2 == false);
+}