《MATLAB信号频谱分析.doc》由会员分享，可在线阅读，更多相关《MATLAB信号频谱分析.doc（20页珍藏版）》请在三一文库上搜索。
1、 题1周期信号频谱分析题目：给定输入信号x(n)为一周期信号，含基波及若干次谐波分量。要求：1、给出信号基波频率及谐波个数 2、基波及各个谐波信号的相位和幅值clearclcA = 0.25 0.5 0.25 ;f = 1 2 4;fi = 0 0.5*pi 0.75*pi;N = 59;fs = 10;t = ( 1 : N ) ./ fs; x_n =A(1).*sin(2*pi*f(1)*t+fi(1)+A(2).*sin(2*pi.*f(2)*t+fi(2)+A(3).* sin(2 pi.*f(3)*t+ fi( 3 ) );save(Signal.mat, x_n )程序如下fs=
2、10;N=60;n=1:N;t=n/fs;x =0.25*sin(2*pi*t)+0.5*sin(2*pi*2*t)+0.25*sin(2*pi*4*t);figure(1)subplot (211); plot(x); grid;y=fft(x,N);mag=abs(y);k=0:length(y)-1;f=fs/N*k;subplot(212);plot(f,mag);xlabel(Frequence(Hz);ylabel(Magnitude);grid;title(对周期信号的一个周期进行频谱分析)figure(2)x=x x x x x ;figure(2)subplot(211);p
3、lot(x);grid;N=length(x)y=fft(x,N);mag=abs(y);k=0:length(y)-1;f=fs/N*ksubplot(212);plot(f,mag);xlabel(Frequence(Hz);ylabel(Magnitude);grid;title(对周期信号的多个周期进行频谱分析)题2数字滤波器设计题目：已知 输入信号x(t)=0.6sin(200t)+sin(400t) +0.3sin(800t) 。要求: 通过MATLAB编程设计一 1)低通滤波器,输出信号x(t)的最低频率成份. 2)高通滤波器,输出信号x(t)的最高频率成份. 3)带通滤波器,输
4、出信号x(t)的中间频率成份. 4)带阻滤波器,输出信号x(t)的最低+最高频率成份. (Rp=1dB,Rs=16dB,Fs=3倍最高频率,截止频率、滤波器设计方法及类型自己选择.)解：编写程序如下：Fs=1200;N=1024;n=0:N;t=n/Fs;x=0.6*sin(200*pi*t)+sin(400*pi*t)+0.3*sin(800*pi*t);y=fft(x,N);mag=abs(y);k=0:length(y)-1;f=Fs/N*k;figure(1)plot(f,mag);xlabel(Frequence(Hz);ylabel(Magnitude);title(N=1024)
5、;grid;Wp=0.12*pi;Ws=0.2*pi;rp=1;rs=16;Fs=1200;%数字滤波器指标wp=Wp*Fs;ws=Ws*Fs;%模拟技术指标n,wn=buttord(wp,ws,rp,rs,s);z,p,k=buttap(n);b,a=zp2tf(z,p,k);b,a=lp2lp(b,a,wn);figure(2)h1=freqs(b,a);freqs(b,a);title(模拟低通的幅频和相频特性图);bz,az=impinvar(b,a,Fs);figure(3)freqz(bz,az,256);title(数字底通的幅频和相频特性图);y1=filter(bz,az,x
6、);y=fft(y1,N);mag=abs(y);k=0:length(y)-1;f=Fs/N*k;figure(10)plot(f,mag);grid;title(输出信号x(t)的最低频率成份)wp=.4*pi;ws=.35*pi;rp=1;rs=16;Fs=1200;%数字滤波器指标wap=2*Fs*tan(wp/2);was=2*Fs*tan(ws/2);n,wn=buttord(wap,was,rp,rs,s);z,p,k=buttap(n);b,a=zp2tf(z,p,k);bt,at=lp2hp(b,a,wap);w=0:pi:4000*pi;figure(4)freqs(bt,
7、at,w);title(模拟高通的幅频和相频特性图);bz,az=bilinear(bt,at,Fs);figure(7)freqz(bz,az,256,1);title(数字高通的幅频和相频特性图);y2=filter(bz,az,x);y=fft(y2,N);mag=abs(y);k=0:length(y)-1;f=Fs/N*k;figure(11)plot(f,mag);grid;title(输出信号x(t)的最高频率成份)wp=0.2*pi 0.35*pi; ws=0.12*pi 0.4*pi;rp=1;rs=30;Fs=1200;%数子滤波器指标wap=2*Fs*tan(wp./2)
8、; was=2*Fs*tan(ws./2);%模拟滤波器指标 n,wn=buttord(wap,was,rp,rs,s);z,p,k=buttap(n);bp,ap=zp2tf(z,p,k);bw=wap(2)-wap(1); w0=sqrt(wap(1)*wap(2);bs,as=lp2bp(bp,ap,w0,bw);w=0:pi:4000*pi;figure(6)freqs(bs,as,w);title(模拟带通的幅频和相频图);bz,az=bilinear(bs,as,Fs);figure(7)freqz(bz,az,256,Fs);title(数字带通的幅频和相频图)y3=filter
9、(bz,az,x);y=fft(y3,N);mag=abs(y);k=0:length(y)-1;f=Fs/N*k;figure(12)plot(f,mag);grid;title(输出信号x(t)的中间频率成份)wp=0.15*pi 0.4*pi;ws=0.2*pi 0.35*pi;rp=1;rs=16;Fs=1200;%数值指标wap=2*Fs*tan(wp./2);was=2*Fs*tan(ws./2);%模拟指标n,wn=buttord(wap,was,rp,rs,s);z,p,k=buttap(n);b,a=zp2tf(z,p,k); bw=wap(2)-wap(1);w0=sqrt
10、(wap(1)*wap(2);bt,at=lp2bs(b,a,w0,bw);figure(8)w=0:pi:4000*pi;freqs(bt,at,w);title(模拟带阻的幅频和相频图)bz1,az1=bilinear(bt,at,Fs);%模拟到数字figure(9)freqz(bz1,az1,256,Fs);title(数字带阻的幅频和相频图)y4=filter(bz1,az1,x);y=fft(y4,N);mag=abs(y);k=0:length(y)-1;f=Fs/N*k;figure(13)plot(f,mag);grid;title(输出信号x(t)的最低+最高频率成份)题3
11、语音信号采集与处理初步题目：1语音信号的采集2语音信号的频谱分析 3设计数字滤波器和画出频率响应 4用滤波器对信号进行滤波 5比较滤波前后语音信号的波形及频谱 6回放和存储语音信号 解：编程如下x,fs,bits=wavread(类似爱情.wav);wavplay(x,fs);X=fft(x,4096);magX=abs(X);angX=angle(X);figure(1)subplot(221);plot(x);title(原始信号波形);subplot(222);plot(X);title(原始信号频谱);subplot(223);plot(magX);title(原始信号幅值);subp
12、lot(224);plot(angX);title(原始信号相位);n=length(x);t=0:1/fs:(n-1)/fs;noise=0.5*sin(100*pi*t)+sin(400*pi*t)+0.5*sin(1600*pi*t);figure(2)plot(t,noise);grid;xlabel(时间 (t);ylabel(幅度（y);title(噪声的时域图);Noise=fft(noise,n);N1=abs(Noise);n1=floor(n/2);f=(0:n1)*fs/n;figure(3);%画出噪声的频谱图plot(f,N1(1:n1+1);grid on; xla
13、bel(频率（f）);ylabel( 幅度（N1）);title(噪声的频谱图);x=x(:,1);x1=noise+x;sound(x1,fs);wavwrite(x1,16,类似爱情1.wav);figure(4) %画出加噪声前后的时域比较图subplot(221);plot(t,x);grid on;xlabel(时间 (t);ylabel(幅度（x); title(加噪声前的时域图);subplot(222);plot(t,x1);grid on;xlabel(时间 (t); ylabel(幅度（Y1); title(加噪声后的时域图);X=fft(x,n); %对加噪声前的采样信号
14、1进行傅立叶变换MagX=abs(X);X1=fft(x1,n); %对加噪声后的采样信号1进行傅立叶变换MagX1=abs(X1);subplot(223); plot(f,MagX(1:n1+1);grid on;xlabel(频率（f）); ylabel( 幅度（MagX）); title(加噪声前的频谱图);subplot(224); plot(f,MagX1(1:n1+1);grid on;xlabel(频率（f）);ylabel( 幅度（MagX1）); title(加噪声后的频谱图);% 构造800Hz的带阻滤波器 %f0_stop1=800; fc=20; Rp=1; Rs=3
15、0; ws_stop1=(f0_stop1-fc)/(fs/2) (f0_stop1+fc)/(fs/2); % 通带的拐角频率wp_stop1=(f0_stop1-5*fc)/(fs/2) (f0_stop1+5*fc)/(fs/2); % 阻带的拐角频率N_stop1,wn_stop1=cheb1ord(wp_stop1,ws_stop1,Rp,Rs); % 求出切比雪夫I型滤波器的阶数N及频率参数wcnum_stop1,den_stop1=cheby1(N_stop1,Rp,wn_stop1,stop); % 求出切比雪夫I型带阻数字滤波器的传递函数模型系数h_stop1,w_stop1
16、=freqz(num_stop1,den_stop1,fs) % 求出离散系统频率响应的数值figure(5);subplot(231); % 画出带阻滤波器的幅频特性图plot(w_stop1*fs/(2*pi),20*log10(abs(h_stop1);xlabel(频率（w_stop1）);ylabel(幅度（h_stop1);title(带阻滤波器在800Hz处的幅频特性图);grid on;B_stop1=filter(num_stop1,den_stop1,x1); %对含噪信号x1进行带阻滤波pause(3);sound(B_stop1,fs);% 读出滤去800Hz噪声后的采
17、样信号wavwrite(B_stop1,fs,16,类似爱情2.wav);subplot(234); % 画出带阻滤波器滤去800Hz噪声后的时域图plot(t,B_stop1); xlabel(时间 (t)）);ylabel(幅度（B_stop1);title(用带阻滤波器滤去800Hz噪声后的时域图);grid on;% % 构造200Hz的带阻滤波器 %f0_stop2=200;ws_stop2=(f0_stop2-fc)/(fs/2) (f0_stop2+fc)/(fs/2); % 通带的拐角频率wp_stop2=(f0_stop2-5*fc)/(fs/2) (f0_stop2+5*f
18、c)/(fs/2); % 阻带的拐角频率N_stop2,wn_stop2=cheb1ord(wp_stop2,ws_stop2,Rp,Rs); % 求出切比雪夫I型滤波器的阶数N及频率参数wcnum_stop2,den_stop2=cheby1(N_stop2,Rp,wn_stop2,stop); % 求出切比雪夫I型带阻数字滤波器的传递函数模型系数h_stop2,w_stop2=freqz(num_stop2,den_stop2,fs) ; % 求出离散系统频率响应的数值subplot(232);plot(w_stop2*fs/(2*pi),20*log10(abs(h_stop2);xla
19、bel(频率（w_stop2）);ylabel(幅度（h_stop2);title(带阻滤波器在200Hz处的幅频特性图);grid on;B_stop2=filter(num_stop2,den_stop2,B_stop1);%对含噪信号Y1进行带阻滤波pause(3);sound(B_stop2,fs); % 读出滤去800Hz和200Hz噪声后的采样信号 wavwrite(B_stop2,fs,16,类似爱情3.wav);subplot(235); % 画出带阻滤波器滤去800Hz和200Hz噪声后的时域图plot(t,B_stop2);xlabel(时间 (t)）);ylabel(幅度
20、（B_stop2);title(带阻滤波器滤去800Hz和200Hz噪声后的时域图);grid on;% 构造50Hz的带阻滤波器 %f0_stop3=50;ws_stop3=(f0_stop3-fc)/(fs/2) ; % 阻带的拐角频率 wp_stop3=(f0_stop3+5*fc)/(fs/2) ; % 通带的拐角频率Wws=0.0075;Wwp=0.0375;N_stop3,wc_stop3=cheb1ord(Wwp,Wws,Rp,Rs,s); % 求出切比雪夫I型滤波器的阶数N及频率参数wcnum_stop3,dem_stop3=cheby1(N_stop3,Rp,wc_stop3
21、,high); % 求出切比雪夫I型滤波器高通数字滤波器的传递函数模型系数h_stop3,w_stop3=freqz(num_stop3,dem_stop3,fs); % 求出离散系统频率响应的数值subplot(233); %画出高通滤波器的幅频特性图plot(w_stop3*fs/(2*pi),(abs(h_stop3);xlabel(频率（w_stop3）);ylabel(幅度（h_stop3);title(高通滤波器在50Hz处的幅频特性图);grid on;B_stop3=filter(num_stop3,dem_stop3,B_stop2);pause(3);sound(B_sto
22、p3,fs); % 读出滤去800Hz和200Hz噪声后的采样信号 wavwrite(B_stop3,fs,16,类似爱情4.wav);subplot(236); % 画出高通滤波器滤去800Hz和200Hz和50Hz噪声后的时域图plot(t,B_stop3);xlabel(时间 (t)）);ylabel(幅度（B_stop3);title(滤去800Hz和200Hz和50Hz噪声后的时域图);grid on;% 用不同的滤波器滤去相应噪声频率后的时域比较 %figure(6);subplot(221);plot(t,x1);grid on; xlabel(时间 (t);ylabel(幅度（
23、x1);title(加噪声后的时域图);subplot(222);plot(t,B_stop1);xlabel( 时间 (t)）);ylabel(幅度（B_stop1);title(用带阻滤波器滤去800Hz噪声后的时域图);grid on;subplot(223);plot(t,B_stop2);xlabel( 时间 (t)）);ylabel(幅度（B_stop2);title( 带阻滤波器滤去800Hz和200Hz噪声后的时域图);grid on;subplot(224);plot(t,B_stop3);xlabel( 时间 (t)）);ylabel(幅度（B_stop3);title(
24、高通滤波器滤去800Hz和200Hz和50Hz噪声后的时域图);grid on; % 用不同的滤波器滤去相应噪声频率后的频域比较 %figure(7);X1=fft(x1,n); %对加噪声后的采样信号1进行傅立叶变换 MagX1=abs(X1); n1=floor(n/2); f=(0:n1)*fs/n; X1_B_stop1=fft(B_stop1,n); %对滤去800Hz噪声后的采样信号进行傅立叶变换 MagX1_B_stop1=abs(X1_B_stop1); X1_B_stop2=fft(B_stop2,n); %对滤去800Hz和200Hz噪声后的采样信号进行傅立叶变换 MagX
25、1_B_stop2=abs( X1_B_stop2); X1_B_stop3=fft(B_stop3,n); %对滤去800Hz和200Hz和50Hz噪声后的采样信号进行傅立叶变换 MagX1_B_stop3=abs( X1_B_stop3); subplot(221);plot(f, MagX1(1:n1+1);grid on; xlabel(频率（f）);ylabel( 幅度（X1）); title(加噪声后的频谱图); subplot(2,2,2);plot(f, MagX1_B_stop1(1:n1+1);grid on; xlabel(频率（f）);ylabel( 幅度（ X1_B_
26、stop1）); title(用带阻滤波器滤去800Hz噪声后的的频谱图); subplot(2,2,3);plot(f, MagX1_B_stop2(1:n1+1);grid on; xlabel(频率（f）);ylabel( 幅度（ X1_B_stop2）); title(用带阻滤波器滤去800Hz和200Hz噪声后的的频谱图); subplot(2,2,4);plot(f, MagX1_B_stop3(1:n1+1);grid on; xlabel(频率（f）);ylabel( 幅度（ X1_B_stop3）); title(用高通滤波器滤去800Hz和200Hz和50Hz噪声后的的频谱
27、图); 程序二：x,fs,bits=wavread(toilet.wav);wavplay(x);pause(3); wavplay(x,11025);pause(3)wavplay(x,44100);X=fft(x,4096);magX=abs(X);angX=angle(X);figure(1)subplot(221);plot(x);title(原始信号波形);subplot(222);plot(X);title(原始信号频谱);subplot(223);plot(magX);title(原始信号幅值);subplot(224);plot(angX);title(原始信号相位);paus
28、e(5)%N阶高通滤波器x=wavread(toilet.wav);sound(x);N=5;wc=0.3;b,a=butter(N,wc,high);X=fft(x);figure(2)subplot(321);plot(x);title(滤波前信号的波形);subplot(322);plot(X);title(滤波前信号的频谱);y=filter(b,a,x);wavwrite(y,toilet1.wav);wavplay(y)Y=fft(y);subplot(323);plot(y);title(IIR滤波后信号的波形);subplot(324);plot(Y);title(IIR滤波后
29、信号的频谱);z=fftfilt(b,x);wavwrite(z,toilet2.wav);wavplay(z)Z=fft(z);subplot(325);plot(z);title(FIR滤波后信号的波形);subplot(326);plot(Z);title(FIR滤波后信号的频谱);pause(5)%N阶低通滤波器x=wavread(toilet.wav);sound(x);X=fft(x);figure(3)subplot(321);plot(x);title(滤波前信号的波形);subplot(322);plot(X);title(滤波前信号的频谱);N=5;wc=0.3;b,a=b
30、utter(N,wc);y=filter(b,a,x);wavwrite(y,toilet3.wav);wavplay(y)Y=fft(y);subplot(323);plot(y);title(IIR滤波后信号的波形);subplot(324);plot(Y);title(IIR滤波后信号的频谱);z=fftfilt(b,x);wavwrite(z,toilet4.wav);wavplay(z)Z=fft(z);subplot(325);plot(z);title(FIR滤波后信号的波形);subplot(326);plot(Z);title(FIR滤波后信号的频谱);pause(5)%2N
31、阶带通滤波器x=wavread(toilet.wav);sound(x);N=5;wc=0.3 ,0.6;b,a=butter(N,wc);X=fft(x);figure(4)subplot(321);plot(x);title(滤波前信号的波形);subplot(322);plot(X);title(滤波前信号的频谱);y=filter(b,a,x);wavwrite(y,toilet5.wav);wavplay(y)Y=fft(y);subplot(323);plot(y);title(IIR滤波后信号的波形);subplot(324);plot(Y);title(IIR滤波后信号的频谱)
32、;z=fftfilt(b,x);wavwrite(z,toilet6.wav);wavplay(z)Z=fft(z);subplot(325);plot(z);title(FIR滤波后信号的波形);subplot(326);plot(Z);title(FIR滤波后信号的频谱);pause(5)%2N阶带阻滤波器x=wavread(toilet.wav);sound(x);N=5;wc=0.2 ,0.7;b,a=butter(N,wc,stop);X=fft(x);figure(5)subplot(321);plot(x);title(滤波前信号的波形);subplot(322);plot(X);title(滤波前信号的频谱);y=filter(b,a,x);wavwrite(y,toilet7.wav);wavplay(y)Y=fft(y);subplot(323);plot(y);title(IIR滤波后信号的波形);subplot(324);plot(Y);title(IIR滤波后信号的频谱);z=fftfilt(b,x);wavwrite(z,toilet8.wav);wavplay(z)Z=fft(z);subplot(325);plot(z);title(FIR滤波后信号的波形);subplot(326);plot(Z);title(FIR滤波后信号的频谱);