Seite 1 von 1
"Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 17:12
von Strawk
Hallo!
Versuche gerade, die Vincenty-Formel in Python zu übersetzen. Bis zum groben Abstand "D" sieht es ganz vielversprechend aus, aber dazu fehlt noch eine Angabe von "w" im Bogenmaß. Wie kann ich einen Wert "w" "im Bogenmaß einsetzen"?
Hier der vorläufige Code:
Code: Alles auswählen
# -*- coding: utf-8 -*-
"""
Created on Fri Apr 26 10:45:34 2019
@author: Karl Kraft
"""
from __future__ import print_function, division
import numpy as np
import math
from math import radians
def vincenty(lat1, lat2, lon1, lon2):
""" calculates the distance in kilometers between two geographical points
"""
'''
lat1 = radians(lat1)
lat2 = radians(lat2)
lon1 = radians(lon1)
lon2 = radians(lon2)
'''
f = 1/298.257223563
# print(f)
a = 6378.137
F = ((lat1+lat2)/2)
G = ((lat1-lat2)/2)
l = ((lon1-lon2)/2)
S = ((math.sin(G))**2) * ((math.cos(l))**2) + ((math.cos(F))**2) * ((math.sin(l))**2)
C = ((math.cos(G))**2) * ((math.cos(l))**2) + ((math.sin(F))**2) * ((math.sin(l))**2)
w = np.arctan(math.sqrt(S/C))
D = 2 * 0.0151440 * a
print(w)
return D
if __name__ == "__main__":
'''
westlicher Ort: Position 1
östlicher Ort: Position 2
'''
'''
# Madrid (1) - Budapest (2)
lat1 = 40.4167754
lat2 = 47.497912
lon1 = -3.7037902
lon2 = 19.040235
'''
'''
# New York (1) - Aachen (2)
lat1 = 40.712784
lat2 = 50.7753455
lon1 = -74.005941
lon2 = 6.0838868
'''
'''
# Berlin (1) - Pretoria (2)
lat1 = 52.520007
lat2 = -25.746111
lon1 = 13.404954
lon2 = 28.188056
'''
# Pretoria (1) - Moskau (2)
lat1 = -25.746111
lat2 = 55.755826
lon1 = 28.188056
lon2 = 37.617300
'''
# Aachen (1) - Mönchengladbach (2)
lat1 = 50.7753455
lat2 = 51.1804572
lon1 = 6.0838868
lon2 = 6.4428041
'''
result = vincenty(lat1, lat2, lon1, lon2)
print("km: %.1f" % (result), "\t")
Grüße
Strawk
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 17:16
von Sirius3
Du mußt nichts einsetzen, da Du ja `w` als Ergebnis bekommst. arctan2 wäre da auch besser.
Und wie immer, Zahlen nicht irgendwo im Code verstecken, sondern als Konstanten am Anfang definieren.
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 17:31
von Strawk
Artan2 braucht aber 2 Werte. Und meinst du D=2*a ?
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 17:38
von Sirius3
Genau, 2 Werte sind besser als einer, weil man durch 0 nicht teilen kann.
Warum benutzt Du numpy und math vermischt?
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 17:41
von Strawk
Und woher soll der zweite Wert kommen? Weil ich es eben nicht besser kann.
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 17:46
von Sirius3
Du hast doch zwei Variablen.
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 17:59
von Strawk
Ich habe arctan von Wurzel aus S/C. Also eine.
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 18:03
von Strawk
Modifizierter Code:
Code: Alles auswählen
# -*- coding: utf-8 -*-
"""
Created on Fri Apr 26 10:45:34 2019
@author: Karl Kraft
"""
from __future__ import print_function, division
import numpy as np
def vincenty(lat1, lat2, lon1, lon2):
""" calculates the distance in kilometers between two geographical points using
the Vincenty-formula
"""
'''
lat1 = radians(lat1)
lat2 = radians(lat2)
lon1 = radians(lon1)
lon2 = radians(lon2)
'''
f = 1/298.257223563
a = 6378.137
F = ((lat1+lat2)/2)
G = ((lat1-lat2)/2)
l = ((lon1-lon2)/2)
S = ((np.sin(G))**2) * ((np.cos(l))**2) + ((np.cos(F))**2) * ((np.sin(l))**2)
C = ((np.cos(G))**2) * ((np.cos(l))**2) + ((np.sin(F))**2) * ((np.sin(l))**2)
w = np.arctan(np.sqrt(S/C))
D = 2 * w * a
T = (np.sqrt(S*C))/w
H1 = (3*T-1) / (2*C)
H2 = (3*T+1) / (2*S)
s = D * ((1+f*H1*(np.sin(F))**2) * (np.cos(G)**2) - f * H2 * (np.cos(F)**2) * (np.sin(G)**2))
return s
if __name__ == "__main__":
'''
westlicher Ort: Position 1
östlicher Ort: Position 2
'''
# Madrid (1) - Budapest (2)
lat1 = 40.4167754
lat2 = 47.497912
lon1 = -3.7037902
lon2 = 19.040235
'''
# New York (1) - Aachen (2)
lat1 = 40.712784
lat2 = 50.7753455
lon1 = -74.005941
lon2 = 6.0838868
'''
'''
# Berlin (1) - Pretoria (2)
lat1 = 52.520007
lat2 = -25.746111
lon1 = 13.404954
lon2 = 28.188056
'''
'''
# Pretoria (1) - Moskau (2)
lat1 = -25.746111
lat2 = 55.755826
lon1 = 28.188056
lon2 = 37.617300
'''
'''
# Aachen (1) - Mönchengladbach (2)
lat1 = 50.7753455
lat2 = 51.1804572
lon1 = 6.0838868
lon2 = 6.4428041
'''
result = vincenty(lat1, lat2, lon1, lon2)
print("km: %.1f" % (result), "\t")
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 18:32
von Sirius3
so:
Code: Alles auswählen
EARTH_RADIUS = 6378.137
EARTH_FLATTENING = 1/298.257223563
def vincenty(lat1, lat2, lon1, lon2):
""" calculates the distance in kilometers between two geographical points using
the Vincenty-formula
"""
lat1 = np.radians(lat1)
lat2 = np.radians(lat2)
lon1 = np.radians(lon1)
lon2 = np.radians(lon2)
f = (lat1 + lat2) / 2
g = (lat1 - lat2) / 2
l = (lon1 - lon2) / 2
sin_g2, cos_g2 = np.sin(g) ** 2, np.cos(g) ** 2
sin_f2, cos_f2 = np.sin(f) ** 2, np.cos(f) ** 2
sin_l2, cos_l2 = np.sin(l) ** 2, np.cos(l) ** 2
s = sin_g2 * cos_l2 + cos_f2 * sin_l2
c = cos_g2 * cos_l2 + sin_f2 * sin_l2
s_sqrt = s ** 0.5
c_sqrt = c ** 0.5
w = np.arctan2(s_sqrt, c_sqrt)
d = 2 * w * EARTH_RADIUS
t = 3 * (s_sqrt * c_sqrt) / w
h1 = (t - 1) / (2*c)
h2 = (t + 1) / (2*s)
return d * ((1 + EARTH_FLATTENING * h1 * sin_f2) * cos_g2
- EARTH_FLATTENING * h2 * cos_f2 * sin_g2)
Re: "Dabei ist w im Bogenmaß einzusetzen"
Verfasst: Sonntag 16. Juni 2019, 19:41
von Strawk
