Python计算机视觉 第8章-图像内容分类

Python计算机视觉 第8章-图像内容分类

8.1 K邻近分类法（KNN）

在分类方法中，最简单且用得最多的一种方法之一就是 KNN（K-Nearest Neighbor ,K邻近分类法），这种算法把要分类的对象（例如一个特征向量）与训练集中已知类标记的所有对象进行对比，并由 k 近邻对指派到哪个类进行投票。这种方法通常分类效果较好，但是也有很多弊端：与 K-means 聚类算法一样，需要预先设定 k 值，k 值的选择会影响分类的性能；此外，这种方法要求将整个训练集存储起来，如果训练集非常大，搜索起来就非常慢。对于大训练集，采取某些装箱形式通常会减少对比的次数 1 从积极的一面来看，这种方法在采用何种距离度量方面是没有限制的；实际上，对于你所能想到的东西它都可以奏效，但这并不意味着对任何东西它的分类性能都很好。另外，这种算法的可并行性也很一般。

KNN 算法步骤

1. 选择 K 值：选择一个整数 ( K )，表示在分类时考虑的邻居数。( K ) 的值通常是一个小的正整数。

2. 计算距离：对于待分类的样本，计算它与训练集中每个样本的距离。常用的距离度量包括欧几里得距离、曼哈顿距离等。

   • 欧几里得距离：
     $$d=\sqrt{\sum _{i=1}^n(x_i-y_i{)}^2}$$
     其中 ( x_i ) 和 ( y_i ) 是样本在第 ( i ) 维的值。

   • 曼哈顿距离：
     $$d=\sum _{i=1}^n|x_i-y_i|$$

3. 找到 K 个最近邻：根据计算出的距离，选择距离待分类样本最近的 ( K ) 个训练样本。

4. 进行投票（分类任务）：对于分类任务，对 ( K ) 个最近邻的类别进行投票，选择出现次数最多的类别作为待分类样本的类别。

5. 回归预测（回归任务）：对于回归任务，计算 ( K ) 个最近邻的标签的平均值或加权平均值，作为待分类样本的预测值。

KNN 算法的优缺点

优点：

• 简单易懂：KNN 是一种非常直观的算法，容易理解和实现。
• 无训练阶段：KNN 不需要显式的训练阶段，适合处理一些小规模数据集。

缺点：

• 计算开销大：每次进行分类时都需要计算距离，计算开销随着数据量增加而显著增加。
• 内存消耗大：需要存储所有训练样本，内存消耗较大。
• 受特征尺度影响：特征的尺度会影响距离计算结果，需要进行特征缩放（例如标准化）。

以下是一段示例代码：

import matplotlib.pyplot as plt
import numpy as np
from sklearn.datasets import load_iris
from sklearn.model_selection import train_test_split
from sklearn.neighbors import KNeighborsClassifier
from sklearn import metrics

# 加载数据
iris = load_iris()
X = iris.data
y = iris.target
feature_names = iris.feature_names
target_names = iris.target_names

# 选择前两个特征以便于可视化
X = X[:, :2]
feature_names = feature_names[:2]

# 划分数据集
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=42)

# 初始化 KNN 分类器
k = 3
knn = KNeighborsClassifier(n_neighbors=k)

# 训练模型
knn.fit(X_train, y_train)

# 进行预测
y_pred = knn.predict(X_test)

# 评估模型
accuracy = metrics.accuracy_score(y_test, y_pred)
print(f'Accuracy: {accuracy:.2f}')

# 绘制决策边界
x_min, x_max = X[:, 0].min() - 1, X[:, 0].max() + 1
y_min, y_max = X[:, 1].min() - 1, X[:, 1].max() + 1
xx, yy = np.meshgrid(np.arange(x_min, x_max, 0.02),
                     np.arange(y_min, y_max, 0.02))
Z = knn.predict(np.c_[xx.ravel(), yy.ravel()])
Z = Z.reshape(xx.shape)

