Unity开发者联盟收藏本站
发帖教程新人必读切换到宽版

 找回密码
 立即注册
Unity开发者联盟»论坛 › 3D Dev Asset 开发资源大全 › 学术研究 Scholarly › matlab信号措置学习（实战代码）
返回列表 发新帖
查看: 210|回复: 0

matlab信号措置学习（实战代码）

[复制链接]
发表于 2024-7-15 18:40 | 显示全部楼层 |阅读模式
没有下载matlab可以打开网页版Octave，很好用Octave Online · Cloud IDE compatible with MATLAB (octave-online.net)
part1
  1. %创建正弦波
  2. %定义信号采样序列。从0s到1s每隔0.001s采样一次，共采样1000次
  3. t=0:0.001:1;
  4. y1=sin(2*pi*50*t);
  5. y2=sin(2*pi*50*t)+2*sin(2*pi*100*t);
  6. y3=y2+0.5*randn(size(t));%加一个噪声信号
  7. %这里的两个50要保持一致，50是取50个采样，最多取t中定义的1000个
  8. subplot(131),plot(t(1:50),y1(1:50)),title('50HZ正弦波')
  9. subplot(132),plot(t(1:50),y2(1:50)),title('50HZ+100HZ正弦波')
  10. subplot(133),plot(t(1:20),y3(1:20)),title('50HZ+100HZ噪声正弦波')
复制代码
  1. %创建周期性方波
  2. t=linspace(0,3*pi)';
  3. x=square(t);%操作函数创建方波信号x
  4. plot(t,x);
  5. xlabel('t/\pi')%t/π,横轴
  6. xlabel('t/pi')%t/pi
  7. %方波是一种非正弦曲线的波形，震荡电路出来的正弦波一般都有谐波分量，方波就是由一系列的谐波分量叠加而成，抱负的方波只有0-1两个值
复制代码
  1. %创建非周期性方波
  2. t=-4:0.001:4;
  3. y=rectpuls(t);%发生非周期方波信号，默认方波宽度为1
  4. z=square(t);%发生周期性的矩形脉冲信号
  5. subplot(131)
  6. plot(t,z,'b','linewidth',4);
  7. axis([-4,4,-1.5,1.5]);
  8. grid on;
  9. xlabel t,ylabel n(t);
  10. subplot(132)
  11. plot(t,y,'r','linewidth',4)
  12. axis([-4,4,-1.5,1.5]);
  13. grid on;
  14. xlabel t,ylabel h(t);
  15. y=2*rectpuls(t,2);%发生指定宽度为2的非周期方波信号。前面的*2是高度*2
  16. subplot(133)
  17. plot(t,y,'linewidth',4)
  18. axis([-4,4,-0.5,2.5]);
  19. grid on;
复制代码
  1. %创建锯齿波
  2. f=50;
  3. T=5*(1/50);%采样时间0.1s，信号周期为5，频率为50HZ
  4. fs=1000;%定义采样率为1000
  5. t=0:1/fs:T-1/fs;
  6. y=sawtooth(2*pi*50*t,0);%周期为2pi，0是向左倾斜的锯齿波
  7. z=sawtooth(2*pi*50*t,1);%周期为2pi，1是向右倾斜的锯齿波
  8. subplot(211);
  9. plot(t,y);
复制代码
  1. %创建指数波
  2. fs=100;
  3. t=0:1/fs:1;
  4. C=1;
  5. a=1;
  6. y=C*exp(a*t);
  7. subplot(211);
  8. plot(t,y);
  9. C=-1;
  10. a=-1;
  11. y=C*exp(a*t);
  12. subplot(212);
  13. plot(t,y,'.');
复制代码
  1. %创建斜坡信号,线性增长信号
  2. fs=10;
  3. t=0:1/fs:1;
  4. x=t;
  5. plot(t,x,'r');
  6. hold on
  7. stem(t,x,'bo');%增状图定名stem
  8. xlabel('时间序列')；
  9. ylabel('x(t)');
  10. title('斜坡序列')
复制代码
  1. %创建插值序号——解决由于信号采样频率降低发生的频率混叠
  2. fs=100;
  3. t=0:1/fs:2-1/fs;
  4. x=sin(2*pi*10*t.^3).^2;%外侧.^2变密,内侧.^3各个幅度变乱犯警则
  5. y=interp(x,4);%输入插值信号函数,将插值信号采样率增加一个因子4
  6. tiledlayout(2,1);
  7. nexttile;
  8. stem(t,x,'filled',markersize,3)
  9. grid on
  10. xlabel('Sample Number')
  11. ylabel('Original')
  12. nexttile;%跟subplot的感化一样
  13. t=1:length(y);
  14. stem(t,y,'filled','markersize',1)
  15. xlabel('Sample Number')
  16. ylabel('Interpolated')
复制代码
  1. %创建添加噪声的正弦波
  2. fs=1000;
  3. t=0:1/fs:1-1/fs;
  4. y1=sin(2*pi*t);
  5. y2=wgn(1000,1,0);%创建1000*1的高斯白噪声y2,噪声样本功率为0
  6. y3=sin(2*pi*t)+wgn(1000,1,0);%创建添加噪声的正弦波
  7. y4=y1+0.5*randn(size(t));%在正弦波上添加正态分布的随机噪声信号，噪声信号的幅值是正弦波的一半
  8. subplot(221),plot(t,y1),title('正弦波信号');
  9. subplot(222),plot(t,y2),title('高斯白噪声信号');
  10. subplot(223),plot(t,y3),title('添加高斯白噪声的正弦波信号');
  11. subplot(224),plot(t,y4),title('添加正态分布噪声的正弦波信号');
复制代码
  1. %创建添加噪声的指数波
  2. fs=1000;
  3. t=0:1/fs:1-1/fs;
  4. x=airy(t*10).*exp(-t,^2);
  5. y1=x+wgn(1000,1,1);%添加高斯白噪声，噪声功率为1
  6. y2=awgn(x,1,'measured');%高斯白噪声信噪比为1，measured指在插手噪声前测定信号强度
  7. subplot(311),plot(t,x),title('指数波信号');
  8. subplot(312),plot(t,y1),title('叠加高斯白噪声的指数波信号');
  9. subplot(313),plot(t,y2),title('添加高斯白噪声的指数波信号');
复制代码
  1. %创建随机波
  2. t=0:100;
  3. N=length(t);
  4. x=rand(1,N);%创建随机信号
  5. y=randn(1,N);%创建正态分布的随机信号
  6. subplot(131);
  7. plot(t,x,'k');
  8. ylabel('x(t)');
  9. subplot(132);
  10. plot(t,y,'r');
  11. ylabel('y(t)');
  12. subplot(133);
  13. stem(t,x,'filled','g');
  14. ylabel('x2(t)');
