1. Home
2. FOURIER directive

# FOURIER directive

Calculates cosine or Fourier transforms of real or complex series.

### Option

`PRINT` = string tokens What to print (`transforms`); default `*`

### Parameters

`SERIES` = variates Real part of each input series Imaginary part of each input series To save real part of each output series To save imaginary part of each output series To save periodogram of each transform

### Description

The Fourier or spectral analysis of time series is described comprehensively by Bloomfield (1976) and Jenkins & Watts (1968). The Fourier transformation of a series calculates the coefficients of the sinusoidal components into which the series can be analysed. There are four types of transformation described below, which are appropriate for different types of symmetry in the series. You may often want the length of the variate holding the supplied series to determine implicitly a natural grid of frequencies at which values of the transform are calculated. Genstat will do this if you have not previously declared the identifier supplied for the transform. Alternatively you may want to determine the transform at a finer grid of frequencies, and you can achieve this by declaring a transform variate that is as long as you require. You can do this only for the two types of Fourier transform that apply to real series.

Series of real numbers are stored in single variates, and series of complex numbers in pairs of variates. You can use the `FOURIER` directive to calculate the cosine transform of the real series { at , t=0…N-1 } stored in a variate `A` by

`FOURIER [PRINT=transform] A`

You calculate the Fourier transform of the complex series { at+ibt , t=0…N-1 } by storing the values at in one variate, `A` say, the corresponding values bt in another, `B` say, and giving the statement:

`FOURIER [PRINT=transform] A; ISERIES=B`

You can restrict the series specified by either the `SERIES` or `ISERIES` parameter to a contiguous set of units. Genstat then applies the transformation only to the restricted series of values. Similarly, you may supply restricted variates with the `TRANSFORM` and `ITRANSFORM` parameters to save the transform: Genstat will then carry out the transformation so as to supply the required number of values. There must be no missing values in the variates in the `SERIES` or `ISERIES` parameters, unless you exclude them by a restriction.

Genstat carries out the Fourier transformation using a fast algorithm which relies on the order of the transformation being highly composite (de Boor 1980). In practice, an appropriate order is a round number such as 300 or 6000, consisting of a digit followed by zeroes. If, however, the order has a large prime factor, the transformation may take much longer. For example, a transformation of order 499 is about 25 times slower than one of order 500. In the descriptions below, therefore, we clearly state the order of each form of the transformation, to illustrate a sensible choice of size.

The cosine transformation of a real series can be used to calculate the spectrum from a set of autocorrelations. Suppose the variate `R` contains the values r0rn, and the variate `F` is to hold the calculated values f0fm of the spectrum. These values correspond to angular frequencies of πj/m; that is, periods of 2m/j, for j=0…m. You apply the transformation by putting

`FOURIER R; TRANSFORM=F`

If `F` has not been declared previously, this statement defines it automatically as a variate with n+1 values (so m=n). If `F` has been declared to have m+1 values, then m must be greater than or equal to n; otherwise Genstat will redeclare `F` to have n+1 values.

The transform is defined when m>n by

fj = r0 + ∑k = 1 … n { 2 × rk × cos(k × π × j / m) }

When m=n the final term in this sum is

rn cos(πj) = rn (-1)j

and it appears without the multiplier 2. The order of the transformation is 2m.

If `R` contains sample autocorrelations, you must multiply it by a variate holding a lag window in order to obtain a smooth spectrum estimate (see Bloomfield 1976, page 166; or Jenkins & Watts 1968, page 243).

The Fourier transformation of a real series can be used to calculate the periodogram of a time series. Suppose the variate `X` of length N contains the supplied series values x0xN1 . The result of the transformation is a set of coefficients a0am of the cosine components and b0bm of the sine components of the series, held in variates `A` and `B`, say. Normally the number of such components is related to the length of the series by taking m=N/2 if N is even or m=(N-1)/2 if N is odd. Then the coefficients correspond to angular frequencies of 2πj/N, which is the same as saying that they correspond to periods N/j for j=0…m. Since by definition b0=0, and bm=0 if N is even, there are N “free” coefficients in `A` and `B` (which you can think of as the real and imaginary parts of a complex transform with values aj+ibj). You can save the periodogram values p0pm in a variate `P`, say: these are the squared amplitudes of the sinusoidal components, and are calculated by Genstat as pj = aj2+bj2.

You obtain the transform by putting

`FOURIER X; TRANSFORM=A; ITRANSFORM=B; PERIODOGRAM=P`

If you want only the periodogram, you can put

`FOURIER X; PERIODOGRAM=P`

If you have not declared `A` previously Genstat defines it automatically, here as a variate of length m+1 where m has the default value defined above. If you have previously declared `A`, it should have length greater than or equal to m+1; otherwise Genstat declares it to have this length. In any case, `B` and `P` should have the same length as `A`, and will be declared (or redeclared) if required.

In the usual case when `A`, `B` or `P` has the default length m+1, the transform is defined by:

aj = ∑t = 0 … N-1 { xt × cos(t × 2π × j / N) } ; j = 0 … m

bj = ∑t = 0 … N-1 { xt × sin(t × 2π × j / N) } ; j = 0 … m

In this case, the order of the transformation is N. If `A`, `B` and `P` have length m′+1 with m′>m, Genstat computes the results at a finer grid of frequencies 2πj/N′, j=0…m′ where N′=2m′. These replace 2πj/N in the above defining sums. The upper limit on the sums remains as N-1, although internally Genstat treats it as N′-1 with the extra values of xNxN′-1 being taken as zero. The order of the transformation is then N′. There are various conventions used for scaling the periodogram with factors 2/m, 1/m or 1/πm. You can apply these by using a `CALCULATE` statement after the transformation. You may also want to apply mean correction to the series before calculating the periodogram.

The Fourier transformation of a complex series is the most general form of the Fourier transformation; the other three types are essentially special cases in which some coefficients are zero or have a symmetric structure. Suppose variates `X` and `Y` contain values x0xN1 and y0yN1, which may be viewed as the real and imaginary parts of the series { xt+iyt, t=0 … N-1 }. The results of the transformation are coefficients a0aN1 and b0bN1 which can be held in variates `A` and `B`, say: these may similarly be considered as parts of complex coefficients at+ibt, t=0 … N-1.

You can do the transformation by putting

`FOURIER SERIES=X; ISERIES=Y; TRANSFORM=A; ITRANSFORM=B`

Both `X` and `Y` must be variates with the same length N. Similarly `A` and `B` must have length N, and if they do not Genstat will declare (or redeclare) them as variates of length N. The order of the transformation is N.

The results are defined by

aj = ∑t = 0 … N-1 { xt × cos(t × 2π × j / N) – yt × sin(t × 2π × j / N) } ; j = 0 … m

bj = ∑t = 0 … N-1 { xt × sin(t × 2π × j / N) + yt × cos(t × 2π × j / N) } ; j = 0 … m

or equivalently in complex form by

(aj + i bj) = ∑t = 0 … N-1 { (xt + i yt) × exp(i t × 2π × j / N) }

The complex transform can be used in cross-spectral analysis.

You can view a Fourier transformation as an orthogonal matrix transformation. Hence its inverse is another Fourier transformation (apart from some simple scaling). You can use this to calculate convolutions. In particular, the correlations of a time series can be obtained by applying the inverse cosine transformation to the periodogram.

The Fourier transform of a conjugate sequence is most easily considered as the reverse of the transformation of a real series, with the roles of the series and the transform interchanged. For the true inverse transformation some simple scaling is also required.

Thus if variates `A` and `B` of length m+1 are supplied containing values a0am and b0bm, which may be viewed as parts of complex coefficients aj+ibj, the result of the transformation is a single real series x0xN1 held in a variate `X` of length N.

`X` can be declared to have length N=2m or N=2m+1 (corresponding to the case N even or odd in the Fourier transformation of a real series). The value of b0 must be zero; also if N=2m, the value of bm must be zero. If either of these conditions is not satisfied, Genstat sets the values of these elements to zero and gives a warning. If `X` has not been declared previously (or has been declared with a length equal to neither 2m nor 2m+1), then it is declared (or redeclared) with a length governed by whether bm is 0: N=2m if bm=0, and N=2m+1 if bm≠0. The value of b0 is checked to be zero as before.

You can obtain the transform using the statement

`FOURIER SERIES=A; ISERIES=B; TRANSFORM=X`

The definition of the transform is, in the case N=2m+1,

xt = a0 + ∑j = 1 … m {2 × (aj × cos(t × 2π × j / N) – bj × sin(t × 2π × j / N)) }

In the case N=2m, the final term in the sum is simply

am cos(tπ) = am (-1)t

and it appears without the multiplier 2. The order of this transformation is N.

Option: `PRINT`.

Parameters: `SERIES`, `ISERIES`, `TRANSFORM`, `ITRANSFORM`, `PERIODOGRAM`.

### Action with `RESTRICT`

You can restrict the series specified by either the `SERIES` or `ISERIES` parameter to a contiguous set of units. Genstat then applies the transformation only to the restricted series of values. Similarly, you may supply restricted variates with the `TRANSFORM` and `ITRANSFORM` parameters to save the transform: Genstat will then carry out the transformation so as to supply the required number of values.

### References

Bloomfield, P. (1976). Fourier Analysis of Time Series: an Introduction. Wiley, New York.

de Boor, C. (1980). FFT as nested multiplication, with a twist. SIAM Journal of Scientific and Statistical Computing, 1, 173-178.

Jenkins, G.M. & Watts, D.G. (1968). Spectral Analysis and its Applications. Holden-Day, San Francisco.

Directive: `CORRELATE`.

Procedures: `DFOURIER`, `MCROSSPECTRUM`, `PERIODTEST`, `REPPERIODOGRAM`, `SMOOTHSPECTRUM`.

Commands for: Time series.

### Example

```" Example FOUR-1: Cosine transformation
Generate a random sequence, then transform it."
CALCULATE Random = URAND(15192; 35)
FOURIER SERIES=Random; TRANSFORM=Cosine1

" Divide the result by 2*(NVAL(Random)-1) to ensure
that retransformed data matches Random."
CALCULATE Cosine1 = Cosine1/2/(NVAL(Random)-1)
FOURIER SERIES=Cosine1; TRANSFORM=Cosine2
PRINT Random,Cosine1,Cosine2
```
Updated on March 8, 2019