PostGIS WKT Raster Beta 0.1 Working Specifications

---

B01c - Being able to import TIFF rasters using raster2pgsql.exe

gdal2wktraster.py - A prototype of the translation utility implemented in Python and with use of ​GDAL and its bindings to Python. Create an SQL commands output to create a table of raster. As input raster file, all GDAL formats are accepted. The script is available as ​gdal2wktraster.py script.

Open Questions:

1. Should we change names of options to strictly follow those used by raster2pgsql and shp2pgsql?
   
   Pierre: Yes! We should follow PostGIS tracks as much as possible in everything we do in this project.
   
   
2. How can I import all the bands from a multiband TIFF?

USAGE:
gdal2wktraster.py -r rasterfile [-r rasterfile] -t [<schema>.]<table> [<options>]

-r <rasterfile> Specifies input raster file. Multiple -r options can be specified for a number of input files or wildcards can be used (? and *). At least one input raster file is REQUIRED.
-t <table> Name of destination table in with or without target schema specified. This parameter is REQUIRED.

OPTIONS:
-s <srid> Set the SRID field. Default is -1.
-b <nbband> Specify the number of band. The number of rasterfile must correspond to this number.
-f <field> Name of target column for raster data. Default column name is rast.
-I Create a GiST index on the raster column.
-o <file> Output file for generated SQL commands. If not specified, stdout is assumed.
-R Simply register the raster in the database. Absolute path to the raster and georeferencing informations are stored instead of the raster actual data.
-F Add a "filename" column containing the original name of the loaded raster file.

(-d|c) Mutually exclusive inserting options:

-d Drops the table, then recreates it and populates it with current raster file data.
-c Creates a new table and populates it, this is the default if you do not specify any options.

-v Switch on excessively verbose mode, useful for debugging.
-h Display this help screen.

---

raster2pgsql.exe - Should be as much similar to ​shp2pgsql.exe as it is possible.

• For Beta 0.1: TIFF only, only the first band, options for tiling.

Open Questions:

1. Pierre: Should there be an option to be able to write directly to the database instead of writing a .sql file by default?
   
   
2. Pierre: If we use GDAL, should we allow more file format than only JPEG and TIFF?
   
   Pierre:I guess this comes for free. We might encounter some problems with some format however. I don't know which one yet.

USAGE:
raster2pgsql [<options>] rasterfile [rasterfile…] [<schema>.]<table>
Create an SQL commands file to create a table of raster. If rasterfile is multiband and –b is not specified, every band are inserted. Multiple band can also be specified using multiple filenames (rasterfile1 is the first band, rasterfile2 the second, etc…). Can process multiple file from a folder. georeference (and pixel size) must exist directly in the files or in a companion World File.

OPTIONS:
-s <srid> Set the SRID field. Default is -1.
-b <nbband> Specify the number of band. The number of rasterfile must correspond to this number.
-P <pixeltypes> Specify the pixels types in which to store each band. Ex. ‘8-bit unsigned integer,16-bit float’. conversion may happens.
-n <nodata values> Specify the nodata value for each bands. Ex. ‘0,0.0’. Default to ‘none’ for each band.
-t <pixels> Divide rasters into <pixels>x<pixels> tiles, one tile per row. Default is to store whole rasters as one row.

(-d|a|b|c|p) Mutually exclusive inserting options:

-d Drops the table, then recreates it and populates it with current raster file data. -a Appends raster file into current table, must be exactly the same pixel size, number of band, nodata value and pixel type. -B Appends raster files as a new bands. When tiled with the –t option, the new band is inserted tiled in the same way as the original band. -c Creates a new table and populates it, this is the default if you do not specify any options. -p Prepare mode, only creates the table.

-r <raster_column> Specify the name of the raster column (mostly useful in append mode).
-D Use postgresql dump format (defaults to sql insert statements).
-I Create a GiST index on the bbox of the raster column.
-? Display this help screen

---

B01d - Being able to get all the properties of a raster (all the “Get” functions).

RT_SRID(raster|geometry) → integer
Return the SRID associated with the raster.

RT_Width(raster) → integer
Return the width of the raster.

RT_Height(raster) → integer
Return the height of the raster.

RT_PixelSizeX(raster) → float64
Return the georeference's X pixel size of the raster. ​See.

RT_PixelSizeY(raster) → float64
Return the georeference's Y pixel size of the raster. ​See.

