znlgis 博客

GIS开发与技术分享

第10章:坐标操作与变换

坐标操作是空间数据处理中的核心能力。Shapely 通过 shapely.coordinates 模块提供了 高效的坐标提取、修改和变换功能。结合 NumPy 的向量化能力和 pyproj 等投影库,可以 实现坐标参考系转换、批量坐标修改等高级操作。本章将全面介绍 Shapely 中的坐标操作函数 及其在实际应用中的使用方法。


get_coordinates:提取坐标

基本用法

get_coordinates(geometry) 将几何对象的坐标提取为 NumPy 数组。

from shapely import get_coordinates
from shapely.geometry import Point, LineString, Polygon
import numpy as np

# 点
point = Point(1, 2)
coords = get_coordinates(point)
print(f"点坐标:\n{coords}")
# [[1. 2.]]

# 线
line = LineString([(0, 0), (1, 1), (2, 0)])
coords = get_coordinates(line)
print(f"线坐标:\n{coords}")
# [[0. 0.]
#  [1. 1.]
#  [2. 0.]]

# 多边形
poly = Polygon([(0, 0), (4, 0), (4, 3), (0, 3)])
coords = get_coordinates(poly)
print(f"多边形坐标:\n{coords}")
print(f"坐标形状: {coords.shape}")

多几何对象

from shapely import get_coordinates
from shapely.geometry import Point, MultiPoint
import numpy as np

# MultiPoint
mpt = MultiPoint([(0, 0), (1, 1), (2, 2)])
coords = get_coordinates(mpt)
print(f"MultiPoint 坐标:\n{coords}")

# 几何数组
geoms = np.array([Point(i, i * 2) for i in range(5)])
coords = get_coordinates(geoms)
print(f"数组坐标:\n{coords}")

include_z 参数

from shapely import get_coordinates
from shapely.geometry import Point

# 3D 点
point_3d = Point(1, 2, 3)

# 默认不包含 Z
coords_2d = get_coordinates(point_3d)
print(f"2D 坐标: {coords_2d}")  # [[1. 2.]]

# 包含 Z
coords_3d = get_coordinates(point_3d, include_z=True)
print(f"3D 坐标: {coords_3d}")  # [[1. 2. 3.]]

复杂几何体

from shapely import get_coordinates
from shapely.geometry import Polygon

# 带洞多边形
poly_with_hole = Polygon(
    [(0, 0), (10, 0), (10, 10), (0, 10)],
    [[(2, 2), (4, 2), (4, 4), (2, 4)]]
)

coords = get_coordinates(poly_with_hole)
print(f"带洞多边形坐标数: {len(coords)}")
print(f"坐标:\n{coords}")

count_coordinates:统计坐标数量

count_coordinates(geometry) 返回几何对象中的坐标点数量。

from shapely import count_coordinates
from shapely.geometry import Point, LineString, Polygon, MultiPolygon
import numpy as np

# 各种几何类型
point = Point(0, 0)
line = LineString([(0, 0), (1, 1), (2, 0)])
poly = Polygon([(0, 0), (1, 0), (1, 1), (0, 1)])

print(f"点坐标数: {count_coordinates(point)}")     # 1
print(f"线坐标数: {count_coordinates(line)}")       # 3
print(f"多边形坐标数: {count_coordinates(poly)}")   # 5 (闭合环,首尾重复)

# 批量统计
geoms = np.array([point, line, poly])
counts = count_coordinates(geoms)
print(f"批量坐标数: {counts}")  # [1, 3, 5]

复杂几何体的坐标计数

from shapely import count_coordinates
from shapely.geometry import Polygon, MultiPolygon

# 带洞多边形
poly_hole = Polygon(
    [(0, 0), (10, 0), (10, 10), (0, 10)],
    [[(2, 2), (4, 2), (4, 4), (2, 4)]]
)
print(f"带洞多边形: {count_coordinates(poly_hole)} 个坐标")

# MultiPolygon
mp = MultiPolygon([
    Polygon([(0, 0), (1, 0), (1, 1), (0, 1)]),
    Polygon([(2, 2), (3, 2), (3, 3), (2, 3)]),
])
print(f"MultiPolygon: {count_coordinates(mp)} 个坐标")

set_coordinates:替换坐标

基本用法

set_coordinates(geometry, coordinates) 用新坐标替换几何对象的坐标,返回新的几何对象。

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Polygon
import numpy as np

