DVIOUT

Transcription

1 第 章 離散フーリエ変換 離散フーリエ変換 これまで 私たちは連続関数に対するフーリエ変換およびフーリエ積分 ( 逆フーリエ変換 ) について学んできました この節では フーリエ変換を離散化した離散フーリエ変換について学びましょう 自然現象 ( 音声 ) などを観測して得られる波 ( 信号値 ; 観測値 ) は 通常 電気信号による連続的な波として観測機器から出力されます しかしながら コンピュータはこの様な連続的な波を直接扱うことができないため AD 変換器 を用いてこの波を一定間隔でサンプリングし コンピュータが直接扱うことのできる離散的な波 ( 数値列 ) に変換して扱います ( 図 を参照のこと ) 前者をアナログ (analog) と呼び 後者をデジタル (didital) と呼びます 観測機器 アナログ信号 AD 変換器 デジタル信号 コンピュータ 図 : アナログからデジタルへ アナログ信号からデジタル信号に変換する装置を A-D 変換器と呼び 逆に デジタル信号からアナログ信号に変換する装置を D-A 変換器と呼びます

2 96 第 章離散フーリエ変換 この離散的な波 ( 数値列 ) から元の連続的な波の性質を調べるために これまで学んだフーリエ変換を離散フーリエ変換に書き換えましょう いま 図 のように観測開始 [ 秒 ] から観測終了 T [ 秒 ] までに観測された連続的な波 x(t) に対して t [ 秒 ] の一定間隔でサンプリングされた 個の数値列 ( 離散的な波 ) を としましょう { x(), x( t), x( t), x( t),, x(( ) t), x(( ) t) } x(t) t t t ( ) t ( ) t T 図 : 離散化された波この離散的な波 ( 数値列 ) を直接フーリエ変換することはできないので 区分求積法 ( リーマン和 ) を使って 近似されたフーリエ変換 (f) を求めます 指数関数 e iπft を ϕ(t) =e iπft とおくと 近似されたフーリエ変換 (f) = Z x(t)e iπft dt = Z T x(t)e iπft dt x()ϕ() t + x( t)ϕ( t) t + x( t)ϕ( t) t + +x(( ) t)ϕ(( ) t) t ( 区分求積法 ) が得られます さらに 関係式 ( t は区間 [,T ] を 等分したもの ) t = T = T [ 秒 ] を代入すると 近似されたフーリエ変換 (f) T n= x(n t)ϕ(n t)

3 離散フーリエ変換 97 が得られます ここで 連続的な波と離散的な波 ( 離散化された波 ) の違いについて比較しておきましょう そもそも フーリエ変換は 連続的な波を三角関数の積み重ねによって表そうとしたものですから いま 連続的な波が 図 のように三角関数の総和で表されている 例えば これを 等分 ( =) して離散的な波 ( 離散化された波 ) として捉え直してみましょう µ πk a k cos t T µ πk b k sin t T k = k = k = k = 図 : 離散化された波のフーリエ変換

4 9 第 章離散フーリエ変換 図 のグラフを注意深く見てみると k =のグラフは 図 のように k =のグラフとして解釈することができます ( 本来 赤線の波を青線の波として解釈することができる ) すなわち 離散的な波 ( 離散化された波 ) では k =,, の波は k =,, の波に吸収され k =,, の三角関数のみで表せることを意味しています 逆に 離散的な波 ( 離散化された波 ) では k =,, の三角関数でしか表せないことになります ( に依存します ) µ πk a k cos t T µ πk b k sin t T k =(k =) k =(k =) k =5(k =) k =6(k =) 図 : 離散化によって吸収される波

