Wie berechnet man den Durchschnitt eines Satzes kreisförmiger Daten?

Ich möchte den Durchschnitt eines Satzes kreisförmiger Daten berechnen. Zum Beispiel könnte ich mehrere Beispiele aus dem Lesen eines Kompasses haben. Das Problem ist natürlich, wie man mit dem Wraparound umgeht. Der gleiche Algorithmus kann für ein Zifferblatt nützlich sein.

Die eigentliche Frage ist komplizierter - was bedeuten Statistiken auf einer Kugel oder in einem algebraischen Raum, der sich "umhüllt", z. B. der additiven Gruppe mod n. Die Antwort ist möglicherweise nicht eindeutig, z. B. könnte der Durchschnitt von 359 Grad und 1 Grad 0 Grad oder 180 Grad sein, aber statistisch gesehen sieht 0 besser aus.

Dies ist ein echtes Programmierproblem für mich und ich versuche, es nicht nur wie ein mathematisches Problem aussehen zu lassen.

Berechnen Sie Einheitsvektoren aus den Winkeln und nehmen Sie den Winkel ihres Durchschnitts.

Diese Frage wird im Buch "Statistics On Spheres", Geoffrey S. Watson, Vorlesungsnotizen der Universität von Arkansas in den mathematischen Wissenschaften, 1983, John Wiley & Sons, Inc., ausführlich untersucht, wie unter http: //catless.ncl erwähnt. ac.uk/Risks/7.44.html#subj4 von Bruce Karsh.

Ein guter Weg, um einen durchschnittlichen Winkel A aus einer Reihe von Winkelmessungen a [i] 0 <= i zu schätzen

                   sum_i_from_1_to_N sin(a[i])
a = arctangent ---------------------------
                   sum_i_from_1_to_N cos(a[i])

Die von starblue angegebene Methode ist rechnerisch äquivalent, aber seine Gründe sind klarer und wahrscheinlich programmatisch effizienter und funktionieren auch im Nullfall gut, also ein großes Lob an ihn.

Das Thema wird nun in Wikipedia und mit anderen Verwendungszwecken, wie z. B. Bruchteilen, ausführlicher untersucht .

Ich sehe das Problem - wenn Sie zum Beispiel einen Winkel von 45 'und einen Winkel von 315' haben, wäre der "natürliche" Durchschnitt 180 ', aber der gewünschte Wert ist tatsächlich 0'.

Ich denke, Starblue ist auf etwas. Berechnen Sie einfach die (x, y) kartesischen Koordinaten für jeden Winkel und addieren Sie die resultierenden Vektoren. Der Winkelversatz des endgültigen Vektors sollte Ihr gewünschtes Ergebnis sein.

x = y = 0
foreach angle {
    x += cos(angle)
    y += sin(angle)
}
average_angle = atan2(y, x)

Ich ignoriere vorerst, dass eine Kompassrichtung im Norden beginnt und im Uhrzeigersinn verläuft, während "normale" kartesische Koordinaten entlang der X-Achse mit Null beginnen und dann gegen den Uhrzeigersinn verlaufen. Die Mathematik sollte unabhängig davon auf die gleiche Weise funktionieren.

FÜR DEN BESONDEREN FALL VON ZWEI WINKELN:

Die Antwort ((a + b) mod 360) / 2 ist FALSCH . Für die Winkel 350 und 2 ist der nächste Punkt 356, nicht 176.

Die Einheitsvektor- und Triggerlösungen können zu teuer sein.

Was ich von einem kleinen Basteln habe, ist:

diff = ( ( a - b + 180 + 360 ) mod 360 ) - 180
angle = (360 + b + ( diff / 2 ) ) mod 360

• 0, 180 -> 90 (zwei Antworten dafür: Diese Gleichung übernimmt die Antwort im Uhrzeigersinn von a)
• 180, 0 -> 270 (siehe oben)
• 180, 1 -> 90,5
• 1, 180 -> 90,5
• 20, 350 -> 5
• 350, 20 -> 5 (alle folgenden Beispiele kehren sich auch richtig um)
• 10, 20 -> 15
• 350, 2 -> 356
• 359, 0 -> 359,5
• 180, 180 -> 180

ackb hat recht, dass diese vektorbasierten Lösungen nicht als echte Durchschnittswerte von Winkeln betrachtet werden können, sondern nur ein Durchschnitt der Einheitsvektor-Gegenstücke sind. Die von ackb vorgeschlagene Lösung scheint jedoch mathematisch nicht zu klingen.

Das Folgende ist eine Lösung, die mathematisch aus dem Ziel der Minimierung (Winkel [i] - avgAngle) ^ 2 (wobei die Differenz gegebenenfalls korrigiert wird) abgeleitet wird, was sie zu einem echten arithmetischen Mittel der Winkel macht.

Zunächst müssen wir uns genau ansehen, in welchen Fällen sich der Unterschied zwischen den Winkeln von dem Unterschied zwischen ihren Gegenstücken mit normaler Anzahl unterscheidet. Betrachten Sie die Winkel x und y. Wenn y> = x - 180 und y <= x + 180, können wir die Differenz (xy) direkt verwenden. Andernfalls müssen wir (y + 360) anstelle von y in der Berechnung verwenden, wenn die erste Bedingung nicht erfüllt ist. Entsprechend müssen wir (y-360) anstelle von y verwenden, wenn die zweite Bedingung nicht erfüllt ist. Da sich die Gleichung der Kurve nur an den Punkten minimiert, an denen sich diese Ungleichungen von wahr zu falsch oder umgekehrt ändern, können wir den gesamten [0,360] -Bereich in eine Reihe von Segmenten unterteilen, die durch diese Punkte getrennt sind. Dann müssen wir nur das Minimum jedes dieser Segmente und dann das Minimum des Minimums jedes Segments finden, das der Durchschnitt ist.

Hier ist ein Bild, das zeigt, wo die Probleme bei der Berechnung von Winkeldifferenzen auftreten. Wenn x im grauen Bereich liegt, liegt ein Problem vor.

Um eine Variable zu minimieren, können wir abhängig von der Kurve die Ableitung dessen nehmen, was wir minimieren möchten, und dann den Wendepunkt finden (an dem die Ableitung = 0 ist).

Hier wenden wir die Idee an, die quadratische Differenz zu minimieren, um die gemeinsame arithmetische Mittelwertformel abzuleiten: sum (a [i]) / n. Die Kurve y = Summe ((a [i] -x) ^ 2) kann auf folgende Weise minimiert werden:

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

Wenden Sie es nun auf Kurven mit unseren angepassten Differenzen an:

b = Teilmenge von a, wobei die richtige (Winkel-) Differenz a [i] -xc = Teilmenge von a, wo die richtige (Winkel-) Differenz (a [i] -360) -x cn = Größe von cd = Teilmenge von a, wo die korrekte (Winkel-) Differenz (a [i] +360) -x dn = Größe von 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

Dies allein reicht nicht aus, um das Minimum zu erreichen, während es für normale Werte mit einer unbegrenzten Menge funktioniert. Das Ergebnis liegt also definitiv im Bereich der Menge und ist daher gültig. Wir brauchen das Minimum innerhalb eines Bereichs (definiert durch das Segment). Wenn das Minimum kleiner als die Untergrenze unseres Segments ist, muss das Minimum dieses Segments an der Untergrenze liegen (da quadratische Kurven nur 1 Wendepunkt haben), und wenn das Minimum größer als die Obergrenze unseres Segments ist, liegt das Minimum des Segments bei obere Grenze. Nachdem wir das Minimum für jedes Segment haben, finden wir einfach dasjenige, das den niedrigsten Wert für das hat, was wir minimieren (Summe ((b [i] -x) ^ 2) + Summe (((c [i] -360)))) ) -b) ^ 2) + Summe (((d [i] +360) -c) ^ 2)).

Hier ist ein Bild der Kurve, das zeigt, wie sie sich an den Punkten ändert, an denen x = (a [i] +180)% 360 ist. Der fragliche Datensatz ist {65,92,230,320,250}.

Hier ist eine Implementierung des Algorithmus in Java, einschließlich einiger Optimierungen, seine Komplexität ist O (nlogn). Sie kann auf O (n) reduziert werden, wenn Sie die vergleichsbasierte Sortierung durch eine nicht vergleichsbasierte Sortierung wie die Radix-Sortierung ersetzen.

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_;
}

Das arithmetische Mittel einer Reihe von Winkeln stimmt möglicherweise nicht mit Ihrer intuitiven Vorstellung überein, wie hoch der Durchschnitt sein sollte. Zum Beispiel ist das arithmetische Mittel der Menge {179,179,0,181,181} 216 (und 144). Die Antwort, an die Sie sofort denken, ist wahrscheinlich 180, es ist jedoch bekannt, dass das arithmetische Mittel stark von Kantenwerten beeinflusst wird. Sie sollten sich auch daran erinnern, dass Winkel keine Vektoren sind, so ansprechend dies auch erscheinen mag, wenn Sie manchmal mit Winkeln arbeiten.

Dieser Algorithmus gilt natürlich auch für alle Größen, die der modularen Arithmetik (mit minimaler Anpassung) entsprechen, wie z. B. die Tageszeit.

Ich möchte auch betonen, dass, obwohl dies im Gegensatz zu den Vektorlösungen ein echter Durchschnitt der Winkel ist, dies nicht unbedingt bedeutet, dass es die Lösung ist, die Sie verwenden sollten, der Durchschnitt der entsprechenden Einheitsvektoren möglicherweise der Wert ist, den Sie tatsächlich haben sollte verwendet werden.

Sie müssen den Durchschnitt genauer definieren . Für den speziellen Fall von zwei Winkeln kann ich mir zwei verschiedene Szenarien vorstellen:

1. Der "wahre" Durchschnitt, dh (a + b) / 2% 360.
2. Der Winkel, der "zwischen" den beiden anderen zeigt, während er im selben Halbkreis bleibt, z. B. für 355 und 5, ist 0 und nicht 180. Dazu müssen Sie prüfen, ob der Unterschied zwischen den beiden Winkeln größer als 180 ist oder nicht. Wenn ja, erhöhen Sie den kleineren Winkel um 360, bevor Sie die obige Formel verwenden.

Ich sehe jedoch nicht, wie die zweite Alternative für den Fall von mehr als zwei Winkeln verallgemeinert werden kann.

Wie bei allen Durchschnittswerten hängt die Antwort von der Wahl der Metrik ab. Für eine gegebene Metrik M ist der Durchschnitt einiger Winkel a_k in [-pi, pi] für k in [1, N] der Winkel a_M, der die Summe der quadratischen Abstände d ^ 2_M (a_M, a_k) minimiert. Für einen gewichteten Mittelwert schließt man einfach die Gewichte w_k in die Summe ein (so dass sum_k w_k = 1). Das ist,

a_M = arg min_x sum_k w_k d ^ 2_M (x, a_k)

Zwei gängige Metrikoptionen sind die Frobenius- und die Riemann-Metrik. Für die Frobenius-Metrik existiert eine direkte Formel, die dem in der Kreisstatistik üblichen Begriff der durchschnittlichen Peilung entspricht. Weitere Informationen finden Sie unter "Mittelwerte und Mittelung in der Rotationsgruppe", Maher Moakher, SIAM Journal on Matrix Analysis and Applications, Band 24, Ausgabe 1, 2002.
http://link.aip.org/link/?SJMAEL/24/1/1

Hier ist eine Funktion für GNU Octave 3.2.4, die die Berechnung durchführt:

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, [email protected]
%     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/'.

Ich möchte eine Methode, die ich verwendet habe, mit einem Mikrocontroller teilen, der keine Gleitkomma- oder Trigonometriefunktionen hatte. Ich musste immer noch 10 Rohlagerwerte "mitteln", um Abweichungen auszugleichen.

