3D数学之-几何图元和几何检测

这段时间满脑子都是向量，射线，三角形,平面，AABB的几何图元和他们之间的空间关系，几何检测分类分解处理组合结果之类的情景。感觉回到的热爱几何和数学的中学时代，希望自己把这样的几何向量化描述，抽象模型和分解组合演绎计算的思想能够不断积累思考下去，变成一种比较强的能力，类似设计模式，实际编程优化重构代码的能力一样。

一、几何图元

1.射线

射线在计算机图形学中占有重要的位置，向量表达式：

p(t) = p0 + td,p0为起点，d = (delta x, delta y)是终点，t属于[0,1]; 或者d是方向单位向量，t是[0,l]。

注意射线起点可能在其它图元的正面或者反面，射线的方向可能正向也 可能反向，射线的长度或许不够长也需要判断。

2.AABB包围盒

AABB包围盒，用 vector3 min, vector3 max来表示实用性较强。

当AABB变换时用变换后的物体重新计算AABB会更加准确，但是直接变换AABB的速度会更快，且该计算可以根据变换过程得到快速计算的方法：

原始变换AABB思路：

1）通过观察向量乘以矩阵乘法每个新的顶点都是来自于原顶点。

2）新的8个顶点需要add函数筛选，遍历每个新的顶点筛选出最小的min.x值，其它也一样，得到新的min,max 新包围盒顶点不是来自于变换前某个顶点，而是从新的8个顶点中add函数筛选取得min.x,min.y,min.z,同样的取得max, min.x来自于原AABB中的一个点，min.y来自于AABB的另一个点或同一点，min.z,max也同理；至于它们min.x要相加那是因为向量乘以矩阵的要求。

AABB快速计算算法：

正负数导致的乘法结果大小不一样，因此前面的1），2）两步可以合并起来得到，一开始就可以根据矩阵上的元素m11，m12,m13等最小的min.x的分量，min.y分量,min.z分量的，矩阵其它行基向量上也一样，于是可以直接变换时候相加得到结果。

这个算法是和变换AABB一样的结果，得到一个至少和原AABB一样大的AABB,或大得多的AABB 。

这个算法难以理解的地方就是整个变换连接起来 ，不要局限在1）步骤了，而要看到2）步骤才能理解；这个想法类似通过分析整个计算过程进行去掉冗余计算，对结果进行提前计算的优化策略。

3.平面

平面用 p.n = d来定义，n是平面法向量，平面的法向量由左手定则决定方向(四指顺时针绕行，拇指指的方向是法向量方向)，平面的正面是从法向量箭头看到法向量起点的面，反面反之。

一般三角形的平面可以直接计算法向量n和代入点得到d值，但是非三角形下的多边形可能会导致不在同一平面上，因向量共线不能求取到正确的n,因此多边形的最佳平面需要最佳平面公式,设多边形有k点：

nx = （z1+z2)(y1-y2) +  （z2+z3)(y2-y3) + ... + （zk+zk+1)(yk-yk+1); // zk+1,yk+1等于z1,y1

ny = （x1+x2)(z1-z2) + （x2+x3)(z2-z3) + ... + （xk+xk+1)(zk-zk+1); // xk+1,zk+1等于x1,z1

nz = （y1+y2)(x1-x2) + （y2+y3)(x2-x3) + ... + （yk+yk+1)(xk-xk+1); // yk+1,xk+1等于y1,z1

d = ( p1.n + p2.n + ... + p k.n) /   k

点q到平面的距离：

a = q.n - d

理解q.n得到该点到平面法向量投影相对原点的距离，平面的法向量是从原点指出的；在平面上的点p,p.n = d是平面上的点投影到平面法向量上的距离。这个距离之差就是点q到平面的垂直距离。

                                       

4.三角形

三角形，需要用三个向量点来定义 vector3 v1, vector3 v2,vector3 v3。三角形网格中因为共用顶点关系可以通过邻接信息去掉冗余的顶点，或者可以通过索引缓存来去掉冗余的顶点。

三角形重要的定理有勾股定理，正弦定理，余弦定理。重要的定义有重心，内心，外心。

3D中的三角形面积为： A = | e1 x e2 | / 2, e1和e2向量的夹角和方向没有限制，任取三角形的两条边向量即可。

