On possède un Modèle Numérique de Terrain sous forme de raster dans une base de données PostgreSQL (PostGIS 2), et on veut le visualiser en 3D sur le Web.

On a utilisé Python et Flask pour exécuter la requête SQL qui va extraire les données du MNT, et transformer le résultat en JSON. Babylon.JS va consommer ce JSON et, via la technologie WebGL, afficher un rendu 3D du terrain dans le navigateur.

Une démonstration du résultat est disponible ici.

Dans cet article, l'accent est mis principalement sur la requête SQL qui va nous permettre d'extraire les données.

Récupération des données

• Requête SQL et PostGIS

Remarque : Cette partie suppose que l'on connaisse un minimum le langage SQL.

Pour alléger les calculs, on va extraire seulement les informations du raster qui intersectent une certaine zone. Il faut donc choisir une zone, de préférence rectangulaire, mais surtout comprise dans le MNT pour espérer avoir un résultat !

Ensuite, la fonction generate_series(), associée à la magie de SQL, va nous permettre de créer une sorte de ''grille'' de points 2D à partir du rectangle (i.e. un échantillonnage, comme pour le drapé de lignes).

Ainsi on a toutes les coordonnées en x et y, il ne nous reste plus qu'à récupérer les valeurs d'altitude de ces points grâce à la fonction ST_Value() de PostGIS.

/!\ On n'oublie pas de faire ceci seulement pour les points du raster qui intersectent les points créés !  (d'où le WHERE ST_Intersects… ). Sinon, la requête interroge le raster pour chaque cellule inutilement…

WITH    
    -- On part d'une zone rectangulaire arbitraire
    buffer AS (
           SELECT ST_MakePolygon('SRID=32632;LINESTRING(
                                  338166 4904645,339302 4904645,339302 4905408,
                                  338166 4905408,338166 4904645)'::geometry) AS geom),


    -- On récupère les extrémités du rectangle
    buffer_extent AS (
           SELECT ST_Xmin(geom) AS xmin, ST_Xmax(geom) AS xmax,
                  ST_Ymin(geom) AS ymin, ST_Ymax(geom) AS ymax
           FROM buffer),


    -- On subdivise les lignes et colonnes du rectangle : par pas de 50m
    columns AS (
        SELECT generate_series(b.xmin::int, b.xmax::int, 50) AS x 
        FROM buffer_extent AS b),
    lines AS (
        SELECT generate_series(b.ymin::int, b.ymax::int, 50) AS y 
        FROM buffer_extent AS b),


    -- On instancie les points de la "grille" à partir des lignes et colonnes
    points2d AS (
        SELECT ST_SetSrid(ST_MakePoint(x, y), 32632) AS geom 
        FROM lines, columns),


    -- On drape la grille de points sur le MNT pour avoir x, y, z
    drape AS (
        SELECT  ST_X(p.geom) AS x, 
                ST_Y(p.geom) AS y, 
                ST_Value(mnt.rast, p.geom) AS z
        FROM mnt, points2d p
        WHERE ST_Intersects(mnt.rast, p.geom)), -- limite aux points de la grille


    average AS (
        SELECT  AVG(z) AS avgz
        FROM drape)


-- On récupère :
SELECT  x, y, z, -- les coordonnées de tous les points
    (xmax+xmin)/2 AS avgx, (ymax+ymin)/2 AS avgy, avgz, -- les 3 moyennes pour le centre
    xmin, xmax, ymin, ymax -- les extrémités 
FROM drape, average, buffer_extent;

• Construction des données

A partir de ces données, on construit 3 choses nécessaires ensuite pour la modélisation en 3D :

• les points (vertices),

• le centre du terrain (center) : utile pour la position de la caméra,

• La résolution, qui comprend le nombre de points sur x (w_sub) et le nombre de points sur y (h_sub).

/!\ Le système de coordonnées de PostGIS est un peu différent de celui de Babylon.js  dans le sens où le premier va associer l'axe z à la hauteur (altitude) alors que pour l'autre, l'altitude sera matérialisée par l'axe des y.

Pour l'exportation des données récupérées, j'ai choisi d'utiliser le format JSON (on notera que le format choisi n'est pour l'instant qu'arbitraire, le but était surtout ici d'afficher les données pour se donner une idée de rendu).

• Limiter la grandeur des données

J'ai remarqué que l'on avait des nombres assez grands (de l'ordre de 10000 ou 100000) mais qu'ils restaient quand même très proches. Donc pour éviter de travailler avec des nombres si grands, j'ai procédé à 2 manipulations.

   - une soustraction des positions minimums en x et en y,

   - une division par 10 de toutes les valeurs x, y et z.

• Flask.jsonify

La transformation des données obtenues de le requête sous forme de JSON se fait grâce à la méthode jsonify de Flask.
On obtient ainsi ce genre de JSON, que l'on va pouvoir exploiter pour initialiser la scene Babylon.js.

{
 "center": [
 56.8, 
 38.15, 
 278.506438678244
 ],

 "resolution": {
 "h_sub": 15.26, 
 "w_sub": 22.72
 },

 "vertices": [
 0.0, 
 273.330004882812, 
 0.0, 
 5.0, 
 269.519995117188, 
 0.0, 
 10.0, 
 269.144995117188,
 ...
 110.0, 
 286.889990234375, 
 75.0
 ]
}

Exploitation des données

L'affichage est réalisé grâce à Babylon.js: une bibliothèque Javascript très puissante qui utilise la technologie WebGL pour faire relativement simplement de la 3D dans le navigateur !

La fonction CreateGround() de la classe Mesh semblait toute appropriée pour afficher ce genre de terrain : On crée un rectangle et on lui donne le nombre de subdivisions.

Ici, nos longueurs et largeurs ne possèdent pas le même nombre de points, donc pas le même nombre de subdivisions ! Il a donc fallu modifier un peu cette fonction pour passer de : 

CreateGround(name, width, height, subdivisions, scene, updatable)

à :

MyCreateGround(name, width, height, width_subdivisions, height_subdivisions, scene, updatable).

Le code des deux fonctions est similaire, mais la seconde fonctionne avec un rectangle qui n'a pas la même résolution en largeur et longueur. Ici un lien vers la nouvelle fonction.

Voici des captures d'écrans de représentations en fils de fer (wireframe) des résultats.

Pas de 100m

Pas de 50m

Pas de 20m