複雑な信号処理プラグインを開発する場合、複数のFIR/IIRフィルタを組み合わせることがよくあります。
そこで,FIR/IIRフィルタの直列/並列接続時の周波数応答をプロットするサンプルコードを作成しました。
回路全体の伝達関数を確認するだけでなく、フィルタ単体や部分的な回路の伝達関数も確認できるように、複数の伝達関数を入力として与えられるように設計しました。
また、カットオフ周波数を周波数応答とともに確認できるように垂直線プロットしています。
なお、直列回路全体の伝達関数は各フィルタの伝達関数の積,並列回路全体の伝達関数は各フィルタの伝達関数の和として計算することができます。
import sys import matplotlib.pyplot as plt import numpy as np import scipy.signal def bode_plot( w, h, sr=None, labels=None, ylim=[-25, 5], to_dB = lambda h : 20 * np.log10(np.abs(h)), xscale="linear", fc=[], ): h_list = [h] if type(h) is not list else h labels = [None] * len(h_list) if type(labels) is not list else labels fc_list = [fc] if type(fc) is not list else fc n_h = len(h_list) alpha = (1/n_h)**0.5 freqs = w/np.pi*(sr/2) if sr is not None else w plt.figure(figsize=(6.0, 8.0)) plt.subplot(2,1,1) for k in range(n_h): plt.plot(freqs, to_dB(h_list[k]), label=labels[k], alpha=alpha) plt.ylabel('Magnitude (dB)') plt.xlabel('Frequency') plt.xscale(xscale) plt.grid() plt.ylim(ylim[0], ylim[1]) for tmp_fc in fc_list: if tmp_fc < freqs[-1]: plt.vlines(tmp_fc, ylim[0], ylim[1], linestyles="dotted", colors ="k", alpha=0.5) if len(labels) > 0: plt.legend(loc='lower right') plt.subplot(2,1,2) for k in range(n_h): plt.plot(freqs, np.angle(h_list[k]) * 180 / np.pi, label=labels[k], alpha=alpha) plt.ylabel('Phase (deg)') plt.xlabel('Frequency') plt.xscale(xscale) plt.grid() plt.ylim(-180, 180) for tmp_fc in fc_list: if tmp_fc < freqs[-1]: plt.vlines(tmp_fc, -180, 180, linestyles="dotted", colors ="k", alpha=0.5) if len(labels) > 0: plt.legend(loc='lower right') plt.tight_layout() plt.show() # パラメタ sr = 44100 fc_iir = sr//10 # カットオフ周波数 fc_fir = sr//20 # カットオフ周波数 order_iir = 7 # フィルタ次数 order_fir = 9 # フィルタ次数 # フィルタ係数算出(バターワースフィルタ) b, a = scipy.signal.butter(order_iir, fc_iir/(sr/2), btype='low') w, h_iir_low = scipy.signal.freqz(b, a) # フィルタ係数算出(FIR) b = scipy.signal.firwin(order_fir, fc_fir/(sr/2), pass_zero=False) a = 1.0 w, h_fir_high = scipy.signal.freqz(b, a) # プロット bode_plot(w, [ h_iir_low, h_fir_high, h_iir_low * h_fir_high, # 直列接続 h_iir_low + h_fir_high, # 並列接続 ], labels=[ "IIR-Low", "FIR-High", "Cascade(IIR-Low, FIR-High)", "Parallel(IIR-Low, FIR-High)" ], sr=sr, fc=[fc_iir, fc_fir], ) np.set_printoptions(precision=3) import code console = code.InteractiveConsole(locals=locals()) console.interact()