2D中的三角形面积, 也可以用坐标系下3边下的有符号面积相加得到，边是由左往右的面积为正数，反之面积为负：

         y+

          |                v2       

          |             /\

          |            /   \

          |      v1  —— v3

          ---------------------------------x+

T =  (y2+y1)(x2-x1)/ 2 + (y3+y2)(x3-x2)/ 2 + (y1+y3)(x1-x3) / 2

   = ( y1(x2 - x3) + y2(x3 - x1) + y3(x1 - x2) ) / 2          // 直接乘出来，合并多项式了下

     = ( (y1 - y3)(x2 - x3) + (y2 - y3)(x3 - x1) ) / 2             // 平移三角形不会改变面积，所有y值都减去了y3

记得v1,v2,v3是顺时针绕序的。

三角形中引入了重心坐标空间的概念，这样在三角形平面上移动就比较好表示了，否则用3D的点或者向量表示会比较复杂。

其实这也是一种分层模块化和转换的思想。

三角形重心坐标是三角形三个顶点向量在三角形表面上贡献的比例。

(b1,b2,b3) <=>b1v1 + b2v2 + b3v3

所以三角形三个顶点的重心分别为：（1,0,0), (0,1,0), (0,0,1)

重心坐标（b1,b2,b3)其实刚好是三角形内部点p划分的三个小三角形的面积和总面积的比：

      y+

          |                v2       

          |             / \

       e3->T3    / p  \e1->T1

          |      v1  —— v3

                       e2->T2

          ---------------------------------x+

         b1 = T1/ T

         b2 = T2/T

         b3 = T3/T

2D中T1，T2，T3，T的面积可以直接求取从而得到三角形的重心坐标。

3D中可以用3D中求取面积的向量叉乘法，也可以通过将三角形投影到面积最大的2D平面上（ 通过法向量的分量可以看出，分量越大的那么三角形投影到该分量上越小，所以投影到垂直于大分量的平面上即得最大投影），从而可以用2D的方法来进行三角形重心坐标的求取。在射线和三角形检查中，先求得射线和三角形平面的交点，然后计算交点的三角形重心坐标，以决定交点是否在三角形内部)。

3D中可以用3D中求取面积的向量叉乘法，向量叉乘会有一个小问题，就是叉乘结果总是正数，当对于三角形重心坐标在三角形外部的，那么不能够正确的表达出来，这个时候将叉乘结果点乘法向量，就会得到有符号的结果（ 因为如果重心坐标在三角形外部，那么与在三角形内部的的两个向量围成的三角形绕序会相反 ）。

三角形重心坐标为(1/3, 1/3,1/3)，计算物理重心时候需要用到。

三角形内心重心坐标为(L1/p, L2/p, L3/p),L是三角形的边长，p是三角形的周长，当需要计算三角形内接圆时候用到。

三角形的外心重心坐标比较复杂些下面给出代码，当需要计算三角形外接圆时候才需要使用到。

// returnRadius是需要返回外接圆半径，  returnG是返回的重心坐标

void caculateOutHeart( Vector3 &v1, Vector3 &v2, Vector3 &v3,float &returnRadius, Vector3 &returnG )

{

 Vector3 e1 = v3 - v2;

 Vector3 e2 = v1 - v3;

 Vector3 e3 = v2 - v1;

 float d1 = -(e2.x * e3.x + e2.y*e3.y + e2.z*e3.z);

 float d2 = -(e3.x*e1.x + e3.y*e1.y + e3.z*e1.z);

 float d3 = -(e1.x*e2.x + e1.y*e2.y + e1.z*e3.z);

 float c1 = d2*d3;

 float c2 = d1*d3;

 float c3 = d1*d2;

 float c = c1 + c2 + c3;

 float cc = 2*c;

 float g1 = ( c2 + c3 ) / cc;

 float g2 = ( c1 + c3 ) / cc;

 float g3 = (c1 + c2 ) / cc;

 returnG = Vector3( g1, g2, g3 );

 float radius = 0.5 * sqrt( (d1+d2) * (d2+d3) * (d3+d1) / c );

 returnRadius = radius;

}

5.复杂多边形

