Geometry转Geography生成无效形状及面积计算异常技术问询

我有两个Polygon/MultiPolygon类型的Geometry形状（包含经纬度点），执行以下语句获取交集：

land_shape.STIntersection(grid_shape) as shape

发现部分交集形状作为Geometry时STIsValid()=1，但转换为Geography（geography::STGeomFromText(shape.STAsText(), 4326)）后STIsValid()≠1。输入的land_shape与grid_shape单独转Geography均有效，排除输入损坏问题。

通过IsValidDetailed()得到报错：

24409: Not valid because some portion of polygon ring (1) lies in the interior of a polygon.
24413: Not valid because of two overlapping edges in curve (1).

另外，即使使用MakeValid()，面积计算仍返回极大值，执行语句：

geography :: STGeomFromText(land_shape.STAsText(), 4326).STIntersection( geography :: STGeomFromText(grid_shape.STAsText(), 4326) ).STArea()

尝试ReorientObject()后问题依旧：

geography :: STGeomFromText(land_shape.STAsText(), 4326).ReorientObject().STIntersection( geography :: STGeomFromText(grid_shape.STAsText(), 4326).ReorientObject() ).STArea()

Geometry与Geography的有效性规则差异 #

• Geometry（平面坐标系）：仅校验平面拓扑规则，只要多边形无明显自交、环闭合，就判定有效；不考虑地球曲率，也不强制环的方向（顺时针/逆时针）。
• Geography（球面坐标系）：严格遵循OGC球面拓扑标准：
  • 多边形环必须符合右手定则（外环顺时针、内环逆时针），否则会被解析成覆盖地球大部分区域的“反多边形”，触发有效性错误。
  • 球面下不允许多边形环自交、边缘重叠，而Geometry的平面交集运算会忽略曲率，可能生成在平面合法但球面拓扑违规的形状。

交集运算的拓扑丢失问题 #

Geometry的交集是基于平面计算的，当涉及大跨度或接近极点的区域时，平面计算的结果转换为球面坐标后，会出现环方向反转、顶点微小偏移导致的自交/重叠，最终触发Geography的有效性校验失败。

1. 直接在Geography层面计算交集 #

跳过Geometry的交集步骤，用预先转换好的Geography对象执行运算，避免平面转球面的拓扑丢失：

-- 先预转Geography
DECLARE @land_geog geography = geography::STGeomFromText(land_shape.STAsText(), 4326);
DECLARE @grid_geog geography = geography::STGeomFromText(grid_shape.STAsText(), 4326);

-- 直接计算Geography交集
SELECT @land_geog.STIntersection(@grid_geog) as intersect_geog;

2. 手动修正多边形环方向 #

如果是环方向错误导致的面积异常：

• 用STArea()判断：若面积接近地球表面积（约5.1亿平方公里），说明环方向反转。
• 手动反转环的顶点顺序，确保多边形指向预期的小区域：

-- 反转外环顶点顺序
SELECT geography::STGeomFromText(
  REPLACE(shape.STAsText(), 'POLYGON((', 'POLYGON((' + REVERSE(STUFF(shape.STAsText(), 1, 10, '')))
, 4326) as corrected_geog;

3. 用STBuffer(0)修复轻微拓扑问题 #

STBuffer(0)能温和消除微小自交、边缘重叠，比MakeValid()更精准：

SELECT @land_geog.STBuffer(0).STIntersection(@grid_geog.STBuffer(0)) as intersect_geog;

4. 拆分MultiPolygon为单个Polygon处理 #

若输入是MultiPolygon，拆分后分别计算交集再合并，降低拓扑复杂度：

SELECT 
  geography::STGeomFromText(
    'MULTIPOLYGON(' + STRING_AGG('(' + p.STAsText() + ')', ',') + ')', 4326
  ) as intersect_geog
FROM (
  SELECT 
    geography::STGeomFromText(land_poly.STAsText(), 4326).STIntersection(geography::STGeomFromText(grid_poly.STAsText(), 4326)) as p
  FROM 
    -- 拆分land_shape的单个Polygon
    (SELECT geometry::STGeomFromText('POLYGON(' + value + ')', 4326) as land_poly 
     FROM STRING_SPLIT(REPLACE(REPLACE(land_shape.STAsText(), 'MULTIPOLYGON((', ''), '))', ''), '),(')) as land_parts,
    -- 拆分grid_shape的单个Polygon
    (SELECT geometry::STGeomFromText('POLYGON(' + value + ')', 4326) as grid_poly 
     FROM STRING_SPLIT(REPLACE(REPLACE(grid_shape.STAsText(), 'MULTIPOLYGON((', ''), '))', ''), '),(')) as grid_parts
  WHERE land_poly.STIntersects(grid_poly) = 1
) as intersect_parts;

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