Richard C. Singleton

An algorithm for computing the mixed radix fast Fourier transform

This paper presents an algorithm for computing the fast Fourier transform, based on a method proposed by Cooley and Tukey. As in their algorithm, the dimension n of the transform is factored (if possible), and n/p elementary transforms of dimension p are computed for each factor p of n . An improved method of computing a transform step corresponding to an odd factor of n is given; with this method, the number of complex multiplications for an elementary transform of dimension p is reduced from (p-1)^{2} to (p-1)^{2}/4 for odd p . The fast Fourier transform, when computed in place, requires a final permutation step to arrange the results in normal order. This algorithm includes an efficient method for permuting the results in place. The algorithm is described mathematically and illustrated by a FORTRAN subroutine.

On computing the fast Fourier transform

and have shown major time savings in using it to compute large transforms on a digital computer. With n a power of two, computing time for this algorithm is proportional to n log2 n, a major improvement over other methods with computing time proportional to n 2. In this paper, the fast Fourier transform algorithm is briefly reviewed and fast difference equation methods for accurately computing the needed trigonometric function values are given. The problem of computing a large Fourier transform on a system with virtual memory is considered, and a solution is proposed. This method has been used to compute complex Fourier transforms of size n = 2 z6 on a computer with 215 words of core storage; this exceeds by a factor of eight the maximum radix two transform size with fixed allocation of this amount of core storage. The method has also been used to compute large mixed radix transforms. A scaling plan for computing the fast Fourier transform with fixed-point arithmetic is also given.

A method for computing the fast Fourier transform with auxiliary memory and limited high-speed storage

A method is given for computing the fast Fourier transform of arbitrarily large size using auxiliary memory files, such as magnetic tape or disk, for data storage. Four data files are used, two in and two out. A multivariate complex Fourier transform ofn=2^{m}data points is computed inmpasses of the data, and the transformed result is permuted to normal order bym - 1additional passes. With buffered input-output, computing can be overlapped with reading and writing of data. Computing time is proportional ton \log_{2} n. The method can be used with as few as three files, but file passing for permutation is reduced by using six or eight files. With eight files, the optimum number for a radix 2 transform, the transform is computed inmpasses without need for additional permutation passes. An ALGOL procedure for computing the complex Fourier transform with four, six, or eight files is listed, and timing and accuracy test results are given. This procedure allows an arbitrary number of variables, each dimension a power of 2.

Spectral analysis of the call of the male killer whale

Poulter previously presented sonagram analyses of underwater sound recordings of various sea mammals, calling attention to the apparent harmonic progression of components of these signals. Other workers have questioned the presence of these harmonics in the original data and have suggested that they may instead be introduced during the analysis. In this paper we use underwater sound recordings of the male killer whale Namu, taken with equipment flat within ±2 dB to beyond 18 342 Hz, the upper limit of the analysis, to show by means of digital spectral analysis methods that a harmonic progression exists. In the signals analyzed, peaks in the estimated spectrum were observed at each integer multiple of the fundamental frequency. In view of the broad frequency response of the recording equipment and the precision of the subsequent digital analysis, we can say with confidence that these harmonics were actually present in the whalers call. For the analysis, portions of an audio recording were converted to digital form. In the conversion, a bandpass filter was used to attenuate power below 40 Hz and above 10 kHz. Digital analysis techniques similar to those proposed by Bingham, Godfrey, and Tukey were then used. For each time span of data to be analyzed, a windowed Fourier transform was first computed, using a fast Fourier transform program. The power spectrum was next computed, as the squared modulus of the windowed transform, and a correction was made for the attenuation of high frequencies during A-to-D conversion. The autocorrelation function was estimated by computing the inverse transform of the power spectrum. A moderate amount of digital smoothing was then applied to the spectrum to reduce irregularities due to noise. The resulting smoothed spectrum is used as an estimate of the power spectral density function.

Generalized Snake-in-the-Box Codes

A snake-in-the-box (SIB) code of order k is defined to be an ordered sequence of binary code words in which adjacent words differ in only one bit, and pairs of code words that are k or more apart in the ordered sequence differ in at least k bit positions. In this paper, constructions for SIB codes of arbitrary order are given, as well as upper bounds on the maximum code sequence length for given order and word size. These codes are potentially useful for binary encoding of analog data. Gray codes are SIB codes of order one, and Kautz has investigated SIB codes of order two. The SIB codes of a given order contain as a subset all SIB codes of higher order.

Algorithms: Algorithm 338: algol procedures for the fast Fourier transform

The following procedures are based on the Cooley-Tukey algorithm [1] for computing the finite Fourier transform of a complex data vector; the dimension of the data vector is assumed here to be a power of two. Procedure <italic>COMPLEXTRANSFORM</italic> computes either the complex Fourier transform or its inverse. Procedure <italic>REALTRANSFORM</italic> computes either the Fourier coefficients of a sequence of real data points or evaluates a Fourier series with given cosine and sine coefficients. The number of arithmetic operations for either procedure is proportional to <italic>n</italic> log<subscrpt>2</subscrpt> <italic>n</italic>, where <italic>n</italic> is the number of data points.

