平面内における２線分どうしの交錯判定と交錯点計算について

やあ子どもたち。夏はいいよね。でも秋もいいよね。秋も深まってきたけど元気でやっているか。今日は平面内における２線分の交錯判定と、交錯点の計算について考えてみるよ。以下考察の中では一応外積とか使ってるから３次元じゃねーかとか一瞬思うかも知れないけど、単にzという数値量の符号を調べてるだけで、x,yしか使ってないわけだから基本的には２次元の話だよ。

●１．線分と線分との交錯判定についてまずは上の左側の図。平面内での２線分、ととの交錯判定は、「片方の線分の乗る直線からみて、もう片方の線分の両端点が、互いに異なる側にあるということを、両線分についてそれぞれ確認」すればよいだけなので外積の計算４回して符号調べるだけ。つまり、線分ABと線分CDの２線分について、以下の２式を同時に満たしていれば交錯、そうでなければ交錯していない。

●２．線分と線分との交錯点の計算について次に交錯点をどう計算するか。２つの線分が互いに交錯するかしないかについては、例えば上述のような交錯判定によりもうわかっているとすると、２線分の乗る２本の直線の交点を一発計算するだけ。つまり２次元の連立方程式を１本解けば求まる。

２次元なので、ｘとｙの成分に書き出して見ると、これは以下のような２ｘ２の行列方程式となる。

交錯することがわかっているなら左辺項の逆行列は必ず存在するので、これをベクトルpについて解いたものがずばり線分ABと、線分CDとの交錯点だ。
●３. プログラムにするとこんな感じ。

bool lc_check_intersection(
                           const float& ax, const float& ay,
                           const float& bx, const float& by,
                           const float& cx, const float& cy,
                           const float& dx, const float& dy
                           )
{
    return
    ((bx-ax)*(cy-by)-(by-ay)*(cx-bx)) * ((ax-bx)*(dy-ay)-(ay-by)*(dx-ax))>=0 &&
    ((dx-cx)*(by-dy)-(dy-cy)*(bx-dx)) * ((cx-dx)*(ay-cy)-(cy-dy)*(ax-cx))>=0;
}
bool lc_calc_intersection(
                          const float& ax, const float& ay,
                          const float& bx, const float& by,
                          const float& cx, const float& cy,
                          const float& dx, const float& dy,
                          float* intersection_point_x,
                          float* intersection_point_y
                          )
{
    const float perp_abx = by-ay;
    const float perp_aby = ax-bx;
    const float perp_cdx = dy-cy;
    const float perp_cdy = cx-dx;
    // ここでdetAの判定をするかどうかは実装者判断ということで。
    float det = perp_abx * perp_cdy - perp_aby * perp_cdx;
    float dot_x = ax*perp_abx + ay*perp_aby;
    float dot_y = cx*perp_cdx + cy*perp_cdy;
    *intersection_point_x = 
      ( perp_cdy * dot_x - perp_aby * dot_y ) / det;
    *intersection_point_y =
      ( -perp_cdx * dot_x + perp_abx * dot_y) / det;
    return 0;
}

一つ目の関数が交錯判定、２つ目の関数が交錯点計算。実際には後者の中から前者の判定を呼ぶ形にしてもいいかも。そしたらもっと短くなる。
さてテスト。冒頭右側の図。同じ色の線分どうしについての交錯テストと交錯した場合は交錯点計算。黄色い線分ペアから徐々に青い線分ペアになってくるとお互い離れちゃってくる感じだから交錯自体が起こらないことがわかる。交錯点には赤い棒を立てた。
ちなみにこんな実装もあるぞ。

void lc_calc_intersection
(const T& a, const T& b, const T& u, const T& v, float* s, float* t,float* det)
{
	// find intersection between
	//  a+su    b+tv
	// by calculating s, t as follows.
	// [s] =    1       [-vy vx    bx-ax]  
	// [ ]  ---------   [               ]
	// [t] = vxuy-uxvy  [-uy ux    by-ay]
	{
		*det = v.x*u.y - u.x*v.y;
		T dba = b - a;
		*s = (-v.y*dba.x + v.x*dba.y) / *det;
		*t = (-u.y*dba.x + u.x*dba.y) / *det;
	}
	return;
}

●４. パフォーマンスとその他ご参考
ここで、「１で交錯判定、２で交錯点計算とロジックを分けなくたって、いきなり２で２線分の乗る直線どうしの交錯を計算してしまって、交錯点の位置で、線分の上で起きている交錯かどうかを調べればよいのではないか」という素朴な疑問が沸いた諸兄もいることと思う。しかし、それは間違いだ。１の交錯判定が独立して存在しているのは、計算量が２に比べて少なくて済むため。実験では少なくとも20％以上の計算コストを節約できた。また、いきなり直線の交錯点計算をする場合は、平行判定をやらなくてはならず、線分領域判定も必要となり、実際には上述２よりもさらに計算コストのかかるロジックとなってしまう。これが１と２をそれぞれ分けて考えた理由だ。
話は変わるが、交錯判定をしたい線分セットが２本ではなく、もっとたくさんの複数線分からなる場合は、全く別の、SweepLineAlgorithmというものがあり、効率も遥かによいようなので、調べてみるのもいいだろう。
ときに黄色から青へのグラデーションって、なんか昔のザナドゥっぽくていいよね笑。では今回は以上の内容だよ。チャオ！