程序设计与计算思维¶

Computer Programming and Computational Thinking

第 5 讲：图像变换 III¶

2025—2026学年度春季学期

清华大学 韩文弢

内容回顾¶

  • Python 语言
    • 元组
    • 函数（参数默认值、关键字参数、Lambda 表达式、函数对象）
  • 变换
    • 线性变换与非线性变换
    • 线性变换与矩阵

准备工作¶

In [1]:
from urllib.request import urlretrieve
import numpy as np
import matplotlib.pyplot as plt
import ipywidgets as widgets

# 设置作图时的中文字体，避免乱码
plt.rcParams['font.sans-serif'] = ['Noto Sans CJK SC']

img_path, _ = urlretrieve('https://pacman.cs.tsinghua.edu.cn/~hanwentao/cpct/05/longcorgi.png')
img = plt.imread(img_path)[:, :, :3]
plt.imshow(img)
Out[1]:
<matplotlib.image.AxesImage at 0x10aa99be0>
No description has been provided for this image

最终效果¶

In [2]:
black = np.array([0, 0, 0])
white = np.array([1, 1, 1])
green = np.array([0, 1, 0])
red = np.array([1, 0, 0])

from math import sin, cos, hypot, atan2

# 恒等
identity = lambda x, y: (x, y)

# 伸缩
scalex = lambda alpha: lambda x, y: (alpha*x, y)
scaley = lambda alpha: lambda x, y: (x, alpha*y)
scale = lambda alpha: lambda x, y: (alpha*x, alpha*y)

# 翻转
swap = lambda x, y: (y, x)
flipy = lambda x, y: (x, -y)

# 旋转
rotate = lambda theta: lambda x, y: (cos(theta)*x + sin(theta)*y, -sin(theta)*x + cos(theta)*y)

# 剪切
shear = lambda alpha: lambda x, y: (x + alpha*y, y)

# 一般线性变换
lin = lambda a, b, c, d: lambda x, y: (a*x + b*y, c*x + d*y)

# 用矩阵表示
def lin_mat(A):
    a = A[0, 0]
    b = A[0, 1]
    c = A[1, 0]
    d = A[1, 1]
    return lambda x, y: lin(a, b, c, d)

# 平移，是仿射变换，但不是线性变换
translate = lambda alpha, beta: lambda x, y: (x + alpha, y + beta)

# 非线性剪切
nonlin_shear = lambda alpha: lambda x, y: (x, y + alpha*x**2)

# 扭曲
warp = lambda alpha: lambda x, y: rotate(alpha*hypot(x, y))(x, y)

# 极坐标转直角坐标
xy = lambda rho, theta: (rho*cos(theta), rho*sin(theta))

# 直角坐标转极坐标
rho_theta = lambda x, y: (hypot(x, y), atan2(y, x))


def transform_image(img, transform_func, pixels=400): 
    """
    对图像应用变换
    
    参数:
    - img: 输入图像 (NumPy 数组)
    - transform_func: 变换函数 (x, y) -> (new_x, new_y)
    - pixels: 输出图像的像素数
    """
    rows, cols = img.shape[:2]
    m = max(cols, rows)
    
    # 创建输出图像
    output = np.zeros((pixels, pixels, img.shape[2]), dtype=img.dtype)
    
    for i in range(pixels):
        for j in range(pixels):
            # 将像素坐标转换为 [-1, 1] 范围
            x, y = scale(2/pixels)(*flipy(*swap(*translate(-pixels/2, -pixels/2)(i, j))))
            
            # 应用变换的逆（因为我们是从输出像素找输入像素）
            new_x, new_y = transform_func(x, y)
            
            # 将坐标转换回像素索引
            src_i, src_j = translate(rows/2, cols/2)(*swap(*flipy(*scale(m/2)(new_x, new_y))))
            src_i = int(src_i)
            src_j = int(src_j)
            
            if 0 <= src_i < rows and 0 <= src_j < cols:
                output[i, j] = img[src_i, src_j]
            else:
                output[i, j] = black  # 黑色背景
    
    return output


def interactive_transform(a=1.0, b=0.0, c=0.0, d=1.0):
    """
    交互式变换演示
    """
    transformed = transform_image(img, lin(a, b, c, d), pixels=400)
    
    plt.figure(figsize=(8, 8))
    plt.imshow(transformed)
    plt.title(f'变换矩阵: a={a:.1f}, b={b:.1f}, c={c:.1f}, d={d:.1f}')
    plt.axis('off')
    plt.show()

# 创建交互式控件
interact = widgets.interactive(
    interactive_transform,
    a=widgets.FloatSlider(value=1.0, min=-1.5, max=1.5, step=0.1, description='a:'),
    b=widgets.FloatSlider(value=0.0, min=-1.5, max=1.5, step=0.1, description='b:'),
    c=widgets.FloatSlider(value=0.0, min=-1.5, max=1.5, step=0.1, description='c:'),
    d=widgets.FloatSlider(value=1.0, min=-1.5, max=1.5, step=0.1, description='d:'),
)

