[postgis-users] ST_Buffer + grid problem

[Ed Linde]

Hi Pierre,
I am trying to compute the "actual distances" in meters before creating the
raster and then re-project it back to SRID 4326.

-- Trying to now get the Top left point of raster
SELECT     ST_X(ST_Transform(ST_SetSRID(ST_Point(8.07734039737749,
57.7505109647578),4326),32632)),
    ST_Y(ST_Transform(ST_SetSRID(ST_Point(8.07734039737749,
57.7505109647578),4326),32632));


-- Width
select  ceil(ST_X(ST_Transform(ST_SetSRID(ST_Point(15.1919565742587,
57.7505109647578),4326),32632)) -
    ST_X(ST_Transform(ST_SetSRID(ST_Point(8.07734039737749,
57.7505109647578),4326),32632)));

-- Height
select  ceil(ST_Y(ST_Transform(ST_SetSRID(ST_Point(8.07734039737749,
57.7505109647578),4326),32632)) -
    ST_Y(ST_Transform(ST_SetSRID(ST_Point(8.07734039737749,
54.4984986588244),4326),32632)));

I have now made it a multiple of 4 so now I will get height = 361932 and
width = 423152

Is this ok? Once I have made the raster, I can then transform it back to
srid 4326 before doing PixelAsPolygons right?

Cheers,
Ed


On Fri, Mar 23, 2012 at 1:55 PM, Pierre Racine
<Pierre.Racine at sbf.ulaval.ca>wrote:

> > Can you please let me know how I can do this? So I did the make empty
> raster ->
> > map algebra -> pixel as polygons approach!
> > Since there is no dependency on data now, is it possible for you to maybe
> > outline the steps you performed at your end to get this working?
> > MapAlgebra needed 2 rasters, so I inputted the same raster twice. So
> please
>
> There is a one raster version of ST_MapAlgebraExpr(). That's the one to
> use....
>
> > --SELECT (md).*
> > --FROM (SELECT ST_MetaData(rast) As md
> > --    FROM (select ST_MakeEmptyRaster( 8, 4, 8.07734039737749,
> > 57.7505109647578, 0.000036,0.000036,0,0, 4326 ) rast) foo) foo2;
> >
> >
> > create table vector_grid as
> > select (ST_PixelAsPolygons(rastfin)).geom ,
> (ST_PixelAsPolygons(rastfin)).val
> > from ( select ST_MapAlgebraExpr(rast, rast, '([rast.x] - 1) * '  ||
> > ST_Width(rast)::text || ' + [rast.y]' ) rastfin from (select
> ST_MakeEmptyRaster( 8,
> > 4, 8.07734039737749, 57.7505109647578, 0.000036,0.000036,0,0, 4326 )
> rast)
> > foo) foo2;
>
> If you want to assign the index values to the raster:
>
> SELECT ST_MapAlgebraExpr(ST_MakeEmptyRaster(8, 4, 8.07734039737749,
> 57.7505109647578, 0.000036,0.000036,0,0, 4326 ), '32BUI', '([rast.x] - 1) *
> 8 + [rast.y]'::text)
>
> If you are happy with the polygons:
>
> SELECT (gvxy).geom, ((gvxy).x - 1) + (gvxy).y cellid
> FROM (SELECT ST_PixelAsPolygons(ST_MakeEmptyRaster(8, 4, 8.07734039737749,
> 57.7505109647578, 0.000036,0.000036,0,0, 4326 )) gvxy) foo
>
> Pierre
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20120323/3c037dff/attachment.html>