图像处理中小波变换滤波

2018-08-30 10:19:31 csdnforyou 阅读数 22756

1. 滤波器介绍

滤波器作为图像处理课程的重要内容,大致可分为两类,空域滤波器和频率域滤波器。本文主要介绍常用的四种滤波器:中值滤波器、均值滤波器、高斯滤波器、双边滤波器,并基于opencv做出实现。空域的滤波器一般可以通过模板对原图像进行卷积。

注意:空域滤波器和频率域滤波器对比

1)空间域指图像本身,空域变换直接对图像中的像素进行操作。

2)图像变换是将图像从空间域变换到某变换域(如 傅立叶变换中的频率域)的数学变换,在变换域 中进行处理,然后通过反变换把处理结果返回到空间域。

3)图像在空域上具有很强的相关性,借助于正交变 换可使在空域的复杂计算转换到频域后得到简化

4)借助于频域特性的分析,将更有利于获得图像的 各种特性和进行特殊处理

 

2、理论知识:

图像的空域滤波无非两种情况,线性滤波和非线性滤波。

      滤波的意思就是对原图像的每个像素周围一定范围内的像素进行运算,运算的范围就称为掩膜。而运算就分两种了,如果运算只是对各像素灰度值进行简单处理(如乘一个权值)最后求和,就称为线性滤波;而如果对像素灰度值的运算比较复杂,而不是最后求和的简单运算,则是非线性滤波;如求一个像素周围3x3范围内最大值、最小值、中值、均值等操作都不是简单的加权,都属于非线性滤波。

常见的线性滤波有:均值滤波、高斯滤波、盒子滤波、拉普拉斯滤波等等,通常线性滤波器之间只是模版系数不同。

非线性滤波利用原始图像跟模版之间的一种逻辑关系得到结果,如最值滤波器,中值滤波器和双边滤波器等。

1、线性滤波

线性滤波器表达公式:,其中均值滤波器和高斯滤波器属于线性滤波器,首先看这两种滤波器

均值滤波器:

模板:

 

从待处理图像首元素开始用模板对原始图像进行卷积,均值滤波直观地理解就是用相邻元素灰度值的平均值代替该元素的灰度值。

高斯滤波器:

高斯滤波一般针对的是高斯噪声,能够很好的抑制图像输入时随机引入的噪声,将像素点跟邻域像素看作是一种高斯分布的关系,它的操作是将图像和一个高斯核进行卷积操作: 

 

模板:通过高斯内核函数产生的

高斯内核函数:

例如3*3的高斯内核模板:

中值滤波:同样是空间域的滤波,主题思想是取相邻像素的点,然后对相邻像素的点进行排序,取中点的灰度值作为该像素点的灰度值。

统计排序滤波器,对椒盐噪声有很好的抑制

详细请参考:数字图像处理之椒盐噪声和中值滤波

中值滤波将窗口函数里面的所有像素进行排序取得中位数来代表该窗口中心的像素值,对椒盐噪声和脉冲噪声的抑制效果特别好,同时又能保留边缘细节,用公式表示是: 

双边滤波(Bilateral filter)也是一种非线性的滤波方法,是结合图像的空间邻近度和像素值相似度的一种折衷处理,同时考虑空域信息和灰度相似性,达到保边去噪的目的,具有简单,非迭代、局部的特点,它比高斯滤波多了一个高斯方差σd,用公式表示就是: 

w(x,y)为加权系数,取决于定义域核和值域核的乘积。

注意:

1)均值模糊无法克服边缘像素信息丢失的缺陷,原因是均值滤波是基于平均权重的。

2)高斯模糊部分克服了该缺陷,但无法完全避免,因为没有考虑像素值的不同。

3)高斯双边模糊-是边缘保留额滤波方法,避免了边缘信息丢失,保留了图像轮廓不变。

3. 实验

filter results

 

结论:从滤波的结果可以看出各种滤波算法对图像的作用非常不同,有些变化非常大,有些甚至跟原图一样。在实际应用时,应根据噪声的特点、期望的图像和边缘特征等来选择合适的滤波器,这样才能发挥图像滤波的最大优点。

 

 

4. C++实现

4.1均值滤波

static void exchange(int& a, int& b)
{    
    int t = 0;
    t = a;
    a = b;
    b = t;
}
 
