通信仿真实验.doc

通信仿真实验学号：10071213 姓名：胡帆目录 1 实验二 滤波器的设计和仿真实现2 2 实验四 随机过程的计算机仿真12 3 实验七 数字基带调制16 4信号的功率谱与带限信道的滤波器等效20 5 实验：基带传输最佳接收机仿真26 6 4psk传输的误码率仿真34 7 4PAM系统的Monte Carlo仿真37 8 信号的功率谱与带限信道的滤波器等效41实验二 滤波器的设计和仿真实现实验日期：2013.10.27实验目的：各种滤波器的设计与仿真实现。实验内容：1.题目一：试设计一个模拟低通滤波器，fp=3500Hz，fs=4500Hz，p=3 dB，s=25dB。分别用巴特沃斯和椭圆滤波器原型，求出其3dB截止频率和滤波器阶数，传递函数，并作出幅频、相频特性曲线。 实验程序1巴特沃斯型：clear;f_p=3500; f_s=4500; R_p=3; R_s=25;n, fn=buttord(f_p,f_s,R_p,R_s, s);Wn=2*pi*fn;b,a=butter(n, Wn, s); f=0:100:10000; s=j*2*pi*f; H_s=polyval(b,s)./polyval(a,s); figure(1);subplot(2,1,1); plot(f, 20*log10(abs(H_s);axis(0 10000 -40 1);xlabel(频率 Hz);ylabel(幅度 dB);subplot(2,1,2); plot(f, angle(H_s); xlabel(频率 Hz);ylabel(相角 rad);figure(2); freqs(b,a);实验结果1： 实验程序2椭圆滤波器型：clear;f_p=3500; f_s=4500; R_p=3; R_s=25; n, fn=ellipord(f_p,f_s,R_p,R_s,s); Wn=2*pi*fn; b,a=ellip(n,R_p,R_s,Wn,s); f=0:100:10000; s=j*2*pi*f; H_s=polyval(b,s)./polyval(a,s); figure(1);subplot(2,1,1); plot(f, 20*log10(abs(H_s); axis(0 10000 -40 1);xlabel(频率 Hz);ylabel(幅度 dB);subplot(2,1,2); plot(f, angle(H_s); xlabel(频率 Hz);ylabel(相角 rad);figure(2); freqs(b,a);实验结果2：2.题目二：试设计一个巴特沃斯型数字低通滤波器，设采样率为8000Hz，fp=2100Hz，fs=2500Hz，p=3dB，s=25dB。并作出幅频、相频特性曲线。实验程序：f_N=8000; f_p=2100; f_s=2500; R_p=3; R_s=25; Ws=f_s/(f_N/2); Wp=f_p/(f_N/2); n, Wn=buttord(Wp,Ws,R_p,R_s); b,a=butter(n, Wn); figure(1);freqz(b,a, 1000, 8000) subplot(2,1,1); axis(0 4000 -30 3)figure(2); f=0:40:4000; z=exp(j*2*pi*f./(f_N); H_z=polyval(b,z)./polyval(a,z); subplot(2,1,1); plot(f, 20*log10(abs(H_z); axis(0 4000 -40 1);xlabel(频率 Hz);ylabel(幅度 dB);subplot(2,1,2); plot(f, angle(H_z); xlabel(频率 Hz);ylabel(相角 rad);实验结果:3.题目三：试设计一个切比雪夫1型高通数字滤波器，采样率为8000Hz，fp=1000Hz，fs=700Hz，p=3dB，s=20dB。并作出幅频、相频特性曲线。实验程序：f_N=8000; f_p=1000; f_s=700; R_p=3; R_s=20; Ws=f_s/(f_N/2); Wp=f_p/(f_N/2); n, Wn=cheb1ord(Wp,Ws,R_p,R_s); b,a=cheby1(n, R_p, Wn, high); freqz(b,a, 1000, 8000)subplot(2,1,1); axis(0 4000 -30 3)实验结果：4.题目四：试设计一个椭圆型带通数字滤波器。设采样率为10000Hz，fp=1000,1500 Hz，fs=600,1900 Hz，p=3dB，s=20dB。并作出幅频、相频特性曲线。实验程序：f_N=10000; f_p=1000, 1500; f_s=600, 1900; R_p=3; R_s=20; Ws=f_s/(f_N/2); Wp=f_p/(f_N/2); n, Wn=ellipord(Wp,Ws,R_p,R_s); b,a=ellip(n, R_p, R_s, Wn); freqz(b,a, 1000, 10000) subplot(2,1,1); axis(0 5000 -30 3)实验结果：5.题目五：试设计一个切比雪夫2型带阻数字滤波器。设采样率为10000Hz，fp=1000,1500 Hz，fs=1200,1300 Hz，p=3dB，s=20dB。并作出幅频、相频特性曲线。实验程序：f_N=10000; f_p=1000, 1500; f_s=1200, 1300; R_p=3; R_s=20; Ws=f_s/(f_N/2); Wp=f_p/(f_N/2); n, Wn=cheb2ord(Wp,Ws,R_p,R_s); b,a=cheby2(n,R_s, Wn, stop); freqz(b,a, 1000, 10000) subplot(2,1,1); axis(0 5000 -35 3)实验结果：6.题目六：在采样率为8000Hz下设计一个在500Hz，1000Hz，1500Hz，2000Hz，.，n*500Hz的地方开槽陷波。陷波带宽(-3dB 处)为60Hz。试设计该滤波器。实验程序：Fs=8000; Ts=1/8000; f0=500; bw=60/(Fs/2); ab=-3 ; n=Fs/f0; num,den = iircomb(n,bw,ab, notch); freqz(num,den, 4000, 8000);axis(0 4000 -30 5);实验结果：7.题目七：用Matlab设计具有下列指标的线性相位FIR低通滤波器：阻带截止频率为0.45和0.8，通带截止频率为0.55和0.7，最大通带衰减为0.15dB，最小阻带衰减为40dB。分别用下面的窗函数来设计滤波器：海明窗、汉宁窗、布莱克曼窗和凯泽窗。对于每种情况，显示其冲激响应系数并画出设计的滤波器增益响应。分析设计结果。实验程序：close allOmegap=0.7*pi;Omegas=0.8*pi;Omegat=18;wp=2*pi*Omegap/Omegat;ws=2*pi*Omegas/Omegat;alphap=0.15;alphas=40;NUM_Hamming=1;NUM_Hann=2;NUM_Blackman=3;NUM_Kaiser=4;c=3.32,3.11,5.56*pi;wc=(ws+wp)/2;delt_w=ws-wp;freq_labels=用海明窗设计的FIR频率响应,用汉宁窗设计的FIR频率响应,用布莱克曼窗设计的FIR频率响应;ht_labels=海明窗FIR冲激响应,汉宁窗FIR冲激响应,布莱克曼窗FIR响应;for filter_kind=NUM_Hamming:NUM_Kaiser,M=ceil(c(filter_kind)/delt_w); N=2*M+1; switch filter_kind, case NUM_Hamming, win=hamming(N); display(海明窗生成的冲激响应系数：, (阶数N=,num2str(N),); figure case NUM_Hann, win=hann(N); display(汉宁窗生成的冲激响应系数：, (阶数N=,num2str(N),); figure case NUM_Blackman, win=blackman(N); display(布莱克曼窗生成的冲激响应系数：,(阶数N=,num2str(N),); figure otherwise disp(error); end n=-M:M; hd=sin(wc*n)./(pi*n); hd(find(n=0)=wc*cos(wc*0)/pi; ht=hd.*win; display(,num2str(ht); subplot(1,2,1); plot(n,ht,.-) title(ht_labels(filter_kind); xlabel(n,FontSize,12); ylabel(ht,fontsize,12); grid on h,w=freqz(ht,1,512); W=w/pi; H=20*log10(abs(h); subplot(1,2,2); hold on title(freq_labels(filter_kind); plot(W,H); xlabel(pi (omega/omega,s); ylabel(增益 (dB); grid onend 凯泽窗程序：Wp=0.7*pi;Ws=0.8*pi;As=40;dw=Ws-Wp;N=ceil(10*pi/dw)+1;beta=0.5842*(As-21)0.4+0.07886*(As-21);b=fir1(N,Wp/pi,low,kaiser(N+1,beta);figure(3)freqz(b,1);title（凯泽窗响应）;实验结果：实验四 随机过程的计算机仿真-实验日期：2013.10.20实验目的：仿真实现各种随机分布的随机数发生器实验内容：1、 题目一：均匀分布随机数的产生用线性同余法，编写Matlab程序，产生均匀分布的随机数。初始种子x(0)自己选择。实验程序1：a=input(Enter multiplier a);c=input(Enter offset c);m=input(Enter modulus m);seed=input(Enter seed );n=1;ix=rem(seed*a+c),m);while (ix=seed) & (nm) disp(Caught in a loop);else text=The period is,num2str(n,15),.; disp(text);end实验结果1：Enter multiplier a241Enter offset c1323Enter modulus m5000Enter seed 1The period is5000.2、 题目二：用反函数法，将均匀分布的随机变量变换为具有单边指数分布的随机变量。编写Matlab程序，产生指数分布的随机数。计算并比较理论pdf和从直方图得到的pdf。指数分布随机变量pdf定义为：，为单位阶跃函数。实验程序2：n=input(输入样点个数 ); a=5; u=rand(1,n); x_exp=-log(2*u)/a; N_sample,x=hist(x_exp,20); subplot(2,1,1)bar(x,N_sample,1) ylabel(样本个数)xlabel(变量x)px=a/2*exp(-a*x); P_hist= N_sample/n; del_x=x(2)-x(1); p_hist=N_sample/n/del_x; subplot(2,1,2)plot(x,px,k,x,p_hist,ok)ylabel(概率密度 )xlabel(变量x)legend(期望pdf,仿真得到pdf,1)实验结果2：N=2000，指数分布3、 题目三：用Matlab编程分别产生标准正态分布、指定均值方差正态分布、瑞利分布、赖斯分布、中心与非中心2分布的随机数，并画出相应的pdf。实验程序3：x=-10:0.01:10y1=normpdf(x,0,1);figureplot(x,y1)xlabel(标准正态分布)ylabel(pdf)y2=normpdf(x,2,1);figureplot(x,y2)xlabel(正态分布)ylabel(pdf)sigma=3s=sigma2y3= x/s.*exp(-x.2/s)%subplot(1,3,3)figureplot(x,y3)xlabel(瑞利分布)ylabel(pdf)N = 100000;K = 0.5;const=1/(2*(K+1);a1=randn(1,N);b2=randn(1,N);a=sqrt(const*(a1+sqrt(2*K).2+a2.2);y4,ai = ksdensity(a);figureplot(ai,y4, bo)xlabel(莱斯分布 )ylabel(pdf)p1 = ncfpdf(x,5,20,10);p = fpdf(x,5,20);figureplot(x,p,-,x,p1,-)xlabel(非中心f分布)ylabel(pdf)实验结果3： 4、 题目四：设输入的随机变量序列X(n)为N=1000独立同分布高斯分布的离散时间序列，均值为0，方差为1，采样间隔0.01s。通过某线性时不变滤波器，输出随机变量序列Y(n)的功率谱密度为：（1） 设计该滤波器（2） 产生随机变量序列Y(n)。实验程序4：fs=100; f=0:1:50; F=2*f/fs; M=10./sqrt(1+(2*pi.*f).2); b,a=yulewalk(20,F,M); plot(F,M) xlabel()ylabel()N=1000; X=randn(1,N); Y=filter(b,a,X) 实验结果4：实验七 数字基带调制实验日期：2013.10.20实验目的：数字通信系统中，基带传输的仿真。实验内容：1.题目一：用MATLAB编程仿真实现四进制脉冲幅度调制(PAM)数字通信系统，并用蒙特卡罗仿真方法计算在信道为加性高斯白噪声时，该系统在不同信噪比下的差错概率。实验程序1：clc;clear all;SNRindB1=0:1:16; SNRindB2=0:0.2:16; for i=1:length(SNRindB1), simu_err_prb(i)=simuPe(SNRindB1(i);end;for i=1:length(SNRindB2), SNR=10(SNRindB2(i)/10); theo_err_prb(i)=(3/2)*Qfunc(sqrt(2/5)*SNR);end;semilogy(SNRindB2,theo_err_prb,SNRindB1,simu_err_prb,o)axis(0 16 0.0001 1)xlabel(SNR in dB)ylabel(Prb of Err)legend(Theoretical,Simulation)function y=gngauss(sgma)y=sgma*randnfunction p=simuPe(snr_in_dB)d=1;A0=-3*d;A1=-d;A2=d;A3=3*d; SNR=10(snr_in_dB/10); sgma=sqrt(5/2)*d2/SNR); N=10000; for i=1:N temp=rand; if (temp1/4) symbol(i)=0; elseif (temp2/4) symbol(i)=1; elseif (temp3/4) symbol(i)=3; else symbol(i)=2; end;end;numoferr=0; for i=1:Nif (symbol(i)=0) r=A0+gngauss(sgma); elseif (symbol(i)=1) r=A1+gngauss(sgma); elseif (symbol(i)=3) r=A2+gngauss(sgma); else r=A3+gngauss(sgma); end; if (abs(r-A0)d | rA0) deci_symbol=0; elseif (abs(r-A1)d) deci_symbol=1; elseif (abs(r-A2)d) deci_symbol=3; else deci_symbol=2; end; if (deci_symbol=symbol(i) numoferr=numoferr+1; end;end;p=numoferr/N;实验结果1：2.题目二：设计FIR根升余弦滤波器，具体指标如下：（1）码片速率为1.28MHz，采样率为4倍码片速率（2）滚降系数0.22，冲激响应序列长度为65实验程序2：N_T=8; R=0.22 Fc=1.28e+6;Fs=4*Fc; Tc=1.0e-6/1.28; Num=rcosfir(R,N_T,4,Tc,sqrt);H,w=freqz(Num,1,1000,whole);H=(H(1:1:501);w=(w(1:1:501);Mag=abs(H);db=20*log10(Mag)/max(Mag);pha=angle(H);plot(w/pi,db);grid;axis(0 1 -60 1);xlabel(归一化角频率 );ylabel(RRC滤波器幅度响应（db）);实验结果2：3.题目三：产生一串（-1.1）均匀分布的随机序列，并对该序列进行脉冲成形滤波。实验程序3：clearclcN_Filter=8; R=0.22Fc=1.28e+6;Rate=4;Fs=Rate*Fc;Tc=1.0e-6/1.28; Delay=8;NSig=16;Sig=1,zeros(1,15);Y = RCOSFLT(Sig, Fc, Fs, fir/sqrt,R, Delay);Xaxis=1:0.25:2*NSig+0.75;figure(1)subplot(211),stem(Sig),grid,Xlabel(k=t/Tc=t/Ts/4),axis(0 16 -0.2 1.2);subplot(212),stem(Xaxis,Y),Xlabel(k=t/Tc=t/Ts/4),axis(0 16 -0.2 0.8),grid实验结果3：信号的功率谱与带限信道的滤波器等效 -实验日期：2013.11.24实验目的：数字调制信号功率谱;使用Remez方法设计滤波器仿真带限信道;使用窗函数法设计滤波器仿真带限信道升余弦滤波器设计;滤波器的频域分析方法实验内容：1. 使用Remez方法设计：实验程序：echo onf_cutoff=2000; f_stopband=3000; fs=50000; f1=2*f_cutoff/fs; f2=2*f_stopband/fs; N=101; F=0 f1 f2 1; M=1 1 0 0; B=remez(N-1,F,M); figure(1);H,W=freqz(B);H_in_dB=20*log10(abs(H);plot(W/(2*pi),H_in_dB);figure(2);plot(W/(2*pi),(180/pi)*unwrap(angle(H);figure(3);plot(zeros(size(0:N-1);hold;stem(0:N-1,B);实验结果：2.使用窗函数法，窗函数为矩形窗，汉宁窗实验程序：echo onLength=101;Fs=50000;W=1000;Ts=1/Fs;n=-(Length-1)/2:(Length-1)/2;t=Ts*n;h=2*W*sinc(2*W*t); N=61;rec_windowed_h=h(Length-N)/2+1:(Length+N)/2);rec_windowed_H,W1=freqz(rec_windowed_h,1);rec_windowed_H_in_dB=20*log10(abs(rec_windowed_H)/abs(rec_windowed_H(1);hanning_window=hanning(N);hanning_windowed_h=h(Length-N)/2+1:(Length+N)/2).*hanning_window.;hanning_windowed_H,W2=freqz(hanning_windowed_h,1);hanning_windowed_H_in_dB=20*log10(abs(hanning_windowed_H)/abs(hanning_windowed_H(1);figure(1);plot(W1/pi, rec_windowed_H_in_dB);hold onplot(W2/pi, hanning_windowed_H_in_dB,g-);legend(rec windowed,hanning windowed)xlabel(f)ylabel(H(f);figure(2);stem(0:N-1,rec_windowed_h,*);hold onstem(0:N-1,hanning_windowed_h,o);legend(rec windowed,hanning windowed)xlabel(n)ylabel(h(n);实验结果：3. 升余弦滤波器设计 实验程序：echo onN=31;T=1;alpha=1/4;n=-(N-1)/2:(N-1)/2; for i=1:length(n), g_T(i)=0; for m=-(N-1)/2:(N-1)/2, g_T(i)=g_T(i)+sqrt(xrc(4*m/(N*T),alpha,T)*exp(j*2*pi*m*n(i)/N); echo off ; end;end;echo on ;g_T=real(g_T) ; n2=0:N-1; G_T,W=freqz(g_T,1);magG_T_in_dB=20*log10(abs(G_T)/max(abs(G_T); g_R=g_T;imp_resp_of_cascade=conv(g_R,g_T); X_T,W = freqz(imp_resp_of_cascade,1);plot(2*W/(2*pi),20*log10(abs(X_T)/max(abs(X_T)xlabel(f)ylabel(GT(f)figure;stem(0:N-1,g_T);xlabel(n)ylabel(GT(n)figure;stem(0:2*(N-1),imp_resp_of_cascade)xlabel(n)ylabel(GT(n)*GR(n)实验结果：实验报告要求：1. 本章所设计的是fir滤波器2. 通带阻带较平坦，过渡带较窄，阻带衰减好3. 汉宁窗下降的比矩形窗多，汉宁窗好4. 平方根升余弦滤波器，是因为在发送端和接收端分别用一个平方根升余弦波器，既能实现升余弦滤波器的作用，也能实现匹配滤波器的实现，从而既能满足奈奎斯特采样定理，又可以提升接收端信噪比，更便于准确接收信号。实验：基带传输最佳接收机仿真实验日期：2013.10.27实验要求：1. SNR=1：13dB2. 三种基带信号仿真误码率在同一图中主程序：echo onSNRindB1=0:1:13;for i1=1:length(SNRindB1),smld_err_prb(i1)=smldPe54(SNRindB1(i1); end;semilogy(SNRindB1, smld_err_prb,*);hold onfor i2=1:length(SNRindB1), smld_err_prb(i2)=smldPe55(SNRindB1(i2);end;semilogy(SNRindB1, smld_err_prb,+);hold onfor i3=1:length(SNRindB1), smld_err_prb(i3)=smldPe56(SNRindB1(i3);end;semilogy(SNRindB1,smld_err_prb,o);xlabel(10log_10(E/N_0)ylabel(Pe);legend(正交信号,反极性信号,开关信号)子程序：子程序smldPe54:function p=smldPe54(snr_in_dB)E=1;SNR=exp(snr_in_dB*log(10)/10); sgma=E/sqrt(2*SNR); N=10000;for i=1:N, temp=rand; if (tempr1), decis=0; else decis=1; end; if (decis=dsource(i), numoferr=numoferr+1; end;end;p=numoferr/N; 子程序smldPe55:function p=smldPe55(snr_in_dB)E=1;SNR=exp(snr_in_dB*log(10)/10); sgma=E/sqrt(2*SNR); N=10000;for i=1:N, temp=rand; if (temp0.5), dsource(i)=0; else dsource(i)=1; endend;numoferr=0;for i=1:N, if (dsource(i)=0), r=-E+gngauss(sgma); else r=E+gngauss(sgma); end; if (r0), decis=0; else decis=1; end; if (decis=dsource(i), numoferr=numoferr+1; end;end;p=numoferr/N; 子程序smldPe56:function p=smldPe56(snr_in_dB)E=1;alpha_opt=1/2;SNR=exp(snr_in_dB*log(10)/10); sgma=E/sqrt(2*SNR); N=10000;for i=1:N, temp=rand; if (temp0.5), dsource(i)=0; else dsource(i)=1; endend;numoferr=0;for i=1:N, if (dsource(i)=0), r=gngauss(sgma); else r=E+gngauss(sgma); end; if (ralpha_opt), decis=0; else decis=1; end; if (decis=dsource(i), numoferr=numoferr+1; end;end;p=numoferr/N;实验结果：散点图：3. =0,0.1,0.2,0.3,0.5,1时接收信号散点图如下：程序1：echo onsigma=0;n0=sigma*randn(100,1);n1=sigma*randn(100,1);n2=sigma*randn(100,1);n3=sigma*randn(100,1);x1=1.+n0;y1=n1;x2=n2;y2=1.+n3;plot(x1,y1,o,x2,y2,*)axis(square)xlabel(sigma=0)实验结果：程序2：echo onsigma=0.1;n0=sigma*randn(100,1);n1=sigma*randn(100,1);n2=sigma*randn(100,1);n3=sigma*randn(100,1);x1=1.+n0;y1=n1;x2=n2;y2=1.+n3;plot(x1,y1,o,x2,y2,*)axis(square)xlabel(sigma=0.1)实验结果：程序3：echo onsigma=0.2;n0=sigma*randn(100,1);n1=sigma*randn(100,1);n2=sigma*randn(100,1);n3=sigma*randn(100,1);x1=1.+n0;y1=n1;x2=n2;y2=1.+n3;plot(x1,y1,o,x2,y2,*)axis(square)xlabel(sigma=0.2)实验结果：程序4：echo onsigma=0.3;n0=sigma*randn(100,1);n1=sigma*randn(100,1);n2=sigma*randn(100,1);n3=sigma*randn(100,1);x1=1.+n0;y1=n1;x2=n2;y2=1.+n3;plot(x1,y1,o,x2,y2,*)axis(square)xlabel(sigma=0.3)实验结果：程序5：echo onsigma=0.5;n0=sigma*randn(100,1);n1=sigma*randn(100,1);n2=sigma*randn(100,1);n3=sigma*randn(100,1);x1=1.+n0;y1=n1;x2=n2;y2=1.+n3;plot(x1,y1,o,x2,y2,*)axis(square)xlabel(sigma=0.5)实验结果：程序6：echo onsigma=1;n0=sigma*randn(100,1);n1=sigma*randn(100,1);n2=sigma*randn(100,1);n3=sigma*randn(100,1);x1=1.+n0;y1=n1;x2=n2;y2=1.+n3;plot(x1,y1,o,x2,y2,*)axis(square)xlabel(sigma=1)实验结果：4psk传输的误码率仿真 -实验日期：2013.11.17一实验目的：1.建立4PSK传输的误码率仿真模型2.多进制符号(及二进制序列)到星座图的映射模块编程3.最大相关判决法的编程4.误码率与误符号率的区别二实验内容：1.画出4psk传输的误码率仿真模型2.画出星座图3.理论误码率：理论误码率程序：SNRindB2=0:0.1:12;for i=1:length(SNRindB2), SNR=exp(SNRindB2(i)*log(10)/10); theo_err_prb(i)=Qfunct(sqrt(2*SNR);end;semilogy(SNRindB2,theo_err_prb);legend(理论误码率 );理论误码率实验结果：4.4psk仿真实验实验完整程序：主程序：echo onSNRindB1=0:1:12;for i=1:length(SNRindB1), pb,ps=cm_sm32(SNRindB1(i); smld_bit_err_prb(i)=pb; smld_symbol_err_prb(i)=ps;end;semilogy