测试信号分析与处理作业实验一二

王锋实验一：利用FFT作快速相关估计一、实验目的a.掌握信号处理的一般方法，了解相关估计在信号分析与处理中的作用。b.熟悉FFT算法程序；熟练掌握用FFT作快速相关估计的算法。c.了解快速相关估计的谱分布的情况。二、实验内容a.读入实验数据[1]。b.编写一利用FFT作相关估计的程序[2]。c.将计算结果表示成图形的形式，给出相关谱的分布情况图。注[1]：实验数据文件名为“Qjt.dat”。实验数据来源：三峡前期工程“覃家沱大桥”实测桥梁振动数据。实验数据采样频率：50Hz。可从数据文件中任意截取几段数据进行分析，数据长度N自定。注[2]：采用Matlab编程。三、算法讨论及分析算法为有偏估计，利用FFT计算相关函数Step1:对原序列补N个零，得新序列x2N(n)Step2:作FFT[x2N(n)]得到X2N(k)2,0,1,......,212kkkNNStep3:取X2N(k)的共轭2()NXk，得221()NXkNStep4:作2201ˆ()()NIFFTXkRmNStep5:调整0ˆ()Rm与ˆ()Rxxm的错位。四、实验结果分析1.该信号可以近似为平稳信号么？可以近似为平稳信号，随机过程的统计特性不随样本的采样时刻而发生变化。取N=8192，分别取间隔m=500,m=700,m=1000,所得到的均值均为0.5366，方差为47369，与时间无关。图1-1自相关函数图(上图表示的R0，下图为调整后的R0)2.该信号是否具有周期性，信噪比如何？loadQjt.dat;%加载数据N=32768;%数据长度i=1:1:N;%提取数据plot(i,Qjt(i));抛去几个极值点，从图1-2可以看出，数据具有一定的周期性，杂音比较少，说明信噪比较高。图1-2数据图3.数据量N对估计结果有何影响？在实验中，用FFT做相关估计101ˆ()()()NmxxnRmxnxnmN，有偏估计均值的表达式为ˆ()()xxxxNmERmRmN有偏估计方差的表达式为22ˆ()1()()()xxxxxxxxRmiRiRimRimN可知，N越大，方差越小，均值越趋于()xxRm，ˆ()xxRm为无偏渐进一致估计，估计效果好。4.采用有偏估计和无偏估计的结果有何不同？哪一种结果更好？无偏估计的均值为101ˆ()()NmxxxxnERmRmNm无偏估计的方差为22ˆ()()()xxxxxxRmiNRiRimNm可以看出均值为无偏估计，但是当m增大时，方差大，估计值分散，综合质量差，有偏估计的方差就没有这样的问题，估计值比较集中，故有偏估计更适合实际应用。五、原程序清单%%计算数据的自相关函数的程序loadQjt.dat;%加载数据N=16384;%数据长度%N=8192;fori=1:1:N%提取数据X(i)=Qjt(i);endnum=2*N;%补零措施fori=N+1:1:numX(i)=0;endXfft=fft(X);%x2N(n)进行FFT变换X2N(k)conXfft=conj(Xfft);%求X2N(k)的共轭R0=ifft((Xfft.*conXfft)/N);%求自相关函数Rxx=0*R0;fori1=1:1:NRxx(i1)=R0(i1+N);Rxx(i1+N)=R0(i1);endn=1:1:num;subplot(2,1,1);stem(n,R0);%绘制R0的图subplot(2,1,2);axisonstem(fftshift(R0));%调整后的R0%%计算数据的均值和均方差的程序loadQjt.dat;%加载数据%N=16384;N=8192;%选取数据量fori=1:1:NX(i)=Qjt(i);%选取数据，数据范围1~8192endold_mean=mean(X)%求数据的均值old_var=var(X)%求数据的方差%选取数据，与原数据的间隔为m=500,数据范围500~8692fori=500:1:(N+500)X(i)=Qjt(i);endnew500_mean=mean(X)%求数据的均值new500_var=var(X)%求数据的方差%选取数据，与原数据的间隔为m=700,数据范围700~8892fori=700:1:(N+700)X(i)=Qjt(i);endnew700_mean=mean(X)%求数据的均值new700_var=var(X)%求数据的方差%选取数据，与原数据的间隔为m=1000,数据范围1000~9192fori=1000:1:(N+1000)X(i)=Qjt(i);endnew1000_mean=mean(X)%求数据的均值new1000_var=var(X)%求数据的方差六、实验后的体会与建议通过该实验，熟练应用了FFT作快速相关估计的算法来解决实际中的问题，更加熟悉了MATLAB相关的操作。建议：课上可以多教授一些关于MATLAB的相关函数的用法，以及常见的误区。实验二：采用Welch法作谱估计一、实验目的a.掌握Welch算法的概念、应用及特点；b.了解谱估计在信号分析中的作用；c.能够利用Welch法对信号作谱估计，对信号的特点加以分析。二、实验内容a.读入实验数据。b.编写一利用Welch法作谱估计的算法程序。c.将计算结果表示成图形的形式，给出信号谱的分布情况图。三、算法讨论及分析四、实验结果及分析图2-1实验结果图1不同的分段点数对Welch谱估计结果有何影响？分段数K增大，每组数据量M减小，则会造成谱估计的偏差增大，曲线平滑；反之,曲线起伏激烈,谱线的尖峰保留。从图2-2中的K值大于图2-3的K值，图2-2中的曲线比图2-3的曲线平滑，但是图2-3的尖峰效果比图2-2好。图2-2K=N/1100时的频谱图图2-3K=N/3000时的频谱图2不同的数据重叠长度Welch对谱估计结果有何影响？如果数据交叠使用，分段数K比平均周期图法增大，则估计方差比平均周期图法小，数据交叠长度越大，分段数K越大，则估计方差会越小，图形会变得平滑，如图2-4所示，上图数据重叠度大，K大，则曲线相对平滑。从式（1）也可以看出，式（1）为谱估计的方差。24sinˆvar1sinxxxMSKM（1）图2-4使用不同的重叠度得到的谱估计图（上：重叠度为97%，下：重叠度为3%）3不同的窗函数对Welch谱估计结果有何影响？使用窗函数，会使Welch谱估计的图形变得平滑，但是，不同的窗函数不同，则其主瓣宽度不同。主瓣宽度越窄，对信号的主瓣频率辨识效果越好，分辨率较高。从图2-5可以看出，使用矩形窗的分辨率高于汉宁窗。图2-5使用不同的窗函数得到的谱估计（上：矩形窗，下：汉宁窗）4周期图法与Welch法的谱估计结果有何不同？周期图法本身的误差就很大，不是无偏估计，而且偏差比较大；welch方法应该会准确一些，它频谱分辨率可能低一些，但估计的方差会小些，所以功率谱看起来平滑一些，另外一个原因是welch方法你是加了窗的。两个估计都是跟准确值有偏差的，只是后者统计方差更小一些，相对更加可靠一点。五、原程序清单loadQjt.dat%加载数据N=length(Qjt);%获得数据的长度L=1800;%确定每段数据的长度K=N/L;%在总数据中共取K段，没有重合%K=2*N/L-1;%在总数据中共取K段，有重合K=fix(K);fornum=1:1:K%使用循环，提取每段数据fori=1:1:LA(i,num)=Qjt((num-1)*L+i);%每段数据占据A阵的一列，没有重合%A(i,num)=Qjt((num-1)*L/2+i);%每段数据占据A阵的一列，有重合endfori2=L+1:1:2048%使用循环，补零A(i2,num)=0;endendfornum=1:1:K%使用循环，对每段数据使用矩形窗截取tempA=A(:,num);%tempw=boxcar(2048);%获取窗函数fori3=1:1:2048%对每段数据使用矩形窗截取A(i3,num)=tempA(i3)*tempw(i3);endendfornum=1:1:K%对截取之后的数据进行FFT变换，并获得傅里叶变换后的数据模tempA=A(:,num);tempFFT=fft(tempA);%对每一段数据进行FFT变换fori3=1:1:2048%获取数据FFT变换后的模Sxi(i3,num)=tempFFT(i3)*conj(tempFFT(i3))/2048;endend%循环结束U=dot(boxcar(2048),boxcar(2048))/2048;%获得修正系数UtempSUM=0*ones(2048,1);fornum=1:1:KtempSUM=Sxi(:,num)+tempSUM;%获得数据FFT变换后的功率谱和endSxx=tempSUM/(K*U);%调整系数subplot(3,1,1);Sxxi=ones(1024,1);forinum=1:1:1024Sxxi(inum)=Sxx(inum);endplot(Sxx);%绘制0~2048点的图[Px,F]=pwelch(Qjt,boxcar(2048),0,2048);%使用pwelch函数计算subplot(3,1,2);plot(Sxxi);subplot(3,1,3);plot(Px);%绘制pwelch所得的图六、实验后的体会和建议通过此次实验，掌握了Welch算法的概念、应用及特点，更加熟悉了MATLAB相关的操作。建议：课上可以多教授一些关于MATLAB的相关函数的用法，以及常见的误区。