static void bubble_sort(int* K, int lenth)
{
    for (int i = 0; i < lenth; i++)
        for (int j = i + 1; j < lenth; j++)
        {
            if (K[i]>K[j])
                exchange(K[i], K[j]);
        }
}
///产生二维的高斯内核
static cv::Mat generate_gassian_kernel(double u, double sigma, cv::Size size)
{
    int width = size.width;
    int height = size.height;
    cv::Mat gassian_kernel(cv::Size(width, height), CV_64FC1);
    double sum = 0;
    double sum_sum = 0;
    for (int i = 0; i < width; i++)
        for (int j = 0; j < height; j++)
        {
            sum = 1.0 / 2.0 / CV_PI / sigma / sigma * exp(-1.0 * ((i - width / 2)*(i - width / 2) + (j - width / 2)*(j - width / 2)) / 2.0 / sigma / sigma);
            sum_sum += sum;
            gassian_kernel.ptr<double>(i)[j] = sum;
        }
    for (int i = 0; i < width; i++)
        for (int j = 0; j < height; j++)
        {
            gassian_kernel.ptr<double>(i)[j] /= sum_sum;
        }
    return gassian_kernel;
}
///均值滤波
void lmt_main_blur(cv::Mat& img_in, cv::Mat& img_out, int kernel_size)
{
    img_out = img_in.clone();
    cv::Mat mat1;
    cv::copyMakeBorder(img_in, mat1, kernel_size, kernel_size, kernel_size, kernel_size, cv::BORDER_REPLICATE);
 
    int cols = mat1.cols;
    int rows = mat1.rows;
    int channels = img_out.channels();
    const uchar* const pt = mat1.ptr<uchar>(0);
    uchar* pt_out = img_out.ptr<uchar>(0);
 
    for (int i = kernel_size; i < rows - kernel_size; i++)
    {
        for (int j = kernel_size; j < cols - kernel_size; j++)
        {
            if (channels == 1)
            {
                long long int sum_pixel = 0;
                for (int m = -1 * kernel_size; m < kernel_size; m++)
                    for (int n = -1 * kernel_size; n < kernel_size; n++)
                    {
                        sum_pixel += pt[(i + m)*cols + (j + n)];
                    }
                img_out.ptr<uchar>(i - kernel_size)[j - kernel_size] = (double)sum_pixel / (kernel_size*kernel_size * 4);
            }
            else if (channels == 3)
            {
                long long int sum_pixel = 0;
                long long int sum_pixel1 = 0;
                long long int sum_pixel2 = 0;
                for (int m = -1 * kernel_size; m < kernel_size; m++)
                    for (int n = -1 * kernel_size; n < kernel_size; n++)
                    {
                        sum_pixel += pt[((i + m)*cols + (j + n))*channels + 0];
                        sum_pixel1 += pt[((i + m)*cols + (j + n))*channels + 1];
                        sum_pixel2 += pt[((i + m)*cols + (j + n))*channels + 2];
                    }
                img_out.ptr<uchar>(i - kernel_size)[(j - kernel_size)*channels + 0] = (double)sum_pixel / (double)(kernel_size*kernel_size * 4);
                img_out.ptr<uchar>(i - kernel_size)[(j - kernel_size)*channels + 1] = (double)sum_pixel1 / (double)(kernel_size*kernel_size * 4);
                img_out.ptr<uchar>(i - kernel_size)[(j - kernel_size)*channels + 2] = (double)sum_pixel2 / (double)(kernel_size*kernel_size * 4);
            }
        }
    }
 
}
///中值滤波
void lmt_median_blur(cv::Mat& img_in, cv::Mat& img_out, int kernel_size)
{
    img_out = img_in.clone();
    cv::Mat mat1;
    cv::copyMakeBorder(img_in, mat1, kernel_size, kernel_size, kernel_size, kernel_size, cv::BORDER_REPLICATE);
 
    int cols = mat1.cols;
    int rows = mat1.rows;
    int channels = img_out.channels();
 
    cv::Mat mat[3];
    cv::Mat mat_out[3];
    cv::split(mat1, mat);
    cv::split(img_out, mat_out);
    for (int k = 0; k < 3; k++)
    {
        const uchar* const pt = mat[k].ptr<uchar>(0);
        uchar* pt_out = mat_out[k].ptr<uchar>(0);
        for (int i = kernel_size; i < rows - kernel_size; i++)
        {
            for (int j = kernel_size; j < cols - kernel_size; j++)
            {
                long long int sum_pixel = 0;
                int* K = new int[kernel_size*kernel_size * 4];
                int ker_num = 0;
                for (int m = -1 * kernel_size; m < kernel_size; m++)
                    for (int n = -1 * kernel_size; n < kernel_size; n++)
                    {
                        K[ker_num] = pt[(i + m)*cols + (j + n)];
                        ker_num++;
                    }
                bubble_sort(K, ker_num);
                mat_out[k].ptr<uchar>(i - kernel_size)[j - kernel_size] = K[ker_num / 2];
            }
        }
    }
    cv::merge(mat_out, 3, img_out);
}
///高斯滤波
void lmt_gaussian_blur(cv::Mat& img_src, cv::Mat& img_dst, cv::Size kernel_size)
{
    img_dst = cv::Mat(cv::Size(img_src.cols, img_src.rows), img_src.type());
    int cols = img_src.cols;
    int rows = img_src.rows;
    int channels = img_src.channels();
    cv::Mat gassian_kernel = generate_gassian_kernel(0, 1, kernel_size);
    int width = kernel_size.width / 2;
    int height = kernel_size.height / 2;
    for (int i = height; i < rows - height; i++)
    {
        for (int j = width; j < cols - width; j++)
        {
            for (int k = 0; k < channels; k++)
            {
                double sum = 0.0;
                for (int m = -height; m <= height; m++)
                {
                    for (int n = -width; n <= width; n++)
                    {
                        sum += (double)(img_src.ptr<uchar>(i + m)[(j + n)*channels + k]) * gassian_kernel.ptr<double>(height + m)[width + n];
                    }
                }
                if (sum > 255.0)
                    sum = 255;
                if (sum < 0.0)
                    sum = 0;
                img_dst.ptr<uchar>(i)[j*channels + k] = (uchar)sum;
            }
        }
    }
 
    
}
///双边滤波
void lmt_bilateral_filter(cv::Mat& img_in, cv::Mat& img_out, const int r, double sigma_d, double sigma_r)
{
    int i, j, m, n, k;
    int nx = img_in.cols, ny = img_in.rows, m_nChannels = img_in.channels();
    const int w_filter = 2 * r + 1; // 滤波器边长  
 
    double gaussian_d_coeff = -0.5 / (sigma_d * sigma_d);
    double gaussian_r_coeff = -0.5 / (sigma_r * sigma_r);
    double  **d_metrix = new double *[w_filter];
    for (int i = 0; i < w_filter; ++i)
        d_metrix[i] = new double[w_filter];
    
    double r_metrix[256];  // similarity weight  
    img_out = cv::Mat(img_in.size(),img_in.type());
    uchar* m_imgData = img_in.ptr<uchar>(0);
    uchar* m_img_outData = img_out.ptr<uchar>(0);
    // copy the original image  
    double* img_tmp = new double[m_nChannels * nx * ny];
    for (i = 0; i < ny; i++)
        for (j = 0; j < nx; j++)
            for (k = 0; k < m_nChannels; k++)
            {
                img_tmp[i * m_nChannels * nx + m_nChannels * j + k] = m_imgData[i * m_nChannels * nx + m_nChannels * j + k];
            }
 
    // compute spatial weight  
    for (i = -r; i <= r; i++)
        for (j = -r; j <= r; j++)
        {
            int x = j + r;
            int y = i + r;
 
            d_metrix[y][x] = exp((i * i + j * j) * gaussian_d_coeff);
        }
 
    // compute similarity weight  
    for (i = 0; i < 256; i++)
    {
        r_metrix[i] = exp(i * i * gaussian_r_coeff);
    }
 
    // bilateral filter  
    for (i = 0; i < ny; i++)
        for (j = 0; j < nx; j++)
        {
            for (k = 0; k < m_nChannels; k++)
            {
                double weight_sum, pixcel_sum;
                weight_sum = pixcel_sum = 0.0;
 
                for (m = -r; m <= r; m++)
                    for (n = -r; n <= r; n++)
                    {
                        if (m*m + n*n > r*r) continue;
 
                        int x_tmp = j + n;
                        int y_tmp = i + m;
 
                        x_tmp = x_tmp < 0 ? 0 : x_tmp;
                        x_tmp = x_tmp > nx - 1 ? nx - 1 : x_tmp;   // 边界处理,replicate  
                        y_tmp = y_tmp < 0 ? 0 : y_tmp;
                        y_tmp = y_tmp > ny - 1 ? ny - 1 : y_tmp;
 
                        int pixcel_dif = (int)abs(img_tmp[y_tmp * m_nChannels * nx + m_nChannels * x_tmp + k] - img_tmp[i * m_nChannels * nx + m_nChannels * j + k]);
                        double weight_tmp = d_metrix[m + r][n + r] * r_metrix[pixcel_dif];  // 复合权重  
 
                        pixcel_sum += img_tmp[y_tmp * m_nChannels * nx + m_nChannels * x_tmp + k] * weight_tmp;
                        weight_sum += weight_tmp;
                    }
 
                pixcel_sum = pixcel_sum / weight_sum;
                m_img_outData[i * m_nChannels * nx + m_nChannels * j + k] = (uchar)pixcel_sum;
 
            } // 一个通道  
 
        } // END ALL LOOP  
    for (i = 0; i < w_filter; i++)
        delete[] d_metrix[i];
    delete[] d_metrix;
}
 

5. Opencv API实现

opencv相关函数简介:

双边滤波函数:bilateralFilter(InputArray src, OutputArray dst, int d, double sigmaColor, double sigmaSpace,int borderType=BORDER_DEFAULT )

   src待滤波图像

   dst滤波后图像

   d滤波器半径

   sigmaColor滤波器值域的sigma

   sigmaSpace滤波器空间域的sigma

   borderType边缘填充方式 BORDER_REPLICATE BORDER_REFLECT BORDER_DEFAULT BORDER_REFLECT_101BORDER_TRANSPARENT BORDER_ISOLATED

 

均值滤波函数:blur(InputArray src, OutputArray dst, Size ksize, Point anchor=Point(-1,-1), intborderType=BORDER_DEFAULT );

   src待滤波图像

   dst滤波后图像

   ksize 均值滤波器的大小

Piont(-1,-1)指中心

   anchor均值滤波器的锚点也就是模板移动点

   borderType边缘填充方式 BORDER_REPLICATE BORDER_REFLECT BORDER_DEFAULT BORDER_REFLECT_101BORDER_TRANSPARENT BORDER_ISOLATED

 

