I have a table of polygons representing travel time isochrones on particular days. For each origin point, there are five isochrone geometries (stored in separate rows). For each origin point, I want to rasterise the five isochrones (a binary NULL or 1), and then combine them into a single raster layer. This raster layer requires a simple map algebra: sum/5, so that each origin will in the end be associated with a single raster layer that has values in [NULL, 0.2, 0.4, 0.6, 0.8, 1] depending on how many of the constituent layers overlap. It is a probability surface.
My data is all stored in Postgres 9.3 (with PostGIS). My problem is that while I want to learn to use PostGIS raster, it seems to have a really steep learning curve, and all of the examples I can find deal with a single raster layer. In the examples, this layer is used as part of a polygon overlay, perhaps averaging the value of the raster for each polygon. I haven't found a replicable example for combining: a) vector --> raster b) map algebra; and c) GROUP BY attribute as per my first paragraph.
I'm OK using GDAL or GRASS if I have to in order to perform this task, but this seems like something PostGIS should be able to handle; it would be convenient to do so given my input data is already PostGIS geometry; and I really want to come to terms with PostGIS raster.
Some sample data structure:
areaid time date isogeom (polygon) 1000 07:15:00 2014-05-05 xxx 1000 07:15:00 2014-05-06 xxy ... 1006 07:15:00 2014-05-05 zzz
I want to rasterise, group by areaid, and then perform the map algebra to come to:
areaid isorast (raster) 1000 aaa 1006 bbb
I have not been successful containing this to PostGIS. My approach has been to convert the vector to raster, dump the rasters to arrays, and perform the combination with numpy arrays via psycopg2, before writing them to a GeoTIFF (to maybe be put back in PostGIS). Not ideal, but do-able.