實驗一：低通采樣定理和內(nèi)插與抽取實現(xiàn)

1、實驗一：低通采樣定理和內(nèi)插與抽取實現(xiàn)一 實驗?zāi)康?. 連續(xù)信號和系統(tǒng)的表示方法，以及坊真方法。2. 用 MATLAB 實現(xiàn)連續(xù)信號采用與重構(gòu)的方法，3. 采樣信號的插值和抽取等重采樣實現(xiàn)方法。4. 用時域采樣信號重構(gòu)連續(xù)時域信號的原理和方法。5. 用 MATLAB 繪圖函數(shù)表示信號的基本方法，實驗數(shù)據(jù)的可視化表示。二 原理1 、時域抽樣定理令連續(xù)信號 xa(t)的傅里葉變換為 Xa(j ), 抽樣脈沖序列 p(t)傅里葉變換為 P(j ),抽樣后的信號 x(t)的傅里葉變換為 X(j )若采用均勻抽樣 , 抽樣周期 Ts, 抽樣頻率為 s=2 fs, 由前面分析可知 :抽樣的過程可以通過抽樣脈

2、沖序列 p(t) 與連續(xù)信號 xa(t)相乘來完成 , 即滿足 : x(t)=xa(t) p(t) ，又周期信號 f(t) 傅里葉變換 為：故可以推得 p(t)的傅里葉變換為 :其中:根據(jù)卷積定理可知 :得到抽樣信號 x(t) 的傅里葉變換為 :其表明 :信號在時域被抽樣后 ,他的頻譜 X(j ) 是連續(xù)信號頻譜 X(j )的形狀 以抽樣頻率 為間隔周期重復(fù)而得到 ,在重復(fù)過程中幅度被 p(t)的傅里葉級數(shù) Pn 加權(quán)。因為Pn只是 n的函數(shù),所以 X(j )在重復(fù)的過程中不會使其形狀發(fā)生變化。假定信號 x(t)的頻譜限制在 - m+ m 的范圍內(nèi), 若以間隔 Ts對 xa(t)進(jìn)行 抽樣,可

3、知抽樣信號 X(t)的頻譜 X(j )是以 s 為周期重復(fù)。顯然 ,若在抽樣的 過程中 s<2 m,則 X(j ) 將發(fā)生頻譜混疊現(xiàn)象 ,只有在抽樣的過程中滿足 s>=2 m 條件,X(j )才不會產(chǎn)生頻譜的混疊 ,接收端完全可以由 x(t) 恢復(fù)原 連續(xù)信號 xa(t),這就是低通信號抽樣定理的核心內(nèi)容。2、信號的重建從頻域看 ,設(shè)信號最高頻率不超過折疊頻率Xa(j )=Xa(j ) | |< s/2Xa(j )=0 | |> s/2 則理想取樣后的頻譜就不會產(chǎn)生混疊 ,故有 :讓取樣信號 x(t) 通過一帶寬等于折疊頻率的理想低通濾波器H(j )=T | |<

4、 s/2H(j )=0 | |> s/2 濾波器只允許通過基帶頻譜 ,即原信號頻譜 ,故:Y(j )=X(j )H(j )=Xa(j )因此在濾波器的輸出得到了恢復(fù)的原模擬信號 :y(t)=xa(t)從時域上看 ,上述理想的低通濾波器的脈沖響應(yīng)為根據(jù)卷積公式可求得理想低通濾波器的輸出為由上式顯然可得則:上式表明只要滿足取樣頻率高于兩倍信號最高頻率 ,連續(xù)時間函數(shù) xa(t)就可 用他的取樣值 xa(nT)來表達(dá)而不損失任何信息， 這時只要把每一個取樣瞬時值與 內(nèi)插函數(shù)式相乘求和即可得出 xa(t),在每一取樣點上 ,由于只有該取樣值所對應(yīng)的 內(nèi)插函數(shù)式不為零 ,所以各個取樣點上的信號值不

