Mathematik ist einfach, sie auf Rechner zu bringen nicht…

Mathe ist ja eigentlich ganz einfach. Die nachfolgende Routine berechnet die kürzeste Distanz zwischen zwei Punkten auf der Erde unter der Berücksichtigung, dass die Erde eine perfekte Kugel ist.
Benutzt wird diese Routine um die Luftlinie zwischen zwei geokodierten Adressen zu ermitteln.

double DistanceBetweenCoordinates(
          double dLatitude1, double dLongitude1, 
          double dLatitude2, double dLongitude2)
{
  // Convert to radient
  dLatitude1 *= M_PI / 180;
  dLongitude1 *= M_PI / 180;
  dLatitude2 *= M_PI / 180;
  dLongitude2 *= M_PI / 180;

  // Get distance
  double dDistance;
  dDistance = sin(dLatitude1) * sin(dLatitude2) + cos(dLatitude1) * 
        cos(dLatitude2) * cos(dLongitude2-dLongitude1);
  // Ohne diese Zeile haben wir ein Problem... dDistance ist evtl. >1
  // dDistance = min(dDistance,1.0);
  dDistance = acos(dDistance) * EARTH_RADIUS;
  return dDistance;
}

Eigentlich ganz einfach. Aber dieser Code hat ein massives Problem, denn dDistance kann vor dem acos durch Rundungsfehler größer als 1 werden.

Für die folgenden Werte passiert etwas Übles:

dLatitude1 0.86208095750908087 
dLongitude1 0.15333764167657613 
dLatitude2 0.86208095750908087 
dLongitude2 0.15333764167657613

Eigentlich müsste glatt 1 herauskommen, denn die Koordinaten sind gleich, aber dDistance wird bei diesen Werten 1.0000000000000002. Der Sinus bzw. Cosinus ergibt eben doch leicht abweichende Werte. Der folgende acos scheitert dann aber und das Resultat wird NaN (Not a Number). Erlaubt sind natürlich nur Werte <=1.0 und >=-1.0. Die Folge ist, dass die Distanz nicht 0.0 ist.

Also muss hier eine kleine Sicherung eingebaut werden, die Aufgrund der Rundungs-Genauigkeitsfehler, das überschreiten von 1 verhindert.
Eigentlich ist Mathe ganz einfach, aber auf einem Rechner ist das alles manchmal ganz anders.

PS: Ich weiß schon warum ich noch nie Fließkommaarithmetik auf PCs mochte.

Ein Gedanke zu „Mathematik ist einfach, sie auf Rechner zu bringen nicht…“

  1. Es gibt noch eine weitere Formel, die weniger anfällig für Rundungsfehler ist, aber zu gleichem Ergebnis führt:

    d = 2*asin(sqrt((sin((lat1-lat2)/2))^2 + cos(lat1)*cos(lat2)*(sin((lon1-lon2)/2))^2))

    Das ergibt die Distanz in Bogenmaß und kann mittels

    d_nm = (pi / (180 * 60)) * d

    in Nautische Meilen oder weiter in andere Entfernungseinheiten konvertiert werden.

Schreibe einen Kommentar

Deine E-Mail-Adresse wird nicht veröffentlicht. Erforderliche Felder sind mit * markiert

I accept that my given data and my IP address is sent to a server in the USA only for the purpose of spam prevention through the Akismet program.More information on Akismet and GDPR.

Diese Website verwendet Akismet, um Spam zu reduzieren. Erfahre mehr darüber, wie deine Kommentardaten verarbeitet werden.