第05章-几何运算与叠加分析

第05章:几何运算与叠加分析

5.1 几何运算概述

几何运算是 GIS 空间分析的核心功能,NetTopologySuite 提供了完整的几何运算支持,包括布尔运算、构造运算和变换运算。

5.1.1 运算分类

分类 运算类型 说明
布尔运算 Union, Intersection, Difference, SymmetricDifference 集合操作
构造运算 Buffer, ConvexHull, Envelope 生成新几何
变换运算 Simplify, Densify, Normalize 几何变换
分析运算 Centroid, InteriorPoint, Boundary 几何属性

5.1.2 准备测试数据

using NetTopologySuite.Geometries;
using NetTopologySuite.IO;

var factory = new GeometryFactory();
var reader = new WKTReader(factory);

// 创建两个重叠的多边形
var polygonA = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});

var polygonB = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(5, 5), new Coordinate(15, 5),
    new Coordinate(15, 15), new Coordinate(5, 15), new Coordinate(5, 5)
});

Console.WriteLine($"多边形A面积: {polygonA.Area}");  // 100
Console.WriteLine($"多边形B面积: {polygonB.Area}");  // 100

5.2 布尔运算

5.2.1 Union - 并集

Union 运算返回两个几何的所有区域:

// 两个多边形的并集
var union = polygonA.Union(polygonB);

Console.WriteLine($"并集类型: {union.GeometryType}");
Console.WriteLine($"并集面积: {union.Area}");       // 175 (100 + 100 - 25重叠)
Console.WriteLine($"并集WKT: {union.AsText()}");

// 多个几何的并集
var polygonC = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(20, 0), new Coordinate(30, 0),
    new Coordinate(30, 10), new Coordinate(20, 10), new Coordinate(20, 0)
});

// 链式并集
var unionAll = polygonA.Union(polygonB).Union(polygonC);
Console.WriteLine($"三个多边形并集面积: {unionAll.Area}");  // 275

// 使用 GeometryCollection 的 Union 方法
var collection = factory.CreateGeometryCollection(new[] { polygonA, polygonB, polygonC });
var collectionUnion = collection.Union();
Console.WriteLine($"集合并集面积: {collectionUnion.Area}");

// 高效的多几何并集(使用 UnaryUnionOp)
using NetTopologySuite.Operation.Union;

var geometries = new List<Geometry> { polygonA, polygonB, polygonC };
var efficientUnion = UnaryUnionOp.Union(geometries);
Console.WriteLine($"高效并集面积: {efficientUnion.Area}");

5.2.2 Intersection - 交集

Intersection 运算返回两个几何的公共区域:

// 两个多边形的交集
var intersection = polygonA.Intersection(polygonB);

Console.WriteLine($"交集类型: {intersection.GeometryType}");  // Polygon
Console.WriteLine($"交集面积: {intersection.Area}");          // 25
Console.WriteLine($"交集WKT: {intersection.AsText()}");

// 不相交几何的交集
var nonIntersecting = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(100, 100), new Coordinate(110, 100),
    new Coordinate(110, 110), new Coordinate(100, 110), new Coordinate(100, 100)
});
var emptyIntersection = polygonA.Intersection(nonIntersecting);
Console.WriteLine($"空交集: {emptyIntersection.IsEmpty}");  // true

// 线与多边形的交集
var line = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(-5, 5), new Coordinate(15, 5)
});
var lineIntersection = polygonA.Intersection(line);
Console.WriteLine($"线与多边形交集类型: {lineIntersection.GeometryType}");  // LineString
Console.WriteLine($"交集长度: {lineIntersection.Length}");  // 10

5.2.3 Difference - 差集

Difference 运算返回在第一个几何中但不在第二个几何中的区域:

// A - B: 在A中但不在B中的部分
var differenceAB = polygonA.Difference(polygonB);
Console.WriteLine($"A-B 面积: {differenceAB.Area}");  // 75

// B - A: 在B中但不在A中的部分
var differenceBA = polygonB.Difference(polygonA);
Console.WriteLine($"B-A 面积: {differenceBA.Area}");  // 75

// 验证关系:A = (A ∩ B) ∪ (A - B)
var reconstructedA = intersection.Union(differenceAB);
Console.WriteLine($"重构A面积: {reconstructedA.Area}");  // 100
Console.WriteLine($"等于原始A: {reconstructedA.EqualsTopologically(polygonA)}");  // true