RT_RotationX(raster) → float64
Return the georeference's X rotation.

RT_RotationY(raster) → float64
Return the georeference's Y rotation. ​See.

RT_UpperLeftX(raster) → float64
Return the georeference's X-coordinate of the center of the upper left pixel. ​See.

RT_UpperLeftY(raster) → float64
Return the georeference's Y-coordinate of the center of the upper left pixel. ​See.

RT_GeoReference(raster) → string
Return the georeference of the raster as a string representing the 6 doubles of an equivalent world file (including the carriage return). ​See.

RT_NumBands(raster) → integer
Return the number of band included in the raster.

RT_BandPixelType(raster, integer) → string
Return the pixel type of the specified 1-based Nth band of raster. Band index is 1-based. The function returns one of the following values:

• 1BB - 1-bit boolean
• 2BUI - 2-bit unsigned integer
• 4BUI - 4-bit unsigned integer
• 8BSI - 8-bit signed integer
• 8BUI - 8-bit unsigned integer
• 16BSI - 16-bit signed integer
• 16BUI - 16-bit unsigned integer
• 32BSI - 32-bit signed integer
• 32BUI - 32-bit unsigned integer
• 16BF - 16-bit float
• 32BF - 32-bit float
• 64BF - 64-bit float

RT_BandNoDataValue(raster, integer) → float32
Return the NoDataValue of the specified 1-based Nth band of raster. Band index is 1-based. The value is always returned as a float32 even if the pixel type is integer.

---

B01e - Being able to intersect vector and raster to produce vector.

RT_GetBBox(raster) → polygon geometry
RT_Envelope(raster|geometry) → polygon geometry
RT_Shape(raster) → polygon geometry
RT_AsPolygon(raster) → polygon geometry set
RT_Intersects(raster|geometry, raster|geometry)
RT_Intersection(raster|geometry, raster|geometry, ‘raster’|’geometry’)→raster/geometry

---

B01f - Being able to return a JPEG, a TIFF or a PNG.

RT_bytea(raster, band) → raster
What is does?

---

Open Question: When exporting a multiband raster to JPEG, TIFF, PNG, SVG or KML, how should we specify the band number in the exporting function.

There is two options to select the band to convert from a multiband raster in all the RT_AsFormat functions.

1. Precede each call with RT_Band() to return a selected band.
   

Pros: This is a general function that can be called before any function that would otherwise require a band parameter.
Cons: This implies creating a temporary raster. This might be more elegant and general but is this too much overhead comparing with having a band parameter?

1. Add a band parameter to each RT_AsFormat function.
   

Pros: Hypothetically less overhead.
Cons: Every functions implying access to a band should then have this parameter when in most case it would be equal to 1. In many cases it makes no sence to have to specify a band parameter since it is the whole raster that we want to export, including all the bands.

Pierre: More I think about it more I think that the first option is the best one…

mloskot: Perhaps there is a compromise in form of two sets of functions: 1) RT_As* which always burn the whole raster (all bands) 2) RT_BandAs* which takes number of band as a parameter and return only this requested band.

---

RT_Band(raster, band) → raster
Return a single band from a multiband raster. If "band" is greater than the value returned by RT_GetNumBands(), the function returns the last band. This function should be used to select a band before converting it to JPEG, TIFF, PNG, SVG or KML with the corresponding function. e.g. RT_AsJPEG(RT_Band(raster, band), quality)

