惯性聚合 高效追踪和阅读你感兴趣的博客、新闻、科技资讯
阅读原文 在惯性聚合中打开

推荐订阅源

H
Heimdal Security Blog
A
Arctic Wolf
K
Kaspersky official blog
V
Vulnerabilities – Threatpost
Threat Intelligence Blog | Flashpoint
Threat Intelligence Blog | Flashpoint
Simon Willison's Weblog
Simon Willison's Weblog
L
LINUX DO - 热门话题
MongoDB | Blog
MongoDB | Blog
T
Threat Research - Cisco Blogs
D
Docker
爱范儿
爱范儿
T
Tenable Blog
C
Check Point Blog
B
Blog
C
Cisco Blogs
Vercel News
Vercel News
The Cloudflare Blog
T
Threatpost
NISL@THU
NISL@THU
T
Tor Project blog
V2EX - 技术
V2EX - 技术
P
Palo Alto Networks Blog
Application and Cybersecurity Blog
Application and Cybersecurity Blog
T
Tailwind CSS Blog
G
GRAHAM CLULEY
P
Privacy & Cybersecurity Law Blog
SecWiki News
SecWiki News
博客园 - 司徒正美
S
Security @ Cisco Blogs
GbyAI
GbyAI
S
Secure Thoughts
Microsoft Security Blog
Microsoft Security Blog
The Register - Security
The Register - Security
Recorded Future
Recorded Future
Cloudbric
Cloudbric
Webroot Blog
Webroot Blog
N
News and Events Feed by Topic
Y
Y Combinator Blog
博客园_首页
T
Troy Hunt's Blog
The Hacker News
The Hacker News
雷峰网
雷峰网
Google DeepMind News
Google DeepMind News
U
Unit 42
AWS News Blog
AWS News Blog
PCI Perspectives
PCI Perspectives
V
Visual Studio Blog
博客园 - 聂微东
有赞技术团队
有赞技术团队
酷 壳 – CoolShell
酷 壳 – CoolShell

博客园 - 我才是银古

第16章:常见问题、排错与最佳实践 第15章:扩展生态、MCAD 与外部集成 第12章:实战案例:机械结构与 3D 打印零件 第14章:构建、测试、调试与贡献流程 第13章:OpenSCAD 源码架构与核心执行流程 第11章:预览、渲染、网格精度与性能优化 第09章:列表推导、递归与算法建模 第08章:参数化零件库与复用设计 第10章:导入导出、命令行与自动化 第06章:CSG 布尔建模方法 第07章:二维图形、拉伸、旋转与投影 第05章:基础几何、坐标系与变换 第04章:参数、变量、函数、模块与作用域 OpenSCAD 教程目录 第03章:OpenSCAD 语言基础 第02章:安装、环境配置与开发工作流 第01章:OpenSCAD 项目全景与学习路线 第02章:源码获取、编译与开发环境配置 第01章:OCCT项目全景与学习路线 第18章:二次开发实战与综合案例 第18章:综合实战案例 第17章:数据交换与协同 第16章:源码架构与二次开发 第15章:插件与自定义工作台开发 第14章:Python脚本宏与自动化 第13章:FEM仿真分析 第12章:CAM数控加工 第11章:SurfaceMesh与逆向工程 第10章:Draft二维绘图与BIM建筑 第09章:工程图TechDraw 第07章:参数化表达式与Spreadsheet 第08章:装配设计Assembly 第06章:Part工作台与几何内核 第05章:PartDesign实体特征建模 第04章:草图Sketcher约束建模 第02章:安装版本与工作环境配置 第03章:界面工作台与基础操作 第01章:项目全景与学习路线 第十二章:插件开发、研究功能与最佳实践 第十章:定时任务与自动化(Cron) 第七章:技能、记忆与自学习闭环 第八章:MCP 集成与上下文文件 第六章:工具系统与终端后端 第五章:模型供应商与配置体系 Hermes Agent 教程目录 第十一章:语音、视觉、浏览器与子代理协作 第四章:CLI/TUI 与会话管理 第十二章:学习路线、实战方案与最佳实践 第十一章:源码结构、开发调试与插件开发 第十章:自动化、远程访问、日志与排障 第九章:Control UI、节点、Canvas 与语音能力 第七章:工具、技能、插件与能力扩展 第八章:安全模型、访问控制与沙箱实践 第六章:Agent 工作区、会话与多智能体路由 第五章:多通道消息接入与聊天平台配置 第四章:配置体系、模型接入与认证管理 第三章:Gateway 架构、协议与运行机制 第二章:安装、环境准备与快速上手 第一章:OpenClaw 项目概览与核心定位 oh-my-openagent 教程目录 09-命令模型回退与配置参考 10-实战案例最佳实践与故障排除 05-工作模式-Ultrawork-Prometheus-Atlas 08-Hooks与MCP系统 06-Category与Skill系统 07-核心工具链 04-智能体全景详解 03-安装与环境配置 02-整体架构与多模型编排机制 01-项目简介与核心理念 01-项目概览与学习路线 02-安装部署与工具适配 03-Skill机制与using-superpowers 05-TDD系统化调试与完成前验证 04-需求澄清方案设计与计划编写 07-并行智能体子智能体与Git-Worktree 第六章:代码审查、反馈处理与分支收尾 08-中国特色Skills与本土团队落地 09-MCP构建工作流执行与自定义Skill 第23章:FreeCAD-Python-API Clipper2 C# 源码解读教程 第19章:PolyTree 多边形树结构 第20章:实际应用与最佳实践 第18章:Minkowski 和与差 第17章:RectClip 矩形裁剪优化 第16章:ClipperOffset 偏移类详解 第15章:填充规则详解 第14章:布尔运算执行流程 第13章:ClipperD 浮点裁剪类 第11章:OutRec 与 OutPt 输出结构 第9章:Active 活动边结构 第10章:Vertex 顶点与 LocalMinima 局部极小值 第12章:Clipper64 裁剪类详解 第7章:高精度运算与128位整数 第8章:ClipperBase 基类详解 第5章:枚举类型与常量定义 第6章:InternalClipper 内部工具类 第2章:核心数据结构 - Point64、PointD 第3章:路径与多边形表示 - Path64、PathD、Paths64、PathsD 第4章:矩形边界 - Rect64、RectD
第10章:坐标操作与变换
我才是银古 · 2026-06-23 · via 博客园 - 我才是银古

第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 数组输入,实现高效的批量处理