数字信号处理
第一章 1.1、产生周期为4的三角脉冲波,其中t是时间数组,T是周期,利用MatLab实现波形,并进行作图 解 >> T=16;t=-2.5*T:0.01:2.5*T; x=sawtooth(t-2*T/4,0.5); plot(t,x) title('三角脉冲波T=4') axis([-20,20,-2,2]) 图 1.2 1.2、单位脉冲序列 在MatLab 中,这一函数可以用zeros 函数实现、并进行作图: 解: >> n=1:50; x=zeros(1,50); x(1)=1; close all; subplot(3,1,1);stem(x);title('张鹏飞制作单位冲击信号序列'); k=-25:25; X=x*(exp(-j*pi/12.5)).^(n'*k); magX=abs(X); subplot(3,1,2);stem(magX);title('张鹏飞制作单位冲击信号的幅度谱'); angX=angle(X); subplot(3,1,3);stem(angX) ; title ('张鹏飞制作单位冲击信号的相位谱') 图 1.3
1.3、产生阶跃序列 解 >> T=5;t=-2*T:3*T; x=1/2+1/2.*sign(t); stem(t,x) title('张鹏飞制作单位阶跃序列') axis([-10,15,-0.2,1.5])
图 1.4 1.4、用DFT计算x(t)=cos(2∏(0.4)t)+ cos(2∏(0.45)t)的频谱 解: >> T0=20;N=32;T=T0/N;A=2*pi*(0.4);B=2*pi*(0.45); t=0:T:T0; x=cos(B*t)+cos(C*t); X=1/N*fft(x,N); f=1/T/N*(-N/2:(N/2-1));%N,f=1/T/N*(-(N-1)/2:(N-1)/2); subplot(2,1,1);stem(f,abs(fftshift(X))); xlabel('Frequency(Hz)');ylabel('magnitude'); subplot(2,1,2);stem(f,angle(fftshift(X))); xlabel('Frequency(Hz)');ylabel('phase'); 图 1.4 设计理论分析方法 连续时间信号的时域分析就是将不同形式的信号波形用不同的时间函数来描述。 任一信号可以在时域对其进行分析和描述,利用傅立叶变换理论也可以对其进行频域分析,以便更好地对信号进行存储,传输和处理,达到提取有用信号的目的。利用傅立叶变换的分析方法,找出与信号时域信号波形相对应的频谱函数的幅度、相位以及能量或功率的分布规律。 设计结论及其分析 如上述例题,连续时间信号通过CTFS分解为不同谐波成分的三角函数或指数函数之和,连续非周期信号分解为频率无限密集的虚指数函数的线性组合,在时域还可以把实信号分解为奇函数与偶函数之和。 DFT变换对在时域、频域都是离散的,可以通过计算机实现数值运算。必须对连续信号进行离散化,进行频谱分析时,参数的选择主要满足时域抽样定理、频域采样定理,依据时域、频域能量分布情况,恰当地选择实域取样间隔T和取样长度L, 设计技术参数 CTFS ICTFS DFT DFS DTFT FFT 和MATLAB各函数参数
思考题: 1、既然可以直接计算CTFT,为什么利用DFT分析连续信号频谱? 答:DFT是对有限长FT的有限点取样,开辟了频域离散化的道路,使数字信号在频域能进行数字运算。 2、若信号持续时间无限,且无解析表达式,如何利用DFT分析其频谱? 答:先对连续信号进行采样,得到采样信号,然后对采样信号进行频谱分析 3、在利用DFT分析连续信号频谱时,会出现哪些误差?如何克服或减少? 答:(1)混叠。在采样前预滤波,滤除高于折叠频率。 (2)栅栏效应。在序列尾部补0。 (3)截断效应 4、如何选择取样频率? 答:为不产生频率混叠失真,应采样频率大于2倍的最高截止频率,即Fs>2Fh. 5、讨论补0对计算结果的影响? 答:改变序列长度N,增加频率采样点数和采样点位置,使原来漏掉的某些频谱分量被检测出来。
第二章 (一) IIR数字滤波器的设计 1.1 已知:设计线性相位带通滤波器,其长度N=15,上下边带截止频率分别为W1= 0.3 π,w2=0.5 π 设计指标为:ωp=0.2 π Rp=0.25dB ωa=0.3 π As=50dB IIR 数字滤波器 解:
>> Window=hamming(16); b=fir1(15,[0.3 0.5],Window); freqz(b,1) wp=0.2*pi;ws=0.3*pi; tr_width=ws-wp; M=ceil(6.6*pi/tr_width)+1; N=0:1:M-1; wc=(ws+wp)/2; hd=ideal_lp(wc,M); w_ham=(boxcar(M))'; h=hd.*w_ham; [db,mag,pha,grd,w]=freqz_m(h,1); delta_w=2*pi/1000; Rp=-(min(db(1:1:wp/delta_w+1))); As=-round(max(db(ws/delta_w+1:1:501))); Close all; subplot(2,2,1);stem(hd);title('理想冲击响应') axis([0 M-1 -0.1 0.3]);ylabel('hd[n]'); subplot(2,2,2);stem(w_ham);title('汉明'); axis([0 M-1 0 1.1]);ylabel('w[n]'); subplot(2,2,3);stem(h);title('实际冲击响应'); axis([0 M-1 -0.1 0.3]);ylabel('h[n]'); subplot(2,2,4);plot(w/pi,db); title('衰减响应'); axis([0 1 -100 10]);ylabel('Decibles');
图 4.1 1.2 现有一个信号:x(n)=1+cos(π*n/4)+ cos(2*π*n/3)设计: IIR低通滤波器,滤除cos(2*π*n/3) 的成分,即想保留的成分为1+cos(π*n/4) 解: >> omegap=pi/4;Apmax=3; omegar=2/3*pi;Armin=50; T=5*pi*10^(-6); wp=[omegap/pi];wr=[omegar/pi]; [N,wc]=buttord(wp,wr,Apmax,Armin,'s'); [b,a]=butter(N,wc,'low'); sysl=tf(b,a) w=linspace(0,2*pi,500); h=freqz(b,a,w); length(h) mag=fftshift(20*log10(abs(h))); n=0:23; x=1+cos(pi/4*n)+cos(2*pi/3*n); X=fft(x,500); X1=fftshift(X); omega=2*pi/500*([0:499]); Y=fftshift(X.*h); title('原信号'); ylabel('幅度');xlabel('频率'); subplot(3,1,1);plot(omega,abs(X1));axis([0,2*pi,0,35]); title('滤波器幅频响应'); ylabel('Gain(dB)');xlabel('Frequency(rad)'); subplot(3,1,2);plot(w,mag);axis([0,2*pi,-100,5]); subplot(3,1,3); plot(omega,abs(Y));axis([0,2*pi,0,25]); 图 4.2 设计理论分析方法 设计IIR数字滤波器的方法有间接法、直接法和计算机辅助法等。设计IIR数字滤波器的任务就是寻求一个因果、物理可实现的系统函数H(z),使它的频响H(ejw)满足所希望得到的低通频域指标,即通带衰减A p,阻带衰减A r,通带截频Wp,阻带截频Wr。采用间接法设计IIR数字滤波器就是给定的指标,先设计一个模拟滤波器,进而通过模拟域与数字域的变换,求得物理可实现的滤波器。有脉冲响应不变法和双线性变换法。 这里常用的是根据频响指标,设计出相应的模拟滤波器H(s),再通过脉冲响应不变法或双线性法转换成数字滤波器H(z)。其过程如图: 双线性不变法 频域指标 设计模拟滤波器 数字滤波器 脉冲响应不变法
设计结论及其分析 脉冲响应不变法由于保持了模拟滤波器的冲击响应在取样时刻的取样值,所以具有较好的时间性。但由于从模拟域变换到数字域的映射关系不是一一对应的关系,因而频谱混叠难以避免,所以实际应用中有它的局限性。 双线性变换法克服了脉冲响应不变法存在频谱混叠的缺点,但在设计过程中一定要注意,由于宽频带压缩到较窄的频率范围,导致了频率畸变。 设计技术参数 A p A r Wp Wr Buttord Butterworth
(二) FIR数字滤波器的设计 4.1、 利用窗函数法设计FIR数字滤波器,指标为: ej2w -∏/4<=w<=∏/4 H(ejw)={ 0 ∏/4<=w<=∏ 解: >> M=600;n=0:M;omegac=2; hd=omegac/pi.*sinc(omegac/pi*(n-M/2)); stem([0:M],hd); h1=hd; H1=fft(h1,512); H1_db=20*log10(abs(H1)); w=0.54-0.46.*cos(2*pi*n./M); stem([0:M],w); h2=hd.*w; H2=fft(h2,512); H2_db=20*log10(abs(H2)); omega=0:pi/255:pi; c=plot(omega,H1_db(1:256),omega,H2_db(1:256),'r'); title('FIR数字低通滤波器') 图 4.3
设计理论分析方法 FTR滤波器具有线性相位的特点,设计核心是求出有限的脉冲响应来逼近给定的频响,常用的设计方法有窗函数发和频率取样法。 设计结论及其分析 通过把理想滤波器的无限长脉冲响应乘以窗函数来产生一个被截短的脉冲响应。 设计技术参数 矩形窗 汉宁窗 汉明窗 思考题: 1、 做出该信号通过滤波器的响应并与做比较两者一致吗?原因是什么? 答:不一致。因为在数字滤波器设计过程中要经过频域变换,若频域响应不变会出现频域混叠,若用双线性法会出现频率关系非线性,都会使设计的滤波器与理想的滤波器存在误差,同时设计过程中取得的阶次,所满足的技术指标的程度都有影响。 2、 IIR滤波器的优缺点是什么? 答:无法实现严格的线性相位,且其反馈有时会造成系统的不稳定。但设计简单一般来说频率选择性越好相位的非线性越严重。 3、 窗口法和频率取样法的优缺点是什么 ? 答:IIR的极点可在单位圆的任何地方,可用较低的阶数获高的选择,有优良的幅度特性,但是相位是非线性的。 4、 脉冲响应不变法和双线性变化法各有什么特点? 答:脉冲响应不变法是数字的单位脉冲响应序列模仿模拟滤波器的冲击响应,让前者正好与后者的采样值相等。由于频域混叠效应,脉冲响应不变法是用于带限的滤波器设计,高通和带阻不予采用脉冲响应不变法,否则要加保护滤波器滤掉高于折叠频率即上的分量。对于带通和低通滤波器,要充分的带限,阻带衰减越大则混叠效应越小。 双线性变化法:可克服第一种方法的缺点——频率失真,可以从稳定的模拟滤波器获得稳定的数字滤波器,其代价是在频率轴上引入了失真。
第三章 随机信号功率谱估计
5.1 某随机信号由两余信号与噪声构成: X(t)=cos(20∏t)+cos(40∏t)+s(t),s(t)为均值为0,方差为1的高斯白噪声。 试分别用周期图法、平均周期图法分析该序列功率谱估计 解: >> Fs=1000; n=0:1/Fs:1; xn=cos(2*pi*20*n)+cos(2*pi*20*n)+randn(size(n)); Pxx=abs(fft(xn,1024)).^2/1001; subplot(2,2,1);plot(n,xn);xlabel('随机信号'); grid on; subplot(2,2,2);plot([0:1023]*Fs/1024,10*log10(Pxx)),xlabel('利用公式'); grid on; window=hamming(1001); [Pxx1,F1]=periodogram(xn,window,1024,Fs); subplot(2,2,3);plot(F1,10*log10(Pxx1)),xlabel('利用函数periodogram'); grid on; noverlap=500; [Pxx2,F2]=psd(xn,1024,Fs,window,noverlap); subplot(2,2,4);plot(F2,10*log10(Pxx2)),xlabel('利用函数psd'); grid on; 图 5.1 设计理论分析方法 数字滤波器的任务就是寻求一个因果、物理可实现的系统函数H(z),使它的频响H(ejw)满足所希望得到的低通频域指标,即通带衰减A p,阻带衰减A r,通带截频Wp,阻带截频Wr。 先找到数字技术指标通带衰减A p,阻带衰减A r,通带截频Wp,阻带截频Wr。然后变换成模拟技术指标,再把带通技术指标变换成低通技术指标,设计出一个模拟低通滤波器后再变换成数字带通滤波器。 设计结论及其分析 采用巴特沃斯设计数字带通滤波器其相位是非线性的,而且其图样随采样点的不同而变化,通带单调下降,使得系统的误差分布不均匀,这种滤波器具有非线性的相频滤波器 思考题 1、利用平均周期图法估计信号的功率谱,选用不同窗函数有什么变换
2、利用经典法AR模型参数发分别进行功率谱估计,结果有什么不同?
3、在AR模型参数发估计功率谱时,使用不同懂得模型阶数将会导致估计的功率谱有什么不同? |
|