Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

I Integrating FFT, python code

  1. Oct 27, 2017 #1
    Hello,
    Thank you for taking time to read my post.

    I have a discrete set of data points that represent an acceleration signal. I want to take the integral of this set of points twice so as to get a function which represents the position over time.
    To accomplish this, I have taken the FFT of the acceleration signal . From there I can extract the complex harmonics and plot the signal with reasonable accuracy. This part I can do and it works great. However I keep running into problems whenever I try to integrate the sinusoids. Can you take a look at my below code and see where I might be going wrong?
    Inputs:
    accel, this is suppose to be the FFT of the signal we want to integrate twice
    w, omega = 2 pi f
    L, how many harmonics you want to use (example, L=10 would use from n=-10,-9.....0...9,10)

    Code (Text):

    #feed a set of FFT points.
    #will alter the coefficients of the real and complex points so as to
    #account for integrating each harmonic twice!
    #ideally, you will be able to input the acceleration and the ouput
    #will be the trajectory of the position
    def int2(accel,w,L):
        sinusoids=[]
        L=int(L)
        period = (2*n.pi)/w
        xaxis = n.linspace(0,period,num=len(accel))
           
        position =[]
        tempMain = 0;
        for i in range(0,(len(accel)-1)):
            for j in range(-L,L):
                a_ = n.real(accel[j])
                b_ = n.imag(accel[j])
                n_ = len(accel)
                if j!=0:
                    j2=j
                    tempMain -= (1/(j2**2))*(a_*n.cos(j*w*xaxis[i])+b_*n.sin(j*w*xaxis[i]))/n_
                else:
                    j2=1
                    tempMain += ((xaxis[i]*xaxis[i])/(2))*(a_+b_)/n_
               
               
                           
            position.append(tempMain)
            tempMain = 0
        #we need to reverse the order of the list. This is bc the
        # algorithim that computes the FFT results in the order
        #being reversed. This is corrected whenever the IFFT is taken.
        position=position[::-1]  
           
        return position
     
     
  2. jcsd
  3. Oct 27, 2017 #2

    Mark44

    Staff: Mentor

    What problems are you having? It would be helpful to know exactly what this code does as opposed to what you expect it to do.
     
  4. Oct 27, 2017 #3
    Hi Mark,

    So the goal of this function is to accept the FFT of a discrete set of data points which represent an acceleration signal a(t) and find the position over time x(t).

    The input is fft(a(t))=a(w).
    I know that for ever harmonic in a(w), that the corresponding time domain signal is An(t) = (real(a(n))*cos(w*n*t) + imag(a(n))*sin(w*n*t))/(total number of points)
    I know that if you integrate cos(nt) you get cos(nt)/n and if you integrate sin(nt) you get sin(nt)/n
    There is a special case, n=0, where this will be a constant(ie n=0 and the sine signal cancel and the cosine goes to 1)

    So I want to integrate everything twice. In this way I will go from acceleration to velocity and then from velocity to position.
    Does that explain things well? Please let me know how I can further elaborate.

    Thank you so much for your feedback :)
     
  5. Oct 28, 2017 #4

    Mark44

    Staff: Mentor

    No to both. ##\int \cos(nt)dt = \frac 1 n \sin(nt)## and ##\int \sin(nt) dt = -\frac 1 n \cos(nt)##, omitting the constants of integration. Could this be why things aren't working for you?
    If you integrate the answers above you get ##-\frac 1 {n^2} \cos(nt)## for cos(nt) integrated twice and ##-\frac 1 {n^2}\sin(nt)## for the second one.
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?
Draft saved Draft deleted