DSP28335实现FFT的工程实现方法
- 1、FFT的工程应用痛点分析
- 2、解决方案
- 3、代码实现
- 4、实际测试与效果分析
- 5、一些新的想法
- 感谢支持，如果这篇文章对你有帮助请您留下宝贵一票，这也将是我持续输出的动力。
1、FFT的工程应用痛点分析
在DSP28335中要实现FFT，不仅需要对FFT的每个步骤熟悉，还需要考虑在现有控制架构下FFT是否影响正常的控制算法运行。在FFT的工程应用的痛点就是它的蝶形运算的循环。要实现工程应用就要解决掉这个魂环对于控制器的影响。那有什么方法可以克服这个问题呢？
2、解决方案
FFT的算法可分解为三步骤，第一步位反转排序，第二步蝶形运算，第三步是构建输出。其中位反转排序是单循环结构，蝶形运算是三级for循环结构，输出构建是单循环结构。为了将所有循环消除，将FFT的每个层级循环嵌入至控制算法的背景循环中。这样FFT的计算实时性将受到了影响，但对于控制算法的影响将降低到了最小。对于控制算法的影响就是每个控制周期多增加了几条执行语句。按照DSP150MHz的执行周期，增加的这几条指令对于控制算法和通讯只增加了2us不到。
3、代码实现
参考"FFT的完整DSP实现“进行了二次开发，代码见。
取消所有for循环的代码如下：
fft.c文件如下：
/*
* xl_fft.c
*
* Created on: 2023-5-12
* Author: xkl
*
* TEST OK
*/
#include "xkfft.h"
#include "math.h"
#include "string.h"
/*
* FFT Algorithm
* === Inputs ===
* x : complex numbers
* N : nodes of FFT. @N should be power of 2, that is 2^(*)
*/
FFT_T sFFT;
void FFtInit(FFT_T *v)
{
memset((void *)&v->In[0], 0, sizeof(FFT_T));
v->SampleEn = 1;
}
int fft(FFT_T *v, int N)
{
TYPE_FFT *x;
static int i,j,l,k,ip;
static int M = 0;
static int le,le2;
static TYPE_FFT_E sR,sI,tR,tI,uR,uI;
x = v->In;
M = (int)(log(N) / log(2));
/*
* bit reversal sorting
*/
switch(v->StepNum)
{
case 0:
l = N / 2;
j = l;
i=1;
v->StepNum = 1;
break;
case 1:
//BitReverse(x,x,N,M);
if(i <= N-2)
{
if (i < j)
{
tR = x[j].real;
tI = x[j].imag;
x[j].real = x[i].real;
x[j].imag = x[i].imag;
x[i].real = tR;
x[i].imag = tI;
}
k = l;
while (k <= j) {
j = j - k;
k = k / 2;
}
j = j + k;
i++;
}
else
{
v->StepNum = 2;
l = 1;
}
break;
case 2:
// FFT 核心开始
/* loop for ceil{log2(N)} */
if(l <= M)
{
le = (int)pow(2,l);
le2 = (int)(le / 2);
uR = 1;
uI = 0;
sR = cos(XK_FFT_PI / le2);
sI = -sin(XK_FFT_PI / le2);
v->StepNum = 3;
// 跳转至第二层循环
j=1;
}
else
{
v->StepNum = 5;
}
break;
case 3:
/* loop for each sub DFT */
if(j<=le2)
{
i=j-1;
// 跳转至第三层循环
v->StepNum = 4;
}
else
{
v->StepNum = 2;
l++;
}
break;
case 4:
/* loop for each butterfly */
if(i<=N-1)
{
ip = i + le2;
tR = x[ip].real * uR - x[ip].imag * uI;
tI = x[ip].real * uI + x[ip].imag * uR;
x[ip].real = x[i].real - tR;
x[ip].imag = x[i].imag - tI;
x[i].real += tR;
x[i].imag += tI;
i+=le;
}
else
{
tR = uR;
uR = tR * sR - uI * sI;
uI = tR * sI + uI *sR;
j++;
v->StepNum = 3;
}
break;
case 5:
MakeOutput(v,N);
v->StepNum = 6;
break;
case 6:
/*谐波含量计算*/
v->StepNum = 0;
v->SampleEn = 1;
break;
default:
break;
}
return 0;
}
/*
* Inverse FFT Algorithm
*/
int ifft(FFT_T *v, int N)
{
int k = 0;
TYPE_FFT *x;
x = v->In;
for (k=0; k<=N-1; k++) {
x[k].imag = -x[k].imag;
}
fft(v, N); /* using FFT */
for (k=0; k<=N-1; k++) {
x[k].real = x[k].real / N;
x[k].imag = -x[k].imag / N;
}
return 0;
}
void MakeInput(FFT_T *v, float N)
{
int i;
for ( i=0;i<N;i++ )
{
v->In[i].real = sin(XK_FFT_PI*2.0*50.5*i/N)+0.5*sin(XK_FFT_PI*2.0*250.5*i/N);
v->In[i].imag = 0.0f;
}
}fft.h文件如下：
/*
* XK_fft.h
*
* Created on: 2023-5-13
* Author: xkl
*/
#ifndef _XK_FFT_H
#define _XK_FFT_H
#define TYPE_FFT_E float /* Type is the same with COMPLEX member */
#ifndef XK_FFT_PI
#define XK_FFT_PI (3.14159265f)
#endif
#ifndef SAMPLE_NODES
#define SAMPLE_NODES (1024)
#endif
typedef struct {
float real;
float imag;
} COMPLEX;
typedef COMPLEX TYPE_FFT; /* Define COMPLEX in Config.h */
typedef struct
{
TYPE_FFT In[SAMPLE_NODES];
TYPE_FFT_E Out[SAMPLE_NODES/2];
unsigned int SampleEn; //采样使能
unsigned int SampleCnt; //采样计数
unsigned int StepNum;
float Thd;
}FFT_T;
extern void FFtInit(FFT_T *v);
extern int fft(FFT_T *v, int N);
extern int ifft(FFT_T *v, int N);
extern void MakeOutput(FFT_T *v,float N);
extern void MakeInput(FFT_T *v, float N);
extern void ThdCacl(FFT_T *v, float N,float FreBase);
extern FFT_T sFFT;
#endif /* XK_FFT_H_ */定时器中断的数据采集代码如下：
/*
* Task6.c
*
* Created on: 2023年5月13日
* Author: Administrator
*/
#include "Task.h"
//-----------------------------------------------
/*
*/
//-----------------------------------------------
#ifndef TASK6_PERIOD
#define TASK6_PERIOD (1/1000.0)
#endif
//-----------------------------------------------
float Task6Time1 = 0;
//-----------------------------------------------
interrupt void Task6(void)
{
long long TimerS;
long long TimerE;
float signal= 0;
//static float TimeDly = 0;
//===========================================
IER |= MINT13; // 中断优先级设定
EINT;
//-------------------------------------------
TimerS = CpuTimer1.RegsAddr->TIM.all;
Tsc = 0.001;
//===========================================
signal = 100*sin(XK_FFT_PI*2.0*50.0*sFFT.SampleCnt*Tsc)+5*sin(XK_FFT_PI*2.0*250.0*sFFT.SampleCnt*Tsc);
if(sFFT.SampleEn == 1)
{
sFFT.In[sFFT.SampleCnt].real = signal;
sFFT.In[sFFT.SampleCnt].imag = 0.0f;
sFFT.SampleCnt++;
if(sFFT.SampleCnt == SAMPLE_NODES)
{
sFFT.SampleEn = 0;
sFFT.SampleCnt = 0;
}
}
//===========================================
TimerE = CpuTimer1.RegsAddr->TIM.all;
if(TimerE > TimerS)
{
Task6Time1 = ((float)(150000 - TimerE + TimerS ))/150.0;
if(Task6Time1 < 0)
{
}
}
else
{
Task6Time1 = ((float)(TimerS - TimerE ))/150.0;
}
}背景循环的执行代码如下：
/************************************************
* File: Main.c
* Created on: 2023年5月13日
* Author: xkl
* Project: DSP_28335_FFT
**********************************************/
// 系统头文件
#include "DSP2833x_Device.h"
#include "DSP2833x_GlobalPrototypes.h"
#include "DSP2833x_Examples.h"
//-----------------------------------------------
// 用户头文件
#include "Task.h"
//-----------------------------------------------
#include "stdlib.h"
#include "stdio.h"
//-----------------------------------------------
long long TimerS = 0;
long long TimerE = 0;
float fftTime1 = 0;
//-----------------------------------------------
// while loop 时间测试
//-----------------------------------------------
void main (void)
{
//===========================================
// 系统初始化
InitSysCtrl();
MemCopy(&RamfuncsLoadStart, &RamfuncsLoadEnd, &RamfuncsRunStart);
InitFlash();
//===========================================
// 初始化 IO
//===========================================
// 禁止中断 初始化 中断向量表
// 禁止全局中断
DINT;
//初始化PIE 控制寄存器
InitPieCtrl();
//禁止 CPU 中断，清除 中断标志
IER = 0x0000;
IFR = 0x0000;
// 初始化中断向量表
InitPieVectTable();
//用户中断函数 映射
EALLOW;
// 定时器 1 函数指初始化
PieVectTable.XINT13 = &Task6;
EDIS;
//------------------------------------------
// 初始化外设
// 初始化定时器
InitCpuTimers();
//-------------------------------------------
// 配置定时器
ConfigCpuTimer(&CpuTimer1, 150, 1000);
//===========================================
// 使能 CPU 中断
IER |= M_INT13; //CpuTimer1 中断使能
//===========================================
// 使能全局中断
EINT;
ERTM;
//===========================================
// 用户 循环
while(1)
{
if(1 == sDataSample.u16SensorCorrectFlag)
{
FFtInit(&sFFT);
sDataSample.u16SensorCorrectFlag = 0;
CpuTimer1Regs.TCR.bit.TSS = 0;
}
TimerS = CpuTimer1.RegsAddr->TIM.all;
if(sFFT.SampleEn == 0)
{
fft(&sFFT,SAMPLE_NODES);
}
TimerE = CpuTimer1.RegsAddr->TIM.all;
fftTime1 = ((float)(TimerS - TimerE ))/150.0;
}
}4、实际测试与效果分析
直接使用中的代码进行测试，总共FFT执行一次需要9470us。使用改进后的代码，每次执行仅需0.2us。从获取结果的时间维度上看，一次性解算完毕需要9470us，而将解算过程分解到每个背景循环中，解算量是不变的，仅仅是背景循环的次数增多。而对于控制器，它最需要的就是控制周期越短越好。所以这个方法实现FFT的工程应用即保证了FFT解算的结果，而且还对实时控制影响降低到了微乎其微的地步。
5、一些新的想法
对于含有控制器的DSP来实现FFT，FFT的结果反映了控制的效果，FFT的结果不是控制的目标。因此，为了保证控制效果的可观测，在不使用谐波测量仪器的前提下，直接使用控制器的一些资源，计算出控制目标的谐波含量还是一个不错的方案。这不仅节省了测量仪器的费用，而且还节省了一些隐形的人力成本。
对于采用FPGA+DSP的控制器架构实现FFT来看，利用FPGA的特性，使用FPGA实现高精度采样，DSP仅完成FFT的解算也是一个不错的方案，这将解决DSP控制中断与高速数据采样中断的冲突，使控制器资源最大化。
本文章为转载内容，我们尊重原作者对文章享有的著作权。如有内容错误或侵权问题，欢迎原作者联系我们进行内容更正或删除文章。
