VC++编程实现对波形数据的频谱分析.docx

VC++编程实现对波形数据的频谱分析.docx

《VC++编程实现对波形数据的频谱分析.docx》由会员分享，可在线阅读，更多相关《VC++编程实现对波形数据的频谱分析.docx（19页珍藏版）》请在冰豆网上搜索。

VC++编程实现对波形数据的频谱分析

VC++编程实现对波形数据的频谱分析

-

　  摘要：

本文介绍了采用离散傅立叶变换（DFT）实现对采样得到的波形数据文件进行频谱分析的一般方法，并且为了提高运算效率、节省中间存储单元，最终采用了"时间抽选奇偶分解"的"库利-图基算法"实现快速离散傅立叶变换，对采样数据进行了高效的频谱分析，并用MicrosoftVisualC++6.0编写实现。

　　关键字：

MicrosoftVisualC++6.0、离散傅立叶变换、快速傅立叶变换、采样

　　一、引言

　　频谱分析是电子工程上一个非常重要的手段，许多计算机辅助电路分析（CAA）类软件都具备这种分析能力，以便电子工程师能清楚的看到某波形的频谱分布情况。

而要对一个输入信号源作频谱分析，将其由时域信号转变为频域信号，就必然要用到傅立叶分析，而无论是在时域还是在频域，都要对连续函数进行积分运算。

很显然，要通过计算机实现此变换必须预先通过抽样将原始的连续数据转变为离散数据，并将计算范围收缩到一个有限区间。

因此在允许一定程度近似的条件下，可以使用"离散傅立叶变换（DFT）"对波形数据进行频谱分析。

　　二、快速傅立叶变换（FFT）算法构成原理

　　要计算一个N点的离散傅立叶变换需要同一个N*N点的W矩阵（关于W矩阵请参阅信号与系统方面的书籍）相运算，随着N值的增大，运算次数显著上升，当点数达到1024时，需要进行复数乘法运算1，048，576次，显然这种算法在实际运用中无法保证当点数较大时的运算速度，无法满足对信号的实时处理。

　　根据W矩阵中W元素的周期性和对称性我们可以将一个N点的DFT运算分解为两组N/2点的DFT运算，然后取和即可，为进一步提高效率，将上述两个矩阵按奇偶顺序逐级分解下去。

当采样点数为2的指数次方M时，可分解为M级子矩阵运算，全部工作量仅为：

　　复数乘法：

M*N/2次

　　复数加法：

N*M次

　　而直接DFT需要的运算量为：

　　复数乘法：

N*N次

　　复数加法：

N*（N-1）次

　　当点数N为几十个点时FFT的优势还不明显，而一旦达到几千、几百个点时优势是十分明显的：

　　N=1024时：

DFT需1048576次运算，FFT仅需5120次运算，改善比204.8。

　　N=2048时：

DFT需4194304次运算，FFT仅需11264次运算，改善比达到372.4。

　　三、"时间抽选奇偶分解快速离散傅立叶变换"的程序实现

　　当采样点数较多时，如变换前和变换后的序列都按自然顺序排列，则中间运算过程会占用大量的中间存储单元，造成效率的低下和存储单元的浪费。

根据FFT的实现原理我们可以对采样序列进行逐次奇偶抽选，打乱以前的次序重新排序，然后按此顺序参加运算，可以实现"即位运算"提高存储单元的利用率。

（一）复数的描述方法

　　进行傅立叶变换时不可避免的要用到复数，而在VC中并没有现成的可用于表示复数的数据类型，可以自己定义一个含有两个成员变量的数据结构来表示复数，这两个成员变量可分别用于表示复数的实部与虚部：

typedefstructtagComplex{

　floatRe;//复数的实部

　floatIm;//复数的虚部

}Complex;

（二）倒序的实现

　　在进行快速傅立叶变换时，可以将输入的时域序列和输出的频域序列都按照自然顺序排列；也可以按照"蝴蝶图"所描述的计算方法对输入的时域序列按奇偶分解后的序列排序而输出的频域序列仍是按自然顺序排列的；还有一中方式是输入的时域序列是不进行抽选的自然序列，而输出的频域序列则是按奇偶分解后的顺序排列的。

