《数字信号处理实验》由会员分享,可在线阅读,更多相关《数字信号处理实验(11页珍藏版)》请在金锄头文库上搜索。
1、实验一实验一自适应滤波器自适应滤波器一、实验目的一、实验目的1、掌握功率谱估计方法2、会用 matlab 对功率谱进行仿真二、实验原理二、实验原理功率谱估计方法有很多种,一般分成两大类,一类是经典谱估计;另一类是现代谱估计。 经典谱估计可以分成两种,一种是BT 法,另一种是周期法;BT 法是先估计自相关函数,然后将相关函数进行傅里叶变换得到功率谱函数。相应公式如下所示:1xx(m ) rNPBTm N |m | 1n 0x*(n ) x (n m )(1 1)(1 2)xx(m )e jwnr周期图法是采用功率谱的另一种定义,但与 BT 法是等价的,相应的功率谱估计如下所示:(ejw)1Pxx
2、N其计算框图如下所示:观测数据x(n)FFTjwnx(n)en0N120nN1(13)取模的平方1/NPxx(ejw)1由于观测数据有限,所以周期图法估计分辨率低,估计误差大。针对经典谱估计的缺点,一般有三种改良方法:平均周期图法、窗函数法和修正的周期图平均法。三、实验要求三、实验要求信号是正弦波加正态零均值白噪声,信噪比为10dB,信号频率为2kHZ,取样频率为 100kHZ。四、实验程序与实验结果四、实验程序与实验结果1用周期图法进行谱估计A、实验程序:%用周期法进行谱估计clear all;N1=128;%数据长度N2=256;N3=512;N4=1024;f=2;%正弦波频率,单位为
3、kHZfs=100;%抽样频率,单位为 kHZn1=0:N1-1;n2=0:N2-1;n3=0:N3-1;n4=0:N4-1;a=sqrt(20);%由信噪比为 10dB 计算正弦信号的幅度2wn1=randn(1,N1);xn1=a*sin(2*pi*f*n1./fs)+wn1;Pxx1=10*log10(abs(fft(xn1).2)/N1);%周期法求功率谱f1=(0:length(Pxx1)-1)/length(Pxx1);wn2=randn(1,N2);xn2=a*sin(2*pi*f*n2./fs)+wn2;Pxx2=10*log10(abs(fft(xn2).2)/N2);f2=
4、(0:length(Pxx2)-1)/length(Pxx2);wn3=randn(1,N3);xn3=a*sin(2*pi*f*n3./fs)+wn3;Pxx3=10*log10(abs(fft(xn3).2)/N3);f3=(0:length(Pxx3)-1)/length(Pxx3);wn4=randn(1,N4);xn4=a*sin(2*pi*f*n4./fs)+wn4;Pxx4=10*log10(abs(fft(xn4).2)/N4);f4=(0:length(Pxx4)-1)/length(Pxx4);subplot(2,2,1);plot(f1,Pxx1);xlabel(频率);
5、ylabel(功率(dB);title(功率谱 Pxx,N=128);subplot(2,2,2);plot(f2,Pxx2);xlabel(频率);ylabel(功率(dB);title(功率谱 Pxx,N=256);subplot(2,2,3);plot(f3,Pxx3);xlabel(频率);ylabel(功率(dB);title(功率谱 Pxx,N=512);subplot(2,2,4);3plot(f4,Pxx4);xlabel(频率);ylabel(功率(dB);title(功率谱 Pxx,N=1024);B、实验仿真结果:2采用汉明窗,分段长度 L=32,用修正的周期图求平均法进
6、行谱估计A:实验程序:clear all;N=512;%数据长度Ns=32;%分段长度f1=2;%正弦波频率,单位为 kHZfs=100;%抽样频率,单位为 kHZn=0:N-1;4a=sqrt(20);%由信噪比为 10dB 计算正弦信号的幅度wn=randn(1,N);xn=a*sin(2*pi*f1*n./fs)+wn;w=hamming(32);%汉明窗Pxx1=abs(fft(w.*xn(1:32),Ns).2)/norm(w)2;Pxx2=abs(fft(w.*xn(33:64),Ns).2)/norm(w)2;Pxx3=abs(fft(w.*xn(65:96),Ns).2)/no
7、rm(w)2;Pxx4=abs(fft(w.*xn(97:128),Ns).2)/norm(w)2;Pxx5=abs(fft(w.*xn(129:160),Ns).2)/norm(w)2;Pxx6=abs(fft(w.*xn(161:192),Ns).2)/norm(w)2;Pxx7=abs(fft(w.*xn(193:224),Ns).2)/norm(w)2;Pxx8=abs(fft(w.*xn(225:256),Ns).2)/norm(w)2;Pxx9=abs(fft(w.*xn(257:288),Ns).2)/norm(w)2;Pxx10=abs(fft(w.*xn(289:320),N
8、s).2)/norm(w)2;Pxx11=abs(fft(w.*xn(321:352),Ns).2)/norm(w)2;Pxx12=abs(fft(w.*xn(353:384),Ns).2)/norm(w)2;Pxx13=abs(fft(w.*xn(385:416),Ns).2)/norm(w)2;Pxx14=abs(fft(w.*xn(417:448),Ns).2)/norm(w)2;Pxx15=abs(fft(w.*xn(449:480),Ns).2)/norm(w)2;Pxx16=abs(fft(w.*xn(481:512),Ns).2)/norm(w)2;Pxx=10*log10(Px
9、x1+Pxx2+Pxx3+Pxx4+Pxx5+Pxx6+Pxx7+Pxx8+Pxx9+Pxx10+Pxx11+Pxx12+Pxx13+Pxx14+Pxx15+Pxx16)/16);5f=(0:length(Pxx)-1)/length(Pxx);plot(f,Pxx);xlabel(频率);ylabel(功率(dB);title(加窗平均周期图法功率谱 Pxx,N=512);grid on;B:实验仿真结果:五参考文献五参考文献:1 丁玉美,阔永红,高新波.数字信号处理-时域离散随机信号处理M.西安:西安电子科技大学出版社,2002.2 万建伟,王玲.信号处理仿真技术M.长沙:国防科技大学出版
10、社,2008.6实验二实验二卡尔曼滤波器的设计卡尔曼滤波器的设计一实验目的一实验目的1.熟悉并掌握卡尔曼滤波、自适应滤波和谱估计的原理。2.可以仿真符合要求的卡尔曼滤波器、 自适应滤波器和各种谱估计方法。3.掌握卡尔曼滤波器的递推公式和仿真方法。4.熟悉 matlab 的用法。二实验原理二实验原理卡尔曼滤波是用状态空间法描述系统的, 由状态方程和测量方程所组成。 卡尔曼滤波用前一个状态的估计值和最近一个观测数据来估计状态变量的当前值,并以状态变量的估计值的形式给出。 其状态方程和量测方程如下所示:xk 1Akxk wkyk Ckxk vk(1 1 )(1 2)其中,k表示时间,输入信号wk是一
11、白噪声,输出信号的观测噪声vk也是一个白噪声,输入信号到状态变量的支路增益等于1,即 B=1;A表示状态变量之间的增益矩阵,可随时间变化,A表示第kk次迭代的取值,C表示状态变量与输出信号之间的增益矩阵,可随时间变化,其信号模型如图 1.1 所示(k用k1代替)。7vkCk+ykWk1+XkZ1Xk1Ak1图 1.1 卡尔曼滤波器的信号模型卡尔曼滤波是采用递推的算法实现的, 其基本思想是先不考虑输入信号wk和观测噪声vk的影响, 得到状态变量和输出信号的估计值,再用输出信号的估计误差加权矫正状态变量的估计值, 使状态变量估计误差的均方值最小。其递推公式如下所示: x ke0.02x k1Hk(
12、yke0.02xk1)1H P(P 1)kkkPke0.04Pk11e0.04Pk(IHk)Pk(112a)(112b)(112c)(112d)假 设 初 始 条 件Ak,Ck,Qk,Rk, yk,xk1,Pk1已 知 , 其 中x0 Ex0,P0 varx0,那么递推流程见图 1.2 所示。Pk1式式(1-5)(1-5)Pk式式(1-4)(1-4)Hk式式(1-3)(1-3)式式(1-6)(1-6)图 1.2卡尔曼滤波递推流程图三实验要求三实验要求一连续平稳的随机信号 x(t),自相关rx() e8xkPk,信号 x(t)为加性噪声所干扰,噪声是白噪声,测量值的离散值y(k)为已知。Matl
13、ab 仿真程序如下:%编卡尔曼滤波递推程序,估计信号 x(t)的波形clear all;clc;Ak=exp(-0.02); %各系数由前面确定;Ck=1; Rk=0.1; p(1)=20; %各初值;Qk=1-exp(-0.04);p1(1)=Ak*p(1)*Ak+Qk; %由 p1 代表 p;x(1)=0; %设信号初值为 0;H(1)=p1(1)*Ck*inv(Ck*p1(1)*Ck+Rk);zk=-3.2,-0.8,-14,-16,-17,-18,-3.3,-2.4,-18,-0.3,-0.4,-0.8,-19,-2.0,-1.2,-11,-14,-0.9,0.8,10,0.2,0.5
14、,-0.5,2.4,-0.5,0.5,-13,0.5,10,-12,0.5,-0.6,-15,-0.7,15,0.5,-0.7,-2.0,-19,-17,-11,-14 %zk为测量出来的离散值;N=length(zk); %要测量的点数;for k=2:Np1(k)=Ak*p(k-1)*Ak+Qk; %未考虑噪声时的均方误差阵;H(k)=p1(k)*Ck*inv(Ck*p1(k)*Ck+Rk); %增益方程;I=eye(size(H(k); %产生和 Hk维数相同的单位矩阵;p(k)=(I-H(k)*Ck)*p1(k); %滤波的均方误差阵;x(k)=Ak*x(k-1)+H(k)*(zk(k
15、)-Ck*Ak*x(k-1); %递推公式;end,x %显示信号 x(k)的数据;9m=1:N;n=m*0.02;plot(n,zk,-r*,n,x,-bo); %便于比较 zk 和 x(k)在同一窗口输出;xlabel(t/s,Fontsize,16);ylabel(z(t),x(t),fontsize,16);title(卡尔曼滤波递推x(t)的估计波形与 z(t)波形,fontsize,16)legend(观测数据 z(t),信号估计值 x(t),2);grid;四实验结果四实验结果五实验小结五实验小结通过卡尔曼滤波估计信号与观测信号比较知,卡尔曼滤波输出的估计信号x(t)与实际观测到的离散值z(t)还是存在一定的误差,卡尔曼滤波是从初始状态1 0就采用递推方法进行滤波, 那么在初值迭代后的一段时间内可能会出现较大的误差,随着迭代进行,各参数逐渐趋于稳定,后面的估计值与观察值的误差就减少了。六参考文献六参考文献1 丁玉美,阔永红,高新波.数字信号处理-时域离散随机信号处理M.西安:西安电子科技大学出版社,2002.2 万建伟,王玲.信号处理仿真技术M.长沙:国防科技大学出版社,2008.1 1