高斯滤波函数:GaussianBlur(InputArray src, OutputArray dst, Size ksize, double sigmaX, double sigmaY=0,int borderType=BORDER_DEFAULT );

   src待滤波图像

   dst滤波后图像

   ksize 高斯滤波器的大小Size(x,y),x和y必须是整数且是奇数。

   sigmaX 高斯滤波器的x方向的滤波器高斯sigma

   sigmaY 高斯滤波器的y方向的滤波器高斯sigma

   borderType边缘填充方式 BORDER_REPLICATE BORDER_REFLECT BORDER_DEFAULT BORDER_REFLECT_101BORDER_TRANSPARENT BORDER_ISOLATED

 

中值滤波函数:medianBlur(InputArray src, OutputArray dst, int ksize );

    src待滤波图像

    dst滤波后图像

    ksize 中值滤波器的大小

函数演示:

void bilateral_filter_show(void)
{
    cv::Mat mat1 = cv::imread("F:\\CVlibrary\\obama.jpg", CV_LOAD_IMAGE_GRAYSCALE); //灰度图加载进来,BGR->HSV 然后取H参数
    if (mat1.empty())
        return;
    cv::imshow("原图像", mat1); 
    cv::Mat src = cv::imread("F:\\CVlibrary\\obama.jpg");
    cv::imshow("原始彩色图像", src);
    std::cout << "channel = " << mat1.channels() << std::endl;
    
    cv::Mat mat3;
    cv::bilateralFilter(src, mat3, 5, 50, 50,cv::BORDER_DEFAULT);
    cv::imshow("opencv给出的双边滤波器", mat3);
    cv::Mat mat4;
    cv::blur(src, mat4, cv::Size(3, 3));
    cv::imshow("均值滤波", mat4);
    cv::Mat mat5;
    cv::GaussianBlur(src, mat5, cv::Size(5, 5), 1,1);
    cv::imshow("高斯滤波器", mat5);
    cv::Mat mat6;
    cv::medianBlur(src, mat6, 3);
    cv::imshow("中值滤波", mat6); 
    cv::Mat mat7;
    lmt_gaussian_blur(src, mat7, cv::Size(5, 5));
    cv::imshow("my gaussian image",mat7);
 
    cv::waitKey(0);
}

高斯、中值、均值、双边滤波的效果

#include "cv.h"
#include "highgui.h"
#include <iostream>

using namespace std;
using namespace cv;

int main(int argc, char* argv[])
{
        Mat src = imread("misaka.jpg");
        Mat dst;

        //参数是按顺序写的

        //高斯滤波
        //src:输入图像
        //dst:输出图像
        //Size(5,5)模板大小,为奇数
        //x方向方差
        //Y方向方差
        GaussianBlur(src,dst,Size(5,5),0,0);
        imwrite("gauss.jpg",dst);
        
        //中值滤波
        //src:输入图像
        //dst::输出图像
        //模板宽度,为奇数
        medianBlur(src,dst,3);
        imwrite("med.jpg",dst);
        
        //均值滤波
        //src:输入图像
        //dst:输出图像
        //模板大小
        //Point(-1,-1):被平滑点位置,为负值取核中心
        blur(src,dst,Size(3,3),Point(-1,-1));
        imwrite("mean.jpg",dst);

        //双边滤波
        //src:输入图像
        //dst:输入图像
        //滤波模板半径
        //颜色空间标准差
        //坐标空间标准差
        bilateralFilter(src,dst,5,10.0,2.0);//这里滤波没什么效果,不明白
        imwrite("bil.jpg",dst);

        waitKey();

        return 0;
}
2018-06-19 01:57:50 hellocsz 阅读数 14833

均值滤波

均值滤波,是图像处理中最常用的手段,从频率域观点来看均值滤波是一种低通滤波器,高频信号将会去掉,因此可以帮助消除图像尖锐噪声,实现图像平滑,模糊等功能。理想的均值滤波是用每个像素和它周围像素计算出来的平均值替换图像中每个像素。采样Kernel数据通常是3X3的矩阵,如下表示:

从左到右从上到下计算图像中的每个像素,最终得到处理后的图像。均值滤波可以加上两个参数,即迭代次数,Kernel数据大小。一个相同的Kernel,但是多次迭代就会效果越来越好。同样,迭代次数相同,Kernel矩阵越大,均值滤波的效果就越明显。


中值滤波

中值滤波也是消除图像噪声最常见的手段之一,特别是消除椒盐噪声,中值滤波的效果要比均值滤波更好。中值滤波是跟均值滤波唯一不同是,不是用均值来替换中心每个像素,而是将周围像素和中心像素排序以后,取中值,一个3X3大小的中值滤波如下:

 


最大最小值滤波

最大最小值滤波是一种比较保守的图像处理手段,与中值滤波类似,首先要排序周围像素和中心像素值,然后将中心像素值与最小和最大像素值比较,如果比最小值小,则替换中心像素为最小值,如果中心像素比最大值大,则替换中心像素为最大值。一个Kernel矩阵为3X3的最大最小值滤波如下:

 


双边滤波

一种同时考虑了像素空间差异与强度差异的滤波器,因此具有保持图像边缘的特性。

先看看高斯滤波器


其中W是权重,i和j是像素索引,K是归一化常量。公式中可以看出,权重只和像素之间的空间距离有关系,无论图像的内容是什么,都有相同的滤波效果。

再来看看双边滤波器,它只是在原有高斯函数的基础上加了一项,如下


其中 I 是像素的强度值,所以在强度差距大的地方(边缘),权重会减小,滤波效应也就变小。总体而言,在像素强度变换不大的区域,双边滤波有类似于高斯滤波的效果,而在图像边缘等强度梯度较大的地方,可以保持梯度


引导滤波

与双边滤波最大的相似之处,就是同样具有保持边缘特性。在引导滤波的定义中,用到了局部线性模型,至于该模型,可以暂时用下图简单的理解


该模型认为,某函数上一点与其邻近部分的点成线性关系,一个复杂的函数就可以用很多局部的线性函数来表示,当需要求该函数上某一点的值时,只需计算所有包含该点的线性函数的值并做平均即可。这种模型,在表示非解析函数上,非常有用。

同理,我们可以认为图像是一个二维函数,而且没法写出解析表达式,因此我们假设该函数的输出与输入在一个二维窗口内满足线性关系,如下


其中,q是输出像素的值,I是输入图像的值,i和k是像素索引,a和b是当窗口中心位于k时该线性函数的系数。其实,输入图像不一定是待滤波的图像本身,也可以是其他图像即引导图像,这也是为何称为引导滤波的原因。对上式两边取梯度,可以得到


即当输入图像I有梯度时,输出q也有类似的梯度,现在可以解释为什么引导滤波有边缘保持特性了。

下一步是求出线性函数的系数,也就是线性回归,即希望拟合函数的输出值与真实值p之间的差距最小,也就是让下式最小


这里p只能是待滤波图像,并不像I那样可以是其他图像。同时,a之前的系数(以后都写为e)用于防止求得的a过大,也是调节滤波器滤波效果的重要参数。通过最小二乘法,我们可以得到


其中,是I在窗口w_k中的平均值,是I在窗口w_k中的方差,是窗口w_k中像素的数量,是待滤波图像p在窗口w_k中的均值。

在计算每个窗口的线性系数时,我们可以发现一个像素会被多个窗口包含,也就是说,每个像素都由多个线性函数所描述。因此,如之前所说,要具体求某一点的输出值时,只需将所有包含该点的线性函数值平均即可,如下


这里,w_k是所有包含像素i的窗口,k是其中心位置。

当把引导滤波用作边缘保持滤波器时,往往有 I = p ,如果e=0,显然a=1, b=0是E(a,b)为最小值的解,从上式可以看出,这时的滤波器没有任何作用,将输入原封不动的输出。如果e>0,在像素强度变化小的区域(或单色区域),有a近似于(或等于)0,而b近似于(或等于),即做了一个加权均值滤波;而在变化大的区域,a近似于1,b近似于0,对图像的滤波效果很弱,有助于保持边缘。而e的作用就是界定什么是变化大,什么是变化小。在窗口大小不变的情况下,随着e的增大,滤波效果越明显。

