[人工智能]【opencv 450 Image Processing】Motion Deblur Filter 运动去模糊滤波器

Goal

在本教程中，您将学习：

运动模糊图像(motion blur image)的 PSF 是什么

如何恢复运动模糊图像restore a motion blur image

Theory

对于退化图像模型理论(degradation image model theory)和 Wiener 滤波器理论(?Wiener filter theory?)，您可以参考教程 Out-of-focus Deblur Filter。 在此页面上，仅考虑线性运动模糊失真(?linear motion blur distortion)。 此页面上的运动模糊图像是真实世界的图像。 模糊是由移动的主体引起的。

What is the PSF of a motion blur image?

线性运动模糊失真的点扩散函数（PSF）是一条线段。 这样的 PSF 由两个参数指定：LEN 是模糊的长度，THETA 是运动的角度。

?

Point spread function of a linear motion blur distortion

线性运动模糊失真的点扩散函数

How to restore a blurred image?

如何恢复模糊的图像？

本页使用维纳滤镜作为还原滤镜，详细内容可以参考教程失焦去模糊滤镜。 为了针对运动模糊情况合成维纳滤波器，需要指定 PSF 的信噪比 (SNR)、LEN 和 THETA。

Source code

您可以在 OpenCV 源代码库的 samples/cpp/tutorial_code/ImgProc/motion_deblur_filter/motion_deblur_filter.cpp 中找到源代码。

/**
* @brief You will learn how to recover an image with motion blur distortion using a Wiener filter
* 您将学习如何使用维纳滤波器恢复具有运动模糊失真的图像
* @author Karpushin Vladislav, karpushin@ngs.ru, https://github.com/VladKarpushin
*/
#include <iostream>
#include "opencv2/imgproc.hpp"
#include "opencv2/imgcodecs.hpp"

using namespace cv;
using namespace std;

void help();
void calcPSF(Mat& outputImg, Size filterSize, int len, double theta); //计算点扩散图像
void fftshift(const Mat& inputImg, Mat& outputImg);//傅里叶之后得到的结果频率范围是0到fs，为了便于进行频率域滤波，也便于观察频谱信息 ，通常将频率范围调整至-fs/2到fs/2，这样就将零频分量（直流分量）迁移到了图像中心，呈现的效果就是中心低频信息，四周外围是高频信息，这个实现我们就称为fftshift。
void filter2DFreq(const Mat& inputImg, Mat& outputImg, const Mat& H);//频域滤波后输出时域图像
void calcWnrFilter(const Mat& input_h_PSF, Mat& output_G, double nsr);//计算维纳滤波器
void edgetaper(const Mat& inputImg, Mat& outputImg, double gamma = 5.0, double beta = 0.2);//使输入图像的边缘逐渐变细，以减少恢复图像中的振铃效应

const String keys =
"{help h usage ? |             | print this message				}"
"{image          |input.png    | input image name				}"
"{LEN            |125          | length of a motion				}"
"{THETA          |0            | angle of a motion in degrees	}"
"{SNR            |700          | signal to noise ratio信噪比	}"
;
// 运动模糊图像恢复算法    包括 PSF 生成、Wiener 滤波器生成和在频域中过滤模糊图像：
int main(int argc, char *argv[])
{
    help();
    CommandLineParser parser(argc, argv, keys);
    if (parser.has("help"))
    {
        parser.printMessage();
        return 0;
    }

    int LEN = parser.get<int>("LEN");
    double THETA = parser.get<double>("THETA");
    int snr = parser.get<int>("SNR");
    string strInFileName = parser.get<String>("image");

    if (!parser.check())
    {
        parser.printErrors();
        return 0;
    }

    Mat imgIn;
    imgIn = imread(strInFileName, IMREAD_GRAYSCALE);
    if (imgIn.empty()) //检查图像是否加载 
    {
        cout << "ERROR : Image cannot be loaded..!!" << endl;
        return -1;
    }

    Mat imgOut;

//! [main]
    // it needs to process even image only 图像的大小为 2*2,4*4 ，8*8 ………………都是可以进行去模糊的
    Rect roi = Rect(0, 0, imgIn.cols & -2, imgIn.rows & -2);//

    //Hw calculation (start)
    Mat Hw, h;
    calcPSF(h, roi.size(), LEN, THETA);
    calcWnrFilter(h, Hw, 1.0 / double(snr));
    //Hw calculation (stop)

    imgIn.convertTo(imgIn, CV_32F);
    edgetaper(imgIn, imgIn);//边缘逐渐变细

    // filtering (start)
    filter2DFreq(imgIn(roi), imgOut, Hw);//频域进行滤波
    // filtering (stop)
//! [main]

    imgOut.convertTo(imgOut, CV_8U);
    normalize(imgOut, imgOut, 0, 255, NORM_MINMAX);
    imwrite("result.jpg", imgOut);
    return 0;
}