poly = Polygon([(0, 0), (1, 0), (1, 1), (0, 1)])
print(f"原始: {list(poly.exterior.coords)}")

# 提取坐标并修改
coords = get_coordinates(poly)
new_coords = coords * 10  # 所有坐标乘以 10
new_poly = set_coordinates(poly, new_coords)
print(f"修改后: {list(new_poly.exterior.coords)}")

坐标偏移

from shapely import get_coordinates, set_coordinates
from shapely.geometry import LineString
import numpy as np

line = LineString([(0, 0), (1, 1), (2, 0), (3, 1)])

# 所有 y 坐标增加 5
coords = get_coordinates(line)
coords[:, 1] += 5
new_line = set_coordinates(line, coords)
print(f"偏移后: {list(new_line.coords)}")

批量坐标替换

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Point
import numpy as np

# 多个点的坐标替换
points = np.array([Point(i, i) for i in range(5)])
coords = get_coordinates(points)

# 应用某种变换
new_coords = coords.copy()
new_coords[:, 0] = coords[:, 0] * 2 + 10
new_coords[:, 1] = coords[:, 1] * 3 - 5

new_points = set_coordinates(points, new_coords)
for p in new_points:
    print(f"  ({p.x:.1f}, {p.y:.1f})")

保持几何结构

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Polygon

# 带洞多边形
poly = Polygon(
    [(0, 0), (10, 0), (10, 10), (0, 10)],
    [[(3, 3), (7, 3), (7, 7), (3, 7)]]
)

# 坐标乘以 2
coords = get_coordinates(poly)
new_poly = set_coordinates(poly, coords * 2)

print(f"原始: 外环顶点={len(poly.exterior.coords)}, 内环数={len(list(poly.interiors))}")
print(f"新的: 外环顶点={len(new_poly.exterior.coords)}, 内环数={len(list(new_poly.interiors))}")
print(f"几何结构保持: {new_poly.is_valid}")

transform:应用变换函数(旧版)

注意shapely.ops.transform 在 Shapely 2.1+ 中已标记为弃用, 推荐使用 shapely.transform 或基于 get_coordinates / set_coordinates 的方式。

基本用法

from shapely.ops import transform
from shapely.geometry import Point, Polygon

# 简单的坐标变换函数
def shift_coords(x, y):
    return x + 10, y + 20

point = Point(1, 2)
shifted = transform(shift_coords, point)
print(f"平移: {shifted}")  # POINT (11 22)

# 多边形变换
poly = Polygon([(0, 0), (1, 0), (1, 1), (0, 1)])
shifted_poly = transform(shift_coords, poly)
print(f"多边形平移: {list(shifted_poly.exterior.coords)}")

复杂变换

import math
from shapely.ops import transform
from shapely.geometry import LineString

# 极坐标变换:将笛卡尔坐标转为极坐标后的点
def cartesian_to_transformed(x, y):
    """自定义坐标变换"""
    new_x = x * math.cos(y * 0.1)
    new_y = x * math.sin(y * 0.1)
    return new_x, new_y

line = LineString([(1, 0), (2, 1), (3, 2), (4, 3)])
transformed = transform(cartesian_to_transformed, line)
print(f"变换后: {list(transformed.coords)}")

transform_coordseq:新版坐标变换(Shapely 2.1+)

shapely.transform 是 Shapely 2.1+ 推荐的坐标变换方式, 它接收一个以 NumPy 数组为输入输出的变换函数。

基本用法

from shapely import transform as shapely_transform
from shapely.geometry import Polygon
import numpy as np

poly = Polygon([(0, 0), (1, 0), (1, 1), (0, 1)])

# 变换函数接收 (N, 2) 或 (N, 3) 的 NumPy 数组
def my_transform(coords):
    """坐标翻倍"""
    return coords * 2

result = shapely_transform(poly, my_transform)
print(f"变换后: {list(result.exterior.coords)}")

向量化变换

from shapely import transform as shapely_transform
from shapely.geometry import LineString
import numpy as np

line = LineString([(0, 0), (1, 0), (2, 1), (3, 1)])

# 非线性变换
def wave_transform(coords):
    """添加波浪效果"""
    result = coords.copy()
    result[:, 1] += np.sin(coords[:, 0] * np.pi) * 0.5
    return result

wavy = shapely_transform(line, wave_transform)
print(f"波浪线坐标: {[f'({x:.2f},{y:.2f})' for x, y in wavy.coords]}")

