czt (Signal Processing Toolbox)

Chirp z-変換

表示

y = czt(x,m,w,a)
y = czt(x)

詳細

y = czt(x,m,w,a) は、信号xのChirp z変換を出力します。Chirp z変換は、wとaで定義される螺旋コンターに沿ったxのz変換です。mは変換の長さを設定するスカラで、wは対象となるz平面の螺旋コンターに沿った点の間の比率であり、また、スカラaは、そのコンター上の出発点を与える複素数です。z平面上のコンター、螺旋、または、"Chirp "は、次式によって与えられます。

z = a*(w.^-(0:m-1))

y = czt(x) は、つぎのデフォルト値を使用します。

m = length(x)

w = exp(j*2*pi/m)

a = 1

これらのデフォルト値を使用して、cztは、単位円周上の等間隔のm点でのxのz変換を出力します。これは、xの離散フーリエ変換、すなわち、fft(x)と等価です。空行列[]により、入力引数にデフォルト値を設定します。

x が行列の場合、czt(x,m,w,a)は、xの各列を変換します。

例題

長さ1013のランダムベクトルxを作成し、cztを用いて、そのDFTを計算します。これは、同じデータ列を使った関数fftより高速です。

randn('state',0);
x = randn(1013,1);
y = czt(x);

cztを使って、フィルタの周波数応答の狭帯域(100～150 Hz)を拡大します。まず、つぎのようにしてフィルタを設計します。

h = fir1(30,125/500,boxcar(31)); % フィルタ

周波数とCZTのパラメータを設定します。

fs = 1000; f1 = 100; f2 = 150;  % ヘルツ単位
m = 1024;
w = exp(-j*2*pi*(f2-f1)/(m*fs));
a = exp(j*2*pi*f1/fs);

フィルタのDFTとCZTを計算します。

y = fft(h,1000);
z = czt(h,m,w,a);

周波数ベクトルを作成し、結果を比較します。

fy = (0:length(y)-1)'*1000/length(y); 
fz = ((0:length(z)-1)'*(f2-f1)/length(z)) + f1;
plot(fy(1:500),abs(y(1:500))); axis([1 500 0 1.2])
title('FFT')
figure
plot(fz,abs(z)); axis([f1 f2 0 1.2])
title('CZT')

アルゴリズム

cztは、設定したChirp コンタ上でz変換を計算する際に、高速コンボリューションを実行するため、データをそれ以上大きい2のベキ乗数のうち最小のデータ数を使います(参考文献[1])。cztは、大きい素数長のデータ列の場合、fftよりはるかに高速に計算を行うことができます。

診断

m、w、または、aがスカラでない場合には、cztはつぎのエラーメッセージを出力します。

Inputs M, W, and A must be scalars.

参考

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

freqz
ディジタルフィルタの周波数応答

`fft`	1次元高速フーリエ変換
`freqz`	ディジタルフィルタの周波数応答

参考文献

[1] Rabiner, L.R., and B. Gold. Theory and Application of Digital Signal Processing. Englewood Cliffs, NJ: Prentice-Hall, 1975. Pgs. 393-399.

csd dct