复制代码
  1. %创建sinc波，sinc(t)=sin(t)/t
  2. x=linspace(-5,5);
  3. y1=sin(pi*x)./(pi*x);
  4. y2=sinc(x);
  5. plot(x,y1,'o',x,y2);
  6. xlabel Time,ylabel Signal;
  7. %在图形中添加图例便于区分曲线
  8. legend('sin signal','sinc signal','Location','NorthWest')
  9. legend boxoff;%图例不显示边框
复制代码
  1. %创建对数扫频信号
  2. t=0:0.001:1;
  3. y=chirp(t,1e-6,1,50,'logarithmic');%初始时刻瞬时频率1e-6，t=1时，瞬时频率为50HZ，扫频方式为logarithmic的对数扫频
  4. plot(t,y);
  5. axis([0.8,1,-1,1]);%调整坐标范围
复制代码
  1. %创建迪利克雷信号diric
  2. x=0:0.001:4*pi;
  3. subplot(211);
  4. plot(x,diric(x,7));%计算x的7级迪利克雷函数，并绘制xinhaobox
  5. axis tight,title('n=7');
  6. subplot(212);
  7. plot(x/pi,diric(x,8));
  8. title('n=8'),xlabel('x/\pi');
复制代码
  1. %创建三角波信号
  2. t=-3:0.1:3; %采样时间6s
  3. y=tripuls(4*t,2*pi,0.2);%采样次数为4t，周期为2pi,斜率为0.2的非周期三角波
  4. z=sawtooth(4*t,0);%周期为2pi的周期性三角波（锯齿波）
  5. subplot(121);
  6. plot(t,y,'r^','linewidth',2);
  7. grid on;
  8. subplot(122);
  9. plot(t,z,'linewidth',2);
  10. grid on;
