椒盐噪声加入及运动模糊后维纳滤波器的实现

####椒盐噪声的加入
代码：

//椒和盐噪声密度均是0.1
void salt(Mat in)
{
	int i, j;
	int n = in.cols*in.rows*0.1;
	for (int k = 0; k < n / 2; k++)
	{
		i = rand() % in.cols;
		j = rand() % in.rows;
		in.at<uchar>(j, i) = 255;
	}
}
void pepper(Mat in)
{
	int i, j;
	int n = in.cols*in.rows*0.1;
	for (int k = 0; k < n / 2; k++)
	{
		i = rand() % in.cols;
		j = rand() % in.rows;
		in.at<uchar>(j, i) = 0;
	}
}
Mat addspnoise(Mat in)
{
	salt(in);
	pepper(in);
	return in;
}

由于是椒盐一起加入的，所以逆谐波滤波出后看起来并不好看
逆谐波滤波Q为-0.5时




逆谐波滤波Q为0.5时

在这里插入图片描述

运动模糊及维纳滤波器的实现：

//a,b,T为运动模糊参数，详见数字图像处理第三版
//运动模糊和维纳放在一个函数中了，以便处理
Mat Winner1(Mat in, float a, float b, float T)
{
	in.convertTo(in, 5);
	int nHeight = in.rows;
	int nWidth = in.cols;
	int P = 2 * nHeight;
	int Q = 2 * nWidth;
	int i, j;
	float U = 0.0;
	float pi_u = 0.0;
	Mat fp = Mat::zeros(P, Q, 5);
	for (int j = 0; j < nHeight; j++)
	{
		for (int i = 0; i < nWidth; i++)
		{
			fp.at<float>(j, i) = in.at<float>(j, i)*pow(-1, j + i);

		}
	}
	Mat planes[] = { Mat_<float>(fp), Mat::zeros(fp.size(), CV_32F) };
	Mat H[] = { Mat::zeros(fp.size(), CV_32F) , Mat::zeros(fp.size(), CV_32F) };
	Mat H_A = Mat::zeros(P, Q, 5);
	Mat out=Mat::zeros(nHeight, nWidth, 5);
	Mat complex;
	merge(planes, 2, complex);
	dft(complex, complex);
	split(complex, planes);
	float P0, P1;
	for (int j = 0; j < P; j++)
	{
		for (int i = 0; i < Q; i++)
		{

			pi_u = CV_PI * ((j - nHeight)*a + (i - nWidth) * b) + 1e-20;
			H[0].at<float>(j, i) = (T / pi_u)*(sin(pi_u))*(cos(pi_u));
			H[1].at<float>(j, i) = -(T / pi_u)*(sin(pi_u))*(sin(pi_u));
			P0 = H[0].at<float>(j, i);
			P1 = H[1].at<float>(j, i);
			H_A.at<float>(j, i) = P0 * P0 + P1 * P1;
			U = H_A.at<float>(j, i) / (H_A.at<float>(j, i) + K);
			planes[0].at<float>(j, i) = U * planes[0].at<float>(j, i);
			planes[1].at<float>(j, i) = U * planes[1].at<float>(j, i);
		}
	}
	Mat complex_cp;
	merge(planes, 2, complex_cp);
	Mat iDft[] = { Mat::zeros(planes[1].size(),CV_32F),Mat::zeros(planes[1].size(),CV_32F) };
	idft(complex_cp, complex_cp);
	split(complex_cp, iDft);

	for (int j = 0; j < nHeight; j++)
	{
		for (int i = 0; i < nWidth; i++)
		{
			out.at<float>(j, i) = iDft[0].at<float>(j, i)*pow(-1, j + i);
		}
	}

	normalize(out, out, 0, 255, CV_MINMAX);
	in.convertTo(in, 0);
	out.convertTo(out, 0);
	return out;
}

输入图像：

运动模糊图像：

Winne滤波递推以及MSE计算，最终K=0.009，MSE=901.999

约束最小二乘y迭代，最终y为