orz

趣味のブログです

座標3点から、円の半径と中心を出すプログラム (python)

前置き

仕事の関係上、図を読み込んで特定の場所に円を描く必要がありました。
そこで、pythonで円の座標と半径を吐き出すプログラムを書いてみました。

プログラムの設計

設計指針ですが、次のHPに書いてあった式を参考にしてみました。
imagingsolution.blog.fc2.com

  • 任意の3点の座標を引数とする
  • (A,B,C)を求めるために, (逆行列) x (x1^2 + y1^2, x2 ... ) を掛け算する
  • 計算した答えの(A,B,C)から中心は(a/2, b/2) である
  • 半径は sqrt(C + (a/2)^2 + (b/2)^2)で求めることができる

コード

def radius(px,py,pz):
  #空の1次元配列用意
  vec = []
  #各座標のx^2 と y^2を計算している
  for i in px,py,pz:
    vec.append(i[0]**2 + i[1]**2)
  #print (vec)
  #空のリストを用意
  matrix = []
  #すべての座標に(0,0,1)を追加して、行列に代入
  for i in px,py,pz:
    i.append(1)
    matrix.append(i)
  #print (matrix)
  #Numpyで扱えるようにarray変換して、逆行列を求める
  matrix = np.array(matrix)
  inv_B = np.linalg.pinv(matrix)
  
  #逆行列と座標のスカラーベクトルをかけて、(a,b,c)を出す
  c = np.dot(inv_B,vec)
  c_sol = []
  #3行3行の行列になっていたので、各行の要素を足して、3行1列のベクトルに変換している
  for x in c:
    c_sol.append(np.sum(x))
  
  #c_solには(2a, 2b, r^2-a^2-b^2)が入っているので、2で割ると中心が座標がでる 
  a,b = c_sol[0]/2, c_sol[1]/2
  print ("円の中心(x,y)は", int(a),",",int(b))
  #
  radius = np.sqrt(c_sol[2] +a**2 + b**2)
  print ("円の半径は",radius)

出力例

p1 = (1,0) , p2 = (0,1) , p3 = (-1,0) として代入計算してみます。
p1 ~ p3は原点を中心とする単位円であることは容易にわかると思います。
検算として、やってみますと、

f:id:orshibuya0926:20210214135335p:plain

あっていますね。
たぶん、どんな円でも行けるかと思いますが、十分な検証はしておりません。