# 绘制图像
plt.contourf(xx, yy, Z, alpha=0.3)
plt.scatter(X[:, 0], X[:, 1], c=y, edgecolor='k', s=20)
plt.xlabel(feature_names[0])
plt.ylabel(feature_names[1])
plt.title(f'KNN Classification (k={k})')
plt.colorbar(label='Class')
plt.show()

实验结果如下图所示：

实验图1 实验结果

8.1.2　用稠密SIFT作为图像特征

我们来看如何对图像进行分类。要对图像进行分类，我们需要一个特征向量来表示一幅图像。在聚类一章我们用平均 RGB 像素值和 PCA 系数作为图像的特征向量；这里我们会介绍另外一种表示形式，即稠密 SIFT 特征向量。

可以使用如下代码：

import sift
from PIL import Image
import numpy as np
import os

def process_image_dsift(imagename, resultname, size=20, steps=10,
                        force_orientation=False, resize=None):
    """ 用密集采样的 SIFT 描述子处理一幅图像，并将结果保存在一个文件中。可选的输入：
        特征的大小 size，位置之间的步长 steps，是否强迫计算描述子的方位 force_orientation
        （False 表示所有的方位都是朝上的），用于调整图像大小的元组 """
    
    # 打开并转换图像为灰度模式
    im = Image.open(imagename).convert('L')
    
    # 根据提供的 resize 参数调整图像大小
    if resize is not None:
        im = im.resize(resize)
    
    # 获取图像尺寸
    m, n = im.size
    
    if imagename[-3:] != 'pgm':
        # 创建一个 pgm 文件
        im.save('tmp.pgm')
        imagename = 'tmp.pgm'
    
    # 创建帧，并保存到临时文件
    scale = size / 3.0
    x, y = np.meshgrid(range(steps, m, steps), range(steps, n, steps))
    xx, yy = x.flatten(), y.flatten()
    frame = np.array([xx, yy, scale * np.ones(xx.shape[0]), np.zeros(xx.shape[0])])
    np.savetxt('tmp.frame', frame.T, fmt='%03.3f')
    
    # 根据是否需要强迫计算描述子的方位，选择相应的命令
    if force_orientation:
        cmmd = f"sift {imagename} --output={resultname} --read-frames=tmp.frame --orientations"
    else:
        cmmd = f"sift {imagename} --output={resultname} --read-frames=tmp.frame"
    
    # 执行命令
    os.system(cmmd)
    print(f'processed {imagename} to {resultname}')

利用类似下面的代码可以计算稠密 SIFT 描述子，并可视化它们的位置：

使用用于定位描述子的局部梯度方向（force_orientation 设置为真），该代码可以在整个图像中计算出稠密 SIFT 特征。图 8-2 显示出了这些位置。

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

# 处理图像并生成 SIFT 描述子
dsift.process_image_dsift('empire.jpg', 'empire.sift', 90, 40, True)

# 读取特征
l, d = sift.read_features_from_file('empire.sift')

# 打开图像
im = np.array(Image.open('empire.jpg'))

# 绘制特征
sift.plot_features(im, l, True)

# 显示图像
plt.show()

图 8-2：在一幅图像上应用稠密 SIFT 描述子的例子

8.1.3　图像分类：手势识别

在这个应用中，我们会用稠密 SIFT 描述子来表示这些手势图像，并建立一个简单的手势识别系统。我们用静态手势（Static Hand Posture）数据库（参见 http://www.idiap.ch/resource/gestures/）中的一些图像进行演示。在该数据库主页上下载数据较小的测试集 test set 16.3Mb

将下载后的所有图像放在一个名为 uniform 的文件夹里，每一类均分两组，并分别放入名为 train 和 test 的两个文件夹中。
用上面的稠密 SIFT 函数对图像进行处理，可以得到所有图像的特征向量。

处理后的结果如下：

图 8-3：6 类简单手势图像的稠密 SIFT 描述子，图像来源于静态手势（Static Hand Posture）数据库

8.2　贝叶斯分类器

另一个简单却有效的分类器是贝叶斯分类器（或称朴素贝叶斯分类器）。贝叶斯分类器是一种基于贝叶斯条件概率定理的概率分类器，它假设特征是彼此独立不相关的（这就是它“朴素”的部分）。贝叶斯分类器可以非常有效地被训练出来，原因在于每一个特征模型都是独立选取的。尽管它们的假设非常简单，但是贝叶斯分类器已经在实际应用中获得显著成效，尤其是对垃圾邮件的过滤。贝叶斯分类器的另一个好处
是，一旦学习了这个模型，就没有必要存储训练数据了，只需存储模型的参数。

该分类器是通过将各个特征的条件概率相乘得到一个类的总概率，然后选取概率最高的那个类构造出来的。

首先让我们看一个使用高斯概率分布模型的贝叶斯分类器基本实现：

import numpy as np

class BayesClassifier(object):
    def __init__(self):
        """ 使用训练数据初始化分类器 """
        self.labels = []  # 类标签
        self.mean = []    # 类均值
        self.var = []     # 类方差
        self.n = 0        # 类别数

    def train(self, data, labels=None):
        """ 在数据 data（n×dim 的数组列表）上训练，标记 labels 是可选的，默认为 0…n-1 """
        if labels is None:
            labels = range(len(data))
        self.labels = labels
        self.n = len(labels)
        for c in data:
            self.mean.append(np.mean(c, axis=0))
            self.var.append(np.var(c, axis=0))

    def classify(self, points):
        """ 通过计算得出的每一类的概率对数据点进行分类，并返回最可能的标记 """
        # 计算每一类的概率
        est_prob = np.array([self.gauss(m, v, points) for m, v in zip(self.mean, self.var)])
        # 获取具有最高概率的索引，该索引会给出类标签
        ndx = est_prob.argmax(axis=0)
        est_labels = np.array([self.labels[n] for n in ndx])
        return est_labels, est_prob

    def gauss(self, mean, var, points):
        """ 计算高斯概率密度函数 """
        coef = 1.0 / np.sqrt(2 * np.pi * var)
        exp = np.exp(-0.5 * ((points - mean) ** 2) / var)
        return coef * exp

该模型每一类都有两个变量，即类均值和协方差。train() 方法获取特征数组列表（每个类对应一个特征数组），并计算每个特征数组的均值和协方差。classify() 方法计算数据点构成的数组的类概率，并选概率最高的那个类，最终返回预测的类标记及概率值，同时需要一个高斯辅助函数：

import numpy as np

def gauss(m, v, x):
    """ 用独立均值 m 和方差 v 评估 d 维高斯分布 """
    if len(x.shape) == 1:
        n, d = 1, x.shape[0]
    else:
        n, d = x.shape

    # 协方差矩阵，减去均值
    S = np.diag(1 / v)
    x = x - m

    # 概率的乘积
    y = np.exp(-0.5 * np.diag(np.dot(x, np.dot(S, x.T))))

    # 归一化并返回
    return y * (2 * np.pi) ** (-d / 2.0) / (np.sqrt(np.prod(v)) + 1e-6)

该函数用来计算单个高斯分布的乘积，返回给定一组模型参数 m 和 v 的概率

将该贝叶斯分类器用于上一节的二维数据，下面的脚本将载入上一节中的二维数据，并训练出一个分类器：

import pickle
import bayes
import imtools
import numpy as np
import matplotlib.pyplot as plt

# 用 Pickle 模块载入二维样本点
with open('points_normal.pkl', 'rb') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

# 训练贝叶斯分类器
bc = bayes.BayesClassifier()
bc.train([class_1, class_2], [1, -1])

# 载入测试数据
with open('points_normal_test.pkl', 'rb') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

# 在某些数据点上进行测试
print(bc.classify(class_1[:10])[0])

# 绘制这些二维数据点及决策边界
def classify(x, y, bc=bc):
    points = np.vstack((x, y))
    return bc.classify(points.T)[0]

imtools.plot_2D_boundary([-6, 6, -6, 6], [class_1, class_2], classify, [1, -1])
plt.show()

该脚本会将前 10 个二维数据点的分类结果打印输出到控制台，输出结果如下：

输出样例图

我们再次用一个辅助函数 classify() 在一个网格上评估该函数来可视化这一分类结果。两个数据集的分类结果如图 8-4 所示；该例中，决策边界是一个椭圆，类似于二维高斯函数的等值线。

图 8-4：用贝叶斯分类器对二维数据进行分类。每个例子中的颜色代表了类标记。正确分类的点用星号表示，误错分类的点用圆点表示，曲线是分类器的决策边界

8.3 支持向量机（SVM）

1. 理论基础

支持向量机的主要目标是找到一个最佳的决策边界（或称为超平面），以将不同类别的数据点分开。其基本思想是：

• 决策边界：在特征空间中，SVM 寻找一个超平面，将不同类别的数据点分开，并且具有最大化类别间隔的特性。
• 最大间隔：SVM 通过最大化超平面与最近数据点之间的间隔来实现分类，这些最近的数据点被称为支持向量。

2. 数学公式

假设我们有一个训练集 ( {(\mathbf{x}i, y_i)}{i=1}^n )，其中 ( \mathbf{x}_i ) 是输入特征向量，( y_i ) 是标签（通常是 +1 或 -1）。我们的目标是找到一个超平面：

$${\mathbf{w}}^T\mathbf{x}+b=0$$

使得间隔最大化。

最大化间隔的数学表达式：

$${\text{max}}_{\mathbf{w},b}\frac{2}{\Vert \mathbf{w}\Vert }$$

约束条件：

$$y_i({\mathbf{w}}^T{\mathbf{x}}_i+b)\ge 1,\forall i$$

其中，( \mathbf{w} ) 是超平面的法向量，( b ) 是偏置项。

为了求解这个优化问题，通常会转化为一个凸优化问题：

$${\text{min}}_{\mathbf{w},b}\frac{1}{2}\Vert \mathbf{w}{\Vert }^2$$

约束条件：

$$y_i({\mathbf{w}}^T{\mathbf{x}}_i+b)\ge 1,\forall i$$

3. 软间隔与核技巧

在实际应用中，数据往往不是线性可分的，因此引入了软间隔（Soft Margin）的概念：

$${\text{min}}_{\mathbf{w},b,{\xi }_i}\frac{1}{2}\Vert \mathbf{w}{\Vert }^2+C\sum _{i=1}^n{\xi }_i$$

约束条件：

$$y_i({\mathbf{w}}^T{\mathbf{x}}_i+b)\ge 1-{\xi }_i$$
$${\xi }_i\ge 0$$

其中，( \xi_i ) 是松弛变量，用于允许一定的分类错误，( C ) 是惩罚因子，控制错误分类的容忍度。

对于非线性可分的数据，SVM 引入了核函数（Kernel Function），将数据映射到高维空间，使其在高维空间中线性可分。常见的核函数包括：

• 线性核：$$K(\mathbf{x},{\mathbf{x}}^{\prime })={\mathbf{x}}^T{\mathbf{x}}^{\prime }$$
• 多项式核：$$K(\mathbf{x},{\mathbf{x}}^{\prime })=({\mathbf{x}}^T{\mathbf{x}}^{\prime }+c{)}^d$$
• 高斯径向基核（RBF 核）：$$K(\mathbf{x},{\mathbf{x}}^{\prime })=\exp \left(-\frac{\Vert \mathbf{x}-{\mathbf{x}}^{\prime }{\Vert }^2}{2{\sigma }^2}\right)$$

8.3.1　使用LibSVM

LibSVM是最好的、使用最广泛的 SVM 实现工具包。LibSVM 为 Python 提供了一个良好的接口（也为其他编程语言提供了接口）。

下面的脚本会载入在前面kNN 范例分类中用到的数据点，并用径向基函数训练一个 SVM 分类器：

import pickle
from svmutil import *

# 用 Pickle 载入二维样本点
with open('points_normal.pkl', 'rb') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

# 转换成列表，便于使用 libSVM
class_1 = list(class_1)
class_2 = list(class_2)
labels = list(labels)
samples = class_1 + class_2  # 连接两个列表

# 创建 SVM 训练数据
prob = svm_problem(labels, samples)

