<p>我不确定,相位变量是在哪里计算的,答案是@Mattijn。</p>
<p>你可以从实数和
交叉谱密度的虚部。</p>
<pre><code>from matplotlib import mlab
# First create power sectral densities for normalization
(ps1, f) = mlab.psd(s1, Fs=1./dt, scale_by_freq=False)
(ps2, f) = mlab.psd(s2, Fs=1./dt, scale_by_freq=False)
plt.plot(f, ps1)
plt.plot(f, ps2)
# Then calculate cross spectral density
(csd, f) = mlab.csd(s1, s2, NFFT=256, Fs=1./dt,sides='default', scale_by_freq=False)
fig = plt.figure()
ax1 = fig.add_subplot(1, 2, 1)
# Normalize cross spectral absolute values by auto power spectral density
ax1.plot(f, np.absolute(csd)**2 / (ps1 * ps2))
ax2 = fig.add_subplot(1, 2, 2)
angle = np.angle(csd, deg=True)
angle[angle<-90] += 360
ax2.plot(f, angle)
# zoom in on frequency with maximum coherence
ax1.set_xlim(9, 11)
ax1.set_ylim(0, 1e-0)
ax1.set_title("Cross spectral density: Coherence")
ax2.set_xlim(9, 11)
ax2.set_ylim(0, 90)
ax2.set_title("Cross spectral density: Phase angle")
plt.show()
fig = plt.figure()
ax = plt.subplot(111)
ax.plot(f, np.real(csd), label='real')
ax.plot(f, np.imag(csd), label='imag')
ax.legend()
plt.show()
</code></pre>
<p>要关联的两个信号的功率谱密度:
<img src="https://i.stack.imgur.com/ozNfl.png" alt="The power spectral density of the two signals to be correlated"/></p>
<p>两个信号的相干性和相位(放大到10赫兹):
<img src="https://i.stack.imgur.com/cVmsf.png" alt="The coherence and the phase of the two signals (zoomed in to 10 Hz)"/></p>
<p>这里是真实的和想象的(!)部分交叉光谱密度:
<img src="https://i.stack.imgur.com/kY1aU.png" alt="real and imaginary part of the cross spectral density"/></p>