濾波器代碼
❶ 設計一個matlab帶通濾波器代碼
你自己整合吧,我沒時間幫你整合,我給你提供一些程序:
絕對正確的代碼:程序1:
fs=22050;
%語音信號采樣頻率為22050
x1=wavread('windows
critical
stop.wav');
%讀取語音信號的數據,賦給變數x1
sound(x1,22050);
%播放語音信號
y1=fft(x1,1024);
%對信號做1024點fft變換
f=fs*(0:511)/1024;
figure(1)
plot(x1)
%做原始語音信號的時域圖形
title('原始語音信號');
xlabel('time
n');
ylabel('fu
n');
figure(2)
freqz(x1)
%繪制原始語音信號的頻率響應圖
title('頻率響應圖')
figure(3)
subplot(2,1,1);
plot(abs(y1(1:512)))
%做原始語音信號的fft頻譜圖
title('原始語音信號fft頻譜')
subplot(2,1,2);
plot(f,abs(y1(1:512)));
title('原始語音信號頻譜')
xlabel('hz');
ylabel('fu');
程序2:
fs=22050;
%語音信號采樣頻率為22050
x1=wavread('windows
critical
stop.wav');
%讀取語音信號的數據,賦給變數x1
t=0:1/22050:(size(x1)-1)/22050;
y1=fft(x1,1024);
%對信號做1024點fft變換
f=fs*(0:511)/1024;
x2=randn(1,length(x1));
%產生一與x長度一致的隨機信號
sound(x2,22050);
figure(1)
plot(x2)
%做原始語音信號的時域圖形
title('高斯隨機雜訊');
xlabel('time
n');
ylabel('fu
n');
randn('state',0);
m=randn(size(x1));
x2=0.1*m+x1;
sound(x2,22050);%播放加雜訊後的語音信號
y2=fft(x2,1024);
figure(2)
plot(t,x2)
title('加噪後的語音信號');
xlabel('time
n');
ylabel('fu
n');
figure(3)
subplot(2,1,1);
plot(f,abs(y2(1:512)));
title('原始語音信號頻譜');
xlabel('hz');
ylabel('fu');
subplot(2,1,2);
plot(f,abs(y2(1:512)));
title('加噪後的語音信號頻譜');
xlabel('hz');
ylabel('fu');
根據以上代碼,你可以修改下面有錯誤的代碼
程序3:雙線性變換法設計butterworth濾波器
fs=22050;
x1=wavread('h:\課程設計2\shuzi.wav');
t=0:1/22050:(size(x1)-1)/22050;
au=0.03;
d=[au*cos(2*pi*5000*t)]';
x2=x1+d;
wp=0.25*pi;
ws=0.3*pi;
rp=1;
rs=15;
fs=22050;
ts=1/fs;
wp1=2/ts*tan(wp/2);
%將模擬指標轉換成數字指標
ws1=2/ts*tan(ws/2);
[n,wn]=buttord(wp1,ws1,rp,rs,'s');
%選擇濾波器的最小階數
[z,p,k]=buttap(n);
%創建butterworth模擬濾波器
[bap,aap]=zp2tf(z,p,k);
[b,a]=lp2lp(bap,aap,wn);
[bz,az]=bilinear(b,a,fs);
%用雙線性變換法實現模擬濾波器到數字濾波器的轉換
[h,w]=freqz(bz,az);
%繪制頻率響應曲線
figure(1)
plot(w*fs/(2*pi),abs(h))
grid
xlabel('頻率/hz')
ylabel('頻率響應幅度')
title('butterworth')
f1=filter(bz,az,x2);
figure(2)
subplot(2,1,1)
plot(t,x2)
%畫出濾波前的時域圖
title('濾波前的時域波形');
subplot(2,1,2)
plot(t,f1);
%畫出濾波後的時域圖
title('濾波後的時域波形');
sound(f1,22050);
%播放濾波後的信號
f0=fft(f1,1024);
f=fs*(0:511)/1024;
figure(3)
y2=fft(x2,1024);
subplot(2,1,1);
plot(f,abs(y2(1:512)));
%畫出濾波前的頻譜圖
title('濾波前的頻譜')
xlabel('hz');
ylabel('fu');
subplot(2,1,2)
f1=plot(f,abs(f0(1:512)));
%畫出濾波後的頻譜圖
title('濾波後的頻譜')
xlabel('hz');
ylabel('fu');
程序4:窗函數法設計濾波器:
fs=22050;
x1=wavread('h:\課程設計2\shuzi.wav');
t=0:1/22050:(size(x1)-1)/22050;
au=0.03;
d=[au*cos(2*pi*5000*t)]';
x2=x1+d;
wp=0.25*pi;
ws=0.3*pi;
wdelta=ws-wp;
n=ceil(6.6*pi/wdelta);
%取整
wn=(0.2+0.3)*pi/2;
b=fir1(n,wn/pi,hamming(n+1));
%選擇窗函數,並歸一化截止頻率
figure(1)
freqz(b,1,512)
f2=filter(bz,az,x2)
figure(2)
subplot(2,1,1)
plot(t,x2)
title('濾波前的時域波形');
subplot(2,1,2)
plot(t,f2);
title('濾波後的時域波形');
sound(f2,22050);
%播放濾波後的語音信號
f0=fft(f2,1024);
f=fs*(0:511)/1024;
figure(3)
y2=fft(x2,1024);
subplot(2,1,1);
plot(f,abs(y2(1:512)));
title('濾波前的頻譜')
xlabel('hz');
ylabel('fu');
subplot(2,1,2)
f2=plot(f,abs(f0(1:512)));
title('濾波後的頻譜')
xlabel('hz');
ylabel('fu');
❷ 帶通濾波器matlab代碼
% 用切比雪夫最佳一致逼近設計線性相位FIR帶通濾波器;
%信號為0.5hz, 0.9hz, 1.1hz和1.5hz的正統信號疊加組成
%通帶為[0.9,1.1]
%頻譜解析度與信號實際長度N成正比
clear all;
f1=0.5;f2=0.9;f3=1.1;f4=1.5;t=0:1203;N=length(t);fs=10;M=512;
x1=sin(2*pi*(f1/fs)*t)+sin(2*pi*(f2/fs)*t)+sin(2*pi*(f3/fs)*t)+sin(2*pi*(f4/fs)*t);
figure(1);
subplot(211);plot(t,x1);title('原信號');
y=fft(x1);
f=(0:1/N:1/2-1/N)*fs;
subplot(212);plot(f,abs(y(1:N/2)));grid;xlabel('hz');%處理前頻譜
wc1=2*f2/fs;wc2=2*f3/fs;wc3=2*f4/fs;%歸一化角頻率,用於下面的f1
f1=[0 wc1-0.05 wc1 wc2 wc2+0.05 1];
A=[0 0 1 1 0 0];%設置帶通或帶阻,1為帶通,0為帶阻
weigh=[1 1 1 ];%設置通帶和阻帶的權重
b=remez(60,f1,A,weigh);%傳函分子
h1=freqz(b,1,M);%幅頻特性
figure(2)
f=(0:1/M:1-1/M)*fs/2;
subplot(211);plot(f,abs(h1));grid;title('帶通');
x2=filter(b,1,x1);
S1=fft(x2);
f=(0:1/N:1/2-1/N)*fs;
subplot(212);plot(f,abs(S1(1:N/2)));grid;xlabel('hz');%處理後頻譜
❸ 數字切比雪夫濾波器MATLAB代碼急求
用Matlab設計切比雪夫濾波器是很簡單的,比如設計I型的
N = 8;%階數
wp=0.2;%歸一化截止頻率
Rp= 1;%通帶波紋
[b a]=cheby1(N,Rp,wp,'low');%設計低通濾波器
freqz(b,a)%頻率響應曲線
❹ 一段matlab低通濾波器程序,求改編成C語言。
這個我剛好做過一個濾波器,事實上對時域信號做FFT,截取一定點數再做逆FFT相當於理想濾波。設計濾波器代碼如下:
f1=100;f2=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('信號濾波後頻域圖');
❺ 求:一個關於FIR帶通濾波器的C語言設計程序 代碼
short h[], short y[])
{
int i, j, sum; for (j = 0; j < 100; j++) {
sum = 0;
for (i = 0; i < 32; i++)
sum += x[i+j] * h[i];
y[j] = sum >> 15;
}
}
2
void fir(short x[], short h[], short y[])
{
int i, j, sum0, sum1;
short x0,x1,h0,h1; for (j = 0; j < 100; j+=2) {
sum0 = 0;
sum1 = 0;
x0 = x[j];
for (i = 0; i < 32; i+=2){
x1 = x[j+i+1];
h0 = h[i];
sum0 += x0 * h0;
sum1 += x1 * h0;
x0 = x[j+i+2];
h1 = h[i+1];
sum0 += x1 * h1;
sum1 += x0 * h1;
}
y[j] = sum0 >> 15;
y[j+1] = sum1 >> 15;
}
}
3
void fir(short x[], short h[], short y[])
{
int i, j, sum0, sum1;
short x0,x1,x2,x3,x4,x5,x6,x7,h0,h1,h2,h3,h4,h5,h6,h7; for (j = 0; j < 100; j+=2) {
sum0 = 0;
sum1 = 0;
x0 = x[j];
for (i = 0; i < 32; i+=8){
x1 = x[j+i+1];
h0 = h[i];
sum0 += x0 * h0;
sum1 += x1 * h0;
x2 = x[j+i+2];
h1 = h[i+1];
sum0 += x1 * h1;
sum1 += x2 * h1;
x3 = x[j+i+3];
h2 = h[i+2];
sum0 += x2 * h2;
sum1 += x3 * h2;
x4 = x[j+i+4];
h3 = h[i+3];
sum0 += x3 * h3;
sum1 += x4 * h3;
x5 = x[j+i+5];
h4 = h[i+4];
sum0 += x4 * h4;
sum1 += x5 * h4;
x6 = x[j+i+6];
h5 = h[i+5];
sum0 += x5 * h5;
sum1 += x6 * h5;
x7 = x[j+i+7];
h6 = h[i+6];
sum0 += x6 * h6;
sum1 += x7 * h6;
x0 = x[j+i+8];
h7 = h[i+7];
sum0 += x7 * h7;
sum1 += x0 * h7;
}
y[j] = sum0 >> 15;
y[j+1] = sum1 >> 15;
}
}
❻ 濾波器matlab代碼
貌似此類專家很少哦~~~無解
❼ MATLAB基於漢寧窗的FIR的低通濾波器的源代碼及注釋
很常見的設計題目
給你一個常式,只需要改一改參數就行了
clear all;
f=[0 0.19 0.2 0.3 0.31 0.59 0.6 0.8 0.81 1];
% 給定頻率軸分點;
m=[0 0 1 1 0 0 1 1 0 0];
% 給定在這些頻率分點上理想的幅頻響應
N1=30;
N2=90;
% 取兩種不同的濾波器長度;
b1=fir2(N1,f,m);
b2=fir2(N2,f,m);
% 得到兩個濾波器;
subplot(311);
stem(b1,'.');grid;
subplot(312);
stem(b2,'.');grid;
M=128;
[h1,w]=freqz(b1,1,M,1);
[h2,w]=freqz(b2,1,M,1);
subplot(313);
plot(w,abs(h1),'b-',w,abs(h2),'g-');grid;
其中,f是歸依化以後的頻率 通過數字濾波器的采樣頻率算出來,根據通帶和阻帶算好f和m就行了
看一看help,這個函數應該有窗函數的選擇 默認情況下是漢明窗
希望能夠幫到你
❽ 幫忙編一個matlab的低通濾波器的程序
clc;clearall;
%歸一化模擬切比雪夫I型低通濾波器的設計
Wp=2*pi*1000;Ws=2*pi*1500;rp=3;rs=30;%設計濾波器的參數
wp=1;ws=Ws/Wp;%頻帶變換得到歸一化濾波器
[N,wc]=cheb1ord(wp,ws,rp,rs,'s');%計算濾波器階數和3dB截止頻率
[z,p,k]=cheb1ap(N,rp);%計算歸一化濾波器的零極點
[b,a]=zp2tf(z,p,k);%計算歸一化濾波器系統函數的系數
w0=0:0.05*pi:2*pi;
[h0,w0]=freqs(b,a,w0);%求頻率響應
figure(1)
plot(w0,20*log10(abs(h0)),'k');
title('歸一化模擬且比雪夫I型低通濾波器');
xlabel('頻率f/Hz');ylabel('幅度/dB');grid;
%一般低通切比雪夫低通濾波器的設計
[B,A]=lp2lp(b,a,Wp);%將歸一化濾波器轉換為一般模擬濾波器
w1=0:2*pi*100:2*pi*30000;
[h1,w1]=freqs(B,A,w1);
figure(2)
plot(w1/(2*pi),20*log10(abs(h1)),'k');
title('一般模擬且比雪夫I型低通濾波器');
xlabel('頻率f/Hz');ylabel('幅度/dB');grid;
%沖激響應不變法設計低通巴特沃斯數字濾波器
Fs=10000;%采樣頻率
Wp1=Wp/Fs;Ws1=Ws/Fs;
rp1=3;rs1=30;%設計濾波器的參數
[N1,Wn]=cheb1ord(Wp,Ws,rp1,rs1,'s')%計算濾波器階數
[b1,a1]=cheby1(N1,rp1,Wn,'s');%樣本AF的系數函數的分子分母系數
[bz,az]=impinvar(b1,a1,Fs);%沖激響應不變法從AF到DF變換
[C1,B1,A1]=dir2par(bz,az)%直接形式轉換為並聯型
w2=[Wp1,Ws1];%數字臨界頻率
[H,f]=freqz(bz,az);%繪制數字濾波器的幅頻特性和相頻特性
[db,mag,pha,grd,f]=freqz_m(bz,az);%擴展函數檢驗濾波器的其他特性
figure(3)
plot(f/pi,db,'k');
title('沖激響應不變法設計低通切比雪夫I型數字濾波器');
xlabel('角頻率w/pi');ylabel('振幅/dB');
axis([0,0.35,-30,5]);grid;
%用設計好的濾波器對信號進濾波處理
figure(4)
f1=500;f2=4000;%輸入信號的頻率
N=100;%數據長度
dt=1/Fs;n=0:N-1;t=n*dt;%采樣間隔和時間序列
x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t);%輸入信號
subplot(2,1,1),plot(t,x),title('輸入信號');
y=filtfilt(bz,az,x);%用濾波器進行濾波處理
y1=filter(bz,az,x);%進行濾波處理
subplot(2,1,2),plot(t,y,t,y1,':'),title('輸出信號');xlabel('時間/s');
legend('filtfilt','filter')%加圖例
freqz_m.m文件
function[db,mag,pha,grd,w]=freqz_m(b,a)
[H,w]=freqz(b,a,1000,'whole');
mag=abs(H);
db=20*log10((mag+eps)/max(mag));
pha=angle(H);
grd=grpdelay(pha);
dir2par.m文件
function[C,B,A]=dir2par(b,a)
M=length(b);N=length(a);
[r1,p1,C]=resiez(b,a);
p=cplxpair(p1,10000000*eps);
I=cplxcomp(p1,p);
r=r1(I);
K=floor(N/2);B=zeros(K,2);A=zeros(K,3);
ifK*2==N;
fori=1:2:N-2
Brow=r(i:1:i+1,:);
Arow=p(i:1:i+1,:);
[Brow,Arow]=resiez(Brow,Arow,[]);
B(fix((i+1)/2),:)=real(Brow);
A(fix((i+1)/2),:)=real(Arow);
end
[Brow,Arow]=resiez(r(N-1),p(N-1),[]);
B(K,:)=[real(Brow)0];A(K,:)=[real(Arow)0];
else
fori=1:2:N-1
Brow=r(i:1:i+1,:);
Arow=p(i:1:i+1,:);
[Brow,Arow]=resiez(Brow,Arow,[]);
B(fix((i+1)/2),:)=real(Brow);
A(fix((i+1)/2),:)=real(Arow);
end
end
cplxcomp.m文件
functionI=cplxcomp(p1,p2)
I=[];
forj=1:1:length(p2)
fori=1:1:length(p1)
if(abs(p1(i)-p2(j))<0.0001)
I=[I,i];
end
end
end
I=I';
❾ 急,幫忙看一段關於濾波器的代碼(2階的低通數字濾波) 簡要把關鍵點說一下就行了!
首先,使用temp=(float)(wave)/16384 將輸入的整數數據轉換成浮點的數據。
使用浮點格式的數據運算後,輸出再使用 wave=(int)(temp*16384);轉換成整數。
先處以16384,是為了防止溢出;運算之後必須還原的~~