5 離散フーリエ変換 99 したがって 離散的な波 ( 離散化された波 ) に含まれる波の最大周波数 f max (= f s /) [Hz] が決定し t [ 秒 ] で [ 回転 ] していることから ( 例では k =の場合 ) 最大の周波数 f max は f max = µf t [Hz] s = t [Hz] となります さらに 離散的な波 ( 離散化された波 ) の基準となる周波数 f [Hz] を f = T [Hz] とおき ( 節参照 ) 式と 式を 式に代入すると 関係式 f max = t = = f [Hz] T (f s = f [Hz]) が得られます 以上より 離散的な波 ( 離散化された波 ) の特徴を捉えるには 周波数, ± f, ± f, ± f,, ± f [Hz] について調べれば十分であることがわかります ( 第 章で述べたように 正の周波数と負の周波数が同時に存在します ) なお 基準となる周波数 f [Hz] を周波数分解能といい f より細かな周波数に分解することはできません また 周波数 f s (= f max )[Hz] をサンプリング周波数といい f s /(=f max )[Hz] 以上の周波数の波は現れません ( 周波数領域の f s [Hz] は 時間領域の T [ 秒 ] に対応しています ) 離散的な波のフーリエ変換に話を戻しましょう 前記の考察から 調べる必要のある周波数 k f [Hz] (k =, ±, ±,, ±/) の近似されたフーリエ変換 (k f) は (k f) = T となります 式と 式より 関係式 = T n= n= x(n t)ϕ(n t) x(n t)e iπ(k f)(n t) f t = が成り立つことに注意すれば 近似されたフーリエ変換 (k f) は 改めて (k f) = T n= π i x(n t)e と書き直すことができます ここで W = e i π µ =cos µ π i sin µ π

6 第 章離散フーリエ変換 とおくと 点 W, W, W,,W は 図 5 のように単位円を 等分した円周上の点となり 関係式 W = を満たします (W は の 乗根となります ) さらに 関係式 W m = W W m = W m m = W が成り立つことから W = W, W = W, W = W, となります ( =) 5 6 (= ) 7 (= ) (= =) 図 5: 等分された単位円周上の点 ( = の場合 ) そこで 負の周波数 k f [Hz] (k =,,,/) の近似されたフーリエ変換 ( k f) を求めると 関係式 ( k f) = T = T = T = T n= n= n= n= = (( k ) f) x(n t)e i πn( k ) = T n= ³ k x(n t) W n T = n= x(n t)w n( k ) ³ x(n t) W k W n T = x(n t)w n( k ) = T n= k x(n t) ³ W n n= k x(n t) ³W n x(n t)e i πn( k )

7 離散フーリエ変換 が得られ k = k とおくと 負の周波数 k f [Hz] (k =,,, /) の近似されたフーリエ変換の値 (k f) は ( f) = (( ) f), ( f) = (( ) f), ( f) = (( ) f), ( (/ ) f) = ((/+) f), ( (/ ) f) = ((/+) f), ( (/) f) = ((/) f) となります したがって k =,,,, の近似されたフーリエ変換 (k f) は ( f), ( f), ( f), ((/ ) f), ((/) f) = ( (/) f), ((/+) f) = ( (/ ) f), ((/+) f) = ( (/ ) f), (( ) f) = ( f), (( ) f) = ( f), (( ) f) = ( f) となります ただし 赤文字は 負の周波数の近似されたフーリエ変換に一致することに注意しましょう また 正の周波数 (/) f の近似されたフーリエ変換は 負の周波数 (/) f の近似されたフーリエ変換として扱うこともできますが 一般的な書籍にならって本テキストでも正の周波数の近似されたフーリエ変換 ((/) f) として扱うことにします 以上をまとめると 離散化された ( 近似された ) フーリエ変換は (k f) = T n= x(n t)w となります さらに 簡潔に記述するために 記号 k = (k f), x n = x(n t) (k =,,,, ) 5 を導入すると 離散化された ( 近似された ) フーリエ変換は k = T n= x n W (k =,,,, )

