【OCR学习笔记】3、OCR图像预处理之Python原生实现算法（中）-阿里云开发者社区

1 简介

图像噪声是指存在于图像数据中不必要的或多余的干扰信息，产生于图像的采集、量化或传输过程，对图像的后处理、分析均会产生极大的影响，因此一种好的去噪方法在去除噪声的同时，还需要保持图像的边界和细节。

2 预处理方法之平滑去噪

OCR的预处理方法可以分为以下几个大类：二值化、平滑去噪和倾斜角检测以及校正，脑图如下：

本章由于篇幅的原因，此文只讲解平滑去噪部分，往后更新会把后续的部分更新结束。

由上图可以看出，平滑去噪主要分为空间滤波、小波阈值滤波和非局部滤波方法以及基于神经网络的滤波方法。

2.1 空间滤波

空间滤波由一个邻域和对该邻域内像素执行的预定义操作组成，滤波器中心遍历输入图像的每个像素点之后就得到了处理后的图像。每经过一个像素点，邻域中心坐标的像素值就替换为预定义操作的计算结果。若在图像像素上执行的是线性操作，则该滤波器称为线性空间滤波器，反之则称为非线性空间滤波器。

2.1.1 线性空间滤波器

常见的线性空间滤波器有平滑线性滤波和高斯滤波器，平滑线性滤波器的输出是包含在滤波器模板邻域内像素的简单平均值，因此也称为均值滤波器。

高斯滤波器
假设构造宽为W、高位H的高斯卷积算子，其中W和H为奇数，锚点位置在，步骤如下：

1、计算高斯矩阵
2、计算高斯矩阵的和
3、高斯矩阵除以其本身的和，即归一化操作，得到高斯卷积算子

高斯卷积算子对应的Python实现如下：

# -*- coding: utf-8 -*-
import numpy as np
from scipy import signal
import math
'''
    得到高斯平滑算子：
    getGaussKernel(sigma,H,W),使用过程中一般H和W一般为奇数，sigma大于零
'''
def getGaussKernel(sigma,H,W):
    '''
        第一步：构建高斯矩阵gaussMatrix
    '''
    gaussMatrix = np.zeros([H,W],np.float32)
    #得到中心点的位置
    cH = (H-1)/2
    cW = (W-1)/2
    #计算1/(2*pi*sigma^2)
    coefficient = 1.0/(2*np.pi*math.pow(sigma,2))
    #
    for r in range(H):
        for c in range(W):
            norm2 = math.pow(r-cH,2) + math.pow(c-cW,2)
            gaussMatrix[r][c] = coefficient*math.exp(-norm2/(2*math.pow(sigma,2)))
    '''
        第二步：计算高斯矩阵的和
    '''
    sumGM = np.sum(gaussMatrix)
    '''
        第三步：归一化，gaussMatrix/sumGM
    '''
    gaussKernel = gaussMatrix/sumGM
    return gaussKernel
#主函数
if __name__ =="__main__":
    gaussKernel = getGaussKernel(2,3,3)
    print(gaussKernel)
    #高斯核gaussKernel是可分解的，可以分解为水平方向和垂直方向的卷积核
    gaussKernel_x = getGaussKernel(2,1,3)
    print(gaussKernel_x)
    gaussKernel_y = getGaussKernel(2,3,1)
    print(gaussKernel_y)
    '''
    水平方向和垂直方向的卷积核进行卷积运算,注意：mode='full',boundary = 'fill',fillvalue=0
    这样的参数设置，得到的结果才和gaussKernel完全相等，否则，边界不相等
    '''
    gaussKernel_xy = signal.convolve2d(gaussKernel_x,gaussKernel_y,mode='full',boundary ='fill',fillvalue=0)
    print(gaussKernel_xy)

均值滤波器
均值滤波，顾名思义就是把图像中的每一个位置的领域的平均值作为该位置的输出值，代码实现与分离的高斯卷积是类似的，只需要将高斯算子替换成均值算子即可。

均值平滑算子对应的Python实现如下：

