Opencv实现的陷波滤波器

Posted phoenixdsg

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了Opencv实现的陷波滤波器相关的知识,希望对你有一定的参考价值。

 在本示例中,共设计了三个函数,分别是巴特沃斯滤波器BLPF()、巴特沃斯陷波滤波器notchFilter_BTW()、高斯陷波滤波器notchFilter_GAUSS()

巴特沃斯陷波滤波器参见书上6.4.5选择性滤波器,高斯陷波滤波器参照6.3.3节中的高斯低通滤波器。

参数说明:

rows——滤波器的行数;
cols——滤波器的列数;
D0——频域截止半径;
n——巴特沃斯的阶数;
cvtype——滤波器的数据类型和通道数,默认为双通道浮点数

 

Mat BLPF(int rows,int cols,float D0,int n=1,int cvtype=CV_32FC2)
{
    Mat filt(rows,cols,cvtype,Scalar::all(0));
    int cx=cols/2,cy=rows/2;
    float D02=D0*D0;
    for(int i=0;i<rows;i++)
    {
        for(int j=0;j<cols;j++)
        {
            int u=cx-j,v=cy-i;//中心坐标
            float Duv2=u*u+v*v;//距离中心半径的平方
            float H0=1-1/(1+pow(Duv2/D02,n));
            int u1=u-130,v1=v-130;
            float D2=u1*u1+v1*v1;
            float H1=1-1/(1+pow(D2/D02,n));
            //.data返回的是uchar*型指针,所以要强制转换成浮点数型
            float* p=(float*)(filt.data+i*filt.step[0]+j*filt.step[1]);

            for(int c=0;c<filt.channels();c++)
            {
                p[c]=H1*H0;
            }

        }
    }
    return filt;
}


Mat notchFilter_BTW(int rows,int cols,std::vector<cv::Point> np,
                float* D,int n=1,int cvtype=CV_32FC2)
{
    Mat filt(rows,cols,cvtype,Scalar::all(0));
    int cx=cols/2,cy=rows/2;
    int numNotch=np.size();
    float* D2=D;
    for(int i=0;i<numNotch;i++)
    {
        D2[i]=D[i]*D[i];
    }
    int uk[numNotch],vk[numNotch];//在画面上的实际陷波坐标点
    int u_k[numNotch],v_k[numNotch];//陷波共轭点
    float Dk[numNotch],D_k[numNotch];//陷波半径r
    float Hk[numNotch],H_k[numNotch];

    for(int i=0;i<rows;i++)
    {
        for(int j=0;j<cols;j++)
        {
            int u=cx-j,v=cy-i;//中心坐标
            for(int s=0;s<numNotch;s++)
            {
                uk[s]=u+np[s].x,vk[s]=v+np[s].y;
                Dk[s]=uk[s]*uk[s]+vk[s]*vk[s];//距离中心半径的平方
                Hk[s]=1-1/(1+pow(Dk[s]/D2[s],n));

                u_k[s]=u-np[s].x,v_k[s]=v-np[s].y;
                D_k[s]=u_k[s]*u_k[s]+v_k[s]*v_k[s];
                H_k[s]=1-1/(1+pow(D_k[s]/D2[s],n));
            }
            //.data返回的是uchar*型指针,所以要强制转换成浮点数型
            float* p=(float*)(filt.data+i*filt.step[0]+j*filt.step[1]);

            for(int c=0;c<filt.channels();c++)
            {
                p[c]=Hk[0]*H_k[0];
                for(int k=1;k<numNotch;k++)
                {
                    p[c]*=Hk[k]*H_k[k];
                }
            }

        }
    }
    return filt;
}
Mat notchFilter_GAUSS(int rows,int cols,std::vector<cv::Point> np,
                float* D,int cvtype=CV_32FC2)
{
    Mat filt(rows,cols,cvtype,Scalar::all(0));
    int cx=cols/2,cy=rows/2;

//    float D02=D0*D0;
    int numNotch=np.size();
    float* D2=D;
    for(int i=0;i<numNotch;i++)
    {
        D2[i]=D[i]*D[i];
    }
    int uk[numNotch],vk[numNotch];//在画面上的实际陷波坐标点
    int u_k[numNotch],v_k[numNotch];//陷波共轭点
    float Dk[numNotch],D_k[numNotch];//陷波半径r
    float Hk[numNotch],H_k[numNotch];

    for(int i=0;i<rows;i++)
    {
        for(int j=0;j<cols;j++)
        {
            int u=cx-j,v=cy-i;//中心坐标
            for(int s=0;s<numNotch;s++)
            {
                uk[s]=u+np[s].x,vk[s]=v+np[s].y;
                Dk[s]=uk[s]*uk[s]+vk[s]*vk[s];//距离中心半径的平方
                Hk[s]=1-exp(-Dk[s]/(D2[s]*2));

                u_k[s]=u-np[s].x,v_k[s]=v-np[s].y;
                D_k[s]=u_k[s]*u_k[s]+v_k[s]*v_k[s];
                H_k[s]=1-exp(-D_k[s]/(D2[s]*2));
            }
            //.data返回的是uchar*型指针,所以要强制转换成浮点数型
            float* p=(float*)(filt.data+i*filt.step[0]+j*filt.step[1]);

            for(int c=0;c<filt.channels();c++)
            {
                p[c]=Hk[0]*H_k[0];
                for(int k=1;k<numNotch;k++)
                {
                    p[c]*=Hk[k]*H_k[k];
                }
            }

        }
    }
    return filt;
}

int main()
{
    Point np[]={Point(130,130),Point(90,130),Point(130,100),Point(90,100)};//输入陷波坐标数组
    vector<Point> vnp(np,np+4);
    float D[4]={5,10,15,20};
//    Mat filt1=notchFilter_BTW(500,600,vnp,D,2);
    Mat filt1=notchFilter_GAUSS(500,600,vnp,D);
    Mat fc1;
    extractChannel(filt1,fc1,0);
    imshow("filter ",fc1);

    waitKey();
    return 0;
}

高斯陷波滤波器的演示实例结果如下:

技术分享图片

下面是巴特沃斯陷波滤波器的演示结果:

技术分享图片

 





以上是关于Opencv实现的陷波滤波器的主要内容,如果未能解决你的问题,请参考以下文章

OpenCV 例程200篇103. 陷波带阻滤波器消除周期噪声干扰

youcans 的 OpenCV 例程 200 篇103. 陷波带阻滤波器消除周期噪声干扰

matlab陷波器的设计

窄带陷波滤波器(Notch filter)

语音处理基于matlab GUI汉宁窗FIR陷波滤波器语音信号加噪去噪含Matlab源码 1711期

理想滤波器为啥不能实现?