新建網(wǎng)頁(yè) 1
球和平面的相交性檢測(cè)
球和平面的靜態(tài)檢測(cè)相對(duì)容易一些,可以用公式12.14來(lái)計(jì)算球心到平面的距離。如果距離小于球半徑,那么它們相交。實(shí)際上還能作一種更靈活的檢測(cè),這種檢測(cè)把相交分為球完全在平面正面,完全在背面,跨平面等三種情況。仔細(xì)分析程序清單13.2:
Listing 13.2: Determining which side of a plane a sphere is on
// Given a sphere and plane, determine which side of the plane
// the sphere is on.
//
// Return values:
//
// < 0 Sphere is completely on the back
// > 0 Sphere is completely on the front
// 0 Sphere straddles plane
int classifySpherePlane(
const Vector3 &planeNormal, // must be normalized
float planeD, // p * planeNormal = planeD
const Vector3 &sphereCenter, // center of sphere
float sphereRadius // radius of sphere
)
{
// Compute distance from center of sphere to the plane
float d = planeNormal * sphereCenter – planeD;
// Completely on the front side?
if (d >= sphereRadius)
return +1;
// Completely on the back side?
if (d <= –sphereRadius)
return –1;
// Sphere intersects the plane
return 0;
}
動(dòng)態(tài)檢測(cè)要稍微復(fù)雜一些。設(shè)平面為靜止的,球作所有的相對(duì)位移。
平面的定義方式一如既往,用標(biāo)準(zhǔn)形式p . n = d,n為單位向量。球由半徑r和初始球心位置c定義。球的位移,由單位向量d指明方向,L代表位移的距離。t從0變化到L,用直線方程c+td計(jì)算球心的運(yùn)動(dòng)軌跡。如圖13.13所示:

不管在平面上的哪一點(diǎn)上發(fā)生相交,在球上的相交點(diǎn)總是固定的,認(rèn)識(shí)到這一點(diǎn)能大大簡(jiǎn)化問(wèn)題。用c-rn來(lái)計(jì)算交點(diǎn),如圖13.14所示:

現(xiàn)在我們知道了球上的相交點(diǎn),就可以利用射線與平面相交性檢測(cè)的方法,替換掉公式13.6中的p0,得到公式13.9:

