[postgis-devel] WKTRaster ST_MapAlgebra

[Pierre Racine]

>In both applications, the map algebra operations are performed on cells that exist in both rasters. If
>the rasters have different extents, then the resulting raster has data only in the region where they
>intersect.  This is also true when rasters have the same extent and contain NODATA values.

That means in:

SELECT ST_MapAlgebra(rast1, rast2, mathExpr) FROM mytable

an implicit ST_Intersects(rast1, rast2) must be performed by the ST_MapAlgebra function between the two raster before any computation is attempted. For sure a query with an explicit Intersects() like this:

SELECT ST_MapAlgebra(rast1, rast2, mathExpr) FROM mytable WHERE ST_Intersects(rast1, rast2)

should be faster since ST_MapAlgebra() is most of the time not even called.

>In both applications, the cell size of the resulting raster is the same as the largest input cell
>size.  So if two rasters are provided as input, one with a cell size of 5, and another with a cell
>size of 9, the resulting raster has a cell size of 9.  In ArcGIS, I was able to confirm that even when
>the grid orientation is different (when comparing 2 different projections, for example), the cells are
>interpolated to the larger cell size.

This let me think that the first task should be to well design ST_Resample() which is planned anyway. This function should be very usefull for any work on MapAlgebra().

>In ArcGIS, when map algebra is performed on two rasters with different projections/coordinate systems,
>the first raster wins.  The second raster is reprojected on the fly to the first coordinate system,
>and the resulting raster has the same coordinate system as the first input raster.  I was not able to
>figure out how to perform map algebra in GRASS with rasters in two different projected coordinate
>systems.

We should not have this problem since we should assume that both raster are in the same coordinate system. Does PostGIS warn us when Intersecting two polygons in two different CS? We should adopt the same warning policy.

>Finally, a question: how about I start fleshing out ST_MapAlgebra (spec) with the following rules?
>* Rasters with different projections cannot be combined.
>* Rasters with different cell sizes will produce an output raster with the largest of the input raster
>cell sizes.

So we have first to resample the one with smallest cell size implicitely (in the function).

>* Rasters with different extents produce an output raster with the extent of the intersection of the
>input rasters.

What append when they have the same CS and pixel size but the intersecting area have different alignment (or rotation)? Does the resample operation imply realignment? Is it possible to resample without realigning (I want to transform to the same pixel size but keep my alignment)?

Pierre