// 多个几何的差集
var multiDiff = polygonA.Difference(polygonB).Difference(polygonC);
Console.WriteLine($"多重差集面积: {multiDiff.Area}");

5.2.4 SymmetricDifference - 对称差集

SymmetricDifference 返回在任一几何中但不在两者交集中的区域:

// 对称差集:(A ∪ B) - (A ∩ B)
var symDiff = polygonA.SymmetricDifference(polygonB);

Console.WriteLine($"对称差集类型: {symDiff.GeometryType}");
Console.WriteLine($"对称差集面积: {symDiff.Area}");  // 150 (75 + 75)

// 验证:对称差集 = (A - B) ∪ (B - A)
var manualSymDiff = differenceAB.Union(differenceBA);
Console.WriteLine($"手动计算对称差集面积: {manualSymDiff.Area}");  // 150

// 对称差集的另一种理解:(A ∪ B) - (A ∩ B)
var unionMinusIntersection = union.Difference(intersection);
Console.WriteLine($"并集减交集面积: {unionMinusIntersection.Area}");  // 150

5.2.5 布尔运算应用示例

public class OverlayAnalysis
{
    private readonly GeometryFactory _factory;
    
    public OverlayAnalysis()
    {
        _factory = new GeometryFactory();
    }
    
    /// <summary>
    /// 计算用地叠加分析
    /// </summary>
    public OverlayResult AnalyzeLandOverlay(Polygon parcel, Polygon buildingArea)
    {
        var result = new OverlayResult();
        
        // 可建设区域
        result.BuildableArea = parcel.Intersection(buildingArea);
        result.BuildableAreaSize = result.BuildableArea.Area;
        
        // 不可建设区域
        result.NonBuildableArea = parcel.Difference(buildingArea);
        result.NonBuildableAreaSize = result.NonBuildableArea.Area;
        
        // 建设覆盖率
        result.CoverageRatio = result.BuildableAreaSize / parcel.Area * 100;
        
        return result;
    }
    
    /// <summary>
    /// 合并相邻地块
    /// </summary>
    public Geometry MergeAdjacentParcels(IEnumerable<Polygon> parcels)
    {
        return UnaryUnionOp.Union(parcels.Cast<Geometry>().ToList());
    }
    
    /// <summary>
    /// 计算道路占地
    /// </summary>
    public RoadOccupancyResult CalculateRoadOccupancy(
        Polygon parcel, 
        LineString road, 
        double roadWidth)
    {
        var result = new RoadOccupancyResult();
        
        // 创建道路缓冲区
        var roadBuffer = road.Buffer(roadWidth / 2);
        
        // 道路占用区域
        result.OccupiedArea = parcel.Intersection(roadBuffer);
        result.OccupiedAreaSize = result.OccupiedArea.Area;
        
        // 剩余地块
        result.RemainingParcel = parcel.Difference(roadBuffer);
        result.RemainingAreaSize = result.RemainingParcel.Area;
        
        // 占用比例
        result.OccupancyRatio = result.OccupiedAreaSize / parcel.Area * 100;
        
        return result;
    }
}

public class OverlayResult
{
    public Geometry BuildableArea { get; set; }
    public double BuildableAreaSize { get; set; }
    public Geometry NonBuildableArea { get; set; }
    public double NonBuildableAreaSize { get; set; }
    public double CoverageRatio { get; set; }
}

public class RoadOccupancyResult
{
    public Geometry OccupiedArea { get; set; }
    public double OccupiedAreaSize { get; set; }
    public Geometry RemainingParcel { get; set; }
    public double RemainingAreaSize { get; set; }
    public double OccupancyRatio { get; set; }
}

5.3 缓冲区运算

5.3.1 基本缓冲区

var factory = new GeometryFactory();

// 点缓冲区
var point = factory.CreatePoint(new Coordinate(0, 0));
var pointBuffer = point.Buffer(10);
Console.WriteLine($"点缓冲区类型: {pointBuffer.GeometryType}");  // Polygon
Console.WriteLine($"点缓冲区面积: {pointBuffer.Area:F2}");       // ≈314.16 (π * 10²)

// 线缓冲区
var line = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(100, 0)
});
var lineBuffer = line.Buffer(5);
Console.WriteLine($"线缓冲区面积: {lineBuffer.Area:F2}");