1. Überprüfen Sie, ob das erste Lager im Bereich von 270 bis 360 oder 0 bis 90 Grad liegt (nördliche zwei Quadranten).
2. Wenn dies der Fall ist, drehen Sie diese und alle nachfolgenden Messwerte um 180 Grad und halten Sie alle Werte im Bereich 0 <= Peilung <360. Andernfalls nehmen Sie die Messwerte so, wie sie kommen.
3. Berechnen Sie nach 10 Messungen den numerischen Durchschnitt unter der Annahme, dass kein Umlauf stattgefunden hat
4. Wenn die 180-Grad-Drehung wirksam war, drehen Sie den berechneten Durchschnitt um 180 Grad, um zu einem "wahren" Lager zurückzukehren.

Es ist nicht ideal; es kann brechen. Ich bin in diesem Fall damit durchgekommen, weil sich das Gerät nur sehr langsam dreht. Ich werde es dort veröffentlichen, falls jemand anderes unter ähnlichen Einschränkungen arbeitet.

Auf Englisch:

1. Erstellen Sie einen zweiten Datensatz mit allen um 180 verschobenen Winkeln.
2. Nehmen Sie die Varianz beider Datensätze.
3. Nehmen Sie den Durchschnitt des Datensatzes mit der kleinsten Varianz.
4. Wenn dieser Durchschnitt aus der verschobenen Menge stammt, verschieben Sie die Antwort erneut um 180.

In Python:

Ein #numpy NX1 Array von Winkeln

if np.var(A) < np.var((A-180)%360):
    average = np.average(A)

else:
    average = (np.average((A-180)%360)+180)%360

Hier ist die vollständige Lösung: (Die Eingabe ist eine Anordnung von Peilungen in Grad (0-360).

public static int getAvarageBearing(int[] arr)
{
    double sunSin = 0;
    double sunCos = 0;
    int counter = 0;

    for (double bearing : arr)
    {
        bearing *= Math.PI/180;

        sunSin += Math.sin(bearing);
        sunCos += Math.cos(bearing);
        counter++; 
    }

    int avBearing = INVALID_ANGLE_VALUE;
    if (counter > 0)
    {
        double bearingInRad = Math.atan2(sunSin/counter, sunCos/counter);
        avBearing = (int) (bearingInRad*180f/Math.PI);
        if (avBearing<0)
            avBearing += 360;
    }

    return avBearing;
}

In Python mit Winkeln zwischen [-180, 180]

def add_angles(a, b):
  return (a + b + 180) % 360 - 180

def average_angles(a, b):
  return add_angles(a, add_angles(-a, b)/2)

Einzelheiten:

Für den Durchschnitt von zwei Winkeln gibt es zwei Durchschnittswerte, die 180 ° voneinander entfernt sind, aber wir möchten vielleicht den näheren Durchschnitt.

Visuell ergibt der Durchschnitt von Blau ( b ) und Grün ( a ) den blaugrünen Punkt:

Winkel 'Wrap Around' (z. B. 355 + 10 = 5), aber Standardarithmetik ignoriert diesen Verzweigungspunkt. Wenn jedoch der Winkel b dem Verzweigungspunkt entgegengesetzt ist, ergibt ( b + g ) / 2 den nächsten Durchschnitt: den blaugrünen Punkt.

Für zwei beliebige Winkel können wir das Problem so drehen, dass einer der Winkel dem Verzweigungspunkt entgegengesetzt ist, eine Standardmittelung durchführen und dann zurückdrehen.

Ich würde den Vektorweg mit komplexen Zahlen gehen. Mein Beispiel ist Python, in das komplexe Zahlen integriert sind:

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)

Beachten Sie, dass Python nicht brauchen eine temporäre neue Liste von Vektoren, die alle der oben in einem Schritt durchgeführt werden kann , bauen; Ich habe mich nur für diesen Weg entschieden, um den Pseudocode zu approximieren, der auch für andere Sprachen gilt.

Hier ist eine vollständige C ++ - Lösung:

#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);
}

