Représentation de cartes topographiques en 3D

Guillaume Duriaud

- Développement d'un petit logiciel de représentation de cartes topographiques en 3D avec Python.
- Exemple d'utilisation de la bibliothèque Vpython.
Réagissez à cet article : Commentez ♪

- Titre : Représentation de cartes topographiques en 3D

- Auteur : Guillaume Duriaud

- Parution : 4 juin 2009

- Mise à jour : 8 juin 2009

- Licence : Copyright ® 2009 Guillaume Duriaud. Aucune reproduction, même partielle, ne peut être faite de ce site et de l'ensemble de son contenu : textes, documents, images, etc. sans l'autorisation expresse de l'auteur. Sinon vous encourez selon la loi jusqu'à trois ans de prison et jusqu'à 300 000 € de dommages et intérêts.

Sommaire

- I. Prérequis

- II. Introduction

- III. Explication du code

- III-1. Importation

- III-2. La fonction divisetab

- III-3. La classe TkScene

- III-3-1. L'initialisation

- III-3-2. La fonction Load

- III-3-3. La fonction Affiche

- III-3-4. La Fonction Refresh

- IV. Téléchargement

I - Prérequis

Voici les versions de bibliothèques utilisées pour développer ce programme :

Langage: Python 2.6.2 ;
Bibliothèques: Numpy 1.3.0, VPython 5.1.

Le programme devrait toutefois fonctionner avec d'autres versions de Python. Pour VPython, il est conseillé tout de même de posséder la dernière, de nombreux changements étant apportés entre les versions apportant son lot d'incompatibilités. Le programme n'a été testé que sur Windows, mais devrait également fonctionner sous Linux ou Mac.

Concernant la bibliothèque VPython, vous pouvez consulter le manuel de référence ou sa traduction.

II - Introduction

Lorsque j'étais en classe préparatoire, j'avais choisi comme sujet de TIPE un sujet en informatique « Reconstitution de Tumulus en 2D et 3D à partir de données métriques » et réalisé en CAML Light. Le travail était assez fastidieux vu qu'il faut à peu près tout gérer soi-même.
Étant depuis passé au Python, je vous propose le programme Python qui reprend ce sujet en étudiant ici uniquement la partie 3D.
Ce programme travaille sur des données topographiques à partir un quadrillage métrique, les numéros de lignes et de colonnes correspondent aux coordonnées X et Y des points topographiques et la valeur de la matrice à la coordonnée en Z. Toutes les coordonnées sont relatives à une référence quelconque et la coordonnée en Z est supposée être une profondeur.
Un niveau de gris (de 0 à 255) est appliqué à la scène en fonction de l'altitude des facettes. Il n'est pas encore possible d'avoir la scène en fil de fer.

Le code d'environ 120 lignes est regroupé dans un unique fichier pytopo3d.py. Le programme n'étant pas très long, il ne devrait poser aucun problème de compréhension.

Pour tester le programme, il suffit d'exécuter le fichier pytopo3d.py (en ayant bien entendu au préalable installé les dépendances) et depuis le menu Fichier de charger une carte. Il y a trois cartes d'exemples (tab2.txt, tab3.txt et tab8.txt).

III - Explication du code

III-1 - Importation

from __future__ import division
from visual import *
import Tkinter
import Tix
import tkFileDialog
import numpy
import copy

La première ligne permet d'utiliser la nouvelle définition de la division qui sera opérationnelle dans Python 3, à savoir que la division / sera flottante.
On importe ensuite les différents modules utiles Tkinter, Visual et NumPy.

III-2 - La fonction divisetab

def divisetab(tab):
rc, cc = tab.shape
newtab = numpy.zeros((2 * rc-1, 2*cc-1))
for i in range(rc):
for j in range(cc):
newtab[2*i, 2*j] = tab[i, j]
if (i<rc-1): newtab[2*i+1, 2*j] =(tab[i, j] + tab[i+1, j]) / 2
if (j<cc-1): newtab[2*i, 2*j+1] =(tab[i, j] + tab[i, j+1]) / 2
if (i<rc-1) and (j<cc-1): newtab[2*i+1, 2*j+1] =(tab[i, j] + tab[i, j+1] + tab[i+1, j] + tab[i+1, j+1]) / 4
return newtab

Cette fonction permet de dédoubler un tableau NumPy représentant la matrice du relevée topographique, en insérant entre chaque élément du tableau un nouvel élément égal à la moyenne de ses deux voisins. Pour le point central, ce sera la moyenne des quatre points les plus proches. Ceci permettra lors de la colorisation de la scène d'améliorer le dégradé des couleurs. Si les pentes ne sont pas trop abruptes, l'effet ne sera toutefois pas réellement visible à l'écran.

III-3 - La classe TkScene

