%?wavread（filename）?讀取一個WAV文件，并返回采樣數據到向量y中，Fs表示采樣頻率?bits表示采樣位數??
Fs=44100;bits=16;
[y?Fs?bits]?=?wavread（‘yinpin.wav‘）;????
%假設聲音文件有兩個聲道，我們只分析第一個聲道，如果要分析第二個聲道可以改為：ft=y（:2）;??
ft?=?y（:1）;??
sigLengthy?=?length（ft）;?%獲取聲音長度??
sigLengthy;
x=ft（1:7938000）;
sigLengthx=length（x）;%??截取三分鐘的音頻點數，180*44100=7938000

%可以使用sound函數來播放聲音??
%sound（x?Fs?bits）????
%定義一個FFT點數
MM=8388608;%2的23次方
%sound（x?Fs?bits）?
%繪制波形圖
?t=（0:sigLengthx-1）/Fs;???
??figure（1）;???
??plot（t?x）;??
??title（‘波形圖‘）
??xlabel（‘時間（s）‘）;??
??ylabel（‘幅度‘）;?
?%幅頻曲線
?xFFT=fft（xMM）;
?XFFT=abs（xFFT）;%取模長
??figure（2）;
?n=0:MM-1;
?f=n*Fs/MM;
??plot（fXFFT）;
??xlim（[0?22050]）;?
??title（‘音頻頻譜‘）
??xlabel（‘頻率（Hz）‘）;??
??ylabel（‘幅度‘）;?
%兩種方法產生高斯白噪聲
%?m=0:?sigLengthx-1;
%?sigma=mean（x.*x）;%噪聲方差
%?noise?=?sqrt（sigma）*randn（[sigLengthx1]）;%噪聲生成
%?AWGNx?=?x+noise;%受寬帶噪聲干擾的音頻信號
??
?%加入?高斯白噪聲
??AWGNx=awgn（x15）;
%???sound（AWGNx?Fs?bits）?
??AWGNxFFT=fft（AWGNxMM）;
?AWGNXFFT=abs（AWGNxFFT）;
?figure（3）;???
?plot（t?AWGNx）;??
?title（‘加入高斯白噪聲之后波形圖‘）
?xlabel（‘時間（s）‘）;??
?ylabel（‘幅度‘）;??
?figure（4）;
?plot（fAWGNXFFT）;
?xlim（[0?22050]）;?
?title（‘加入高斯白噪聲之后音頻頻譜‘）
?xlabel（‘頻率（Hz）‘）;???
?ylabel（‘幅度‘）;


%加入高頻噪聲
%設置余弦10000Hz?噪聲函數
?fn=10000;
?noise=0.05*sin（2*pi*fn/Fs.*（1:sigLengthx）‘）;
?highnoise=noise+x;
?highnoiseFFT=fft（highnoiseMM）;
?HighnoiseFFT=abs（highnoiseFFT）;
%??
?figure（5）;???
?plot（t?highnoise）;??
?title（‘加入高頻噪聲之后波形圖‘）
?xlabel（‘時間（s）‘）;??
?ylabel（‘幅度‘）;
??
?figure（6）;
?plot（fHighnoiseFFT）;
?xlim（[0?22050]）;?
?title（‘加入高頻噪聲之后音頻頻譜‘）
?xlabel（‘頻率（Hz）‘）;??
?ylabel（‘幅度‘）;
%??sound（?highnoise?Fs?bits）?


%設置低通濾波器
%設置濾波器的性能參數?通帶邊界為4000Hz??阻帶邊界4500Hz???
%通帶波紋0.01?阻帶波紋0.01
c=firpmord（[4000?4500][10][0.01?0.001]Fs‘cell‘）;
b=firpm（c{:}）‘;
Bfft=fft（bMM）;
len