课程思政

1 导入栅格数据

1.1 建立postgis_raster拓展

CREATE EXTENSION IF NOT EXISTS postgis_raster SCHEMA public;

如果上面的拓展建立成功，输入命令后

SELECT postgis_full_version();

应该得到类似下面的输出结果

POSTGIS="3.2.0 3.2.0" [EXTENSION] PGSQL="140" GEOS="3.10.1-CAPI-1.16.0" PROJ="7.2.1" GDAL="GDAL 3.3.3, released 2021/10/25" LIBXML="2.9.9" LIBJSON="0.12" LIBPROTOBUF="1.2.1" WAGYU="0.5.0 (Internal)" TOPOLOGY RASTER

启用对数据库外栅格的支持并使用所有驱动程序

ALTER DATABASE postgres SET postgis.enable_outdb_rasters = true;
ALTER DATABASE postgres SET postgis.gdal_enabled_drivers = 'ENABLE_ALL';

1.2 创建模式

创建ch07模式

CREATE SCHEMA ch07;

1.3 导入数据

利用raster2pgsql导入栅格数据

raster2pgsql -C D:\postgisdata\ch07\data\rs.tif ch07.china2 | psql -h localhost -U postgres -p 5432 -d postgres

在QGIS中可视化结果

2 利用PostGIS函数创建栅格数据

2.2 创建表

创建表

CREATE TABLE ch07.bag_o_rasters(
    rid serial primary key, rast_name text, rast raster
);

2.3 利用ST_AsRaster将几何类型转换为Raster

您可以使用 ST_AsRaster函数从几何图形转换栅格。当有以下需求时：

• 将几何图形输出到图像格式以供查看。

• 在栅格上叠加几何图形以突出特定区域或合并感兴趣的边界、道路或点。

• 存储有关您可以使用栅格分析函数族查询的区域的数值统计信息。

2.3.1 根据几何数据生成单波段的栅格数据

INSERT INTO ch07.bag_o_rasters(rast_name, rast)  
WITH a1 AS (
	SELECT ST_Buffer(
		ST_GeomFromText(
			'LINESTRING(
				448252 5414206,448289 5414317,448293 5414330,
				448324 5414417,448351 5414495
			)', 
			32631), 
		10
	) As geom
)   -- <1>
SELECT 'disprop road', ST_AsRaster(geom,50,500) FROM a1 -- <2>
UNION ALL
SELECT 'disprop road', ST_AsRaster(geom,50,500) FROM a1
	UNION ALL
	SELECT 'proport fixed w road',
	    ST_AsRaster( geom, 200, 
	        ( ST_YMax(geom) - ST_YMin(geom) ) *  200 /
	        ( ST_XMax(geom) - ST_XMin(geom)  )::integer)
	FROM a1; -- <3>

2.3.2 由坐标定位的几何

WITH 
    r AS 
        (SELECT 
            ST_MakeEmptyRaster(
                500,500,445000,5415000,2,-2,0,0,32631
            ) As rast
        ), -- <1>
    g AS 
        (SELECT ST_Buffer(
            ST_GeomFromText(
                'LINESTRING(
                    448252 5414206,448289 5414317,448293 5414330,
                    448324 5414417,448351 5414495)',
                32631), 
            10) As geom
        ) -- <2>
INSERT INTO ch07.bag_o_rasters(rast_name, rast)     
SELECT 'canvas aligned road', ST_AsRaster(geom,rast,'8BUI'::text) -- <3>
FROM r CROSS JOIN g;

2.4 利用raster2pgsql加载栅格数据

在cmd命令行中将栅格图片导入到pg数据库表中

raster2pgsql -e -R -a D:/postgisdata/ch07/data/adbadge_tall.png ch07.bag_o_rasters | psql -U postgres -d postgres -h localhost -p 5432

更新名称

UPDATE ch07.bag_o_rasters SET rast_name = 'Reo' WHERE rast_name IS NULL;

2.5 从头开始构建栅格数据

INSERT INTO ch07.bag_o_rasters(rast_name, rast) 
SELECT 
    'Raster 1 band scratch', 
    ST_AddBand( -- <1>
        ST_MakeEmptyRaster( -- <2>
            500,500,445000,5415000,2,-2,0,0,32631
        ),
        '8BUI'::text,255,0 -- <3>
    ) As rast;

设置像素值

UPDATE ch07.bag_o_rasters AS b
SET rast = ST_SetValue(rast,1,10,20,146)
WHERE b.rast_name = 'Raster 1 band scratch';

2.6 构建热力图

WITH heatmap As ( -- <1>
    SELECT array_agg( -- <2>
        (ST_Buffer(  
            ST_Translate(
                ST_SetSRID(
                    ST_Point(445500,5414500), 32631
                ),
                -500 + i * 150,
                -200 + 160 * i
            ),
            i * 50),
            50 + i * 15.0
        )::geomval 
    ) As gvals
    FROM generate_series(-3,4) As i
)
INSERT INTO ch07.bag_o_rasters(rast_name, rast)  
SELECT 
    'Raster 1 band heatmap', 
    ST_SetValues(rast,1, heatmap.gvals) As rast -- <3>
FROM ch07.bag_o_rasters As b CROSS JOIN heatmap
WHERE b.rast_name = 'Raster 1 band scratch';