Dissémination de données géoréférencées – qgis-server et OpenLayers – Partie 2/2

4. À la pointe de la technologie : OpenLayers 3

Dans la course aux fonctionnalités, nous ne cessons de casser ce qui marche pour mettre à jour les nouvelles fonctionnalités. Pourquoi donc casser ce bel exemple OpenLayers2 pour le remplacer par OpenLayers3 ? La déformation excessive des cartes aux latitudes élevées est un défaut bien connu de la projection de Mercator utilisée par défaut par tous les environnements cartographiques de la toile (projection de code 900913 – l’écriture leet de Google – puisque la norme choisie n’a pas été validée par l’instance qu’est l’EPSG). Étant conscients des déficiences de la projection de Mercator, nous voudrions pouvoir sélectionner un mode de projection localement tangent à la région considérée – à savoir Universal Transverse Mercator (UTM). Et nous en arrivons à l’excuse pour passer de OpenLayers 2 à 3 : la version actuelle d’OpenLayers ne permet pas la projection au vol des images, et impose donc une sortie déformée, particulièrement gênante lorsque nous nous éloignons de l’équateur.

Nous reprenons donc toutes nos investigations, mais cette fois dans le contexte d’OpenLayers 3, avec la nécessité pour chaque nouvelle couche de préciser le mode de projection d’entrée (toutes nos données sont stockées en coordonnées sphériques, donc WGS84 de code EPSG:4326. OpenLayers ne connaît pas tous les modes de projection de sortie : nous devons les définir selon nos besoins (voir figure 5).

Fig. 5: Gauche : le passage à une projection locale WGS84/UTM33N résout la déformation observée auparavant dans OpenLayers 2 et son Google Mercator. Cependant, cette projection n’a pas prétention de s’étendre à tout le globe, et réduire le grossissement met en évidence la déformation sur les régions adjacentes, ici le nord Canada et la Sibérie (droite). Cet exemple est disponible sur jmfriedt.sequanux.org/reproj.html.

 

Heureusement, proj4 se charge de cette opération pour nous, et http://spatialreference.org/ref/epsg/wgs-84-utm-zone-33n/ nous informe de la façon de définir WGS84/UTM33N (la bande qui couvre la Norvège – nous aurions choisi UTM31N pour la France). Ainsi, le mode de projection qui va nous intéresser se définit par :

Et nous pouvons faire appel aux outils de conversion de la forme :

Ici la projection de la carte respecte la norme EPSG (WGS84/UTM33N s’appelle EPSG:32633) que nous avons défini par proj4 auparavant, tandis que les données en entrée (ici issues de la fonction WM(T)S du serveur) sont au format WGS84 en coordonnées sphériques (aussi nommée EPSG:4326).

5. Ajout automatique de toutes les couches WMS d’un projet qgis-server dans OpenLayers2

Une dernière question que nous pouvons nous poser, dans le contexte d’un travail collaboratif où chaque utilisateur est susceptible d’ajouter ses propres couches, tient en la capacité à automatiser la génération d’une page Web contenant toutes les couches fournies par un serveur QGis. Cela signifie donc récupérer la description XML des capacités du serveur, découper cette description pour extraire le nom des couches disponibles, et finalement générer automatique le script JavaScript affichant toutes ces couches. De telles fonctionnalités sont fournies (paquet php-xml dans la distribution Debian) par simplexml. Ainsi, un premier script qui se contente de lister toutes les couches disponibles en recherchant d’abord les fonctionnalités WMS du serveur, puis en étudiant [5] la liste des propriétés (FeatureTypeList) et finalement les propriétés de chacun de ces objets (FeatureType) pour en extraire le nom (FeatureTypeList->FeatureType->Name) ressemble à :

Ayant compris ce mécanisme, il devient simple de générer les champs décrivant chaque couche WMS affichée, en lui attribuant le même nom de couche que celui utilisé dans QGis. Le résultat est :

Ce qui donne la figure 6. Il serait d’une part très fastidieux d’ajouter à la main toutes les couches de ce projet, mais surtout dans cet exemple l’affichage s’adapte automatiquement aux nouvelles couches ajoutées par les divers contributeurs au projet.

Fig. 6: Ajout automatique de toutes les couches fournies par un serveur WMS sur dans un contexte d’OpenLayers2. Le script d’exemple /opt/qgis-server/plugins/HelloServer qui affiche le logo de QGis a été volontairement laissé en place pour distinguer le serveur en production (qgis.sequanux.org) du serveur de tests (127.0.0.1).

 

6. Ajout d’images orthorectifiées produites par microdrone

Nous avions présenté le flux de traitement permettant de générer une série d’images orthorectifiées et de modèles numériques d’élévation associés avec le souci de comparer ces jeux de données entre eux, mais sans prétention de les positionner dans un référentiel absolu partagé par convention entre plusieurs utilisateurs. Il s’avère [6] que la solution de simplement translater un jeu de données par rapport à l’autre ne suffit plus lorsque la zone considérée s’étend sur plusieurs hectares. Nous avons observé que dans ce cas, le positionnement par GPS (mono-fréquence, telle que fournie sur microdrone DJI Phantom3 Professional) n’est exact qu’à une dizaine de mètres, alors que nous visons un positionnement au décimètre près pour des comparaisons de vols successifs. Il s’avère que vouloir corriger les défauts d’échelle et de position du modèle numérique d’élévation en lui appliquant les corrections permettant de superposer les images orthorectifiées est une approche peu judicieuse, et qu’il vaut mieux insérer dans le flot de traitement l’étape d’exploitation des points de contrôle au sol pour corriger les défauts du modèle de caméra et de leur position (http://forum-micmac.forumprod.com/campari-residuals-differ-from-the-point-cloud-t881.html). Pour notre part, les points de contrôle au sol (GCP) sont acquis a posteriori dans Google Earth – les coordonnées sphériques fournies avec 6 décimales sont précises à 11 cm à l’équateur : 7 points aisément reconnaissables dans les jeux de données sont d’une part identifiée dans Google Earth (bandes de parking, barrière, coin à la base de bâtiments) avec des coordonnées converties de coordonnées sphériques (WGS84) vers un référentiel projeté localement tangent à la Terre (WGS84/UTM31N), et d’autre part leur position (en pixel) identifiée sur les images acquises par drone. Ce couple de fichiers (coordonnées dans l’espace v.s coordonnées sur les photographies) alimente l’outil Campari de Micmac. Après traitement dans ces conditions, nous observons une exactitude de position par rapport aux couches vectorielles de Openstreetmaps mais surtout par rapport à une image aérienne de l’IGN (BD ORTHO, convertie de Lambert93 à WGS84/UTM31N puisque nous avons vu que OpenLayers 2 ne sait pas le faire à la volée) meilleure qu’une vingtaine de centimètres – erreur attribuable à notre manque d’assiduité à pointer les points de contrôle. Ce jeu de données est disponible sur http://qgis.sequanux.org/cgi-bin/project/femto/qgis_mapserv.fcgi.

Fig. 7: Résultat de l’insertion de trois images orthorectifiées, en comparaison de données de référence que sont les images aériennes de la BD ORTHO de l’IGN et les données vectorielles de Opensteetmaps. En haut à gauche : deux images orthorectifiées acquises à 30 minutes d’intervalle positionnées au moyen de 7 GCPs. En haut à droite : le même jeu de données, positionné uniquement par les positions GPS du drone au moment de la prise de vue. Le trait rouge, indiquant la différence de position, est long de 9,8 m. En bas à gauche : superposition d’une image produite par drone (moitié gauche) avec une image IGN (moitié droite). En bas à droite : comparaison de deux orthophotos produites à un mois d’intervalle.

Conclusion

Nous avons proposé un environnement de travail permettant de disséminer des informations géoréférencées soit vers des utilisateurs de logiciels dédiés de gestion de systèmes d’informations géographiques (SIG) ou vers les API Web associées (OpenLayers). Est-ce que ce mode de dissémination est réellement utilisé en pratique ? L’institut polaire norvégien (NPI) propose ses informations dans ce format tel que décrit sur http://geodata.npolar.no/#basemap-services. Ainsi, sous QGis, configurer l’onglet WM(T)S vers http://geodata.npolar.no/arcgis/rest/services/Basisdata/NP_Ortofoto_Svalbard_WMTS_25833/MapServer/WMTS? donne accès aux informations bitmap que sont les photographies aériennes de très haute résolution (jusquà 16 cm/pixel !) tandis que l’url http://geodata.npolar.no/arcgis/services/CryoClim/glaciers/MapServer/WmsServer? dans l’onglet WFS donne accès aux informations vectorielles que sont les limites des glaciers extraites des images satellitaires SPOT.

Le lecteur est encouragé à reproduire ces expériences sur ses propres zones d’intérêt. Compte tenu de la résolution médiocre de Landsat, le résultat est peu convaincant pour la vallée de Chamonix, mais les glaciers groenlandais (Illulisat près de la baie de Disco) ou la banquise antarctique sont parfaitement appropriés à cette démonstration, voir la mer d’Aral sur http://www.nasa.gov/mission_pages/landsat/news/40th-top10-aralsea.html (un peu pénible à assembler car de superficie supérieure à ce que le site de l’USGS permet d’analyser) ou la mer morte sur http://earthobservatory.nasa.gov/IOTD/view.php?id=77592 et pour laquelle l’extension des bassins d’exploitation du sel est particulièrement visible. Le lecteur saura-t-il calculer la surface ainsi affectée ? Indice : tracer un polygone avec l’outil New Shapefile puis dans la table des attributs, créer un nouveau champ pour chaque polygone avec la variable $area.

Une autre voie d’étude qui semble intéressante dans le cadre de la diffusion sur le Web de données géoréférencées est qgis2web, greffon de QGis actuellement inexploitable sous Debian testing/sid compte tenu d’une dépendance cassée. Finalement, étant donné que toutes les informations pertinentes à une utilisation en surface sont disponibles sur OpenStreetMaps, ne serait-il pas temps de commencer à considérer la cartographie des services sous terrains (eau, gaz, électricité), qui même s’ils sont propriété de leurs exploitants respectifs, sont devenus des services nécessaires au quotidien d’un habitant d’Europe occidentale actuel : chaque kilomètre de route contient plusieurs dizaines de kilomètres de services sous terrains [7] qui ne demandent qu’à être cartographiés à l’occasion des ouvertures de routes.

Jean-Michel FRIEDT
[Institut FEMTO-ST, dpt. Temps-fréquence, Besançon]

Émile Carry
[Institut FEMTO-ST, dpt. Temps-fréquence, Besançon]

Références et Notes

[1] PIERROT-DESEILLIGNY M. et FRIEDT J.-M, « La photogrammétrie pour tous : MicMac pour la reconstruction 3D de scènes géoréférencées à partir de photographies numériques », tutorial à FOSS4G-fr 2016, décrit sur jmfriedt.free.fr/foss4g_2016
[2] FRIEDT J.-M. Friedt, TOLLE F., et BERNARD É., « Utilisation de Micmac pour la génération de modèle numérique d’élévation par traitement d’images acquises par microdrone », GNU/Linux Magazine France 191, pp.48–57 (Mars 2016)
[3] ERLE S., GIBSON R., et WALSCH J., « Building the geospatial web », dans « Mapping Hacks – Tips & Tools for Electronic Cartography », O’Reilly (2005)

[4] GIBSON R. et ERLE S., « Google Maps Hacks – Tips & Tools for Geographic Searching and Remixing », O’Reilly (2006)
[5] Les divers champs se déduisent de l’analyse de la sortie de wget -O – « http://127.0.0.1/cgi-bin/mon_projet/qgis_mapserv.fcgi?SERVICE=WFS&VERSION=1.0.0&REQUEST=GetCapabilities »
[6] LISEIN J., PINEUX N., PIERROT-DESEILLIGNY, DEGRÉ A., et LEJEUNE P.,
« Détection de l’érosion dans un bassin versant agricole par comparaison d’images multidates acquises par drone », Colloque Drones et moyens légers aéroportés d’observation, 26/06/2014 (Montpellier), disponible sur http://orbi.ulg.ac.be/handle/2268/171616
[7]
Chaque kilomètre de route à Hong-Kong recouvre 47 km de services enterrés,
http://www.uti.hk/media/attachments/2nd_ICUMAS_full_paper.pdf ou http://www.scmp.com/news/hong-kong/article/1647088/small-army-utility-specialists-keeps-hongkongers-safe-pipes-cables, ou en d’autres termes il existe 47 câbles, tuyaux, fibres optiques et autres modes de transmission de fluides et d’informations sous chaque route.

Remerciements

J.-P Culas (CM-Drones, Besançon) a effectué les vols au-dessus du parking de FEMTO-ST. Les membres du forum Micmac (http://forum-micmac.forumprod.com/) ont une fois de plus partagé leurs connaissances pour corriger les défaillances dans les traitements que nous proposions initialement.

Retrouvez cet article (et bien d’autres) dans GNU/Linux Magazine n°200, disponible sur la boutique et sur la plateforme de lecture en ligne Connect !