Skip to Content

Nous adressons toutes nos pensées à la famille de notre ami Jérôme !

http://www.forumsig.org/showthread.php/43488-Disparition-de-Phoenix

OGR: que la force soit avec les formats virtuels

Niveau Intermédiaire
Logiciels utilisés OGR
Plateforme Windows | Mac | Linux | FreeBSD

Beaucoup de Sigistes connaissent ou du moins ont entendu parler du format virtuel de GDAL (GDAL Virtual Format Tutorial ) pour les fichiers de type raster et de ses multiples possibilités (The power of GDAL virtual formats ou Image Mosaicking with GDAL), mais il me semble que peu d'entre eux connaissent ce même format pour les couches vectorielles.

Principes:

Les caractéristiques sont exposées dans OGR Virtual Format. Tout comme celui pour les rasters, il s'agit d'un simple fichier au format XML. Ce fichier VRT (.vrt) se traite comme n'importe quel  fichier vectoriel (fichier shapefile, etc.):

<OGRVRTDataSource>
    <OGRVRTLayer name="points">
        <SrcDataSource>/chemin/points.csv</SrcDataSource>
        <GeometryType>wkbPoint</GeometryType>
        <LayerSRS>EPSG:31370</LayerSRS>
        <GeometryField encoding="PointFromColumns" x="x" y="y"/>
    </OGRVRTLayer>
</OGRVRTDataSource>

Dans le cas présent, il s'agit d'un simple fichier CSV qui pourrait être ouvert avec un logiciel SIG comme QGIS plus simplement. Mais, il est possible d'aller plus loin comme:

  • la création directe d'un fichier 3D, non reconnu come tel par QGIS, mais bien par des logiciels comme GRASS GIS (toujours à partir d'un simple fichier CSV):
<OGRVRTDataSource>
    <OGRVRTLayer name="points">
        <SrcDataSource>/chemin/points3D.csv</SrcDataSource>
        <GeometryType>wkbPoint25D</GeometryType>
        <LayerSRS>EPSG:31370</LayerSRS>
        <GeometryField encoding="PointFromColumns" x="x" y="y" z="z"/>
    </OGRVRTLayer>
</OGRVRTDataSource>
  • le changement de système de projection de la couche (depuis la version 1.10 de GDAL):
<OGRVRTDataSource>
    <OGRVRTWarpedLayer>
        <OGRVRTLayer name="test">
            <SrcDataSource>/chemin/test.shp</SrcDataSource>
        </OGRVRTLayer>
        <TargetSRS>EPSG:4326</TargetSRS>
    </OGRVRTWarpedLayer>
</OGRVRTDataSource>
  • la redéfinition des noms et des formats des colonnes (depuis la version 1.7 de GDAL);
  • le découpage d'une couche en fonction d'un polygone (depuis la version 1.7 de GDAL, voir plus loin);
  • la possibilité d'utiliser des indexes spatiaux;
  • la possibilité de charger plusieurs couches à la fois (comme avec GDAL);
  • la possibilité de découper une couche en fonction des attributs;
  • ....

De plus, tous les types de support peuvent être traités: fichiers textes, fichiers shapefiles, bases de données, base de données spatiale, avec ODBC, etc.

avec une base SQLite non spatiale:

<OGRVRTDataSource>
    <OGRVRTLayer name="Formations">
        <SrcDataSource>/chemin/test.sqlite</SrcDataSource> 
        <SrcLayer>Form</SrcLayer> 
        <GeometryType>wkbPoint</GeometryType>
            <LayerSRS>EPSG:31370</LayerSRS>
        <GeometryField encoding="PointFromColumns" x="lng" y="lat"/> 
    </OGRVRTLayer>
</OGRVRTDataSource>

avec un fichier dbf:

<OGRVRTDataSource>
    <OGRVRTLayer name="Formations">
        <SrcDataSource>/chemin/test.dbf</SrcDataSource> 
        <GeometryType>wkbPoint</GeometryType>
        <LayerSRS>EPSG:31370</LayerSRS>
        <GeometryField encoding="PointFromColumns" x="lng" y="lat"/> 
    </OGRVRTLayer>
</OGRVRTDataSource>

