【youcans 的 OpenCV 例程200篇】228. 特征描述之 extendLBP 改进算子

特征通常是针对于图像中的某个目标而言的。针对目标所在区域的特征描述符（Region descriptors），称为区域特征描述子。

局部二值模式（LBP，Local binary patterns）是一种用来描述图像局部纹理特征的算子，它具有旋转不变性和灰度不变性的优点。 LBP 特征计算简单、效果较好，在计算机视觉领域得到了广泛的应用。

4.2.2 改进的 LBP 纹理特征描述子

基本 LBP 算子

基本的 LBP 算子定义在 3×3 的窗口内，以窗口中心像素为阈值，与相邻的 8 个像素的灰度值比较，大于阈值则标记为 1，否则标记为 0。从右上角开始顺时针旋转，排列 8 个 0/1标记值，得到一个 8 位二进制数，就是窗口中心像素点的 LBP 值。
$LBP_{P,R} (x_c,y_c) = \sum_{p=0}^{P-1} S(g_p-g_c)*2^p\\ S(g_p-g_c) =$

{\begin{cases} 1, g_{p} \geq g_{c} \\ 0, g_{p} < g_{c} \end{cases}

$\begin{cases} 1, \quad g_p \ge g_c\\ 0, \quad g_p \lt g_c \end{cases}$

L B P_{P, R} (x_{c}, y_{c}) = p = 0 \sum P - 1 S (g_{p} - g_{c}) * 2^{p} S (g_{p} - g_{c}) = {1, g_{p} \geq g_{c} 0, g_{p} < g_{c}

基本的 LBP 纹理特征描述子只覆盖了一个固定半径范围内的小区域。这种特征描述方法是随尺度变化的，当图像尺度变化时 LBP 特征编码也会发生变化，因此在大尺寸图像时就不能准确提取到所需的纹理特征，不能反映所描述的纹理信息。

圆形可变半径模式

为了满足尺度、灰度和旋转不变性的要求，Ojala 等对 LBP 算子进行了改进，将 3×3 邻域扩展到任意邻域，并用圆形邻域代替了方形邻域。改进算子允许在半径为 R 的圆形邻域内有 P 个采样点，称为扩展 LBP 算子（Extended LBP，Circular LBP）。

每个采样点的值可以通过下式计算：
$x_p = x_c + R*cos(2 \pi p /P) \\ y_p = y_c - R*sin(2 \pi p /P) \\$
其中 $x_c,y_c)$ 为邻域中心， $x_p,y_p)$ 为采样点 p。如果采样点 p 的坐标不是整数位置，则用双线性插值方法估计其灰度值：

[\begin{matrix} f (0, 0), f (0, 1) \\ f (1, 0), f (1, 1) \end{matrix}]

$\begin{bmatrix} f(0,0), f(0,1)\\ f(1,0), f(1,1) \end{bmatrix}$

[\begin{matrix} 1 - y \\ y \end{matrix}]

$\begin{bmatrix} 1-y\\ y \end{bmatrix}$

f (x, y) = [1 - x, x] [f (0, 0), f (0, 1) f (1, 0), f (1, 1)] [1 - y y]

旋转不变模式

LBP 算子是灰度不变的，但不是旋转不变的。图像旋转后的 LBP 值是不同的，从而影响识别精度。。

Maenpaa等提出具有旋转不变性的 LBP 算子，不断旋转圆形邻域得到一系列 LBP 值，将最小的 LBP 值作为该邻域的 LBP 值，从而具有旋转不变性。

等价模式（Uniform Pattern）

一个 LBP 特征大量不同的二进制组合，且随邻域集内采样点数的增加而以指数形式增长，如半径为 R 的圆形区域内含有 P 个采样点的 LBP 算子具有 $2^P$ 种不同模式。模式种类数量太多不利于纹理的提取、分类、识别及存取，因此需要研究用少量模式数量来表示图像的纹理特征。

Ojala提出了“等价模式（Uniform Pattern）”，把某个最多有两次 0/1 跳变的二进制组合所对应的 LBP 定义为一个等价模式类，其它所有的模式都定义为混合模式类。这种方法将模式数量由 $2^P$ 种减少为 $P (P - 1) + 2$ 种，显著减少了特征向量的维数，可以减少高频噪声的影响。

此外，还有各种改进的 LBP 方法，例如：

TLBP，中心像素与周围所有像素比较，而不是选择 P 个采样点
DLBP，考察四个方向的灰度变化，每个方向用 2bits 编码
MLBP，用采样点像素的平均值代替中心像素进行比较处理
MB-LBP，将图像级联分块，以小区域代替像素单位进行处理
VLBP，对于动态图像序列，考虑前 p 帧图像和后 p 帧图像的 LBP 特征
RGB-LBP，对彩色图像的 RGB 颜色分量分别计算 LBP 后再进行连接

