Great, Thank you so much for this. However, I do not seem to have the ST_MapAlgebraExpr(rast, rast.....) function. Was this implemented in rev 8001?
Andreas 2011/11/24 Tom van Tilburg <[email protected]> > Andreas, > > We did approx. the same thing for non-quadrate polygons. Perhaps it might > be useful: > > Step 1: Make a raster from every polygon, based on the grid specifications > of the elevation raster. Here is also the solution: the raster cells will > only be created for cells that have their midpoint *inside* your geometry. > ST_AsRaster(a.geom, b.rast, '<pixtype>') > > Step 2: Overlay the elevation raster with the raster you just created and > keep only the values of the elevation raster > ST_MapAlgebraExpr( > <raster from geometry> > ,b.rast > ,'rast2' -- <-- keep only raster 2 value > , '<pixtype>','INTERSECTION','0','0',0 > ) > > Step 3: get the mean from the statistics on the resulting raster > (ST_SummaryStats( > (ST_Union( -- < --- we did a UNION because we occasionaly had vectors > crossing tiled rasters > <overlay raster> > )).rast > )).mean As avg_height > > That did the trick. Complete script is below. > > I suspect your method of doing a ST_Intersection for every pix. makes it > slower because it creates a geometry first that you do not really need. > > Cheers, > Tom > > ---------- > FULL SCRIPT > > > SELECT > a.gid As id, > (ST_SummaryStats( > (ST_Union( > ST_MapAlgebraExpr( > ST_AsRaster(a.geom, b.rast, '32BF') > ,b.rast > ,'rast2', '32BF','INTERSECTION','0','0',0 > ) > )).rast > )).mean As avg_height > > FROM > polygons.grid a LEFT JOIN > rasters.elev b > ON ST_Intersects(a.geom, b.rast) > GROUP BY a.gid > > > > On 24-11-2011 11:14, Andreas Forø Tollefsen wrote: > > Hi, > > I am trying to calculate the average pixel value in a elevation raster > inside quadrate polygons. > However, I am not getting the correct values from my query: > > SELECT gid, AVG(((foo.geomval).val)) as avgmnt > FROM (SELECT p.gid, ST_Intersection(p.cell, r.rast) AS geomval FROM > mountain r, priogrid_land p WHERE ST_Intersects(p.cell, r.rast, ) AND p.gid > =186124) AS foo > GROUP BY gid ORDER BY gid; > > The problem here is that the ST_Intersects(geom, rast) takes into > consideration the pixels that is outside, but touches the border of the > quadrate polygons. Then, the average values for each quadrate polygon is > affected by pixels inside other polygons. This will potentially lead to a > flawed result. > So what I want is to be able to calculate the average value for the pixels > INSIDE the polygon excluding those outside. > > How can i restrict the AVG pixel value to be calculated only for pixels > that is inside the polygon, and not the pixels that touch the outside of > the border? > > Thanks! > > Best, > Andreas > > > > > _______________________________________________ > postgis-users mailing > [email protected]http://postgis.refractions.net/mailman/listinfo/postgis-users > > > > _______________________________________________ > postgis-users mailing list > [email protected] > http://postgis.refractions.net/mailman/listinfo/postgis-users > >
_______________________________________________ postgis-users mailing list [email protected] http://postgis.refractions.net/mailman/listinfo/postgis-users
