Commit 7ee39614 authored by Maik Froechtenicht's avatar Maik Froechtenicht

Fix for vtkTriangle::TrianglesIntersect() Issue #17092

parent ddda0845
Pipeline #99186 running with stage
......@@ -984,7 +984,7 @@ double Determinant(double a[3], double b[3], double c[3], double d[3])
c[0]-d[0], c[1]-d[1], c[2]-d[2] );
}
static const double eps = 256 * std::numeric_limits<double>::epsilon();
static const double eps = 256*std::numeric_limits<double>::epsilon();
// The orientation values are chosen so that any combination of 3 will produce
// a unique value.
......@@ -1000,9 +1000,9 @@ int Orientation(const double p1[2], const double p2[2], const double p3[2])
// Return 4 if the path connecting p1,p2,p3 is counterclockwise.
// Return 2 if the path connecting p1,p2,p3 is clockwise.
// Return 1 if the points are colinear.
double v1[2] = { p2[0] - p1[0], p2[1] - p1[1] };
double v2[2] = { p3[0] - p1[0], p3[1] - p1[1] };
double signedArea = v1[0] * v2[1] - v1[1] * v2[0];
double v1[2] = { p2[0]-p1[0], p2[1]-p1[1] };
double v2[2] = { p3[0]-p1[0], p3[1]-p1[1] };
double signedArea = v1[0]*v2[1] - v1[1]*v2[0];
if ( std::abs( signedArea ) < eps )
{
return Colinear;
......@@ -1013,8 +1013,8 @@ int Orientation(const double p1[2], const double p2[2], const double p3[2])
int CoplanarTrianglesIntersect(double p1[2], double q1[2], double r1[2],
double p2[2], double q2[2], double r2[2])
{
// Determine whether or not triangle T1 = (p1, q1, r1) intersects triangle
// T2 = (p2, q2, r2), assumming that they are coplanar. This method is adapted
// Determine whether or not triangle T1 = (p1,q1,r1) intersects triangle
// T2 = (p2,q2,r2), assumming that they are coplanar. This method is adapted
// from Olivier Devillers, Philippe Guigue. Faster Triangle-Triangle
// Intersection Tests. RR-4488, IN-RIA. 2002. <inria-00072100>
......@@ -1022,11 +1022,11 @@ int CoplanarTrianglesIntersect(double p1[2], double q1[2], double r1[2],
// counterclockwise.
if ( Orientation( p1, q1, r1 ) == Clockwise )
{
std::swap(q1, r1);
std::swap(q1,r1);
}
if ( Orientation( p2, q2, r2 ) == Clockwise )
{
std::swap(q2, r2);
std::swap(q2,r2);
}
// Next, we compute the orientation of p1 w.r.t. the edges that comprise T2
......@@ -1041,12 +1041,12 @@ int CoplanarTrianglesIntersect(double p1[2], double q1[2], double r1[2],
// then p1 lies on an edge of T2.
int sumOfSigns = p1Orientation[0] + p1Orientation[1] + p1Orientation[2];
static const int Three_CounterClockwise = 3 * Counterclockwise;
static const int Two_Colinear_One_Clockwise = 2 * Colinear + Clockwise;
static const int Two_Colinear_One_Counterclockwise = (2 * Colinear +
static const int Three_CounterClockwise = 3*Counterclockwise;
static const int Two_Colinear_One_Clockwise = 2*Colinear + Clockwise;
static const int Two_Colinear_One_Counterclockwise = (2*Colinear +
Counterclockwise);
static const int One_Colinear_Two_Counterclockwise = (Colinear +
2 * Counterclockwise);
2*Counterclockwise);
if (sumOfSigns == Three_CounterClockwise || // condition 1
sumOfSigns == Two_Colinear_One_Clockwise || // condition 2
......@@ -1069,7 +1069,7 @@ int CoplanarTrianglesIntersect(double p1[2], double q1[2], double r1[2],
for (index = 0; index < 3; index++)
{
if (p1Orientation[index] == Counterclockwise &&
p1Orientation[(index + 2) % 3] == Clockwise)
p1Orientation[(index+2)%3] == Clockwise)
{
break;
}
......@@ -1080,13 +1080,13 @@ int CoplanarTrianglesIntersect(double p1[2], double q1[2], double r1[2],
return 0;
}
double* T2[3] = {p2, q2, r2};
double* T2[3] = {p2,q2,r2};
p2 = T2[index];
q2 = T2[(index + 1) % 3];
r2 = T2[(index + 2) % 3];
q2 = T2[(index+1)%3];
r2 = T2[(index+2)%3];
// First decision tree (p1 belongs to region R1)
if (p1Orientation[(index + 1) % 3] == Counterclockwise)
if (p1Orientation[(index+1)%3] == Counterclockwise)
{
if (Orientation( r2, p2, q1 ) != Clockwise) // Test I
{
......@@ -1174,7 +1174,7 @@ int CoplanarTrianglesIntersect(double p1[2], double q1[2], double r1[2],
}
else
{
// The paper has an error here.
// The paper has an error here.
// Paper: if (Orientation( r2, p2, r1 ) == Clockwise) // Test V.a
// Fix 1: if (Orientation( p1, p2, r1 ) == Clockwise) // Test V.a
// Fix 2:
......@@ -1191,7 +1191,7 @@ int CoplanarTrianglesIntersect(double p1[2], double q1[2], double r1[2],
}
else
{
// The paper has an error here: q1 is in Region R25 when (p1, q2, q1) is
// the paper has an error here: q1 is in Region R25 when (p1,q2,q1) is
// clockwise.
if (Orientation( p1, q2, q1 ) != Counterclockwise) // Test III.b
{
......@@ -1261,15 +1261,14 @@ int CoplanarTrianglesIntersect(double p1[2], double q1[2], double r1[2],
}
// Determine whether or not triangle (p1, q1, r1) intersects triangle
// (p2, q2, r2). This method is adapted from Olivier Devillers, Philippe Guigue.
// Faster Triangle-Triangle Intersection Tests. RR-4488, IN-RIA. 2002.
// <inria-00072100>
// Determine whether or not triangle (p1,q1,r1) intersects triangle (p2,q2,r2).
// This method is adapted from Olivier Devillers, Philippe Guigue. Faster
// Triangle-Triangle Intersection Tests. RR-4488, IN-RIA. 2002. <inria-00072100>
int vtkTriangle::TrianglesIntersect(double p1[3], double q1[3], double r1[3],
double p2[3], double q2[3], double r2[3])
{
// Triangle T1 = (p1, q1, r1) and lies in plane Pi1
// Triangle T2 = (p2, q2, r2) and lies in plane Pi2
// Triangle T1 = (p1,q1,r1) and lies in plane Pi1
// Triangle T2 = (p2,q2,r2) and lies in plane Pi2
// First, we determine whether T1 intersects Pi2
double det1[3] = { Determinant( p2, q2, r2, p1 ),
......@@ -1282,7 +1281,7 @@ int vtkTriangle::TrianglesIntersect(double p1[3], double q1[3], double r1[3],
// The triangles are coplanar. We pick the Cartesian principal plane that
// maximizes their projected area and perform the query in 2-D.
double v1[3], v2[3];
for (int i = 0; i < 3; i++)
for (int i=0;i<3;i++)
{
v1[i] = q1[i] - p1[i];
v2[i] = r1[i] - p1[i];
......@@ -1291,7 +1290,7 @@ int vtkTriangle::TrianglesIntersect(double p1[3], double q1[3], double r1[3],
vtkMath::Cross( v1, v2, normal );
int index = 0;
for (int i = 1; i < 3; i++)
for (int i=1;i<3;i++)
{
if (std::abs( normal[index] ) < std::abs( normal[i] ))
{
......@@ -1304,7 +1303,7 @@ int vtkTriangle::TrianglesIntersect(double p1[3], double q1[3], double r1[3],
return CoplanarTrianglesIntersect(&p1[1], &q1[1], &r1[1],
&p2[1], &q2[1], &r2[1]);
}
else if (index == 1)
if (index == 1)
{
double p1_[2] = { p1[0], p1[2] };
double q1_[2] = { q1[0], q1[2] };
......@@ -1322,13 +1321,13 @@ int vtkTriangle::TrianglesIntersect(double p1[3], double q1[3], double r1[3],
}
bool degenerate = false;
double* points[3] = {p1, q1, r1};
for (int i = 0; i < 3; i++)
double* points[3] = {p1,q1,r1};
for (int i=0;i<3;i++)
{
if (std::abs( det1[i] ) < eps)
{
degenerate = true;
if (vtkTriangle::PointInTriangle(points[i], p2, q2, r2, eps))
if (vtkTriangle::PointInTriangle(points[i], p2,q2,r2, eps))
{
return 1;
}
......@@ -1374,8 +1373,7 @@ int vtkTriangle::TrianglesIntersect(double p1[3], double q1[3], double r1[3],
int index1;
for (index1 = 0; index1 < 3; index1++)
{
int sumOfSigns = ( det1[(index1 + 1) % 3] > 0. ) +
( det1[(index1 + 2) % 3] > 0. );
int sumOfSigns = ( det1[(index1+1)%3] > 0. ) + ( det1[(index1+2)%3] > 0. );
if (sumOfSigns != 1)
{
break;
......@@ -1383,17 +1381,16 @@ int vtkTriangle::TrianglesIntersect(double p1[3], double q1[3], double r1[3],
}
assert(index1 >= 0 && index1 < 3);
double* T1[3] = {p1, q1, r1};
double* T1[3] = {p1,q1,r1};
p1 = T1[index1];
q1 = T1[(index1 + 1) % 3];
r1 = T1[(index1 + 2) % 3];
q1 = T1[(index1+1)%3];
r1 = T1[(index1+2)%3];
bool swap1 = (det1[index1] < -eps);
int index2;
for (index2 = 0; index2 < 3; index2++)
{
int sumOfSigns = ( det2[(index2 + 1) % 3] > 0. ) +
( det2[(index2 + 2) % 3] > 0. );
int sumOfSigns = ( det2[(index2+1)%3] > 0. ) + ( det2[(index2+2)%3] > 0. );
if (sumOfSigns != 1)
{
break;
......@@ -1401,26 +1398,26 @@ int vtkTriangle::TrianglesIntersect(double p1[3], double q1[3], double r1[3],
}
assert(index2 >= 0 && index2 < 3);
double* T2[3] = {p2, q2, r2};
double* T2[3] = {p2,q2,r2};
p2 = T2[index2];
q2 = T2[(index2 + 1) % 3];
r2 = T2[(index2 + 2) % 3];
q2 = T2[(index2+1)%3];
r2 = T2[(index2+2)%3];
bool swap2 = (det2[index2] < -eps);
if (swap1)
{
std::swap(q2, r2);
std::swap(q2,r2);
}
if (swap2)
{
std::swap(q1, r1);
std::swap(q1,r1);
}
// The final step is to determine whether or not the line segments formed by
// the intersection of T1 and Pi2 and the intersection of T2 and Pi1 overlap.
// This is done by checking the following predicate:
// Determinant(p1, q1, p2, q2) <= 0. ^ Determinant(p1, r1, r2, p2) <= 0.
// Determinant(p1,q1,p2,q2) <= 0. ^ Determinant(p1,r1,r2,p2) <= 0.
if ((Determinant( p1, q1, p2, q2 ) <= 0.) &&
(Determinant( p1, r1, r2, p2 ) <= 0.))
{
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment