GPX Daten auswerten für Abbiegehinweise

mit matplotlib, NumPy, pandas, SciPy, SymPy und weiteren mathematischen Programmbibliotheken.
Antworten
sroth
User
Beiträge: 5
Registriert: Donnerstag 6. Juli 2023, 13:28

Hallo zusammen,

bisher habe ich es geschafft zu ermitteln, in welchem Winkel 3 GPX-Koordinaten zueinander liegen.
Das möchte ich nutzen um zu berechnen an welchem Punkt (koordinate) des Tracks ein Abbiegehinweis
eingefügt werden soll.

Ich habe leider keine Idee, wie ich die Richtung (rechts/links) ermitteln soll. Vielleicht hat hier jemand einen Geistesblitz dazu?

Winkelberechnung:

Code: Alles auswählen

import math


# Beispielkoordinaten
lat1 = 51.37246756
lat2 = 51.37253485
lat3 = 51.37140207
lon1 = 6.5132431
lon2 = 6.51523736
lon3 = 6.51518346
lat12 = (lat1 + lat2)/2*0.01745
lat23 = (lat2 + lat3)/2*0.01745
lat13 = (lat1 + lat3)/2*0.01745
dx12 = 111.3 * math.cos(lat12) * (lon1 - lon2)
dy12 = 111.3 * (lat1 - lat2)
distance12 = math.sqrt(dx12 * dx12 + dy12 * dy12)
print(distance12)
dx23 = 111.3 * math.cos(lat23) * (lon2 - lon3)
dy23 = 111.3 * (lat2 - lat3)
distance13 = math.sqrt(dx23 * dx23 + dy23 * dy23)
print(distance13)
dx13 = 111.3 * math.cos(lat13) * (lon1 - lon3)
dy13 = 111.3 * (lat1 - lat3)
distance23 = math.sqrt(dx13 * dx13 + dy13 * dy13)
print(distance13)
b = distance12
a = distance23
c = distance13
winkelalpha = math.acos((b**2 + c**2 - a**2) / (2 * c * b))
winkelbetha = math.acos((a**2 + c**2 - b**2) / (2 * a * c))
winkelgamma = math.acos((a**2 + b**2 - c**2) / (2 * a * b))

print('alpha', winkelalpha / 0.01745)
print('betha', winkelbetha / 0.01745)
print('gamma',winkelgamma / 0.01745)
Vielleicht seid Ihr ja talentierter in mathe?

Vielen Dank vorab
Sven
__deets__
User
Beiträge: 14545
Registriert: Mittwoch 14. Oktober 2015, 14:29

Rechts/Links ist relativ zur Bewegung davor. Es muss also eine Gerade geben, in die man sich vorher bewegt hat. Die hat einen richtungsvektor. Den dreht man um 90 Grad (links oder rechts ist egal, muss nur immer gleich sein), als Normale. und das skalarprodukt deines Vektors in die neue Richtung mit der Normalen ist dann positiv oder negativ abhängig davon ob’s in die gleiche oder entgegengesetzte Richtung zeigt.
sroth
User
Beiträge: 5
Registriert: Donnerstag 6. Juli 2023, 13:28

Hier ist der berichtigte Code:

Code: Alles auswählen

import math

#  <wpt lat="51.37245634" lon="6.50727828">
#  <wpt lat="51.37261336" lon="6.52132793">
#  <wpt lat="51.36827275" lon="6.52157946">
#  <wpt lat="51.37567510" lon="6.52391508">

# Beispielkoordinaten
lat1 = 51.37245634
lat2 = 51.37261336
latrechts = 51.36827275
latlinks = 51.37567510
lon1 = 6.50727828
lon2 = 6.52132793
lonrechts = 6.5215794
lonlinks = 6.52391508
lat12 = (lat1 + lat2)/2*0.01745
dx12 = 111.3 * math.cos(lat12) * (lon1 - lon2)
dy12 = 111.3 * (lat1 - lat2)
distance12 = math.sqrt(dx12 * dx12 + dy12 * dy12)
b = distance12
print('Strecke 1-2 = ',b,'km')

