电话按键音的识别

电话按键音的识别
————————————————————————————————作者：————————————————————————————————日期：
一、 实验目的
1. 本实验的内容基于对电话通信系统中拨号音的识别仿真实验。

主要涉及到电话
拨号音识别的基本原理和识别的主要方法。

利用matlab 软件以及fft 算法实现对电话通信系统中拨号音的识别。

并进一步利用matlab 中的图形用户界面制作简单直观的模拟界面，使其对拨号音的识别有个基本的了解。

2. 能够利用矩阵不同的基频合成0-9不同按键的拨号音，并能够对不同的拨号音
加以正确的识别。

进一步画出简单的图形，从而实现对电话拨号音系统的简单的实验仿真。

3.了解学习DTMF （双音多频）相关知识，知道双音多频的信号是用两个特定的单音频率信号的组合来代表数字或功能。

二、实验仪器
安装有MATLAB 软件的计算机一台。

三、实验原理
1. DTMF （双音多频） （1） 基本概念
电话拨号有两种，脉冲和音频，所谓音频也成双音多频（DTMF ）信号的拨号方式即是电话拨号时每按一个键，有两个音频频率叠加成一个双音频信号，十六个按键由八个音频频率区分。

（2） 编制规则
具体DTMF 编制规则如表1所示
表1 DTMF 键盘
高群/Hz
低群/Hz 1209 3 697 1 2 3 A 770 4 5 6 B 852
7
8
9
C
941 * 0 # D
双音多频的拨号键盘是4×4的矩阵，每一行代表一个低频，每一列代表一个高频。

用户每按一个键就发送一个高频和低频的正弦信号组合。

2. 电话按键音识别
MATLAB 中audioread 函数可以实现多种音频信号的读取，调用形式为
[y,fs]=audioread()
其中y 为所读取的音频数据； 而fs 为采样频率；
调用的为指定载入的音频文件名称。

用fft 频谱分析公式
/2/w T fs f fs π=Ω=Ω= (1)
2/w N k π=⋅ (2)
由公式（1）和（2）可得：
(3)
四．实验内容
本实验要实现电话拨号音(DTMF)信号的检测的识别，可以通过直接计算傅里叶变换时的输入信号组成的频率。

这里采用FFT 算法对信号进行解码分析。

首先对接收到的数字信号作FFT 分析。

计算出其频幅谱，继而得到功率谱，组成输入信号的频率必定对应功率谱的峰值。

对于连续的双音多频(DTMF)信号。

需要把有效的数字拨号信号从静音间隔中分割提取出来，然后再用FFT 算法对信号进行解码分析，得出电话拨号码。

五、实验过程与分析
1. 读取音频文件
利用MATLAB 中audioread 读取
[x,fs]=audioread('电话按键音.m4a'); 2. 声道提取
电话按键音是双声道，本实验提取声道1进行分析： x=x(:,1)； 3. 端点提取
该按键音为一段连续的11位数字拨号音，分析时需要对它们分别进行处理，而这11位数字间的时间间隔一般情况是没有规律的，无法直接用for 循环控制截取，而要是手动一段一段地进行截取，换一个音频信号就得再次重新截取，因此考虑到可能存在噪音的情况，基于短时能量和短时平均过零率来进行端点检测，实现端点的提取，对于不同的音频信号，只需要修改部分可变参数便可以适用。

4. 分析
根据音频实际情况设置对应参数：
framelen= floor(fs*40/1000);%floor 向下取整，帧长 frameinc= floor(fs*10/1000);%帧移 % 进行分帧，每帧长framelength ，
voice_min_len =15;% 最短语音长度150ms
%语音段的最短长度，若语音段长度小于此值（如果语音段中的静音帧数未超过此值，），则认为其为一段噪音
unvoice_min_len = 5;%结束段最小持续50ms n3=0.05;%过零率下限参数0.05 检测结果见图1
结果： 将11个检测到的号码音频左右端点提取出来
对信号进行fft 变换后，得到的频谱是呈现对称的，在此截取前部分进行观察，得到频谱图如图2所示。