多边形的凹凸检测，可以用所有多边形边向量之间点乘计算cos theta，得到所有外角  theta的和，如果等于2pi那么是凸多边形，如果是凹的那么会大于2pi。或者向量之间取向量(Vi-1 - Vi)和向量(Vi+1 - Vi)计算内角和，如果是凹多边形，那么会出现 (Vi-1 - Vi)和(Vi+1 - Vi)是外角的情况，且该外角会比内角小，所以如果内角和等于(n-2)*pi,那么是凸多边形，如果小于(n-2)*pi，那么是凹多边形。

对多边形的分解，需要用扇形分解即可，一般少用到。

二、几何检测

几何检测中比较重要的是射线和平面的相交检测，射线和AABB的相交性检测，射线和三角形的相交性检测，射线和球的相交性检测；球和平面的相交性检测， AABB和平面的相交检测，两个AABB间静态和动态的检测。

还有更多的检测，具体的算法见AABB类。还有更多的检测，如果需要还要进行分析检测。

因为算法比较多和繁琐，下面总结下自己理解这些算法过程中和研究这些算法中得到的思想方法。

1.抽象建模为向量的思维

对研究的对象向量化，平面，射线，AABB,球都用向量表示，研究的空间关系，距离为向量的模，夹角是向量的cos thta，法向量为sin theta。

2.多运动简化为单一运动，多轴上的变换简化为单一轴的变换最后合并简化问题

如果是两个对象是动态类型的检测，那么将运动都假设放置到一个物体上，一个假设为静止的。向量的加法和减法刚好可以做到。

   射线，三角形间的检测，需要考虑不同情形，投影到面上，AABB间检测需分开每个轴考虑，最后对时间轴求交集。

3.利用向量构建关系，引入辅助中间量，分解合并转化未知量，代入方程解决问题

利用向量的点乘，叉乘，和三角关系勾股正弦余弦，三角形相似，三角函数构建关系。根据几何意义引入中间辅助线辅助量，分解合并式子得到未知量。如果不行那么强解三角形也是可以的。得到未知量，将未知量代入物体向量方程得到结果。

有些问题，可以提前结束，比如计算比例，那么双方都不用除2是可以的，有些可以比较除法，避免平方根运算得到结果。

4.分类讨论和注意精度

对于关系式的结果，注意分类讨论，射线原点在平面正面反面；射线方向是否和需要检查的对象之间同向或者反向；射线的长度不够长；几何体内部或者几何体相交了；注意前提条件和运算的整个过程；在平面的正面反面，共线平行，共面的几何关系等情况。

还有因为浮点数精度问题，很多计算结果在合理的误差范围内，都是可以接受的。

具体的几何检测方面用到的众多几何图元和几何图元间的相交检测关系，见 AABB几何检测类的代码。

三、AABB几何检测类

遇到的问题：

考虑前提：前提条件需要注意，比如n.x=0，导致一开始就可以对x方向上做判断； 一些方向一开始就是反着的；一些点乘结果就是负数的。

考虑连续变换：AABB利用观察变换过程得到快速的算法，这个算法是由重新筛选和筛选后重新add组成的，所以需要注意前后的变换。

考虑原点方向和长度，分解和组合： 需要考虑到射线原点在背面(正面反面)，射线方向反了(向量方向和移动方向），射线的长度短了（没有图元间距离长)。 平行的状况，方向相反了，方位错位了会导致不3D空间中不想交。

// AABB3.h

#ifndef _AABB3_H_INCLUDE_

#define _AABB3_H_INCLUDE_

#ifndef _VECTOR3_H_INCLUDE_

#include "Vector3.h"

#endif

inline float mymax(float a,float b){ return a > b? a : b;}

inline float mymin(float a,float b){ return a < b? a : b;}

inline void myswap(float &a,float &b){float temp = a;a = b;b = temp;}

class Matrix4x3;//Matrix4_3

//AABB3是一个轴对齐的3D 矩形AABB边界框类，所谓轴对齐，就是AABB的min顶点处是坐标轴的惯性坐标系原点

class AABB3{

public:

 Vector3 min;

 Vector3 max;

 //查询各种参数

 Vector3 size() const {return max - min;}

 float xSize() const {return max.x - min.x;}

 float ySize() const {return max.y - min.y;}

