射线三角形交 C++

Question

我正在测试射线是否与三角形相交。目前，我使用以下代码测试三角形和从原点到三角形中点的射线之间是否存在相交：

Ray<float> *ray = new Ray<float>(Vec3<float>(0), chosenTriangle->GetTriangleMidpoint()); 

这里是我用来存储向量运算的Vec3对象:

template<typename T>
class Vec3
{
public:
    T x, y, z;
    Vec3() : x(T(0)), y(T(0)), z(T(0)) { }
    Vec3(T xx) : x(xx), y(xx), z(xx) { }

    Vec3(T xx, T yy, T zz) : x(xx), y(yy), z(zz) {}
    Vec3& normalize()
    {
        T nor2 = length2();
        if (nor2 > 0) {
            T invNor = 1 / sqrt(nor2);
            x *= invNor, y *= invNor, z *= invNor;
        }
        return *this;
    }

    Vec3<T> operator * (const T &f) const { return Vec3<T>(x * f, y * f, z * f); }
    Vec3<T> operator * (const Vec3<T> &v) const { return Vec3<T>(x * v.x, y * v.y, z * v.z); }
    T dot(const Vec3<T> &v) const { return x * v.x + y * v.y + z * v.z; }
    Vec3<T> operator - (const Vec3<T> &v) const { return Vec3<T>(x - v.x, y - v.y, z - v.z); }
    Vec3<T> operator + (const Vec3<T> &v) const { return Vec3<T>(x + v.x, y + v.y, z + v.z); }
    bool operator == (const Vec3<T> &v) { return x == v.x && y == v.y && z == v.z; }
    Vec3<T> operator - () const { return Vec3<T>(-x, -y, -z); }
    T length2() const { return x * x + y * y + z * z; }
    T length() const { return sqrt(length2()); }
    Vec3<T> CrossProduct(Vec3<T> other)
    { 
        return Vec3<T>(y*other.z - other.y*z, x*other.z - z*other.x, x*other.y - y*other.x); 
    }
    friend std::ostream & operator << (std::ostream &os, const Vec3<T> &v)
    {
        os << "[" << v.x << " " << v.y << " " << v.z << "]";
        return os;
    }

选定的三角形和光线具有以下值，其中vertA，vertB和vertC是三角形的顶点，并且在表示三角形的对象中找到。

计算指定射线和交点之间是否存在相交的代码如下。此代码在三角形对象方法中找到，其中vertA，vertB和vertC是全局变量。

bool CheckRayIntersection(Vec3<T> &o, Vec3<T> &d)
{
    Vec3<T> e1 = vertB - vertA;
    Vec3<T> e2 = vertC - vertA;
    Vec3<T> p = d.CrossProduct(e2);
    T a = e1.dot(p);

    if(a == 0)
        return false;

    float f = 1.0f/a;

    Vec3<T> s = o - vertA;
    T u = f * s.dot(p);
    if(u < 0.0f || u > 1.0f)
        return false;

    Vec3<T> q = s.CrossProduct(e1);
    T v = f * d.dot(q);

    if(v < 0.0f || u+v > 1.0f)
        return false;

    T t = f * e2.dot(q);

    return (t >= 0);
    
}

我仍然从函数得到了一个错误的返回值，但我认为它应该返回true，因为通过三角形中点的向量应该在三角形中点处相交。有人可以告诉我我的代码哪里有问题吗？或者错误的返回false实际上是正确的吗？

- Adrian De Barro

---

@SigTerm 没有给出完整的代码，他们只是三角形的坐标，将会更新。 - Adrian De Barro

看起来你正在实现Möller-Trumbore算法。我的直觉告诉我，你的一个Vec3方法可能有误，因为你的实现具有我所提到的那个方法的“形状”。 - Rerito

在进行叉积运算之前，您应该对e1和e2进行归一化处理。 - MichaelCMS

@Rerito：好知道。本来想给他写一个“未优化”的逐步碰撞检测，但因为太长而放弃了。 - SigTerm

@AdrianDeBarro：真的很离题，但是看看这个光线追踪教程：http://www.flipcode.com/archives/Raytracing_Topics_Techniques-Part_1_Introduction.shtml。在那里实现了光线/三角形相交（它也有代码）。你至少可以检查一下你的结果是否与另一个实现一致。 - MichaelCMS

Answer 1

有了你的数据，我通过对光线方向进行归一化（这是代码中唯一明显的更改）来获得一致的结果。

以下是代码实现（我参考了论文，但它并不是非常优化）：

struct quickVect
{

