程序设计与计算思维¶

Computer Programming and Computational Thinking

第 3 讲：图像变换¶

2025—2026学年度春季学期

清华大学 韩文弢

内容回顾¶

  • Python 语言
    • 对象、类型、值
    • 内置类型
      • 简单类型：整数、浮点数、布尔型、空值
      • 序列类型：字符串、列表
    • 第三方库的类型：np.ndarray
    • 赋值
    • 分支：if 语句
    • 循环：while 语句、for 语句
    • 函数：定义和调用
  • 图像
    • 用矩阵（或张量）来表示
    • 通过运算进行变换

准备工作¶

In [1]:
import numpy as np
import matplotlib.pyplot as plt

from urllib.request import urlretrieve

philip_file, _ = urlretrieve("https://pacman.cs.tsinghua.edu.cn/~hanwentao/cpct/01/philip.png")
philip = plt.imread(philip_file)
plt.imshow(philip)
Out[1]:
<matplotlib.image.AxesImage at 0x10e621160>
No description has been provided for this image

采样¶

为了将图像像素化（pixelate），需要对它进行采样操作（sampling）。

由于像素化的图像所需的像素较少，可以先进行降采样（downsampling），每隔若干行或若干列选取一个像素。

In [2]:
r = 20
downsample_philip = philip[::r, ::r]
plt.imshow(downsample_philip)
Out[2]:
<matplotlib.image.AxesImage at 0x10e6bf390>
No description has been provided for this image

降采样之后要让图像恢复原来的大小，需要进行升采样（upsampling），把每个像素扩展成一个正方形。

升采样可以用克罗内克积（Kronecker product）实现。

In [3]:
upsample_philip = np.kron(downsample_philip, np.full((r, r, 1), 1))
plt.imshow(upsample_philip)
Out[3]:
<matplotlib.image.AxesImage at 0x10e74ad50>
No description has been provided for this image

克罗内克积¶

数学上，克罗内克积是两个任意大小的矩阵间的运算，表示为 $\otimes$。简单地说，就是将前一个矩阵的每个元素乘上后一个完整的矩阵。

定义：如果 $A$ 是一个 $m \times n$ 的矩阵，而 $B$ 是一个 $p \times q$ 的矩阵，克罗内克积 $A \otimes B$ 是一个 $mp \times nq$ 的分块矩阵。

$$ A \otimes B = \begin{bmatrix} a_{11} B & \cdots & a_{1n}B \\ \vdots & \ddots & \vdots \\ a_{m1} B & \cdots & a_{mn} B \end{bmatrix} $$

进一步展开就是：

$$ A \otimes B = \begin{bmatrix} a_{11} b_{11} & a_{11} b_{12} & \cdots & a_{11} b_{1q} & \cdots & \cdots & a_{1n} b_{11} & a_{1n} b_{12} & \cdots & a_{1n} b_{1q} \\ a_{11} b_{21} & a_{11} b_{22} & \cdots & a_{11} b_{2q} & \cdots & \cdots & a_{1n} b_{21} & a_{1n} b_{22} & \cdots & a_{1n} b_{2q} \\ \vdots & \vdots & \ddots & \vdots & & & \vdots & \vdots & \ddots & \vdots \\ a_{11} b_{p1} & a_{11} b_{p2} & \cdots & a_{11} b_{pq} & \cdots & \cdots & a_{1n} b_{p1} & a_{1n} b_{p2} & \cdots & a_{1n} b_{pq} \\ \vdots & \vdots & & \vdots & \ddots & & \vdots & \vdots & & \vdots \\ \vdots & \vdots & & \vdots & & \ddots & \vdots & \vdots & & \vdots \\ a_{m1} b_{11} & a_{m1} b_{12} & \cdots & a_{m1} b_{1q} & \cdots & \cdots & a_{mn} b_{11} & a_{mn} b_{12} & \cdots & a_{mn} b_{1q} \\ a_{m1} b_{21} & a_{m1} b_{22} & \cdots & a_{m1} b_{2q} & \cdots & \cdots & a_{mn} b_{21} & a_{mn} b_{22} & \cdots & a_{mn} b_{2q} \\ \vdots & \vdots & \ddots & \vdots & & & \vdots & \vdots & \ddots & \vdots \\ a_{m1} b_{p1} & a_{m1} b_{p2} & \cdots & a_{m1} b_{pq} & \cdots & \cdots & a_{mn} b_{p1} & a_{mn} b_{p2} & \cdots & a_{mn} b_{pq} \end{bmatrix} $$

例如：

$$ \begin{bmatrix} 1 & 2 \\ 3 & 1 \\ \end{bmatrix} \otimes \begin{bmatrix} 0 & 3 \\ 2 & 1 \\ \end{bmatrix} = \begin{bmatrix} 1\cdot 0 & 1\cdot 3 & 2\cdot 0 & 2\cdot 3 \\ 1\cdot 2 & 1\cdot 1 & 2\cdot 2 & 2\cdot 1 \\ 3\cdot 0 & 3\cdot 3 & 1\cdot 0 & 1\cdot 3 \\ 3\cdot 2 & 3\cdot 1 & 1\cdot 2 & 1\cdot 1 \\ \end{bmatrix} = \begin{bmatrix} 0 & 3 & 0 & 6 \\ 2 & 1 & 4 & 2 \\ 0 & 9 & 0 & 3 \\ 6 & 3 & 2 & 1 \end{bmatrix} $$