# -*- coding: utf-8 -*-
import sys
import numpy as np
from scipy import signal
import cv2
# 均值平滑
def meanBlur(image,H,W,_boundary='fill',_fillvalue=0):
    # H、W均不为零
    if H==0 or W==0:
        print 'W or H is not zero'
        return image
    # 没有对均值平滑算子进行分离
    #meanKernel = 1.0/(H*W)*np.ones([H,W],np.float32)
    #result = signal.convolve2d(image,meanKernel,mode='same',boundary = _boundary,fillvalue=_fillvalue)
    # 卷积后进行数据类型转换,得到均值平滑的结果
    #result = result.astype(np.uint8)
    #return result
    # 因为均值算子是可分离的卷积核，根据卷积运算的结合律
    # 可以先进行水平方向的卷积，
    # 再进行垂直方向的卷积
    # 首先水平方向的均值平滑
    meanKernel_x = 1.0/W*np.ones([1,W],np.float32)
    i_conv_mk_x = signal.convolve2d(image,meanKernel_x,mode='same',boundary = _boundary,fillvalue=_fillvalue)
    # 然后对得到的水平卷积的结果再进行垂直方向的卷积
    meanKernel_y = 1.0/H*np.ones([H,1],np.float32)
    i_conv_xy = signal.convolve2d(i_conv_mk_x,meanKernel_y,mode='same',boundary = _boundary,fillvalue=_fillvalue)
    i_conv_xy = np.round(i_conv_xy)
    # 卷积后的结果进行数据类型转换，得到均值平滑的结果
    result = i_conv_xy.astype(np.uint8)
    return result
#主函数:示例
if __name__ =="__main__":
    if len(sys.argv)>1:
        image = cv2.imread(sys.argv[1],cv2.CV_LOAD_IMAGE_GRAYSCALE)
    else:
        print("Usge:python meanBlur.py imageFile")
    # 均值滤波卷积核的宽高均设为 2*halfWinSize+1
    halfWinSize = 1
    MAX_HALFWINSIZE = 20
    cv2.namedWindow("meanBlur",1)
    # 回调函数，均值滤波
    def callback_meanBlur(_halfWinSize):
        result = meanBlur(image,2*_halfWinSize+1,2*_halfWinSize+1,_boundary='symm',_fillvalue=0)
        cv2.imshow("meanBlur",result)
    callback_meanBlur(halfWinSize)
    cv2.createTrackbar("winSize/2","meanBlur",halfWinSize,MAX_HALFWINSIZE,callback_meanBlur)
    latexImage = meanBlur(image,29,29,'symm')
    cv2.imwrite("latexImage.png",latexImage)
    cv2.waitKey(0)
    cv2.destroyAllWindows()

2.1.2 非线性空间滤波器

中值滤波
是一种非线性数字滤波器。与上述滤波不同，其首先是对邻域内的值进行排序，中值作为输出。该滤波器消除椒盐噪声的效果比平滑线性空间滤波器更好。
中值滤波对应的Python实现如下：

# -*- coding: utf-8 -*-
import sys
import numpy as np
import cv2
# 中值滤波 
def medianBlur(image,winSize):
    # 图像的宽高
    rows,cols = image.shape
    # 窗口的宽高，均为奇数
    winH,winW = winSize
    halfWinH = (winH-1)/2
    halfWinW = (winW-1)/2
    # 中值滤波后的输出图像
    medianBlurImage = np.zeros(image.shape,image.dtype)
    for r in range(rows):
        for c in range(cols):
            # 判断边界
            rTop = 0 if r-halfWinH < 0 else r-halfWinH
            rBottom = rows-1 if r+halfWinH > rows-1 else r+halfWinH
            cLeft = 0 if c-halfWinW < 0 else c-halfWinW
            cRight = cols-1 if c+halfWinW > cols-1 else c+halfWinW
            # 取中值的区域
            region = image[rTop:rBottom+1,cLeft:cRight+1]
            # 求中值
            medianBlurImage[r][c] = np.median(region)
    return medianBlurImage
# 主函数
if __name__ =="__main__":
    if len(sys.argv)>1:
        image = cv2.imread(sys.argv[1],cv2.CV_LOAD_IMAGE_GRAYSCALE)
    else:
        print("Usge:python medianBlur.py imageFile")
    # 显示原图
    cv2.imshow("image",image)
    # 中值滤波
    medianBlurImage = medianBlur(image,(3,3))
    # 显示中值滤波后的结果
    cv2.imshow("medianBlurImage",medianBlurImage)
    cv2.imwrite("medianBlurImage.jpg",medianBlurImage)
    cv2.waitKey(0)
    cv2.destroyAllWindows()

