Profil d'itinéraire depuis un fichier GPX en python
Depuis une trace GPX, dressez un profil annoté automatiquement des villes traversées
Introduction
Amoureux du cyclisme, et du sport en général, bienvenue !
Cet article traite du rendu d'une trace GPX en profil d'itinéraire, comme ceux que l'on voit sur les grands tours par exemple.
La trace GPX que nous allons prendre en exemple dans cet article est celle d'une cyclosportive française : "la Pyrénéenne", édition 2023. Il s'agit d'une course qui comprend l'ascension du Col du Tourmalet.
La trace GPX utilisée dans cet article est disponible sur Kaggle ici.
Voici le rendu sur lequel nous allons aboutir à la fin de cet article :
Fichiers GPX
Tout d’abord, voyons ce qu’est un fichier GPX et comment il est structuré. GPX est un fichier XML utilisé pour échanger des données GPS entre applications. Un fichier GPX comprend les coordonnées GPS, ainsi que l'altitude d'une multitude de points (appelés waypoints) formant ainsi un itinéraire. Étant un format de sérialisation de données basé sur XML, un fichier GPX est structuré par des balises.
Pour simplifier les choses, disons qu'un fichier GPX a une structure à plusieurs niveaux. Il peut inclure une liste d'itinéraires (marqués par une balise rte). Pour chacun d'eux, le fichier peut présenter une liste de traces (tag trk), puis une liste de segments reflétant la réception satellite (tag rkseg) puis enfin, le maillage le plus fin, une liste de waypoints caractérisés par, entre autres, la latitude , la longitude et l'altitude (balise trkpt).
Voici un extrait du fichier gpx que nous utiliserons dans cet article :
<trk> <name>pyreneenne2023_V2</name> <trkseg> <trkpt lat="43.06821" lon="0.1469556"> <ele>545.000000</ele> </trkpt> <trkpt lat="43.068223" lon="0.14695"> <ele>545.000000</ele> </trkpt> ...
Un fichier GPX peut être bien plus complexe que celui-ci. Heureusement pour nous, il existe des bibliothèques qui nous permettront de les décomposer facilement. Dans notre cas, nous utiliserons gpxpy.
Import des librairies et des données
Dans le cadre de cet article, nous utiliserons les librairies pandas, numpy, certaines fonctions mathématiques pour calculer les distances et enfin mathplotlib afin de tracer notre profil.
import numpy as np # linear algebra
import pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)
from math import sin, cos, sqrt, atan2, radians
import os
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches
Comme dit précédemment, notre "parser" sera gpxpy. Nous utiliserons également le package Nominatim de la bibliothèque geopy pour résoudre les coordonnées GPS en emplacements physiques.
import gpxpy
import gpxpy.gpx
from geopy.geocoders import Nominatim
Importons notre fichier gpx et décomposons-le.
try:
gpx_file = open(os.path.join("/kaggle/input/gpx-trace", "pyreneenne2023.gpx"), "r")
gpx = gpxpy.parse(gpx_file)
except IOError:
print("Couldn't open the input GPX file. Ensure it's in the 'input' dir.")
raise()
Fonctions utiles autour des coordonnées GPS
Nous n'allons pas détailler leur utilisation maintenant, cependant notez que nous aurons besoin de 3 choses :
- déterminer un emplacement physique à partir de ses coordonnées GPS,
- à l'inverse, déterminer les coordonnées GPS d'un lieu physique,
- enfin, calculez la distance entre 2 points GPS.
def getLocationFromCoords(geolocator, latitude, longitude):
"""
Converting GPS coordinates to location
Parameter
----------
geolocator : Nominatim API instance
latitude : -
longitude : GPS coordinates in WGS84 projection
Return
----------
dictionnary of location features
"""
coord = str(latitude)+","+str(longitude)
location = geolocator.reverse(coord)
return location.raw
def getCoordsFromLocation(geolocator, cityName):
"""
Converting location to GPS coordinates
Parameter
----------
geolocator : a Nominatim API instance
cityName : name of location (string)
Return
----------
Location object
"""
coords = geolocator.geocode(cityName)
return coords
def getDistanceBetweenCoords(lat1, long1, lat2, long2):
"""
Computing distance in km between 2 points
Parameter
----------
lat1, long1 : latitude and longitude of first point
lat2, long2 : latitude and longitude of second point
Return
----------
Distance in Km
"""
earthRadius = 6371
c = 0
if all((lat1, long1, lat2, long2)):
latFrom = radians(lat1)
lonFrom = radians(long1)
latTo = radians(lat2)
lonTo = radians(long2)
dlon = lonTo - lonFrom
dlat = latTo - latFrom
a = sin(dlat / 2)**2 + cos(latFrom) * cos(latTo) * sin(dlon / 2)**2
c = 2 * atan2(sqrt(a), sqrt(1 - a))
return earthRadius * c
Segments
Sauvegardons à présent nos waypoints dans un dataframe et, au passage, calculons quelques données supplémentaires : kilomètres cumulés et segments.
L'itinéraire sur lequel nous travaillons fait environ 102 kilomètres de long. Cependant, il y a des milliers de waypoints dans le fichier gpx. C'est pourquoi, afin de gagner un temps de calcul précieux sur les tâches ultérieures, nous allons "découper" la piste par segment en fonction de ce que nous appellerons une unité de segment. Comme unité, nous choisirons 1 km.
trackDf = [] #final dataframe initialization
segUnit = 1 #segment unit, here 1 km
km = 0 #cumulated km initialization
lastLat = None #previously read latitude, None at the moment
lastLon = None #previously read longitude, None at the moment
#We iterate on gpx structure then ...
for track in gpx.tracks:
for segment in track.segments:
for point in segment.points:
# ... we compute the distance between current waypoint and last one
distFromLastStep = getDistanceBetweenCoords(lastLat, lastLon, point.latitude, point.longitude)
km = km + distFromLastStep #cumulated km is updated
# ... current segment determination
segment = km // segUnit
trackDf.append([segment, km, point.latitude, point.longitude, point.elevation])
# ... current coords become previous ones
lastLat = point.latitude
lastLon = point.longitude
trackDf = pd.DataFrame(trackDf)
trackDf.columns = ['segment', 'km', 'latitude', 'longitude', 'elevation']
trackDf.head(5)
Pourcentages d'inclinaison
Nous avons introduit précédemment l'unité de segment mais n'avons pas expliqué son utilité. Comme nous l'avons dit, il existe des milliers de waypoints et nous souhaitons tous les conserver afin de tracer une courbe bien détaillée. Cependant, nous souhaitons annoter notre profil d'itinéraire avec des informations de pente et pour cela, nous calculerons le pourcentage d'inclinaison du début à la fin d'un segment.
Calculons les poucentages d'inclinaison pour tous les segments dans un nouveau dataframe : slopesDf
slopesDf = [] # dataframe initialization
# For each segment we calculated previously ...
for segment in trackDf['segment'].unique():
#get first step of current segment then its elevation
firstSt = trackDf.loc[trackDf[trackDf['segment'].eq(segment)]['km'].idxmin()]
fromElv = firstSt['elevation']
#get last step of current segment then its elevation
lastSt = trackDf.loc[trackDf[trackDf['segment'].eq(segment)]['km'].idxmax()]
toElv = lastSt['elevation']
#determine current segment lenght (in km)
segmentLength = trackDf[trackDf['segment'].eq(segment)]['km'].max() - trackDf[trackDf['segment'].eq(segment)]['km'].min()
#calculate slope of current segment(in %)
segmentSlope = (toElv - fromElv) / (segmentLength * 1000) * 100
slopesDf.append([segment, segmentSlope])
slopesDf = pd.DataFrame(slopesDf)
slopesDf.columns = ['segment', 'slope']
slopesDf.head(5)
Nous disposons désormais d'un dataframe représentant notre trace, découpé par segment, et d'un autre composé des pentes associées à chacune d'elles. Comme vous pouvez le deviner, nous allons fusionner ces deux dataframes.
trackDf = pd.merge(trackDf, slopesDf, on='segment')
trackDf.head()
Emplacements physiques
Dans cette partie, nous allons traiter la dernière étape de ce que nous pouvons appeler l'enrichissement des données. En effet, nous souhaitons ajouter les noms des principaux lieux traversés par l'itinéraire.
Tout d’abord, prenons un exemple pour comprendre clairement ce dont nous allons discuter ensuite. Affichons le nom du lieu ayant pour latitude 42.911039 et pour longitude 0.177621. Ce waypoint est repris de l'itinéraire sur lequel nous travaillons :
# Nomatim API instanciation
geolocator = Nominatim(user_agent="MyApp")
# coords to location conversion
getLocationFromCoords(geolocator, 42.911039, 0.177621)
{'place_id': 103784147,
'licence': 'Data © OpenStreetMap ...,
'osm_type': 'way',
'osm_id': 153833856,
'lat': '42.9110386270338',
'lon': '0.17762102379162775',
'class': 'highway',
'type': 'secondary',
'place_rank': 26,
'importance': 0.10000999999999993,
'addresstype': 'road',
'name': 'Avenue du Tourmalet',
'display_name': 'Avenue du Tourmalet,
La Mongie, Bagnères-de-Bigorre,
Hautes-Pyrénées, Occitanie, France
métropolitaine, 65200, France',
'address': {'road': 'Avenue du Tourmalet',
'hamlet': 'La Mongie',
'village': 'Bagnères-de-Bigorre',
'municipality': 'Bagnères-de-Bigorre',
'county': 'Hautes-Pyrénées',
'ISO3166-2-lvl6': 'FR-65',
'state': 'Occitanie',
'ISO3166-2-lvl4': 'FR-OCC',
'region': 'France métropolitaine',
'postcode': '65200',
'country': 'France',
'country_code': 'fr'},
'boundingbox': ['42.9100848',
'42.9110690',
'0.1698093',
'0.1795048']}
Comme vous le remarquez, un emplacement possède de nombreuses propriétés. Celles qui nous intéresseront sont « type », « hamlet », « village » et « municipality ». Pour information, ce waypoint fait référence à la station bien connue, La Mongie, située sur les pentes du Tourmalet.
Un hameau est dans la zone administrative d'un village qui est, lui-même, dans la zone administrative d'une commune. Si la règle était aussi claire qu’elle n'y paraît, tout irait bien. Il faut simplement, dans ce cas, conserver le niveau le plus fin. Malheureusement, ce n'est pas le cas. Regardons un contre-exemple avec un autre waypoint de notre itinéraire :
getLocationFromCoords(geolocator, 42.919790, 0.194096)
{'place_id': 103719106,
'licence': 'Data © OpenStreetMap ...,
'osm_type': 'way',
'osm_id': 153833872,
'lat': '42.91980478971772',
'lon': '0.19409219130305133',
'class': 'highway',
'type': 'secondary',
'place_rank': 26,
'importance': 0.10000999999999993,
'addresstype': 'road',
'name': 'D 918',
'display_name': 'D 918,
Bagnères-de-Bigorre, Hautes-Pyrénées,
Occitanie, France métropolitaine,
65200, France',
'address': {'road': 'D 918',
'village': 'Bagnères-de-Bigorre',
'municipality': 'Bagnères-de-Bigorre',
'county': 'Hautes-Pyrénées',
'ISO3166-2-lvl6': 'FR-65',
'state': 'Occitanie',
'ISO3166-2-lvl4': 'FR-OCC',
'region': 'France métropolitaine',
'postcode': '65200',
'country': 'France',
'country_code': 'fr'},
'boundingbox': ['42.9197099',
'42.9205736',
'0.1934308',
'0.1963245']}
Comme on peut le constater, il n'y a pas de hameau ce qui fait que le niveau le plus fin est le village de "Bagnères-de-Bigorre". Le problème est que Bagnères-de-Bigorre est une ville assez éloignée de ce point, comme on peut le voir ci-dessous :
coords = getCoordsFromLocation(geolocator, 'Bagnères-de-Bigorre')
getDistanceBetweenCoords(42.919790, 0.194096, coords.latitude, coords.longitude)
16.576867273238797
Bagnères-de-Bigorre est à 16Km du waypoint que nous venons de tester. En effet, dans des regions peu denses, les découpages administratifs sont parfois assez étendus et on ne peut pas, simplement, associer à un waypoint le lieu découvert. C'est pourquoi, afin de ne conserver que les informations pertinentes, nous allons appliquer ces quelques règles :
- Le type de lieu peut prendre de nombreuses valeurs différentes : primary, secondary, zoo, house, fuel, unclassified, etc ... Pour garantir que nous nous concentrerons sur les lieux principaux, nous conserverons uniquement les lieux notés primary et secondary,
- Pour chaque emplacement restant, nous vérifierons si la commune rattachée se trouve à moins de 1,5km. De cette façon, nous éviterons d'afficher une ville lorsque l'itinéraire traverse un lieu sans nom éloigné d'elle.
Veuillez noter que le code suivant aurait pu être implémenté dans la boucle précédente. Toutefois, par souci de clarté, il a été isolé.
placesDf = [] # dataframe initialization
lastPlace = '' # last visited place name, currently empty
placeGroup = 0 # in order to group waypoints by name places, we
# will manage a group id
# For each segment ...
for segment in trackDf['segment'].unique():
# ... we get first step of current segment, its latitude, longitude and elevation
firstSt = trackDf.loc[trackDf[trackDf['segment'].eq(segment)]['km'].idxmin()]
fromElv = firstSt['elevation']
fromLat = firstSt['latitude']
fromLong = firstSt['longitude']
# ... we determine the location
location = getLocationFromCoords(geolocator, fromLat, fromLong)
# ... by default, place name is the municipality
placeName = location['address']['municipality']
# ... if current location is categoried as primary or secondary
if location['type'] in ('primary', 'secondary'):
# ... place name is overwritten by hamlet or, if not, village name
if "hamlet" in location['address']:
placeName = location['address']['hamlet']
elif "village" in location['address']:
placeName = location['address']['village']
# ... as we said before, we calculate reals coords of current municipality, then
# the distance to current waypoint
realCoords = getCoordsFromLocation(geolocator, placeName+' '+location['address']['municipality'])
realDist = getDistanceBetweenCoords(fromLat, fromLong, realCoords.latitude, realCoords.longitude)
# ... finally, we keep location if found place name is within 1.5km from current waypoint
if realDist < 1.5:
if placeName != lastPlace:
placeGroup += 1
lastPlace = placeName
placesDf.append([segment, placeName, fromElv, realDist, placeGroup])
placesDf = pd.DataFrame(placesDf)
placesDf.columns = ['segment', 'place', 'elevation', 'distreal', 'group']
placesDf.head(5)
Nous en avons presque terminé. Vous remarquerez peut-être, dans l'extrait ci-dessus, que certaines villes apparaissent plus d'une fois. Cela peut se produire lorsque l'itinéraire traverse une zone étendue. C'est pourquoi nous avons ajouté un groupe. Il suffit de supprimer les lignes en double sur la base du groupe et de conserver la plus proche uniquement.
placesDf = placesDf.sort_values(['group', 'distreal']).drop_duplicates(subset=['group'], keep='first')
Voici le résultat :
print(placesDf.to_markdown())
Traçage du profil
Nous disposons désormais de tous les éléments nécessaires pour tracer notre profil d’itinéraire.
Tout d'abord, en fonction des pourcentages de pente, nous souhaitons afficher une couleur différente. Pour ce faire, nous allons créer quelques tableaux :
# slopes rules divisions
slopesTable = [lambda x: x < 2,
lambda x: (x >= 2) & (x < 4),
lambda x: (x >= 4) & (x < 5),
lambda x: (x >= 5) & (x < 8),
lambda x: x >= 8,]
# slopes associated colors
slopesColor = ['palegreen', 'yellow', 'orange', 'orangered', 'maroon']
# slopes legend
slopesDescr = ['inf 2%', '2 ~ 4%', '4 ~ 5%', '5 ~ 8%', 'sup 8%']
Ensuite, comme vous le savez peut-être, l’altitude peut parfois apparaître imprécise. Pour atténuer cet effet, nous allons lisser légèrement nos données. Créons une fonction à cet effet :
def smoothProfile(signal,L=10):
res = np.copy(signal)
for i in range (1,len(signal)-1):
L_g = min(i,L)
L_d = min(len(signal)-i-1,L)
Li=min(L_g,L_d)
res[i]=np.sum(signal[i-Li:i+Li+1])/(2*Li+1)
return res
Enfin, il ne nous reste plus qu'à tracer le profil. Voici quelques points que nous allons ajouter pour améliorer le résultat :
- Hormis le départ et l'arrivée, nous n'afficherons pas deux noms de lieux consécutifs dans un même périmètre de 4 km. Si nous le faisions, il pourrait y avoir en effet des chevauchements.
- Nous afficherons en arrière-plan une seconde courbe grisée pour simuler un effet 3D. Bien sûr, nous aurions pu implémenter une projection 3D avec mathplotlib mais ce n'est pas le résultat que nous recherchions, du moins pour ce profil.
legend = [] # legend initialization
style = dict(size=10, color='gray') # style used for annotations
# We first define the canvas, then axis (right and top ones will be hidden)
fig, ax = plt.subplots(figsize=(12,4))
plt.xlabel("Kilometers")
plt.ylabel("Elevation")
ax.spines[['right', 'top']].set_visible(False)
# Rendering of backhround gray curve (for a 3D effect)
# This curve is shifted 1km further and 30m above
plt.fill_between(trackDf['km']+1, smoothProfile(trackDf['elevation'])+30, color='gray', zorder=0)
# For each slope category we defined, profile is filled and legend is updated
for i in range(0, len(slopesTable)):
plt.fill_between(trackDf['km'], smoothProfile(trackDf['elevation']), where = (slopesTable[i])(trackDf['slope']), color=slopesColor[i], zorder=1)
legend.append(mpatches.Patch(color=slopesColor[i], label=slopesDescr[i]))
# We want annotations of place names to be over profile
annotationsAnchor = trackDf['elevation'].max() * 1.1;
lastSegment = 0 #previously visited segment
# Place names annotation (every 4km, except for start and finish)
for index, row in placesDf.iterrows():
if row['segment'] > lastSegment + 4 or row['segment']==0 or row['segment']==placesDf['segment'].max():
ax.annotate(row['place'],
xy=(row['segment'], row['elevation']),
xytext=(row['segment'], annotationsAnchor),
arrowprops=dict(arrowstyle="-", color='lightgray'),
horizontalalignment = 'center',
rotation=90, **style)
lastSegment = row['segment']
# Shrink current axis by 20%
box = ax.get_position()
ax.set_position([box.x0, box.y0, box.width * 0.8, box.height])
# Put the legend to the right of the current axis
ax.legend(handles=legend, loc='center left', bbox_to_anchor=(1, 0.5))
plt.show()
Crédits
Miniature issue de vecstock sur Freepik
Retrouvez dans la rubrique "Nos datasets" toutes les données dont vous aurez besoin pour tester et pratiquer !