lat2rechts = (lat2 + latrechts)/2*0.01745
dx23 = 111.3 * math.cos(latrechts) * (lon2 - lonrechts)
dy23 = 111.3 * (lat2 - latrechts)
distance23 = math.sqrt(dx23 * dx23 + dy23 * dy23)
a = distance23
print('Strecke 2-3 = ',a,'km')

lat1rechts = (lat1 + latrechts)/2*0.01745
dx13 = 111.3 * math.cos(lat1rechts) * (lon1 - lonrechts)
dy13 = 111.3 * (lat1 - latrechts)
distance13 = math.sqrt(dx13 * dx13 + dy13 * dy13)
c = distance13
print('Strecke 1-3 = ',c,'km')


winkelalpha = math.acos((b**2 + c**2 - a**2) / (2 * c * b))
winkelbetha = math.acos((a**2 + c**2 - b**2) / (2 * a * c))
winkelgamma = math.acos((a**2 + b**2 - c**2) / (2 * a * b))

print('rechts abbiegen:')
print('alpha', winkelalpha / 0.01745)
print('betha', winkelbetha / 0.01745)
print('gamma',winkelgamma / 0.01745)


lat12 = (lat1 + lat2)/2*0.01745
dx12 = 111.3 * math.cos(lat12) * (lon1 - lon2)
dy12 = 111.3 * (lat1 - lat2)
distance12 = math.sqrt(dx12 * dx12 + dy12 * dy12)
b2 = distance12
print('Strecke 1-2 = ',b2,'km')

lat2links = (lat2 + latlinks)/2*0.01745
dx23 = 111.3 * math.cos(latlinks) * (lon2 - lonlinks)
dy23 = 111.3 * (lat2 - latlinks)
distance23 = math.sqrt(dx23 * dx23 + dy23 * dy23)
a2 = distance23
print('Strecke 2-4 = ',a2,'km')

lat1links = (lat1 + latlinks)/2*0.01745
dx13 = 111.3 * math.cos(lat1links) * (lon1 - lonlinks)
dy13 = 111.3 * (lat1 - latlinks)
distance13 = math.sqrt(dx13 * dx13 + dy13 * dy13)
c2 = distance13
print('Strecke 1-4 = ',c2,'km')

winkelalpha = math.acos((b2**2 + c2**2 - a2**2) / (2 * c2 * b2))
winkelbetha = math.acos((a2**2 + c2**2 - b2**2) / (2 * a2 * c2))
winkelgamma = math.acos((a2**2 + b2**2 - c2**2) / (2 * a2 * b2))

print('links abbiegen:')
print('alpha', winkelalpha / 0.01745)
print('betha', winkelbetha / 0.01745)
print('gamma',winkelgamma / 0.01745)
__deets__ hat geschrieben: Donnerstag 6. Juli 2023, 18:47 Rechts/Links ist relativ zur Bewegung davor. Es muss also eine Gerade geben, in die man sich vorher bewegt hat. Die hat einen richtungsvektor. Den dreht man um 90 Grad (links oder rechts ist egal, muss nur immer gleich sein), als Normale. und das skalarprodukt deines Vektors in die neue Richtung mit der Normalen ist dann positiv oder negativ abhängig davon ob’s in die gleiche oder entgegengesetzte Richtung zeigt.
Kannst Du evtl. etwas konkreter werden wie ich das angehen kann. Das ist mir etwas zu abstrakt. Sorry
sroth
User
Beiträge: 5
Registriert: Donnerstag 6. Juli 2023, 13:28

Habe den Testcode nochmal überarbeitet, es ist aber scheinbar ein Fehler vorhanden? Vielleich könnt Ihr helfen?

Code: Alles auswählen

import math

