Changes between Version 66 and Version 67 of WKTRasterTutorial01


Ignore:
Timestamp:
Jun 17, 2010, 1:27:53 PM (14 years ago)
Author:
pracine
Comment:

Legend:

Unmodified
Added
Removed
Modified
  • WKTRasterTutorial01

    v66 v67  
    9696}}}
    9797
    98 We downloaded 13 zipped files. If you have 7-Zip installed, you can simply select all the zip files, right-click on them and select "7-Zip->Extract files..." 7-Zip will extract every TIF file in the same folder. Each file is 6000 pixels wide by 6000 pixels high. You can quickly preview the images if you have a good images viewer like !IrfanView installed. Efficient raster/vector analysis with PostGIS WKT Raster requires raster files to be split into tiles when loaded in the database. Fortunately gdal2wktraster.py has an option to tile the images the size we like. We will tile them into 100 pixels x 100 pixels resulting in 6000/100 * 6000/100 * 13 = 46800 tiles. gdal2wktraster.py also accepts wildcard so we will be able to load our 13 big rasters in one unique command.
     98We downloaded 13 zipped files. If you have 7-Zip installed, you can simply select all the zip files, right-click on them and select "7-Zip->Extract files..." 7-Zip will extract every TIF file in the same folder. Each file is 6000 pixels wide by 6000 pixels high. You can quickly preview the images if you have a good images viewer like !IrfanView installed. Efficient raster/vector analysis with PostGIS WKT Raster requires raster files to be split into tiles when loaded in the database. Fortunately gdal2wktraster.py has an option to tile the images the size we like. We will tile them into 50 pixels x 50 pixels resulting in 6000/50 * 6000/50 * 13 = 187200 tiles. gdal2wktraster.py also accepts wildcard so we will be able to load our 13 big rasters in one unique command.
    9999
    100100The command to load the raster files in PostGIS WKT Raster looks like:
    101101
    102102{{{
    103     >"C:/Program Files/PostgreSQL/8.4/bin/gdal2wktraster.py" -r C:\Temp\TutData\SRTM\tif\*.tif -t srtm_tiled -s 4326 -k 100x100 -I > C:\Temp\TutData\SRTM\srtm.sql
     103    >"C:/Program Files/PostgreSQL/8.4/bin/gdal2wktraster.py" -r C:\Temp\TutData\SRTM\tif\*.tif -t srtm_tiled -s 4326 -k 50x50 -I > C:\Temp\TutData\SRTM\srtm.sql
    104104}}}
    105105
     
    110110The -s option, identical to shp2pgsql.exe one, is required to specify the spatial reference system ID. In this case the raster are in "WGS 84" having the SRID number 4326. Unlike some GIS software, PostGIS does not support on the fly reprojection so that we cannot do operations on table stored with different spatial reference systems. As we could see, the caribou point layer and the SRTM images are in different spatial reference systems. The points are in "NAD 83/Quebec Lambert" and SRTM are in "WGS 84". We will have to deal with this problem later.
    111111
    112 The -k option specify the size of the tiles we want to load in PostGIS. Every input raster will be split into 100x100 tiles. This dimension is a good compromise allowing efficient raster/vector analysis. It is better if the size of the tiles is a divider of the size of each raster. Otherwise the last colomns and rows of tiles of each raster will be filled with nodata values. This might have an impact on performance but not on the result since WKT Raster analysis functions ignore nodata values.
     112The -k option specify the size of the tiles we want to load in PostGIS. Every input raster will be split into 50x50 pixels tiles. This dimension is a good compromise allowing efficient raster/vector analysis. It is better if the size of the tiles is a divider of the size of each raster. Otherwise the last colomns and rows of tiles of each raster will be filled with nodata values. This might have an impact on performance but not on the result since WKT Raster analysis functions ignore nodata values.
    113113
    114114The -I option tells the loader to create a spatial index on the raster tile table. The index is very important as it allow PostGIS WKT Raster to restrict his computing efforts only to the tiles involved in a spatial operation. In this tutorial for example, the intersection operations will be performed only on the tiles that actually intersects with the caribou points and it is much faster to search for those tiles if they are spatially indexed than try them one after the other sequentially in the raster table.
     
    323323}}}
    324324       
    325 This query takes about 30 minutes to complete. It can greatly be optimized by invoquing the st_intersection only once with the help of a subquery producing record values that we split later in the main query like this:
     325This query takes about 13 minutes to complete. It can greatly be optimized by invoquing the st_intersection only once with the help of a subquery producing record values that we split later in the main query like this:
    326326
    327327{{{
     
    338338}}}
    339339       
    340 This version takes only 17 minutes. Almost half the time spent by the preceding one... This little increase in complexity certainly worths the gain in performance. The size of the tiles has also a major effect here. If we restart the whole process with 60 pixels x 60 pixels tiles (60 is a divider of 6000) the intersection operation takes 9 minutes. 200 pixels x 200 pixels: 47 minutes... 30 pixels x 30 pixels: 5 minutes. You choose...
     340This version takes only 7.4 minutes. Almost half the time spent by the preceding one... This little increase in complexity certainly worths the gain in performance. The size of the tiles has also a major effect here. If we restart the whole process with 100 pixels x 100 pixels tiles the intersection operation takes 9 minutes. 200 pixels x 200 pixels: 47 minutes... 30 pixels x 30 pixels: 5 minutes. You choose...
    341341
    342342Once again you can visualise the result of the query by doing this query in OpenJUMP (you might have to increase RAM allowed to OpenJUMP by changing "-Xmx256M" to "-Xmx1024M" in "OpenJUMPInstallPath/bin/openjump.bat"):