## Maxima – 130 – Funzioni numeriche – 2 Continuo da qui, copio dal Reference Manual, PDF scaricabile da qui, sono a p.384.

`fft (x)`
Computes the complex fast Fourier transform. `x` is a list or array (named or unnamed) which contains the data to transform. The number of elements must be a power of 2. The elements must be literal numbers (integers, rationals, floats, or bigfloats) or symbolic constants, or expressions `a + b*%i` where `a` and `b` are literal numbers or symbolic constants.

`fft` returns a new object of the same type as `x`, which is not modified. Results are always computed as floats or expressions `a + b*%i` where `a` and `b` are floats. If bigfloat precision is needed the function `bf_fft` can be used instead as a drop-in replacement of `fft` that is slower, but supports bfloats. In addition if it is known that the input consists of only real values (no imaginary parts), `real_fft` can be used which is potentially faster.

The discrete Fourier transform is defined as follows. Let `y` be the output of the transform. Then for `k` from `0` through `n - 1`,
`y[k] = (1/n) sum(x[j] exp(+2 %i %pi j k / n), j, 0, n - 1)`

As there are various sign and normalization conventions possible, this definition of the transform may differ from that used by other mathematical software.

When the data `x` are real, real coefficients `a` and `b` can be computed such that
`x[j] = sum(a[k]*cos(2*%pi*j*k/n)+b[k]*sin(2*%pi*j*k/n), k, 0, n/2)`
with
`a = realpart (y)`
`b = 0`
and, for `k` from `1` through `n/2 - 1`,
`a[k] = realpart (y[k] + y[n - k])`
`b[k] = imagpart (y[n - k] - y[k])`
and
`a[n/2] = realpart (y[n/2])`
`b[n/2] = 0`

`load(fft)` loads this function.

See also `inverse_fft` (inverse transform), `recttopolar`, and `polartorect`. See `real_fft` for FFTs of a real-valued input, and `bf_fft` and `bf_real_fft` for operations on bigfloat values.

``````(%i1) load ("fft") \$

(%i2) fpprintprec : 4 \$

(%i3) L : [1, 2, 3, 4, -1, -2, -3, -4] \$

(%i4) L1 : fft (L);
(%o4) [0.0, 1.81 %i - 0.1035, 0.0, 0.3106 %i + 0.6035, 0.0,
0.6035 - 0.3106 %i, 0.0, (- 1.81 %i) - 0.1035]
(%i5) L2 : inverse_fft (L1);
(%o5) [1.0, 4.44E-16 %i + 2.0, 3.0, 4.0 - 4.44E-16 %i, - 1.0,
(- 4.44E-16 %i) - 2.0, - 3.0, 4.44E-16 %i - 4.0]
(%i6) lmax (abs (L2 - L));
(%o6)                              4.44E-16``````

Complex data.

``````(%i1) load ("fft") \$

(%i2) fpprintprec : 4 \$

(%i3) L : [1, 1 + %i, 1 - %i, -1, -1, 1 - %i, 1 + %i, 1] \$

(%i4) L1 : fft (L);
(%o4) [0.5, 0.5, 0.25 %i - 0.25, (- 0.3535 %i) - 0.3535, 0.0, 0.5,
(- 0.25 %i) - 0.25, 0.3535 %i + 0.3535]
(%i5) L2 : inverse_fft (L1);
(%o5) [1.0, 1.0 %i + 1.0, 1.0 - 1.0 %i, - 1.0, - 1.0, 1.0 - 1.0 %i,
1.0 %i + 1.0, 1.0]
(%i6) lmax (abs (L2 - L));
(%o6)                                 0.0``````

Computation of sine and cosine coefficients.

``````(%i1) load ("fft") \$

(%i2) fpprintprec : 4 \$

(%i3) L : [1, 2, 3, 4, 5, 6, 7, 8] \$

(%i4) n : length (L) \$

(%i5) x : make_array (any, n) \$

(%i6) fillarray (x, L) \$

(%i7) y : fft (x) \$

(%i8) a : make_array (any, n/2 + 1) \$

(%i9) b : make_array (any, n/2 + 1) \$

(%i10) a : realpart (y) \$

(%i11) b : 0 \$

(%i12) for k : 1 thru n/2 - 1 do
(a[k] : realpart (y[k] + y[n - k]),
b[k] : imagpart (y[n - k] - y[k]));
(%o12)                               done
(%i13) a[n/2] : y[n/2] \$

(%i14) b[n/2] : 0 \$

(%i15) listarray (a);
(%o15)                 [4.5, - 1.0, - 1.0, - 1.0, - 0.5]
(%i16) listarray (b);
(%o16)                    [0, 2.414, 1.0, 0.4142, 0]
(%i17) f(j) := sum (a[k]*cos(2*%pi*j*k/n) + b[k]*sin(2*%pi*j*k/n), k, 0, n/2) \$

(%i18) makelist (float (f (j)), j, 0, n - 1);
(%o18)             [1.0, 8.0, 7.0, 6.0, 5.0, 4.0, 3.0, 2.0]``````

`real_fft (x)`
Computes the fast Fourier transform of a real-valued sequence `x`. This is equivalent to performing `fft(x)`, except that only the first `N/2+1` results are returned, where `N` is the length of `x`. `N` must be power of two.

No check is made that `x` contains only real values.

The symmetry properites of the Fourier transform of real sequences to reduce he complexity. In particular the first and last output values of real_fft are purely real.

For larger sequencyes, `real_fft` may be computed more quickly than `fft`.

Since the output length is short, the normal `inverse_fft` cannot be directly used. Use `inverse_real_fft` to compute the inverse.

`inverse_real_fft (y)`
Computes the inverse Fourier transfrom of `y`, which must have a length of `N/2+1` where `N` is a power of two. That is, the input `x` is expected to be the output of `real_fft`.

No check is made to ensure that the input has the correct format. (The first and last
elements must be purely real.)

`bf_inverse_fft (y)`
Computes the inverse complex fast Fourier transform. This is the bigfloat version of `inverse_fft` that converts the input to bigfloats and returns a bigfloat result.

`bf_fft (y)`
Computes the forward complex fast Fourier transform. This is the bigfloat version of `fft` that converts the input to bigfloats and returns a bigfloat result.

`bf_real_fft (x)`
Computes the forward fast Fourier transform of a real-valued input returning a bigfloat result. This is the bigfloat version of `real_fft`.

`bf_inverse_real_fft (y)`
Computes the inverse fast Fourier transfrom with a real-valued bigfloat output. This is the bigfloat version of `inverse_real_fft`.

Posta un commento o usa questo indirizzo per il trackback.