A Modification des raster
AA Modifier la résolution
La fonction ST_Resize prend un raster en argument et retourne le même raster ré-échantilloné. La géographie totale du raster ne change pas mais le nombre de points change – donc leur résolution – pour correspondre aux dimensions données en argument.
Si pour X ou Y le nombre de pixels en sortie n’est pas un multiple de celui en entrée, un algorithme est utilisé pour interpoler la valeur des pixels « qui tombent entre deux ». Les algorithmes disponibles sont par ordre de rapidité (inverse de la qualité) : NearestNeighbor, Bilinear, Cubic, CubicSpline, Lanczos.
Exemple avec notre modèle numérique de terrain (MNT) importé dans le dernier module :
SELECT rid,
ST_Resize
(rast, 0.05, 0.05) rast FROM mnt_25m
Dans cet exemple, on ré-échantillonne le raster en ne gardant sur chaque dimension que 5 % des points (soit 0,0025 % des points au total).
Afin de visualiser le résultat avec QGIS, on crée une table à partir de cette requête :
CREATE TABLE mnt_resize AS
SELECT rid,
ST_Resize
(rast, 0.05, 0.05) rast FROM mnt_25m;
Sans oublier d’ajouter les contraintes (sinon QGIS n’affichera pas la couche) :
SELECT AddRasterConstraints(’mnt_resize’, ’rast’);
Au lieu de donner des pourcentages à ST_Resize (en nombre réels), on peut lieu donner des dimensions absolues. Étant donné que les tuiles font ici 100 pixels de côté, on aurait pu écrire en nombre entiers (attention : 5 est un nombre entier, 5.0 est un nombre réel) :
SELECT rid,
ST_Resize
(rast, 5, 5) rast FROM mnt_25m
Cependant, parce que notre raster de taille 401x401 compte 5 tuiles sur chaque dimension (4 tuiles de 100 pixels et une de 1 pixel), la dernière tuile plus petite se « densifie » et monte en résolution pour avoir le même nombre de pixels que les autres (5x5), au lieu de rester en 5x1 ou 1x5. La converture raster devient alors hétérogène en résolution, ce qui n’est pas un problème en soi tant qu’on le prend en compte dans les opérations qu’on souhaite mener ensuite.
AB Reprojection
Comme pour les géométries, la fonction ST_Transform reprojette le raster en entrée dans le référentiel spatial (SRID) donné. À cet effet, PostGIS ré-échantillonne les pixels par application de la fonction de transformation du référentiel source vers celui de destination.
CREATE TABLE mnt_4326 AS
SELECT rid,
ST_Transform
(rast, 4326) rast FROM mnt_25m;
SELECT AddRasterConstraints(’mnt_4326’, ’rast’);
On observe des perturbations en noir, à cause de la division en tuiles : chaque tuile a été reprojetée indépendemment. La solution consiste à donner en argument à ST_Transform un raster de référence au lieu du SRID. Pour notre exemple, on utilisera simplement ST_Union pour agréger les tuiles en un seul raster avant de reprojeter. En conséquence la table mnt_4326_2 aura un raster monolithique non-tuilé :
CREATE TABLE mnt_4326_2 AS
SELECT 1 as rid,
ST_Transform
(
ST_Union
(rast), 4326) rast FROM mnt_25m;
SELECT AddRasterConstraints(’mnt_4326_2’, ’rast’);
AC Extraction de région ou polygone
La fonction ST_Clip permet de couper un raster sur un polygone donné.
CREATE TABLE mnt_clip1 AS
SELECT 1 as rid,
ST_Clip
(
ST_Union
(rast, 1), 1,
ST_Buffer
(
ST_Centroid
(
ST_Envelope
(
ST_Union
(rast))), 3000),
false) as rast
FROM mnt_25m;
SELECT AddRasterConstraints(’mnt_clip1’, ’rast’);
On travaille sur la couverture raster dans sa globalité donc on aggrège avec ST_Union . L’argument suivant 1 est le numéro de canal à utiliser : notre MNT n’en a qu’un de toute façon.
À partir de l’enveloppe ( ST_Envelope ), on calcule le centroïde qui est donc le centre de notre MNT. ST_Buffer transforme le point en un disque (qui est physiquement un polygone), lequel est passé à ST_Clip comme forme pour le découpage.
Le dernier argument, false , indique à ST_Clip de ne pas couper à la plus petite enveloppe possible contenant toutes les valeurs qui ne sont pas « no data ».
On observe cependant que du noir entoure le disque : c’est à cause d’une mauvaise gestion des « no data ». À l’aide de la fonction d’interrogation de QGIS, on s’aperçoit que ST_Clip a mis des valeur zéro pour toute la zone coupée, qui n’est pas considérée comme un « no data ».
Quelques ajustements sont à faire :
CREATE TABLE mnt_clip1 AS
SELECT 1 AS rid,
ST_SetBandNoDataValue
(
ST_Clip
(
ST_Union
(rast), 1,
ST_Buffer
(
ST_Centroid
(
ST_Envelope
(
ST_Union
(rast))), 3000),
0, false),
0) AS rast
FROM mnt_25m;
SELECT AddRasterConstraints(’mnt_clip1’, ’rast’);
