FFT.cpp

00001 // --------------------------------------------------------------------------------------------------------
00002 //
00003 //  September 2018
00004 //  Author: Valeria Toffoli, Rohan Gurav 
00005 // --------------------------------------------------------------------------------------------------------
00006 //
00007 //  FFT.h
00008 //
00009 // --------------------------------------------------------------------------------------------------------
00010 // 
00011 //  This library provides all the functions necessary to obtain the FFT form an analog ADXL connected with   
00012 //  EV-COG-AD3029 or EV-COG-AD4050 Board. 
00013 //
00014 //  Permission is hereby granted, free of charge, to any person obtaining
00015 //  a copy of this software and associated documentation files (the
00016 //  "Software"), to deal in the Software without restriction, including
00017 //  without limitation the rights to use, copy, modify, merge, publish,
00018 //  distribute, sublicense, and/or sell copies of the Software, and to
00019 //  permit persons to whom the Software is furnished to do so, subject to
00020 //  the following conditions:
00021 //
00022 //  The above copyright notice and this permission notice shall be
00023 //  included in all copies or substantial portions of the Software.
00024 //
00025 //  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
00026 //  EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
00027 //  MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
00028 //  NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
00029 //  LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
00030 //  OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
00031 //  WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
00032 //
00033 // --------------------------------------------------------------------------------------------------------
00034 
00035 
00036 
00037 #include <stdint.h>
00038 #include <math.h>
00039 #include <complex>
00040 #include "mbed.h"
00041 #include "FFT.h"
00042 
00043 
00044 int FFT::log2(int N)    //funzione per calcolare il logaritmo in base 2 di un intero
00045 {
00046   int k = N, i = 0;
00047   while(k) {
00048     k >>= 1;
00049     i++;
00050   }
00051   return i - 1;
00052 }
00053 
00054 int FFT::check(int n)    //usato per controllare se il numero di componenti del vettore di input è una potenza di 2
00055 {
00056   return n > 0 && (n & (n - 1)) == 0;
00057 }
00058 
00059 int FFT::reverse(int N, int n)    //calcola il reverse number di ogni intero n rispetto al numero massimo N
00060 {
00061   int j, p = 0;
00062   for(j = 1; j <= log2(N); j++) {
00063     if(n & (1 << (log2(N) - j)))
00064       p |= 1 << (j - 1);
00065   }
00066   return p;
00067 }
00068 
00069 void FFT::order(complex<double>* f1, int N)     //dispone gli elementi del vettore ordinandoli per reverse order
00070 {
00071   complex<double> f2[N];
00072   for(int i = 0; i < N; i++)
00073     f2[i] = f1[reverse(N, i)];
00074   for(int j = 0; j < N; j++)
00075     f1[j] = f2[j];
00076 }
00077 
00078 void FFT::transform(complex<double>* f, int N)     //calcola il vettore trasformato
00079 {
00080   order(f, N);    //dapprima lo ordina col reverse order
00081   complex<double> W[N / 2]; //vettore degli zeri dell'unità.
00082                             //Prima N/2-1 ma genera errore con ciclo for successivo
00083                            //in quanto prova a copiare in una zona non allocata "W[N/2-1]"
00084   W[1] = polar(1., -2. * M_PI / N);
00085   W[0] = 1;
00086   for(int i = 2; i < N / 2; i++)
00087     W[i] = pow(W[1], i);
00088   int n = 1;
00089   int a = N / 2;
00090   for(int j = 0; j < log2(N); j++) {
00091     for(int i = 0; i < N; i++) {
00092       if(!(i & n)) {
00093         /*ad ogni step di raddoppiamento di n, vengono utilizzati gli indici */
00094         /*'i' presi alternativamente a gruppetti di n, una volta si e una no.*/
00095         complex<double> temp = f[i];
00096         complex<double> Temp = W[(i * a) % (n * a)] * f[i + n];
00097         f[i] = temp + Temp;
00098         f[i + n] = temp - Temp;
00099       }
00100     }
00101     n *= 2;
00102     a = a / 2;
00103   }
00104 }
00105 
00106 void FFT::fourier(complex<double>* f, int N, double d)
00107 {
00108   transform(f, N);
00109   for(int i = 0; i < N; i++)
00110     f[i] *= d; //moltiplica il vettore per il passo in modo da avere il vettore trasformato effettivo
00111 }