任意の3点を通る球体の中心の求め方

回答(6)です。

下のほうで求めた中心座標の検証に使った、CAD(Rhinoceros)での求めかたの例も書いておきます。

1)3点を描く。
http://mcnc.hp.infoseek.co.jp/cgi-bin/img-box/img20100217062901.jpg

2)3点と中心とする球を描く。
http://mcnc.hp.infoseek.co.jp/cgi-bin/img-box/img20100217062916.jpg

3)球の交線の交点が答え。
http://mcnc.hp.infoseek.co.jp/cgi-bin/img-box/img20100217062932.jpg

簡単に求められて、確認も目視で分かりやすいですが、CADが高価だった一昔前でしたらとても贅沢な手段ですね。

Excelのソルバーは使ったことがなかったので、私も勉強になりました。
ついでに、任意の4点を通る球の中心と半径も、同じ方法で求められることを確認しました。
きっかけをいただき、ありがとうございました。

お礼

2010/03/10 13:17

返答遅くなりました。回答何度もありがとうございます。
ソルバーまで検証して頂きありがとうございます。
私は2Ｄ-ＣＡＤしかなく3Ｄ検証が出来なく困っていました。

質問者

回答(5)です。

>ソルバーを使ってもなぜかうまくいきません。

ソルバーも試してみました。

質問で紹介された参考URLの内容そのままではダメで、
制約条件を設定したらうまくいきました。

仮に、中心に3点の最小値を入れて計算させると、
http://mcnc.hp.infoseek.co.jp/cgi-bin/img-box/img20100216221739.jpg

中心座標1が求まります。
http://mcnc.hp.infoseek.co.jp/cgi-bin/img-box/img20100216221804.jpg

同様に、中心に3点の最大値を入れて計算させると、中心座標2が求まります。
さすがExcel、すばらしいですね。

回答(3)です。回答(1)さんへ、

>青軸上ならどこでもかける
>
>任意の3点を通る最小の球体の中心　って縛りをつけないと

球の半径rを指定するので2つ以内に決まりますよ。
こんな感じ。
http://mcnc.hp.infoseek.co.jp/cgi-bin/img-box/img20100216073946.jpg

もう、殆ど答えは出ていますが、中心点Ｄ（a,b,c)で半径ｒの円を考えて、
点Ａを通る場合と点Ｂを通る場合と、点Ｃを通る場合とで式をたてると、
求まります。

お礼

2010/03/10 13:15

返答遅くなりました。回答ありがとうございます。
計算するとなると大変ですね。

質問者

下記URLで、回答(2)さんの内容からもうちょっと先のところまでお話がされています。

(BIGLOBEなんでも相談室)
http://soudan1.biglobe.ne.jp/qa195295.html

Perlで確認してみました。
-----
# 3点
my($x1,$y1,$z1)=(-51,-21,-11);
my($x2,$y2,$z2)=( 62, 12, 22);
my($x3,$y3,$z3)=(-73, 83, 33);
# 半径
my $r1=100;

my $tm01=$x1**2-$x2**2+$y1**2-$y2**2+$z1**2-$z2**2;
my $tm02=$x1**2-$x3**2+$y1**2-$y3**2+$z1**2-$z3**2;
my $tm11=-2*($x1-$x2)*($z1-$z3)+2*($x1-$x3)*($z1-$z2);
my $tm12=-2*($y1-$y2)*($z1-$z3)+2*($y1-$y3)*($z1-$z2);
my $tm13=$tm01*($z1-$z3)-$tm02*($z1-$z2);
my $tm21=-2*($x1-$x2)*($y1-$y3)+2*($x1-$x3)*($y1-$y2);
my $tm22=-2*($z1-$z2)*($y1-$y3)+2*($z1-$z3)*($y1-$y2);
my $tm23=$tm01*($y1-$y3)-$tm02*($y1-$y2);
my $tma=1+$tm11**2/$tm12**2+$tm21**2/$tm22**2;
my $tmb=-2*$x1+2*($y1+$tm13/$tm12)*$tm11/$tm12+2*($z1+$tm23/$tm22)*$tm21/$tm22;
my $tmc=$x1**2+($y1+$tm13/$tm12)**2+($z1+$tm23/$tm22)**2-$r1**2;
my $xq1=(-$tmb+sqrt($tmb**2-4*$tma*$tmc))/2/$tma;
my $xq2=(-$tmb-sqrt($tmb**2-4*$tma*$tmc))/2/$tma;
my $yq1=-$tm13/$tm12-$tm11/$tm12*$xq1;
my $yq2=-$tm13/$tm12-$tm11/$tm12*$xq2;
my $zq1=-$tm23/$tm22-$tm21/$tm22*$xq1;
my $zq2=-$tm23/$tm22-$tm21/$tm22*$xq2;

printf("xq1=%9.3f yq1=%9.3f yq1=%9.3f\n",$xq1,$yq1,$zq1);
printf("xq2=%9.3f yq2=%9.3f yq2=%9.3f\n",$xq2,$yq2,$zq2);
-----

出た答え。
-----
xq1= -2.410 yq1= 63.052 yq1= -34.965
xq2= -20.432 yq2= 11.190 yq2= 78.607
-----

CADで作図したものと一致したので大丈夫と思います。

お疲れさまでした。

すみません、間違いがありました。
後ろのほうです。

+++++ 誤 +++++
printf("xq1=%9.3f yq1=%9.3f yq1=%9.3f\n",$xq1,$yq1,$zq1);
printf("xq2=%9.3f yq2=%9.3f yq2=%9.3f\n",$xq2,$yq2,$zq2);
-----

出た答え。
-----
xq1= -2.410 yq1= 63.052 yq1= -34.965
xq2= -20.432 yq2= 11.190 yq2= 78.607
-----
++++++++++++++

+++++ 正 +++++
printf("xq1=%9.3f yq1=%9.3f zq1=%9.3f\n",$xq1,$yq1,$zq1);
printf("xq2=%9.3f yq2=%9.3f zq2=%9.3f\n",$xq2,$yq2,$zq2);
-----

出た答え。
-----
xq1= -2.410 yq1= 63.052 zq1= -34.965
xq2= -20.432 yq2= 11.190 zq2= 78.607
-----
++++++++++++++