程序运行最终结果为：
The telephone number is:
/f fs k N =*
图1 端点检测结果图 图2 电话按键音fft 频谱图 5. 可靠性检验
自己录制一段拨号音，检验程序的可靠性，这里为了 为了验证程序的正确、可靠性，又分别录制了我自己手机号以及1234567890两段不同格式不同号码，在有噪声条件下的按键音，检验过程如下 （1）(.mp3格式) 参数设置同上
图3 端点检测结果图 图4 电话按键音的fft 频谱图 程序运行结果：
The telephone number is:
图5 端点检测结果图 图6 电话按键音的fft 频谱图
（2）1234567890(.ogg 格式)
该音频录制时，有的按键音持续时间较短，需要对检测初始参数进行修改
0.51 1.5
2 2.5
3 3.54x 10
5
-1
01
s p e e c h
样本点
100
200
300
400500
600
700
800
100200
e n e r g y
帧数
100200300
400500600700800
5z c r
帧数
1234567x 10
5
-1
01
s p e e c h
样本点
200
400
600
8001000
1200
1400
1600
200400
e n e r g y
帧数
200400600
8001000120014001600
50
z c r
帧数
0.51 1.52
2.53
3.54
4.55
5.5x 10
5
-1
01
s p e e c h
样本点
200400
60080010001200
100e n e r g y
帧数
200
400
600800
1000
1200
2040
z c r
帧数
6008001000120014001600
2
46
x 10
-4
fft 频谱图
频率/Hz 幅度
6008001000120014001600
2
46x 10
-4
fft 频谱图
频率/Hz 幅度
6008001000120014001600
1
234x 10
-4
fft 频谱图
频率/Hz 幅度
6008001000120014001600
1
23
x 10
-4
fft 频谱图
频率/Hz 幅度
6008001000120014001600
1
23
x 10-4
fft 频谱图频率/Hz 幅度
6008001000120014001600
0.5
11.5
2x 10
-4
fft 频谱图频率/Hz 幅度
6008001000120014001600
0.5
11.52
2.5x 10-4
fft 频谱图频率/Hz
幅度
6008001000120014001600
1
23x 10-4
fft 频谱图频率/Hz
幅度
6008001000120014001600
0.5
11.52
x 10-4
fft 频谱图频率/Hz
幅度
6008001000120014001600
0.5
11.52
2.5x 10-4
fft 频谱图频率/Hz
幅度
修改参数：
voice_min_len =6;% 最短语音长度60ms
n3=0.02;%过零率下限参数0.02
程序运行结果：
The telephone number is: 1234567890
经验证，该程序可行。

六、小结
1.原本想着，进行端点检测时考虑了噪声的影响进行滤波或许就不必要了，但是将滤波器加进程序后发现，对于不同的音频信号，程序有更好的适应性，端点检测结果的精确度居然也有所提高，不用再反复修改设置参数，然后就将滤波器加入程中了。

2.实验过程中，首先是手动截取11段数字音频验证号码可检测出来后，在进行端点检测的构思，最后考虑滤波。

3.通过本次实验，发现自己在MATLAB函数的应用与程序设计中仍有很大的进步空间。

4.有效的端点检测技术不仅能在语音识别系统中减少数据的采集量，节约处理时间，还能排除无声段或噪声段的干扰。

对于端点检测的程序，其原理经查阅资料明白一些，但不多，查阅到语音端点检测参考程序之后，对其进行了步骤增加、减少等操作，使其适合本实验，考虑到一些参数对于不同的音频需要做不同的改变，便将其提取到主程序中，方便修改。

5.在引用的端点检测代码原理不很了解情况下，结合本实验对其进行部分修改比较困难。

6.对于待检测信号的第一个数字滤波效果不好，未达到理想的效果。

7.先设计滤波器对信号进行滤波，滤去部分噪声，再设计考虑噪声的端点检测函数，相当于信号进行了双重滤波处理，是结果能够更加准确。

