scipy에서 주파수 분석(FFT)의 시간 변화 시각화(스펙트로그램)

scipy0.16부터의 기능으로, 취득한 신호에 FFT 걸었지만 시간 변화의 가시화를 용이하게 쓸 수 있게 되었다. 임의의 데이터로 성음 분석 같은 것이 가능.
내 Mac 환경에서 anaconda에서 파이썬을 넣고 scipy 버전이 오래되었기 때문에,
conda update scipy

하고 Scipy를 새로 만들었습니다.
시간 변화를 쫓는 것이 아닌 경우는 scipy의 fftpack을 사용한다.

파이썬 사용 환경

나는 anaconda에서 구축, python2.7,scipy0.16,Mac OSX

비정상 주파수 분석

scipy.fftpack의 함수의 FFT에서는 정상 신호의 신호의 시각화는 가능하지만, 비정상 신호의 시간 방향의 주파수 변화를 시각화하기 어렵다. scipy.signal의 spectrogram을 사용하면 FFT 한 결과의 시간 변화를 시각화 할 수 있습니다.

예를 들어, 자신의 수중 데이터로 해 본다.
여기서는 fft.py라는 스크립트와 같은 디렉토리에 data.CSV라는 1열 시간, 2열 데이터를 둔 CSV 파일을 넣었다.
데이터는 10kHz로 취득한 것.

코드

fft.py

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

plt.close('all')

input_file = u"data.CSV"
(time, data) = np.loadtxt(input_file,unpack=True, delimiter=",",usecols = (0,1))

fs = 10000.0 # サンプリング周波数
f,t,Sxx = signal.spectrogram(data, fs, nperseg=512)

plt.figure()
plt.pcolormesh(t,f,Sxx,vmax=1e-6)
plt.xlim([0,18])
plt.xlabel(u"時間 [sec]")
plt.ylabel(u"周波数 [Hz]")
# plt.colorbar()
plt.show()

signal.spectrogram의 nperseg는 주파수의 분할수 여기를 바꾸면 분할수가 바뀐다. 기본 256.
주파수의 입법에 따라서는 컬러 맵이 보기 어려워지므로, matplolib의 pcolormesh의 vmax/vmin의 값을 지정해 주면 좋다.

참고