音声を可視化したい(2) - 理想のユーザ・インターフェイスを求めて

2011-11-01 - 理想のユーザ・インターフェイスを求めての続き。

波のデータをスペクトル解析したい。どの周波数成分がドミナントになっているかとか。まずは、フーリエ変換のライブラリを試す。(参考サイト：高速フーリエ変換（FFT） - 人工知能に関する断創録)

テストとして、周期1秒のサインカーブをフーリエ変換する。使用するのは、numpy.fftライブラリ。

Discrete Fourier Transform (numpy.fft) — NumPy v1.15 Manual

numpy.fft.fftだと、離散データ $a_{m}$ のフーリエ係数 $A_{k}$ は、
$A_k = \sum_{m=0}^{n-1}a_m \exp(-2\pi imk/n)$ where $k=0,...,n-1$
と表せる。一方で、逆変換は、
$a_m=\frac{1}{n}\sum_{k=0}^{n-1}A_k \exp(2\pi imk/n)$ where $m=0,...,n-1$
となる。

$y=\sin(2\pi t)$ をフーリエ変換し、周波数ごとの振幅をプロットするのが下のコード。

import numpy as np
import matplotlib.pyplot as plt

lmax = 2.0
dsc = 0.01

t = np.arange(0, lmax, dsc)　#Time=2秒までのデータを200分割
y = np.sin(2*np.pi*t)

f = np.fft.fftfreq(t.size, dsc)　#t.size=200
s = np.fft.fft(y)

am = np.sqrt(s.real**2 + s.imag**2) #フーリエ係数の振幅を求める

plt.subplot(211)　#元の波形をプロット
plt.plot(t, y, 'ro', markersize=4)
plt.xlabel('Time [sec]')

plt.subplot(212)　#周波数ごとの振幅
plt.plot(f, am, 'bo', markersize=6)
plt.xlabel('Frequency [Hz]')
plt.axis([0, 10, 0, 150]) #正の周波数成分のみプロットする

plt.show()

実行結果。