GRASS GIS pas à pas pour les débutants: 10- les modules d'interpolation spatiale: splines régularisées avec tension (RST) et splines bilinéaires et bicubiques avec régularisation de Tykhonov

	Niveau	Débutant
	Logiciels utilisés	GRASS GIS
	Plateforme	Windows | Mac | Linux | FreeBSD

Petit rappel, c'est quoi une spline ?

« Dans le domaine mathématique de l'analyse numérique, une spline est une fonction définie par morceaux par des polynômes.
Dans les problèmes d'interpolation, la méthode des splines est souvent préférée à l'interpolation polynomiale, car on obtient des résultats similaires en se servant de polynômes ayant des degrés inférieurs » dans Wikipedia fr.wikipedia.org/wiki/Spline

La couche test sur laquelle les modules seront appliqués ((et les sondages, pour les géologues) a été présentée dans la partie 7 de même que les rappels:

n'oubliez pas de lire l'aide (onglet Manual) pour comprendre les options (onglets Optionals, Settings ou Parameters) de chaque interface présentée.

  = couches vectorielles points et lignes

  = couches rasters « points» et « lignes»

= couches vectorielles points et lignes 3D

Splines régularisées avec tension (RST)

C’est une méthode d'interpolation qui cherche une surface lissée la moins déformée possible passant par tous les points d'origine. Elle a été développée théoriquement par Helena Mitasova et son équipe et portée en FORTRAN pour GRASS dès 1993.

La meilleure analogie est l'étirement d'une feuille en caoutchouc de manière à ce qu'elle passe par tous les points traités.  Le premier paramètre, la tension, fixe la rigidité et l'élasticité de la feuille. Le deuxième paramètre contrôle l'« adoucissement » de la feuille, c'est-à-dire sa tendance à passer par les points analysés.

Ces deux aspects ont été séparés dans les implémentations postérieures dans d'autres logiciels SIG comme ArcGIS qui sépare les splines régularisées des splines de tension (help.arcgis.com/fr/arcgisdesktop/10.0/help/index.html#//009z00000078000000) ou dans SAGA GIS (sourceforge.net/apps/trac/saga-gis/wiki/grid_spline)

C'est sans doute le module le plus avancé de GRASS GIS du fait de sa longue histoire. Avant de commencer, autant vous prévenir tout de suite, c'est aussi un des modules les plus compliqués à utiliser du fait des très nombreuses options disponibles. Leurs rôles sont très difficiles à appréhender de prime abord (allez voir grass.fbk.eu/gdp/html_grass64/v.surf.rst.html pour en avoir une idée, il m'a fallu un temps fou pour en maîtriser certains). Une bonne introduction est donnée dans « Processing digital terrain models by regularized spline with tension: tuning interpolation parameters for different input datasets » de Cebecauer, T., Hofierka, J., Suri, M., dans Proceedings of the Open source GIS - GRASS users conference 2002 - Trento, Italy, 11-13 September 2002 disponible à www.ing.unitn.it/~grass/conferences/GRASS2002/proceedings/proceedings/pdfs/Cebecauer_Tomas.pdf

Dans le cas idéal, la surface interpolée passera donc par les points d'origine, mais, du fait de l'adoucissement et/ou de la tension, elle pourra s'étendre au-delà de la plage des valeurs originales. (voir grass.osgeo.org/wiki/RST_Spline_Surfaces)

En entrée:

Commandes:

C'est la commande v.surf.rst, spéciale, car elle permet, si on le désire, de générer plusieurs surfaces en même temps (surface, pente, aspect, courbure, etc.).

Interface:

Ce module possède de nombreux paramètres en option...

Résultat:

Le paramètre Tension a été expliqué, Smoothing, comme son nom l'indique, contrôle l'adoucissement de la surface, c'est-à-dire sa capacité à passer par tous les points, ce qui est le cas s’il est égal ou proche de 0.

• paramètres Smoothing (adoucissement) et Tension (rigidité, élasticité) mis ou proches de 0. La surface passe bien par tous les points dans le premier cas, la surface est une « plaque », ne passant par aucun point dans le deuxième cas

Comment varient ces paramètres ?

• variations du paramètre Tension seul avec un adoucissement (Smoothing) constant:

• variations du paramètre Smoothing (adoucissement) avec une tension (Tension) constante:

On peut observer que l'interpolation de la surface résultante arrive vite à un équilibre: peu de variations entre les valeurs de Smoothing 5 à 25 par exemple. Si on le désire, GRASS fournit des outils pour estimer cet équilibre.

Il est aussi possible de fixer le nombre maximum de points utilisés pour le calcul d'un segment (segmax):

• variations:

Et bien d'autres comme theta qui permet de fixer une direction préférentielle de traitement (dans le cas des couches géologiques en fonction de leur direction illustrées dans Grass Gis et Paraview: la modélisation géologique 3D pour les démunis, ou les enthousiastes...). Pour vous donnez une idée de la précision qu'il est possible d'obtenir, voici les paramètres utilisés pour construire une surface géologique qui passe exactement par les points rouges dans la figure ci-dessous, je la contrains par:

• le paramètre Smoothing  (= 0.01, pour passer exactement par tous les points);
• une Tension élevée (de 100 à 200, pour éviter les effets d'« escaliers »);
• theta, la direction préférentielle, qui est ici la direction de l'axe des plis (en bleu);
• et le paramètre scalex=0.5 qui signifie que je donne deux fois plus d'importance à l'axe x (élongation) qu'à l'axe y car la stucture est allongée dans cette direction;

Dans le cas de très nombreux points, comme les données LIDAR, il est possible de segmenter le processus de calcul pour éviter des temps trop longs et des fichiers temporaires trop lourds.

Cerise sur le gâteau, il est aussi possible de sérier le traitement par lots en fonction de la variable traitée (ici les altitudes):

Et je vous fais grâce des autres surfaces obtenues destinées à contrôler le résultat...

Conclusions

J'espère vous avoir montré que c'est la méthode la plus avancée dans GRASS GIS. L'aide du module affirme qu'il est généralement possible de se contenter du résultat obtenu avec les options par défaut, mais, une fois les paramètres compris, on a toujours envie d'aller plus loin... C'est la combinaison de toutes les options qu'il est le plus difficile de maitriser.

Son grand défaut est le lissage: la surface obtenue est plus lisse que celles obtenues par les autres méthodes. En effet, en adoucissant les gradients d’altitude, l’écart entre les points est réduit et les ondulations de relief peuvent être gommées (pas de falaises ici...).

Cette méthode est préconisée au sein de l’IGN Fr pour des cas précis:

« la modélisation qui permet la meilleure représentation des zones sommitales est l’interpolation RST. Comme l’interpolation RST cherche à obtenir un aspect lissant sur l’ensemble de la surface, les sommets et les zones de talwegs se trouvent être un peu corrigés du manque d’information. Par cette interpolation, ils obtiennent une forme générale assez convexe sur les zones sommitales et concaves sur les zones de talwegs qui pallie, dans une certaine mesure, l’erreur de sous-estimation de ces zones. C’est l’interpolation qui présente le biais le moins important» dans Rousseaux, 2006, B.9 Caractérisation d’erreurs sur un modèle numérique de terrain en fonction de zones morphologiques, Bulletin d’information scientifique et technique de l’IGN n° 75 (2006/1), recherche.ign.fr/doc/BI75/B09_75_ROUSSEAUX.pdf).

Le traitement des couches vectorielles lignes sera vu dans la partie GRASS GIS pas à pas pour les débutants: 11- les modules d'interpolation spatiale: traitement des courbes de niveau.

Pour les spécialistes, l'interpolation par splines est équivalente à un Krigeage avec une covariance généralisée fixée mais ici, il n'y a aucune hypothèse préalable à vérifier, seul le résultat compte.

Pour les géologues

v.surf.rst permet aussi très facilement de traiter les surfaces des couches recoupées par les sondages et d'illustrer l'effet du paramètre smoothing (0.01 à gauche, 5 à droite à tension constante):

Splines bicubiques ou bilinéaires avec régularisation de Tykhonov (B-splines)

Les B-splines sont une généralisation des courbes de Bézier que nous utilisons régulièrement dans des logiciels comme Inkscape ou Adobe Illustrator. En mathématiques, les interpolations bilinéaires et bicubiques sont des extensions des interpolations linéaires ou cubiques pour interpoler les points de données sur une grille régulière à deux dimensions. L'interpolation bilinéaire se base sur des polynômes d'ordre 1, la bicubique, sur des polynomes d'ordre 3 (alors que La méthode du plus proche voisin, vue dans la partie 8, est une interpolation polynomiale d'ordre 0).

En simplifiant outrageusement, les techniques consistent à trouver les fonctions splines passant par tous les points traités dans chaque maille d'une grille. L'algorithme par splines bicubiques est l'algorithme d'interpolation numérique le plus précis mais il nécessite un paramètre supplémentaire: la taille d'une maille de la grille. Leurs principaux défauts sont un adoucissement des surfaces résultantes. Le paramètre de régularisation de Tykhonov permet donc d'agir pour contrôler cet adoucissement: avec une faible valeur, la surface interpolée passera par ou près des points originaux, avec une plus grande, la surface sera plus lissée.

En entrée:

Commande:

v.surf.bspline

Interface:

Les paramètres à prendre en compte sont plus simples que ceux du module précédent: longueur, largeur de chaque maille de traitement, algorithme d'interpolation (bicubique ou bilinéaire) et paramètres de régularisation de Tykhonov.

Résultat:

Le résultat peut être soit une couche matricielle (ici), soit une couche vectorielle, mais pas les deux en même temps.

Comme toujours, pour les comprendre, faisons varier les paramètres :

• longueur x largeur d'une maille de traitement:

 

• Paramètre de Tychonov

Si on le désire, le module peut calculer les paramètres qui permettent d'affiner le paramètre de Tychonov.

Conclusions

Dans GRASS, les avantages et les défauts sont sensiblement les mêmes que ceux de l'algorithme RST.

Tous les traitements ont été effectués sur Mac OS X avec GRASS GIS 6.4.1 ou  6.4.2  et Inkscape pour les figures.

Site officiel : GRASS GIS
Autres Liens : GRASS GIS pas à pas pour les débutants: 7 - les modules d'interpolation spatiale: introduction, passage d'une couche vectorielle à une couche matricielle (v.to.rast) et d'une couche vectorielle 2D à une d'une couche vectorielle 3D (v.to.3D, v.drape)
Autres Liens : GRASS GIS pas à pas pour les débutants: 8 - les modules d'interpolation spatiale: partitionnement de l'espace, Delaunay (TIN), polygones de Voronoï, algorithmes des plus proches voisins
Autres Liens : GRASS GIS pas à pas pour les débutants: 9 - les modules d'interpolation spatiale: algorithmes avec pondération, inverse de la distance (IDW) et Krigeage géostatistique
Autres Liens : GRASS GIS pas à pas pour les débutants: 11- les modules d'interpolation spatiale: traitement des isolignes (courbes de niveau,etc.) et conclusions générales

licence Creative Commons Paternité-Pas d'Utilisation Commerciale-Partage des Conditions Initiales à l'Identique Commerciale 2.0 France