MSP430 FFT算法实现

2019-07-23 15:30发布

*****************main   programe********************/
#include  <msp430x14x.h>

#include       <math.h>     
#include       <stdio.h>     
#include       <stdlib.h>     
#include      "typedef.h "
      

float       result[8];   
struct     compx   s[8];      
int       Num=8;
const   float   pp=3.141592653589793;

extern void FFT(struct   compx *s,int num);


double mpr[8],mpi[8],mfr[8],mfi[8];
int n,k,l,il;
extern void kbfft();
extern void FFT_dingdian();
main()
{

    //   ????
    int   i;
    FILE   *fin,*fout;
    if((fin=fopen( "indata.txt ", "r "))==NULL)
    {  
        printf( "can 't   open   infile ");
        exit(0);
    };
    for(i=0;i <Num;i++)
    {   
        //fscanf(fin, "%lf ",&s.real);
        s.real=cos(2*3.1415926*(30)*(i)/80.0);//
        s.imag=(double)0;
        mpr=s.real;
        mpi=0;        
        
    };
    fclose(fin);
    //////////////
      
    //FFT(s,Num);
   
    //n=8;k=3;l=0;il=0;
    //kbfft(mpr,mpi,n,k,mfr,mfi,l,il);

   
    FFT_dingdian();
      
    //????
    if((fout=fopen( "fftresult.txt ", "w "))==NULL)
    {
        printf( "can 't   open   outfile ");
        exit(0);
    };
   
    for(i=0;i <Num;i++)
    {
        // result=sqrt(pow(s.real,2)+pow(s.imag,2));
        fprintf(fout, "%lf+(%lf)i   ",s.real,s.imag);
    };
   
    fclose(fout);

}

//---定点---
#include  <msp430x14x.h>
#include "math.h"
#define NN 8
//************************余弦表**********************************
const  float  COS_tab[91]=
                  {1.0000,0.9998,0.9993,0.9986,0.9975,0.9962,0.9945,0.9925,0.9903,
                   0.9877,0.9848,0.9816,0.9781,0.9744,0.9703,0.9659,0.9614,0.9563,
                   0.9511,0.9455,0.9397,0.9336,0.9272,0.9205,0.9135,0.9063,0.8988,
                   0.8910,0.8829,0.8746,0.8660,0.8572,0.8480,0.8387,0.8290,0.8192,
                   0.8090,0.7986,0.7880,0.7771,0.7760,0.7547,0.7431,0.7314,0.7193,
                   0.7071,0.6947,0.6820,0.6691,0.6561,0.6428,0.6293,0.6157,0.6018,
                   0.5878,0.5736,0.5592,0.5446,0.5299,0.5150,0.5000,0.4848,0.4695,
                   0.4540,0.4384,0.4226,0.4067,0.3907,0.3746,0.3584,0.3420,0.3256,
                   0.3090,0.2924,0.2757,0.2588,0.2419,0.2250,0.2079,0.1908,0.1736,
                   0.1564,0.1392,0.1219,0.1045,0.0872,0.0698,0.0523,0.0349,0.0174,
                   0.0000};

char Re(char M,char k)
{
  char m=1;
  int n;
  while(M>0){m=m*2;M--;}//m==pow(2,M)
  n=360 * k /m;
  return n;
  
}

void FFT_dingdian( void )
{
    signed char flag1=1,flag2=1;
    char  i,j,k,R,I;
    char p,p1;
    //float i5,i6,i7,i8;   
    char   M,L;
    int    a,b,x=1,y=1,z=1;
    int    N1;
   
    struct comp
    {
      float re;
      float im;
    }aa[NN]={{1,0},{0,0},{-1,0},{0,0},{1,0},{0,0},{-1,0},{0,0}};//,{9,0},{10,0},
             //{11,0},{12,0},{13,0},{14,0},{15,0},{16,0}};//
   
    struct comp bb,bb1,bb2;
   

    L=0;
    N1=NN;
    while(N1>0)
    {N1=N1/2;L++;}// L=log(NN)/log(2);
    L=L-1;
   
   //输入变址
   for(i=0,j=0;i<NN;i++)
    {
      if(i>=j){k=NN/2;}
      else
      {
        k=NN/2;
        bb=aa;
        aa=aa[j];
        aa[j]=bb;
        
      }
      while(k<=j){j=j-k;k=k/2;
                  if(k==j)break;}
      
       j=j+k;  
    }  
  //变址结束
   

   
  //由下开始FFT计算     
   
    for(M=1;M<=L;M++)
    {
         x=M-1;
         a=1;
         while(x>0){a=a*2;x--;}//a==2的M-1次方
      for(k=0;k<a;k++)//a==2的M-1次方,即每一级有2的M-1次方个不同的旋转因子,步进值是1
                                      //---同时也是蝶形结两结点间距
        {
             y=M;
             b=1;
             while(y>0){b=b*2;y--;}//b==2的M次方
             //p=0;
          for(p=k;p<NN-1;p=p+b)//b==2的M次方----相同旋转因子的蝶形结间距
              {
               
               
               p1=p+a;//a==2的M-1次方
               
               R=Re(M,k);//----------------------=360*k/(2^M)--k最大为2^M-1
               I=Re(M,k);//Re(M,k)==Im(M,k)
           
               if(R>=0&&R<=90){flag1=1;}
               else if(R>90&&R<=180){ R=180-R;flag1=-1;}
               else if(R>180&&R<=270){ R=R-180;flag1=-1;}
               else if(R>270&&R<=360){ R=360-R;flag1=1;}
               
               
               if(I>=0&&I<=90){I=90-I;flag2=1;}   
               else if(I>90&&I<=180){I=I-90;flag2=1;}
               else if(I>180&&I<=270){I=270-I;flag2=-1;}
               else if(I>270&&I<=360){I=I-270;flag2=-1;}
           
              bb.re=aa[p].re;
              bb.im=aa[p].im;
              
           
              bb1.re=aa[p1].re*(flag1)*COS_tab[R];
              bb1.im=aa[p1].im*(flag2)*COS_tab[I];
              bb2.re=aa[p1].re*(flag2)*COS_tab[I];
              bb2.im=aa[p1].im*(flag1)*COS_tab[R];
           aa[p].re=bb.re+bb1.re  + bb1.im;
           aa[p].im=bb.im+bb2.im - bb2.re ;
           aa[p1].re=bb.re-bb1.re- bb1.im;
           aa[p1].im=bb.im-bb2.im + bb2.re ;
           }     
         }
      }   
   
  //FFT计算结束
   return;   
}
友情提示: 此问题已得到解决,问题已经关闭,关闭后问题禁止继续编辑,回答。
该问题目前已经被作者或者管理员关闭, 无法添加新回复
8条回答
cehuafan
1楼-- · 2019-07-24 18:16
 精彩回答 2  元偷偷看……
cehuafan
2楼-- · 2019-07-24 20:46
现在很流行频谱的设计

一周热门 更多>