这三种方式各有优点，第一种对输入、输出不需要进一步排序，但由于自然排序不符合"蝴蝶图"运算规律，会占用大量中间存储单元。

而后两种则无须中间存储单元，但需要倒一次序。

权衡利弊，当采样点较多时还是采用后两种方式好，多一次倒序运算对现在的高性能计算机而言并不是什么负担。

下面代码用于对原始采样序列的时间抽选奇偶分解工作，其中A、N分别表示指向采样序列复数数组的指针和序列的长度。

intNV2=N/2;

intNM1=N-1;

intI,J,K=0;

ComplexT;//用于中介的复数变量T

I=J=1;

while（I<=NM1）

{

　if（I

　{

　　T=A[J-1];//将A[J-1]的内容和A[I-1]的内容互换，借助于中间变量T

　　A[J-1]=A[I-1];

　　A[I-1]=T;

　}

　K=NV2;

　while（K

　{

　　J-=K;

　　K/=2;

　}

　J+=K;

　I++;

}

　　（三）时域信号的频谱分析

　　首先要将从外设输入或采集的时域波形数据经抽样量化后，通过CFile类的Open（……）、Read（……）等成员函数将其读取到缓存中，并将其转化为复变量存放于复变量数组A中，同时需要验证以下数据量的长度是否为2的整数次幂，如若不是则必须用0来补齐，否则无法用"蝴蝶图"进行分解运算。

下面代码用于完成对原始采样时域序列的快速傅立叶变换，A、M分别表示指向原始采样数据数组的指针和序列长度的2的整数次幂：

……

ComplexU,W,T;

intLE,LE1,I,J,IP;

intN=（int）pow（2,M）;

//在此采用的是时间抽选奇偶分解方式，所以在参加运算前首先要对时间序列进行倒序

ReverseOrder（A,N）;

intL=1;

while（L<=M）

{

　LE=（int）pow（2,L）;

　LE1=LE/2;

　U.Re=1.0f;

　U.Im=0.0f;

　W.Re=（float）cos（PI/（1.0*LE1））;//计算W算子的值

　W.Im=（float）-1.0*sin（PI/（1.0*LE1））;

　if（abs（W.Re）<1.0e-12）

　　W.Re=0.0f;

　if（abs（W.Im）<1.0e-12）

　　W.Im=0.0f;

　　J=1;

　　while（J<=LE1）

　　{

　　　I=J;

　　　while（I<=N）

　　　{

　　　　IP=I+LE1;

　　　　T.Re=（float）A[IP-1].Re*U.Re-A[IP-1].Im*U.Im;//计算复数运算A*U

　　　　T.Im=（float）A[IP-1].Re*U.Im+A[IP-1].Im*U.Re;

　　　　A[IP-1].Re=（float）A[I-1].Re-T.Re;//计算复数运算A-T

　　　　A[IP-1].Im=（float）A[I-1].Im-T.Im;

　　　　A[I-1].Re+=T.Re;//计算复数运算A+T

　　　　A[I-1].Im+=T.Im;

　　　　I+=LE;

　　　}

　　　floattemp=U.Re;

　　　U.Re=（float）U.Re*W.Re-U.Im*W.Im;//计算复数运算U*W

　　　U.Im=（float）temp*W.Im+U.Im*W.Re;

　　　J++;

　　}

　　L++;

　}

……

　　上述代码执行完毕时，原先存放着时域数值的复变量数组内存放的就是经过分析后的频域值了，对此数据可以通过绘图将频域波形直观的显示出来，也可以将其存成数据文件，以备进一步使用。

　　四、测试及运算结果分析

　　编译运行程序，打开一三角脉冲的数据文件，并将分析结果保存，该三角脉冲幅度为1，持续时间2毫秒，采样时抽样时间间隔是20微秒，延拓周期（数据记录长度）为10毫秒，采样点数目500点，取2的整数次幂512个样点。

下附该三角脉冲频谱的计算结果及误差分析：

