说明
介绍在PotGIS中将点数据转换为栅格数据。
关键字: raster
、point
、PostGIS
环境准备
- Postgresql版本:
PostgreSQL 14.0, 64-bit
- PostGIS版本:
POSTGIS="3.3.2"
- QGIS版本:
3.28.3-Firenze
基本步骤
一、数据准备
测试数据中有一张点数据表,坐标系3857
。
CREATE TABLE IF NOT EXISTS public.test_point
(
geom geometry(Point,3857),
id bigint,
h numeric,--指定要转换为栅格值的列,如高程、温度、风速等
CONSTRAINT test_point_pkey PRIMARY KEY (id)
)
TABLESPACE pg_default;
ALTER TABLE public.test_point
OWNER to postgres;
二、二维点数据转换为三维点
DROP TABLE IF EXISTS test_point_3d;
CREATE TABLE test_point_3d
as select h,ST_Force3D(ST_Transform(geom,3857),tp.h) as geom
from test_point tp;
三、将三维点数据转换为栅格
WITH inputs AS
(SELECT 16::INTEGER AS pixelsize,'linear:radius:100' AS algorithm,
ST_Collect(geom) AS geom, ST_Expand(ST_Collect(geom), 10) AS ext
FROM test_point_3d ),sizes AS
(SELECT ceil((ST_XMax(ext) - ST_XMin(ext))/pixelsize)::integer AS width,
ceil((ST_YMax(ext) - ST_YMin(ext))/pixelsize)::integer AS height,
ST_XMin(ext) AS upperleftx,ST_YMax(ext) AS upperlefty
FROM inputs )SELECT 1 AS rid,
ST_InterpolateRaster( geom,
algorithm,
ST_SetSRID(ST_AddBand(ST_MakeEmptyRaster(width,
height,upperleftx,upperlefty,pixelsize),
'32BF'), ST_SRID(geom)) ) AS rast FROM sizes, inputs;
其中ST_InterpolateRaster 参见
:
https://postgis.net/docs/RT_ST_InterpolateRaster.html
其参数algorithm
参见:
https://gdal.org/programs/gdal_grid.html#interpolation-algorithms