第07章 - 几何处理与空间分析

7.1 GeometryUtil概述

7.1.1 功能分类

GeometryUtil 是OGU4Net的几何处理核心类，提供以下功能：

分类	功能	方法示例
格式转换	WKT/GeoJSON/Geometry互转	Wkt2Geometry, Geometry2Geojson
空间关系	判断几何关系	Intersects, Contains, Within
空间分析	几何运算	Buffer, Intersection, Union
属性计算	测量和属性	Area, Length, Centroid
拓扑验证	有效性检查	IsValid, IsSimple
几何操作	简化、凸包等	Simplify, ConvexHull

7.1.2 使用前提

using OpenGIS.Utils.Geometry;
using OpenGIS.Utils.Configuration;

// GeometryUtil依赖GDAL，确保已初始化
GdalConfiguration.ConfigureGdal();

7.2 格式转换

7.2.1 WKT与Geometry互转

using OpenGIS.Utils.Geometry;
using OgrGeometry = OSGeo.OGR.Geometry;

// WKT → Geometry
string wkt = "POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))";
OgrGeometry geom = GeometryUtil.Wkt2Geometry(wkt);

// Geometry → WKT
string wktBack = GeometryUtil.Geometry2Wkt(geom);
Console.WriteLine(wktBack);
// 输出: POLYGON ((0 0,10 0,10 10,0 10,0 0))

7.2.2 WKT与GeoJSON互转

// WKT → GeoJSON
string wkt = "POINT (116.404 39.915)";
string geojson = GeometryUtil.Wkt2Geojson(wkt);
Console.WriteLine(geojson);
// 输出: { "type": "Point", "coordinates": [ 116.404, 39.915 ] }

// 注意：GeoJSON → WKT 当前不支持
// GDAL/OGR不支持直接解析GeoJSON字符串
// 请使用WKT格式，或通过GdalReader从文件加载GeoJSON

7.2.3 Geometry → GeoJSON

// Geometry对象转GeoJSON
var geom = GeometryUtil.Wkt2Geometry("LINESTRING (0 0, 1 1, 2 0)");
string geojson = GeometryUtil.Geometry2Geojson(geom);
Console.WriteLine(geojson);
// 输出: { "type": "LineString", "coordinates": [ [ 0.0, 0.0 ], [ 1.0, 1.0 ], [ 2.0, 0.0 ] ] }

7.3 空间关系判断

7.3.1 关系类型

OGU4Net支持OGC定义的空间关系：

关系	方法	说明
相交	Intersects	A和B有公共部分
包含	Contains	A完全包含B
被包含	Within	A完全在B内部
相接	Touches	A和B仅边界接触
交叉	Crosses	A和B部分交叉
重叠	Overlaps	A和B部分重叠
不相交	Disjoint	A和B完全分离

7.3.2 使用Geometry对象

var polyA = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
var polyB = GeometryUtil.Wkt2Geometry("POLYGON ((5 5, 15 5, 15 15, 5 15, 5 5))");
var point = GeometryUtil.Wkt2Geometry("POINT (5 5)");

// 相交判断
bool intersects = GeometryUtil.Intersects(polyA, polyB);
Console.WriteLine($"A与B相交: {intersects}");  // true

// 包含判断
bool contains = GeometryUtil.Contains(polyA, point);
Console.WriteLine($"A包含点: {contains}");  // true

// 在内部判断
bool within = GeometryUtil.Within(point, polyA);
Console.WriteLine($"点在A内: {within}");  // true

// 相接判断
var touchPoly = GeometryUtil.Wkt2Geometry("POLYGON ((10 0, 20 0, 20 10, 10 10, 10 0))");
bool touches = GeometryUtil.Touches(polyA, touchPoly);
Console.WriteLine($"A与C相接: {touches}");  // true

// 不相交判断
var farPoly = GeometryUtil.Wkt2Geometry("POLYGON ((100 100, 110 100, 110 110, 100 110, 100 100))");
bool disjoint = GeometryUtil.Disjoint(polyA, farPoly);
Console.WriteLine($"A与D不相交: {disjoint}");  // true