Algorithms: Algorithm 339: an algol procedure for the fast Fourier transform with arbitrary factors

The following procedures are based on the Cooley-Tukey algor i thm [1] for comput ing the finite Fourier t r ans fo rm of a complex da ta vector; the dimension of the da t a vector is assumed here to be a power of two. Procedure COMPLEXTRANSFORM computes ei ther the complex Fourier t ransform or its inverse. Procedure REALTRANSFORM computes ei ther the Fourier coefficients of a sequence of real da ta points or eva lua tes a Fourier series wi th given cosine and sine coefficients. The number of ar i thmet ic operat ions for ei ther procedure is proport ional to n logs n, where n is the number of da ta points. Procedures FFT2, REVFFT2, REORDER, and REAL TRAN are building blocks, and are used in the two complete procedures ment ioned above. The fas t t r ans fo rm can be computed in a number of different ways, and these bui lding block procedures were wri t ten so as to make practical the comput ing of large t rans forms on a system wi th vir tual memory. Using a method proposed by Singleton [2], d a t a is accessed in sub-sequences of consecutive a r ray elements , and as m u ch comput ing as possible is done in one section of the d a t a before moving on to another. Procedure FFT2 computes the Fourier t r ans fo rm of da ta in normal order, giving a resu l t in reverse b ina ry order. Procedure REVFFT2 computes the Fourier t r ans fo rm of da ta in reverse b ina ry order and leaves the resul t in normal b inary order. Procedure REORDER permutes a complex vector f rom b inary to reverse b ina ry order or f rom reverse b inary to b inary order; this procedure also permutes real da ta in prepara t ion for efficient use of the complex Fourier t ransform. Procedures FFT2, REVFFT2, and REORDER m a y also be used to compute mul t iva r i a te Fourier t ransforms. The procedure R E A L T R A N is used to unscramble and combine the complex t rans forms of the even and odd numbered e lements of a sequence of real d a t a points . This procedure is not restr ic ted to powers of two and can be used whenever the number of da t a points is even.

A short bibliography on the fast Fourier transform

166 L. E. Alsop and A. A. Nowroozi, “Fast Fourier analysis,” J. Geophys. Res., vol. 71, pp. 5482-5483, November 15, 1966. €3. Andrews, “A high-speed algorithm for the computer generation of Fourier transforms,” IEEE Trans. Computers (Short Notes), vol. C-17, pp. 373.375, April 1968. J. S . Bailey, “A fast Fourier transform without multiplications,” Proc. Symp. on Computer Processing in Communications, vol. 19, MKI Symposia Ser. New York: Polytechnic Press, 1969. V. Benignus, “Estimation of the coherence spectrum and its confidence interval using the fast Fourier transform,” this issue, pp. 145-150. G. D. Bergland, “The fast Fourier transform recursive equations for arbitrary length records,” Math. Computation, vol. 21, pp, 236-238, April 1967. -9 “A fast Fourier transform algorithm using base eight iterations,” Math. Computation, vol. 22, pp. 275-279, April 1968. -, “A fast Fourier transform algorithm for realvalued series,” Commun. A C M , vol. 11, pp. 703--710, October 1968. -, “A radix-eight fast Fourier transform subroutine for real-valued series,” this issue, pp. 138144. -, “A guided tour of the fast Fourier transform,” IEEE Spectrum (to be published). “Fast Fourier transform hardware implementations. I. An overview. 11. A survey,’’ this issue,

Load-Sharing Core Switches Based on Block Designs

Designs for load-sharing zero-noise core switches have been proposed by Constantine, Marcus, and Chien. Blachman class has proposed a core memory wiring plan which with modification can be converted to a load-sharing zero-noise switch. An examination of these switch plans shows that they have a common relationship to a class of mathematical structures known to mathematicians and statisticians as balanced incomplete block designs. This relationship is formulated, and it is then shown that all balanced incomplete block designs lead to load-sharing zero-noise switches. Three methods of forming the winding matrix for a switch are given, and expressions for the load-sharing factor, set bias, and reset bias in terms of the balanced incomplete block design parameters are derived for each switch type. Similarly, partially balanced incomplete block designs are shown to lead to low-noise load-sharing switches. Switch operation under fault conditions is briefly discussed. Most of the known load-sharing core switch types can be viewed as based on either balanced or partially balanced incomplete block designs. A review of the available block designs indicates that a number of new switches can be based on these designs. A modification of a distributed memory model proposed by C. Rosen is discussed. With wiring plans based on block designs, it appears possible to construct very-large-capacity memory units which are relatively insensitive to wiring errors.

Real-time spectral analysis on a small general-purpose computer

There is a growing need for methods of quickly estimating the changing frequency content of a nonstationary signal. In this paper we describe a method for doing spectral analysis in real time on a small general-purpose digital computer, and discuss some of the theoretical and practical problems of developing similar systems for other computers.