　　频率（Hz）FFT求得　　　　X（f）　　　　误差

　　0.00　　1.00006E-03　1.00000E-03　6.10352E-08

　　100.00　9.67593E-04　9.67531E-04　6.14332E-08

　　200.00　8.75203E-04　8.75150E-04　6.25092E-08

　　300.00　7.36904E-04　7.36849E-04　6.39413E-08

　　400.00　5.72852E-04　5.72787E-04　6.52926E-08

　　500.00　4.05351E-04　4.05285E-04　6.61362E-08

　　600.00　2.54638E-04　2.54572E-04　6.61847E-08

　　700.00　1.35403E-04　1.35338E-04　6.53870E-08

　　800.00　5.47602E-05　5.46963E-05　6.39612E-08

　　900.00　1.20072E-05　1.19448E-05　6.23453E-08

　　1000.00　6.10719E-08　1.17757E-32　6.53870E-08

　　1100.00　8.05672E-06　7.99613E-06　6.05985E-08

　　1200.00　2.43706E-05　2.43095E-05　6.11450E-08

　　1300.00　3.93026E-05　3.92400E-05　6.25965E-08

　　1400.00　4.68226E-05　4.67581E-05　6.45128E-08

　　1500.00　4.50979E-05　4.50316E-05　6.62543E-08

　　1600.00　3.58664E-05　3.57992E-05　6.71930E-08

　　1700.00　2.30135E-05　2.29466E-05　6.69399E-08

　　1800.00　1.08697E-05　1.08042E-05　6.55073E-08

　　1900.00　2.74348E-06　2.68014E-05　6.33390E-08

　　2000.00　6.11826E-08　1.17757E-32　6.11826E-08

　　2100.00　2.25379E-06　2.19395E-06　5.98376E-08

　　2200.00　7.29243E-06　7.23256E-06　5.98625E-08

　　2300.00　1.25974E-05　1.25360E-05　6.13467E-08

　　2400.00　1.59746E-05　1.59107E-05　6.38421E-08

　　2500.00　1.62779E-05　1.62114E-05　6.64915E-08

　　2600.00　1.36254E-05　1.35571E-05　6.83226E-08

　　2700.00　9.16539E-06　9.09679E-06　6.86075E-08

　　2800.00　4.53216E-06　4.46500E-06　6.71550E-08

　　2900.00　1.21487E-06　1.15945E-06　6.44190E-08

　　注：

在此，FFT运算结果都倍乘了系数10毫秒（0.01秒）。

　　在分析结果中产生了误差，是由于待分析的连续时间信号不具备离散性或周期性，也可能有无限长度。

为了适应FFT方法的需要，对波形进行了抽样和截断，这样再用程序分析采样数据必然会引入误差，从分析结果可以看出，频率越高，误差波动也越大，此分析结果产生的误差在允许范围之内，是一个可以满意的近似。

实践证明，本程序的算法是正确可靠的。

　　小结：

　　DFT尤其是FFT的应用已遍及各个科学领域，"DFT的应用"与"FFT的应用"几乎成为同义语。

通过本文介绍和程序示例可以清楚的看到FFT方法在直接处理离散信号数据的作用，而且也可以很好的用于对连续时间信号分析的逼近。

本程序在Windows2000Professional下、由MicrosoftVisualC++6.0编译通过

用VB实现数字波形显示程序减小字体增大字体作者：

佚名来源：

本站整理发布时间：

2009-07-1518:

03:

17

id=126063>

　　摘要：

本文详细介绍了在VB集成环境下数字波形高速显示的方法，同时对双通道波形显示和数字滤波方法也进行了介绍。

　　关键词：

数字；波形；显示；滤波

　　1前言：

　　随着计算机技术及电子技术的发展，数字采集技术在检测领域的应用越来越广泛，检测速度越来越高，检测的数据量越来越大，特别是在无损检测领域，将检测数据通过计算机处理后绘制出波形，并实时显示，对及时发现伤损、分析伤损具有重要意义。

　　2波形显示

　　检测数据通常是离散的数据，将离散的数据绘制出波形，可通过在两点间连接线段的方法实现。

