图像平滑opencv c++实现

图像噪声

噪声图像是指在图像中存在随机干扰或不规则变化的图像。以下是几种常见的噪声图像类型:

  1. 高斯噪声图像:高斯噪声是指在图像中每个像素值上添加的服从高斯分布的随机噪声。它是最常见的图像噪声类型之一。

    //对原始图像添加高斯噪声
    Mat addGaussianNoise(Mat image, double mean, double stddev)
    {    //image 原始图像 mean  stddev
    	// 创建随机数生成器
    	default_random_engine generator;
    	normal_distribution<double> distribution(mean, stddev);
    
    	// 添加高斯噪声
    	Mat noise_image = image.clone();
    	for (int i = 0; i < noise_image.rows; i++)
    	{
    		for (int j = 0; j < noise_image.cols; j++)
    		{
    			double noise = distribution(generator);
    			int pixel = noise_image.at<uchar>(i, j) + noise;
    			if (pixel < 0)
    			{
    				pixel = 0;
    			}
    			else if (pixel > 255)
    			{
    				pixel = 255;
    			}
    			noise_image.at<uchar>(i, j) = pixel;
    		}
    	}
    
    	return noise_image;
    }
    
  2. 椒盐噪声图像:椒盐噪声是指在图像中随机出现亮白或暗黑的像素点,使得这些像素点的值变得与周围像素不一致。它模拟了图像传感器或信号传输中的随机错误。

    //对原始图像添加椒盐噪声
    void addSaltAndPepperNoise(cv::Mat &src, cv::Mat&dst, double noise_ratio, int border_size) {
    	//src原始图像,dst结果图像,noise_ratio噪声比例,border_size边界距离
    	dst = src.clone();
    	int rows = src.rows - 2 * border_size;
    	int cols = src.cols - 2 * border_size;
    	int num_noise_pixels = static_cast<int>(noise_ratio * rows * cols);
    
    	for (int i = 0; i < num_noise_pixels; ++i) {
    		//随机点位
    		int row = rand() % rows + border_size;
    		int col = rand() % cols + border_size;
    
    		if (rand() % 2 == 0) {
    			dst.at<uchar>(row, col) = 0; // 椒噪声
    		}
    		else {
    			dst.at<uchar>(row, col) = 255; // 盐噪声
    		}
    	}
    }
    
  3. 瑞利噪声图像:瑞利噪声是指在图像中存在强度不均匀的干扰,常见于无线通信或医学图像中。它的幅度遵循瑞利分布。

  4. 指数噪声图像:指数噪声是指在图像中存在幅度随机变化的干扰,它的幅度遵循指数分布。

  5. 拉普拉斯噪声图像:拉普拉斯噪声是指在图像中存在有尖峰和宽峰的干扰,它的幅度遵循拉普拉斯分布。

  6. 泊松噪声图像:泊松噪声是指在图像中存在随机出现的光子计数的波动干扰,它的幅度遵循泊松分布。泊松噪声在低光条件下的图像中常见。

这些是常见的噪声类型,它们可以单独存在,也可以以组合形式出现在图像中。此外,还有其他一些特定于特定应用领域的噪声类型,例如相机传感器噪声、压缩噪声等。

图像平滑

图像平滑是指在图像中减少噪声或过滤掉不必要的细节,使图像更加平滑和清晰的过程。以下是几种常用的图像平滑方法:

  1. 均值滤波(Mean Filter):将每个像素的值替换为其邻域像素的平均值。这种方法简单有效,但可能会导致图像细节的模糊。

  2. 中值滤波(Median Filter):将每个像素的值替换为其邻域像素值的中值。中值滤波能够有效地去除椒盐噪声等离群点,对于保留边缘细节效果较好。

  3. 高斯滤波(Gaussian Filter):通过将每个像素与其邻域像素进行加权平均,使用高斯函数作为权重。高斯滤波器对于平滑图像并保持边缘信息非常有效。

  4. 双边滤波(Bilateral Filter):在计算加权平均值时考虑像素的空间距离和像素值之间的相似性。这种滤波方法可以在保持边缘细节的同时进行噪声抑制。

  5. 邻近平均滤波(Neighborhood Averaging Filter):将每个像素的值替换为其邻域像素的加权平均值,其中较远的像素具有较小的权重。该方法可用于平滑图像并减少噪声。