// 多边形缓冲区
var polygon = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});
var polygonBuffer = polygon.Buffer(5);
Console.WriteLine($"多边形缓冲区面积: {polygonBuffer.Area:F2}");

// 负缓冲区(内缩)
var innerBuffer = polygon.Buffer(-2);
Console.WriteLine($"内缩缓冲区面积: {innerBuffer.Area}");  // 36 (6×6)

5.3.2 缓冲区样式

using NetTopologySuite.Operation.Buffer;

var factory = new GeometryFactory();
var line = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(50, 0), new Coordinate(50, 50)
});

// 圆形端点(默认)
var roundBuffer = line.Buffer(10);
Console.WriteLine($"圆形端点缓冲区: {roundBuffer.NumPoints} 个点");

// 平头端点
var bufferParams = new BufferParameters
{
    EndCapStyle = EndCapStyle.Flat
};
var flatBuffer = BufferOp.Buffer(line, 10, bufferParams);
Console.WriteLine($"平头端点缓冲区: {flatBuffer.NumPoints} 个点");

// 方头端点
bufferParams.EndCapStyle = EndCapStyle.Square;
var squareBuffer = BufferOp.Buffer(line, 10, bufferParams);
Console.WriteLine($"方头端点缓冲区: {squareBuffer.NumPoints} 个点");

// 设置圆弧分段数(控制平滑度)
bufferParams.QuadrantSegments = 8;  // 默认是8
bufferParams.EndCapStyle = EndCapStyle.Round;
var smoothBuffer = BufferOp.Buffer(line, 10, bufferParams);
Console.WriteLine($"高精度缓冲区: {smoothBuffer.NumPoints} 个点");

// 低精度缓冲区(更快,更少点)
bufferParams.QuadrantSegments = 2;
var roughBuffer = BufferOp.Buffer(line, 10, bufferParams);
Console.WriteLine($"低精度缓冲区: {roughBuffer.NumPoints} 个点");

5.3.3 连接样式

var factory = new GeometryFactory();
var angledLine = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(50, 50), new Coordinate(100, 0)
});

// 圆角连接(默认)
var roundJoin = new BufferParameters { JoinStyle = JoinStyle.Round };
var roundBuffer = BufferOp.Buffer(angledLine, 10, roundJoin);

// 斜角连接
var bevelJoin = new BufferParameters { JoinStyle = JoinStyle.Bevel };
var bevelBuffer = BufferOp.Buffer(angledLine, 10, bevelJoin);

// 尖角连接
var mitreJoin = new BufferParameters 
{ 
    JoinStyle = JoinStyle.Mitre,
    MitreLimit = 5.0  // 限制尖角的最大长度
};
var mitreBuffer = BufferOp.Buffer(angledLine, 10, mitreJoin);

Console.WriteLine($"圆角连接点数: {roundBuffer.NumPoints}");
Console.WriteLine($"斜角连接点数: {bevelBuffer.NumPoints}");
Console.WriteLine($"尖角连接点数: {mitreBuffer.NumPoints}");

5.3.4 单边缓冲区

using NetTopologySuite.Operation.Buffer;

var factory = new GeometryFactory();
var line = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(100, 0)
});

// 左侧缓冲区
var leftParams = new BufferParameters { IsSingleSided = true };
var leftBuffer = BufferOp.Buffer(line, 10, leftParams);

// 右侧缓冲区(使用负距离)
var rightBuffer = BufferOp.Buffer(line, -10, leftParams);

Console.WriteLine($"左侧缓冲区: {leftBuffer.AsText()}");
Console.WriteLine($"右侧缓冲区: {rightBuffer.AsText()}");

5.3.5 可变宽度缓冲区

using NetTopologySuite.Operation.Buffer;

var factory = new GeometryFactory();

// 创建线
var line = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0),
    new Coordinate(50, 0),
    new Coordinate(100, 0)
});

// 可变宽度缓冲区(需要使用 OffsetCurve)
// 注意:NetTopologySuite 2.x 支持偏移曲线
var offsetCurve = line.Buffer(5);  // 基本缓冲区

// 模拟可变宽度:分段处理
var segment1 = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(50, 0)
});
var segment2 = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(50, 0), new Coordinate(100, 0)
});

var buffer1 = segment1.Buffer(5);   // 宽度5
var buffer2 = segment2.Buffer(10);  // 宽度10

var variableBuffer = buffer1.Union(buffer2);
Console.WriteLine($"可变宽度缓冲区面积: {variableBuffer.Area}");

