Fórmules de Vincenty

Les fórmules de Vincenty formen un algorisme molt eficient per al càlcul de la distància entre dos punts de la superfície d'un el·lipsoide de revolució. Són utilitzades àmpliament en geodèsia per calcular distàncies sobre la superfície de la Terra atès que requereix un nombre petit d'operacions tot i donar una precisió de 0.5mm (0,000015 "), molt millor que el mètode tradicional de la fórmula del haversine que s'utilitzava anteriorment en trigonometria esfèrica. L'algorisme va ser publicat per Thaddeus Vincenty el 1975.

Utilitza un mètode iteratiu.

Algorisme

El següent algorisme expressa les fórmules d'una forma senzilla de calcular:

a, b = semieixos major i menor del el·lipsoide
φ1, φ2 = latitud geodèsica
L = diferència en longitud
f = aplatament del el·lipsoide (a-b)/a
U1 = atan((1−f).tanφ1) (U is ‘latitud reduïda’) 
U2 = atan((1−f).tanφ2)
λ = L, λ′ = 2π
do while abs(λ−λ′) > 10−12 (implica un error < 0.06mm)
{
sinσ = √[ (cosU2.sinλ)² + (cosU1.sinU2 − sinU1.cosU2.cosλ)² ]
cosσ = sinU1.sinU2 + cosU1.cosU2.cosλ
σ = atan2(sinσ, cosσ)
sinα = cosU1.cosU2.sinλ / sinσ
cos²α = 1 − sin²α (trig identity; §6) 
cos2σm = cosσ − 2.sinU1.sinU2/cos²α
C = f/16.cos²α.[4+f.(4−3.cos²α)]
λ′ = λ
λ = L + (1−C).f.sinα.{σ+C.sinσ.[cos2σm+C.cosσ.(−1+2.cos²2σm)]} (emprar cos2σm=0 si se està al llarg de l'equador)
}
u² = cos²α.(a²−b²)/b²
A = 1+u²/16384.{4096+u².[−768+u².(320−175.u²)]}
B = u²/1024.{256+u².[−128+u².(74−47.u²)]}
Δσ = B.sinσ.{cos2σm+B/4.[cosσ.(−1+2.cos²2σm) − B/6.cos2σm.(−3+4.sin²σ).(−3+4.cos²2σm)]}
s = b.A.(σ−Δσ)
α1 = atan2(cosU2.sinλ, cosU1.sinU2 − sinU1.cosU2.cosλ)
α2 = atan2(cosU1.sinλ, −sinU1.cosU2 + cosU1.sinU2.cosλ)

On:
*s és la distància (mateixes unitats que a, b)
*α1 és l'azimut inicial
*α2 és l'azimut final (en direcció p1→p2)

La fórmula pot no tenir solució per a dos punts gairebé antipodals. Limitar el nombre d'iteracions per evitar aquest cas.

L'el·lipsoide més estès és el WGS84, per al qual:

a = 6.378.137 m (± 2 m)
b = 6 356 752.3142 m
f = 1/298.257223563

Exemple per una prova

Cas de prova extret de Geoscience Austràlia, usant WGS-84:

Flinders Peak 37° 57'03 72.030" S, 144° 25'29 52.440" E
Buninyong 37° 39'10 15.610" S, 143° 55'35 38.390" E
Resultat:
s 54.972,271m
α1 306° 52'05 37"
α2 127° 10'25 07"(≡ 307 · 10'25 07" p1 → p2)

Vegeu també

• Algorisme de navegació

Referències

• Direct and Inverse Solutions of Geodesics on the ellipse with application of nested equations . Publicació original de Vincenty.

Enllaços externs

• Codi font disponible (llicència LGPL) per Chris Veness