2直線の最短距離となる座標を算出したいです。

11 Ansichten (letzte 30 Tage)
ryo tanaka
ryo tanaka am 7 Okt. 2019
Kommentiert: ryo tanaka am 17 Okt. 2019
2直線の最短距離座標を図を参考に算出しました。
係数s、tはsymsを用いて定義しました。
しかし、データ数が627912個もあります。
matlabを走らせると計算時間が12時間ほどかかりました。
もっと簡単な算出方法、短時間で計算ができるようなスクリプトを作りたいのですが、中々上手くいきません。
良い方法があれば教えてほしいです。
無題.jpg

Antworten (4)

Kazuya
Kazuya am 8 Okt. 2019
627912点のデータか何のデータか分かりませんが、2直線の式をまず求めて
で距離を計算するのが計算負荷が小さそうですがどうでしょ。
  6 Kommentare
Kazuya
Kazuya am 8 Okt. 2019
画像も見えて状況がつかめてきました。Yoshio さんのいう通り matlabFunction を使う方法がよいのかも?
係数s、t を syms で定義して処理したコードも見せて頂くことは可能ですか?考えるベースにできればと思いまして。例えば上で頂いたデータに対して適用した形とか。。
ryo tanaka
ryo tanaka am 10 Okt. 2019
回答ありがとうございます。
返信が遅くなり、大変申し訳ありませんでした。
係数s、t を syms で定義して処理したコードは添付画像通りに算出していますので、、
コードを添付した方がよろしいでしょうか?

Melden Sie sich an, um zu kommentieren.


Yoshio
Yoshio am 8 Okt. 2019
もし、最短距離だけを求めるのであれば、
から、,, とすればベクトル化した計算式ができそうです。
  6 Kommentare
Yoshio
Yoshio am 10 Okt. 2019
Bearbeitet: Yoshio am 10 Okt. 2019
確認ありがとうございます。O = [25, 0, 0]で良かったです。
何回か投稿を編集されているようで、私が見たときは最短距離だけ出ればよいかと思ったのですが、座標も必要なのですね。
その場合、t,s を求める計算をv1, v2, P12等で記述できれば高速化可能と思います。
ryo tanaka
ryo tanaka am 10 Okt. 2019
何度か投稿の編集を行っていましたので、伝えきれていなかったと思います。
失礼いたしました。
なるほど、t、sの記述方法を変えることで計算時間の短縮ができそうですね。

Melden Sie sich an, um zu kommentieren.


Yoshio
Yoshio am 8 Okt. 2019
添付の図が見えませんので、良く分からない所がありますが、シンボリックが解をmatlabFunctionに変換できるので、一度数式で得られた解を数値計算向けの関数にして利用するのはどうでしょうか? 一度試してみてください。
  1 Kommentar
ryo tanaka
ryo tanaka am 8 Okt. 2019
回答ありがとうございます。添付した図が見れず申し訳ありません。修正しました。
matlabFunctionは知りませんでした。一度試してみようと思います。ありがとうございます。

Melden Sie sich an, um zu kommentieren.


Yoshio
Yoshio am 10 Okt. 2019
Bearbeitet: Yoshio am 10 Okt. 2019
直線モデルを以下のように置く
とすると、最短距離の条件から
よって
これをについて解くと
ここで
上記に基づき、ベクトル化を考慮してプログラムすると
% B,O,C,A = [x座標,y座標,z座標]
B=[39.8125,-0.9625,-198;39.9875,-1.8375,-198;39.9875,-0.9625,-198;39.9875,-0.7875,-198;40.1625,-7.2625,-198];
A=[-41.2125,-0.4375,-198;-41.2125,-0.2625,-198;-41.0375,-6.2125,-198;-41.0375,-0.9625,-198;-41.0375,-0.7875,-198];
C=[25,0,0];
O=[-25,0,0];
A = A';
B = B';
C = C';
O = O';
% A = rand(3,972);
% B = rand(3,646);
AA = repmat(A,1,size(B,2));
BB = reshape(repmat(B,size(A,2),1),3,size(B,2)*size(A,2));
v1 = AA-O;
v2 = BB-C;
P12 = O-C;
v1xv2 = cross(v1,v2);
d = abs(P12'*v1xv2./vecnorm(v1xv2));
tic
x1 = O;
x2 = C;
x12 = x1-x2;
v1v1 = sum(v1.*v1);
v1v2 = sum(v1.*v2);
v2v2 = sum(v2.*v2);
D = -v1v1.*v2v2 + v1v2.^2;
b = -[sum(v1.*x12);sum(v2.*x12)];
t1 = (-v2v2.*b(1,:)+v1v2.*b(2,:))./D;
t2 = (-v1v2.*b(1,:)+v1v1.*b(2,:))./D;
P1 = x1+t1.*v1;
P2 = x2+t2.*v2;
d0 = vecnorm(P1-P2);
toc
max(abs(d0-d))
min(abs(D))
972×646=627912の組だと約0.2秒でした
dとd0の差はepsの10倍程なので、合っているかと思います。なお、2つの直線が全く平行の場合は、D=0となり、解は不定となるので、この場合は別途考慮する必要がありますね。
  8 Kommentare
Yoshio
Yoshio am 16 Okt. 2019
[1点目のAO直線に対するBC直線1~end、2点目のAO直線に対するBC直線1~end、・・・、end点目のAO直線に対するBC直線1~end]
という順番でならんでいる
という記述の意味がよくわからないのですが、MATLABの良さは、アイディアを簡単に試せる点にありますので、サンプルデータを作って、確認していただくのが良いかと思います。
ryo tanaka
ryo tanaka am 17 Okt. 2019
質問の仕方が悪かったです。すみません。
簡単なサンプルデータで確認してみます。

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu 数学 finden Sie in Help Center und File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!