《现代信号处理大作业题目+答案.pdf》由会员分享，可在线阅读，更多相关《现代信号处理大作业题目+答案.pdf（17页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、研究生“现代信号处理”课程大型作业研究生“现代信号处理”课程大型作业（以下四个题目任选三题做）（以下四个题目任选三题做）1.请用多层感知器（MLP）神经网络误差反向传播（BP）算法实现异或问题（输入为，并画出学习曲线。其中，非线性函数X X 0 0;0 1;1 0;11T，要求可以判别输出为 0 或 1）采用 S 型 Logistic 函数。2.试用奇阶互补法设计两带滤波器组(高、低通互补)，进而实现四带滤波器组；并画出其频响。滤波器设计参数为：Fp1.7KHz,Fr2.3KHz,Fs8KHz,Armin70dB。3.根据现代数字信号处理（姚天任等，华中理工大学出版社，2001）第四章附录提供
2、的数据(pp.352-353)，试用如下方法估计其功率谱，并画出不同参数情况下的功率谱曲线：1）Levinson 算法2）Burg 算法3）ARMA 模型法4）MUSIC 算法4.图 1 为均衡带限信号所引起失真的横向或格型自适应均衡器(其中横向 FIR 系统长M=11),系统输入是取值为1 的随机序列x(n)，其均值为零；参考信号d(n)x(n 7)；信道具有脉冲响应：2(n2)1)n 1,2,31cos(h(n)2W0其它式中W用来控制信道的幅度失真(W=24,如取 W=2.9,3.1,3.3,3.5 等)，且信道受到均2值为零、方差v 0.001(相当于信噪比为 30dB)的高斯白噪声v
3、(n)的干扰。试比较基于下列几种算法的自适应均衡器在不同信道失真、不同噪声干扰下的收敛情况(对应于每一种情况,在同一坐标下画出其学习曲线):1）横向/格-梯型结构 LMS 算法2）横向/格-梯型结构 RLS 算法并分析其结果。图 1横向或格-梯型自适应均衡器参考文献参考文献1 姚天任,孙洪.现代数字信号处理M.武汉:华中理工大学出版社,20012 杨绿溪.现代数字信号处理M.北京:科学出版社,20073 S.K.Mitra.孙洪等译.数字信号处理基于计算机的方法(第三版)M.北京:电子工业出版社,20064S.Haykin,郑宝玉等译.自适应滤波器原理(第四版)M.北京:电子工业出版社,200
4、35J.G.Proakis,C.M.Rader,F.Y.Ling,etc.Algorithms for Statistical Signal Processing M.Beijing:Tsinghua University Press,2003一、请用多层感知器（请用多层感知器（MLPMLP）神经网络误差反向传播（）神经网络误差反向传播（BPBP）算法实现异或问题）算法实现异或问题（输入为（输入为X X0 0;0 1;1 0;11T，要求可以判别输出为，要求可以判别输出为 0 0 或或 1 1），并画出学习曲线。其，并画出学习曲线。其中，非线性函中，非线性函 数采用数采用 S S 型型 Log
5、isticLogistic 函数。函数。1、原理：反向传播(BP)算法：(1)、多层感知器的中间隐层不直接与外界连接，其误差无法估计。(2)、反向传播算法：从后向前（反向）逐层“传播”输出层的误差，以间接算出隐层误差。分两个阶段：正向过程：从输入层经隐层逐层正向计算各单元的输出反向过程：由输出层误差逐层反向计算隐层各单元的误差，并用此误差修正前层的权值。2、流程图：开始选择初始值前向计算求所有神经元的输出计算输出层j3、程序：%使用了 3 层结构，第二层隐藏层4 个单元。2，3 层都使用 Logisitic 函数。%训练 xor 数据。functionmlp()f=fopen(XOR.txt)
6、;A=fscanf(f,%g,3 inf);A=A;p=A(1:2,:);%训练输入数据t=A(3,:);%desire outtrain_num,input_scale=size(p);%规模fclose(f);accumulate_error=zeros(1,3001);alpha=0.5;%学习率threshold=0.005;%收敛条件 e2 thresholdwd1=0;wd2=0;bd1=0;bd2=0;circle_time=0;hidden_unitnum=4;%隐藏层的单元数w1=rand(hidden_unitnum,2);%4 个神经元，每个神经元接受 2 个输入w2=r
7、and(1,hidden_unitnum);%一个神经元，每个神经元接受 4 个输入b1=rand(hidden_unitnum,1);b2=rand(1,1);while 1temp=0;circle_time=circle_time+1;for i=1:train_num%前向传播a0=double(p(i,:);%第 i 行数据n1=w1*a0+b1;a1=Logistic(n1);%第一个的输出n2=w2*a1+b2;a2=Logistic(n2);%第二个的输出a=a2;%后向传播敏感性e=t(i,:)-a;accumulate_error(circle_time)=temp+abs
8、(e)2;temp=accumulate_error(circle_time);s2=F(a2)*e;%输出层 delta 值s1=F(a1)*w2*s2;%隐层 delta 值%修改权值wd1=alpha.*s1*a0;wd2=alpha.*s2*a1;w1=w1+wd1;w2=w2+wd2;bd1=alpha.*s1;bd2=alpha.*s2;b1=b1+bd1;b2=b2+bd2;end;%end of forif accumulate_error(circle_time)3001break;end;%end of ifend;%end of whileplot(accumulate_
9、error,m);grid;xlabel(学习次数)ylabel(误差)disp(计算误差=,num2str(accumulate_error(circle_time);disp(迭代次数=,num2str(circle_time);%测试a0=double(0 0);n1=w1*a0+b1;a1=Logistic(n1);n2=w2*a1+b2;a2=Logistic(n2);a=a2;disp(0 0=,num2str(a);a0=double(0 1);n1=w1*a0+b1;%thena1=Logistic(n1);n2=w2*a1+b2;a2=Logistic(n2);a=a2;di
10、sp(0 1=,num2str(a);a0=double(1 0);n1=w1*a0+b1;a1=Logistic(n1);n2=w2*a1+b2;a2=Logistic(n2);a=a2;disp(1 0=,num2str(a);a0=double(1 1);n1=w1*a0+b1;a1=Logistic(n1);n2=w2*a1+b2;a2=Logistic(n2);a=a2;disp(1 1=,num2str(a);m=0;%-function a=Logistic(n)a=1./(1+exp(-n);%-function result=F(a)r,c=size(a);result=ze
11、ros(r,r);for i=1:rresult(i,i)=(1-a(i)*a(i);end;4、实验结果：计算误差=0.0049993迭代次数=27060 0=0.0231820 1=0.9631101 0=0.9653901 1=0.0433745、学习曲线图：图 1.MLP二、试用用奇阶互补法设计两带滤波器组二、试用用奇阶互补法设计两带滤波器组(高、低通互补高、低通互补)，进而实现四带，进而实现四带滤波器组；并画出其频响。滤波器设计参数为：滤波器组；并画出其频响。滤波器设计参数为：F Fp p1.71.7KHzKHz,F Fr r2.32.3KHzKHz,F Fs s8 8KHzKHz,
12、A Arminrmin70dB70dB。1、设计步骤：(1)对Fp、Fr进行预畸Fp);FsFrr tg();Fsp tg(2)计算cp*r，判断c是否等于1，即该互补滤波器是否为互补镜像滤波器(3)计算相关系数1 12k1(10k pr;0.1Apmin1)(100.1Ar min1);k1 k2;1 1kq0;21k5913q q0 2q015q0150q0;N lg(k12/16)/lgq;i;(N为 奇数)u 1i;(N为 偶数)22qi12m0(1)m1mqm(m1)sin(21 2(1)mqm(2m 1)u)N;2mcos(u)NN N N2;N1 N2;42vi(1 ki2)(1
13、i2/k);i2v2i;i 1,N21 22i2v2i1;i 1,N121 2i1i(4)互补镜像滤波器的数字实现Ai2i2i;Bi;2i2iAi Z2Bi Z21H1(Z);i 1,N1H2(Z)Z;i 1,N2221 A Z1 B Ziiii1HL(Z)H1(Z)H2(Z);22、程序：function filter2()Fp=1700;Fr=2300;Fs=8000;Wp=tan(pi*Fp/Fs);Wr=tan(pi*Fr/Fs);Wc=sqrt(Wp*Wr);k=Wp/Wr;k1=sqrt(sqrt(1-k2);q0=0.5*(1-k1)/(1+k1);q=q0+2*q05+15*q
14、09+150*q013;N=11;N2=fix(N/4);M=fix(N/2);N1=M-N2;for jj=1:Ma=0;for m=0:5a=a+(-1)m*q(m*(m+1)*sin(2*m+1)*pi*jj/N);%N is odd,u=jendab=0;for m=1:5b=b+(-1)m*q(m2)*cos(2*m*pi*jj/N);endbW(jj)=2*q0.25*a/(1+2*b);V(jj)=sqrt(1-k*W(jj)2)*(1-W(jj)2/k);endfor i=1:N1alpha(i)=2*V(2*i-1)/(1+W(2*i-1)2);endfor i=1:N2be
15、ta(i)=2*V(2*i)/(1+W(2*i)2);endfor i=1:N1a(i)=(1-alpha(i)*Wc+Wc2)/(1+alpha(i)*Wc+Wc2);endfor i=1:N2b(i)=(1-beta(i)*Wc+Wc2)/(1+beta(i)*Wc+Wc2);endw=0:0.0001:0.5;LP=zeros(size(w);HP=zeros(size(w);for n=1:length(w)z=exp(j*w(n)*2*pi);H1=1;for i=1:N1H1=H1*(a(i)+z(-2)/(1+a(i)*z(-2);endH2=1/z;for i=1:N2H2=H
16、2*(b(i)+z(-2)/(1+b(i)*z(-2);endLP(n)=abs(H1+H2)/2);HP(n)=abs(H1-H2)/2);endplot(w,LP,k,w,HP,m);%hold on;xlabel(数字频率);ylabel(幅度);3、实验结果：图 2.两带滤波器4、四带滤波器组程序：function filterfourFp=1700;Fr=2300;Fs=8000;Wp=tan(pi*Fp/Fs);Wr=tan(pi*Fr/Fs);Wc=sqrt(Wp*Wr);k=Wp/Wr;k1=sqrt(sqrt(1-k2);q0=0.5*(1-k1)/(1+k1);q=q0+2
17、*q05+15*q09+150*q013;N=11;N2=fix(N/4);M=fix(N/2);N1=M-N2;for jj=1:Ma=0;for m=0:5a=a+(-1)m*q(m*(m+1)*sin(2*m+1)*pi*jj/N);%N is odd,u=jendb=0;for m=1:5b=b+(-1)m*q(m2)*cos(2*m*pi*jj/N);endW(jj)=2*q0.25*a/(1+2*b);V(jj)=sqrt(1-k*W(jj)2)*(1-W(jj)2/k);endfor i=1:N1alpha(i)=2*V(2*i-1)/(1+W(2*i-1)2);endfor i
18、=1:N2beta(i)=2*V(2*i)/(1+W(2*i)2);endfor i=1:N1a(i)=(1-alpha(i)*Wc+Wc2)/(1+alpha(i)*Wc+Wc2);endfor i=1:N2b(i)=(1-beta(i)*Wc+Wc2)/(1+beta(i)*Wc+Wc2);endw=0:0.0001:0.5;LLP=zeros(size(w);LHP=zeros(size(w);HLP=zeros(size(w);HHP=zeros(size(w);for n=1:length(w)z=exp(j*w(n)*2*pi);H1=1;for i=1:N1H1=H1*(a(i)
19、+z(-2)/(1+a(i)*z(-2);endH21=1;for i=1:N1H21=H21*(a(i)+z(-4)/(1+a(i)*z(-4);endH2=1/z;for i=1:N2H2=H2*(b(i)+z(-2)/(1+b(i)*z(-2);endH22=1/(z2);for i=1:N2H22=H22*(b(i)+z(-4)/(1+b(i)*z(-4);endLP=(H1+H2)/2);HP=(H1-H2)/2);LLP(n)=abs(H21+H22)/2*LP);LHP(n)=abs(H21-H22)/2*LP);HHP(n)=abs(H21+H22)/2*HP);HLP(n)=
20、abs(H21-H22)/2*HP);endplot(w,LLP,k,w,LHP,m,w,HLP,g,w,HHP,b);xlabel(数字频率);ylabel(幅度);5、实验结果：图 3.四带滤波器组三、根据现代数字信号处理根据现代数字信号处理（姚天任等，华中理工大学出版社，（姚天任等，华中理工大学出版社，20012001）第）第四章附录提供的数据四章附录提供的数据(pp.352-353)(pp.352-353)，试用如下方法估计其功率谱，并画出不同，试用如下方法估计其功率谱，并画出不同参数情况下的功率谱曲线：参数情况下的功率谱曲线：1 1）LevinsonLevinson 算法算法2 2）
21、BurgBurg 算法算法3 3）ARMAARMA 模型法模型法4 4）MUSICMUSIC 算法算法1、Levinson 算法分析：(1)、由输入数据估计自相关函数，一种渐近无偏估计（称之为取样自相关函数）的公式为：(m)1RxxNN1 mn0 x(n)x(m n),*km N 1(2)、根据估计所得到的自相关函数，用下面的迭代公式估算 AR 模型参数：R(0)ak,iR*(i)2ki1ak,0 0Dkak,iR(k 1i),i0kk1Dk2k2k21(1k1)k2*ak1,i ak,ik1ak,k1i,i 1,2,kak1,k1 k1(3)、对于 AR（p）模型，按以上述递推公式迭代计算直
22、到k 1 p时为止。将算出的模型参数代入下式即可得到功率谱估计值：(e)Sxxj2p1ap,ie jwii1p2程序：function sigma2,a=levinson(signal_source,p)%阶数由 p 确定N=length(signal_source);%确定自相关函数for m=0:N-1R(m+1)=sum(conj(signal_source(1:N-m).*signal_source(m+1:N)/N;end%设置迭代初值a1=-R(2)/R(1);sigma2=(1-abs(a1)2)*R(1);gamma=-a1;%开始迭代for k=2:psigma2(k)=R(
23、1)+sum(a1.*conj(R(2:k);D=R(k+1)+sum(a1.*R(k:-1:2);gamma(k)=D/sigma2(k);sigma2(k)=(1-abs(gamma(k)2)*sigma2(k);a1=a1-gamma(k)*conj(fliplr(a1),-gamma(k);enda=1 a1;%计算功率谱估计值sigma2=real(sigma2);p=15;%p 分别为 15、30、45、60sigma2,a=Levinson(signal_in_complex,p);%计算功率谱f1=linspace(-0.5,0.5,512);%从-0.5 到 0.5 生成 5
24、12 个等间隔数据for k=1:512S1(k)=10*log10(sigma2(end)/(abs(1+sum(a(2:end).*exp(-j*2*pi*f1(k)*1:p)2);%公式(2.3.7)并以 dB 表示end；实验结果：图 4.Levinson 算法2、Burg 算法分析：(1)、设输入数据序列为x(n)，0 n N 1，对前后向预测误差之和求偏导，得反射系数*2ek1(n)ek1(n1)nkN1k(enkN1k1(n)ek1(n1)22前后向预测误差递推公式：ek(n)1 e(n)*kkk ek1(n)1ek1(n1)ak,i ak 1,ikak 1,k i,i 1,2,
25、3,.,k,ak 1,0 1(2)、重复以上步骤直至k=p，根据迭代得到的AR 模型参数计算功率谱，计算功率谱的公式同上面算法。程序：function sigma2,a=burg(signal_source,p)N=length(signal_source);ef=signal_source;eb=signal_source;sigma2=sum(signal_source*signal_source)/N;a=;for k=1:pefp=ef(k+1:end);ebp=eb(k:end-1);gamma(k)=2*efp*ebp/(efp*efp+ebp*ebp);sigma2(k+1)=(
26、1-abs(gamma(k)2)*sigma2(k);ef(k+1:end)=efp-gamma(k)*ebp;eb(k+1:end)=ebp-gamma(k)*efp;a=a-gamma(k)*conj(fliplr(a),-gamma(k);end;a=1 a;sigma2=real(sigma2);实验结果：图 5.Burg 算法3、ARMA 算法分析：(1)、用 x(n)通过 A(z)，得到 y(n)。(2)、用一无穷阶的 AR 模型近似 MA 模型。用 Burg 算法可得到此近似 AR 模型的参数以及激励白噪声的功率。一般此 AR 模型的阶数应大于 MA 模型阶数的两倍，以获得较好的
27、近似效果。(3)、可以证明，将上一步求出的近似 AR 模型参数视为时间序列，则 MA 模型就可视为一线性预测滤波器，按最小均方误差准则就可以求出 MA 模型参数，方法同样可用 Burg 算法。这样，ARMA 模型的参数就全部估计出来了，根据以下公式即可算出功率谱：(ej)2Sxxp1bie jwii1p2q21aie jwii1程序：function a,b,sigma2=arma(signal_source,p,q)N=length(signal_source);M=N-1;r=xcorr(signal_source,unbiased);for k=1:pR(:,k)=(r(N+q-k+1:
28、N+M-k).;enda1=(-pinv(R)*(r(N+q+1:N+M).).;a=1 a1;Y=filter(1 a1,1,signal_source);pp=4*q;sigma,a1=burg(Y,pp);sigma2=sigma(2:end);sigma,a2=burg(a1(2:end),q);b=a2;实验结果：图 6.ARMA 算法4、MUSIC 算法分析：(1)构造自相关阵R(1)R(N 1)R(0)R(0)R(N 2)R(1)R R R(N 1)R(N 2)R(0)自相关函数可用有偏估计代替。i i，i=1,2,N。(2)计算 R 的特征向量v v(3)估计 R 的维数 M，
29、方法有 AIC 和 MDL 法。(4)根据剩余特征向量计算 MUSIC 谱。PMUSIC(f)1iM 1e eNH Hi iv v2程序：function S=music(X,p,ii)N=length(X);r=xcorr(X,biased);clear Rfor k=1:NR(:,k)=(r(N-(k-1):2*N-k).;endV,D=eig(R);f=linspace(-0.5,0.5,512);S0=fft(V(:,p+1:end),512);for k=1:512S(k)=10*log10(1/(S0(k,:)*S0(k,:);endS=S(257:512)S(1:256);subplot(2,2,ii);plot(f,S,b);title(MUSIC:,int2str(p),维);xlabel(归一化频率),ylabel(功率谱/dB);hold on;实验结果：图 7.MUSIC 算法