MATLAB完成拨号音的识别

记录一下使用MATLAB完成拨号音的识别的课程大作业的经历。

先看几张测试效果图(外放设备：iPhone 11，录音设备：AirPods Pro，录音环境：安静的家中)


最后累计做了大约20次测试左右，发现识别效果在100%，能认为这个算法还是比较ok的。

二话不说，先上核心代码！

load('recorder_filter.mat');%load滤波器设置，带通600-1600
bohao = [697,1209;697,1336;697,1477;770,1209;770,1336;770,1477;852 1209;852 1336;852 1477;941 1336];%从1-9-0的频率
[recorder,fs]=audioread('生成音频文件.wav');%录音文件
N = length(recorder);%音频文件的长度

上面这一段代码的功能大致有以下几个：

• 1、把之前用MATLAB设计好的保存为mat文件的带通滤波器load进来准备等下使用；
• 2、预先明确拨号音的频率，用一个二维数组保存
• 3、读取之前录音下来的文件，获取recorder(数据)和fs(频率)、N(长度)

% 对录音进行滤波
recorder = filter(recorder_filter,1,recorder);

% 求短时能量
wlen=200; inc=80;          % 给出帧长和帧移
win=hanning(wlen);         % 给出海宁窗
X=enframe(recorder,win,inc)';     % 分帧
fn=size(X,2);              % 求出帧数
time=(0:N-1)/fs;           % 计算出信号的时间刻度
for i=1 : fn
    u=X(:,i);              % 取出一帧
    u2=u.*u;               % 求出能量
    En(i)=sum(u2);         % 对一帧累加求和
end
frameTime=frame2time(fn,wlen,inc,fs);   % 求出每帧对应的时间
figure(2);
subplot 211; plot(time,recorder); % 画出时间波形 
title('MORSE语音波形');
ylabel('幅值'); xlabel(['时间/s' 10 '(a)']);
subplot 212; plot(frameTime,En)     % 画出短时能量图
title('短时能量');
ylabel('幅值'); xlabel(['时间/s' 10 '(b)']);

这里其实我也不太懂，我只是了解了一下基本的“短时能量”相关的知识”。我推荐如果想具体了解的可以看一下这本书《MATLAB在语音信号分析和合成中的应用》，上面讲的比较详细。

% 通过对短时能量求反,利用findpeaks函数找波谷
En_reverse = [];
En_reverse = max(En)*3 - En;%取反
[minv,minl]=findpeaks(En_reverse,'minpeakdistance',100);%找波谷，波谷间隔为100帧
hold on;plot(frameTime(minl),En(minl),'o','color','r');hold off;

比如原本信号读取进来之后经过上面的短时能量分析之后是这个样子的plot

然后我们用findpeaks这个函数之后得到的plot就是下面这个样子的

可以看到这个时候的短时能量就被我们成功的分成了有能量的和没有能量的了。关于findpeaks这个函数可以具体的去看MATLAB的官方的documentation。
既然完成了上面的分割了，那接下来就进行简单的FFT分析就行了，分析了之后的结果进行一个相似度计算就能知道分割后的语音片段包含什么频率的拨号音了。

% 对音频进行分割，得到单独的单一拨号的片段
En(En<0.00001) = 0;
target = En(minl);
point = [];
for i =1:length(minl)
    point(i) = find(time==frameTime(minl(i)));
end
% 进行FFt运算得到频率
figure(3);
fs_result = [];
for i = 1:length(point)-1
    temp = recorder(point(i)+1:point(i+1));
    temp = filter(recorder_filter,1,temp);
    subplot(2,1,1);plot((point(i)+1:point(i+1))/fs,temp);
    [P1,f] = fft_recorder(temp,fs);
    [pk1,lc1] = findpeaks(P1,'SortStr','descend','NPeaks',2);
    fs_result = [fs_result;f(lc1)];
    subplot(2,1,2);plot(f,P1);axis([0 2000 0 inf]);
end
% 对fs_result进行排序，运算找到最接近的频率
for i = 1:length(fs_result)
   fs_result(i,:) = sort(fs_result(i,:));
end
err = [];
number = [];
for i =1:length(fs_result)
    for j = 1:length(bohao)
        err = [err (fs_result(i,1)-bohao(j,1))^2+(fs_result(i,2)-bohao(j,2))^2];
    end
    if min(err) < 10
        temp = find(err == min(err));
        if temp == 10
            temp = 0;
        end
        number = [number temp];
        err = [];
    end
end
disp(number);

最后的最后还是老规矩，放上百度云。
链接：https://pan.baidu.com/s/1pY7ug2uI06lbAAYUn3EkVQ
提取码：unl4

2021.12.15更新百度云链接

链接：https://pan.baidu.com/s/1-gNFwaWShtVVsbi5yr37Fg
提取码：q2mx