Python: Matplotlib の cohere 関数でコヒーレンスを推定

以前 に Matplotlib の psd, csd 関数で簡単に周波数応答が求められることを紹介しました。
Matplotlib にはコヒーレンスを推定する cohere 関数もあるので、簡単にコヒーレンスも求められます。

実際に周波数応答を見る場合にはコヒーレンスも求め、信号にノイズ等の影響が生じていないかを判定します。
詳しくは 小野測器さん のページ等で説明されています。

cohere(x, y)でコヒーレンス \(C_{xy}\) を推定します。関数の内部では Welch 法で各信号の PSD, CSD を計算しています。

# coding: utf-8
import numpy as np
from scipy import signal
import matplotlib.pyplot as plt
from matplotlib import mlab

n = 1024*4
fs = 512
dt = 1/fs

b = signal.firwin(30, 0.2, window="boxcar")
x = np.random.randn(n)
y = signal.filtfilt(b, 1, x)

Pyx, f = mlab.csd(y, x, NFFT=1024, Fs=fs,
		  window=mlab.window_hanning, noverlap=1024/2)
Pxx, f = mlab.psd(x, NFFT=1024, Fs=fs,
		  window=mlab.window_hanning, noverlap=1024/2)
FRF = Pyx/Pxx

Cyx, f = mlab.cohere(y, x, NFFT=1024, Fs=fs,
		  window=mlab.window_hanning, noverlap=1024/2)

plt.figure()
# plt.plot(f, 10*np.log10(np.abs(FRF)))
plt.plot(f, Cyx)
plt.xlim(0, fs/2)
plt.xlabel("Frequency[Hz]")
plt.show()

15071001.png

コメント

Comments powered by Disqus
書籍更新情報
2016-10-21
Pythonによる科学技術計算 基礎編
PDF版の販売を開始しました。
販売ページはこちら

2016-09-09
Pythonによる科学技術計算 基礎編
1.2版への更新が可能になりました。
サポートページはこちら
電子書籍
Pythonによる科学技術計算 基礎編
Kindle ストア、Leanpubで販売中です
Pythonによる科学技術計算 基礎編
PDF版の販売はこちら
同人誌
技術書典(2016/6/25)
Emacs/org-modeのPDF作成術
電子版をBOOTHで販売中です
Emacs/org-modeのPDF作成術
Share