脳波データを演算ソフトのScilabでFFTしてみるが。。変

前回、OpenBCIの脳波データはいろいろ分析するため非同期（リアルタイムには行わない）でやる方針にした。で、分析ツールとしてはMatlabが有名ですが、めちゃくちゃ高いので、同等のことができるフリーソフトであるScilabを使ってみた。まだ全然使いこなしていませんが、統計とか微積とかいろいろできる、数値演算ソフトのようです。
簡単なスクリプトでCSVを読み込んだり、FFTしたりグラフ化できるそうなので、前回サンプリングした脳波をFFTしてみました。
Scilabのアプリを起動してj簡単なコマンドを実行したところ。

CSVを取り込んで、グラフ化するスクリプト
(1秒分で250サンプリングなので、データ数は250としました)



x=linspace(1,250,250)           // x軸用の配列を作る

y=csvRead("C:\home\xxx\conv_pick.txt")  // CSVを読み込む

plot(x,y)                       // x,yのデータを使ってグラフ描画

f:id:chakoku:20150523181238p:image:w250:left

１ｃｈ分の脳波データ（1秒間、250サンプル）をグラフ化したところ。フィルタ処理をまだ入れていないので、60HzのACノイズや基準電位？のドリフトが発生しています。

ひとまずこれをそのままFFTしてみた。以下が続くスクリプト



yy=fft(y)                       // サンプルデータをFFT

yy_abs=abs(yy)                  // FFTデータを絶対値化

plot(x,yy_abs)                  // FFTデータを使ってグラフ描画

f:id:chakoku:20150523181239p:image:w250:left

FFTした結果。。１の部分が飛びぬけて大きいので、自分の理解では、サンプルデータ内部で1分周（１Hz）の成分が非常に大きいということか。あるいは、、これはDC成分なのか？？　FFT不勉強でちと分からず。いずれにせよ、波形が単純増加しているので、本来の周波成分は0付近に消えてしまったということか。。。

f:id:chakoku:20150523193029p:image:w250:left

それならばということで、、極短なピークになっている周波数成分を手で０に修正して再度描画させた。すると左のようなスペクトルが得られた。中央付近（X軸のメモリで表現すると120あたり、単純に考えると120Hzと思えるが、毎秒240サンプリングなので、2.4で割る必要がある。(←本当か！？)）でピークになってるのは、５０Hz〜６０Hzに該当すると思われ（自信ない。。）これはAC電源の成分が出ていると思われる。X軸のメモリで言うと、40(17Hz)あたりから20(8Hz)にかけて上昇していると理解しました。にわかの理解なので自信ありませんが。。
※2.4で割るのは間違っている。多分。今回はサンプリングレートとサンプル数が同じなので周波数変換時に補正は不要※だったら、本当に120Hzの信号か！？あるいは、、120Hzはサンプリング240Hzで取得できる最大周波数なのでノイズなのか？？

■まとめ
脳波のサンプリングデータをFFTするには、まずフィルタでDC成分やAC電源のノイズを除去すべき。次は、、ScilabでDC成分除去やバンドパスフィルタ（ノッチフィルタ）を適用する予定。。

■補足
Scilabを理解するため、入門書を買った（ジュンク堂等でもあまり置いてないので、Amazonで半分あてずっぽうで注文）。届いたのを読んでみると、文法等の説明は1/5ぐらいで、あとは行列演算とか微積、統計での使い方であった。文法の説明がほとんどないのだけど、Scilabってあまり手続き的な使い方はしないということか。。データに対して関数による演算がメインだから？