#  <wpt lat="51.37245634" lon="6.50727828">
#  <wpt lat="51.37261336" lon="6.52132793">
#  <wpt lat="51.36827275" lon="6.52157946">
#  <wpt lat="51.37567510" lon="6.52391508">
# lon = Längengrad
# lat = Breitengrad

# Beispielkoordinaten
lat1 = 51.37245634
lat2 = 51.37261336
latrechts = 51.36827275
latlinks = 51.37567510
lon1 = 6.50727828
lon2 = 6.52132793
lonrechts = 6.5215794
lonlinks = 6.52391508
lat12 = (lat1 + lat2)/2*0.01745
dx12 = 111.3 * math.cos(lat12) * (lon1 - lon2)
dy12 = 111.3 * (lat1 - lat2)
distance12 = math.sqrt(dx12 * dx12 + dy12 * dy12)
b = distance12
print('Strecke 1-2 = ',b,'km')

# A = (lat1|lon1)
# C = (lat2|lon2]
# B = (latrechts|lonrechts)
# vac = ((lat2-lat1)|(lon2-lon1))
# vab = ((latrechts-lat1)|(lonrechts-lon1))

sp_rechts = (lat2-lat1) * (latrechts-lat1) + (lon2-lon1) * (lonrechts-lon1)
print('rechts abbiegen',sp_rechts)

sp_links = (lat2-lat1) * (latlinks-lat1) + (lon2-lon1) * (lonlinks-lon1)
print('links abbiegen',sp_rechts)

lat2rechts = (lat2 + latrechts)/2*0.01745
dx23 = 111.3 * math.cos(latrechts) * (lon2 - lonrechts)
dy23 = 111.3 * (lat2 - latrechts)
distance23 = math.sqrt(dx23 * dx23 + dy23 * dy23)
a = distance23
print('Strecke 2-3 = ',a,'km')

lat1rechts = (lat1 + latrechts)/2*0.01745
dx13 = 111.3 * math.cos(lat1rechts) * (lon1 - lonrechts)
dy13 = 111.3 * (lat1 - latrechts)
distance13 = math.sqrt(dx13 * dx13 + dy13 * dy13)
c = distance13
print('Strecke 1-3 = ',c,'km')


winkelalpha = math.acos((b**2 + c**2 - a**2) / (2 * c * b))
winkelbetha = math.acos((a**2 + c**2 - b**2) / (2 * a * c))
winkelgamma = math.acos((a**2 + b**2 - c**2) / (2 * a * b))

print('rechts abbiegen:')
print('alpha', winkelalpha / 0.01745)
print('betha', winkelbetha / 0.01745)
print('gamma',winkelgamma / 0.01745)


lat12 = (lat1 + lat2)/2*0.01745
dx12 = 111.3 * math.cos(lat12) * (lon1 - lon2)
dy12 = 111.3 * (lat1 - lat2)
distance12 = math.sqrt(dx12 * dx12 + dy12 * dy12)
b2 = distance12
print('Strecke 1-2 = ',b2,'km')

lat2links = (lat2 + latlinks)/2*0.01745
dx23 = 111.3 * math.cos(latlinks) * (lon2 - lonlinks)
dy23 = 111.3 * (lat2 - latlinks)
distance23 = math.sqrt(dx23 * dx23 + dy23 * dy23)
a2 = distance23
print('Strecke 2-4 = ',a2,'km')

lat1links = (lat1 + latlinks)/2*0.01745
dx13 = 111.3 * math.cos(lat1links) * (lon1 - lonlinks)
dy13 = 111.3 * (lat1 - latlinks)
distance13 = math.sqrt(dx13 * dx13 + dy13 * dy13)
c2 = distance13
print('Strecke 1-4 = ',c2,'km')

winkelalpha = math.acos((b2**2 + c2**2 - a2**2) / (2 * c2 * b2))
winkelbetha = math.acos((a2**2 + c2**2 - b2**2) / (2 * a2 * c2))
winkelgamma = math.acos((a2**2 + b2**2 - c2**2) / (2 * a2 * b2))

