Fourier series coefficients

Calculate the Fourier coefficients of the series expansion of a function, and the amplitude and phase spectra. The script contains some theory and 3 different methods to calculate the coefficients.

USAGE
fourier_coeff(fun,t0,T)
fourier_coeff(fun,t0,T,M)
fourier_coeff(fun,t0,T,M,N)
fourier_coeff(fun,t0,T,M,N,method)
fourier_coeff(fun,t0,T,M,N,method,res)
fourier_coeff(fun,t0,T,M,N,method,res,num_P)

INPUT
- FUN: character string representing the function with "t" as the independent variable (e.g. '10*cos(2*pi*3*t-pi/4)'). Dot-arithmetic operators must be used (.* ./ .\ .^). FUN must be defined for [T0,T0+T]
- T0: initial "t" for the definition of FUN
- T: period of the function
- M: number of frequencies (default: 5)
- N: number of data points per period (default: 100)
- METHOD: 1 (least-squares), 2 (integrals [default]) or 3 (FFT)
- RES: 1 (plot results) or 0 (do not plot [default])
- NUM_P: number of periods the function will be sampled at (default: 1)

OUTPUT
- FREQ: frequencies
- COEFF: Fourier series coefficients in the form in the form [a0 a1 ... aM b1 ... bM], where f(t) = a0 + \sum_{m=1}^M { a_m cos(2pi m t/T) + b_m sin(2pi m t/T) }
So the corresponding frequencies are: 0, 1/T, 2/T, ..., M/T
- APSPEC: the first column contains the amplitude spectrum, and the second column the phase spectrum
- If RES=1:
Figure the original function and the Fourier series expansion, and another with the amplitude and phase spectra