基于+MATLAB+的时域信号采样及频谱分析

一、课程设计题目：基于MATLAB的时域信号采样及频谱分析二、基本要求：①掌握数字信号处理的基本概念、基本理论和基本方法；②学会MATLAB的使用，掌握MATLAB的程序设计方法；③学会用MATLAB对信号进行分析和处理；④信号的各参数需由键盘输入，输入不同参数即可得不同的x(t)和x(n)；⑤撰写课程设计论文，用数字信号处理基本理论分析结果。三、设计方法与步骤：①画出连续时间信号的时域波形及其幅频特性曲线，其中幅度因子A＝444.128，衰减因子a＝222.144，模拟角频率)()sin()(0tutAetxatΩ=−0Ω＝222.144；②对信号进行采样，得到采样序列，其中T＝)(tx500),()sin()(0≤Ω=−nnunTAenxanTsf1为采样间隔，通过改变采样频率可改变T，画出采样频率分别为200Hz，500Hz，1000Hz时的采样序列波形；③对不同采样频率下的采样序列进行频谱分析，绘制其幅频和相频曲线，对比各频率下采样序列和的幅频曲线有无差别，如有差别说明原因。)(nx)(tx④设系统单位抽样响应为)()(5nRnh=，求解当输入为时的系统响应，画出,,的时域波形及幅频特性曲线，并利用结果验证卷积定理的正确性（此内容将参数设置为A＝1，a＝0.4，)(nx)(ny)(nx)(nh)(ny0Ω＝2.0734，T＝1）。⑤用FFT对信号,,进行谱分析，观察与④中结果有无差别。)(nx)(nh)(ny⑥由采样序列恢复出连续时间信号，画出其时域波形，对比与原)(nx)(1tx)(1tx连续时间信号的时域波形，计算并记录两者最大误差。)(tx四、详细程序及仿真波形分析：1、连续时间信号x(t)及其200Hz/500Hz/1000Hz频率抽样信号函数x(n)%绘制信号x(n)的幅度谱和相位谱n=0:50;%定义序列的长度是50A=input('请输入A的值A:');%设置信号的有关参数a=input('请输入a的值a:');w0=input('请输入w0的值w0:');T1=0.005;T2=0.002;T3=0.001;T0=0.001;x=A*exp(-a*n*T0).*sin(w0*n*T0);%pi是MATLAB定义的π，信号乘可采用“.*”y1=A*exp(-a*n*T1).*sin(w0*n*T1);%pi是MATLAB定义的π，信号乘可采用“.*”y2=A*exp(-a*n*T2).*sin(w0*n*T2);%pi是MATLAB定义的π，信号乘可采用“.*”y3=A*exp(-a*n*T3).*sin(w0*n*T3);%pi是MATLAB定义的π，信号乘可采用“.*”closeall%清除已经绘制的x(n)图形subplot(2,1,1);stem(n,x),gridon%绘制x(n)的图形title('离散时间信号')subplot(2,1,2);plot(n,x),gridontitle('连续时间信号')figure(2)subplot(3,1,1);stem(n,y1),gridontitle('200Hz理想采样信号序列');%设置结果图形的标题subplot(3,1,2);stem(n,y2),gridontitle('500Hz连续时间信号')subplot(3,1,3);stem(n,y3),gridontitle('1000Hz连续时间信号')k=-25:25;W=(pi/12.5)*k;w=W/pi;Y1=y1*(exp(-j*pi/12.5)).^(n'*k);figure(3)subplot(2,1,1);plot(w,abs(Y1));grid,xlabel('w'),ylabel('幅度');title('200Hz理想采样信号序列的幅度谱');axis([-2201000]);subplot(2,1,2);plot(w,angle(Y1));grid,xlabel('w'),ylabel('幅角');title('200Hz理想采样信号序列的相位谱')Y2=y2*(exp(-j*pi/12.5)).^(n'*k);figure(4)subplot(2,1,1);plot(w,abs(Y2));grid,xlabel('w'),ylabel('幅度');title('500Hz理想采样信号序列的幅度谱');axis([-2201000]);subplot(2,1,2);plot(w,angle(Y2));grid,xlabel('w'),ylabel('幅角');title('500Hz理想采样信号序列的相位谱')Y3=y3*(exp(-j*pi/12.5)).^(n'*k);figure(5)subplot(2,1,1);plot(w,abs(Y3));grid,xlabel('w'),ylabel('幅度');title('1000Hz理想采样信号序列的幅度谱');axis([-2201000]);subplot(2,1,2);plot(w,angle(Y3));grid,xlabel('w'),ylabel('幅角');title('1000Hz理想采样信号序列的相位谱')分析：采样频率为1000Hz时没有失真，500Hz时有横线，产生失真，200Hz时横线加长，失真增大。说明采样频率越大失真越小。2、设系统单位抽样响应为，求解当输入为时的系统响应，画出,,的时域波形及幅频特性曲线，并利用结果验证卷积定理的正确性（此内容将参数设置为A＝1，a＝0.4，＝2.0734，T＝1）。)()(5nRnh=)(nx)(ny)(nx)(nh)(ny0Ωn=1:50;%定义序列的长度是50hb=zeros(1,50);%注意：MATLAB中数组下标从1开始hb(1)=1;hb(2)=1;hb(3)=1;hb(4)=1;hb(5)=1;closeall;subplot(3,1,1);stem(hb);title('系统hb[n]');m=1:50;T=1;%定义序列的长度是和采样率A=1;a=0.4;T=1;w0=2.0734;x=A*exp(-a*m*T).*sin(w0*m*T);%pi是MATLAB定义的π，信号乘可采用“.*”subplot(3,1,2);stem(x);title('输入信号x[n]');y=conv(x,hb);subplot(3,1,3);stem(y);title('输出信号y[n]');figure(2)subplot(3,1,1);plot(n,hb),gridontitle('矩形序列时域波形');subplot(3,1,2);plot(m,x),gridontitle('输入信号x[n]时域波形');subplot(3,1,3);plot(m,y),gridontitle('输出信号y[n]时域波形');分析：在数字信号处理中经常要进行卷积运算，MATLAB中有一个内部函数conv可以计算两个有限长序列的卷积，该函数计算的两个序列都是从n=0开始3、用FFT对信号,,进行谱分析，观察与④中结果有无差别。)(nx)(nh)(nyn=1:50;%定义序列的长度是50hb=zeros(1,50);%注意：MATLAB中数组下标从1开始hb(1)=1;hb(2)=1;hb(3)=1;hb(4)=1;hb(5)=1;closeall;subplot(3,1,1);m=1:50;T=1;%定义序列的长度是和采样率A=1;a=0.4;T=1;w0=2.0734;x=A*exp(-a*m*T).*sin(w0*m*T);%pi是MATLAB定义的π，信号乘可采用“.*”y=conv(x,hb);subplot(3,1,1);plot(n,abs(fft(hb)))title('h(n)的FFT')subplot(3,1,2);plot(n,abs(fft(x)))title('x(n)的FFT')subplot(3,1,3);plot(abs(fft(y)))title('y(n)的FFT')分析：MATLAB中，计算矢量x的DFT及其逆变换的函数分别为fft和ifft，这两个函数采用了混合算法，当N为质数时，采用的是原始的DFT算法。函数是用机器语言编写的，执行速度混快。N点的FFT调用形式为fft(x,N)。如果x的长度小于N，则补零使其成为N点序列；如果省略N点，即以fft(x)形式调用，则按矢量x的长度进行计算；如果x表示一个矩阵，则调用后计算出每列的N点的FFT。4、由采样序列恢复出连续时间信号，画出其时域波形，对比与原连续时间信号的时域波形，计算并记录两者最大误差。)(nx)(1tx)(1tx)(txA=input('pleaseinputtheA:');%设置信号的有关参数a=input('pleaseinputthea:');W0=input('pleaseinputtheW0:');fs=input('pleaseinputthefs:');%采样率n=0:49;%定义序列的长度T=1/fs;t0=10/a;Dt=1/(5*a);t=0:Dt:t0;xa=A*exp(-a*t).*sin(W0*t);K1=50;k1=0:1:K1;W1max=2*pi*500;W1=W1max*k1/K1;w1=W1/pi;Xa=xa*exp(-j*t'*W1);x=A*exp(-a*n*T).*sin(W0*n*T);figure(1);subplot(4,1,1);plot(t*1000,xa);title('连续时间信号x(t)');axis([0t0*10000200]);grid,xlabel('t:毫秒'),ylabel('x(t)');subplot(4,1,2);plot(w1,abs(Xa));title('连续时间信号频谱Xa(w1)');subplot(4,1,3);stem(x)%绘制x(n)图形grid,xlabel('n'),ylabel('x(n)');title('采样序列x(n)');x1=spline(n*T,x,t);grid,xlabel('t:毫秒'),ylabel('x(t)');subplot(4,1,4);plot(t*1000,x1);axis([0t0*10000200]);title('由x(n)恢复x1(t)');grid,xlabel('t:毫秒'),ylabel('x1(t)');errror=max(abs(x1-xa));k2=-25:25;W2=(pi/12.5)*k2;w2=W2/pi;X=x*(exp(-j*pi/12.5)).^(n'*k2);%序列的傅里叶变换函数figure(2);subplot(2,1,1);plot(w2,abs(X));grid,xlabel('w2'),ylabel('幅度');title('输入信号幅度谱')axis([-2201000]);subplot(2,1,2);plot(w2,angle(X));grid,xlabel('w2'),ylabel('幅角');title('输入信号相位谱');axis([-22-55]);分析：恢复曲线与原信号曲线相同，说明恢复误差很小，如果采样频率减小，误差增大，采样频率增大，则恢复误差更小。采样频率应遵循乃奎斯特定理。五、调试总结：Axis函数设置图形坐标。Conv函数进行卷积运算。Stem用还绘制离散图谱，plot用来绘制连续的函数图形。fft进行fft运算。subplot(3,1,1);plot(n,abs(fft(hb)))title('h(n)的FFT')subplot(3,1,2);plot(n,abs(fft(x)))title('x(n)的FFT')subplot(3,1,3);plot(abs(fft(y)))title('y(n)的FFT')最后一行plot(abs(fft(y)))如果写成plot(n，abs(fft(y)))会出现长度不够，没有图形输出。六、设计总结：在