%处理思路及相关matlab代码
%适用于声音的相关检测
%第一步，提取音频
[X,Fs,Wmod,Fidx]=readwav('C:\Users\Emannel\Desktop\matlab\车声.wav');%C:\Users\Emannel\Desktop\matlab\车声.wav
%采样频率Fs
X=X(:,1); %提取一个声道X
X=X-mean(X); % 消除直流分量
X=X/max(abs(X)); % 幅值归一
N=length(X); %获取信号数组长度N
time=(0:N-1)/Fs; %得到每个点对应的时间time
% for i=1:N
% if abs(X(i))<0.05
% X(i)=0;
% end
% end
plot(time,X); %时间信号图线
ylim([-1,1]);
grid on;
%分帧加窗
wlen=1024; %窗长wlen
inc=512; %帧增量inc
y=enframe(X,hamming(wlen),inc)'; %数据按列排布
Fn=size(y,2); %计算帧数fn
frametime=frame2time(Fn,wlen,inc,Fs); %计算每帧对应的时间frametime
%判断含噪无话段
IS=zeros(Fn,1); %定义无话段位置一维变量
Ru=zeros(Fn,1);
for k=1:Fn %计算自相关函数
oneframe=y(:,k);
ru=xcorr(oneframe);
Ru(k)=max(ru);
% oneframe=abs(oneframe);
% if mean(oneframe) > 0.008
% IS(k)=1;
%end
end
Rum=multimidfilter(Ru,10); %中值滤波处理
Rum=Rum/max(Rum);
NIS=50;
thredth=max(Rum(1:NIS)); % 计算阈值
T1=1.1*thredth;
T2=1.3*thredth;
%第二步，端点检测
%开始端点检测
[voiceseg,vsl,SF,NF]=vad_param1D(Rum,T1,T2); % 自相关函数的端点检测
% 作图
figure
subplot 211; plot(time,X,'k');
title('语音波形');
ylabel('幅值'); axis([0 max(time) -1 1]);
% subplot 312; plot(time,signal,'k');
% title(['加噪语音波形(信噪比' num2str(SNR) 'dB)']);
% ylabel('幅值'); axis([0 max(time) -1 1]);
subplot 212; plot(frametime,Rum,'k');
title('短时自相关函数'); axis([0 max(time) 0 1.2]);
xlabel('时间/s'); ylabel('幅值');
line([0,frametime(Fn)], [T1 T1], 'color','k','LineStyle','--');
line([0,frametime(Fn)], [T2 T2], 'color','k','LineStyle','-');
for k=1 : vsl % 标出语音端点
nx1=voiceseg(k).begin; nx2=voiceseg(k).end;
fprintf('%4d %4d %4d\n',k,nx1,nx2);
subplot 211;
line([frametime(nx1) frametime(nx1)],[-1 1],'color','k','LineStyle','-');
line([frametime(nx2) frametime(nx2)],[-1 1],'color','k','LineStyle','--');
end
%语音减噪，LMS自适应法（效果不是很好）
% X=X-mean(X); % 消除直流分量
% X=X/max(abs(X)); % 幅值归一
% M=32;
% mu=0.01;
% %r2=wgn(N,1,0.5); % 产生高斯白噪声
% r2=rand(size(X));
%
% SNR=5;
% b=fir1(31,0.5); % 设计FIR滤波器,代替H
% r21=filter(b,1,r2); % FIR滤波
% [r1,r22]=add_noisedata(X,r21,Fs,Fs,SNR);% 产生带噪语音，信噪比为SNR
%
% b=fir1(31,0.5); % 设计FIR滤波器,代替H
% h = adaptfilt.lms(M,mu); % LMS滤波
% [y,e] = filter(h,r2,X);
% output=e; % LMS滤波输出
% snr=SNR_singlech(X,output); % 计算滤波后的信噪比
% fprintf('snr=%5.4f\n',snr);
% %sound(r1,Fs); % 从声卡发声比较
% %pause(25)
% sound(output,Fs);
% % 作图
% figure
% subplot 211; plot(time,X,'k'); ylabel('幅值')
% ylim([-1 1 ]); title('原始语音信号');
% subplot 212; plot(time,output,'k');
% ylim([-1 1 ]); title('LMS滤波输出语音信号');
% xlabel('时间/s'); ylabel('幅值')
%谱减法（效果一般，带有音乐噪声）
% X=X-mean(X); % 消除直流分量
% X=X/max(abs(X)); % 幅值归一
%
% IS=1; % 设置前导无话段长度
% SNR=5; % 设置信噪比SNR
% signal=Gnoisegen(X,SNR); % 叠加噪声
% snr1=SNR_singlech(X,signal); % 计算初始信噪比
% overlap=wlen-inc; % 求重叠区长度
% NIS=fix((IS*Fs-wlen)/inc +1); % 求前导无话段帧数
% a=4; b=0.001; % 设置参数a和b
% output=simplesubspec(signal,wlen,inc,NIS,a,b);% 谱减
% sound(output);
% snr2=SNR_singlech(X,output); % 计算谱减后的信噪比
% fprintf('snr2=%5.4f \n',snr2);
%
% % 作图
% figure
% plot(time,output,'k');grid;%hold on;
% title('谱减后波形'); ylabel('幅值'); xlabel('时间/s');
%Boll谱减法
%X=X-mean(X); % 消除直流分量
%X=X/max(abs(X)); % 幅值归一
SNR=5; % 设置信噪比
signal=Gnoisegen(X,SNR); % 叠加噪声
IS=1; % 设置IS
output=SSBoll79(X,Fs,IS); % 调用SSBoll79函数做谱减
ol=length(output); % 把output补到与x等长
if ol<N
output=[output; zeros(N-ol,1)];
end
snr=SNR_singlech(X,output); % 计算谱减后的信噪比
fprintf('snr=%5.4f \n',snr);
sound(output);
% 作图
figure;
plot(time,output,'k');grid on; ylim([-1 1]);
title('谱减后波形'); ylabel('幅值'); xlabel('时间/s');
%互相关函数的求解
%原有的信号上添加噪声，并截取其中一段（时间延迟）作为互相关的另一项
for k=1:Fn %计算自相关函数
oneframe=y(:,k);
ru=xcorr(oneframe);
Ru(k)=max(ru);
% oneframe=abs(oneframe);
% if mean(oneframe) > 0.008
% IS(k)=1;
%end
end
% x=x';
% y=[zeros(1,50) x(1:250-50)];
% x1=awgn(x,5,'measured');%加入SNR为零的正态分布白噪声随机序列,'measured'表测量信号噪声.x1=awgn(sx,0,'measured','db' )
% y1=awgn(y,5,'measured');%第二种改法y1=awgn(sy,0,'measured','db' );
% rxx=xcorr(x1,x1,'biased');%求信号x1自相关函数,度应是1024
% ryy=xcorr(y1,y1,'biased');%求信号y1自相关函数
% rxy=xcorr(y1,x1,'biased');%求信号x1,y1互相关函数[rxx,lags]=xcorr(x1,x1,'biased')
% sxx=fft(rxx,2*N-1);%求x1功率谱密度
% syy=fft(ryy,2*N-1);%求y1功率谱密度
% sxy=fft(rxy,2*N-1);%求互相关谱密度,
% cxy=(sxy.*sxy)/(sxx.*syy); %最大似然加权
% Hs=cxy./((abs(sxy)).*(1-cxy)); %最大似然加权
% sxyf=Hs.*sxy;
% TDOA=real(ifft(sxyf,2*N-1));
% TDOA=TDOA/sum(abs(TDOA));
% t1=-249:0.001:249;
% final=interp1(n1,TDOA,t1,'spline');%插值
% c=max(abs(final));
% k=(find(abs(final)==c)-249000)/1000
% m=k/fs %换算出信号的延迟时间
% plot(n1,abs(TDOA)) ;
% strmax=sprintf('(%4.4f,%4.4f)',k,c);
% text(k+5.2,c,strmax)
% ylabel('幅度');
% xlabel('采样点');
% title('广义互相关(ml)');
% grid on
评论0