fft (MATLAB Function 参考文献)

1-次元高速フーリエ変換

表示

Y = fft(X)
Y = fft(X,n)
Y = fft(X,[],dim)
Y = fft(X,n,dim)

定義

関数 X = fft(x) と x = ifft(X) は、長さ N のベクトルに対するつぎのように表せる変換と逆変換を実現します。

ここで、

は、単位円上での N 番目の根です。

詳細

Y = fft(X) は、高速フーリエ変換(FFT)アルゴリズムを使って、ベクトル X の離散フーリエ変換(DFT)を出力します。

X が行列の場合、fft は、行列の各列のフーリエ変換を出力します。

X が、多次元配列の場合、fft は、最初に1でない次元をもつものを取り扱います。

Y = fft(X,n) は、n 点のDFT を出力します。X の長さが、nより短い場合、X は、nの長さになるようにゼロを付加します。また、X の長さが、n よりも長い場合、X は、nより長い分を打ち切ります。X が行列の場合、列の長さが、上と同様な方法で調整されます。

Y = fft(X,[],dim) と Y = fft(X,n,dim) は、次元 dim に沿って、FFT 演算を適用します。

例題

フーリエ変換の一般的な使い方は、ノイズを含む時間領域の記録から、信号の周波数成分を検出することです。1000 Hz でサンプリングされたデータを考えます。信号は、 50 Hz と 120 Hz で、それに平均値がゼロのノイズを重ねています。

t = 0:0.001:0.6;
x = sin(2*pi*50*t)+sin(2*pi*120*t);
y = x + 2*randn(size(t));
plot(y(1:50))
title('Signal Corrupted with Zero-Mean Random Noise')
xlabel('time (seconds)')

オリジナル信号を調べて周波数成分を識別することは困難です。周波数領域に変換するには、ノイズを含む信号 y の離散フーリエ変換を、512個の降即フーリエ変換(FFT)から求めてください。

Y = fft(y,512);

種々の周波数でのパワーの尺度になるパワースペクトルは、つぎのように計算します。

Pyy = Y.* conj(Y) / 512;

周波数に対する大きさの応答を表す軸上に最初の257点をグラフ化します。

f = 1000*(0:256)/512;
plot(f,Pyy(1:257))
title('Frequency content of y')
xlabel('frequency (Hz)')

これは、DC 成分からナイキスト周波数までの y の周波数成分を表現します(信号は、強いピークを示します)。

アルゴリズム

FFT関数(fft, fft2, fftn, ifft, ifft2, ifftn) は、FFTW [3],[4]と呼ばれるライブラリをベースにしています。N が構成要素(that is, when N = N1N2)のとき、N 点の DFT を計算するには、FFTW ライブラリは、Cooley-Tukey アルゴリズム [1]を使って、問題を分割し、まず、サイズ N2 の N1 の変換を計算し、その後で、サイズ N1 の N2 の変換を計算を計算します。分割は、問題が、それぞれのマシンの固有の固定サイズ"codelets"の中の一つを使って解くことのできるまで、N1- と N2-点の DFT を共に繰り返し適用します。つぎつぎの繰り返しの codelets は、Cooley-Tukey の修正版[5]、素因数分解アルゴリズム [6]とsplit-radix アルゴリズム [2]を組み合わせたいくつかのアルゴリズムを使用します。N の特別な分解は、階層的に選択されます。

N が素数の場合、FFTW ライブラリは、まず、N 点問題をRader アルゴリズム [7]を使って、3つの(N-1) 点問題に分割します。上述した Cooley-Tukey 分割を使って、(N-1) 点のDFT を計算します。

N に対して、実数入力の DFTs は、複素数入力の DFT の計算時間が約半分になります。しかし、N が大きな素数の場合、スピードの違いはありません。

fft に対する実行時間は、変換を行う長さに依存します。長さが2のベキ乗の場合、最も高速です。また、小さな素数のみの組み合わせからなるものがつぎに高速になります。そして、大きい素数に分解される場合には、徐々に時間を要するようになります。

参考

dftmtx, filter, and freqz(Signal Processing Toolbox)

fft2, fftn, fftshift, ifft

参考文献

Cooley, J. W. and J. W. Tukey, "An アルゴリズム for the Machine Computation of the Complex Fourier Series," Mathematics of Computation, Vol. 19, April 1965, pp. 297-301.

Duhamel, P. and M. Vetterli, "Fast Fourier Transforms: A Tutorial Review and a State of the Art," Signal Processing, Vol. 19, April 1990, pp. 259-299.

FFTW (http://www.fftw.org)

Frigo, M. and S. G. Johnson, "FFTW: An Adaptive Software Architecture for the FFT," Proceedings of the International Conference on Acoustics, Speech, and Signal Processing, Vol. 3, 1998, pp. 1381-1384.

Oppenheim, A. V. and R. W. Schafer, Discrete-Time Signal Processing, Prentice-Hall, 1989, p. 611.

Oppenheim, A. V. and R. W. Schafer, Discrete-Time Signal Processing, Prentice-Hall, 1989, p. 619.

Rader, C. M., "Discrete Fourier Transforms when the Number of Data Samples Is Prime," Proceedings of the IEEE, Vol. 56, June 1968, pp. 1107-1108.

feval fft2