我有一个几何列存储2D点一个巨大的表(gps_points)。我试图做到的是运行一个查询输出像
id | freq
-------------
1 | 365
2 | 1092
3 | 97
...
其中,“id”为我的总边界框和“频率”中的小矩形的唯一标识符是落在该特定的矩形内的点的数量。
所以,我已经定义了一个PostGIS的表:
create table sub_rects (
id int,
geom geometry)
我然后运行一个脚本外,我在那里生成1000×1000这样的矩形和创建它们的多边形,所以我得到一百万行是这样的:
insert into sub_rects values(1,ST_GeomFromText('POLYGON((1.1 1.2, 1.1 1.4, 1.5 1.4, 1.5 1.2, 1.1 1.2))'));
当然除了每个多边形本身得到一组新的坐标,以在1000×1000网格在我的GPS数据的边框坐标匹配其实际位置,以及ID获取每个元组更新。
然后我生成空间索引并在此表中的主键索引。
最后,我可以运行这个表,我的原始数据表(gps_points)与
select id, count(*) from sub_rects r join gps_points g on r.geom && g.geom group by id;
这给了我我试图输出。问题是,它需要永远载入所有的小多边形,而我想生成具有不同数量的矩形或运行在具有不同的基本坐标数据集的地图每一次,我不得不放弃sub_rects并产生和重新加载它。
是否有这样做的更好的办法?我不需要图形输出。我只需要生成数据。不必产生支撑台(sub_rects)从外部将是非常好的,我怀疑有办法更少的计算完成同样的事情的昂贵的方法。我更希望不安装任何额外的软件。
ETA:按照意见中的要求,这里是查询计划(我家的机器上,所以更小的数据集和其他表名,但同样的计划):
gisdb=# explain analyse select g.id id, count(*) from gridrect g join broadcast b on g.geom && b.wkb_geometry group by g.id;
QUERY PLAN
-------------------------------------------------------------------------------------------------------------------------------------------------------------------
GroupAggregate (cost=0.57..177993.58 rows=10101 width=12) (actual time=14.740..3528.600 rows=1962 loops=1)
Group Key: g.id
-> Nested Loop (cost=0.57..144786.36 rows=6621242 width=4) (actual time=13.948..3050.741 rows=1366376 loops=1)
-> Index Scan using gridrect_id_idx on gridrect g (cost=0.29..485.30 rows=10201 width=124) (actual time=0.079..6.582 rows=10201 loops=1)
-> Index Scan using broadcast_wkb_geometry_geom_idx on broadcast b (cost=0.29..12.78 rows=137 width=32) (actual time=0.011..0.217 rows=134 loops=10201)
Index Cond: (g.geom && wkb_geometry)
Planning time: 0.591 ms
Execution time: 3529.320 ms
(8 rows)
和2:
按照在回答我建议修改的建议有此代码:
(SELECT row_number() OVER (ORDER BY geom) id, geom
FROM (SELECT st_geomfromtext(
concat('Polygon((', x || ' ' || y, ',', x + xstep || ' ' || y, ',', x + xstep || ' ' || y + ystep,
',', x || ' ' || y + ystep, ',', x || ' ' || y, '))')) geom
FROM (SELECT x, y
FROM (SELECT generate_series(xmin, xmin + xdelta, xstep) x) x,
(SELECT generate_series(ymin, ymin + ydelta, ystep) y) y) foo) bar);
其中XMIN,YMIN,xdelta,ydelta,XSTEP和ystep全部由外部脚本计算,但可能只是如果包裹上面的函数调用中也被计算为一个Postgres功能的一部分。从生成该临时表,然后执行该查询是两个数量级,比我在做什么,最初更快。
两件事情。首先创建SQL级(从pgAdmin的举例)表。
create table polygons as
select st_geomfromtext(concat('Polygon((',x||' '||y,',',x||'
'||y+0.2,',',x+0.4||' '||y+0.2,',',x+0.4||' '||y,',',x||' '||y,'))')) geom
FROM (select generate_series(0,199.9,0.2) x) x,
(select generate_series(0,199.9,0.4) y) y
创建索引
创建使用要旨(的geom)多边形索引;
然后使用您的查询或这一个。检查哪一个会更快,你的情况
select id, count(*)
from sub_rects r
join gps_points g on st_dwithin(r.geom, p.geom, 0)
由ID组;
下面是从生成边界框的网格的例子:
https://gis.stackexchange.com/questions/16374/how-to-create-a-regular-polygon-grid-in-postgis
要生成密度数据,首先尝试所有的数据创建一个临时表,然后得到计数。在我的经验,下面的比全部梳理到一个查询有所加快:
create temp table rect_points as
select r.id as rect_id, p.id as point_id
from sub_rects r, gps_points p
where p.geom && r.geom;
create index idx on rect_points (rect_id);
select rect_id, count(*) from rect_points group by rect_id;