RT_AsJPEG(raster, quality) → JPEG as "bytea"
Return the raster as a JPEG encoded as a PostgreSQL bytea. By default quality is set to 75, but this option can be used to select other values. Values must be in the range 10-100. Low values result in higher compression ratios, but poorer image quality. Values above 95 are not meaningfully better quality but can but substantially larger. (copied from ​http://www.gdal.org/frmt_jpeg.html)

---

Open Question: Is JPEG export limited to raster having 8 bit unsigned integer pixeltype (8BUI)?

​See how GDAL do it. It converts only 8 bits rasters. Should we do the same?

Otherwise, how do we convert other types to 8BUI? e.g. 16BUI or 8BSI?

Pierre: It might be more simple to ignore pixeltypes other than 8BUI but it would be very convenient to have a way to quickly export elevation data for example as a JPEG. It would be nice to have an elegant solution to this. Maybe something inspired from MapServer.

Proposition one (Pierre): RT_AsJPEG could simply (optionally when the pixeltype is not 8BUI) map the RT_Maximum() and RT_Minimum() value to 0-255. RT_Maximum() and RT_Minimum() are not in the spec yet but this could be on nice usage of it. They will imply caching the min and max when importing and editing. Both function should ignore the NoDataValues. They could also be two parameters passed to RT_AsJPEG(raster, quality, min, max).

Proposition two: There could also be just one parameter (string) defining a mapping method:

• Method "None": No mapping. This is possible only for 8BUI.

• Method "MaxMinValue": Use the Max and the Min cached in the raster. e.g. for 16BSI (min, max) → (-2033, 2456) → (round((-2033 - -2033)/(2456 - -2033)*255), round((2456 - -2033)/(2456 - -2033)*255)) → (0, 255).
  

This is equivalent to RT_AsJPEG(raster, quality, RT_Minimum(rast), RT_Maximum(rast))

• Method "MaxMinType": Use the Max and the Min allowed by the type. e.g. for 16BSI (min, max) → (-2033, 2456) → (round((-2033 - -32768)/(32767 - -32768)*255), round((2456 - -32768)/(32767 - -32768)*255)) → (120, 137)
  

This would be equivalent to RT_AsJPEG(raster, quality, RT_BandPixelTypeMin(rast), RT_BandPixelTypeMax(rast)). Both functions (RT_BandPixelTypeMin & BandPixelTypeMax) are not yet planned and I could not find an SQL query that returns the equivalent range for a type. ​One possible solution.

mloskot: ATM, I have no thoughts on this issue.

---

Open Question: Is JPEG export limited to raster having 1 or 3 bands?

​See how GDAL do it. It converts only 1 or 3 band rasters. Should we do the same? In this case 1 band rasters would be exported as a greyscale JPEG having R G and B identical and 3 band rasters would be interpreted as R, G and B.

Pierre: I think the answer should be yes. I don't see how we could have a 2 band raster fit into RGB.

mloskot: I agree, the answer should be yes.

---

RT_AsTIFF(raster, compression) → TIFF as "bytea"
Return the raster as a JPEG encoded as a PostgreSQL bytea. If raster is a multiband raster and no band were selected with RT_Band() every band are written to the resulting TIFF.

compression=[JPEG/LZW/PACKBITS/DEFLATE/CCITTRLE/CCITTFAX3/CCITTFAX4/NONE]: Set the type of compression to use. None is the default. The CCITT compression should only be used with 1bit (NBITS=1) data. JPEG should only be used with Byte data. When using JPEG add a number specifying the quality. 75 is the default. e.g. RT_AsTIFF(raster, "JPEG60") (copied from ​http://www.gdal.org/frmt_gtiff.html)

---

Open Question: What if we want to export only the first two band of a three band layer?

Maybe we need a RT_RasterFromBands(band1,band2,etc…) to reconstitute a multiband raster from multiple sources (having the same width, height, pixelsize, etc…)

mloskot: or RT_RasterFromBands(bands) where bands is ARRAY[int]. For instance, RT_RasterFromBands(ARRAY[1,3]) will burn new raster from 1 and 3 bands of input raster.

---

RT_AsPNG(raster, band) → PNG as "bytea"

---

RASTER_COLUMNS Metadata Table

The following metadata table provides external applications and possibly internal tools the ability to discover tables with raster data, and to learn various information about those datasets.

Column Name	Type	Constraints	Comments
r_table_catalog	character varying(256)	NOT NULL
r_table_schema	character varying(256)	NOT NULL
r_table_name	character varying(256)	NOT NULL
r_column	character varying(256)	NOT NULL
srid	integer	NOT NULL	Use 0 for unknown SRID
pixel_types	ARRAY=VARCHAR=	NOT NULL	an array of pixeltypes, one per band. The band count is implicit in the size of this array.
out_db	boolean	NOT NULL	false if internal tiles, true if tiles are references to files outside the database
regular_blocking	boolean	NOT NULL	false by default, true if all blocks are equal sized, abutted and non-overlapping, started at top left origin (see below)
nodata_values	ARRAY=DOUBLE=		an array of nodata values, one per band. The entry may be NULL to indicate no nodata values.
pixelsize_x	double		width of a pixel in geounits (per SRID)
pixelsize_y	double		height of a pixel in geounits (per SRID)
blocksize_x	integer		the width of a block in pixels (NULL if irregular)
blocksize_y	integer		the height of a block in pixels (NULL if irregular)
extent	GEOMETRY		a polygon geometry containing all raster tiles, or NULL if predefined bounds are not known. For "regular_blocking" cases this geometry will be a simple rectangle. In other cases it might be an irregular polygon.

If the regular_blocking field is true a number of restrictions are placed on the raster column that is defined:

1. All tiles must have the same size (blocksize_x and blocksize_y).
2. All tiles must be non-overlapping, and appear on regular block grid.
3. The top left block must start at the top left corner of the extent.
4. The right most column, and bottom row of blocks may have portions that extend beyond the raster extent. These areas will be assumed to be nodata and not part of the described raster.
5. The extent field must be a simple rectangle (non-rotated).

It is permissable to for regular_blocking rasters to omit some blocks/tiles (sparse rasters) in which case the missing tiles should be assumed to be all nodata or zero valued.

Specification Comments following the Code Sprint

This section is dedicated to discussion that took place during the ​Toronto Code Sprint 2009.

1. Bug tracking - Launch separate Google project or join PostGIS bug tracker, as a module? Answer: We will track our bug in the PostGIS Google Code project. PostGIS plan on moving (back) to Trac soon.
   
   
2. Wiki: is there a chance to have a better wiki for documentation, schedule, specs, road map, users inputs? Actual PostGIS wiki limitations: security, no images, no tables (or sophisticated page structure) (For this Leo and Regina are drafting up a PSC and then once we have a PSC we'll vote to locate on OSGEO. Main issue is Refractions server is kind of old)
   
   
   1. Is Google Code Wiki well working?
      
      (remember export issues with Google code - see gvSIG - Cuba)
      
      
3. Windows build. What is the plan?
   
   
   1. First step is to review Windows/Visual C++ native build configuration for PostGIS. ML: I will try to investigate and solve outstanding issues like PGXS makefiles generation, soon.
      
      
   2. ML: Prefer native build toolsets on supported platform (Windows/Visual C++, GCC/Linux, etc.). On Windows, avoid Cygwin, MinGW, wherever possible. Issue directly related to point a).
      
      
   3. Is cmake a possibility?
      
      
4. GiST indexing: is wrapping PostGIS Gist function in plpgslq sufficient? From Paul comments it seems actually that cast are sufficients.
   
   
5. API: Any other SQL raster function you would like to see?
   
   
   1. ML: RT_AsJPEG - What JPEG implementation to use?
      
      
   2. ML: RT_AsTIFF - Do we aim to depend on libgeotiff?
      
      
   3. ML: Is there any function that would require/fit GDAL dependency?
      
      
   4. Pierre: We should maybe investigate the cost of using GDAL and have a general RT_AsImage(format, params) instead of using the individual tiff, jpeg and png libraries. Those will anyway also have a cost I guess… And we will probably have to link with GDAL for AsPolygon(rast) and AsRaster(geom) anyway.
      
      
6. Implementing raster → vector conversion for intersections(). How to do it? Possible using only plpgsql? So we don't have to link with liblwgeos and we can have our own release schedule.
   • polygonize with GDAL or GRASS, that generates some internal structure for those libraries, and you turn that into a WKT form of polygon. That can be passed to PostGIS, and you don't need to link to liblwgeom.
     
     
7. Best practice for band referencing: having a band parameter (index?) or always use RT_Band()
   
   -e.g. SELECT RT_AsTIFF(RT_Band(raster,2),60) FROM coverandtemp
   
   SELECT RT_AsTIFF(RT_Band(raster,2),RT_Band(raster,1),60) FROM coverandtemp
   or
   SELECT RT_AsTIFF(raster,2,60) FROM coverandtemp
   
   Pierre: Having a band parameter (or a list of band) for each function trying to access a band would be preferable to avoid useless memcopy of band data.
   
   
8. Architecture & design issues:
   
   
   • How to handle no NODATA? Pierre from Frank input: We'll need a HasNoDataValue flag for each band.
     
     
9. Documentation
   
   
   • Setup a structure (likely on a Wiki), sections for developers and for users
     
     
   • Update existing documents like RFCs - they are out of date.
     
     
   • Describe key terms: canonical format, serialized format, transfer format (WKB), in-db raster, out-db raster, etc. This should help people to grasp the WKT Raster key elements quickly, and help them to decide which document/RFC to read in details.
     
     
10. Development Roadmap
    
    
    • Review current schedule.
    • Add note on a website where users should (or should not) expect 1 release - clear statement.