　　2.1用Line方法显示波形

　　VB提供了Line画直线方法，可在窗体上增加一个图片框控件，适当设置图片的大小和背景颜色，用Line方法将离散数据按检测顺序连接成线段，即可将波形显示在图片框中。

但该方法显示波形速度较慢，不适合高速显示的应用。

　　2.2WindowsAPI函数显示波形

　　在VB中两点间连线的另一种方法是用WindowsAPI函数，Win32API提供了以下两个函数，联合使用可实现波形的快速显示，经过测试，显示速度比使用Line方法快70%以上。

　　LineTo函数：

　　函数功能：

画出由数组定义的点连接的一系列线段。

　　函数原型：

BOOLLineTo（HDChdc,intnXEnd,intnYEnd）；

　　参数：

　　hdc:

设备环境句柄。

　　nXEnd:

定义线段终点的X坐标。

　　nYEnd:

定义线段终点的Y坐标。

　　返回值：

若函数调用成功，则返回非0值；若函数调用失败，则返回值为0。

　　MoveToEx函数：

　　函数功能：

将当前位置更新为指定的点，并有选择的返回原先的位置。

　　函数原型：

BOOLMoveToEx（HDChdc,intX,intY,LPPOINTlpPoint）；

　　参数：

　　hdc:

设备环境句柄。

　　X:

定义新位置的X坐标（逻辑坐标）。

　　Y:

定义新位置的Y坐标（逻辑坐标）。

　　lpPoint:

指向一个POINT结构，结构中存放原先的位置。

若此参数为NULL,则不返回原先的位置

　　返回值：

若函数调用成功，则返回非0值；若函数调用失败，则返回值为0。

　　在连接线段时，首先将检测数据放入一个数组中，用MoveToEx函数定位画线的起始点坐标，然后用LineTo函数画出起始点至下一个点之间的线段，再用MoveToEx将画线的起始点定位到下一个点，继续用LineTo函数画线，如此循环，即可将离散点连接成波形。

　　例：

　　zz=MoveToEx（Picture1.hdc,i,Mwave（i）,LpPoint1）

　　zz=LineTo（Picture1.hdc,i,Mwave（i+1））

实时波形显示界面

   通常计算机需要接收外部实时发送的数据并用十分形象的方式显示出来。

例如柱状图、饼图等等。

本应用程序则采用波形的形式显示，并接将之设计为可以接收多路数据的波形显示界面。

本应用程序在VS2008环境下调试通过，源码下载连接如下

一.程序界面

点击该图放大

上图中显示了两路波形即三角波和正弦波，当然这两路波形是由程序计算出来的并不是从外部接收的。

实际工作中则可以配合串口通讯设备接收它发来的数据并显示出来。

二.波形控件类介绍

本程序的实现主要依赖于那个波形显示控件。

从下载连接那里可以下载该源码，里面的文件中，2DPushGraph.h和2DPushGraph.cpp即是该控件的类的定义文件和实现文件。

分析之后得到该类的一些信息。

其成员函数包括：

程序代码：

COLORREFm_crTextColor;//标签文字颜色

   COLORREFm_crBGColor;//背景颜色

   COLORREFm_crGridColor;//栅格颜色

   CString  m_strMaxLabel;//最大值之处的标签

   CString  m_strMinLabel;//最小值之处的标签

   bool   m_bShowMinMax;//最小最大值显示

   bool   m_bShowGrid;//是否显示栅格

   bool   m_bStylesModified;//是否被修改

   int   m_nMoveOffset;//偏移

   int   m_nMaxCoords;//最大缓冲

   int   m_nMaxPeek;//显示数据的最大值

   int   m_nMinPeek;//显示数据的最小值

   int   m_nGridSize;//栅格间距

   int   m_nPeekOffset;//峰值偏移

其实现方法包括：

程序代码：

