信号的频谱分析及DSP实现

《现代信号处理课程设计》课程设计报告设计题目 信号的频谱分析及DSP实现目录第1章绪论 1设计背景 1设计要求 3设计思路简介 3第2章系统开发平台与环境 3CCS开发环境 3SEED-DEC2812开发实验箱 33...............................4FFT算法实现 4频谱分析系统设计 8源程序代码 9第4章系统仿真 14仿真设置 14仿真图 16第6章总结 20参考文献 201章绪论1.1设计背景1.1设计背景随着数字技术与计算机技术的发展，数字信号处理(DSP)技术已深入到各个学科领域。近些年来，数字信号处理技术同数字计算器、大规模集成电路等，有了突飞猛进的发展。在数字信号处理中，离散傅里叶变换(Discrete．TimeFourierTransform，DFTDFT（FFTDFTFFTFFTFFTFFT数字信号处理器（DSP）是一种可编程的高性能处理器。它不仅是一种适用于数数字信号处理器（DSP）是一种可编程的高性能处理器。它不仅是一种适用于数字信号处理，而且在图像处理、语音处理、通信等领域得到广泛的应用。DSP字信号处理，而且在图像处理、语音处理、通信等领域得到广泛的应用。DSP理器中集成有高速的乘法硬件，能快速的进行大量的乘法加法运算。理器中集成有高速的乘法硬件，能快速的进行大量的乘法加法运算。11965(FFT)算法，FFT36.4%；影响，在加汉宁（Hanning）窗时，只进行幅值恢复时的最大幅值误差仍高达9070差较大的问题。1975差较大的问题。1975JohnC．Burges确测量问题。1983确测量问题。1983ThomasGrandkeHanning高了离散高次谐波参数的分析精度。1993年，丁康和谢明提出了三点卷积法幅高了离散高次谐波参数的分析精度。1993年，丁康和谢明提出了三点卷积法幅的一些问题。1994年，谢明、丁康等提出和发展了比例频谱校正方法，使内差的一些问题。1994年，谢明、丁康等提出和发展了比例频谱校正方法，使内差行了深入系统的分析和研究。1996年，余佳兵，史铁林等提出了采用复调制细行了深入系统的分析和研究。1996年，余佳兵，史铁林等提出了采用复调制细校正以解决密集频率成分的离散频谱的校正问题。1997年，谢明、丁康等分析校正以解决密集频率成分的离散频谱的校正问题。1997年，谢明、丁康等分析识别和校正两个己发生干涉的密集频率成分的自动判定和校正的方法。1998年识别和校正两个己发生干涉的密集频率成分的自动判定和校正的方法。1998年NN/22确的频率。1999年，丁康、谢明等提出了对连续时域信号分前后两段作傅里叶确的频率。1999年，丁康、谢明等提出了对连续时域信号分前后两段作傅里叶进行频率和相位校正。2001年，徐培民、杨积东、闻邦椿提出了自动识别和修进行频率和相位校正。2001年，徐培民、杨积东、闻邦椿提出了自动识别和修1-61.2CCCSFFT仪对常用信号进行频谱分析。1.322.1CCSCCS式信号处理程序的编制和测试，它能够加速开发进程，提高工作效率。TMS320C54xCCS（CodeComposerStudio）FFT22.1CCSCCS式信号处理程序的编制和测试，它能够加速开发进程，提高工作效率。CCS提供了基本的代码生成工具，它们具有一系列的调试、分析能力。CCS支持如下图1所示的开发周期的所有阶段。3设计

编程和编译创建工程文件、编写源代码和配置文件

调试语法检查、探测点设置和日志保存等

分析实时调试、统计和跟踪图12.2SEED-DEC2812开发实验箱SEED-DECxxxxDSPDSPDSPSEED-DEC2812嵌入式DSP开发板原理框图如图2所示：图23FFT算法实现快速傅里叶变换（FFT）是离散傅里叶变换的快速算法，它是根据离散傅里叶变换的奇、偶、虚、实等特性，对离散傅里叶变换的算法进行改进获得的，它对离散傅里叶变换并没有新的发现。有限长序列x(n)及其频域表示X(k)可由以下离散傅立叶变换得出4（8）

x(k)=DFT[x(n)]=

Nn0

x(n)WnkN

0kN1（9）

x(n)=IDFT[X(k)]=j2nk

1N1X(kWnkN Nnkn0

0kN1其中WnkeN

N 。式（8）称为离散傅立叶正变换，式（9）称为离散傅立叶逆变换，x(n)与X(k)构成了离散傅立叶变换对。N次复数乘法和N-1全部X(k)(0kN1N2次复数乘法和N(N-1)次复数加法。实现一次复数乘法需要四次实数乘法和两次实数加法，一次复数加法需要两次实数加法，因此直接计算全部X(k)共需要4N2次实数乘法和2N(2N-1)次实数加法。当N减少计算离散傅里叶变换运算量的问题变得至关重要。为减少运算量，提高运算速度，就必须改进算法。计算DFT过程中需要完成程中的运算，从而减少计算DFT所需的时间。如前所述，N点的DFT的复乘次数等于N2。显然，把N点的DFT分解为几个较短的DFT，可是乘法的次数大大减少。另外，旋转因子Wm具有明显的周期性和N对称性，其周期为：WmlNeN其对称性表现为：

j2(mlN) j2N N

WNWmWNm或NmWmN N N NFFT算法就是不断地把长序列的DFT分解成几个短序列的DFT，并利用Wm的N周期性和对称性来减少DFT的运算次数。Wnk具有以下固有特性：NWnk的周期性：WnkWnN

Wn(kN）N N N N5Wnk的对称性：Wnk（W

Wn(Nk)N N N NWnk的可约性：WnW ,W WnN N N/n N Nn另外，WN/2(kN

Wk。N N N利用Wnk的上述特性，将x(n)或X(k)序列按一定规律分解成短序列进行运N算，这样可以避免大量的重复运算，提高计算DFT的运算速度。算法形式有很多种，但基本上可以分为两大类，即按时间抽取(DecimationInTime，DIT)FFT算法和按频率抽取(DecimationInFrequency，DIF)FFT算法。如果序列x(n)的长度N2M，其中M是整数(如果不满足此条件，可以人为地增补零值点来达到)，在时域上按奇偶抽取分解成短序列的DFT，使最小DFT运算单元为2FFT运算中最小DFT运算单元称为基为基-2时间抽取FFT(DIT-FFT)算法]将x(n)按n为奇偶分解成两个子序列，当n为偶数时，令n=2r；当n为奇数时，令n=2r+l；可得到（10）

x(2r)x1

(r),x(2r1)x2

N(r),r0,..., 12则其DFT可写成X(k)

N2

x(2r)W2rkN

N2

x(2r1)W(2r1)kNr0 r0N2r0

x(r)W2rkN

N2r0

x(r)W(2r1)kNN2r0

N12x(2r)Wrk x(2rrkN/2 N2r0（11）

X(k)Wk1 N

(k)2X(k)X1

(k均分别是N/2x1

(n)x2

(n)的DFT，而且r与k的取值满足0,1,…,N/2-1。而X(k)是一个N点的DFT，因此式（11）只计算了X(k)的前N/2的值。由DFT和Wnk的性质可得到X(k)的后N/2的值为：N6N N kN NX(k )X(k )W 2X(k )

(k)Wk

(k) （12）2 1

N 2 2 1 N 2式和式表明，只要计算出两个N/2DFTx1

(kx2

(k，经NX(k)N2MN/22M2DFTXm(p)Xm(q)

m1

p

m1

(q)为输出数据，W

nk为旋转N因子，则对于基-2DIT-FFT算法，蝶形运算的基本公式为{Xm1X

(p)(q)X

(p)m(p)X

(qm km m N其图形表示如图3所示，称Xm(P)为上结点，Xm(q)为下结点。图3时间抽取蝶形运算单元对于一个8点的N=8的基-2DIT-FFT流图，如图4所示。图4N=8DIT-FFT运算流图7根据上述算法原理及运算流图，可以得出基-2DIT-FFT的基本特点，特点如下。N2MMN/2个蝶形运算单元，总共所N N需蝶形运算个数为

M log N。2 2 2运算量估计：每个蝶形运算需要一次复数乘法和两次复数加(减)法，N点FFTN共需要

log2

NNlog2

N次复数加(减)法。实际上有些蝶形运算不需要做复乘。这同一组存储器中的同一位置，不需要另辟存储单元，直到最后输出。2FFT输出的X(k)X(n)是按的倒序存入存储单元，即为倒序输入，正序输出。这种顺序看起来相当杂乱，然而它是有规律的，即位码倒序规则。8点FFT的三次迭代运算可以看出WkN代中，只有一种类型的蝶形运算系数，即W0，参加蝶形运算的两个数据点间隔8为l；在第二级迭代中，有两种类型的蝶形运算系数，分别是W0和W2，参加蝶8 8形运算的两个数据点间隔为2；在第三级迭代中，有四种类型的蝶形运算系数，分别是W0，W1，W2，W3，参加蝶形运算的两个数据点间隔4。可见，每次迭8 8 8 8最多，参加蝶形运算的两个数据点的间隔也最大，为N/2。频谱分析系统设计硬件方案系统的硬件实现方案如图5所示：8待待检测信号信号隔离电源、复位时钟上位机信号滤波 AD转换TMS320F2812USB电平转换图5信号前向处理电路设计信号的前向处理主要进行了隔离和低通滤波。隔离采用精密隔离放大器IS0124来实现，它通过封装在器件内部的2pF的隔离电容来实现输入级和输出级信号的可使信号受损的噪声通过信号的传播，将每个信号分离到一个干净的信号子系IS0124AD对信号进行模拟滤波。模拟滤波采用连续时间滤波器MAX274来完成。AD转换设计F2812内部带有12位的AD转换模块，考虑到信号采集和分析的精度，本系统外扩了14位的AD转换芯片MAXl320。其独立的采样／保持(T／H)电路为每一通道提供了同时采样，可以保持输入信号间的相对相位信息。软件设计思路系统上电，进人main()函数后，首先进行系统的初始化，包括系统时钟初始化、XINTFSCI定时器初始化等，初始化完成后，启动定时器，使能中断，接着进入主循环，系统开始正常运转。开始开始波形发生FFT初始化工作变量计算步长按照编码逆序排列输入序列调用波形发生子程序产生波形(3个正弦波)用标准C的sin计算当前波形值9(128点)用蝶形算法计算FFT计算功率谱结束计算功率谱63.3源程序代码#include"DSP281x_Device.h"//DSP281xHeaderfileIncludeFile#include"DSP281x_Examples.h"//DSP281xExamplesIncludeFile#include"DSP281x_Device.h"//DSP281xHeaderfileIncludeFile#include"DSP281x_Examples.h"//DSP281xExamplesIncludeFile#include"f2812a.h"#include"math.h"#definePI3.1415926#defineSAMPLENUMBER128//Perruptvoidadc_isr(void);voidInitForFFT();voidMakeWave();//voidFFT(floatdataR[SAMPLENUMBER],floatdataI[SAMPLENUMBER]);//Globalvariablesusedinthisexample:Uint16LoopCount;Uint16ConversionCount;Uint16Voltage1[1024];Uint16Voltage2[1024];inti;intINPUT[SAMPLENUMBER],DATA[SAMPLENUMBER];floatfWaveR[SAMPLENUMBER],fWaveI[SAMPLENUMBER],w[SAMPLENUMBER];floatsin_tab[SAMPLENUMBER],cos_tab[SAMPLENUMBER];main(){InitSysCtrl();DINT;//¹ØÖжÏInitPieCtrl();IER=0x0000;IFR=0x0000;InitPieVectTable();//Interruptsthatareusedinthisexamplearere-mappedto10//ISRfunctionsfoundwithinthisfile.EALLOW; //ThisisneededtowritetoEALLOWprotectedPieVectTable.ADCINT=&adc_isr;EDIS;AdcRegs.ADCTRL1.bit.RESET=1; //ResettheADCasm("RPT#10||NOP");AdcRegs.ADCTRL3.all=0x00C8;AdcRegs.ADCTRL3.bit.ADCBGRFDN=0x3;AdcRegs.ADCTRL3.bit.ADCPWDN=1; //PoweruprestofADC//EnableADCINTinPIEPieCtrlRegs.PIEIER1.bit.INTx6=1;IER|=M_INT1;//EnableCPUInterrupt1EINT; //EnableGlobalinterruptINTMERTM; //EnableGlobalrealtimeinterruptLoopCount=0;ConversionCount=0;//ConfigureADCAdcRegs.ADCMAXCONV.all=0x0001; //Setup2conv'sonAdcRegs.ADCCHSELSEQ1.bit.CONV00=0x0;//SetupADCINA3as1stSEQ1//conv.AdcRegs.ADCCHSELSEQ1.bit.CONV01=0x1;//SetupADCINA2as2ndSEQ1//conv.AdcRegs.ADCTRL2.bit.EVA_SOC_SEQ1=1; //EnableEVASOCtostart//SEQ1AdcRegs.ADCTRL2.bit.INT_ENA_SEQ1=1; //EnableSEQ1interrupt//(everyEOS)//ConfigureEVA//AssumesEVAClockisalreadyenabledinInitSysCtrl();EvaRegs.T1CMPR=0x0080; //SetupT1compareEvaRegs.T1PR=0x10; //SetupperiodregisterEvaRegs.GPTCONA.bit.T1TOADC=1; //EnableEVASOCinEVAEvaRegs.T1CON.all=0x1042; //Enabletimer1//(upcountmode)//WaitforADCinterruptwhile(1){LoopCount++;}}11interruptvoid adc_isr(void){Voltage1[ConversionCount]=AdcRegs.ADCRESULT0>>4;Voltage2[ConversionCount]=AdcRegs.ADCRESULT1>>4;//If40conversionshavebeenlogged,startoverif(ConversionCount==1023){ConversionCount=0;}elseConversionCount++;//ReinitializefornextADCsequenceAdcRegs.ADCTRL2.bit.RST_SEQ1=1; //ResetSEQ1AdcRegs.ADCST.bit.INT_SEQ1_CLR=1; //ClearINTSEQ1bitPieCtrlRegs.PIEACK.all=PIEACK_GROUP1; //Acknowledgeinterrupt//toPIEInitForFFT();MakeWave();for(i=0;i<SAMPLENUMBER;i++){fWaveR[i]=INPUT[i];fWaveI[i]=0.0f;w[i]=0.0f;}FFT(fWaveR,fWaveI);for(i=0;i<SAMPLENUMBER;i++){DATA[i]=w[i];}while(1);//breakpointreturn;}voidInitForFFT(){inti;for(i=0;i<SAMPLENUMBER;i++)12{sin_tab[i]=sin(PI*2*i/SAMPLENUMBER);cos_tab[i]=cos(PI*2*i/SAMPLENUMBER);}}voidFFT(floatdataR[SAMPLENUMBER],floatdataI[SAMPLENUMBER]){intx0,x1,x2,x3,x4,x5,x6,xx;inti,j,k,b,p,L;floatTR,TI,temp;/**********followingcodeinvertsequence************/for(i=0;i<SAMPLENUMBER;i++){x0=x1=x2=x3=x4=x5=x6=0;x0=i&0x01;x1=(i/2)&0x01;x2=(i/4)&0x01;x3=(i/8)&0x01;x4=(i/16)&0x01;x5=(i/32)&0x01;x6=(i/64)&0x01;xx=x0*64+x1*32+x2*16+x3*8+x4*4+x5*2+x6;dataI[xx]=dataR[i];}for(i=0;i<SAMPLENUMBER;i++){dataR[i]=dataI[i];dataI[i]=0;}/**************followingcodeFFT*******************/for(L=1;L<=7;L++){/*for(1)*/b=1;i=L-1;while(i>0){b=b*2;i--;}/*b=2^(L-1)*/for(j=0;j<=b-1;j++)/*for(2)*/{p=1;i=7-L;while(i>0)/*p=pow(2,7-L)*j;*/{p=p*2;i--;}p=p*j;for(k=j;k<128;k=k+2*b)/*for(3)*/{TR=dataR[k];TI=dataI[k];temp=dataR[k+b];13dataR[k]=dataR[k]+dataR[k+b]*cos_tab[p]+dataI[k+b]*sin_tab[p];dataI[k]=dataI[k]-dataR[k+b]*sin_tab[p]+dataI[k+b]*cos_tab[p];dataR[k+b]=TR-dataR[k+b]*cos_tab[p]-dataI[k+b]*sin_tab[p];dataI[k+b]=TI+temp*sin_tab[p]-dataI[k+b]*cos_tab[p];}/*ENDfor(3)*/}/*ENDfor(2)*/}/*ENDfor(1)*/for(i=0;i<SAMPLENUMBER/2;i++