(1)
本日から，情報処理応用の一つとして「信号処理」を扱う． 信号処理とは，実験などで得られた何らかの波形データを（数学に基づいた）方法で加工し，その特徴を抽出したり，雑音を軽減したり，必要な情報を抽出したりする手法の総称である． 広い意味では，2次元データである画像などの解析（画像処理と呼ばれ，顔認識やモザイク処理，色調変換や画像圧縮など）も信号処理と呼ばれるが，ここでは簡単に1次元データの処理を考えよう．
ある信号波形データがあったとして，この波形から何か情報を抽出したいとする．
たとえば，音声信号から内容を推定したり，楽器や自動車の振動波形からその減衰特性や共振特性を把握したい，各種気象データから長期の気候変動を知りたい場合などがこれにあたる．
実際にはこれらの波形は単純な正弦波などには程遠く，複雑な信号波形であるため，そのまま観察していても得られる情報は少なく，よく分からない場合が多い．
そこで，その信号にどのような周波数（振動数※）成分がどのくらいの大きさで含まれているかを解析することによって，どのような振動成分が支配的であるかを調べると，その性質を把握することが出来る．
そのための解析手法としてフーリエ変換が用いられる．
まず，フーリエ変換（Fourier Transform)がどのような変換で，どのような性質かを理解しよう．
フーリエ変換とは，任意の関数を三角関数の和（積分）で表す変換であり，もともとはフランスの数学者 Joseph Fourier により，固体の熱伝導方程式を解析するために提案された方法である．
連続な関数 f (t) に対して，フーリエ変換は以下のように定義される．
(1)
ω は角周波数であり，周波数 f とは ω = 2πf なる関係がある．
（指数部が純虚数の指数関数は周期関数となる）
さて，計測器などでコンピュータに信号を取り組んだ場合（サンプリングと呼ぶ），必ず離散的，かつ，有限区間のデータの集合となる．
そこで，フーリエ変換の式を離散データに対して計算できるようにすると，以下のようになる．
ここで，下のデータをC言語の配列風に f[t] とする．
(2)
と表される． N はデータの個数である． これが離散フーリエ変換の基本形である．
右辺の周期関数の部分は，ω が 2π，4π，6π と増加するにつれ，周波数が増加する．
即ち，離散フーリエ変換は，いろいろな周波数の正弦波と，元の波形 f との積和（積の和，つまり内積）を計算することと等価である．
このように離散的なデータに対してのフーリエ変換を離散フーリエ変換(Discrete Fourier Transform, DFT)と呼び，コンピュータで計算できるのはこのDFTである．
| ※周波数と振動数 どちらも単位はHzで似たようなものと考えられそうだが，以下のように区別される． 周波数は，電波，磁界，音などの信号を対象とした表現に用いられる． また，「共振周波数」などの言葉の通り，積極的に発振させる意図を持った利用がなされる． 振動数は，機械の構造物や機構の応答に用いられる表現で利用される． また，原子や分子の運動に関しても振動数が用いられる． 「固有振動数」という言葉の通り，対象とする物理系の持つ特性を考える場合に用いられる． |
ところが，DFTはデータ数が大きくなると著しく計算負荷が高く（計算時間がN2に比例）なるため，小さなデータ数しか計算できない問題があった．
1965，James CooleyとJohn Tukeyにより，DFTの計算を劇的に高速化するアルゴリズムが開発された．
これが高速フーリエ変換(Fast Fourier Transform, FFT)である．
高速フーリエ変換は，データ数Nが小さな素数の積（最も良いのは2の累乗）の場合，Nが大きくなっても計算時間がさほど増加せず（NlogNに比例），画像データなどの巨大なデータも一瞬で変換できるため，実用性が高まった．
フーリエ変換では，変換結果は一般的に複素数になる．
変換した関数の絶対値の2乗（＝パワー）を縦軸に，横軸に周波数（または振動数）をとりグラフに表したものを，「パワースペクトル」と言う．
また，縦軸を振幅（の絶対値）としたものを，「振幅スペクトル」と呼び，いずれも周波数解析の重要な方法となっている．
以下に，パワースペクトル例を示す．
スペクトルの一例．横軸が周波数，縦軸が振幅（またはその二乗であるパワー）を表す．
どの周波数の振動が，どのくらいの割合で含まれるかを表している．
(出典)http://www.mathworks.co.jp/jp/help/matlab/ref/fft.html
リンク先は「MATLAB」という数学・技術演算を行うための商用ソフトウェアのリファレンスであり，信号波形のデータから上図のようなグラフを生成することができる．
情報処理教室にはMATLABがインストールされている端末があるので，試してみよう．
そのほか，Windowsメディアプレイヤーのソフトなどにもこれを利用した視覚効果がある．
フーリエ変換以前に，もともと連続的な信号を計測機器などにより一定時間間隔で取得し，離散的なデータ列に変換することを「サンプリング」と呼び，日本語では「標本化」と言う．
一般的には，注目する量の変化よりも速い速度＝細かい時間刻みで，その時の値を順次計測する．
この時間刻みのことをサンプリング間隔と呼び，単位は時間(sec, min等)となる．
たとえば，「一日の気温の変化」をサンプリングするには1時間おきの計測で十分だが，「エアコンの出力制御」では1分，あるいはもう少し短い時間間隔で室内の温度をサンプリングする必要がある．
サンプリング間隔の逆数を「サンプリング周波数」と呼び，単位は Hz である．
計測データは一定のサンプリング時間間隔おきに配列に格納されていく． 時間軸上のデータは基本的に一次元配列データであるが，画像データなどは2次元配列データとなる．
多くの物理変数は連続であるが，計算機で扱うために「離散化」して，コンピュータに格納できる形式（一般的には配列）にする．
どのような時間間隔（サンプリング周波数）でサンプリングしたかが分からないと，後で周波数の算出ができなくなる．
サンプリングする時間の間隔を「サンプリング周期」と呼び Δt とすれば，データ配列の 0 番目には f(0)，1番目には f(Δt)，2番目にはf(2Δt)・・・，といった計測データが順次入る．
通常，サンプリング周期は一定とする．
ある関数を別の関数に何らかの方法で変形することを「変換」と呼ぶ．
キーワード：フーリエ変換，離散フーリエ変換，高速フーリエ変換の違い．
今週の実習では，フーリエ変換のプログラミングの準備として，波形生成などいくつかの関数を用意する．
FFTによる解析のため，以下に挙げる関数を準備せよ．
ヘッダファイルおよび実装ファイル名はそれぞれ
153R000000-07-dft.h ,
153R000000-07-dft.cpp
とする．
main関数を含むファイルは153R000000-07-1.cpp とする.
（000000の部分は学生番号）
各関数は，main関数から呼び出し，正しく動作するか確認すること．
後半のダミーの波形データ生成関数については，得られた波形の形状をグラフにして確認すること．
グラフの提出は不要．
（画面に出力したものExcelに貼り付けるか，ファイルにcsv形式で出力しExcelで描画．）
以下，特に指定が無い限り
x は信号データの格納された実数型配列，
n は配列の要素数（データ個数）
である．
どんな n が渡されてきても，うまく動く関数を作ろう．
n=100程度で動作チェックすること．
void clear_array(double x[], int n); // 要素数 n の配列 x に，すべて 0 を代入する（初期化）．void square_wave(double x[], int n); // 要素数 n の配列 x に，任意の周波数（例えば 5 ）の矩形波を書き込む．振幅は任意，例えば -1 から 1 など．void sin_waves(double x[], int n, double f1, double f2, double f3);
#include <math.h>
const double pi = 2.0*asin(1.0);
for(i=0; i<n; i++) {
x[i] = ( 2.0 * sin( 2*pi*i/n * f1 )
+ 0.5 * cos( 2*pi*i/n * f2 )
+ 0.25 * sin( 2*pi*i/n * f3 );
}
| square_wave関数，矩形波 | sin_waves関数，正弦波の和 |
 n=50, 周波数 5 の矩形波の例 |
 正弦波，余弦波の和の例 |
離散フーリエ変換(DFT)と高速フーリエ変換(FFT)アルゴリズムについて文献等を調査し，簡単に説明せよ．
用いた文献名を明記すること．
テキストファイル(.txt)またはWordファイル(.docx)にて提出．
ファイル名は153R000000-07-2.txt
または
153R000000-07-2.docx とする．
課題1で生成した波形データにたいして，本ページ内の離散フーリエ変換の式(2)，及び先週までに作成した複素数計算関数を用いて，F[0], F[2π], F[4π], ...などの値を計算してみよ．
計算結果はテキストファイル 153R000000-07-3.txt に貼り付けて提出せよ．