在滤波效果上,引导滤波和双边滤波差不多,在一些细节上,引导滤波较好。引导滤波最大的优势在于,可以写出时间复杂度与窗口大小无关的算法,因此在使用大窗口处理图片时,其效率更高。


[python] view plain copy
  1. def guidedfilter(I,p,r,eps):  
  2.     '''''I:引导图图; 
  3.     p:输入图(p=I); 
  4.     r :半径: 
  5.     eps:regulation  
  6.     f:为窗口半径为r的均值滤波器; 
  7.     corr:相关; 
  8.     var:方差; 
  9.     cov:协方差 
  10.     '''  
  11.     height,width = I.reshape()  
  12.     m_I = cv2.boxFilter(I,-1,(r,r))  #f_mean(I) 均值滤波blur和盒式滤波一样  
  13.     m_p = cv2.boxFilter(p,-1,(r,r))  #f_mean(p)  
  14.   
  15.   
  16.     m_II = cv2.boxFilter(I*I,-1,(r,r)) #f_mean(I.*I)  
  17.     m_Ip = cv2.boxFilter(I * p, -1, (r, r))  #f_mean(I.*p)  
  18.   
  19.     var_I = m_II-m_I*m_I   #求方差:corr_I -mean_I.*mean_I  
  20.     cov_Ip = m_Ip - m_I * m_p  #协方差: #cov_Ip-mean_I.*mean_p  
  21.   
  22.     a = cov_Ip/(var_I+eps)  #cov_Ip./(var_I+eps)  
  23.     b = m_p-a*m_I   #mean_p -a.*mean_I  
  24.     m_a = cv2.boxFilter(a,-1,(r,r))  #mean_a  
  25.     m_b = cv2.boxFilter(b,-1,(r,r))  #mean_b  
  26.     return m_a*I+m_b  
2019-05-01 15:51:58 qq_30815237 阅读数 13365

      小波指的是一种能量在时域非常集中的波,它的能量有限,都集中在某一点附近,而且积分的值为零,这说明它与傅里叶波一样是正交波。

    图像的傅里叶变换是将图像信号分解为各种不同频率的正弦波。同样,小波变换是将图像信号分解为由原始小波位移和缩放之后的一组小波。小波在图像处理里被称为图像显微镜,原因在于它的多分辨率分解能力可以将图片信息一层一层分解剥离开来。剥离的手段就是通过低通和高通滤波器

     小波变换可以和傅里叶变换结合起来理解。傅里叶变换是用一系列不同频率的正余弦函数去分解原函数,变换后得到是原函数在正余弦不同频率下的系数。小波变换使用一系列的不同尺度的小波去分解原函数,变换后得到的是原函数在不同尺度小波下的系数。不同的小波通过平移与尺度变换分解,平移是为了得到原函数的时间特性,尺度变换是为了得到原函数的频率特性。

小波变换步骤:

1.把小波w(t)和原函数f(t)的开始部分进行比较,计算系数C。系数C表示该部分函数与小波的相似程度。

2.把小波向右移k单位,得到小波w(t-k),重复1。重复该步骤直至函数f结束.

3.扩展小波w(t),得到小波w(t/2),重复步骤1,2.

4.不断扩展小波,重复1,2,3.

haar小波:

我这里使用的haar小波,缩放函数是[1 1],小波函数是[1 -1]。是最简单的小波了。

图像二维离散小波变换 :

     图像的二维离散小波分解和重构过程如下图所示,分解过程可描述为:首先对图像的每一行进行 1D-DWT,获得原始图像在水平方向上的低频分量 L 和高频分量 H,然后对变换所得数据的每一列进行 1D-DWT,获得原始图像在水平和垂直方向上的低频分量 LL、水平方向上的低频和垂直方向上的高频 LH、水平方向上的高频和垂直方向上的低频 HL 以及水平和垂直方向上的的高频分量 HH。

    重构过程可描述为:首先对变换结果的每一列进行以为离散小波逆变换,再对变换所得数据的每一行进行一维离散小波逆变换,即可获得重构图像。由上述过程可以看出,图像的小波分解是一个将信号按照低频和有向高频进行分离的过程,分解过程中还可以根据需要对得到的 LL 分量进行进一步的小波分解,直至达到要求。

对于二维图像Haar变换不再从一个方向进行滤波,而是从水平和竖直两个方向进行低通和高通滤波(水平和竖直先后不影响),用图像表述如图所示:图中a表示原图,图b表示经过一级小波变换的结果,h1 表示水平反向的细节,v1 表示竖直方向的细节,c1表示对角线方向的细节,b表示下2采样的图像。图c中表示继续进行Haar小波变换。一级Haar小波变换实际效果如图3所示

                       

matlab实例

小波去噪实现步骤:

(1)二维信号的小波分解。选择一个小波和小波分解的层次N,然后计算信号s到第N层的分解。
(2)对高频系数进行阈值量化。对于从1~N的每一层,选择一个阈值,并对这一层的高频系数进行软阈值量化处理。
(3)二维小波重构。根据小波分解的第N层的低频系数和经过修改的从第一层到第N的各层高频系数,计算二维信号的小波重构

Matlab函数介绍

1、dwt2是二维单尺度小波变换,其可以通过指定小波或者分解滤波器进行二维单尺度小波分解。DWT2的一种语法格式:[cA,cH,cV,cD]=dwt2(X,'wname');也就是说DWT2只能对某个输入矩阵X进行一次分解。

[cA1, cH1, cV1, cD1] = dwt2(I_noise, 'haar');
figure
subplot(221), imshow(cA1, []);
subplot(222), imshow(cH1, []);
subplot(223), imshow(cV1, []);
subplot(224), imshow(cD1, []);

 

可以看出,第一张图是图像的近似,相当于图像的低频部分,而其它三张图是图像的轮廓,也就是水平,垂直和对角三个方向的细节。是图像的高频部分。至此,各变量的维数如下所示。

                                       

2、wavedec2函数

     该函数用于对多尺度二维小波进行分解,其常用调用格式:[C,S] = wavedec2(X,N,'wname'):

  • X:要进行小波分解的图像; 
  • N :指定分解的层数; 
  • wname:指定用什么小波基进行分解。 
  • 输出: 
  • c:为各层分解系数; 
  • s: 各层分解系数长度,也就是大小。

用小波函数wname对信号X在尺度N上的二维分解,其中N为大于1的正整数。可以对输入矩阵X进行N次分解。C代表分解系数的组合,是一个向量:   
                                   
     C的大小为 [1,img_height×img_width];A(N)是图像第N层的近似表示,尺度最小,在金字塔中就是每层的下采样的图像,而H、V、D分别表示图像的水平高频分量,垂直高频分量,对角高频分量。正如我们在金字塔概念中所了解的,在第N-1层下采样到N层,N层的图像维度(尺度)是变小了,也就意味着在下采样过程中丢失了信息,而这些丢失的信息实质是高频信息,那么这些信息在小波分解中可以通过HVD这些高频分量来保存。 
    这里贴上小波分解之后的结果图,直观地感受一下。这里对原始图像进行三层小波分解。红框a表示的就是近似图像。
 

    需要指出的是,每一次的小波分解都是在近似图像上进行分解。S 是储存各层分解系数长度的,即第一行是A(N)的长度,第二行是H(N)|V(N)|D(N)|的长度,第三行是 H(N-1)|V(N-1)|D(N-1)的长度,倒数第二行是H(1)|V(1)|D(1)长度,最后一行是原始图像img的长度(大小)。 这里原始图像是512×512,并进行了3层的小波分解。对应的s内容如下图: 
  
                                       