display(interact)

interactive(children=(FloatSlider(value=1.0, description='a:', max=1.5, min=-1.5), FloatSlider(value=0.0, desc…

对象变换与坐标变换¶

  • 对象变换是指通过改变对象的坐标来实现变换，如点 $(x, y) \to (2x, 2y)$ 是进行了 2 倍的放大
  • 坐标变换是指通过改变坐标轴来实现变换，例如 $x$ 和 $y$ 坐标轴都缩小一半，点的坐标也会发生上述变换

这两种变换的效果可能相同（点的坐标发生同样的变化），但是含义不同。

为了做想要的变换时有一致的效果，我们将图像的坐标系变换到 $[-1, 1]^2$ 的区域里再进行变换。

image.png

完整的过程¶

完整的变换过程如下所示，注意我们是从变换后的图像中取每个像素，然后计算对应变换前的图像中哪个像素：

image.png

逆变换¶

如果 $f$ 是从 2 维向量到 2 维向量的函数，我们定义 $f$ 的逆，记为 $f^{-1}$，具有“撤销” $f$ 效果的性质，即

$$f(f^{-1}(v)) = v$$

和

$$f^{-1}(f(v)) = v$$

这个等式可能对所有 $v$ 成立，也可能对某些 $v$ 成立。

示例：放大和缩小¶

In [3]:
# 验证缩放变换的逆
v = np.random.rand(2)

# scale(2) 的逆是 scale(0.5)
T = scale(0.5)
T_inv = scale(2)

# 应用 T 然后应用 T_inv
result = T_inv(*T(*v))

print(f'原始向量: {v}')
print(f'变换后: {result}')
print(f'相等？ {np.allclose(v, result)}')

原始向量: [0.82284975 0.77955379]
变换后: (np.float64(0.822849746691098), np.float64(0.7795537925491557))
相等？ True
In [4]:
# 验证旋转变换的逆
v = np.random.rand(2)

# rotate(θ) 的逆是 rotate(-θ)
theta = np.radians(30)
T = rotate(theta)
T_inv = rotate(-theta)

result = T_inv(*T(*v))

print(f'原始向量: {v}')
print(f'变换后: {result}')
print(f'相等？ {np.allclose(v, result)}')

原始向量: [0.11806597 0.27107506]
变换后: (np.float64(0.1180659687326121), np.float64(0.2710750580133623))
相等？ True

求解逆变换¶

逆变换到底是在做什么？

再次考虑缩放。假设将输入向量 $v$ 缩放 2 得到输出向量 $u$：

$$u = 2 v$$

现在假设想进行反向操作。如果给定 $u$，如何找到 $v$？

在这个特定情况下，有 $v = \frac{1}{2} u$。

考虑更加一般的情况，如果有一个线性变换，可以写成

$$u = A \, v$$

其中 $A$ 是矩阵。

如果给定 $u$ 并想反向找到 $v$，需要求解线性方程组。

通常可以（什么时候不可以？）解这个方程组找到一个新矩阵 $B$ 使得

$$v = B \, u$$

即 $B$ 撤销了 $A$ 的效果。那么，

$$v = (B \, A) v$$

因此 $B A$ 必须是单位矩阵，称 $B$ 为 $A$ 的矩阵逆，并写成

$$B = A^{-1}$$

对于 $2 \times 2$ 的矩阵，可以写出矩阵逆的公式。对于更大规模的矩阵，通常需要计算机运行算法来找到逆。

In [5]:
# 验证矩阵逆
A = np.random.randn(2, 2)
v = np.random.rand(2)

# 使用矩阵逆
A_inv = np.linalg.inv(A)

# 应用 A 然后应用 A_inv
result = A_inv @ (A @ v)

print(f'原始向量: {v}')
print(f'变换后: {result}')
print(f'相等？ {np.allclose(v, result)}')

原始向量: [0.83213606 0.74802034]
变换后: [0.83213606 0.74802034]
相等？ True
In [6]:
# 验证 (A*B)^(-1) = B^(-1) * A^(-1)
A = np.random.randn(2, 2)
B = np.random.randn(2, 2)

left = np.linalg.inv(A @ B)
right = np.linalg.inv(B) @ np.linalg.inv(A)

print(f'(A*B)^(-1):\n{left}')
print(f'B^(-1) * A^(-1):\n{right}')
print(f'相等？ {np.allclose(left, right)}')

(A*B)^(-1):
[[-0.81035038 -2.08689412]
 [ 0.52604683  0.31146011]]
B^(-1) * A^(-1):
[[-0.81035038 -2.08689412]
 [ 0.52604683  0.31146011]]
相等？ True

2x2 矩阵逆的公式¶

如果 $A = \begin{pmatrix} a & b \\ c & d \end{pmatrix}$，则

$$A^{-1} = \frac{1}{ad-bc} \begin{pmatrix} d & -b \\ -c & a \end{pmatrix}$$

In [7]:
# 手动计算 2x2 矩阵逆
def inv_2x2_manual(A):
    a, b = A[0, 0], A[0, 1]
    c, d = A[1, 0], A[1, 1]
    det = a * d - b * c
    if abs(det) < 1e-10:
        raise ValueError("矩阵不可逆")
    return np.array([[d, -b], [-c, a]]) / det

A = np.random.randn(2, 2)
print(f'NumPy 逆:\n{np.linalg.inv(A)}')
print(f'手动计算逆:\n{inv_2x2_manual(A)}')
print(f'相等？ {np.allclose(np.linalg.inv(A), inv_2x2_manual(A))}')

NumPy 逆:
[[ 0.72936789 -1.0798709 ]
 [-0.11181915  2.1705302 ]]
手动计算逆:
[[ 0.72936789 -1.0798709 ]
 [-0.11181915  2.1705302 ]]
相等？ True

非线性变换的逆变换¶

如果有一个非线性变换 $T$，能求逆吗？换句话说，如果 $u = T(v)$，能否解这个方程找到用 $v$ 关于 $u$ 的函数表示？

通常这是一个困难的问题，有时可以解析地做到，但通常不能。

尽管如此，有一些数值方法有时可以解这些方程，例如牛顿法。在 scipy 中可以使用 fsolve 函数来进行数值求解。

In [8]:
from scipy.optimize import fsolve

# 使用数值方法求非线性变换的逆
def inverse_numerical(f, y, x0=(0.0, 0.0)):
    """
    数值求解非线性变换的逆
    f: 变换函数 (x, y) -> (u, v)
    y: 目标输出 (u, v)
    x0: 初始猜测
    """
    def equations(x):
        u, v = f(x[0], x[1])
        return [u - y[0], v - y[1]]
    
    solution = fsolve(equations, x0)
    return solution[0], solution[1]

# 测试：非线性剪切的逆
alpha = 0.5
original_point = (0.5, 0.5)
transformed = nonlin_shear(alpha)(*original_point)
# 创建逆函数
inverse_nonlin_shear = lambda alpha: lambda x, y: inverse_numerical(nonlin_shear(alpha), (x, y))
recovered = inverse_nonlin_shear(alpha)(*transformed)
# 等价于 recovered = inverse_numerical(nonlin_shear(alpha), transformed)

print(f'原始点: {original_point}')
print(f'变换后: {transformed}')
print(f'恢复的点: {recovered}')

原始点: (0.5, 0.5)
变换后: (0.5, 0.625)
恢复的点: (np.float64(0.5), np.float64(0.5))

图像变换的完整流程¶

下面实现完整的图像变换（反向生成）：

In [9]:
def transform_image(img, transform_func, pixels=400): 
    """
    对图像应用变换
    
    参数:
    - img: 输入图像 (NumPy 数组)
    - transform_func: 变换函数 (x, y) -> (new_x, new_y)
    - pixels: 输出图像的像素数
    """
    rows, cols = img.shape[:2]
    m = max(cols, rows)
    
    # 创建输出图像
    output = np.zeros((pixels, pixels, img.shape[2]), dtype=img.dtype)
    
    for i in range(pixels):
        for j in range(pixels):
            # 将像素坐标转换为 [-1, 1] 范围
            x, y = scale(2/pixels)(*flipy(*swap(*translate(-pixels/2, -pixels/2)(i, j))))
            
            # 应用变换的逆（因为我们是从输出像素找输入像素）
            new_x, new_y = transform_func(x, y)
            
            # 将坐标转换回像素索引
            src_i, src_j = translate(rows/2, cols/2)(*swap(*flipy(*scale(m/2)(new_x, new_y))))
            src_i = int(src_i)
            src_j = int(src_j)
            
            if 0 <= src_i < rows and 0 <= src_j < cols:
                output[i, j] = img[src_i, src_j]
            else:
                output[i, j] = black  # 黑色背景
    
    return output
In [10]:
# 示例：应用旋转变换
theta = np.radians(30)
rotated_img = transform_image(img, rotate(-theta))

fig, axes = plt.subplots(1, 2, figsize=(12, 6))
axes[0].imshow(img)
axes[0].set_title('原始图像')
axes[0].axis('off')
axes[1].imshow(rotated_img)
axes[1].set_title(f'旋转 {np.degrees(theta):.1f}°')
axes[1].axis('off')
Out[10]:
(np.float64(-0.5), np.float64(399.5), np.float64(399.5), np.float64(-0.5))
No description has been provided for this image
In [11]:
# 示例：应用剪切变换
alpha = 0.5
sheared_img = transform_image(img, shear(-alpha))

fig, axes = plt.subplots(1, 2, figsize=(12, 6))
axes[0].imshow(img)
axes[0].set_title('原始图像')
axes[0].axis('off')
axes[1].imshow(sheared_img)
axes[1].set_title(f'剪切 (α={alpha})')
axes[1].axis('off')
Out[11]:
(np.float64(-0.5), np.float64(399.5), np.float64(399.5), np.float64(-0.5))