fftfunct.c

来自「dspic开发控制程序,有助开发pic单片机程序」· C语言 代码 · 共 117 行

#include <p30f6014.h>
#include <dsp.h>
#include "common.h"

#define FFT_SIZE 256
#define LOG2N 8       		// LOG2(FFT_SIZE)


//外部函数，来自AsmUtils.s
extern void CalcLogMagnitude(fractcomplex* pIn, int* pOut,
										 unsigned int Size); 
extern void CalcLpIirArray( unsigned int Size, fractional* pIn, int* pIIROut, unsigned int LPK); 

static unsigned int page;
static fractcomplex* offset;

/* FFT输入和输出值须存放在 Y space*/
static fractcomplex FftBuf[FFT_SIZE] __attribute__ 
								( (section (".ybss"), aligned(FFT_SIZE*4) ) );

//FFT Twiddle 系数表
//  Wreal[k] = 32767 * cos( 2*PI*k/FFT_SIZE )
//  Wimag[k] = -32767 * sin( 2*PI*k/FFT_SIZE )
//  0 <= k < FFT_SIZE/2
const fractcomplex TWIDDLE_FACTORS[FFT_SIZE/2] = 
{
	 32767,      0,  32757,   -804,  32728,  -1608,  32678,  -2410, 
	 32609,  -3212,  32521,  -4011,  32412,  -4808,  32285,  -5602, 
	 32137,  -6393,  31971,  -7179,  31785,  -7962,  31580,  -8739, 
	 31356,  -9512,  31113, -10278,  30852, -11039,  30571, -11793, 
	 30273, -12539,  29956, -13279,  29621, -14010,  29268, -14732, 
	 28898, -15446,  28510, -16151,  28105, -16846,  27683, -17530, 
	 27245, -18204,  26790, -18867,  26319, -19519,  25832, -20159, 
	 25329, -20787,  24811, -21402,  24279, -22005,  23731, -22594, 
	 23170, -23170,  22594, -23731,  22005, -24279,  21403, -24811, 
	 20787, -25329,  20159, -25832,  19519, -26319,  18868, -26790, 
	 18204, -27245,  17530, -27683,  16846, -28105,  16151, -28510, 
	 15446, -28898,  14732, -29268,  14010, -29621,  13279, -29956, 
	 12539, -30273,  11793, -30571,  11039, -30852,  10278, -31113, 
	  9512, -31356,   8739, -31580,   7962, -31785,   7179, -31971, 
	  6393, -32137,   5602, -32285,   4808, -32412,   4011, -32521, 
	  3212, -32609,   2411, -32678,   1608, -32728,    804, -32757, 
	     0, -32767,   -804, -32757,  -1608, -32728,  -2410, -32678, 
	 -3212, -32609,  -4011, -32521,  -4808, -32412,  -5602, -32285, 
	 -6392, -32137,  -7179, -31971,  -7962, -31785,  -8739, -31580, 
	 -9512, -31356, -10278, -31113, -11039, -30852, -11793, -30571, 
	-12539, -30273, -13278, -29956, -14010, -29621, -14732, -29268, 
	-15446, -28898, -16151, -28510, -16846, -28105, -17530, -27683, 
	-18204, -27245, -18867, -26790, -19519, -26319, -20159, -25832, 
	-20787, -25329, -21402, -24812, -22005, -24279, -22594, -23731, 
	-23170, -23170, -23731, -22594, -24279, -22005, -24811, -21403, 
	-25329, -20787, -25832, -20159, -26319, -19519, -26790, -18868, 
	-27245, -18204, -27683, -17530, -28105, -16846, -28510, -16151, 
	-28898, -15446, -29268, -14732, -29621, -14010, -29956, -13279, 
	-30273, -12539, -30571, -11793, -30852, -11039, -31113, -10278, 
	-31356,  -9512, -31580,  -8739, -31785,  -7962, -31971,  -7179, 
	-32137,  -6393, -32285,  -5602, -32412,  -4808, -32521,  -4011, 
	-32609,  -3212, -32678,  -2411, -32728,  -1608, -32757,  -804
};

