Programming Computer Vision with Python （学习笔记二）

摘要：首先介绍跟图像处理显示有关两个库和，然后介绍增强图像对比度的实现原理。直方图均衡化就是为了达到这个目的，均衡化后的图像，像素落在每个灰度级上的个数是相等的。

首先介绍跟图像处理、显示有关两个库：NumPy和Matplotlib，然后介绍增强图像对比度的实现原理。

NumPy是Python用于科学计算的基础库，提供了一些很有用的概念，如：N维数组对象，可用于表示向量、矩阵、图像数据等，另外还包含了线性代数及其运算函数。NumPy的数组对象在本书示例中会被大量使用，它可以作诸如矩阵乘法、变换、向量乘法和正态化等运算，我们通过这些运算来实现图像对齐、图像分类、图像扭转等。
这是一个基础库，通常不需要额外安装。

N维数组在NumPy中对应的数据类型是ndarry，有时使用别名array（即numpy.array）。但要注意的是，它与Python的内置类型array是两回事，不要混淆，Python内置array类型只处理一维数组，其功能远不及ndarray。ndarray中的所有元素的存储类型是一样的，下面对ndarray一些重要的属性进行说明：

ndarray.ndim
数组维度

ndarray.shape
对于一个n×m矩阵，shape返回元组(n,m)

ndarray.size
数组的所有元素个数

ndarray.dtype
数组元素的数据类型

ndarray.itemsize
数据中每个元素的类型长度（单位byte）

ndarray.data
包含数组所有元素的buffer，通常我们只是使用数组下标来获取元素的值

构造
用Python的数组表示来构造ndarray，很直观：

In [3]: import numpy as np

In [5]: a = np.array([[0,1,2], 
                      [3,4,5]]) 

In [6]: a.shape
Out[6]: (2, 3)

In [7]: a.ndim
Out[7]: 2

In [8]: a.dtype.name
Out[8]: "int64"

In [9]: a.itemsize
Out[9]: 8

In [10]: a.size
Out[10]: 6

In [11]: type(a)
Out[11]: numpy.ndarray

构建dnarray时可以指定元素的类型：

In [12]: b = np.array([0,1,2],dtype=int16)

In [13]: b.itemsize
Out[13]: 2

我们最常用的是想把一幅图像转为np.array表示，而PIL的Image类可以处理大部分的图像格式，所以从Image转为np.array很有用，如：

from PIL import Image
import numpy as np
im = np.array(Image.open("Selection_001.png"))

注：Image对象之所以能直接转为ndarray类型，是因为Image类实现了ndarray的data和shape等接口。

其它一些有用的构造方法：

np.zeros( (n, m) ) 构建n乘m数组，其中元素初始化为0

np.ones( (n, m) ) 同上，但元素初始化为1

np.empty( (n, m) ) 同上，但元素不作初始化

np.arange([start,] stop[, step,], dtype=None) 构建1维数组，元素的值从start到stop，增加步长为step

In [75]: np.arange(5)
Out[75]: array([0, 1, 2, 3, 4])

In [76]: np.arange(5, 10)
Out[76]: array([5, 6, 7, 8, 9])

In [77]: np.arange(5, 10, 2)
Out[77]: array([5, 7, 9])

np.linspace( start, stop, item_count ) 构建1维数组，元素从start到stop，元素个数为item_count，所以元素的增加步长是自动计算的： (to - from) / （item_count - 1)

In [63]: np.linspace(5,10,2)
Out[63]: array([  5.,  10.])

In [64]: np.linspace(5,10,3)
Out[64]: array([  5. ,   7.5,  10. ])

In [65]: np.linspace(5,10,4)
Out[65]: array([  5.        ,   6.66666667,   8.33333333,  10.        ])

In [66]: np.linspace(5,10,5)
Out[66]: array([  5.  ,   6.25,   7.5 ,   8.75,  10.  ])

基本运算
两个数组的+-<>*运算，作用于两个数组相对应位置的元素，结果是一个新数组：

In [22]: a
Out[22]: 
array([[1, 2, 3],
       [4, 5, 6]])

In [23]: b
Out[23]: 
array([[ 1.,  1.,  1.],
       [ 1.,  1.,  1.]])

In [24]: a + b
Out[24]: 
array([[ 2.,  3.,  4.],
       [ 5.,  6.,  7.]])

In [25]: a - b
Out[25]: 
array([[ 0.,  1.,  2.],
       [ 3.,  4.,  5.]])

In [26]: a < b
Out[26]: 
array([[False, False, False],
       [False, False, False]], dtype=bool)

In [30]: c
Out[30]: 
array([[1, 1, 1],
       [2, 2, 2]])

In [31]: a * c
Out[31]: 
array([[ 1,  2,  3],
       [ 8, 10, 12]])

