应用Butterworth滤波器去除正弦噪声

Question

我试图消除这幅图像中的正弦噪声：

​

​

以下是它的DFT谱(在应用了对数和任意强度标度之后)：

​

​

我已经有了一个巴特沃斯过滤器来应用于这个图像。它会击倒中频的峰值。加载后，我会小心地将它从0.255扩展到0..1.0。这是过滤器：

​

​

结果并不好：

​

​

我的问题：

• 为什么图像中仍然存在大量的噪声？
• 为什么结果比原来的图像更暗？过滤器显然没有触及直流项，所以我希望平均强度是相同的。
• 为什么过滤器只去除一些峰值？它来自于一本教科书，所以我倾向于相信它是正确的，但是光谱中还有其他的峰值--它们也是噪音的一部分吗？我试着用同心滤镜把它们移除，但效果不佳，使图像变暗，无法识别。

我从冈萨雷斯和伍兹的“数字图像处理”一书中摘取了这张照片(裁剪)和滤镜。在它们的例子中，周期噪声通过滤波完全消除，图像的平均强度保持不变。

我的下载图像和过滤器，DFT，过滤，IDFT的源代码如下：

import cv

def unshift_crop(comp, width, height):
    result = cv.CreateImage((width, height), cv.IPL_DEPTH_8U, 1)
    for x in range(height):
        for y in range(width):
            real, _, _, _ = cv.Get2D(comp, x, y)
            real = int(real) * ((-1)**(x+y))
            cv.Set2D(result, x, y, cv.Scalar(real))
    return result

def load_filter(fname):
    loaded = cv.LoadImage(fname, cv.CV_LOAD_IMAGE_GRAYSCALE)
    flt = cv.CreateImage(cv.GetSize(loaded), cv.IPL_DEPTH_32F, 2)
    width, height = cv.GetSize(loaded)
    for i in range(width*height):
        px, _, _, _ = cv.Get1D(loaded, i)
        #cv.Set1D(flt, i, cv.Scalar(px/255.0, 0))
        cv.Set1D(flt, i, cv.Scalar(px/255.0, px/255.0))
    return flt

if __name__ == '__main__':
    import sys
    fname, filt_name, ofname = sys.argv[1:]
    img = cv.LoadImage(fname, cv.CV_LOAD_IMAGE_GRAYSCALE)
    width, height = cv.GetSize(img)
    src = cv.CreateImage((width*2, height*2), cv.IPL_DEPTH_32F, 2)
    dst = cv.CreateImage((width*2, height*2), cv.IPL_DEPTH_32F, 2)
    cv.SetZero(src)
    for x in range(height):
        for y in range(width):
            px, _, _, _ = cv.Get2D(img, x, y)
            px = float(px) * ((-1) ** (x+y))
            cv.Set2D(src, x, y, cv.Scalar(px, 0))
    cv.DFT(src, dst, cv.CV_DXT_FORWARD)
    flt = load_filter(filt_name)
    cv.Mul(dst, flt, src)
    cv.DFT(src, dst, cv.CV_DXT_INV_SCALE)
    result = unshift_crop(dst, width, height)
    cv.SaveImage(ofname, result)

编辑

在原始源中存在一个错误，其中过滤器假想组件被加载为零。这就是导致结果图像看起来比实际更暗的原因。我已经解决了这个问题，并评论了相关的内容。

使用@0x69提供的固定源和过滤器(是的，我知道它不是真正的巴特沃斯过滤器，但在这个阶段我很乐意尝试任何东西)，结果是：

​

​

比我一开始要做的要好，但还是没有我所期望的那么好。有人能打败这个吗？我怀疑，在更多的缺口，以消除剩余的高峰，可能会有一些好处。

编辑2

我已经联系过作者了。这是他们的回应：

问题是实验中使用的图像是浮点，而书中显示的图像(以及下载中提供的原始图像)是8位。这是印刷等所需的。 为了复制实验，你必须从无噪音图像开始，然后添加你自己的噪音。

Answer 1

我试过使用这样的修改过滤器：

我得到的是这个->

我不能完全解释结果，但我最好的猜测是，正弦噪声通过与主图像信号相互作用而产生第二，第三，.谐波噪声波结果也很不理想，似乎还有一些噪声谐波残留在这里.顺便说一句，谢谢你提出有趣的问题。

编辑：

我第二次尝试改进过滤器。过滤器：

过滤结果：

看来这一次没有明显的正弦波噪声图案。

Answer 2

我记得几年前在一个图像处理课程中玩过这个图像，我得到了和你一样的结果。

我不知道教科书的作者是如何得到他们在书中显示的图像的，但是他们一定做了更多的事情，然后应用了巴特沃斯过滤器。正如你提到的，有更多的山峰，所以他们有可能使用更多的巴特沃斯过滤器来去除这些。

然而，对我来说，图像的平均值并没有改变。你试过计算这两幅图像的平均值并进行比较吗？这可能是因为显示时的缩放会导致图像变暗。

Answer 3

我觉得你把事情弄得有点复杂。如果你愿意的话就在Matlab上做吧。

它给你很好的结果。

%  Question: Filtration in Frequency Domain

im = imread('applo_noisy.tif');
FT = fft2(double(im));
FT1 = fftshift(FT);%finding spectrum
%imtool(abs(FT1),[]);

m = size(im,1);
n = size(im,2);

t = 0:pi/15:2*pi;
xc=(m+150)/2; % point around which we filter image
yc=(n-150)/2;

r=200;   
r1 = 40;

xcc = r*cos(t)+xc;
ycc =  r*sin(t)+yc;
xcc1 = r1*cos(t)+xc;
ycc1 =  r1*sin(t)+yc;

mask = poly2mask(double(xcc),double(ycc), m,n); % Convert region-of-interest polygon to mask
mask1 = poly2mask(double(xcc1),double(ycc1), m,n); % generating mask for filtering

% a=51;
% b=5;
% a(b) = 0; % 51 0 0 0 0

mask(mask1)=0;

FT2=FT1;
FT2(mask)=0;%cropping area or bandreject filtering

output = ifft2(ifftshift(FT2));
imtool(output,[]);