Previous tutorial: Removing out-of-focus blur filter
Next tutorial: Anisotropic image segmentation using gradient structure tensors
| Original author | Karpushin Vladislav |
|---|---|
| Compatibility | OpenCV >= 3.0 |
Goals
In this tutorial you will learn
- What is PSF for motion blur images
- How to restore motion blur images
Theory
For information on degenerate image model theory and Wiener filtering theory, see the tutorial “Defocus Blur Filter“. This page only considers linear motion blur distortion. The motion blurred images on this page are real-world images. Blur is caused by moving subjects.
What is PSF for motion blurred images?
The point spread function (PSF) of linear motion blur distortion is a line segment. This PSF is specified by two parameters: LEN is the length of the blur and THETA is the angle of motion.
Point spread function of linear motion blur deformation
How to restore blurred images?
This page uses the Wiener filter as the restoration filter. For details, please refer to the tutorial “Defocus blur filter“. In order to synthesize a Wiener filter in the case of motion blur, it is necessary to specify the signal-to-noise ratio (SNR) of the PSF, LEN em> andTHETA.
Source code
You can find the source code in the OpenCV source code repository at samples/cpp/tutorial_code/ImgProc/motion_deblur_filter/motion_deblur_filter.cpp.
#include <iostream>
#include "opencv2/imgproc.hpp"
#include "opencv2/imgcodecs.hpp"
using namespace cv;
using namespace std;
void help();
void calcPSF(Mat & amp; outputImg, Size filterSize, int len, double theta);
void fftshift(const Mat & amp; inputImg, Mat & amp; outputImg);
void filter2DFreq(const Mat & amp; inputImg, Mat & amp; outputImg, const Mat & amp; H);
void calcWnrFilter(const Mat & amp; input_h_PSF, Mat & amp; output_G, double nsr);
void edgetaper(const Mat & amp; inputImg, Mat & amp; 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}"
;
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()) //Check whether the image has been loaded
{<!-- -->
cout << "ERROR : Image cannot be loaded..!!" << endl;
return -1;
}
Mat imgOut;
//Only need to process even-numbered images
Rect roi = Rect(0, 0, imgIn.cols & amp; -2, imgIn.rows & amp; -2);
//Calculate Hw (start)
Mat Hw,h;
calcPSF(h, roi.size(), LEN, THETA);
calcWnrFilter(h, Hw, 1.0 / double(snr));
//Calculate Hw (stop)
imgIn.convertTo(imgIn, CV_32F);
edgetaper(imgIn, imgIn);
// Filter (start)
filter2DFreq(imgIn(roi), imgOut, Hw);
// filter (stop)
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;
}
void calcPSF(Mat & amp; 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];
}
void fftshift(const Mat & amp; inputImg, Mat & amp; 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);
}
void filter2DFreq(const Mat & amp; inputImg, Mat & amp; outputImg, const Mat & amp; 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];
}
void calcWnrFilter(const Mat & amp; input_h_PSF, Mat & amp; output_G, double nsr)
{<!-- -->
Math_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);
denom + = nsr;
divide(planes[0], denom, output_G);
}
void edgetaper(const Mat & amp; inputImg, Mat & amp; 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);
}
Description
Motion blur image restoration algorithms include PSF generation, Wiener filter generation and frequency domain blur image filtering:
// Only need to process even-numbered images Rect roi = Rect(0, 0, imgIn.cols & amp; -2, imgIn.rows & amp; -2); //Calculate Hw (start) Mat Hw,h; calcPSF(h, roi.size(), LEN, THETA); calcWnrFilter(h, Hw, 1.0 / double(snr)); //Calculate Hw (stop) imgIn.convertTo(imgIn, CV_32F); edgetaper(imgIn, imgIn); // Filter (start) filter2DFreq(imgIn(roi), imgOut, Hw); // filter (stop)
The function calcPSF() forms a PSF from the input parameters LEN and THETA (unit: degrees):
void calcPSF(Mat & amp; 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];
}
The function edgetaper() tapers the edges of the input image to reduce ringing effects in the restored image:
void edgetaper(const Mat & amp; inputImg, Mat & amp; 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);
}
The functions calcWnrFilter(), fftshift(), and filter2DFreq() implement image filtering by specifying the frequency domain PSF. These functions are copied from the Defocus Filter tutorial.
Results
The image below is a real image with motion blur distortion. The license plates of both vehicles could not be read. The red mark indicates the license plate location.
Motion blurred images. Unable to read license plate
Below is the restoration result of the black car license plate. The calculated results are LEN = 125, THETA = 0, SNR = 700.
Restored image of black car license plate
Below is the restoration result of a white car license plate. The calculated results are LEN = 78, THETA = 15, SNR = 300.
Restored image of white car license plate
Values for signal-to-noise ratio, LEN, and THETA were selected manually to obtain the best visual results. The THETA parameter is consistent with the direction of travel of the car, while the LEN parameter depends on the speed of the car. Although the effect is not perfect, it at least gives us an idea of the content of the image. After some hard work, the car license plate can now be read.
Notes
The parameters LEN and THETA are the most important. LEN and THETA should be adjusted first, then SNR< /em>.
You can also find a quick video demonstration of license plate recovery methods on YouTube.