void help()
{
    cout << "2018-08-14" << endl;
    cout << "Motion_deblur_v2" << endl;
    cout << "You will learn how to recover an image with motion blur distortion using a Wiener filter" << endl;
}

//! [calcPSF] 根据输入参数 LEN 和 THETA（以度为单位）形成一个 PSF     线性运动模糊失真的点扩散函数（PSF）是一条线段
void calcPSF(Mat& outputImg, Size filterSize, int len, double theta)
{
    Mat h(filterSize, CV_32F, Scalar(0));
    Point point(filterSize.width / 2, filterSize.height / 2);//中心点  椭圆圆心
    ellipse(h, point, Size(0, cvRound(float(len) / 2.0)), 90.0 - theta, 0, 360, Scalar(255), FILLED);//椭圆
    Scalar summa = sum(h);
    outputImg = h / summa[0];//归一化
}
//! [calcPSF]

//! [fftshift] fft变换后进行频谱搬移 https://blog.csdn.net/zhaitianbao/article/details/117955380
void fftshift(const Mat& inputImg, Mat& outputImg)
{
    outputImg = inputImg.clone();
    int cx = outputImg.cols / 2;
    int cy = outputImg.rows / 2;
    Mat q0(outputImg, Rect(0, 0, cx, cy));
    Mat q1(outputImg, Rect(cx, 0, cx, cy));
    Mat q2(outputImg, Rect(0, cy, cx, cy));
    Mat q3(outputImg, Rect(cx, cy, cx, cy));
    Mat tmp;
    q0.copyTo(tmp);
    q3.copyTo(q0);
    tmp.copyTo(q3);
    q1.copyTo(tmp);
    q2.copyTo(q1);
    tmp.copyTo(q2);
}
//! [fftshift]

//! [filter2DFreq] 在频域中过滤模糊图像，再转换到时域
void filter2DFreq(const Mat& inputImg, Mat& outputImg, const Mat& H)
{
    Mat planes[2] = { Mat_<float>(inputImg.clone()), Mat::zeros(inputImg.size(), CV_32F) };
    Mat complexI;
    merge(planes, 2, complexI);
    dft(complexI, complexI, DFT_SCALE);

    Mat planesH[2] = { Mat_<float>(H.clone()), Mat::zeros(H.size(), CV_32F) };
    Mat complexH;
    merge(planesH, 2, complexH);
    Mat complexIH;
    mulSpectrums(complexI, complexH, complexIH, 0);//频域相乘进行滤波

    idft(complexIH, complexIH);//变换到时域
    split(complexIH, planes);
    outputImg = planes[0];//时域输出滤波后图像
}
//! [filter2DFreq]

//! [calcWnrFilter] 计算维纳滤波器
void calcWnrFilter(const Mat& input_h_PSF, Mat& output_G, double nsr)
{
    Mat h_PSF_shifted;
    fftshift(input_h_PSF, h_PSF_shifted);//点扩散函数 先进行搬移
    Mat planes[2] = { Mat_<float>(h_PSF_shifted.clone()), Mat::zeros(h_PSF_shifted.size(), CV_32F) };
    Mat complexI;
    merge(planes, 2, complexI);
    dft(complexI, complexI);//傅里叶变换
    split(complexI, planes);
    Mat denom;
    pow(abs(planes[0]), 2, denom);//|H|^2
    denom += nsr;//|H|^2+1/snr
    divide(planes[0], denom, output_G);//H/(|H|^2+1/snr)  ->output_G
}
//! [calcWnrFilter]

//! [edgetaper]使输入图像的边缘逐渐变细，以减少恢复图像中的振铃效应
//通过edgetaper在反卷积之前调用该函数，可以减少噪声放大和沿图像边界的振铃效应。 图像恢复对噪声功率参数的敏感性降低
void edgetaper(const Mat& inputImg, Mat& outputImg, double gamma, double beta)//对图像边缘进行模糊处理
{
    int Nx = inputImg.cols;
    int Ny = inputImg.rows;
    Mat w1(1, Nx, CV_32F, Scalar(0));
    Mat w2(Ny, 1, CV_32F, Scalar(0));

    float* p1 = w1.ptr<float>(0);
    float* p2 = w2.ptr<float>(0);
    float dx = float(2.0 * CV_PI / Nx);
    float x = float(-CV_PI);
    for (int i = 0; i < Nx; i++)
    {
        p1[i] = float(0.5 * (tanh((x + gamma / 2) / beta) - tanh((x - gamma / 2) / beta)));
        x += dx;
    }
    float dy = float(2.0 * CV_PI / Ny);
    float y = float(-CV_PI);
    for (int i = 0; i < Ny; i++)
    {
        p2[i] = float(0.5 * (tanh((y + gamma / 2) / beta) - tanh((y - gamma / 2) / beta)));
        y += dy;
    }
    Mat w = w2 * w1;
    multiply(inputImg, w, outputImg);
}
//! [edgetaper]

