Fehler in Berechnungen

Christian Gast #1722050
Lesenswert?
• ▲
▼ 
Hallo,

ich habe eine Funktion zum Umwandeln von GPS-Koordinaten ins UTM-Format 
geschrieben. Als ich sie zuerst auf dem PC in VB ausprobiert hatte, 
kamen exakte Ergebnisse raus.
Jetzt hab ich sie aber in C für einen ATmega32 geschrieben und der 
Rechtswert lag um ca. 1,5m daneben, was nicht so schlimm wäre. Beim 
Hochwert stimmten aber nur die ersten drei Stellen, was einem Fehler von 
über 1000m entspricht. Ich habe sowohl in VB als auch in C 
Gleitkommazahlen einfacher Genauigkeit verwendet.
Beim Hochwert ist schon die Meridianbogenlänge so ungenau.
Woran kann dieser große Fehler liegen? sind die float-Variablen zu 
ungenau?

Hier der Code:

1
  const float A  = 6378137.0;    // grosse WGS84-Halbachse

2
  const float B  = 6356752.3;    // kleine WGS84-Halbachse

3
  const float K0 = 0.9996;    // UTM Streckfaktor

4
  const float L0 = 0.1570796;    // UTM-Zentralmeridian fuer Zone 32 (Bogenmass)

5
  const float e2 = 0.0067395;    // (pow(A, 2) - pow(B, 2)) / pow(A, 2); 

6
  const float n  = 0.0016792;    // (A - B) / (A + B);

7


8
  float Lon = 7 * (M_PI / 180);   // Geographische Laenge als Bogenmass

9
  float Lat = 49 * (M_PI / 180);  // Geographische Breite als Bogenmass

10
  float G;            // Meridianbogenlaenge

11
  float Hw;            // Hochwert

12
  float Rw;            // Rechtswert

13


14
  // Parameter zur Berechnung der UTM-Koordinaten

15
  float nu2 = (e2 / (1 - e2)) * pow(cos(Lat), 2);

16
  float t   = tan(Lat);

17
  float N   = A / (sqrt(1 - e2 * pow(sin(Lat), 2)));

18
  float l   = Lon - L0;

19
  

20
  // Meridianbogenlaenge berechnen

21
  G = (A / (1 + n)) * ((1 + pow(n, 2) / 4) * Lat - (3 / 2) * (n - pow(n, 3) / 8) * sin(2 * Lat) + (15 / 16) * pow(n, 2) * sin(4 * Lat) - (35 / 48) * pow(n, 3) * sin(6 * Lat));

22
  

23
  // Hochwert berechnen

24
  Hw = K0 * (G + ((N * t) / 2) * pow(l, 2) * pow(cos(Lat), 2) + ((N * t) / 24) * (5 - pow(t, 2) + 9 * nu2) * pow(l, 4) * pow(cos(Lat), 4) + ((N * t) / 720) * (61 - 58 * pow(t, 2) + pow(t, 4) + 270 * nu2 - 330 * pow(t, 2) * nu2) * pow(l, 6) * pow(cos(Lat), 6));

25
  

26
  // Rechtswert berechnen

27
  Rw = 500000 + K0 * (N * l * cos(Lat) + (N / 6) * (1 - pow(t, 2) + nu2) * pow(l, 3) * pow(cos(Lat), 3) + (N / 120) * (5 - 18 * pow(t, 2) + pow(t, 4) + 14 * nu2 - 58 * pow(t, 2) * nu2) * pow(l, 5) * pow(cos(Lat), 5) + (N / 5040) * (61 - 479 * pow(t, 2) + 179 * pow(t, 4)) * pow(l, 7) * pow(cos(Lat), 7));




Vielen Dank für eure Hilfe!

Grüße, CHristian
Antwort Zitieren Bearbeiten Thread verschieben Thread sperren Anmeldepflicht aktivieren Anpinnen Thread löschen Thread mit anderem zusammenführen Beitrag melden
Rolf Magnus Gast #1722149
Lesenswert?
• ▲
▼ 
> Woran kann dieser große Fehler liegen? sind die float-Variablen zu
> ungenau?

Höchstwahrscheinlich. Kann es sein, daß in deinem VB-Code die 
Berechnungen mit doppelter Genauigkeit durchgeführt werden und damit 
auch die Zwischenwerte genauer sind? Dann wird erst beim Schreiben in 
die Variable die Genauigkeit reduziert. Beim AVR dagegen sind die 
Berechnungen immer mit einfacher Genauigkeit. Sind die Ergebnisse des 
C-Code genauer, wenn du ihn auf dem PC ausführst?
Antwort Zitieren Bearbeiten Löschen Beitrag melden
Karl H. #1722287
Lesenswert?
• ▲
▼ 
Christian schrieb:

> Woran kann dieser große Fehler liegen? sind die float-Variablen zu
> ungenau?

