読者です 読者をやめる 読者になる 読者になる

線分を扱うプログラムはベクトルで解く

プログラミング

最近、別々の人から同じ質問をされたのでエントリーにまとめることにします。
その質問というのが「線分ABと点Pが与えられたとき、AB上でもっともPに近い点を求めるには?」というもの。
垂線をおろして交点を求めるだけの簡単なプログラムのように思えて、これはちょっと工夫が必要です。

誰にでも思いつくナイーブな解法

垂線をおろして交点を求めればいいわけです。もし交点がなければ線分の端点AかBのどちらかが「最も近い点」になるはず。
実際に手順を書いてみましょう。

  1. ABの傾きaを求める。
  2. 垂線の傾きは -1/a。ただしABが垂直なら垂線の傾きは0。垂線の傾きをbとする。
  3. 点Pを通り傾きがbとなる直線を求める。【一次方程式を解く】
  4. ABを直線の式で表し、垂線との交点を求める。【連立一次方程式を解く】
  5. 交点の座標がAとBの間にあるなら、交点が求める点。
  6. 交点の座標がAの外側ならAが求める点、Bの外側ならBが求める点。

なんと連立一次方程式を解く必要があります。これは面倒くさいですね。
しかもこの方式、平面図形の場合にしか使えません。空間図形になるとそもそも「垂線の傾き」という概念がないですから。
ここで役に立つのが、高校で習ったベクトルです。

ベクトルを使ったスマートな解法

ベクトルxの長さのことを|x|と表現します。まずは考え方。

  • 線分ABを、基点Aとベクトルaで表す。AからPまでのベクトルをbとする。
  • aとbのなす角をθとすると内積abは、|a||b|cosθ。これを|a|で割ると、線分AB上に下ろした交点までのAからの距離を意味します。
  • |b|cosθが 0 以下なら A が求める点、|a|以上ならBが求める点、0より大きく|a|未満だったら、A + a・(|b|cosθ/|a|)が求める点です。

これを実際にプログラムにしてみましょう。コツとして、ベクトルの長さ|a|を求めるには平方根計算が必要になるので、|a|2のまま扱った方が効率がいいという点に留意します(|a|2を出すだけなら内積a・aを求めるだけです)。

  1. ベクトルaはB-A。ベクトルbはP-A。
  2. 内積iはab。交点までの距離比率rはi/|a|2
  3. rが0以下ならAを返す、1以上ならBを返す、それ以外ならA + r・aを返す。

C言語で書き下ろします*1。二次元ベクトルを表すvector2Dという型を用意して使っていますが、オブジェクト指向言語なら内積を求めるメソッド・引き算をするメソッドを用意することで三次元など一般的な次元にも対応可能。

typedef struct{
	double x;
	double y;
} vector2D;

vector2D nearest( vector2D A, vector2D B, vector2D P ){
	vector2D a, b;
	double r;
	
	a.x = B.x - A.x;
	a.y = B.y - A.y;
	b.x = P.x - A.x;
	b.y = P.y - A.y;
	
	r = (a.x*b.x + a.y*b.y) / (a.x*a.x + a.y*a.y);
	
	if( r<= 0 ){
		return A;
	}else if( r>=1 ){
		return B;
	}else{
		vector2D result;
		result.x = A.x + r*a.x;
		result.y = A.y + r*a.y;
		return result;
	}
}

まとめ

線分上で最も近い点を求める必要があるなんてそんな頻繁にはないことです。だからこの解法を覚えておく必要はありません。
大事なのは

  • 図形のプログラムはベクトルですんなり書ける場合がある。
  • 高校では、かなり実用的なことを習っている。

と覚えておくことです。

*1:手元にコンパイラがないので、書き下ろしただけです。見たところ間違いないけど、エラーがあったら教えてください