下载文档原格式
Matlab平台下的通信系统建模与仿真2.1 Matlab平台简介及MATLAB的特点2.1.1 Matlab平台简介MATLAB软件的最初版本是由Cleve Moler博士等组成的名为MathWorks 公司在1992年推向市场的，软件名称MATLAB为矩阵（matrix）和实验室（laboratory）两个英文单词的前三个字母的组合（即矩阵实验室）。
所以，MATLAB非常适合矩阵运算,这也是MATLAB软件的最大特点之一。
2.1.2 MATLAB的特点（1）语言简洁，使用方便，运算符丰富，库函数极其丰富。
利用其丰富的库函数避开繁杂的子程序编程任务。
（2）语法限制不严格，程序设计自由。
例如，在MATLAB里，用户无需对矩阵预定义就可以使用。
（3）MATLAB既具有结构化的控制语句（如for循环、while循环、break 语句和if语句），又有面向对象编程的特性。
（4）程序的可移植性很好，基本上不做修改就可以在各种型号的计算机和操作系统上运行。
（5）MATLAB的图形功能强大。
数据的可视化非常简单,同时MATLAB 还具有较强的编辑图形界面的能力。
（6）MATLAB的缺点是软件成本高；其次MATLAB的程序不用编译等预处理，也不生成可执行文件，程序为解释执行，和其他高级程序相比，程序的执行速度较慢。
（7）MATLAB由核心部分和各种可选的工具箱两部分组成。
核心部分中包含数百个核心内部函数。
MATLAB的另一重大特色是功能强劲的工具箱。
其工具箱又可分为两类：功能性工具箱和学科性工具箱。
功能性工具箱能用于多种学科，主要用来扩充其符号计算功能、图示建模仿真功能、文字处理功能以及与硬件实时交互功能。
而学科性工具箱是专业性比较强的，如control、toolbox、signal processing toolbox、communication toolbox等。
（8）源程序的开放性。
2.2 Matlab平台下的系统建模与仿真方法在上一章我们介绍了OFDM系统的基本原理，阐明了如果在通信系统中一如OFDM可以提高系统的抗多径干扰能力，从而提高系统的稳定性和可靠性。
3.5 MATLAB 绘图实训3.5.1 实训目的1.学会MA TLAB 绘图的基本知识；2.掌握MA TLAB 子图绘制、图形注释、图形编辑等基本方法；3.学会通过MA TLAB 绘图解决一些实际问题；4.练习二维、三维绘图的多种绘图方式，了解图形的修饰方法；5.学会制作简单的MATLAB 动画。
图3-46 炮弹发射示意图3.5.2 实训内容1. 炮弹发射问题（1）炮弹发射的基础知识炮弹以角度α射出的行程是时间的函数，可以分解为水平距离)(t x 和垂直距离)(t y 。
)cos()(0αtv t x = %水平方向的行程；205.0)sin()(gt tv t y -=α %垂直方向的行程； 其中，0v 是初速度；g 是重力加速度，为9.82m/s ；t 是时间。
（2）炮弹发射程序举例：分析以下程序以及图3-47各个图形的实际意义。
a=pi/4; v0=300; g=9.8;t=0:0.01:50; x=t*v0*cos(a);y=t*v0*sin(a)-0.5*g*t.^2;subplot(221);plot(t,x);grid;title(‘时间-水平位移曲线'); subplot(222);plot(t,y);grid;title(‘时间-垂直位移曲线');subplot(223);plot(x,y);grid;title(‘水平位移-垂直位移曲线'); subplot(224);plot(y,x);grid;title(‘垂直位移-水平位移曲线');图3-4745角发射曲线 （3）编程解决炮弹发射问题①假设在水平地面上以垂直于水平面的角度向上发射炮弹，即发射角90=α，假设初速度分别为[310，290，270]m/s,试绘制时间-垂直位移曲线，编程求取最高射程；绘图要求：◆ 标题设为“炮弹垂直发射问题”； ◆ 在图上通过添加文本的方式表明初速度； ◆ 在x 轴标注“时间”； ◆ 在y 轴上标注“垂直距离”； ◆ 添加网格线； ◆ 将310m/s 的曲线改为线粗为2的红色实线； ◆ 将290m/s 的曲线改为线粗为3的绿色点划线； ◆ 将270m/s 的曲线改为线粗为2的蓝色长点划线；a=pi/2; v1=310; g=9.8;t=0:0.01:50; x1=t*v1*cos(a);y1=t*v1*sin(a)-0.5*g*t.^2;plot(t,y1);grid; title('炮弹垂直发射问题'); xlabel('时间'); ylabel('垂直距离'); hold on; v2=290;x2=t*v2*cos(a);y2=t*v2*sin(a)-0.5*g*t.^2; plot(t,y2); v3=270;x3=t*v3*cos(a);y3=t*v3*sin(a)-0.5*g*t.^2; plot(t,y3);zgsc=[max(y1); max(y2); max(y3)] %三次发射的最高射程 运行结果如下：zgsc =1.0e+003 * 4.9031 4.29083.7194最高射程分别为:4903.1米,4290.8米,3719.4米。
% ch3example1A.mclear;f_p=2400; f_s=5000; 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'); % 计算H(s)f=0:100:10000; % 计算频率点和频率范围s=j*2*pi*f; % s=jw=j*2*pi*fH_s=polyval(b,s)./polyval(a,s); % 计算相应频率点处H(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); % 也可用指令freqs直接画出H(s)的频率响应曲线。
% ch3example1B.mclear;f_p=2400; f_s=5000; 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'); % 计算H(s)f=0:100:10000; % 计算频率点和频率范围s=j*2*pi*f; % s=jw=j*2*pi*fH_s=polyval(b,s)./polyval(a,s); % 计算相应频率点处H(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); % 也可用指令freqs直接画出H(s)的频率响应曲线。
补充内容：模拟调制系统的MATLAB 仿真1。
抽样定理为了用实验的手段对连续信号分析,需要先对信号进行抽样(时间上的离散化），把连续数据转变为离散数据分析。
抽样（时间离散化）是模拟信号数字化的第一步.Nyquist 抽样定律：要无失真地恢复出抽样前的信号，要求抽样频率要大于等于两倍基带信号带宽.抽样定理建立了模拟信号和离散信号之间的关系,在Matlab 中对模拟信号的实验仿真都是通过先抽样，转变成离散信号，然后用该离散信号近似替代原来的模拟信号进行分析的。
【例1】用图形表示DSB 调制波形)4cos()2cos(t t y ππ= 及其包络线。
clf％％计算抽样时间间隔fh=1;%%调制信号带宽（Hz ）fs=100*fh ；％%一般选取的抽样频率要远大于基带信号频率，即抽样时间间隔要尽可能短. ts=1/fs ；％%根据抽样时间间隔进行抽样，并计算出信号和包络t=(0：ts:pi/2)'；%抽样时间间隔要足够小,要满足抽样定理.envelop=cos （2＊pi ＊t)；%%DSB 信号包络y=cos(2＊pi ＊t).＊cos(4*pi ＊t ）;％已调信号%画出已调信号包络线plot （t ，envelop ，'r ：','LineWidth',3）；hold onplot(t ，-envelop,'r:’,’LineWidth'，3）;%画出已调信号波形plot(t ，y,’b'，'LineWidth ’,3）;axis(［0,pi/2,-1,1］)%hold off%xlabel （'t ’);％写出图例【例2】用图形表示DSB 调制波形)6cos()2cos(t t y ππ= 及其包络线。
clf％％计算抽样时间间隔fh=1;％%调制信号带宽(Hz ）fs=100＊fh;%抽样时间间隔要足够小，要满足抽样定理。
matlab通信仿真实例在Matlab中进行通信系统的仿真，可以涉及到多种不同的通信技术和协议，包括调制解调、信道编码、多址接入等。
以下以OFDM系统为例，介绍Matlab 中通信仿真的实例。
OFDM（正交频分复用）是一种常用于现代通信系统中的技术，它将高速数据流分割成多个较低速的子流，并将每个子流分配到不同的子载波上。
优点是能够抵抗多径效应和频率选择性衰落，并提供高数据速率。
首先，我们需要创建一个包含OFDM系统参数的结构体。
例如：ofdmParam.M = 16; % 子载波数量ofdmParam.K = 4; % 用于混合多路复用的用户数量ofdmParam.N = ofdmParam.M * ofdmParam.K; % 总子载波数量ofdmParam.CP = 16; % 循环前缀长度接下来，我们可以生成用于OFDM仿真的数据流。
例如，我们可以使用随机整数生成器生成一系列整数，并将其转换为复数形式的调制符号：data = randi([0, ofdmParam.M-1], 1, ofdmParam.N);dataMod = qammod(data, ofdmParam.M);然后，我们可以创建一个包含OFDM信号的函数。
在OFDM系统中，生成的数据符号将分配到不同的子载波上，然后在时域中通过插入循环前缀进行叠加：function[ofdmSignal] = createOFDMSignal(dataMod, ofdmParam) ofdmSignal = [];for k = 0:ofdmParam.K-1% 提取相应的数据符号，并进行IFFTofdmData =ifft(dataMod(k*ofdmParam.M+1:(k+1)*ofdmParam.M));% 添加循环前缀ofdmDataWithCP = [ofdmData(end-ofdmParam.CP+1:end), ofdmData];% 将OFDM符号添加到OFDM信号中ofdmSignal = [ofdmSignal, ofdmDataWithCP];endend将OFDM信号传输到信道中，我们可以使用加性高斯白噪声（AWGN）信道模型来模拟实际通信环境：EbNo = 10; % 信噪比snr =10*log10(ofdmParam.N*ofdmParam.M/(ofdmParam.N*ofdmParam.M+1 )*(10^(EbNo/10)));ofdmSignalNoisy = awgn(ofdmSignal, snr, 'measured');最后，我们可以对接收到的OFDM信号进行解调和信号恢复。
实验一2：for i=1:50,j=1:50;A(i,j)=i+j-1;endC=flipud(A),B=fliplr(A),A(1:10,1:10)=zeros;D=A3：clear;a=255*rand(50);disp(a),a(find(a<128))=0;a(find(a>128))=255;w=a;disp(w)实验二：1：（1）clear;n=0:23;N=12;y1=sin(2*4*pi*n/N);y2=sin(2*5*pi*n/N);y3=sin(2*7*pi*n/N);y4=sin(2*10*pi*n/N);plot(n,y1,'+r-',n,y2,'-.',n,y3,'k-',n,y4,'*g-');axis([0 25 -1 2 ]);xlabel('横轴n');ylabel('纵轴y');text(0.8,0.9,'y1');text(19.6,0.9,'y2');text(14,0.9,'y3');text(5,0.9,'y4');title('plot tu xing');legend('y1=sin(2*4*pi*n/N)','y2=sin(2*5*pi*n/N)','y3=sin(2*7*pi*n/N)','y4=sin(2*10*pi*n/N)')（2）clear;n=0:23;N=12;y1=sin(2*4*pi*n/N);y2=sin(2*5*pi*n/N);y3=sin(2*7*pi*n/N);y4=sin(2*10*pi*n/N);stem(n,y1,'+r-');hold on;stem(n,y2,'b-');hold on;stem(n,y3,'k-');hold on;stem(n,y4,'*g-');axis([0,25,-2,4]);xlabel('横轴n');ylabel('纵轴y');text(0.8,1.2,'y1');text(3,1.4,'y2');text(14,1.2,'y3');text(10.9,1.3,'y4');title('n from 0 to 23');legend('y1=sin(2*4*pi*n/N)','y2=sin(2*5*pi*n/N)','y3=sin(2*7*pi*n/N)','y4=sin(2*10*pi*n/N)') 2：clear all;N=input('input N=');T=1;fs=1/T;N_sample=128;dt=1/N_sample;t=0:dt:10*T-dt;n=-N:N;Fn=(-1/2)*sinc(n/2).*exp(-j*n*pi/2);Fn(N+1)=0;ft=zeros(1,length(t));for m=-N:N;ft=ft+ Fn(m+N+1)*exp(j*m*pi*2*fs*t);endplot(t,ft);legend('时域曲线ft=ft+ Fn(m+N+1)*exp(j*m*pi*2*fs*t)',4)3：N_sample=128;dt=1/N_sample;t=0:dt:1-dt;st=[ones(1,N_sample/2),-ones(1,N_sample/2)];subplot(311);plot(t,st);axis([0 1 -2 2]);xlabel('t');ylabel('s(t）'); subplot(312);[f,sf]=T2F(t,st);plot(f,abs(sf));hold on;axis([-10 10 0 1]);xlabel('f');ylabel('|S(f)|');sff=(-1/2)*exp(-j*f*pi/2).*sinc(f/2);plot(f,abs(sff));[t,st]=F2T(f,sf);subplot(313);axis([0 1 -2 2]);xlabel('t');ylabel('恢复的s(t)');plot(t,st)子函数：（1）：function[t,st]=F2T(f,sf)df=f(2)-f(1);Fmx=f(end)-f(1)+df;dt=1/Fmx;N=length(sf);T=dt*N;t=0:dt:T-dt;sff=fftshift(sf);st=Fmx*ifft(sff);（2）：function[f,sf]=T2F(t,st)dt=t(2)-t(1);T=t(end);df=1/T;N=length(st);f=-N/2*df:df:N/2*df-df;sf=fft(st);sf=T/N*fftshift(sf);2图：M=4 M=5M=7 M=10M=4 M=5M=7 M=10N=10 N=10002468100246810实验三：1.（2）clear;p=2;N=5000;u=randn(5,4);c=sqrt(p);u=u*c;power_u=var(u);subplot(211)plot(u(1:20));grid on;ylabel('u(n)');xlabel('n');subplot(212)hist(u,50);grid on;ylabel('histogram of u(n)');2：clear all;N=500;disounce=(sign(rand(1,N)-0.5+eps)+1)/2; n=8;temp1=ones(1,n);temp0=zeros(1,n);new_disounce=[];for i=1:length(disounce)if disounce(i)==0new_disounce=[new_disounce temp0];elsenew_disounce=[new_disounce temp1];endendT=0.01;t=0:T/n:T/n*(length(new_disounce)-1);subplot(211);plot(t,new_disounce)axis([min(t)-0.1,max(t)+0.1,min(new_disounce)-0.1,max(new_disounce)+0.1]) y=new_disounce;Nf=length(t);Y=fft(y,Nf);Pyy=abs(Y).^2/Nf;dt=t(2)-t(1);df=1/dt;f=df*(0:(Nf-1)/2)/Nf;subplot(212);plot(f,10*log10(Pyy(1:((Nf-1)/2+1))));title('功率谱密度');xlabel('频率(HZ)');grid3：clear all;SNRindB1=0:1:12;SNRindB2=0:0.01:12;for i=1:length(SNRindB1)smld_err_prb(i)= exa0308_Fun_singlepe(SNRindB1(i));endsemilogy(SNRindB1,smld_err_prb,'r*');for i=1:length(SNRindB2)SNR=10^(SNRindB2(i)/10);theo_err_prb(i)=(1/2)*erfc(sqrt(SNR/2));endhold on;semilogy(SNRindB2,theo_err_prb);grid on子函数：function [p]=exa0308_Fun_singlepe(snr_in_dB)E=1;SNR=10^(snr_in_dB/10);sgma=E/sqrt(SNR*2);N=10000;dsource=(sign(rand(1,N)-0.5+eps)+1)/2;numoferr=0;for i=1:Nif (dsource(i)==0)r=-E+sgma*randn;elser=E+sgma*randn;endif(r<0.5*E)decis=0;elsedecis=1;endif (decis~=dsource(i)) numoferr=numoferr+1; end end p=numoferr/N;:实验3的图 1：>> x=2*rand(5,4) x =0.7375 1.9921 0.9260 0.5703 1.3593 1.8842 0.2846 1.5663 0.0472 1.6797 0.9457 1.9128 0.8051 0.0293 0.9237 0.82740.8319 1.4452 1.7069 0.359802468101214161820-4-2024u (n )nu =-1.2765 1.1282 0.4296 -1.4450 -1.3768 0.7019 1.3848 -3.8243 -1.9091 -1.3563 1.7813 -2.9299 -0.8431 0.1011 0.4232 0.9075 0.4459 0.3491 -1.6545 0.89682：0.20.40.60.811.21.41.61.820.51050100150200250300350400-400-200200功率谱密度频率(HZ)3：0246810121010101010100。
matlab无线通信仿真代码
本文分享了一份基于MATLAB平台的无线通信仿真代码。
该代码包含了常用的调制解调算法、信道模型、编码解码等模块，通过调用这些模块可以实现多种无线通信系统的仿真。
具体来说，该代码可以用于仿真以下几种无线通信系统：
1. 传统的模拟调制系统，包括AM、FM、PM等模式的调制解调。
2. 数字调制系统，包括PSK、QAM、FSK等模式的调制解调。
3. 多径衰落信道模型，包括AWGN信道、瑞利衰落信道、莱斯衰落信道等。
4. 信道编码解码，包括卷积码、Turbo码、LDPC码等。
除此之外，该代码还包含了多个演示程序，可以直观地展示无线通信系统的仿真结果。
这些演示程序涉及到信道容量、误码率、信噪比等指标，可以帮助用户深入理解无线通信系统的性能。
总之，该MATLAB无线通信仿真代码简单易用，功能齐全，可用于学术研究、教学实验、工程应用等多个领域。
- 1 -。
基于Matlab（Simulink）《通信原理》实验仿真（模拟部分）基于Matlab(Simulink)《通信原理》实验仿真(模拟部分)摘要模拟通信在通信系统中的使用非常广泛，而MATLAB(Simulink)是用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，它可以解决包括信号和图像处理、通讯等众多应用领域中的问题。
利用MATLAB集成环境下的M文件和Simulink工具箱可以完成通讯系统设计与仿真，本文主要是利用MATLAB集成环境下的M文件，编写程序来实现AM、FM、VSB调制与解调过程，并分别绘制出其信号波形。
再通过Simulink工具箱对模拟通信系统进行建模仿真。
Simulation of communication in a communication system is very extensive, and the use of MATLAB (Simulink) is used to algorithm development, data visualization, data analysis and numerical calculation of the senior technical calculation language and interactive environment, it can solve the including signal and image processing, communicationetc many applications in question. MATLAB integration environment Mfiles and Simulink tool box can complete communication system design and simulation, this paper is mainly use of MATLAB integration environment, programming of the M files to achieve AM, FM, VSB modulation and demodulation process, and separately plot its signal waveform. Again through Simulink communication system toolbox of simulation modeling simulation .关键词模拟信号;AM;FM;VSB;调制解调;MATLAB(Simulink)目录:第一章绪论第二章理论与方法2.1 matlab简介2.2 Simulink简介2.3 通信原理概述第三章设计方案3.1用 MATLAB的M文件进行模拟调制3.1.1 基于MATLAB的AM信号调制与解调3.1.2 基于MATLAB的FM信号调制与解调3.1.3 基于MATLAB的VSB信号调制与解调3.2用Simulink对模拟通信系统进行建模仿真3.2.1频分复用和超外差接收机的仿真模型3.2.2调频立体声接收机模型第四章小结参考文献致谢附录1:程序清单第一章绪论调制在通信系统中的作用至关重要。
% ch3example1A.mclear;f_p=2400; f_s=5000; 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'); % 计算H(s)f=0:100:10000; % 计算频率点和频率范围s=j*2*pi*f; % s=jw=j*2*pi*fH_s=polyval(b,s)./polyval(a,s); % 计算相应频率点处H(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); % 也可用指令freqs直接画出H(s)的频率响应曲线。
% ch3example1B.mclear;f_p=2400; f_s=5000; 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'); % 计算H(s)f=0:100:10000; % 计算频率点和频率范围s=j*2*pi*f; % s=jw=j*2*pi*fH_s=polyval(b,s)./polyval(a,s); % 计算相应频率点处H(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); % 也可用指令freqs直接画出H(s)的频率响应曲线。
% ch3example2A.mf_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); % 计算H(z)figure(1);freqz(b,a, 1000, 8000) % 作出H(z)的幅频相频图, freqz(b,a, 计算点数, 采样率)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); % 计算相应频率点处H(s)的值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');% ch3example3A.mf_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'); % 计算H(z)freqz(b,a, 1000, 8000) % 作出H(z)的幅频相频图, freqz(b,a, 计算点数, 采样率) subplot(2,1,1); axis([0 4000 -30 3])% ch3example4A.mf_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); % 计算H(z)freqz(b,a, 1000, 10000) % 作出H(z)的幅频相频图, freqz(b,a, 计算点数, 采样率) subplot(2,1,1); axis([0 5000 -30 3])% ch3example5A.mf_N=10000; %采样率f_p=[1000, 1500]; f_s=[1200, 1300]; R_p=3; R_s=30; %设计要求指标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'); %计算带阻H(z)系数freqz(b,a, 1000, 10000) %作出H(z)的幅频相频图, freqz(b,a, 计算点数, 采样率) subplot(2,1,1); axis([0 5000 -35 3])% ch3example6A.mFs=8000; Ts=1/8000; % 采样率f0=500; % 梳状滤波器开槽基频率bw=60/(Fs/2); % 归一化开槽带宽ab=-3 ; % 计算开槽带宽位置处的衰减分贝值n=Fs/f0; % 计算滤波器阶数[num,den] = iircomb(n,bw,ab, 'notch'); % 计算H(z)freqz(num,den, 4000, 8000); % 作出H(z)的幅频相频图axis([0 4000 -30 5]);% ch3example10.msim('ch3example10.mdl');bar([-10:1/5:10-1/5],histdata);axis([-10 10 0 0.05]);% 作出归一化统计直方图% ch3example10prg.msim('ch3example10.mdl');bar([-10:1/5:10-1/5],histdata*5);axis([-10 10 0 0.25]);% 作出归一化统计直方图% ch3example11prg1.mt=-1: 0.01 :10 ; % 计算时间范围f_t=exp(-t).*(t>0);subplot(3,1,1);plot(t,f_t); % 时域波形axis([-1 10 -0.1 1.1]); xlabel('time (sec)');w=-40: 0.1 : 40; % 计算角频率范围F_w=1./(1+j*w); % 频谱理论结果subplot(3,1,2);plot(w, abs(F_w)); % 频域幅度谱axis([-40 40 0 1.1]);xlabel('freq (rad/s)');subplot(3,1,3);plot(w, angle(F_w)); % 频域相位谱axis([-40 40 -pi/2 pi/2]);xlabel('freq (rad/s)');% ch3example11prg2.mw_m=40; % 截断频率T=pi/w_m; % 采样间隔L=5;t=0: T :L ; % 时域截断x_t=exp(-t).*(t>0); %信号序列N=length(x_t); %序列长度(点数)%-----------------------------------X_k=fft(x_t); % FFT计算%-----------------------------------w0=2*pi/(N*T); % 离散频率间隔kw=2*pi/(N*T) .*[0: N-1]; % 离散频率样点X_kw=T.*X_k; % 乘以T得到连续傅利叶变换频谱的样值%-----------------------------------subplot(2,1,1);plot(kw, abs(X_kw),'.','MarkerSize', 14); % 作出数值计算的幅度谱点axis([-40 90 -0.1 1.1]);xlabel('freq (rad/s)');hold on; % 保持前面作图曲线不被擦除w=-40: 0.1 : 40;X_w=1./(1+j*w); % 理论计算频谱表达式plot(w, abs(X_w)); % 作图对比subplot(2,1,2);plot(kw, angle(X_kw),'.','MarkerSize', 14); % 作出数值计算的相位谱点hold on;plot(w, angle(X_w)); % 频域相位谱axis([-40 90 -pi pi]);xlabel('freq (rad/s)');% ch3example12prg1.mDf=5; % 频率间隔f_s=2*500; % 采样率N=f_s/Df; % 序列点数t=0:1./f_s:(N-1)./f_s; % 计算时间段freq=0:Df:(N-1)*Df; % 计算频率段f_t=2*sin(2*pi*100*t)+cos(2*pi*180*t); % 信号F_f=1/f_s*fft(f_t,N); % 用FFT计算频谱plot(freq-f_s/2,abs(fftshift(F_f))); % 将零频率移动到FFT中心xlabel('频率Hz'); ylabel('幅度谱'); % 并作出幅度频谱% ch3example14prg1.mfs=200; % 采样率Delta_f = 1; % 频率分辨率T = 1/fs; % 时间分辨率L=1/ Delta_f; % 时域截取长度N= floor(fs/Delta_f)+1; % 计算截断信号的采样点数t=0:T:L; % 截取时间段和采样时间点freq=0: Delta_f :fs; % 分析的频率范围和频率分辨率f_t=(sin(2*pi*50*t)+0.7*sin(2*pi*75*t))';% 在截取范围内的分析的信号时域波形f_t_rectwin= rectwin(N).*f_t./sqrt(sum(abs(rectwin(N).^2))./N);% 矩形窗（功率归一化）f_t_hamming= hamming(N) .*f_t./sqrt(sum(abs(hamming(N).^2))./N);% 海明窗（功率归一化）f_t_hann = hann(N) .*f_t./sqrt(sum(abs(hann(N).^2))./N);% 汉宁窗（功率归一化）F_w_rectwin =T.* fft(f_t_rectwin, N);% 进行N点FFT，并乘以采样时间间隔T得到频谱F_w_hamming =T.* fft(f_t_hamming, N); % 加海明窗的频谱F_w_hann =T.* fft(f_t_hann, N); % 加汉宁窗的频谱figure(1); subplot(2,2,1);plot(t,f_t);title('Original Signal');subplot(2,2,2);plot(t, f_t_rectwin);title('Rectwin Windowing');subplot(2,2,3);plot(t, f_t_hamming);title('hamming Windowing');subplot(2,2,4);plot(t, f_t_hann);title('hanning Windowing');figure(2);subplot(3,1,1);plot(freq, 10*log10(abs(F_w_rectwin)));title('Rectwin Windowing Spectrum');ylabel('幅度dB');axis([0,200,-50,0]);grid on;subplot(3,1,2);plot(freq, 10*log10(abs(F_w_hamming)));title('hamming Windowing Spectrum');ylabel('幅度dB');axis([0,200,-50,0]);grid on;subplot(3,1,3);plot(freq, 10*log10(abs(F_w_hann)));title('hanning Windowing Spectrum');ylabel('幅度dB');axis([0,200,-50,0]);grid on;p_original_signal=var(f_t) % 计算原始信号的平均功率p_hannwindowed=sum((abs(F_w_hann/T)).^2)/(N^2) % 计算加窗后的信号功率% ch3example16prg1.mfs=500; % 采样率Df=1; % 频率分辨率N=floor(fs/Df)+1; % 计算的序列点数t=0:1/fs:(N-1)/fs; % 截取信号的时间段F=0:Df:fs; % 功率谱估计的频率分辨率和范围xk=sin(2*pi*50*t)+2*sin(2*pi*130*t)+randn(1,length(t));% 截取时间段上的离散信号样点序列Pxx=abs(fft(xk)).^2/(N^2); % 功率谱估计Pav_timedomaim=sum(xk.^2)/N % 在时域计算信号功率Pav_freqdomain=sum(Pxx) % 通过功率谱计算信号功率plot(F,10*log10(Pxx));xlabel('freq Hz');ylabel('PSD dB') % 作出功率谱密度图% ch3example16prg2.mfs=500; % 采样率Df=1; % 频率分辨率N=floor(fs/Df)+1; % 计算的序列点数t=0:1/fs:(N-1)/fs; % 截取信号的时间段F=0:Df:fs; % 功率谱估计的频率分辨率和范围xk=sin(2*pi*50*t)+2*sin(2*pi*130*t)+randn(1,length(t));% 截取时间段上的离散信号样点序列Pxx=(abs(fft(xk(1:167))).^2+...abs(fft(xk(168:334))).^2+...abs(fft(xk(335:501))).^2)/3/((N/3)^2);Pav_timedomaim=sum(xk.^2)/N % 在时域计算信号功率Pav_freqdomain=sum(Pxx) % 通过功率谱计算信号功率plot(0:3:fs,10*log10(Pxx));xlabel('freq Hz');ylabel('PSD dB');%作出功率谱密度图% ch3example16prg3.mfs=500; % 采样率Df=1; % 频率分辨率N=floor(fs/Df)+1; % 计算的序列点数t=0:1/fs:(N-1)/fs; % 截取信号的时间段F=0:Df:fs; % 功率谱估计的频率分辨率和范围xk=sin(2*pi*50*t)+2*sin(2*pi*130*t)+randn(1,length(t));% 截取时间段上的离散信号样点序列Pxx=(abs(fft(xk(1:167))).^2+...abs(fft(xk(83:249))).^2+...abs(fft(xk(168:334))).^2+...abs(fft(xk(250:416))).^2+...abs(fft(xk(335:501))).^2)/5/((N/3)^2);Pav_timedomaim=sum(xk.^2)/N % 在时域计算信号功率Pav_freqdomain=sum(Pxx) % 通过功率谱计算信号功率plot(0:3:fs,10*log10(Pxx));xlabel('freq Hz');ylabel('PSD dB');%作出功率谱密度图% ch3example16prg4.mfs=500; % 采样率Df=1; % 频率分辨率N=floor(fs/Df)+1; % 计算的序列点数t=0:1/fs:(N-1)/fs; % 截取信号的时间段F=0:Df:fs; % 功率谱估计的频率分辨率和范围xk=sin(2*pi*50*t)+2*sin(2*pi*130*t)+randn(1,length(t));% 截取时间段上的离散信号样点序列w=hamming(167)';%海明窗w=w*sqrt(167/sum(w.*w));% 使得海明窗与矩形窗等能量，即加窗后不对信号功率产生影响Pxx=(abs(fft(w.*xk(1:167))).^2+...abs(fft(w.*xk(83:249))).^2+...abs(fft(w.*xk(168:334))).^2+...abs(fft(w.*xk(250:416))).^2+...abs(fft(w.*xk(335:501))).^2)/5/((N/3)^2);Pav_timedomaim=sum(xk.^2)/N % 在时域计算信号功率Pav_freqdomain=sum(Pxx) % 通过功率谱计算信号功率plot(0:3:fs,10*log10(Pxx));xlabel('freq Hz');ylabel('PSD dB');%作出功率谱密度图% ch3example16prg5.mfs=500; % 采样率Df=1; % 频率分辨率N=floor(fs/Df)+1; % 计算的序列点数t=0:1/fs:(N-1)/fs; % 截取信号的时间段F=0:Df:fs; % 功率谱估计的频率分辨率和范围xk=sin(2*pi*50*t)+2*sin(2*pi*130*t)+randn(1,length(t));% 截取时间段上的离散信号样点序列[Pxx,F] = psd(xk,512,500,hamming(256),128);plot(F,10*log10(Pxx/(512/2))); xlabel('freq Hz');ylabel('PSD dB')Pav=sum(Pxx/(512/2)) % 通过功率谱计算信号功率% ch3example23prg1.mkc=150e3; % Hz/V VCO控制灵敏度omega_n=2*pi*15e3/2.5; % PLL自然角频率K=2*pi*(0.5*kc); % 估算环路增益zeta=1; % 临界阻尼tau_1=K/((omega_n).^2);tau_2=2*zeta/omega_n-1/K;freq=0:10:100e3; % 计算频率范围0到100KHzs=j*2*pi*freq;G_s=(1+tau_2*s)./(1+tau_1*s); % 环路滤波器传递函数figure(1);semilogx(freq,(abs(G_s)));% 作出环路滤波器的频率响应title('环路滤波器幅频响应');xlabel('Hz');ylabel('|G(s)|');grid on;b=[tau_2,1]; % 环路滤波器分子系数向量a=[tau_1,1]; % 环路滤波器分母系数向量H_s=(G_s*K./s)./(1+G_s*K./s);figure(2);semilogx(freq,20*log10(abs(H_s)));% 作出闭环频率响应title('PLL线性相位模型闭环频率响应');xlabel('Hz');ylabel('20log|H(s)|(dB)'); grid on;。
