Great, Thank you so much for this.<div>However, I do not seem to have the ST_MapAlgebraExpr(rast, rast.....) function. Was this implemented in rev 8001?</div><div><br></div><div>Andreas<br><br><div class="gmail_quote">2011/11/24 Tom van Tilburg <span dir="ltr"><<a href="mailto:tom.van.tilburg@gmail.com">tom.van.tilburg@gmail.com</a>></span><br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
  
    
  
  <div bgcolor="#FFFFFF" text="#000000">
    Andreas,<br>
    <br>
    We did approx. the same thing for non-quadrate polygons. Perhaps it
    might be useful:<br>
    <br>
    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.<br>
        ST_AsRaster(a.geom, b.rast, '<pixtype>')<br>
    <br>
    Step 2: Overlay the elevation raster with the raster you just
    created and keep only the values of the elevation raster<br>
    ST_MapAlgebraExpr(<br>
                <raster from geometry><br>
                ,b.rast<br>
                ,'rast2' -- <-- keep only raster 2 value<br>
                , '<pixtype>','INTERSECTION','0','0',0<br>
            )<br>
    <br>
    Step 3: get the mean from the statistics on the resulting raster<br>
    (ST_SummaryStats(<br>
        (ST_Union( -- < --- we did a UNION because we occasionaly had
    vectors crossing tiled rasters<br>
            <overlay raster><br>
        )).rast<br>
    )).mean As avg_height  <br>
    <br>
    That did the trick. Complete script is below.<br>
    <br>
    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.<br>
    <br>
    Cheers,<br>
     Tom<br>
    <br>
    ----------<br>
    FULL SCRIPT<br>
    <br>
    <br>
    SELECT <br>
    a.gid As id,<br>
    (ST_SummaryStats(<br>
        (ST_Union(<br>
            ST_MapAlgebraExpr(<br>
                ST_AsRaster(a.geom, b.rast, '32BF')<br>
                ,b.rast<br>
                ,'rast2', '32BF','INTERSECTION','0','0',0<br>
            )<br>
        )).rast<br>
    )).mean As avg_height<br>
    <br>
    FROM <br>
    polygons.grid a LEFT JOIN <br>
    rasters.elev b<br>
        ON ST_Intersects(a.geom, b.rast)<br>
    GROUP BY a.gid<div><div></div><div class="h5"><br>
    <br>
    <br>
    On 24-11-2011 11:14, Andreas Forø Tollefsen wrote:
    </div></div><blockquote type="cite"><div><div></div><div class="h5">Hi,
      <div><br>
      </div>
      <div>I am trying to calculate the average pixel value in a
        elevation raster inside quadrate polygons.</div>
      <div>However, I am not getting the correct values from my query:</div>
      <div><br>
      </div>
      <div>
        <div>SELECT gid, AVG(((foo.geomval).val)) as avgmnt </div>
        <div>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 </div>
        <div>GROUP BY gid ORDER BY gid;</div>
      </div>
      <div><br>
      </div>
      <div>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.</div>
      <div>So what I want is to be able to calculate the average value
        for the pixels INSIDE the polygon excluding those outside.</div>
      <div><br>
      </div>
      <div>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?</div>
      <div><br>
      </div>
      <div>Thanks!</div>
      <div><br>
      </div>
      <div>Best,</div>
      <div>Andreas</div>
      <div><br>
      </div>
      <div><br>
      </div>
      <br>
      <fieldset></fieldset>
      <br>
      </div></div><pre>_______________________________________________
postgis-users mailing list
<a href="mailto:postgis-users@postgis.refractions.net" target="_blank">postgis-users@postgis.refractions.net</a>
<a href="http://postgis.refractions.net/mailman/listinfo/postgis-users" target="_blank">http://postgis.refractions.net/mailman/listinfo/postgis-users</a>
</pre>
    </blockquote>
    <br>
  </div>

<br>_______________________________________________<br>
postgis-users mailing list<br>
<a href="mailto:postgis-users@postgis.refractions.net">postgis-users@postgis.refractions.net</a><br>
<a href="http://postgis.refractions.net/mailman/listinfo/postgis-users" target="_blank">http://postgis.refractions.net/mailman/listinfo/postgis-users</a><br>
<br></blockquote></div><br></div>