Es nimmt die Winkel in Form eines Doppelvektors an und gibt den Durchschnitt einfach als Doppel zurück. Die Winkel müssen in Grad angegeben werden, und natürlich ist der Durchschnitt auch in Grad.

Basierend auf Alnitaks Antwort habe ich eine Java-Methode zur Berechnung des Durchschnitts mehrerer Winkel geschrieben:

Wenn Ihre Winkel im Bogenmaß sind:

public static double averageAngleRadians(double... angles) {
    double x = 0;
    double y = 0;
    for (double a : angles) {
        x += Math.cos(a);
        y += Math.sin(a);
    }

    return Math.atan2(y, x);
}

Wenn Ihre Winkel in Grad sind:

public static double averageAngleDegrees(double... angles) {
    double x = 0;
    double y = 0;
    for (double a : angles) {
        x += Math.cos(Math.toRadians(a));
        y += Math.sin(Math.toRadians(a));
    }

    return Math.toDegrees(Math.atan2(y, x));
}

Hier ist eine Idee: Erstellen Sie den Durchschnitt iterativ, indem Sie immer den Durchschnitt der Winkel berechnen, die am nächsten beieinander liegen, und dabei ein Gewicht beibehalten.

Eine andere Idee: Finden Sie die größte Lücke zwischen den angegebenen Winkeln. Suchen Sie den Punkt, der ihn halbiert, und wählen Sie dann den entgegengesetzten Punkt auf dem Kreis als Referenznullpunkt aus, um den Durchschnitt daraus zu berechnen.

Stellen wir diese Winkel mit Punkten auf dem Kreisumfang dar.

Können wir annehmen, dass alle diese Punkte auf dieselbe Hälfte des Kreises fallen? (Andernfalls gibt es keine offensichtliche Möglichkeit, den "Durchschnittswinkel" zu definieren. Stellen Sie sich zwei Punkte auf dem Durchmesser vor, z. B. 0 ° und 180 ° - ist der Durchschnitt 90 ° oder 270 °? Was passiert, wenn wir 3 oder mehr haben? Punkte gleichmäßig verteilen?)

Mit dieser Annahme wählen wir einen beliebigen Punkt auf diesem Halbkreis als "Ursprung" und messen den gegebenen Satz von Winkeln in Bezug auf diesen Ursprung (nennen wir dies den "relativen Winkel"). Beachten Sie, dass der relative Winkel einen absoluten Wert von weniger als 180 Grad hat. Nehmen Sie schließlich den Mittelwert dieser relativen Winkel, um den gewünschten Durchschnittswinkel zu erhalten (natürlich relativ zu unserem Ursprung).

Es gibt keine einzige "richtige Antwort". Ich empfehle, das Buch KV Mardia und PE Jupp, "Directional Statistics" (Wiley, 1999), für eine gründliche Analyse zu lesen.

(Ich möchte nur meinen Standpunkt aus der Schätzungstheorie oder der statistischen Inferenz teilen.)

Nimbles Versuch besteht darin, die MMSE ^ -Schätzung einer Reihe von Winkeln zu erhalten, aber es ist eine der Möglichkeiten, eine "gemittelte" Richtung zu finden. Man kann auch eine MMAE ^ -Schätzung oder eine andere Schätzung als "gemittelte" Richtung finden, und dies hängt von Ihrem metrischen quantifizierenden Richtungsfehler ab. oder allgemeiner in der Schätzungstheorie die Definition der Kostenfunktion.

^ MMSE / MMAE entspricht dem minimalen mittleren quadratischen / absoluten Fehler.

ackb sagte: "Der durchschnittliche Winkel phi_avg sollte die Eigenschaft haben, dass sum_i | phi_avg-phi_i | ^ 2 minimal wird ... sie mitteln etwas, aber keine Winkel."

