opencv中人脸识别算法的基本原理（二）

opencv中人脸识别算法的基本原理（二）

上一篇博客介绍了opencv自带的人脸识别方法中的Eigenfaces和Fisherfaces，本文主要介绍最后一种LBPH方法的原理和过程。

LBP算子介绍

在介绍LBPH之前先要了解LBP算子的基本原理。LBP是Local Binary Pattern的简称，即局部二值模式。它是一种有效的纹理描述算子，具有旋转不变性和灰度不变性等显著特点。

基本的LBP算子

如图所示，最早的LBP算子定义为在3 * 3的窗口，以窗口中心像素点为阈值，将相邻8个像素的灰度值与中心像素比较，若周围像素值大于中心值则该像素位置标记为1，否则为0。这样每个8邻域就会产生一个8位的数，再按照其位置赋予不同权重求和一个整数就得到了该窗口的LBP值。

圆形LBP算子（扩展LBP算子）

opencv内使用的圆形LBP算子。基本的LBP算子最大缺陷在于它只覆盖了一个固定半径范围内的小区域，这显然不能满足不同尺寸和频率纹理的要求。为了适应不同尺度纹理特征，Ojala等对LBP算子进行了改进，将3 * 3邻域扩展到了任意领域，并用圆形代替了正方形，改进后的LBP算子允许半径为R的圆形邻域内有任意多像素点。

常见的圆形LBP算子如图所示，由于使用圆形来作为窗口，不能保证所有像素点都落在圆内，因此需要使用插值算法对没有完全落在像素位置点的点进行计算灰度值，opencv使用的就是双线性插值算法。以中心点C(Xc，Yc)为圆心，R为半径作为LBP算子的窗口，C周围等距离分布着P个点，则邻域内某点A(Xa,Ya)的坐标可以表示为

Xa = Xc + R * cos(2 * pi * i)/P;
Ya = Yc - R * sin(2 * pi * i)/P;
i = 0,1,2...P;

LBP等价模式

一个LBP算子可以产生多种二进制模式随着邻域内采样点增加，二进制模式种类急剧增加，如3 * 3邻域内有8个采样点，即有2^8=256中模式；5 * 5邻域内有20个采样点，就有2^20=1048576种二进制模式。这么多的二进制模式对纹理的提取或者识别都提出了挑战。为了解决这个问题，等价模式（Uniform Pattern）被提出用来对LBP算子的模式种类进行降维。等价模式的定义为：当某个局部二进制模式所对应的循环二进制数从0到1或从1到0最多只有2次跳变。不符合这个定义的都归为混合模式。00000000、00000111、10001111都是等价模式。

等价模式大大减少了二进制模式种类同时不会丢失信息。模式的数量可以由原来的2^n减少为n*(n-1) + 2种，n表示了邻域内的采样点数目。等价模式代表了图像的边缘等关键模式，如下图所示，白点表示1，黑点表示0.等价模式类占总模式的绝大多数。

旋转不变的LBP算子

上面介绍的LBP算子是灰度不变的，但是是旋转改变的。Maenpaa等人LBP算子进行了拓展，将圆形邻域进行旋转得到的LBP值，并取其最小值作为LBP值，使得新的LBP算子具有旋转不变性。旋转不变的LBP值求解过程如下图所示，最终得到的LBP值为15。

LBP算子的特点

一般分析纹理的方法包括统计分析法和结构分析法。前者是基于图像灰度值的分布和相互关系，来找出反应这些关系的特征。后者是分析纹理的结构，从中获取结构特征。任何纹理都同时包含了统计特征和结构特征，因此单一的某种分析方法都不能取得满意的结果。LBP方法可以看成结合了两种方法，一方面LBP具有原始纹理和布局规则，因此具有结构特点；同时LBP又可以看出图像经过非线性滤波之后的统计。因此，LBP方法能广泛应用于各种纹理图像识别。

opencv中的LBPH

opencv中的LBPH人脸识别的大致过程是：

1. 将训练集中人脸图调整为一维矩阵，并将调整后的人脸图像全部合成一个矩阵A；
2. 选则检测窗口对矩阵A用LBP算子处理：
   
   1. 计算检测窗口（圆形）上每一点的像素值；
   2. 对没有完全落在像素位置点的点进行双线性插值；
   3. 若检测窗口的像素值大于中心点像素值进行对比，则标记为1，否则标记为0。进而得到全部LBP矩阵；
3. 根据LBP矩阵计算每个检测窗口的空间直方图，最终得到一个一维矩阵，即训练集的特征向量；
4. 当输入测试人脸图片识别时，也需要计算其LBP矩阵然后根据空间直方图获取特征向量；
5. 然后通过比较测试人脸的特征向量和训练集的特征向量，识别结果就是训练集中与测试人脸“最近”的。