 float zSize() const {return max.z - min.z;}

 Vector3 center() const{return (min + max) * 0.5f;}

 //提取八个顶点中的一个

 Vector3 corner(int i) const;

 //清空矩形边界框

 void empty();

 //向矩形边界框中添加点

 void add(const Vector3 &p);

 //向矩形边界框中添加AABB

 void add(const AABB3 &box);

 //变换矩形边界框(用矩阵来变换)，计算新的AABB，非常实用

 void setToTransformedBox(const AABB3 &box,const Matrix4x3 &m);//Matrix4_3

 //包含相交性测试

 bool isEmpty() const;

 //返回true,如果矩形边界框包含改点

 bool contains(const Vector3 &p) const;

 //返回矩形边界框上的最近点，是伏帖在AABB表面，而不是8个顶点中的一个

 Vector3 closestPointTo(const Vector3 &p) const;//closest最靠近

 //返回true,如果和球相交

 bool intersectsSphere(const Vector3 &center,float radius) const;//intersect相交

 //和参数射线的相交性测试,如果不相交则返回值大于1的1e30f，非常实用

 float rayIntersect(const Vector3 &rayOrg,const Vector3 &rayDelta,Vector3 *returnNormal = 0) const;

 //判断矩形边框在平面的那一面

 int classifyPlane(const Vector3 &n,float d) const;

 //AABB和平面的动态相交性测试，非常实用

 float intersectPlane(const Vector3 &n,float planeD,const Vector3 &dir) const;

};

//检测两个AABB的相交性，如果相交返回true,还可以返回相交部分的AABB

bool intersectAABBs(const AABB3 &box1,const AABB3 &box2,AABB3 *boxIntersect = 0);

//返回运动AABB和静止AABB相交时的参数点，如果不相交则返回值大于1，非常实用

float intersectMovingAABB(const AABB3 &stationaryBox,const AABB3 &movingBox,const Vector3 &d);

#endif

// AABB3.cpp

#include <assert.h>

#include <stdlib.h>

#include "AABB3.h"

#include "Matrix4x3.h" //Matrix4_3.h

//#include "CommonStuff.h"

/*

    6 -------7

   /|       /|

  / |   / |        +y

 2-------- 3 |         |   +z

 | |  |  |         |  /

 |  |4----|--|5        | /

 | /      | /          |/

 0 --------1           o------------- +x

 min点是0点

 max点是7点

*/

Vector3 AABB3::corner(int i) const{ //i为0~7传递进来即可

 assert(i >= 0);

 assert(i <= 7);

 return Vector3(

  (i & 1) ? max.x : min.x,

  (i & 2) ? max.y : min.y,

  (i & 4) ? max.z : min.z

  );

}

//将其赋值为极大和极小值，以清空矩形边界框

void AABB3::empty(){

 const float kBigNumber = 1e37f; // 1e37 == 1*(10的37次方),若为1e-37 == 1*（10的-37次方).

 min.x = min.y = min.z = kBigNumber;//若设置为0,0,0则对于运算中用到min,max作为分母则导致问题;或屏幕上显示为黑点

 max.x = max.y = max.z = -kBigNumber;

}

//向矩形边界框中加一个点,确保该点在框内，如果不在则扩张矩形边界框以包含该点

void AABB3::add(const Vector3 &p){

 if (p.x < min.x) min.x = p.x;

 if (p.x > max.x) max.x = p.x;

 if (p.y < min.y) min.y = p.y;

 if (p.y > max.y) max.y = p.y;

 if (p.z < min.z) min.z = p.z;

 if (p.z > max.z) max.z = p.z;

}

//向矩形边界框中添加AABB,确保要添加的AABB在类AABB内，否则扩张类的AABB

void AABB3::add(const AABB3 &box){

 if (box.min.x < min.x) min.x = box.min.x;

 if (box.max.x > max.x) max.x = box.max.x;

 if (box.min.y < min.y) min.y = box.min.y;

 if (box.max.y > max.y) max.y = box.max.y;

 if (box.min.z < min.z) min.z = box.min.z;

 if (box.max.z > max.z) max.z = box.max.z;

}

//变换矩形边界框并计算新的AABB

//记住，这将得到一个至少和原AABB一样大的AABB,或大得多的AABB

