La gestione dei raster Tiff con GDAL e MapServer non funziona sempre alla perfezione come si vorrebbe. Si deve fare sempre attenzione al tipo di raster con cui si ha a che fare, potrebbe trattarsi di Float32, Int16 oppure Byte.
Con gdalinfo è possibile ottenere l'informazione, ecco alcuni esempi di output su raster costituiti da una singola banda:
Band 1 Block=2160x1 Type=Float32, ColorInterp=Gray NoData Value=-9999
Band 1 Block=2160x1 Type=Int16, ColorInterp=Gray NoData Value=-9999
Band 1 Block=64x64 Type=Byte, ColorInterp=Gray
La conversione con gdal_translate - ad esempio dal formato Arc/Info ASCII Grid - potrebbe non indovinare il formato corretto, con il risultato di troncare alcuni valori.
Ecco come forzare l'utilizzo del tipo adeguato:
gdal_translate -of GTiff -ot Float32 minrainfal.asc minrainfal.tif gdal_translate -of GTiff -ot Int16 maxrainfal.asc maxrainfal.tif
Per ottimizzare la classificazione di un raster Float32 o Int16 con MapServer, conviene aggiungere le direttive per la gestione dei raster tramite LAYER.PROCESSING, in particolare SCALE e SCALE_BUCKETS. Ecco un esempio in cui i valori compresi tra -55.0 e 45.0 vengono classificati in 24 classi:
LAYER
NAME "temperature"
PROCESSING "SCALE=-55.00,45.00"
PROCESSING "SCALE_BUCKETS=24"
CLASS
NAME "< -50"
EXPRESSION ([pixel] < -50.83)
STYLE
COLOR 0 53 255
END
END
CLASS
NAME "-50 ÷ -46"
EXPRESSION ([pixel] >= -50.83 and [pixel] < -46.67)
STYLE
COLOR 0 106 255
END
END
...
Vedere in proposito MapServer Raster Data. Tali direttive non servono (anzi creano problemi) con raster a 8 bit.
Per verificare l'elaborazione del raster si può utilizzare shp2img, in questo modo:
shp2img -m /path/to/file.map -all_debug 10 -o /path/to/raster ... msDrawRasterGDAL_16BitClassification(GEBCO): scaling to 8913 buckets from range=-5600.5,3312.5. ...
Alcune caratteristiche di un file TIFF possono incidere sulle performance di MapServer. Su file molto grossi ad esempio la memorizzazione in blocchi (256×256) piuttosto che in righe di pixel può migliorare i tempi di accesso. Anche la compressione può incidere negativamente sulle prestazioni.
Ecco un comando che converte un file TIFF in un'altro file TIFF, l'organizzazione sarà a blocchi di 256×256 pixel e nessuna compressione. Le informazioni geografiche saranno memorizzate in un world file (.tfw) esterno:
gdal_translate -co TILED=YES -co TFW=YES -co COMPRESS=NONE src.tif dst.tif
Con gdalinfo si potrà vedere il blocksize:
Band 1 Block=256x256 Type=Byte, ColorInterp=Palette
Quei mattacchioni della NASA fanno scaricare il BlueMarble in formato PNG, non georiferito!
Bisognerebbe creare un world file secondo la struttura documentata da Georeference with World Files. Il problema è sapere i parametri! Per fortuna abbiamo un GeoTiff del TrueMarble alla stessa risoluzione e nello stesso sistema di riferimento, con gdalinfo vediamo i parametri del GeoTiff:
gdalinfo TrueMarble.2km.21600x10800.tif ... Origin = (-180.000000000000000,90.000000000000000) Pixel Size = (0.016666666666667,-0.016666666666667) ...
Un trucco per creare il world file a partire da un GeoTiff usando i tool gdal è il seguente:
gdal_translate -co "TFW=YES" TrueMarble.2km.21600x10800.tif temp.tif mv temp.tfw TrueMarble.2km.21600x10800.tfw rm temp.tif
Creiamo il world file con lo stesso nome del file .png, ma estensione .wld oppure .tfw (un world file con tale estensione viene trovato sia da MapServer che dalla libreria GDAL). Le coordinate del centro del pixel in alto a sinistra sono le cordinate dell'origine riportate da gdalinfo, meno la metà di un pixel:
0.0166666667 0.0000000000 0.0000000000 -0.0166666667 -179.9916666667 89.9916666667
Peccato che la velocità di accesso non sia accettabile, il manuale di GDAL recita: PNG files are linearly compressed, so random reading of large PNG files can be very inefficient.
Una soluzione non ottimale è quella di convertire in Tiff non compresso con ImageMagick (il world file da usare rimane esattamente lo stesso):
convert world_topo_bathy.png -compress none world_topo_bathy.tif
La soluzione migliore è creare un GeoTiff completo, che non ha bisogno del world file perché incorpora tutte le informazioni geografiche necessarie. A tale scopo si può usare il comando:
gdal_translate -a_srs EPSG:4326 world_topo_bathy.png world_topo_bathy.tif
Nel GeoTiff risultante vengono incorporate sie le informazioni relative al sistema di riferimento che quelle contenute nel world file.
Se il file è troppo grande conviene creare le overview, dette anche piramidi, altrimenti uno zoom alla massima estensione obbliga a leggere e ridimensionare al volo tutto il raster.
In genere si generano delle copie dell'immagine con risoluzione 1/2, 1/4, 1/8, … dell'originale. Le overview possono essere contenute all'interno dello stesso GeoTiff:
gdaladdo TrueMarble.8km.5400x2700.tif 2 4 8
Se invece vogliamo tenere le overview in un file separato (che avrà estensione .tif.ovr) ed utilizzare un algoritmo di ricampionamento migliore (gauss, ma anche average è migliore del nearest predefinito):
gdaladdo -ro -r gauss TrueMarble.8km.5400x2700.tif 2 4 8
Con le versioni di GDAL < 1.6.0 si deve usare una sintassi diversa, il file esterno con le overview avrà formato Erdas Image ed estensione .aux:
gdaladdo --config USE_RRD YES TrueMarble.8km.5400x2700.tif 2 4 8
MapServer e la libreria Gdal accedono automaticamente alle overview anche se si trovano nel file separato .aux.
Anche se nelle tile sono state aggiunte le overview, ai fattori di scala maggiori (zoom minore) è necessario accedere a tutti i file (tile) della copertura per generare la mappa. Questo fa decadere di molto le prestazioni, è quindi opportuno creare una versione ridotta della copertura, contenuta in un singolo file.
Tale file conviene che sia TIFF memorizzato a blocchi ed eventualmente deve comprendere le overview opportune. Ecco come usare gdal_merge per ottenere il risultato:
gdal_merge.py -o output.tif -of GTiff -co TILED=YES -co TFW=YES -ps 0.9765625 0.9765625 2001/*.jpg
La dimensione del pixel (in metri) viene calcolata dividento la larghezza della copertura per il numero di pixel desiderati. Ad esempio se si ha una superficie larga 16 km e si vuole creare un'immagine larga 16384 pixel:
(670510 - 654510) / 16384 = 0.9765625