【问题标题】:Convolving image with kernel in Fourier domain傅里叶域中带核的卷积图像
【发布时间】:2019-07-19 13:07:37
【问题描述】:

我在我的图像和卷​​积核周围使用零填充,将它们转换为傅里叶域,然后将它们反相以获得卷积图像,请参见下面的代码。然而,结果是错误的。我期待一个模糊的图像,但输出是四个移位的季度。为什么输出错误,如何修复代码?

输入图片:

卷积结果:

from PIL import Image,ImageDraw,ImageOps,ImageFilter
import numpy as np 
from scipy import fftpack
from copy import deepcopy
import imageio
## STEP 1 ##
im1=Image.open("pika.jpeg")
im1=ImageOps.grayscale(im1)
im1.show()
print("s",im1.size)
## working on this image array
im_W=np.array(im1).T
print("before",im_W.shape)
if(im_W.shape[0]%2==0):
im_W=np.pad(im_W, ((1,0),(0,0)), 'constant')
if(im_W.shape[1]%2==0):
im_W=np.pad(im_W, ((0,0),(1,0)), 'constant')
print("after",im_W.shape)
Boxblur=np.array([[1/9,1/9,1/9],[1/9,1/9,1/9],[1/9,1/9,1/9]])
dim=Boxblur.shape[0]

##padding before frequency domain multipication
pad_size=(Boxblur.shape[0]-1)/2
pad_size=int(pad_size)
##padded the image(starts here)

p_im=np.pad(im_W, ((pad_size,pad_size),(pad_size,pad_size)), 'constant')
t_b=(p_im.shape[0]-dim)/2
l_r=(p_im.shape[1]-dim)/2
t_b=int(t_b)
l_r=int(l_r)

##padded the image(ends here)

## padded the kernel(starts here)
k_im=np.pad(Boxblur, ((t_b,t_b),(l_r,l_r)), 'constant')
print("hjhj",k_im)
print("kernel",k_im.shape)

##fourier transforms image and kernel
fft_im = fftpack.fftshift(fftpack.fft2(p_im))
fft_k  = fftpack.fftshift(fftpack.fft2(k_im))
con_in_f=fft_im*fft_k
ifft2 = abs(fftpack.ifft2(fftpack.ifftshift(con_in_f)))
convolved=(np.log(abs(ifft2))* 255 / np.amax(np.log(abs(ifft2)))).astype(np.uint8)
final=Image.fromarray(convolved.T)
final.show()
u=im1.filter(ImageFilter.Kernel((3,3), [1/9,1/9,1/9,1/9,1/9,1/9,1/9,1/9,1/9], scale=None, offset=0))
u.show()

【问题讨论】:

  • 查看您的代码,您想删除您在那里执行的所有其他 fft 班次(您只需要我指出的那个),它们是多余的。我还看到一个log,我不知道为什么会出现,它会改变你的颜色。
  • 当我们零填充时,内核应该在中心还是左上角,我很困惑如何将空间域中的内核具体转换为频域。跨度>

标签: python image-processing scipy fft convolution


【解决方案1】:

离散傅里叶变换 (DFT) 以及扩展的 FFT(计算 DFT)在两者的第一个元素(对于图像,左上角像素)中具有原点输入和输出。这就是我们经常在输出上使用fftshift函数的原因,以便将原点转移到我们更熟悉的位置(图像中间)。

这意味着我们需要将一个 3x3 均匀加权模糊内核转换为如下所示,然后再将其传递给 FFT 函数:

1/9  1/9  0  0  ... 0  1/9
1/9  1/9  0  0  ... 0  1/9
  0    0  0  0  ... 0    0
...  ...               ...
  0    0  0  0  ... 0    0
1/9  1/9  0  0  ... 0  1/9

即内核的中间位于图像的左上角,中间上方和左侧的像素环绕并出现在图像的右端和下端。

我们可以使用ifftshift 函数来做到这一点,该函数在填充后应用于内核。填充内核时,我们需要注意原点(内核的中间)位于内核映像k_im 内的位置k_im.shape // 2(整数除法)。最初原点位于[3,3]//2 == [1,1]。通常,我们匹配大小的图像大小是偶数,例如[256,256]。原点位于[256,256]//2 == [128,128]。这意味着我们需要在左侧和右侧(以及底部和顶部)填充不同的数量。我们需要小心计算这个填充:

sz = img.shape  # the sizes we're matching
kernel = np.ones((3,3)) / 9
sz = (sz[0] - kernel.shape[0], sz[1] - kernel.shape[1])  # total amount of padding
kernel = np.pad(kernel, (((sz[0]+1)//2, sz[0]//2), ((sz[1]+1)//2, sz[1]//2)), 'constant')
kernel = fftpack.ifftshift(kernel)

请注意,输入图像img 不需要填充(但如果您想强制使用 FFT 更便宜的大小,您可以这样做)。也不需要在乘法之前将fftshift应用于FFT的结果,然后立即反转这个移位,这些移位是多余的。仅当您想显示傅里叶域图像时才应使用fftshift。最后,对过滤后的图像应用对数缩放是错误的。

生成的代码是(我使用 pyplot 进行显示,根本没有使用 PIL):

import numpy as np
from scipy import misc
from scipy import fftpack
import matplotlib.pyplot as plt

img = misc.face()[:,:,0]

kernel = np.ones((3,3)) / 9
sz = (img.shape[0] - kernel.shape[0], img.shape[1] - kernel.shape[1])  # total amount of padding
kernel = np.pad(kernel, (((sz[0]+1)//2, sz[0]//2), ((sz[1]+1)//2, sz[1]//2)), 'constant')
kernel = fftpack.ifftshift(kernel)

filtered = np.real(fftpack.ifft2(fftpack.fft2(img) * fftpack.fft2(kernel)))
plt.imshow(filtered, vmin=0, vmax=255)
plt.show()

请注意,我取的是逆 FFT 的实部。虚部应仅包含非常接近零的值,这是计算中舍入误差的结果。取绝对值虽然很常见,但却是不正确的。例如,您可能希望对包含负值的图像应用过滤器,或应用产生负值的过滤器。在这里取绝对值会产生人工制品。如果逆 FFT 的输出包含显着不同于零的虚值,则填充过滤内核的方式存在错误。

还要注意这里的内核很小,因此模糊效果也很小。为了更好地看到模糊效果,请制作更大的内核,例如np.ones((7,7)) / 49

【讨论】:

  • 这是完美的,这对所有类型的内核都适用吗,比如 Sobel 内核或任何内核?
  • @ShivamThakur:是的,会的。我在这里忘记说的一件事:通过 FFT 进行的卷积仅对较大的内核有效。对于像这里的 3x3 内核或 Sobel 内核这样的小内核,最好不要使用 FFT。
  • 我做了 Sobel 和其他 3x3 内核,似乎我得到了预期的完美结果,我在傅里叶域中进行卷积以加快处理速度。
  • @ShivamThakur:正如我所说,对于 Sobel 或其他 3x3 内核,这是非常低效的。只需使用 scipy.signal.convolve2d 之类的东西。给你的代码计时,假设它更快,不要以这种方式实现卷积!
  • 我不允许使用除 scipys fft 之外的任何库,在卷积之前我们是否必须水平和垂直翻转内核?
猜你喜欢
  • 2015-12-17
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
相关资源
最近更新 更多