西电数字信号处理大作业

西电-数字信号处理大作业数字信号处理上机大作业实验一：信号、系统及系统响应(1)简述实验目的及实验原理。1.实验目的熟悉连续信号经理想采样前后的频谱变化关系，加深对时域采样定理的理解。熟悉时域离散系统的时域特性。利用卷积 方法 快递客服问题件处理详细方法山木方法pdf计算方法pdf华与华方法下载八字理论方法下载 观察分析系统的时域特性。掌握序列傅里叶变换的计算机实现方法，利用序列的傅里叶变换对连续信号、离散信号及系统响应进行频域分析。2.实验原理与方法时域采样。LTI系统的输入输出关系。按实验步骤附上实验过程中的信号序列、系统单位脉冲响应及系统响应序列的时域和幅频特性曲线，并对所得结果进行分析和解释。Matlab源程序如下：A=1;T1=1/1000;T2=1/300;T3=1/200;a=25*pi;w0=30*pi;n=0:99;x1=A*exp(-a*n*T1).*sin(w0*n*T1);x2=A*exp(-a*n*T2).*sin(w0*n*T2);x3=A*exp(-a*n*T3).*sin(w0*n*T3);m=linspace(-pi,pi,10000);X1=x1*exp(-j*n'*m);%n'与m构造矩阵，xi向量与矩阵每一列相乘对应元素相加，构成DTFT后的矩阵X2=x2*exp(-j*n'*m);X3=x3*exp(-j*n'*m);figure(1);subplot(3,2,1)plot(m/pi,abs(X1));xlabel('\omega/π');ylabel('|H(e^j^\omega)|');title('采样频率为1000Hz时的幅度谱');subplot(3,2,3)plot(m/pi,abs(X2));xlabel('\omega/π');ylabel('|H(e^j^\omega)|');title('采样频率为300Hz时的幅度谱');subplot(3,2,5)plot(m/pi,abs(X3));xlabel('\omega/π');ylabel('|H(e^j^\omega)|');title('采样频率为200Hz时的幅度谱');subplot(3,2,2)plot(n,abs(x1));xlabel('n');ylabel('x1(t)');title('采样频率为1000Hz时的时域波形');subplot(3,2,4)plot(n,abs(x2));xlabel('n');ylabel('x2(t)');title('采样频率为300Hz时的时域波形');subplot(3,2,6)plot(n,abs(x3));xlabel('n');ylabel('x3(t)');title('采样频率为200Hz时的时域波形');波形图如下：②时域离散信号、系统和系统响应分析。Matlab源程序如下:xb=[1];xc=ones(1,10);ha=ones(1,10);hb=[1,2.5,2.5,1];y=conv(xb,hb);n1=0:length(y)-1;n2=0:length(hb)-1;figure(1)subplot(2,1,1);stem(n1,y,'filled');xlabel('n');ylabel('y(n)');title('y(n)的时域响应');subplot(2,1,2);stem(n1,hb,'filled');xlabel('n');ylabel('hb(n)');title('hb(n)的时域相应');w=linspace(-pi,pi,10000);Y=y*exp(-j*n1'*w);Hb=hb*exp(-j*n2'*w);figure(2)subplot(2,2,1);plot(w/pi,abs(Y));xlabel('\omega/π');ylabel('幅度');title('DTFT[y(n)]的幅度');subplot(2,2,2);plot(w/pi,angle(Y));xlabel('\omega/π');ylabel('相位');title('DTFT[y(n)]的相位');subplot(2,2,3);plot(w/pi,abs(Hb));xlabel('\omega/π');ylabel('幅度');title('DTFT[Hb(n)]的幅度');subplot(2,2,4);plot(w/pi,angle(Hb));xlabel('\omega/π');ylabel('相位');title('DTFT[Hb(n)]的相位');z=conv(xc,ha);n3=0:length(z)-1;Z=z*exp(-j*n3'*w);figure(3);subplot(2,1,1);plot(w/pi,abs(Z));xlabel('\omega/π');ylabel('幅度');title('DTFT[z(n)]的幅度');subplot(2,1,2);plot(w/pi,angle(Z));xlabel('\omega/π');ylabel('相位');title('DTFT[z(n)]的相位');波形图如下：系统ha(n)对信号xc(n)的响应特性：③由 题 快递公司问题件快递公司问题件货款处理关于圆的周长面积重点题型关于解方程组的题及答案关于南海问题 目②中图2即可验证卷积定理。(3)实验中的主要结论。时域采样定理只有当采样频率fs.max大于信号中最高频率fmax的2倍时(fs.max>2fmax)，即满足奈奎斯特定律的时候，采样之后的数字信号才能够不发生混叠，保存原有信号，不失真。卷积定理：时域卷积等于频域相乘。任何函数和单位脉冲函数卷积得到的都是它本身。当所取得N不同时，卷积出来的结果也不同。(4)思考题在分析理想采样序列特性的实验中，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同？它们所对应的模拟频率是否相同？为什么？答：数字频率度量不相同，但他们所对应的模拟频率相同。由w=Ω*Ts得，采样间隔变化时模拟频率对应的数字频率会有相应的变化，故其度量会有所变化。在卷积定理验证的实验中，如果选用不同的频域采样点数M值，例如，选M=10和M=20，分别做序列的傅里叶变换，求得的结果有无差异？答：有差异，所得到的结果点数不同实验二：用FFT作谱分析(1)简述实验目的及实验原理。实验目的进一步加深DFT算法原理和根本性质的理解(因为FFT只是DFT的一种快速算法，所以FFT的运算结果必然满足DFT的根本性质)。熟悉FFT算法原理和FFT子程序的应用。学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应用FFT。(2)结合实验中所得给定典型序列幅频特性曲线，与理论结果比拟，并分析说明误差产生的原因以及用FFT作谱分析时有关参数的选择方法。1)Matlab源程序如下N=64;n=0:999;fs=50;T=1/fs;x1=ones(1,4);x2=[1,2,3,4,4,3,2,1];x3=[4,3,2,1,1,2,3,4];x4=cos(0.25*pi*n);x5=sin(0.125*pi*n);x6=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);X10=fft(x1,N);X20=fft(x2,N);X30=fft(x3,N);X40=fft(x4,N);X50=fft(x5,N);X60=fft(x6,N);X1=fftshift(X10);X2=fftshift(X20);X3=fftshift(X30);X4=fftshift(X40);X5=fftshift(X50);X6=fftshift(X60);k1=0:31;k2=32:63;w1=2*pi*fs*k1/N;w2=(2*pi*k2/N-2*pi)*fs;w=[w2w1];figure(1)subplot(2,1,1);n1=0:length(x1)-1;stem(n1,x1);xlabel('n');ylabel('x1(n)');title('x1〔n〕的时域波形');subplot(2,1,2);stem(w,abs(X1));xlabel('w');ylabel('|X1(k)|');title('DFT(x1)幅频特性')figure(2)subplot(2,1,1);n2=0:length(x2)-1;stem(n2,x2);xlabel('n');ylabel('x2(n)');title('x2〔n〕的时域波形');subplot(2,1,2);stem(w,abs(X2));xlabel('w');ylabel('|X2(k)|');title('DFT(x2)幅频特性')figure(3)subplot(2,1,1);n3=0:length(x3)-1;stem(n3,x3);xlabel('n');ylabel('x3(n)');title('x3〔n〕的时域波形');subplot(2,1,2);stem(w,abs(X3));xlabel('w');ylabel('|X3(k)|');title('DFT(x3)幅频特性')figure(4)stem(w,abs(X4));xlabel('w');ylabel('|X4(k)|');title('DFT(x4)幅频特性')figure(5)stem(w,abs(X5));xlabel('w');ylabel('|X5(k)|');title('DFT(x5)幅频特性')figure(6)stem(w,abs(X6));xlabel('w');ylabel('|X6(k)|');title('DFT(x6)幅频特性')波形图如下：x(n)=x4(n)+x5(n)8点离散傅立叶变换Matlab源程序如下：N=8;n=0:7;x4=cos(0.25*pi*n);x5=sin(0.125*pi*n);xn=x4+x5;X10=fft(xn,N);X1=fftshift(X10);figure(1)stem(abs(X1));axisauto;xlabel('k');ylabel('|Xn(k)|');title('8点DFT[x(n)]幅频特性')波形图如下：16点离散傅立叶变换Matlab源程序如下：M=16;n=0:15;x4=cos(0.25*pi*n);x5=sin(0.125*pi*n);xn=x4+x5;X20=fft(xn,M);X2=fftshift(X20);stem(abs(X2));axisauto;xlabel('k');ylabel('|Xn(k)|');title('16点DFT[x(n)]幅频特性')波形图如下：x(n)=x4(n)+jx5(n)8点离散傅立叶变换Matlab源程序如下：N=8;n=0:7;x4=cos(0.25*pi*n*T);x5=sin(0.125*pi*n*T);xn=x4+j*x5;X10=fft(xn,N);X1=fftshift(X10);figure(1)stem(abs(X1));xlabel('w');ylabel('|Xn(k)|');title('8点DFT[x(n)]幅频特性')波形图如下：16点离散傅立叶变换Matlab源程序如下：M=16;n=0:15;x4=cos(0.25*pi*n);x5=sin(0.125*pi*n);xn=x4+j*x5;X20=fft(xn,M);X2=fftshift(X20);figure(1)stem(abs(X2));axisauto;xlabel('k');ylabel('|Xn(k)|');title('16点DFT[x(n)]幅频特性')波形图如下：(3)总结实验所得主要结论。序列的N点离散傅里叶变换是其z变换在单位圆上的N点等间隔采样，也是其离散时间傅里叶变换在[0,2π]区间上的N点等间隔采样。DFT的变换区间长度N不同，表示对X(z)的单位圆上，或对其离散时间傅里叶变换在[0,2π]区间上的采样间隔和采样点数不同，因而DFT的变换结果不同。信号的DFT变换分为共轭对称局部合共个反对称局部，即信号的实部对应离散傅里叶变换的共轭对称局部，虚部对应离散信号的共轭反对称局部。(4)思考题在N=8时，x2(n)和x3(n)的幅频特性会相同吗?为什么?N=16呢?答：N=8时幅频特性一样，N=16时幅频特性不一样。如果周期信号的周期预先不知道，如何用FFT进行谱分析?答：设一个定长的m值,先取2m,看2m/m的误差是否大，如大的话再取4m,看4m/2m的误差是否大，如不大，4m(4倍的m值)那么可近似原来点的谱分析。实验三：用窗函数法设计FIR数字滤波器(1)简述实验目的及实验原理。实验目的(1)掌握用窗函数法设计FIR数字滤波器的原理和方法。(2)熟悉线性相位FIR数字滤波器特性。(3)了解各种窗函数对滤波特性的影响。(2)实验结果矩型窗、升余弦窗、改良升余弦窗和二阶升余弦窗的窗函数及频谱特Matlab源程序如下：N=31;wdrec=boxcar(N);[hrec,w]=freqz(wdrec,1);figure(1)subplot(1,2,1)stem(wdrec);xlabel('n');ylabel('\omega(n)');title('矩形窗函数')subplot(1,2,2);plot(w/pi,10*log(abs(hrec)/abs(hrec(1))));xlabel('ω/Π');ylabel('W(\omega)/dB')title('矩形窗函数幅度特性')wdhan=hanning(N);[hhan,w]=freqz(wdhan,1);figure(2)subplot(1,2,1)stem(wdhan);xlabel('n');ylabel('\omega(n)');title('汉宁窗函数')subplot(1,2,2);plot(w/pi,10*log(abs(hhan)/abs(hhan(1))));xlabel('ω/Π');ylabel('W(\omega)/dB')title('汉宁窗函数幅度特性')wdham=hamming(N);[hham,w]=freqz(wdham,1);figure(3)subplot(1,2,1)stem(wdham);xlabel('n');ylabel('\omega(n)');title('汉明窗函数')subplot(1,2,2);plot(w/pi,10*log(abs(hham)/abs(hham(1))));xlabel('ω/Π');ylabel('W(\omega)/dB')title('汉明窗函数幅度特性')wdblack=blackman(N);[hblack,w]=freqz(wdblack,1);figure(4)subplot(1,2,1)stem(wdblack);xlabel('n');ylabel('\omega(n)');title('布莱克曼窗函数')subplot(1,2,2);plot(w/pi,10*log(abs(hblack)/abs(hblack(1))));xlabel('ω/Π');ylabel('W(\omega)/dB')title('布莱克曼窗函数幅度特性')波形图如下：利用不同窗函数与hd(n)构成系统，其单位脉冲相应及频谱特性Matlab源程序如下：N=31;a=(N-1)/2;fp=120;fs=150;Fs=1000;%以上数据自由选取wc=(fp+fs)/Fs;h1=fir1(N-1,wc,boxcar(N));h2=fir1(N-1,wc,hanning(N));h3=fir1(N-1,wc,hamming(N));h4=fir1(N-1,wc,blackman(N));[H1,w]=freqz(h1,1);[H2,w]=freqz(h2,1);[H3,w]=freqz(h3,1);[H4,w]=freqz(h4,1);figure(1)subplot(1,2,1)stem(h1);xlabel('n');ylabel('h1(n)');title('矩形窗数字滤波器单位脉冲响应')subplot(1,2,2);plot(w/pi,10*log(abs(H1)/abs(H1(1))));xlabel('ω/Π');ylabel('W(\omega)/dB')title('矩形窗数字滤波器幅频特性')figure(2)subplot(1,2,1)stem(h2);xlabel('n');ylabel('h2(n)');title('汉宁窗数字滤波器单位脉冲响应')subplot(1,2,2);plot(w/pi,10*log(abs(H2)/abs(H2(1))));xlabel('ω/Π');ylabel('W(\omega)/dB')title('汉宁窗数字滤波器幅频特性')figure(3)subplot(1,2,1)stem(h3);xlabel('n');ylabel('h3(n)');title('汉明窗数字滤波器单位脉冲响应')subplot(1,2,2);plot(w/pi,10*log(abs(H3)/abs(H3(1))));xlabel('ω/Π');ylabel('W(\omega)/dB')title('汉明窗数字滤波器幅频特性')figure(4)subplot(1,2,1)stem(h4);xlabel('n');ylabel('h4(n)');title('布莱克曼窗数字滤波器单位脉冲响应')subplot(1,2,2);plot(w/pi,10*log(abs(H4)/abs(H4(1))));xlabel('ω/Π');ylabel('W(\omega)/dB')title('布莱克曼窗数字滤波器幅频特性')波形图如下：(3)总结用窗函数法设计FIR滤波器的主要特点。FIR滤波器目前常用的设计方法有窗函数法和频率采样法,窗函数法是从时域进行设计,而频率采样法是从频域进行设计.窗函数法由于简单、物理意义清晰,因而得到了较为广泛的应用.窗函数法设计的根本思想是：首先根据技术指标要求,选取适宜的阶数N和窗函数的类型w(n),使其幅频特性逼近理想滤波器幅频特性.其次,因为理想滤波器的hd(n)是无限长的,所以需要对hd(n)进行截断,数学上称这种方法为窗函数法.简而言之,用窗函数法设计FIR滤波器是在时域进行的,先用傅里叶变换求出理想滤波器单位抽样相应hd(n),然后加时间窗w(n)对其进行截断,以求得FIR滤波器的单位抽样响应h(n).(3)思考题(1)如果给定通带截止频率和阻带截止频率以及阻带最小衰减，如何用窗函数法设计线性相位低通滤波器?写出设计步骤。答：①根据所需设计的数字滤波器类型，确定线性相位数字滤波器类型，②选择适宜的窗函数；③确定理想低通数字滤波器的频率响应函数；④计算理想低通数字滤波器的单位脉冲响应；⑤加窗得到设计结果；如果要求用窗函数法设计带通滤波器，且给定上、下边带截止频率为ω1和ω2，试求理想带通的单位脉冲响应hd(n)。答：hd(n)=1/2π∫w1w0e-jwaejnwdw+1/2π∫w2w3e-jwaejnwdw〔其中w0=-w0-wc,w1=-w0+wc,w2=w0-wc,w3=w0+wc〕计算整理后可得：hd(n)=2/((n-a)*π)*sin[(n-a)wc]*cos[(n-a)w0]=2wc/π*sa[(n-a)wc]*cos[(n-a)w0]