PostGIS实现栅格数据地理配准【ST_SetGeoReference】等六大函数

全部文章内容请转公众号【时空智能计算】及交流群867520770，原创不易，求关注支持，更多开源GIS相关知识技能分享，免费提供学习问答交流。

本文全面详细介绍PostGIS处理栅格数据的地理配准相关函数，包括栅格数据的坐标系、起点、像素尺度、旋转倾斜等参数设置，ST_SetGeoReference、ST_SetRotation、ST_SetSRID、ST_SetScale、ST_SetSkew、ST_SetUpperLeft等六个核心函数。

一、ST_SetGeoReference

ST_SetGeoReference 是 PostGIS 中用于配置栅格数据 6 个核心地理配准参数的函数，决定了栅格在地理空间中的位置、像素尺度、旋转、倾斜等属性。

(一)核心原理

— 变体1：文本格式（支持GDAL/ESRI格式）
raster ST_SetGeoReference(raster rast, text georefcoords, text format=GDAL);

— 变体2：数值参数（精准指定6个参数，2.1.0+新增）
raster ST_SetGeoReference(raster rast, double precision upperleftx, double precision upperlefty, double precision scalex, double precision scaley, double precision skewx, double precision skewy);

栅格的地理配准本质是通过 6 个参数建立 “像素坐标” 与 “地理坐标” 的映射关系，这 6 个参数是：

ST_SetGeoReference 的核心作用是批量设置这 6 个参数，替代逐个调用ST_SetUpperLeftX、ST_SetScaleX等单参数函数，提升效率。

(二)关键参数

两种格式的参数顺序完全相同，但upperleftx/upperlefty的计算逻辑不同（这是最容易错误的点）：

示例对比：

若传入文本’10 0 0 -10 5.1 -4.9’（ESRI 格式）：

• scalex=10，scaley=-10
• 实际 upperleftx = 5.1 – 10×0.5 = 0.1
• 实际 upperlefty = -4.9 – (-10)×0.5 = 0.1
• 最终效果与 GDAL 格式’10 0 0 -10 0.1 0.1’一致（示例中因浮点精度显示为 0.0999999999999996）。

(三)使用场景

当栅格的空间位置偏移时，通过该函数批量调整 6 个参数，修正地理坐标映射：

— 校正遥感影像的地理配准偏差（GDAL格式）
UPDATE satellite_img
SET rast = ST_SetGeoReference(rast, '10 0 0 -10 100.5 200.8', 'GDAL')
WHERE rid = 1;

将不同来源、不同尺度的栅格统一为相同的地理配准参数，便于拼接 / 叠加分析：

— 批量设置所有栅格为相同的地理配准（数值参数变体）
UPDATE raster_collection
SET rast = ST_SetGeoReference(rast, 100, 200, 5, –5, 0, 0)
WHERE srid = 4326;

从 GDAL/ESRI 格式的栅格文件导入数据后，按源文件的配准参数还原空间位置：

— 从ESRI格式的文本参数还原地理配准
SELECT ST_SetGeoReference(rast, '5 0 0 -5 105.2 305.8', 'ESRI') AS registered_rast
FROM imported_raster;

(四)注意事项

若栅格包含数据库外存储的波段（如外部文件），修改地理参考可能导致无法访问外部数据。解决方案：修改前先通过ST_AddBand将外部波段转为数据库内存储，或确认外部文件的地理参考与新参数兼容。

• 不要混淆 GDAL/ESRI 格式的坐标逻辑，尤其是 ESRI 格式的upperleftx/upperlefty是像素中心坐标，需提前转换；
• 优先使用数值参数变体，避免文本格式的解析错误，且参数含义更直观。

• 大部分场景下 skewx/skewy 为 0（栅格无旋转），仅在处理倾斜 / 旋转栅格时设置；
• 非 0 倾斜系数会导致栅格像素不再是轴对齐，可能增加后续空间分析的计算开销。

• 示例中 rid=2 的结果出现0.0999999999999996而非精确 0.1，是浮点运算的正常现象，可通过ROUND函数修正显示：

SELECT rid, ROUND((ST_Metadata(rast)).upperleftx, 1) AS upperleftx FROM foo;

(五)最佳实践

数值参数无需解析文本，避免格式错误，且参数顺序明确（upperleftx→upperlefty→scalex→scaley→skewx→skewy）：

— 推荐：数值参数（直观、无解析错误）
ST_SetGeoReference(rast, 0.1, 0.1, 10, –10, 0, 0)

— 不推荐：文本格式（需记忆参数顺序，易出错）
ST_SetGeoReference(rast, '10 0 0 -10 0.1 0.1', 'GDAL')

修改后通过ST_Metadata/ST_GeoReference验证参数是否正确：

— 验证地理配准参数
SELECT
ST_GeoReference(rast) AS georef_text, — 返回文本格式的6参数
(ST_Metadata(rast)).scalex,
(ST_Metadata(rast)).upperleftx
FROM raster_table;

地理配准后需设置正确的 SRID，否则栅格的地理坐标无空间参考意义：

— 先设置地理配准，再设置SRID（WGS84）
UPDATE raster_table
SET rast = ST_SetSRID(ST_SetGeoReference(rast, 0.1, 0.1, 10, –10, 0, 0), 4326);

二、ST_SetRotation

1. 函数签名

raster ST_SetRotation(raster rast, float8 rotation);

2. 函数描述

均匀旋转栅格，旋转以弧度为单位。

三、ST_SetSRID

1. 函数签名

raster ST_SetSRID(raster rast, integer srid);

2. 函数描述

将栅格上的 SRID 设置为特定整数值。

注意：此函数不会以任何方式转换栅格 – 它只是设置定义其当前所在坐标参考系统的空间参考的元数据。对于以后的转换很有用。

四、ST_SetScale

1. 函数签名

raster ST_SetScale(raster rast, float8 xy);
raster ST_SetScale(raster rast, float8 x, float8 y);

2. 函数描述

以坐标参考系为单位设置像素的 X 和 Y 大小。 数字单位/像素宽度/高度。 如果只传入一个单位，则假设 X 和 Y 是相同的数字。

注意：ST_SetScale 与ST_Rescale 的不同之处在于 ST_SetScale 不会对栅格重新采样以匹配栅格范围。 它仅更改栅格的元数据（或地理参考）以纠正最初错误指定的缩放比例。 ST_Rescale 会生成具有不同宽度和高度的栅格，计算结果以适合输入栅格的地理范围。 ST_SetScale 不修改栅格的宽度或高度。

五、ST_SetSkew

1. 函数签名

raster ST_SetSkew(raster rast, float8 skewxy);

raster ST_SetSkew(raster rast, float8 skewx, float8 skewy);

2. 函数描述

设置地理参考 X 和 Y 倾斜（或旋转参数）。 如果只传入一个，则将X和Y设置为相同的值。

六、ST_SetUpperLeft

1. 函数签名

raster ST_SetUpperLeft(raster rast, double precision x, double precision y);

2. 函数描述

将栅格像素左上角的值设置为投影的 X 和 Y 坐标。