8 第 章離散フーリエ変換 となります この式を離散フーリエ変換 (Discrete Fourier Transform; DFT) と呼びます ( 信号処理用に最適化されています ) なお 離散フーリエ変換の離散周波数スペクトル密度 T (f) は T (f) T (k f) = k T = n= x n W (k =,,,, ) となります 理解を深めるために 例として =の場合を見てみましょう 具体的に離散フーリエ変換を書き下してみると = T x W + x W + x W + x W, = T x W + x W + x W + x W, = T x W + x W + x 6 + x W, = T x W + x W + x W x W = となり 関係式 =(=W ) を用いると = T x W + x W + x W + x W, = T x W + x W + x W + x W, = T x W + x W + x W + x W, = T x W + x W + x W + x W = となります また 行列の形式 = T で表現することもできます さらに x x x x =, = i, =, = i を代入すれば 具体的に離散フーリエ変換 = T i i i i x x x x

9 離散フーリエ変換 を求めることができます もう つの例として =の場合を見てみましょう 具体的に離散フーリエ変換を書き下してみると = T x W + x W + x W + x W + x W + x 5 W + x 6 W + x 7 W, = T x W + x W + x W + x W + x + x 5 W 5 + x 6 W x 7 W, = T x W + x W + x + x W 6 + x + x 5 W + x 6 W + x 7 W, = T x W + x W + x W 6 + x W 9 + x W + x 5 W 5 + x 6 W + x 7 W, = T x W + x + x + x W + x W + x 5 W + x 6 W + x 7 W, 5 = T x W + x W 5 + x W + x W 5 + x W + x 5 W 5 + x 6 W 5 + x 7 W =, 6 = T x W + x W 6 + x W + x W + x W + x 5 W + x 6 W 6 + x 7 W =, 7 = T x W + x W 7 + x W + x W + x + x 5 W 5 + x x 7 W = となり 関係式 =(=W ) を用いると = T x W + x W + x W + x W + x W + x 5 W + x 6 W + x 7 W, = T x W + x W + x W + x W + x + x 5 W 5 + x 6 W x 7 W, = T x W + x W + x + x W 6 + x W + x 5 W + x x 7 W, = T x W + x W + x W 6 + x W + x + x 5 W 7 + x 6 W 5 + x 7 W, = T x W + x + x W + x + x W + x 5 + x 6 W + x 7 W, 5 = T x W + x W 5 + x W + x W 7 + x + x 5 W + x 6 W 6 + x 7 W =, 6 = T x W + x W 6 + x + x W + x W + x 5 W 6 + x 6 + x 7 W =, 7 = T x W + x W 7 + x W 6 + x W 5 + x + x 5 W + x 6 W + x 7 W = となります また =の場合と同様に 行列の形式 6 W = T で表現することもできます さらに W =, W = i, = i, W = i, =, W 5 = + i, 6 = i, W 7 = + i x x x x x x 5 x 6 x 7

10 第 章離散フーリエ変換 を代入すれば 具体的に離散フーリエ変換 i i i + i i + i = T i i i i i i i + i i + i i i + i i i i i i i i 7 + i i + i i i i を求めることができます 例として 連続的な波 x(t) =t ( 区間 [, ]; T =[ 秒 ]) を一定間隔 ( =) でサンプリングして得られた離散化された波 ( 離散的な波 ; 数値列 ; データ ; リスト ) ½ ¾,,,,, 5, 6, 7 の離散フーリエ変換について考察しておきましょう ( 図 6 参照 ) x x x x x x 5 x 6 x 7 図 6: 離散化された波の例 (x(t) =t) 行列の形式 i i i + i i + i = i i i i i i i + i i + i i i + i i i i i i i i 7 + i i + i i i i より 離散フーリエ変換 k (k =,,, 7) を計算すると = 7, = +i( + ), = +i, = +i( + ), =, 5 = = +i( ), 6 = = i, 7 = = +i( ) / / / / / 5/ 6/ 7/

