2

現在、次の問題に悩まされています。

私は軌道データ（つまり、経度と緯度のデータ）を持っています。これを補間して線形フィッティングを見つけます（matlabでpolyfitとpolyvalを使用）。

私がやりたいのは、x 軸 (経度軸) が最適な線上にあるように軸を回転させることです。したがって、データはこの (回転された) 軸上にあるはずです。

私が試したのは、適合線の傾きから回転行列を評価することです (1 年生の多項式の式の m y=mx+q) 。

[cos(m) -sin(m);sin(m) cos(m)]

そして、元のデータにこの行列を掛けます...無駄に！

データが期待される x 軸上ではなく、中央にあるプロットを取得し続けます。

私は何が欠けていますか？

助けてくれてありがとう！

よろしくお願いします、

ウィンターミュート

4

2 に答える 2

1

いくつかのこと：

  1. 線形関数 y=mx+b がある場合、その線の角度は でありatan(m)、 ではありませんm。これらは小さなm', but very different for largem` の場合とほぼ同じです。

  2. 2+ 次の polyfit の線形成分は、1 次の polyfit の線形成分とは異なります。1 回は作業レベルで、もう 1 回は 1 次近似で、データを 2 回近似する必要があります。

  3. 勾配mが与えられた場合、三角関数を使用するよりも回転行列を計算するより良い方法があります (例: cos(atan(m)))。ジオメトリを実行するときは常に三角関数を避け、線形代数演算に置き換えるようにしています。これは通常より高速であり、特異点に関する問題が少なくなります。以下のコードを参照してください。

  4. この方法は、一部の軌道で問題を引き起こす可能性があります。たとえば、南北の軌道を考えてみましょう。しかし、それはより長い議論です。

説明した方法と上記の注意事項を使用して、これを実装するサンプル コードを次に示します。

%Setup some sample data
long = linspace(1.12020, 1.2023, 1000);
lat  = sin ( (long-min(long)) / (max(long)-min(long))*2*pi  )*0.0001 + linspace(.2, .31, 1000);

%Perform polynomial fit
p = polyfit(long, lat, 4);

%Perform linear fit to identify rotation
pLinear = polyfit(long, lat, 1);
m = pLinear(1);  %Assign a common variable for slope
angle = atan(m);

%Setup and apply rotation
%    Compute rotation metrix using trig functions
rotationMatrix = [cos(angle) sin(angle); -sin(angle) cos(angle)];

%    Compute same rotation metrix without trig
a = sqrt(m^2/(1+m^2));  %a, b are the solution to the system:    
b = sqrt(1/(1+m^2));    %    {a^2+b^2 = 1}, {m=a/b}
%                       %That is, the point (b,a) is on the unit
%                       circle, on a line with slope m
rotationMatrix = [b a; -a b];  %This matrix rotates the point (b,a) to (1,0)

%    Generally you rotate data after removing the mean value
longLatRotated = rotationMatrix * [long(:)-mean(long) lat(:)-mean(lat)]';

%Plot to confirm
figure(2937623)
clf

subplot(211)
hold on
plot(long, lat, '.')
plot(long, polyval(p, long), 'k-')
axis tight
title('Initial data')
xlabel('Longitude')
ylabel('Latitude')

subplot(212)
hold on;
plot(longLatRotated(1,:), longLatRotated(2,:),'.b-');
axis tight
title('Rotated data')
xlabel('Rotated x axis')
ylabel('Rotated y axis')
于 2013-03-11T18:40:39.583 に答える
0

回転行列で探している角度は、水平に対する線の角度です。これは、次の理由から、勾配の逆正接として見つけることができます。

tan(\theta) = Opposite/Adjacent = Rise/Run = slope

そのためt = atan(m)、線を回転して水平に戻したいことに注意して、回転行列を次のように定義します。

R = [cos(-t) sin(-t)
     sin(-t) cos(-t)]

これで、ポイントを回転させることができますR

于 2013-03-11T18:37:55.047 に答える