射線和三角形的相交性檢測(cè)
在圖形學(xué)和計(jì)算幾何中射線與三角形的相交性檢測(cè)是非常重要的。因?yàn)槿狈ι渚€和復(fù)雜物體間相交性檢測(cè)的方法,我們通常用三角網(wǎng)格代表(或至少是近似代表)物體表面,再作射線和三角網(wǎng)格的相交性檢測(cè)。第一步是計(jì)算射線和包含該三角形的平面的交點(diǎn),第二步是通過(guò)計(jì)算交點(diǎn)的重心坐標(biāo),來(lái)判斷它是否在三角形中。
為了使測(cè)試效率盡可能高,要使用下列技巧:
(1)在檢測(cè)中盡可能地返回負(fù)值(沒(méi)有相交),這稱(chēng)作"提前結(jié)束(early out)"。
(2)盡可能地延遲昂貴的數(shù)學(xué)運(yùn)算,如除法。有兩個(gè)原因:第一,如果并不需要昂貴運(yùn)算的結(jié)果,比如說(shuō)遇到了提前結(jié)束的情況,那么執(zhí)行這些運(yùn)算的時(shí)間就白白浪費(fèi)了。第二,它給了編譯器更多的空間以利用現(xiàn)代處理器的指令管道的優(yōu)點(diǎn)。在準(zhǔn)備除法運(yùn)算時(shí),它產(chǎn)生執(zhí)行其他測(cè)試的代碼(可能導(dǎo)致提前結(jié)束)。所以,在執(zhí)行期間,如果確實(shí)需要除法運(yùn)算的結(jié)果,該結(jié)果可能已經(jīng)被計(jì)算出來(lái),或至少已經(jīng)部分被計(jì)算出來(lái)了。
(3)只檢測(cè)與三角形正面的相交,這幾乎可以節(jié)省一半的檢測(cè)時(shí)間。
// Ray-triangle intersection test.
//
// Algorithm from Didier Badouel, Graphics Gems I, pp 390-393
float rayTriangleIntersect(
const Vector3 &rayOrg, // origin of the ray
const Vector3 &rayDelta, // ray length and direction
const Vector3 &p0, // triangle vertices
const Vector3 &p1, // .
const Vector3 &p2, // .
float minT) // closest intersection found so far. (Start with 1.0)
{
// We'll return this huge number if no intersection is detected
const float kNoIntersection = 1e30f;
// Compute clockwise edge vectors.
Vector3 e1 = p1 – p0;
Vector3 e2 = p2 – p1;
// Compute surface normal. (Unnormalized)
Vector3 n = crossProduct(e1, e2);
// Compute gradient, which tells us how steep of an angle
// we are approaching the *front* side of the triangle
float dot = n * rayDelta;
// Check for a ray that is parallel to the triangle or not pointing toward the front face
// of the triangle.
//
// Note that this also will reject degenerate triangles and rays as well. We code this in a
// very particular way so that NANs will bail here. (This does not behave the same as
// "dot >= 0.0f" when NANs are involved.)
if (!(dot < 0.0f))
return kNoIntersection;
// Compute d value for the plane equation. We will use the plane equation with d on the right side:
//
// Ax + By + Cz = d
float d = n * p0;
// Compute parametric point of intersection with the plane containing the triangle, checking at the
// earliest possible stages for trivial rejection.
float t = d – n * rayOrg;
// Is ray origin on the backside of the polygon? Again, we phrase the check so that NANs will bail.
if (!(t <= 0.0f))
return kNoIntersection;
// Closer intersection already found? (Or does ray not reach the plane?)
//
// since dot < 0:
//
// t/dot > minT
//
// is the same as
//
// t < dot * minT
//
// (And then we invert it for NAN checking
)
if (!(t >= dot * minT))
return kNoIntersection;
// OK, ray intersects the plane. Compute actual parametric point of intersection.
t /= dot;
assert(t >= 0.0f);
assert(t <= minT);
// Compute 3D point of intersection
Vector3 p = rayOrg + rayDelta * t;
// Find dominant axis to select which plane
// to project onto, and compute u's and v's
float u0, u1, u2;
float v0, v1, v2;
if (fabs(n.x) > fabs(n.y))
{
if (fabs(n.x) > fabs(n.z))
{
u0 = p.y – p0.y;
u1 = p1.y – p0.y;
u2 = p2.y – p0.y;
v0 = p.z – p0.z;
v1 = p1.z – p0.z;
v2 = p2.z – p0.z;
}
else
{
u0 = p.x – p0.x;
u1 = p1.x – p0.x;
u2 = p2.x – p0.x;
v0 = p.y – p0.y;
v1 = p1.y – p0.y;
v2 = p2.y – p0.y;
}
}
else
{
if (fabs(n.y) > fabs(n.z))
{
u0 = p.x – p0.x;
u1 = p1.x – p0.x;
u2 = p2.x – p0.x;
v0 = p.z – p0.z;
v1 = p1.z – p0.z;
v2 = p2.z – p0.z;
}
else
{
u0 = p.x – p0.x;
u1 = p1.x – p0.x;
u2 = p2.x – p0.x;
v0 = p.y – p0.y;
v1 = p1.y – p0.y;
v2 = p2.y – p0.y;
}
}
// Compute denominator, check for invalid.
float temp = u1 * v2 – v1 * u2;
if (!(temp != 0.0f))
return kNoIntersection;
temp = 1.0f / temp;
// Compute barycentric coords, checking for out-of-range at each step
float alpha = (u0 * v2 – v0 * u2) * temp;
if (!(alpha >= 0.0f))
return kNoIntersection;
float beta = (u1 * v0 – v1 * u0) * temp;
if (!(beta >= 0.0f))
return kNoIntersection;
float gamma = 1.0f - alpha - beta;
if (!(gamma >= 0.0f))
return kNoIntersection;
// Return parametric point of intersection
return t;
}
還有一個(gè)能優(yōu)化昂貴計(jì)算的策略沒(méi)體現(xiàn)在上述代碼中:即預(yù)先計(jì)算結(jié)果。如果像多邊形向量這樣的值預(yù)先被計(jì)算出來(lái)的話(huà),就可以采用更加優(yōu)化的策略。
射線和AABB的相交性檢測(cè)
檢測(cè)AABB和射線的相交性非常重要,因?yàn)楦鶕?jù)檢測(cè)的結(jié)果可以避免對(duì)更復(fù)雜物體的測(cè)試。(例如,我們要檢測(cè)射線與多個(gè)由三角網(wǎng)格組成的物體的相交性,可以先計(jì)算射線和三角網(wǎng)格的AABB的相交性。有時(shí)候可以一次就排除整個(gè)物體,而不必去檢測(cè)這個(gè)物體的所有三角形。)
Woo提出一種方法,先判斷矩形邊界框的哪個(gè)面會(huì)相交,再檢測(cè)射線與包含這個(gè)面的平面的相交性。如果交點(diǎn)在盒子中,那么射線與矩形邊界框相交,否則不存在相交。
cAABB3中rayIntersect()就是用Woo的技術(shù)來(lái)實(shí)現(xiàn)的。
//---------------------------------------------------------------------------
// Parametric intersection with a ray. Returns parametric point
// of intsersection in range 0
1 or a really big number (>1) if no
// intersection.
//
// From "Fast Ray-Box Intersection," by Woo in Graphics Gems I, page 395.
//
// See 12.9.11
//---------------------------------------------------------------------------
float AABB3::rayIntersect(const Vector3& rayOrg, // origin of the ray
const Vector3& rayDelta, // length and direction of the ray
Vector3* returnNormal) const // optionally, the normal is returned
{
// We'll return this huge number if no intersection
const float kNoIntersection = 1e30f;
// Check for point inside box, trivial reject, and determine parametric distance to each front face.
bool inside = true;
float xt, xn;
if (rayOrg.x < min.x)
{
xt = min.x - rayOrg.x;
if (xt > rayDelta.x)
return kNoIntersection;
xt /= rayDelta.x;
inside = false;
xn = -1.0f;
}
else if (rayOrg.x > max.x)
{
xt = max.x - rayOrg.x;
if (xt < rayDelta.x)
return kNoIntersection;
xt /= rayDelta.x;
inside = false;
xn = 1.0f;
}
else
xt = -1.0f;
float yt, yn;
if (rayOrg.y < min.y)
{
yt = min.y - rayOrg.y;
if (yt > rayDelta.y)
return kNoIntersection;
yt /= rayDelta.y;
inside = false;
yn = -1.0f;
}
else if (rayOrg.y > max.y)
{
yt = max.y - rayOrg.y;
if (yt < rayDelta.y)
return kNoIntersection;
yt /= rayDelta.y;
inside = false;
yn = 1.0f;
}
else
yt = -1.0f;
float zt, zn;
if (rayOrg.z < min.z)
{
zt = min.z - rayOrg.z;
if (zt > rayDelta.z)
return kNoIntersection;
zt /= rayDelta.z;
inside = false;
zn = -1.0f;
}
else if (rayOrg.z > max.z)
{
zt = max.z - rayOrg.z;
if (zt < rayDelta.z)
return kNoIntersection;
zt /= rayDelta.z;
inside = false;
zn = 1.0f;
}
else
zt = -1.0f;
// Inside box?
if (inside)
{
if (returnNormal != NULL)
{
*returnNormal = -rayDelta;
returnNormal->normalize();
}
return 0.0f;
}
// Select farthest plane - this is
// the plane of intersection.
int which = 0;
float t = xt;
if (yt > t)
{
which = 1;
t = yt;
}
if (zt > t)
{
which = 2;
t = zt;
}
switch (which)
{
case 0: // intersect with yz plane
{
float y = rayOrg.y + rayDelta.y * t;
if (y < min.y || y > max.y)
return kNoIntersection;
float z = rayOrg.z + rayDelta.z * t;
if (z < min.z || z > max.z)
return kNoIntersection;
if (returnNormal != NULL)
{
returnNormal->x = xn;
returnNormal->y = 0.0f;
returnNormal->z = 0.0f;
}
}
break;
case 1: // intersect with xz plane
{
float x = rayOrg.x + rayDelta.x * t;
if (x < min.x || x > max.x)
return kNoIntersection;
float z = rayOrg.z + rayDelta.z * t;
if (z < min.z || z > max.z)
return kNoIntersection;
if (returnNormal != NULL)
{
returnNormal->x = 0.0f;
returnNormal->y = yn;
returnNormal->z = 0.0f;
}
}
break;
case 2: // intersect with xy plane
{
float x = rayOrg.x + rayDelta.x * t;
if (x < min.x || x > max.x)
return kNoIntersection;
float y = rayOrg.y + rayDelta.y * t;
if (y < min.y || y > max.y)
return kNoIntersection;
if (returnNormal != NULL)
{
returnNormal->x = 0.0f;
returnNormal->y = 0.0f;
returnNormal->z = zn;
}
}
break;
}
// Return parametric point of intersection
return t;
}