为何相同dacf数据的FFT（快速傅里叶变换）结果不一致？

我有一个名为yclip.txt的单列数据文件，程序通过ps.acf_dacf()函数计算该数据自相关函数（ACF）的梯度（即dacf），再对dacf执行FFT得到ADP。但出现异常：直接对计算得到的dacf做FFT时，ADP的零频分量无预期峰值；将dacf保存为txt文件后重新加载再执行FFT，ADP结果却正常。两种情况下dacf的标准差一致，尝试过扁平化数组、切换numpy/scipy的FFT函数，差异依然存在。

直接计算dacf并执行FFT的代码（结果异常） #

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from scipy.fft import fft, fftfreq
import pulsar as ps

# 加载数据
data = pd.read_csv('yclip.txt', header = None)
data = np.array(data)
acf, dacf = ps.acf_dacf(data)
n = len(dacf)

# 对dacf做傅里叶变换
adp = np.fft.fft(dacf)
adp = (np.abs(adp[:n//2]))**2

adp = adp/np.max(adp)                       # 归一化

adpfreq = fftfreq(n, 0.65536)[:n//2]        # 生成频率轴

# 绘图
fig, ax = plt.subplots(1, 2)
ax[0].plot(dacf)
ax[0].set_xlabel("Lags")
ax[0].set_ylabel("d/dt(ACF)")
ax[0].set_title("d/dt(ACF)")
ax[1].plot(adpfreq, adp)
ax[1].set_xlabel("Frequency (kHz)")
ax[1].set_ylabel("Power") 
ax[1].set_title("ADP")
plt.suptitle("直接计算dacf后做FFT")
plt.show()

该代码生成的ADP无零频峰值，结果不符合预期。

保存再加载dacf后执行FFT的代码（结果正常） #

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from scipy.fft import fft, fftfreq
import pulsar as ps

# 从文件加载dacf
dacf = pd.read_csv('dacf.txt', header = None)
dacf = np.array(dacf)
n = len(dacf)

# 对dacf做傅里叶变换
adp = np.fft.fft(dacf)
adp = (np.abs(adp[:n//2]))**2
adpy = adp/np.max(adp)                           # 归一化

adpfreq = np.fft.fftfreq(n, 0.65536)[:n//2]      # 生成频率轴

# 绘图
fig, ax = plt.subplots(1, 2)
ax[0].plot(dacf)
ax[0].set_xlabel("Lags")
ax[0].set_ylabel("d/dt(ACF)")
ax[0].set_title("d/dt(ACF)")
ax[1].plot(adpfreq, adp)
ax[1].set_xlabel("Frequency (kHz)")
ax[1].set_ylabel("Power") 
ax[1].set_title("ADP")
plt.suptitle("保存再加载dacf后做FFT")
plt.show()

该代码生成的ADP有正常的零频峰值。

ps.acf_dacf()函数定义 #

def acf_dacf(y):

  if np.any(y != y[0]):
    yacf = acf(y, nlags = len(y))
    dyacf = np.gradient(yacf)
  else:
    yacf = np.ones(len(y))
    dyacf = np.zeros(len(y))

  return yacf, dyacf

原因分析 #

核心问题出在数组维度的隐性差异：

1. 直接计算得到的dacf是二维数组（输入的data是二维数组，np.gradient处理二维数组时会保留其维度，形状为(N, 1)）；而保存再加载后的dacf，虽然np.array(dacf)返回的也是二维，但在绘图和FFT处理时，numpy会自动将单维度的二维数组视为一维处理？不，实际是：numpy的fft函数对二维数组会按列分别计算FFT，而你取adp[:n//2]时，得到的是二维数组的前半部分，后续的平方、归一化操作都会保留二维结构，最终绘图时虽然看起来和一维数组的图一致，但FFT的计算逻辑已经完全不同——二维数组的FFT结果包含了对列维度的计算，导致频谱异常。
2. 验证方法：在第一种代码中，FFT前添加dacf = dacf.flatten()将二维数组转为一维，再执行后续操作，结果会和第二种代码完全一致。
3. 补充说明：标准差一致是因为numpy计算标准差时会自动展平数组，但FFT对维度的敏感性远高于标准差计算，因此维度差异不会体现在标准差上，却会彻底改变FFT结果。

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