Affichage de l'altitude DEM-MNT

[mirelo]

Evidemment, puisqu'il a été écrit seulement pour les fichiers 3" (c'est une ébauche !).
Je n'ai pas cherché beaucoup sur le site que tu m'avais indiqué, pensant que les fichiers 1" étaient seulement pour l'Amérique.
Je vais corriger ça sans tarder.

[jcr83]

Ah, OK ! Ce n'était pas précisé dans ton message, j'ai donc utilisé les fichiers au pas de 1" qui étaient déjà installés sur mon PC.

[mirelo]

Et pour moi, les fichiers 3" étaient une évidence !
La nouvelle version 2Altihgt1.py

2Altihgt1.zip

(1.15 KiB) Downloaded 447 times

accepte les 2 valeurs, mais comme les fichiers portent le même nom, attention à la confusion.
J'ai laissé l'image "brute", 1 pixel pour 1 m en altitude, 1 pixel pour 1 "pas" à l'horizontale, ce qui donne des images éventuellement très aplaties, mais sans altérations.

Juste une petite remarque :
Pour apprécier la qualité de l'interpolation (c'est quand même le but premier de ce script) une grille à 3" serait, me semble-t-il, plus parlante.

[jcr83]

Cette nouvelle version fonctionne bien avec les fichiers au pas de 1". Je vais coder l'algorithme en C. Si je te fournis le code source, as-tu les outils nécessaires pour le compiler ? Si oui, cela permettra de comparer les résultats avec ceux de ton programme en Python.

[mirelo]

N'étant pas vraiment familier du C, je préfère te proposer un jeu d'essai sur le fichier N43E006.hgt au pas de 1",
avec tous les résultats intermédiaires.
lat, lon : 43.53184 6.60626
X-X0 = 0.535999999999
Zi0-1 [290 , 291 , 292 , 293 ]
A, B, C,Z : 0.0, 0.0, 1.0, 291.53599999999915
Zi0 [301 , 302 , 303 , 303 ]
A, B, C,Z : -0.5, 0.5, 1.0, 302.60265267199907
Zi0+1 [304 , 305 , 305 , 306 ]
A, B, C,Z : 1.0, -1.5, 0.5, 304.99104665600021
Zi0+2 [306 , 307 , 308 , 309 ]
A, B, C,Z : 0.0, 0.0, 1.0, 307.53599999999915

Y-Y0 = 0.624000000009
Zi [291.53599999999915, 302.60265267199907, 304.99104665600021, 307.53599999999915]
A, B, C,Z : 4.417409024, -8.756538368, 6.727523328, 304.464341972
Bon codage

[mirelo]

Après m'être penché un peu plus attentivement sur le problème des "effets de bord", dus à la présence de zones sans données (No-data areas), je voudrais compléter mon algorithme de la manière suivante :
Pour chaque étape du calcul :
Si c'est un des 2 points du segment central qui manque, on renvoie l'altitude fictive –32768, soit :
si Zi0 == –32768 : Zi = –32768
si Zi0+1 == –32768 : Zi = –32768
Si c'est un des points extérieurs qui manque (ou les 2), on le( ou les) remplace par son (ou leur) "ombre", située sur la droite [Zi0 Zi0+1], soit :
si Zi0-1 == -32768 : Zi0-1 = 2 * Zi0 – Zi1
si Zi0+2 == -32768 : Zi0+2 = 2 * Zi1 – Zi0
On utilise ainsi le maximum d'information disponible pour faire les calculs en mode dégradé, intermédiaire entre les interpolations bi-cubique et bi-linéaire.