双边滤波(Bilateral Filter)
不仅如高斯滤波一样考虑到像素间的欧氏距离，还关注到像素范围域中的辐射差异（如像素与中心间的相似程度、颜色强度和深度距离等），这两个权重域分别称为：空间域（Spatial Domain S）和像素范围域（Range Domain R）。权重的计算公式如下：

其中，和是平滑参数，和分别是像素和的强度函数。计算权重之后，对结果做归一化，公式如下：

其中，是去噪后的像素强度。

双边滤波对应的Python实现如下：

# -*- coding: utf-8 -*-
import sys
import numpy as np
import cv2
import math
# 基于空间距离的权重模板 ( 和计算高斯算子的过程是一样的 )
def getClosenessWeight(sigma_g,H,W):
    r,c = np.mgrid[0:H:1, 0:W:1]
    r -= (H-1)/2
    c -= (W-1)/2
    closeWeight = np.exp(-0.5*(np.power(r,2)+np.power(c,2))/math.pow(sigma_g,2))
    return closeWeight
# BilateralFiltering 双边滤波，返回的数据类型为浮点型
def bfltGray(I,H,W,sigma_g,sigma_d):
    # 构建空间距离的权重模板
    closenessWeight = getClosenessWeight(sigma_g,H,W)
    # 模板的中心点位置
    cH = (H -1)/2
    cW = (W -1)/2
    # 图像矩阵的行数和列数
    rows,cols = I.shape
    #双边滤波后的结果
    bfltGrayImage = np.zeros(I.shape,np.float32)
    for r in range(rows):
        for c in range(cols):
            pixel = I[r][c]
            # 判断边界
            rTop = 0 if r-cH < 0 else r-cH
            rBottom = rows-1 if r+cH > rows-1 else r+cH
            cLeft = 0 if c-cW < 0 else c-cW
            cRight = cols-1 if c+cW > cols-1 else c+cW
            # 权重模板作用的区域
            region = I[rTop:rBottom+1,cLeft:cRight+1]
            # 构建灰度值相似性的权重因子
            similarityWeightTemp = np.exp(-0.5*np.power(region -pixel,2.0)/math.pow(sigma_d,2))        
            closenessWeightTemp = closenessWeight[rTop-r+cH:rBottom-r+cH+1,cLeft-c+cW:cRight-c+cW+1]
            # 两个权重模板相乘
            weightTemp = similarityWeightTemp*closenessWeightTemp
            weightTemp = weightTemp/np.sum(weightTemp)
            bfltGrayImage[r][c] = np.sum(region*weightTemp)
    return bfltGrayImage
# 主函数
if __name__ =="__main__":
    if len(sys.argv)>1:
        image = cv2.imread(sys.argv[1],cv2.CV_LOAD_IMAGE_GRAYSCALE)
    else:
        print("Usge:python BFilter.py imageFile")
    # 显示原图
    cv2.imshow("image",image)
    cv2.imwrite("image.png",image)
    # 将灰度值归一化
    image = image/255.0
    # 双边滤波
    bfltImage = bfltGray(image,33,33,10,0.8)
    # 显示双边滤波的结果
    cv2.imshow("BilateralFiltering",bfltImage)
    bfltImage = bfltImage*255.0
    bfltImage = np.round(bfltImage)
    bfltImage = bfltImage.astype(np.uint8)
    # 保存双边滤波的结果
    cv2.imwrite("BilateralFiltering.png",bfltImage)
    cv2.waitKey(0)
    cv2.destroyAllWindows()

2.2 小波阈值去噪

对于二维图像信号，可分别在水平和垂直的方向使用滤波器，从而实现二维小波多分辨率分解。如下图为经过小波三层多分辨率分解之后子图像的划分。

其中，LL子带是在两个方向利用低通小波滤波器产生的小波系数，表示图像的近似；HL子带是在行方向利用低通小波滤波器后，再用高通小波滤波器在列方向卷积产生的系数，表示图像水平方向的奇异特性；LH子带是在行方向利用高通小波滤波器后，再用低通小波滤波器在列方向卷积产生的系数，表示图像垂直方向的奇异特性；HH子带是在两个方向利用高通小波滤波器产生的小波系数，表示图像的对角边缘特性。通常情况下，噪声部分包含在HL、LH和HH中，只要对这三部分作相应的小波系数处理，最后进行重构即可达到消噪的目的。