图像合成（线性组合）¶

数学中的一个重要概念是线性组合。这个概念结合了

  • 缩放一个对象
  • 组合多个对象

通过线性组合，可以实现将多个对象的缩放版本组合在一起。

下面尝试对图像进行线性组合操作。首先载入新的图像。

In [4]:
corgis_file, _ = urlretrieve("https://pacman.cs.tsinghua.edu.cn/~hanwentao/cpct/03/corgis.png")
corgis = plt.imread(corgis_file)[:, :, :3]  # 去掉 alpha 通道（透明度）
plt.imshow(corgis)
Out[4]:
<matplotlib.image.AxesImage at 0x10e7f1450>
No description has been provided for this image

下面对图像中每个像素的颜色进行缩放。注意过渡缩放（颜色分量大于 1）会导致饱和。

In [5]:
c = 0.5
plt.imshow(c * corgis)
Out[5]:
<matplotlib.image.AxesImage at 0x10e85b390>
No description has been provided for this image

交互式结果¶

可以使用 ipywidgets 来创建可以交互的结果。

In [6]:
from ipywidgets import interact

def scale_color(c):
    plt.imshow((c * corgis).clip(0, 1))
    
# 进行交互式操作，最后的分号用于忽略本行代码的输出
interact(scale_color, c=(0.0, 3.0, 0.01));

