Schneider P., Eberly D.H. Geometric Tools for Computer Graphics
Подождите немного. Документ загружается.
614 Chapter 11 Intersection in 3D
if (max1 < min0 || max0 < min1)
return false;
}
}
} else {
// polygons are parallel (coplanar, C0.normal did not separate)
// test edge normals for C0
for (i = 0; i < C0.M; i++) {
D = Cross(C0.normal, C0.E(i));
ComputeInterval(C0, D, min0, max0);
ComputeInterval(C1, D, min1, max1);
if (max1 < min0 || max0 < min1)
return false;
}
// test edge normals for C1
for (i1 = 0; i1 < 3; i1++) {
D = Cross(C1.normal, C1.E(i));
ComputeInterval(C0, D, min0, max0);
ComputeInterval(C1, D, min1, max1);
if (max1 < min0 || max0 < min1)
return false;
}
}
return true;
}
In the presence of a floating-point arithmetic system, the comparison to the zero
vector of the cross product of polygon normals should be modified to use relative
error. The selected threshold
epsilon is a threshold on the value of the squared
sine of the angle between the two normal vectors, based on the identity n
0
×n
1
=
n
0
n
1
sin(θ),whereθ is the angle between the two vectors. The choice of epsilon
is at the discretion of the application.
// comparison when the normals are unit length
Point N0xN1 = Cross(C0.normal, C1.normal);
float N0xN1SqrLen = Dot(N0xN1, N0xN1);
if (N0xN1SqrLen >= epsilon) {
// polygons are (effectively) not parallel
} else {
// polygons are (effectively) parallel
}
11.11 The Method of Separating Axes 615
// comparison when the normals are not unit length
Point N0xN1 = Cross(C0.normal, C1.normal);
float N0xN1SqrLen = Dot(N0xN1, N0xN1);
float N0SqrLen = Dot(C0.normal, C0.normal);
float N1SqrLen = Dot(C1.normal, C1.normal);
if (N0xN1SqrLen >= epsilon * N0SqrLen * N1SqrLen) {
// polygons are (effectively) not parallel
} else {
// polygons are (effectively) parallel
}
11.11.2 Separation of Moving Convex Polyhedra
The structure of this algorithm is similar to that of the two-dimensional problem. See
Section 7.7 for the details. The pseudocode for testing for intersection of two convex
polyhedra of positive volume is listed below.
bool TestIntersection(ConvexPolyhedron C0, Point W0, ConvexPolyhedron C1,
Point W1, float tmax, float& tfirst, float& tlast)
{
W = W1 - W0; // process as if C0 stationary, C1 moving
tfirst = 0; tlast = INFINITY;
// test faces of C0 for separation
for (i = 0; i < C0.L; i++) {
D = C0.F(i).normal;
ComputeInterval(C0, D, min0, max0);
ComputeInterval(C1, D, min1, max1);
speed = Dot(D, W);
if (NoIntersect(tmax, speed, min0, max0, min1, max1, tfirst, tlast))
return false;
}
// test faces of C1 for separation
for (j = 0; j < C1.L; j++) {
D = C1.F(j).normal;
ComputeInterval(C0, D, min0, max0);
ComputeInterval(C1, D, min1, max1);
speed = Dot(D, W);
if (NoIntersect(tmax, speed, min0, max0, min1, max1,
tfirst, tlast))
return false;
}
616 Chapter 11 Intersection in 3D
// test cross products of pairs of edges
for (i = 0; i < C0.M; i++) {
for (j = 0; j < C1.M; j++) {
D = Cross(C0.E(i), C1.E(j));
ComputeInterval(C0, D, min0, max0);
ComputeInterval(C1, D, min1, max1);
speed = Dot(D, W);
if (NoIntersect(tmax, speed, min0, max0, min1, max1, tfirst, tlast) )
return false;
}
}
return true;
}
The function NoIntersect is exactly the one used in the two-dimensional prob-
lem. In the case of two convex polyhedra that are planar polygons, the full pseudocode
is not provided here but can be obtained by replacing each block of the form
if (max1 < min0 || max0 < min1)
return false;
by a block of the form
speed = Dot(D, W);
if (NoIntersect(tmax, speed, min0, max0, min1, max1, tfirst, tlast))
return false;
11.11.3 Intersection Set for Stationary Convex Polyhedra
The find-intersection query for two stationary convex polyhedra is a special example
of Boolean operations on polyhedra. Section 13.6 provides a general discussion for
computing Boolean operations, in particular for computing the intersection of two
polyhedra.
11.11.4 Contact Set for Moving Convex Polyhedra
Given two moving convex objects C
0
and C
1
, initially not intersecting, with velocities
w
0
and w
1
,ifT>0 is the first time of contact, the sets C
0
+T w
0
={X + T w
0
: X ∈
C
0
}and C
1
+T w
1
={X +T w
1
: X ∈C
1
}are just touching with no interpenetration.
As indicated earlier for convex polyhedra, the contact is one of face-face, face-edge,
face-vertex, edge-edge, edge-vertex, or vertex-vertex. The analysis is slightly more
11.11 The Method of Separating Axes 617
complicated than that of the 2D setting, but the ideas are the same—the relative ori-
entation of the convex polyhedra to each other must be known to properly compute
the contact set. We recommend reading Section 7.7 first to understand those ideas
before finishing this section.
The
TestIntersection function can be modified to keep track of which vertices,
edges, or faces are projected to the end points of the projection interval. At the
first time of contact, this information is used to determine how the two objects are
oriented with respect to each other. If the contact is vertex-vertex, vertex-edge, or
vertex-face, then the contact point is a single point, a vertex. If the contact is edge-
edge, the contact is typically a single point, but can be an entire line segment. If the
contact is edge-face, the contact set is a line segment. Finally, if the contact is face-
face, the intersection set is a convex polygon. This is the most complicated scenario
and requires a two-dimensional convex polygon intersector. Each end point of the
projection interval is generated by either a vertex, an edge, or a face. Similar to
the implementation for the two-dimensional problem, a two-character label can be
associated with each polyhedron to indicate the projection type. The single character
labels are
V for a vertex projection, E for an edge projection, and F for a face projection.
The nine two-character labels are
VV, VE, VF, EV, EE, EF, FV, FE, and FF. The first letter
corresponds to the minimum of the interval, and the second letter corresponds to
the maximum of the interval. A convenient data structure for storing the labels, the
projection interval, and the indices of the polyhedron components that project to the
interval end points is exactly the one used in the two-dimensional problem:
Configuration
{
float min, max;
int index[2];
char type[2];
}
The projection interval is [min, max]. As an example, if the projection type is VF,
index[0] is the index of the vertex that projects to the minimum and index[1] is the
index of the face that projects to the maximum.
Two configuration objects are delcared,
Cfg0 for polyhedron C
0
and Cfg1 for
polyhedron C
1
. In the two-dimensional problem, the block of code for each poten-
tial separating direction
d had the
ComputeInterval calls replaced by ProjectNormal
and ProjectGeneral. The function ProjectNormal knows that an edge normal is being
used and only the other extreme projection must be calculated. The function
Pro-
jectGeneral
determines the projection type at both interval extremes. The analogous
functions must be implemented for the three-dimensional problem, but they apply
to the separation tests involving face normals. The separation tests for cross products
of pairs of edges require
ProjectGeneral to replace both calls to ComputeGeneral since
it is not guaranteed in these cases that only the edges project to interval extremes.
618 Chapter 11 Intersection in 3D
The pseudocode for the projection of a polyhedron onto a normal line for one of
its own faces is listed below. The code was simpler in the two-dimensional problem
because of the linear ordering of the vertices and edges. When the current minimum
projection occurs twice, at that point you know that the value is the true minimum
and that an edge must project to it. In the three-dimensional case, there is no linear
ordering, so the code is more complicated. If the current minimum occurs three
times, then the vertices must be part of a face that is an extreme face along the current
direction under consideration. At this point you can return safely knowing that no
other vertices need to be processed. The complication is that once you have found
a third vertex that projects to the current minimum, you must determine which
face contains those vertices. Also possible is that the current minimum occurs twice
and is generated by an edge perpendicular to the current direction, but that edge is
not necessarily an extreme one with respect to the direction. An early exit from the
function when the current minimum occurs twice is not possible.
void ProjectNormal(ConvexPolyhedron C, Point D, int faceindex,
Configuration Cfg)
{
// store the vertex indices that map to current minimum
int minquantity, minindex[3];
// project the face
minquantity = 1;
minindex[0] = C.IndexOf(C.F(faceindex).V(0));
Cfg.min = Cfg.max = Dot(D,C.V(minindex[0]));
Cfg.index[1] = faceindex;
Cfg.type[1] = ’F’;
for (i = 0; i < C.N; i++) {
value = Dot(D, C.V(i));
if (value < Cfg.min) {
minquantity = 1;
minindex[0] = i;
Cfg.min = value;
Cfg.index[0] = i;
Cfg.type[0] = ‘V’;
} else if (value == Cfg.min && Cfg.min < Cfg.max) {
minindex[minquantity++] = i;
if (minquantity == 2) {
if (C.ExistsEdge(minindex[0], minindex[1])) {
// edge is parallel to initial face
Cfg.index[0] = C.GetEdgeIndexFromVertices(minindex);
Cfg.type[0] = ‘E’;
}
11.11 The Method of Separating Axes 619
// else: two nonconnected vertices project to current
// minimum, the first vertex is kept as the current extreme
} else if ( minquantity ==3){
// Face is parallel to initial face. This face must project
// to the minimum and no other vertices can do so. No need
// to further process vertices.
Cfg.index[0] = C.GetFaceIndexFromVertices(minindex);
Cfg.type[0] = ’F’;
return;
}
}
}
}
The pseudocode for general projection of a polyhedron onto the line is listed
below.
void ProjectGeneral(ConvexPolyhedron C, Point D, Configuration Cfg)
{
Cfg.min = Cfg.max = Dot(D, C.V(0));
Cfg.index[0] = Cfg.index[1] = 0;
for (i = 1; i < C.N; i++) {
value = Dot(D, C.V(i));
if (value < Cfg.min) {
Cfg.min = value;
Cfg.index[0] = i;
} else if (value > Cfg.max) {
Cfg.max = value;
Cfg.index[1] = i;
}
}
Cfg.type[0] = Cfg.type[1] = ‘V’;
for(i=0;i<2;i++) {
for each face F sharing C.V(Cfg.index[i]) {
if (F.normal parallel to D) {
Cfg.index[i] = C.GetIndexOfFace(F);
Cfg.type[i] = ’F’;
break;
}
}
if (Cfg.type[i] != ’F’) {
620 Chapter 11 Intersection in 3D
for each edge E sharing C.V(Cfg.index[i]) {
if (E perpendicular to D) {
Cfg.index[i] = C.GetIndexOfEdge(E);
Cfg.type[i] = ‘E’;
break;
}
}
}
}
}
The NoIntersect function that was modified in two dimensions to accept con-
figuration objects instead of projection intervals is used exactly as is for the three-
dimensional problem. With all such modifications,
TestIntersection has the equiv-
alent formulation:
bool TestIntersection(ConvexPolyhedron C0, Point W0, ConvexPolyhedron C1, Point W1,
float tmax, float& tfirst, float& tlast)
{
W = W1 - W0; // process as if C0 stationary, C1 moving
tfirst = 0; tlast = INFINITY;
// test faces of C0 for separation
for (i = 0; i < C0.L; i++) {
D = C0.F(i).normal;
ProjectNormal(C0, D, i, Cfg0);
ProjectGeneral(C1, D, Cfg1);
speed = Dot(D, W);
if (NoIntersect(tmax, speed, Cfg0, Cfg1, Curr0, Curr1, side, tfirst, tlast) )
return false;
}
// test faces of C1 for separation
for (j = 0; j < C1.L; j++) {
D = C1.F(j).normal;
ProjectNormal(C1, D, j, Cfg1);
ProjectGeneral(C0, D, Cfg0);
speed = Dot(D, W);
if (NoIntersect(tmax, speed, Cfg0, Cfg1, Curr0, Curr1, side, tfirst, tlast))
return false;
}
// test cross products of pairs of edges
for (i = 0; i < C0.M; i++) {
11.11 The Method of Separating Axes 621
for (j = 0; j < C1.M; j++) {
D = Cross(C0.E(i), C1.E(j));
ProjectGeneral(C0, D, Cfg0);
ProjectGeneral(C1, D, Cfg1);
speed = Dot(D, W);
if (NoIntersect(tmax, speed, Cfg0, Cfg1, Curr0, Curr1, side, tfirst, t last))
return false;
}
}
return true;
}
The FindIntersection pseudocode has exactly the same implementation as
TestIntersection, but with one additional block of code (after all the loops) that is
reached if there will be an intersection. When the polyhedra intersect at time T , they
are effectively moved with their respective velocities, and the contact set is calcu-
lated. The pseudocode is shown below. The intersection is a convex polyhedron and
is returned in the last argument of the function. If the intersection set is nonempty,
the return value is
true. Otherwise, the original moving convex polyhedra do not
intersect, and the function returns
false.
bool FindIntersection(ConvexPolyhedron C0, Point W0, ConvexPolyhedron C1, Point W1,
float tmax, float& tfirst, float& tlast, ConvexPolyhedron& I)
{
W = W1 - W0; // process as if C0 stationary, C1 moving
tfirst = 0; tlast = INFINITY;
// test faces of C0 for separation
for (i = 0; i < C0.L; i++) {
D = C0.F(i).normal;
ProjectNormal(C0, D, i,Cfg0);
ProjectGeneral(C1, D, Cfg1);
speed = Dot(D, W);
if (NoIntersect(tmax, speed, Cfg0, Cfg1, Curr0, Curr1, side, tfirst, tlast))
return false;
}
// test faces of C1 for separation
for (j = 0; j < C1.L; j++) {
D = C1.F(j).normal;
ProjectNormal(C1, D, j, Cfg1);
ProjectGeneral(C0, D, Cfg0);
speed = Dot(D, W);
622 Chapter 11 Intersection in 3D
if (NoIntersect(tmax, speed, Cfg0, Cfg1, Curr0, Curr1, side, tfirst, tlast))
return false;
}
// test cross products of pairs of edges
for (i = 0; i < C0.M; i++) {
for (j = 0; j < C1.M; j++) {
D = Cross(C0.E(i), C1.E(j));
ProjectGeneral(C0, D, Cfg0);
ProjectGeneral(C1, D, Cfg1);
speed = Dot(D, W);
if (NoIntersect(tmax, speed, Cfg0, Cfg1, Curr0, Curr1, side, tfirst, tlast))
return false;
}
}
// compute the contact set
GetIntersection(C0, W0, C1, W1, Curr0, Curr1, side, tfirst, I);
return true;
}
The intersection calculator pseudocode is shown below.
void GetIntersection(ConvexPolyhedron C0, Point W0, ConvexPolyhedron C1,
Point W1, Configuration Cfg0, Configuration Cfg1, int side, float tfirst,
ConvexPolyhedron I)
{
if (side == 1) {
// C0-max meets C1-min
if (Cfg0.type[1] == ‘V’) {
// vertex-{vertex/edge/face} intersection
I.InsertVertex(C0.V(Cfg0.index[1]) + tfirst * W0);
} else if (Cfg1.type[0] == ‘V’) {
// {vertex/edge/face}-vertex intersection
I.InsertVertex(C1.V(Cfg1.index[0]) + tfirst * W1);
} else if (Cfg0.type[1] == ‘E’) {
Segment E0Moved = C0.E(Cfg0.index[1]) + tfirst * W0;
if (Cfg1.type[0] == ’E’) {
Segment E1Moved = C1.E(Cfg1.index[0]) + tfirst * W1;
FindSegmentIntersection(E0Moved, E1Moved, I);
} else {
ConvexPolygon F1Moved = C1.F(Cfg1.index[0]) + tfirst * W1;
FindSegmentPolygonIntersection(E0Moved, F1Moved, I);
}
11.11 The Method of Separating Axes 623
} else if (Cfg1.type[0] == ‘E’) {
Segment E1Moved = C1.E(Cfg1.index[0]) + tfirst * W1;
if (Cfg0.type[1] == ’E’) {
Segment E0Moved = C0.E(Cfg0.index[1]) + tfirst * W0;
FindSegmentIntersection(E1Moved, E0Moved, I);
} else {
ConvexPolygon F0Moved = C0.F(Cfg0.index[1]) + tfirst * W0;
FindSegmentPolygonIntersection(E1Moved, F0Moved, I);
}
} else {
// Cfg0.type[1] == ‘F‘ && Cfg1.type[0] == ‘F’
// face-face intersection
ConvexPolygon F0Moved = C0.F(Cfg0.index[1]) + tfirst * W0;
ConvexPolygon F1Moved = C1.F(Cfg1.index[0]) + tfirst * W1;
FindPolygonIntersection(F0Moved, F1Moved, I);
}
} else if ( side == -1 ) {
// C1-max meets C0-min
if (Cfg1.type[1] == ‘V’) {
// vertex-{vertex/edge/face} intersection
I.InsertVertex(C1.V(Cfg1.index[1]) + tfirst * W1);
} else if (Cfg0.type[0] == ‘V’) {
// {vertex/edge/face}-vertex intersection
I.InsertVertex(C0.V(Cfg0.index[0]) + tfirst * W0);
} else if (Cfg1.type[1] == ‘E’) {
Segment E1Moved = C1.E(Cfg1.index[1]) + tfirst * W1;
if (Cfg0.type[0] == ’E’) {
Segment E0Moved = C0.E(Cfg0.index[0]) + tfirst * W0;
FindSegmentIntersection(E1Moved, E0Moved, I);
} else {
ConvexPolygon F0Moved = C0.F(Cfg0.index[0]) + tfirst * W0;
FindSegmentPolygonIntersection(E1Moved, F0Moved, I);
}
} else if (Cfg0.type[0] == ‘E’) {
Segment E0Moved = C0.E(Cfg0.index[0]) + tfirst * W0;
if (Cfg1.type[1] == ’E’) {
Segment E1Moved = C1.E(Cfg1.index[1]) + tfirst * W1;
FindSegmentIntersection(E0Moved, E1Moved, I);
} else {
ConvexPolygon F1Moved = C1.F(Cfg1.index[1]) + tfirst * W1;
FindSegmentPolygonIntersection(E0Moved, F1Moved, I);
}
} else {
// Cfg1.type[1] == ‘F’ && Cfg0.type[0] == ‘F’