znlgis 博客

GIS开发与技术分享

第9章:仿射变换

仿射变换是几何学和计算机图形学中的基本操作,包括平移、旋转、缩放、错切等。Shapely 通过 shapely.affinity 模块提供了完整的 2D/3D 仿射变换支持。本章将详细介绍仿射变换的 数学原理和 Shapely 中各种变换函数的使用方法。


仿射变换概述

什么是仿射变换

仿射变换是一种保持”直线性”和”平行性”的几何变换。即:

  • 变换前是直线的,变换后仍是直线
  • 变换前平行的线段,变换后仍然平行

常见的仿射变换包括:

变换类型 说明
平移(Translation) 沿指定方向移动
旋转(Rotation) 绕某点旋转指定角度
缩放(Scale) 按比例放大或缩小
错切(Skew) 沿某轴倾斜
镜像(Reflection) 沿某轴翻转

2D 仿射变换矩阵

2D 仿射变换可以用 3×3 矩阵表示:

| a  b  xoff |   | x |   | a*x + b*y + xoff |
| d  e  yoff | × | y | = | d*x + e*y + yoff |
| 0  0  1    |   | 1 |   | 1                |

在 Shapely 中,2D 变换矩阵用 6 个元素表示:[a, b, d, e, xoff, yoff]

3D 仿射变换矩阵

3D 变换使用 12 个元素:[a, b, c, d, e, f, g, h, i, xoff, yoff, zoff]

对应的变换:

x' = a*x + b*y + c*z + xoff
y' = d*x + e*y + f*z + yoff
z' = g*x + h*y + i*z + zoff

affine_transform:通用仿射变换

基本用法

affine_transform(geometry, matrix) 对几何对象应用任意仿射变换。

from shapely.affinity import affine_transform
from shapely.geometry import Point, LineString, Polygon

# 单位矩阵(不做变换)
identity = [1, 0, 0, 1, 0, 0]
point = Point(1, 2)
result = affine_transform(point, identity)
print(f"单位变换: {result}")  # POINT (1 2)

# 平移矩阵:x+10, y+20
translate_matrix = [1, 0, 0, 1, 10, 20]
result = affine_transform(point, translate_matrix)
print(f"平移: {result}")  # POINT (11 22)

2D 缩放矩阵

from shapely.affinity import affine_transform
from shapely.geometry import Polygon

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

# x 方向放大 2 倍,y 方向放大 3 倍
scale_matrix = [2, 0, 0, 3, 0, 0]
scaled = affine_transform(poly, scale_matrix)
print(f"原始面积: {poly.area}")     # 1.0
print(f"缩放后面积: {scaled.area}")  # 6.0
print(f"缩放后坐标: {list(scaled.exterior.coords)}")

旋转矩阵

import math
from shapely.affinity import affine_transform
from shapely.geometry import Point

# 绕原点旋转 90 度
angle = math.radians(90)
cos_a = math.cos(angle)
sin_a = math.sin(angle)
rotation_matrix = [cos_a, -sin_a, sin_a, cos_a, 0, 0]

point = Point(1, 0)
rotated = affine_transform(point, rotation_matrix)
print(f"旋转90°: ({rotated.x:.6f}, {rotated.y:.6f})")  # (0, 1)

组合变换矩阵

import math
from shapely.affinity import affine_transform
from shapely.geometry import Polygon

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

# 先缩放2倍,再旋转45度,再平移(10, 5)
# 注意:矩阵乘法顺序是从右到左

# 缩放
s = 2
# 旋转
angle = math.radians(45)
c, s_val = math.cos(angle), math.sin(angle)

# 组合矩阵:T × R × S
# 其中 S = [s,0,0,s,0,0], R = [c,-s,s,c,0,0], T = [1,0,0,1,10,5]
a = s * c
b = s * (-s_val)
d = s * s_val
e = s * c
combined = [a, b, d, e, 10, 5]

result = affine_transform(poly, combined)
print(f"组合变换结果: {result}")
print(f"面积: {result.area:.2f}")  # 面积放大4倍

3D 变换

from shapely.affinity import affine_transform
from shapely.geometry import Point