批量变换

from shapely import transform as shapely_transform
from shapely.geometry import Point
import numpy as np

# 批量几何体变换
points = np.array([Point(i, i * 2) for i in range(5)])

def offset_transform(coords):
    result = coords.copy()
    result[:, 0] += 100
    result[:, 1] += 200
    return result

new_points = shapely_transform(points, offset_transform)
for p in new_points:
    print(f"  ({p.x:.0f}, {p.y:.0f})")

自定义坐标变换函数

旋转变换

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Polygon
import numpy as np

def rotate_coords(geometry, angle_deg, cx=0, cy=0):
    """绕指定点旋转几何体"""
    angle = np.radians(angle_deg)
    cos_a = np.cos(angle)
    sin_a = np.sin(angle)

    coords = get_coordinates(geometry)
    # 平移到原点
    dx = coords[:, 0] - cx
    dy = coords[:, 1] - cy
    # 旋转
    new_x = dx * cos_a - dy * sin_a + cx
    new_y = dx * sin_a + dy * cos_a + cy

    new_coords = np.column_stack([new_x, new_y])
    return set_coordinates(geometry, new_coords)


poly = Polygon([(0, 0), (2, 0), (2, 1), (0, 1)])
rotated = rotate_coords(poly, 45, cx=1, cy=0.5)
print(f"自定义旋转: {[f'({x:.2f},{y:.2f})' for x, y in rotated.exterior.coords]}")

极坐标转笛卡尔坐标

from shapely import get_coordinates, set_coordinates
from shapely.geometry import LineString
import numpy as np

def polar_to_cartesian(geometry):
    """将极坐标 (r, theta) 转换为笛卡尔坐标 (x, y)"""
    coords = get_coordinates(geometry)
    r = coords[:, 0]
    theta = coords[:, 1]
    x = r * np.cos(theta)
    y = r * np.sin(theta)
    return set_coordinates(geometry, np.column_stack([x, y]))


# 创建极坐标下的螺旋线
angles = np.linspace(0, 4 * np.pi, 100)
radii = np.linspace(1, 5, 100)
spiral_polar = LineString(zip(radii, angles))

spiral_cartesian = polar_to_cartesian(spiral_polar)
print(f"螺旋线顶点数: {len(spiral_cartesian.coords)}")
coords = get_coordinates(spiral_cartesian)
print(f"X 范围: [{coords[:, 0].min():.2f}, {coords[:, 0].max():.2f}]")
print(f"Y 范围: [{coords[:, 1].min():.2f}, {coords[:, 1].max():.2f}]")

坐标抖动(添加噪声)

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Polygon
import numpy as np

def add_noise(geometry, std=0.1, seed=42):
    """给坐标添加随机噪声"""
    rng = np.random.default_rng(seed)
    coords = get_coordinates(geometry)
    noise = rng.normal(0, std, coords.shape)
    return set_coordinates(geometry, coords + noise)


poly = Polygon([(0, 0), (10, 0), (10, 10), (0, 10)])
noisy = add_noise(poly, std=0.3)
print(f"原始面积: {poly.area:.2f}")
print(f"抖动后面积: {noisy.area:.2f}")
print(f"坐标: {[f'({x:.2f},{y:.2f})' for x, y in noisy.exterior.coords]}")

与 pyproj 集成实现坐标参考系转换

WGS84 → Web Mercator 示例

# 注意:需要安装 pyproj: pip install pyproj
# from pyproj import Transformer
# from shapely.ops import transform
# from shapely.geometry import Point, Polygon

# 创建坐标转换器
# transformer = Transformer.from_crs("EPSG:4326", "EPSG:3857", always_xy=True)

# WGS84 坐标(经度, 纬度)
# beijing = Point(116.4, 39.9)

# 转换到 Web Mercator
# beijing_mercator = transform(transformer.transform, beijing)
# print(f"WGS84: {beijing}")
# print(f"Web Mercator: {beijing_mercator}")

# 以下是不依赖 pyproj 的近似 Web Mercator 变换示例
from shapely import get_coordinates, set_coordinates
from shapely.geometry import Point
import numpy as np
import math

def wgs84_to_web_mercator_approx(geometry):
    """WGS84 到 Web Mercator 的近似转换"""
    coords = get_coordinates(geometry)
    lon = coords[:, 0]
    lat = coords[:, 1]

    x = lon * 20037508.34 / 180.0
    y = np.log(np.tan((90 + lat) * np.pi / 360.0)) / (np.pi / 180.0)
    y = y * 20037508.34 / 180.0

    return set_coordinates(geometry, np.column_stack([x, y]))

beijing = Point(116.4, 39.9)
beijing_mercator = wgs84_to_web_mercator_approx(beijing)
print(f"WGS84: ({beijing.x}, {beijing.y})")
print(f"Web Mercator 近似: ({beijing_mercator.x:.2f}, {beijing_mercator.y:.2f})")

使用 pyproj 进行投影转换的完整示例

# 完整的 pyproj 集成示例(需要安装 pyproj)
#
# from pyproj import Transformer
# from shapely.ops import transform
# from shapely.geometry import Polygon
#
# # WGS84 (EPSG:4326) → UTM Zone 50N (EPSG:32650)
# transformer = Transformer.from_crs(
#     "EPSG:4326", "EPSG:32650", always_xy=True
# )
#
# # 北京附近的区域(经纬度)
# area_wgs84 = Polygon([
#     (116.3, 39.8), (116.5, 39.8),
#     (116.5, 40.0), (116.3, 40.0),
# ])
#
# # 转换到 UTM 投影
# area_utm = transform(transformer.transform, area_wgs84)
#
# # 现在可以进行准确的面积和距离计算
# print(f"WGS84 面积(度²): {area_wgs84.area:.6f}")
# print(f"UTM 面积(m²): {area_utm.area:.0f}")
# print(f"UTM 面积(km²): {area_utm.area / 1e6:.2f}")

Shapely 2.x + pyproj 推荐方式

# Shapely 2.x 推荐使用 get_coordinates / set_coordinates 配合 pyproj
#
# from pyproj import Transformer
# from shapely import get_coordinates, set_coordinates
# from shapely.geometry import Polygon
# import numpy as np
#
# transformer = Transformer.from_crs("EPSG:4326", "EPSG:32650", always_xy=True)
#
# poly = Polygon([(116.3, 39.8), (116.5, 39.8), (116.5, 40.0), (116.3, 40.0)])
#
# # 提取坐标
# coords = get_coordinates(poly)
#
# # 批量转换
# x_new, y_new = transformer.transform(coords[:, 0], coords[:, 1])
# new_coords = np.column_stack([x_new, y_new])
#
# # 设置新坐标
# poly_utm = set_coordinates(poly, new_coords)
# print(f"投影后面积: {poly_utm.area:.0f} m²")

批量坐标操作

处理几何数组

from shapely import get_coordinates, set_coordinates, count_coordinates
from shapely.geometry import Point, LineString, Polygon
import numpy as np

# 不同类型的几何体
geoms = np.array([
    Point(1, 2),
    LineString([(0, 0), (3, 4)]),
    Polygon([(0, 0), (1, 0), (1, 1), (0, 1)]),
])

# 统计坐标数量
counts = count_coordinates(geoms)
print(f"各几何体坐标数: {counts}")
print(f"总坐标数: {counts.sum()}")

# 提取所有坐标
all_coords = get_coordinates(geoms)
print(f"所有坐标形状: {all_coords.shape}")
print(f"坐标:\n{all_coords}")

统一坐标偏移

from shapely import get_coordinates, set_coordinates
from shapely import box
import numpy as np

# 创建多个矩形
geoms = np.array([box(i, i, i + 1, i + 1) for i in range(5)])

# 将所有几何体的坐标平移
coords = get_coordinates(geoms)
coords[:, 0] += 100  # x 偏移 100
coords[:, 1] += 200  # y 偏移 200
shifted_geoms = set_coordinates(geoms, coords)

# 验证
from shapely import centroid
for i, g in enumerate(shifted_geoms):
    c = centroid(g)
    print(f"矩形{i}: 质心=({c.x:.1f}, {c.y:.1f})")

坐标精度控制

set_precision / get_precision

set_precision 可以将几何对象的坐标精度设置为指定的网格大小。

from shapely import set_precision, get_precision
from shapely.geometry import Point, Polygon

# 设置精度
point = Point(1.123456789, 2.987654321)
snapped = set_precision(point, grid_size=0.01)
print(f"原始: {point}")
print(f"精度 0.01: {snapped}")

snapped2 = set_precision(point, grid_size=0.1)
print(f"精度 0.1: {snapped2}")

