图像旋转和频域缩放？答案

【问题标题】：Image rotation and scaling the frequency domain?图像旋转和频域缩放？
【发布时间】：2012-12-06 13:02:08
【问题描述】：

我正在编写一些代码来使用相位相关来恢复测试图像相对于模板的旋转、缩放和平移，a la Reddy & Chatterji 1996。我对原始测试图像进行 FFT 以找到比例因子和旋转角度，但随后我需要 旋转和缩放 测试图像的 FFT 以获得平移。 p>

现在我可以在空间域中应用旋转和缩放然后进行 FFT，但这似乎有点低效 - 是否有可能获得旋转/缩放图像的傅立叶系数直接频域？

编辑 1： 好的，我按照 user1816548 的建议玩了一下。对于 90o 的倍数的角度，我可以得到模糊合理的旋转，尽管图像的极性发生了奇怪的变化。不是 90o 的倍数的角度会给我带来非常滑稽的结果。

编辑 2： 我对图像应用了零填充，并且在旋转它时包裹了 FFT 的边缘。我很确定我正在围绕 FFT 的直流分量旋转，但对于不是 90o 倍数的角度，我仍然会得到奇怪的结果。

示例输出：

可执行的 Numpy/Scipy 代码：

import numpy as np
from scipy.misc import lena
from scipy.ndimage.interpolation import rotate,zoom
from scipy.fftpack import fft2,ifft2,fftshift,ifftshift
from matplotlib.pyplot import subplots,cm

def testFourierRotation(angle):

    M = lena()
    newshape = [2*dim for dim in M.shape]
    M = procrustes(M,newshape)

    # rotate, then take the FFT
    rM = rotate(M,angle,reshape=False)
    FrM = fftshift(fft2(rM))

    # take the FFT, then rotate
    FM = fftshift(fft2(M))
    rFM = rotatecomplex(FM,angle,reshape=False)
    IrFM = ifft2(ifftshift(rFM))

    fig,[[ax1,ax2,ax3],[ax4,ax5,ax6]] = subplots(2,3)

    ax1.imshow(M,interpolation='nearest',cmap=cm.gray)
    ax1.set_title('Original')
    ax2.imshow(rM,interpolation='nearest',cmap=cm.gray)
    ax2.set_title('Rotated in spatial domain')
    ax3.imshow(abs(IrFM),interpolation='nearest',cmap=cm.gray)
    ax3.set_title('Rotated in Fourier domain')
    ax4.imshow(np.log(abs(FM)),interpolation='nearest',cmap=cm.gray)
    ax4.set_title('FFT')
    ax5.imshow(np.log(abs(FrM)),interpolation='nearest',cmap=cm.gray)
    ax5.set_title('FFT of spatially rotated image')
    ax6.imshow(np.log(abs(rFM)),interpolation='nearest',cmap=cm.gray)
    ax6.set_title('Rotated FFT')
    fig.tight_layout()

    pass

def rotatecomplex(a,angle,reshape=True):
    r = rotate(a.real,angle,reshape=reshape,mode='wrap')
    i = rotate(a.imag,angle,reshape=reshape,mode='wrap')
    return r+1j*i

def procrustes(a,target,padval=0):
    b = np.ones(target,a.dtype)*padval
    aind = [slice(None,None)]*a.ndim
    bind = [slice(None,None)]*a.ndim
    for dd in xrange(a.ndim):
        if a.shape[dd] > target[dd]:
            diff = (a.shape[dd]-target[dd])/2.
            aind[dd] = slice(np.floor(diff),a.shape[dd]-np.ceil(diff))
        elif a.shape[dd] < target[dd]:
            diff = (target[dd]-a.shape[dd])/2.
            bind[dd] = slice(np.floor(diff),target[dd]-np.ceil(diff))
    b[bind] = a[aind]
    return b

【问题讨论】：

我其实也对此感兴趣。如果您查看空间旋转图像的 FFT，您会看到从水平轴上的 400 到垂直中心线的尖峰。但这些线不在旋转的 FFT 中。另外，您只绘制幅度，相位图如何？你能把这些也发一下吗？

标签： math image-processing numpy fft transformation

【解决方案1】：

我不确定这是否已经解决，但我相信我已经解决了您关于第三个图中观察到的效果的问题：

您观察到的这种奇怪效果是由于您实际计算 FFT 的来源。本质上，FFT 从数组的第一个像素 M[0][0] 开始。但是，您围绕 M[size/2+1,size/2+1] 定义旋转，这是正确的方法，但错误:)。傅里叶域是从M[0][0] 计算出来的！如果您现在在傅立叶域中旋转，您将围绕M[0][0] 而不是围绕M[size/2+1,size/2+1] 旋转。我无法完全解释这里到底发生了什么，但你也会得到我以前得到的相同效果。为了在傅里叶域中旋转原始图像，您必须首先将 2D fftShift 应用于原始图像 M，然后计算 FFT、旋转、IFFT，然后应用 ifftShift。这样你的图像旋转中心和傅里叶域的中心就会同步。

AFAI 记得我们还在两个单独的数组中旋转实部和虚部，然后将它们合并。我们还在复数上测试了各种插值算法，但效果不大：)。它在我们的包裹中pytom。

但是，这可能会损失很多，但另外两个移位并不是很快，除非您指定一些时髦的数组索引算法。

【讨论】：

【解决方案2】：

好吧，旋转和缩放的图像会导致旋转和缩放（反比例）傅立叶变换。

还要注意，旋转和缩放在像素数上都是线性的，而 FFT 是 O(w*logw*h*logh)，所以最终它实际上并没有那么昂贵。

【讨论】：

您好，感谢您的回复。我可能最终只会在空间域中进行转换，但我想了解它在傅立叶域中的工作方式 - 请参阅我的编辑。
确保您在正确的位置进行旋转 - 根据您的约定，它是中心或左上角。基本上，图像的最亮点是您应该旋转的地方。您可以通过拍摄一些相当常规的图像并转换图像及其 FT 并查看它们之间的关系来轻松检查。
嗯，还是有点难受——我肯定是在围绕直流分量旋转 FFT。这可能与我在对其应用旋转时如何填充 FFT 以处理边缘有关吗？
您是在完全围绕 DCT 旋转还是大约旋转？因为数组的长度是均匀的（例如 8），所以如果它不围绕 DCT 像素的中心旋转，那么围绕 DCT 的旋转将是错误的。我认为您获得的结果与您围绕 DCT 像素的 corner 旋转所获得的结果差不多......
好的，我尝试围绕 DCT 像素的精确中心旋转（fftshift 将其放置在 ceil((dim-1)/2.) 以获取偶数和奇数尺寸）。即使 FFT 中最亮的像素在旋转后保持在同一位置，但当我进行 IFFT 时，我仍然会得到奇怪的“万花筒”效果。

【解决方案3】：

我意识到这很晚了，但我只是想在这里回答这个问题，因为我回顾了我关于移位不变性的基本知识。问题是您在旋转之前扩展了傅立叶空间（如考虑混叠）。查看旋转图像的 FT：轴向尖峰（别名）出现在傅里叶旋转的 IFT 中没有的边缘。

您应该旋转然后处理锯齿。因为您正在考虑混叠（以周期 = 像素数循环傅立叶空间），然后通过旋转来消除这种努力，所以您会导致混叠出现在最终图像中。本质上，您是在分散傅立叶别名，因此将图像空间别名放在一起。

90 度旋转可以顺利进行，因为没有锯齿； k空间的角完美匹配。

【讨论】：