# 设置 SVM 参数（使用 RBF 核函数）
param = svm_parameter('-t 2')

# 在数据上训练 SVM
m = svm_train(prob, param)

# 在训练数据上分类效果如何？
res = svm_predict(labels, samples, m)

我们用与前面一样的方法载入数据集，但是这次需要把数组转成列表，因为LibSVM 不支持数组对象作为输入。这里，我们用 Python 的内建函数 map() 进行转换，map() 函数中用到了对角一个元素都会进行转换的 list() 函数。紧接着我们创建了一个 svm_problem 对象，并为其设置了一些参数。调用 svm_train() 求解该优化问题用以确定模型参数，然后就可以用该模型进行预测了。最后一行调用 svm_predict()，用求得的模型 m 对训练数据分类，并显示出在训练数据中分类的正确率，打印输出结果如下：

样例输出结果

结果表明该分类器完全分开了训练数据，400 个数据点全部分类正确。

载入其他数据集，并对该分类器进行测试，将数据转成列表，结果如图8-5所示：

图 8-5：用支持向量机 SVM 对二维数据进行分类。在这两幅图中，我们用不同颜色标识类标记。正确分类的点用星号表示，错误分类的点用圆点表示，曲线是分类器的决策边界

8.3.2　再论手势识别

在多类手势识别问题上使用 LibSVM 相当直观。LibSVM 可以自动处理多个类，我们只需要对数据进行格式化，使输入和输出匹配 LibSVM 的要求。

下面的代码会载入训练数据测试数据：

features = list(features)
test_features = list(test_features)

# 为标记创建转换函数
transl = {}
for i, c in enumerate(classnames):
    transl[c], transl[i] = i, c

# 创建 SVM 训练数据
prob = svm_problem(convert_labels(labels, transl), features)
param = svm_parameter('-t 0')  # 使用线性核函数

# 在数据上训练 SVM
m = svm_train(prob, param)

# 在训练数据上分类效果如何
res = svm_predict(convert_labels(labels, transl), features, m)

# 测试 SVM
res = svm_predict(convert_labels(test_labels, transl), test_features, m)[0]
res = convert_labels(res, transl)

与之前一样，我们调用 map() 函数将数组转成列表；因为 LibSVM 不能处理字符串标记，所以这些标记也需要转换。字典 transl 会包含一个在字符串和整数标记间的变换。你可以试着在控制台上打印该变换，看看其对应变换关系。参数 -t 0 设置分类器是线性分类器，决策边界在 10 000 维特征原空间中是一个超平面。

现在，对标记进行比较：

acc = sum(1.0*(res==test_labels)) / len(test_labels)
print 'Accuracy:', acc
print_confusion(res,test_labels,classnames)

用线性核函数得出的分类结果如下：

分类结果

用 PCA 将维数降低到 50，分类正确率变为：

可以看出，当特征向量维数降低到原空间数据维数的 1/200 时，结果并不差。

8.4　光学字符识别

作为一个多类问题实例，让我们来理解数独图像。OCR（Optical Character Recognition，光学字符识别）是一个理解手写或机写文本图像的处理过程。一个常见的例子是通过扫描文件来提取文本，例如书信中的邮政编码或者谷歌图书（http://books.google.com/）里图书馆卷的页数。这里我们看一个简单的在打印的数独图形中识别数字的OCR 问题。数独是一种数字逻辑游戏，规则是用数字 1-9 填满 9×9 的网格，使每一行每一列和每个 3×3 的子网格包含数字1-9。

8.4.1　训练分类器

对于这种分类问题，我们有 10 个类：数字 1…9，以及一些什么也没有的单元格。我们给定什么也没有的单元格的类标号是 0，这样所有类标记就是 0…9。我们会用已经剪切好的数独单元格数据集来训练一个 10 类的分类器文件 sudoku_images.zip中有两个文件夹“ocr data”和“sudokus”，后者包含了不同条件下的数独图像集，我们稍后讲解。

8.4.2　选取特征