III-3-1 - L'initialisation

class TkScene3D:
def __init__(self):
self.root = Tix.Tk()
self.root.title('Options')

mainmenu = Tkinter.Menu(self.root)
menuFile = Tkinter.Menu(mainmenu, tearoff = 0)
menuFile.add_command(label="Charger une carte", command = self.Load)
mainmenu.add_cascade(label = "Fichier", menu = menuFile)
self.root.config(menu = mainmenu)

self.lblLissage = Tkinter.Label(self.root, text = "Lissage")
self.lblLissage.grid(column = 0, row = 0, sticky = Tkinter.W)
self.cbbLissage = Tix.ComboBox(self.root, editable = 1, dropdown = 1)
self.cbbLissage.grid(column = 1, row = 0, sticky = Tkinter.W)
self.cbbLissage.insert(0, "Aucun")
self.cbbLissage.insert(1, "Un peu")
self.cbbLissage.insert(2, "Beaucoup")
self.cbbLissage.pick(0)
self.lblZ = Tkinter.Label(self.root, text = "Echelle Z")
self.lblZ.grid(column = 0, row = 1, sticky = Tkinter.W)
self.eZ = Tkinter.Entry(self.root)
self.eZ.insert(Tkinter.END, "1")
self.eZ.grid(column = 1, row = 1, sticky = Tkinter.W)
self.btnRefresh = Tkinter.Button(self.root, text = "Affiche", command = self.Refresh)
self.btnRefresh.grid(column = 0, row = 4, columnspan = 2)
self.f = frame()
self.model = None
self.tab = {"X":1, "Y": 1, "Z": 1, "Array":numpy.array([])}

self.eZ.bind('<Return>', self.Refresh)

L'interface graphique est composée d'une scène VPython et d'une petite fenêtre Tkinter où quelques options permettront de modifier cette dernière.
La fenêtre Tkinter contient un menu menuFile à partir duquel vous chargerez une carte.
La ComboBox cbbLissage définit le nombre de divisions du tableau représentant la carte, ce qui jouera sur la qualité du dégradé des couleurs. Pour l'option « Aucun », on travaille avec le tableau brut. Avec l'option « Un peu », on fait subir une division au tableau et avec l'option « Beaucoup », on lui fait subir deux divisions successives.
l'Entry eZ permet de modifier l'échelle en Z. Au chargement de la carte, on se place en coordonnées normées. En modifiant la valeur de eZ, vous pouvez ainsi mettre en évidence plus nettement les reliefs, ce qui est notamment primordial lorsqu'on travaille avec des relevés topographiques de tumulus, où le terrain est quasiment plat.
Le Button btnRefresh permet alors de rafraîchir la scène.

L'attribut tab est un dictionnaire qui contient les paramètres X, Y et Z de la carte à savoir la distance entre deux points selon l'axe X, la distance entre deux points selon l'axe Y et l'échelle des valeurs en Z, le tout étant exprimé dans la même unité.

La scène va utiliser plusieurs attributs :
model est un objet VPython qui va contenir un ensemble de petites surfaces triangulaires ;
f, une frame, est une structure dans laquelle on place model.

III-3-2 - La fonction Load

def Load(self):
strfile = tkFileDialog.askopenfilename(filetypes = [("All", "*"),("Fichiers txt","*.txt")])
a = open(strfile)
fl = map(float, a.readline().split(";"))
self.tab["X"] = fl[0]
self.tab["Y"] = fl[1]
self.tab["Z"] = fl[2]
j = 0
ltab = []
for i in a:
ltab.append([])
l = map(int, i.split(";"))
ltab[j].extend(l)
j += 1
self.tab["Array"] = numpy.array(ltab)
self.eZ.delete(0, Tkinter.END)
self.eZ.insert(Tkinter.END, str(self.tab["Z"]))
self.cbbLissage.pick(0)
self.Refresh()

Cette fonction permet de charger une carte. La structure d'une carte est simple:
La première ligne contient les valeurs X, Y, Z précédemment définies, séparées par un point-virgule.
Ensuite sont entrées les valeurs en Z des données topographiques de la carte selon un quadrillage métrique, chaque valeur étant séparée par un point-virgule avec passage à la ligne selon le quadrillage. Ne pas oublier que Z est considérée comme une profondeur. À chaque chargement d'une nouvelle carte, on repasse en coordonnées normées.

III-3-3 - La fonction Affiche

def Affiche(self, tab):
rc, cc = tab.shape
pos = []
divi = 2 ** (int(self.cbbLissage.subwidget('listbox').curselection()[0])-1)
eZ = 1 / float(self.eZ.get())
eX = float(self.tab["X"])
eY = float(self.tab["Y"])

