3D SAT Problem

I wasn't finding a separating axis at all . Here's the solution:

if ( !TestFaces(_psiShape0, _psiShape1) ) {
	return false;
}
if ( !TestFaces( _psiShape1, _psiShape0 ) ) {
	return false;
}
if ( !TestEdges( _psiShape0, _psiShape1 ) ) {
	return false;
}
return true;

You also want to check for parallel edges in the edge test. A random normal can generate a false separation or penetration.

E.g. you could use: |e1 x e2| = |e1| * |e2| * sin( alpha ) *and* sin( alpha ) ~= alpha for alpha << 1

You also want to check for parallel edges in the edge test. A random normal can generate a false separation or penetration.

E.g. you could use: |e1 x e2| = |e1| * |e2| * sin( alpha ) *and* sin( alpha ) ~= alpha for alpha << 1

Thanks for pointing that out.

I'm skipping the edge-edge test if the cross product between them is small than 0.00001f.

Now (little off-topic), I'm facing with the problem of detecting the side planes of the reference face. The problem is that I'm not using an good structure (such the half-edge) to get the connected face planes.

The incident face is just the face of the second object that is most on the negative direction of the reference face right?

Do you have any recommendation?

I'm skipping the edge-edge test if the cross product between them is small than 0.00001f.

If you don't normalize the edge vectors you need to take their length into account. Otherwise this would just be an absolute tolerance. Here is a good discussion how to use relative and absolute tolerances together: http://realtimecollisiondetection.net/pubs/Tolerances/

Now (little off-topic), I'm facing with the problem of detecting the side planes of the reference face.

You don't clip against the neighboring planes. For every edge of the reference face you build a 'side' plane with the reference face normal. E.g. say you have a plane class which takes a normal and point on the plane:

// Build plane (v are the vertices and n is the normal of the reference face)

RnPlane ClipPlane( RnCross( v[ 2 ] - v[ 1 ], n ), v[ 1 ] );

ClipPlane.Normalize();

// Clip (e.g. using Sutherland-Hodgeman)

aPolyhedron.Clip( ClipPlane );

The incident face is just the face of the second object that is most on the negative direction of the reference face right?

Yes, essentially the most anti-parallel face! In practice just the face with the smallest (negative) dot product with the reference face normal.

Do you have any recommendation?

Did you see my presentations on the SAT? You can have a look at my 2013 and 2015 presentations here:

http://www.valvesoftware.com/company/publications.html

Oh, I see. However, I'm confused about the vertices connections. For example, I've a working Sutherland-Hodgman clipping, but it assumes that the vertices are wind ordered in order to get correct edges connections (A to B, A to C, etc.). Another example is that if the edges are not correctly ordered, the side planes creation won't generate correct planes (which I've checked).

For testing purposes, I've manually created an OBB which the edges orders are CCW. In the picture below, the edges are 0-3, 3-2, etc.



Nice GDC slides, BTW.

You mean A to B, B to C, C to D and then X to A I assume. Yes, that is correct, but all you need is an integer array for the edges and then store an offset into that array per face. At the first index you store the number of vertices for that face and the following entries contain the vertex indices for that face. In your example above the top face followed by the front face would be indexed like this { 4, 1, 2, 6, 5, 4, 0, 3, 2, 1, ... }

Then in code you have an array of face planes and a parallel array of offsets into the topology list. In code then this could look something like this:

for ( all faces i )
{
    int offset = mTopologyOffsets[ i ];

    int vertexCount = mFaceTopology[ offset ];
    for ( all vertices k = 1 as long as k <= vertexCount )
    {
        int vertexIndex = mFaces[ offset + k ]
        Vector3 vertex = mFaceTopology[ vertexIndex ];
    }

} 

This is easy and can get you going for a simple box and convex polyhedra.

Oh. Thanks!

I'll try it.

Do you compute all edges during pre-processing step or compute its normals on the fly? I still think I'll need to implement QuickHull and half-edges since it's the most viable structure to acess adjacency information.

- You compute the vertices, local face planes and topology during pre-processing. You then compute the edge cross products at runtime.

- Edges are just indices into the vertex array.

HTH, I am not sure if I understand what the question is.

