■3球面の交点を求める計算式(付:2円の交点)

 先に作成した
   「GPS簡易シミュレーションプログラム」
においては、3球面の交点として受信機(カーナビやケータイ)の位置を算出しています。ここでは3個の球面の交点の計算式を掲げます。また参考までに2円の交点の計算式についても付記します。

・3球面の交点
・Pi(xi, yi, zi)[i=1,2,3] を中心とする半径ri[i=1,2,3]の球の交点P(x,y,z)は次の3元連立方程式を
 解くことで求められます。
  (x1-x)2 + (y1-y)2 + (z1-z)2 = r12 ・・・(式1.1)
  (x2-x)2 + (y2-y)2 + (z2-z)2 = r22 ・・・(式1.2)
   (x3-x)2 + (y3-y)2 + (z3-z)2 = r32 ・・・(式1.3)
・(式1.2)および(式1.3)から(式1.1)を引いて整理すると
   (x2-x1)x + (y2-y1)y + (z2-z1)z = -(A2-A1)/2 ・・・(式2.1)
   (x3-x1)x + (y3-y1)y + (z3-z1)z = -(A3-A1)/2 ・・・(式2.2)
    ここで、
     Ai = ri2 - xi2 - yi2 - zi2  [i=1,2,3] 
  また、
   xij = xi-xj, yij = yi-yj, zij = zi-zj  [i=1,2,3]
      Ai1 = -(Ai-A1)/2                       [i=2,3]
  とおくと(式2.1)、(式2.2)は
   x21x + y21y + z21z = A21          ・・・(式3.1)
   x31x + y31y + z31z = A31          ・・・(式3.2)
・(式3.1)、(式3.2)より
   x = [(A21y31 - A31y21) + (y21z31 - y31z21)z] / D ・・・(式4.1)
   y = [(A31x21 - A21x31) + (x31z21 - x21z31)z] / D ・・・(式4.2)
  あるいは、
   x = B0 + B1z                   ・・・(式5.1)
   y = C0 + C1z                   ・・・(式5.2)
  ここで、
     D = x21y31 - y21x31   (≠0 とする、注1)
     B0 = (A21y31 - A31y21)/D, B1 = (y21z31 - y31z21)/D
     C0 = (A31x21 - A21x31)/D, C1 = (x31z21 - x21z31)/D
・(式5.1)、(式5.2)を(式1.1)に代入して整理すると、zに関する次の2次方程式が得られます。
   Ez2 + 2Fz + G = 0            ・・・(式6)
  ここで、
     E = B12 + C12 + 1
     F = B1(B0-x1) + C1(C0-y1) - z1
     G = (B0-x1)2 + (C0-y1)2 + z12 - r12
・この2次方程式を解くと((F2 - EG)1/2≧0 として)
   z = [ -F ± (F2 - EG)1/2 ] / E ・・・(式7)
 これを(式5.1)、(式5.2)に代入すると x、y が得られます。

(注1)D = 0の場合は別の手順で計算する必要があります(詳細省略)。
(注2)判別式:(F2 - EG) が負のとき、交点なし。
・2円の交点
・Pi(xi, yi)[i=1,2,3] を中心とする半径ri[i=1,2,3]の円の交点P(x,y)は次の2元連立方程式を
 解くことで求められます。
  (x1-x)2 + (y1-y)2 = r12 ・・・(式a1.1)
  (x2-x)2 + (y2-y)2 = r22 ・・・(式a1.2)
・(式a1.2)から(式a1.1)を引いて整理すると
   (x2-x1)x + (y2-y1)y = -(A2-A1)/2 ・・・(式a2)
    ここで、
     Ai = ri2 - xi2 - yi2   [i=1,2] 
  また、
   xij = xi-xj, yij = yi-yj  [i=1,2,3]
      A21 = -(A2-A1)/2
  とおくと(式a2)は
   x21x + y21y = A21    ・・・(式a3)
・(式a3)より(x21≠0 とする、注3)
   x = (A21 - y21y) / x21 ・・・(式a4)
  あるいは、
   x = B0 + B1y      ・・・(式a5)
  ここで、
     B0 = A21/x21
     B1 = -y21/x21
・(式a5)を(式a1.1)に代入して整理すると、yに関する次の2次方程式が得られます。
   Ey2 + 2Fy + G = 0       ・・・(式a6)
  ここで、
     E = B12 + 1
     F = B1(B0-x1) - y1
     G = (B0-x1)2 + y12 - r12
・この2次方程式を解くと((F2 - EG)1/2≧0 として)
   y = [ -F ± (F2 - EG)1/2 ] / E ・・・(式a7)
 これを(式a5)に代入すると x が得られます。

(注3)x21 = 0の場合は別の手順で計算する必要があります(詳細省略)
3球面の交点を求める簡単なアプリ (D = 0 の場合の処理方法も記載)

ホーム