我们首先要确定选取怎样的特征向量来表示每一个单元格里的图像。有很多不错的选择；这里我们将会用某些简单而有效的特征。输入一个图像，下面的函数将返回一个拉成一组数组后的灰度值特征向量：

def compute_feature(im):
 """ 对一个 ocr 图像块返回一个特征向量 """
 # 调整大小，并去除边界
 norm_im = imresize(im,(30,30))
 norm_im = norm_im[3:-3,3:-3]
 return norm_im.flatten()

compute_feature() 函数用到 imtools 模块中的尺寸调整函数 imresize()，来减少特征向量的长度。我们还修剪掉了大约 10% 的边界像素，因为这些修剪掉的部分通常是网格线的边缘部分，如图 8-6 所示。

图 8-6：用于训练 10 类数独 OCR 分类器的训练样本图像

现在我们用下面的函数来读取训练数据：

def load_ocr_data(path):
    """返回路径中所有图像的标记及 OCR 特征"""
    # 对以 .jpg 为后缀的所有文件创建一个列表
    imlist = [os.path.join(path, f) for f in os.listdir(path) if f.endswith('.jpg')]
    # 创建标记
    labels = [int(imfile.split('/')[-1][0]) for imfile in imlist]
    # 从图像中创建特征
    features = []
    for imname in imlist:
        im = array(Image.open(imname).convert('L'))
        features.append(compute_feature(im))
    return array(features), labels

上述代码将每一个 JPEG 文件的文件名中的第一个字母提取出来做类标记，并且这些标记被作为整型数据存储在 lables 列表里；用上面的函数计算出的特征向量存储在一个数组里。

8.4.3　多类支持向量机

在得到了训练数据之后，我们接下来要学习一个分类器，这里将使用多类支持向量机。代码和上一节中的代码类似：

from svmutil import *

# 训练数据
features, labels = load_ocr_data('training/')
# 测试数据
test_features, test_labels = load_ocr_data('testing/')

# 训练一个线性 SVM 分类器
features = map(list, features)
test_features = map(list, test_features)
prob = svm_problem(labels, features)
param = svm_parameter('-t 0')
m = svm_train(prob, param)

# 在训练数据上分类效果如何
res = svm_predict(labels, features, m)

# 在测试集上表现如何
res = svm_predict(test_labels, test_features, m)

该代码会训练出一个线性 SVM 分类器，并在测试集上对该分类器的性能进行测试，你可以通过调用最后两个 svm_predict() 函数得到以下输出结果：

样例输出结果

训练集中的 1409 张图像在 10 类中都被完美地分准确了，在测试集上识别性能也在 99% 左右。

8.4.4　提取单元格并识别字符

有了识别单元格内容的分类器后，下一步就是自动地找到这些单元格。一旦我们解决了这个问题，就可以对单元格进行裁剪，并把裁剪后的单元格传给分类器。我们假设数独图像是已经对齐的，其水平和垂直网格线平行于图像的边，如图 8-7 所示。在这些条件下，我们可以对图像进行阈值化处理，并在水平和垂直方向上分别对像素值求和。由于这些经阈值处理的边界值为 1，而其他部分值为 0，所以这些边界处会给出很强的响应，可以告诉我们从何处进行裁剪。

下面函数接受一幅灰度图像和一个方向，返回该方向上的 10 条边界：

from scipy.ndimage import measurements

def find_sudoku_edges(im, axis=0):
    """ 对一幅对齐后的数独图像查找单元格的边界 """
    # 阈值处理，处理后对每行（列）相加求和
    trim = 1 * (im < 128)
    s = trim.sum(axis=axis)
    
    # 寻找连通域
    s_labels, s_nbr = measurements.label(s > (0.5 * max(s)))
    
    # 计算各连通域的质心
    m = measurements.center_of_mass(s, s_labels, range(1, s_nbr + 1))
    
    # 对质心取整，质心即为相线条所在位置
    x = [int(x[0]) for x in m]
    
    # 只要检测到 4 条粗线，便在这 4 条粗线之间添加直线
    if len(x) == 4:
        dx = diff(x)
        x = [
            x[0], x[0] + dx[0] / 3, x[0] + 2 * dx[0] / 3,
            x[1], x[1] + dx[1] / 3, x[1] + 2 * dx[1] / 3,
            x[2], x[2] + dx[2] / 3, x[2] + 2 * dx[2] / 3, x[3]
        ]
        if len(x) == 10:
            return x
        else:
            raise RuntimeError('Edges not detected.')

