Scipy的welch函数和MATLAB的pwelch函数结果不一致

7

我在使用Python中的scipy.signal方法welch(https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.signal.welch.html)时遇到了困难。该方法用于估计时间信号的频谱,但是与MATLAB的pwelch方法相比,即使使用相同的参数(窗口大小、重叠等),输出结果也完全不同。以下是我在每种语言中的代码,并且输入文件和输出文件在此链接中:

https://www.dropbox.com/s/2ch36phbbmjfhqg/inputs_outputs.zip?dl=0

输入是一个二维数组,其中行代表时间步长,每列是一个信号片段。输出中的列是从输入中相应列的频谱。

Python:

import numpy as np
from scipy.signal import welch, get_window
input = np.genfromtxt('python_input.csv', delimiter=',')
fs = 128

window = get_window('hamming', fs*1)
ff,yy = welch(input, fs=fs, window = window, noverlap = fs/2, nfft=fs*2, 
axis=0, scaling="density", detrend=False)
np.savetxt("python_spectrum.csv", 10*np.log10(yy), delimiter=",")

MATLAB:

input       = csvread('matlab_input.csv');
fs          = 128
win         = hamming(fs);
[pxx,f]     = pwelch(input ,win,[],[],fs,'psd');
csvwrite('matlab_spectrum.csv',pxx);

我怀疑scipy是问题所在,因为它的输出不符合我之前使用的滤波器(4阶butterworth带通滤波器,0.3到35 Hz,使用filtfilt)的反映 - 然而,MATLAB的输出却可以:

使用MATLAB中的imagesc显示每个方法的输出

这里还有一些元素差异的图表:

元素差异(y轴应该为0-64!)(第三个图表中我排除了最极端的值)

我尝试了一个简单的正弦波,在两种编程语言中都运行良好 - 所以我完全迷失了。
1个回答

4

你的Python和MATLAB文件排序方式不同,这导致了错误。虽然数组的形状相同,但在Python中,值是按行排序的,在MATLAB中则是按列排序。

你可以通过重新调整输入数组并进行转置来解决此问题。这将使得值的排序方式与MATLAB相同:

input = input.reshape((input.shape[1], input.shape[0])).T

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接