Python图像处理（一）基本操作

（仅个人学习摘抄）

脚本命名千万千万不要用python里面的库还有什么文件名一样，会出错！！！

1. PIL

注意：打开图片的时候，地址斜杠再Python中有其他含义，所以要避免歧义

更多函数使用：https://effbot.org/imagingbook/image.htm

 

 

2. Matplotlib

2.1 绘制图像、点和线

输出的图片，y 轴的原点是从上面往下的增加的，和其他的坐标系显示有点不一样。因为在 PyLab 库中，左上角是坐标原点。

图像中，点标记和线默认的颜色是蓝色。

不显示坐标轴：axis('off')

matplotlib_1.py

 

2.2 图像轮廓和直方图

绘制图像的轮廓（或者其他二位函数的等轮廓线），首先需要将图像灰度化：

from PIL import Image
from pylab import *

# 读取图像到数组
im = array(Image.open(r'D:\test\pic\test.jpg').convert('L'))

# 新建一个图像
figure()
# 不使用颜色信息
gray()
# 在原点的左上角显示轮廓图像
contour(im, origin='image')
axis('equal')
axis('off')

# 直方图
figure()
hist(im.flatten(),128)
show()

 

　　图像的直方图用来表征图像像素值的分布情况。用一定数目的小区间（bin）来指定表征像素值的范围，每个小区间会得到落入该小区间表示范围的像素数目。图像的直方图可以使用 hist() 函数绘制：

　　　　figure()

　　　　hist(im.flatten(),128)

　　　　show()

hist() 函数的第二个参数指定小区间的数目。因为 hist() 只接受一维数组作为输入，所以在绘制直方图之前，必须先对图像进行压平处理。flatten() 将任意数组按照行优先准则转换成一维数组。

 

2.3 交互式标注

通过 ginput() 函数实现交互式标注

# 交互式标注
from PIL import Image
from pylab import *

im = array(Image.open(r'D:\test\pic\test.jpg'))
imshow(im)
print('Please click 3 plots')
x = ginput(3)
print('you click:', x)
show()

在绘图区域单击三次，将会显示点击的坐标。

 

 

3. NumPy

3.1 图像数组表示

　　调用 array() 方法将图像转换成 NumPy 的数组对象，NumPy 中的数组十多维的，可以用来表示向量、矩阵和图像。

　　每行的第一个元组表示图像数组的大小（行、列、颜色通道），紧接着的字符串表示数组元素的数据类型。因为图像通常被编码成无符号八位整数（uint8），所以在第一种情况下，对图像进行灰度化处理，并且在创建数组时使用额外的参数 “f”；该参数将数据类型转换为浮点型。注意：由于灰度图像没有颜色信息，所以在形状元组中，它只有两个数值。

　　数组中的元素可以使用下标访问。位于坐标 i、j，以及颜色通道 k 的像素值可以如下访问：

　　　　value = im[i, j, k]

　　多个数组元素可以使用数组切片方式访问。切片方式返回的是以指定间隔下标访问该数组的元素值。

如果仅使用一个下标，则该下标为行下标。

 

3.2 灰度变换

array() 变换的相反操作用 PIL 的 fromarray() 函数完成：

　　　　pil_im = Image.fromarray(im)

如果之前一些操作将“uint8”数据类型转换为其他数据类型，那么在创建 PIL 图像之前，需要将数据类型转换回来：

　　　　pil_im = Image.fromarray(uint8(im))

如果不确定输入的数据类型，安全起见，先转换回来。注意，NumPy 总是将数组数据类型转换成能够表示数据的“最低”数据类型。对浮点数做乘积或除法操作会使整数类型的数组变成浮点型。

 

3.3 图像缩放

 

3.4 直方图均衡化

　　直方图均衡化是指将一幅图像的灰度直方图变平，使变换后的图像中每个灰度值的分布概率都相同。直方图均衡化是对图像灰度值进行归一化的非常好的方法，并且可以增强图像的对比度。

　　直方图均衡化的变换函数是图像中像素值的累积分布函数（cumulative distribution function，简称为 cdf，将像素值的范围映射到目标范围的归一化操作）。

　　interp() 函数有两个输入参数，一个是灰度图像，一个是直方图中使用小区间的数目。函数返回直方图均衡化后的图像，和用来做像素值映射的累积分布函数。cdf[-1] 是为了归一化到 0~1 范围。

　　直方图均衡化后图像的对比度增强了，原先图像灰色区域的细节变得清晰。

调用其他 .py 文件中的函数，如上例：

①　import imtools

　　imtools.histeq(im)

②　from imtools import histeq

　　histeq(im)

 

3.5 图像平均

　　图像平均是减少噪声的一种简单的方式。

　　该函数可以直接跳过不能打开的图像。还可以使用 mean() 函数计算平均图像。mean() 函数需要将所有的图像堆积到一个数组中，比较占内存。

 

3.6 图像的主成分分析（PCA）

