数字信号处理西电

数字信号处理上机第一次实验

实验一：

设给定模拟信号()1000t

a x t e

-=，的单位是ms 。

(1) 利用MATLAB 绘制出其时域波形和频谱图（傅里叶变换），估计其等效带宽（忽略谱分

量降低到峰值的3%以下的频谱）。 (2) 用两个不同的采样频率对给定的进行采样。

○1

。

○2 。

比较两种采样率下的信号频谱，并解释。 实验一MATLAB 程序： （1）

○

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',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

t ()a x t ()()15000s a f x t x n =以样本秒采样得到。()()11j x n X e ω画出及其频谱()()1

1000s a f x t x n =以样本得到。()()

11j x n X e ω画出及其频谱

○2 clc;

fs=1000;

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',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

实验三：

设，，编写MATLAB 程序，计算： （１） ５点圆周卷积； （２） ６点圆周卷积； （３） 线性卷积；

（４） 画出的，和时间轴对齐。 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 on

title('五点圆周卷积y1(n)');

xlabel('n'),ylabel('y1(n)');axis([0 6 0 15]) subplot(312)

(){}11,2,2x n =(){}21,2,3,4x n =()1y n ()2y n ()3y n ()1y n ()2y n ()3y n

grid on

title('六点圆周卷积y2(n)');

xlabel('n'),ylabel('y2(n)');axis([0 6 0 15])

subplot(313)

stem(y3);

grid on

title('线性卷积y3(n)');

xlabel('n'),ylabel('y3(n)');axis([0 6 0 15])

给定因果系统： （１） 求系统函数并画出零极点示意图。 （２） 画出系统的幅频特性和相频特性。

（３） 求脉冲响应并画序列图。

提示：在ＭＡＴＬＡＢ中，zplane(b,a) 函数可画零极点图；Freqz(b,a,N)可给出范围内均匀间隔的点频率响应的复振幅；Impz(b,a,N)可求的逆变换（即脉冲响应）。

clc

a = [1,0]

b = [1,-0.9] figure(1) zplane(b,a);

title('零极点分布图') w=[-3*pi:0.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)');

()()()0.91y n y n x n =-+()H z ()j H e

ω

()?ω()h n []

0,πN ()H z

数字信号处理第二次实验

1. 给定模拟信号()()()2sin 45cos 8x t t t ππ=+，对其进行采样，用DFT （FFT ）进行信号频谱分析。

（１） 确定最小采样频率和最小采样点数。

（２） 若以()0.010:1t n n N ==-秒进行采样，至少需要取多少采样点？

（３） 用DFT 的点数50,100N =画出信号的N 点DFT 的幅度谱，讨论幅度谱结果。 （４） N 分别为64N =和60N =，能否分辨出信号的所有频率分量。 （５） 在（３）和（４）的条件下做补0 FFT ，分析结果。 （６） 在不满足最小采样点数的情况下做补0DFT ，观察是否可以分辨出两个频率分量。

程序如下： clear

close all clc

%(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('以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*n4*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,'.');

title('N=60幅度谱')

xlabel('n')

ylabel('X4')

%(5)补0DFT

N5=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,'.');

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. 设雷达发射线性调频信号()()2exp 2h t j t πμ=，13510μ=?，采样率9210s f =?，采样点数

20000N =。回波信号()()()12s t h t h t ττ=-+-，6110τ-=，62 1.110τ-=?。

（１） 画出()h t 的频谱。

（２） 利用DFT 的时延性质产生()s t ，比较直接在时域产生和在频域产生（再变换到时域）

的结果是否相同。

（３） 匹配滤波的结果是()()()y t s t h t *=*-，（“*”表示线性卷积）。分别用直接线性卷积

和DFT 的卷积定理求解()y t 。比较二者结果，并记录两种方法的运行时间（用tic ，toc 指令）。

（４） 画出()y t 的频谱。

程序如下： 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)频谱') 第二问

时域构造回波信号，时延通过补零实现 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)))))%将回波信号转换到频域并将零频搬至频谱中央

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 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('时域法与频域法回波之差')

第三问

tic

y_shiyu=conv(s_pinyu,conj(fliplr(ht)));%%%时域匹配滤波

toc

y_shiyu_quchu=y_shiyu(1,N:end);%%%%去暂态点

figure;

plot([0:2*N-1],abs(y_shiyu_quchu))

xlabel('距离单元')

ylabel('幅度')

title('时域匹配滤波')

tic

y_pinyu=ifft(fft(s_pinyu).*conj(P));频域匹配滤波

toc

figure;

plot([0:2*N-1],abs(y_pinyu))

xlabel('距离单元')

ylabel('幅度')

title('频域匹配滤波')

第四问

figure

plot([-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）用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带阻滤波器所有四种类型的最低阶数。指标如下：12kHz的采样率，2.1kHz和4.5kHz

的通带边界频率，2.7kHz和3.9kHz的阻带边界频率，

0.6dB的通带波纹，45dB的最小阻带衰减。并在同一张

图中画出每种滤波器的频率响应。

用到的函数：butter，buttord，cheb2ord，chebl1，cheby2，