这些方法可以单独使用,也可以组合使用,根据具体的图像和应用需求选择适合的方法。此外,还有其他一些图像平滑方法,如双向滤波、总变差滤波等,可以根据具体情况进行选择和尝试。

中值滤波器

a)原始中值滤波器

01 原理

中值滤波是基于排序统计理论的一种能有效抑制噪声的非线性信号处理技术,中值滤波的基本原理是把数字图像或数字序列中一点的值用该点的一个邻域中各点值的中值代替,让周围的像素值接近的真实值,从而消除孤立的噪声点。方法是用某种结构的二维滑动模板,将板内像素按照像素值的大小进行排序,生成单调上升(或下降)的为二维数据序列。

02代码实现

//中值滤波器
void MedianFlitering(Mat  src, Mat & dst, int ksize)
{
	int R = ksize / 2;//窗口半径

	dst = src.clone();
	//对图像边界进行处理
	cv::Mat  paddedSrc;
	copyMakeBorder(src, paddedSrc, R, R, R, R, BORDER_REFLECT);//对图像边界添加宽度为R个像素的镜像边框

	for (int y = R; y < dst.rows - R; y++)
	{
		for (int x = R; x < dst.cols - R; x++)
		{
			vector<uchar> values;
			for (int j = -R; j <= R; j++) {
				for (int i = -R; i <= R; i++) {
					values.push_back(paddedSrc.at<uchar>(y + j, x + i));
				}
			}
			sort(values.begin(), values.end());//对像素值进行排序
			dst.at<uchar>(y, x) = values[values.size() / 2];
		}
	}
}

b)直方图加速中值滤波器

01原理

直方图加速中值滤波器的基本原理是利用直方图统计和查找表的思想,通过预先计算滤波窗口内像素值的直方图,并根据直方图信息快速找到中值。

02代码实现

//用直方图对中值滤波器加速
void FastMedianFilter(const Mat& src, Mat& dst, int ksize) {


	int R = ksize / 2;
	int window_size = ksize * ksize;
	int median_position = window_size / 2;//中间值位置

	dst = src.clone();
	//对图像边界进行处理
	Mat paddedSrc;
	copyMakeBorder(src, paddedSrc, R, R, R, R, BORDER_REPLICATE);//BORDER_REPLICATE复制边界填充

	for (int y = R; y < dst.rows - R; y++) {
		for (int x = R; x < dst.cols - R; x++) {
			int histogram[256] = { 0 };//设置直方图
			for (int j = -R; j <= R; j++) {
				for (int i = -R; i <= R; i++) {
					histogram[src.at<uchar>(y + j, x + i)]++;//统计在该灰度值点的个数
				}
			}
			int count = 0;
			for (int i = 0; i < 256; i++) {
				count += histogram[i];//对窗口内点进行计数
				if (count > median_position) {
					//对边缘噪点进行处理
					dst.at<uchar>(y - R, x - R) = i;//忽略填充边界
					break;
				}
			}
		}
	}
}

c)实验

(1)不同窗口下三种实现方式速度对比(对椒盐噪声)

3*3(second)5*5(second)7*7(second)
opencv中值滤波器0.00026470.00136740.0060189
中值滤波0.06468330.1244980.179576
直方图加速中值滤波器0.02005180.01925340.0214629
  • 窗口越小速度越快

中值滤波器实现中,复杂度主要来自两部分:

  1. 对于每个像素,遍历一个大小为 ksize*ksize 的窗口并收集值。如果 ksizeK,那么这部分的复杂度是 O(K^2)

  2. 对于收集的每个像素值,对它们进行排序。在一般情况下,排序算法的时间复杂度是 O(n log n)。在这个例子中, n 是窗口大小,也就是 K^2

对于每个像素,总的复杂度是 O(K^2 log K^2)。又因为对图像中的每个像素都执行这个操作,所以,如果图像的大小是 M*NM 行,N 列),那么整个中值滤波器的时间复杂度就是 O(M*N*K^2 log K^2)

直方图方法,可以将中值滤波器的复杂度降低到 O(M*N*K^2)

(2)讨论分析中值滤波在什么样的噪声图像下效果比较好

对比噪声类型:椒盐噪声,高斯噪声

在这里插入图片描述

处理椒盐噪声的效果较好

双边滤波器

双边高斯滤波器

01原理

双边滤波的核函数是空间域核与像素范围域核的综合结果:在图像的平坦区域,像素值变化很小,对应的像素范围域权重接近于1,此时空间域权重起主要作用,相当于进行高斯模糊;在图像的边缘区域,像素值变化很大,像素范围域权重变大,从而保持了边缘的信息。

img

02代码实现

//双边高斯滤波器
void bilateralFilterCustom(const Mat& src, Mat& dst, int d, double sigmaColor, double sigmaSpace) {
	int radius = d / 2;  // 计算滤波器半径
	double colorCoeff = -0.5 / (sigmaColor * sigmaColor);  // 颜色权重系数
	double spaceCoeff = -0.5 / (sigmaSpace * sigmaSpace);  // 空间权重系数

	Mat srcPadded;  // 定义扩展后的源图像
	copyMakeBorder(src, srcPadded, radius, radius, radius, radius, BORDER_REFLECT);  // 对源图像进行扩展

	dst = Mat::zeros(src.size(), src.type());  // 初始化输出图像

	// 遍历扩展后的源图像
	for (int y = radius; y < srcPadded.rows - radius; y++) {
		for (int x = radius; x < srcPadded.cols - radius; x++) {
			double sum = 0.0;  // 初始化滤波器加权和
			double wSum = 0.0;  // 初始化权重和
			// 遍历滤波器窗口
			for (int j = -radius; j <= radius; j++) {
				for (int i = -radius; i <= radius; i++) {
					double spaceDist2 = i * i + j * j;  // 计算空间距离的平方
					double colorDist2 = pow(srcPadded.at<uchar>(y + j, x + i) - srcPadded.at<uchar>(y, x), 2);  // 计算颜色差的平方
					double w = exp(spaceDist2 * spaceCoeff) * exp(colorDist2 * colorCoeff);  // 计算权重
					sum += srcPadded.at<uchar>(y + j, x + i) * w;  // 累加加权和
					wSum += w;  // 累加权重和
				}
			}
			dst.at<uchar>(y - radius, x - radius) = round(sum / wSum);  // 计算输出像素值
		}
	}
}

直方图加速的双边滤波

01 原理

02具体实现步骤

根据给定的源图像,使用直方图加速的双边滤波器的自定义实现步骤如下:

  1. 计算空间权重:创建一个大小为d x d的空间权重矩阵space_weight,遍历每个像素位置(i, j),计算空间距离((i - d / 2)^2 + (j - d / 2)^2)的负指数值,并将其赋给对应位置的空间权重值。

  2. 创建颜色权重查找表:创建一个大小为256 x 1的颜色权重矩阵color_weight,遍历每个灰度值i,计算其平方的负指数值,并将其赋给对应位置的颜色权重值。

  3. 对源图像进行边界扩展:创建一个扩展后的源图像padded_src,使用copyMakeBorder函数将原图像src的边界进行镜像填充,填充宽度为d/2

  4. 初始化目标图像:创建一个与源图像大小相同的目标图像dst,并将其初始化为全零图像。

  5. 应用直方图加速双边滤波器:对于每个像素位置(y, x)在源图像中,遍历滤波器窗口。对于窗口内的每个位置(dy, dx),计算当前像素与中心像素的灰度差color_diff,获取对应的空间权重和颜色权重。累加每个窗口位置的加权像素值到sum中,并累加权重到wsum中。

  6. 计算输出像素值:通过除以权重和wsum,将sum除以wsum得到输出像素值,并使用static_cast将其转换为uchar类型。将计算得到的像素值赋给目标图像dst的相应位置。