---- Sie quantifizieren Fehler im quadratischen Sinne und es ist einer der häufigsten Wege, jedoch nicht der einzige. Die Antwort, die von den meisten Menschen hier bevorzugt wird (dh Summe der Einheitsvektoren und Ermittlung des Ergebniswinkels), ist tatsächlich eine der vernünftigen Lösungen. Es ist (kann bewiesen werden) der ML-Schätzer, der als die "gemittelte" Richtung dient, die wir wollen, wenn die Richtungen der Vektoren als von-Mises-Verteilung modelliert werden. Diese Verteilung ist nicht ausgefallen und nur eine periodisch abgetastete Verteilung von einem 2D-Guassianer. Siehe Gl. (2.179) in Bishops Buch "Mustererkennung und maschinelles Lernen". Auch hier ist es keineswegs die einzig beste, um die "durchschnittliche" Richtung darzustellen, aber es ist durchaus vernünftig, sowohl eine gute theoretische Begründung als auch eine einfache Implementierung zu haben.

Nimble sagte: "ackb hat Recht, dass diese vektorbasierten Lösungen nicht als echte Durchschnittswerte von Winkeln betrachtet werden können, sondern nur ein Durchschnitt der Gegenstücke zu Einheitsvektoren sind."

----das ist nicht wahr. Die "Einheitsvektor-Gegenstücke" enthüllen die Information über die Richtung eines Vektors. Der Winkel ist eine Größe ohne Berücksichtigung der Länge des Vektors, und der Einheitsvektor ist etwas mit zusätzlichen Informationen, dass die Länge 1 ist. Sie können Ihren "Einheits" -Vektor so definieren, dass er die Länge 2 hat, es spielt keine Rolle.

Hier ist eine vollständig arithmetische Lösung, bei der gleitende Durchschnitte verwendet werden und die Werte normalisiert werden. Es ist schnell und liefert korrekte Antworten, wenn sich alle Winkel auf einer Seite des Kreises befinden (innerhalb von 180 ° zueinander).

Es ist mathematisch äquivalent zum Addieren des Versatzes, der die Werte in den Bereich (0, 180) verschiebt, den Mittelwert berechnet und dann den Versatz subtrahiert.

Die Kommentare beschreiben, welchen Bereich ein bestimmter Wert zu einem bestimmten Zeitpunkt annehmen kann

// 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;
}

Nun, ich bin sehr spät zur Party, dachte aber, ich würde meine 2 Cent hinzufügen, da ich keine endgültige Antwort finden konnte. Am Ende habe ich die folgende Java-Version der Mitsuta-Methode implementiert, die hoffentlich eine einfache und robuste Lösung bietet. Zumal die Standardabweichung sowohl eine Maßstreuung liefert als auch, wenn sd == 90, anzeigt, dass die Eingabewinkel zu einem mehrdeutigen Mittelwert führen.

EDIT: Eigentlich habe ich festgestellt, dass meine ursprüngliche Implementierung noch weiter vereinfacht werden kann, und zwar besorgniserregend einfach, wenn man die Konversation und Trigonometrie berücksichtigt, die in den anderen Antworten vor sich gehen.

/**
 * The Mitsuta method
 *
 * @param angles Angles from 0 - 360
 * @return double array containing
 * 0 - mean
 * 1 - sd: a measure of angular dispersion, in the range [0..360], similar to standard deviation.
 * Note if sd == 90 then the mean can also be its inverse, i.e. 360 == 0, 300 == 60.
 */
public static double[] getAngleStatsMitsuta(double... angles) {
    double sum = 0;
    double sumsq = 0;
    for (double angle : angles) {
        if (angle >= 180) {
            angle -= 360;
        }
        sum += angle;
        sumsq += angle * angle;
    }

    double mean = sum / angles.length;
    return new double[]{mean <= 0 ? 360 + mean: mean, Math.sqrt(sumsq / angles.length - (mean * mean))};
}

... und für alle (Java) Geeks da draußen können Sie den obigen Ansatz verwenden, um den mittleren Winkel in einer Linie zu erhalten.

Arrays.stream(angles).map(angle -> angle<180 ? angle: (angle-360)).sum() / angles.length;

Alnitak hat die richtige Lösung. Die Lösung von Nick Fortescue ist funktional dieselbe.

