%*****參考文獻：EMD-based?Filtering?（EMDF）?of?Low-Frequency?Noise?for?Speech
clear?all;
close?all;
clc;
fs=16000;
frame_Len=512;
Shift_Len=128;
fid?=?fopen?（‘SI923.WAV‘?‘r‘）;%?Open?file?for?reading?
if?fid?==-?1
errordlg?（?‘Specify?_valid?file?name!‘?‘ERROR:?Invalid?file?name‘）;?
else?
fseek?（fid?1024?-1）;
Input_data0?=?fread?（fid?inf?‘int16‘）;
fclose?（fid）;?
end;?
?
?figure;
?plot（Input_data0‘y‘）;
?title（‘input?Data??‘）?;

Input_data=Input_data0-mean（Input_data0）;
inputnoise=100*（rand（length（Input_data0）1）-0.5）;
Input_data=Input_data+?inputnoise;

%?figure;
%?plot（inputnoise.^2）;
figure;
plot（Input_data）;
title（‘Data?of?Adding?Noise‘）?;
%?d0=round（Input_data）;
%?d_t=int16（d0）;
%?wavplay（d_tfs）;
%?wavwrite（d_tfs16‘includenoise.wav‘）;
%?
%********************傳統(tǒng)信噪比?***************
%??SNR=10*log10（（sum（Input_data0.^2）/sum（Input_data.^2）））
%***********************************************

Input_data=Input_data‘;
DataLen=length（Input_data）;????????????????%輸入信號長度
frame_Num=fix（DataLen/Shift_Len）;????????????????%計算需要處理的幀數??
Win_hm=hamming（frame_Len）‘;????????????%漢明窗


%*******************P_speech_absence****************************
%*********參考文獻：Noise?Spectrum?Estimation?in?Adverse?Environments:Improved?Minima
%Controlled?Recursive?Averaging********
?%定義部分
%??guji_noise=[];
Noise_power（1:frame_Len）=0;??%?噪聲估計
Sf_first=zeros（1frame_Len）;????????%?Sf
Smin_first=zeros（1frame_Len）;???????????%第一次更新最小功率頻譜值Smin
S_first=zeros（1frame_Len）;??????????????%第一次回歸平滑?S
Smin_second=zeros（1frame_Len）;?????%第二次更新最小功率譜??Smin^
S_second=zeros（1frame_Len）;????????????%第二次回歸平滑?S^
Sf_second=zeros（1frame_Len）;??????%Sf^
I=zeros（1frame_Len）;??????????????????%初步估計語音存在概率??I（kl）
D_window_first=1e20*ones（120?frame_Len）;???%?求S最小值搜索窗?D
D_window_second=1e20*ones（120?frame_Len）;%?求S^最小值搜索窗?D

?frame_Data（1:frame_Len）=Input_data（1:frame_Len）.*Win_hm;??%第一幀
?Fre_Data（1:frame_Len）=fft（frame_Dataframe_Len）;??%變換到頻域
?Data_Amp=abs（Fre_Data）;
?
?Noise_power=Data_Amp.^2;???%初始化Noise_power
?%******************用前10幀估算噪聲，計算輸入分段信噪比*************
?Noise_power0=0;
?SEGSNR_in=0;
?SEGSNR_out=0;
?for?j=1:5
??frame_Data（1:frame_Len）=Input_data（（j-1）*Shift_Len+1:（j-1）*Shift_Len+frame_Len）;
?Noise_power0=Noise_power0+?sum（frame_Data.^2）;
?end?
?Noise_power0=?Noise_power0/5;?
?%******************************************************************
%??guji_noise=[guji_noiseNoise_power];
%??figure;
%??plot（?Noise_power）;
?Sf_first=Data_Amp（1）^2;??????%計算Sf（k0）
??for?k=2:frame_Len-1
?????Sf_first（k）=0.25*Data_Amp（k-1）^2+0.5*Data_Amp（k）^2+0.25*Data_Amp（k+1）^2;?
??end
??Sf_first（frame_Len）=Data_Amp（frame_Len）^2;

%初始化
Smin_first=Sf_first;????????%?Smin（k0）?=?Sf（k0）
S_first=Sf_first;?????????????%?S（k0）=Sf（k0）
Smin_second=Sf_first;???????????%Smin（k0）^=Sf（k0）
S_second=Sf_first;????????????%S（k0）^=Sf（k0）
Sf_second=Sf_first;??????????%Sf（k0）^=Sf（k0）
D_window_first（1:）=Sf_first;?%搜索窗內第一個值?為Sf（k

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----
?????目錄???????????0??2016-09-26?22:09??IMCRA\
?????文件????????8394??2013-08-26?11:36??IMCRA\IMCRA.m
?????文件???????95848??2013-01-15?11:23??IMCRA\SI923.WAV