數(shù)字信號(hào)處理實(shí)驗(yàn)報(bào)告二與三

實(shí)驗(yàn)二用FFT進(jìn)行譜分析一.實(shí)驗(yàn)?zāi)康模?進(jìn)一步加深DFT算法原理和根本性質(zhì)的理解〔因?yàn)镕FT只是DFT的一種快速算法，所以FFT的運(yùn)算結(jié)果必然滿足DFT的根本性質(zhì)〕。熟悉FFT程序結(jié)構(gòu)及編程方法。2熟悉應(yīng)用FFT對(duì)確定信號(hào)進(jìn)行譜分析方法，熟悉FFT算法原理和FFT子程序的應(yīng)用。3學(xué)習(xí)用FFT對(duì)連續(xù)信號(hào)和時(shí)域離散信號(hào)進(jìn)行譜分析的方法，了解可能出現(xiàn)的分析誤差及其原因，以便在實(shí)際中正確應(yīng)有FFT。二.實(shí)驗(yàn)內(nèi)容：用matlab編程產(chǎn)生并畫出信號(hào)x1(n)、x2(n)、x3(n)、x4(n)、x5(n)。圖1圖1用matlab編制FFT函數(shù)對(duì)上述信號(hào)進(jìn)行頻譜分析，并畫出上述信號(hào)譜圖。三.實(shí)驗(yàn)結(jié)果〔1〕1.%Thisprogrammistogeneratesignalx1(n)=R4(n).k=-6:6;x=[zeros(1,6),ones(1,4),zeros(1,3)];stem(k,x);〔信號(hào)圖如圖1〕title('圖1');2.n=-5:1:10;x=(n+1).*(n>=0&n<=3)+(8-n).*(n>=4&n<=7)+0;stem(n,x)；title('圖2');3.n=-5:10;x=(4-n).*(n>=0&n<=3)+(n-3).*(n>=4&n<=7);stem(n,x);title('圖3');4.n=-10:10;x=cos(pi/4*n);stem(n,x);title('圖4');5.n=-10:10;x=sin(pi/8*n);stem(n,x);title('圖5');實(shí)驗(yàn)結(jié)果〔2〕:FFT算法functiony=myditfft(x)%y=myditfft(x)%本程序?qū)斎胄蛄衳實(shí)現(xiàn)DIT-FFT基2算法，點(diǎn)數(shù)取大于等于x長度的2的冪次%x為給定時(shí)間序列%y為x的離散傅立葉變換m=nextpow2(x);N=2^m;%求x的長度對(duì)應(yīng)的2的最低冪次miflength(x)<N;%假設(shè)x的長度不是2的冪，補(bǔ)零到2的整數(shù)冪x=[x,zeros(1,N-length(x))];endnxd=bin2dec(fliplr(dec2bin([1:N]-1,m)))+1;%求1:2^m數(shù)列的倒序y=x(nxd);%將x倒序排列作為y的初始值formm=1:m;%將DFT作m次基2分解,從左到右，對(duì)每次分解作DFT運(yùn)算Nmr=2^mm;u=1;%旋轉(zhuǎn)因子u初始化為WN^0=1WN=exp(-i*2*pi/Nmr);%本次分解的根本DFT因子WN=exp(-i*2*pi/Nmr)forj=1:Nmr/2;%本次跨越間隔內(nèi)的各次蝶形運(yùn)算fork=j:Nmr:N;%本次蝶形運(yùn)算的跨越間隔為Nmr=2^mmkp=k+Nmr/2;%確定蝶形運(yùn)算的對(duì)應(yīng)單元下標(biāo)t=y(kp)*u;%蝶形運(yùn)算的乘積項(xiàng)y(kp)=y(k)-t;%蝶形運(yùn)算y(k)=y(k)+t;%蝶形運(yùn)算endu=u*WN;%修改旋轉(zhuǎn)因子,多乘一個(gè)根本DFT因子WNendend1.k=-6:6;x=[zeros(1,6),ones(1,4),zeros(1,3)];y=myditfft(x);k=-6:9;stem(k,y);xlabel('m');ylabel('X[M]');title('FFT圖');2.n=-5:1:10;x=(n+1).*(n>=0&n<=3)+(8-n).*(n>=4&n<=7)+0;y=myditfft(x);stem(n,y);xlabel('n');ylabel('X[M]');title('FFT圖');3.n=-5:10;x=(4-n).*(n>=0&n<=3)+(n-3).*(n>=4&n<=7);y=myditfft(x);stem(n,y);xlabel('n');ylabel('X[M]');title('FFT3');4.n=-10:10;x=cos(pi/4*n);y=myditfft(x);n=-10:21;stem(n,y);xlabel('n');ylabel('X[M]');title('FFT4');5.n=-10:10;x=sin(pi/8*n);y=myditfft(x);n=-10:21;stem(n,y);xlabel('n');ylabel('X[M]');title('FFT5');四.簡要答復(fù)以下問題：①在N=8時(shí)，x2(n)和x3(n)的幅頻特性會(huì)相同嗎？為什么？N=16呢？答：不相同。當(dāng)N=8時(shí)，序列x1(n)和x2(n)中相同的元素值對(duì)應(yīng)的n值是不同的，所乘的旋轉(zhuǎn)因子的值也不同，因而得到的最終結(jié)果也是不同的。同理，N=16時(shí)，所得的幅頻特性也是不同的。②FFT在什么條件下也可以用來分析周期信號(hào)序列的頻譜？如果正弦信號(hào)系統(tǒng)sin(2πf0k),f0=0.1Hz，用16點(diǎn)FFT來做DFT運(yùn)算，得到的頻譜是信號(hào)本身的真實(shí)譜嗎？為什么？答：由于FFT算法對(duì)序列長度的要求是N=2^M，M為正整數(shù)。所以，當(dāng)周期信號(hào)序列一個(gè)周期的長度滿足N=2^M〔M為正整數(shù)〕的條件時(shí)，F(xiàn)FT可以用來分析周期信號(hào)的頻譜。不是真實(shí)的頻譜。因?yàn)樾蛄械闹芷贜=10不是2的整數(shù)次冪，所以不是真實(shí)的。實(shí)驗(yàn)三：用雙線性變換法設(shè)計(jì)IIR數(shù)字濾波器一.實(shí)驗(yàn)?zāi)康模?.掌握用雙線性變換法設(shè)計(jì)IIRDF的原理及具體設(shè)計(jì)方法，熟悉用雙線性變換法設(shè)計(jì)IIRDF的計(jì)算機(jī)編程。2.觀察用雙線性變換法設(shè)計(jì)的DF的頻響特性，了解雙線性變換法的特點(diǎn)。3.熟悉用雙線性變換法設(shè)計(jì)BW和CB型DF的全過程。二，實(shí)驗(yàn)原理：為了克服沖激響應(yīng)不變法產(chǎn)生的頻率混疊現(xiàn)象，這是從S平面到Z平面的標(biāo)準(zhǔn)變換z＝esT的多值對(duì)應(yīng)關(guān)系導(dǎo)致的，為了克服這一缺點(diǎn)，產(chǎn)生了雙線性變換法。下式就是模擬濾波器和經(jīng)采樣后的數(shù)字濾波器之間的變換關(guān)系：三.實(shí)驗(yàn)內(nèi)容：1.讀懂所給參考程序，熟悉程序的整體結(jié)構(gòu)和功能。2.設(shè)計(jì)一個(gè)CBI型低通DF，通帶截頻fp=3000Hz，衰耗滿足Apmax=3dB，阻帶截頻fS=3400Hz,衰耗ASmin=31dB,取樣頻率fs=8000Hz。畫出其模擬濾波器及數(shù)字濾波器頻響特性曲線。程序：FS=8000;Fl=3000;Fh=3400;%通帶、阻帶截止頻率Rp=3;Rs=31;wp=Fl*2*pi;%臨界頻率采用角頻率表示ws=Fh*2*pi;%臨界頻率采用角頻率表示wp1=wp/FS;%求數(shù)字頻率ws1=ws/FS;%求數(shù)字頻率OmegaP=2*FS*tan(wp1/2);%頻率預(yù)畸OmegaS=2*FS*tan(ws1/2);%頻率預(yù)畸%選擇濾波器的最小階數(shù)[n,Wn]=cheb1ord(OmegaP,OmegaS,Rp,Rs,'s');%此處是代入經(jīng)預(yù)畸變后獲得的歸一化模擬頻率參數(shù)[bt,at]=cheby1(n,Rp,Wn,'s');%設(shè)計(jì)一個(gè)n階的契比雪夫I型模擬濾波器subplot(2,1,1);[h,w]=freqz(bt,at);plot(w*FS/(2*pi),abs(h));grid;xlabel('頻率/Hz');ylabel('幅值');title('模擬濾波器頻響特性');[bz,az]=bilinear(bt,at,FS);%雙線性變換為數(shù)字濾波器subplot(2,1,2);[H,W]=freqz(bz,az);%求解數(shù)字濾波器的頻率響應(yīng)plot(W*FS/(2*pi),abs(H));grid;xlabel('頻率/Hz');ylabel('幅值');title('數(shù)字濾波器頻響特性');3、設(shè)計(jì)一個(gè)BW型低通DF，滿足：通帶截頻fp=100Hz，衰耗滿足Apmax=3dB，阻帶截頻fS=400Hz,衰耗ASmin=15dB,取樣頻率fs=2000Hz。畫出其模擬濾波器及數(shù)字濾波器頻響特性曲線。FS=2000;Fl=100;Fh=400;%通帶、阻帶截止頻率Rp=3;Rs=15;wp=Fl*2*pi;%臨界頻率采用角頻率表示ws=Fh*2*pi;%臨界頻率采用角頻率表示wp1=wp/FS;%求數(shù)字頻率ws1=ws/FS;%求數(shù)字頻率OmegaP=2*FS*tan(wp1/2);%頻率預(yù)畸OmegaS=2*FS*tan(ws1/2);%頻率預(yù)畸%選擇濾波器的最小階數(shù)[n,Wn]=buttord(OmegaP,OmegaS,Rp,Rs,'s');%此處是代入經(jīng)預(yù)畸變后獲得的歸一化模擬頻率參數(shù)[bt,at]=butter(n,Wn,'s');%設(shè)計(jì)一個(gè)n階的契比雪夫I型模擬濾波器subplot(2,1,1);[h,w]=freqz(bt,at);plot(w*FS/(2*pi),abs(h));grid;xlabel('頻率/Hz');ylabel('幅值');title('模擬濾波器頻響特性');[bz,az]=bilinear(bt,at,FS);%雙線性變換為數(shù)字濾波器[H,W]=freqz(bz,az);%求解數(shù)字濾波器的頻率響應(yīng)subplot(2,1,2);plot(W*FS/(2*pi),abs(H));grid;xlabel('頻率/Hz');ylabel('幅值');title('數(shù)字濾波器頻響特性');設(shè)計(jì)一個(gè)BW型高通DF，滿足：通帶截頻fp=400Hz，衰耗滿足Apmax=3dB，阻帶截頻fS=350Hz,衰耗ASmin=15dB,取樣頻率fs=1000Hz。畫出其模擬濾波器及數(shù)字濾波器頻響特性曲線。fp=400;fs=350;ft=1000;rp=3;rs=15;wp=fp/(ft/2);ws=fs/(ft/2);%利用Nyquist頻率進(jìn)行歸一化[n,wc]=buttord(wp,ws,rp,rs,'s');%求模擬濾波器的最小階數(shù)和截止頻率[b,a]=butter(n,wc,'high','s');%設(shè)計(jì)高通模擬濾波器系數(shù)b，a%Omega1=linspace(0,1000,5);%[h,w]=freqs(b,a,Omega1);[h,w]=freqs(b,a);subplot(2,1,1);plot(w*ft/(2*pi),abs(h));grid;title('模擬濾波器');xlabel('頻率/Hz');ylabel('幅值');[N,Wc]=buttord(wp,ws,rp,rs);%求數(shù)字濾波器的最小階數(shù)和截止頻率[B,A]=butter(n,wc,'high');%設(shè)計(jì)高通數(shù)字濾波器系數(shù)b，a[H,W]=freqz(B,A);%繪出頻率響應(yīng)曲線subplot(2,1,2);plot(W*ft/(2*pi),abs(H));grid;title('數(shù)字濾波器');xlabel('頻率/Hz');ylabel('幅值');設(shè)計(jì)一個(gè)CB型帶通DF，滿足：通帶邊界頻率為100Hz～500Hz，通帶衰耗小于3dB，過渡帶寬20Hz，阻帶衰耗大于15dB，取樣頻率fs=2000Hz。畫出其模擬濾波器及數(shù)字濾波器頻響特性曲線。fp=[100,500];fs=[80,520];ft=2000;rp=3;rs=15;wp=fp/(ft/2);ws=fs/(ft/2);[n,wn]=cheb2ord(wp,ws,rp,rs,'s');[b,a]=cheby2(n,rs,wn,'s');[h,w]=freqs(b,a);subplot(2,1,1);plot(w*ft/(2*pi),abs(h));grid;title('模擬濾波器');xlabel('頻率/Hz');ylabel('幅值');[N,Wn]=cheb2ord(wp,ws,rp,rs);[B,A]=cheby2(N,r