interactive(children=(FloatSlider(value=1.5, description='c', max=3.0, step=0.01), Output()), _dom_classes=('w…

下面尝试进行多幅图像的组合，需要另外一幅图像，例如可以通过上下颠倒来得到。

In [7]:
upsidedown_corgis = corgis[::-1, :]
plt.imshow(upsidedown_corgis)
Out[7]:
<matplotlib.image.AxesImage at 0x10eb05090>
No description has been provided for this image

对两幅图像进行组合，看看有什么效果。

In [8]:
plt.imshow(0.5 * upsidedown_corgis + 0.5 * corgis)
Out[8]:
<matplotlib.image.AxesImage at 0x10eb63250>
No description has been provided for this image

凸组合¶

如果线性组合中所有系数都非负且和为 1，就说这是一个凸组合。

例如，用 α 和 1-α 作为两个相加为 1 的系数，用不同的 α 值来缩放两幅图像，赋予它们不同的权重。

In [9]:
def mix_image(alpha):
    plt.imshow(alpha * upsidedown_corgis + (1 - alpha) * corgis)

interact(mix_image, alpha=(0.0, 1.0, 0.01));

interactive(children=(FloatSlider(value=0.5, description='alpha', max=1.0, step=0.01), Output()), _dom_classes…

Photoshop 中的操作¶

在使用 Photoshop 处理图像时，会发现滤镜（Photoshop 滤镜参考）很有用。例如：

  1. 模糊
  2. 锐化
  3. 风格化 → 查找边缘
  4. 像素化
  5. 扭曲

其中一些变换（如模糊、锐化、查找边缘）是通过计算卷积（convolution）来实现的。卷积非常高效，在当前的人工智能（机器学习）中经常出现，特别是在图像识别场景中。

图像滤镜（卷积）¶

卷积的概念¶

在做图像变换或检测等操作时，为了考虑像素点周围的情况，需要一种操作能够去以每个像素为中心，对其周围的像素进行某些运算，并汇总成一个数值。卷积是可以实现这一目标的操作，用来描述具体操作的矩阵被称为卷积核（kernel）。

2D_Convolution_Animation.gif

卷积的实现¶

卷积可以通过多重循环来实现。

In [10]:
def my_convolve(image, kernel):
    # 计算卷积后的图像大小（不考虑边界情况），创建结果数组
    ih, iw, ic = image.shape
    kh, kw = kernel.shape
    oh, ow, oc = ih - kh + 1, iw - kw + 1, ic
    out = np.zeros((oh, ow, oc))

    # 卷积计算
    for y in range(oh):  # 图像的每行
        for x in range(ow):  # 图像的每列
            for c in range(ic):  # 图像的每个通道
                for i in range(kh):  # 卷积核的每行
                    for j in range(kw):  # 卷积核的每列
                        out[y, x, c] += image[y + i, x + j, c] * kernel[i, j]
                        
    return out

卷积的计算复杂度¶

卷积的计算复杂度主要取决于乘法的次数，也就是（图像中的像素数）×（核中的单元格数）。

思考：从复杂度的角度来看，为什么小核比大核好？

硬件加速器¶

一些重要的计算可以通过使用专用硬件来加速，最常见加速器是 GPU。

图形处理单元（GPU，Graphical Processing Units）：

  • 最初是作为图像渲染器设计的
  • 在其他非常规律的计算中也相当快速
  • 卷积由于其规律的结构，是一种非常适合 GPU 的操作

恒等¶

恒等返回原始图像，卷积核是

$$ \begin{bmatrix} 0 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 0 \\ \end{bmatrix} $$

In [11]:
identity = np.array([
    0, 0, 0,
    0, 1, 0,
    0, 0, 0
]).reshape(3, 3)

corgis0 = my_convolve(corgis, identity)
plt.imshow(corgis0)
Out[11]:
<matplotlib.image.AxesImage at 0x10ec6fb10>
No description has been provided for this image

结果正确，但是速度有点慢。原因是 Python 语言动态特性导致效率的损失，不适用于大量计算的场景。

在 Jupyter 中，可以用 %timeit 命令来测量运行时间：

In [12]:
%timeit -n1 -r1 my_convolve(corgis, identity)

5.31 s ± 0 ns per loop (mean ± std. dev. of 1 run, 1 loop each)

使用 SciPy 中的卷积实现¶

可以使用 SciPy 软件包中的卷积实现来获得更快的运行速度。

In [13]:
from scipy.ndimage import convolve

corgis0 = convolve(corgis, identity, axes=(0, 1))  # 需要明确指定做卷积的维度
plt.imshow(corgis0)
Out[13]:
<matplotlib.image.AxesImage at 0x10f593890>
No description has been provided for this image

结果一致，速度也很快。SciPy 中的卷积是用 C 语言实现的，然后给出了 Python 调用的接口。

In [14]:
%timeit convolve(corgis, identity, axes=(0, 1))

1.67 ms ± 8.89 μs per loop (mean ± std. dev. of 7 runs, 1,000 loops each)

边界探测¶

边界探测需要考虑四周的情况，卷积核是

$$ \begin{bmatrix} 0 & -1 & 0 \\ -1 & \ \ 4 & -1 \\ 0 & -1 & 0 \end{bmatrix} $$

In [15]:
edge_detect = np.array([
    0, -1, 0,
    -1, 4, -1,
    0, -1, 0
]).reshape(3, 3)

#corgis1 = convolve(corgis, edge_detect, axes=(0, 1))
#plt.imshow(corgis1)

corgis1 = np.abs(convolve(corgis, edge_detect, axes=(0, 1))).sum(axis=2)
plt.imshow(corgis1, cmap='gray')
Out[15]:
<matplotlib.image.AxesImage at 0x10f9b9f90>
No description has been provided for this image

锐化¶

在原始图像的基础上增强边界可以实现锐化效果，卷积核是

$$ \begin{bmatrix} 0 & -1 & 0 \\ -1 & \ \ 5 & -1 \\ 0 & -1 & 0 \\ \end{bmatrix} $$

In [16]:
sharpen = identity + edge_detect
corgis2 = convolve(corgis, sharpen, axes=(0, 1)).clip(0, 1)
plt.imshow(corgis2)
Out[16]:
<matplotlib.image.AxesImage at 0x10fa13ed0>
No description has been provided for this image

模糊¶

跟周围做平均可以达到模糊的效果，卷积核是

$$ \frac{1}{9} \begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \\ \end{bmatrix} $$

In [17]:
box_blur = np.full((3, 3), 1) / 9
corgis3 = convolve(corgis, box_blur, axes=(0, 1))
plt.imshow(corgis3)
Out[17]:
<matplotlib.image.AxesImage at 0x10fac9e50>
No description has been provided for this image

思考：上述操作为什么不需要 clip？

水平方向变化探测¶

对左右两边的像素求差值可以探测水平方向上的变化，卷积核是

$$ \begin{bmatrix} -1 & 0 & 1 \\ -1 & 0 & 1 \\ -1 & 0 & 1 \\ \end{bmatrix} $$

In [18]:
dx = np.array([
    -1, 0, 1,
    -1, 0, 1,
    -1, 0, 1
]).reshape(3, 3) / 2

corgis4 = np.abs(convolve(corgis, dx, axes=(0, 1))).sum(axis=2)
plt.imshow(corgis4, cmap='gray')
Out[18]:
<matplotlib.image.AxesImage at 0x10fb7c050>
No description has been provided for this image

竖直方向变化探测¶

对上下两边的像素求差值可以探测竖直方向上的变化，卷积核是

$$ \begin{bmatrix} -1 & -1 & -1 \\ 0 & 0 & 0 \\ 1 & 1 & 1 \\ \end{bmatrix} $$

In [19]:
dy = dx.T
corgis5 = np.abs(convolve(corgis, dy, axes=(0, 1))).sum(axis=2)
plt.imshow(corgis5, cmap='gray')
Out[19]:
<matplotlib.image.AxesImage at 0x10fbda210>
No description has been provided for this image

卷积的滤镜应用小结¶

卷积可以用来做图像的变换和检测。

思考：卷积核中所有元素之和与对应的操作是变换还是检测有什么关系？

卷积的应用：卷积神经网络¶

在图像识别等任务中，先要识别各个局部，然后再看更全局的信息。这时就可以用卷积来提取局部信息，称为卷积神经网络（Convolutional Neural Networks，CNN）。

image.png

本讲小结¶

更多的图像变换方法：

  • 采样
  • 线性组合
  • 卷积