S表示每一层分解结果的维数,如果进行n层小波分解,S 的大小是(n+1)*2,最后一行表示的是原始图像的size。

                              

3、wdcbm2函数

     [thr,nkeep] = wdcbm2(c,s,alpha,m) 返回与level相关的阈值thr和要保持的系数数NKEEP, 函数用于去噪或压缩。使用基于Birge-Massart策略的小波系数选择规则获得thr。通常,alpha= 1.5用于压缩,alpha= 3用于去噪。使用wdcbm2选择各层的独立阈值。

     [C,S]是要由wavedec2函数得到的进行去噪或压缩的图像的小波分解结构,level j = size(S,1)-2.

    THR是3*j的矩阵,THR(:,j)包含对于level j情况下,水平,对角线和垂直三个方向的阈值。 NKEEP是长度为j的向量,NKEEP(j)包含要保持在级别j情况下系数的数量。

j,M和ALPHA定义策略:

  1. 在j + 1级(和更粗略的级别),一切都保留。
  2. 对于从1到j的级别i,n_i最大系数保持为n_i = M /(j + 2-i)^ ALPHA。
  3. M的默认值是M = prod(S(1,:))最粗近似系数的数量。

4、wdencmp函数

    [XC,CXC,LXC,PERF0,PERFL2] = wdencmp('lvd',C,S,'wname',N,THR,SORH) 
      函数wdencmp用于一维或二维信号的消噪或压缩。wname是所用的小波函数,gbl(global的缩写)表示每层都采用同一个阈值进行处理,lvd表示每层用不同的阈值进行处理,N表示小波分解的层数,THR为阈值3*N,SORH表示选择软阈值还是硬阈值(分别取为’s’和’h’),XC是消噪或压缩后的信号,[CXC,LXC]是XC的小波分解结构,PERF0和PERFL2是恢复和压缩L^2的范数百分比, 是用百分制表明降噪或压缩所保留的能量成分。

5、wrcoef2函数

   wrcoef2函数是用来重建一幅图像的系数,其实就是根据小波分解之后的系数c来重建其对应的图像。重建好的图像的尺度与原始图像一致。即无论你要重构哪个层的系数,最终它的维度都是和原始图像的尺度一致。其调用形式如下: 
        X = wrcoef2(‘type’,c,s,’wname’,N) 

  • type :指定要进行重构的小波系数,如a –近似图像 ;h – 水平高频分量;v – 垂直高频分量;d–对角高 
  • 频分量; 
  • c: 是小波分解函数wavedec2分解的小波系数; 
  • s: 是wavedec2分解形成的尺度; 
  • wname :指定小波基; 
  • N :指定重构的小波系数所在的层。 默认重构最大层的系数,N = size(S,1)-2。N所指的层数是如何表示的?比如将图像小波分解成3层,那么N = 3是代表256×256那一层,还是64×64那一层?N=3 代表的是64×64那一层

wrcoef2 的过程就相当于 appcoef2 或者 detcoef2 (抽取系数)后再进行 upcoef2(重构)。

clear;
close all;

file = 'lena_gray_512.tif';
img  = imread(file);
img = double(img);
% 对图像进行3层的小波分解
N = 3; % 设置分解层数
[c,s] = wavedec2(img,N,'db1');

% 对各层的近似图像a进行重构
a1 = wrcoef2('a',c,s,'db1',1);
a2 = wrcoef2('a',c,s,'db1',2);
a3 = wrcoef2('a',c,s,'db1',3);

6、appcoef2 函数

appcoef2适用于2维图像,其主要是为了提取小波分解中形成的近似图像,即低频分量。 
 A = appcoef2(c,s,’wname’,N)  

  • c:小波分解的小波系数 
  • s:小波分解的对应尺度 
  • wname :指定小波基 
  • N :指定小波系数所在的层数 

7、detcoef2

函数detcoef2 用来对二维离散小波变换的高频部分系数进行提取。 其调用形式为: 
D = detcoef2(O,c,s,N)  

  • O:指定提取哪个高频分量,取值分别为:’h’ –水平高频 or ‘v’ – 垂直高频 or ‘d’ – 对角高频; 
  • c:小波系数矩阵; 
  • s:尺度矩阵;

8、wthcoef2函数

该函数用于对二维信号的小波系数阈值进行处理,常用调用格式:
    NC = wthcoef2('type',C,S,N,T,SORH):返回经过小波分解结构[C,S]进行处理后的新的小波分解向量NC,[NC,S]即构成一个新的小波分解结构。N是一个包含高频尺度的向量,T是相应的阈值,且N和T长度须相等。返回'type'(水平、垂直、对角线)方向的小波分解向量NC。参数SORH用来对阈值方式进行选择,当SORH = 's'时,为软阈值,当SORH = 'h'时,为硬阈值。

9、重构函数 waverec2

    waverec2函数是wavedec2的反函数,返回的结果X就是原始图像。其基于小波分解结构[c,s]对矩阵X进行多级小波重构,其中[c,s]是wavedec2函数的返回值。其调用格式如下: 
      X = waverec2(c,s,’wname’) )  

  • c: 系数矩阵 
  • s: 尺度矩阵 
  • wname : 指定小波基 

值得注意的是,X = waverec2(c,s,’wname’) 相当于 X = appcoef2(c,s,’wname’,0)。

如何进行小波分解:

假设{x1,x2}是一个由两个元素组成的信号,定义这两个元素的平准和细节为:

a = (x1+x2)/2 ;d=(x1-x2)/2

变换实例如下:

    

from:https://blog.csdn.net/qq_39936376/article/details/80809770

from:http://blog.sina.com.cn/s/blog_84024a4a0101fn02.html

from:https://blog.csdn.net/Chaolei3/article/details/80940459

2016-05-02 16:18:14 liaozali9022 阅读数 16382

一、前已完成任务情况 

    1、概况

设计题目:基于正交变换与自适应滤波的图像去噪算法

设计目的:设计一种基于正交变换域自适应滤波器的的图像去噪算法,在消除图像噪声的同时尽可能地保留图像固有的信息。

提取出三个关键词:正交变换、自适应滤波、图像去噪

matlab设计流程:

    与单纯运用某种自适应算法相比,基于小波分解的自适应滤波算法在收敛速度和稳定性上都有了很大的提高

 

2、小波变换的基本理论

 

 

 

       示意:尺度越大,采用越大的时间窗,尺度越小,采用越短的时间窗,即尺度与频率成反比。

在时频两域都具有表征信号局部特征的能力,其在低频部分具有较高的频率分辨

率和较低的时间分辨率,在高频部分具有较高的时间分辨率和较低的频率分辨率,

所以被誉为分析信号的数字显微镜

 

       小波变换与傅里叶变换对比的优点:短时傅立叶分析是把一个短时间的窗函数加在信号上,再对这一部分加窗的信号做傅立叶变换,当然,假定了在短时间内信号是平稳的。然后移动窗函数,对整个时间轴上的信号做短时傅立叶分析。但是,又有问题出现了,短时傅立叶分析一旦确定了窗函数,窗函数的形状便不再改变,于是短时傅立叶分析对于信号的时频分辨率也随之确定。

而我们希望稍微智能一点,即对于变化剧烈的信号,时间分辨率能高一点,频率分辨率可以降低,而对于变化平坦的信号,频率分辨率要高一些,时间分辨率可以低点。于是小波分析便产生了。