7.3.3 使用WKT

// 基于WKT的空间关系判断（更便捷）
string wktA = "POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))";
string wktB = "POINT (5 5)";

bool intersects = GeometryUtil.IntersectsWkt(wktA, wktB);
Console.WriteLine($"相交: {intersects}");  // true

bool contains = GeometryUtil.ContainsWkt(wktA, wktB);
Console.WriteLine($"包含: {contains}");  // true

7.3.4 实际应用：点在面内查询

// 查询某个点落在哪个行政区
var point = "POINT (116.404 39.915)";

var districts = OguLayerUtil.ReadLayer(DataFormatType.SHP, "districts.shp");

foreach (var feature in districts.Features)
{
    if (GeometryUtil.ContainsWkt(feature.Wkt, point))
    {
        var name = feature.GetValue("Name");
        Console.WriteLine($"点位于: {name}");
        break;
    }
}

7.4 空间分析

7.4.1 缓冲区分析

// 基于Geometry对象
var geom = GeometryUtil.Wkt2Geometry("POINT (0 0)");
var buffered = GeometryUtil.Buffer(geom, 10.0);
Console.WriteLine(GeometryUtil.Geometry2Wkt(buffered));

// 基于WKT（更便捷）
string wkt = "LINESTRING (0 0, 10 0)";
string bufferedWkt = GeometryUtil.BufferWkt(wkt, 5.0);
Console.WriteLine(bufferedWkt);

// 缓冲区应用：查找范围内的要素
var centerPoint = "POINT (116.404 39.915)";
var searchArea = GeometryUtil.BufferWkt(centerPoint, 0.01);  // 约1公里

var pois = OguLayerUtil.ReadLayer(
    DataFormatType.SHP,
    "pois.shp",
    spatialFilterWkt: searchArea
);
Console.WriteLine($"范围内有 {pois.GetFeatureCount()} 个POI");

7.4.2 交集运算

// 两个多边形求交
var polyA = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
var polyB = GeometryUtil.Wkt2Geometry("POLYGON ((5 5, 15 5, 15 15, 5 15, 5 5))");

var intersection = GeometryUtil.Intersection(polyA, polyB);
Console.WriteLine(GeometryUtil.Geometry2Wkt(intersection));
// 输出: POLYGON ((5 5,10 5,10 10,5 10,5 5))

// 基于WKT
string wktA = "POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))";
string wktB = "POLYGON ((5 5, 15 5, 15 15, 5 15, 5 5))";
string intersectionWkt = GeometryUtil.IntersectionWkt(wktA, wktB);

7.4.3 并集运算

// 两个几何合并
var geomA = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 5 0, 5 5, 0 5, 0 0))");
var geomB = GeometryUtil.Wkt2Geometry("POLYGON ((5 0, 10 0, 10 5, 5 5, 5 0))");

var union = GeometryUtil.Union(geomA, geomB);
Console.WriteLine(GeometryUtil.Geometry2Wkt(union));

// 多个几何合并
var geometries = new List<OgrGeometry>
{
    GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 5 0, 5 5, 0 5, 0 0))"),
    GeometryUtil.Wkt2Geometry("POLYGON ((5 0, 10 0, 10 5, 5 5, 5 0))"),
    GeometryUtil.Wkt2Geometry("POLYGON ((0 5, 5 5, 5 10, 0 10, 0 5))")
};

var merged = GeometryUtil.Union(geometries);

// 基于WKT的多几何合并
var wktList = new[]
{
    "POLYGON ((0 0, 5 0, 5 5, 0 5, 0 0))",
    "POLYGON ((5 0, 10 0, 10 5, 5 5, 5 0))"
};
string mergedWkt = GeometryUtil.UnionWkt(wktList);

7.4.4 差集运算

// A与B的差集（A中不属于B的部分）
var polyA = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
var polyB = GeometryUtil.Wkt2Geometry("POLYGON ((5 5, 15 5, 15 15, 5 15, 5 5))");

var difference = GeometryUtil.Difference(polyA, polyB);
Console.WriteLine(GeometryUtil.Geometry2Wkt(difference));