　　PCA（Principal Component Analysis，主成分分析）是一个非常有用的降维技巧。它可以在使用尽可能少维数的前提下，尽量多的保持训练数据的信息。一幅 100 X 100 像素的小灰度图像，也有 10000维，可以看成 10000 维空间中的一个点。PCA 产生的投影矩阵可以被视为将原始坐标变换到现有的坐标系，坐标系中的各个坐标按照重要性递减排列。

　　为了对图像数据进行 PCA 变换，图像需要转换成一维向量表示，可以用 NumPy 类库中的 flatten() 进行变换。

　　将变平的图像堆积起来，可以得到一个矩阵，矩阵的一行表示一幅图像。在计算主方向之前，所有的行图像按照平均图像进行了中心化。我们通常使用 SVD（Singular Value Decomposition，奇异值分解）方法来计算主成分；但当矩阵维数很大时，SVD 的计算非常慢，此时通常不用 SVD 分解。

PCA 操作代码：

from PIL import Image
from numpy import *

def pca(x):
    ''' 主成分分析：
        输入：矩阵 X，其中该矩阵中存储训练数据，每一行为一条训练数据
        返回：投影矩阵（按照维度的重要性排序）、方差和均值'''

    # 获取维数
    num_data,dim = X.shape

    # 数据中心化
    mean_X = X.mean(axis=0)
    X = X - mean_X

if dim>num_data:
    # PCA-使用紧致技巧
    M = dot(X,X.T)            # 协方差矩阵
    e,EV = linalg.eigh(M)    # 特征值和特征向量
    tmp = dot(X.T,EV).T     # 这就是紧致技巧
    V = tmp[::-1]            # 由于最后的特征向量是我们所需要的，所以需要将其逆转
    S = aqrt(e)[::-1]        # 由于特征值是按照递增顺序排列的，所以需要将其逆转
    for i in range(V.shape[1]):
        V[:i] /= S
    else:
        # PCA- 使用 SVD 方法
        U,S,V = linalg.svd(X)
        V = V[:num_data]    # 仅仅返回前 num_data 维数据才合理

# 返回投影矩阵、方差和均值
return V,S,mean_X

　　该函数首先通过减去每一维的均值将数据中心化，然后计算协方差矩阵对应最大特征值的特征向量，此时可以使用简明的技巧或者 SVD 分解。这里使用 range() 函数，该函数的输入参数为一个整数 n，函数返回整数 0~(n-1) 的一个列表。也可以使用 arange() 函数返回一个数组，xrange() 函数返回一个产生器（可能会提升速度）。

　　如果数据个数小于向量的维数，不用 SVD 分解，而是计算维数更小的协方差矩阵 XXT 的特征向量。通过计算对应前 k（k 是降维后的维数）最大特征值的特征向量，可以使上面的 PCA 操作更快。

　　对图像进行 PCA 变换，图像的名称保存在列表 imlist 中，调用 pca.py 文件。

　　图像从一维表示转换成二维图像，使用 reshape() 函数。

 

3.7 pickle 模块

　　pickle 模块可以接受几乎所有的 Python 对象，并将其转换成字符串表示，该过程叫做封装（pickling）。从字符串表示中重构该对象，称为拆封（unpickling）。这些字符串表示可以方便的存储和传输。

例如保存上面的图像的平均图像和主成分：

在上例中，许多对象可以保存在同一个文件中。pickle 模块中有很多不同的协议可以生成 ,pkl 文件；如果不确定的话，最好以二进制文件的形式读取和写入。在其他 Python 会话中载入数据，只需要使用 load() 方法：

　　注意，载入对象的顺序必须和先前保存的一样。Python 中有一个用 C 语言写的优化版本，叫做 cpickle 模块，该模块和标准 pickle 模块完全兼容。

　　使用 with 语句处理文件的读写操作，可以自动打开和关闭文件（即使在文件打开时发生错误）。

 

 

4. SciPy

4.1 图像模糊

　　图像的高斯模糊是非常经典的图像卷积例子。本质上，图形模糊就是将（灰度）图像 I 和一个高斯核进行卷积操作：

　　　　　　　　I_σ = I*G_σ

其中 * 表示卷积操作；G_σ 是标准差为 σ 的二维高斯核，定义为：

　　　　　　　　

高斯模糊通常是图像处理操作的一部分，比如图像插值操作，兴趣点计算等。

　　SciPy 有用来做滤波操作的 scipy.ndimage.filters 模块。该模块使用快速一维分离的方式来计算卷积。

代码：scipy_1.py 

 

 

4.2 图像导数

　　图像强度变化情况可以用灰度图像 I 的 x 和 y 方向导数 I_x 和 I_y 描述。

　　图像导数大多数可以通过卷积简单的实现：

　　　　　　I_x=I*D_x 和 I_y = I*D_y

对于 Dx 和 Dy，通常选择 Prewitt 滤波器：

　　　　

Sobel 滤波器：

　　　　

 

from PIL import Image
from numpy import *
from pylab import *
from scipy.ndimage import filters

im = array(Image.open(r'D:\test\pic\test.jpg').convert('L'))

# Sobel 导数滤波器
imx = zeros(im.shape)
filters.sobel(im,1,imx)

imy = zeros(im.shape)
filters.sobel(im,0,imy)

