Matlab实现音频降噪

可用于音频信号处理的学习，兼容Octave！

%*****参考文献：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‘）;
% 
%********************传统信噪比 ***************
%  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