# 3D 点的平移
point_3d = Point(1, 2, 3)
# 3D 矩阵:[a,b,c, d,e,f, g,h,i, xoff,yoff,zoff]
matrix_3d = [1, 0, 0, 0, 1, 0, 0, 0, 1, 10, 20, 30]
result = affine_transform(point_3d, matrix_3d)
print(f"3D 平移: {result}")  # POINT Z (11 22 33)

translate:平移变换

基本用法

translate(geometry, xoff=0.0, yoff=0.0, zoff=0.0) 将几何对象沿指定方向平移。

from shapely.affinity import translate
from shapely.geometry import Point, Polygon, LineString

# 点的平移
point = Point(0, 0)
moved = translate(point, xoff=5, yoff=3)
print(f"平移后: {moved}")  # POINT (5 3)

# 多边形的平移
poly = Polygon([(0, 0), (2, 0), (2, 2), (0, 2)])
moved_poly = translate(poly, xoff=10, yoff=20)
print(f"平移后坐标: {list(moved_poly.exterior.coords)}")

# 线的平移
line = LineString([(0, 0), (1, 1), (2, 0)])
moved_line = translate(line, xoff=-1, yoff=5)
print(f"平移后线: {list(moved_line.coords)}")

3D 平移

from shapely.affinity import translate
from shapely.geometry import Point

point_3d = Point(1, 2, 3)
moved = translate(point_3d, xoff=10, yoff=20, zoff=30)
print(f"3D 平移: {moved}")

批量平移示例

from shapely.affinity import translate
from shapely.geometry import Polygon

# 创建网格排列的多边形
base = Polygon([(0, 0), (1, 0), (1, 1), (0, 1)])
grid = []
for i in range(5):
    for j in range(5):
        grid.append(translate(base, xoff=i * 1.5, yoff=j * 1.5))

print(f"创建了 {len(grid)} 个多边形")
print(f"第一个: {list(grid[0].exterior.coords)}")
print(f"最后一个: {list(grid[-1].exterior.coords)}")

rotate:旋转变换

基本用法

rotate(geometry, angle, origin='center', use_radians=False) 将几何对象绕指定点旋转。

from shapely.affinity import rotate
from shapely.geometry import Polygon

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

# 绕几何中心旋转 45 度
rotated = rotate(poly, 45)
print(f"旋转 45°: {rotated}")
print(f"面积不变: {rotated.area:.2f}")  # 面积不变

origin 参数

origin 参数指定旋转中心:

from shapely.affinity import rotate
from shapely.geometry import Polygon, Point

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

# 绕几何中心旋转(默认)
r1 = rotate(poly, 90, origin='center')
print(f"绕中心: {[f'({x:.0f},{y:.0f})' for x, y in r1.exterior.coords]}")

# 绕质心旋转
r2 = rotate(poly, 90, origin='centroid')
print(f"绕质心: {[f'({x:.0f},{y:.0f})' for x, y in r2.exterior.coords]}")

# 绕原点旋转
r3 = rotate(poly, 90, origin=Point(0, 0))
print(f"绕原点: {[f'({x:.0f},{y:.0f})' for x, y in r3.exterior.coords]}")

# 绕指定点旋转
r4 = rotate(poly, 90, origin=Point(1, 1))
print(f"绕(1,1): {[f'({x:.0f},{y:.0f})' for x, y in r4.exterior.coords]}")

角度单位

import math
from shapely.affinity import rotate
from shapely.geometry import Point

point = Point(1, 0)

# 使用度(默认)
r_deg = rotate(point, 90, origin=Point(0, 0))
print(f"90度旋转: ({r_deg.x:.6f}, {r_deg.y:.6f})")

# 使用弧度
r_rad = rotate(point, math.pi / 2, origin=Point(0, 0), use_radians=True)
print(f"π/2弧度旋转: ({r_rad.x:.6f}, {r_rad.y:.6f})")

旋转方向

Shapely 中正角度表示逆时针旋转,负角度表示顺时针旋转:

from shapely.affinity import rotate
from shapely.geometry import Point

