[postgis-users] Raster pixel value

Andreas Forø Tollefsen andreasft at gmail.com
Thu Nov 24 04:20:25 PST 2011


I know. RTFM before asking dumb questions. Set the ST_Summary(rast, false)
so it does not exclude_nodata_value? Right?

Andreas

2011/11/24 Andreas Forø Tollefsen <andreasft at gmail.com>

> Great. I updated to latest rev and added the st_union.sql to my function
> list. Now it works.
> However, I need to ask one more question.
> When running the script, it calculates the average pixel value inside the
> polygon. It does not take into account the NULL DATA values.
> Is there any way i can make the script calculate the average pixel value
> but count NULL DATA values as 0?
>
> For instance, one pixel in a cell belonging to the Comoros have a value
> (0.95). The rest is NULL. Hence, it seems a bit wrong that the whole
> polygon get 0.95 when almost none of its area have any elevation?
>
> Any ideas?
>
> Thank you so much for this info Tom!
>
> Andreas
>
> 2011/11/24 Tom van Tilburg <tom.van.tilburg at gmail.com>
>
>>  Yes, I think so.
>> Current windows build is *r8221
>> *
>> Also, I think ST_Union(raster) (which I used in the example) is not
>> included in this version yet.
>> You would have to download a prototype from.
>>
>> http://trac.osgeo.org/postgis/browser/trunk/raster/scripts/plpgsql/st_union.sql
>>
>> Chrs,
>>  Tom
>>
>>
>> On 24-11-2011 12:07, Andreas Forø Tollefsen wrote:
>>
>> 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 <tom.van.tilburg at gmail.com>
>>
>>>  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 listpostgis-users at postgis.refractions.nethttp://postgis.refractions.net/mailman/listinfo/postgis-users
>>>
>>>
>>>
>>> _______________________________________________
>>> postgis-users mailing list
>>> postgis-users at postgis.refractions.net
>>> http://postgis.refractions.net/mailman/listinfo/postgis-users
>>>
>>>
>>
>>
>> _______________________________________________
>> postgis-users mailing listpostgis-users at postgis.refractions.nethttp://postgis.refractions.net/mailman/listinfo/postgis-users
>>
>>
>>
>> _______________________________________________
>> 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/20111124/982fa467/attachment.html>


More information about the postgis-users mailing list