avec un fichier comprenant une géométrie au format WKT (voir GDAL and OGR Virtual Format, reading *.csv file and creating shapefile of lines):

<OGRVRTDataSource>
       <OGRVRTLayer name="Linestest">
           <SrcDataSource SEPARATOR="SEMICOLON">/chemin/wkt.csv</SrcDataSource>
           <SrcLayer>wkt</SrcLayer>
           <GeometryType>wkbLineString</GeometryType>
           <GeometryField encoding="WKT" field="lignes"/>
       </OGRVRTLayer>
</OGRVRTDataSource>

avec des données provenant de ArcSDE (voir Connecting to ArcSDE databases)

<OGRVRTDataSource>
       <OGRVRTLayer name="layer name in qgis"> 
           <SrcDataSource>SDE:servername,5151,,user,pwd,SCHEMA.TABLENAME,SDE.DEFAULT</SrcDataSource> 
           <SrcLayer>SCHEMA.TABLENAME</SrcLayer> 
           <GeometryType>wkbMultiPolygon</GeometryType> 
           <LayerSRS>+proj=tmerc +lat_0=39.66666666666666 +lon_0=-8.131906111111112 +k=1 +x_0=200180.598 +y_0=299913.01 +ellps=intl +units=m +towgs84=-223.237,110.193,36.649 +no_defs</LayerSRS>
       </OGRVRTLayer> 
</OGRVRTDataSource>

Je ne continuerai pas avec les possibilités basiques du format pour me focaliser sur la partie la plus intéressante, les traitements avec des requêtes SQL.

Requêtes SQL

Depuis la version 1.10 de GDAL, les commandes sont de 2 types:

OGR SQL

En introduction, commençons par un cas simple qui pourrait être réglé par le logiciel SIG mais que nous traiterons ici avec un fichier virtuel. Il s'agit d'une carte géologique dont je veux extraire le tracé d'une formation (BUR ici) particulière. J'utiliserai QGIS pour visualiser les résultats.

Le fichier vrt est alors le suivant:

<OGRVRTDataSource>
    <OGRVRTLayer name="BUR Burnot">
        <SrcDataSource>/chemin/varform.shp</SrcDataSource>
        <SrcSQL>SELECT * from varform WHERE "VAR_FORM" = "BUR"</SrcSQL>
    </OGRVRTLayer>
</OGRVRTDataSource>

Et le résultat est:

En pratique, tous les exemples illustrés dans SQL dans OGR peuvent être utilisés.

  • découper une zone d'un shapefile en fonction d'un polygone;
<OGRVRTDataSource>
    <OGRVRTLayer name="clip">
        <SrcDataSource>/chemin/varform.shp</SrcDataSource>
        <SrcRegion clip="true">POLYGON((170504.41570972 116282.08448781,170576.34282934 113909.22348016,175203.46747074 114069.66549297,175169.19847681 116087.90004614,170504.41570972 116282.08448781))</SrcRegion>
    </OGRVRTLayer>
</OGRVRTDataSource>

Résultat:

  •  il y a aussi moyen d'effectuer des opérations géométriques simples comme l'union entre 2 couches:

Le fichier est alors le suivant (union des BUR entre les 2 couches):

<OGRVRTDataSource>
	<OGRVRTUnionLayer name="union_couche_BUR">
		<OGRVRTLayer name="couche1">
			<SrcDataSource>/chemin/biesm.shp</SrcDataSource>
			<SrcSQL>SELECT * from biesm WHERE "VAR_FORM" = "BUR"</SrcSQL>
		</OGRVRTLayer>
		<OGRVRTLayer name="couche2">
			<SrcDataSource>/chemin/goz.shp</SrcDataSource>
			<SrcSQL>SELECT * from goz WHERE "VAR_FORM" = "BUR"</SrcSQL>
		</OGRVRTLayer>
	</OGRVRTUnionLayer>
</OGRVRTDataSource>

Résultat:

  • quelques autres fonctions sont aussi disponibles comme la sélection en fonction des surfaces/aires: 
<OGRVRTDataSource>
    <OGRVRTLayer name="surfaces">
        <SrcDataSource>/chemin/varform.shp</SrcDataSource>
        <SrcSQL>SELECT * FROM varform  WHERE OGR_GEOM_AREA > 10000000</SrcSQL>
    </OGRVRTLayer>
