Matthias,

This is a different approach but it might be useful. This is untested, but hopefully it will give you the idea:

select
    floor((st_x(a.the_geom)-b.minx)/1000)::integer as col,
    floor((st_y(a.the_geom)-b.miny)/1000)::integer as row,
    count(*) as cnt
  from
    (
        select floor(min(st_x(the_geom)))::integer as minx,
               floor(min(st_y(the_geom)))::integer as miny
          from points ) as b,
    points a
  group by col, row
  order by cnt desc;

With the col, row, values and the count you can create the polygons after the fact.

-Steve W

On 10/17/2011 6:07 PM, Matthias Ludwig wrote:
Hello,

for a statistical analysis I need to create a polygon grid to count the points 
inside. The analysis area measures about 15x15km. At the moment I use the code 
below. It's functioning but is slow for high point counts(>5000-100000 points) 
and small polygon sizes (5-20m). Because the points are clustered and the big area 
size it's not possible to insert all polygons and delete those which are empty. 
Presently one polygon is inserted for every point in the first step. In the second 
step dublicated polygons where deleted. How can I insert one polygon in a specific 
place?

CREATE TABLE point_fishnet (gid serial not null primary key, area text, count 
integer default 0);
SELECT addgeometrycolumn('point_fishnet','the_geom', 900913, 'POLYGON', 2);

INSERT INTO point_fishnet(the_geom, area)
        (
                SELECT ST_Translate(bb.the_geom, x_koord, y_koord) AS the_geom,
                        'germany'
                FROM    generate_series(
                                (
                                        SELECT  
floor(min(ST_X(the_geom)))::integer
                                        FROM    points
                                ),
                                (
                                        SELECT  
ceil(max(ST_X(the_geom)))::integer
                                        FROM    points                  
                                ),
                                1000    --Offset = Polygonsize
                        ) AS x_koord,
                        generate_series(
                                (
                                        SELECT  
floor(min(ST_Y(the_geom)))::integer
                                        FROM    points
                                ),
                                (
                                        SELECT  
ceil(max(ST_Y(the_geom)))::integer
                                        FROM    points
                                ),
                                1000    --Offset = Polygonsize
                        ) AS y_koord,
                        (
                                SELECT ST_MakeEnvelope(0, 0, 1000, 1000, 
900913) AS the_geom
                        ) AS bb,
                        (
                                SELECT  the_geom
                                FROM    points
                        ) AS point
                WHERE   ST_Intersects(
                                ST_Translate(bb.the_geom, x_koord, y_koord),
                                point.the_geom
                        )
        );

-- delete dublicates
DELETE FROM point_fishnet
WHERE gid NOT IN
(
        SELECT min(gid)
        FROM point_fishnet
        GROUP BY the_geom HAVING count(*)>= 1
);
_______________________________________________
postgis-users mailing list
[email protected]
http://postgis.refractions.net/mailman/listinfo/postgis-users

_______________________________________________
postgis-users mailing list
[email protected]
http://postgis.refractions.net/mailman/listinfo/postgis-users

Reply via email to