Als erstes ersetzt du mal alle pow Aufrufe
Das ist doch Unsinn pow(n,2) über Logarithmen rechnen zu lassen, wenn es 
ein simples n*n auch tut.

>   G = (A / (1 + n)) * ((1 + pow(n, 2) / 4) * Lat - (3 / 2)

Aha.
Dir ist bewusst, dass 3 / 2 ein Ergebnis von 1 liefert und keineswegs 
1.5?

http://www.mikrocontroller.net/articles/FAQ#Datentypen_in_Operationen

Ditto für noch viele andere Divisionen mit konstanten Zahlen, die in 
deinem Code vorkommen.

Bei deinen Berechnungen wirst du allerdings nie ein einigermassen 
genaues Ergebnis erhalten. 6 signifikante Stellen sind nun mal nicht 
viel. Aber ersetz erst mal die offensichtlichen Fehler, dann sollte es 
schon genauer werden.

Und wenn du dann auch noch gleiche Berechnungen vorziehst und in 
Hilfsvariablen ablegst, dann kannst du sogar die Laufzeit noch ein wenig 
drücken. sin(Lat) und cos(Lat) sind zb solche Kandidaten.
Antwort Zitieren Bearbeiten Löschen Beitrag melden
Christian Gast #1722401
Lesenswert?
• ▲
▼ 
Hallo,

jetzt hab ich mein genaues Ergebnis... ich hatte wirklich übersehen, 
dass die ganzzahligen Divisionen ja gar kein float ergeben.

Ab welcher Potenz lohnt es sich eigentlich, pow anstatt von 
Multiplizieren zu verwenden?

Vielen Dank und Grüße,
Christian
Antwort Zitieren Bearbeiten Löschen Beitrag melden
Karl H. #1722419
Lesenswert?
• ▲
▼ 
Christian schrieb:

> Ab welcher Potenz lohnt es sich eigentlich, pow anstatt von
> Multiplizieren zu verwenden?

Musst du ausprobieren.
Du kannst allerdings davon ausgehen, dass pow diese Berechnung hier 
durchführt

      exp( exponent * ln( mantisse ) )

eine andere Chance hat pow gar nicht um nicht-ganzzahlige Exponenten zu 
berechnen.

weder ln noch exp sind aber 'billige' Funktionen. In ihnen stecken 
selbst wieder eine ganze Menge Multiplikationen.

Benutze pow wenn du musst!
Praktisch alle Potenzen die in technischen Formeln vorkommen, können 
einfacher berechnet werden, indem man sich selbst eine Potenz dafür 
macht. In technischen Formenl kommen häufig Quadrate, 3te Potenzen und 
ganz selten 4te Potenzen vor. Bei all diesen bist du mit Ausschreiben 
(oder Makro oder Inline-Funktion) besser bedient (sowohl was Genauigkeit 
als auch Geschwindigkeit angeht) als mit pow.
Antwort Zitieren Bearbeiten Löschen Beitrag melden
Andreas F. #1722442
Lesenswert?
• ▲
▼ 
Karl heinz Buchegger schrieb:
> Bei all diesen bist du mit Ausschreiben
> (oder Makro oder Inline-Funktion) besser bedient (sowohl was Genauigkeit
> als auch Geschwindigkeit angeht) als mit pow.

Und wie so oft kommt es auf den Compiler an ;-)

Gerade ausprobiert, gcc 4.3.2 für x86 mit -O2 macht tatsächlich aus 
"pow(a, 2.0)" eine einfache Multiplikation:

1
   3:   dd 45 08                fldl   0x8(%ebp)

2
    return pow(a, 2.0);

3
   6:   d8 c8                   fmul   %st(0),%st


Bei 3.0 ruft er aber die Library-Funktion auf.

Andreas
Antwort Zitieren Bearbeiten Löschen Beitrag melden
Thread beobachten

Antwort schreiben

Hinweis: der ursprüngliche Beitrag ist mehr als 6 Monate alt. Bitte hier nur auf die ursprüngliche Frage antworten, für neue Fragen einen neuen Beitrag erstellen.

Formatierung
  • [c]C-Code[/c]
  • [code]Code in anderen Sprachen, ASCII-Zeichnungen[/code]
  • [math]Formel in LaTeX-Syntax[/math]
  • [[Titel]] – Link zu einem Artikel
  • Verweis auf einen anderen Beitrag: Rechtsklick auf dessen Beitrags-ID (#8078628), „Adresse kopieren“, in den Text einfügen

Mehr zur Formatierung

Längeren Sourcecode als Anhang, nicht in den Text. JPG nur für Fotos und Scans; Zeichnungen und Screenshots als PNG oder GIF. Siehe Bildformate.

Dein Name und deine E-Mail-Adresse

Schon registriert? Anmelden

Wird nicht angezeigt. Wir schicken dir einen Link; sobald du damit bestätigst, erscheint dein Beitrag.

Mit dem Abschicken legst du einen kostenlosen Account an und erkennst die Nutzungsbedingungen und die Datenschutzerklärung an.