J'ai par ailleurs réalisé un petit programme Srtmign0.py
http://dl.free.fr/rIeo8grbN,
toujours en Python, destiné à boucher les trous de la grille SRTM.hgt (3") sur le territoire français, en calculant les valeurs manquantes à partir de la grille MNT250-93 de l'IGN disponible ici :
http://professionnels.ign.fr/DISPLAY/000/527/330/5273304/France_250_ASC_L93.zip
Cette grille est complète, mais au pas de 250 m, et s'arrête aux frontières. Il peut donc rester des trous, mais hors de France.
Bien sûr, ce ne sont que des valeurs estimées, mais elles m'ont permis cette fois d'atteindre l'Aiguille du Midi (avec des erreurs raisonnables !).

Pour cela, on place dans un même dossier : Srtmign0.py, la grille MNT250_L93_FRANCE.ASC et les fichiers 3".hgt,
puis on répond Lat et Long correspondant au fichier à traiter et on n'a plus qu'à attendre un peu pour voir le résultat.

[jcr83]

J'ai reprogrammé l'algorithme en C. En voici une illustration :

interpolation bicubique

bicubic.jpg (23.36 KiB) Viewed 7726 times

En bleu, le résultat de l'interpolation bilinéaire, et en rouge celle de l'interpolation bicubique. On voit bien l'avantage de cette dernière.
Je l'intégrerai donc dans la prochaine version de ttMaps. Merci à Mirelo pour les formules.

[mirelo]

Je n'avais pas fait de simulation pour l'effet de bord (mea culpa ) et je me rends compte en regardant celle-ci qu'on peut avoir un meilleur algorithme pour le cas des données manquantes :

si Zi0-1 == -32768 :
Zi0-1 = Zi0+2 - 3 * (Zi0+1 – Zi0)
ou directement (ce qui conduit au même résultat) :
Ai = 0
Bi = (Zi0+2 + Zi0) / 2 – Zi0+1
Ci = 2 * Zi0+1 – (3 * Zi0 + Zi0+2) / 2

si Zi0+2 == -32768 :
Zi0+2 = Zi0-1 + 3 * ( Zi0+1 – Zi0)
ou directement (ce qui conduit au même résultat) :
Ai = 0
Bi = (Zi0+1 + Zi0-1) / 2 – Zi0
Ci = (Zi0+1 - Zi0-1) / 2

et enfin (comme précédemment) :
si Zi0-1 et Zi0+2 == -32768 :
Zi0-1 = 2 * Zi0 – Zi0+1
Zi0+2 = 2 * Zi0+1 – Zi0
ou directement (ce qui conduit au même résultat) :
Ai = 0
Bi = 0
Ci = Zi0+1 – Zi0

On fait ainsi une interpolation parabolique à partir des 3 points connus, ce qui revient à supposer que la variation de pente est uniforme et non plus contrainte à tendre arbitrairement vers la corde.
L'interpolation linéaire n'est alors utilisée que dans le cas exceptionnel où seuls 2 points sont connus et là, on ne peut vraiment rien faire de plus !
La différence sautera aux yeux, même si elle n'est pas très importante en terme d'altitude.

Interpolation.jpg (4.57 KiB) Viewed 7698 times

Cette fois, je pense vraiment qu'on ne peut pas faire mieux

[marcuslecactus]

je déterre ce vieux sujet...
dimanche dernier, l'altitude enregistrée dans la trace a perdu le nord : données manquantes et farfelues, la montagne, les sapins et les nuages n'ont pas aidées...
j'exploite ensuite la trace avec un tableur, du coup mon profil était tout moche.
j'ai bien utilisé http://www.gpsvisualizer.com pour récupérer des données altimétriques correctes (enfin plausibles, au vu du profil ainsi obtenu), mais cette application remplace juste la série des altitudes GPS par celle des DEM (parmi 3 sources).
je voudrais juste ajouter une colonne "altitude DEM" en choisissant ma source de DEM (par exemple http://www.viewfinderpanoramas.org/dem3.html vu que je suis en plein dedans). Et ainsi comparer la précision des résultats.
Malheureusement 2Altihgt1.py n'est plus disponible en téléchargement, je ne sais si ce script aurait pu m'aider.
Peut-être y a t'il d'autres utilitaires qui peuvent m'aider ?
merci

[mirelo]

@marcuslecactus : J'ai ressorti 2Altihgt1.py pour le remettre sur le forum, bien que je doute un peu qu'il corresponde à tes attentes. Il s'agissait en effet d'un programme d'essai destiné à mettre en œuvre mon algorithme d'interpolation bicubique simplement le long d'un segment de droite reliant 2 points déterminés.
A la même époque, un autre sujet http://ttmaps.free.fr/forums/viewtopic.php?f=21&t=135, qui avait été relégué dans le chapitre des sujets divers (Divers < Autres sujets de discussion), parlait du profil des itinéraires. A cette occasion, j'avais modifié le programme précédent afin de lui faire enchaîner le calcul pour tous les segments d'un itinéraire. Il s'agit de 2Profiliti0c.py qui avait également disparu des écrans radar et que tu pourrais maintenant utiliser, après avoir transformé ta trace en itinéraire.