1

我正在尝试将 MATLAB 程序转换为 python。我在设置交叉功率谱密度函数并使用 Matlab 获得匹配结果时遇到问题。

MATLAB代码中使用的函数编写如下：

[Pxy,f] = cpsd(x,y,M,round(M/2),M,fs);

在我的代码中可用的文档中，我读到：M = 128（FFT 点数）和 fs = 25.0（采样频率 [Hz]）。x 和 y 是加速度数据的行向量 1x751。

使用的函数有六个参数，所以我假设这[pxy,f] = cpsd(x,y,window,noverlap,f,fs)是程序员打算从 MATLAB 库中调用的函数，它是文档中唯一可用的具有六个参数的函数（请参见此处）

此函数返回 f 中指定的频率处的交叉功率谱密度估计值。 （令我烦恼的是 f 没有定义为频率，但变量 M 被传递到那里，它是 FFT 点的数量，但让我们假设这不是一个错误）。

现在，我想用scipy.signal.csd来转换这个函数，但是有两个问题：

  1. 窗口在 MatLab 中定义为整数，但 scipy 的 csd 只允许窗口作为元组、字符串或类似数组的对象；
  2. 在 scipy 的 csd 中，没有一个参数允许返回特定频率的交叉功率谱密度估计。

对于数字 1，我定义了一个窗口如下： window = hamming(M, sym=False) 我选择汉明窗口，因为它是在 MATLAB 的 csd 中将窗口作为整数传递时指定的默认窗口（“如果窗口是整数，则 cpsd 将 x 和 y 分为长度窗口的片段和每个片段的窗口都具有该长度的汉明窗口。”）并且鉴于我正在进行光谱分析，因此并没有使其对称，因此使用周期性窗口是有意义的。

对于 2 号，我没有解决方案。

这是我在 python 代码中设置的函数：

M = 128
fs = 25.0
window = hamming(M, sym=False)
noverlap = np.round(M/2)
f, fxx = signal.csd(x,y,fs=fs,window=window,noverlap=noverlap

结果在 Pxy（交叉功率谱密度）方面不匹配，但在频率方面是完美的。这些是 matlab 结果中的第一个元素：

1.3590e-05
3.4354e-05
4.5282e-05
6.2549e-05
5.7965e-05
4.9697e-05
5.5413e-05

虽然这是我从 Python 中得到的：

2.04688576e-06
3.37540142e-05
4.51821900e-05
6.19997501e-05
5.78926181e-05
5.00058106e-05
5.53681683e-05

我尝试在 matlplotlib （在此处记录）中使用简单的交叉谱密度函数，如下所示：

fxx, f = mlab.csd(x,y,NFFT=M,Fs=fs,noverlap=noverlap)

而且我获得了更多的匹配结果，但仍然不完美。

1.18939608e-05
3.45206717e-05
4.56902859e-05
6.45083475e-05
5.73594952e-05
5.01539145e-05
5.34534933e-05

目标不是消除转换中可能出现的数值误差，而是使用匹配输入来操作交叉功率谱密度

有人可以帮忙吗？非常感谢提前！！！

4

0 回答 0