最终,函数返回输出图像dst

03代码实现

//用直方图对双边滤波器进行加速
void fastBilateralFilter(const Mat& src, Mat& dst, int d, double sigma_color, double sigma_space) {
	//计算空间权重
	Mat space_weight(d, d, CV_32F);
	for (int i = 0; i < d; i++) {
		for (int j = 0; j < d; j++) {
			float g = static_cast<float>(exp(-((i - d / 2) * (i - d / 2) + (j - d / 2) * (j - d / 2)) / (2 * sigma_space * sigma_space)));
			space_weight.at<float>(i, j) = g;
		}
	}

	// 创建颜色权重查找表
	const int L = 256;
	Mat color_weight(L, 1, CV_32F);//生成颜色直方图,颜色差距最大为255
	for (int i = 0; i < L; i++) {
		float g = static_cast<float>(exp(-(i * i) / (2 * sigma_color * sigma_color)));
		color_weight.at<float>(i, 0) = g;
	}

	// 对源图像进行边界扩展
	Mat padded_src;
	copyMakeBorder(src, padded_src, d / 2, d / 2, d / 2, d / 2, BORDER_REFLECT);

	// 初始化目标图像
	dst = Mat::zeros(src.size(), src.type());

	// 应用直方图加速双边滤波器
	for (int y = 0; y < src.rows; y++) {
		for (int x = 0; x < src.cols; x++) {
			float sum = 0;
			float wsum = 0;
			for (int dy = -d / 2; dy <= d / 2; dy++) {
				for (int dx = -d / 2; dx <= d / 2; dx++) {
					int color_diff = abs(padded_src.at<uchar>(y + dy + d / 2, x + dx + d / 2) - padded_src.at<uchar>(y + d / 2, x + d / 2));//计算颜色差距color_diff,一定会落在颜色直方图内 
					float weight = space_weight.at<float>(dy + d / 2, dx + d / 2) * color_weight.at<float>(color_diff, 0);
					sum += weight * padded_src.at<uchar>(y + dy + d / 2, x + dx + d / 2);
					wsum += weight;
				}
			}
			dst.at<uchar>(y, x) = static_cast<uchar>(sum / wsum);
		}
	}
}

c)问题

  • 双边高斯滤波和传统高斯滤波对比
    在这里插入图片描述

  • 三种实现方式速度对比

opencv双边高斯滤波器0.014423
双边高斯滤波器1.51627
直方图加速双边高斯滤波器0.0336481

在这里插入图片描述

在不同噪声下进行比较

噪声类型:椒盐噪声、高斯噪声

在这里插入图片描述

维纳滤波器

01原理

维纳滤波器,也被称为最小均方误差滤波器,是一种线性滤波器,用于恢复受到噪声影响的图像。其基本的工作原理是通过最小化输出误差的平方均值来估计原始图像。