5.4 构造运算

5.4.1 ConvexHull - 凸包

凸包是包含所有点的最小凸多边形:

var factory = new GeometryFactory();

// 从多点计算凸包
var points = factory.CreateMultiPoint(new Point[]
{
    factory.CreatePoint(new Coordinate(0, 0)),
    factory.CreatePoint(new Coordinate(10, 0)),
    factory.CreatePoint(new Coordinate(5, 5)),
    factory.CreatePoint(new Coordinate(10, 10)),
    factory.CreatePoint(new Coordinate(0, 10)),
    factory.CreatePoint(new Coordinate(3, 3))  // 内部点
});

var convexHull = points.ConvexHull();
Console.WriteLine($"凸包类型: {convexHull.GeometryType}");
Console.WriteLine($"凸包点数: {convexHull.NumPoints}");  // 5(不包含内部点)
Console.WriteLine($"凸包: {convexHull.AsText()}");

// 从凹多边形计算凸包
var concavePolygon = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(5, 5),  // 凹点
    new Coordinate(0, 10), new Coordinate(0, 0)
});

var polygonHull = concavePolygon.ConvexHull();
Console.WriteLine($"原多边形面积: {concavePolygon.Area}");
Console.WriteLine($"凸包面积: {polygonHull.Area}");

// 凸包用于简化边界
var complexLine = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(5, 10),
    new Coordinate(10, 0), new Coordinate(15, 10),
    new Coordinate(20, 0)
});
var lineHull = complexLine.ConvexHull();
Console.WriteLine($"线的凸包: {lineHull.AsText()}");

5.4.2 Envelope - 外包矩形

var factory = new GeometryFactory();

var polygon = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(1, 2), new Coordinate(10, 3),
    new Coordinate(8, 15), new Coordinate(3, 10), new Coordinate(1, 2)
});

// 获取外包矩形
var envelope = polygon.Envelope;
Console.WriteLine($"外包矩形: {envelope.AsText()}");

// 获取外包矩形信息
var envelopeInternal = polygon.EnvelopeInternal;
Console.WriteLine($"最小X: {envelopeInternal.MinX}");
Console.WriteLine($"最大X: {envelopeInternal.MaxX}");
Console.WriteLine($"最小Y: {envelopeInternal.MinY}");
Console.WriteLine($"最大Y: {envelopeInternal.MaxY}");
Console.WriteLine($"宽度: {envelopeInternal.Width}");
Console.WriteLine($"高度: {envelopeInternal.Height}");
Console.WriteLine($"面积: {envelopeInternal.Area}");
Console.WriteLine($"中心: ({envelopeInternal.Centre.X}, {envelopeInternal.Centre.Y})");

// 扩展外包矩形
var expandedEnvelope = new Envelope(envelopeInternal);
expandedEnvelope.ExpandBy(5);  // 各方向扩展5
var expandedRect = factory.ToGeometry(expandedEnvelope);
Console.WriteLine($"扩展后外包矩形: {expandedRect.AsText()}");

5.4.3 Centroid 和 InteriorPoint

var factory = new GeometryFactory();

// 简单多边形的质心
var square = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});
var centroid = square.Centroid;
Console.WriteLine($"正方形质心: {centroid.AsText()}");  // POINT (5 5)

// 凹多边形的质心(可能在多边形外部)
var concave = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(100, 0),
    new Coordinate(100, 100), new Coordinate(50, 10),  // 深凹
    new Coordinate(0, 100), new Coordinate(0, 0)
});
var concaveCentroid = concave.Centroid;
var concaveInterior = concave.InteriorPoint;

Console.WriteLine($"凹多边形质心: {concaveCentroid.AsText()}");
Console.WriteLine($"质心在多边形内: {concave.Contains(concaveCentroid)}");

Console.WriteLine($"凹多边形内部点: {concaveInterior.AsText()}");
Console.WriteLine($"内部点在多边形内: {concave.Contains(concaveInterior)}");  // 总是 true

// 带洞多边形
var shell = factory.CreateLinearRing(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(100, 0),
    new Coordinate(100, 100), new Coordinate(0, 100), new Coordinate(0, 0)
});
var hole = factory.CreateLinearRing(new Coordinate[]
{
    new Coordinate(40, 40), new Coordinate(60, 40),
    new Coordinate(60, 60), new Coordinate(40, 60), new Coordinate(40, 40)
});
var donut = factory.CreatePolygon(shell, new[] { hole });

