2次ベジェ曲線と任意点の距離を求める

icon 項目のみ表示/展開表示の切り替え

2次ベジェ曲線

2次ベジェ曲線と任意点との距離を求める方法を記述しています。
サンプルとしてマウスカーソルで任意の点を指定すると曲線に一番近い点を求め直線で結びます。
また、曲線の媒介変数の値と座標及び距離を表示することができます。
2次ベジェ曲線は、2点の始終点と1点の制御点より曲線を表現する方法で以下の式で表すことができます。
x(t)=(1-t)^2 \times x0+ 2(1-t) t \times x1+t^2 \times x2
y(t)=(1-t)^2 \times y0+ 2(1-t) t \times y1+t^2 \times y2
0 \leqq t \leqq 1
以下の展開公式を使用して展開しtについて整理し方程式に変形します。

(a-b)^2=a^2-2ab+b^2

x(t)=(1-2t+t^2) x_0 + 2 t \cdot x_1 -2t^2 x_1 + t^2 x_2
=x_0-2t \cdot x_0+t^2 x_0 + 2 t \cdot x_1 -2t^2 x_1 + t^2 x_2
= t^2 x_0 -2t^2 x_1 + t^2 x_2-2t \cdot x_0 +2 t \cdot x_1  + x_0
=(x_0 -2 x_1 +x_2) t^2  +(-2 x_0 +2 x_1) t + x_0

下図の丸い円上でマウスを左ボタンを押しながらカーソルを移動させると丸い円が移動します。左ボタンを離すと確定されます。
マウスカーソルの位置により任意の点を指定すると曲線への最短距離を求め直線で表します。
マウスカーソルを移動させると最短距離の場合の媒介変数の値、曲線上の座標、距離を表示します。
マウスがなくタッチ操作ができる環境では、下図下のラジオボタンで移動させる点を選択しその後図上でタップすると点を一度移動させることができます。
任意点の位置によっては最少の距離が2つ発生する場合がありますが、本ページでは媒介変数の値が小さい方のみ表示しています。





選択

点との距離

点(X0,Y0)と点(X1,Y1)との距離は、以下の式で表せます。
\sqrt{ (X_1-X_0)^2+(Y_1-Y_0)^2} 任意の点(Xp,Yp)とベジェ曲線x(t),y(t)との距離は、任意の点とベジェ曲線の距離が一番小さくなるtを求める必要があります。
計算を楽にするため任意点の座標を0,0になるようにベジェ曲線を移動させまます。
任意点の距離の2乗は、以下の式で表せます。
d=x(t)^2+y(t)^2
=\bigl( (x_0 -2 x_1 +x_2) t^2  +(-2 x_0 +2 x_1) t + x_0 \bigr) ^2+\bigl( (y_0 -2 y_1 +y_2) t^2  +(-2 y_0 +2 y_1) t + y_0 \bigr) ^2
各係数をa~hに置き換えると
a=x_0 -2 x_1 +x_2 , b=-2 x_0 +2 x_1 ,c=x_0
d=y_0 -2 y_1 +y_2 , e=-2 y_0 +2 y_1 ,f= y_0
r(t)=\bigl( a t^2+ b t+c \bigr) ^2+\bigl( d t^2+ e t+f \bigr) ^2
r(t)=\bigl( a^2 t^4 + 2a  b t^3 + b^2 t^2 + 2 a  c t^2 + 2 b  c t + c^2 \bigr ) + \bigl( d^2 t^4 + 2d  e t^3 + e^2 t^2 + 2 d  f t^2 + 2 e  f t + f^2 \bigr )
r(t)=a^2 t^4 + d^2 t^4 +2a  b t^3 + 2d  e t^3+ b^2 t^2 + 2 a  c t^2 + e^2 t^2++ 2 d  f t^2+ 2 b  c t      + 2 e  f t + c^2+ f^2
r(t)=(a^2 +d^2 )t^4 +  (2a  b+2d  e  )t^3 +  (b^2+ 2ac+e^2 +2 d  f ) t^2 + (2 b  c + 2 e  f  )t  + c^2+ f^2

微分

r'(t)=4(a^2 +d^2 )t^3 +  3(2a  b+2d  e  )t^2 +  2(b^2+ 2ac+e^2 +2 d  f ) t + (2 b  c + 2 e  f  )
微分値が0となる場合のtが最短距離の候補となります。
3次式ですので最大3個の解が得られます。
微分式の各項の係数をg~で置き換えます。
g=4(a^2+d^2) , h=6ab+6de , k=2b^2+4ac+2e^2+4df , m=2bc+2ef
gt^3+ht^2+kt+m=0
上の式のtを解の公式で求めます。