七、实验程序
（1）主程序：
[x,fs]=audioread('电话按键音.m4a');%待检测数据
%[x,fs]=audioread('10月11日.mp3');%检验号码
%[x,fs]=audioread('11_11_1.ogg');%检验号码1234567890
x=x(:,1);%提取声道1
x=lv_bo_1(x,fs);
sound(x,fs);%播放滤波后的音频
%% 调整各参数
framelen= floor(fs*40/1000);%floor向下取整，帧长
frameinc= floor(fs*10/1000);%帧移
% 进行分帧，每帧长framelength，
voice_min_len =15;% 最短语音长度60ms ,150ms
%语音段的最短长度，若语音段长度小于此值（如果语音段中的静音帧数未超过此值，），则认为其为一段噪音
unvoice_min_len = 5;%结束段最小持续50ms
n3=0.02;%过零率下限参数0.02，0.05
%% end
[fan_wei]=vad(framelen,frameinc,x,fs,voice_min_len,unvoice_min_len,n3);%引用时注意修改函数vad中的最短语音长度
number=' ';
figure
for i=1:length(fan_wei)
y=x(fan_wei(i,1):fan_wei(i,2));
n=0:length(y)-1;%建立一个信号等长的序列;
%=================频域图像==============================
N=length(y);%取信号矩阵的长度
FFT=fft(y,N);%N点傅里叶变换
mag=abs(FFT)/(N/2);%还原真实幅值
f=n*fs/N;%频率序列
%
f=f(1:fix(N/2));mag=mag(1:fix(N/2));%
subplot(3,4,i)%图2：频谱图
%取l/2作图
plot(f,mag);axis([600 1700 min(mag) max(mag)]);
title('fft频谱图');%标题：
xlabel('频率/Hz');%标注横坐标
ylabel('幅度');%标注纵坐标
grid on;%打开网格线
%数字判断
a1=(697+770)/2;a2=(770+852)/2;a3=(852+941)/2;
a4=(1209+1336)/2;a5=(1336+1477)/2;a6=(1477+1633)/2;
if(f(mag==max(mag(f<=1000)))<a1)
if(f(mag==max(mag(f>1100)))<a4)
tel=1;
elseif(f(mag==max(mag(f>1100)))<a5)
tel=2;
elseif(f(mag==max(mag(f>1100)))<a6)
tel=3;
else
tel='A';
end
elseif(f(mag==max(mag(f<=1000)))<a2)
if(f(mag==max(mag(f>1100)))<a4)
tel=4;
elseif(f(mag==max(mag(f>1100)))<a5)
tel=5;
elseif(f(mag==max(mag(f>1100)))<a6)
tel=6;
else
tel='B';
end
elseif(f(mag==max(mag(f<=1000)))<a3)
if(f(mag==max(mag(f>1100)))<a4)
tel=7;
elseif(f(mag==max(mag(f>1100)))<a5)
tel=8;
elseif(f(mag==max(mag(f>1100)))<a6)
tel=9;
else
tel='C';
end
else
if(f(mag==max(mag(f>1100)))<a4)
tel='*';
elseif(f(mag==max(mag(f>1100)))<a5)
tel=0;
elseif(f(mag==max(mag(f>1100)))<a6)
tel='#';
else
tel='D';
end
end
number=[number,num2str(tel)];
end
disp(['The telephone number is:',number])
（2）滤波器
function [f1] = lv_bo_1( x,fs )
%% 双线性变换法设计Butterworth滤波器
t=0:1/fs:(size(x)-1)/fs;
Au=0.03;
d=(Au*cos(2*pi*10000*t))';
x2=x+d;
%数字滤波器指标fp=3Hz,fr=4Hz,Rp=1db,As=20db,fs=10Hz wp=0.25*pi;%通带边界频率（归一化）：wp=fp*2*pi/fs
ws=0.3*pi;%阻带边界频率（归一化）：ws=fr*2*pi/fs
Rp=1;%通带波纹
Rs=15;
%将模拟指标转换成数字指标
Fs=fs;
Ts=1/Fs;%时域最小间隔，即时域分辨率
wp1=2/Ts*tan(wp/2);%模拟低通原型滤波器通带频率
ws1=2/Ts*tan(ws/2); %模拟低通原型滤波器阻带频率
[N,Wn]=buttord(wp1,ws1,Rp,Rs,'s');%选择滤波器的最小阶数[Z,P,K]=buttap(N); %创建butterworth模拟滤波器
[Bap,Aap]=zp2tf(Z,P,K);
[b,a]=lp2lp(Bap,Aap,Wn);
%双线性变换得到分子和分母的系数向量bz,az
[bz,az]=bilinear(b,a,Fs);%用双线性变换法实现模拟滤波器到数字滤波器的转换f1=filter(bz,az,x2);%进行滤波处理
（3）端点检测
function [fan_wei]=vad(a1,a2,x,fs,n1,n2,n3)
%归一化
x=double(x);
x=x/max(abs(x));
framelen= a1;
frameinc= a2;
y=enframe(x,hanning(framelen),frameinc);%分帧
%计算短时间能量
amp=sum(abs(y),2);
%amp = sum(abs(enframe(filter([1 -0.9375], 1, x), framelen, frameinc)).^2, 2);
%计算过零率
tmp1=enframe(x(1:length(x)-1),framelen,frameinc);
tmp2=enframe(x(2:length(x)),framelen,frameinc);
signs=(tmp1.*tmp2)<0; diffs=(tmp1-tmp2)>n3;
zcr=sum(signs.*diffs,2);zcr=[zcr;zcr(end)];
zcr_yu=0.2*mean(zcr);yuzhi=0.2*mean(amp);%
N=length(amp);
for n=1:N
if amp(n)<yuzhi||(zcr(n)<zcr_yu)
continue;
end
kaitou=n;
break;
end
for n=N:-1:1
if amp(n)<yuzhi||(zcr(n)<zcr_yu)
continue;
end
jiewei=n;
break;
end
noise=[amp(1:kaitou-1);amp(jiewei+1:end)];
noise_mean=mean(noise);noise_var = std(noise);
speech_mean=mean(amp(kaitou:jiewei));%noise_mean
yuzhi1= 0.3*speech_mean; %
yuzhi2= max(0.3*speech_mean , (noise_mean + noise_var)*1.3);%
noise_zcr=[zcr(1:kaitou-1);zcr(jiewei+1:end)];
noise_zcr_mean = mean(noise_zcr);
noise_zcr_std = std(noise_zcr);
speech_zcr_mean=mean(zcr(kaitou:jiewei));
zcr_yu1 = 0.3*speech_zcr_mean;
zcr_yu2 = max(0.3*speech_zcr_mean , (noise_zcr_mean+noise_zcr_std)*0.3); st = [];en = [];
bstart_state = 0;bend_state = 0;
segment = 0;unvoice = 0;
voice_min_len = n1;% 最短语音长度70ms
unvoice_min_len = n2;%结束段最小持续50ms
st_candicate = 0;en_candicate = 0;
%% 开始端点检测
for i = 1:N
if((amp(i)>=yuzhi2&&zcr(i)>=zcr_yu1)&&~bstart_state) %find start
bstart_state=1; %
if(~st_candicate)
st_candicate=i;
end
segment =segment + 1;
elseif((amp(i)>=yuzhi2||zcr(i)>=zcr_yu1)&&bstart_state)%
if(unvoice>=unvoice_min_len) %
st = [st; st_candicate];
en_candicate = en_candicate + unvoice_min_len - 1;
en = [en; en_candicate];bstart_state = 0; bend_state = 1;
unvoice = 0; segment = 0; st_candicate = 0;
else
unvoice = 0;bend_state = 0;
segment = segment + 1;
end
elseif((amp(i)<yuzhi2&&zcr(i)<zcr_yu1)&&bstart_state) %
if segment>=voice_min_len %
unvoice = unvoice + 1;
if ~bend_state %
en_candicate = i;
end
bend_state = 1; %
else %
bstart_state = 0; segment = 0;
bend_state = 0;unvoice = 0;
st_candicate = 0;en_candicate = 0;
end
elseif((amp(i)>=yuzhi2||zcr(i)>=zcr_yu1)&&~bstart_state) %
if ~ st_candicate
st_candicate = i;
end
else%
st_candicate = 0;
continue;
end
end
if(unvoice>=unvoice_min_len )%
st = [st; st_candicate];
en = [en; (en_candicate + unvoice_min_len -1)];
segment = 0;
end
if( segment >= voice_min_len) %
st = [st; st_candicate]; en = [en; N];
end
figure(1);
subplot(3,1,1),plot(x); %原始语音波形
axis([1,length(x),-1,1]);ylabel('speech');xlabel('样本点');
for k=1:length(st) line([st(k)*frameinc,st(k)*frameinc],[-1,1],'linestyle',':','color','blu
e','LineWidth',2);
line([en(k)*frameinc,en(k)*frameinc],[-1,1],'linestyle',':','color','red ','LineWidth',2);
end
subplot(3,1,2),plot(amp);%原始语音能量
axis([1,length(amp),0,max(amp)]),ylabel('energy');,xlabel('帧数');
line([1,N],[yuzhi1,yuzhi1],'color','yellow','LineWidth',2);%由语音能量
line([1,N],[yuzhi2,yuzhi2],'color','red','LineWidth',2);%由噪声平均能量
和语音能量比较而得
for k=1:length(st)
line([st(k),st(k)],[min(amp),max(amp)],'linestyle',':','color','blue','L ineWidth',2);
line([en(k),en(k)],[min(amp),max(amp)],'linestyle',':','color','red','Li neWidth',2);
end
subplot(3,1,3),plot(zcr);%原始语音过零率
axis([1,length(zcr),0,max(zcr)]),ylabel('zcr');xlabel('帧数');
line([1,N],[zcr_yu1,zcr_yu1],'color','yellow','LineWidth',2);%由语音能量line([1,N],[zcr_yu2,zcr_yu2],'color','red','LineWidth',2);%语音加噪音
for k=1:length(st)
line([st(k),st(k)],[min(zcr),max(zcr)],'linestyle',':','color','blue','L ineWidth',2);
line([en(k),en(k)],[min(zcr),max(zcr)],'linestyle',':','color','red','Li neWidth',2);
end
fan_wei=[st,en]*frameinc;。