fft1024.c

C6000DSP的C语言程序示例 一共有3个 内容不错

#include "i_cmplx.h"     	/* definition of the complex type */  
#include "twiddle1024.h" 	/* quantised and scaled Twiddle factors */
#define LL 1024          	/* Maximum length of FFT */ 

/* fft radix-2 funtion using Decimation In Frequency */

//#pragma CODE_SECTION(fft, "mycode"); // Makes the program run from internal memory

void fft(COMPLEX *Y, int N) /* FFT(input sample array, # of points)     */
	{
   	int temp1R, temp1I, temp2R,temp2I;  /* 32 bits temporary storage for */
   									   /* intermediate results			*/					
   	short tempR, tempI, c, s;      /* 16 bits temporary storages for  */
   	                               /* variables */
   	int TwFStep,  /* Step between twiddle factors */
   		TwFIndex, /* Index of twiddle factors */
   		BLStep,    /* Step for incrementing butterfly index */
   		BLdiff,   /* Difference between upper and lower butterfly legs */
   		upperIdx,
   		lowerIdx, /* upper and lower indexes of buterfly leg */ 
   		i, j, k;  /* loop control variables */
   		
   	BLdiff=N;
   	TwFStep=1;
   	for(k=N;k>1;k=(k>>1)) /* Do Log(base 2)(N) Stages */
   		{
      	BLStep=BLdiff;
      	BLdiff=BLdiff>>1; 
      	TwFIndex=0;
      	for(j=0;j<BLdiff;j++)/* Nbr of twiddle factors to use=BLDiff  */
      		{
      		c=w[TwFIndex].real;
      		s=w[TwFIndex].imag;
      		TwFIndex=TwFIndex+TwFStep;                 
      		/* Now do N/BLStep butterflies */  
      		for(upperIdx=j;upperIdx<N;upperIdx+=BLStep)
         		{                              
/* Calculations inside this loop avoid overflow by shifting left once
   the result of every adittion/substration and by shifting left 15 
   places the result of every multiplication. Double precision temporary
   results (32-bit) are used in order to avoid losing information because
   of overflow. Final DFT result is scaled by N (number of points), i.e.,
   2^(Nbr of stages) =2^(log(base 2) N) = N								*/
             		
            	lowerIdx=upperIdx+BLdiff;
            	temp1R     = (Y[upperIdx].real - Y[lowerIdx].real)>>1;
            	temp2R     = (Y[upperIdx].real + Y[lowerIdx].real)>>1;
            	Y[upperIdx].real  =  (short) temp2R;
            	temp1I     = (Y[upperIdx].imag - Y[lowerIdx].imag)>>1;
            	temp2I     = (Y[upperIdx].imag + Y[lowerIdx].imag)>>1;
            	Y[upperIdx].imag  =  (short) temp2I;
            	temp2R     = (c*temp1R - s*temp1I)>>15;
            	Y[lowerIdx].real  = (short) temp2R;
            	temp2I     = (c*temp1I + s*temp1R)>>15;
            	Y[lowerIdx].imag  =  (short) temp2I;
            	}
         	}
         	TwFStep = TwFStep<<1; /* update separation of twiddle factors)*/
       	}

/* bit reversal for resequencing data */

	j=0;
   for (i=1;i<(N-1);i++)
   	{
      k=N/2;
      while (k<=j)
      	{
         j = j-k;
         k=k/2;
         }
      j=j+k;
      if (i<j)
      	{
         tempR=Y[j].real;
         tempI=Y[j].imag;
         Y[j].real=Y[i].real;
         Y[j].imag=Y[i].imag;
         Y[i].real=tempR;
         Y[i].imag=tempI;
         }
      }
  	return;
	}