var donutCentroid = donut.Centroid;
var donutInterior = donut.InteriorPoint;

Console.WriteLine($"环形质心: {donutCentroid.AsText()}");
Console.WriteLine($"质心在环形内: {donut.Contains(donutCentroid)}");  // false(在洞里)

Console.WriteLine($"环形内部点: {donutInterior.AsText()}");
Console.WriteLine($"内部点在环形内: {donut.Contains(donutInterior)}");  // true

5.4.4 Boundary - 边界

var factory = new GeometryFactory();

// 点的边界(空)
var point = factory.CreatePoint(new Coordinate(5, 5));
Console.WriteLine($"点的边界: {point.Boundary.AsText()}");  // GEOMETRYCOLLECTION EMPTY

// 线的边界(端点)
var line = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 10)
});
Console.WriteLine($"线的边界: {line.Boundary.AsText()}");  // MULTIPOINT ((0 0), (10 10))

// 闭合线的边界(空)
var closedLine = factory.CreateLinearRing(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 0)
});
Console.WriteLine($"闭合线的边界: {closedLine.Boundary.AsText()}");  // MULTIPOINT EMPTY

// 多边形的边界(外环和内环)
var shell = factory.CreateLinearRing(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(20, 0),
    new Coordinate(20, 20), new Coordinate(0, 20), new Coordinate(0, 0)
});
var hole = factory.CreateLinearRing(new Coordinate[]
{
    new Coordinate(5, 5), new Coordinate(5, 15),
    new Coordinate(15, 15), new Coordinate(15, 5), new Coordinate(5, 5)
});
var polygon = factory.CreatePolygon(shell, new[] { hole });

Console.WriteLine($"带洞多边形边界类型: {polygon.Boundary.GeometryType}");
Console.WriteLine($"带洞多边形边界: {polygon.Boundary.AsText()}");

5.5 几何简化

5.5.1 Douglas-Peucker 简化

using NetTopologySuite.Simplify;

var factory = new GeometryFactory();

// 创建复杂线
var complexLine = factory.CreateLineString(new Coordinate[]
{
    new Coordinate(0, 0),
    new Coordinate(1, 0.1),
    new Coordinate(2, 0),
    new Coordinate(3, 0.1),
    new Coordinate(4, 0),
    new Coordinate(5, 0.1),
    new Coordinate(6, 0),
    new Coordinate(7, 0.1),
    new Coordinate(8, 0),
    new Coordinate(9, 0.1),
    new Coordinate(10, 0)
});

Console.WriteLine($"原始线点数: {complexLine.NumPoints}");

// Douglas-Peucker 简化
var dpSimplifier = new DouglasPeuckerSimplifier(complexLine);
dpSimplifier.DistanceTolerance = 0.2;  // 容差
var simplified = dpSimplifier.GetResultGeometry();

Console.WriteLine($"简化后点数: {simplified.NumPoints}");
Console.WriteLine($"简化后: {simplified.AsText()}");

// 使用便捷方法
var simplified2 = complexLine.Simplify(0.2);
Console.WriteLine($"便捷方法简化后点数: {simplified2.NumPoints}");

5.5.2 拓扑保持简化

using NetTopologySuite.Simplify;

var factory = new GeometryFactory();

// 创建复杂多边形
var complexPolygon = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0),
    new Coordinate(5, 0.1),
    new Coordinate(10, 0),
    new Coordinate(10.1, 5),
    new Coordinate(10, 10),
    new Coordinate(5, 10.1),
    new Coordinate(0, 10),
    new Coordinate(-0.1, 5),
    new Coordinate(0, 0)
});

Console.WriteLine($"原始多边形点数: {complexPolygon.NumPoints}");
Console.WriteLine($"原始多边形有效: {complexPolygon.IsValid}");

// 拓扑保持简化
var tpSimplifier = new TopologyPreservingSimplifier(complexPolygon);
tpSimplifier.DistanceTolerance = 0.5;
var tpSimplified = tpSimplifier.GetResultGeometry();

Console.WriteLine($"拓扑保持简化后点数: {tpSimplified.NumPoints}");
Console.WriteLine($"简化后有效: {tpSimplified.IsValid}");

// 对比:Douglas-Peucker 简化
var dpSimplified = complexPolygon.Simplify(0.5);
Console.WriteLine($"DP简化后点数: {dpSimplified.NumPoints}");
Console.WriteLine($"DP简化后有效: {dpSimplified.IsValid}");

