比较两个频谱

Question

我试图比较两个频谱，但我对几个点感到困惑。

一个设备以40 Hz采样，另一个以100 Hz采样，因此我不确定是否需要考虑这一点。无论如何，我已经从这两个设备产生了频谱，现在我希望比较它们。我如何在每个点上做相关性，以便我在每个点上得到皮尔逊相关性。当然，我知道如何做一个整体，但我希望看到强相关性的点和不太强的点？

Answer 1

如果您正在计算功率谱密度P(f)，那么如何对原始信号x(t)进行采样并不重要。你可以直接和定量地比较这两种光谱。为了确保你已经计算了谱密度，你可以显式地检查帕赛沃斯定理：

$ \int P(f) df = \int x(t)^2 dt $

当然，你必须考虑哪些频率是实际评估的，记住，快速傅立叶变换给出的频率是f= 1/T，直到或低于奈奎斯特频率f_ny =1/(2dT)，这取决于x(t)中的样本数量是偶数还是奇数。

以下是psd的python示例代码

def psd(x,dt=1.):

"""Computes one-sided power spectral density of x.

PSD estimated via abs**2 of Fourier transform of x
Takes care of even or odd number of elements in x:
- if x is even both f=0 and Nyquist freq. appear once
- if x is odd f=0 appears once and Nyquist freq. does not appear

Note that there are no tapers applied: This may lead to leakage!

Parseval's theorem (Variance of time series equal to integral over PSD) holds and can be checked via
print ( np.var(x), sum(Px*f[1]) )
Accordingly, the etsimated PSD is independent of time series length

Author/date: M. von Papen / 16.03.2017
"""

N = np.size(x)
xf = np.fft.fft(x)
Px = abs(xf)**2./N*dt
f  = np.arange(N/2+1)/(N*dt)

if np.mod(N,2) == 0:
    Px[1:N/2] = 2.*Px[1:N/2]
else:
    Px[1:N/2+1] = 2.*Px[1:N/2+1]

# Take one-sided spectrum
Px = Px[0:N/2+1]

return Px, f