</OGRVRTDataSource>

Résultat:

SQLite SQL dialect

Ce nouveau format, apparu dans la version 1.10 de GDAL, permet d'utiliser le SQL utilisé par SQLite (SQL As Understood By SQLite), mais aussi toutes les fonctions spatiales de SpatiaLite (SpatiaLite 4.0.0: SQL functions reference list ).

A titre de démonstration, nous allons appliquer plusieurs traitements à une même couche (simplification, zone tampon, enveloppe convexe (Convex hull):

<OGRVRTDataSource>
    <OGRVRTLayer name="BUR">
        <SrcDataSource>/chemin/bur.shp</SrcDataSource>
    </OGRVRTLayer>
    <OGRVRTLayer name="BUR_simplifié">
        <SrcDataSource>/chemin/bur.shp</SrcDataSource>
        <SrcSQL dialect="sqlite">SELECT Simplify(geometry,10) from bur</SrcSQL>
    </OGRVRTLayer>
    <OGRVRTLayer name="BUR_tanpon">
        <SrcDataSource>/chemin/bur.shp</SrcDataSource>
        <SrcSQL dialect="sqlite">SELECT Buffer(geometry,200) from bur</SrcSQL>
    </OGRVRTLayer>
    <OGRVRTLayer name="BUR_convexhull">
        <SrcDataSource>/chemin/bur.shp</SrcDataSource>
        <SrcSQL dialect="sqlite">SELECT ConvexHull(geometry) from bur</SrcSQL>
    </OGRVRTLayer>
</OGRVRTDataSource>

L'ouverture de cette couche vrt dans QGIS conduit au dialogue suivant:

Et en sélectionnant toutes les couches, le résultat est:

Attention, tout n'est pas possible. Ainsi si je veux extraire les limites de polygones (ExteriorRing()), cela ne marchera pas s'il y a plusieurs polygones dans une couche:

<OGRVRTDataSource>
    <OGRVRTLayer name="multipolygones">
        <SrcDataSource>/chemin/bur.shp</SrcDataSource>
        <SrcSQL dialect="sqlite">SELECT ExteriorRing(Geometry) from bur</SrcSQL>
    </OGRVRTLayer>
    <OGRVRTLayer name="polygone">
        <SrcDataSource>/chemin/bur1.shp</SrcDataSource>
        <SrcSQL dialect="sqlite">SELECT ExteriorRing(Geometry) from bur1</SrcSQL>
    </OGRVRTLayer>
</OGRVRTDataSource>

 

Seul le deuxième résultat donnera une couche valide. L'autre ne fournira qu'une table non spatiale:

Ceci veut dire que l'on peut parfaitement se servir de ces formats virtuels pour créer des tables non spatiales à sa guise (pour les jointures, etc.).

À quoi cela peut-il bien servir ?

Tout ça, c'est intéressant, me direz-vous, mais à quoi cela peut-il servir en pratique puisque tout peut être fait par un logiciel SIG ?

Cela permet de « bidouiller » ce que l'on veut à partir de couches vectorielles avant traitement dans un logiciel SIG ou pour les traitements sans logiciel SIG:

Avec GDAL/OGR seul:

  • transformation directe en un autre format:
ogr2ogr -f KML -t_srs EPSG:4326 test.kml test.vrt
gdal_grid -zfield "z" -a invdist:power=2.0:smoothing=1.0 -txe 85000 89000 -tye 894000 890000 *-outsize 400 400 -of GTiff -ot Float64 -l dem test.vrt test.tiff

 Avec Python:

Avec d'autres programmes:

Conclusions

Dans la série « Comment traiter les données SIG sans logiciel SIG », ce format remporte certainement la palme.

Tous les traitements ont été réalisés sur Mac OS X avec la version 1.10 de GDAL de KyngChaos

Site officiel : GDAL/OGR


Creative Commons License
licence Creative Commons Paternité-Pas d'Utilisation Commerciale-Pas de Modification 2.0 France

Commentaires

Poster un nouveau commentaire

Le contenu de ce champ sera maintenu privé et ne sera pas affiché publiquement.