frequency - 比较两个频谱

Question

我正在尝试比较两个频谱，但我对很多点感到困惑。

一个设备以 40 Hz 采样，另一个以 100 Hz 采样，所以我不确定是否需要考虑到这一点。无论如何，我已经从这两种设备中产生了频谱，现在我想比较这些。如何在每个点上进行相关性，以便在每个点上获得皮尔逊相关性。我当然知道如何做一个整体的，但我想看到强相关性的点和不那么强的点？

Answer 1

如果您正在计算功率谱密度 P(f)，那么原始信号 x(t) 的采样方式并不重要。您可以直接定量地比较两个光谱。为确保您已计算出光谱密度，您可以明确检查 Parsevals 定理：

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

当然，您必须考虑实际评估哪些频率 请记住 fft 为您提供从 f = 1/T 直到或低于奈奎斯特频率 f_ny = 1/(2 dt) 的频率，具体取决于 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