Archives du mot-clé géométries

Transformations de coordonnées avec GDAL

Cette petite note technique traite d’une des activités les plus fréquentes quand on travaille avec des données géographiques: la transformation de coordonnées entre différentes projections. L’incontournable PostGIS est pourvu de la fonction st_transform(the_geom,srid), qui passe n’importe quelle géométrie dans le système de coordonnées voulu. Or, il n’est pas toujours possible pour des contraintes techniques ou autres soucis de performances de solliciter PostGIS.

La bibliothèque GDAL est un outil très adapté dans les cas où vous voulez manipuler des données géographiques en C++ ou Python. Celle ci a la mérite d’être bien développée, même s’il n’y a pas vraiment une grande communauté et beaucoup de tutoriel sur ses fonctionnalités. Pour preuve, le site de GDAL n’est autre qu’une documentation générée par Doxygen.

Il n’est pourtant pas très compliqué, en C++ ou Python, de transformer des coordonnées avec GDAL, supposant qu’on utilise le langage C++:

  • Installer la bibliothèque  GDAL, sur des systèmes de type Debian ou Ubuntu:
sudo  apt-get install libgdal1-dev
  • Il est nécessaire de faire les inclusions suivantes:
#include <gdal/ogr_geometry.h>
#include <gdal/ogr_spatialref.h>
#include <gdal/ogr_srs_api.h>s
  • Pour effectuer des transformations d’un système à l’autre, on déclare d’abord deux objets représentant des systèmes de coordonnées, avec la classe OGRSpatialReference, puis l’objet destiné à transformer des coordonnées (classe OGRCoordinateTransformation). Typiquement, on utilise le code EPSG, pour représenter un système de coordonnées. Ici, on instancie la projection 4326 (latitude et longitude) et la projection 25572, qu’on appelle aussi Lambert II, et qui a l’avantage en France de représenter fidèlement les distances en mètres. On cherche alors à transformer les coordonnées de Lambert II vers longitude/latitude
OGRSpatialReference refLonLat;
refLonLat.importFromEPSG(4326);
OGRSpatialReference refLambert;
refLambert.importFromEPSG(27572);
OGRCoordinateTransformation * lambertToLonLat = OGRCreateCoordinateTransformation(&refLambert,&refLonLat);
  • Lorsqu’on veut convertir un vecteur de deux coordonnées, avec les variables x et y, on utilise l’objet de transformation avec la méthode Transform
double x = ...;
double y = ...;
lambertToLonLat->Transform(1,&x,&y);

A noter, les coordonnées sont converties sur place, c’est à dire que les variables (ou tables) qui lui sont données sont modifiées. Si le premier argument vaut 1, on peut lui donner de simples valeurs de type double, ou bien des tableaux de double dans le cas contraire…

Pour convertir des coordonnées en Python, un code Python équivalent serait:

try:
from osgeo import ogr
except:
import ogr
try:
from osgeo import osr
except:
import osr
refLambert = osr.SpatialReference()
refLambert.ImportFromEPSG(27572)
refLonLat = osr.SpatialReference()
refLonLat.ImportFromEPSG(4326)
coordTrans = osr.CoordinateTransformation(refLambert,refLonLat)
point = ogr.Geometry(ogr.wkbPoint)
x = ...
y = ...
point.AddPoint(x,y)
point.Transform(coordTrans)