postgis_raster

pig repo add pgdg -u          # 添加 PGDG 仓库并更新缓存

CREATE EXTENSION postgis_raster CASCADE;  -- 依赖: postgis

CREATE EXTENSION postgis_raster;

# 将 GeoTIFF 以 100x100 瓦片导入，创建空间索引，使用 COPY
raster2pgsql -s 4326 -t 100x100 -I -C -M elevation.tif public.dem | psql mydb

# 追加到已有表
raster2pgsql -s 4326 -t 100x100 -a more_data.tif public.dem | psql mydb

-- 获取栅格尺寸和像素大小
SELECT rid,
    ST_Width(rast) AS width,
    ST_Height(rast) AS height,
    ST_ScaleX(rast) AS pixel_size_x,
    ST_ScaleY(rast) AS pixel_size_y,
    ST_NumBands(rast) AS bands,
    ST_SRID(rast) AS srid
FROM dem LIMIT 5;

-- 获取指定点的值
SELECT ST_Value(rast, ST_SetSRID(ST_MakePoint(-73.985, 40.748), 4326)) AS elevation
FROM dem
WHERE ST_Intersects(rast, ST_SetSRID(ST_MakePoint(-73.985, 40.748), 4326));

-- 获取列/行位置的值（波段 1）
SELECT ST_Value(rast, 1, 10, 20) FROM dem WHERE rid = 1;

SELECT (ST_SummaryStats(rast)).*
FROM dem WHERE rid = 1;
-- 返回：count、sum、mean、stddev、min、max

-- 将栅格裁剪到多边形边界
SELECT ST_Clip(rast, geom) AS clipped_rast
FROM dem, boundaries
WHERE ST_Intersects(rast, geom);

-- 单栅格地图代数：高程分类
SELECT ST_MapAlgebra(rast, 1, NULL,
    'CASE WHEN [rast] > 100 THEN 1 WHEN [rast] > 50 THEN 2 ELSE 3 END') AS classified
FROM dem;

-- 双栅格地图代数：两个 DEM 的差值
SELECT ST_MapAlgebra(a.rast, 1, b.rast, 1, '[rast1] - [rast2]') AS diff
FROM dem_old a, dem_new b
WHERE ST_Intersects(a.rast, b.rast);

-- 将栅格像素转换为矢量点
SELECT (ST_PixelAsPoints(rast)).*
FROM dem WHERE rid = 1;

-- 将栅格转换为多边形（每个唯一值一个）
SELECT (ST_DumpAsPolygons(rast)).*
FROM dem WHERE rid = 1;

-- 栅格与矢量相交并获取值
SELECT p.name, ST_Value(d.rast, p.geom) AS elevation
FROM dem d, points p
WHERE ST_Intersects(d.rast, p.geom);

-- 重采样为不同的像素大小
SELECT ST_Rescale(rast, 0.001, -0.001) FROM dem;

-- 重投影到不同的 SRID
SELECT ST_Transform(rast, 3857) FROM dem;

-- 导出为 GeoTIFF（二进制）
SELECT ST_AsTIFF(rast) FROM dem WHERE rid = 1;

-- 导出为 PNG
SELECT ST_AsPNG(rast) FROM dem WHERE rid = 1;

-- 导出为 JPEG
SELECT ST_AsJPEG(rast) FROM dem WHERE rid = 1;