Quickhull and HE are indeed a very good solution for this problem. But the one I showed here is fine as well and easy to construct for e.g. boxes. The half-edge structure becomes really useful if you start pruning edge pairs that don't build a Minkowski face. But this is an optimization and I would handle this at the end when everything is working correctly.

Thanks.

I'm using supporting points as you described in your slide to get a SA (sep. axis), which is working correctly. But I have some (quite big) questions.

Here:



#1 The distance is signed, so:

CVector3 vSupportToPointOnPlane = vSupport - (p3WorldFacePlane.m_vNormal * p3WorldFacePlane.m_fDist);
float fSeparation = p3WorldFacePlane.m_vNormal.Dot(vSupportToPointOnPlane);

correct?

#2 When using object 1 as the reference face, in the case that the minimum penetration is found by an axis on object 2, I need to flip the sign of the normal (second query in this case). Do I need to reiterate through all faces of object 1 to get the correct index based on the flipped sep. axis direction (which comes from query 2). Eg. checking the face with largest dot product in the opposite direction of the normal coming from object 2 to object1?

This is quite confusing since the object 1 is always the reference face, and using that as the SA looks correct (but is not the actual axis when its not pointing towards object 2).

#3 I'm using your approach to get the correct side planes of the reference face, but I'm getting an wrong behaviour which I think is related to the distance of the side plane to the origin. Here's how I'm doing:

CVector3 vA = _mTransform0 * _pcpHull0->pvVerts[ui32FaceVertexIndexA];
CVector3 vB = _mTransform0 * _pcpHull0->pvVerts[ui32FaceVertexIndexB];

CPlane3& p3SidePlane = p3ReferenceSidePlanes[I];

p3SidePlane.m_vNormal = (vB - vA).Cross(vReference);
p3SidePlane.m_fDist = p3SidePlane.m_vNormal.Dot(vA);
p3SidePlane.Normalize();

