Python函数频谱图像异常：2-5号函数频谱黑屏及log10除零警告

核心问题原因 #

1. 除零警告与黑色频谱：计算频谱时直接对DFT结果取log10，而DFT中大量元素为0，log10(0)会产生-inf值。Matplotlib显示时，这些-inf会被当作最小值，导致绝大多数区域显示为黑色，仅少数非零点因数值差异无法被肉眼识别。
2. 频谱可视化位置问题：默认的DFT结果低频分量位于图像角落，未居中显示，即使修复log问题，也难以直观观察频谱特征。

修复步骤 #

1. 避免log10计算中的零值 #

在取绝对值后添加一个极小的偏移量（如1e-10），防止出现log10(0)的情况：

magnitude_spectrum = 20 * np.log10(np.abs(dft_signal) + 1e-10)

2. 频谱中心化（fftshift） #

使用np.fft.fftshift将低频分量移到图像中心，便于观察频谱分布：

shifted_spectrum = np.fft.fftshift(magnitude_spectrum)

3. 统一可视化参数 #

为imshow指定cmap='gray'保证灰度显示，必要时可通过vmin和vmax调整显示范围，增强对比度。

修改后的完整代码 #

import numpy as np
import matplotlib.pyplot as plt

# Dimensiunile imaginii
M, N = 256, 256

# Generare coordonate n1, n2
n1 = np.arange(0, M)
n2 = np.arange(0, N)
n1, n2 = np.meshgrid(n1, n2)

# Definirea celor cinci funcții date
def signal_function_1(n1, n2):
    return np.sin(2 * np.pi * n1 + 3 * np.pi * n2)

def signal_function_2(n1, n2):
    return np.sin(4 * np.pi * n1) + np.cos(6 * np.pi * n2)

def signal_function_3(n1, n2):
    Y = np.zeros_like(n1, dtype=complex)
    Y[0, 5] = Y[0, N-5] = 1
    return np.real(np.fft.ifft2(Y))

def signal_function_4(n1, n2):
    Y = np.zeros_like(n1, dtype=complex)
    Y[5, 0] = Y[N-5, 0] = 1
    return np.real(np.fft.ifft2(Y))

def signal_function_5(n1, n2):
    Y = np.zeros_like(n1, dtype=complex)
    Y[5, 5] = Y[N-5, N-5] = 1
    return np.real(np.fft.ifft2(Y))

# Calculul semnalelor
signal_1 = signal_function_1(n1, n2)
signal_2 = signal_function_2(n1, n2)
signal_3 = signal_function_3(n1, n2)
signal_4 = signal_function_4(n1, n2)
signal_5 = signal_function_5(n1, n2)

# Calculul DFT bidimensional
dft_signal_1 = np.fft.fft2(signal_1)
dft_signal_2 = np.fft.fft2(signal_2)
dft_signal_3 = np.fft.fft2(signal_3)
dft_signal_4 = np.fft.fft2(signal_4)
dft_signal_5 = np.fft.fft2(signal_5)

# 修复频谱计算：添加偏移量避免log10(0)，并中心化频谱
def compute_spectrum(dft_signal):
    mag = np.abs(dft_signal)
    # 添加极小值防止log10(0)
    mag_spectrum = 20 * np.log10(mag + 1e-10)
    # 中心化频谱
    return np.fft.fftshift(mag_spectrum)

magnitude_spectrum_1 = compute_spectrum(dft_signal_1)
magnitude_spectrum_2 = compute_spectrum(dft_signal_2)
magnitude_spectrum_3 = compute_spectrum(dft_signal_3)
magnitude_spectrum_4 = compute_spectrum(dft_signal_4)
magnitude_spectrum_5 = compute_spectrum(dft_signal_5)

# Vizualizare pentru prima funcție
plt.figure(figsize=(12, 8))

plt.subplot(221), plt.imshow(signal_1, cmap='gray'), plt.title('Imagine 1')
plt.subplot(222), plt.imshow(magnitude_spectrum_1, cmap='gray'), plt.title('Spectru 1 (Centrat)')

# Vizualizare pentru a doua funcție
plt.subplot(223), plt.imshow(signal_2, cmap='gray'), plt.title('Imagine 2')
plt.subplot(224), plt.imshow(magnitude_spectrum_2, cmap='gray'), plt.title('Spectru 2 (Centrat)')

plt.tight_layout()
plt.show()

plt.figure(figsize=(12, 8))
# Vizualizare pentru a treia funcție
plt.subplot(221), plt.imshow(signal_3, cmap='gray'), plt.title('Imagine 3')
plt.subplot(222), plt.imshow(magnitude_spectrum_3, cmap='gray'), plt.title('Spectru 3 (Centrat)')

plt.tight_layout()
plt.show()

# Vizualizare pentru a patra funcție
plt.figure(figsize=(12, 8))

plt.subplot(221), plt.imshow(signal_4, cmap='gray'), plt.title('Imagine 4')
plt.subplot(222), plt.imshow(magnitude_spectrum_4, cmap='gray'), plt.title('Spectru 4 (Centrat)')

# Vizualizare pentru a cincea funcție
plt.subplot(223), plt.imshow(signal_5, cmap='gray'), plt.title('Imagine 5')
plt.subplot(224), plt.imshow(magnitude_spectrum_5, cmap='gray'), plt.title('Spectru 5 (Centrat)')

plt.tight_layout()
plt.show()

额外说明 #

• 对于signal3-5，它们是从频域特定点逆变换得到的信号，再次FFT后理论上会还原出频域的点。修复后可以清晰看到这些点在中心化频谱中的位置。
• 添加plt.tight_layout()可以自动调整子图间距，避免标题重叠。

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