Computer Programs in Biomedicine 7 (1977) 267-272 © Elsevier/North-Holland Biomedical Press
A PORTABLE INTEGER FFT IN FORTRAN Donald M. MONRO
Engineering in Medicine Laboratory, Imperial College, London SW7 2BT, England
A radix 2 integer fast Fourier transform is organised for compatibility among small computer systems by simulating binary fraction multiplications and bit reversal as functions in basic FORTRAN. After installation the efficiency can be improved by machine code replacement of some parts. Fast Fourier transformation
Portable
1. Introduction
is the transformed sequence. The inverse procedure is
In the acquisition and analysis of biological signals the usual minicomputer installation presents s o m e combination of restrictions in memory size and speed of floating point arithmetic. The Fourier series, in the form of a fast Fourier transform (FFT), must be the most widely used standard library facility but is certain to be inconvenienced by these limitations. Therefore an integer FFT is a logical complement to the standard repertoire of FFT programs although it is customarily regarded as highly machine dependent. It is demonstrated here that it can be generalised and implemented at a fundamental level of FORTRAN which renders it portable even among restricted FORTRAN compilers, In preparing an integer FFT it is usual to employ a certain amount of machine code, but in this program these parts are simulated in FORTRAN. The implementation has all the advantages of an integer FFT except that the speed is below the ultimate. However, simple replacement of these FORTRAN simulations will improve the speed greatly,
N- l
lk = ~
Jn ej 27rkn/N,
k = 0, 1 .... , N - 1 .
n =0 (2) The particular style of FFT adopted is a version of the Sande-Tukey algorithm decomposed as described by Gentleman and Sande [1 ]. A number of efficient implementations have been reported for floating point arithmetic; those of Monro [ 2 - 4 ] cover a range of applications with compatible programs. This FFT is arranged as radix 2, giving a compact program although the greatest speeds are normally achieved using higher or mixed radices. The coding avoids FORTRAN features which some compilers do not support. The method of circumventing the variable dimension conforms to common practice (but not the standard for FORTRAN) in assigning the dimension 1 to arrays IRE and IMG in the subroutine. The generation of sine and cosine values is done with the recursion and correction described by Singleton [5]. The calling sequence for the program is described in the comments in the listing.
2. Implementation The FFT is calculated according to the definition N- 1
1 ~ Ik e-i2rrkn/N
Jn=~k=O
3. Binary fractions and scaling '
n = O , 1 ..... N - 1
'
The integer arithmetic of FORTRAN can produce results which overflow. However, by treating the integers as fixed point fractions overflow on multiplication
(1) for a data array {Ik} o f N members. The array {Jn} 267
268
D.M. Monro, A portable integer FFT in FORTRAN ]i[
A PORTABLE INTEGER FFT IN FORTRAN Donald M. MONRO
Engineering in Medicine Laboratory, Imperial College, London SW7 2BT, England
A radix 2 integer fast Fourier transform is organised for compatibility among small computer systems by simulating binary fraction multiplications and bit reversal as functions in basic FORTRAN. After installation the efficiency can be improved by machine code replacement of some parts. Fast Fourier transformation
Portable
1. Introduction
is the transformed sequence. The inverse procedure is
In the acquisition and analysis of biological signals the usual minicomputer installation presents s o m e combination of restrictions in memory size and speed of floating point arithmetic. The Fourier series, in the form of a fast Fourier transform (FFT), must be the most widely used standard library facility but is certain to be inconvenienced by these limitations. Therefore an integer FFT is a logical complement to the standard repertoire of FFT programs although it is customarily regarded as highly machine dependent. It is demonstrated here that it can be generalised and implemented at a fundamental level of FORTRAN which renders it portable even among restricted FORTRAN compilers, In preparing an integer FFT it is usual to employ a certain amount of machine code, but in this program these parts are simulated in FORTRAN. The implementation has all the advantages of an integer FFT except that the speed is below the ultimate. However, simple replacement of these FORTRAN simulations will improve the speed greatly,
N- l
lk = ~
Jn ej 27rkn/N,
k = 0, 1 .... , N - 1 .
n =0 (2) The particular style of FFT adopted is a version of the Sande-Tukey algorithm decomposed as described by Gentleman and Sande [1 ]. A number of efficient implementations have been reported for floating point arithmetic; those of Monro [ 2 - 4 ] cover a range of applications with compatible programs. This FFT is arranged as radix 2, giving a compact program although the greatest speeds are normally achieved using higher or mixed radices. The coding avoids FORTRAN features which some compilers do not support. The method of circumventing the variable dimension conforms to common practice (but not the standard for FORTRAN) in assigning the dimension 1 to arrays IRE and IMG in the subroutine. The generation of sine and cosine values is done with the recursion and correction described by Singleton [5]. The calling sequence for the program is described in the comments in the listing.
2. Implementation The FFT is calculated according to the definition N- 1
1 ~ Ik e-i2rrkn/N
Jn=~k=O
3. Binary fractions and scaling '
n = O , 1 ..... N - 1
'
The integer arithmetic of FORTRAN can produce results which overflow. However, by treating the integers as fixed point fractions overflow on multiplication
(1) for a data array {Ik} o f N members. The array {Jn} 267
268
D.M. Monro, A portable integer FFT in FORTRAN ]i[