在图像拍摄过程中由于各种原因会造成图像退化,图像退化模型如下
g ( x , y ) = h ( x , y ) ⋆ f ( x , y ) + η ( x , y ) g(x,y)=h(x,y)⋆f(x,y)+η(x,y) g(x,y)=h(x,y)f(x,y)+η(x,y)
傅里叶变换后
G ( u , v ) = H ( u , v ) F ( u , v ) + N ( u , v ) G(u,v)=H(u,v)F(u,v)+N(u,v) G(u,v)=H(u,v)F(u,v)+N(u,v)
如果退化函数和加性噪声都考虑,空域滤波器无法解决图像退化问题,逆滤波效果因为噪声的存在会变得非常差,这个时候就需要用到维纳滤波,(维纳滤波的推导写在结论中)维纳滤波公式如下:
F ( u , v ) = [ H ( u , v ) 1 ∣ H ( u , v ) ∣ 2 + S η ( u , v ) / S f ( u , v ) ∣ H ( u , v ) ∣ 2 ] G ( u , v ) F^(u,v)=[ H(u,v)1∣H(u,v)∣ 2 +S η(u,v)/S f (u,v)∣H(u,v)∣ 2 ]G(u,v) F(u,v)=[H(u,v)1H(u,v)2+Sη(u,v)/Sf(u,v)H(u,v)2]G(u,v)
其中$ Sη(u,v)=∣N(u,v)∣2$ 为噪声的功率谱,这个我们可以通过用户输入的方差构造一个噪声图像 N ( u , v ) N(u,v) N(u,v)并计算功率谱;

S f ( u , v ) = ∣ F ( u , v ) ∣ 2 Sf(u,v)=∣F(u,v)∣2 Sf(u,v)=∣F(u,v)2为输入图像的功率谱,

退化函数进行了简化,将退化函数置为1,因此维纳滤波公式简化为:
F ( u , v ) = [ S f ( u , v ) S f ( u , v ) + S η ( u , v ) ] G ( u , v ) F ^ ( u , v ) = [ S f ( u , v ) S f ( u , v ) + S η ( u , v ) ] G ( u , v ) F ( u , v ) = [ S f ( u , v ) + S η ( u , v ) S f ( u , v ) ] G ( u , v ) F ^ ( u , v ) = [ S f ( u , v ) S f ( u , v ) + S η ( u , v ) ] G ( u , v ) \hat{F}(u, v)=\left[\frac{S_{f}(u, v)}{S_{f}(u, v)+S_{\eta}(u, v)}\right] G(u, v) F ^ (u,v)=[ S f (u,v)+S η (u,v) S f (u,v) ]G(u,v) F(u,v)=[Sf(u,v)Sf(u,v)+Sη(u,v)]G(u,v)F^(u,v)=[Sf(u,v)+Sη(u,v)Sf(u,v)]G(u,v)F(u,v)=[Sf(u,v)+Sη(u,v)Sf(u,v)]G(u,v)

维纳滤波器主要作用

  1. 噪声抑制: 维纳滤波器可以用于抑制图像中的噪声,例如高斯噪声、白噪声等。由于维纳滤波器可以根据噪声的功率谱和信号的功率谱动态调整滤波参数,因此能够更好地保留图像的细节信息。
  2. 图像去模糊: 当图像受到模糊影响时,例如运动模糊、聚焦模糊等,维纳滤波器可以用于恢复清晰的图像。维纳滤波器的去模糊效果在知道模糊函数(也称为点扩散函数)的情况下尤为有效。
  3. 图像复原: 在一些应用中,例如医学成像、遥感成像等,我们可能需要从受到噪声和模糊影响的图像中恢复出原始的、清晰的图像。这时,维纳滤波器是一种常用的图像复原工具。