数组A与B的乘积：A.dot(B)或np.dot(A, B)。
对+=和*=等运算符产生的结果，直接修改调用数组自身，而不是返回新数组。
其它一些有用的运算操作：np.sin, np.cos, np.exp（指数), np.sqrt（开方）等。

下标访问

In [45]: a
Out[45]: 
array([[ 0,  1,  2,  3],
       [10, 11, 12, 13],
       [20, 21, 22, 23],
       [30, 31, 32, 33],
       [40, 41, 42, 43]])

In [46]: a[2,3]  #访问行下标为2，列下标为3的元素
Out[46]: 23

In [47]: a[0:5, 1] #访问行下标从0到5（不含），列下标为1的元素
Out[47]: array([ 1, 11, 21, 31, 41])

In [50]: a[:, 1] #访问所有行，但列下标为1的元素
Out[50]: array([ 1, 11, 21, 31, 41])

In [51]: a[1:3] #访问行下标从1到3（不含）的元素
Out[51]: 
array([[10, 11, 12, 13],
       [20, 21, 22, 23]])

In [52]: a[-1] #访问最后一行
Out[52]: array([40, 41, 42, 43])

In [2]: x = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
In [3]: x[1:7:2]  #指定起始、结束（不含）以及步长
Out[3]: array([1, 3, 5])

变形
展开为一维数组：

In [53]: a = np.array([[1,2],[3,4]])  #2乘2数组

In [54]: a
Out[54]: 
array([[1, 2],
       [3, 4]])

In [57]: b = a.ravel()  #展开为1维数组，返回新数组

In [58]: b
Out[58]: array([1, 2, 3, 4])

In [59]: b.reshape(2, 2) #变形为2乘2数组，返回新数组
Out[59]: 
array([[1, 2],
       [3, 4]])

In [60]: b.resize(2, 2)  #变形为2乘2数组，直接修改本身

In [61]: b
Out[61]: 
array([[1, 2],
       [3, 4]])

有了以上的了解，我们来看看实际的应用例子。先读取一张图片，把它转为ndarray类型，再看其数组属性：

In [88]: from PIL import Image
In [89]: import numpy as np

In [91]: im = np.array(Image.open("Selection_001.png"))  #用PIL.Image读取图像，并转为ndarray数组
In [92]: print im.shape, im.dtype
(240, 568, 3) uint8  #表示图像数据240行，568列，颜色通道数3，以uint8类型存储

In [93]: im_l = np.array(Image.open("Selection_001.png").convert("L"))  #转为灰度图像
In [94]: print im_l.shape, im_l.dtype
(240, 568) uint8  #灰度图像没有颜色通道信息

矩阵
class numpy.matrix(data, dtype=None, copy=True)
从data中构造一个矩阵对象，data可以是ndarray也可以是字符串，若data为ndarray，则copy表示是否复制data来构造。

In [4]: np.matrix("1 2; 3 4")
Out[4]: 
matrix([[1, 2],
        [3, 4]])

In [5]: np.matrix([[1, 2], [3, 4]])
Out[5]: 
matrix([[1, 2],
        [3, 4]])

还可以使用以下两个函数来构造矩阵：
numpy.mat(data, dtype=None)或numpy.asmatrix(data, dtype=None)，两个只是名字不一样，都相当于numpy.matrix(data, copy=False)。

矩阵类提供了一些矩阵运算的方便的接口，如：
getT：返回转置矩阵
getI: 返回可逆矩阵的逆矩阵
getH：返回共轭转置矩阵
getA：返回矩阵的ndarray

Matplotlib是一个用于科学计算及制图方面的强大的开源库，支持很多常见的图形图表，如：

虽然Matplotlib功能很强大，我们可能只是用到它很少的一些接口，比如画图像的轮廓和灰度图像的柱状图。
安装Matplotlib

sudo apt-get install python-matplotlib

pylab和pyplot
为简化画图工作，Matplotlib的pyplot模块提供了与MATLAB相似的接口，并且可以跟IPython配合使用。
需要注意的是，书中的代码示例使用的是Matplotlib.pylab这个模块：

from PIL import Image
from pylab import *
im = array(Image.open("empire.jpg"))  #读图并转为ndarray
imshow(im)

根据Matplotlib官网上的pyplot和pylab的关系说明得知：使用pylab只是为了import时方便起见，import pylab相当于import了pyplot和numpy模块中大部分的接口，虽然有些例子还这样用，但已经不被推荐使用，而是推荐使用pyplot。另外，pyplot模块内置了状态机，它能自动生成必要的图例和坐标轴等信息，可以简化画图代码。

对图像进行灰度变换的目的是为了：

改善画质，使图像更加清晰

有选择地突出图像中感兴趣的特征或抑制图像中某些不需要的特征，使图像与视觉响应特性相匹配

改变图像的直方图分布，增加图像对比度

