python - 显示某些函数频谱的正确图像时出错

Question

我有一个 Python 代码，可以为五个不同的函数生成图像及其光谱，但存在一个问题，即函数 2 到 5 为光谱生成黑色图像，而不是所需的输出。

问题：

RuntimeWarning: divide by zero encountered in log10
  magnitude_spectrum_2 = 20 * np.log10(np.abs(dft_signal_2))

signal_3-5 相同。

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)


magnitude_spectrum_1 = 20 * np.log10(np.abs(dft_signal_1))
magnitude_spectrum_2 = 20 * np.log10(np.abs(dft_signal_2))
magnitude_spectrum_3 = 20 * np.log10(np.abs(dft_signal_3))
magnitude_spectrum_4 = 20 * np.log10(np.abs(dft_signal_4))
magnitude_spectrum_5 = 20 * np.log10(np.abs(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), plt.title('Spectru 1')

# 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), plt.title('Spectru 2')

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), plt.title('Spectru 3')

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), plt.title('Spectru 4')

# 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), plt.title('Spectru 5')

plt.show()

我尝试解决这个问题，但没有成功。

Answer 1

为了避免被零除的错误，可以只计算非零样本的对数：

magnitude_spectrum_1 = 20 * np.log10(np.abs(dft_signal_1), where=dft_signal_1 > 0)
magnitude_spectrum_2 = 20 * np.log10(np.abs(dft_signal_2), where=dft_signal_2 > 0)
magnitude_spectrum_3 = 20 * np.log10(np.abs(dft_signal_3), where=dft_signal_3 > 0)
magnitude_spectrum_4 = 20 * np.log10(np.abs(dft_signal_4), where=dft_signal_4 > 0)
magnitude_spectrum_5 = 20 * np.log10(np.abs(dft_signal_5), where=dft_signal_5 > 0)