//观察向量乘以矩阵计算结果表达式,且从得到的新8个顶点中计算新的最大点最小点，见add函数，

// 它不是来自于变换前某个顶点而是从新的8个顶点中add筛选取得min.x,min.y,min.z,同样的取得max,

// min.x来自于原AABB中的一个点，min.y来自于AABB的另一个点或同一点，min.z,max也同理；至于它们要相加那是因为向量乘以矩阵的要求。

// 通过观察向量乘以矩阵乘法，和得到新的8个顶点后的add函数筛选，正负数导致的乘法结果大小，避免了先计算八个顶点的变换后的值，再从八个顶点计算新的AABB.

void AABB3::setToTransformedBox(const AABB3 &box, const Matrix4x3 &m){ //box是需要变换的AABB,本类为结果AABB

 if(box.isEmpty()){

  empty();

  return;

 }

 //从平移部分开始,AABB和轴对齐的，所以Matrix平移部分就是AABB的平移部分，接着在惯性坐标系下计算新的顶点

 min = max = getTranslation(m);

 //依次检查矩阵的九个元素，计算新的AABB

 //m的第一行

 if (m.m11 > 0.0f){

  min.x += m.m11 * box.min.x; max.x += m.m11 * box.max.x;

 }else{

        min.x += m.m11 * box.max.x; max.x += m.m11 * box.min.x;

 }

 if (m.m12 > 0.0f){

  min.y += m.m12 * box.min.x; max.y += m.m12 * box.max.x;

 }else{

  min.y += m.m12 * box.max.x; max.y += m.m12 * box.min.x;

 }

 if(m.m13 > 0.0f){

  min.z += m.m13 * box.min.x; max.z += m.m13 * box.max.x;

 }else{

  min.z += m.m13 * box.max.x; max.z += m.m13 * box.min.x;

 }

 //m的第二行

 if(m.m21 > 0.0f){

  min.x += m.m21 * box.min.y; max.x += m.m21 * box.max.y;

 }else{

  min.x += m.m21 * box.max.y; max.x += m.m21 * box.min.y;

 }

 if(m.m22 > 0.0f){

  min.y += m.m22 * box.min.y; max.y += m.m22 * box.max.y;

 }else{

  min.y += m.m22 * box.max.y; max.y += m.m22 * box.min.y;

 }

 if(m.m23 > 0.0f){

  min.z += m.m23 * box.min.y; max.z += m.m23 * box.max.y;

 }else{

  min.z += m.m23 * box.max.y; max.z += m.m23 * box.min.y;

 }

 //m的第三行

 if (m.m31 > 0.0f){

  min.x += m.m31 * box.min.z; max.x += m.m31 * box.max.z;

 }else{

  min.x += m.m31 * box.max.z; max.x += m.m31 * box.min.z;

 }

 if (m.m32 > 0.0f){

  min.y += m.m32 * box.min.z; max.y += m.m32 * box.max.z;

 }else{

  min.y += m.m32 * box.max.z; max.y += m.m32 * box.min.z;

 }

 if (m.m33 > 0.0f){

  min.z += m.m33 * box.min.z; max.z += m.m33 * box.max.z;

 }else{

  min.z += m.m33 * box.max.z; max.z += m.m33 * box.min.z;

 }

}

bool AABB3::isEmpty() const{ //因为empty函数设置矩阵空时候是反过来了；包括计算中导致矩阵反过来了，那么也是空

 return (min.x > max.x) || (min.y > max.y) || (min.z > max.z);

}

bool AABB3::contains(const Vector3 &p) const {

 return

  (p.x >= min.x && p.x <= max.x) &&  //取址对象 > 算术+左右移 > 关系+位运算 >逻辑 > 赋值

  (p.y >= min.y && p.y <= max.y) &&

  (p.z >= min.z && p.z <= max.z);

}

Vector3 AABB3::closestPointTo(const Vector3 &p) const {

 Vector3 r;

 if (p.x < min.x){

  r.x = min.x;

 }else if(p.x > max.x){

  r.x = max.x;

 }else{//这样就不会落在8个顶点上，实现在外部的点拉近到AABB表面

  r.x = p.x;

 }

 if(p.y < min.y){

  r.y = min.y;

 }else if(p.y > max.y){

  r.y = max.y;

 }else{

  r.y = p.y ;

 }

 if(p.z < min.z){

  r.z = min.z;

 }

 else if(p.z > max.z){

  r.z = max.z;

 }else{

  r.z = p.z;

 }

 return r;

}

