深入理解地理信息系统（GIS）中的面矢量：使用Python的shapely库操作Polygon对象

ytkz2024-04-182024-09-09

在地理信息系统（GIS）中，矢量数据是最基本的数据类型之一。它们用于表示地理空间中的对象，包括它们的形状、位置和属性。矢量数据的主要类型包括点（Point）、线（LineString）和面（Polygon）。在这个讨论中，我们将重点关注面矢量的数据结构和操作。

面矢量（Polygon）是由一个外部轮廓和零个或多个内部轮廓（空洞）组成的。外部轮廓定义了多边形的形状和大小，而内部轮廓则定义了多边形内部的空洞。理解这种数据结构是理解和操作地理数据的关键，无论是从文件中读取数据，还是创建自己的多边形对象。

在本文中，我们将首先介绍如何使用Python的shapely库创建和操作Polygon对象。然后，我们将介绍如何使用fiona库从shp文件中读取矢量数据。最后，我们将介绍多面（MultiPolygon）的概念，这是由多个面（Polygon）组成的几何对象，常用于表示一组不连续的地理区域，例如群岛。

让我们开始吧！

怎么读取shp文件

首先我们学习一个知识点，如果将shp文件读取为Polygon对象。我们用到了Python的fiona库

import fiona
with fiona.open(polygonfile, 'r') as f:
        polygons = [shape(feature['geometry']) for feature in f]

polygons的属性如下所示：

polygons是一个列表，第一个元素即是第一个面矢量（也可能是其他的shp类型）。

面矢量（Polygon）的数据结构

这里，我们构造数据，一步一步地深入学习面矢量。

在Python的shapely库中，一个Polygon对象可以被创建。

一个Polygon对象的数据结构如下：

外部轮廓（Exterior）：这是一个坐标点的列表，表示多边形的外部边界。这些坐标点定义了多边形的形状和大小。在创建Polygon对象时，外部轮廓是必需的。
内部轮廓（Interiors）：这是一个坐标点列表的列表，表示多边形内部的空洞。每个内部轮廓都是一个坐标点的列表，表示一个空洞的边界。在创建Polygon对象时，内部轮廓是可选的。

from shapely.geometry import Polygon

# 多边形的外部轮廓
exterior = [(0, 0), (10, 0), (10, 10), (0, 10)]
# 创建多边形
polygon = Polygon(exterior)

对面矢量可视化

def plot_polygon(polygon1):
        import numpy as np
        import matplotlib.pyplot as plt
    # 提取多边形的外部和内部坐标
        ext = np.array(polygon1.exterior.coords.xy)
        int_coords = [np.array(interior.coords.xy) for interior in polygon1.interiors]

        # 创建图形和轴
        fig, ax = plt.subplots()

        # 绘制外部多边形
        ax.fill(*ext, alpha=0.5, fc='r', ec='none')

        # 绘制内部多边形（空洞）
        for interior in int_coords:
            ax.fill(*interior, alpha=1, fc='w', ec='none')

        # 设置坐标轴比例，并显示图形
        ax.set_aspect('equal', 'box')
        plt.show()
from shapely.geometry import Polygon

# 多边形的外部轮廓
exterior = [(0, 0), (10, 0), (10, 10), (0, 10)]
# 创建多边形
polygon = Polygon(exterior)
plot_polygon(polygon)

单个面的可视化如下：

from shapely.geometry import Polygon

# 多边形的外部轮廓
exterior = [(0, 0), (10, 0), (10, 10), (0, 10)]

# 多边形的内部空洞
interior = [(4, 4), (6, 4), (6, 6), (4, 6)]

# 创建多边形
polygon = Polygon(exterior, [interior])

带空洞的面矢量可视化如下：

那带两个空洞呢？

from shapely.geometry import Polygon

# 多边形的外部轮廓
exterior = [(0, 0), (10, 0), (10, 10), (0, 10)]

# 多边形的内部空洞
interior1 = [(4, 4), (6, 4), (6, 6), (4, 6)]
interior2 = [(1, 1), (1, 2), (2, 2), (2, 1)] # 第二个空洞
# 创建多边形
polygon = Polygon(exterior, [interior1,interior2])

带一个空洞的面矢量可视化如下：

那带三个呢？

from shapely.geometry import Polygon

# 多边形的外部轮廓
exterior = [(0, 0), (10, 0), (10, 10), (0, 10)]

# 多边形的内部空洞
interior1 = [(4, 4), (6, 4), (6, 6), (4, 6)]
interior2 = [(1, 1), (1, 2), (2, 2), (2, 1)] # 第二个空洞
interior3 = [(7, 7), (7, 8), (8, 8), (8, 7)] # 第三个空洞
# 创建多边形
polygon = Polygon(exterior, [interior1,interior2,interior3])

带三个空洞的面矢量可视化如下：

在地理信息系统（GIS）中，”多面”（MultiPolygon）是由多个面（Polygon）组成的几何对象。每个面可以有自己的外部边界和一个或多个内部边界（空洞）。多面常常用于表示一组不连续的地理区域，例如群岛。

在Python的shapely库中，我们可以使用MultiPolygon类来创建一个多面对象。下面是一个例子：

from shapely.geometry import Polygon, MultiPolygon

# 创建两个多边形
polygon1 = Polygon([(0, 0), (5, 0), (5, 5), (0, 5)])
polygon2 = Polygon([(6, 6), (11, 6), (11, 11), (6, 11)])

# 创建一个多面
multipolygon = MultiPolygon([polygon1, polygon2])