Explanation

运动模糊图像恢复算法包括 PSF 生成、Wiener 滤波器生成和在频域中过滤模糊图像：

    // it needs to process even image only
    Rect roi = Rect(0, 0, imgIn.cols & -2, imgIn.rows & -2);
    //Hw calculation (start)
    Mat Hw, h;
    calcPSF(h, roi.size(), LEN, THETA);
    calcWnrFilter(h, Hw, 1.0 / double(snr));
    //Hw calculation (stop)
    imgIn.convertTo(imgIn, CV_32F);
    edgetaper(imgIn, imgIn);
    // filtering (start)
    filter2DFreq(imgIn(roi), imgOut, Hw);
    // filtering (stop)

函数 calcPSF() 根据输入参数 LEN 和 THETA（以度为单位）形成一个 PSF：

void calcPSF(Mat& outputImg, Size filterSize, int len, double theta)
{
    Mat h(filterSize, CV_32F, Scalar(0));
    Point point(filterSize.width / 2, filterSize.height / 2);
    ellipse(h, point, Size(0, cvRound(float(len) / 2.0)), 90.0 - theta, 0, 360, Scalar(255), FILLED);
    Scalar summa = sum(h);
    outputImg = h / summa[0];
}

函数 edgetaper() 使输入图像的边缘逐渐变细，以减少恢复图像中的振铃效应：振铃效应（ringing artifacts） - ostartech - 博客园 (cnblogs.com)

void edgetaper(const Mat& inputImg, Mat& outputImg, double gamma, double beta)
{
    int Nx = inputImg.cols;
    int Ny = inputImg.rows;
    Mat w1(1, Nx, CV_32F, Scalar(0));
    Mat w2(Ny, 1, CV_32F, Scalar(0));
    float* p1 = w1.ptr<float>(0);
    float* p2 = w2.ptr<float>(0);
    float dx = float(2.0 * CV_PI / Nx);
    float x = float(-CV_PI);
    for (int i = 0; i < Nx; i++)
    {
        p1[i] = float(0.5 * (tanh((x + gamma / 2) / beta) - tanh((x - gamma / 2) / beta)));
        x += dx;
    }
    float dy = float(2.0 * CV_PI / Ny);
    float y = float(-CV_PI);
    for (int i = 0; i < Ny; i++)
    {
        p2[i] = float(0.5 * (tanh((y + gamma / 2) / beta) - tanh((y - gamma / 2) / beta)));
        y += dy;
    }
    Mat w = w2 * w1;
    multiply(inputImg, w, outputImg);
}

函数 calcWnrFilter()、fftshift() 和 filter2DFreq() 实现了在频域中通过指定的 PSF 对图像进行过滤。 这些函数是从教程 Out-of-focus Deblur Filter 中复制而来的。

Result

下面您可以看到具有运动模糊失真的真实世界图像。 两辆车的车牌都是不可读的。 红色标记显示汽车的车牌位置。

Motion blur image. The license plates are not readable

运动模糊图像。 车牌不可读

下面你可以看到黑色车牌的恢复结果。 结果是使用 LEN = 125、THETA = 0、SNR = 700 计算得出的。

?The restored image of the black car license plate

黑色车牌的还原图

?

下面你可以看到白色车牌的恢复结果。 结果是使用 LEN = 78、THETA = 15、SNR = 300 计算得出的。

SNR、LEN 和 THETA 的值是手动选择的，以提供最佳的视觉效果。 THETA 参数与小车的运动方向一致，LEN 参数取决于小车的运动速度。 结果并不完美，但至少它给了我们图像内容的暗示。 经过一番努力，现在可以读取汽车牌照。

笔记

参数 LEN 和 THETA 是最重要的。 您应该先调整 LEN 和 THETA，然后再调整 SNR。

您还可以在 YouTube https://youtu.be/xSrE0hdhb4o 上找到车牌恢复方法的快速视频演示。

?