point = Point(1, 0)
origin = Point(0, 0)

# 逆时针 90°
ccw = rotate(point, 90, origin=origin)
print(f"逆时针90°: ({ccw.x:.1f}, {ccw.y:.1f})")  # (0, 1)

# 顺时针 90°
cw = rotate(point, -90, origin=origin)
print(f"顺时针90°: ({cw.x:.1f}, {cw.y:.1f})")    # (0, -1)

scale:缩放变换

基本用法

scale(geometry, xfact=1.0, yfact=1.0, zfact=1.0, origin='center') 按比例缩放。

from shapely.affinity import scale
from shapely.geometry import Polygon

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

# 均匀放大 2 倍
scaled = scale(poly, xfact=2, yfact=2)
print(f"放大2倍面积: {scaled.area}")  # 16.0

# 非均匀缩放
stretched = scale(poly, xfact=3, yfact=1)
print(f"x方向拉伸面积: {stretched.area}")  # 12.0

origin 参数

from shapely.affinity import scale
from shapely.geometry import Polygon, Point

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

# 绕中心缩放(默认)
s1 = scale(poly, xfact=2, yfact=2, origin='center')
print(f"绕中心: {list(s1.exterior.coords)}")

# 绕原点缩放
s2 = scale(poly, xfact=2, yfact=2, origin=Point(0, 0))
print(f"绕原点: {list(s2.exterior.coords)}")

# 绕质心缩放
s3 = scale(poly, xfact=2, yfact=2, origin='centroid')
print(f"绕质心: {list(s3.exterior.coords)}")

镜像操作

通过负缩放因子实现镜像:

from shapely.affinity import scale
from shapely.geometry import Polygon, Point

# 三角形
tri = Polygon([(0, 0), (3, 0), (1, 2)])

# X 轴镜像(上下翻转)
mirror_x = scale(tri, xfact=1, yfact=-1, origin=Point(0, 0))
print(f"X轴镜像: {list(mirror_x.exterior.coords)}")

# Y 轴镜像(左右翻转)
mirror_y = scale(tri, xfact=-1, yfact=1, origin=Point(0, 0))
print(f"Y轴镜像: {list(mirror_y.exterior.coords)}")

# 中心点镜像
mirror_center = scale(tri, xfact=-1, yfact=-1, origin=Point(0, 0))
print(f"中心镜像: {list(mirror_center.exterior.coords)}")

缩小操作

from shapely.affinity import scale
from shapely.geometry import Polygon

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

# 缩小到 50%
small = scale(poly, xfact=0.5, yfact=0.5)
print(f"原始面积: {poly.area}")
print(f"缩小后面积: {small.area}")
print(f"缩小后坐标: {list(small.exterior.coords)}")

skew:错切 / 倾斜变换

基本用法

skew(geometry, xs=0.0, ys=0.0, origin='center', use_radians=False) 对几何对象进行错切变换。

  • xs:沿 x 轴方向的错切角度
  • ys:沿 y 轴方向的错切角度
from shapely.affinity import skew
from shapely.geometry import Polygon

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

# X 方向错切 30°
skewed_x = skew(poly, xs=30)
print(f"X错切: {list(skewed_x.exterior.coords)}")
print(f"面积不变: {skewed_x.area:.2f}")  # 面积不变

# Y 方向错切 15°
skewed_y = skew(poly, ys=15)
print(f"Y错切: {list(skewed_y.exterior.coords)}")

组合错切

from shapely.affinity import skew
from shapely.geometry import Polygon

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

# 同时在两个方向错切
skewed = skew(poly, xs=20, ys=10)
print(f"组合错切面积: {skewed.area:.2f}")  # 面积不变
print(f"坐标: {[f'({x:.2f},{y:.2f})' for x, y in skewed.exterior.coords]}")

origin 参数

from shapely.affinity import skew
from shapely.geometry import Polygon, Point

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

# 绕中心错切(默认)
s1 = skew(poly, xs=30, origin='center')

# 绕原点错切
s2 = skew(poly, xs=30, origin=Point(0, 0))

# 绕指定点错切
s3 = skew(poly, xs=30, origin=Point(1, 1))

for i, s in enumerate([s1, s2, s3], 1):
    coords = [f'({x:.2f},{y:.2f})' for x, y in s.exterior.coords]
    print(f"方案{i}: {coords}")

使用弧度

import math
from shapely.affinity import skew
from shapely.geometry import Polygon

poly = Polygon([(0, 0), (2, 0), (2, 2), (0, 2)])
skewed = skew(poly, xs=math.pi / 6, use_radians=True)  # 30度
print(f"弧度错切: {list(skewed.exterior.coords)}")

变换组合

链式变换

多个变换可以依次应用,实现复杂的几何变换:

from shapely.affinity import translate, rotate, scale
from shapely.geometry import Polygon, Point

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

# 步骤1: 绕原点旋转 45°
step1 = rotate(poly, 45, origin=Point(0, 0))
print(f"旋转后: {[f'({x:.2f},{y:.2f})' for x, y in step1.exterior.coords]}")

# 步骤2: 放大 2 倍
step2 = scale(step1, xfact=2, yfact=2, origin=Point(0, 0))
print(f"缩放后: {[f'({x:.2f},{y:.2f})' for x, y in step2.exterior.coords]}")

# 步骤3: 平移到 (10, 5)
step3 = translate(step2, xoff=10, yoff=5)
print(f"平移后: {[f'({x:.2f},{y:.2f})' for x, y in step3.exterior.coords]}")

# 面积验证:应该是原来的 4 倍
print(f"原始面积: {poly.area}, 最终面积: {step3.area:.2f}")

通用变换函数

from shapely.affinity import affine_transform
import math

def compose_transforms(*matrices):
    """组合多个 2D 仿射变换矩阵(从右到左应用)"""
    # 将 6 元素列表转为 3x3 矩阵
    def to_3x3(m):
        return [
            [m[0], m[1], m[4]],
            [m[2], m[3], m[5]],
            [0, 0, 1],
        ]

    def from_3x3(m):
        return [m[0][0], m[0][1], m[1][0], m[1][1], m[0][2], m[1][2]]

    def mat_mul(a, b):
        result = [[0] * 3 for _ in range(3)]
        for i in range(3):
            for j in range(3):
                for k in range(3):
                    result[i][j] += a[i][k] * b[k][j]
        return result

    result = to_3x3(matrices[-1])
    for m in reversed(matrices[:-1]):
        result = mat_mul(to_3x3(m), result)
    return from_3x3(result)


# 使用组合变换
from shapely.geometry import Polygon

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

# 定义单独变换矩阵
scale_m = [2, 0, 0, 2, 0, 0]  # 放大2倍
angle = math.radians(45)
rotate_m = [math.cos(angle), -math.sin(angle),
            math.sin(angle), math.cos(angle), 0, 0]
translate_m = [1, 0, 0, 1, 10, 5]  # 平移

# 组合为一个矩阵
combined = compose_transforms(translate_m, rotate_m, scale_m)
result = affine_transform(poly, combined)
print(f"组合变换结果面积: {result.area:.2f}")

实际应用场景

场景 1:坐标系转换辅助

from shapely.affinity import translate, scale
from shapely.geometry import Polygon

# 将局部坐标转换到地图坐标
local_building = Polygon([(0, 0), (20, 0), (20, 15), (0, 15)])

# 已知:原点对应 (500000, 4000000),比例尺 1:1
map_building = translate(local_building, xoff=500000, yoff=4000000)
print(f"地图坐标: {list(map_building.exterior.coords)[:2]}...")

场景 2:建筑物旋转对齐

from shapely.affinity import rotate
from shapely.geometry import Polygon, Point

# 建筑物轮廓(倾斜放置)
building = Polygon([
    (100, 100), (120, 110), (115, 130), (95, 120)
])

# 将建筑物旋转对齐到正北方向
# 假设需要旋转 -26.57 度
centroid = building.centroid
aligned = rotate(building, -26.57, origin=centroid)
print(f"对齐前: {[f'({x:.0f},{y:.0f})' for x, y in building.exterior.coords]}")
print(f"对齐后: {[f'({x:.1f},{y:.1f})' for x, y in aligned.exterior.coords]}")