复制代码
<hr/>part2
啦啦啦啦部门加上跑的图啦
  1. %熟悉二维图形的绘制的一些东西
  2. %显示4*4图形分割
  3. t1=(0:11)/11*pi;
  4. t2=(0:400)/400*pi;%数字越大，越密集
  5. t3=(0:50)/50*pi;
  6. y1=cos(t1).*cos(5*t1);
  7. y2=cos(t2).*cos(5*t2);
  8. y3=cos(t3).*cos(5*t3);
  9. subplot(221),plot(t1,y1,&#39;r.&#39;);
  10. title(&#39;(1)点过少的离散图形&#39;);
  11. subplot(222),plot(t1,y1,t1,y1,&#39;r.&#39;);
  12. title(&#39;(2)点过少的持续图形&#39;);
  13. subplot(223),plot(t2,y2,&#39;r.&#39;);
  14. title(&#39;(3)点过多的离散图形&#39;);
  15. subplot(224),plot(t3,y3);
  16. title(&#39;(4)实线画持续图形&#39;);
复制代码
  1. %图窗布局的应用
  2. x=linspace(-pi,pi);
  3. y=sin(x);
  4. tiledlayout(2,2);%2*2共四个图块
  5. nexttile;%激活第一个图块
  6. plot(x);
  7. nexttile;%激活第一个图块
  8. plot(x,y);
  9. nexttile([1,2]);%创建第三个图块，该图块占据一行两列的坐标区
  10. plot(x,y);
复制代码
  1. %保持命令的应用
  2. x=linspace(-pi,pi);
  3. y1=sin(x).*exp(x);
  4. plot(x,y1);
  5. hold on;%使得第二条曲线叠加在第一条曲线上面，而不是直接覆盖
  6. y2=cos(x).*sin(x);
  7. plot(x,y2);
  8. hold off %封锁图像保持命令，此时绘制第三条曲线就会覆盖掉前两个
  9. y3=2*sin(3*x);
  10. plot(x,y3)
复制代码
  1. %调整坐标系
  2. x=0:pi/100:2*pi;
  3. y=cos(x).^5;
  4. plot(x,y);
  5. axis([0,pi,-2,2]);
复制代码
  1. %绘制三角函数
  2. x=linspace(0,10*pi,100);
  3. plot(x,cos(x).^2+sin(x));
  4. title(&#39;三角函数之和&#39;);
  5. xlabel x坐标,ylabel y坐标
复制代码
  1. %绘制函数图形并标注
  2. x=linspace(-5,5);
  3. y=x.^3+exp(x);
  4. plot(x,y);
  5. xt=[-3,3];%设置标注位置
  6. yt=[-24,46];%设置标注位置
  7. str={&#39;loccal min&#39;,&#39;local max&#39;};
  8. text(xt,yt,str);%添加标注
复制代码
  1. %绘制向量图形，并标出函数名
  2. plot(1:10);
  3. gtext(&#39;My Plot&#39;,&#39;Color&#39;,&#39;red&#39;,&#39;FontSize&#39;,14);
复制代码
  1. %添加绘图注释
  2. x=0:0.1:5;
  3. y1=exp(0.5*x).*cos(2*x);
  4. y2=cos(x).^5;
  5. plot(x,y1,&#39;r-&#39;,x,y2,&#39;b-.&#39;);
  6. title(&#39;函数曲线&#39;); %添加名字
  7. legend(&#39;函数1&#39;,&#39;函数2&#39;); %添加图例
  8. xlabel 自变量x,ylabel 函数值y;
  9. grid on;%显示网格
复制代码
<hr/>part3
%傅里叶变换与谱分析！！！本节有很多matlab&#39;中的函数没有在octave中实现！以下代码建议在matlab中跑！
按照本身的理解，大部门逆变换都相当于对变换的撤销动作

  1. %信号变换的本质大将信号从时域转化为频域
  2. fs=100;
  3. N=1000;
  4. n=0:N-1;
  5. t=n/fs;
  6. x=diric(t,7);%随时间变化的迪利克雷信号
  7. subplot(311);
  8. plot(t,x);
  9. title(&#39;原始函数&#39;);
  10. subplot(312);
  11. y=fft(x);%快速傅里叶变换
  12. plot(t,y);
  13. title(&#39;傅里叶变换&#39;);
  14. subplot(313);
  15. z=ifft(y);
  16. plot(t,z);
  17. title(&#39;傅里叶逆变换&#39;);%发现z函数和x函数一样
复制代码

  1. %计算周期方波信号与非周期方波信号的频域变换
  2. fs=100;
  3. t=0:1/fs:2;
  4. x=[square(2*pi*t);rectpuls(2*pi*t)];%定义函数矩阵,两个元素分袂是周期方波信号与非周期方波信号
  5. subplot(311);
  6. plot(t,x,&#39;linewidth&#39;,4);
  7. axis([0 2 -1.5 1.5]);
  8. title(&#39;原始信号&#39;);
  9. subplot(312);
  10. y=fft2(x);%二维快速傅里叶变换
  11. plot(t,y);
  12. title(&#39;傅里叶变换后信号&#39;);
  13. axis([0 2 -11 20]);
  14. subplot(313);
  15. z=ifft2(y);
  16. plot(t,z);
  17. title(&#39;二维傅里叶逆变换&#39;);%一维傅里叶变换是正弦信号的叠加，而二维傅里叶变换是正弦平面波的叠加。
  18. axis([0 2 -1.5 1.5]);
复制代码


在matlab中，颠末fft变换后，数据的频率范围是从[0,fs]摆列的。而一般，我们在画图或者讨论的时候，是从[-fs/2,fs/2]的范围进行分析。因此，需要将颠末fft变换后的图像的[fs/2,fs]部门移动到[-fs/2,0]这个范围内。而fftshift就是完成这个功能。凡是，如果想得到所见的中间是0频的图像，颠末fft变换后，都要再颠末fftshift。


  1. %计算余弦波噪声信号的傅里叶变换
  2. fs=100;%信号采样率100HZ
  3. t=0:1/fs:1-1/fs;%信号采样时间序列
  4. f=1;%设置信号频率为1
  5. x=cos(2*pi*f*t)+rand(1,100);%在余弦信号中添加随机噪声
  6. subplot(311);
  7. plot(t,x);
  8. title(&#39;原始信号&#39;);
  9. subplot(312);
  10. y=fftshift(x);%对信号x进行0频分量转移计算
  11. plot(t,y);
  12. title(&#39;傅里叶变换&#39;)
  13. subplot(313);
  14. z=ifftshift(y);%对信号x进行0频分量转移逆运算
  15. plot(t,z);
  16. title(&#39;傅里叶逆变换&#39;);
复制代码
  1. %(N维傅里叶变换)对信号进行fft变换
  2. fs=100;%信号采样率100HZ
  3. f1=100;f2=200;%定义两个信号的频率
  4. t=0:1/fs:1;%信号采样时间序列
  5. x1=sin(2*pi*100*t)+cos(2*pi*200*t);%原始x
  6. subplot(311);
  7. plot(x1);
  8. title(&#39;原始信号&#39;);
  9. x=sin(2*pi*100*t)+cos(2*pi*200*t)+rand(size(t));%x收到噪声污染
  10. subplot(312),plot(x(1:50));%绘制时域内的信号，噪声污染后很难看出正弦波的成分
  11. title(&#39;噪声污染信号&#39;);
  12. Y=fft(x,512);%对x进行512点的傅里叶变换
  13. f=fs*(0:256)/512;
  14. subplot(313),plot(f,Y(1:257));%绘制频域内信号，可以明显看出信号中100HZ和200HZ的两个频率分量
  15. title(&#39;加噪声信号频域变换&#39;);
复制代码

%离散傅里叶变换
  1. N=20;%定义信号采样点数
  2. n=0:N-1;%定义信号采样时间序列
  3. x=exp(n);%定义原始指数信号
  4. D=dftmtx(N);%求DFT矩阵
  5. X=x*D;
  6. subplot(211);
  7. stem(n,x),title(&#39;原始信号x&#39;);
  8. subplot(212);
  9. stem(abs(X));title(&#39;离散傅里叶变换X&#39;);
复制代码
  1. %复倒谱分析 matlab跑
  2. %复倒谱是指一个信号的傅里叶变换对数的傅里叶反变换
  3. fs=1000;
  4. t=0:1/fs:1-1/fs;
  5. x=square(20*pi*t);%创建方波信号
  6. x=vco(x,[0.1 0.5]*fs,fs);%对方波信号进行调频，第二个参数是调频范围！！！！！！！！注意此处VCO功能属于Octave Forge的信号包，但尚未实现
  7. y=fft(x);%快速傅里叶变换
  8. z=cceps(x);%复倒谱分析
  9. ax1=subplot(311),plot(t,x);
  10. title(&#39;原始信号&#39;);
  11. ax2=subplot(312);
  12. plot(t,y);
  13. title(&#39;原始信号进行fft&#39;);
  14. ax3=subplot(313);
  15. plot(t,z);
  16. title(&#39;原始信号复倒谱&#39;);
复制代码
  1. %傅里叶同步压缩变换  在matlab中跑(m跑)
  2. fs=1000;
  3. t=0:1/fs:1;
  4. x=square(20*pi*t)+0.2*randn(size(t));%受污染的方波信号
  5. subplot(211),plot(t,x);
  6. title(&#39;原始信号&#39;);
  7. subplot(212);
  8. fsst(x,fs,&#39;yaxis&#39;);%计算并绘制信号的傅里叶同步压缩变换
  9. title(&#39;原始信号傅里叶同步压缩clc变换&#39;);
复制代码
  1. %对信号进行傅里叶同步压缩变换和逆变换 m跑
  2. fs=1000;
  3. t=0:1/fs:3;
  4. x=chirp(t,30,2,5).*exp(-(2*t-3).^2)+2;%创建高斯调制啁啾信号，信号直流值为2，初始啁啾频率为30HZ，2s后衰减到5Hz
  5. subplot(311),plot(t,x);%绘制时域内原始信号
  6. xlabel(&#39;Time(s)&#39;);
  7. title(&#39;原始信号x&#39;);
  8. subplot(312);
  9. [y,f]=fsst(x,fs); %计算信号的傅里叶同步压缩变换
  10. fsst(x,fs); %绘制信号的傅里叶同步压缩变换
  11. xlabel(&#39;Time(s)&#39;);
  12. title(&#39;原始信号傅里叶同步压缩变换&#39;);
  13. z=ifsst(y);%进行傅里叶同步压缩逆变换，反改变换重建信号
  14. t=(0:length(y)-1)/fs; %定义信号变换后的采样时间序列
  15. subplot(313),plot(t,z);%绘制重构信号
  16. xlabel(&#39;Time(s)&#39;);
  17. title(&#39;重构信号&#39;);
  18. %成果是第一个和第三个信号一致
复制代码
  1. %计算信号dft变换（使用goertzel算法）
  2. %goertzel算法目的是从给定的采样中求出某一特定频率信号的能量，用于有效性评价
  3. fs = 1000;
  4. f1=10e3;f2=1.24e3;
  5. t=0:1/fs:1;
  6. x= cos(2*pi*t*f1)+cos(2*pi*t*f2)+randn(size(t));%受污染的余弦波信号
  7. subplot(2,1,1),plot(t,x);
  8. title(&#39;原始信号x&#39;)
  9. N = (length(x)+1)/2; %计算采样点数
  10. f = (fs/2)/N*(0:N-1); %计算频率
  11. indxs = find(f>1.2e3 & f<1.3e3);%限制频率范围在1.2~1.3kHZ之间
  12. y= goertzel(x,indxs);%在指定频率范围内对信号使用二阶goertzel算法计算DFT
  13. subplot(2,1,2),plot(t,y);
  14. title(&#39;原始信号DFT变换&#39;);
复制代码
%短时傅里叶变换
  1. %短时傅里叶变换用于确按时变信号其局部区域正弦波的频率与相位；可以理解为对一段长信号，截取每一段时间的短信号做fft，将得到的频谱图时间沿时间轴摆列，即可得到时频的云图。用于信号截短的东西叫 窗函数（宽度相当于时间长度），窗越小，时域特性越明显，频域特性越不明显。
  2. fs = 1000;  
  3. t = -4:1/fs:4;
  4. x = sinc(t); %辛格信号
  5. x = vco(x,[0.1 0.4]*fs,fs); %对sinc信号进行调频，调频范围是[0.1 0.4]
  6. %进行短时傅里叶变换，使用长度为256的Kaiser窗口，形状参数为5，指定重叠长度为220个样本，DFT长度为512个点
  7. stft(x,fs,&#39;Window&#39;,kaiser(256,5),&#39;OverlapLength&#39;,220,&#39;FFTLength&#39;,512);
  8. view(-45,65); %更改视图，将STFT显示为瀑布图
  9. colormap jet; %颜色设置为JET
复制代码




瀑布图
  1. %对信号执行短时傅里叶变换和逆变换
  2. fs = 1000;
  3. t = 0:1/fs:2-1/fs;
  4. f0 = 100;%初始瞬时频率
  5. f1 = 500;%参考时间瞬时频率
  6. x =chirp(t,f0,1,f1,&#39;quadratic&#39;,[],&#39;concave&#39;);%先行扫频余弦信号，扫频方式为quadratic二次扫频，[]暗示忽略初始相位，谱图形状为凸抛物线
  7. win = hamming(100,&#39;periodic&#39;); %设计一个长度为100的周期性hamming窗
  8. noverlap = 80; %样本重叠数为80
  9. subplot(3,1,1)
  10. plot(t,x)
  11. xlabel(&#39;时间t&#39;);ylabel(&#39;原始信号x(t)&#39;);
  12. title(&#39;原始信号&#39;);
  13. subplot(3,1,2)
  14. [y,f,t] = stft(x,fs,&#39;Window&#39;, win,&#39;OverlapLength&#39;,noverlap);%使用重叠样本数80计算信号的短时傅里叶变换
  15. plot(f,y) ;
  16. xlabel(&#39;频率f&#39;);ylabel(&#39;变换信号y(t)&#39;);
  17. title(&#39;STFT变换&#39;);
  18. subplot(3,1,3)
  19. [iy,t] = istft(y,fs,&#39;Window&#39;, win,&#39;OverlapLength&#39;,noverlap);
  20. plot(t,iy) %使用重叠样本数80计算信号的短时傅里叶逆变换
  21. xlabel(&#39;时间t&#39;);ylabel(&#39;变换信号iy(t)&#39;);
  22. title(&#39;ISTFT变换&#39;);
复制代码
  1. %绘制短时傅里叶变换谱图spectrogram
  2. t = -2:1/1e3:2;%操作线性间隔值函数定义信号采样时间序列，信号持续时间4s
  3. f = 10;
  4. y = sin(2*pi*f*t.^3);
  5. spectrogram(y);
复制代码
  1. %对比 窗函数与长度 对信号的短时傅里叶变换谱图 的影响
  2. %加窗本色是用一个所谓的窗函数与原始的时域信号作乘积的过程（当然加窗也可以在频域进行，但时域更为遍及），使得相乘后的信号似乎更好地满足傅立叶变换的周期性要求。
  3. fs = 3000;
  4. t = 0:1/fs: 1-1/fs;
  5. x= chirp(t,100,1,400,&#39;quadratic&#39;) + chirp(t,350,1,50);%创建叠加的先行扫频余弦信号
  6. %创建四个长度分歧的taylorwin窗口
  7. window1= taylorwin(256);       
  8. window2= taylorwin(512);
  9. window3= taylorwin(1024);
  10. window4= taylorwin(2056);
  11. %创建四个长度分歧的矩形窗口
  12. window11= rectwin(256);
  13. window21= rectwin(512);
  14. window31= rectwin(1024);
  15. window41= rectwin(2056);
  16. tiledlayout(5,2)%五行两列的分布图布局
  17. nexttile;
  18. spectrogram(x,32);%将信号分为32个样本段，计算并绘制信号的短时傅里叶变换谱图
  19. nexttile;
  20. spectrogram(x,window1,&#39;yaxis&#39;);%操作长度为256泰勒窗计算并绘制信号的短时傅里叶变换谱图，指定频率轴绘制在y轴
  21. nexttile;
  22. spectrogram(x,window1);
  23. nexttile;
  24. spectrogram(x,window2);
  25. nexttile;
  26. spectrogram(x,window3);
  27. nexttile;
  28. spectrogram(x,window4);
  29. nexttile;
  30. spectrogram(x,window11);
  31. nexttile;
  32. spectrogram(x,window21);
  33. nexttile;
  34. spectrogram(x,window31);
  35. nexttile;
  36. spectrogram(x,window41);
  37. %结论：对于短时傅里叶变换，最重要的是窗口长度的拔取。当频域刻度和平移步长足够密时，增加的只是生成图像的大小，但是物理层面的分辩率没有改变。改变物理层面分辩率的是窗口长度
  38. %窗口长度大：频率清晰，时间模糊
  39. %窗口长度小：频率模糊，时间清晰
复制代码
  1. %对比重叠样本数对短时傅里叶变换谱图的影响
  2. fs = 1000;
  3. t = 0:1/fs: 1-1/fs;
  4. f = 2;
  5. x= exp(1j*pi*sin(2*pi*f*t)*32);
  6. noverlap1=10; %定义分歧的信号重叠样本数
  7. noverlap2=20;
  8. tiledlayout(3,1)
  9. nexttile;
  10. spectrogram(x,32,noverlap1,64,&#39;centered&#39;,&#39;yaxis&#39;);%分成32个样本段绘制信号的短时傅里叶变换中心双边谱图，重叠样本数10（重叠样本数必需小于样本段数32，频率显示在y轴）
  11. nexttile;
  12. spectrogram(x,32,noverlap2,64,&#39;centered&#39;,&#39;yaxis&#39;); %重叠样本数20，对比重叠样本数对频谱的影响
  13. nexttile;
  14. spectrogram(x,32,noverlap1,64,&#39;centered&#39;,&#39;reassigned&#39;,&#39;yaxis&#39;);%计算并绘制信号的短时傅里叶变换中心双边谱图，将能量集中
复制代码
  1. %绘制短时傅里叶变换单边与双边谱图
  2. fs = 1000;
  3. t = 0:1/fs: 1-1/fs;
  4. f = 2;
  5. x= chirp(t,100,1,400,&#39;quadratic&#39;);%扫频余弦信号
  6. Nsample=32; %样本分割段数
  7. noverlap=10; %信号重叠样本数
  8. nfft=64; %DFT点数
  9. tiledlayout(3,1)
  10. nexttile;
  11. spectrogram(x,Nsample,noverlap,64,&#39;onesided&#39;,&#39;yaxis&#39;); %onesided将信号分成32个样本段，绘制信号的短时傅里叶变换单边谱图，频率显示在y轴
  12. nexttile;
  13. spectrogram(x,Nsample,noverlap,64,&#39;centered&#39;,&#39;yaxis&#39;);%centered将信号分成32个样本段，绘制信号的短时傅里叶变换中心双边谱图
  14. nexttile;
  15. spectrogram(x,Nsample,noverlap,64,&#39;twosided&#39;,&#39;yaxis&#39;); %twosided估计信号从头分配的双边谱图
复制代码
  1. %创建对称凸二次啁啾信号
  2. fs=1000;
  3. t = -2:1/fs:2;
  4. f0 = 100;
  5. f1 = 200;
  6. x = chirp(t,f0,1,f1,&#39;quadratic&#39;,[],&#39;convex&#39;);
  7. y = vco(x,[0.1  0.8]*fs,fs); %调频
  8. tiledlayout(2,2)
  9. nexttile;
  10. pspectrum(x,t,&#39;spectrogram&#39;,&#39;TimeResolution&#39;,0.1, ...
  11.     &#39;OverlapPercent&#39;,99,&#39;Leakage&#39;,0.85)%计算并绘制啁啾信号的谱图。指定相邻段之间99%的重叠和0.85的频谱泄露。
  12. title(&#39;频谱图&#39;)
  13. nexttile;
  14. stft(x,fs); %计算并绘制信号的STFT
  15. title(&#39;短时傅里叶变换&#39;)
  16. nexttile;
  17. spectrogram(x,100,80,1024,fs,&#39;yaxis&#39;) %计算并绘制啁啾的短时傅里叶变换谱图。将信号分成100个样本段，指定相邻段之间重叠的80个样本、DFT长度为1024个样本
  18. title(&#39;短时傅里叶变换谱图&#39;)
  19. nexttile;
  20. xspectrogram(x,y,100 ,80,1024,fs,&#39;yaxis&#39;) %计算并绘制啁啾的交叉功率谱图。将信号分成100个样本段，指定相邻段之间重叠的80个样本、DFT长度为1024个样本
  21. title(&#39;短时傅里叶变换交叉功率谱图&#39;)
复制代码

part4 频谱估计
  1. %采样点数对正弦波信号的傅里叶变换的影响
  2. fs = 100;
  3. N1=100;
  4. N2=200;
  5. n1 = 0:N1-1;
  6. n2 = 0:N2-1;
  7. t1=n1/fs;
  8. t2=n2/fs;
  9. x1 = sin(2*pi*4*t1);%两个正弦波信号
  10. x2 = sin(2*pi*4*t2);
  11. subplot(3,1,1)
  12. plot(t1,x1,&#39;*&#39;,t2,x2,&#39;r&#39;)
  13. xlabel(&#39;时间/s&#39;);
  14. legend(&#39;N1=100&#39;,&#39;N2=200&#39;);
  15. title(&#39;原始函数&#39;)
  16. grid on;
  17. subplot(3,1,2)
  18. y1=fft(x1,N1);
  19. y2=fft(x2,N2);
  20. plot(t1,y1,t2,y2)
  21. xlabel(&#39;时间/s&#39;);
  22. title(&#39;傅里叶变换&#39;)
  23. legend(&#39;N1=100&#39;,&#39;N2=200&#39;);
  24. grid on;
  25. subplot(3,1,3)
  26. z1=abs(y1);%计算信号1傅里叶变换后的振幅
  27. z2=abs(y2);%计算信号2傅里叶变换后的振幅
  28. f1=(0:N1-1)*fs/N1;%定义频率序列
  29. f2=(0:N2-1)*fs/N2;
  30. plot(f1(1:N1/2),z1(1:N1/2)*2/N1,f2(1:N2/2),z2(1:N2/2)*2/N2);
  31. xlabel(&#39;频率/Hz&#39;);ylabel(&#39;振幅&#39;);
  32. title(&#39;傅里叶幅值变换&#39;)
  33. legend(&#39;N1=100&#39;,&#39;N2=200&#39;);
  34. grid on;
复制代码
  1. %绘制正弦波形的单边和双边频谱
  2. f0 = 100;
  3. fs = 5000;
  4. n=1:1:1000;
  5. N = length(n);
  6. t =n/fs;
  7. x = sin(2*pi*f0*t);
  8. subplot(3,1,1)
  9. plot(n,x)
  10. title(&#39;原始信号&#39;)
  11. y =fft(x);
  12. y2=abs(y/N); %计算信号双边谱的幅度！！！！公式记得
  13. subplot(3,1,2)
  14. plot(n,y2)
  15. xlabel(&#39;f(Hz)&#39;)
  16. ylabel(&#39;|P2(f)|&#39;)
  17. title(&#39;FFT信号双边频谱&#39;)
  18. y1=y2(1:N/2+1);
  19. %实际上，双边谱是由单边谱/2所得频谱，然后正负对称频率叠加到一起，而零频率不需要/2
  20. y1(2:end-1)=2*y1(2:end-1); %计算信号单边谱幅度！！！！公式记得
  21. f = fs*(0:N/2)/N; %计算信号频谱幅度所对应的频率
  22. subplot(3,1,3)
  23. plot(f,y1)
  24. xlabel(&#39;f(Hz)&#39;)
  25. ylabel(&#39;|P1(f)|&#39;) %信号单边频谱
  26. title(&#39;FFT信号单边频谱&#39;)
复制代码
  1. %绘制正弦信号的相位谱
  2. fs=1000; %定义采样序列为1000Hz
  3. N=1024; %定义采样点数
  4. n=0:N-1; %定义采样点数序列
  5. t=n/fs; %定义采样时间序列
  6. y=sin(2*pi*100*t); %定义正弦信号
  7. Y=fft(y,N); %计算FFT信号
  8. A=abs(Y); %求FFT信号的幅值
  9. f=n*fs/N; %计算信号频率
  10. subplot(311)
  11. %奈奎斯特频率：为了防止信号混叠，需要定义的最小采样频率，此频率称为奈奎斯特频率。奈奎斯特频率是离散信号系统采样频率的一般f_nyquist=fs/2.按照对称性，用fft做信号频谱分析，只需要考察0~f_nyquist范围内的幅频特性即可f_nyquist=fs*(0:N/2)/N
  12. plot(f(1:N/2),A(1:N/2)); %绘制奈奎斯特频率相位谱
  13. xlabel(&#39;频率/hz&#39;),ylabel(&#39;幅值&#39;),title(&#39;幅值谱&#39;);
  14. grid on;
  15. subplot(312)
  16. ph=2*angle(Y(1:N/2)); %求信号的相位
  17. plot(f(1:N/2),ph(1:N/2));  %绘制频率相位谱
  18. xlabel(&#39;频率/hz&#39;),ylabel(&#39;相位&#39;),title(&#39;相位谱&#39;);
  19. grid on;
  20. subplot(313)
  21. ph1=ph*180/pi; %将信号相位由弧度转化为角度
  22. plot(f(1:N/2),ph1(1:N/2)); %绘制频率相位谱
  23. xlabel(&#39;频率/hz&#39;),ylabel(&#39;相位&#39;),title(&#39;相位谱&#39;);
  24. grid on;
复制代码
  1. %创建余弦波FFT相位谱
  2. fs = 100;
  3. t = 0:1/fs: 4-1/fs;
  4. f1=10;f2=30;f3=40;
  5. x=cos(2*pi*f1*t)+2*cos(2*pi*f2*t-pi/4)+4*cos(2*pi*f3*t+pi/2); %创建3个余弦波组成的信号。第一个余弦波相位是0，幅值是1；第二个余弦波相位为-pi/4，幅值是2；第三个余弦波相位为pi/2，幅值是4
  6. x= awgn(x,10,&#39;measured&#39;); %创建添加高斯白噪声的余弦波，高斯白噪声功率值为10
  7. subplot(221)
  8. lx = length(x); %获取傅里叶变换后信号的长度
  9. f = (-lx/2:lx/2-1)/lx*fs;
  10. stem(f,abs(x)),title(&#39;原始信号序列图&#39;)%计算信号的幅值，绘制频率-幅值函数
  11. xlabel &#39;Frequency (Hz)&#39;
  12. ylabel &#39;|x|&#39;
  13. grid
  14. y = fft(x);
  15. subplot(222)
  16. stem(f,abs(y)),title(&#39;FFT信号序列图&#39;)%计算FFT信号的幅值，绘制频率-幅值函数
  17. xlabel &#39;Frequency (Hz)&#39;
  18. ylabel &#39;|y|&#39;
  19. grid
  20. subplot(223)
  21. z = fftshift(y);%将傅里叶变换后信号中的零频分量移到频谱中心
  22. stem(f,abs(z)),title(&#39;FFT频谱序列图&#39;) ;%计算信号的傅里叶变换，将变换幅值绘制为频率函数
  23. xlabel &#39;Frequency (Hz)&#39;
  24. ylabel &#39;|z|&#39;
  25. grid
  26. subplot(224)
  27. tol=1e-6; %定义小幅值变换值
  28. z(abs(z)<tol)= 0;  %删除小幅值变换值
  29. theta = angle(z);  %计算变换的相位
  30. stem(f,theta/pi),title(&#39;FFT相位频谱序列图&#39;)  %绘制频率-相位序列图
  31. xlabel &#39;Frequency (Hz)&#39;
  32. ylabel &#39;Phase / \pi&#39;
  33. grid
复制代码
  1. %光谱泄露对正弦信号频谱的影响
  2. fs = 100;
  3. t = 0:1/fs:2-1/fs;
  4. x=sin(2*pi*20*t);
  5. subplot(311)
  6. pspectrum(x,t,&#39;Leakage&#39;,1) %将泄漏增加到最大值1，解决了紧密间隔的调子，但掩盖了附近的微弱调子
  7. subplot(312)
  8. pspectrum(x,t,&#39;Leakage&#39;,0) %将泄漏降低到最小值0，以牺牲光谱分辩率为代价，将泄漏降低到最小
  9. subplot(313)
  10. pspectrum(x,t,&#39;Leakage&#39;,0.5) %将泄露显示默认值0.5
复制代码
  1. %信号功率谱单边谱与双边谱（三角波）
  2. f = 50;
  3. T = 10*(1/f);%10个信号周期
  4. fs = 1000;
  5. t = 0:1/fs:T-1/fs;
  6. x = sawtooth(2*pi*f*t,1/2);%发生基频为50Hz的三角波的10个周期
  7. subplot(221)
  8. plot(t,x)
  9. title(&#39;三角波信号&#39;)
  10. subplot(222)
  11. pspectrum(x,fs,&#39;Leakage&#39;,0.91) %默认情况下绘制功率谱的单边谱，泄漏量0.91
  12. subplot(223)
  13. pspectrum(x,fs,&#39;Leakage&#39;,1,&#39;TwoSided&#39;,true) %TwoSided绘制双边光谱，泄漏量为1
  14. subplot(224)
  15. [p,f,t] = pspectrum(x,fs,&#39;spectrogram&#39;); %默认情况下绘制单边谱
  16. waterfall(f,t,p&#39;); %瀑布图形式显示光谱
  17. xlabel(&#39;Frequency (Hz)&#39;)
  18. ylabel(&#39;Time (seconds)&#39;)
  19. view([60 45])
  20. colormap hot %颜色设置为JEF
复制代码
  1. %创建对称凸二次啁啾波频谱
  2. fs=1e3;
  3. t = -2:1/fs:2;
  4. f0 = 100;
  5. f1 = 200;
  6. x = chirp(t,f0,1,f1,&#39;quadratic&#39;,[],&#39;convex&#39;);%余弦扫频信号
  7. subplot(221)
  8. pspectrum(x,t,&#39;spectrogram&#39;,&#39;TimeResolution&#39;,0.1, ...
  9.     &#39;OverlapPercent&#39;,99,&#39;Leakage&#39;,0.85) %计算绘制啁啾谱图，将信号分割成分段，时间分辩率为0.1s，相邻段之间重叠为99%，频谱泄漏量为0.85
  10. subplot(222)
  11. pspectrum(x,fs,&#39;spectrogram&#39;, ...
  12.     &#39;FrequencyLimits&#39;,[100 300],&#39;TimeResolution&#39;,1)%计算信号谱图，频带设置在100~300，谱图时间分辩率为1
  13. subplot(223)
  14. pspectrum(x,fs,&#39;FrequencyLimits&#39;,[100 300],&#39;Leakage&#39;,1) %计算信号功率谱，频带设置在100~300，频谱泄露量为1
  15. subplot(224)
  16. pspectrum(x,fs,&#39;persistence&#39;, ...
  17.     &#39;FrequencyLimits&#39;,[100 300],&#39;TimeResolution&#39;,1) %计算信号的持续谱，频带设置在100~300，谱图时间分辩率为1
复制代码
  1. %绘制信号谱图
  2. fs = 1000;
  3. t=0:1/fs:1;
  4. x= exp(2j*pi*100*cos(2*pi*2*t)) + randn(size(t))/100;%受噪声污染的指数波信号
  5. subplot(2,1,1)
  6. [sp,fp,tp] = pspectrum(x,fs,&#39;spectrogram&#39;,&#39;FrequencyResolution&#39;,10); %指定频率分辩率为10Hz，计算光谱图[光谱，频谱频率，时间值]
  7. mesh(tp,fp,sp) %三位网格图命令绘制光谱图
  8. xlabel(&#39;Time (s)&#39;)
  9. ylabel(&#39;Frequency (Hz)&#39;)
  10. subplot(2,1,2)
  11. waterfall(fp,tp,sp&#39;);%绘制光谱瀑布图
  12. xlabel(&#39;Frequency (Hz)&#39;)
  13. ylabel(&#39;Time (seconds)&#39;)
  14. view([30 45])%调整视图视角
复制代码
  1. %创建信号的谱图和重分配谱图
  2. fs = 1000;
  3. t=0:1/fs:1;
  4. x= (1+0.5*cos(2*pi*5*t)).*cos(2*pi*50*t+5*sin(2*pi*10*t));%余弦波信号
  5. x=hilbert(x);%对信号进行希尔伯特变换，返回离散时间解析信号
  6. subplot(2,2,1)
  7. pspectrum(x,fs,&#39;spectrogram&#39;); %计算绘制光谱图
  8. subplot(2,2,2)
  9. pspectrum(x,fs,&#39;spectrogram&#39;,&#39;FrequencyResolution&#39;,10,&#39;Reassign&#39;,true) %指定频率分辩率带宽10Hz，计算绘制从头分配的光谱图
  10. subplot(2,2,3)
  11. pspectrum(x,fs,&#39;spectrogram&#39;,&#39;TimeResolution&#39;,0.2) %用0.2s的时间分辩率计算光谱图
  12. subplot(2,2,4)
  13. pspectrum(x,fs,&#39;spectrogram&#39;,&#39;TimeResolution&#39;,0.2,&#39;Reassign&#39;,true)%用0.2s的时间分辩率计算绘制从头分配的光谱图
复制代码

%功率谱几个概念：
%1.信号以波的形式存在着，会发生功率。单元频带的信号功率称为功率谱，它显示在必然区域中信号功率随着频率变化的分布情况。
%2.对于确定性的信号，出格长短周期的确定性信号，常用能量谱描述。对随机信号，无法用确定的时间函数来暗示，也就无法用频谱暗示，往往用功率谱描述其频率特性
%3.功率谱是功率谱密度函数的简称，是单元频率的信号功率（单元频率内的信号能量）。暗示了信号功率随着频率的变化情况；

%能量谱几个概念
%1.对于任意的时间信号x(t),这个信号可以是任意随时间变化的物理量。在对信号进行能量分析时，不加区分的将其视为 1欧姆的电阻上的电流。这个单元电阻的能量属性，就是这个信号的能量属性
%2.如果信号是能量信号，通过傅里叶变换，就很容易分手分歧频域分量所对应的能量。

%功率谱密度PSD的几个概念
%1.定义了信号或时间序列的功率如何随频率分布
%2.cpsd命令用于计算信号的交叉功率谱密度
  1. %创建信号交叉功率谱
  2. fs = 1000;
  3. t = 0:1/fs: 1-1/fs;
  4. f0 = 200;
  5. x = sin(2*pi*f0*t).^2;
  6. subplot(3,1,1)
  7. plot(t,x) %绘制时域图
  8. xlabel(&#39;时间序列t&#39;);ylabel(&#39;原始信号x(t)&#39;);
  9. title(&#39;原始信号&#39;);
  10. subplot(3,1,2)
  11. pspectrum(x,t,&#39;Leakage&#39;,1) %绘制泄漏增为最大值的品偶
  12. title(&#39;信号频谱&#39;);
  13. subplot(3,1,3)
  14. r = wgn(1000,1,0);%创建高斯白噪声
  15. cpsd(x,r,triang(500),128,2048,fs); %绘制信号三角窗分割的交叉功率谱，重叠样本数128，DFT点数2048，采样频率fs
  16. title(&#39;信号CPSD&#39;);
复制代码
  1. %分析余弦信号的交叉功率谱
  2. fs = 1000;
  3. t = 0:1/fs:2-1/fs;
  4. x = cos(2*pi*t*100)+0.25*randn(size(t));
  5. tau = pi/100; %定义信号相位
  6. y = cos(2*pi*100*(t-tau))+0.25*randn(size(t));
  7. cpsd(x,y,[],[],[],fs)  %计算并绘制交叉功率谱密度的幅值
复制代码

part 5 脉冲信号分析

%脉冲信号时指瞬间俄然变化，感化时间极短的电压或电流，可以是周期性反复的，也可以长短周期或单次的。脉冲信号是一种离散信号（y轴不持续），形状多种多样。
  1. %创建三角脉冲序列
  2. t = 0:1/1e3:60;%定义信号采样率为1000Hz，信号采样时间为60s
  3. f=0.05;%脉冲信号频率为0.05Hz
  4. d = [0:2:60;sin(2*pi*f*(0:2:60))]&#39;;%信号偏移量
  5. y = pulstran(t,d,&#39;tripuls&#39;);%通过信号偏移生成非周期三角波脉冲序列
  6. plot(t,y,&#39;*&#39;);
  7. xlabel(&#39;Time (s) &#39;), ylabel(&#39;Waveform&#39;)
复制代码
  1. %创建方波脉冲序列
  2. t = 0:1/1e3:60;
  3. d = 0:1/4:1;%偏移量
  4. y = pulstran(t,d,&#39;rectpuls&#39;,0.1);%生成非周期性的矩形脉冲信号
  5. subplot(121)
  6. plot(t,y)
  7. xlabel(&#39;t&#39;)
  8. ylabel(&#39;h(t)&#39;)
  9. grid on
  10. axis([0 1 -0.1 1.1])
  11. t1=0:1/1E3:1;
  12. d1=0:1/3:1;
  13. x = square(t1,0.5);%生成方波信号
  14. y1=pulstran(t1,d1,x,1e3);%生成方波脉冲信号序列，采样率1000Hz
  15. subplot(122)
  16. plot(t1,y1)
  17. xlabel(&#39;t&#39;)
  18. ylabel(&#39;h(t)&#39;)
  19. grid on
  20. axis([0 1 -5 2])
复制代码

%最常见的脉冲波就是方波。脉冲信号可以用来暗示信息，也可以用来作为载波。
%抱负脉冲序列，是理论上没有外来信号干扰的抱负状态合成的
  1. %绘制抱负采样信号序列
  2. n=0:50;%定义采样时间序列长度是50
  3. C=444.128;%定义指数波常量C
  4. a=50*sqrt(2.0)*pi;%定义指数波常量a
  5. T=0.001;%定义采样率f0为1000，采样时间间隔T为采样率的倒数
  6. w0=50*sqrt(2.0)*pi;%定义正弦波角频率
  7. x=C*exp(-a*n*T).*sin(w0*n*T);%采样信号函数
  8. subplot(211),plot(n,x,&#39;g*&#39;);%绿色*绘制
  9. title(&#39;采样信号序列&#39;);
  10. subplot(212),stem(n,x, &#39;LineStyle&#39;,&#39;-.&#39;,...
  11.      &#39;MarkerFaceColor&#39;,&#39;red&#39;,...
  12.      &#39;MarkerEdgeColor&#39;,&#39;yellow&#39;);%绘制x（n）的火柴杆图像，线形为-.，标识表记标帜填充颜色为红色，轮廓颜色为黄色
  13. title(&#39;抱负采样信号序列&#39;);
复制代码

%单元脉冲信号，又称单元采样序列。特点是，n=0时取值1，n=其他值时取值为0
%单元冲激信号是在t=0时刻，取值无穷大，在整个区间内对时间的积分为1，意义是 强度为1
%单元阶跃信号是在n<0时，值为0；n>=0时，值为1（单元阶跃是单元脉冲的累加）
  1. %创建单元脉冲信号
  2. %单元脉冲信号可以暗示为h=[1,zeros(1,n-1)]:zero(m)是生成m阶全零矩阵
  3. N=10;%定义信号序列的采样点数10
  4. x=[1,zeros(1,N-1)];%创建单元脉冲信号序列
  5. stem(x);%输入针线图命令stem，绘制单元脉冲信号序列
  6. hold on
  7. x=[1,zeros(1,N-1)];%创建单元脉冲信号序列
  8. x(1)=0;%定义信号中冲击值，n=1时取值为0
  9. t=linspace(-10,0,10);
  10. stem(t,x);%绘制单元脉冲信号序列
  11. title(&#39;单元脉冲信号序列&#39;);
  12. axis([-10 10 -0.5 1.5])
复制代码

part 5 脉冲信号的分析

脉冲信号是指 瞬间俄然变化，感化时间极短的电压或电流，可以是周期反复的，也可以长短周期单次的，是一种离散信号。pulstran命令用来生成离散信号
  1. %创建三角脉冲序列
  2. t = 0:1/1e3:60;%采样时间
  3. f=0.05;
  4. d = [0:2:60;sin(2*pi*f*(0:2:60))]&#39;;%偏移量
  5. y = pulstran(t,d,&#39;tripuls&#39;);%生成序列
  6. plot(t,y,&#39;*&#39;);
  7. xlabel(&#39;Time (s) &#39;), ylabel(&#39;Wavefomat
复制代码
  1. %创建方波脉冲信号
  2. t = 0:1/1e3:60;
  3. d = 0:1/4:1;  %偏移量
  4. y = pulstran(t,d,&#39;rectpuls&#39;,0.1);%采样率为0.1，非周期性的矩形脉冲信号
  5. subplot(121)
  6. plot(t,y)%脉冲信号时序图
  7. xlabel(&#39;t&#39;)
  8. ylabel(&#39;h(t)&#39;)
  9. grid on
  10. axis([0 1 -0.1 1.1])
  12. t1=0:1/1E3:1;
  13. d1=0:1/3:1;  
  14. x = square(t1,0.5);
  15. y1=pulstran(t1,d1,x,1e3);%1e3 暗示采样率为1000，方波脉冲信号序列
  16. subplot(122)   
  17. plot(t1,y1)
  18. xlabel(&#39;t&#39;)
  19. ylabel(&#39;h(t)&#39;)
  20. grid on
  21. axis([0 1 -5 2])
复制代码

实际上 合成的序列中还会存在噪声干扰，抱负脉冲序列用stem函数生成
  1. %绘制抱负采样信号序列
  2. n=0:50;%定义采样时间序列的长度是50
  3. C=444.128;%指数波常量
  4. a=50*sqrt(2.0)*pi;%指数波常量
  5. T=0.001;%采样时间为0.001，采样频率为1000
  6. w0=50*sqrt(2.0)*pi;%角频率
  7. x=C*exp(-a*n*T).*sin(w0*n*T);%采样信号函数
  8. subplot(211),plot(n,x,&#39;g*&#39;);
  9. title(&#39;采样信号序列&#39;);
  10. subplot(212),stem(n,x, &#39;LineStyle&#39;,&#39;-.&#39;,...
  11.      &#39;MarkerFaceColor&#39;,&#39;red&#39;,...
  12.      &#39;MarkerEdgeColor&#39;,&#39;yellow&#39;);
  13. title(&#39;抱负采样信号序列&#39;);
复制代码

本帖子中包含更多资源

您需要 登录 才可以下载或查看，没有账号？立即注册

×
回复

使用道具 举报

返回列表 发新帖
懒得打字嘛，点击右侧快捷回复 【右侧内容，后台自定义】
高级模式
B Color Image Link Quote Code Smilies
您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

小黑屋|手机版|Unity开发者联盟 ( 粤ICP备20003399号 )

GMT+8, 2024-11-11 05:43 , Processed in 0.104234 second(s), 27 queries .

Powered by Discuz! X3.5 Licensed

© 2001-2024 Discuz! Team.

快速回复 返回顶部 返回列表