fft library for mbed
Dependents: 4180_Tuner mbed_capstone 4180_EditThis_copy 4180_EditThis_copy_Demo_Test
Fork of FFT by
Diff: FFT.cpp
- Revision:
- 0:e3af07c00c13
diff -r 000000000000 -r e3af07c00c13 FFT.cpp --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/FFT.cpp Fri Feb 11 00:35:40 2011 +0000 @@ -0,0 +1,132 @@ +/* + @file FFT.cpp + @version: 1.0 + @author: Suky + @web www.micros-designs.com.ar + @date 10/02/11 +*/ +#include "FFT.h" + +// Extracted from Numerical Recipes in C +void vFFT(float data[], unsigned int nn){ +/*Replaces data[1..2*nn] by its discrete Fourier transform, if isign is input as 1; or replaces +data[1..2*nn] by nn times its inverse discrete Fourier transform, if isign is input as -1. +data is a complex array of length nn or, equivalently, a real array of length 2*nn. nn MUST +be an integer power of 2 (this is not checked for!).*/ + unsigned int n,mmax,m,j,istep,i; + double wtemp,wr,wpr,wpi,wi,theta; + float tempr,tempi; + + #define SWAP(a,b) tempr=(a);(a)=(b);(b)=tempr + + n=nn << 1; + j=1; + for (i=1;i<n;i+=2) { + if(j>i){ + SWAP(data[j],data[i]); + SWAP(data[j+1],data[i+1]); + } + m=n >> 1; + while (m >= 2 &&j>m){ + j-=m; + m >>= 1; + } + j+=m; + } + + mmax=2; + while (n > mmax) { + istep=mmax << 1; + theta=(6.28318530717959/mmax); + wtemp=sin(0.5*theta); + wpr = -2.0*wtemp*wtemp; + wpi=sin(theta); + wr=1.0; + wi=0.0; + for (m=1;m<mmax;m+=2) { + for (i=m;i<=n;i+=istep) { + j=i+mmax; + tempr=wr*data[j]-wi*data[j+1]; + tempi=wr*data[j+1]+wi*data[j]; + data[j]=data[i]-tempr; + data[j+1]=data[i+1]-tempi; + data[i] += tempr; + data[i+1] += tempi; + } + wr=(wtemp=wr)*wpr-wi*wpi+wr; + wi=wi*wpr+wtemp*wpi+wi; + } + mmax=istep; + } +} + +// Extracted from Numerical Recipes in C +void vRealFFT(float data[], unsigned int n){ +/*Calculates the Fourier transform of a set of n real-valued data points. Replaces this data (which +is stored in array data[1..n]) by the positive frequency half of its complex Fourier transform. +The real-valued rst and last components of the complex transform are returned as elements +data[1] and data[2], respectively. n must be a power of 2. This routine also calculates the +inverse transform of a complex data array if it is the transform of real data. (Result in this case +must be multiplied by 2/n.)*/ + unsigned long i,i1,i2,i3,i4,np3; + float c1=0.5,c2,h1r,h1i,h2r,h2i; + double wr,wi,wpr,wpi,wtemp,theta; + theta=3.141592653589793/(double) (n>>1); + + c2 = -0.5; + vFFT(data,n>>1); + wtemp=sin(0.5*theta); + wpr = -2.0*wtemp*wtemp; + wpi=sin(theta); + wr=1.0+wpr; + wi=wpi; + np3=n+3; + for (i=2;i<=(n>>2);i++) { + i4=1+(i3=np3-(i2=1+(i1=i+i-1))); + h1r=c1*(data[i1]+data[i3]); + h1i=c1*(data[i2]-data[i4]); + h2r = -c2*(data[i2]+data[i4]); + h2i=c2*(data[i1]-data[i3]); + data[i1]=h1r+wr*h2r-wi*h2i; + data[i2]=h1i+wr*h2i+wi*h2r; + data[i3]=h1r-wr*h2r+wi*h2i; + data[i4] = -h1i+wr*h2i+wi*h2r; + wr=(wtemp=wr)*wpr-wi*wpi+wr; + wi=wi*wpr+wtemp*wpi+wi; + } + data[1] = (h1r=data[1])+data[2]; + data[2] = h1r-data[2]; + +} + + +void vCalPowerf(float Input[],float Power[], unsigned int n){ + unsigned char k,j; + + for(k=0,j=0;k<n;k++,j+=2){ + Power[k]=sqrt(Input[j]*Input[j]+Input[j+1]*Input[j+1]); + } +} + +void vCalPowerInt(float Input[],unsigned char Power[], unsigned int n){ + unsigned char k,j; + + for(k=0,j=0;k<n;k++,j+=2){ + Power[k]=sqrt(Input[j]*Input[j]+Input[j+1]*Input[j+1]); + } +} + +void vCalPowerLog(float Input[],unsigned char Power[], unsigned int n){ + unsigned char k,j; + float Temp; + + for(k=0,j=0;k<n;k++,j+=2){ + if((Input[j]!=0) && (Input[j+1]!=0)){ + Temp=sqrt(Input[j]*Input[j]+Input[j+1]*Input[j+1]); + Power[k]=10.0*log10(Temp); + }else{ + Power[k]=0; + } + } + +}