12.6
フーリエ変換と信号処理
音声信号などを周波数スペクトルに変換して、そのスペクトルを操作することに よってより複雑な信号処理が可能となる。ここでは、離散信号を周波数変換する離散 フーリエ変換と、それを高速に計算する FFT(Fast Fourier Transform) について紹介 する。 12.6.1 フーリエ級数 • 対象とする信号は、周期 T の信号 g(t) g(t) = 1 2a0+ ∞ X m=1 (amcos(2π T mt) + bmsin( 2π T mt)) (1) ( am= T2 RT 0 g(t) cos(2πT mt)dt m = 0, 1, 2, . . . bm= T2 RT 0 g(t) sin(2πT mt)dt m = 1, 2, 3 . . . (2) • g(t) は、周期 T の周期関数なので、積分区間は 0 < t < T ではなく、t0< t < t0+ T でよい。 • am, bm: 周波数m T の周波数成分を表す。周波数が 0, 1 T, 2 T, . . . の飛び飛びの 値のところに成分をもち、線 スペクトルと言われる。 12.6.2 複素フーリエ級数 • フーリエ級数を複素数 ei2πTmtを用いて指数形式に表現にしたもの。 • オイラーの式より、 cos(2π T mt) = 1 2 ³ ei2πTmt+ e−i2πTmt ´ (3) sin(2π T mt) = 1 2i ³ ei2πTmt− e−i2πTmt ´ (4) これらを式 (1) のPの中の項に代入すると amcos(2π T mt) + bmsin( 2π T mt) = am× (3) + bm× (4) = 1 2(am− ibm) e i2π T mt+1 2(am+ ibm) e −i2π T mt = Cmei 2π T mt+ Bme−i2πTmt ただし、 Cm= 12(am− ibm) Bm= 12(am+ ibm) また、 Cm= ¯Bm (複素共役) である。よって、(1) は、次のように書きなおされる。 g(t) = a0 2 + ∞ X m=1 ³ Cmei 2π Tmt+ Bme−i2πTmt ´ (5) 以下では、この式をさらに整理してみる。まず、 C0≡ a0 2 = Cme i2π Tmt|m=0 (6)
とし、さらに、式 (5) の右辺の Bme−i 2π Tmtで m → −m とすると、 ∞ X m=1 Bme−i 2π Tmt= −∞ X m=−1 B−mei 2π Tmt そこで、B−mを改めて Cm (m < 0) で表すと、 = −∞ X m=−1 Cmei 2π T mt (7) 式 (6) と (7) を、式 (5) に代入すると、 g(t) = C0+ ∞ X m=1 Cmei 2π Tmt+ −∞ X m=−1 Cmei 2π Tmt = ∞ X m=−∞ Cmei 2π T mt (8) • ここで、Cm = ¯Bmであったから、次の共役の関係が成り立つ。Bm = C−m より、 Cm= ¯C−m (9) • フーリエ係数 Cmをまとめておこう。 Cm= 1 2(am− ibm) (10) に、式 (2) を代入して、 Cm = 1 T Z T 0 g(t)(cos(2π T mt) − i sin( 2π T mt))dt = 1 T Z T 0 g(t)e−i2π T mtdt (11) なお、この式は複素積分の形をしているが、計算するには、 Cm = 1 T Z T 0 g(t) cos(2π T mt)dt − i 1 T Z T 0 g(t) sin(2π T mt))dt として、二つの実積分に分けて計算すればよい。 Cm C-m Cm eimt C-m e-imt 図 1: Cmと ¯Cm φm b /2m a /2m 図 2: フーリエ係数 • Cmと am, bmとの関係 式 (10) の am, bmは、式 (1) のフーリエ級数の係数であった (図 2) ( Cmの実部 = Re(Cm) = 12am Cmの虚部 = Im(Cm) = −12bm
式 (1) を次の形に変形すると g(t) = 1 2a0+ ∞ X m=1 q am2+ bm2sin µ 2π T mt + φm ¶ ここで、 Ã sin φm= p bm am2+ bm2 , cos φm= p am am2+ bm2 ! よって、 ( |Cm| = p am2+ bm2 ⇔ 周波数mTの振幅の絶対値 arg Cm= φm ⇔ 位相差 12.6.3 フーリエ変換 • 複素フーリエ級数の式 (11) を、ωm= 2mπT とおいて、さらに積分区間を −T2 だ けずらすと(周期関数なので積分の値は変わらない) Cm= 1 T Z T 2 −T 2 g(t)e−iωmtdt ここで、G(ωm) T = Cmとおくと、 G(ωm) = Z T 2 −T 2 g(t)e−iωmtdt となる。この式で T → ∞ (2πT → 0) とすると、次式に収束する。 G(ω) = Z ∞ −∞ g(t)e−iωtdt これを g のフーリエ変換という。 • フーリエ級数は周期信号を変換したが、フーリエ変換は非周期信号を変換でき、 周波数領域 −∞ < ω < ∞ で連続なスペクトルをもつ。
12.6.4 離散フーリエ変換 (DFT, Discrete Foureir Transform)
• 名前は離散フーリエ変換だが、その実は、離散信号に対するフーリエ級数。た だし、フーリエ級数は周期信号を対象とする。そこで非周期信号に対して周期 を仮定して離散時間フーリエ級数を計算したものが離散フーリエ変換であると 考えることができる。 • 図 3 に示すように、離散的に連続する信号から、長さ T の区間の信号を切り出 してきて、それを繰り返しつないで、周期 T の関数と見なす (基本周波数 T1、 角周波数 2πT)。 • 時間軸領域の離散化Ã サンプリング間隔 ∆ (サンプリング周波数 Sr= 1/∆) 分割数 N T = N ∆
g 0 t T T : 図 3: 非周期信号から長さ T の区間を切り出し、周期 T の信号と考える。 0 ∆ ∆ k (N-1)∆ N =T∆ g(k ) = g∆ k 図 4: フーリエ係数を求める積分を矩形領域の和として近似 • 図 4 に示すように、式 (11) の積分を長方形の面積の和で表すと、 Cm= 1 T N −1X k=0 ³ g(k∆)e−i2πT mk∆· ∆ ´ ここで、g(k∆) = gk(第 k 番目のサンプル値) とおき、また ∆ を P の外に括り 出す。 =∆ T ÃN −1 X k=0 gke−i 2π Tmk∆ ! そして、∆/T = 1/N の関係を使うと、 Cm= 1 N N −1X k=0 gke−i 2π Nmk m ∈ (−∞, ∞) (12) このような計算を離散フーリエ変換 (DFT) といい、DFT 数 Gmを Gm= N Cm で定義する。 Gm≡ N Cm= N −1X k=0 gke−i 2π Nmk m ∈ (−∞, ∞) (13) • 帯域制限 (m の制限) Gmは、周波数mT のスペクトル成分を表す。サンプリング周波数 Srが標本化定 理に従って選ばれているとすると、g(t) は、Sr=∆1 の 12 以下の周波数成分し か含まないので、 − 1 2∆≤ m T ≤ 1 2∆ この辺々に T = N ∆ を掛けると、 −N 2 ≤ m ≤ N 2
式 (11) では、m ∈ (−∞, ∞) であったが、上式の範囲外の m の Gmは 0 とし てよい。従って、DFT 係数は、 Gm= N −1X k=0 gke−i 2π Nmk −N 2 ≤ m ≤ N 2 (14) • 式 (14) は、一見複雑な式だが、 WN = e−i 2π N とおくと Gm= N −1X k=0 gkWNmk • 対称性 Gmの共役を取ってみると、 Gm = N −1X k=0 gkWN km = N −1X k=0 gk(e−i 2π N)km = N −1X k=0 gk(WN−1)km = N −1X k=0 gkWN−km= G−m よって、Gm= G−mとなり、このことから、次の関係が成り立つ。 – Gm= am+ bmi、G−m= a−m+ b−mi と表すと、 am = a−m (実部は正負で線対称) bm = −b−m (虚部は原点に対して点対称) – |G−m| = |Gm| つまり、図 5 に示すように、N 個の信号を DFT 変換すると、0 を中心にして、 左右に N2 個のスペクトルが対称に現れる。 |G| 0 N/2 N/2 図 5: DFT スペクトルの対称性 • Gmの周期性を用いた公式 (m が負の成分をシフト) (WN)N = (e−i 2π N)N = e−i2π= 1 の関係を使うと、(WN)n+N = (WN)n(WN)N = (WN)nとなる。一般に、(WN)n = (WN)mod(n,N ) これを使って Gm+N を計算してみると、 Gm+N = X k gkWN(N +m)k= Gm
よって、−N2 ≤ m < 0 で Gmを求める代わりに、それに N を加えたN2 ≤ m < N で Gmを求めてもよい. Gm= N −1X k=0 gkWNkm m = 0, . . . , N − 1 ただし、m ≥ N2 の Gmは、元の Gm−N (m − N は負) に相当する. Gm−N = Gm (m ≥ N 2) • まとめ Gm = N −1X k=0 gkWNkm (m = 0, 1, . . . , N − 1) (15) gk = 1 N N −1X m=0 GmWN−km (k = 0, 1, . . . , N − 1) (16)
この第 2 式は、DFT の逆変換で逆離散フーリエ変換 (IDFT, Inverse Discrete Fourier Transform) と呼ばれる。gk= g(k∆) であるから、式 (8) で t = k∆ と すると、 gk = g(k∆) = ∞ X m=−∞ Cmei 2π Tmk∆= N −1X m=0 1 NGme i2π Nmk= 1 N N −1X m=0 GmWN−km ここで、Pは Gmの定義される範囲だけで加えればよく、また、Cm= N1Gm と ∆/T = 1/N の関係を用いた。 • 8 点 DFT の計算例 N 個のデータ点を使った DFT を”N 点 DFT” という。N = 8 の時の G0, G1と G2を具体的に計算してみよう。 G0 = 7 X k=0 gkW80·k= 7 X k=0 gkW80= 7 X k=0 gk= g0+ g1+ g2+ · · · + g7 G1 = 7 X k=0 gkW81·k= 7 X k=0 gkW8k = g0W0 8 + g1W81+ g2W82+ g3W83+ g4W84+ g5W85+ g6W86+ g7W87 = g0+ g1e−iπ
4 + g2e−iπ2 + g3e−i3π4 + g4e−iπ+ g5e−i5π4 + g6e−i3π2 + g7e−i7π4 = g0+ g1( √ 2 2 − i √ 2 2 ) + g2(−i) + g3(− √ 2 2 − i √ 2 2 ) +g4+ g5(− √ 2 2 + i √ 2 2 ) + g6i + g7( √ 2 2 + i √ 2 2 ) G2 = 7 X k=0 gkW82·k = g0W82·0+ g1W82·1+ g2W82·2+ g3W82·3+ g4W82·4+ g5W82·5+ g6W82·6+ g7W82·7 = g0W80+ g1W82+ g2W84+ g3W86+ g4W88+ g5W810+ g6W812+ g7W814 = · · ·
• 具体例を図 6 に示す。N2 < m < N の Gmは Gm−N (Gm= Gm−N = ¯GN −m) を表しているので、図 6 では N = 8 より、4 < m < 8 に対して G5 = G−3 = G¯3 G6 = G−2 = G¯2 G7 = G−1 = G¯1 G4を中心として、G1, G2, G3と G5, G6, G7の、実部は線対称、虚部は点対称 となる。 1 0 1 2 3 4 5 6 7 -3 -2 -1 [Hz] 1 0 1 2 3 4 5 6 7 -3 -2 -1 [Hz] 0 1 2 3 4 5 6 7 -3 -2 -1 [Hz] 1 0 1 2 3 4 5 6 7 -3 -2 -1 [Hz] 1 0 1 2 3 4 5 6 7 -3 -2 -1 [Hz] 1 0 1 2 3 4 5 6 7 -3 -2 -1 [Hz] 1 0 1 2 3 4 5 6 7 -3 -2 -1 [Hz] 1 0 1 2 3 4 5 6 7 -3 -2 -1 [Hz] 1 図 6: 正弦波およびその重ね合わせの DFT の例 (f0= 1Hz; サンプリング周波数 8Hz の場合)
12.6.5 FFT (Fast Fourier Transform)
• DFT を高速に計算する方法が FFT である。
N = 4 の場合を考えてみる。 G0 = 3 X k=0 gkW40k = g0W40+ g1W40+ g2W40+ g3W40 G1 = 3 X k=0 gkW41k = g1W40+ g1W41+ g2W42+ g3W43 G2 = 3 X k=0 gkW42k = g1W40+ g1W42+ g2W44+ g3W46 G3 = 3 X k=0 gkW43k = g1W40+ g1W43+ g2W46+ g3W49 これを行列で表すと、 G0 G1 G2 G3 = W0 4 W40 W40 W40 W0 4 W41 W42 W43 W0 4 W42 W44 W46 W0 4 W43 W46 W49 g0 g1 g2 g3 (17) 一つの Gmを計算するのに、乗算が 4 回と加算が 3 回必要なことがわかる。よっ て、すべての Gmを計算するのには、乗算が 4 × 4 回と加算が 4 × 3 回必要 である。よって、一般に N 個のデータの DFT には、乗算が N2回と加算が N (N − 1) 回必要なことがわかるだろう。 FFT は、この計算回数を、それぞれ (N/2)(log2(N − 1)) と N log2N までに減 らす計算方法であり、非常に高速な計算方法である。 • 基数 2 の周波数間引き型 FFT データの数 N を 2 のべき乗とする。例えば、N = 22 = 4 や N = 210 = 1024 である。 • N = 2 の場合 G0= g0W20+ g1W20= g0+ g1 G1= g0W20+ g1W21= g0− g1 これを行列で表すと、 " G0 G1 # = " W0 2 W20 W0 2 W21 # " g0 g1 # ここで、2 × 2 行列 W2を次のように定義する。 W2≡ " W0 2 W20 W0 2 W21 # = " 1 1 1 −1 # これを図で表すと、図 7 のようになる。 g0 g1 -1 G0 G1 図 7: 2 点 FFT の計算の流れ
• 次に、N = 4 の場合を考える (式 (17))。まず、添え字 m = 0, 1, 2, 3 を次のよ うに変換し、Gmを並び替える。 1. m を 2 進数で表す。 2. ビットの並びを逆転する。 3. 10 進数に戻し、m0とする。 m 2 進数 逆転 m0 0 00 00 0 1 01 10 2 2 10 01 1 3 11 11 3 式 (17) の Giを m0の順に並べ替える。 G0 G2 G1 G3 = W0 4 W40 W40 W40 W0 4 W42 W44 W46 W0 4 W41 W42 W43 W0 4 W43 W46 W49 g0 g1 g2 g3 このように、周波数領域の Gmを並び替えることから、周波数間引き型 FFT と呼ばれる。時間領域の並び替えを行う方法もあり、時間間引き型 FFT と呼 ばれる。 • この 4 × 4 行列を W4とする。 W4≡ W0 4 W40 W40 W40 W0 4 W42 W44 W46 W0 4 W41 W42 W43 W0 4 W43 W46 W49 これを4つのブロックに分けて考える。 W4≡ W0 4 W40 W40 W40 W0 4 W42 W44 W46 W0 4 W41 W42 W43 W0 4 W43 W46 W49 ≡ " A B C D # A, B, C, D について調べて見る。 A = " W0 4 W40 W0 4 W42 # B = " W0 4 W40 W4 4 W46 # = " W0 4 W40 W0 4 W42 # = A ここで、WN2k = (e−i2π N)2k= (e−i 2π N/2)k = Wk N/2の関係を使うと、W40= W20、 W2 4 = W21であるから、 A = B = " W0 2 W20 W0 2 W21 # = W2
一方、 C = " W0 4 W41 W0 4 W43 # = " W0 4 W40 W0 4 W42 # " W0 4 0 0 W1 4 # = " W0 2 W20 W0 2 W21 # " W0 4 0 0 W1 4 # = W2V2 ただし、 V2≡ " W0 4 0 0 W1 4 # D = " W2 4 W43 W6 4 W49 # = " W2 4 W43 W2 4 W45 # = W2 4 " W0 4 W41 W0 4 W43 # = W2 4C = −1 · W2V2 以上をまとめると W4 = " W2 W2 W2V2 −W2V2 # = " W2 O2 O2 W2 # " I2 O2 O2 V2 # " I2 I2 I2 −I2 # ここで、O2は 2 × 2 の零行列で、I2は 2 × 2 の単位行列。 よって、 G0 G2 G1 G3 = " W2 O2 O2 W2 # " I2 O2 O2 V2 # " I2 I2 I2 −I2 # g0 g1 g2 g3 さて、この式の右辺の右から三つの項の掛け算を g0 0 g0 1 g0 2 g0 3 = " I2 O2 O2 V2 # " I2 I2 I2 −I2 # g0 g1 g2 g3 と表すと、 G0 G2 G1 G3 = " W2 O2 O2 W2 # g0 0 g0 1 g0 2 g0 3 となり、これは次の二つの 2 点 FFT を計算していることに等しいことが分か るだろう。 " G0 G2 # = W2 " g0 0 g0 1 # " G1 G3 # = W2 " g0 2 g0 3 # • 行列の要素を書き出せば、 G0 G2 G1 G3 = " W2 W2 # 1 1 W0 4 W1 4 1 0 1 0 0 1 0 1 1 0 −1 0 0 1 0 −1 g0 g1 g2 g3
g0 g1 g3 g2 W2 -1 -1 W 0 4 W14 G0 G2 G3 G1 W2 図 8: 4 点 FFT の計算の流れ。W2は 2 点 FFT。 ここで、行列の要素がブランクのところは 0 である。これを図で表すと、図 8 のようになる。 • N = 23= 8 の場合の 8 点 FFT は、同様にして 4 点 FFT に帰着できる。まず、 添え字を並び替える。操作は 4 点 FFT と同じ。 m 2 進数 逆順 m0 0 000 000 0 1 001 100 4 2 010 010 2 3 011 110 6 4 100 001 1 5 101 101 5 6 110 011 3 7 111 111 7 m0の順に G iを並べて、 G0 G4 G2 G6 G1 G5 G3 G7 = " W4 O4 O4 W4 # " I4 O4 O4 V4 # " I4 I4 I4 −I4 # g0 g1 g2 g3 g4 g5 g6 g7 = · W4 O4 O4 W4 ¸ 1 1 1 1
O
4O
4 W0 8 W1 8 W2 8 W3 8 1 0 0 0 1 0 0 0 0 1 0 0 0 1 0 0 0 0 1 0 0 0 1 0 0 0 0 1 0 0 0 1 1 0 0 0 −1 0 0 0 0 1 0 0 0 −1 0 0 0 0 1 0 0 0 −1 0 0 0 0 1 0 0 0 −1 g0 g1 g2 g3 g4 g5 g6 g7 これを図で表すと、図 9 のようになる。 • FFT の計算量 計算の流れ図で見ると、FFT の計算が図 10 に示す演算の組み合わせになって いることが分かるだろう。これをバタフライ演算という。 N 点の FFT をやるための、バタフライ演算の段数は、log2N となる。各段で N/2 回のバタフライ演算が行われる。また、各バタフライ演算では、2 回の加g0 g1 g3 g2 g5 g4 g6 g7 W4 W4 -1 -1 -1 -1 W08 W18 W28 W38 G0 G4 G6 G2 G5 G1 G3 G7 図 9: 8 点 FFT の計算の流れ。W4は 4 点 FFT。 -1 W i N 図 10: バタフライ演算 算が行われるので、加算の回数は、 log2N × N 2 × 2 = N log2N 一方乗算は、WNi を掛けるところでだけ行われる。最後の 2 点 FFT のステー ジでは乗算がないので、乗算が行われるステージの数は、log2N − 1 となる。 それらの各ステージで N/2 回の乗算が行われるので、乗算の回数は (log2N − 1) × N 2 となる。
12.6.6 離散コサイン変換 (DCT, Discrete Cosine Transform)
• 離散コサイン変換 (DCT) は DFT の変形と考えることができる。DCT はすべ ての計算を実数で行うので、DFT より簡単といえる。 • DCT は、図 11 のように、N 個のデータ N 個のデータ g0, g1, . . . , gN −1を鏡像 に反転して並べた 2N 個のデータに対して離散フーリエ変換を適用する。つま り、g0, g1, . . . , gN −1を反転して gN −1, gN −2, . . . , g0をつくり、それらをつない で次のような 2N 個のデータからなる信号を作る。 gN −1, gN −2, . . . , g0, g0, g1, . . . , gN −1 この左半分の添え字を先頭から −N, −N + 1, . . . , −1, 0, 1, . . . , N − 1 と付け替 え、改めて次のデータ g0kを定義する。 g−N0 , g0−N +1, g−N +20 , . . . , g−10 , g00, g01, . . . , g0N −1 g0 kを式 (15) と (16) を用いて離散フーリエ変換し G0mを計算すると、次式を 得る。 G0m = N −1X k=−N gk0W2Nmk (18) gk0 = 1 2N N −1X m=−N G0mW2N−km (19)
図 11: 離散コサイン変換、DCT まず、式 (18) を調べてみる。Pを分解して G0m= N −1X k=0 g0kW2Nmk+ −1 X k=−N g0kW2Nmk 右辺第 1 項は、k = 0, . . . , N − 1 では、gk0 = gkなので、 右辺第 1 項 = N −1X k=0 gkW2Nmk 右辺第 2 項は、k0= −k − 1 とおき、k0= 0, . . . , N − 1 で g0 −k0−1= gk0である ことを用いると、 右辺第 2 項 = 0 X k0=N −1 g0 −k0−1Wm(−k 0−1) 2N = 0 X k0=N −1 gk0Wm(−k 0−1) 2N 総和の上限と下限を入れ替えて、k0を k に置き換えると、 = N −1X k=0 gkW2Nm(−k−1) 以上より G0mは次のようになり、さらに変形を加えると、 G0 m = N −1X k=0 gk(W2Nmk+ W2Nm(−k−1)) = N −1X k=0 gkW− 1 2m 2N (W m(k+1 2) 2N + W m(−k−1 2) 2N ) = N −1X k=0 gkW− 1 2m 2N (e−i 2π 2Nm(k+12)+ e−i2N2π−m(k+12)) = W−12m 2N N −1X k=0 gk2 cos(π Nm(k + 1 2)) ここで、次式のように Gmを定義する。 Gm≡ 1 2W 1 2m 2N G0m (20) これを用いて、 Gm = N −1X k=0 gkcosmπ(2k + 1) 2N (21)
このように、Gmが実数となり、元の(折り返す前の)N 個のデータから、複 素数を使わずにコサイン関数の和として DCT 係数 Gmを表現できる。 式 (21) では、k は 0 から N − 1 まで動くのに対して、m は、−N から N − 1 まで動くが、次のような性質がある。まず m = −N に対して、 G0 −N = W −1 2(−N ) 2N N −1X k=0 gk2 cos(π N(−N )(k + 1 2)) = WN2 2N N −1X k=0 gk2 cos(−π(k + 1 2)) = 0 よって、G0−N= G−N = 0 となる。また、 G0 −m = W −1 2(−m) 2N N −1X k=0 gk2 cos(π N(−m)(k + 1 2)) = W−12m 2N N −1X k=0 gk2 cos(π N(−m)(k + 1 2)) = G0 m 次に、上の関係式も用いて、式 (19) で与えられた逆変換を計算してみる。式 (21) のように、順変換では gkしか使わないので、式 (19) も gkだけについて (つまり g0 kの内 k = 0, 1, . . . , N − 1 についてだけ) 計算すればよい。 gk= 1 2N N −1X m=−N G0 mW2N−km この式の Σ の項を四つに分解する。 N −1X m=−N G0 mW2N−km = G0−NW2NkN + G00W2N0 + N −1X m=1 G0 mW2N−km+ −1 X m=−N G0 mW2N−km = G0 −NW2NkN + G00+ N −1X m=1 G0 mW2N−km+ −1 X m=−N +1 G0 mW2N−km 右辺第 1 項は、上で示したように 0 で、第 2 項は、式 (20) より、 G0 0= 2G0 第 3 項は、(G0mを Gmに書き直すために) (2W− m 2 2N )(12W m 2 2N) = 1 を G0mに乗じ て、次のように変形する。 N −1X m=1 G0 mW2N−km = N −1X m=1 (2W−m2 2N )( 1 2W m 2 2N)G0mW2N−km = N −1X m=1 (2W−m2 2N )GmW2N−km = N −1X m=1 2GmW−m(k+ 1 2) 2N 右辺第 4 項は、m0 = −m とおいて、第 3 項と同様にして、さらに、G0 −m= G0m を用いると、
1 X m0=N −1 G0−m0Wkm 0 2N = N −1X m0=1 G0 mWkm 0 2N = N −1X m=1 (2Wm2 2N)( 1 2W m 2 2N)G0mW2Nkm = N −1X m=1 (2Wm2 2N) 1 2W m 2 2NG0mW2Nkm = N −1X m=1 2GmW m 2 2NW2Nkm = N −1X m=1 2GmWm(k+ 1 2) 2N 最後で、Gmが実数であることから Gm= Gmを用いた。以上をまとめて、 gk = 1 2N(2G0+ N −1X m=1 2Gm(W−m(k+ 1 2) 2N + W m(k+1 2) 2N ) = 2 N( 1 2G0+ N −1X m=1 Gmcosmπ(2k + 1) 2N ) • このように、元の N 個のデータを対称に並べ替えて 2N 個のデータにして離散 フーリエ変換を適用すると、次のような、実数でコサイン関数からなる式を導 くことができる。 Gm= N −1X k=0 gkcosmπ(2k + 1) 2N (m = 0, 1, . . . , N − 1) これを DCT という。DCT 係数を使って、元の信号 gkは、次の逆離散コサイ ン変換で表される。 gk = 2 N( 1 2G0+ N −1X m=1 Gmcosmπ(2k + 1) 2N ) (k = 0, 1, . . . , N − 1) なお、Gmと gk の係数の取り方にいくつかの流儀がある。DCT については、 画像の節でさらに説明する。 12.6.7 フーリエ変換の音声信号への応用 • FFT/DCT フィルタ ディジタルフィルターによる処理は時間領域の中で行われたが、FFT では、信 号を時間領域から周波数領域へ変換することができるので、音声の周波数スペ クトルを直接処理することができるようになる。つまり、この周波数スペクト ルを目的に合わせて変形し、それを IDFT で音声信号に戻すことによってフィ ルタリングなどを行うことができる。DCT も利用される。 • 知覚符号化による圧縮 人間の聴覚には、最小可聴限界特性と聴覚マスキング効果という特性があり、 これを利用して高能率符号化を行う。最小可聴限界特性とは、人間はあるレベ ルより小さい音を聞き取ることができないことで、この特性が周波数に依存し ている。聴覚マスキング効果とは、ある周波数の大きな音がある場合、その周 辺の周波数の音は隠れて聞き取れなくなることをいう。
そこで、入力されたオーディオ信号を、所定のサンプル数毎にブロック化し、 それを周波数変換して周波数スペクトルを求め、それを元に周波数の区間毎に 聴覚のマスキングを計算し、予め設定された周波数帯域ごとの許容する量子化 雑音に基づいて、フーリエ係数を量子化し圧縮する。これによって、MPEG (Moving Picture Coding Experts Group) などの方式では、高い圧縮率を得て いる。