[postgis-users] Raster pixel value
Tom van Tilburg
tom.van.tilburg at gmail.com
Thu Nov 24 06:42:11 PST 2011
Andreas,
Would it be possible to share code and/or sample data?
With me St_SummaryStats seems to work correct, even when I set all my
raster values to 1.
Chrs,
Tom
On 24-11-2011 13:42, Andreas Forø Tollefsen wrote:
> I double checked the values that I get from this query and they do not
> seem to be correct.
> For polygon cells that have only pixels with the value 1 inside only
> have a mean of 0.85 (see image).
> The correct value should be a mean of 1. After some visual controls I
> am not sure if these numbers are correct.
> There are many polygon cells with only pixel value 1, but no polygon
> cell with a higher mean than .87.
>
> Any idea why this could be? Both the raster and polygons are in srid 4326.
>
> Andreas
>
>
>
> 2011/11/24 Andreas Forø Tollefsen <andreasft at gmail.com
> <mailto:andreasft at gmail.com>>
>
> 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
> <mailto: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
> <mailto: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
>> <mailto: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 list
>>> postgis-users at postgis.refractions.net <mailto:postgis-users at postgis.refractions.net>
>>> http://postgis.refractions.net/mailman/listinfo/postgis-users
>>
>>
>> _______________________________________________
>> postgis-users mailing list
>> postgis-users at postgis.refractions.net
>> <mailto:postgis-users at postgis.refractions.net>
>> http://postgis.refractions.net/mailman/listinfo/postgis-users
>>
>>
>>
>>
>> _______________________________________________
>> postgis-users mailing list
>> postgis-users at postgis.refractions.net <mailto:postgis-users at postgis.refractions.net>
>> http://postgis.refractions.net/mailman/listinfo/postgis-users
>
>
> _______________________________________________
> postgis-users mailing list
> postgis-users at postgis.refractions.net
> <mailto:postgis-users at postgis.refractions.net>
> http://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/56c9e005/attachment.html>
More information about the postgis-users
mailing list