   float x,y,z;
   float l;
 };

#define DOT(v1,v2) (v1.x*v2.x + v1.y*v2.y+v1.z*v2.z)
#define CROSS(rez,v1,v2) \
rez.x  = v1.y*v2.z - v1.z*v2.y; \
rez.y  = v1.z*v2.x - v1.x*v2.z; \
rez.z  = v1.x*v2.y - v1.y*v2.x;

#define SUB(rez,v1,v2) \
rez.x = v1.x-v2.x; \
rez.y = v1.y-v2.y; \
rez.z = v1.z-v2.z;


#define LENGTH(v) (sqrtf(v.x* v.x + v.y*v.y + v.z*v.z))

#define NORMALIZE(v) \
v.l = LENGTH(v); \
v.x = v.x / v.l; \
v.y = v.y / v.l; \
v.z = v.z / v.l;

#define EPSILON 0.000001f

//#define TEST_CULL

bool testIntersection(quickVect& v1, quickVect& v2, quickVect& v3, quickVect& orig,quickVect& dir)
{
 quickVect e1,e2,pvec,qvec,tvec;

 SUB(e1,v2,v1);
 SUB(e2,v3,v1);

 CROSS(pvec,dir,e2);

 NORMALIZE(dir);
 //NORMALIZE(pvec);
 float det = DOT(pvec,e1);
#ifdef TEST_CULL
if (det <EPSILON)
{

    return false;
}
SUB(tvec,orig,v1);
float u = DOT(tvec,pvec);
if (u < 0.0 || u > det)
{

    return false;
}
CROSS(qvec,tvec,e1);
float v = DOT(dir,qvec);
if (v < 0.0f || v + u > det)
{

    return false;
}
#else
 if (det < EPSILON && det > -EPSILON )
 {

     return false;
 }

 float invDet = 1.0f / det;
 SUB(tvec,orig,v1);
// NORMALIZE(tvec);
 float u = invDet * DOT(tvec,pvec);
 if (u <0.0f || u > 1.0f)
 {

     return false;
 }
 CROSS(qvec,tvec,e1);
// NORMALIZE(qvec);
 float v = invDet* DOT(qvec,dir);
 if (v < 0.0f || u+v > 1.0f)
 {

     return false;
 }
#endif
 return true;
}

- MichaelCMS

Answer 2

将MichaelCMS的答案直接翻译为glm可用的代码。

// must normalize direction of ray
bool intersectRayTri(Tri& tri, glm::vec3 o, glm::vec3 n) {
    glm::vec3 e1, e2, pvec, qvec, tvec;

    e1 = tri.v2 - tri.v1;
    e2 = tri.v3 - tri.v1;
    pvec = glm::cross(n, e2);
    
    n = glm::normalize(n);
    //NORMALIZE(pvec);
    float det = glm::dot(pvec, e1);

    if (det != 0)
    {
        float invDet = 1.0f / det;
        tvec = o - tri.v1;
        // NORMALIZE(tvec);
        float u = invDet * glm::dot(tvec, pvec);
        if (u < 0.0f || u > 1.0f)
        {

            return false;
        }
        qvec = glm::cross(tvec, e1);
        // NORMALIZE(qvec);
        float v = invDet * glm::dot(qvec, n);
        if (v < 0.0f || u + v > 1.0f)
        {

            return false;
        }
    }
    else return false;
    return true; // det != 0 and all tests for false intersection fail
}

- paulytools