// 对称差集（两者不重叠的部分）
var symDiff = GeometryUtil.SymDifference(polyA, polyB);
Console.WriteLine(GeometryUtil.Geometry2Wkt(symDiff));

7.4.5 叠加分析应用

// 应用示例：计算建设用地与保护区的冲突区域

var buildingLand = OguLayerUtil.ReadLayer(DataFormatType.SHP, "building_land.shp");
var protectedArea = OguLayerUtil.ReadLayer(DataFormatType.SHP, "protected_area.shp");

// 合并所有保护区
var protectedGeoms = protectedArea.Features
    .Where(f => !string.IsNullOrEmpty(f.Wkt))
    .Select(f => GeometryUtil.Wkt2Geometry(f.Wkt))
    .ToList();
var mergedProtected = GeometryUtil.Union(protectedGeoms);

// 检查每个建设用地地块与保护区的交集
var conflicts = new List<(OguFeature, double)>();

foreach (var feature in buildingLand.Features)
{
    if (string.IsNullOrEmpty(feature.Wkt)) continue;
    
    var landGeom = GeometryUtil.Wkt2Geometry(feature.Wkt);
    
    if (GeometryUtil.Intersects(landGeom, mergedProtected))
    {
        var intersect = GeometryUtil.Intersection(landGeom, mergedProtected);
        var conflictArea = GeometryUtil.Area(intersect);
        
        if (conflictArea > 0)
        {
            conflicts.Add((feature, conflictArea));
        }
    }
}

Console.WriteLine($"发现 {conflicts.Count} 处冲突:");
foreach (var (feature, area) in conflicts)
{
    Console.WriteLine($"  地块 {feature.Fid}: 冲突面积 {area:F2}");
}

7.5 属性计算

7.5.1 面积计算

// 计算面积
var geom = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
double area = GeometryUtil.Area(geom);
Console.WriteLine($"面积: {area}");  // 100

// 基于WKT
double areaWkt = GeometryUtil.AreaWkt("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");

// 注意：面积单位取决于坐标系
// 地理坐标系(4326)：面积单位是平方度（无实际意义）
// 投影坐标系：面积单位是平方米

7.5.2 长度计算

// 计算长度（线的长度或面的周长）
var line = GeometryUtil.Wkt2Geometry("LINESTRING (0 0, 3 4)");
double length = GeometryUtil.Length(line);
Console.WriteLine($"线长: {length}");  // 5

var poly = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
double perimeter = GeometryUtil.Length(poly);
Console.WriteLine($"周长: {perimeter}");  // 40

// 基于WKT
double lengthWkt = GeometryUtil.LengthWkt("LINESTRING (0 0, 10 0)");

7.5.3 质心计算

// 计算质心
var poly = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
var centroid = GeometryUtil.Centroid(poly);
Console.WriteLine(GeometryUtil.Geometry2Wkt(centroid));  // POINT (5 5)