首先对图像进行阈值化处理，对灰度值小于 128 的暗区域赋值为 1，否则为 0；然后在特定的方向上（如 axis=0 或 1）对这些经阈值处理后的像素相加求和。Scipy.ndimage 包含 measurements 模块，该模块在二进制或标记数组中对于计数及测量区域是非常有用的。首先，labels() 找出二进制数组中相连接的部件；该二进制数组是通过求和后取中值并进行阈值化处理得到的。然后，center_of_mass() 函数计算每个独立组件的质心。你可能得到 4 个或 10 个点，这主要依赖于数独平面造型设计（所有的线条是等粗细的或子网格线条比其他的粗）。在 4 个点的情况下，会以一定的间隔插入 6 条直线。如果最后的结果没有 10 条线，则会抛出一个异常。

sudokus 文件夹里包含不同难易程度的数独图像，每幅图像都对应一个包含数独真实值的文件，我们可以用它来检查识别结果。有一些图像已经和图像的边框对齐，从中挑选一幅图像，用以检查图像裁剪及分类的性能：

找到边界后，从每一个单元格提取出 crops。将裁剪出来的这些单元格传给同一特征提取函数，并将提取出来的特征作为训练数据保存在一个数组中。通过 loadtxt()读取数独图像的真实标记，用 svm_predict() 函数对这些特征向量进行分类，在控制台上打印出的结果应该为：

样例输出结果

这里使用的只是其中较简单的图像

如果用一个 9×9 的子图绘制这些经裁剪后的单元格，它们应该和图 8-7（右图）类似

图 8-7：一个检测并裁剪这些数独网格区域的例子：一幅数独网格图像 ( 左 )；9×9 裁剪后的图像，每个独立单元都会被送到 OCR 分类器中 ( 右 )

8.4.5　图像校正

如果你对上面分类器的性能还算满意，那么下一个挑战便是如何将它应用于那些没有对齐的图像。这里我们将用一种简单的图像校正方法来结束本章数独图像识别的例子，使用该校正方法的前提是网格的 4 个角点都已经被检测到或者手工做过标记。图 8-8（左）中是一幅在进行识别时受角度影响剧烈的图像。

一个单应矩阵可以像上面的例子那样映射网格以使边缘能够对齐，我们这里所要做的就是估算该变换矩阵。下面的例子手工标记 4 个角点，然后将图像变换为一个1000×1000 大小的方形图像：

from scipy import ndimage
import homography

imname = 'sudoku8.jpg'
im = array(Image.open(imname).convert('L'))

# 标记角点
figure()
imsshow(im)
gray()
x = ginput(4)

# 左上角、右上角、右下角、左下角
fp = array([array([p[1], p[0], 1]) for p in x]).T
tp = array([[0, 0, 1], [0, 1000, 1], [1000, 1000, 1], [1000, 0, 1]]).T

# 估算单应矩阵
H = homography.H_from_points(tp, fp)

# 辅助函数，用于进行几何变换
def warpfcn(x):
    x = array([x[0], x[1], 1])
    xt = dot(H, x)
    xt = xt / xt[2]
    return xt[0], xt[1]

# 用全透视变换对图像进行变换
im_g = ndimage.geometric_transform(im, warpfcn, (1000, 1000))

这里用到了 scipy.ndinmage 模块中一个更加普遍的变换函数 geometric_transform()，该函数获取一个 2D 到 2D 的映射，映射为另一个二维来取代变化矩阵，所以我们需要一个辅助函数（该例中用一个三角形的分段仿射变换），变换后的图像如图 8-8 中右图所示

图 8-8：用全透视变换对一幅图像进行校正的例子。四个角点被手工标记的数独原图（左），变换为 1000×1000 大小的方形图（右）