[x,fs11]=audioread('I:\旧数据\xinhao\ke she\第五\shi miao.mp3');
fs=1000;
%sound(x,fs11); %声音
%x=x(:,1); %单声道
Fs=fs; %频率
N=length(x); %音乐的长度
n=0:N-1;
f=n*fs/N; %频率序列
tt=(N-1)/fs;
ts=1/fs;
t=0:ts:tt;
X=fft(x,N);
figure(1);
subplot(211);plot(t,x);title('音频时域波形');
xlabel('t');ylabel('x');
subplot(212);plot(f,X);title('音频频域波形');
xlabel('f/hz');ylabel('X');
%设计噪声
nt1=cos(2*pi*300*n/fs); %第一种噪声
%nt1=nt.';
nt2=2*randn(1,N)-1; %第二种噪声
nt=nt1+nt2;
nt=nt.';
%fp=23000;f_s=18000;
fp=400;f_s=390;
wp=fp*2/fs; ws=f_s*2/fs
Ap=3;As=18;
[N1,wn]=buttord(wp,ws,Ap,As);
[b,a]=butter(N1,wn,'high');
[h,w]=freqz(b,a);
ntt=filter(b,a,nt);
%ntt=ntt.';
Ntt=fft(ntt,N);
figure(2);
plot(f,Ntt);title('噪声频域波形');
xlabel('f/hz');ylabel('Ntt');
%加噪过程
x_ntt=x + ntt; %音频+噪声
%x_ntt=x_ntt.';
%sound(x_ntt,fs11);
x_ntt_fft=fft(x_ntt,N);
figure(3);
subplot(211);plot(t,x_ntt);title('噪声时域波形');
xlabel('t');ylabel('x-noise');
subplot(212);plot(f,x_ntt_fft);title('加噪音频频域波形');
xlabel('f/hz');ylabel('x-noise-fft');
%设置滤波器
Ap=3;As=18;
fp1=300;f_s1=350;
%fp1=15000;f_s1=17000; %观察加噪音频频带和音频频带得到滤波器阻带通带的截止频率
wp1=fp1*2/fs; ws1=f_s1*2/fs
Ap=3;As=18;
%[Nn,w2]=cheb1ord(wp1,ws1,Ap,As); %切比雪夫滤波器
%[b1,a1]=cheby1(Nn,Ap,w2);
%[Nn,w2]=ellipord(wp1,ws1,Ap,As); %椭圆滤波器
%[b1,a1]=ellip(Nn,Ap,As,w2);
[Nn,w2]=buttord(wp1,ws1,Ap,As); %巴特沃斯滤波器
[b1,a1]=butter(Nn,w2);
[H1,w3]=freqz(b1,a1);
figure(4);
plot(w/pi,abs(H1));title('滤波器幅频图');
xlabel('w/pi');ylabel('H1');
%进行滤波
x2=filter(b1,a1,x_ntt);
N2=length(x2);
n2=0:N2;
f2=n*fs/N2;
X2=fft(x2,N2);
sound(x2,fs11);
figure(5);
subplot(211);plot(t,x2);title('滤波后的音频时域图');
xlabel('f/hz');ylabel('x');
subplot(212);plot(f2,X2);title('滤波后的音频频域图');
xlabel('f/hz');ylabel('X');
|