Hi Kim, Thanks for your answer. However, we want raster as an output, since we want to be able to use the summarystats function. Please elaborate how you think ST_PixelAsPolygons should solve out issue? Thanks.
Andreas 2013/6/18 Kim Bisgaard <[email protected]> > Hi, > > Try to use 'ST_PixelAsPolygons(ST_Clip(n.rast, p.cell))' > instead of 'ST_Intersects(n.rast, p.cell)' > > Regards, > Kim > > > On 2013-06-18 11:03, Andreas Forø Tollefsen wrote: > > Hi, > > We are working on a raster summarystats script to calculate various > statistics for the pixels within fishnet polygons. > > Our raster cell size is 0.0083333333333... x 0.0083333333333... degrees > while our quadrat polygons are 0.5 x 0.5 decimal degrees. > This should give us 60x60 raster pixels within each of our polygons. > ArcGIS zonal statistics returns a pixel count of 3600 in addition to other > statistics. > However, PostGIS returns 3721 pixel count. > > We do not really understand why, but it seems that our query includes > some pixels that are outside of the polygon, but still touches the vertices > of the polygon and are therefore included in the calculation. > Are there any way of modifying our script to return the same result as > ArcGIS? > Thanks! > > Andreas > > script: > > /* This query makes one raster for each PRIO-GRID cell. Clip and union > is the procedure. */ > INSERT INTO nightlightsprio (gid, "year", rast) > (SELECT gid, "year", ST_Union(raster) as rast > FROM > (SELECT p.gid, n."year", ST_Clip(n.rast, p.cell) as raster > FROM nightlights n, priogridyear p > WHERE ST_Intersects(n.rast, p.cell) > AND n."year" = p."year" > ) > as priorast > GROUP BY gid, "year"); > > > /* Default BandNoDataValue is 0. Raster value 0 means no light, not no > data. Setting to NULL. This produces correct results. */ > UPDATE nightlightsprio2 SET rast = ST_SetBandNoDataValue(rast, 1, NULL); > > > ALTER TABLE nightlightsprio2 ADD COLUMN nightlights_sum double precision, > ADD COLUMN nightlights_mean double precision, > ADD COLUMN nightlights_sd double precision, > ADD COLUMN nightlights_min double precision, > ADD COLUMN nightlights_max double precision, > ADD COLUMN nightlights_count integer; > > UPDATE nightlightsprio2 SET nightlights_sum = > (ST_SummaryStats(rast)).sum; > UPDATE nightlightsprio2 SET nightlights_mean = > (ST_SummaryStats(rast)).mean; > UPDATE nightlightsprio2 SET nightlights_sd = > (ST_SummaryStats(rast)).stddev; > UPDATE nightlightsprio2 SET nightlights_min = (ST_SummaryStats(rast)).min; > UPDATE nightlightsprio2 SET nightlights_max = (ST_SummaryStats(rast)).max; > UPDATE nightlightsprio2 SET nightlights_count = > (ST_SummaryStats(rast)).count; > > > > _______________________________________________ > postgis-users mailing > [email protected]http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users > > > -- > Kim Bisgaard > > Application Development Division Phone: +45 3915 7562 (direct) > Danish Meteorological Institute Fax: +45 3915 7460 (division) > > > _______________________________________________ > postgis-users mailing list > [email protected] > http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users > >
_______________________________________________ postgis-users mailing list [email protected] http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users