02具体步骤

  1. 创建几个Mat对象,包括pad(经过零填充的原始图像),cpx(双通道频域图像),dst(滤波后的图像)。

  2. 获取最佳的傅里叶变换尺寸,使其为2的指数。

  3. 使用copyMakeBorder函数对原始图像进行零填充,获得填充后的图像pad

  4. 创建一个临时图像tmpR,将参考图像ref调整为与pad相同的尺寸。

  5. 调用GetSpectrum函数,传入临时图像tmpR,获取参考图像的频谱refSpectrum

  6. 创建一个临时图像tmpN,生成服从标准正态分布的噪声图像。

  7. 调用GetSpectrum函数,传入临时图像tmpN,获取噪声图像的频谱noiseSpectrum

  8. 创建一个双通道Mat数组planes,其中第一个通道为pad的浮点型副本,第二个通道初始化为与pad相同尺寸的全零图像。

  9. 使用merge函数将两个通道合并为一个复数频域图像cpx

  10. 使用dft函数对复数频域图像cpx进行离散傅里叶变换。

  11. 使用split函数将复数频域图像cpx分离为实部和虚部图像,并存储在planes数组中。

  12. 计算维纳滤波因子,

  13. 将参考图像的频谱refSpectrum除以参考图像频谱与噪声图像频谱之和,得到factor

  14. 使用multiply函数将实部和虚部图像与维纳滤波因子factor相乘,得到滤波后的频域图像。

  15. 使用merge函数将实部和虚部图像重新合并为一个复数频域图像cpx

  16. 使用idft函数对复数频域图像cpx进行反离散傅里叶变换,并使用DFT_SCALEDFT_REAL_OUTPUT标志进行缩放和提取实部。

  17. 使用Rect函数裁剪滤波后的图像dst,以去除填充的部分,使其恢复到原始大小。

  18. 将滤波后的图像dst转换为8位无符号整数类型,并将其作为函数的返回值。

03代码实现

/维纳滤波
Mat GetSpectrum(const Mat &src)//求功率谱
{
	Mat dst, cpx;//dst功率谱,cpx
	Mat planes[] = { Mat_<float>(src), Mat::zeros(src.size(), CV_32F) };//分离为两个通道
	merge(planes, 2, cpx);//合并通道
	dft(cpx, cpx);
	split(cpx, planes);//分离通道
	magnitude(planes[0], planes[1], dst);//计算频域幅度值
	//频谱就是频域幅度图的平方
	multiply(dst, dst, dst);
	return dst;
}
Mat WienerFilter(const Mat &src, const Mat &ref, int stddev)
{

	Mat pad, cpx, dst;//pad是原图像0填充后的图像,cpx是双通道频域图,dst是滤波后的图像

	//获取傅里叶变化最佳图片尺寸,为2的指数
	int m = getOptimalDFTSize(src.rows);
	int n = getOptimalDFTSize(src.cols);

	//对原始图片用0进行填充获得最佳尺寸图片
	copyMakeBorder(src, pad, 0, m - src.rows, 0, n - src.cols, BORDER_CONSTANT, Scalar::all(0));//BORDER_CONSTANT添加恒定的像素值

	//获得参考图片频谱
	Mat tmpR(pad.rows, pad.cols, CV_8U);
	resize(ref, tmpR, tmpR.size());
	Mat refSpectrum = GetSpectrum(tmpR);

	//获得噪声频谱
	Mat tmpN(pad.rows, pad.cols, CV_32F);
	randn(tmpN, Scalar::all(0), Scalar::all(stddev));//生成正态分布随机的矩阵,高斯噪声图像
	Mat noiseSpectrum = GetSpectrum(tmpN);

	//对src进行傅里叶变换
	Mat planes[] = { Mat_<float>(pad), Mat::zeros(pad.size(), CV_32F) };
	merge(planes, 2, cpx);
	dft(cpx, cpx);
	split(cpx, planes);

	//维纳滤波因子,这里退化函数设置为1,简化运算
	Mat factor = refSpectrum / (refSpectrum + noiseSpectrum);
	multiply(planes[0], factor, planes[0]);//实部计算
	multiply(planes[1], factor, planes[1]);//虚部计算

	//重新合并实部planes[0]和虚部planes[1]
	merge(planes, 2, cpx);

	//进行反傅里叶变换
	idft(cpx, dst, DFT_SCALE | DFT_REAL_OUTPUT);
	// 裁剪到原始大小
	dst = dst(Rect(0, 0, src.cols, src.rows));
	dst.convertTo(dst, CV_8UC1);
	return dst;
}

问题

(1)对模糊图像进行复原
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-PDpPASXT-1689073627300)(C:\Users\张翔宇\AppData\Roaming\Typora\typora-user-images\image-20230603213233215.png)]

  • 1
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值