PostGIS栅格子集表查询返回异常及优化创建方法咨询

背景 #

• 使用raster2pgsql创建存储栅格数据的表grib_data
• 创建表turb_mod，计划存储grib_data中值超过特定阈值的点子集
• 创建表的查询语句：

WITH turb AS (SELECT rid, rast, (ST_PixelAsPoints(rast)).val AS val
FROM grib_data
)
SELECT rid, rast INTO turb_mod
FROM turb WHERE val > 0.5;

• 创建时返回“SELECT 53”，预期turb_mod包含53行数据

问题 #

执行SELECT (ST_PixelAsPoints(rast)).x AS x FROM turb_mod;时，返回了原表的所有记录，而非预期的53条

疑问解答 #

1. 为何查询未返回53条记录？ #

问题出在turb_mod的创建逻辑上：

• ST_PixelAsPoints(rast)会把单个栅格拆分成多条像素点记录，CTE里的turb是将每个栅格的所有像素点展开，每条记录对应一个像素点+完整的原始栅格
• WHERE val > 0.5筛选出了53个符合条件的像素点对应的rid和rast，但这里的rast是整个原始栅格，不是单个像素点
• 因此turb_mod里的53行记录实际是53个完整原始栅格（可能存在重复，同一栅格可能有多个符合条件的像素点），再次执行ST_PixelAsPoints(rast)时，会把每个完整栅格的所有像素点拆分出来，最终返回原表所有记录

2. 更优的创建子集表方法 #

如果需要的是符合阈值的像素点（用于空间聚类），应直接存储这些像素点的几何信息，而非完整栅格。推荐两种方式：

方式一：直接提取像素点字段

CREATE TABLE turb_mod AS
SELECT
  rid,
  (ST_PixelAsPoints(rast)).geom AS geom,
  (ST_PixelAsPoints(rast)).val AS val
FROM grib_data
WHERE (ST_PixelAsPoints(rast)).val > 0.5;

方式二：用LATERAL JOIN更清晰处理像素点展开

CREATE TABLE turb_mod AS
SELECT
  g.rid,
  p.geom,
  p.val
FROM grib_data g
LATERAL ST_PixelAsPoints(g.rast) p
WHERE p.val > 0.5;

这两种方式创建的turb_mod，每条记录对应一个符合条件的像素点，包含空间几何（geom）和值（val），后续执行空间聚类（如ST_ClusterDBSCAN）时可直接基于geom字段操作，效率更高。

内容的提问来源于stack exchange，提问作者swepilot