2.c

FFT汇编

// q.cpp : Defines the entry point for the console application.
//

//#include "stdafx.h"


/* FileName:computation.c
 * Description:对AD采样结果进行FFT计算，并根据估算的频率设定通带
 * Programmer:luckuis
 * Data:23/04/08 9:31
*/

/******函数头文件声明*********/
#include <msp430x16x.h>
//#include "dl430fn.h"
#include "stdio.h"
#include "math.h"
//#include "computation.h"
//#include "define.h"
#include <math.h>

#define FFT_SAMPLE_NUM 128
#define FFT_CACULATE_NUM 128
#define FFT_MAX_OF_GRADE 7
#define FFT_TOTAL_OF_NUM FFT_CACULATE_NUM
#define FS 100.0

float fftFreq=0.0;
float temp10=0;

__no_init int fft_num[FFT_SAMPLE_NUM];
__no_init int tableVOfFFT[FFT_CACULATE_NUM];  //FFT计算结果的虚部
__no_init float tableROfFFT_1[FFT_CACULATE_NUM/2];//功率谱计算。
//------------------------级循环单元----------------------------
//级一层使用的内存单元
//------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
const int maxOfGrade=FFT_MAX_OF_GRADE;                          //FFT计算过程中总共需要计算的级数
const int totalOfNum=FFT_TOTAL_OF_NUM;                        //参与FFT计算的数据点数
int gradeOverFlowFlag;                        //级运算结束后，根据这位判断同级中最大的数是否超出Q13表示范围，如果超过了就将这一级的运算结果整体右移2位
int counterOfFFTShiftTimes;                            //整个FFT计算完成后，总共移位的次数。
int tempDataOfOverFlowJudge;                 //蝶形因子计算完成后，需要判断结果是否超出Q13，判断过程中的变量，就是放在该临时变量中。
int gradeOfW;                                //旋转因子级系数，初始值为数据总点数，第一级开始就除以2，以后每开始一级都除以2
int numOfGrade;                             //级数，128点就分为7级  
int maxOfW;                                  //最大旋转因子系数，在同一级计算中，是根据旋转因子系数进行循环的，所有的旋转因子计算完成后，该级就计算完成了
//int stepOfBf;                                //同级运算中，对于一个蝶形因子，低位运算单元与高位运算单元的距离
int stepBetweenW;                            //同级运算中，对于具有相同旋转因子的蝶形单元将被放在一起计算，该内存中存放在这些蝶形单元之间跨距
//------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
int factorOfW;                               //存放旋转因子系数，是根据旋转因子系数进行循环的，所有的旋转因子计算完成后，该级就计算完成了
int shiftOfLowBf;                            //存放蝶形因子低位偏移量，然后蝶形因子就根据该偏移量找到对应的内存单元中存放的蝶形因子系数，在蝶形运算过程中，同级相同蝶形因子计算过程中，一个蝶形单元运算完成后，该偏移量将加上该级的蝶形运算单元跨距，开始新的同旋转单元蝶形因子计算。
int shiftOfW;                                //存放着旋转因子偏移量，根据该偏移量查找旋转因子
//------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
int shiftOfHighBf; 
int L1VCOS;       
int L1VSIN;        
int L1RCOS;        
int L1RSIN; 

int startAdress;                //第2级循环中，用来设置3级循环的蝶形单元的起始地址
int currentW;                   //表示第2级循环中，当前的旋转因子

int temp1;
int temp2;
int temp11;
int temp20=1;
unsigned int k[4]={0,0,0,0};

int trueStartAdress;

extern unsigned long asmFFT();

void findMaxes();






void main(void)
{
  WDTCTL = WDTPW + WDTHOLD; 
  unsigned int i=0;

  for(i=0;i<FFT_TOTAL_OF_NUM;i++)
  {
    temp10=2.0*4096.0*sin(2*3.1415926*33*(i/FS));
    fft_num[i]=(int)temp10;
    tableVOfFFT[i]=0.0;
  }
  
 
  asmFFT();
  
  unsigned long  yrm=0;
      unsigned long  ytm=0;
      unsigned int tempHig=0;
      unsigned int tempLow=0;
      
          /*功率谱计算*/
      for (i=0; i<totalOfNum; i+=2)
      {
        //yrm=fft_num[i];
        
        //调用乘法器实现乘法运算
        //实部虚部平方结果都用Q13 表示，现将平方后的数据也用Q13表示，因此需要就将平方后的数据右移13位。
        MPYS=fft_num[i];
        OP2=fft_num[i]; 
        tempHig=RESHI;
        tempHig=tempHig<<3;
        tempLow=RESLO;
        tempLow=tempLow>>13;  
        tempHig=tempLow+tempHig; 
        //yrm=(RESHI<<16)+RESLO;
        yrm=tempHig;

        //调用乘法器实现乘法运算
        MPYS=tableVOfFFT[i];
        OP2=tableVOfFFT[i]; 
        
        tempHig=RESHI;
        tempHig=tempHig<<3;
        tempLow=RESLO;
        tempLow=tempLow>>13;
        tempHig=tempLow+tempHig; 
        
        ytm=tempHig; 
        
        fft_num[i]=yrm+ytm; 
      	//tableROfFFT_1[i]=(tableROfFFT[i]*tableROfFFT[i]+tableVOfFFT[i]*tableVOfFFT[i]);
      }
      
      findMaxes();
      __no_operation();
}


void findMaxes()
{
  float fs=FS;
  int i=0;
  float max=0;
  float Nf;
  float tempfloat=0;  
  
  for (i=0; i<totalOfNum; i+=2)  //找出功率谱的最大值的序号
  {
     tempfloat=fft_num[i];
     if (tempfloat>max)
     {
	max=tempfloat;
	k[0]=i;
     }
  } 
  
  unsigned int exchange=0;
  unsigned int flag=0;
  for(i=0;i<maxOfGrade;i++)   //最大值的码位倒置
  {
    if ((k[0]&temp20)>0)
      flag=1;
    else
      flag=0;
    exchange=(exchange<<1)+flag;
    temp20=temp20<<1;
  }
  k[1]=exchange+1;          //确定最大值前后的序号
  k[2]=exchange-1;
  k[3]=exchange;
  
  temp20=1;
  exchange=0;
  for(i=0;i<maxOfGrade;i++)  //最大值后的序号，码位倒置
  {
    if ((k[1]&temp20)>0)
      flag=1;
    else
      flag=0;
    exchange=(exchange<<1)+flag;  
    temp20=temp20<<1;
  }  
  k[1]=exchange;
  
  temp20=1;
  exchange=0;
  for(i=0;i<maxOfGrade;i++)  //最大值前的序号，码位倒置
  {
    if ((k[2]&temp20)>0)
      flag=1;
    else
      flag=0;    
    exchange=(exchange<<1)+flag;
    temp20=temp20<<1;
  }  
  k[2]=exchange;
  
   /*重心校正法校正   */
      if(fft_num[k[1]]<fft_num[k[2]])
        Nf = -sqrt(fft_num[k[2]] ) / ( sqrt(fft_num[k[0]] ) + sqrt(fft_num[k[2]]) ) + k[3];     
      else      
        Nf = sqrt(fft_num[k[1]]) / ( sqrt(fft_num[k[0]] ) + sqrt(fft_num[k[1]]) ) + k[3];  

  fftFreq=Nf*fs/(float)totalOfNum;
  
}