[SciPy-User] Help with fast fourier transform

[Chris Michalski]

Looks like your problem is one of expectations.  

You did an FFT on real data - specifically a "pulse" that is 100 point long surrounded by 450 "zeros" on either side.  Since the FFT function computes a complex FFT - you got negative and positive frequencies.  As the FFT samples come out, they represent "DC" (zero freq) through just shy of fs/2 (fs = sample rate), then the output wraps to -fs/2 which then advances to just shy of  "0".

Your original output is correct - it just needs to be reordered.  The "desired" output you referenced (the picture) takes the last "N/2" complex points from the FFT function and moves then to the first N/s points.   Once the data is reordered, before you plot you have to set the frequency resolution of the FFT.  In a complex FFT output, the first sample in the reordered output would correspond to -fs/2 and the uppermost frequency would be (about fs/2)  The complication is that there is a zero frequency (DC offset) term.  Since the FFT returns an even number of samples (specifically a power of two), the convention I remember is that there are more "negative frequency" samples than positive frequency samples.

By the way.  There is a much simpler way to get a pulse, specifically:

import numpy as np

f = np.zeros(1000, dtype='d')
f[450:550] = 1.0

If "F" is the FFT of "f" (the time function), then the reordering can be accomplished with

# Create new array to hold the reordered output (holding complex data)
size = len(F)
FF = np.array(size, dtype='D')

FF[0:size/2] = F[size/2:]
FF[size/2+1:]  = F[0:size/2]

Getting the hang of thinking of operations as vectored code rather than loops can take some practice - but well worth the results.  

chris


On Apr 29, 2010, at 9:17 AM, Oscar Gerardo Lazo Arjona wrote:

> Hello! I'm new to this mailing list and to numpy in general.
> I need to calculate fft for my optics class. But I'm having trouble.
> 
> This is the code I'm using (commented):
> 
> import pylab
> import numpy as np
> 
> def g(x):
>   if x>450 and x<550:
>       return 1
>   else:
>       return 0
>      f=[g(x) for x in range(0,1000)]
> f=np.array(f)
> 
> #this funtion can be ploted as
> #http://i41.tinypic.com/k0shnk.png
> #which is a step function "centered" at 500
> 
> 
> F=np.fft.fft(f)
> #when calculate the fft of f i get an array of complex numbers
> #whose absolute value can be ploted as
> #http://i40.tinypic.com/1zcecxx.png
> 
> #But that is not the desired output.
> #Instead of that i expect something that can be ploted like this
> #http://i39.tinypic.com/1zmch0g.png
> 
> 
> #what i think must be happening because my function f
> #has an offset of 500 (it's supposed to be centereed at 500)
> 
> #So i think it all reduces to two options:
> #somehow telling fft to consider the origin at 500
> #(indicate thetell the offset)
> 
> #Or make fft accept a list of points like
> #[[-2,0],[-1,1],[0,1],[1,1],[2,0]]
> #so that it can know the position of the step relative to the origin
> 
> Please help!
> 
> Oscar
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user