PythonでFFTやスペクトログラムからバンド計算をする方法

  • このエントリーをはてなブックマークに追加

Pythonを使うと簡単にFFTやスペクトログラムを計算できますが、より深い分析を行う場合は任意の帯域(バンド)毎に統計量や技術計算をかける場合があります。ここでは、Pythonで得られたスペクトルからバンド計算をする方法を紹介します。

こんにちは。wat(@watlablog)です。ここではFFTやスペクトログラムからバンド毎の計算を行うTipsを紹介します

バンド計算の概要

FFTとスペクトログラムのおさらい

これからFFTやスペクトログラムからバンド計算を行う方法について記載しますが、そもそも「FFTとは?」「スペクトログラムって?」という方は是非以下の記事をまず読んでみて下さい。

PythonでFFT実装!SciPyのフーリエ変換まとめ

Pythonで音のSTFT計算を自作!スペクトログラム表示する方法

バンド計算とは?

バンド計算を説明する前に、まずは上記FFTから得られる情報をみてみましょう。

時間波形に対しFFTをかけると、以下のような周波数波形(ここでは振幅を例にしています)を得ます。この波形を見ると、まず目に入ってくるのは図中に示したピークの情報です。

FFTから得られる情報の例

ピークの値を定量的に求める場合は「PythonでFFT波形から任意個数のピークを自動検出する方法」で紹介した方法があります。

しかし、単純にピークの値を知りたいわけではなく、○[Hz]-○[Hz]の範囲の統計量(平均・分散・標準偏差等)を知りたい時もよくあります。

バンドの情報の例

この周波数領域の任意範囲というのがバンドBand帯域)と呼ばれるものです。統計量の他に音や振動の場合は全バンドのオーバーオール値や部分バンドのパーシャルオーバーオール値といった振幅レベルの指標を求めるというのも、製品開発の中では頻繁に行われています。

この記事ではバンド毎の計算を行う事をバンド計算と呼んでいますが、特に一般用語ではないと思われるので使用の際はご注意下さい。

ちなみに、スペクトログラムに対して任意のバンドを設定して統計量を計算すると、下図のような時間によって統計量が変化する結果を得ます。

今回はこれらのバンド計算を書いていきます!

Pythonでバンド計算をするコード

FFT波形の場合

まずは基本としてFFT波形の場合です。FFT波形は1Dの配列なので、Numpyのスライスを用いて周波数範囲を指定して計算に使用するデータ(extracted)を抽出しています。

この時、周波数[Hz]を指標(何番目のデータか)に変換する計算を入れている事に注意して下さい。

以下のコードでは平均と標準偏差を計算してコンソールに表示、さらにグラフ上に抽出した結果と計算結果を描画しています。

基本は「PythonでFFT実装!SciPyのフーリエ変換まとめ」の1フレームフーリエ変換と同じですので、是非こちらの記事もご覧下さい。

意図した部分が抽出できており、統計量を計算する事ができました。ちなみにここで算出している標準偏差は不偏標準偏差です。Numpyデフォルトだと母標準偏差のようですが、Pandasデフォルトに合わせました。

以下の図は10[Hz]から20[Hz]までのデータで計算をしています。

FFT波形からバンド計算を行った結果の例

スペクトログラムの場合

次はスペクトログラムの例です。スペクトログラムは2Dデータですが、やっている事は先ほどの1Dデータと変わらずNumpyのスライスをしているだけです。

calc_spec_band()関数がスペクトログラムからバンド計算をする部分です。

プロット部分等少しいじっていますが、基本のコードは「PythonでFFT波形から任意個数のピークを自動検出する方法」とほぼ同じなので、詳細はこちらを参考にして下さい。

以下が結果です。左のグラフは500[Hz]から1500[Hz]の平均、不偏標準偏差、ピーク値をプロットしています。

スペクトログラムならではの時間変化を確認する事ができました。

スペクトログラムからバンド計算を行った結果の例

まとめ

本記事ではFFTとスペクトログラムのおさらいを行い、バンド計算の概要を説明しました。

そして1Dと2Dのデータに対してバンド計算を行うPythonコードを紹介しました。

結果、単純にピークを検出する以上の情報をそれぞれの結果から得ることができました。

今回は短いTipsでしたが、是非Pythonによる信号処理をお楽しみ下さい。

バンド計算は簡単だけど有用な情報が沢山です!
Twitterでも関連情報をつぶやいているので、wat(@watlablog)のフォローお待ちしています!

SNSでもご購読できます。

コメント

  1. hiroko より:

    図が分かりやすくいつもありがとうございます!!バンド毎の分析で参考にさせていただいています。
    例えば、10[Hz]から20[Hz]の積分値や、全体に占める10[Hz]から20[Hz]の積分値の割合はどのように計算すればよいでしょうか。

    1. wat より:

      いつもご訪問ありがとうございます。
      例えばスペクトログラムのコードをそのまま使えば、
      extracted = fft_array[i][fl_i:fh_i]
      で10[Hz]から20[Hz]を抽出することは可能です。

      その後、means[i] = np.mean(extracted)等の代わりに、
      積分の演算を行えば可能と思います。
      単純な離散データの合算であればnp.sum()を使えば良いと思いますが、
      ある程度正確な積分を行う場合はscipy.integrateでSimpson則やRomberg則による積分を行うのも良いかも知れません。

      全体の積分は周波数バンドを抽出する前のfft_array[i]を使って計算すれば良いと思いますので、
      forループの外で計算して最後に比を取れば目的の演算ができると思います。

      積分演算は周波数波形の演算でよく使われるオーバーオールという計算につながりますので、
      そのうちこのブログでも記事にしたいですが、まずは上記方針でご検討頂ければと思います!

コメントを残す

*