このページの内容は最新ではありません。最新版の英語を参照するには、ここをクリックします。

離散フーリエ変換

離散フーリエ変換 (DFT) は、デジタル信号処理の基本となるツールです。製品の基礎をなしているのは高速フーリエ変換 (FFT) で、短い実行時間で DFT を計算します。多くのツールボックス関数 (Z 領域周波数応答、スペクトル解析とケプストラム解析、および一部のフィルター設計関数やフィルター実装関数を含む) には、FFT が組み込まれています。

MATLAB® 環境には、関数 fft と ifft があり、それぞれ離散フーリエ変換と逆離散フーリエ変換の計算に使用できます。入力シーケンス x とこのシーケンスから変換した X (単位円周上における等間隔の周波数での離散時間フーリエ変換) に対し、この 2 つの関数によって次の関係が実装されます。

$$X(k+1)=\sum _{n=0}^{N-1}x(n+1)W_N^{kn}$$

と

$$x(n+1)=\frac{1}{N}\sum _{k=0}^{N-1}X(k+1)W_N^{-kn}.$$

これらの方程式では、MATLAB のベクトル インデックス方式のため、級数の添字は 0 ではなく 1 から始まります。また、次の関係があります。

$$W_N=e^{-j2\pi /N}.$$

メモ: MATLAB では、慣例により、関数 fft に負の j が使用されます。これは工学上の慣例であり、物理学や純粋数学では一般に正符号の j が使用されます。

関数 fft では、単一の入力引数 x を使用して、入力ベクトルまたは入力行列の DFT が計算されます。x がベクトルの場合、関数 fft ではベクトルの DFT が計算されます。x が方形配列の場合は、関数 fft では配列の各列の DFT が計算されます。

例として、次のような時間ベクトルと信号を作成しましょう。

t = 0:1/100:10-1/100;                     % Time vector
x = sin(2*pi*15*t) + sin(2*pi*40*t);      % Signal

信号の DFT および変換後のシーケンスの振幅と位相を計算します。小さい振幅の変換の値をゼロに設定して位相を計算する際の丸め誤差を減らします。

y = fft(x);                               % Compute DFT of x
m = abs(y);                               % Magnitude
y(m<1e-6) = 0;
p = unwrap(angle(y));                     % Phase

この振幅と位相を角度でプロットするには、次のコマンドを入力します。

f = (0:length(y)-1)*100/length(y);        % Frequency vector

subplot(2,1,1)
plot(f,m)
title('Magnitude')
ax = gca;
ax.XTick = [15 40 60 85];

subplot(2,1,2)
plot(f,p*180/pi)
title('Phase')
ax = gca;
ax.XTick = [15 40 60 85];

関数 fft の 2 番目の引数は、変換する点の数 n を指定し、DFT の長さを表しています。

n = 512;
y = fft(x,n);
m = abs(y);
p = unwrap(angle(y));
f = (0:length(y)-1)*100/length(y);

subplot(2,1,1)
plot(f,m)
title('Magnitude')
ax = gca;
ax.XTick = [15 40 60 85];

subplot(2,1,2)
plot(f,p*180/pi)
title('Phase')
ax = gca;
ax.XTick = [15 40 60 85];

この場合、関数 fft によって入力シーケンスの長さが n より短い場合にはゼロが付加され、n より長い場合にはシーケンスが切り取られます。n が指定されていない場合は、入力シーケンスの長さが既定値として使用されます。関数 fft の実行時間は、実行する DFT の長さ n に左右されます。このアルゴリズムの詳細については、fft のリファレンス ページを参照してください。

メモ: 結果として返される FFT 振幅は A*n/2 で、A は元の振幅、n は FFT データ点数です。これは、FFT データ点の数がデータ サンプルの数以上である場合にのみ真になります。FFT データ点の数が少ない場合は、FFT 振幅は元の振幅より上記の分だけ小さくなります。

逆離散フーリエ変換関数ifftも入力シーケンスを受け入れることができ、オプションで変換する点の数も指定できます。次の例を見てください。元のシーケンス x と、再構成されたシーケンスは丸め誤差の範囲内で一致しています。

t = 0:1/255:1;
x = sin(2*pi*120*t);
y = real(ifft(fft(x)));

figure
plot(t,x-y)

また、このツールボックスには、2 次元 FFT とその逆変換に使用する関数の fft2 および ifft2 が含まれています。これらの関数は、2 次元の信号やイメージ処理に便利です。Goertzel 関数は DFT を計算するもう 1 つのアルゴリズムで、ツールボックスに含まれています。この関数は、長い信号の部分に対する DFT の計算に有効です。

場合によっては、関数 fft または関数 fft2 の出力を再配列して、ゼロ周波数成分がシーケンスの中央になるようにすると便利なことがあります。関数 fftshift により、ゼロ周波数成分をベクトルまたは行列の中央に移動することができます。

参考

fft | fft2 | fftshift | goertzel | ifft | ifft2