Misplaced Pages

Chirp Z-transform

Article snapshot taken from Wikipedia with creative commons attribution-sharealike license. Give it a read and then ask your questions in the chat. We can research this topic together.
(Redirected from Bluestein's FFT algorithm)
This article needs additional citations for verification. Please help improve this article by adding citations to reliable sources. Unsourced material may be challenged and removed.
Find sources: "Chirp Z-transform" – news · newspapers · books · scholar · JSTOR (January 2013) (Learn how and when to remove this message)

The chirp Z-transform (CZT) is a generalization of the discrete Fourier transform (DFT). While the DFT samples the Z plane at uniformly-spaced points along the unit circle, the chirp Z-transform samples along spiral arcs in the Z-plane, corresponding to straight lines in the S plane. The DFT, real DFT, and zoom DFT can be calculated as special cases of the CZT.

Specifically, the chirp Z transform calculates the Z transform at a finite number of points zk along a logarithmic spiral contour, defined as:

X k = ∑ n = 0 N − 1 x ( n ) z k − n {\displaystyle X_{k}=\sum _{n=0}^{N-1}x(n)z_{k}^{-n}}
z k = A ⋅ W − k , k = 0 , 1 , … , M − 1 {\displaystyle z_{k}=A\cdot W^{-k},k=0,1,\dots ,M-1}

where A is the complex starting point, W is the complex ratio between points, and M is the number of points to calculate.

Like the DFT, the chirp Z-transform can be computed in O(n log n) operations where n = max ( M , N ) n=\max(M,N) . An O(N log N) algorithm for the inverse chirp Z-transform (ICZT) was described in 2003, and in 2019.

Bluestein's algorithm

Bluestein's algorithm expresses the CZT as a convolution and implements it efficiently using FFT/IFFT.

As the DFT is a special case of the CZT, this allows the efficient calculation of discrete Fourier transform (DFT) of arbitrary sizes, including prime sizes. (The other algorithm for FFTs of prime sizes, Rader's algorithm, also works by rewriting the DFT as a convolution.) It was conceived in 1968 by Leo Bluestein. Bluestein's algorithm can be used to compute more general transforms than the DFT, based on the (unilateral) z-transform (Rabiner et al., 1969).

Recall that the DFT is defined by the formula

X k = ∑ n = 0 N − 1 x n e − 2 π i N n k k = 0 , … , N − 1. {\displaystyle X_{k}=\sum _{n=0}^{N-1}x_{n}e^{-{\frac {2\pi i}{N}}nk}\qquad k=0,\dots ,N-1.}

If we replace the product nk in the exponent by the identity

n k = − ( k − n ) 2 2 + n 2 2 + k 2 2 {\displaystyle nk={\frac {-(k-n)^{2}}{2}}+{\frac {n^{2}}{2}}+{\frac {k^{2}}{2}}}

we thus obtain:

X k = e − π i N k 2 ∑ n = 0 N − 1 ( x n e − π i N n 2 ) e π i N ( k − n ) 2 k = 0 , … , N − 1. {\displaystyle X_{k}=e^{-{\frac {\pi i}{N}}k^{2}}\sum _{n=0}^{N-1}\left(x_{n}e^{-{\frac {\pi i}{N}}n^{2}}\right)e^{{\frac {\pi i}{N}}(k-n)^{2}}\qquad k=0,\dots ,N-1.}

This summation is precisely a convolution of the two sequences an and bn defined by:

a n = x n e − π i N n 2 {\displaystyle a_{n}=x_{n}e^{-{\frac {\pi i}{N}}n^{2}}}
b n = e π i N n 2 , {\displaystyle b_{n}=e^{{\frac {\pi i}{N}}n^{2}},}

with the output of the convolution multiplied by N phase factors bk. That is:

X k = b k ∗ ( ∑ n = 0 N − 1 a n b k − n ) k = 0 , … , N − 1. {\displaystyle X_{k}=b_{k}^{*}\left(\sum _{n=0}^{N-1}a_{n}b_{k-n}\right)\qquad k=0,\dots ,N-1.}

This convolution, in turn, can be performed with a pair of FFTs (plus the pre-computed FFT of complex chirp bn) via the convolution theorem. The key point is that these FFTs are not of the same length N: such a convolution can be computed exactly from FFTs only by zero-padding it to a length greater than or equal to 2N–1. In particular, one can pad to a power of two or some other highly composite size, for which the FFT can be efficiently performed by e.g. the Cooley–Tukey algorithm in O(N log N) time. Thus, Bluestein's algorithm provides an O(N log N) way to compute prime-size DFTs, albeit several times slower than the Cooley–Tukey algorithm for composite sizes.

The use of zero-padding for the convolution in Bluestein's algorithm deserves some additional comment. Suppose we zero-pad to a length M ≥ 2N–1. This means that an is extended to an array An of length M, where An = an for 0 ≤ n < N and An = 0 otherwise—the usual meaning of "zero-padding". However, because of the bk–n term in the convolution, both positive and negative values of n are required for bn (noting that b–n = bn). The periodic boundaries implied by the DFT of the zero-padded array mean that –n is equivalent to M–n. Thus, bn is extended to an array Bn of length M, where B0 = b0, Bn = BM–n = bn for 0 < n < N, and Bn = 0 otherwise. A and B are then FFTed, multiplied pointwise, and inverse FFTed to obtain the convolution of a and b, according to the usual convolution theorem.

Let us also be more precise about what type of convolution is required in Bluestein's algorithm for the DFT. If the sequence bn were periodic in n with period N, then it would be a cyclic convolution of length N, and the zero-padding would be for computational convenience only. However, this is not generally the case:

b n + N = e π i N ( n + N ) 2 = b n [ e π i N ( 2 N n + N 2 ) ] = ( − 1 ) N b n . {\displaystyle b_{n+N}=e^{{\frac {\pi i}{N}}(n+N)^{2}}=b_{n}\left=(-1)^{N}b_{n}.}

Therefore, for N even the convolution is cyclic, but in this case N is composite and one would normally use a more efficient FFT algorithm such as Cooley–Tukey. For N odd, however, then bn is antiperiodic and we technically have a negacyclic convolution of length N. Such distinctions disappear when one zero-pads an to a length of at least 2N−1 as described above, however. It is perhaps easiest, therefore, to think of it as a subset of the outputs of a simple linear convolution (i.e. no conceptual "extensions" of the data, periodic or otherwise).

z-transforms

Bluestein's algorithm can also be used to compute a more general transform based on the (unilateral) z-transform (Rabiner et al., 1969). In particular, it can compute any transform of the form:

X k = ∑ n = 0 N − 1 x n z n k k = 0 , … , M − 1 , {\displaystyle X_{k}=\sum _{n=0}^{N-1}x_{n}z^{nk}\qquad k=0,\dots ,M-1,}

for an arbitrary complex number z and for differing numbers N and M of inputs and outputs. Given Bluestein's algorithm, such a transform can be used, for example, to obtain a more finely spaced interpolation of some portion of the spectrum (although the frequency resolution is still limited by the total sampling time, similar to a Zoom FFT), enhance arbitrary poles in transfer-function analyses, etc.

The algorithm was dubbed the chirp z-transform algorithm because, for the Fourier-transform case (|z| = 1), the sequence bn from above is a complex sinusoid of linearly increasing frequency, which is called a (linear) chirp in radar systems.

See also

References

  1. ^ A study of the Chirp Z-transform and its applications - Shilling, Steve Alan
  2. "Chirp Z-transform - MATLAB czt". www.mathworks.com. Retrieved 2016-09-22.
  3. Martin, Grant D. (November 2005). "Chirp Z-Transform Spectral Zoom Optimization with MATLAB®" (PDF).
  4. Bostan, Alin (2003). Algorithmique efficace pour des operations de base en Calcul formel (PDF) (PhD). Ecole polytechnique.
  5. Bostan, Alin; Schost, Éric (2005). "Polynomial evaluation and interpolation on special sets of points". Journal of Complexity. 21 (4): 420–446. doi:10.1016/j.jco.2004.09.009.
  6. Engineers Solve 50-Year-Old Puzzle in Signal Processing – Inverse Chirp Z-Transform, By IOWA STATE UNIVERSITY OCTOBER 10, 2019
  7. ^ Bluestein, L. (1970-12-01). "A linear filtering approach to the computation of discrete Fourier transform". IEEE Transactions on Audio and Electroacoustics. 18 (4): 451–455. doi:10.1109/TAU.1970.1162132. ISSN 0018-9278.
  8. "Bluestein's FFT Algorithm". DSPRelated.com.

General

External links

Categories: