[postgis-users] RE : RE : RE : get unique values from raster intersecting apolygon

Pierre Racine Pierre.Racine at sbf.ulaval.ca
Mon May 28 13:20:35 PDT 2012 

If what you want is the mean of all the pixel centroid falling inside the polygon you do a query like this:

SELECT (ST_SummaryStatsAgg(ST_Clip(rast, geom))).*
FROM clip_build, nclheat
WHERE ST_Intersects(rast, geom)

If what you want is the area weighted mean of all the raster values inside the polygon you do a query like this:

SELECT gt.id,  (aws).weightedmean
FROM (SELECT ST_AreaWeightedSummaryStats(gv) aws
             FROM (SELECT ST_Intersection(rt.rast, gt.geom) gv
                          FROM rasttable rt, geomtable gt
                          WHERE ST_Intersects(rt.rast, gt.geom)
                        ) foo
             GROUP BY gt.id
           ) foo2

In the first case each pixel value is taken into account only once. In the second case only the proportion of the pixel intersecting each polygon is taken into account for each polygon. The first case should be faster than the second case.

See this older post for more details:

http://postgis.refractions.net/pipermail/postgis-users/2012-April/033831.html

Pierre

> -----Original Message-----
> From: postgis-users-bounces at postgis.refractions.net [mailto:postgis-users-
> bounces at postgis.refractions.net] On Behalf Of Francois Hugues
> Sent: Monday, May 28, 2012 3:31 PM
> To: PostGIS Users Discussion
> Subject: [postgis-users] RE : RE : RE : get unique values from raster intersecting
> apolygon
> 
> Ok ! I don't know if I'll really manage to do better but I'll try one last time before
> stopping to disturb everybody.
> 
> I want all the values intersecting the polygons but I don't want a pixel to be
> returned twice or more (if my raster was totally covered I would have as many
> couple gid/value as number of pixels in my raster).
> 
> I attached one more picture to this message. The numbers are the values
> associated to each pixel of the raster. Regarding this picture, I would like a result
> like this (I didn't report all values especially those at the borders of both polygon,
> values 15 and 20 because I don't know which polygon cover the most the pixel),
> the interesting ones are from a-13 to b-23
> 
> gid|val
>  A |2
>  A |4
>  A |5
>  A |6
>  A |8
>  A |9
>  A |10
>  A |11
>  A |12
>  A |13
>  A |17
>  A |18
>  B |16
>  B |19
>  B |22
>  B |23
>  B |24
>  B |25
>  B |26
>  B |27
>  B |30
>  B |31
>  B |32
>  B |33
> 
> Hugues.
> 
> 
> -------- Message d'origine--------
> De: postgis-users-bounces at postgis.refractions.net de la part de Pierre Racine
> Date: lun. 28/05/2012 20:26
> À: PostGIS Users Discussion
> Objet : Re: [postgis-users] RE :  RE :  get unique values from raster intersecting
> apolygon
> 
> What is not clear in your description is: How do you want the value associated
> to each polygon to be computed from the many pixel the polygon is intersecting
> with? Do you want:
> 
> -the raster value corresponding to one point inside the polygon
> -the mean of all the pixel centroid falling inside the polygon
> -the area weighted mean of all the raster values inside the polygon
> -the value most represented inside the polygon
> -etc...
> 
> There are many ways to compute this value.
> 
> Pierre
> 
> > -----Original Message-----
> > From: postgis-users-bounces at postgis.refractions.net [mailto:postgis-users-
> > bounces at postgis.refractions.net] On Behalf Of Francois Hugues
> > Sent: Monday, May 28, 2012 12:40 PM
> > To: PostGIS Users Discussion
> > Subject: [postgis-users] RE : RE : get unique values from raster intersecting
> > apolygon
> >
> > Thanks but i think the way I made my picture was not really clear and I
> > understand now how you read it !
> >
> > I already know the existence of st_value for points. Here, the red point on my
> > picture are the values of the raster and not a table of points : if it is build from
> a
> > DEM, you have the elevation measure based on a regular pattern and finally
> you
> > get a grid and in the centre of each cell the value measured (my red points)
> >
> > For my problem I have two table : one made of polygons, the other a raster
> > made of tiles. I want to get value where polygons intersect raster but a celle
> > must not appear twice. My raster will be totally covered by polygons and at
> the
> > end I want a gid|val table with as many lines as my raster has values.
> >
> > I hope it is clearer now but, for sure, my english does not help !
> >
> > Thanks,
> >
> > Hugues.
> >
> >
> > -------- Message d'origine--------
> > De: postgis-users-bounces at postgis.refractions.net de la part de Pierre Racine
> > Date: lun. 28/05/2012 16:56
> > À: PostGIS Users Discussion
> > Objet : Re: [postgis-users] RE :  get unique values from raster intersecting
> > apolygon
> >
> > I this case the right query looks like this:
> >
> > SELECT gid, ST_Value(rast, geom) val
> > FROM yourraster table, yourpointtable
> > WHERE ST_Intersects(rast, geom)
> >
> > Make sure to tile and index your raster properly (smaller is better). You can do
> > this at import time.
> >
> > The result for cell outside any polygon will be "null".
> >
> > Pierre
> >
> > > -----Original Message-----
> > > From: postgis-users-bounces at postgis.refractions.net [mailto:postgis-users-
> > > bounces at postgis.refractions.net] On Behalf Of Francois Hugues
> > > Sent: Monday, May 28, 2012 9:06 AM
> > > To: PostGIS Users Discussion
> > > Subject: [postgis-users] RE : get unique values from raster intersecting
> > apolygon
> > >
> > > I want one value per red point (which are raster data).
> > >
> > > The gid of one polygon will appear as many times as necessary regarding
> how
> > it
> > > intersects the raster.
> > >
> > > Hugues.
> > >
> > >
> > > -------- Message d'origine--------
> > > De: postgis-users-bounces at postgis.refractions.net de la part de Pierre
> Racine
> > > Date: lun. 28/05/2012 14:24
> > > À: PostGIS Users Discussion
> > > Objet : Re: [postgis-users] get unique values from raster intersecting
> apolygon
> > >
> > > You want one value per polygon or one value per red point? It is not clear to
> > me
> > > from your drawing.
> > >
> > > Pierre
> > >
> > > > -----Original Message-----
> > > > From: postgis-users-bounces at postgis.refractions.net [mailto:postgis-
> users-
> > > > bounces at postgis.refractions.net] On Behalf Of Francois Hugues
> > > > Sent: Saturday, May 26, 2012 5:29 PM
> > > > To: postgis-users at postgis.refractions.net
> > > > Subject: [postgis-users] get unique values from raster intersecting a
> polygon
> > > >
> > > > Hello list,
> > > >
> > > > I have to get values from raster where they intersect polygons from a layer
> > > with
> > > > numerous polygons which have shared boundaries, but I don't want to
> have
> > > > duplicate values from raster. The final result must be something like : gid |
> > val.
> > > >
> > > > To explain my case, I made a picture.
> > > >
> > > > I think ST_intersection could help but I don't know what will be the exact
> > > result.
> > > > In my example, we have a grid (a DEM for example), each red point is the
> > value
> > > > in the centre of each cell and two polygons A and B. I think cells totally
> > within
> > > > polygons are not a problem, but in green I have highlighted some
> ambiguous
> > > > cases :
> > > > - 1 and 2 are cells shared by both polygons but not in the same proportions
> > and
> > > I
> > > > would like to get value 1 for polygon A but not for B and value 2 for
> polygon
> > B
> > > > but not for A
> > > > - What will be the result for cell number 3 and for the other green cell
> > without
> > > > number (and for all the ones which are not totally within a polygon) ? I
> think
> > > the
> > > > value is the same for the whole cell (st_dumpaspolygon) and will be
> returned
> > > > whatever the part of the cell covered by the polygon, but if someone could
> > > > confirm.
> > > >
> > > > If st_intersection use st_DumpAsPolygon, could it be better to work
> directly
> > > with
> > > > this polygonized table using its own gist index and dealing with shared
> areas
> > of
> > > > cells between different polygons ?
> > > >
> > > > Hugues.
> > > >
> > > >
> > > >
> > > >
> > >
> > > _______________________________________________
> > > postgis-users mailing list
> > > 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
> 
> _______________________________________________
> postgis-users mailing list
> postgis-users at postgis.refractions.net
> http://postgis.refractions.net/mailman/listinfo/postgis-users
>

More information about the postgis-users mailing list