5.5.3 VW 简化算法

using NetTopologySuite.Simplify;

var factory = new GeometryFactory();

// 创建复杂线
var complexLine = factory.CreateLineString(Enumerable
    .Range(0, 100)
    .Select(i => new Coordinate(i, Math.Sin(i * 0.1) * 5))
    .ToArray());

Console.WriteLine($"原始线点数: {complexLine.NumPoints}");

// VW (Visvalingam-Whyatt) 简化
var vwSimplifier = new VWSimplifier(complexLine);
vwSimplifier.DistanceTolerance = 1.0;
var vwSimplified = vwSimplifier.GetResultGeometry();

Console.WriteLine($"VW简化后点数: {vwSimplified.NumPoints}");

5.6 几何验证与修复

5.6.1 验证几何

using NetTopologySuite.Operation.Valid;

var factory = new GeometryFactory();

// 创建有效多边形
var validPolygon = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});

Console.WriteLine($"有效多边形 IsValid: {validPolygon.IsValid}");

// 创建自相交多边形(无效)
var selfIntersecting = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 10),
    new Coordinate(10, 0), new Coordinate(0, 10), new Coordinate(0, 0)
});

Console.WriteLine($"自相交多边形 IsValid: {selfIntersecting.IsValid}");

// 获取详细验证信息
var validator = new IsValidOp(selfIntersecting);
var validationError = validator.ValidationError;

if (validationError != null)
{
    Console.WriteLine($"错误类型: {validationError.ErrorType}");
    Console.WriteLine($"错误位置: ({validationError.Coordinate.X}, {validationError.Coordinate.Y})");
    Console.WriteLine($"错误消息: {validationError.Message}");
}

5.6.2 修复几何

using NetTopologySuite.Geometries.Utilities;

var factory = new GeometryFactory();

// 自相交多边形
var invalid = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 10),
    new Coordinate(10, 0), new Coordinate(0, 10), new Coordinate(0, 0)
});

Console.WriteLine($"修复前有效: {invalid.IsValid}");

// 方法1:使用 Buffer(0) 修复
var fixed1 = invalid.Buffer(0);
Console.WriteLine($"Buffer(0)修复后有效: {fixed1.IsValid}");
Console.WriteLine($"Buffer(0)修复后类型: {fixed1.GeometryType}");

// 方法2:使用 GeometryFixer(NTS 2.5+)
var fixer = new GeometryFixer(invalid);
var fixed2 = fixer.GetResult();
Console.WriteLine($"GeometryFixer修复后有效: {fixed2.IsValid}");

// 方法3:自定义修复(针对特定问题)
public static Geometry FixGeometry(Geometry geometry)
{
    if (geometry.IsValid) return geometry;
    
    // 尝试 Buffer(0) 修复
    var fixed = geometry.Buffer(0);
    if (fixed.IsValid && !fixed.IsEmpty)
    {
        return fixed;
    }
    
    // 如果是多边形,尝试使用 GeometryFixer
    if (geometry is Polygon || geometry is MultiPolygon)
    {
        var fixer = new GeometryFixer(geometry);
        fixed = fixer.GetResult();
        if (fixed.IsValid && !fixed.IsEmpty)
        {
            return fixed;
        }
    }
    
    throw new InvalidOperationException("无法修复几何");
}

5.6.3 常见无效几何类型

var factory = new GeometryFactory();

// 1. 自相交多边形
var selfIntersect = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 10),
    new Coordinate(10, 0), new Coordinate(0, 10), new Coordinate(0, 0)
});
Console.WriteLine($"自相交: {new IsValidOp(selfIntersect).ValidationError?.ErrorType}");

// 2. 孔在外环外部
var shell = factory.CreateLinearRing(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});
var holeOutside = factory.CreateLinearRing(new Coordinate[]
{
    new Coordinate(20, 0), new Coordinate(30, 0),
    new Coordinate(30, 10), new Coordinate(20, 10), new Coordinate(20, 0)
});
var invalidHole = factory.CreatePolygon(shell, new[] { holeOutside });
Console.WriteLine($"孔在外部: {new IsValidOp(invalidHole).ValidationError?.ErrorType}");

// 3. 重复点
var duplicatePoints = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0), new Coordinate(10, 0),  // 重复
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});
// 注意:重复点通常不会导致几何无效,但可能影响某些操作