In the case of the top face of the box (the image I've posted), I have A and B 1-5, 5-6, 6-2, 2-1, etc.

I think the problem is not related to the winding order of the face edges (I've tried switching between CCW and CW but no succeed).

Maybe I'm missing something on the Sutherland-Hodgman algorithm. I couldn't find any implementation that addresses all issues with the algorithm.

But when tested using a single plane to clip arbitrary points, and I got a clipped polygon with no problem.

Do you have any suggestion?

Here's the clip routine on the case you want to check:

static void SutherlandHodgman(
		const CVector3* _pvPolyVerts, unsigned int _ui32TotalPolyVerts,
		const CPlane3* _p3Planes, unsigned int _ui32TotalPlanes,
		std::vector<CVector3>& _vClippedPoly) {

		//Set initial verts.
		_vClippedPoly.swap(std::vector<CVector3>(_ui32TotalPolyVerts));
		for (unsigned int I = 0; I < _ui32TotalPolyVerts; ++I) {
			_vClippedPoly[I] = _pvPolyVerts[I];
		}

		std::vector<CVector3>* pvOutput = &_vClippedPoly;

		for ( unsigned int I = 0; I < _ui32TotalPlanes; ++I ) {
			const CPlane3& p3Plane = _p3Planes[I];
			
			std::vector<CVector3> vInput(*pvOutput);

			pvOutput->clear();

			if (!vInput.size()) { break; }

			CVector3 vA = vInput[vInput.size() - 1];

			for (unsigned int J = 0; J < vInput.size(); ++J) {
				CVector3 vB = vInput[J];

				float fDA = p3Plane.m_vNormal.Dot(vA) + p3Plane.m_fDist;
				float fDB = p3Plane.m_vNormal.Dot(vB) + p3Plane.m_fDist;

				bool bInsideA = fDA > 0.0f;
				bool bInsideB = fDB > 0.0f;

				CVector3 vP; //Intersection point between a line and a plane.
				if ( fDA == fDB ) {
					vP = ( vA + vB ) * 0.5f; //The mid point seems to be the most reasonable one.
				}
				else {
					float fS = fDA / (fDA - fDB);
					vP = vA + (vB - vA) * fS;
				}

				if (bInsideB) {
					if (!bInsideA) {
						pvOutput->push_back(vP);
					}
					pvOutput->push_back(vB);
				}
				else {
					if (bInsideA) {
						pvOutput->push_back(vP);
					}
				}
				vA = vB;
			}
		}

		_vClippedPoly.swap(*pvOutput);
	}

I've made this snap-shot:



The side planes are correct, assuming that the red face is the reference face and the blue one is the incident. The 4 points on pink are the side planes displaced by its offset (which is correct).

#1:

A plane is defined by a normal and a point on the plane. The plane equation has the form: n*x - n*p = n*x - d = 0. You can now plug in any point into the plane equation to qualify whether the point is in front of, behind or on the plane. If the normal has unit length (called Hesse-Form) the result is the signed distance. So yes, the distance is signed, but it is simply computed like this:

Separation = Dot( Normal, Support ) - PlaneOffset // Where PlaneOffset is the d above or m_fDist in your code.

Make sure that the plane offset has the correct sign. Sometimes people store a plane as: a*x + b*y + b*z + d = 0 to make it more more SIMD friendly. in that case d = -n*p (note the minus!)

#2:

These are two different things if I understand you correctly. For the separating axis you don't need to do anything special, but for the contact normal you might need to flip it if it comes from object 2 as by definition we assume that it points from object 1 to object 2.

I hope I understand your question correctly, otherwise provide a simple example what your problem is.

#3:

It is difficult to check your code since I cannot step through it. 'Real-Time Collision Detection' has a correct implementation of the algorithm which is pretty easy to understand. I copy&paste my implementation here. Maybe you can compare it and find your mistake yourself.

A nice thing about collision detection is that you can visualize everything. So in your case I would draw the face normals and make sure that always point away from the cube. Make sure all edge directions are correct and build a CW or CCW loop around each face. Finally draw the side plane normals and make sure they point to the outside as well. After you have checked that all data is correct I would start investigating your clipping implementation. Otherwise you are just chasing waterfalls.

//--------------------------------------------------------------------------------------------------
void rnClipPolygon( rnPolygon& Out, const rnPolygon& Polygon, const rnPlane& Plane )
	{
        RN_ASSERT( Out.Empty() ); 
	RN_ASSERT( Polygon.Size() >= 3 );
	
	rnVertex Vertex1 = Polygon.Back();
	float Distance1 = rnDistance( Plane, Vertex1.Position );
	
	for ( int Index = 0; Index < Polygon.Size(); ++Index )
		{
		rnVertex Vertex2 = Polygon[ Index ];
		float Distance2 = rnDistance( Plane, Vertex2.Position );
	
		if ( Distance1 <= 0.0f && Distance2 <= 0.0f )
			{
			// Both vertices are behind the plane - keep vertex2
			Out.PushBack( Vertex2 );
			}
		else if ( Distance1 <= 0.0f && Distance2 > 0.0f )
			{
			// Vertex1 is behind the plane, vertex2 is in front -> intersection point
			float Fraction = Distance1 / ( Distance1 - Distance2 );
			rnVector3 Position = Vertex1.Position + Fraction * ( Vertex2.Position - Vertex1.Position );
	
			// Keep intersection point 
			rnVertex Vertex;
			Vertex.Position = Position;
			Out.PushBack( Vertex );
			}
		else if ( Distance2 <= 0.0f && Distance1 > 0 )
			{
			// Vertex2 is behind the plane, vertex1 is in front -> intersection point
			float Fraction = Distance1 / ( Distance1 - Distance2 );
			rnVector3 Position = Vertex1.Position + Fraction * ( Vertex2.Position - Vertex1.Position );
	
			// Keep intersection point 
			rnVertex Vertex;
			Vertex.Position = Position;
			Out.PushBack( Vertex );
	
			// And also keep vertex2
			Out.PushBack( Vertex2 );
			}
	
		// Keep vertex2 as starting vertex for next edge
		Vertex1 = Vertex2;
		Distance1 = Distance2;
		}
	}
 

#2:
This are two things. For the separating axis you don't need to do anything special, but for the contact normal you might need to flip it if it comes from object 2 as by definition we assume that it points from object 1 to object 2.

I hope I understand your question correctly, otherwise provide a simple example what your problem is.

Yes, perfectly answered ^^.

#3:
It is difficult to check your code since I cannot step through it. 'Real-Time Collision Detection' has a correct implementation of the algorithm which is pretty easy to understand. I copy&paste my implementation here. Maybe you can compare it and find your mistake yourself. A nice thing about collision detection is that you can visualize everything. So in your case I would draw the face normals and make sure that alway point away from the cube. Make sure all edge directions are correct and build CW or CCW loop around each face. Finally draw the side plane normals and make sure they point outside as well. After you have checked that all data is correct I would start investigating your clipping implementation. Otherwise you are just chasing waterfalls.

The face normals are pointing in the correct direction. On the image the reference face is the red one, and the incident the blue one, which is actually -X and +X, respectively. I've ommited just for claricity purposes.

Thanks for sharing your code, I would never find that the RTCD book had the implementation . So many good topics on that book.

Are the side planes normals pointing towards the center of the reference object (enclosing it) or outwards? In my implementation the clipping routine will discard the points if the normals are outwards. Sorry, it's because it didn't became explicitly on the slide (just in the figure) when you say "keep all contacts below the reference face and after that project on it".

Fixed. I was generating unecessary polygon points.

The side planes normals are pointing outward for me. It really only depends on what your clipping code expects. The clipping code I posted expects outward pointing normals.

I've fixed this part (just updated the post) .

Now it's working as it should .

As you said before, the contact normal flipping is just a matter of contact resolution direction. The only change, in order to clip the correct reference face, is invert the reference object:

if ( bIsFaceContactA && bIsFaceContactB ) {
		if (fqQuery0.fSeparation < fqQuery1.fSeparation) {
			CreateFaceContact(fqQuery0, fqQuery1, mTransform0, pcpHull0, mTransform1, pcpHull1);
		}
		else {
			CreateFaceContact(fqQuery1, fqQuery0, mTransform1, pcpHull1, mTransform0, pcpHull0);
		}

		//Flip contact normal if necessary here.
	}
	else {
		CreateEdgeContact(eqQuery, mTransform0, pcpHull0, mTransform1, pcpHull1);
	}

am I correct?

Yup, the only remaining problem is now that you favor face over edge contacts and avoid feature flip-flops between faces. I don't know what bIsFaceContactA and bIsFaceContactB are doing, but a simple way to enforce this would be something along these lines (assuming MKS for the absolute tolerance):

const float kLinearSlop = 0.005f;
const float kRelEdgeTolerance = 0.90f;
const float kRelFaceTolerance = 0.98f;
const float kAbsTolerance = 0.5f * kLinearSlop;

if ( EdgeQuery.Separation > kRelEdgeTolerance * rnMax( FaceQuery1.Separation, FaceQuery2.Separation ) + kAbsTolerance )
	{
	rnBuildEdgeContact( Out, Transform1, Hull1, Transform2, Hull2, EdgeQuery );
	}
else
	{
	if ( FaceQuery2.Separation > kRelFaceTolerance * FaceQuery1.Separation + kAbsTolerance )
		{
		// Face contact (2)
		bool FlipNormal = true;
		rnBuildFaceContact( Out, Transform2, Hull2, Transform1, Hull1, FaceQuery2, FlipNormal );
		}
	else
		{
		// Face contact (1)
		bool FlipNormal = false;
		rnBuildFaceContact( Out, Transform1, Hull1, Transform2, Hull2, FaceQuery1, FlipNormal );
		}
	}

For the logic part here remember that at this point all separations are negative (since we didn't find a separating axis). So greater means closer to zero and max( s1, s2 ) returns the one closer to zero as well.

Just enfatizing that all edge combinations must be supporting edges (a face in the Minkowski Diff.) in order to get the correct distance of the possible separating axis. Basically, on the brute-force way, this is checking if all edges of each objects are the closest ones in respect to the edge cross-product normal (similar to a supporting point). This blows up the hole SAT performance if not optimized.

Thanks Dirk for all your answers .

Yeah, you cannot just check the edges of the minimizing faces. I see people do this, but it is wrong. Also just checking the leaving edges from the vertices of the minimizing faces is not enough. There are configurations where this can fail as well.

Again, for performance the trick is to *not* call the SAT at all and utilize temporal coherence. If you find a separating axis check it in the next frame first. Chances are high it is still a separating axis. If you find two touching features inspect the relative orientation in the next frame and try to build a contact directly from these features if the relative configuration didn't change within some reasonable tolerance.