snapped3 = set_precision(point, grid_size=1.0)
print(f"精度 1.0: {snapped3}")

多边形精度控制

from shapely import set_precision
from shapely.geometry import Polygon

poly = Polygon([
    (0.00001, 0.00002),
    (1.99998, 0.00003),
    (1.99999, 2.00001),
    (0.00002, 1.99999),
])

# 捕捉到 0.001 网格
precise = set_precision(poly, grid_size=0.001)
print(f"精度控制后: {list(precise.exterior.coords)}")

获取当前精度

from shapely import get_precision, set_precision
from shapely.geometry import Point

point = Point(1.5, 2.5)
print(f"默认精度: {get_precision(point)}")  # 0.0(表示全精度)

snapped = set_precision(point, grid_size=0.1)
print(f"设置后精度: {get_precision(snapped)}")  # 0.1

批量精度控制

from shapely import set_precision
from shapely.geometry import Point
import numpy as np

points = np.array([
    Point(1.123, 2.456),
    Point(3.789, 4.012),
    Point(5.345, 6.678),
])

# 批量设置精度
precise_points = set_precision(points, grid_size=0.1)
for p in precise_points:
    print(f"  ({p.x}, {p.y})")

Z 坐标和 M 坐标的处理

Z 坐标

from shapely import get_coordinates
from shapely.geometry import Point, LineString

# 3D 几何
point_3d = Point(1, 2, 3)
line_3d = LineString([(0, 0, 0), (1, 1, 10), (2, 0, 20)])

# 提取包含 Z 的坐标
coords_z = get_coordinates(line_3d, include_z=True)
print(f"3D 坐标:\n{coords_z}")

# 修改 Z 值
import numpy as np
coords_z[:, 2] *= 2
from shapely import set_coordinates
new_line = set_coordinates(line_3d, coords_z)
print(f"Z 值翻倍后: {list(new_line.coords)}")

检查是否有 Z 坐标

from shapely import has_z
from shapely.geometry import Point
import numpy as np

p2d = Point(1, 2)
p3d = Point(1, 2, 3)

print(f"2D 点有Z: {has_z(p2d)}")  # False
print(f"3D 点有Z: {has_z(p3d)}")  # True

# 批量检查
geoms = np.array([p2d, p3d])
print(f"批量检查: {has_z(geoms)}")

添加或移除 Z 坐标

from shapely import get_coordinates, set_coordinates, force_2d, force_3d
from shapely.geometry import LineString
import numpy as np

# 2D 线
line_2d = LineString([(0, 0), (1, 1), (2, 0)])

# 添加 Z 坐标
line_3d = force_3d(line_2d, z=10.0)
coords = get_coordinates(line_3d, include_z=True)
print(f"添加Z后: {list(line_3d.coords)}")

# 移除 Z 坐标
line_flat = force_2d(line_3d)
print(f"移除Z后: {list(line_flat.coords)}")

实际应用

应用 1:坐标偏移校正

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Polygon
import numpy as np

# 假设由于 GPS 偏移,所有坐标需要校正
parcels = [
    Polygon([(100.001, 200.002), (110.003, 200.001),
             (110.002, 210.003), (100.000, 210.001)]),
    Polygon([(110.003, 200.001), (120.002, 200.003),
             (120.001, 210.002), (110.002, 210.003)]),
]

# 校正偏移量
dx, dy = -0.002, -0.001

geoms = np.array(parcels)
coords = get_coordinates(geoms)
coords[:, 0] += dx
coords[:, 1] += dy
corrected = set_coordinates(geoms, coords)

for i, g in enumerate(corrected):
    print(f"地块{i} 校正后: {[f'({x:.3f},{y:.3f})' for x, y in g.exterior.coords]}")

应用 2:投影转换(纯数学实现)

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Polygon
import numpy as np

def simple_mercator(geometry):
    """简化的 Mercator 投影"""
    R = 6378137.0  # 地球半径(米)
    coords = get_coordinates(geometry)
    lon_rad = np.radians(coords[:, 0])
    lat_rad = np.radians(coords[:, 1])

    x = R * lon_rad
    y = R * np.log(np.tan(np.pi / 4 + lat_rad / 2))

    return set_coordinates(geometry, np.column_stack([x, y]))

# 经纬度多边形
area_geo = Polygon([
    (116.0, 39.5), (117.0, 39.5),
    (117.0, 40.5), (116.0, 40.5),
])

