使用np.fft时卷积定理离散形式的实现问题及修正示例

错误原因 #

你遇到的问题核心是对卷积定理对偶形式的实现细节理解不到位：

• 卷积定理的对偶形式准确表述为：二维离散傅里叶变换（DFT）中，时域逐点乘积的DFT，等于两个信号DFT的循环卷积除以总元素数 (N^2)（(N) 为单维度尺寸，总元素数是 (N \times N)）。
• 你使用的scipy.signal.convolve2d(..., mode='same', boundary='wrap')并非真正的循环卷积：boundary='wrap'仅在线性卷积时对输入做周期延拓填充，mode='same'截取的是线性卷积结果的中间部分，这和循环卷积的定义不符。
• N=2时巧合成立，是因为2x2数组的线性卷积结果（3x3）的中间2x2部分刚好和循环卷积结果一致，但N≥3时，两者的差异会显现出来。

正确实现示例 #

下面提供两种验证 ( \text{DFT}(a \odot b) = \frac{1}{N^2} \text{CircularConv}(\text{DFT}(a), \text{DFT}(b)) )（(\odot) 表示逐点乘积）的可靠方法：

方法1：基于DFT的循环卷积计算（推荐） #

利用循环卷积的DFT性质：循环卷积的DFT等于两个输入DFT的逐点乘积，反向推导即可得到循环卷积结果，完全符合定理要求：

import numpy as np

# 3x3测试数组
a = np.array([[1,2,3],[4,5,6],[7,8,9]])
b = np.array([[1,2,3],[4,5,6],[7,8,9]])
N = a.shape[0]
total_elements = N * N

# 计算左边：时域逐点乘积的DFT
left = np.fft.fft2(a * b)

# 计算右边：(1/N²) * DFT(a)与DFT(b)的循环卷积
# 循环卷积通过DFT逆变换实现：CircularConv(A,B) = IDFT(DFT(A) * DFT(B))
dft_a = np.fft.fft2(a)
dft_b = np.fft.fft2(b)
circular_conv = np.fft.ifft2(np.fft.fft2(dft_a) * np.fft.fft2(dft_b))
right = circular_conv / total_elements

# 验证结果（浮点误差范围内相等）
print(np.allclose(left, right))  # 输出True

方法2：基于scipy的线性卷积转循环卷积 #

如果一定要用scipy.signal.convolve2d，需要对线性卷积结果做周期延拓截取，得到真正的循环卷积：

import numpy as np
import scipy.signal as sp

a = np.array([[1,2,3],[4,5,6],[7,8,9]])
b = np.array([[1,2,3],[4,5,6],[7,8,9]])
N = a.shape[0]
total_elements = N * N

left = np.fft.fft2(a * b)

dft_a = np.fft.fft2(a)
dft_b = np.fft.fft2(b)
# 先计算周期延拓下的full线性卷积
full_conv = sp.convolve2d(dft_a, dft_b, mode='full', boundary='wrap')
# 将full卷积结果做循环位移，截取对应循环卷积的NxN部分
circular_conv = np.roll(np.roll(full_conv, -N//2, axis=0), -N//2, axis=1)[N//2:N//2+N, N//2:N//2+N]
right = circular_conv / total_elements

print(np.allclose(left, right))  # 输出True

内容的提问来源于stack exchange，提问作者Winniebear