/* 汉宁窗(Hanning)数据表*/
static fractional HanningWindow[FFT_SIZE] = 
{
     0,    19,    44,    78,   123,   177,   241,   314,   398,   490,
   593,   705,   826,   957,  1097,  1247,  1405,  1572,  1749,  1934,
  2128,  2330,  2541,  2761,  2988,  3224,  3467,  3718,  3977,  4244,
  4517,  4798,  5086,  5381,  5682,  5989,  6303,  6623,  6949,  7281,
  7618,  7960,  8308,  8660,  9017,  9378,  9744, 10113, 10487, 10864,
 11244, 11627, 12013, 12402, 12793, 13187, 13582, 13979, 14377, 14777,
 15178, 15579, 15981, 16383, 16785, 17187, 17588, 17989, 18389, 18787,
 19184, 19579, 19973, 20364, 20753, 21139, 21522, 21902, 22279, 22653,
 23022, 23388, 23749, 24106, 24458, 24806, 25148, 25485, 25817, 26143,
 26463, 26777, 27084, 27385, 27680, 27968, 28249, 28522, 28789, 29048,
 29299, 29542, 29778, 30005, 30225, 30436, 30638, 30832, 31017, 31194,
 31361, 31519, 31669, 31809, 31940, 32061, 32173, 32276, 32368, 32452,
 32525, 32589, 32643, 32688, 32722, 32747, 32762, 32767, 32767, 32762,
 32747, 32722, 32688, 32643, 32589, 32525, 32452, 32368, 32276, 32173,
 32061, 31940, 31809, 31669, 31519, 31361, 31194, 31017, 30832, 30638,
 30436, 30225, 30005, 29778, 29542, 29299, 29048, 28789, 28522, 28249,
 27968, 27680, 27385, 27084, 26777, 26463, 26143, 25817, 25485, 25148,
 24806, 24458, 24106, 23749, 23388, 23022, 22653, 22279, 21902, 21522,
 21139, 20753, 20364, 19973, 19579, 19184, 18787, 18389, 17989, 17588,
 17187, 16785, 16383, 15981, 15579, 15178, 14777, 14377, 13979, 13582,
 13187, 12793, 12402, 12013, 11627, 11244, 10864, 10487, 10113,  9744,
  9378,  9017,  8660,  8308,  7960,  7618,  7281,  6949,  6623,  6303,
  5989,  5682,  5381,  5086,  4798,  4517,  4244,  3977,  3718,  3467,
  3224,  2988,  2761,  2541,  2330,  2128,  1934,  1749,  1572,  1405,
  1247,  1097,   957,   826,   705,   593,   490,   398,   314,   241,
   177,   123,    78,    44,    19,     0
};

//函数原型
void CalcFFT( fractional* pIn, int* pOut);

//计算输入值的FFT子程序
void CalcFFT( fractional* pIn, int* pOut)
{
int i;
int tmp;
	VectorWindow(FFT_SIZE, pIn, pIn, HanningWindow);	//加汉宁窗
	for (i=0; i<FFT_SIZE; i++)
	{
		FftBuf[i].real = pIn[i];	//输入的数据是实数
		FftBuf[i].imag = 0;			  
	}
	asm("mov #psvpage(_TWIDDLE_FACTORS),w0\nmov w0,_page");
	asm("mov #psvoffset(_TWIDDLE_FACTORS),w0\nmov w0,_offset");
	FFTComplexIP( LOG2N, FftBuf, offset, page );
	BitReverseComplex(LOG2N, FftBuf);
	CalcLogMagnitude( FftBuf, pIn, FFT_SIZE/2);	//创建Log2 幅值
#if 0	//如果要平滑显示
	CalcLpIirArray( FFT_SIZE/2, pIn, pOut, 8192); 
#else
	VectorCopy(FFT_SIZE/2, pOut, pIn);
#endif
}