print('links abbiegen:')
print('alpha', winkelalpha / 0.01745)
print('betha', winkelbetha / 0.01745)
print('gamma',winkelgamma / 0.01745)
Ergebnis:

Code: Alles auswählen

Strecke 1-2 =  0.9765254879647761 km
rechts abbiegen 0.0002002688233061789
links abbiegen 0.0002002688233061789
Strecke 2-3 =  0.4832748508440217 km
Strecke 1-3 =  1.0975589715292424 km
rechts abbiegen:
alpha 26.12427535426178
betha 62.831400814729896
gamma 91.07828678744349
Strecke 1-2 =  0.9765254879647761 km
Strecke 2-4 =  0.3640149651888166 km
Strecke 1-4 =  1.210355537101945 km
links abbiegen:
alpha 14.746188523070716
betha 43.063892555040646
gamma 122.22388187832375
__deets__
User
Beiträge: 14545
Registriert: Mittwoch 14. Oktober 2015, 14:29

Das ist alles viel zu kompliziert. Mit Vektorrechnung wird das ganz simpel:

Gegeben drei Punkte A, B, C, die in der Reihenfolge abgelaufen werden. Die Frage ist, ob von B nach C gegenueber dem Teilstueck on A nach B ein Richtungswechsel stattgefunden hat, und ob das rechts oder links war.

Wenn wir die Notation P1P2 einfuehren fuer einen Vektor der von P2 nach P1 zeigt, dann ist also die Frage, ob Vektor CB von Vektor BA signifikant abweicht in der Richtung, und in welche Richtung.

Aus der Vektorrechnung ist weiterhin bekannt, dass das Skalarprodukt zweier normalisierter Vektoren (also Vektoren mit Laenge 1) den Winkel zwischen diesen beiden als Cosinus ergibt.

Normalisieren ist einfach nur das dividieren beider Komponenten durch die Laenge

Code: Alles auswählen

l = sqrt(x_BA**2 + y_BA**2)
nBA = x_BA / l, y_BA / l
x_BA und y_BA stehen hierbei fuer die jeweilige Komponente von BA.

Notieren wir nBA als normalisierten Vektor von BA, dann ist also einfach der Winkel zwischen den beiden aufeinanderfolgenden Schritten

Code: Alles auswählen

w = arcos(nBA*nCB)
Wenn wir jetzt wissen wollen, ob wir links oder rechts gelaufen sind, koennen wir ausnutzen, dass das Skalarprodukt zwischen zwei Vektoren

- positiv ist, wenn sie in die gleiche Richtung zeigen
- 0, wenn sie rechtwinklig zueinander sind
- negativ, wenn sie voneinander wegzeigen.

Um also zu wissen, ob CB links oder rechts von BA abzweigt, bilden wir zuerst eine Normale von BA, die nach links zeigt, lBA. Die berechnet sich einfach als

Code: Alles auswählen

x_lBA, y_lBA = -y_BA, x_BA
Kurze Gegenprobe: der Vektor 0, 1 zeigt im Koordinatensystem nach oben. Obige Formel angewandt ist dann -1, 0, zeigt also nach links. Der Vektor 1, 0 zeigt nach rechts. Obige Formel angewandt liefert 0, 1. Zeigt also nach oben, 90 grad links gedreht von rechts. Weiter Proben kannst du ja mal selbst durchspielen.

Damit stimmt also unsere Berechnung, und

Code: Alles auswählen

links = lBA * CB > 0
rechts = lBA * CB < 0
Achtung, es kann vorkommen, dass links & rechts False sind, denn dann geht CB in gerader Verlaengerung von BA!

Fang das an so umzusetzen, und zeig, was du daraus gemacht hast.
sroth
User
Beiträge: 5
Registriert: Donnerstag 6. Juli 2023, 13:28