// 4. 未闭合环
// LinearRing 会自动闭合,但手动创建可能出问题

5.7 几何变换

5.7.1 仿射变换

using NetTopologySuite.Geometries.Utilities;

var factory = new GeometryFactory();

var polygon = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});

// 创建仿射变换
var transform = new AffineTransformation();

// 平移
transform.Translate(5, 5);
var translated = transform.Transform(polygon);
Console.WriteLine($"平移后: {translated.AsText()}");

// 旋转(围绕原点)
transform = new AffineTransformation();
transform.Rotate(Math.PI / 4);  // 45度
var rotated = transform.Transform(polygon);
Console.WriteLine($"旋转后: {rotated.AsText()}");

// 围绕中心点旋转
var centroid = polygon.Centroid;
transform = new AffineTransformation();
transform.Translate(-centroid.X, -centroid.Y);  // 移到原点
transform.Rotate(Math.PI / 4);                   // 旋转
transform.Translate(centroid.X, centroid.Y);    // 移回
var rotatedAroundCenter = transform.Transform(polygon);
Console.WriteLine($"围绕中心旋转后: {rotatedAroundCenter.AsText()}");

// 缩放
transform = new AffineTransformation();
transform.Scale(2, 2);  // 放大2倍
var scaled = transform.Transform(polygon);
Console.WriteLine($"缩放后面积: {scaled.Area}");  // 400

// 组合变换
transform = new AffineTransformation();
transform.Scale(2, 2);
transform.Rotate(Math.PI / 4);
transform.Translate(100, 100);
var combined = transform.Transform(polygon);
Console.WriteLine($"组合变换后: {combined.AsText()}");

5.7.2 坐标变换

using NetTopologySuite.Geometries.Utilities;

// 自定义坐标变换
public class FlipYTransform : ICoordinateSequenceFilter
{
    public bool Done => false;
    public bool GeometryChanged => true;
    
    public void Filter(CoordinateSequence seq, int i)
    {
        var y = seq.GetY(i);
        seq.SetY(i, -y);  // Y轴翻转
    }
}

var factory = new GeometryFactory();
var polygon = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});

// 创建副本并应用变换
var flipped = (Polygon)polygon.Copy();
flipped.Apply(new FlipYTransform());

Console.WriteLine($"原始: {polygon.AsText()}");
Console.WriteLine($"Y翻转后: {flipped.AsText()}");

5.7.3 Normalize - 标准化

var factory = new GeometryFactory();

// 创建两个拓扑相同但坐标顺序不同的多边形
var poly1 = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(0, 0), new Coordinate(10, 0),
    new Coordinate(10, 10), new Coordinate(0, 10), new Coordinate(0, 0)
});

var poly2 = factory.CreatePolygon(new Coordinate[]
{
    new Coordinate(10, 0), new Coordinate(10, 10),
    new Coordinate(0, 10), new Coordinate(0, 0), new Coordinate(10, 0)
});

Console.WriteLine($"Equals: {poly1.Equals(poly2)}");  // true(拓扑相等)
Console.WriteLine($"EqualsExact: {poly1.EqualsExact(poly2)}");  // false

// 标准化后比较
var norm1 = (Geometry)poly1.Copy();
norm1.Normalize();

var norm2 = (Geometry)poly2.Copy();
norm2.Normalize();

Console.WriteLine($"标准化后EqualsExact: {norm1.EqualsExact(norm2)}");  // true
Console.WriteLine($"EqualsNormalized: {poly1.EqualsNormalized(poly2)}");  // true

5.8 本章小结

本章详细介绍了 NetTopologySuite 的几何运算与叠加分析:

  1. 布尔运算:Union、Intersection、Difference、SymmetricDifference
  2. 缓冲区运算:基本缓冲区、缓冲区样式、单边缓冲区
  3. 构造运算:ConvexHull、Envelope、Centroid、InteriorPoint、Boundary
  4. 几何简化:Douglas-Peucker、拓扑保持简化、VW 简化
  5. 几何验证与修复:IsValidOp、GeometryFixer
  6. 几何变换:仿射变换、坐标变换、标准化

5.9 下一步

下一章我们将学习空间分析算法,包括:

  • 线性参考
  • 最近点计算
  • 多边形分解
  • 三角剖分

相关资源

posted @ 2025-12-29 10:22  我才是银古  阅读(2)  评论(0)    收藏  举报