Für den Sonderfall wo

(Summe (x_Komponente) = 0,0 && Summe (y_Komponente) = 0,0) // zB 2 Winkel von 10 und 190 Grad pro Stück.

Verwenden Sie 0,0 Grad als Summe

Rechnerisch müssen Sie diesen Fall testen, da atan2 (0., 0.) undefiniert ist und einen Fehler erzeugt.

Der durchschnittliche Winkel phi_avg sollte die Eigenschaft haben, dass sum_i | phi_avg-phi_i | ^ 2 minimal wird, wobei der Unterschied in [-Pi, Pi) liegen muss (weil er möglicherweise kürzer ist, um umgekehrt zu sein!). Dies wird leicht erreicht, indem alle Eingabewerte auf [0, 2Pi] normalisiert werden, ein laufender Durchschnitt phi_run beibehalten wird und die Normalisierung | phi_i-phi_run | gewählt wird zu [-Pi, Pi) (durch Addition oder Subtraktion von 2Pi). Die meisten der oben genannten Vorschläge machen etwas anderes, was nicht der Fall ist diese minimale Eigenschaft hat, dh sie mitteln etwas , aber keine Winkel.

Ich habe das Problem mit Hilfe der Antwort von @David_Hanak gelöst. Wie er sagt:

Der Winkel, der "zwischen" den beiden anderen zeigt, während er im selben Halbkreis bleibt, z. B. für 355 und 5, ist 0 und nicht 180. Dazu müssen Sie prüfen, ob der Unterschied zwischen den beiden Winkeln größer als 180 ist oder nicht. Wenn ja, erhöhen Sie den kleineren Winkel um 360, bevor Sie die obige Formel verwenden.

Also habe ich den Durchschnitt aller Winkel berechnet. Und dann erhöhen Sie alle Winkel, die kleiner als dieser sind, um 360. Berechnen Sie dann den Durchschnitt neu, indem Sie sie alle addieren und durch ihre Länge dividieren.

        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;

Funktioniert perfekt.

Python-Funktion:

from math import sin,cos,atan2,pi
import numpy as np
def meanangle(angles,weights=0,setting='degrees'):
    '''computes the mean angle'''
    if weights==0:
         weights=np.ones(len(angles))
    sumsin=0
    sumcos=0
    if setting=='degrees':
        angles=np.array(angles)*pi/180
    for i in range(len(angles)):
        sumsin+=weights[i]/sum(weights)*sin(angles[i])
        sumcos+=weights[i]/sum(weights)*cos(angles[i])
    average=atan2(sumsin,sumcos)
    if setting=='degrees':
        average=average*180/pi
    return average

Sie können diese Funktion in Matlab verwenden:

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 

Eine Lösung und eine kleine Erklärung finden Sie unter folgendem Link für JEDE Programmiersprache: https://rosettacode.org/wiki/Averages/Mean_angle

Zum Beispiel C ++ - Lösung :

#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;
}

Ausgabe:

Mean Angle for 1st set : -0.000000 degrees
Mean Angle for 2nd set : -90.000000 degrees
Mean Angle for 3rd set : 20.000000 degrees

Oder Matlab-Lösung :

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

Während die Antwort von starblue den Winkel des durchschnittlichen Einheitsvektors angibt, ist es möglich, das Konzept des arithmetischen Mittels auf Winkel zu erweitern, wenn Sie akzeptieren, dass es möglicherweise mehr als eine Antwort im Bereich von 0 bis 2 * pi (oder 0 ° bis) gibt 360 °). Beispielsweise kann der Durchschnitt von 0 ° und 180 ° entweder 90 ° oder 270 ° betragen.

Das arithmetische Mittel hat die Eigenschaft, der Einzelwert mit der minimalen Summe der quadratischen Abstände zu den Eingabewerten zu sein. Der Abstand entlang des Einheitskreises zwischen zwei Einheitsvektoren kann leicht als inverser Kosinus ihres Punktprodukts berechnet werden. Wenn wir einen Einheitsvektor wählen, indem wir die Summe des quadratischen inversen Kosinus des Punktprodukts unseres Vektors und jedes Eingabeeinheitsvektors minimieren, haben wir einen äquivalenten Durchschnitt. Denken Sie auch hier daran, dass es in Ausnahmefällen zwei oder mehr Mindestwerte geben kann.

