matlab通信系统仿真实验报告

matlab通信系统仿真实验报告 MATLAB通信系统仿真实验报告 专业 姓名 学号 指导教师 实验学时 实验时间 实验地点 实验一、MATLAB的基本使用与数学运算 目的:学习MATLAB的基本操作，实现简单的数学运算程序。 内容: 1-1 要求在闭区间[0，2π]上产生具有10个等间距采样点的一维数组。试用两种不同的指令实现。 运行代码:x=[0:2*pi/9:2*pi] 运行结果: 1-2 用M文件建立大矩阵x x=[ 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3.1 3.2 3.3 3.4 3.5 3.6 3.7 3.8 3.9] 代码:x=[ 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 3.1 3.2 3.3 3.4 3.5 3.6 3.7 3.8 3.9] m_mat 运行结果: 1-3已知A=[5，6;7，8]，B=[9，10;11,12]，试用MATLAB分别计算 A+B,A*B,A.*B,A^3，A.^3，A/B,A\B. 代码:A=[5 6;7 8] B=[9 10;11 12] x1=A+B X2=A-B X3=A*B X4=A.*B X5=A^3 X6=A.^3 X7=A/B X8=A\B 运行结果: 1-4任意建立矩阵A，然后找出在[10,20]区间的元素位置。 程序代码及运行结果: 代码:A=[12 52 22 14 17;11 10 24 03 0;55 23 15 86 5 ] c=A>=10&A<=20 运行结果: 1-5 总结:实验过程中，因为对软件太过生疏遇到了些许困难，不过最后通过查与同学交流都解决了。例如第二题中，将文件保存在了D盘，而导致频频出错，最后发现必须保存在MATLAB文件之下才可以。第四题中，逻辑语言运用到了ij，也出现问题，虽然自己纠正了问题，却也不明白错在哪了，在老师的讲解下知道位置定位上不能用ij而应该用具体的整数。总之第一节实验收获颇多。 实验二、MATLAB程序的编写 目的:掌握顺序结构、选择结构、循环结构程序设计方法。学会编写函数。 内容: 2-1编写程序，建立向量N=[1,2,3,4,5]，然后利用向量N产生下列向量; (1)2,4,6,8,10 (2)1/2，1，3/2，2，5/2 (3)1，1/2，1/3，1/4，1/5 (4)1，1/4，1/9，1/16，1/25 代码:N=[1,2,3,4,5] X1=N*2 X2=N/2 X3=1./N X4=X3*X3 运行结果: 2-2从键盘输入一个三位整数，将他反向输出，如输入为639，输出为936.输入一个百分制成绩，要求输出成绩等级A,B,C,D,E。其中90~100分为A，80~89分为B，70~79分为C，60~69分为D，60分以下为E。 要求: (1)分别用if语句 代码:clear m=input('请输入一个三位数:') m1=fix(m/100); m2=rem(fix(m/10),10); m3=rem(m,10); n=m1+m2*10+m3*100; disp(n); (2)clear; Mark=input('请输入成绩:'); Rank=cell(1,5); S=struct('Marks',Mark, 'Rank',Rank); for i=1:10; a{i}=89+i; b{i}=79+i; c{i}=69+i; d{i}=59+i; e{i}=0+i; q{i}=9+i; g{i}=19+i; h{i}=29+i; m{i}=39+i; n{i}=49+i; end; for i=1:5; switch S(i).Marks case 100 S(i).Rank='A'; case a S(i).Rank='A'; case b S(i).Rank='B'; case c S(i).Rank='C'; case d S(i).Rank='D'; case e S(i).Rank='E'; case q S(i).Rank='E'; case g S(i).Rank='E'; case h S(i).Rank='E'; case m S(i).Rank='E'; case n S(i).Rank='E'; otherwise S(i).Rank='成绩输入错误'; end end disp([num2str(S(i).Marks),blanks(3),S(i).Rank]);disp(''); 运行结果: —3输入20个两位随机数，求其中的最大数最小数。要求分别用循环结构和调2 用MATLAB的max函数、min函数实现。 (1)a=fix(rand(1,20)*100) ma=max(a) mi=min(a) 运行结果: (2)a=fix(rand(1,20)*100); for i=1:20; max=a(1); min=a(1); if maxa(i); min=a(i); end end max min 运行结果: 2-6写出下列程序输出结果 (1) s=0; a=[12,13,14;15,16,17;18,19,20;21,22,23]; for k=a for j=1:4 if rem(k(j),2)~=0 s=s+k(j); end end end s 运行结果: (2) global x x=1:2:5; y=2:2:6; sub(y); x y (3) function fun=sub(z) global x z=3*x; x=x+z; 运行结果: 总结:第二次实验，对软件的使用比较熟练了，但还是遇到了些许问题。在运算符号的使用中，应当注意“.*”的使用，在最初因为不太会运用遇到了些困难，后来通过同学讨论和翻阅课本找到了答案。2—2中的第二种方法是按照课本例题改编的，有些啰嗦，不多至少是结果正确。还有2—6中刚开始没能正常输出，在老师的指导下知道(2)(3)是一起使用，算是运用到了函数调用。好在最后所有题目都得到了满意的结果。 实验三、MATLAB图形处理 目的:能够根据数据绘制各种形状的二、三维图形。 3-1绘制曲线y=x^3+x+1,x的取值范围为[-5,5] 代码: x=-5:0.01:5 y=x.^3+x+1 plot(x,y) 运行结果: 3-4有一组测量数据满足y=exp(-a*t),t的变化范围为0~10，用不同的线性和标 记点画出a=0.1,a=0.2和a=0.5 三种情况下的曲线。 代码: t=0:0.1:10; y1=exp(-0.1*t); y2=exp(-0.2*t); y3=exp(-0.5*t); title('t from 0 to 10'); plot(t,y1,t,y2,t,y3); xlabel('Variable t'); ylabel('Variable y'); text(0.8,1.5,'曲线y1=exp^{-0.1t}'); text(2.5,1.1,'曲线y1=exp^{-0.2t}'); text(0.8,1.5,'曲线y1=exp^{-0.5t}'); legend('y1','y2','y3') 运行结果: 3-7绘制饼图，x=[66 49 71 56 38],并将第五个切块分离出来。 代码: x=[66 49 71 56 38]; subplot(1,2,1); pie(x); subplot(1,2,2); pie(x,[0,0,0,0,1]); 运行结果: 总结:这次实验，比较有成就感，并没有遇到什么太复杂的困难，但是软件操作上出现了写麻烦，一不小心将软件页面的各个功能窗口关上了，颇费周折终于找到了那些功能窗口，但是整个页面都有些混乱。好在还是将题目做了出来，图出现的时候感觉特别有成就感。真的说明一件事情，英语学不好很麻烦啊。 实验四、MATLAB仿真模拟调制 目的:能用MATLAB仿真调幅信号和调角信号。 5-1 用在区间[0，2]内的信号 m(t)=t 0<=t<=1;m(t)=-t+2 1<=t<=2; 以DSB-AM方式调制一个载波频率为25HZ、幅度为1的载波产生已调信号u(t)。写一个Matlab的M文件，并用该文件作下面的题: (1) 画出已调信号; (2) 求已调信号的功率; (3) 求已调信号的振频谱，并与消息信号m(t)的频谱作比较。 程序代码: dt=0.01; %时间采样间隔 fc=25; T=1; N=floor(T/dt); t1=[0:N]*dt; t2=t1+1; %t=[t1 t2]; mt1=t1; %信源 mt2=-t2+2; %DSB-AM modulation dsb1=mt1.*cos(2*pi*fc*t1); dsb2=mt2.*cos(2*pi*fc*t2); subplot(2,2,1); plot(t1,dsb1);hold on;plot(t2,dsb2); pwr1=mt1.^2; pwr2=mt2.^2; subplot(2,2,2); plot(t1,pwr1);hold on;plot(t2,pwr2); [mtf1,mtfft1]=FFT_SHIFT(t1,mt1); [mtf2,mtfft2]=FFT_SHIFT(t2,mt2); subplot(2,2,3); plot(mtf1,abs(mtfft1));hold on;plot(mtf2,abs(mtfft2)); 运行结果: 5-2 设AM调整时，输入信号为没(t)=0.2sin1000pi*t+0.5cos1000exp2 *pi*t,A=1,载波中心频率fc=10khz (1) 用MATLAB画出AM信号的波形及其频谱 程序代码: 1、function [f,sf]=FFT_SHIFT(t,st) df=t(2)-t(1); T=t(end); df=1/T; N=length(t); f=[-N/2:N/2-1]*df; sf=fft(st); sf=fftshift(sf); 2、dt=0.00001; %时间采样间隔 fm1=500;fm2=500*1.414; %信源频率 fc=10000; %载波中心频率 T=0.01; N=floor(T/dt); t=[0:N-1]*dt; mt=0.2*sin(2*pi*fm1*t)+0.5*cos(2*pi*fm2*t); %信源 %AM modulation A=1; am=(A+mt).*cos(2*pi*fc*t); [f,AMf]=FFT_SHIFT(t,am); subplot(311); plot(t,mt); subplot(312); plot(t,am); subplot(313); plot(f,AMf); 运行结果: 5-3 设FM调制时，调频器的输入信号为一个周期性的锯齿波，锯齿波的一个周期为信号g(t)=t0<=t<1,g(t)=0其他，FM的中心频率fc=100hz,Kfm=10hz,试做 (1) 画出调频后的信号波形及其振幅谱 (2) 若接收端采用鉴频器进行解调，且AWGN信道的功率密度谱为N0/2，试画出当解 调器输入信噪比0dB,10Db,20dB时的解调输出信号，并与原信号进行比较。 程序代码: 1、function [f,sf]=FFT_SHIFT(t,st) df=t(2)-t(1); T=t(end); df=1/T; N=length(t); f=[-N/2:N/2-1]*df; sf=fft(st); 2.dt=0.001; %时间采样间隔 fc=100; T=1; N=floor(T/dt); t=[0:N]*dt; kf=10; mt=t; mti=t.^2/2; fmt=cos(2*pi*fc*t+2*pi*kf*mti); figure(1); subplot(2,1,1); plot(t,fmt);hold on;plot(t,mt,'r'); [f,ft]=FFT_SHIFT(t,fmt); subplot(2,1,2); plot(f,abs(ft)); s=1/2;% 调制信号功率是A^2/2; sn=100;db=s/(10^(sn/10));%求白噪声的方差。 noise1=sqrt(db)*randn(size(t));%产生高斯白噪声。 fmt1=fmt+noise1; figure(2); subplot(2,1,1); plot(t,fmt1); N=length(fmt); dfmt=zeros(1,N); for k=1:N-1 %已调信号微分 dfmt(k)=(fmt1(k+1)-fmt1(k))/dt; end envlp=abs(hilbert(dfmt));%求瞬时幅度，即包络 subplot(2,1,2); plot(t,envlp); 运行结果: 总结:这一次实验遇到的问题是花费了两节课解决的。第一个题目很快的完成后，第二个跟第三个题目都遇到了问题。第二题不出图，第三题图不正确，反复检查修改了好多遍程序还是不对，一直提示[f,ft]=FFT_SHIFT(t,fmt)程序这句话有问题，就在打算放弃的时候，终于在老师的提醒下知道，是因为没有调用函数，讲课本后面附录中的函数体加上之后，图就顺利的出来了。说明对程序的理解还不够到位，需要继续学习和努力。 实验五、MATLAB仿真模拟信号的数字传输 目的:能用MATLAB仿真函数的抽样、量化过程，掌握信号编码方法。 6-1设低通信号s(t)=sin2*pi*t+0.5cos4*pi*t (1) 画出该低通信号的波形; (2) 画出抽样速率为fs=4Hz的抽样序列; (3) 从抽样序列恢复出原始信号; (4) 当抽样速率fs=2hz时，画出恢复出的抽样信号。 运行代码:1. clear all;close all; dt=0.01; t=0:dt:10; xt=sin(2*pi*t)+0.5*cos(4*pi*t); [f,xf]=FFT_SHIFT(t,xt); fs=4; sdt=1/fs; t1=0:sdt:10; st=sin(2*pi*t1)+0.5*cos(4*pi*t1); [f1,sf]=FFT_SHIFT(t1,st); t2=-50:dt:50; gt=sinc(fs*t2); stt=INSERT0(st,sdt/dt); xt_t=conv(stt,gt); figure(1) subplot(3,1,1); 原始信号 '); plot(t,xt);title(' subplot(3,1,2); stem(t1,st);title('抽样信号'); subplot(3,1,3); t3=-50:dt:60+sdt-dt; plot(t3,xt_t);title('抽样信号恢复'); axis([0 10 -1 1]) 2. function[out]=INSERT0(d,M) N=length(d); out=zeros(1,M*N); for i=0:N-1 out(i*M+1)=d(i+1); end; 3. function [f,sf]=FFT_SHIFT(t,st) df=t(2)-t(1); T=t(end); df=1/T; N=length(t); f=[-N/2:N/2-1]*df; sf=fft(st); 运行结果: 6-2用一个均匀量化器对零均值、单位方差的高斯源进行量化，这个量化器在区间[-10，10]内均匀量化。假定量化电平设在各量化区域的中间点，求出并画出量化电平数为N=3.4.5.6.7.8.9.10时，量化产生的均方失真作为量化电平数N的函数。 t=[0:0.01:1]; s=normrnd(0,1,size(t)); err=zeros(1,8); for i=1:8 [err(i),x_qtz]=myquantizer(s,i+2,-10.10); end N=[3,4,5,6,7,8,9,10]; Plot(n,err); function [err,x_qtz]=myquantizer(x,n,l,h) xmax=max(l,h); x_qtz=x/xmax;delta=2/n; q=delat*[0:n-1]-(n-1)/2*delat; for i=1:n index=find((q(i)-delat/2<=x_qtz)&(x_qtz<=q(i)+delta/2)); x_qtz(index)=q(i)*ones(1,length(index));end x_qtz=x_qtz*xmax; err=sum((x-x_qtz).^2)/length(x); function [sqnr,x_qtz,code]=UniPcm(x,n) xmax=max(abs(x)); x_qtz=x/xmax; b_qtz=x_qtz; delta=2/n; q=delat*[0:n-1]-(n-1)/2*delta; for i=1:n; index=delat=find((q(i)-delta/2<=x_qtz)&x_qtz<=q(i)+delta/2)); x_qtz(index)=(q(i)*ones(1,length(index)); b_qtz(find(x_qtz==q(i)))=(i-1)*ones(1,length(find(x_qtz==q(i)))); end x_qtz=x_qtz*xmax; nu=ceil(log2(n)); code=zeros(length(x),nu); for i=1:length(x) for j=nu:-1:0 if(fix(b_qtz(i)/(2^j))==1) code(i,nu-j)=1; b_qtz(i)=b_qtz(i)-2^j; end end end 6-6试编写A律13折线近似法的编码与解码程序。 function code=myAcode(x) iw=[0 16 32 64 128 256 512 1024] 实验六、MATLAB仿真数字信号的基带传输 目的:能够绘制常用码型，码型功率谱和眼图。 7-1 设二进制符号序列为1011010010，以矩形脉冲为例，利用Matlab分别画出 相应的单极性不归零、双极性不归零、单极性归零和双极性归零波形。 程序代码: 1、单极性不归零 function y=snrz(x); t0=200; t=0:1/t0:length(x); for i=1:length(x) if x(i)==1 for j=1:t0 y=((i-1)*t0+j)=1; end else for j=1:t0 y=((i-1)*t0+j)=0; end end end y=[y,x(i)]; plot(t,y); title(‘1 0 1 1 0 1 0 0 1 0’); axis([0,i,-0.1,1.1]); 2、双极性不归零:将程序1中的y=((i-1)*t0+j)=0改为y=((i-1)*t0+j)=-1 3、单极性归零 function y=srz(x); t0=200; t=0:1/t0:length(x); for i=1:length(x) if x(i)==1 for j=1:t0/2 y=((2*i-2)*t0/2+j)=1; y=((2*i-1)*t0/2+j)=0; end else for j=1:t0 y=((i-1)*t0+j)=0; end end end y=[y,x(i)]; plot(t,y); title(‘1 0 1 1 0 1 0 0 1 0’); grid on; axis([0,i,-0.1,1.1]); 4、双极性归零:将程序3中的for j=1:t0 y=((i-1)*t0+j)=0改为 for j=1:t0/2 y=((2*i-2)*t0/2+j)=-1; y=((2*i-1)*t0/2+j)=0; 运行结果: 7-4 假设随机二进制序列“10110001”，“1”码对应的基带波形为升余弦波形，持续的时间为Ts;”0”码对应的基带波形与“1”码相反。画出基带波形以及眼图。 程序代码: 1、Ts=1;N=15; eye_num=6; a=1;N_data=1000; dt=Ts/N, t=-3*Ts:dt:3*Ts; d=sign(randn(1,N_data)); dd=INSERT0(d,N); ht=sinc(t/Ts).*(cos(a*pi*t/Ts))./(1-4*a^)*t.^2/Ts^2+eps); st=conv(dd,ht); tt=-3*Ts:dt:(N_data+3)*N*dt-dt; subplot(2 1 1);plot(tt,st);axis([0 20 -1.2 1.2]); xlabel(‘t/Ts’);ylable(‘基带信号’);subplot(2 1 2); ss=zero(1,eye_num*N);ttt=0:dt:eye_num*N*dt-dt; for k=3:50 ss=st(k*N+1:(k+eye_num)*N); drawnow; plot(ttt,ss);hold on;end; xlabl(‘t/Ts’);ylable(‘基带信号眼图’); 2. function[out]=INSERT0(d,M) N=length(d); out=zeros(1,M*N); for i=0:N-1 out(i*M+1)=d(i+1); end; 运行结果: