[postgis-users] Rasters in PostGIS

[Roger André]

I had a project where I needed to read the pixel values from a large mosaic
for several hundred points.  The fastest method I was able to come up with
for the query was to write a Python script that did the following:

- open the raster with gdal
- use the point coordinates as input parameters to a function that
calculates the correct raster row and column based on the dimensions of the
raster, and its pixel dimensions

The IO for this tool was crude, as I used a delimited text file to contain
all of the points of interest, but it worked.

Before that, I had also "vectorized" a fairly coarse resolution raster and
stored the pixel geometries as polygons in PostGIS.  It worked very well for
doing relational queries against other data sets, but was murderously slow
for the renderer to style and display.

In both of these examples, I have not come up with a good way to extract a
specific sub-region of the data based on pixel values and return a raster.
Pulling out bounding boxes via the "gdal_translate -projwin" function is
trivial though, so maybe if all you want to do is retrieve an area that
surrounds a point, or contains a group of points, you could do that?

A combination of these methods (geometry in PostGIS, raster on disk) might
be a good solution.  I've attached the Python scripts I used for the
vectorizing, and the raster lookup, in case they might be useful starting
points.

Best of luck,

Roger
--

On Fri, Apr 17, 2009 at 8:44 AM, yoda2nd <yoda2nd at niacc.org> wrote:

> Actually I saw that WKT Raster was avalible, but I was unsure about it with
> it only in the beta stages and there does not appear to be much in the way
> of documentation yet.  I did check out the read me file in the SVN
> repository, and it seems to conflict with the documentation.
> http://trac.osgeo.org/postgis/wiki/WKTRaster/Documentation01
> https://svn.osgeo.org/postgis/spike/wktraster/README
>
> The read me file says that it will not work with the current version of
> PostGIS, but the online documentation says it will.  Not sure which one is
> correct.
>
>
>
> ------------------------------------------------------------------------
>
> I've done this before using a similar approach you have and using GDAL to
> read the elevations from the GeoTiff and dump them into an elevations table
> with point geoms.  If I had to do it again, now that we have WKT Raster in
> the works, I'd probably try that if you don't mind compiling yourself.
>
>
> http://trac.osgeo.org/postgis/wiki/WKTRaster
>
> Mateusz has a nice test drive example of WKTRaster
> http://mateusz.loskot.net/2009/03/30/crunching-wkt-rasters/
>
>
> Hope that helps,
> Regina
>
>
> _______________________________________________
> postgis-users mailing list
> postgis-users at postgis.refractions.net
> http://postgis.refractions.net/mailman/listinfo/postgis-users
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20090417/5c03b878/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: img2shp_v2.py
Type: text/x-python
Size: 2479 bytes
Desc: not available
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20090417/5c03b878/attachment.py>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: raster_intersect.py
Type: text/x-python
Size: 2352 bytes
Desc: not available
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20090417/5c03b878/attachment-0001.py>