Dieses Konzept könnte auf eine beliebige Anzahl von Dimensionen erweitert werden, da der Abstand entlang der Einheitskugel genauso berechnet werden kann wie der Abstand entlang des Einheitskreises - der inverse Kosinus des Punktprodukts zweier Einheitsvektoren.

Für Kreise könnten wir diesen Durchschnitt auf verschiedene Arten lösen, aber ich schlage den folgenden O (n ^ 2) -Algorithmus vor (Winkel sind im Bogenmaß und ich vermeide die Berechnung der Einheitsvektoren):

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

Wenn alle Winkel innerhalb von 180 ° zueinander liegen, könnten wir einen einfacheren O (n) + O (Sortier) -Algorithmus verwenden (wiederum unter Verwendung des Bogenmaßes und unter Vermeidung der Verwendung von Einheitsvektoren):

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

Um Grad zu verwenden, ersetzen Sie pi einfach durch 180. Wenn Sie mehr Dimensionen verwenden möchten, müssen Sie höchstwahrscheinlich eine iterative Methode verwenden, um den Durchschnitt zu ermitteln.

Das Problem ist sehr einfach. 1. Stellen Sie sicher, dass alle Winkel zwischen -180 und 180 Grad liegen. 2. a Addiere alle nicht negativen Winkel, nimm ihren Durchschnitt und zähle, wie viele. 2. b. Füge alle negativen Winkel hinzu, nimm ihren Durchschnitt und zähle, wie viele. 3. Nehmen Sie die Differenz von pos_average minus neg_average. Wenn die Differenz größer als 180 ist, ändern Sie die Differenz in 360 minus Differenz. Ansonsten ändern Sie einfach das Vorzeichen der Differenz. Beachten Sie, dass der Unterschied immer nicht negativ ist. Der Average_Angle entspricht dem pos_average plus Differenz multipliziert mit dem "Gewicht", der negativen Anzahl geteilt durch die Summe aus negativer und positiver Anzahl

Hier ist ein Java-Code für durchschnittliche Winkel, ich denke, es ist ziemlich robust.

public static double getAverageAngle(List<Double> angles)
{
    // r = right (0 to 180 degrees)

    // l = left (180 to 360 degrees)

    double rTotal = 0;
    double lTotal = 0;
    double rCtr = 0;
    double lCtr = 0;

    for (Double angle : angles)
    {
        double norm = normalize(angle);
        if (norm >= 180)
        {
            lTotal += norm;
            lCtr++;
        } else
        {
            rTotal += norm;
            rCtr++;
        }
    }

    double rAvg = rTotal / Math.max(rCtr, 1.0);
    double lAvg = lTotal / Math.max(lCtr, 1.0);

    if (rAvg > lAvg + 180)
    {
        lAvg += 360;
    }
    if (lAvg > rAvg + 180)
    {
        rAvg += 360;
    }

    double rPortion = rAvg * (rCtr / (rCtr + lCtr));
    double lPortion = lAvg * (lCtr / (lCtr + rCtr));
    return normalize(rPortion + lPortion);
}

public static double normalize(double angle)
{
    double result = angle;
    if (angle >= 360)
    {
        result = angle % 360;
    }
    if (angle < 0)
    {
        result = 360 + (angle % 360);
    }
    return result;
}

Ich habe eine andere Methode als @Starblue, die "richtige" Antworten auf einige der oben angegebenen Winkel gibt. Beispielsweise:

• angle_avg ([350,10]) = 0
• angle_avg ([- 90,90,40]) = 13,333
• angle_avg ([350,2]) = 356

Es wird eine Summe über die Unterschiede zwischen aufeinanderfolgenden Winkeln verwendet. Der Code (in 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