小波分析同样对信号在一短时间内做加窗(小波函数)分析,只是这个窗,即小波函数既可以在时间轴上移动,又可以伸缩,伸缩便意味着时频分辨率的改变。

       小波变换的物理表现:分解的结果是产生长度减半的两个部分,一个是经低通滤波器产生的原始信号的平滑部分,另一个则是经高通滤波器产生的原始信号细节部分。

    不同j所确定的频带是独立的,随j变化相互独立的频带覆盖了整个频率轴,分辨率j反应了频带的位置和带宽把信号分解到一系列相互独立的频带上分辨率j-1时的近似信号=f(x)在分辨率为j时的近似部分+细节部分

再对低频部分进行相似运算

 

3、小波变换在matlab中的使用

[C, S] = wavedec2(x, n, wname); % 对图像进行小波分解

小波变换的层数选择:通常分解层数过多,并且对所有的各层小波空间的系数都进行阈值处理会造成信号的信息丢失严重,去噪后的信噪比反而下降,同时导致运算量增大,使处理变慢.

分解层数过少则去噪效果不理想,信噪比提高不多,但不会出现信噪比下降的情况.

小波函数的选取:按常理选取了sym4

 

[C,S]的含义C=[A(n)|H(n)|D(n)|………H(1)|V(1)|D(1)]

h ,v ,d 分别反映水平、垂直、对角线方向

S(1,:)=size of approximation coefficients(n)

下图的解释更直观一点:

 

 

 

 

 

 

 

 

部分程序运行:原始图像→加噪图像(同时加入高斯噪声和椒盐噪声)→小波分解后的图像显示

.M=imread('detfinger1.png'); 

%读取MATLAB中的名为detfinger的图像  分辨率:512×512

subplot(1,2,1) 

imshow(M);           %显示原始图像 

title('原始图像') 

 

P1=imnoise(M,'gaussian',0.02);     

%加入高斯躁声,方差为0.02的高斯噪声

 

P2=imnoise(P1,'salt &pepper',0.02);

%同时加入高斯噪声和椒盐躁声 

subplot(1,2,2) 

imshow(P2);%加入椒盐躁声后显示图像 

title('加入高斯椒盐躁声后');

%%对加入高斯噪声的图像进行小波分解

P2=double(P2);

[CA,CH,CV,CD]= dwt2(P2,'db1 ');

y =[CA,CH;CV,CD];

y=uint8(y);

subplot(3,1,3); 

imshow(y); 

title('小波分解');

4、滤波器的特性

我的思路:先了解维纳滤波器的原理再进一步改进,理解自适应滤波器,进行自适应滤波器的设计,再进行自适应滤波器的组合设计,从而实现去噪功能

 

维纳滤波器:根据平稳随机信号的全部过去和当前的观察数据来估计信号的当前值,在最小均方误差的条件下得到系统的传递函数,参数是固定的,适用于平稳随机信号。

卡尔曼滤波器:根据当前时刻数据的观测值和前一时刻对该一时刻的预测值进行递推数据。它自动调节本身的冲击响应特性(自动调节数字滤波器的系数),以适应信号变化的特性,从而达到最优化滤波。它的参数是时变的,适用于非平稳随机信号。

这两种滤波器最优滤波的条件:噪声的统计特性先验已知

但实际应用中,常常无法得到这些统计特征的先验知识(统计特性是随时间变化的),因此实现不了最优滤波。

自适应滤波器:无法得到一些统计特性的先验知识(统计特性随时间变化)

而自适应滤波器在输入信号统计特性未知或者统计特性变化时,能够自动调节自身的参数,使其按照某种准则达到最优滤波。由于自适应滤波器具有这种特性,自提出以来,在实际工程中的众多领域得到了广泛应用。

 

5、如何设计和建立自适应滤波器

   自适应滤波器工作原理是:系统能够按照某种算法自动调节权系数,使其实际输出和期望输出的均方误差达到最小值。

   自适应滤波器的结构:有FIRIIR两种。由于IIR滤波器存在稳定性的问题,因此一般采用FIR滤波器。由于FIR滤波器横向结构的算法具有容易实现和计算量少等优点,在对收敛速度不是很快的场合,多采用FIR作为自适应滤波器结构。

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

自适应滤波算法:

自适应滤波器的核心部分就是自适应滤波算法,基于不同的准则,最常用的自适应滤波算法是LMS算法和RLS算法,LMS算法以最小均方误差为准则,RLS算法(递归最小二乘算法)以最小误差平方和为准则。

LMS算法由于其具有计算量小,稳定性好,易于实现等优点,被广泛使用。

LMS算法中均方误差表示为:

使均方误差达到最小值时,得到最佳权系数w*,应满足下列方程:在实际中,LMS 迭代算法是以最速下降法为原则进行的,用公式表示为:

式中,是自适应步长  ( n) 为迭代的梯度。

