基于Welch周期图的湍流密度谱积分无法还原流向速度波动问题

我尝试使用Welch周期图计算湍流密度谱，结果看似正常，但通过谱积分无法得到预期的流向速度波动值。以下是我的代码：

sr = 833.3  # Sampling rate i.e. samples per second (1/dt with dt=0.0012)
Nysquid = sr/2

# Perform Welch's periodogram with Hann window and 50% overlap
seg_len = 4
segment = int(seg_len * sr)
myhann = signal.get_window('hann', segment)

# obtain Power density with Hann window and 50% overlap
myparams2 = dict(fs = sr, nperseg = segment, window = myhann, 
                 noverlap = segment/2, scaling = 'density', return_onesided=True)

# xdat is an array composed of the streamwise velocity at a given
# location and for many times.
freqn, xn = signal.welch(x = xdat, **myparams2)

我通过实测数据计算速度波动的代码如下：

U = xdat
U_std = np.std(U)
print("uu=", 0.5 * U_std**2) 

理论上我应该能通过积分谱还原该值，但未能成功。请问使用能量谱还原速度波动的正确表达式是什么？

核心问题出在窗口能量修正和功率谱密度的积分逻辑上，以下是正确的处理方式：

1. 明确功率谱密度的积分关系 #

当scaling='density'时，Welch返回的是功率谱密度(PSD)，单位为(速度)^2 / Hz。对单侧谱（return_onesided=True）积分的结果等于信号的均方值（即U_std²），而你需要的0.5*U_std²是湍流动能波动，因此积分后取一半即可。

2. 窗口能量修正 #

Hann窗口会引入能量损失，因为窗口系数的平方和不等于样本数。你有两种修正方式：

方式一：窗口初始化时自动归一化 #

调用signal.get_window时添加fftbins=True参数，让窗口满足FFT能量归一化要求：

myhann = signal.get_window('hann', segment, fftbins=True)

方式二：手动修正谱值 #

计算窗口的能量修正因子，对得到的PSD进行修正：

window_correction = 1 / (np.sum(myhann**2) / segment)
xn_corrected = xn * window_correction

3. 正确的积分代码 #

完成修正后，通过积分还原速度波动：

# 计算频率间隔
df = freqn[1] - freqn[0]
# 积分PSD得到均方值
uu_from_spectrum = np.trapz(xn_corrected, freqn)  # 或用 np.sum(xn_corrected * df)
# 输出湍流动能波动值
print("uu from spectrum=", 0.5 * uu_from_spectrum)

关键验证 #

修正后，积分得到的0.5*uu_from_spectrum应该和你直接计算的0.5*U_std²基本一致，误差来源于Welch方法的统计特性和窗口的微小影响。

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