for i in range(rc-1):
for j in range(cc-1):
pos.append((i*eX/divi, j*eY/divi, tab[i, j] / eZ))
pos.append(((i*eX/divi), (j+1)*eY/divi, tab[i, j+1] / eZ))
pos.append(((i+0.5)*eX/divi, (j+0.5)*eY/divi, (tab[i, j]+ tab[i, j+1]+tab[i+1, j] + tab[i+1, j+1]) / (4*eZ)))

pos.append((i*eX/divi, (j+1)*eY/divi,     tab[i, j+1] / eZ))
pos.append(((i+1)*eX/divi, (j+1)*eY/divi, tab[i+1, j+1] / eZ))
pos.append(((i+0.5)*eX/divi, (j+0.5)*eY/divi, (tab[i, j]+ tab[i, j+1]+tab[i+1, j] + tab[i+1, j+1]) / (4*eZ)))

pos.append(((i+1)*eX/divi, j*eY/divi, tab[i+1, j] / eZ))
pos.append((i*eX/divi, j*eY/divi, tab[i, j] / eZ))
pos.append(((i+0.5)*eX/divi, (j+0.5)*eY/divi, (tab[i, j]+ tab[i, j+1]+tab[i+1, j] + tab[i+1, j+1]) / (4*eZ)))

pos.append(((i+1)*eX/divi, (j+1)*eY/divi, tab[i+1, j+1] / eZ))
pos.append(((i+1)*eX/divi, j*eY/divi, tab[i+1, j] / eZ))
pos.append(((i+0.5)*eX/divi, (j+0.5)*eY/divi, (tab[i, j]+ tab[i, j+1]+tab[i+1, j] + tab[i+1, j+1]) / (4*eZ)))

try: self.model.frame.visible = 0
except Exception, err: print err
self.f = frame()
self.model = faces( pos = pos, frame = self.f)

mini = tab.min() / eZ
maxi = tab.max() / eZ
print mini, maxi
mini, maxi = min(mini, maxi), max(mini, maxi)
for i in range(len(pos)):
self.model.color[i] = (1 - (self.model.pos[i,2] - mini) / (maxi-mini),
1 - (self.model.pos[i,2] - mini) / (maxi-mini),
1 - (self.model.pos[i,2] - mini) / (maxi-mini))

scene.center = (cc*eX/(2*divi), rc*eY/(2*divi), (maxi+mini)/2)
scene.forward = (0, 1, 1)
scene.up = (0, 0, -1)
scene.lights = []
scene.ambient = 1

La partie concernant Visual Python est entièrement contenue dans cette fonction. L'instruction int(self.cbbLissage.subwidget('listbox').curselection()[0]) permet de récupérer l'index sélectionné dans la ComboBox.
divi est valeur du nombre de divisions subit par le tableau et sera utile pour conserver l'échelle entre les différentes coordonnées. Il en sera de même pour les paramètres eX, eY et eZ.

On va ensuite remplir l'attribut model. Cet objet est une face de VPython et se construit de la manière suivante: il prend comme paramètre le paramètre pos, une liste de sommets, la taille de la liste étant donc un multiple de 3, chaque groupe de trois représentant une surface triangulaire. On peut également lui indiquer un paramètre frame qui fera que tous les sommets seront inscrits dans une même structure, ce qui sera utile quand on voudra remettre à nu la scène avant de réafficher une nouvelle carte ou la même, mais en ayant changer les paramètres d'affichage. La normale de chaque surface triangulaire (ce qui permet de définir le côté de la face qui est éclairé) est calculée avec l'ordre des sommets. Il donc important de respecter toujours le même sens dans l'écriture des sommets.

Au niveau de la scène, on choisit une lumière ambiante non directionnelle pure afin d'avoir la même vue, quel que soit l'endroit d'où l'on voit la scène. Comme on travaille par défaut avec des profondeurs (plutôt que des altitudes), le vecteur up pointant vers le haut de l'espace est (0, 0, -1). Si vous travaillez avec des altitudes, il suffira de modifier le fichier de la carte pour que la valeur de Z soit négative.

III-3-4 - La Fonction Refresh

def Refresh(self):
tabdivi = copy.deepcopy(self.tab["Array"])
for i in range(int(self.cbbLissage.subwidget('listbox').curselection()[0])):
tabdivi = divisetab(tabdivi)
self.Affiche(tabdivi)

Cette fonction lance l'affichage ou le réaffichage de la scène si les paramètres ont été modifiés de la carte. La première chose est d'effectuer la division de la matrice puis de lancer la fonction Affiche sur cette dernière.

IV - Téléchargement

Version

Date

Taille

Mode FTP

Mode HTTP de secours

1.0.0.0.0

2009.06.08

3.59 Ko

pytopo3d.zip

pytopo3d.zip