这里来估计均方误差 MSE 的梯度,并以此梯度估计值代替最速下降法中理论情况下的梯度真值,LMS( 最小均方算法进行梯度估计的方法是以误差信号的每一次迭代的瞬间时平方值代替其均方值,并以此来估计梯度。

在自适应滤波器的设计过程中,滤波器阶数 N和步长μ 等参数的选取对仿真结果的影响总结如下.

(1)自适应算法中步长 μ 的选择:μ 的大小影响算法的收敛速度。  其中,为矩阵的最大特征值μ 取值小,收敛速度慢,计算工作量大,但滤波性能较好。μ 取值大,收敛速度快,计算工作量小,滤波性能较差。但 μ 大到一定值时,收敛速度变化不明显,且取值过大,会造成计算溢出.

(2)自适应滤波器阶数 的影响: N 不可以任意选取,需要根据经验加上实际仿真验证比较才能最终确定.

当阶数 取值大时,迭代次数增加,收敛速度变快。但当阶数 大到一定程度,收敛速度变化不明显,且可能引起系数迭代过程不收敛.

 

 

自适应滤波器原理:

输入信号x(n)

、参考信号d(n)、误差信号e(n)

e(n)y(n)d(n)的误差信号,根据最小均方误差算法调节自适应滤波器算法的参数,来优化滤波器结构,从而实现尽可能地保留图像固有的信息。

LMS算法最核心的思想是用平方误差代替均方误差,基本的LMS算法为:

w(k,n+ 1)=w(k,n)+ 2μe(n)x(n-k)

其中,w(k,n)w(k,n+ 1)分别为迭代前后的系数值;nn+ 1为前后两个时刻;k= 0,1,

… ,N- 1,N为滤波器的阶数;μ为收敛因子;e(n)=d(n)-y(n)=d(n)-xT(n)w(n)=d(n)-wT(n)x(n)为误差信号;x(n-k)为输入信号;y(n)=xT(n)w(n)为输出信号

6、matlab仿真

 

进行小波变换→  进行滤波后→  再进行图像重构→ 再小波反变换进行图像重构→先用维纳滤波器再用自适应滤波器替换

 

 

 

%%对加入高斯噪声的图像进行小波分解

J=double(J);

[CA,CH,CV,CD]= dwt2(J,'db1 ');

% [c,s]=wavedec2(j,1,'haar');

y =[CA,CH;CV,CD];

y=uint8(y);

subplot(2,2,3); 

imshow(y);

title('小波分解');

 

%xx=wiener2(CA,[3 3]); %对加噪图像进行二维自适应维纳滤波   %滤波器窗口大小

yy=wiener2(CH,[3 3]); %对加噪图像进行二维自适应维纳滤波   %滤波器窗口大小

zz=wiener2(CV,[3 3]); %对加噪图像进行二维自适应维纳滤波   %滤波器窗口大小

qq=wiener2(CD,[3 3]); %对加噪图像进行二维自适应维纳滤波   %滤波器窗口大小

 

%%%图像重构    

XX = idwt2(CA,yy,zz,qq,'db1');    

 

XX=uint8(XX);

subplot(2,2,4); 

imshow(XX); 

title('重构以后的图像');

 

                                                     

 

 

 

 

 

 

 

6自适应滤波器的组合设计:

    但是当信号源中混有多种干扰噪声时,单独一个滤波器不能达到抵消多种干扰噪声的要求。针对信号源中混有的多种干扰噪声问题,提出了自适应滤波器的组合设计问题。

分为线性自适应滤波器和非线性自适应滤波器,非线性自适应滤波器具有更强的信号处理功能,但是计算比较复杂,实际用的更多的是线性自适应滤波器。

 

总结:Mallat算法和LMS算法结合。

虽然算法简单,运算量小,易于实现,但由于它的收敛速度对输入信号的自相关函数矩阵特征值的分布敏感如果分布太散,即最大值与最小值差异太大,收敛速度就会很慢,因此利用小波的时一频局部特性,就可以减小自适应滤波器输入向量自相关阵特征值的分散程度,大大增加了算法的收敛步长,提高了算法的收敛速度和稳定性。

 

(不选)      基于小波分解重构的自适应滤波算法比仅进行小波分解的自适应滤波算法滤波精度要稍高一些,但由于需要信号重构,算法的复杂度随之增加,因此该算法实时性相对稍差。

2018-09-02 23:27:36 Eastmount 阅读数 23615

该系列文章是讲解Python OpenCV图像处理知识,前期主要讲解图像入门、OpenCV基础用法,中期讲解图像处理的各种算法,包括图像锐化算子、图像增强技术、图像分割等,后期结合深度学习研究图像识别、图像分类应用。希望文章对您有所帮助,如果有不足之处,还请海涵~

该系列在github所有源代码:https://github.com/eastmountyxz/ImageProcessing-Python
PS:请求帮忙点个Star,哈哈,第一次使用Github,以后会分享更多代码,一起加油。

同时推荐作者的C++图像系列知识:
[数字图像处理] 一.MFC详解显示BMP格式图片
[数字图像处理] 二.MFC单文档分割窗口显示图片
[数字图像处理] 三.MFC实现图像灰度、采样和量化功能详解
[数字图像处理] 四.MFC对话框绘制灰度直方图
[数字图像处理] 五.MFC图像点运算之灰度线性变化、灰度非线性变化、阈值化和均衡化处理详解
[数字图像处理] 六.MFC空间几何变换之图像平移、镜像、旋转、缩放详解
[数字图像处理] 七.MFC图像增强之图像普通平滑、高斯平滑、Laplacian、Sobel、Prewitt锐化详解

前文参考:
[Python图像处理] 一.图像处理基础知识及OpenCV入门函数
[Python图像处理] 二.OpenCV+Numpy库读取与修改像素
[Python图像处理] 三.获取图像属性、兴趣ROI区域及通道处理

本篇文章主要讲解Python调用OpenCV实现图像平滑,包括四个算法:均值滤波、方框滤波、高斯滤波和中值滤波。全文均是基础知识,希望对您有所帮助。知识点如下:
1.图像平滑
2.均值滤波
3.方框滤波
4.高斯滤波
5.中值滤波

PS:本文介绍图像平滑,想让大家先看看图像处理的效果,后面还会补充一些基础知识供大家学习。文章参考自己的博客及网易云课堂李大洋老师的讲解,强烈推荐大家学习。

PSS:2019年1~2月作者参加了CSDN2018年博客评选,希望您能投出宝贵的一票。我是59号,Eastmount,杨秀璋。投票地址:https://bss.csdn.net/m/topic/blog_star2018/index

五年来写了314篇博客,12个专栏,是真的热爱分享,热爱CSDN这个平台,也想帮助更多的人,专栏包括Python、数据挖掘、网络爬虫、图像处理、C#、Android等。现在也当了两年老师,更是觉得有义务教好每一个学生,让贵州学子好好写点代码,学点技术,"师者,传到授业解惑也",提前祝大家新年快乐。2019我们携手共进,为爱而生。

一.图像平滑

1.图像增强
图像增强是对图像进行处理,使其比原始图像更适合于特定的应用,它需要与实际应用相结合。对于图像的某些特征如边缘、轮廓、对比度等,图像增强是进行强调或锐化,以便于显示、观察或进一步分析与处理。图像增强的方法是因应用不同而不同的,研究内容包括:(参考课件和左飞的《数字图像处理》)

2.图像平滑
图像平滑是一种区域增强的算法,平滑算法有邻域平均法、中指滤波、边界保持类滤波等。在图像产生、传输和复制过程中,常常会因为多方面原因而被噪声干扰或出现数据丢失,降低了图像的质量(某一像素,如果它与周围像素点相比有明显的不同,则该点被噪声所感染)。这就需要对图像进行一定的增强处理以减小这些缺陷带来的影响。
简单平滑-邻域平均法

3.邻域平均法
图像简单平滑是指通过邻域简单平均对图像进行平滑处理的方法,用这种方法在一定程度上消除原始图像中的噪声、降低原始图像对比度的作用。它利用卷积运算对图像邻域的像素灰度进行平均,从而达到减小图像中噪声影响、降低图像对比度的目的。
但邻域平均值主要缺点是在降低噪声的同时使图像变得模糊,特别在边缘和细节处,而且邻域越大,在去噪能力增强的同时模糊程度越严重。

首先给出为图像增加噪声的代码。

# -*- coding:utf-8 -*-
import cv2
import numpy as np

#读取图片
img = cv2.imread("test.jpg", cv2.IMREAD_UNCHANGED)
rows, cols, chn = img.shape

#加噪声
for i in range(5000):    
    x = np.random.randint(0, rows) 
    y = np.random.randint(0, cols)    
    img[x,y,:] = 255

cv2.imshow("noise", img)
           
#等待显示
cv2.waitKey(0)
cv2.destroyAllWindows()

输出结果如下所示:



二.均值滤波

1.原理
均值滤波是指任意一点的像素值,都是周围N*M个像素值的均值。例如下图中,红色点的像素值为蓝色背景区域像素值之和除25。

其中红色区域的像素值均值滤波处理过程为: ((197+25+106+156+159)+ (149+40+107+5+71)+ (163+198+**226**+223+156)+ (222+37+68+193+157)+ (42+72+250+41+75)) / 25

其中5*5的矩阵称为核,针对原始图像内的像素点,采用核进行处理,得到结果图像。

提取1/25可以将核转换为如下形式:

2.代码
Python调用OpenCV实现均值滤波的核心函数如下:
result = cv2.blur(原始图像,核大小)
其中,核大小是以(宽度,高度)表示的元祖形式。常见的形式包括:核大小(3,3)和(5,5)。
K=19[111111111](1) K=\frac{1}{9}\left[ \begin{matrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{matrix} \right] \tag{1}
K=125[1111111111111111111111111](2) K=\frac{1}{25}\left[ \begin{matrix} 1 & 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 1 \\ \end{matrix} \right] \tag{2}

代码如下所示:

#encoding:utf-8
import cv2  
import numpy as np  
import matplotlib.pyplot as plt
 
#读取图片
img = cv2.imread('test01.png')
source = cv2.cvtColor(img,cv2.COLOR_BGR2RGB)
 
#均值滤波
result = cv2.blur(source, (5,5))
 
#显示图形
titles = ['Source Image', 'Blur Image']  
images = [source, result]  
for i in xrange(2):  
   plt.subplot(1,2,i+1), plt.imshow(images[i], 'gray')  
   plt.title(titles[i])  
   plt.xticks([]),plt.yticks([])  
plt.show()  

输出结果如下图所示:

核设置为(10,10)和(20,20)会让图像变得更加模糊。

如果设置为(1,1)处理结果就是原图,核中每个权重值相同,称为均值。



三.方框滤波

方框滤波和均值滤波核基本一致,区别是需不需要均一化处理。OpenCV调用boxFilter()函数实现方框滤波。函数如下:
result = cv2.boxFilter(原始图像, 目标图像深度, 核大小, normalize属性)
其中,目标图像深度是int类型,通常用“-1”表示与原始图像一致;核大小主要包括(3,3)和(5,5),如下所示。

K=19[111111111](3) K=\frac{1}{9}\left[ \begin{matrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{matrix} \right] \tag{3}
K=125[1111111111111111111111111](4) K=\frac{1}{25}\left[ \begin{matrix} 1 & 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 1 \\ \end{matrix} \right] \tag{4}

normalize属性表示是否对目标图像进行归一化处理。当normalize为true时需要执行均值化处理,当normalize为false时,不进行均值化处理,实际上为求周围各像素的和,很容易发生溢出,溢出时均为白色,对应像素值为255。

在图像简单平滑中,算法利用卷积模板逐一处理图像中每个像素,这一过程可以形象地比作对原始图像的像素一一进行过滤整理,在图像处理中把邻域像素逐一处理的算法过程称为滤波器。平滑线性滤波器的工作原理是利用模板对邻域内像素灰度进行加权平均,也称为均值滤波器。

代码如下所示:

#encoding:utf-8
import cv2  
import numpy as np  
import matplotlib.pyplot as plt
 
#读取图片
img = cv2.imread('test01.png')
source = cv2.cvtColor(img,cv2.COLOR_BGR2RGB)
 
#方框滤波
result = cv2.boxFilter(source, -1, (5,5), normalize=1)
 
#显示图形
titles = ['Source Image', 'BoxFilter Image']  
images = [source, result]  
for i in xrange(2):  
   plt.subplot(1,2,i+1), plt.imshow(images[i], 'gray')  
   plt.title(titles[i])  
   plt.xticks([]),plt.yticks([])  
plt.show()  

代码中使用5*5的核,normalize=1表示进行归一化处理,此时与均值滤波相同,输出结果如下图所示:

下面是图像左上角处理前后的像素结果:

print(source[0:3, 0:3, 0])
#[[115 180 106]
# [ 83 152  72]
# [ 55  58  55]]
print(result[0:3, 0:3, 0])
#[[92 90 78]
# [92 89 77]
# [82 80 72]]

如果省略参数normalize,则默认是进行归一化处理。如果normalize=0则不进行归一化处理,像素值为周围像素之和,图像更多为白色。

#encoding:utf-8
import cv2  
import numpy as np  
import matplotlib.pyplot as plt
 
#读取图片
img = cv2.imread('test01.png')
source = cv2.cvtColor(img,cv2.COLOR_BGR2RGB)
 
#方框滤波
result = cv2.boxFilter(source, -1, (5,5), normalize=0)
 
#显示图形
titles = ['Source Image', 'BoxFilter Image']  
images = [source, result]  
for i in xrange(2):  
   plt.subplot(1,2,i+1), plt.imshow(images[i], 'gray')  
   plt.title(titles[i])  
   plt.xticks([]),plt.yticks([])  
plt.show() 

输出结果如下图所示:

上图很多像素为白色,因为图像求和结果几乎都是255。如果设置的是2*2矩阵,只取四个像素结果要好些。
result = cv2.boxFilter(source, -1, (2,2), normalize=0)



四.高斯滤波

为了克服简单局部平均法的弊端(图像模糊),目前已提出许多保持边缘、细节的局部平滑算法。它们的出发点都集中在如何选择邻域的大小、形状和方向、参数加平均及邻域各店的权重系数等。
图像高斯平滑也是邻域平均的思想对图像进行平滑的一种方法,在图像高斯平滑中,对图像进行平均时,不同位置的像素被赋予了不同的权重。高斯平滑与简单平滑不同,它在对邻域内像素进行平均时,给予不同位置的像素不同的权值,下图的所示的 3 * 3 和 5 * 5 领域的高斯模板。

高斯滤波让临近的像素具有更高的重要度,对周围像素计算加权平均值,较近的像素具有较大的权重值。如下图所示,中心位置权重最高为0.4。

Python中OpenCV主要调用GaussianBlur函数,如下:
dst = cv2.GaussianBlur(src, ksize, sigmaX)
其中,src表示原始图像,ksize表示核大小,sigmaX表示X方向方差。注意,核大小(N, N)必须是奇数,X方向方差主要控制权重。

K(3,3)=[0.050.10.050.10.40.10.050.10.05](5) K(3,3)=\left[ \begin{matrix} 0.05 & 0.1 & 0.05 \\ 0.1 & 0.4 & 0.1 \\ 0.05 & 0.1 & 0.05 \end{matrix} \right] \tag{5}
K(5,5)=[1121113431248421343111211](4) K(5,5)=\left[ \begin{matrix} 1 & 1 & 2 & 1 & 1 \\ 1 & 3 & 4 & 3 & 1 \\ 2 & 4 & 8 & 4 & 2 \\ 1 & 3 & 4 & 3 & 1 \\ 1 & 1 & 2 & 1 & 1 \\ \end{matrix} \right] \tag{4}

代码如下:

#encoding:utf-8
import cv2  
import numpy as np  
import matplotlib.pyplot as plt
 
#读取图片
img = cv2.imread('test01.png')
source = cv2.cvtColor(img,cv2.COLOR_BGR2RGB)
 
#高斯滤波
result = cv2.GaussianBlur(source, (3,3), 0)

#显示图形
titles = ['Source Image', 'GaussianBlur Image']  
images = [source, result]  
for i in xrange(2):  
   plt.subplot(1,2,i+1), plt.imshow(images[i], 'gray')  
   plt.title(titles[i])  
   plt.xticks([]),plt.yticks([])  
plt.show()  

输出结果如下所示:

如果使用15*15的核,则图形将更加模糊。



五.中值滤波

1.概念
在使用邻域平均法去噪的同时也使得边界变得模糊。而中值滤波是非线性的图像处理方法,在去噪的同时可以兼顾到边界信息的保留。选一个含有奇数点的窗口W,将这个窗口在图像上扫描,把窗口中所含的像素点按灰度级的升或降序排列,取位于中间的灰度值来代替该点的灰度值。 例如选择滤波的窗口如下图,是一个一维的窗口,待处理像素的灰度取这个模板中灰度的中值,滤波过程如下:

如下图所示,将临近像素按照大小排列,取排序像素中位于中间位置的值作为中值滤波的像素值。

2.代码
OpenCV主要调用medianBlur()函数实现中值滤波。图像平滑里中值滤波的效果最好。
dst = cv2.medianBlur(src, ksize)
其中,src表示源文件,ksize表示核大小。核必须是大于1的奇数,如3、5、7等。

代码如下所示:

#encoding:utf-8
import cv2  
import numpy as np  
import matplotlib.pyplot as plt
 
#读取图片
img = cv2.imread('test01.png')
 
#高斯滤波
result = cv2.medianBlur(img, 3)

#显示图像
cv2.imshow("source img", img)
cv2.imshow("medianBlur", result)

#等待显示
cv2.waitKey(0)
cv2.destroyAllWindows()

输出结果如下图所示:

常用的窗口还有方形、十字形、圆形和环形。不同形状的窗口产生不同的滤波效果,方形和圆形窗口适合外轮廓线较长的物体图像,而十字形窗口对有尖顶角状的图像效果好。中值滤波对于消除孤立点和线段的干扰十分有用,尤其是对于二进噪声,但对消除高斯噪声的影响效果不佳。对于一些细节较多的复杂图像,可以多次使用不同的中值滤波。

希望文章对大家有所帮助,如果有错误或不足之处,还请海涵。
(By:Eastmount 2018-09-01 早8点 https://blog.csdn.net/Eastmount/)