= '''GDAL Driver for PostGIS Raster Working Specifications''' =
{{{
#!div  style='background-color: #F4F4F4; padding: 5px; border: 1px solid gray; float: right; margin-left: 5px; width: 260px; font-size: small;' > 

'''Quick Links'''

 * [wiki:WKTRaster PostGIS Raster Home Page]
 * [wiki:WKTRaster/PlanningAndFunding Planning & Funding]

 * [wiki:WKTRaster/SpecificationWorking02 Working Specifications for PostGIS 2.0]

 * [wiki:WKTRaster/SpecificationFinal01 Old Final Specifications for Beta 0.1.6]

}}}

== '''Current status of the driver (February 2011)''' ==

The driver is:
 * Able to read in-db evenly blocked rasters (all blocks with same size)

 * Able to read in-db one-row-rasters:
   * If the table really has more than one row: using -where clause in connection string
   * If the table has more than one row: the table must have been marked as "regularly blocked table", with -k in loader

 * Able to generate two kind of raster object based on two modes:
   - ONE_RASTER_PER_ROW ('mode = 1' in connection string, or nothing): The default mode. Each table row is considered as an independent raster. If the table required has more than one row, and no -where clause has been specified in connection string, all the table rows will be considered as reported as Subdatasets. Unless you specify the other working mode
   - ONE_RASTER_PER_TABLE ('mode = 2' in connection string): Each table is considered as a raster coverage, and each row is a raster tile. 

 * Too slow (Reads the entire table metadata for constructing the [http://www.gdal.org/classGDALDataset.html GDALDataset object], and needs one server round per [http://www.gdal.org/classGDALRasterBand.html GDALRasterBand::IReadBlock] call)

The driver is not:
 * Able to read out-db rasters (developed, but not tested, and with known bugs)
 * Able to create new rasters
 * Able to manage all the [http://trac.osgeo.org/postgis/raw-attachment/wiki/WKTRaster/Documentation01/WKTRasterArrangements.gif PostGIS Raster arrangements]
 * Able to provide a color interpretation for bands


== '''Design principles''' ==

=== Topic: The basis ===
The main class of a GDAL driver is [http://www.gdal.org/classGDALDataset.html GDALDataset]: A set of associated raster bands. So, 1 GDALDataset must be able to contain:
  * An untiled image stored in a raster table's row.
  * A tiled image stored in a raster table (regular or irregular, rectangular or not, with or without missing tiles, with or without overlapping between tiles)
  * A raster object coverage from the rasterization of a vector coverage stored in a raster table (regular or irregular, rectangular or not, with or without missing tiles, with or without overlapping between tiles)

In the first case, 1 GDALDataset = 1 PostGIS Raster object. In the other two cases, 1 GDALDataset = Several PostGIS Raster objects. For this reason, the GDAL PostGIS Raster driver has '''2 working modes''': ONE_RASTER_PER_TABLE, ONE_RASTER_PER_ROW.

However, currently the driver only deals with continuous tiled raster layers, when all the raster tiles are the same size, snap to the same grid and do not overlap (the ideal case).

'''Open question''': Are 2 working modes enough to manage all the raster arrangements?

 Pierre: I think yes. We have to distinguish "want we want to produce" from "what we have to deal with". The two modes answer "want we want to produce" and the different table arrangement are "what we have to deal with". 

 From a GDAL user point of view I know there is a bunch of raster rows in the DB and there is only two things I want to do: extract those rasters rows one by one creating one raster per row or treat them all as a single raster and blend them all together. Furthermore, I want to be able to SELECT those rows using a WHERE statement. If I want a single raster from the db, I have to build my WHERE clause accordingly. There is no need for an extra mode for this. Beside, I don't want to know or have to know what is the raster table arrangement. I expect the driver is able to deal with them all.

 Then, the driver has to deal with all the possible arrangement of those selected rows in both mode (overlap, gaps, missing tiles, etc...). You tried to enumerate the posssible arrangement above but I think there is only two cases: the tiles are regularly tiled or they are not, whatever the number of tile there is (1 or more). To me the irregular case is a generalization of the first one. 

 If, and only if, you can optimize the regularly tiled case, then you write is as an exception. The problem is to make sure the table is REALLY regularly tiled without relying on the user knowledge. Just the introduction of the -a option to raster2pgsql.py allowing to append tiles to an existing table make the "regularly blocked" flag untrustable. If really we want to maintain this flag we will have to create something like a ST_ValidateRegularBlocking aggregate function.

 
----
=== Topic: Constructing the GDALDataset object ===

To construct a GDALDataset object, the driver must:
  * Open the dataset (create db connection)
  * Read some data about the dataset (metadata): srid, georeference information, projection information, raster data size, band information (number of bands, pixel size, color interpretation, if present), any other driver-specific dataset related information (i.e.: in our case, schema and table name)
  * Construct the structure for raster bands, with instances of [http://www.gdal.org/classGDALRasterBand.html GDALRasterBand] class. You need to provide some basic information: data type (pixel size), block size (GDAL contains a concept of the natural block size of rasters so that applications can organized data access efficiently for some file formats) and color interpretation (if any).

The metadata must be read from the raster table, using SQL functions like ST_Extent (used for raster data extent), ST_Metadata (used for general raster metadata) or functions like ST_SRID, ST_Width, ST_Height, etc. When your GDALDataset matches only one raster row (a raster tile) this is not a problem. But when your GDALDataset matches a whole raster table (ONE_RASTER_PER_TABLE mode), you have 2 options:
  * Call the functions over the whole table and filter the result (i.e.: select distinct st_srid(rast) from raster_table, select distinct st_metadata(rast) from raster table). It can be a really slow operation, but you can check if all tiles are like expected (for example: if they are the same size, if they share the same srid, if they overlap or not, etc)

 '''Pierre comment''': I think the driver should not try to detect bad raster arrangement with SQL queries. It should just get what it needs from the DB and burn rasters tiles as they come. If it detected that the arrangement is regular then burn them as regular. If they are not, then burn them accordingly.

  * Call the functions limiting the output to one result. Fast operation, but may be incorrect

 '''Pierre comment''': What might be incorrect? There should not be different srid, pixeltype, or pixelsize in the same table. We have have to warn that we do not support this bad arrangement yet.

Currently, the driver takes the first (and slow) option. That caused performance problems (see ticket #497)

'''Open question''': How to fetch the information needed to construct the GDALDataset? Pay attention to the fact that '''you are not asking for raster data yet'''. You only need metadata, for constructing the basic GDALDataset object.

 Pierre: So just knowing the raster extent, the pixelsize and the pixeltype is not sufficient? You could do this with a quick query SELECT ST_Extent(rast), min(ST_BandPixelSize(rast, band)), min(ST_BandPixelType(rast, band)) assuming all tiles have the same pixelsize and pixeltype.

----
=== Topic: !Reading/Writing raster data ===

Once constructed the basic structure (GDALDataset object and related GDALRasterBand objects), you need to choose the strategy for raster data reading/writing:

  * ''Natural'' block oriented r/w: The driver reads/writes data in equal sized blocks. The potentially more efficient way of r/w data. Really, the natural block size for this dataset is chosen during GDALRasterBand creation. So, '''it's driver's responsibility to provide the desired value for block size'''. To use this method, your driver must provide an implementation of [http://www.gdal.org/classGDALRasterBand.html#09e1d83971ddff0b43deffd54ef25eef IReadBlock].

  * Region oriented r/w: The driver reads/writes arbitrary regions of data. It's a potentially less efficient method, because you have to take care of '''data type translation''' if the data type of the buffer is different than that of the GDALRasterBand. You also must takes care of '''image decimation / replication''' if the buffer size (nBufXSize x nBufYSize) is different than the size of the region being accessed (nXSize x nYSize). To use this method, your driver must provide an implementation of [http://www.gdal.org/classGDALRasterBand.html#5497e8d29e743ee9177202cb3f61c3c7 IRasterIO].

 '''Pierre question:''' In which case "the data type of the buffer is different than that of the GDALRasterBand"?

 '''Pierre question:''' In which case "the buffer size (nBufXSize x nBufYSize) is different than the size of the region being accessed (nXSize x nYSize)"?

 '''Pierre question:''' Do we agree that "Region oriented" is a general case of "''Natural'' block oriented" and that IReadBlock can be implemented as a wrapper around IRasterIO?

Clearly, there's no best method for reading/writing data in our case. In the ideal case of regulary blocked rasters, with no overlapping and same grid for all tiles, the block oriented r/w is the more appropiate strategy. But in the rest of the cases, a more general r/w method must be provided.

Currently, the natural block oriented r/w method is the one implemented for the driver. This is a limitation for 2 reasons:
  * Obviously, it only fits one raster arrangement
  * Each !ReadBlock call forces a new server round, constructing a Box and getting the raster row that contains it. This can be really slow, in case of huge raster coverages (question raised in ticket #497 too). 

 '''Pierre question:''' Why does each IReadBlock force a new server round and IRasterIO not?

'''Open question''': How to get the needed metadata in case of ONE_RASTER_PER_TABLE arrangement. As argued in ticket #497, executing ''ST_Extent'' or ''ST_Metadata'' without limits over a big table can be a really heavy process.

 '''Pierre question:''' How long take a ST_Extent query on a 1 000 000 tiles indexed table?

'''Open question''': What should be the general r/w algorithm? 
jorgearevalo: I think the strategy ''read as much data as you can'' should be the right one, to minimize server rounds. This is: construct a query that, using ''ST_Intersects'', fetches as much rows as possible. This query would be executed in ''IRasterIO'' method. But I don't know how to choose the geographic limits for the query (how much data is ''as much data as you can''?)

 '''Pierre question:''' What do you mean by "I don't know how to choose the geographic limits"? Can't you deduce them from the parameter passed to the function and the metadata of the rasters?

----