一連の循環データの平均を計算したいと思います。たとえば、コンパスの読み取りからいくつかのサンプルがあるかもしれません。もちろん問題は、ラップアラウンドにどのように対処するかです。同じアルゴリズムが文字盤にも役立つ場合があります。
実際の質問はもっと複雑です-統計は球または「ラップアラウンド」する代数空間、たとえば加法群modnで何を意味しますか。答えは一意ではない場合があります。たとえば、359度と1度の平均は0度または180度である可能性がありますが、統計的には0の方が適切です。
これは私にとって本当のプログラミングの問題であり、数学の問題のように見えないようにしようとしています。
角度から単位ベクトルを計算し、それらの平均の角度を取得します。
問題はわかります。たとえば、45 フィートの角度と 315 フィートの角度がある場合、「自然な」平均は 180 フィートですが、必要な値は実際には 0 フィートです。
スターブルーは何かに取り組んでいると思います。各角度の (x, y) デカルト座標を計算し、結果のベクトルを足し合わせるだけです。最終的なベクトルの角度オフセットは、必要な結果になるはずです。
x = y = 0
foreach angle {
x += cos(angle)
y += sin(angle)
}
average_angle = atan2(y, x)
コンパスの方位は北から始まり、時計回りに進むのに対し、「通常の」デカルト座標は X 軸に沿ってゼロから始まり、反時計回りに進むことを今のところ無視しています。数学は関係なく同じように機能するはずです。
2 つの角度の特別な場合:
答え((a + b) mod 360 ) / 2は間違っています。角度 350 と 2 の場合、最も近い点は 176 ではなく 356 です。
単位ベクトルと三角関数の解は高すぎる可能性があります。
少しいじって得たものは次のとおりです。
diff = ( ( a - b + 180 + 360 ) mod 360 ) - 180
angle = (360 + b + ( diff / 2 ) ) mod 360
ackb の言うとおり、これらのベクトル ベースのソリューションは角度の真の平均とは見なされず、対応する単位ベクトルの平均にすぎません。ただし、ackb の提案する解決策は、数学的には正しくないようです。
以下は、(angle[i] - avgAngle)^2 (必要に応じて差が修正される) を最小化するという目標から数学的に導き出されたソリューションであり、角度の真の算術平均になります。
最初に、角度間の差が通常の数の対応物の差と異なる場合を正確に調べる必要があります。角度 x と y を考えます。y >= x - 180 かつ y <= x + 180 の場合、差 (xy) を直接使用できます。それ以外の場合、最初の条件が満たされない場合は、計算で y の代わりに (y+360) を使用する必要があります。同様に、2 番目の条件が満たされない場合は、y の代わりに (y-360) を使用する必要があります。最小化する曲線の方程式は、これらの不等式が true から false に、またはその逆に変化する点でのみ変化するため、[0,360) の範囲全体をこれらの点で区切られた一連のセグメントに分割できます。次に、これらの各セグメントの最小値と、各セグメントの最小値の最小値 (平均) を見つけるだけです。
これは、角度差の計算で問題が発生する場所を示す画像です。x が灰色の領域にある場合、問題が発生します。
変数を最小化するには、曲線に応じて、最小化したいものの導関数を取得してから、転換点 (導関数 = 0 の場所) を見つけることができます。
ここでは、差の二乗を最小化するという考え方を適用して、一般的な算術平均式 sum(a[i])/n を導き出します。曲線 y = sum((a[i]-x)^2) は、次の方法で最小化できます。
y = sum((a[i]-x)^2)
= sum(a[i]^2 - 2*a[i]*x + x^2)
= sum(a[i]^2) - 2*x*sum(a[i]) + n*x^2
dy\dx = -2*sum(a[i]) + 2*n*x
for dy/dx = 0:
-2*sum(a[i]) + 2*n*x = 0
-> n*x = sum(a[i])
-> x = sum(a[i])/n
調整された差異を使用して曲線に適用します。
b = 正しい (角度) 差 a[i]-x c = 正しい (角度) 差 (a[i]-360)-x である a のサブセット cn = c のサイズ d = a のサブセット正しい (角度) 差 (a[i]+360)-x dn = d のサイズ
y = sum((b[i]-x)^2) + sum(((c[i]-360)-b)^2) + sum(((d[i]+360)-c)^2)
= sum(b[i]^2 - 2*b[i]*x + x^2)
+ sum((c[i]-360)^2 - 2*(c[i]-360)*x + x^2)
+ sum((d[i]+360)^2 - 2*(d[i]+360)*x + x^2)
= sum(b[i]^2) - 2*x*sum(b[i])
+ sum((c[i]-360)^2) - 2*x*(sum(c[i]) - 360*cn)
+ sum((d[i]+360)^2) - 2*x*(sum(d[i]) + 360*dn)
+ n*x^2
= sum(b[i]^2) + sum((c[i]-360)^2) + sum((d[i]+360)^2)
- 2*x*(sum(b[i]) + sum(c[i]) + sum(d[i]))
- 2*x*(360*dn - 360*cn)
+ n*x^2
= sum(b[i]^2) + sum((c[i]-360)^2) + sum((d[i]+360)^2)
- 2*x*sum(x[i])
- 2*x*360*(dn - cn)
+ n*x^2
dy/dx = 2*n*x - 2*sum(x[i]) - 2*360*(dn - cn)
for dy/dx = 0:
2*n*x - 2*sum(x[i]) - 2*360*(dn - cn) = 0
n*x = sum(x[i]) + 360*(dn - cn)
x = (sum(x[i]) + 360*(dn - cn))/n
これだけでは最小値を取得するのに十分ではありませんが、制限のないセットを持つ通常の値では機能するため、結果は間違いなくセットの範囲内にあり、したがって有効です。範囲内の最小値が必要です (セグメントによって定義されます)。最小値がセグメントの下限よりも小さい場合、そのセグメントの最小値は下限にある必要があり (二次曲線には転換点が 1 つしかないため)、最小値がセグメントの上限よりも大きい場合、セグメントの最小値は上界。各セグメントの最小値が得られたら、最小化する値が最も低いセグメントを見つけます (sum((b[i]-x)^2) + sum(((c[i]-360 )-b)^2) + sum(((d[i]+360)-c)^2))。
これは曲線のイメージで、x=(a[i]+180)%360 の点で曲線がどのように変化するかを示しています。問題のデータセットは {65,92,230,320,250} です。
これは Java でのアルゴリズムの実装で、いくつかの最適化を含みます。その複雑さは O(nlogn) です。比較ベースのソートを基数ソートなどの非比較ベースのソートに置き換えると、O(n) に減らすことができます。
static double varnc(double _mean, int _n, double _sumX, double _sumSqrX)
{
return _mean*(_n*_mean - 2*_sumX) + _sumSqrX;
}
//with lower correction
static double varlc(double _mean, int _n, double _sumX, double _sumSqrX, int _nc, double _sumC)
{
return _mean*(_n*_mean - 2*_sumX) + _sumSqrX
+ 2*360*_sumC + _nc*(-2*360*_mean + 360*360);
}
//with upper correction
static double varuc(double _mean, int _n, double _sumX, double _sumSqrX, int _nc, double _sumC)
{
return _mean*(_n*_mean - 2*_sumX) + _sumSqrX
- 2*360*_sumC + _nc*(2*360*_mean + 360*360);
}
static double[] averageAngles(double[] _angles)
{
double sumAngles;
double sumSqrAngles;
double[] lowerAngles;
double[] upperAngles;
{
List<Double> lowerAngles_ = new LinkedList<Double>();
List<Double> upperAngles_ = new LinkedList<Double>();
sumAngles = 0;
sumSqrAngles = 0;
for(double angle : _angles)
{
sumAngles += angle;
sumSqrAngles += angle*angle;
if(angle < 180)
lowerAngles_.add(angle);
else if(angle > 180)
upperAngles_.add(angle);
}
Collections.sort(lowerAngles_);
Collections.sort(upperAngles_,Collections.reverseOrder());
lowerAngles = new double[lowerAngles_.size()];
Iterator<Double> lowerAnglesIter = lowerAngles_.iterator();
for(int i = 0; i < lowerAngles_.size(); i++)
lowerAngles[i] = lowerAnglesIter.next();
upperAngles = new double[upperAngles_.size()];
Iterator<Double> upperAnglesIter = upperAngles_.iterator();
for(int i = 0; i < upperAngles_.size(); i++)
upperAngles[i] = upperAnglesIter.next();
}
List<Double> averageAngles = new LinkedList<Double>();
averageAngles.add(180d);
double variance = varnc(180,_angles.length,sumAngles,sumSqrAngles);
double lowerBound = 180;
double sumLC = 0;
for(int i = 0; i < lowerAngles.length; i++)
{
//get average for a segment based on minimum
double testAverageAngle = (sumAngles + 360*i)/_angles.length;
//minimum is outside segment range (therefore not directly relevant)
//since it is greater than lowerAngles[i], the minimum for the segment
//must lie on the boundary lowerAngles[i]
if(testAverageAngle > lowerAngles[i]+180)
testAverageAngle = lowerAngles[i];
if(testAverageAngle > lowerBound)
{
double testVariance = varlc(testAverageAngle,_angles.length,sumAngles,sumSqrAngles,i,sumLC);
if(testVariance < variance)
{
averageAngles.clear();
averageAngles.add(testAverageAngle);
variance = testVariance;
}
else if(testVariance == variance)
averageAngles.add(testAverageAngle);
}
lowerBound = lowerAngles[i];
sumLC += lowerAngles[i];
}
//Test last segment
{
//get average for a segment based on minimum
double testAverageAngle = (sumAngles + 360*lowerAngles.length)/_angles.length;
//minimum is inside segment range
//we will test average 0 (360) later
if(testAverageAngle < 360 && testAverageAngle > lowerBound)
{
double testVariance = varlc(testAverageAngle,_angles.length,sumAngles,sumSqrAngles,lowerAngles.length,sumLC);
if(testVariance < variance)
{
averageAngles.clear();
averageAngles.add(testAverageAngle);
variance = testVariance;
}
else if(testVariance == variance)
averageAngles.add(testAverageAngle);
}
}
double upperBound = 180;
double sumUC = 0;
for(int i = 0; i < upperAngles.length; i++)
{
//get average for a segment based on minimum
double testAverageAngle = (sumAngles - 360*i)/_angles.length;
//minimum is outside segment range (therefore not directly relevant)
//since it is greater than lowerAngles[i], the minimum for the segment
//must lie on the boundary lowerAngles[i]
if(testAverageAngle < upperAngles[i]-180)
testAverageAngle = upperAngles[i];
if(testAverageAngle < upperBound)
{
double testVariance = varuc(testAverageAngle,_angles.length,sumAngles,sumSqrAngles,i,sumUC);
if(testVariance < variance)
{
averageAngles.clear();
averageAngles.add(testAverageAngle);
variance = testVariance;
}
else if(testVariance == variance)
averageAngles.add(testAverageAngle);
}
upperBound = upperAngles[i];
sumUC += upperBound;
}
//Test last segment
{
//get average for a segment based on minimum
double testAverageAngle = (sumAngles - 360*upperAngles.length)/_angles.length;
//minimum is inside segment range
//we test average 0 (360) now
if(testAverageAngle < 0)
testAverageAngle = 0;
if(testAverageAngle < upperBound)
{
double testVariance = varuc(testAverageAngle,_angles.length,sumAngles,sumSqrAngles,upperAngles.length,sumUC);
if(testVariance < variance)
{
averageAngles.clear();
averageAngles.add(testAverageAngle);
variance = testVariance;
}
else if(testVariance == variance)
averageAngles.add(testAverageAngle);
}
}
double[] averageAngles_ = new double[averageAngles.size()];
Iterator<Double> averageAnglesIter = averageAngles.iterator();
for(int i = 0; i < averageAngles_.length; i++)
averageAngles_[i] = averageAnglesIter.next();
return averageAngles_;
}
一連の角度の算術平均は、平均がどうあるべきかについての直感的な考えと一致しない場合があります。たとえば、セット {179,179,0,181,181} の算術平均は 216 (および 144) です。すぐに思いつく答えはおそらく 180 ですが、算術平均はエッジの値に大きく影響されることはよく知られています。また、角度はベクトルではないことも覚えておく必要があります。
もちろん、このアルゴリズムは、時刻など、剰余演算 (最小限の調整で) に従うすべての量にも適用されます。
また、ベクトル ソリューションとは異なり、これが角度の真の平均であっても、必ずしもそれが使用すべきソリューションであることを意味するわけではなく、対応する単位ベクトルの平均が実際の値である可能性があることも強調したいと思います。使用する必要があります。
平均をより正確に定義する必要があります。2 つの角度の特定のケースでは、2 つの異なるシナリオを考えることができます。
ただし、2 つ以上の角度の場合に 2 番目の選択肢を一般化する方法がわかりません。
浮動小数点または三角関数の機能を持たないマイクロコントローラーで使用した方法を共有したいと思います。変動を滑らかにするために、10 個の生の方位測定値を「平均化」する必要がありました。
それは理想的ではありません。それは壊れる可能性があります。この場合、デバイスは非常にゆっくりとしか回転しないため、うまくいきました。他の誰かが同様の制限の下で働いていることに気付いた場合に備えて、そこに掲載します。
すべての平均と同様に、答えはメトリックの選択によって異なります。所定のメトリック M の場合、[1,N] の k に対する [-pi,pi] のいくつかの角度 a_k の平均は、距離の二乗和 d^2_M(a_M,a_k) を最小化する角度 a_M です。加重平均の場合、単純に合計に重み w_k を含めます (sum_k w_k = 1 となるように)。あれは、
a_M = arg min_x sum_k w_k d^2_M(x,a_k)
メトリックの 2 つの一般的な選択肢は、フロベニウス メトリックとリーマン メトリックです。フロベニウス メトリックには、循環統計における平均方位の通常の概念に対応する直接式が存在します。詳細については、Maher Moakher 著、SIAM Journal on Matrix Analysis and Applications、第 24 巻、第 1 号、2002 年、「回転のグループにおける平均と平均化」を参照してください。
http://link.aip.org/link/?SJMAEL/24/1/1
以下は、計算を行う GNU Octave 3.2.4 の関数です。
function ma=meanangleoct(a,w,hp,ntype)
% ma=meanangleoct(a,w,hp,ntype) returns the average of angles a
% given weights w and half-period hp using norm type ntype
% Ref: "Means and Averaging in the Group of Rotations",
% Maher Moakher, SIAM Journal on Matrix Analysis and Applications,
% Volume 24, Issue 1, 2002.
if (nargin<1) | (nargin>4), help meanangleoct, return, end
if isempty(a), error('no measurement angles'), end
la=length(a); sa=size(a);
if prod(sa)~=la, error('a must be a vector'); end
if (nargin<4) || isempty(ntype), ntype='F'; end
if ~sum(ntype==['F' 'R']), error('ntype must be F or R'), end
if (nargin<3) || isempty(hp), hp=pi; end
if (nargin<2) || isempty(w), w=1/la+0*a; end
lw=length(w); sw=size(w);
if prod(sw)~=lw, error('w must be a vector'); end
if lw~=la, error('length of w must equal length of a'), end
if sum(w)~=1, warning('resumming weights to unity'), w=w/sum(w); end
a=a(:); % make column vector
w=w(:); % make column vector
a=mod(a+hp,2*hp)-hp; % reduce to central period
a=a/hp*pi; % scale to half period pi
z=exp(i*a); % U(1) elements
% % NOTA BENE:
% % fminbnd can get hung up near the boundaries.
% % If that happens, shift the input angles a
% % forward by one half period, then shift the
% % resulting mean ma back by one half period.
% X=fminbnd(@meritfcn,-pi,pi,[],z,w,ntype);
% % seems to work better
x0=imag(log(sum(w.*z)));
X=fminbnd(@meritfcn,x0-pi,x0+pi,[],z,w,ntype);
% X=real(X); % truncate some roundoff
X=mod(X+pi,2*pi)-pi; % reduce to central period
ma=X*hp/pi; % scale to half period hp
return
%%%%%%
function d2=meritfcn(x,z,w,ntype)
x=exp(i*x);
if ntype=='F'
y=x-z;
else % ntype=='R'
y=log(x'*z);
end
d2=y'*diag(w)*y;
return
%%%%%%
% % test script
% %
% % NOTA BENE: meanangleoct(a,[],[],'R') will equal mean(a)
% % when all abs(a-b) < pi/2 for some value b
% %
% na=3, a=sort(mod(randn(1,na)+1,2)-1)*pi;
% da=diff([a a(1)+2*pi]); [mda,ndx]=min(da);
% a=circshift(a,[0 2-ndx]) % so that diff(a(2:3)) is smallest
% A=exp(i*a), B1=expm(a(1)*[0 -1; 1 0]),
% B2=expm(a(2)*[0 -1; 1 0]), B3=expm(a(3)*[0 -1; 1 0]),
% masimpl=[angle(mean(exp(i*a))) mean(a)]
% Bsum=B1+B2+B3; BmeanF=Bsum/sqrt(det(Bsum));
% % this expression for BmeanR should be correct for ordering of a above
% BmeanR=B1*(B1'*B2*(B2'*B3)^(1/2))^(2/3);
% mamtrx=real([[0 1]*logm(BmeanF)*[1 0]' [0 1]*logm(BmeanR)*[1 0]'])
% manorm=[meanangleoct(a,[],[],'F') meanangleoct(a,[],[],'R')]
% polar(a,1+0*a,'b*'), axis square, hold on
% polar(manorm(1),1,'rs'), polar(manorm(2),1,'gd'), hold off
% Meanangleoct Version 1.0
% Copyright (C) 2011 Alphawave Research, robjohnson@alphawaveresearch.com
% Released under GNU GPLv3 -- see file COPYING for more info.
%
% Meanangle is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or (at
% your option) any later version.
%
% Meanangle is distributed in the hope that it will be useful, but
% WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
% General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program. If not, see `http://www.gnu.org/licenses/'.
英語で：
パイソンでは：
角度の #numpy NX1 配列
if np.var(A) < np.var((A-180)%360):
average = np.average(A)
else:
average = (np.average((A-180)%360)+180)%360
複素数を使用してベクトルを使用します。私の例は、組み込みの複素数を持つ Python にあります。
import cmath # complex math
def average_angle(list_of_angles):
# make a new list of vectors
vectors= [cmath.rect(1, angle) # length 1 for each vector
for angle in list_of_angles]
vector_sum= sum(vectors)
# no need to average, we don't care for the modulus
return cmath.phase(vector_sum)
Python はベクトルの一時的な新しいリストを作成する必要がないことに注意してください。上記のすべてを 1 つのステップで実行できます。この方法を選択したのは、他の言語にも適用可能な疑似コードを近似するためです。
完全な C++ ソリューションは次のとおりです。
#include <vector>
#include <cmath>
double dAngleAvg(const vector<double>& angles) {
auto avgSin = double{ 0.0 };
auto avgCos = double{ 0.0 };
static const auto conv = double{ 0.01745329251994 }; // PI / 180
static const auto i_conv = double{ 57.2957795130823 }; // 180 / PI
for (const auto& theta : angles) {
avgSin += sin(theta*conv);
avgCos += cos(theta*conv);
}
avgSin /= (double)angles.size();
avgCos /= (double)angles.size();
auto ret = double{ 90.0 - atan2(avgCos, avgSin) * i_conv };
if (ret<0.0) ret += 360.0;
return fmod(ret, 360.0);
}
double のベクトルの形式で角度を取り、単純に double として平均を返します。角度は度単位でなければならず、もちろん平均も度単位です。
唯一の「正解」はありません。徹底的な分析のために、KV Mardia と PE Jupp の著書「Directional Statistics」(Wiley、1999 年) を読むことをお勧めします。
（推定理論または統計的推論からの私の見解を共有したいだけです）
Nimble の試行は、一連の角度の MMSE^ 推定を取得することですが、「平均化された」方向を見つけることは選択肢の 1 つです。MMAE^ 推定値、または「平均化された」方向である他の推定値を見つけることもできます。これは、方向の誤差を定量化するメトリックに依存します。またはより一般的に推定理論では、コスト関数の定義。
^ MMSE/MMAE は、最小平均二乗/絶対誤差に対応します。
ackb は、「平均角度 phi_avg には、sum_i|phi_avg-phi_i|^2 が最小になるというプロパティが必要です...それらは何かを平均しますが、角度は平均しません」
---- 平均二乗の意味でエラーを定量化します。これは最も一般的な方法の 1 つですが、唯一の方法ではありません。ここでほとんどの人が好む答え (つまり、単位ベクトルの合計と結果の角度を取得する) は、実際には合理的な解決策の 1 つです。ベクトルの方向がフォン ミーゼス分布としてモデル化されている場合、必要な「平均化された」方向として機能するのは (証明できる) ML 推定量です。この分布は凝ったものではなく、2D Guassian から定期的にサンプリングされた分布にすぎません。式を参照してください。(2.179) ビショップの著書「パターン認識と機械学習」。繰り返しますが、「平均的な」方向性を表すのにこれが唯一の最良のものではありませんが、理論的な正当性と単純な実装の両方を備えた非常に合理的なものです。
Nimble は次のように述べています。
- - 本当じゃない。「対応する単位ベクトル」は、ベクトルの方向の情報を明らかにします。角度はベクトルの長さを考慮しない量であり、単位ベクトルは長さが 1 であるという追加情報を含むものです。「単位」ベクトルを長さ 2 に定義できますが、実際には問題ではありません。
これは、移動平均を使用し、値を正規化するように注意する完全な算術ソリューションです。すべての角度が円の片側にある場合 (互いに 180° 以内) は、高速で正しい答えが得られます。
これは、値を範囲 (0, 180) にシフトするオフセットを加算し、平均を計算してからオフセットを減算することと数学的に同等です。
コメントは、特定の値がいつでも取ることができる範囲を説明します
// angles have to be in the range [0, 360) and within 180° of each other.
// n >= 1
// returns the circular average of the angles int the range [0, 360).
double meanAngle(double* angles, int n)
{
double average = angles[0];
for (int i = 1; i<n; i++)
{
// average: (0, 360)
double diff = angles[i]-average;
// diff: (-540, 540)
if (diff < -180)
diff += 360;
else if (diff >= 180)
diff -= 360;
// diff: (-180, 180)
average += diff/(i+1);
// average: (-180, 540)
if (average < 0)
average += 360;
else if (average >= 360)
average -= 360;
// average: (0, 360)
}
return average;
}
アイデアは次のとおりです。重みを維持しながら、互いに最も近い角度の平均を常に計算することにより、平均を繰り返し作成します。
別のアイデア: 与えられた角度の間の最大のギャップを見つけます。2 等分する点を見つけて、円の反対側の点を参照ゼロとして選択し、そこから平均を計算します。
これらの角度を円周上の点で表しましょう。
これらすべての点が円の同じ半分にあると仮定できますか? (そうでなければ、「平均角度」を明確に定義する方法はありません。直径上の 2 点、たとえば 0 度と 180 度を考えてみてください --- 平均は 90 度または 270 度ですか? 3 つ以上の場合はどうなりますか?ポイントを均等に分散しますか?)
この仮定の下で、その半円上の任意の点を「原点」として選び、この原点に対する与えられた角度のセットを測定します (これを「相対角度」と呼びます)。相対角度の絶対値は厳密に 180 度未満であることに注意してください。最後に、これらの相対角度の平均をとって、目的の平均角度を取得します (もちろん、原点に対して)。
Alnitakには適切なソリューションがあります。NickFortescueのソリューションは機能的に同じです。
特別な場合のために
（sum（x_component）= 0.0 && sum（y_component）= 0.0）//たとえば、10度と190度の2つの角度。
合計として0.0度を使用します
atan2（0。、0。）は未定義であり、エラーが生成されるため、計算上、このケースをテストする必要があります。
この関数は Matlab で使用できます。
function retVal=DegreeAngleMean(x)
len=length(x);
sum1=0;
sum2=0;
count1=0;
count2=0;
for i=1:len
if x(i)<180
sum1=sum1+x(i);
count1=count1+1;
else
sum2=sum2+x(i);
count2=count2+1;
end
end
if (count1>0)
k1=sum1/count1;
end
if (count2>0)
k2=sum2/count2;
end
if count1>0 && count2>0
if(k2-k1 >= 180)
retVal = ((sum1+sum2)-count2*360)/len;
else
retVal = (sum1+sum2)/len;
end
elseif count1>0
retVal = k1;
else
retVal = k2;
end
@David_Hanak からの回答の助けを借りて問題を解決しました。彼が述べているように：
同じ半円に留まりながら他の 2 つの角度を指す角度。たとえば、355 と 5 の場合、これは 180 ではなく 0 になります。これを行うには、2 つの角度の差が 180 より大きいかどうかを確認する必要があります。か否か。その場合は、上記の式を使用する前に、小さい方の角度を 360 ずつ増やしてください。
だから私がしたことは、すべての角度の平均を計算することでした. 次に、これよりも小さいすべての角度を 360 ずつ増やします。次に、それらをすべて加算して長さで割ることにより、平均を再計算します。
float angleY = 0f;
int count = eulerAngles.Count;
for (byte i = 0; i < count; i++)
angleY += eulerAngles[i].y;
float averageAngle = angleY / count;
angleY = 0f;
for (byte i = 0; i < count; i++)
{
float angle = eulerAngles[i].y;
if (angle < averageAngle)
angle += 360f;
angleY += angle;
}
angleY = angleY / count;
完璧に動作します。
任意のプログラミング言語について、次のリンクで解決策と簡単な説明を確認できます: https://rosettacode.org/wiki/Averages/Mean_angle
たとえば、C++ ソリューション:
#include<math.h>
#include<stdio.h>
double
meanAngle (double *angles, int size)
{
double y_part = 0, x_part = 0;
int i;
for (i = 0; i < size; i++)
{
x_part += cos (angles[i] * M_PI / 180);
y_part += sin (angles[i] * M_PI / 180);
}
return atan2 (y_part / size, x_part / size) * 180 / M_PI;
}
int
main ()
{
double angleSet1[] = { 350, 10 };
double angleSet2[] = { 90, 180, 270, 360};
double angleSet3[] = { 10, 20, 30};
printf ("\nMean Angle for 1st set : %lf degrees", meanAngle (angleSet1, 2));
printf ("\nMean Angle for 2nd set : %lf degrees", meanAngle (angleSet2, 4));
printf ("\nMean Angle for 3rd set : %lf degrees\n", meanAngle (angleSet3, 3));
return 0;
}
出力：
Mean Angle for 1st set : -0.000000 degrees
Mean Angle for 2nd set : -90.000000 degrees
Mean Angle for 3rd set : 20.000000 degrees
またはMatlabソリューション：
function u = mean_angle(phi)
u = angle(mean(exp(i*pi*phi/180)))*180/pi;
end
mean_angle([350, 10])
ans = -2.7452e-14
mean_angle([90, 180, 270, 360])
ans = -90
mean_angle([10, 20, 30])
ans = 20.000
starblue の回答は平均単位ベクトルの角度を示していますが、0 から 2*pi (または 0° から360°）。例えば、0°と 180°の平均は、90°または 270°のいずれかになります。
算術平均には、入力値までの距離の二乗和が最小の単一の値であるという特性があります。2 つの単位ベクトル間の単位円に沿った距離は、内積の逆余弦として簡単に計算できます。ベクトルと各入力単位ベクトルの内積の二乗逆余弦の和を最小化することによって単位ベクトルを選択すると、同等の平均が得られます。繰り返しになりますが、例外的なケースでは、2 つ以上の最小値が存在する可能性があることに注意してください。
単位球に沿った距離は、単位円に沿った距離とまったく同じ方法で計算できるため、この概念は任意の数の次元に拡張できます。これは、2 つの単位ベクトルの内積の逆余弦です。
円については、さまざまな方法でこの平均を解くことができますが、次の O(n^2) アルゴリズムを提案します (角度はラジアンで、単位ベクトルの計算は避けます)。
var bestAverage = -1
double minimumSquareDistance
for each a1 in input
var sumA = 0;
for each a2 in input
var a = (a2 - a1) mod (2*pi) + a1
sumA += a
end for
var averageHere = sumA / input.count
var sumSqDistHere = 0
for each a2 in input
var dist = (a2 - averageHere + pi) mod (2*pi) - pi // keep within range of -pi to pi
sumSqDistHere += dist * dist
end for
if (bestAverage < 0 OR sumSqDistHere < minimumSquareDistance) // for exceptional cases, sumSqDistHere may be equal to minimumSquareDistance at least once. In these cases we will only find one of the averages
minimumSquareDistance = sumSqDistHere
bestAverage = averageHere
end if
end for
return bestAverage
すべての角度が互いに 180° 以内にある場合は、より単純な O(n)+O(sort) アルゴリズムを使用できます (ここでもラジアンを使用し、単位ベクトルの使用を避けます)。
sort(input)
var largestGapEnd = input[0]
var largestGapSize = (input[0] - input[input.count-1]) mod (2*pi)
for (int i = 1; i < input.count; ++i)
var gapSize = (input[i] - input[i - 1]) mod (2*pi)
if (largestGapEnd < 0 OR gapSize > largestGapSize)
largestGapSize = gapSize
largestGapEnd = input[i]
end if
end for
double sum = 0
for each angle in input
var a2 = (angle - largestGapEnd) mod (2*pi) + largestGapEnd
sum += a2
end for
return sum / input.count
度数を使用するには、円周率を 180 に置き換えるだけです。より多くの次元を使用する場合は、反復法を使用して平均を求める必要があります。
上記の角度のいくつかに「正しい」答えを与える@Starblueとは異なる方法があります。例えば：
連続する角度間の差の合計を使用します。コード (Matlab):
function [avg] = angle_avg(angles)
last = angles(1);
sum = angles(1);
for i=2:length(angles)
diff = mod(angles(i)-angles(i-1)+ 180,360)-180
last = last + diff;
sum = sum + last;
end
avg = mod(sum/length(angles), 360);
end