5、變。內(nèi)容1. 連續(xù)時間信號線性濾波實現(xiàn)設(shè)計一個 Butterworth 模擬帶通濾波器，設(shè)計指標(biāo)為：通帶頻率： 1000- 2000Hz，兩側(cè)過渡帶寬 500Hz，通帶波紋 1dB，阻帶衰減 100dB。假設(shè)一個信號， 其中 f1=100Hz,f2=1500Hz,f3=2900Hz 。信號的采樣頻率為 10000Hz。試將原信號 與通過該濾波器的模擬信號進(jìn)行比較。參考程序如下：wp=1000 2000*2*pi;ws=500 2500*2*pi;Rp=1;Rs=100; %濾波器設(shè)計參數(shù)，對于給定 Hz 應(yīng)乘以 2N,Wn=buttord(wp,ws,Rp,Rs,'s'); %

6、求得濾波器的最小階數(shù)和截止頻率w=linspace(1,3000,1000)*2*pi;%設(shè)置繪制頻率響應(yīng)的頻率點b,a=butter(N,Wn,'s');%設(shè)計模擬 Butterworth 濾波器H=freqs(b,a,w);%計算給定頻率點的復(fù)數(shù)頻率響應(yīng)magH=abs(H);phaH=unwrap(angle(H); % 計算幅頻響應(yīng)和相頻響應(yīng) plot(w/(2*pi),20*log10(magH);%以頻率為橫坐標(biāo)繪制幅頻響應(yīng)xlabel(' 頻率 /Hz');ylabel(' 振幅 /dB'); title('Butterwo

7、rth 模擬帶通濾波器 ');hold on;plot(1000 1000,ylim,'r');plot(2000 2000,ylim,'r');% 繪帶邊界 grid on figure(2) dt=1/10000;%模擬信號采樣間隔f1=100;f2=1500;f3=2900;% 輸入信號的三個頻率成分 t=0:dt:0.04;%給定模擬時間段x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); % 輸入信號 H=tf(b,a); %濾波器在 MATLAB系統(tǒng)中的表示 y,t1=lsim(H,

8、x,t);%模擬輸出subplot(2,1,1),plot(t,x),title(' 輸入信號 ') % 繪出輸入信號 subplot(2,1,2),plot(t1,y)%繪制輸出信號title(' 輸出信號 '),xlabel(' 時間 /s')2. 理想矩形濾波器的時域表示clear; clc; fh=100; for I=1:4 k=I; fs=k*2*fh; N=10*k; n=-N:N; dt=1/fs; T=N*dt; t=-T:dt:T;%h=2*fh/fs*sinc(n/k);%h=sinc(n/k); subplot(2,2,I

9、) plot(t,h); hold on ； stem(t,h); hold on; plot(t,zeros(length(t),'linewidth',3); title('fs/2fh=k,k=',num2str(k),'fontsize',28); axis('off') endfs/2fh=k,k=2fs/2fh=k,k=4fs/2fh=k,k=1 fs/2fh=k,k=33. 連續(xù)時間信號的采樣和重建 、分別用 150HZ及 300HZ對信號采樣源信號為：fa=5*sin(2*pi*40*t1)+1.8*sin(4*p

10、i*40*t1)+0.8*sin(5*pi*40*t1) ，用 150Hz 的頻率 對 f(t) 進(jìn)行采樣，其采樣圖如圖 1 所示；用 300Hz 的頻率對 f(t) 進(jìn)行采樣，其采樣圖如圖 2 所示。程序如下： fs1=150;t1=-0.1:1/fs1:0.1; fa=5*sin(2*pi*40*t1)+1.8*sin(4*pi*40*t1)+0.8*sin(5*pi*40*t1); figure(1);plot(t1,fa),xlabel('fs1=150Hz 時 ,fa 采樣時域圖 '); hold off;fs2=300;t2=-0.1:1/fs2:0.1;fb=5*

11、sin(2*pi*40*t2)+1.8*sin(4*pi*40*t2)+0.8*sin(5*pi*40*t2);時,fb 采樣時域圖 ');0.1figure(2);plot(t2,fb),xlabel('fs2=300Hz-8-0.1-0.08-0.06-0.04-0.020 0.020.040.060.08fs1=150Hz時 ,fa 采 樣 時 域 圖8-0.1 -0.08 -0.06 -0.04 -0.020 0.02 0.04 0.06 0.08 0.1fs2=300Hz時 ,fb 采 樣 時 域 圖圖1 150HZ采樣頻率對信號采樣圖圖 2 300HZ采樣頻率對信號

12、采樣圖 、對信號進(jìn)行快速離散傅里葉變換將兩個采樣信號進(jìn)行快速離散傅里葉變換 (FFT), 用 150Hz的頻率對 f(t) 進(jìn)行 采樣，其采樣后快速傅立葉變換頻譜圖圖 3 所示；用 300Hz 的頻率對 f(t) 進(jìn)行 采樣，其采樣后快速傅立葉變換頻譜圖圖 4 所示。程序如下：f=40;fs=150;N=300;k=0:N-1;t=-0.1:1/fs:0.1;w1=150*k/N;fa=5*sin(2*pi*f*t)+1.8*sin(4*pi*f*t)+0.8*sin(5*pi*f*t);xfa=fft(fa,N);xf1=abs(xfa);時,fa 經(jīng)fft 后頻譜圖 . 單位： Hz

13、9;);figure(1);plot(w1,xf1),xlabel('fs=150Hz f=40;fs=300;N=300;k=0:N-1;t=-0.1:1/fs:0.1;w2=300*k/N fb=5*sin(2*pi*f*t)+1.8*sin(4*pi*f*t)+0.8*sin(5*pi*f*t);xfb=fft(fb,N);xf2=abs(xfb);figure(2);plot(w2,xf2),xlabel('fs=300Hz 時,fb 經(jīng)fft 后頻譜圖 .單位： Hz ');80fs=150Hz 時 ,fa 經(jīng) fft 后 頻 譜 圖 . 單 位 ： Hz圖

14、3 150HZ采樣后經(jīng) FFT后頻譜圖圖4 300HZ采樣后經(jīng) FFT后頻譜圖 、信號的重建我們可以通過利用內(nèi)插法把原信號從采樣信號中恢復(fù)出來，觀察信號在滿 足怎樣的采樣條件下能夠恢復(fù)為原信號，圖 5和圖 6分別為恢復(fù)后的原信號 程序如下：Wm=180*pi;Wc=Wm;fs1=150;Ws=2*pi*fs1; n=-800:800;nTs1=n/fs1; fa=5.1*sin(2*pi*40*nTs1)+1.8*sin(4*pi*40*nTs1)+0.8*sin(5*pi*40*nTs1); Dt=1/fs1;t1=-0.1:Dt:0.1;fa1=fa/fs1*Wc/pi*sinc(Wc/p

15、i)*(ones(length(nTs1),1)*t1-nTs1'*ones(1,length( t1);figure(1);plot(t1,fa1);axis(-0.1 0.1 -8 8); xlabel('fs=150Hz 時 ,fa 利用內(nèi)插由樣本重建原信號圖 .');Wm=180*pi;Wc=Wm;fs2=300;Ws=2*pi*fs2; n=-800:800;nTs2=n/fs2;fb=5.1*sin(2*pi*40*nTs2)+1.8*sin(4*pi*40*nTs2)+0.8*sin(5*pi*40*nTs2); Dt=1/fs2;t1=-0.1:Dt:0

16、.1;fb1=fb/fs2*Wc/pi*sinc(Wc/pi)*(ones(length(nTs2),1)*t1-nTs2'*ones(1,length( t1);figure(2);plot(t1,fb1);axis(-0.1 0.1 -8 8); xlabel('fs=300Hz 時 ,fb 利用內(nèi)插由樣本重建原信號圖 .'); grid;8圖5 150HZ 采樣后的信號的重建信號8圖6 300HZ 采樣后的信號的重建信號4. 采樣信號的抽取和插值t = 0:.00025:1;% Time vectorx = sin(2*pi*30*t) + sin(2*pi*60

17、*t); figure(1) subplot(211) stem(t(1:120),x(1:120);hold on y = decimate(x,4);% View the original and decimated signals: stem(x(1:120), axis(0 120 -2 2) % Original signaltitle( 'Original Signal' )subplot(212)stem(y(1:30)title( 'Decimated Signal' )% Decimated signalDecimated Signal2%Fi

18、le_C3:UpSampl.m%該程序仿真通過零階保持內(nèi)插對信號進(jìn)行上抽樣，分析上抽樣對信號頻譜的影響clear allclcf0=0.06;% 信號數(shù)字頻率N=256;dt=1;% 抽樣時間t=0:N-1*dt;sig=sin(2*pi*f0*t);I=4; % 內(nèi)插因子N1=N*I;Addsig(1:N1)=0;for k=1:Nfor i=1:IdLIAddsig(I*k-I+i)=sig(k); % 對原始數(shù)字信號補零內(nèi)插， 抽樣時間變?yōu)?endendfigure(1)subplot(2,1,1)t2=0:N1-1;stem(t,sig);title('Original'

19、;)xlabel('Sample time /n')axis(0 19 -1.2 1.2)subplot(2,1,2)stem(Addsig);axis(0 80 -1.2 1.2)xlabel('Sample time /n')title('Signal output of Zero-order hold interpolator')Am=fft(sig);AddAm=fft(Addsig);f1=0:2/N:2-1/N;f2=0:2/N1:2-1/N1;figure(2)subplot(2,1,1) plot(f1,20*log10(abs(

20、Am);axis(0 1 -20 60);ylabel('Amplitude(dB)')xlabel('omega/ pi')title('Amplitude Spectrum')subplot(2,1,2)plot(f2,20*log10(abs(AddAm)axis(0 1 -20 60);ylabel('Amplitude(dB)')xlabel('omega/ pi')title('Amplitude Spectrum of Zero-order hold interpolation signal&

21、#39;)b,a=butter(20,1/I);% 設(shè)計截止頻率為 pi/l，階數(shù)為 10 的低通巴特沃思濾波器y=filter(b,a,Addsig);% 對補零后的信號進(jìn)行低通濾波，完成上抽樣過程UpAm=fft(y);figure(3) subplot(2,1,1) plot(y);axis(0 1023 -1.2 1.2) title('output signal out of Interpolator') xlabel('Sample time /n') subplot(2,1,2)plot(f2,20*log10(abs(UpAm) axis(0 1

22、 -30 50);ylabel('Amplitude(dB)') xlabel('omega/ pi')%File_C3:UpSampl.m%該程序仿真通過零階保持內(nèi)插對信號進(jìn)行上抽樣，分析上抽樣對信號頻譜的影響 clear allclcf0=0.06; % 信號數(shù)字頻率N=256;dt=1; % 抽樣時間t=0:N-1*dt; sig=sin(2*pi*f0*t);I=4;% 內(nèi)插因子N1=N*I;Addsig(1:N1)=0;for k=1:Nfor i=1:IdLIAddsig(I*k-I+i)=sig(k); % 對原始數(shù)字信號補零內(nèi)插， 抽樣時間變?yōu)?

23、endendfigure(1) subplot(2,1,1) t2=0:N1-1; stem(t,sig); title('Original') xlabel('Sample time /n') axis(0 19 -1.2 1.2) subplot(2,1,2) stem(Addsig);axis(0 80 -1.2 1.2)xlabel('Sample time /n')title('Signal output of Zero-order hold interpolator') Am=fft(sig);AddAm=fft(Addsig); f1=0:2/N:2-1/N;f2=0:2/N1:2-1/N1;figure(2) subplot(2,1,1) plot(f1,20*log10(abs(Am);axis(0 1 -20 60); ylabel('Amplitude(dB)') xlabel('omega/ pi') title('Amplitude