レイと三角形

ray-triangle.png

重心座標で三角形を表すと

(1)
\begin{equation} T(u,v,w) = uA+vB+wC \end{equation}

この(u,v,w)が重心座標です。
u+v+w=1です。
三角形の内側にある点の条件は、

(2)
\begin{align} 0 \leq u,v,w \leq 1 \end{align}

です。
三角形を定義するためにu+v+w=1と決まっているので、こんな風に式を変形できます。
u=1-u-v

(3)
\begin{equation} T(v,w)=A+v(B-A)+w(C-A) \end{equation}

線分を定義する二つの点P,Qがあるとします。

(4)
\begin{align} R(t) = P + t(Q-P), 0 \leq t \leq 1 \end{align}

この線分の式と、三角形の式を組み合わせると
線分と三角形の交点の条件は

(5)
\begin{align} T(v,w)=A+v(B-A)+w(C-A)=P+t(Q-P)\newline (P-Q)t+(B-A)v+(C-A)w=P-A \end{align}

知りたい変数,t,v,wを左辺に持っていく
この式を行列で表現すると

(6)
\begin{align} [(P-Q),(B-A),(C-A)] \begin{vmatrix} t\\ v\\ w \end{vmatrix} = [(P-A)] \end{align}

この式はクラメルの公式により、こんな風に解けます

(7)
\begin{align} t=\frac{det[(P-A),(B-A),(C-A)]}{det[(P-Q),(B-A),(C-A)]}\newline v=\frac{det[(P-A),(P-A),(C-A)]}{det[(P-Q),(B-A),(C-A)]}\newline w=\frac{det[(P-Q),(B-A),(P-A)]}{det[(P-Q),(B-A),(C-A)]} \end{align}

三重積より、

(8)
\begin{align} det[\vec{a} \vec{b} \vec{c}] = \vec{a} \cdot (\vec{b} \times \vec{c}) \end{align}

なのでこういう風になります

(9)
\begin{align} t=(P-A) \cdot \vec{n}/d \newline v=(C-A) \cdot \vec{e}/d \newline w = -(B-A) \cdot \vec{e}/d \end{align}

ここで
$\vec{n}=(B-A)\times (C-A)$つまり、nは三角形の法線です。
$d=(P-Q) \cdot \vec{n}$….線分と三角形の法線のなす角度を示しています。
$\vec{e} = (P-Q) \times (P-A)$

ソースコード

//pqは線分を示す2つの点
//abcは三角形を示す3つの点
//出力として、線分と三角形の交点を三角形の重心座標で表現します。
int IntersectSegmentTriangle(Point p,Point q,Point a,Point b,Point c,
 float &u,float &v,float &w,float &t){
   Vector ab=b-a;
   Vector ac=c-a;
   Vector qp=p-q;
   //三角形の法線を計算
   Vector n=Cross(ab,ac);
   //分母のdを計算 d<=0ならば、線分は三角形と平行か、三角形から逸れている
   float d=Dot(qp,n);
   if(d<=0.0f) return 0;
   //tの値を求めます
   Vector ap=p-a;
   t=Dot(ap,n);
   if(t<0.0f) return 0;
   if(t>d) return 0;//線分ではなく、線やレイなど無限に広がる線のテストの場合はここは省きます。
   //重心座標で、交点が三角形の内側に入っているかテストします
   Vector e=cross(qp,ap);
   v=Dot(ac,e);
   if(v<0.0f||v>d)return 0;
   w=-Dot(ab,e);
   if(w<0.0f||v+w>d) return 0;
 
   float ood=1.0f/d;
   t*=ood;
   v*=ood;
   w*=ood;
   u = 1.0f-v-w;
   return 1;
}

サポートサイト Wikidot.com