最简单的灰度变换就是反转颜色，示例：

In [88]: from PIL import Image
In [89]: import numpy as np
In [90]: import matplotlib.pyplot as plt

In [97]: im = np.array(Image.open("cover.png").convert("L"))
In [98]: plt.gray()  #不加的话，显示出来的图像会有颜色
In [100]: plt.imshow(im)
In [102]: plt.show()

In [103]: im2 = 255 - im
In [104]: plt.imshow(im2)
In [105]: plt.show()

反转前：

反转后：

灰度变换的一个很有用的例子就是直方图均衡化，这里的直方图指图像的灰度直方图，因为我们要示例的是灰度图像，每个像素用8bit表示，值从0到255，共有256个灰度级。但通常的图像像素值，都没有完全占用这256个级别，很多像素的灰度值集中在一起，这样导致灰度之间的变化不明显，如果我们把图像的灰度级按比例拉伸到256级，可以使得像素灰度级差距增大，从而使图像看起来更清晰，对比度更强一些。直方图均衡化就是为了达到这个目的，均衡化后的图像，像素落在每个灰度级上的个数是相等的。而且原图像的第i个灰度累积和（即落在[0,i]区间所有像素个数）与均衡化后的第i个灰度累积和相等，即原图像累积和按0到255的比例进行变换。所以下面将使用累积分布函数（cumulative distribution function，简称cdf）。

直方图数据的统计将借助numpy.histogram函数来获得：

numpy.histogram(a, bins=10, range=None, normed=False, weights=None, density=None)

传入数组及直方图的柱的数目（柱也可由X轴点的系列指定），统计落在各个柱区间的元素的个数。

参数：
a: 数组，需要扁平化
bins： bin指的是直方图中的“柱”，取值对应X轴上的区间[x,y)，此参数可选，传入int表示等宽柱的数量，也支持非等宽柱的设置
range：(float, float），可选，指定柱的最低和最高值
normed：bool，可选，NumPy1.6弃用，建议使用density参数
density：bool，可选，False表示函数返回的是落在每个柱区间的元素的数量，若为True，函数返回的是由`概率密度分布函数`对每个柱计算出来的值

返回值：
hist：ndarray，如density参数所说
bin_edges：柱的边界数组，length(hist) + 1，即X轴上柱之间的分割点形成的数组

示例：
In [8]: a = np.array([0,1,2,3,4])
In [9]: np.histogram(a, 5)
Out[9]: (array([1, 1, 1, 1, 1]), #a中落在以下各个区间的元素的个数
array([ 0. ,  0.8,  1.6,  2.4,  3.2,  4. ])) #柱的边界（区间），自动均分

In [10]: np.histogram(a, 5, density=True)
Out[10]: (array([ 0.25,  0.25,  0.25,  0.25,  0.25]), #概率密度分布
 array([ 0. ,  0.8,  1.6,  2.4,  3.2,  4. ]))

而累积和的计算需要用到numpy.cumsum函数：

numpy.cumsum(a, axis=None, dtype=None, out=None)
示例：
In [21]: a = np.array([1,2,3,4,5,6])
In [22]: np.cumsum(a)
Out[22]: array([ 1,  3,  6, 10, 15, 21])

现在来写一个函数实现直方图均衡化：

import numpy as np
from PIL import Image
import matplotlib.pyplot as plt

def histeq(im,nbr_bins=256):
    imhist,bins = np.histogram(im.flatten(),nbr_bins,density=True) #对每个元素求概率密度
    cdf = imhist.cumsum() #对概率密度数组求累积和
    cdf = 255 * cdf / cdf[-1] #累积和变换到0-255区间
    im2 = np.interp(im.flatten(),bins[:-1],cdf) #线性插值
    return im2.reshape(im.shape), cdf #还原图像维度

im = np.array(Image.open("hist-sample.jpg").convert("L"))
im2,cdf = histeq(im)

plt.gray()
plt.subplot(221) #2行2列，第1个图
plt.imshow(im)
plt.subplot(222) #2行2列，第2个图
plt.hist([x for x in im.flatten() if x < 250], 128)
plt.subplot(223)
plt.imshow(im2)
plt.subplot(224)
plt.hist([x for x in im2.flatten() if x < 250], 128)
plt.show()

效果对比如下，上面的是原图及直方图，下面的是均衡化后的图及直方图：

明显看出，均衡化后的图对比度要更强一些。

多图像平均法是一个用于降噪和美化图片的简单方法。假设多张图像具有相同尺寸，一个计算方法就是把所有图像的数据相加起来再除以图像数目从而得到图像的平均值。这个操作使用ndarray的+=和/=运算符就可以完成。
另一个实现的方法就是使用numpy.mean()函数，放在后面再讲。

下一个笔记内容讲图像的主成分分析(PCA)。