// 基于WKT
string centroidWkt = GeometryUtil.CentroidWkt("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
Console.WriteLine(centroidWkt);  // POINT (5 5)

7.5.4 内部点

// 获取保证在几何内部的点（对于凹多边形更准确）
var concavePoly = GeometryUtil.Wkt2Geometry(
    "POLYGON ((0 0, 10 0, 10 10, 5 5, 0 10, 0 0))");

var interiorPoint = GeometryUtil.InteriorPoint(concavePoly);
Console.WriteLine(GeometryUtil.Geometry2Wkt(interiorPoint));

7.5.5 距离计算

// 计算两个几何之间的距离
var pointA = GeometryUtil.Wkt2Geometry("POINT (0 0)");
var pointB = GeometryUtil.Wkt2Geometry("POINT (3 4)");

double distance = GeometryUtil.Distance(pointA, pointB);
Console.WriteLine($"距离: {distance}");  // 5

// 判断是否在指定距离内
bool isNear = GeometryUtil.IsWithinDistance(pointA, pointB, 10.0);
Console.WriteLine($"在10单位内: {isNear}");  // true

7.6 几何操作

7.6.1 边界和外包矩形

// 获取边界
var poly = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
var boundary = GeometryUtil.Boundary(poly);
Console.WriteLine(GeometryUtil.Geometry2Wkt(boundary));
// 输出: LINESTRING (0 0,10 0,10 10,0 10,0 0)

// 获取外包矩形
var envelope = GeometryUtil.Envelope(poly);
Console.WriteLine(GeometryUtil.Geometry2Wkt(envelope));

7.6.2 凸包

// 计算凸包
var points = GeometryUtil.Wkt2Geometry(
    "MULTIPOINT ((0 0), (10 0), (5 5), (10 10), (0 10), (5 3))");
var convexHull = GeometryUtil.ConvexHull(points);
Console.WriteLine(GeometryUtil.Geometry2Wkt(convexHull));

7.6.3 简化

// 简化几何（减少顶点数）
var complexLine = GeometryUtil.Wkt2Geometry(
    "LINESTRING (0 0, 1 0.1, 2 0, 3 0.1, 4 0, 5 0.1)");

var simplified = GeometryUtil.Simplify(complexLine, 0.2);
Console.WriteLine(GeometryUtil.Geometry2Wkt(simplified));

// 基于WKT
string simplifiedWkt = GeometryUtil.SimplifyWkt(
    "LINESTRING (0 0, 1 0.1, 2 0, 3 0.1, 4 0)", 0.2);

7.6.4 密化

// 密化几何（增加顶点）
var sparseLine = GeometryUtil.Wkt2Geometry("LINESTRING (0 0, 10 0)");
var densified = GeometryUtil.Densify(sparseLine, 1.0);  // 每隔1个单位加点
Console.WriteLine(GeometryUtil.Geometry2Wkt(densified));

7.7 拓扑验证

7.7.1 有效性验证

using OpenGIS.Utils.Engine.Model;

// 验证几何有效性
var validPoly = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
var result1 = GeometryUtil.IsValid(validPoly);
Console.WriteLine($"有效: {result1.IsValid}");  // true

// 自相交的无效多边形
var invalidPoly = GeometryUtil.Wkt2Geometry(
    "POLYGON ((0 0, 10 10, 0 10, 10 0, 0 0))");
var result2 = GeometryUtil.IsValid(invalidPoly);
Console.WriteLine($"有效: {result2.IsValid}");  // false
Console.WriteLine($"错误: {result2.ErrorMessage}");

7.7.2 简单性检查

// 检查几何是否简单（无自相交）
var simpleLine = GeometryUtil.Wkt2Geometry("LINESTRING (0 0, 10 10)");
var simpleResult = GeometryUtil.IsSimple(simpleLine);
Console.WriteLine($"简单: {simpleResult.IsSimple}");  // true

var selfIntersect = GeometryUtil.Wkt2Geometry("LINESTRING (0 0, 10 10, 0 10, 10 0)");
var complexResult = GeometryUtil.IsSimple(selfIntersect);
Console.WriteLine($"简单: {complexResult.IsSimple}");  // false
Console.WriteLine($"原因: {complexResult.Reason}");

7.7.3 拓扑验证应用

// 验证图层中所有几何的有效性
var layer = OguLayerUtil.ReadLayer(DataFormatType.SHP, "polygons.shp");

var invalidFeatures = new List<(int Fid, string Error)>();

foreach (var feature in layer.Features)
{
    if (string.IsNullOrEmpty(feature.Wkt)) continue;
    
    try
    {
        var geom = GeometryUtil.Wkt2Geometry(feature.Wkt);
        var result = GeometryUtil.IsValid(geom);
        
        if (!result.IsValid)
        {
            invalidFeatures.Add((feature.Fid, result.ErrorMessage ?? "Unknown error"));
        }
    }
    catch (Exception ex)
    {
        invalidFeatures.Add((feature.Fid, ex.Message));
    }
}

if (invalidFeatures.Any())
{
    Console.WriteLine($"发现 {invalidFeatures.Count} 个无效几何:");
    foreach (var (fid, error) in invalidFeatures)
    {
        Console.WriteLine($"  FID {fid}: {error}");
    }
}
else
{
    Console.WriteLine("所有几何都是有效的");
}

7.8 几何类型判断

7.8.1 获取几何类型

using OpenGIS.Utils.Engine.Enums;

var point = GeometryUtil.Wkt2Geometry("POINT (0 0)");
var type = GeometryUtil.GetGeometryType(point);
Console.WriteLine($"类型: {type}");  // POINT

var line = GeometryUtil.Wkt2Geometry("LINESTRING (0 0, 1 1)");
type = GeometryUtil.GetGeometryType(line);
Console.WriteLine($"类型: {type}");  // LINESTRING

var poly = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");
type = GeometryUtil.GetGeometryType(poly);
Console.WriteLine($"类型: {type}");  // POLYGON

7.8.2 几何维度和点数

var geom = GeometryUtil.Wkt2Geometry("POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0))");

// 维度（0=点，1=线，2=面）
int dim = GeometryUtil.Dimension(geom);
Console.WriteLine($"维度: {dim}");  // 2

// 点数
int numPoints = GeometryUtil.NumPoints(geom);
Console.WriteLine($"点数: {numPoints}");  // 5

// 是否为空
bool isEmpty = GeometryUtil.IsEmpty(geom);
Console.WriteLine($"是否为空: {isEmpty}");  // false

7.9 几何比较

7.9.1 精确比较

var geomA = GeometryUtil.Wkt2Geometry("POINT (0 0)");
var geomB = GeometryUtil.Wkt2Geometry("POINT (0 0)");
var geomC = GeometryUtil.Wkt2Geometry("POINT (0.0001 0)");

// 精确比较
bool exact = GeometryUtil.EqualsExact(geomA, geomB);
Console.WriteLine($"A精确等于B: {exact}");  // true

// 带容差的比较
bool tolerance = GeometryUtil.EqualsExactTolerance(geomA, geomC, 0.001);
Console.WriteLine($"A约等于C(容差0.001): {tolerance}");  // true

// 拓扑相等
bool topo = GeometryUtil.EqualsTopo(geomA, geomB);
Console.WriteLine($"A拓扑等于B: {topo}");  // true

7.10 小结

本章详细介绍了GeometryUtil提供的几何处理功能：

格式转换：WKT/GeoJSON/Geometry互转
空间关系：Intersects, Contains, Within等判断
空间分析：Buffer, Intersection, Union, Difference
属性计算：Area, Length, Centroid, Distance
几何操作：Simplify, ConvexHull, Envelope
拓扑验证：IsValid, IsSimple

这些功能是GIS应用开发的核心基础，掌握它们可以实现丰富的空间分析功能。

posted @ 2025-12-03 15:51 我才是银古阅读(0) 评论(0) 收藏举报

刷新页面返回顶部

梦开始的地方

慢慢来会比较快！

第07章 - 几何处理与空间分析

第07章 - 几何处理与空间分析

7.1 GeometryUtil概述

7.1.1 功能分类

7.1.2 使用前提

7.2 格式转换

7.2.1 WKT与Geometry互转

7.2.2 WKT与GeoJSON互转

7.2.3 Geometry → GeoJSON

7.3 空间关系判断

7.3.1 关系类型

7.3.2 使用Geometry对象

7.3.3 使用WKT

7.3.4 实际应用：点在面内查询

7.4 空间分析

7.4.1 缓冲区分析

7.4.2 交集运算

7.4.3 并集运算

7.4.4 差集运算

7.4.5 叠加分析应用

7.5 属性计算

7.5.1 面积计算

7.5.2 长度计算

7.5.3 质心计算

7.5.4 内部点

7.5.5 距离计算

7.6 几何操作

7.6.1 边界和外包矩形

7.6.2 凸包

7.6.3 简化

7.6.4 密化

7.7 拓扑验证

7.7.1 有效性验证

7.7.2 简单性检查

7.7.3 拓扑验证应用

7.8 几何类型判断

7.8.1 获取几何类型

7.8.2 几何维度和点数

7.9 几何比较

7.9.1 精确比较

7.10 小结

公告