关于用MATLAB设计对信号进行频谱分析和滤波处理的程序

完整的程序

%写上标题

%设计低通滤波器：

[N,Wc]=buttord()

%估算得到Butterworth低通滤波器的最小阶数N和3dB截止频率Wc

[a,b]=butter(N,Wc) %设计Butterworth低通滤波器

[h,f]=freqz() %求数字低通滤波器的频率响应

figure(2) % 打开窗口2

subplot(221) %图形显示分割窗口

plot(f,abs(h)) %绘制Butterworth低通滤波器的幅频响应图

title(巴氏低通滤波器'')

grid %绘制带网格的图像

sf=filter(a,b,s) %叠加函数S经过低通滤波器以后的新函数

subplot(222)

plot(t,sf) %绘制叠加函数S经过低通滤波器以后的时域图形

xlabel('时间 (seconds)')

ylabel('时间按幅度')

SF=fft(sf,256) %对叠加函数S经过低通滤波器以后的新函数进行256点的基—2快速傅立叶变换

w= %新信号角频率

subplot(223)

plot())%绘制叠加函数S经过低通滤波器以后的频谱图

title('低通滤波后的频谱图')

%设计高通滤波器

[N,Wc]=buttord()

%估算得到Butterworth高通滤波器的最小阶数N和3dB截止频率Wc

[a,b]=butter(N,Wc,'high')%设计Butterworth高通滤波器

[h,f]=freqz()%求数字高通滤波器的频率响应

figure(3)

subplot(221)

plot())%绘制Butterworth高通滤波器的幅频响应图

title('巴氏高通滤波器')

grid%绘制带网格的图像

sf=filter()%叠加函数S经过高通滤波器以后的新函数

subplot(222)

plot(t,sf)%绘制叠加函数S经过高通滤波器以后的时域图形

xlabel('Time(seconds)')

ylabel('Time waveform')

w%新信号角频率

subplot(223)

plot())%绘制叠加函数S经过高通滤波器以后的频谱图

title('高通滤波后的频谱图')

%设计带通滤波器

[N,Wc]=buttord([)

%估算得到Butterworth带通滤波器的最小阶数N和3dB截止频率Wc

[a,b]=butter(N,Wc)%设计Butterworth带通滤波器

[h,f]=freqz()%求数字带通滤波器的频率响应

figure(4)

subplot(221)

plot(f,abs(h))%绘制Butterworth带通滤波器的幅频响应图

title('butter bandpass filter')

grid%绘制带网格的图像

sf=filter(a,b,s)%叠加函数S经过带通滤波器以后的新函数

subplot(222)

plot(t,sf)%绘制叠加函数S经过带通滤波器以后的时域图形

xlabel('Time(seconds)')

ylabel('Time waveform')

SF=fft()%对叠加函数S经过带通滤波器以后的新函数进行256点的基—2快速傅立叶变换

w=( %新信号角频率

subplot(223)

plot('))%绘制叠加函数S经过带通滤波器以后的频谱图

title('带通滤波后的频谱图')

1.for循环处出错2.拟合函数写错

下面是修改过的

fn=1600%采样频率

fp=[80,150]%要求保留的频段

fs=[60,180]

Rp=3 Rs=30

Wp=fp/(fn/2)%计算归一化角频率

Ws=fs/(fn/2)

t=0:1/fn:1

x1=sin(2*pi*300*t) %这里X1是噪声 需要通过带阻滤波滤除

x2=sin(2*pi*120*t)%这里X2是信号，

x=x1+x2

[n,Wn]=buttord(Wp,Ws,Rp,Rs)%计算阶数和截止频率

[b,a]=butter(n,Wn)%计算H(z)分子、分母多项式系数

y2=filtfilt(b,a,x)

for i=1:fn

if (y2(i)<0 &y2(i+1)>0)

x3=[i,i+1,i+2]

y3=[y2(i),y2(i+1),y2(i+2)]

p=polyfit(x3,y3,1)

end

end

这个我刚好做过一个滤波器，事实上对时域信号做FFT，截取一定点数再做逆FFT相当于理想滤波。设计滤波器代码如下：

f1=100f2=200%待滤波正弦信号频率

fs=2000%采样频率

m=(0.3*f1)/(fs/2)%定义过度带宽

M=round(8/m)%定义窗函数的长度

N=M-1%定义滤波器的阶数

b=fir1(N,f2/fs)%使用fir1函数设计滤波器

%输入的参数分别是滤波器的阶数和截止频率

figure(1)

[h,f]=freqz(b,1,512)%滤波器的幅频特性图

%[H,W]=freqz(B,A,N)当N是一个整数时函数返回N点的频率向量和幅频响应向量

plot(f*fs/(2*pi),20*log10(abs(h)))%参数分别是频率与幅值

xlabel('频率/赫兹')ylabel('增益/分贝')title('滤波器的增益响应')

figure(2)

subplot(211)

t=0:1/fs:0.5%定义时间范围和步长

s=sin(2*pi*f1*t)+sin(2*pi*f2*t)%滤波前信号

plot(t,s)%滤波前的信号图像

xlabel('时间/秒')ylabel('幅度')title('信号滤波前时域图')

subplot(212)

Fs=fft(s,512)%将信号变换到频域

AFs=abs(Fs)%信号频域图的幅值

f=(0:255)*fs/512%频率采样

plot(f,AFs(1:256))%滤波前的信号频域图

xlabel('频率/赫兹')ylabel('幅度')title('信号滤波前频域图')

figure(3)

sf=filter(b,1,s)%使用filter函数对信号进行滤波

%参数分别为滤波器系统函数的分子和分母多项式系数向量和待滤波信号输入

subplot(211)

plot(t,sf)%滤波后的信号图像

xlabel('时间/秒')ylabel('幅度')title('信号滤波后时域图')

axis([0.2 0.5 -2 2])%限定图像坐标范围

subplot(212)

Fsf=fft(sf,512)%滤波后的信号频域图

AFsf=abs(Fsf)%信号频域图的幅值

f=(0:255)*fs/512%频率采样

plot(f,AFsf(1:256))%滤波后的信号频域图

xlabel('频率/赫兹')ylabel('幅度')title('信号滤波后频域图')

---