例程 14.8：特征描述之 extendLBP 改进算子

    # 14.8 特征描述之 extendLBP 改进算子
    def basicLBP(gray):
        height, width = gray.shape
        dst = np.zeros((height, width), np.uint8)
        kernelFlatten = np.array([1, 2, 4, 128, 0, 8, 64, 32, 16])  # 从左上角开始顺时针旋转
        for h in range(1, height-1):
            for w in range(1, width-1):
                LBPFlatten = (gray[h-1:h+2, w-1:w+2] >= gray[h, w]).flatten()  # 展平为一维向量, (9,)
                dst[h, w] = np.vdot(LBPFlatten, kernelFlatten)  # 一维向量的内积
        return dst

    # extend LBP，在半径为 R 的圆形邻域内有 N 个采样点
    def extendLBP(gray, r=3, n=8):
        height, width = gray.shape
        ww = np.empty((n, 4), np.float)  # (8,4)
        p = np.empty((n, 4), np.int)  # [x1, y1, x2, y2]
        for k in range(n):  # 双线性插值估计坐标偏移量和权值
            # 计算坐标偏移量 rx，ry
            rx = r * np.cos(2.0 * np.pi * k / n)
            ry = -(r * np.sin(2.0 * np.pi * k / n))
            # 对采样点分别进行上下取整
            x1, y1 = int(np.floor(rx)), int(np.floor(ry))
            x2, y2 = int(np.ceil(rx)), int(np.ceil(ry))
            # 将坐标偏移量映射到 0-1
            tx = rx - x1
            ty = ry - y1
            # 计算插值的权重
            ww[k, 0] = (1 - tx) * (1 - ty)
            ww[k, 1] = tx * (1 - ty)
            ww[k, 2] = (1 - tx) * ty
            ww[k, 3] = tx * ty
            p[k, 0], p[k, 1], p[k, 2], p[k, 3] = x1, y1, x2, y2

        dst = np.zeros((height-2*r, width-2*r), np.uint8)
        for h in range(r, height-r):
            for w in range(r, width-r):
                center = gray[h, w]  # 中心像素点的灰度值
                for k in range(n):
                    # 双线性插值估计采样点 k 的灰度值
                    # neighbor = gray[i+y1,j+x1]*w1 + gray[i+y2,j+x1]*w2 + gray[i+y1,j+x2]*w3 + gray[i+y2,j+x2]*w4
                    x1, y1, x2, y2 = p[k,0], p[k,1], p[k,2], p[k,3]
                    gInterp = np.array([gray[h+y1,w+x1], gray[h+y2,w+x1], gray[h+y1,w+x2], gray[h+y2,w+x2]])
                    wFlatten = ww[k,:]
                    grayNeighbor = np.vdot(gInterp, wFlatten)  # 一维向量的内积
                    # 由 N 个采样点与中心像素点的灰度值比较，构造 LBP 特征编码
                    dst[h-r, w-r] |= (grayNeighbor > center) << (np.uint8)(n-k-1)
        return dst

    # 特征描述之 extendLBP 改进算子
    img = cv2.imread("../images/fabric1.png", flags=1)
    gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)  # 灰度图像

    # 1) skimage 特征检测
    from skimage.feature import local_binary_pattern
    timeBegin = cv2.getTickCount()
    lbpSKimage = local_binary_pattern(gray, 8, 1)
    timeEnd = cv2.getTickCount()
    time = (timeEnd-timeBegin)/cv2.getTickFrequency()
    print("1) skimage.feature 封装:", round(time, 4))

    timeBegin = cv2.getTickCount()
    imgLBP1 = basicLBP(gray)  # 从右上角开始顺时针旋转
    timeEnd = cv2.getTickCount()
    time = (timeEnd-timeBegin)/cv2.getTickFrequency()
    print("2) basicLBP:", round(time, 4))

    timeBegin = cv2.getTickCount()
    r1, n1 = 3, 8
    imgLBP2 = extendLBP(gray, r1, n1)
    timeEnd = cv2.getTickCount()
    time = (timeEnd-timeBegin)/cv2.getTickFrequency()
    print("3) extendLBP(r={},n={}):{}".format(r1, n1, round(time, 4)))

    timeBegin = cv2.getTickCount()
    r2, n2 = 5, 8
    imgLBP3 = extendLBP(gray, r2, n2)
    timeEnd = cv2.getTickCount()
    time = (timeEnd-timeBegin)/cv2.getTickFrequency()
    print("4) extendLBP(r={},n={}):{}".format(r1, n1, round(time, 4)))

    plt.figure(figsize=(9, 6))
    plt.subplot(231), plt.axis('off'), plt.title("origin")
    plt.imshow(cv2.cvtColor(img, cv2.COLOR_BGR2RGB))
    plt.subplot(232), plt.axis('off'), plt.title("gray")
    plt.imshow(gray, 'gray')
    plt.subplot(233), plt.axis('off'), plt.title("LBP(skimage)")
    plt.imshow(lbpSKimage, 'gray')
    plt.subplot(234), plt.axis('off'), plt.title("basic LBP")
    plt.imshow(imgLBP1, 'gray')
    plt.subplot(235), plt.title("extend LBP (r={},n={})".format(r1,n1))
    plt.imshow(imgLBP2, 'gray'), plt.axis('off')
    plt.subplot(236), plt.title("extend LBP (r={},n={})".format(r2,n2))
    plt.imshow(imgLBP3, 'gray'), plt.axis('off')
    plt.tight_layout()
    plt.show()