[SciPy-user] Fourier series

Anne Archibald aarchiba@physics.mcgill...
Wed Oct 22 05:44:01 CDT 2008

2008/10/22 Nils Wagner <nwagner@iam.uni-stuttgart.de>:

> Is there a function in scipy to compute the Fourier
> coefficients
> a_0, a_1, b_1, a_2, b_2 of a periodic function f(t)=f(t+T)
> http://en.wikipedia.org/wiki/Fourier_series
> An example would be appreciated.

Not exactly. The usual way to do such a thing would be to evaluate the
function on a grid and use a fast fourier transform:

n = 1024

xs = np.arange(n)*T/float(n)
ys = f(xs)
ft = np.fft.rfft(ys)

You will need to adjust the normalizations somewhat. You will also
need to choose a large enough n to capture all features of interest -
I recommend running it with at least 2 and 4 times the number of
Fourier coefficients you need and comparing. (The extra ones are to
avoid missing features that fall between grid points.)

This method assumes your computational cost would be dominated by the
Fourier transform. If instead your computational cost is dominated by
the function evaluations, you'll have to come up with something more
clever that does adaptive sampling to find all points of interest in
the function. You might try something like using scipy.integrate.quad
and storing the results of all evaluations. You'd then use these
points to evaluate your Fourier transform.


More information about the SciPy-user mailing list