场景 3:创建对称图案

from shapely.affinity import rotate, translate, scale
from shapely import union_all
from shapely.geometry import Polygon, Point

# 基础形状
petal = Polygon([(0, 0), (1, 0.3), (2, 0), (1, -0.3)])

# 通过旋转创建花朵图案
petals = []
for angle in range(0, 360, 60):
    rotated = rotate(petal, angle, origin=Point(0, 0))
    petals.append(rotated)

flower = union_all(petals)
print(f"花朵图案面积: {flower.area:.2f}")
print(f"花瓣数: {len(petals)}")

场景 4:空间校正

from shapely.affinity import affine_transform
from shapely.geometry import Polygon

# 假设有一组控制点对应关系
# 原始坐标 -> 校正后坐标
# 通过最小二乘法或其他方法求得仿射变换参数
# 这里直接给出一个校正矩阵
correction_matrix = [1.0001, 0.0002, -0.0001, 0.9999, 0.5, -0.3]

# 对一批图形应用校正
parcels = [
    Polygon([(100, 200), (150, 200), (150, 250), (100, 250)]),
    Polygon([(200, 300), (280, 300), (280, 380), (200, 380)]),
]

for i, parcel in enumerate(parcels):
    corrected = affine_transform(parcel, correction_matrix)
    print(f"地块{i} 校正前质心: ({parcel.centroid.x:.1f}, {parcel.centroid.y:.1f})")
    print(f"地块{i} 校正后质心: ({corrected.centroid.x:.1f}, {corrected.centroid.y:.1f})")

变换的数学验证

旋转矩阵验证

import math
from shapely.affinity import rotate
from shapely.geometry import Point

# 验证旋转的正确性
point = Point(3, 0)
origin = Point(0, 0)

for angle in [30, 45, 60, 90, 180, 270]:
    rotated = rotate(point, angle, origin=origin)
    rad = math.radians(angle)
    expected_x = 3 * math.cos(rad)
    expected_y = 3 * math.sin(rad)
    print(f"{angle:3d}°: 实际=({rotated.x:.4f}, {rotated.y:.4f}), "
          f"理论=({expected_x:.4f}, {expected_y:.4f})")

面积守恒验证

from shapely.affinity import rotate, skew, translate
from shapely.geometry import Polygon

poly = Polygon([(0, 0), (3, 0), (3, 2), (0, 2)])
original_area = poly.area

# 旋转保持面积
rotated = rotate(poly, 37)
print(f"旋转后面积: {rotated.area:.6f} (原始: {original_area})")

# 平移保持面积
translated = translate(poly, 100, 200)
print(f"平移后面积: {translated.area:.6f}")

# 错切保持面积
skewed = skew(poly, xs=25)
print(f"错切后面积: {skewed.area:.6f}")

缩放面积变化

from shapely.affinity import scale
from shapely.geometry import Polygon, Point

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

for xf, yf in [(2, 1), (1, 3), (2, 3), (0.5, 0.5)]:
    scaled = scale(poly, xfact=xf, yfact=yf, origin=Point(0, 0))
    ratio = scaled.area / poly.area
    expected = abs(xf * yf)
    print(f"xfact={xf}, yfact={yf}: "
          f"面积比={ratio:.2f}, 理论值={expected:.2f}")

本章小结

本章介绍了 Shapely 中的仿射变换功能:

函数 说明 保持面积
affine_transform 通用仿射变换 取决于矩阵
translate 平移
rotate 旋转
scale 缩放 ✗(按因子变化)
skew 错切

关键要点:

  1. 2D 变换使用 6 元素矩阵 [a, b, d, e, xoff, yoff]
  2. origin 参数控制变换中心:'center''centroid'Point 对象
  3. rotate 默认使用角度制,use_radians=True 切换到弧度
  4. 负缩放因子实现镜像操作
  5. 多个变换可链式组合,也可通过矩阵乘法合并为单一变换

下一章将介绍坐标操作与变换,学习如何提取、修改坐标并与外部投影库集成。