//球与AABB的相交性检查，如果最近点到球心距离小于半径则返回真

bool AABB3::intersectsSphere(const Vector3 &center,float radius) const {

 Vector3 closestPoint = closestPointTo(center);

 return distanceSquared(center,closestPoint) < radius*radius;

}

//射线和AABB相交性的检测，来自woo的算法:先检测相交的平面；再求射线和平面相交。

//函数的返回值是射线的t值在[0,1]之间，如果大于1则不相交

float AABB3::rayIntersect(const Vector3 &rayOrg,

        const Vector3 &rayDelta,

        Vector3 * returnNormal) //返回相交平面的法向量

        const {

//检测射线和AABB相交于那个平面的选择原理:将射线rayDelta分解为x,y,z方向向量，最迟到达AABB区域内部的分量向量起着决定作用.

//最迟到达AABB区域的向量在总的射线rayDelta中所消耗的比率也越大；计算后比率最大的和该分量轴为法向量的AABB平面相交。

 const float kNoIntersection = 1e30f;

 bool inside = true;

 // xt是x方向上相交的t比率，xn是x方向上需要返回射线于AABB相交表面的法向量的x分量表面的正面对着入射射线

 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;

 }

   //是否在矩形边界框内

 if (inside){

  if (returnNormal != NULL){

   *returnNormal = -rayDelta;

   returnNormal->normalize();

  }

  return 0.0f;

 }

 //选择最远的平面(也就是最迟到达的平面),发生相交的地方

 int which = 0;

 float t = xt;

 if (yt > t){

  which = 1;

  t = yt;

 }

 if (zt > t){

  which = 2;

  t = zt;

 }

 switch (which){

 case 0: //x轴最迟到达(法向量在x轴)，在yz平面上相交

  {

     float y = rayOrg.y + rayDelta.y*t;

     if (y<min.y || y>max.y) return kNoIntersection;//y轴上面在最迟到达交点时候，向上或向下穿出去了，不在AABB区域了

     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: //y轴最迟到达(法向量在y轴)，在xz平面上相交

  {

   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: //z轴最迟到达(法向量在z轴)，在xy平面上相交

  {

   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 t;

}

// 射线与AABB的检测，该射线不用考虑射线方向向量的大小，假设方向向量是无限长的，这个判断在用鼠标在3D世界中拾取常用到
// 避免参数方程，代入球面方程判断解，只能用球形包围盒的模式
bool AABB3::RayIntersectOneAxisOK(const float minValue, const float maxValue, const float fRayOrg,  const float fRayDelta)
{
    // 后面也可以做该判断，但是这里提前检查下，如果在包围盒里面就直接返回，因为在区域内的概率也是很大的
    if( fRayOrg >= minValue && fRayOrg <= maxValue )
    {
        return true;
    }
    // nx是相交平面的平面法向量分量值
    bool bIsOK = false;
    float nAxis = -1.0f;
    if( fRayOrg > maxValue )
    {
        nAxis = 1.0f;
    }

    if( nAxis < 0)
    {
        float fXMax = fRayOrg + (maxValue - fRayOrg) * fRayDelta;
        // fXMin不用考虑了，因为一般情况下最小的都不会相交，就算相交通过fXMax可以同样得到判断
        if( fXMax >= minValue && fXMax <= maxValue)
        {
            bIsOK = true;
        }
    }
    else
    {
        float fXMin = fRayOrg + (fRayOrg - minValue) * fRayDelta;
        if( fXMin >= minValue && fXMin <= maxValue)
        {
            bIsOK = true;
        }
    }

    return bIsOK;
}

bool    AABB3::rayIntersectLimitless(const Vector3 &rayOrg, Vector3 rayDelta)
{
    // 因为射线方向向量没有固定长度的，那么先要单位化
    rayDelta.normalize();
    return ( RayIntersectOneAxisOK(min.x, max.x, rayOrg.x, rayDelta.x)
        && RayIntersectOneAxisOK(min.y, max.y, rayOrg.y, rayDelta.y)
        && RayIntersectOneAxisOK(min.z, max.z, rayOrg.z, rayDelta.z));
}

//判断静止的平面和AABB的相交性，如果返回值大于0则在平面正面，小于0在背面，等于0相交

// 此为用p.n = d,任意三个分量可以决定一个顶点，所以根据n.x,n.y,n.z正负求取即可；至于minD需要相加是因为p.n=minD

// 本身向量乘法需要相加。

int AABB3::classifyPlane(const Vector3 &n, float d) const{

 float minD,maxD;

 if(n.x > 0.0f){

  minD = min.x * n.x;

  maxD = max.x * n.x;

 }else{

  minD = max.x * n.x;

  maxD = min.x * n.x;

 }

 if(n.y > 0.0f){

  minD += min.y * n.y;

  maxD += max.y * n.y;

 }else{

  minD += max.y * n.y;

  maxD += min.y * n.y;

 }

 if(n.z > 0.0f){

  minD += min.z * n.z;

  maxD += max.z * n.z;

 }else{

  minD += max.z * n.z;

  maxD += min.z * n.z;

 }

 // AABB在平面的正面

 if(minD >= d){

  return 1;

 // AABB在平面的反面

 }else if(maxD <= d){

  return -1;

 }else{

 // 相交了

  return 0;

 }

}

//动态AABB和平面的检测:AABB可以分解为射线和可以求得最近点和最远点的顶点集。

//1.只考虑正面检测，则向量点乘小于0为正面运动,平面的法向量是起点在原点方，指向远离坐标原点位置的；

//2.求得最近点和最远点距离平面的距离做初步判断；

//3.求得最近点射线到平面的距离

float AABB3::intersectPlane(

       const Vector3 &n,

       float planeD,

       const Vector3 &dir // AABB移动的方向和大小dir属于[0,l]

       ) const {

 assert(fabs(n*n - 1.0f) < 0.01f);

 assert(fabs(dir*dir - 1.0f) < 0.01f);

 const float kNoIntersection = 1e30f;

 float dot =  n * dir;

 if (dot >= 0.0f) // dir要向着平面正面移动，如果是逆着平面正面移动或者垂直平面法向量移动都是没有运动相交的，这时可以用静态检查下是否相交。

  return kNoIntersection;

 // 计算AABB到平面的最近点和最远点沿着n方向距离原点的距离

 float minD, maxD;

 if(n.x > 0.0f){

  minD = n.x * min.x;maxD = n.x * max.x;

 }else{

  minD = n.x * max.x; maxD = n.x * min.x;

 }

 if (n.y > 0.0f){

  minD += n.y * min.y; maxD += n.y * max.y;

 }else{

  minD += n.y * max.y; maxD += n.y * min.y;

 }

 if (n.z > 0.0f){

  minD += n.z * min.z; maxD += n.z * max.z;

 }else{

  minD += n.z * max.z; maxD += n.z * min.z;

 }

 // 在平面背面了，且移动方向是平面法向量相反的方向，所以不会相交的了

 if(maxD < planeD)

  return kNoIntersection;

 // maxD >= planeD 且 minD < planeD的情况，也就是开始就相交了

 if(minD < planeD)

  return 0.0f;

 // planeD - minD是平面到AABB垂直距离的负数，dot是AABB到在n方向上总投影的距离的负数，除法t刚好得到倍数t。

 float t = (planeD - minD) / dot;

 return t;

}

//

//全局非成员函数代码

//

//检测两个AABB是否相交，如果相交则返回true,也可以返回相交部分

bool intersectAABBs(const AABB3 &box1,

     const AABB3 &box2,

     AABB3 *boxIntersect /* = 0 */){

 // 存在交点的情况下才相交

 if (box1.min.x > box2.max.x) return false;

 if (box1.min.y > box2.max.y) return false;

 if (box1.min.z > box2.max.z) return false;

 if (box2.min.x > box1.max.x) return false;

 if (box2.min.y > box1.max.y) return false;

 if (box2.min.z > box1.max.z) return false;

 

 // 相交部分的AABB

 if (boxIntersect != NULL){

  boxIntersect->min.x = mymax(box1.min.x,box2.min.x);

  boxIntersect->min.y = mymax(box1.min.y,box2.min.y);

  boxIntersect->min.z = mymax(box1.min.z,box2.min.z);

  boxIntersect->max.x = mymin(box1.max.x,box2.max.x);

  boxIntersect->max.y = mymin(box1.max.y,box2.max.y);

  boxIntersect->max.z = mymin(box1.max.z,box2.max.z);

 }

 return true;

}

//两个动态AABB间的检测

//前提是：两个AABB都是和坐标轴对齐的，如果不对齐则用其它方法。

//检测方法：

//1.合成：将两个的运动合成了一个的运动，一个静止。

//2.分解：将一个矢量运动分解了x,y,z轴上的三个独立的运动。

//3.分解相交公式:te = (Smin - Mmax(0)) / di ;di是某个轴上的分量

//                tl = (Smax - Mmin(0)) / di;

//4.分量合成：将三个[te,tl]的分量合成，得到交集最小的te,就是相交的时刻。

float intersectMovingAABB(const AABB3 &stationaryBox,

        const AABB3 &movingBox,

         const Vector3 &d) {

 const float kNoIntersection = 1e30f;

 // 进入的时刻

 float tEnter = 0.0f;

 // 离开的时刻

 float tLeave = 1.0f;

 //提前检测：

 //1.如果某一维度上面di为零，则开始不相交就返回

 //2.如果两个维度上相交为空，则也返回

 if(d.x == 0.0f){

  if(stationaryBox.min.x > movingBox.max.x ||

   movingBox.min.x > stationaryBox.max.x)

   return kNoIntersection;

   }else{

    float oneOverD = 1.0f / d.x;

    float xEnter = (stationaryBox.min.x - movingBox.max.x) * oneOverD;

    float xLeave = (stationaryBox.max.x - movingBox.min.x) * oneOverD;

    if(xEnter > xLeave){//可鞥是从左往右的移动，AABB对齐特性，故要求交换

     myswap(xEnter,xLeave);

    }

    // enter取得大的

    if (xEnter > tEnter) tEnter = xEnter;

    // leave取得小的

    if (xLeave < tLeave) tLeave = xLeave;

    // xEnter是要求大于0，但是如果xLeave取得负数，也就是oneOverD为负数，背离StaticBox移动的那么将不会相交

    if (tEnter > tLeave){

     return kNoIntersection;

    }

 }

 if(d.y == 0.0f){

  if((stationaryBox.min.y > movingBox.max.y)||

    (stationaryBox.max.y < movingBox.min.y)){

    return kNoIntersection;

  }

 }else{

  float oneOverD = 1.0f / d.y;

  float yEnter = (stationaryBox.min.y - movingBox.max.y) * oneOverD;

  float yLeave = (stationaryBox.max.y - movingBox.min.y) * oneOverD;

  if (yEnter > yLeave){

   myswap(yEnter,yLeave);

  }

  if(yEnter > tEnter) tEnter = yEnter;

  if(yLeave < tLeave) tLeave = yLeave;

  // 如果yEnter很晚进入，tLeave还是之前小的;或者yLeave很早离开，tEnter还是原来的那么会

  // 产生tEnter > tLeave这样不会相交;

  if (tEnter > tLeave){

   return kNoIntersection;

  }

 }

 if (d.z == 0.0f){

  if((stationaryBox.min.z > movingBox.max.z)||

     (stationaryBox.max.z < movingBox.min.z)){

      return kNoIntersection;

  }

 }else{

  float oneOverD = 1.0f / d.z;

  float zEnter = (stationaryBox.min.z - movingBox.max.z) * oneOverD;

  float zLeave = (stationaryBox.max.z - movingBox.min.z) * oneOverD;

  if(zEnter > zLeave){

   myswap(zEnter,zLeave);

  }

  if(zEnter > tEnter) tEnter = zEnter;

  if(zLeave < tLeave) tLeave = zLeave;

  // 如果zEnter很晚进入，tLeave还是之前小的;或者zLeave很早离开，tEnter还是原来的那么会

  // 产生tEnter > tLeave这样不会相交;

  if(zEnter > zLeave) {

   return kNoIntersection;

  }

 }

 // 返回相交的时间

 return tEnter;

}