area_proj = simple_mercator(area_geo)
print(f"投影前面积(度²): {area_geo.area:.4f}")
print(f"投影后面积(m²): {area_proj.area:.0f}")
print(f"投影后面积(km²): {area_proj.area / 1e6:.2f}")

应用 3:坐标系统一

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Point, LineString, Polygon
import numpy as np

def unify_coordinates(geometries, target_precision=0.001):
    """统一几何体的坐标精度"""
    geoms = np.array(geometries)
    coords = get_coordinates(geoms)

    # 四舍五入到目标精度
    coords = np.round(coords / target_precision) * target_precision

    return set_coordinates(geoms, coords)


geoms = [
    Point(1.00001, 2.00002),
    LineString([(0.00003, 0.00004), (1.99997, 1.99998)]),
    Polygon([(0.00001, 0.00001), (3.00002, 0.00003),
             (3.00001, 3.00002), (0.00002, 3.00001)]),
]

unified = unify_coordinates(geoms)
for g in unified:
    print(f"  {g.geom_type}: {g}")

应用 4:坐标范围归一化

from shapely import get_coordinates, set_coordinates, total_bounds
from shapely.geometry import Polygon
import numpy as np

# 多个不同位置和大小的多边形
polys = np.array([
    Polygon([(10, 20), (30, 20), (30, 40), (10, 40)]),
    Polygon([(50, 60), (70, 60), (70, 80), (50, 80)]),
    Polygon([(25, 35), (45, 35), (45, 55), (25, 55)]),
])

# 获取总范围
minx, miny, maxx, maxy = total_bounds(polys)
print(f"坐标范围: x=[{minx}, {maxx}], y=[{miny}, {maxy}]")

# 归一化到 [0, 1] 范围
coords = get_coordinates(polys)
coords[:, 0] = (coords[:, 0] - minx) / (maxx - minx)
coords[:, 1] = (coords[:, 1] - miny) / (maxy - miny)
normalized = set_coordinates(polys, coords)

minx2, miny2, maxx2, maxy2 = total_bounds(normalized)
print(f"归一化后: x=[{minx2:.2f}, {maxx2:.2f}], y=[{miny2:.2f}, {maxy2:.2f}]")

坐标操作性能提示

向量化 vs 循环

from shapely import get_coordinates, set_coordinates
from shapely.geometry import Point
import numpy as np
import time

# 创建大量点
n = 10000
points = np.array([Point(i * 0.01, i * 0.02) for i in range(n)])

# 方法1:向量化(推荐)
start = time.time()
coords = get_coordinates(points)
coords[:, 0] += 100
coords[:, 1] += 200
result1 = set_coordinates(points, coords)
t1 = time.time() - start

# 方法2:逐个变换(慢)
from shapely.affinity import translate
start = time.time()
result2 = [translate(p, xoff=100, yoff=200) for p in points]
t2 = time.time() - start

print(f"向量化: {t1:.4f}s")
print(f"逐个变换: {t2:.4f}s")
print(f"加速比: {t2 / t1:.1f}x")

内存效率

from shapely import get_coordinates
from shapely.geometry import Point
import numpy as np

# get_coordinates 返回的是连续内存数组
points = np.array([Point(i, i) for i in range(1000)])
coords = get_coordinates(points)

print(f"数组形状: {coords.shape}")
print(f"数据类型: {coords.dtype}")
print(f"内存大小: {coords.nbytes} bytes")
print(f"内存连续: {coords.flags['C_CONTIGUOUS']}")

本章小结

本章全面介绍了 Shapely 的坐标操作功能:

函数 说明
get_coordinates 提取坐标为 NumPy 数组
count_coordinates 统计坐标点数量
set_coordinates 替换坐标(生成新几何体)
shapely.ops.transform 旧版坐标变换(已弃用)
shapely.transform 新版坐标变换(推荐)
set_precision 设置坐标精度
get_precision 获取坐标精度
force_2d / force_3d Z 坐标处理

关键要点:

  1. 使用向量化操作get_coordinates + NumPy 操作 + set_coordinates 是最高效的方式
  2. 与 pyproj 集成:通过 get_coordinates/set_coordinatesshapely.ops.transform 实现投影转换
  3. 精度控制set_precision 可以统一坐标精度,避免浮点问题
  4. Z 坐标支持include_z 参数控制是否处理 Z 坐标
  5. 批量操作:所有函数都支持 NumPy 数组输入,实现高效的批量处理