小波去噪的基本思路：

1）二维图像的小波分解：选择一个小波和小波分解的层次N，计算信号s到第N层的分解。

2）对高频系数进行阈值量化：对1～N的每一层选择阈值，并对该层的高频系数进行软阈值量化处理。

3）二维小波重构：根据小波分解第N层的低频系数和经过修改的第一到第N层的各层高频系数，计算二维信号的小波重构。

2.3 非局部方法

2.3.1 NL-Means

A.Buades等人于2005年首次将Non-Local理论应用于图像去噪，提出了基于图像所有像素的非局部均值算法(NL-means)，他们对去噪方法Dh做了如下定义：

其中，v是带噪声的图像，n是滤波参数，通常取决于噪声的标准偏差。理想情况下，比v更加平滑，为白噪声。给定离散噪声图像，对于像素i，估计值NL[v]|（i）为图像中所有像素的加权平均值：

其中，权重取决于像素i、j间的相似性，且。而所谓的相似性，则进一步取决于强度灰度级向量和的相似性。其中，表示以像素k为中心固定大小的正方形邻域，该相似度被定义为加权欧式距离的递减函数：，a是高斯核的标准差。与相近的灰阶邻域像素在平均值上有较大的权重：

其中，是归一化常数，参数h表示影响过滤程度。通过控制指数函数的衰减，可以降低权重。

2.3.2 BM3D

BM3D(Block-Matching 3D Filtering)是当前效果最好的算法之一。它与NL-means算法具有相似之处在于同样运用到了非局部块匹配的思想，但其复杂度远高于NL-Means。

BM3D共包含两大步骤，每一步里又包括相似块分组、协同滤波和聚合。找到相似块后，不同于NL-means的均值处理，BM3D会对其做域转换操作，利用协同过滤以降低相似块自身的噪声，并在聚合阶段对相似块进行加权处理。

BM3D算法具体步骤：

第1步：基础估计

1)相似块分组：首先在噪声图像中选择参照块，在参照块周围适当大小区域内进行搜索，寻找若干个差异度最小的块，与参照块一起整合成一个三维的矩阵。距离函数的公式如下：

其中，表示以像素R为中心的参照块，为搜索块，为第一步的块大小，为L2距离，Y'为硬阈值操作，为归一化后的二维线性操作。根据距离找到的相似块集合可用如下公式表示：

其中，为判断是否相似的超参数。

2)协同滤波:先使用归一化的3D线性变换降低相似块中的噪声，再使用对应的反变换得到处理后的相似块：

其中，为归一化的3D线性变换，为对应的反变换，为硬阈值操作，为处理后的相似块。这一步可以有效避免NL-means均值操作导致的相似块信息冗余和引入参照块自身包含的噪声这两大问题。

3)聚合:对上一步处理后得到的相似块进行加权平均值操作得到目标块的像素值：

其中，为权重。

第2步：最终估计

1)相似块分组:对步骤一初步处理后的图像，重新计算L2距离，得到相似块集合：

其中，为判断是否为相似块的超参数。

2)协同过滤:对第一步得到的相似块集合做域变换后，使用维纳收缩系数加权，再经过反变换得到新的集合：

其中，为维纳收缩系数。为相似块集合。

3)聚合:与步骤1的聚合操作相似，最终结果与2)一致：

BM3D算法对应的pybm3d库函数使用如下：

import numpy as np
import skimage.data
from skimage.measure import compare_psnr
import pybm3d
noise_std_dev = 40
img = skimage.data.astronaut()
noise = np.random.normal(scale=noise_std_dev, size=img.shape).astype(img.dtype)
noisy_img = img + noise
out = pybm3d.bm3d.bm3d(noisy_img, noise_std_dev)
noise_psnr = compare_psnr(img, noisy_img)
out_psnr = compare_psnr(img, out)
print("PSNR of noisy image: ", noise_psnr)
print("PSNR of reconstructed image: ", out_psnr)