C3：音声信号処理

C3：音声信号処理
担当：中村和晃（情報通信工学部門）
TA：長澤優佑（馬場口研究室）
1
概略・目的
本実験では，録音した音声信号に対し，コンピュータを使った簡単な加工を行うことにより，
ディジタル信号処理の基礎を修得する．なお，処理を施すためのプログラムは C++を使用して作
成する（C++は C 言語の完全上位互換であるため，C++に不慣れであれば，C 言語の構文や関数
をそのまま用いても良い）．
2
信号処理
信号処理の概略を以下の図に示す．
図 1: 信号処理の概略
I. 音声をマイクにより電気信号に変換する． II. 連続的に変化する電気信号をディジタル信号（離散的な数値の系列）に変換する．具体的に
は，サンプリング，量子化といった操作を行う．2.1，2.2 節参照．
III. 得られたディジタル信号を加工し，新たなディジタル信号を作る．
IV. ディジタル信号を連続的に変化する電気信号に変換する．
V. 電気信号をスピーカーにより音声に変換する．
1
ディジタル信号処理の本質的な部分は III. の操作が担っている．図 1 に示した通り，II. の操作
により得られるディジタル信号は単なる数値の並びである．III. の操作は，得られた数値列をも
とに新たな数値列を生成する処理であると言える．入力された数値列 x[0], x[1], x[2], · · · から新
しい数値列 y[0], y[1], y[2], · · · をどう作るか，これがディジタル信号処理のキーポイントである．
2.1
サンプリング
コンピュータは，連続的に変化する信号をそのまま扱うことができない．そこで，そのような
信号を扱う際には，一定の間隔で信号の値をサンプルし，離散的なデータとして扱うことが必要
となる．この操作をサンプリングと呼ぶ．
図 2: サンプリング
図 2 において，あるデータから次のデータまでの時間間隔 Ts をサンプリング周期と呼び，fs = T1s
をサンプリング周波数と呼ぶ．つまり，サンプリング周波数とは，1 秒間にデータをサンプルする
回数を表す（単位は [Hz]）．一般的な音楽データのサンプリング周波数は 44.1 [kHz] であり，す
なわち，1 秒間に 44100 回サンプルされた信号が使用されている．この実験では，それより少な
い 8000 [Hz] でサンプリングされた信号を扱う．
2.2
量子化
コンピュータは，ある時刻にサンプルした値をメモリに蓄えるとき，その値を有限個の 0 と 1
からなるデータとして記録する．この操作を量子化という．コンピュータは，ひとつの値を蓄え
るのに有限個の 0 と 1 しか使用できないため，もともとの値を正確に表現することはできない．
従って，何らかの近似を行ってデータを蓄えることになる．データの値を表現するのに使える 0
と 1 の個数を量子化ビット数と言うが，この数が多いほど，より正確にデータを表現できる．例
えば，量子化ビット数が 2 ビットの場合，値を 00，01，10，11 の 4 段階でしか表現できないが，
16 ビットの場合，216 = 65536 段階で表すことができる．
一般的な音楽データでは 16 ビットの量子化が採用されている．この実験でも同様に 16 ビット
で量子化された信号を扱う． 2
実験課題
3
以下の実験課題を実施すること．なお，
「時間が余っている人向け」としたものは早く終わった
人用の暇つぶし課題であるので，必ずしも実施しなくて良い（特に解説もしない）．
実験課題 1
(i) WavePad を使用して自分の声の録音を試みよ．録音する音声はどのようなものでも構わ
ない．録音方法については付録 1 を参照せよ．
• 本実験で使用するノートパソコンにはマイクが内蔵されているため，外部マイクを接
続する必要はない．
• 周りの人の声が混ざらないよう，適宜録音タイミングをずらすこと．
• 録音した音声ファイルの保存先はデスクトップ上にある C3 フォルダとせよ．なお，こ
の “C3” は実際にはショートカットであり，実体は以下の場所にある．
C:\cygwin64\home\student\C3
(ii) 録音した音声に対し，amp.cpp，phase_inverter.cpp，rev.cpp の 3 つを適用し，処理
結果を確認せよ．また，それぞれのソースファイルを開いて実装方法を調べ，以降のプ
ログラミングの参考とせよ．
• プログラムのコンパイルと実行は Cygwin Terminal 上で行う．デスクトップ上のショー
トカットをダブルクリック後
cd C3
と入力し，作業用フォルダに移動せよ．
• コンパイル方法は次の通り（ここでは amp.cpp を例とする）．
./C3gcc.sh amp.cpp
上記のコマンドを実行すると，特にコンパイルエラーがなければ，同名の実行ファイ
ル（例えば amp.exe）が生成される1 ．
• プログラムの実行方法は次の通り（ここでは amp.exe を例とする）．
./amp input.wav output.wav
ここで input.wav は入力ファイル名，output.wav は出力ファイル名である．
実験課題 2
入力 x[n] に対し，
y[2n] = x[n],
y[2n+1] = 0
(1)
として出力 y[n] を生成する処理をアップサンプリングと呼ぶ．
1
C3gcc.sh は本実験のために用意したシェルスクリプトであり，一般的な UNIX コマンドではないので，本実験以
外では使用できない．興味のある人はシェルスクリプトの入門書を参照．
3
(i) このアップサンプリング処理を実行するプログラムを実装せよ（up.cpp 内の関数 up を
完成させればよい）．また，実装したプログラムに実験課題 1(i) で録音した音声を入力
し，結果を確認せよ．特に，出力信号のスペクトログラムを入力信号のスペクトログラ
ムと比較することにより，出力信号の高周波域で雑音が生じていることを確認せよ．ス
ペクトロラムの表示方法については付録 2 を参照．
(ii) 高周波域で雑音が生じる理由について理論的に考察せよ（講義時間中に解説する）．
(iii) ［時間が余っている人向け］レポート課題 6 の処理を実装し，結果を確認してみよ．
(iv) ［時間が余っている人向け］式 (1) を拡張して
y[Kn] = x[n],
y[Kn+k] = 0
(k = 1, 2, · · · , K −1)
(2)
とした場合，どのような結果となるか．実際に実装し，得られた結果について考察を加
えよ．ただし K ≥ 3 とする．
実験課題 3
入力 x[n] に対し，
1
(x[n] + x[n−1])
2
として出力 y[n] を生成するシステムを平均化フィルタと呼ぶ．
y[n] =
(3)
(i) この平均化フィルタを実行するプログラムを実装せよ（mean.cpp 内の関数 mean を完成
させればよい）．また，実験課題 2(i) で生成したアップサンプリング後の信号（高周波
域に雑音を含むもの）に対し平均化フィルタを適用し，その効果を確認せよ．
(ii) 平均化フィルタの性質について理論的に考察せよ（講義時間中に解説する）．
(iii) ［時間が余っている人向け］式 (3) を拡張して
D−1
1 ∑
y[n] =
x[n−d]
D
(4)
d=0
とした場合，どのような結果となるか．実際に実装し，得られた結果について考察を加
えよ．ただし D ≥ 3 とする．
実験課題 4
入力 x[n] に対し，
1
(x[n] − x[n−1])
2
として出力 y[n] を生成するシステムを差分フィルタと呼ぶ．
y[n] =
(5)
(i) この差分フィルタを実行するプログラムを実装せよ（diff.cpp 内の関数 diff を完成さ
せればよい）．また，実験課題 2(i) で生成したアップサンプリング後の信号（高周波域
に雑音を含むもの）に対し差分フィルタを適用し，その効果を確認せよ．
4
(ii) 差分フィルタの性質について理論的に考察せよ（講義時間中に解説する）．
(iii) ［時間が余っている人向け］レポート課題 7 の処理を実装し，結果を確認してみよ．ま
た，その結果について考察を加えよ．
(iv) ［時間が余っている人向け］レポート課題 7 の処理において，D が 3 以上の奇数とした
場合はどうなるか．実際に実装し，得られた結果について考察を加えよ．
実験課題 5
サンプリング周波数 fs = 8000 [Hz] の下で f [Hz] の連続 cos 波をサンプリングすることにより離
散 cos 波を生成することを考える．
(i) 上述の離散 cos 波を生成するプログラムを実装せよ（cos_wave.cpp 内の関数 cos_wave
を完成させればよい）．
• f [Hz] の連続 cos 波を式で書くと
y(t) = A cos(2πf t)
(6)
となる（A は振幅）．この信号に対し，まず 0 番目のデータ y[0] を時刻 t = 0 でサンプ
ルし，以降，Ts 秒ごとに次のデータをサンプルすることにより離散信号 y[n] を生成す
るものとする．このとき，n 番目のデータ y[n] は時刻 t = nTs でサンプルされるため，
その値は次のようになる．
y[n] = A cos(2πf nTs )
• 上記のサンプリング間隔 Ts とは，すなわちサンプリング周期のことであり，これはサ
ンプリング周波数 fs の逆数に等しい．
(ii) 上記 (i) のプログラムを用いて f = 1000k (k = 0, 1, 2, · · · , 12) [Hz] の離散コサイン波
を生成し，どのような音に聞こえるか（何 [Hz] の音に聞こえるか）を確認せよ．
• この実験課題では
./cos_wave 1000 output.wav
のようにしてプログラムを実行する．ここで，
「1000」の部分は生成する cos 波の周波
数（単位は [Hz]）であり，output.wav は出力ファイル名である．
(iii) 上記 (ii) のような現象が起きる理由について考察せよ（講義時間中に解説する）．
(iv) ［時間が余っている人向け］式 (6) を
∑
y[n] =
Ak cos(2πfk nTs )
k
のように拡張すると周波数の異なる複数の cos 波の合成波を得ることができる（Ak , fk
はそれぞれ k 個目の cos 波の振幅と周波数）．このような合成波はどのような音に聞こ
えるか，実際に実装して確かめてみよ．
5
実験課題 6
カラオケなどで使われるエコーは，次式に示すように，入力信号を少しずつ遅らせた信号を足し
合わせることによって実現できる.
y[n] = a0 x[n] + a1 x[n−d] + a2 x[n−2d]
(7)
(i) 上式の処理を実装し，d, a0 , a1 , a2 を適切に調節することによって，エコー効果を実現
するシステムを作成せよ（echo.cpp 内の関数 echo を完成させればよい）.
(ii) エコーのかかっている音声信号からエコー効果を取り除くシステムをエコーキャンセラ
という．エコーキャンセラはエコーシステムの逆システムである．このことを踏まえ，エ
コーキャンセラを表現する式を伝達関数の観点から求めよ（講義時間中に解説する）．
(iii) 上記 (ii) で求めた式に基づいてエコーキャンセラを実装せよ（echo_canceler.cpp 内の
関数 echo_canceler を完成させればよい）．また，実装したプログラムに上記 (i) の結
果を入力し，元の信号が復元されることを確認せよ．
(iv) ［時間が余っている人向け］式 (7) では足し合わせる遅延信号の数は有限であったが，
y[n] = x[n] + ax[n−d] +
∞
∑
bk x[n−kd]
k=2
のように，無限個の遅延信号を足し合わせることを考える．このようなシステムの実現
法を考え，実際に実装してみよ．
実験課題 7
長さ N の離散信号 x[n]（この x[n] は実信号に限らない）に対する N 点離散フーリエ変換（Discrete
Fourier Transform; DFT）は
N−1
∑
X[k] =
x[n]WNnk
(8)
n=0
で与えられる（0 ≤ k < N ）．ただし WN = e−j N である．
2π
(i) 上式の直接計算により X[k] (0 ≤ k < N ) を求めるプログラムを実装せよ（DFT.cpp
内の関数 DFT を完成させればよい）．また，実装したプログラムに sample4096.wav,
sample8192.wav, sample16384.wav, sample32768.wav を入力したときの実行時間を
それぞれ測定せよ．
• この実験課題の出力は音声データではないので，出力ファイルの形式は.wav ではなく
.csv または.txt とすること．
• プログラムの実行時間は常に同じではなく，他の処理や OS の状態によってある程度
の範囲で変動する．DFT の直接計算では元々の実行時間が長いため変動の影響は少な
いが，次に述べる FFT では，元々の実行時間が短い分，変動の影響が相対的に大きく
なる．このため，少なくとも FFT においては，実行時間の計測を複数回実行し，その
平均値や中央値を採用することが望ましい．
6
(ii) 付録 3-6 にて紹介している通り，DFT の効率的な計算法として高速フーリエ変換（Fast
Fourier Transform; FFT）がある．FFT により X[k] (0 ≤ k < N ) を求めるプログラム
を実装せよ（FFT.cpp 内の関数 FFT を完成させればよい）．また，実装したプログラム
に sample4096.wav, sample8192.wav, sample16384.wav, sample32768.wav を入力し
たときの実行時間をそれぞれ測定するとともに，処理結果が実験課題 7 のときと一致す
ることを確認せよ．
(iii) ［時間が余っている人向け］レポート課題 8 で述べているような計算方法を考案し，実
際に実装してみよ．
実験課題 8
ある閾値 θ1 , θ2 (0 ≤ θ1 < θ2 ≤ π) に対し振幅特性が


 0 (0 ≤ Ω < θ1 )