11 となります 波の特徴を調べるために離散周波数スペクトル密度 離散フーリエ変換 5 T (k f) = k T = k = k (k =,,,, 7) を求め 表にまとめると表 のようになります なお T =より 周波数分解能 ( 基準となる周波数 ) f が f = =[Hz] になることと 青文字 (,,, ) が正の周波数の離 T 散周波数スペクトル密度 および 赤文字 ( 5 =, 6 =, 7 = ) が負の周波数の離散周波数スペクトル密度となることに注意しましょう k 周波数 k f [Hz] Re T (k f) Im T (k f) 表 : 離散周波数スペクトル密度 ここで 連続的な波の周波数スペクトル密度 ( フーリエ変換 ) T (f) = f T = (f) = (f) = Z x(t) e iπft dt = Z = +e i(πf) + i (πf) e i(πf) (πf) t e iπft dt と離散化された波の離散周波数スペクトル密度 T (k f) を比較するために両方のグラフを重ねて描くと 図 7 のようになります 周波数スペクトル密度 T (f) のグラフに対して離散周波数スペクトル密度 T (k f) の各点がずれている理由は この節の前半で述べたように 離散化することで k =, ±5, ±6, ±7, の周波数 k f の離散周波数スペクトル密度 T (k f) が k =, ±, ±, ±, + の周波数 k f の離散周波数スペクトル密度 T (k f) に吸収されているためです したがって このような誤差を少なくするためにはサンプリング数 を大きくとる必要があります 図 の =の場合 =の場合に比べて k =, ±9, ±, ±, の周波数 k f の離散周波数スペクトル密度 T (k f) が k =, ±, ±, ±, ±, ±5, ±6, ±7, + の周波数 k f の離散周波数スペクトル密度 T (k f) に分散されて吸収されるため誤差が少なくなっていることが確認できます このように 波を離散フーリエ変換するには十分なサンプリング数が必要であることがわかります ただし サンプリング数を増やすとそれだけ多くの計算が必要になることにも注意しなければなりません 実際に離散フーリエ変換で解析や加工を行なう場合は 計算量に注意しながら 対象となる波の特徴を十分に捉えられるだけの最低限のサンプリング数を確保する必要があります

12 第 章離散フーリエ変換 Re T (f) と Re T (k f) Im T (f) と Im T (k f) 図 7: 周波数スペクトル密度の比較 ( =) Re T (f) と Re T (k f) Im T (f) と Im T (k f) 図 : 周波数スペクトル密度の比較 ( =) さらに 得られた離散周波数スペクトル密度 T (k f) から元の連続的な波 x(t) の性質を調べ 元の連続的な波 x(t) を三角関数で近似した波 ( 関数 ) を求めると =の場合 x(t) a + / k= µ µ µ πk πk a k cos t + b k sin t T T = a + ak cos(π(k f)t)+b k sin(π(k f)t) k= = a + a cos πt + a cos πt + a cos 6πt + a cos πt +b sin πt + b sin πt + b sin 6πt + b sin πt = 7 cos πt cos πt cos 6πt cos πt + sin πt sin πt sin 6πt + sin πt

13 離散フーリエ変換 7 となります ( 第 節の 波のフーリエ変換 を参照 ) なお 各三角関数のスペクトル ( 振幅 ) は a k =Re T (k f)+re T ( k f), (k =,,,, ) b k = Im T (k f)+im T ( k f) を計算することで得られます ( ただし T ( / f) =で計算 ) このように 三角関数で近似された波を求めることで 区間 [,T ] における元の連続的な波の性質をとらえることができます 特に 元の連続的な波が三角関数で構成されている場合 元の連続的な波そのものを知ることができます ( 例題, 例題, 例題 ) この例の締めくくりとして 前述で得た三角関数で近似された波と元の連続的な波のグラフを描いて比較しておきます ( 図 9) = = 図 9: 三角関数で近似された波 必ず Im T () = および Im T (/ f) =となります また 離散化されたスペクトル ( 数値列 ) に現れない離散周波数スペクトル密度の値は となります なぜそうなるか それぞれ理由を考えてみてください a =Re T ( f)+re T ( f) =Re T () + Re T () = Re T (), a / =Re T (/ f)+re T ( / f) =Re T (/ f)+=re T (/ f), b / = Im T (/ f)+im T ( / f) =+= 三角関数で近似された波は 全区間 [, ] に拡張され 区間 [,T ] を積分区間とする周期 T を持つフーリエ級数 ( 周期関数 ) となります

