数字信号处理西电.docx

1、数字信号处理西电数字信号处理西电数字信号处理上机第一次实验实验一：设给定模拟信号，的单位是ms。(1) 利用MATLAB绘制出其时域波形和频谱图（傅里叶变换），估计其等效带宽（忽略谱分量降低到峰值的3%以下的频谱）。(2) 用两个不同的采样频率对给定的进行采样。比较两种采样率下的信号频谱，并解释。实验一MATLAB程序：（1） clc; fs=5000; ts=1/fs; N=1000; t=(-N:N)*ts; s=exp(-abs(t); plot(t,s,linewidth,1.5) xlabel(时间) ylabel(幅度) set(gca,fontweight,b,fontsize,

2、12) SPL=N*100; figure sp=fftshift(fft(s,SPL); sp=sp/max(sp)*100; freqb=-fs/2:fs/SPL:fs/2-fs/SPL; plot(freqb,abs(sp) xlabel(频率) ylabel(频谱幅度)set(gca,fontweight,b,fontsize,12)yy=abs(abs(sp)-3);aa,freqind=min(yy);(freqind-SPL/2)*fs/SPLclc; fs=1000; ts=1/fs; N=1000; t=(-N:N)*ts; s=exp(-abs(t); plot(t,s,l

3、inewidth,1.5) xlabel(时间) ylabel(幅度) set(gca,fontweight,b,fontsize,12) SPL=N*100; figure sp=fftshift(fft(s,SPL); sp=sp/max(sp)*100; freqb=-fs/2:fs/SPL:fs/2-fs/SPL; plot(freqb,abs(sp) xlabel(频率) ylabel(频谱幅度)set(gca,fontweight,b,fontsize,12)yy=abs(abs(sp)-3);aa,freqind=min(yy);(freqind-SPL/2)*fs/SPL实验三

4、：设，编写MATLAB程序，计算：（1） 点圆周卷积；（2） 点圆周卷积；（3） 线性卷积；（4） 画出的，和时间轴对齐。a = 1,2,2;b = 1,2,3,4;y1 = cconv(a,b,5)y2 = cconv(a,b,6)y3 = conv(a,b)figure(1);subplot(311)stem(y1);grid ontitle(五点圆周卷积y1(n);xlabel(n),ylabel(y1(n);axis(0 6 0 15)subplot(312)stem(y2);grid ontitle(六点圆周卷积y2(n);xlabel(n),ylabel(y2(n);axis(0

5、6 0 15)subplot(313)stem(y3);grid ontitle(线性卷积y3(n);xlabel(n),ylabel(y3(n);axis(0 6 0 15)实验四：给定因果系统：（1） 求系统函数并画出零极点示意图。（2） 画出系统的幅频特性和相频特性。（3） 求脉冲响应并画序列图。提示：在中，zplane(b,a) 函数可画零极点图；Freqz(b,a,N)可给出范围内均匀间隔的点频率响应的复振幅；Impz(b,a,N)可求的逆变换（即脉冲响应）。clca = 1,0b = 1,-0.9figure(1)zplane(b,a);title(零极点分布图)w=-3*pi:0

6、.01:3*pi;h,phi=freqz(b,a,w);figure(2);subplot(2,1,1);plot(w, abs(h);grid on;title(幅频特性);xlabel(f/Hz),ylabel(H(f); subplot(2,1,2);plot(w, phi);grid on;title(相频特性);xlabel(f/Hz),ylabel(W(f);数字信号处理第二次实验1. 给定模拟信号，对其进行采样，用DFT（FFT）进行信号频谱分析。（1） 确定最小采样频率和最小采样点数。（2） 若以秒进行采样，至少需要取多少采样点？（3） 用DFT的点数画出信号的点DFT的幅度谱

7、，讨论幅度谱结果。（4） 分别为和，能否分辨出信号的所有频率分量。（5） 在（）和（）的条件下做补0 FFT，分析结果。（6） 在不满足最小采样点数的情况下做补0DFT，观察是否可以分辨出两个频率分量。程序如下：clearclose allclc%(1)确定最小采样频率和最小采样点数w1=4*pi;w2=8*pi;f1=w1/(2*pi);f2=w2/(2*pi);disp(最小采样频率:)fs1=2*max(f1,f2);disp(fs1);f=f2-f1;disp(最小采样点数:)N=ceil(fs1/f);disp(N);%（2）t=0.01ns采样T=0.01;fs2=1/T;disp

8、(以t=0.01ns采样,最少采样点数为：)N0=fs2/f;disp(N0);%（3）（4）N=50,100,64,60时的幅度谱w1=4*pi;w2=8*pi;f1=w1/(2*pi);f2=w2/(2*pi);N1=50;N2=100;N3=64;N4=60;n1=0:N1-1;n2=0:N2-1;n3=0:N3-1;n4=0:N4-1;x1=2*cos(w1*n1*T)+5*cos(w2*n1*T);x2=2*cos(w1*n2*T)+5*cos(w2*n2*T);x3=2*cos(w1*n3*T)+5*cos(w2*n3*T);x4=2*cos(w1*n4*T)+5*cos(w2*n

9、4*T);X1=abs(fft(x1,N1);X2=abs(fft(x2,N2);X3=abs(fft(x3,N3);X4=abs(fft(x4,N4);figure(1)subplot(2,2,1);stem(n1,X1,.);title(N=50幅度谱)xlabel(n)ylabel(X1)subplot(2,2,2);stem(n2,X2,.);title(N=100幅度谱)xlabel(n)ylabel(X2)subplot(2,2,3);stem(n3,X3,.);title(N=64幅度谱)xlabel(n)ylabel(X3)subplot(2,2,4);stem(n4,X4,.

10、);title(N=60幅度谱)xlabel(n)ylabel(X4)%(5)补0DFTN5=200;n5=0:N5-1;X5=abs(fft(x1,N5);X6=abs(fft(x2,N5);X7=abs(fft(x3,N5);X8=abs(fft(x4,N5);figure(2)subplot(2,2,1);stem(n5,X5,.);title(补0后N=50幅度谱)xlabel(n)ylabel(X5)subplot(2,2,2);stem(n5,X6,.);title(补0后N=100幅度谱)xlabel(n)ylabel(X6)subplot(2,2,3);stem(n5,X7,.

11、);title(补0后N=64幅度谱)xlabel(n)ylabel(X7)subplot(2,2,4);stem(n5,X8,.);title(补0后N=60幅度谱)xlabel(n)ylabel(X8)%（6）N=2时不满足最小采样点数N6=2;n9=0:N6-1;x9=2*cos(w1*n9*T)+5*cos(w2*n9*T);X9=abs(fft(x9,N5);figure(3)stem(n5,X9,.);xlabel(n)ylabel(X9)title(N=2时补0后的幅度谱)运行结果：2. 设雷达发射线性调频信号，采样率，采样点数。回波信号，。（1） 画出的频谱。（2） 利用DFT

12、的时延性质产生，比较直接在时域产生和在频域产生（再变换到时域）的结果是否相同。（3） 匹配滤波的结果是，（“”表示线性卷积）。分别用直接线性卷积和DFT的卷积定理求解。比较二者结果，并记录两种方法的运行时间（用tic，toc指令）。（4） 画出的频谱。程序如下：figure;plot(-0.5:1/(fft_num):0.5-1/(fft_num),. fftshift(20*log10(abs(fft(ht,fft_num)%将线性调频信号转换到频域并将零频搬至频谱中央axis(-0.5 0.5 10 50)xlabel(归一化频率)ylabel(幅度/dB)title(h(t)频谱)第二问

13、时域构造回波信号，时延通过补零实现s_shiyu=zeros(1,shiyan1*fs),ht,zeros(1,N-shiyan1*fs)+zeros(1,shiyan2*fs),ht,zeros(1,N-shiyan2*fs);figure;plot(0:2*N-1,abs(s_shiyu)axis(0 2*N-1 0 2.5)xlabel(距离单元)ylabel(幅度)title(时域法s(t)figure;plot(-0.5:1/(fft_num):0.5-1/(fft_num),. fftshift(20*log10(abs(fft(s_shiyu,fft_num)%将回波信号转换到频

14、域并将零频搬至频谱中央axis(-0.5 0.5 10 60)xlabel(归一化频率)ylabel(幅度/dB)title(时域法s(t)频域)频域构造回波信号，时延通过DFT时延性质产生L=2*N;P=fft(ht,L);P_shiyan1=P.*exp(-j*2*pi*fs*0:L-1/L*shiyan1);%目标1频谱P_shiyan2=P.*exp(-j*2*pi*fs*0:L-1/L*shiyan2);%目标2频谱s_pinyu=ifft(P_shiyan1)+ifft(P_shiyan2);figure;plot(0:2*N-1,abs(s_pinyu)axis(0 2*N-1

15、0 2.5)xlabel(距离单元)ylabel(幅度)title(频域法s(t)figure;plot(-0.5:1/(fft_num):0.5-1/(fft_num),. fftshift(20*log10(abs(fft(s_pinyu,fft_num)axis(-0.5 0.5 10 50)axis(-0.5 0.5 10 60)xlabel(归一化频率)ylabel(幅度/dB)title(频域法s(t)频域)figure;plot(0:2*N-1,abs(s_shiyu-s_pinyu)axis(0 2*N-1 0 1)xlabel(距离单元)ylabel(幅度)title(时域法

16、与频域法回波之差)第三问 ticy_shiyu=conv(s_pinyu,conj(fliplr(ht);%时域匹配滤波tocy_shiyu_quchu=y_shiyu(1,N:end);%去暂态点figure;plot(0:2*N-1,abs(y_shiyu_quchu)xlabel(距离单元)ylabel(幅度)title(时域匹配滤波)ticy_pinyu=ifft(fft(s_pinyu).*conj(P);频域匹配滤波tocfigure;plot(0:2*N-1,abs(y_pinyu)xlabel(距离单元)ylabel(幅度)title(频域匹配滤波)第四问 figureplot

17、(-0.5:1/L:0.5-1/L,fftshift(20*log10(abs(fft(y_pinyu)xlabel(归一化频率)ylabel(幅度/dB)title(y(t)频域)运行结果：Elapsed time is 0.235671 seconds.Elapsed time is 0.001764 seconds. 数字信号处理上机第三次实验1. IR滤波器设计（1） 用matlab确定一个数字IIR低通滤波器所有四种类型的最低阶数。指标如下：40kHz的采样率，4kHz的通带边界频率，8kHz的阻带边界频率，0.5dB的通带波纹，40dB的最小阻带衰减。并在同一张图中画出每种w。（2

18、） 用matlab确定一个数字IIR高通滤波器所有四种类型的最低阶数。指标如下：3500Hz的采样率，1050Hz的通带边界频率，600Hz的阻带边界频率，1dB的通带波纹，50dB的最小阻带衰减。并在同一张图中画出每种滤波器的频率响应。（3） 用matlab确定一个数字IIR带通滤波器所有四种类型的最低阶数。指标如下：7kHz的采样率，1.4kHz和2.1kHz的通带边界频率，1.05kHz和2.45kHz的阻带边界频率，0.4dB的通带波纹，50dB的最小阻带衰减。并在同一张图中画出每种滤波器的频率响应。（4） 用matlab确定一个数字IIR带阻滤波器所有四种类型的最低阶数。指标如下：1

19、2kHz的采样率，2.1kHz和4.5kHz的通带边界频率，2.7kHz和3.9kHz的阻带边界频率，0.6dB的通带波纹，45dB的最小阻带衰减。并在同一张图中画出每种滤波器的频率响应。用到的函数： butter，buttord，cheb2ord，chebl1，cheby2，ellip，ellipord.程序如下：(1) clcclear allclose allfc=40;fp=4;fs=8;rp=0.5;rs=40;wp=2*pi*fp/fc;ws=2*pi*fs/fc;disp(In buttord)n,wc=buttord(wp,ws,rp,rs,s)b,a=butter(n,wc,

20、low,s);w=0:0.001:6;h,w=freqs(b,a,w);h=20*log10(abs(h);plot(w,h,r-)disp(In cheb1ord)n,wpo=cheb1ord(wp,ws,rp,rs,s)b,a=cheby1(n,rp,wpo,low,s);h,w=freqs(b,a,w);h=20*log10(abs(h);hold onplot(w,h,b-)disp(In cheb2ord)n,wso=cheb2ord(wp,ws,rp,rs,s)b,a=cheby2(n,rs,wso,low,s);h,w=freqs(b,a,w);h=20*log10(abs(h)

21、;hold onplot(w,h,k-)disp(In ellipord)n,wc=ellipord(wp,ws,rp,rs,s)b,a=ellip(n,rp,rs,wc,low,s);h,w=freqs(b,a,w);h=20*log10(abs(h);hold onplot(w,h,m-)title(滤波器的频率响应)legend(巴特沃斯,切比雪夫型,切比雪夫型,椭圆)grid onxlabel(w)ylabel(h)运行结果：巴特沃斯n = 9wc = 0.7533切比雪夫型n = 5wpo = 0.6283切比雪夫型n = 5wso = 1.2069椭圆n = 4wc = 0.628

22、3(2) clcclear allclose allfc=3500;fp=1050;fs=600;rp=1;rs=50;wp=2*pi*fp/fc;ws=2*pi*fs/fc;disp(In buttord)n,wc=buttord(wp,ws,rp,rs,s)b,a=butter(n,wc,high,s);w=0:0.001:6;h,w=freqs(b,a,w);h=20*log10(abs(h);plot(w,h,r-)disp(In cheb1ord)n,wpo=cheb1ord(wp,ws,rp,rs,s)b,a=cheby1(n,rp,wpo,high,s);h,w=freqs(b,

23、a,w);h=20*log10(abs(h);hold onplot(w,h,b-)disp(In cheb2ord)n,wso=cheb2ord(wp,ws,rp,rs,s)b,a=cheby2(n,rs,wso,high,s);h,w=freqs(b,a,w);h=20*log10(abs(h);hold onplot(w,h,k-)disp(In ellipord)n,wc=ellipord(wp,ws,rp,rs,s)b,a=ellip(n,rp,rs,wc,high,s);h,w=freqs(b,a,w);h=20*log10(abs(h);hold onplot(w,h,m-)ti

24、tle(滤波器的频率响应)legend(巴特沃斯,切比雪夫型,切比雪夫型,椭圆)axis(0,6,-100,0)grid onxlabel(w)ylabel(h)运行结果：巴特沃斯n = 12wc = 1.7402切比雪夫型n = 7wpo = 1.8850切比雪夫型n = 7wso = 1.2049椭圆n = 5wc = 1.8850（3）clcclear allclose allfc=7;fp=1.4,2.1;fs=1.05,2.45;rp=0.4;rs=50;wp=2*pi*fp/fc;ws=2*pi*fs/fc;disp(In buttord)n,wc=buttord(wp,ws,rp

25、,rs,s)b,a=butter(n,wc ,s);w=0:0.001:6;h,w=freqs(b,a,w);h=20*log10(abs(h);plot(w,h,r-)disp(In cheb1ord)n,wpo=cheb1ord(wp,ws,rp,rs,s)b,a=cheby1(n,rp,wpo ,s);h,w=freqs(b,a,w);h=20*log10(abs(h);hold onplot(w,h,b-)disp(In cheb2ord)n,wso=cheb2ord(wp,ws,rp,rs,s)b,a=cheby2(n,rs,wso ,s);h,w=freqs(b,a,w);h=20

26、*log10(abs(h);hold onplot(w,h,k-)disp(In ellipord)n,wc=ellipord(wp,ws,rp,rs,s)b,a=ellip(n,rp,rs,wc ,s);h,w=freqs(b,a,w);h=20*log10(abs(h);hold onplot(w,h,m-)title(滤波器的频率响应)legend(巴特沃斯,切比雪夫型,切比雪夫型,椭圆)axis(0,6,-100,20)grid onxlabel(w)ylabel(h)运行结果：巴特沃斯n = 12wc = 1.2305 1.9250切比雪夫型n = 7wpo = 1.2566 1.8

27、850切比雪夫型n = 7wso = 1.1050 2.1437椭圆n = 5wc = 1.2566 1.8850（4）clcclear allclose allfc=12;fp=2.1,4.5;fs=2.7,3.9;rp=0.5;rs=40;wp=2*pi*fp/fc;ws=2*pi*fs/fc;disp(In buttord)n,wc=buttord(wp,ws,rp,rs,s)b,a=butter(n,wc,stop,s);w=0:0.001:6;h,w=freqs(b,a,w);h=20*log10(abs(h);plot(w,h,r-)disp(In cheb1ord)n,wpo=cheb1ord(wp,ws,rp,rs,s)b,a=cheby1(n,rp,wpo,stop,s);h,w=freqs(b,a,w);h=20*log10(abs(h);