|X(Ω)| =
1 (θ1 ≤ Ω ≤ θ2 )

 0 (θ < Ω ≤ π)
2
となるようなフィルタをバンドパスフィルタ（Band-pass Filter; BPF）と呼ぶ．逆に，


 1 (0 ≤ Ω < θ1 )
|X(Ω)| =
0 (θ1 ≤ Ω ≤ θ2 )

 1 (θ < Ω ≤ π)
2
となるようなフィルタをバンドストップフィルタ（Band-stop Filter; BSF）と呼ぶ．FFT を利用
してこれらのフィルタを実装し（bandpass.cpp 内の関数 bandpass および bandstop.cpp 内の関
数 bandstop を完成させればよい），その効果を確認せよ2 ．
• 上記の θ1 , θ2 は正規化角周波数の領域における閾値であり，人間にとって必ずしも直感
的ではないので，実際には θ1 , θ2 を周波数の領域に変換した α1 , α2 (0 ≤ α1 < α2 ≤ f2s )
を入力するものとする．
• N 点 DFT により計算された X[k] は，周波数の領域では
リエ係数に相当する．
2
k
N fs
[Hz] の成分に対するフー
通常，これらのフィルタを離散フーリエ変換に基づいて実現することはせず，平均化フィルタ等と同様に時間領域
で処理を実現する．この実験課題の内容は，FFT を使用した一種の「遊び」である．
7
4
レポート
形式・内容：
• 以下のレポート課題に答えること．基本課題への回答は必須とする．発展課題には必
ずしも回答しなくて良いが，回答した場合は加点する．
• PDF もしくは Microsoft Word 形式で作成したファイルをメールで提出すること．メー
ルタイトル（subject）と本文（body）の双方に学籍番号と氏名を忘れずに書くこと．
• 各実験課題で作成したプログラムのソースファイルも併せて提出すること．ただし「時
間が余っている人向け」としたものについては提出しなくてよい．
締め切り： 2 週間後の木曜日 17:00 まで
提出先： [email protected]
【基本課題】レポート課題 1
以下の音声データについて，各々のスペクトログラムを図示し，その内容について考察せよ．
(a) 実験課題 1 で録音した音声データ
(b) 音声データ (a) に実験課題 2 のアップサンプリング処理を施したデータ
(c) 音声データ (b) に実験課題 3 の平均化フィルタを適用したデータ
(d) 音声データ (b) に実験課題 4 の差分フィルタを適用したデータ
スペクトログラムとは何を表現したものであるか，データ (b) のスペクトログラム中にはどのよ
う成分が現れているか，それらの各成分はそれぞれ何に相当するか，どのデータとどのデータを
比較することで何が分かるのか，などを（必要に応じて数式などを交えて）論じ，論理的に考察
すること．
【基本課題】レポート課題 2
実験課題 5(ii) で生成した離散 cos 波について，入力した周波数 f を横軸，実際に聞こえた音の
周波数 fˆ を縦軸にとり，f と fˆ の関係をグラフとして図示せよ．また，そのグラフのような結果
が得られる理由について，数式を交えて論じよ．その際，数式中に現れる変数が具体的にどのよ
うな値となるか，各 k ごとに述べること．
【基本課題】レポート課題 3
実験課題 6 で実装したエコーシステムについて，処理結果の音声波形（スペクトログラムでは
ない）を図示するとともに，その際に用いたパラメータ d, a0 , a1 , a2 の値を示せ．また，そのパ
ラメータ値が妥当と考えられる理由について論じよ．エコーとはどのような現象か，式 (7) にお
ける x[n], x[n−d], x[n−2d] はそれぞれ何に相当するか，パラメータ d, a0 , a1 , a2 はそれぞれど
のような意味を持つか，などを踏まえて説明すること．
8
【基本課題】レポート課題 4
N 点 DFT を直接計算する場合と FFT により計算する場合のそれぞれについて，必要となる複
素乗算の回数を述べよ（答えだけでなく理由も述べること）．その回数（理論値）と実験課題 7 で
計測した実行時間を比較し，結果について考察せよ．
【基本課題】レポート課題 5
本実験の感想を述べよ．
【発展課題】レポート課題 6
実験課題 2 では
y[2n] = x[n],
y[2n+1] = 0
としてアップサンプリングを行った．これを
y[2n] = x[n],
y[2n+1] =
x[n] + x[n+1]
2
(9)
のように変更すると，高周波域で生じる雑音を低減できる．この理由を周波数特性の観点から説
明せよ．
［ヒント］
入力 x[n]，出力 y[n] の DTFT をそれぞれ X(Ω), Y (Ω) とする．実験課題 2 にて証明した通り，
式 (1) のアップサンプリングでは
Y (Ω) = X(2Ω)
という関係式が成り立つ．ここで，仮に式 (9) のアップサンプリングにおいて X(Ω) と Y (Ω) の
間の関係式が
Y (Ω) = H(Ω)X(2Ω)
のようになったとする．このとき，上式の意味するところは，x[n] を式 (1) の方法でアップサン
プリングして，その後に周波数特性が H(Ω) であるようなフィルタをかけると，式 (9) の方法で
アップサンプリングした場合と同じ出力が得られる，ということである．実験課題 3 の内容から，
振幅特性 |H(Ω)| が低周波域で大きな値をとり高周波域で小さな値をとるような関数であれば，そ
の H(Ω) により表されるフィルタはローパスフィルタとなり，式 (1) の方法でアップサンプリン
グしたときに生じる高周波域の成分が減衰され，雑音低減の効果が生まれることになる．以上の
議論から，
• 式 (9) のアップサンプリングでは X(Ω) と Y (Ω) の関係が Y (Ω) = H(Ω)X(2Ω) のような形
になること
• その H(Ω) により表されるフィルタがローパスフィルタとしての振幅特性を持つこと
を示せば，式 (9) により高周波域の雑音が低減できることを説明できたことになる．
9
【発展課題】レポート課題 7
実験課題 4 で導入した差分フィルタ
y[n] =
1
(x[n] − x[n−1])
2
を拡張し，
y[n] =
D−1
1 ∑
x[n−d](−1)d
D
(10)
d=0
という式で表されるフィルタを考える（ただし D は 4 以上の偶数とする）．このフィルタは基本
的にはハイパスフィルタとして機能するが，特定の周波数の成分を完全に遮断してしまう効果も
ある．この理由を振幅特性の観点から説明せよ．また，遮断される周波数 fc (D) は具体的にどの
ような式となるか，答えよ．ただし，この問題においてサンプリング周波数は fs = 8000 [Hz] で
あるとする．
［ヒント］
式 (10) で表されるフィルタの周波数特性を H(Ω) とする．実験課題 3, 4 の内容から，このフィ
ルタが特定の周波数成分を遮断するということは，|H(Ω)| = 0 を満たす Ω が存在するというこ
とである．そのような Ω を Ωc (D) とおくことにする．正規化角周波数 Ω と周波数 f の関係は付
録 3-1 より
fs
f =
Ω
2π
であるから，fc (D) は
fc (D) =
fs
Ωc (D)
2π
で与えられることになる．
では次に，H(Ω) が具体的にはどのような式となるかを考える．以下，入力 x[n]，出力 y[n] の
離散時間フーリエ変換をそれぞれ X(Ω), Y (Ω) とおく．DTFT の線形性および時間シフトに関す
る性質を踏まえ，式 (10) の両辺に対しその DTFT を計算すると，
Y (Ω) =
D−1
D−1
}
1 ∑{
X(Ω) ∑ ( −jΩ )d
−e
X(Ω)e−jdΩ (−1)d =
D
D
d=0
d=0
となる．従って，
D−1
Y (Ω)
1 ∑ ( −jΩ )d
H(Ω) =
=
−e
X(Ω)
D
d=0
という式が得られる．
10
【発展課題】レポート課題 8
実験課題 7 で実装した FFT では，N 点 DFT の計算に際し N2 log N 回の複素乗算が必要とな
る．しかし，入力 x[n] (n = 0, · · · , N −1) が実信号，すなわち任意の n について x[n] が実数で
ある場合に限り，複素乗算の回数をおおよそ N4 log N 回まで抑えられる方法がある．その方法を
述べよ．また，なぜその方法により複素乗算の回数が抑えられるのかを説明せよ．なお，この問
題において N は 2 のべき乗であるとする．
［ヒント］
実信号 x[n] の離散フーリエ変換を X[k] (k = 0, · · · , N −1) とすると，
X[N −k] =
=
=
=
N−1
∑
n=0
N−1
∑
n=0
N−1
∑
n=0
N−1
∑
n(N−k)
x[n]WN
x[n]WNnN WN−nk
x[n]WN−nk
x[n] WNnk
・
・
・ x[n] は実数なので x[n] = x[n]
n=0
=
N−1
∑
x[n]WNnk = X[k]
n=0
であることから，X[N −k] は X[k] の複素共役となることが分かる．ここで，複素共役を求める
のに複素乗算は必要ないことに注意．
11
付録 1: WavePad による音声の録音
本実験では，音声の録音・再生や波形データ・スペクトル確認のためのソフトウェアとして
WavePad（お試し版）を用いる．お試し版の WavePad はフリーウェアであり，下記のサイトか
らダウンロードすることで自宅の PC 等でも使用可能である．ただし，使用目的は非営利のもの
に限られるので注意されたい．
http://www.nch.com.au/wavepad/jp/
WavePad による音声の録音方法は次の通りである．まず，マイクが PC に接続されていること
を確認したのち，WavePad を起動する（ただし，本実験ではマイク内臓のノート PC を使用する
ため，事前にマイクを接続する必要はない）．起動すると，図 3 に示すような画面が立ち上がる．
図 3: WavePad の画面
次に，図 3 の A. の部分をクリックし，サンプリング周波数などの設定を変更する．本実験では，
サンプリング周波数として 8000 [Hz] を採用するので，図 4 の画面でサンプルレートを「8000」に
合わせる．また，チャンネルは「モノラル」に合わせる．以上の後，
「デフォルト設定」をクリッ
クしてこれらの設定を反映させる．録音は，図 3 の B. の部分をクリックした瞬間から開始され，
C. の部分をクリックした時点で終了となる．終了後，D. の部分をクリックして録音した音声を
ファイルに保存する．このとき，図 5 に示すような「波形設定」という画面が表示されるが，特
に気にせず「OK」を押せばよい．
なお，本実験では，音声データの保存先はデスクトップ以下の C3 フォルダとすること．また，
保存形式としては「.wav」を選択すること．
12
図 4: WavePad サンプリング周波数等設定画面
図 5: WavePad セーブ画面
付録 2： WavePad によるスペクトログラムの表示
スペクトログラムとは，各時刻において信号にどの周波数の成分がどの程度含まれているかを
画像で表したものである．WavePad では，図 6 に示す手順によりスペクトログラムを表示できる．
図 6: WavePad スペクトログラム
図 6 のスペクトログラムにおいて，横軸は時間軸を，縦軸は周波数の軸を表す．色は各時刻に
おける各周波数成分の強さに対応しており，白に近いほど対応する成分が強いことを示している．
13
付録 3: 信号処理の理論
信号処理の理論について，この実験の基礎となる内容を簡単にまとめておくので，必要に応じ
て参照されたい．
付録 3-1: 周波数・角周波数・正規化角周波数
信号の解析とは，その信号がどのような成分から構成されているかを調べることに他ならない．
そのためには信号を複数の「成分」に分解することが必要であり，分解のためのツール（道具）が
フーリエ変換およびその亜種（離散時間フーリエ変換や離散フーリエ変換など）である．これら
の変換は一般に数式を用いて説明されるが，数式中の文字が表す周波数，角周波数，正規化角周
波数などの概念を正しく把握できていなければ，その式が表す内容を理解することは難しい．そ
こでまず，これらの概念について説明する．
なお，この実験指導書では，周波数を f で，角周波数を小文字の ω で，正規化角周波数を大文
字の Ω で表すことにしている．
周波数
周期信号に関して，1 秒間に何周期分の波が現れるかを求めたものを周波数と呼ぶ．周期 21 [秒]
の周期信号では 1 秒間に 2 周期分の波が現れる．周期 13 [秒] の周期信号では 1 秒間に 3 周期分の
波が現れる．一般化すると，周期 T [秒] の周期信号では 1 秒間に T1 周期分の波が現れる．従って，
その周波数 f は周期 T に対して
1
(11)
f = [Hz]
T
となる．
角周波数
周期信号の 1 周期分を 1 回転分の角度 2π [rad] に対応させ，信号が 1 秒間に何 rad 回転するかを
求めたものを角周波数と呼ぶ．周波数 1 [Hz] の信号は 1 秒間に 1 回転するので，その角周波数は
2π [rad/秒] である．周波数 2 [Hz] の信号は 1 秒間に 2 回転するので，その角周波数は 2π × 2 = 4π
[rad/秒] である．一般化すると，周波数 f [Hz] の周期信号の角周波数 ω は
ω = 2πf [rad/秒]
(12)
となる．
正規化角周波数
離散信号，すなわちディジタル信号においては，時間の単位として「秒」よりも「サンプル」を
用いた方が議論がしやすい．これを踏まえ，
「1 秒間に何 rad 回転するか」ではなく「1 サンプル
の間に何 rad 回転するか」を求めたものを正規化角周波数と呼ぶ3 ．サンプリング周波数が fs [Hz]
の場合，1 秒間に fs 個のサンプルが等間隔に存在するので，1 サンプルあたりの回転角度（正規
3
連続信号では，
「サンプル」というものがそもそも存在しないため，正規化角周波数の概念も存在しない．
14
化角周波数）は，1 秒あたりの回転角度（角周波数）の
角周波数 ω は
ω
Ω=
[rad]
fs
1
fs
となる．従って，正規化角周波数 Ω と
(13)
という関係で結ばれる．サンプリング周波数 fs がサンプリング周期 Ts の逆数であることに注意
すると，上式は
Ω = ωTs
(14)
と書くこともできる．式 (12) および式 (13) より，正規化角周波数 Ω と周波数 f は
Ω=
2πf
fs
⇐⇒ f =
fs
Ω
2π
(15)
という関係で結ばれることが分かる．
付録 3-2: フーリエ変換
任意の連続信号 x(t) は，様々な角周波数の sin 波と cos 波の重ね合わせにより
∫ ∞
{a(ω) cos(ωt) + b(ω) sin(ωt)}dω
x(t) =
(16)
0
として表現できる4 ．フーリエ変換 (Fourier Transform) とは，x(t) から上式の係数 a(ω), b(ω)
を求める変換である．これに対し，式 (16) 自身は，係数 a(ω), b(ω) から元の信号 x(t) を復元す
る方法を示しており，逆フーリエ変換 (Inverse Fourier Transform) と呼ばれる．一般的な教
科書では，フーリエ変換は
∫
∞
X(ω) =
x(t)e−jωt dt
(17)
−∞
と定義され，その値は複素数となる．これに対し式 (16) の a(ω), b(ω) は実数であり，両者が同じ
ものであるとは想像しにくいかもしれない．また，一般的な教科書では逆フーリエ変換は
∫ ∞
1
x(t) =
X(ω)ejωt dω
(18)
2π −∞
と定義され，式 (16) とは全く異なるもののように思えるかもしれない．以下では，これらが本質
的に同じものであることを説明する．
高校数学で学習した内容を思い返してもらうと，三角関数では合成公式
√
a cos θ + b sin θ =
a2 + b2 cos(θ + β)
が成り立つ．ただし β は
−b
sin β = √
,
a2 + b 2
a
cos β = √
2
a + b2
を満たす角度である．この合成公式を用いると，式 (16) は
∫ ∞
x(t) =
A(ω) cos(ωt + β(ω))dω
(19)
0
4
厳密には関数 x(t) が有界かつ不連続な点を含まない場合に限られるが，実用上は任意と考えてほぼ問題ない．
15
のように変形でき，cos 波のみを重ね合わせた形として捉えることもできる．この形における
√
A(ω) = {a(ω)}2 + {b(ω)}2 を角周波数 ω の成分の振幅と呼び，β(ω) を同じく角周波数 ω の成
分の位相5 と呼ぶ．振幅および位相は実数値として定義されるため，本来，これらの計算に複素数
は必須ではない．しかし，式 (19) を実際に計算するのは容易ではない．複素数は，式をより簡潔
な形に直して計算を容易化するために導入される．従って，フーリエ変換が複素数領域で定義さ
れていることに本質的な意味はない6 ．フーリエ変換の理解に際しては，この点をしっかりと意識
しておくことが肝要である．
さて，計算を容易化するため，a(ω) を実部に，−b(ω) を虚部に持つような複素数 X(ω) を考え
る．すなわち
X(ω) = a(ω) − jb(ω)
とおく．X(ω) をこのように定義すると，|X(ω)| = A(ω) かつ ∠X(ω) = β(ω) となることに注意
されたい．ここで，オイラーの公式
ejθ = cos θ + j sin θ
に従えば，
1
j
cos θ = (ejθ + e−jθ ), sin θ = − (ejθ − e−jθ )
2
2
であるから，式 (16) は次のように変形できる．
}
∫ ∞{
ejωt + e−jωt
ejωt − e−jωt
a(ω)
− jb(ω)
dω
x(t) =
2
2
0
}
∫ ∞{
a(ω) − jb(ω) jωt a(ω) + jb(ω) −jωt
=
e +
e
dω
2
2
0
∫
∫
1 ∞
1 ∞
jωt
=
X(ω)e dω +
X(ω)e−jωt dω
2 0
2 0
(20)
上記の X(ω) は本来 ω ≥ 0 の範囲で定義されるものであるが，X(−ω) = X(ω) としてその定義域
を拡張する．これにより式 (20) は更に次のように変形できる．
∫
∫
1 ∞
1 ∞
jωt
X(ω)e dω +
X(−ω)e−jωt dω
x(t) =
2 0
2 0
∫
∫
1 ∞
1 0
′
=
X(ω)ejωt dω +
X(ω ′ )ejω t dω ′
2 0
2 −∞
∫ ∞
1
=
X(ω)ejωt dω
(21)
2 −∞
上式は，定数倍の違いを除いて，一般的な逆フーリエ変換の定義式 (18) と等価であり，かつ，こ
ちらの方が式 (16) や式 (19) よりも簡潔な式となっていることが分かる．
以上の議論から分かる通り，フーリエ変換は信号 x(t) を様々な cos 波の重ね合わせとして表現
したときの各成分の振幅と位相を求める変換であり，実数領域での変換として理解することが可
能なものである．複素数を導入し，その領域で変換式が定義されるのはあくまで計算上の都合で
あり，複素フーリエ係数 X(ω) はフーリエ変換の本質ではない．重要なのは振幅 |X(ω)| や位相
∠X(ω) などの実数値である．また，成分として重要なのは ω ≥ 0 の範囲に対応するもののみであ
る．ω < 0 の範囲に対応する成分は，逆フーリエ変換の結果が常に実数となるようにするための
補助成分と解釈すれば良い．
5
6
正確には「初期位相」と呼ぶべきか．
これは実信号のみを対象とする場合の話である．複素信号のフーリエ変換を考える場合はこの限りではない．
16
付録 3-3: 離散時間フーリエ変換 (DTFT)
フーリエ変換は連続信号がどのような成分から構成されているかを解析するためのツールであ
る．これに対し，離散信号，すなわちディジタル信号がどのような成分から構成されているかを
解析するツールが離散時間フーリエ変換 (Discrete Time Fourier Transform; DTFT) であ
る．一般に DTFT は
∞
∑
X(Ω) =
x[n]e−jnΩ
(22)
n=−∞
と定義される．フーリエ変換の定義式 (17) と比較すると，積分が総和に変わっている点，変数が
角周波数 ω から正規化角周波数 Ω に変わっている点に大きな違いがある．しかし，DTFT とフー
リエ変換は実は同じものである．以下では，このことを離散信号に「むりやり」フーリエ変換を
適用することにより確認する．
以下，連続信号 x(t) をサンプリング周期 Ts でサンプルして離散信号 y[n] を生成する場合につい
て考える．具体的には，時刻 t = 0 で 0 番目のサンプルを取り y[0] = x(0) とする．時刻 t = Ts で 1
番目のサンプルを取り y[1] = x(Ts ) とする．時刻 t = 2Ts で 2 番目のサンプルを取り y[2] = x(2Ts )
とする．一般化すると，時刻 t = nTs で n 番目のサンプルをとることにより
y[n] = x(nTs )
(23)
として y[n] を生成する（n < 0 の場合も含む）．
ここで，時刻 t = nTs (n = 0, ±1, ±2, · · · ) においてのみ非ゼロの定数値を取り，それ以外の
任意の時刻において 0 を取るような関数 γ(t) を考える．このとき，上記のサンプリング処理は，
x(t) に γ(t) を乗じる処理であるとの見方ができる（図 7 参照）7 ．これにより得られる連続信号
x(t)γ(t) は，離散信号 y[n] を連続信号へと「むりやり」拡張したものとみなせる．この x(t)γ(t)
にフーリエ変換を適用することを考える．
図 7: デルタ関数の無限列により x(t) を「むりやり」離散化
7
図 7 のようにして γ(t) をデルタ関数の無限列で表すと，x(t)γ(t) の t = nTs における値は実際には +∞ となる．
これに対し，y[n] の値は全て有限である．従って，x(t)γ(t) と y[n] は厳密には異なる．もし x(t)γ(t) の値が t = nTs
においても有限であるとすると，x(t)γ(t) の積分値は積分範囲に関わらず 0 となるため，本文の議論は成立しなくな
る．要するに，連続信号を「むりやり」離散化する際にデルタ関数を導入し「有限」を「無限」にすり替えているとこ
ろがこの議論のトリックである．
17
上記のような関数 γ(t) は，デルタ関数の無限列を用いて
γ(t) =
∞
∑
δ(t − nTs )
n=−∞
と表すことができる．これを用いると，x(t)γ(t) のフーリエ変換は，角周波数 ω を用いて
∫ ∞
F[x(t)γ(t)] =
x(t)γ(t)e−jωt dt
−∞
]
∫ ∞[ ∑
∞
=
x(t)δ(t − nTs ) e−jωt dt
−∞
=
n=−∞
∞ [∫
∑
∞
−∞
n=−∞
{x(t)e
−jωt
]
}δ(t − nTs )dt
(24)
と書ける8 ．ここで，任意の連続関数 q(t) に対し
∫ ∞
q(t)δ(t − t0 ) = q(t0 )
−∞
であることに注意すると，式 (24) は更に
∞
∑
F[x(t)γ(t)] =
=
n=−∞
∞
∑
x(nTs )e−jω(nTs )
x(nTs )e−jn(ωTs )
n=−∞
のように展開できる．これに式 (14) の関係式と式 (23) のサンプリングの式を代入すると，
∞
∑
F[x(t)γ(t)] =
y[n]e−jnΩ
n=−∞
となる．これは離散信号 y[n] の DTFT と完全に一致する9 ．
以上の議論から，DTFT は離散信号を「むりやり」連続とみなして強引にフーリエ変換を行う
ことと等価であり，両者は実は同じものである．従って，DTFT はフーリエ変換と全く同様の性
質を持つ．例えば，DTFT もフーリエ変換と同じく実数の範囲内で理解できるものであり，複素
数の導入は計算を容易化するための工夫に過ぎない．従って，y[n] の DTFT を Y (Ω) としたとき，
複素数 Y (Ω) 自体は本質的な意味を持たない．重要なのは |Y (Ω)| や ∠Y (Ω) といった実数値の方
である．絶対値 |Y (Ω)| は，信号 y[n] に含まれる正規化角周波数 Ω の成分の振幅を表す．同様に，
∠Y (Ω) は，正規化角周波数 Ω の成分の位相を表す．
fs
余談であるが，式 (15) の関係式から，正規化角周波数 Ω の成分とは，周波数 2π
Ω の成分に他
ならない．従って，fs = 8000 [Hz] を採用している本実験においては，
8000
0 = 0 [Hz]
2π
8000
Ω=π
⇐⇒ f =
π = 4000 [Hz]
2π
8000
Ω = 2π
⇐⇒ f =
2π = 8000 [Hz]
2π
という対応関係に従って DTFT の結果（各周波数成分の振幅・位相）を理解すればよい．
Ω=0
⇐⇒
f=
8
厳密には無限和と積分の交換可能性を議論する必要があるが，ここでは立ち入らない．
本文でも触れたが，フーリエ変換が角周波数 ω を変数として定義されるのに対し，DTFT は正規化角周波数 Ω を
変数として定義されることに注意．
9
18
付録 3-4: DTFT の性質
DTFT の性質について，この実験に関係するものを幾つか紹介する10 ．
線形性
DTFT は線形性を持つ．線形性とは，信号 w[n] を別の信号 x[n], y[n] の線形和として w[n] =
αx[n] + βy[n] のように与えたとき，それらの DTFT X(Ω), Y (Ω), W (Ω) の間にも同じ関係
W (Ω) = αX(Ω) + βY (Ω)
が成り立つことを言う．数式による証明は次の通りである．
W (Ω) =
=
∞
∑
n=−∞
∞
∑
w[n]e−jnΩ
(αx[n] + βy[n])e−jnΩ
n=−∞
∞
∑
= α
x[n]e−jnΩ + β
n=−∞
∞
∑
y[n]e−jnΩ = αX(Ω) + βY (Ω)
n=−∞
時間シフト
信号 x[n] に d サンプル分の遅延を加えた信号を q[n] = x[n−d] とおく．x[n], q[n] の DTFT を
それぞれ X(Ω), Q(Ω) とすると，Q(Ω) = X(Ω)e−jdΩ となる．すなわち，x[n−d] の DTFT は
X(Ω)e−jdΩ である．これは，時間領域で d サンプル分の遅延を与えることが，周波数領域で e−jdΩ
を乗じることと等価である，ということを意味する．数式による証明は次の通りである．
Q(Ω) =
=
=
=
∞
∑
n=−∞
∞
∑
q[n]e−jnΩ
x[n−d]e−jnΩ
n=−∞
∞
∑
x[m]e−j(m+d)Ω
m=−∞
[ ∞
∑
]
x[m]e
−jmΩ
e−jdΩ = X(Ω)e−jdΩ
m=−∞
周期性
DTFT は必ず周期 2π で周期的となる（ここで言う「周期 2π 」とは正規化角周波数 Ω の観点か
ら見た場合の周期である）．このことを以下に示す．
DTFT の定義式
∞
∑
X(Ω) =
x[n]e−jnΩ
n=−∞
10
これらの性質は DTFT が収束することを前提として成り立つものである．
19
に従って X(Ω + 2π) を計算すると，
X(Ω + 2π) =
=
=
∞
∑
x[n]e−jn(Ω+2π)
n=−∞
∞
∑
n=−∞
∞
∑
x[n]e−jnΩ ej(−2πn)
x[n]e−jnΩ = X(Ω) n=−∞
となる（n は整数なので ej(−2πn) = 1）．上式から，X(Ω) が周期 2π で周期的であることが自明
に導かれる．
さて，式 (13) の関係式より Ω = fωs であるから，X(Ω) が周期 2π で周期的ということは，角周
波数 ω の観点から見れば，X(ω) = X(Ω)|Ω= ω は周期 2πfs で周期的ということになる．さらに，
fs
式 (14) の関係式より Ω =
2πf
fs
であるから，周波数 f の観点から見れば，X(f ) = X(Ω)|Ω= 2πf は
fs
周期 fs で周期的となる．これは，サンプリング周波数 fs が大きければ大きいほど，すなわち密に
サンプリングすればするほど周波数スペクトルの周期は長くなるということを意味する．fs → ∞
の極限では，離散信号 x[n] は元の連続信号に一致するが，このとき周波数スペクトルの周期は無
限大となり，周期性が失われる．図 8 に示すように，周波数スペクトルの周期が短いほど，元の
連続信号には存在していなかったはずの「偽の成分」が多数現れる．この偽の成分はサンプリン
グにより生じる「歪み」のようなものと捉えることができ，サンプリング周波数が小さいほどこ
の歪みが多くなるものと解釈できる．
図 8: サンプリングにより生じる周波数スペクトルの「歪み」
20
対称性
本実験では処理の対象として実信号のみを考えている．その限りにおいて，DTFT は対称性を
持つ．具体的に言うと，実信号 x[n] の DTFT を X(Ω) としたとき，X(−Ω) は X(Ω) の複素共役
となり，従って振幅 |X(Ω)| は Ω = 0 の軸に関して，位相 ∠X(Ω) は原点 (0, 0) に関して，それぞ
れ対称となる11 ．このことを以下に示す．
x[n] が実信号のとき，任意の n について x[n] = x[n] が成り立つので，
X(−Ω) =
∞
∑
x[n]ejnΩ
n=−∞
=
∞
∑
x[n] ejnΩ
n=−∞
=
∞
∑
x[n]e−jnΩ = X(Ω)
n=−∞
となる．上式から自明に |X(Ω)| = |X(−Ω)| かつ ∠X(Ω) = −∠X(−Ω) が導かれる．
以上の事実と前項で述べた周期性を合わせて考えると，次式に示すように，Ω = 0 の場合だけ
でなく Ω = kπ （k は任意の整数）の場合に関する対称性も示すことができる．
X(kπ + Ω) = X(−kπ − Ω)
= X(−kπ − Ω + 2kπ)
= X(kπ − Ω)
上式は，X(kπ + Ω) が X(kπ − Ω) の複素共役となることを示しており，従って |X(kπ + Ω)| =
|X(kπ − Ω)| かつ ∠X(kπ + Ω) = −∠X(kπ − Ω) が成り立つ．すなわち，振幅スペクトル |X(Ω)|
は Ω = kπ の軸に関して対称（線対称）となり，位相スペクトル ∠X(Ω) は点 (kπ, 0) に関して対
称（点対称）となる．
以上に述べた対称性と周期性の議論から，DTFT を用いて実信号の解析を行う場合，周波数ス
ペクトル X(Ω) のうち 0 ≤ Ω ≤ π の範囲に注目すれば十分と言える．
付録 3-5: 離散フーリエ変換 (DFT)
離散信号を理論的に解析するためのツールとしては DTFT があれば十分である．しかし，コン
ピュータによる実験的な解析を行いたい場合，DTFT では不十分である．DTFT により得られる
周波数スペクトルは連続的であるが，連続的なものはコンピュータでは扱えないからである．そ
こで，サンプリングにより連続信号を離散化したのと同様に，周波数スペクトルも離散化しよう
という発想から生まれたツールが離散フーリエ変換 (Discrete Fourier Transform; DFT) で
ある．DFT は一般に
N
−1
∑
2π
(25)
X[k] =
x[n]e−j N nk
n=0
11
そもそも，付録 3-2 で述べたように，負の角周波数に対応するフーリエ係数は，逆フーリエ変換の結果が実数とな
ることを保証するために X(−ω) = X(ω) として導入したものである．これは離散信号を対象とする DTFT でも同様
であるので，実信号の DTFT に対称性があるのは当たり前と言えば当たり前である．
21
と定義される．回転因子 WN = e−j N を導入し
2π
N
−1
∑
X[k] =
x[n]WNnk
(26)
n=0
と表記することもある．DFT の詳細な解説は市販の教科書や「ディジタル信号処理」の講義に譲
ることとし，ここでは DTFT と DFT の関係について論じることにする．
上述の通り，DFT の導入はコンピュータによる信号の実験的解析が目的であるから，その解析
対象を長さが有限の信号に限ることにする．無限長の信号はコンピュータでは扱えないからであ
る．さて，有限長の信号は，それを無限に繰り返すことで周期的な無限長信号に拡張することが
できる．この周期的な無限長信号に DTFT を適用すると，その周波数スペクトルは離散的となる．
つまり，ある意味これだけで DFT が実現できたことになる．このことを以下に示す．
周期的な離散信号 x[n] を考え，その周期を N とおく（すなわち x[n] = x[n + N ]）．このとき，
x[n] の DTFT を X(Ω) とすると，X(Ω) は
X(Ω) =
=
∞
∑
x[n]e−jnΩ
n=−∞
∞
∑
x[n + N ]e
−jnΩ −jN Ω jN Ω
e
e
= e
∞
∑
jN Ω
n=−∞
x[n + N ]e−j(n+N )Ω
n=−∞
のように計算できる．ここで，n + N = n′ として無限和部分の変数を n から n′ に変換すると，
X(Ω) = e
jN Ω
∞
∑
′
x[n′ ]e−jn Ω = ejN Ω X(Ω)
n′ =−∞
となる．従って
(
)
1 − ejN Ω X(Ω) = 0
jN Ω ̸= 1
が成り立つ．N Ω が 2π の定数倍でないとき，すなわち Ω ̸= 2π
N k（k は任意の整数）のとき，e
2π
であることから X(Ω) = 0 となる．これは，Ω = N k となる位置でしか X(Ω) は値を持たないと
いうことを意味し，すなわち X(Ω) は離散的となる．
X(Ω) が Ω = 2π
N k となる位置でしか値を持たないのであれば，周波数スペクトルの定義域を
2π
Ω = N k が満たされる位置のみに限定し，定義式を
X[k] = X(Ω)|Ω= 2π k =
N
∞
∑
x[n]e−jn N k
2π
(27)
n=−∞
のように変更しても問題はない．ここで，上式の無限和部分について，n = lN から n = (l+1)N−1
までの総和（l は任意の整数）は，
∑
(l+1)N−1
x[n]e
k
−jn 2π
N
=
N−1
∑
x[n+lN ]e−j(n+lN ) N k
2π
n=0
n=lN
=
=
N−1
∑
n=0
N−1
∑
n=0
22
x[n]e−jn N k e−j(2π)lk
2π
x[n]e−jn N k
2π
となり n = 0 から n = N −1 までの総和に一致する（lk は整数なので e−j(2π)lk = 1）．この事実か
ら，式 (27) は同じものを延々と足し続けており計算に無駄が生じていることが分かる．そこで，
総和の計算範囲を 1 周期分，n = 0 から n = N −1 までに限定してしまうことにする12 ．すると，
X[k] =
N−1
∑
x[n]e−j N nk
2π
n=0
となる．これは DFT の定義式 (25) と完全に一致する．
以上のことから，DFT は，有限長の離散信号を繰り返しにより周期的な無限長信号に拡張し，
その無限長信号に DTFT を適用することとほぼ等価であることが分かる．従って，DTFT の諸性
質は DFT でも全く同様に成り立つ．例えば，付録 3-4 で述べたように，DTFT には対称性があ
り X(π + Ω) は X(π − Ω) の複素共役となるが，これと同様の性質が DFT にも存在する．前述の
N
関係式 Ω = 2π
N k から，DTFT における Ω = π は DFT における k = 2 に相当するので，DFT の
係数 X[k] は k = N2 に関して対称となる．すなわち，任意の整数 k ′ について，N が偶数であれば
[
]
[
]
[
]
[
]
X N2 + k ′ は X N2 − k ′ の複素共役となり，N が奇数であれば X N 2+1 + k ′ は X N 2−1 − k ′
の複素共役となる．
付録 3-6: 高速フーリエ変換 (FFT)
DFT を式 (25) に従って計算すると非常に時間がかかる．特に信号長が長い場合に顕著であり，
式 (25) による DFT の直接計算は実用的ではない．
DFT を効率よく計算する方法として，高速フーリエ変換 (Fast Fourier Transform; FFT)
という手法が知られている．FFT にはいくつかの種類があるが，そのうちの一つである「時間間
引き FFT」と呼ばれる方法を以下に示す．その他の方法については，市販の教科書やディジタル
信号処理の講義を参照されたい．なお，FFT は，DFT すなわち離散フーリエ変換を高速に計算
する方法であり，通常のフーリエ変換を高速に計算する方法ではない．間違えて覚えてしまう人
もいるので注意してほしい．
信号 x[n] からその DFT 係数 X[k] を計算することを考える．ただし，x[n] の信号長を N とし，
N は 2 のべき乗であるとする．このとき，N は偶数であるから，自然数 m を用いて N = 2m と
おける．ここで，x[n] の偶数番目のみを取り出した長さ m の信号を xe [n] とおく．同様に，x[n]
の奇数番目のみを取り出した長さ m の信号を xo [n] とおく．式で表すと，それぞれ
xe [n] = x[2n],
xo [n] = x[2n+1]
である．さらに，xe [n], xo [n] の m 点 DFT をそれぞれ Xe [k], Xo [k] とおく（0 ≤ k < m）．この
とき，X[k] は
X[k] =
=
=
m−1
∑
n=0
m−1
∑
n=0
m−1
∑
x[2n]WN2nk +
m−1
∑
(2n+1)k
x[2n+1]WN
n=0
2nk
xe [n]W2m
+
WNk
nk
xe [n]Wm
+ WNk
n=0
m−1
∑
2nk
xo [n]W2m
n=0
m−1
∑
nk
xo [n]Wm
n=0
12
本当の理由は「計算に無駄が生じているから」ではなく「同じものを無限回足している式 (27) の定義では X[k] が
全て +∞ となり都合が悪いから」である．ここでは付録 3-3 のときとは逆に「無限」を「有限」にすり替えている．
23
と計算できる．k < m のとき，
Xe [k] =
m−1
∑
nk
xe [n]Wm
かつ
Xo [k] =
n=0
m−1
∑
nk
xo [n]Wm
n=0
であることから，
X[k] = Xe [k] + WNk Xo [k]
(28)
となることが分かる．一方，k ≥ m のとき，k ′ = k−m < m とおくと，
X[k] = X[k ′ + m] =
=
m−1
∑
n=0
m−1
∑
′
′
nk
nm
xe [n]Wm
Wm
+ WNk WNm
m−1
∑
′
nk
nm
xo [n]Wm
Wm
n=0
m−1
∑
′
′
nk
xe [n]Wm
− WNk
n=0
nk
xo [n]Wm
′
n=0
′
= Xe [k ] −
′
WNk Xo [k ′ ]
(29)
となることが分かる．式 (28)(29) より，N 点 DFT X[k] は，
X[k] = Xe [k] + WNk Xo [k]
X[k + m] = Xe [k] − WNk Xo [k]
(30)
として 2 つの m 点 DFT に分解できる（ただし 0 ≤ k < m）．これは N = 2m 点 DFT を 1 つ計
算する場合よりも短い時間で計算できる．分解後の m 点 DFT の計算に対しても上と同様の式変
形を再帰的に適用することで，更に効率的な計算が可能となる．
付録 3-7: 逆離散フーリエ変換 (IDFT) とその計算法
長さ N の離散信号 x[n] について，その DFT を X[k] (0 ≤ k < N ) とすると，x[n] は X[k] から
N−1
1 ∑
X[k]WN−nk
x[n] =
N
(31)
k=0
として復元できる（WN は式 (26) で述べた回転因子である）．これを逆離散フーリエ変換 (Inverse
Discrete Fourier Transform; IDFT) と呼ぶ．式 (31) の両辺について，その複素共役をとると，
x[n] =
=
N−1
1 ∑
X[k]WN−nk
N
1
N
k=0
N−1
∑
X[k]WNkn =
k=0
N−1
∑
X[k]
k=0
N
WNkn
X[k]
となる．上式の右辺は， N (k = 0, 1, · · · , N −1) を複素離散信号とみなして DFT（逆変換で
はなく正変換）を行ったものに等しい．従って，FFT を実装したプログラムがあれば，以下の手
順で高速に IDFT を計算できる．
1. X[k] の複素共役 X[k] を計算し，これをさらに
2.
X[k]
N
1
N
倍することにより
X[k]
N
を求める．
(k = 0, 1, · · · , N −1) を複素離散信号とみなし FFT（正変換）を適用する．
3. ステップ 2. の結果は x[n] に一致するので，その複素共役をとり x[n] を得る．
24
付録 3-8: z 変換
ここまで，ディジタル信号を解析するためのツールとして DTFT を取り上げ，これについて
長々と述べてきた．しかし，そもそもの問題として，式 (22) で定義される DTFT は無限和であ
り，必ずしも有限の値に収束するとは限らない．例えば x[n] = 2n のように，n の増大に伴って信
号の値が発散するようなケースでは，DTFT は収束せず，周波数解析を行うこと自体に意味がな
くなる．そこで，式 (22) の e−jnΩ に e−rn という項を乗じ，収束性を高めてやろうという発想が
出てくる．これにより，先の x[n] = 2n のようなケースでも，r に十分大きな値を設定してやれば
2e−r < 1 となるため，x[n]e−rn = (2e−r )n は発散しなくなる．
さて，e−rn を導入することにより式 (22) は
X(r, Ω) =
∞
∑
x[n]e
−jnΩ −rn
e
n=−∞
=
∞
∑
(
)−n
x[n] er+jΩ
(32)
n=−∞
のように書き直せる．ここで，er+jΩ は絶対値 er ，偏角 Ω の複素数であるため，任意の複素数に
ついて，それに対応する (r, Ω) が一意に定まる．つまり，er+jΩ は複素変数 z で置き換えること
ができる．この置き換えにより，式 (32) は
X(z) =
∞
∑
x[n]z −n
(33)
n=−∞
と書ける．式 (33) により定義される変換を z 変換という．上述の導出方法を考慮すれば自明であ
るが，r = 0 のとき，すなわち z = ejΩ のとき，z 変換は DTFT に一致する．これは，DTFT が
z 変換の特殊なケースであることを示している．
付録 3-5 で述べた通り，ディジタル信号を解析するツールとしては DTFT があれば十分である．
それでも z 変換を導入することの利点は，DTFT における e−jnΩ という項が z 変換では z −n とい
う項に置き換わっており，より簡単な形の式となっているため，計算の利便性が高いことが挙げ
られる13 ．より計算が簡単な z 変換を用いてフィルタやシステムを解析し，最後に z = ejΩ を代
入することで周波数特性や振幅特性を得るというやり方は常套手段と呼べるものである．加えて，
z がどのような値であれば式 (33) が収束するのかを考えることで，DTFT が収束するか否かを判
定することも可能である．これも z 変換を導入することの利点の一つである．
以下，z 変換に関連する話題のうち本実験に関係するものについて紹介する（z 変換の詳細につ
いては市販の教科書やディジタル信号処理の講義を参照）．
伝達関数
任意の信号 x[n] を入力として受け取り，何らかの処理を加えたのち処理後の信号 y[n] を出力
する系を一般に「システム」という．課題 3, 4 で実装するような「フィルタ」もシステムであ
る．さて，システムへの入力 x[n] とシステムからの出力 y[n] に対し，各々の z 変換をそれぞれ
X(z), Y (z) とする．このとき，
Y (z)
H(z) =
X(z)
をそのシステムの伝達関数と呼ぶ．伝達関数は，その解析を通してシステムの様々な特性を知る
ことができるため，非常に重要である．
13
付録 3-2 を思い出してほしい．フーリエ変換に複素数を導入したのも計算の利便性が理由であった．
25
周波数特性／周波数応答
伝達関数 H(z) に z = ejΩ を代入し，
H(Ω) = H(z)|z=ejΩ
としたものを，システムの周波数特性または周波数応答という．直感的に分かる通り，x[n], y[n]
Y (Ω)
の DTFT をそれぞれ X(Ω), Y (Ω) とすると，H(Ω) は X(Ω) に一致する．すなわち，周波数特性
は，出力信号の DTFT を入力信号の DTFT で割ったものと同じである．
振幅特性・位相特性
周波数特性 H(Ω) の絶対値 |H(Ω)| をシステムの振幅特性という．また，周波数特性 H(Ω) の
偏角 ∠H(Ω) をシステムの位相特性という．振幅特性は，正規化角周波数 Ω の成分がシステムに
よって何倍に増幅／減衰するかを示すものである．一方，位相特性は，正規化角周波数 Ω の成分
に対しシステムがどの程度の遅延を与えるかを示すものである．
逆システム
図 9 のように，システム L1 の出力 y[n] を受け取って別の信号 w[n] を生成するシステム L2 を
考える．L2 により y[n] から原信号 x[n] が復元されるとき，すなわち w[n] = x[n] のとき，L2 を
L1 の逆システムという．
図 9: 逆システムの例
図 9 において，x[n], y[n], w[n] の z 変換をそれぞれ X(z), Y (z), W (z) とおく．また，システ
ム L1 , L2 の伝達関数がそれぞれ H1 (z), H2 (z) であるとする．前述の通り，伝達関数とは出力の
z 変換を入力の z 変換で割ったものであるから，システム L1 について
Y (z) = H1 (z)X(z)
(34)
W (z) = H2 (z)Y (z)
(35)
が，システム L2 について
が成立する．式 (34) を式 (35) に代入することにより
W (z) = H2 (z)H1 (z)X(z)
(36)
が得られるが，L2 が L1 の逆システムである場合，w[n] = x[n] より W (z) = X(z) となるため，
式 (36) は
1
1 = H2 (z)H1 (z) ⇐⇒ H2 (z) =
H1 (z)
となる．上式は，逆システムの伝達関数が，元のシステムの伝達関数に対し，その逆数で与えら
れることを示している．
26
付録 3-8: z 変換の性質
z 変換の性質について，この実験に関係するものを幾つか紹介する14 ．
線形性
z 変換は DTFT の拡張とも言えるため，DTFT と同様に線形性が成り立つ．すなわち，信号
w[n] を別の信号 x[n], y[n] の線形和として w[n] = αx[n] + βy[n] のように与えたとき，それらの
z 変換 X(z), Y (z), W (z) の間にも同じ関係
W (z) = αX(z) + βY (z)
が成り立つ．数式による証明は次の通りである．
W (z) =
=
∞
∑
w[n]z −n
n=−∞
∞
∑
(αx[n] + βy[n])z −n
n=−∞
∞
∑
= α
x[n]z −n + β
n=−∞
∞
∑
y[n]z −n = αX(z) + βY (z)
n=−∞
時間シフト
時間シフトに関して，z 変換にも DTFT と同様の性質がある．
信号 x[n] に d サンプル分の遅延を加えた信号を q[n] = x[n−d] とおく．x[n], q[n] の z 変換をそ
れぞれ X(z), Q(z) とすると，Q(z) = X(z)z −d となる．すなわち，x[n−d] の z 変換は X(z)z −d
である．これは，時間領域において d サンプル分の遅延を与えることが，z 変換領域において z −d
を乗じることと等価であるということを意味する．数式による証明は次の通りである．
Q(z) =
=
=
=
∞
∑
n=−∞
∞
∑
q[n]z −n
x[n−d]z −n
n=−∞
∞
∑
x[m]z −(m+d)
m=−∞
[ ∞
∑
]
x[m]z −m z −d = X(z)z −d
m=−∞
14
これらの性質は z 変換が収束することを前提として成り立つものである．なお，z 変換の収束性に関しては「収束
半径」という概念が存在するが，その解説は市販の教科書やディジタル信号処理の講義に譲る．
27

Download Report