14 第 章離散フーリエ変換 例題 ある波 x(t) を [ 秒 ]( 区間 [, ]) にわたって観測し 一定間隔でサンプリングしたところ 離散化された波 ( 離散的な波 ; 数値列 ; データ ; リスト ) {,,, } を得た 離散フーリエ変換 ( 離散周波数スペクトル密度 ) を使って ある波 x(t) を求めなさい 解答例 T =より 周波数分解能 f = = =5 [Hz] である さらに =よ T り 離散周波数スペクトル密度を求めると Ã! T (k f) = k = T x n W = x n (k =,,, ) T T T = = n= n= i i i i x x x x = 6 6 = となる = となることに注意すれば 各周波数における離散周波数スペクトル密度は以下の通りである k 周波数 k f [Hz] Re T (k f) 5 5 Im T (k f) 表より ある波は x(t) は 角速度 ω =πf = ± π [rad/ 秒 ]( 周波数 f = ±5 [Hz]) の cos 波形の波であることがわかる なお cos 波形の波の振幅は となる 以上より ある波 x(t) は Re T (5) + Re T ( 5) = 5+5 = x(t) =cos π t 5 5 である * (( k) f) =(( k) f) より = ( f) =(( ) f) =(( ) f) = * 表が作成できれば 節と同様に 波の特徴を調べることができます * 答えが正しいか グラフを描いて確認しましょう

15 離散フーリエ変換 9 例題 ある波 x(t) を [ 秒 ]( 区間 [, ]) にわたって観測し 一定間隔でサンプリングしたところ 離散化された波 ( 離散的な波 ; 数値列 ; データ ; リスト ) {,,, } を得た 離散フーリエ変換 ( 離散周波数スペクトル密度 ) を使って ある波 x(t) を求めなさい 解答例 T =より 周波数分解能 f = = =[Hz] である さらに =より T 離散周波数スペクトル密度を求めると Ã! T (k f) = k = T x n W = x n (k =,,, ) T T T = = n= i i i i n= x x x x = i +i = 5 i5 5+i5 となる = となることに注意すれば 各周波数における離散周波数スペクトル密度は以下の通りである k 周波数 k f [Hz] Re T (k f) 5 5 Im T (k f) 5 5 表より ある波は x(t) は 角速度 ω =πf = ±π [rad/ 秒 ]( 周波数 f = ± [Hz]) の cos 波形の波および sin 波形の波を合わせた波であることがわかる なお cos 波形の波および sin 波形の波の振幅は Re T () + Re T ( ) = 5+5 =, となる 以上より ある波 x(t) は である Im T () + Im T ( ) = ( 5) + 5 = x(t) =cosπt +sinπt

16 第 章離散フーリエ変換 例題 ある波 x(t) を [ 秒 ]( 区間 [, ]) にわたって観測し 一定間隔でサンプリングしたところ 離散化された波 ( 離散的な波 ; 数値列 ; データ ; リスト ) {,,, } を得た 離散フーリエ変換 ( 離散周波数スペクトル密度 ) を使って ある波 x(t) を求めなさい 解答例 T =より 周波数分解能 f = = =5 [Hz] である さらに =よ T り 離散周波数スペクトル密度を求めると Ã! T (k f) = k = T x n W = x n (k =,,, ) T T T = = n= n= i i i i x x x x = = となる = となることに注意すれば 各周波数における離散周波数スペクトル密度は以下の通りである k 周波数 k f [Hz] 5 5 Re T (k f) Im T (k f) 表より ある波は x(t) は 角速度 ω =πf =π [rad/ 秒 ]( 周波数 f =[Hz]) の cos 波形の波 および定数項 ( 直流 ) を合わせた波であることがわかる なお cos 波形の波の振幅は となる 以上より ある波 x(t) は である Re T () + Re T ( )=+= x(t) =+cosπt * k = / の場合 負の周波数 (/) f [Hz] の波は なので 正の周波数 (/) f [Hz] の波のみで扱います ちなみに 離散化された波では 表以外の波は全て です * 定数項はフーリエ級数の a に対応します