3次方程式の解(Cubic equation)
3次方程式の解の導出参照
ax^3+bx^2+cx+d=0\ (a \neq 0)
\displaystyle x_n=e \sqrt[3]{\frac{-2b^3+9abc-27a^2d}{54a^3}+ \frac{ \sqrt{ 3(- b^2c^2 +27a^2d^2 -18 abcd +4b^3d +4ac^3 )}}{18a^2}}
\displaystyle +f \sqrt[3]{\frac{-2b^3+9abc-27a^2d}{54a^3}- \frac{ \sqrt{ 3(- b^2c^2 +27a^2d^2 -18 abcd +4b^3d +4ac^3 )}}{18a^2}} -\frac{b}{3a}
n=0\ \ e=1,f=1
\displaystyle n=1\ \ e=\frac{-1 + \sqrt{3} i}{2},f=\frac{-1 - \sqrt{3} i}{2}
\displaystyle n=2\ \ e=\frac{-1 - \sqrt{3} i}{2},f=\frac{-1 + \sqrt{3} i}{2}

解の中で実数で0~1があればその中で一番距離が近いものを用います。
更に、0と1で距離が近い方を採用します。

ソースコード

本ページの基本的な部分のソースコードを以下に示す。
以下のソースのダウンロード(quadratic_bezier2.zip)
<!DOCTYPE html>
<html lang="ja">
<head>
  <meta charset="utf-8">
</head>
<body>
<span id="text2"></span><br />
<span id="text1"></span><br />
<span id="text3"></span><br />
<svg onmousemove="move()" xmlns="http://www.w3.org/2000/svg" viewBox="0 0 320 240" width="320" height="240">
  <g id="svg1">
    <rect x="0" y="0" width="320" height="240" stroke="blue" fill="none"/>
    <path id="dline1" d="M200,150L180,94" fill="none" stroke="black"/>
    <path id="b1" d="M50,60Q100,150 250,60" fill="none" stroke="black"/>
    <circle id="c1" cx="50" cy="60" r="8" fill="red" stroke="black" onmousedown="down1()" onmouseup="up1()"/>
    <circle id="c2" cx="100" cy="150" r="8" fill="red" stroke="black" onmousedown="down2()" onmouseup="up2()"/>
    <circle id="c3" cx="250" cy="60" r="8" fill="red" stroke="black" onmousedown="down3()" onmouseup="up3()"/>
    <circle id="cb" cx="200" cy="150" r="8" fill="green" stroke="black" onmousedown="down4()" onmouseup="up4()"/>
  </g>
</svg><br />
<form name="test1Form">
    <fieldset id="id1">
      <legend>選択</legend>
      <label><input name="name1" onclick="selB()" type="radio" value="1" />始点</label>
      <label><input name="name1" onclick="selB()" type="radio" value="2" />制御点</label>
      <label><input name="name1" onclick="selB()" type="radio" value="3" />終点</label>
      <label><input name="name1" onclick="selB()" type="radio" value="4" />任意点</label>
      <label><input name="name1" onclick="selB()" type="radio" checked="true" value="5" />フリー</label>
    </fieldset>
  </form>
<br />

<script type="text/javascript">
  var demo_run=false;
// ベジェ曲線の諸元
  var xyn=[	[ 50,60 ],
    		[ 100,150 ],
    		[250,60 ]];
  var pxy=[ 200,150 ];
// 3次方程式の解を求める
// ax^3+bx^2+cx+d=0のxを求める
  function cubic_equ(a,b,c,d){
    var p=-b*b/3/a/a+c/a;
    var q=2*b*b*b/27/a/a/a-b*c/3/a/a+d/a;
    var r2=81*q*q+12*p*p*p;
    var r_real=0;
    var r_img=0;
    if(r2<0){
      r_img=Math.sqrt(-r2);
    }else{
      r_real=Math.sqrt(r2);
    }
    var u3_real;
    var u3_img
    if(r_img){
      u3_real=-9*q/18;
      u3_img=r_img/18;
      v3_real=-9*q/18;
      v3_img=-r_img/18;
    }else{
      u3_real=(-9*q+r_real)/18;
      u3_img=0;
      v3_real=(-9*q-r_real)/18;
      v3_img=0;
    }
    var u_real;
    var u_img
    if(u3_img ){
      var z=Math.sqrt(u3_real*u3_real+u3_img*u3_img);
      var t=Math.atan2(u3_img,u3_real);
      z=Math.pow(z,1/3);
      t=t/3;
      u_real=z*Math.cos(t);
      u_img=z*Math.sin(t);
    }else{
      if(u3_real<0)
        u_real=-Math.pow(-u3_real,1/3);
      else
        u_real=Math.pow(u3_real,1/3);
      u_img=0;
    }
    var v_real;
    var v_img
    if(v3_img){
      var z=Math.sqrt(v3_real*v3_real+v3_img*v3_img);
      var t=Math.atan2(v3_img,v3_real);
      z=Math.pow(z,1/3);
      t=t/3;
      v_real=z*Math.cos(t);
      v_img=z*Math.sin(t);
    }else{
      if(v3_real<0)
        v_real=-Math.pow(-v3_real,1/3);
      else
        v_real=Math.pow(v3_real,1/3);
      v_img=0;
    }
    var omega1_real=-0.5;
    var omega1_img=Math.sqrt(3)/2;
    var omega2_real=-0.5;
    var omega2_img=-Math.sqrt(3)/2;
    var y0_real,y0_img;
    var y1_real,y1_img;
    var y2_real,y2_img;
    y0_real=u_real+v_real;
    y0_img=u_img+v_img;
    y1_real=omega1_real*u_real-omega1_img*u_img + omega2_real*v_real-omega2_img*v_img;
    y1_img=omega1_img*u_real+omega1_real*u_img + omega2_img*v_real+omega1_real*v_img
    y2_real=omega2_real*u_real-omega2_img*u_img + omega1_real*v_real-omega1_img*v_img;
    y2_img=omega2_img*u_real+omega2_real*u_img + omega1_img*v_real+omega1_real*v_img

    var x0_real,x0_img;
    var x1_real,x1_img;
    var x2_real,x2_img;
    x0_real=y0_real-b/(3*a);
    x0_img=y0_img;
    x1_real=y1_real-b/(3*a);
    x1_img=y1_img;
    x2_real=y2_real-b/(3*a);
    x2_img=y2_img;
    var xn=new Array(6);
    xn[0]=x0_real;
    xn[1]=x0_img;
    xn[2]=x1_real;
    xn[3]=x1_img;
    xn[4]=x2_real;
    xn[5]=x2_img;
    return xn;
  }
  // パラメータxynで媒介変数tのときの曲線上の点の座標を求める
  function bezier2(xyn,t){
    var xy=new Array(2);
    xy[0]=(1-t)*(1-t)*xyn[0][0]+2*(1-t)*t*xyn[1][0]+t*t*xyn[2][0];
    xy[1]=(1-t)*(1-t)*xyn[0][1]+2*(1-t)*t*xyn[1][1]+t*t*xyn[2][1];
    return xy;
  }

  function selB(){
    var rbs=document.forms["test1Form"]["name1"];
    var rb=event.target;
    var o=document.getElementById('id2');
    mode=eval(rb.value) | 0x80;
  }
  // マウスカーソルが移動しているときに呼び出される
  function move2(x,y){
    if(demo_run==true){
      var x0=xyn[0][0]-x;
      var y0=xyn[0][1]-y;
      var x1=xyn[1][0]-x;
      var y1=xyn[1][1]-y;
      var x2=xyn[2][0]-x;
      var y2=xyn[2][1]-y;
      var a=x0-2*x1+x2;
      var b=-2*x0+2*x1;
      var c=x0;
      var d=y0-2*y1+y2;
      var e=-2*y0+2*y1;
      var f=y0;
      var g=4*(a*a+d*d);
      var h=6*a*b+6*d*e;
      var k=2*b*b+4*a*c+2*e*e+4*d*f;
      var m=2*b*c+2*e*f;
      var tn=cubic_equ(g,h,k,m);
      var n;
      var tdn=0;
      var dmin=0;
      var i=0;
      var xy;
      for(n=0;n<3;n++){
        if(Math.abs(tn[n*2+1])<1e-14){ // 実数解の場合
          if(0<=tn[n*2] && tn[n*2]<=1){
            xy=bezier2(xyn,tn[n*2]);
            var dx=xy[0]-x;
            var dy=xy[1]-y;
            var d=dx*dx+dy*dy;
            if(i==0){
              dmin=d;
              tdn=tn[n*2];
            }else{
              if(d<dmin){
                dmin=d;
                tdn=tn[n*2];
              }
            }
            ++i;
          }
        }
      }
      if(i==0){ // 有効な実数解がない場合
        xy=bezier2(xyn,0);
        var dx=xy[0]-x;
        var dy=xy[1]-y;
        var d=dx*dx+dy*dy;
        dmin=d;
        tdn=0;
      }else{
        xy=bezier2(xyn,0);
        var dx=xy[0]-x;
        var dy=xy[1]-y;
        var d=dx*dx+dy*dy;
        if(dmin>d){
          dmin=d;
          tdn=0;
        }
      }
      xy=bezier2(xyn,1);
      dx=xy[0]-x;
      dy=xy[1]-y;
      d=dx*dx+dy*dy;
      if(dmin>d){
        dmin=d;
        tdn=1;
      }
      var dline = document.getElementById('dline1');
      xy=bezier2(xyn,tdn);
      dline.setAttribute("d","M"+x+","+y+"L"+xy[0]+","+xy[1]);
      var text1 = document.getElementById('text1');
      dmin=Math.sqrt(dmin);
      text1.innerText=x+","+y+" t="+tdn.toFixed(3)+" "+xy[0].toFixed(1)+","+xy[1].toFixed(1)+" d="+dmin.toFixed(1);
    }
    var c = document.getElementById('cb');
    c.setAttribute("cx",x);
    c.setAttribute("cy",y);
  }

  var mode=0;

  function down1(){
    mode=1;
  }
  function up1(){
    mode=0;
  }
  function down2(){
    mode=2;
  }
  function up2(){
    mode=0;
  }
  function down3(){
    mode=3;
  }
  function up3(){
    mode=0;
  }
  function down4(){
    mode=4;
  }
  function up4(){
    mode=0;
  }

  // マウスカーソルが移動しているときに呼び出される
  function move(){
    var s = document.getElementById('svg1');
    var r = s.getBoundingClientRect();
    var x=Math.round(event.clientX-r.left);
    var y=Math.round(event.clientY-r.top);
    event.preventDefault();
    if((mode & 0xf)==1){
      xyn[0][0]=x;
      xyn[0][1]=y;
      var b = document.getElementById('b1');
      b.setAttribute("d","M"+xyn[0][0]+","+xyn[0][1]+"Q"+xyn[1][0]+","+xyn[1][1]+" "+xyn[2][0]+","+xyn[2][1]);
      var c = document.getElementById('c1');
      c.setAttribute("cx",x);
      c.setAttribute("cy",y);
      move2(pxy[0],pxy[1]);
    }
    if((mode & 0xf)==2){
      xyn[1][0]=x;
      xyn[1][1]=y;
      var b = document.getElementById('b1');
      b.setAttribute("d","M"+xyn[0][0]+","+xyn[0][1]+"Q"+xyn[1][0]+","+xyn[1][1]+" "+xyn[2][0]+","+xyn[2][1]);
      var c = document.getElementById('c2');
      c.setAttribute("cx",x);
      c.setAttribute("cy",y);
      move2(pxy[0],pxy[1]);
    }
    if((mode & 0xf)==3){
      xyn[2][0]=x;
      xyn[2][1]=y;
      var b = document.getElementById('b1');
      b.setAttribute("d","M"+xyn[0][0]+","+xyn[0][1]+"Q"+xyn[1][0]+","+xyn[1][1]+" "+xyn[2][0]+","+xyn[2][1]);
      var c = document.getElementById('c3');
      c.setAttribute("cx",x);
      c.setAttribute("cy",y);
      move2(pxy[0],pxy[1]);
    }
    if((mode & 0xf)==4){
      pxy[0]=x;
      pxy[1]=y;
      move2(x,y);
    }
    if(mode & 0x80){
      mode=0;
      var rbs=document.forms["test1Form"]["name1"];
      rbs[4].checked=true;
    }
    var text2 = document.getElementById('text2');
    text2.innerHTML="始点   "+xyn[0][0]+","+xyn[0][1]+"<br />制御点 "+xyn[1][0]+","+xyn[1][1]+"<br />終点   "+xyn[2][0]+","+xyn[2][1];
  }

  window.onload = function(){
// ページ読み込み時に実行したい処理
    demo_run=true;
  }
</script>
</body>
</html>