本文原創著作權歸 百度
laizhishen 所有，如有轉載，請按如下方式顯式標明原創作者及原文出處，以示尊重！！
作者：laizhishen
原文：http://hi.baidu.com/laizhishen/blog/item/3d206d209cca9c54ac34de46.htmlFrustum裁剪是CLOD中很重要的一個演算法，很多文章都是一句話就過去，或者直接給出代碼。但是數學推導很少給出，本文章的目的就是解釋大家看這些代碼中的疑問。透視投影是將相機空間中的點從視錐體(view frustum)變換到規則觀察體(Canonical View Volume，CCV)中，即是世界空間的視錐體(view frustum)中的任何一個點，如果經過投影矩陣變換後，它必定規則觀察體(Canonical View Volume，CCV)中，也即是在(-1,-1,0) ~ (1,1,1)之間的值1所示這裡 我們定義P0=(-1,-1,0) P1=(1,-1,0)P2=(1,-1,1)P3=(-1,-1,1)P4=(-1,1,0)P5=(1,1,0)P6=(1,1,1)P7=(-1,1,1)我們通過這8個點構建6個面，通過3點共面，假設面是Pos0，pos1，pos2構成，u=Pos1-Pos0，v=Pos2-Pos0,那麼法向量n=u×v。則d=-(n×Pos0)。得到每個平面的平面公式n，d，從而得到A，B，C，D，(n的xn，yn，zn，d，就是a，b，c，d)。這裡要注意構造面試後P點的順序，d3D是按順時針來構造面。法線都是由frustum裡到外。Near： (P0,P4,P5) n=(0,0,-1),d=0 0x+0y-1z+0=0 Far: (P2,P6,P7) n=(0,0,1),d=-1 0x+0y+1z-1=0 Left: (P0,P3,P7) n=(-1,0,0),d=-1 -1x+0y+0z-1=0 Right: (P1,P5,P6) n=(1,0,0),d=-1 1x+0y+0z-1=0 Top: (P4,P7,P6) n=(0,1,0),d=-1 0x+1y+0z-1=0 Bottom: (P0,P1,P2) n=(0,-1,0),d=-1 0x-1y+0z-1=0 我們假設這六個平面中某個平面上存在一個點P’(x’,y’,c’,1)，那它的平面方程為A’x’ + B’y’+ C’z’+ D’= 0，在進行投影變換之前的座標為P(x,y,z,1)，在世界空間中的平面方程Ax + By + Cz + D = 0，如果P‘點在CVV中，那點P(x,y,z,1)必定在view Frustum中。原理一1.推理由介紹我們可以得到下面的等式 |A’| |A| (x’,y’,c’,1) × |B’| = 0 (x,y,z,1) × |B| = 0|C’| |C| |D’| |D| 同時(x’,y’,c’,1) × Tvproj= (x,y,z,1)注意Tviewproj矩陣=Tview×Tproj(即是攝像機矩陣和投影矩陣相乘。是將全局座標轉換到視平面的變化矩陣)結合著3個等式我們可以得到|A| |A’|Tviewproj × |B| = |B’| 等式1|C| |C’||D| |D’|通過等式1，我們可以求出6個面在全局座標系下的平面方程Ax + By + Cz + D = 02代碼 void CFrustum::InitFrustum1(const D3DXMATRIX& aoViewMatrix,const D3DXMATRIX& aoProjMatrix){ D3DXMATRIX loComboMatrix;D3DXMATRIX loInvComboMatrix;D3DXMatrixMultiply(&loComboMatrix,&aoViewMatrix,&aoProjMatrix);// 求得view * proj的逆矩陣.D3DXMatrixInverse(&loInvComboMatrix, NULL, &loComboMatrix );// 如果經過投影矩陣，所有的三維全局座標的點都變為(-1,-1,0) ~ (1,1,1)之間的值.// 將同次空間的臨界值填入moCVVPos.moCVVPos[0].x = -1.0f; moCVVPos[0].y = -1.0f; moCVVPos[0].z = 0.0f;moCVVPos[1].x = 1.0f; moCVVPos[1].y = -1.0f; moCVVPos[1].z = 0.0f;moCVVPos[2].x = 1.0f; moCVVPos[2].y = -1.0f; moCVVPos[2].z = 1.0f;moCVVPos[3].x = -1.0f; moCVVPos[3].y = -1.0f; moCVVPos[3].z = 1.0f;moCVVPos[4].x = -1.0f; moCVVPos[4].y = 1.0f; moCVVPos[4].z = 0.0f;moCVVPos[5].x = 1.0f; moCVVPos[5].y = 1.0f; moCVVPos[5].z = 0.0f;moCVVPos[6].x = 1.0f; moCVVPos[6].y = 1.0f; moCVVPos[6].z = 1.0f;moCVVPos[7].x = -1.0f; moCVVPos[7].y = 1.0f; moCVVPos[7].z = 1.0f;//將cvv座標轉換回全局座標for(int i = 0; i < 8; i++ )D3DXVec3TransformCoord( &moCVVPos[i], &moCVVPos[i], &loInvComboMatrix );// 通過得到的全局座標製作平截頭體平面.// 向量由平截頭體內部指向外部的平面.D3DXPlaneFromPoints(&moFrustumPlane[FRUSTUM_PLANE_FRONT], moCVVPos , moCVVPos+4, moCVVPos+5); // 近平面(near)D3DXPlaneFromPoints(&moFrustumPlane[FRUSTUM_PLANE_BACK], moCVVPos+2, moCVVPos+6, moCVVPos+7); // 遠平面(far)D3DXPlaneFromPoints(&moFrustumPlane[FRUSTUM_PLANE_LEFT], moCVVPos , moCVVPos+3, moCVVPos+7); // 左平面(left)D3DXPlaneFromPoints(&moFrustumPlane[FRUSTUM_PLANE_RIGHT], moCVVPos+1, moCVVPos+5, moCVVPos+6); // 右平面(right)D3DXPlaneFromPoints(&moFrustumPlane[FRUSTUM_PLANE_TOP], moCVVPos+4, moCVVPos+7, moCVVPos+6); // 上平面(top)D3DXPlaneFromPoints(&moFrustumPlane[FRUSTUM_PLANE_BOTTOM], moCVVPos , moCVVPos+1, moCVVPos+2); // 下平面(bottom)}上面的代碼很容易理解，完全按原理來實現，但是當我們每次調用這個函數，我們每次都要求Tvproj矩陣的逆矩陣，同時還要重新用點來構造每個frustum平面，並且調用的是D3DXPlaneFromPoints函數。我們還有更好的方法嗎？下面我們接著最佳化我們的代碼，最佳化之前我們接著分析原理二1.推理我們知道在P(x,y,z,1) 經過Tviewproj矩陣變換後得到點P’(x’,y’,z’,w’),這個點在前面的推導過程中,保證是在frustum的某個平面上的.這個時候我們知道如果對P的每個分量都除以w’就可以把P’歸一化到一個長方體的空間，即x’/w’, y’/w’在[-1,1]區間，z’/w’在[0,1]區間，所以如果有一個投影轉換後的點P1’,它的x1’,y1’在[-w’,w’]區間,z1’在[0,w’]區間，這個點肯定就在視錐體內.|a11,a12,a13,a14|假設Tviewproj = [v0,v1,v2,v3] = |a21,a22,a23,a24||a31,a32,a33,a34||a41,a42,a43,a44|其中其中V0,v1,v2,v3是四個列向量.P ’= P×Tviewproj = (P?v0, P?v1, P?v2, P?v3 ) = (x’,y’,z’,w’) 根據上面的x的範圍我們有: -w’ <= x’ 就是 -p’?v3 <= p’?v0 就是 p’?v0 + p’?v3 <= 0 得到 P’?(v0+v3) <= 0 把列向量換成矩陣的元素有 (x,y,z,1).(m_11 + m_14, m12 + m24, m13 + m34, m14 + m44 ) <= 0 就是 (m_11+m_14)*x + (m_12 + m_24)*y + (m_13 + m_34)*z + (m_14 + m_44)*1 <=0 簡單地看這是一個 A*x + B*y + C*z + D*1 <= 0 描述了一個半空間，就是平面A*x + B*y + C*z + D*1 = 0右邊的空間所以我們知道是錐體的左裁減面為(m_11+m_14)*x + (m_12 + m_24)*y + (m_13 + m_34)*z + (m_14 + m_44)*1= 0. 相應地可以計算出 6個裁減面 Left= (v0 + v3).Right =(v3 - v0) Bottom =(v3 + v1) Top = (v3 - v1 ) Near = (v2) Far = (v3 - v2) 裁減的時候把點帶入公式Ax+By+Cz+Dw看大與0還是小與0就可以知道在平面的裡面還是外面，在實際計算中需要對A，B，C，D進行歸一化。因為我們在處理座標的時候，點幾乎都是經過歸一化的。2.代碼void CFrustum::InitFrustum(const D3DXMATRIX& aoViewMatrix,const D3DXMATRIX& aoProjMatrix){ D3DXMATRIX loComboMatrix;D3DXMatrixMultiply(&loComboMatrix,&aoViewMatrix,&aoProjMatrix);// calculate the planes// NearD3DXPLANE* lpPlane = &moFrustumPlane[FRUSTUM_PLANE_NEAR];lpPlane->a = loComboMatrix._14 + loComboMatrix._13; lpPlane->b = loComboMatrix._24 + loComboMatrix._23;lpPlane->c = loComboMatrix._34 + loComboMatrix._33;lpPlane->d = loComboMatrix._44 + loComboMatrix._43;D3DXPlaneNormalize(lpPlane,lpPlane);// FarlpPlane = &moFrustumPlane[FRUSTUM_PLANE_FAR];lpPlane->a = loComboMatrix._14 - loComboMatrix._13; lpPlane->b = loComboMatrix._24 - loComboMatrix._23;lpPlane->c = loComboMatrix._34 - loComboMatrix._33;lpPlane->d = loComboMatrix._44 - loComboMatrix._43;D3DXPlaneNormalize(lpPlane,lpPlane);//LeftlpPlane = &moFrustumPlane[FRUSTUM_PLANE_LEFT];lpPlane->a = loComboMatrix._14 + loComboMatrix._11; // LeftlpPlane->b = loComboMatrix._24 + loComboMatrix._21;lpPlane->c = loComboMatrix._34 + loComboMatrix._31;lpPlane->d = loComboMatrix._44 + loComboMatrix._41;D3DXPlaneNormalize(lpPlane,lpPlane);// RightlpPlane = &moFrustumPlane[FRUSTUM_PLANE_RIGHT];lpPlane->a = loComboMatrix._14 - loComboMatrix._11; lpPlane->b = loComboMatrix._24 - loComboMatrix._21;lpPlane->c = loComboMatrix._34 - loComboMatrix._31;lpPlane->d = loComboMatrix._44 - loComboMatrix._41;D3DXPlaneNormalize(lpPlane,lpPlane);// ToplpPlane = &moFrustumPlane[FRUSTUM_PLANE_TOP];lpPlane->a = loComboMatrix._14 - loComboMatrix._12; lpPlane->b = loComboMatrix._24 - loComboMatrix._22;lpPlane->c = loComboMatrix._34 - loComboMatrix._32;lpPlane->d = loComboMatrix._44 - loComboMatrix._42;D3DXPlaneNormalize(lpPlane,lpPlane);// BottomlpPlane = &moFrustumPlane[FRUSTUM_PLANE_BOTTOM];lpPlane->a = loComboMatrix._14 + loComboMatrix._12; // BottomlpPlane->b = loComboMatrix._24 + loComboMatrix._22;lpPlane->c = loComboMatrix._34 + loComboMatrix._32;lpPlane->d = loComboMatrix._44 + loComboMatrix._42;D3DXPlaneNormalize(lpPlane,lpPlane);}第2種方法在沒有看見之間，先看見的是代碼http://www.racer.nl/reference/vfc_markmorley.htm ，自己想了好久都沒考慮清楚。最後看見解釋，鄙視自己的數學邏輯推理。3d看來困難的就在于思考方式，同樣的代碼處理，第二種方式明顯快很多。判斷裁剪我們已經有了裁剪體的方程，當我們需要裁剪一個頂點的時候，這六個方程已經足夠了。但是我們要判斷一個地區的可見度時，我們進行一些額外的計算。2所示，一個物體和投影體的關係大致可以分為：包圍、被包圍、相交和相離四種情況。圖中最大的淺藍色的矩形包圍了整個投影體。深綠色的小矩形則完全被投影體包圍。淺綠色的矩形和投影體相交。這三種情況下物體都是可以被看到的。剩下紅色的矩形則和投影體相離、只有它完全不可見圖2當處理節點的可見度的時候，由於節點的不規則性。我們還需要引入包圍體的概念。所謂的包圍體，就是用一個比較簡單的幾何體去度量另外一個比較複雜的幾何體，讓它剛好能包圍另外一個幾何體。比較合適的包圍體外形有矩形、正方形和球體。其中球體處理最為簡單，但是近似度也最差。我們為每一個節點都建立一個包圍體，只要測試這個包圍體，我們就可以決定一個節點的可見度，由於包圍體肯定大於這個節點，因此我們可以保證不會有任何可見的節點被裁剪在投影體之外BOOL CFrustum::PointInFrustum( float afX, float afY, float afZ ){// A*x+B*y+C*z+D = 0 is in plane,for(int i = 0; i < FRUSTUM_PLANE_COUNT; i++ ){D3DXPLANE &loPlane = moFrustumPlane[i];//減少函數調用//if(D3DXPlaneDotCoord(&moFrustumPlane[i], &D3DXVECTOR3(afX, afY, afZ)) < 0.0f)if(loPlane.a * afX + loPlane.b * afY + loPlane.c * afZ + loPlane.d < 0.0f)return FALSE;}return TRUE;}BOOL CFrustum::SphereInFrustum( float afX, float afY, float afZ, float AfRadius ){// A*x+B*y+C*z+D = -radius is in plane,for(int i = 0; i < FRUSTUM_PLANE_COUNT; i++ ){D3DXPLANE &loPlane = moFrustumPlane[i];//減少函數調用，直接用公式運算//if(D3DXPlaneDotCoord(&moFrustumPlane[i], &D3DXVECTOR3(afX, afY, afZ)) < -AfRadius)if( loPlane.a * afX + loPlane.b * afY + loPlane.c * afZ + loPlane.d <= -AfRadius ){return false;}}return TRUE;}BOOL CFrustum::CubeInFrustum( float afX, float afY, float afZ, float aiSize,BOOL & abIsCompletelyContained ){float lfAlfaX = afX + aiSize;float lfDeltaX = afX - aiSize;float lfAlfaY = afY + aiSize;float lfDeltaY = afY - aiSize;float lfAlfaZ = afZ + aiSize;float lfDeltaZ = afZ - aiSize;DWORD ldwNumPointInFrustum = 0;for(int i = 0; i < FRUSTUM_PLANE_COUNT; i++ ){int j = 8;BOOL lbIsInAllPlanes = TRUE;D3DXPLANE &loPlane = moFrustumPlane[i];//if(D3DXPlaneDotCoord(&loPlane[i], &D3DXVECTOR3(lfDeltaX,lfDeltaY, lfDeltaZ)) < 0.0f)if(loPlane.a * lfDeltaX + loPlane.b * lfDeltaY + loPlane.c * lfDeltaZ + loPlane.d < 0.0f){lbIsInAllPlanes = FALSE;j--;}//if(D3DXPlaneDotCoord(&loPlane[i], &D3DXVECTOR3(lfAlfaX,lfDeltaY, lfDeltaZ)) < 0.0f)if(loPlane.a * lfAlfaX + loPlane.b * lfDeltaY + loPlane.c * lfDeltaZ + loPlane.d < 0.0f){lbIsInAllPlanes = FALSE;j--;}//if(D3DXPlaneDotCoord(&loPlane[i], &D3DXVECTOR3(lfDeltaX,lfAlfaY, lfDeltaZ)) < 0.0f)if(loPlane.a * lfDeltaX + loPlane.b * lfAlfaY + loPlane.c * lfDeltaZ + loPlane.d < 0.0f){lbIsInAllPlanes = FALSE;j--;}//if(D3DXPlaneDotCoord(&loPlane[i], &D3DXVECTOR3(lfAlfaX,lfAlfaY, lfAlfaZ)) < 0.0f)if(loPlane.a * lfAlfaX + loPlane.b * lfAlfaY + loPlane.c * lfDeltaZ + loPlane.d < 0.0f){lbIsInAllPlanes = FALSE;j--;}//if(D3DXPlaneDotCoord(&loPlane[i], &D3DXVECTOR3(lfDeltaX,lfDeltaY, lfAlfaZ)) < 0.0f)if(loPlane.a * lfDeltaX + loPlane.b * lfDeltaY + loPlane.c * lfAlfaZ + loPlane.d < 0.0f){lbIsInAllPlanes = FALSE;j--;}//if(D3DXPlaneDotCoord(&loPlane[i], &D3DXVECTOR3(lfAlfaX,lfDeltaY, lfAlfaZ)) < 0.0f)if(loPlane.a * lfAlfaX + loPlane.b * lfDeltaY + loPlane.c * lfAlfaZ + loPlane.d < 0.0f){lbIsInAllPlanes = FALSE;j--;}//if(D3DXPlaneDotCoord(&loPlane[i], &D3DXVECTOR3(lfDeltaX,lfDeltaY, lfAlfaZ)) < 0.0f)if(loPlane.a * lfDeltaX + loPlane.b * lfDeltaY + loPlane.c * lfAlfaZ + loPlane.d < 0.0f){lbIsInAllPlanes = FALSE;j--;}//if(D3DXPlaneDotCoord(&loPlane[i], &D3DXVECTOR3(lfAlfaX,lfAlfaY, lfAlfaZ)) < 0.0f)if(loPlane.a * lfAlfaX + loPlane.b * lfAlfaY + loPlane.c * lfAlfaZ + loPlane.d < 0.0f){lbIsInAllPlanes = FALSE;j--;}// if none contained, return FALSE.if(0 == j)return FALSE;// update counter if they were all in front of plane.if(lbIsInAllPlanes)++ldwNumPointInFrustum;}abIsCompletelyContained = (BOOL)(ldwNumPointInFrustum == FRUSTUM_PLANE_COUNT);return TRUE;}上面這些原理都來至於網路，一些細節是我在理解過程中自己整理，沒有什麼特別之處，主要是歸類而已。在學習3D中一定要將原理吃透。這樣當積累到一定的時候，你自己就會發現3D也不過如此，但關鍵自己一定要勤於動手和自己推理一遍。