boolCreateFromStatic（UINTnStaticID,CWnd*pParent）;//该方法用于创建控件实例，通常在对话框初始化函数中调用，nStaticID为控件ID，pPatent为创建在那个窗口中的ID。

   LPCTSTR   GetLabelForMax（）const;//获取最大值处的标签

   LPCTSTR   GetLabelForMin（）const;//获取最大值处的标签

   COLORREF  GetBGColor（）const;//获得背景色   

   COLORREF  GetGridColor（）const;//获得栅格颜色   

   voidSetBGColor（COLORREFcrColor）;

   voidSetGridColor（COLORREFcrColor）;

   voidSetTextColor（COLORREFcrColor）;

   COLORREF  GetLineColor（UINTuiLineID）;//该方法用于获得ID号为uiLineID的波形的线条颜色，因为要在波形控件中画线必须先调用boolAddLine（UINTuiLineID,COLORREFcrColor）;创建一个波形，ID号为uiLineID。

   boolSetLineColor（COLORREFcrColor,UINTuiLineID）;

   COLORREF  GetTextColor（）const;//获得标签文本的颜色

   int  GetGridSize（）  const;//获得栅格间隔

   voidSetGridSize（unsignedshortusWidthAndHeight）;

   int  GetMaxPeek（）   const;//获得最大值纵坐标

   int  GetMinPeek（）   const;//获得最小值纵坐标   

   voidSetPeekRange（intnMin,intnMax）;

   voidSetMaxPeek（intnMax）;

   voidSetMinPeek（intnMin）;

   unsignedshortGetInterval（）const;   //获得横坐标间距

   voidSetInterval（unsignedshortusInterval）;

   voidSetLabelForMax（LPCTSTRlpszLabel）;

   voidSetLabelForMin（LPCTSTRlpszLabel）;

   boolAddLine（UINTuiLineID,COLORREFcrColor）;//该方法用于在该控件上创建一个波形。

   voidRemoveLine（UINTuiLineID）;//清除ID号为uiLineID的波形

   boolPush（intnMagnitude,UINTuiLineID）;//在ID号为uiLineID号的波形上添加一个数据点，nMagnitude为幅度

   voidShowAsBar（UINTuiLineID,boolbAsBar）;//柱状图显示

   voidUpdate（）;//更新

三.应用编程

运用该控件在VS2008环境下的编程步骤如下：

1.建立一个对话框的MFC工程，在对话框上按照上图所示的界面布置控件。

其中波形控件那里布置一个PictureControl控件将其ModalFrame和Type均属性设置为true，其他均设置为False。

注意给PictureControl取的ID!

后面编程将会用到。

2.将波形控件类的定义文件和实现文件拷贝至你的工程目录下。

但这实际上并没有将该类真正添加到你的工程下，需手动添加类。

常规操作，不详述。

3.在对话框的定义和实现文件中分别添加如下代码：

程序代码：

#include"2DPushGraph.h"

4.在对话框定义文件中（我给的供下载的例程中的是DataRealTimeDlg.h这个文件）中定义一个该控件类的变量：

程序代码：

private:

        C2DPushGraphm_PushGraph;

5.在对话框的实现文件中（我给的供下载的例程中的是DataRealTimeDlg.cpp这个文件）的对话框初始化函数中添加如下代码：

程序代码：

   m_PushGraph.CreateFromStatic（IDC_REALCTRL,this）;   //这个IDC_REALCTRL即是那个PictureControl控件的ID号。

   m_PushGraph.ModifyStyle（0,WS_THICKFRAME）;   //设置风格

6.现在不妨试着运行以下该程序，应该可以观察到那个控件显示出来了吧。

在添加以下代码：

程序代码：

   m_PushGraph.AddLine（m_sin,  RGB（255,255,255））;

   m_PushGraph.AddLine（m_tra,  RGB（255,0,0））

这两行代码分别添加了一个正弦波形，ID号为m_sin，一个三角波形，ID号为m_tra。

不过运行之后并没有数据点绘制出来。

那是应为还没有调用boolPush（intnMagnitude,UINTuiLineID）函数添加数据点。

试着添加几个数据点再运行即可观察到波形。

7.还有一些控件的响应代码看看那个下载的例程吧。

四.结束语

   成功了吧，AnyProblem，Contactmeplease！

控件设计]Nobi'sStatusChart-野