opencv中的LBPHFaceRecognizer类实现了实现人脸识别。其中void train(InputArrayOfArrays _in_src, InputArray _in_labels)函数和int predict(InputArray src) const函数都会通过Mat elbp(InputArray src, int radius, int neighbors)函数使用LBP算子：

//计算扩展的LBP
static Mat elbp(InputArray src, int radius, int neighbors) 
{
    Mat dst;
    elbp(src, dst, radius, neighbors);
    return dst;
}

而接着通过void elbp(InputArray src, OutputArray dst, int radius, int neighbors)函数根据输入图像进行分类：

static void elbp(InputArray src, OutputArray dst, int radius, int neighbors)
{
    int type = src.type();
    switch (type) {
    case CV_8SC1:   elbp_<char>(src,dst, radius, neighbors); break;
    case CV_8UC1:   elbp_<unsigned char>(src, dst, radius, neighbors); break;
    case CV_16SC1:  elbp_<short>(src,dst, radius, neighbors); break;
    case CV_16UC1:  elbp_<unsigned short>(src,dst, radius, neighbors); break;
    case CV_32SC1:  elbp_<int>(src,dst, radius, neighbors); break;
    case CV_32FC1:  elbp_<float>(src,dst, radius, neighbors); break;
    case CV_64FC1:  elbp_<double>(src,dst, radius, neighbors); break;
    default:
        String error_msg = format("Using Original Local Binary Patterns for feature extraction only works on single-channel images (given %d). Please pass the image data as a grayscale image!", type);
        CV_Error(Error::StsNotImplemented, error_msg);
        break;
    }
}

elbp_函数实现了扩展的LBP：

template <typename _Tp> static
inline void elbp_(InputArray _src, OutputArray _dst, int radius, int neighbors) {
    //获取源矩阵
    Mat src = _src.getMat();
    // 为结果分配空间，注意结果比源矩阵外围少了两圈
    _dst.create(src.rows-2*radius, src.cols-2*radius, CV_32SC1);
    Mat dst = _dst.getMat();
    dst.setTo(0);
    for(int n=0; n<neighbors; n++) {
        // 采样点
        float x = static_cast<float>(radius * cos(2.0*CV_PI*n/static_cast<float>(neighbors)));
        float y = static_cast<float>(-radius * sin(2.0*CV_PI*n/static_cast<float>(neighbors)));
        // 
        int fx = static_cast<int>(floor(x));    //向下采样
        int fy = static_cast<int>(floor(y));
        int cx = static_cast<int>(ceil(x));     //向上采样
        int cy = static_cast<int>(ceil(y));
        // 小数部分
        float ty = y - fy;
        float tx = x - fx;
        // 为双线性插值设置比重，离哪个坐标近比重大
        float w1 = (1 - tx) * (1 - ty);
        float w2 =      tx  * (1 - ty);
        float w3 = (1 - tx) *      ty;
        float w4 =      tx  *      ty;
        // 进行插值
        for(int i=radius; i < src.rows-radius;i++) {
            for(int j=radius;j < src.cols-radius;j++) {
                // 计算插值
                float t = static_cast<float>(w1*src.at<_Tp>(i+fy,j+fx) + w2*src.at<_Tp>(i+fy,j+cx) + w3*src.at<_Tp>(i+cy,j+fx) + w4*src.at<_Tp>(i+cy,j+cx));
                // floating point precision, so check some machine-dependent epsilon
                dst.at<int>(i-radius,j-radius) += ((t > src.at<_Tp>(i,j)) || (std::abs(t-src.at<_Tp>(i,j)) < std::numeric_limits<float>::epsilon())) << n;
            }
        }
    }
}

通过计算好的LBP矩阵，计算其空间直方图分布，得到图像的特征向量：

static Mat spatial_histogram(InputArray _src, int numPatterns,
                             int grid_x, int grid_y, bool /*normed*/)
{
    Mat src = _src.getMat();
    // calculate LBP patch size
    int width = src.cols/grid_x;
    int height = src.rows/grid_y;
    // allocate memory for the spatial histogram
    Mat result = Mat::zeros(grid_x * grid_y, numPatterns, CV_32FC1);
    // return matrix with zeros if no data was given
    if(src.empty())
        return result.reshape(1,1);
    // initial result_row
    int resultRowIdx = 0;
    // iterate through grid
    for(int i = 0; i < grid_y; i++) {
        for(int j = 0; j < grid_x; j++) {
            Mat src_cell = Mat(src, Range(i*height,(i+1)*height), Range(j*width,(j+1)*width));
            Mat cell_hist = histc(src_cell, 0, (numPatterns-1), true);
            // copy to the result matrix
            Mat result_row = result.row(resultRowIdx);
            cell_hist.reshape(1,1).convertTo(result_row, CV_32FC1);
            // increase row count in result matrix
            resultRowIdx++;
        }
    }
    // return result as reshaped feature vector
    return result.reshape(1,1);
}