magnitude = sqrt(imx**2+imy**2)

 

　　Sobel() 函数第二个参数用来表示选择 x 或者 y 方向导数（1：x 方向；0：y 方向），第三个参数保存输出的变量。正导数显示为亮的像素，负导数显示为暗的像素，灰色区域表示导数的值接近于零。

缺点：滤波器的尺寸需要随着图像分辨率的变化而变化。为了在噪声方面更加稳健，以及在任意尺度上计算导数，使用高斯导数滤波器。

　　　　　　　　I_x=I*G_σx 和 I_y = I*G_σy

其中，G_σx 和 G_σy 表示 G_σ 在 x 和 y 方向上的导数，G_σ 为标准差为 σ 的高斯函数。

　　该函数的第三个参数指定对每个方向计算哪种类型的导数，第二个参数为使用的标准差。

 

4.3 形态学：对象计数

　　形态学（或数学形态学）是度量和分析基本形状的图像处理方法的基本框架与集合。形态学通常用于处理二值图像，但是也能够用于灰度图像。二值图像是指图像的每个像素只能取两个值，通常是 0 和 1。二值图像通常是，在计算物体的数目，或者度量其大小时，对一幅图像进行阈值化后的结果。

 

4.4 一些有用的 SciPy 模块

1、读写 .mat 文件

读取 Matlab 的 .mat 文件，使用 scipy.io 模块：

　　　　data = scipy.io.loadmat('test.mat')

data 对象包含一个字典，字典中的键对应于保存在原始 ,mat 文件中的变量名。由于这些变量是数组格式的，因此可以很方便的保存到 ,mat 文件中。只需创建一个字典（其中包含你想要保存的所有变量），然后使用 savemat() 函数：

　　　　data = {}

　　　　data['X'] = x

　　　　scipy.io.savemat('test.mat',data)

上面保存的是数组 x，所以当读到 Matlab 中时，变量的名字仍为 x。

2、以图象形式保存数组

imsave() 函数可以从 scipy.misc 模块中载入。要将数组 im 保存到文件中，可以使用下面的命令：

　　　　from scipy.misc import imsave

　　　　imsave('test.jpg',im)

scipy.misc 模块同样包含了著名的 Lena 测试图像：

　　　　lena = scipy.misc.lena()

返回一个 512X512 的灰度图像数组。

 

 

5. 图像去噪

　　ROF（Rudin-Osher-Fatemi）去噪模型。ROF 模型具有很好的性质：处理后的图像更平滑，同时保持图像边缘和结构信息。

ROF 模型去噪：

 

from numpy import *

def denoise(im,U_init,tolerance=0.1,tau=0.125,tv_weight=100):
    ''' 使用 A.Chanbolle(2005)在公式(11)中的计算步骤实现 Rudin-Osher-Fatemi(ROF)去噪模型
        输入：含有噪声的输入图像（灰度图像）、U的初始值、YV正则项权值、步长、停业条件
        输出：去噪和去除纹理后的图像、纹理残留'''

    m,n = im.shape        # 噪声图像的大小

    # 初始化
    U = U_init
    Px = im         # 对偶域的 x 分量
    Py = im         # 对偶域的 y 分量
    error = 1

    while(error > tolerance):
        Uold = U

        # 原始变量的梯度
        GradUx = roll(U,-1,axis=1)-U         # 变量 U 梯度的 x 分量
        GradUy = roll(U,-1,axis=0)-U         # 变量 U 梯度的 y 分量

        # 更新对偶变量
        PxNew = Px + (tau/tv_weight)*GradUx
        PyNew = Py + (tau/tv_weight)*GradUy
        NormNew = maximum(1,sqrt(PxNew**2+PyNew**2))

        Px = PxNew/NormNew        # 更新 x 分量（对偶）
        Py = PyNew/NormNew        # 更新 y 分量（对偶）

        # 更新原始变量
        RxPx = roll(Px,1,axis=1)        # 对 x 分量进行向右 x 轴平移
        RyPy = roll(Py,1,axis=0)        # 对 y 分量进行向右 y 轴平移

        DivP = (Px-RxPx)+(Py-RyPy)        # 对偶域的散度
        U = im + tv_weight*DivP            # 更新原始变量

        # 更新误差
        error = linalg.norm(U-Uold)/sqrt(n*m)

    return U,im-U            # 去噪后的图像和纹理残留

 

利用一个合成图像调用该函数：

from PIL import Image
from numpy import *
from numpy import random
from scipy.ndimage import filters
import rof

# 使用噪声创建合成图像
im = zeros((500,500))
im[100:400,100:400] = 128
im[200:300,200:300] = 255
im = im + 30*random.standard_normal((500,500))

G = filters.gaussian_filter(im,10)
U,T = rof.denoise(im,im)

im = Image.fromarray(im)
im.show()
G = Image.fromarray(G)
G.show()
U = Image.fromarray(U)
U.show()

ROF 模型去噪后的图像保留了边缘和图像的结构信息，同时模糊了“噪声”。

 

 

学习书目《Python计算机视觉编程》