为什么我的 2D fft 中会出现多行零？

Question

我正在尝试复制论文中的结果。

“沿恒定纬度（东西）和经度（南北）截面的空间和时间二维傅里叶变换（2D-FT）用于表征40度以南模拟通量变化的光谱。 " - Lenton 等人（2006 年） 公布的数字显示“2D-FT 方差的对数”。

我试图创建一个由相似数据的季节性周期以及噪声组成的数组。我将噪声定义为原始数组减去信号数组。

这是我用来绘制信号阵列在纬度上平均的 2D-FT 的代码：

import numpy as np

from numpy import ma
from matplotlib import pyplot as plt
from Scientific.IO.NetCDF import NetCDFFile 

### input directory 
indir = '/home/nicholas/data/'

### get the flux data which is in
### [time(5day ave for 10 years),latitude,longitude]
nc = NetCDFFile(indir + 'CFLX_2000_2009.nc','r')
cflux_southern_ocean  = nc.variables['Cflx'][:,10:50,:]
cflux_southern_ocean  = ma.masked_values(cflux_southern_ocean,1e+20) # mask land
nc.close()

cflux = cflux_southern_ocean*1e08 # change units of data from mmol/m^2/s

### create an array that consists of the seasonal signal fro each pixel
year_stack = np.split(cflux, 10, axis=0)
year_stack = np.array(year_stack)
signal_array = np.tile(np.mean(year_stack, axis=0), (10, 1, 1))
signal_array = ma.masked_where(signal_array > 1e20, signal_array) # need to mask
### average the array over latitude(or longitude)
signal_time_lon = ma.mean(signal_array, axis=1)

### do a 2D Fourier Transform of the time/space image
ft = np.fft.fft2(signal_time_lon)
mgft = np.abs(ft)   
ps = mgft**2 
log_ps = np.log(mgft)
log_mgft= np.log(mgft)

ft 的每第二行完全由零组成。为什么是这样？ 是否可以在信号中添加一个随机的小数字以避免这种情况。

signal_time_lon = signal_time_lon + np.random.randint(0,9,size=(730, 182))*1e-05

编辑：添加图像并阐明含义

rfft2 的输出看起来仍然是一个复数数组。使用 fftshift 将图像的边缘移动到中心；无论如何，我仍然有一个功率谱。我希望我得到零行的原因是我为每个像素重新创建了时间序列。 ft[0, 0] 像素包含信号的平均值。所以 ft[1, 0] 对应于一个正弦曲线，在起始图像的行中的整个信号上都有一个周期。

这是使用以下代码的起始图像：

plt.pcolormesh(signal_time_lon); plt.colorbar(); plt.axis('tight')

这是使用以下代码的结果：

ft = np.fft.rfft2(signal_time_lon)
mgft = np.abs(ft)   
ps = mgft**2                    
log_ps = np.log1p(mgft)
plt.pcolormesh(log_ps); plt.colorbar(); plt.axis('tight')

图像中可能不清楚，但只有每隔一行包含完全零。每十个像素 (log_ps[10, 0]) 就是一个高值。其他像素（log_ps[2, 0]、log_ps[4, 0] 等）的值非常低。

Answer 1

考虑以下示例：

In [59]: from scipy import absolute, fft

In [60]: absolute(fft([1,2,3,4]))
Out[60]: array([ 10.        ,   2.82842712,   2.        ,   2.82842712])

In [61]: absolute(fft([1,2,3,4, 1,2,3,4]))
Out[61]: 
array([ 20.        ,   0.        ,   5.65685425,   0.        ,
         4.        ,   0.        ,   5.65685425,   0.        ])

In [62]: absolute(fft([1,2,3,4, 1,2,3,4, 1,2,3,4]))
Out[62]: 
array([ 30.        ,   0.        ,   0.        ,   8.48528137,
         0.        ,   0.        ,   6.        ,   0.        ,
         0.        ,   8.48528137,   0.        ,   0.        ])

如果X[k] = fft(x) 和Y[k] = fft([x x])，则Y[2k] = 2*X[k] 对应k in {0, 1, ..., N-1}，否则为零。

因此，我会调查您的signal_time_lon 是如何平铺的。这可能就是问题所在。