[raster] Optimize ST_Intersection to return ST_Value when intersecting with points

[pracine]

It might be much quicker to provide the value associated with a pixel intersecting with a point geometry when using ST_Intersection() by accessing the pixel value with ST_Value instead of vectorizing the whole tile. We have first to fix st_value(rast, geom) so it behave exactly as st_intersect on edge of raster and pixels.

This, obviously, do not work when the intersecting geometry is a line or a polygon. In this case we should try to vectorize only the intersecting part of the raster tile.

[pracine]

Adressed in r5687

[pracine]

Actually no. it was not adressed. To do...

[pracine]

For this we need ST_Value to behave like ST_Intersection by returning a SET OF 2 or 4 values when the requested coordinated fall between two pixels or four pixels.

Actually it is more the opposite: The ST_Value variant taking a point geometry should be "downtimized" (only 3 occurrence in Google !!) as a wrapper around ST_Intersection so that when a point located the intersection of 2 or 4 pixels is passed, the function returns 2 or 4 values since the right behavior is the one of ST_Intersection in this case. For sure this would make this variant of ST_Value much slower since ST_Intersection has to polygonize the tile first which ST_Value do not do right now (it converts the point to a x and a y raster coordinate and returns the value at this coordinate).

To make it to be fast, we must then first optimize ST_Intersection by being able to polygonize a subregion of the tile (the one intersecting with the geometry). This would generally optimize ST_Intersection for any geometry. In the case of ST_Value this region would reduce to a maximum of 4 pixels and hence would be quite fast.

It is already possible to pass a mask to the GDAL polygonize function. That mean we could construct a mask with the area we want to polygonize. Would this be fast? For sure it would be more efficient to be able to just pass x1,y2, x2,y2 to the function to delimitate this area. To see with GDAL people.

[robe]

Pierre,

Actually I rather liked the way ST_Value worked returning only one value. Its convenient in many cases. If you are going to have it return multiple -- then call it ST_values and please leave ST_Value alone. Like for example when I'm grabbing an elevation as a rule of thumb -- I really only want one value. If you must consider all pixels then average them.

[pracine]

I'm still thinking anout it and consult our big PostGIS raster lidar freak user here (Etienne) to try to find the best options possible. At a first glance:

-There should be an option saying to ST_Value should return nodata value as NULL or nodata value. NULL by default.

-There should be an option to specify which one of the four pixel value to return when the point geometry falls on the edge of 2 or more pixels.

-ST_Value should continue to return only one value; use ST_Intersection if you want a vector-like behavior with pixel edges (i.e. possibly many values).

[pracine]

The options in the second case would be 'UL', 'UR', 'LL', 'LR' for upper left, upper right, lower left, lower right. The default would be 'LL' like it is right now (I have to verify).

[pracine]

I will close this one has the problem has somehow changed.

A better optimization is to NOT vectorize the whole tile whatever kind of geometry is passed. This is ticket #1583.

As for the optionnal parameter to ST_Value, I will open a new ticket (#2116).