Vielen Dank, das ist eine Grundlage.
Benutzeravatar
__blackjack__
User
Beiträge: 13572
Registriert: Samstag 2. Juni 2018, 10:21
Wohnort: 127.0.0.1
Kontaktdaten:

Ich würde noch ein zweites Achtung hintersetzen: Die Eigenschaften von Gleitkommazahlen können dazu führen das „geradeaus“ trotzdem ganz leicht von 0 abweicht. Aber Du willst ja sowieso einen *deutlichen* Winkel haben für eine Richtungsänderungsansage. Denn bei beispielsweise 3° Abweichung zum vorherigen Teilstück sagt man ja nicht „biegen sie jetzt bitte Links ab“ an. 😉
„Debugging is twice as hard as writing the code in the first place. Therefore, if you write the code as cleverly as possible, you are, by definition, not smart enough to debug it.“ — Brian W. Kernighan
sroth
User
Beiträge: 5
Registriert: Donnerstag 6. Juli 2023, 13:28

So, ich habe einmal etwas rumgebastelt.
Leider sehe ich jetzt noch nicht was falsch ist, da ich zu viele Bäume im Weg stehen habe...

Was ich habe:

Code: Alles auswählen

import math
import numpy

#  <wpt lat="51.37245634" lon="6.50727828">
#  <wpt lat="51.37261336" lon="6.52132793">
#  <wpt lat="51.36827275" lon="6.52157946">
#  <wpt lat="51.37567510" lon="6.52391508">
# lon = Längengrad
# lat = Breitengrad

# Beispielkoordinaten
lat1 = 51.37245634
lat2 = 51.37261336
latrechts = 51.36827275
latlinks = 51.37567510
lon1 = 6.50727828
lon2 = 6.52132793
lonrechts = 6.5215794
lonlinks = 6.52391508

x_BA = lat2 - lat1
y_BA = lon2 - lon1
l_AB = math.sqrt(x_BA**2 + y_BA**2)
print(l_AB)
nBA = x_BA / l_AB, y_BA / l_AB
print(nBA)

x_CA = (latlinks - lat2)
y_CA = (lonlinks - lon2)
l_AC = math.sqrt(x_CA**2 + y_CA**2)
print(l_AC)
nCA = x_CA / l_AC, y_CA / l_AC
print(nCA)

x_DA = (latrechts - lat2)
y_DA = (lonrechts - lon2)
l_AD = math.sqrt(x_DA**2 + y_DA**2)
print(l_AD)
nDA = x_DA / l_AD, y_DA / l_AD
print(nDA)


winkel_links = math.acos(numpy.dot(nBA, nCA))
print('Winkel links: ', winkel_links)
Ergebnis:

Code: Alles auswählen

0.014050527406574075
(0.01117538121235962, 0.9999375534775952)
0.004008440713057717
(0.7638231968919771, 0.6454255370604108)
0.0043478882613285814
(-0.9983260238324291, 0.05783727292094942)
Winkel links:  0.8580404461634757
Winkel rechts :  1.5241023733494934
Ich glaube ich habe das emfohlene gemacht??? Oder nicht? Diese Art Mathe ist nicht mein Horizont. Damit habe ich zum ersten Mal zu tun!
__deets__
User
Beiträge: 14545
Registriert: Mittwoch 14. Oktober 2015, 14:29

Nee. Das ist Quatsch mit „Winkel links“. Das gibt es nicht. Es gibt einen Winkel, weil der cosinus symmetrisch um die x-Achse ist, und daher genau nicht die Richtung unterschieden werden kann.

CA ist auch Unfug. Es gibt zwei interessante Vektoren - BA (erster Schritt), und CB (zweiter Schritt).

Warum sind da plötzlich 4 Punkte definiert?

Nimm dir ein Stück Papier. Mal dir 3 Punke darauf. Nimm deren Koordinaten & rechne das durch. Wenn du den Weg da mal nachvollzogen hast, Portier das für 3 Punkte nach Python.
Antworten