数字信号处理综合MATLAB设计 双音多频拨号系统_数字信号处理实验六双音多频拨号

一、实验目的
1. 了解数字信号处理当今应用的基本情况。
2. 对该课程做系统地总结。
3．将所学知识运用到实践中，能够学以致用。

二、实验原理
双音多频 DTMF（ Dual Tone Multi-Frequency ）信号，是用两个特定的单音频率信号的组合来代表数字或功能。在 DTMF 电话机中有 16 个按键，其中 10 个数字键 0 —9 ， 6 个功能键 * 、 # 、 A 、 B 、 C 、 D 。其中 12 个按键是我们比较熟悉的按键另外由第 4 列确定的按键作为保留，作为功能键留为今后他用。 根据 CCITT 建议，国际上采用 697Hz 、770Hz 、852Hz 、94lHz 低频群及 1209Hz 、1336Hz 、1477Hz 、1633Hz高频群。对应的 DTMF 信号用sin(2πf₁t) + sin(2πf₂t) 表示。从低频群和高频群任意各抽出一种频率进行组合，共有 16 种组合，代表 16 种不同的数字键或功能，每个按键唯一地由一组行频和列频组成，如表1 所示。[10]

（1）DTMF信号的产生合成
        一般用数字方法产生DTMF信号：规定用8kHz对DTMF信号进行采样，采样后得到时域离散信号为：x(n) = sin(2πf1n / 8000) + sin(2πf2n / 8000)（式4-1）
        有计算法和查表法两种方法可形成上述序列，计算法运行速度慢，查表法速度快但要占用一定的存储空间。因为采样频率是8kHz，既要求每125ms输出一个样本，得到的序列再送到D/A转换器。输出经过平滑滤波便是连续的DTMF信号， DTMF信号通过电话线路再送到交换机。
（2）DTMF信号的检测
        要实现电话拨号音（ DTMF ）信号的检测识别，可以先将收到的连续信号经过 A/D 变换，变成数字信号再进行检测。可用 DFT(FFT)对双音多频信号进行频谱分析，由信号的幅度谱判断信号的两个频率，最后确定对应的数字和符号。当检测的频率数目较少时，用滤波器组实现更适合。FFT 是 DFT 的快速算法，但当要计算的频率点数目远小于 DFT 的区间长度时，用 FFT 快速算法的效果并不明显，而且还要占用很多内存，因此不如直接用 DFT 合适。戈泽尔（Goertzel）算法是计算 DFT 的一种线性滤波方法，用该算法实现软件实验。对于戈泽尔算法建议同学们自己查阅相关资料，掌握算法的基本思想。
（3） 检测 DTMF 信号的DFT 参数选择
        用 DFT 检测模拟信号 DTMF 信号所含的两个音频频率，即为用 DFT 对模拟信号进行频谱分析的问题。根据数字信号处理所学知识，要确定三个参数：采样频率 Fs ，DFT 的变换点数 N ，需要对信号的观察时间的长度Tp  。这三个参数不能随意选取，要根据对信号频谱分析要求确定。对信号频谱分析有三个要求：
        1．频谱分析的分辨率。观察要检测到的 8 个频率，相邻间隔最小的是第一和第二个频率，间隔是 73Hz，即要求 DFT 至少能分辨相隔 73Hz 的两个信号。所以对信号的观察时间应为Tp  = 1/ F = 1/ 73 = 13.7ms 。考虑到可靠性，应留有富裕量，要求按键时间在 40ms 以上。
        2. 频谱分析的频率范围  要检测信号的频率范围是 697—1633Hz，但考虑到存在语音干扰，除检测这 8 个频率外，还要检测他们的二次倍频的幅度大小。这样频谱分析的频率范围是 697-3266Hz。按照采样定理，则有最小采样频率为Fs min   = 6.35kHz ，系统规定 Fs  = 8kHz ，这样可算出对信号的最少采样点数为Nmin   = Tp min Fs  ≈ 110 。
        3. 检测频率的准确性  这是一个用 DFT 检测正弦波频率是否准确的问题。序列的 DFT 是序列频域函数在0 − 2π  区间的等间隔采样，如果是周期序列，截取周期序列的整个周期进行 DFT，其采样点刚好在周期信号的频率上，DFT 的幅度最大处准确的是信号频率。分析 DTMF 信号不可能经过采样得到周期序列，因此存在一个检测频率的准确性问题。

        DFT 的频率采样点频率为wk  = 2πk / N (k = 0,1,2,..., N − 1) ，相应的在模拟域的采样点频率为 fk   = Fsk / N (k = 0,1,2,..., N − 1) ，要选择一个合适的 N 值使得用该公式计算出的 fk  能接近要检测的实际频率。或者将实际的频率值带入能得到的 k值接近整数值。经过分析研究，发现 N=205 是最好的。随意将实际频率带入上面的公式可得到最接近整数的 K 值。例如：

 通过一系列的分析，可确定： Fs  = 8kHz, N = 205,Tp  ≥ 40ms
        工具箱函数 Goerztel 的调用格式为：Xgk=goertzel(xn,K+1)
其中 xn 是被变换的时域序列，用于 DTMF 信号检测时就是 DTMF 信号的 205 个采样值。K 是要求计算的 DFT 频点的序号向量。用 N 表示 xn 的长度，则要求 0 ≤ K ≤ N − 1。由以上的分析可知只计算 DTMF 信号 8 个基频时K=[18,20,22,24,31,34,38,42]
        Xgk 是变换结果向量，存放的是由 K 指定的频率点的 DFT[x(n)]的值[10]。
三、 实验内容
        要求输入 6 位或 8 位电话号码，利用戈泽尔算法求出各数字信号的频谱并绘出频谱图形，再检测输出电话号码。
参考流程图：

 四、实验报告要求
1.简述实验目的及原理。明白具体的问题解决的思路是什么。
2.整理好经过运行并证明是正确的实验程序并加上注释，绘出相应的图形。

f1=697;
f2=770;
f3=852;
f4=941;
F1=1209;
F2=1336;
F3=1477;
F4=1633;
N=205;
tmin=40/1000;
fs=8000;
dt=1/fs;
N1=floor(tmin/dt);
if N<=N1
t=[0:N-1]*dt;
else
error('The time is not enough!');
end
k1=sin(2*pi*f1*t)+sin(2*pi*F1*t);
k2=sin(2*pi*f1*t)+sin(2*pi*F2*t);
k3=sin(2*pi*f1*t)+sin(2*pi*F3*t);
ka=sin(2*pi*f1*t)+sin(2*pi*F4*t);
k4=sin(2*pi*f2*t)+sin(2*pi*F1*t);
k5=sin(2*pi*f2*t)+sin(2*pi*F2*t);
k6=sin(2*pi*f2*t)+sin(2*pi*F3*t);
kb=sin(2*pi*f2*t)+sin(2*pi*F4*t);
k7=sin(2*pi*f3*t)+sin(2*pi*F1*t);
k8=sin(2*pi*f3*t)+sin(2*pi*F2*t);
k9=sin(2*pi*f3*t)+sin(2*pi*F3*t);
kc=sin(2*pi*f3*t)+sin(2*pi*F4*t);
km=sin(2*pi*f4*t)+sin(2*pi*F1*t);
k0=sin(2*pi*f4*t)+sin(2*pi*F2*t);
kj=sin(2*pi*f4*t)+sin(2*pi*F3*t);
kd=sin(2*pi*f4*t)+sin(2*pi*F4*t);
key=['1','2','3','a';'4','5','6','b';'7','8','9','c';'*','0','#','d'];
k=[18,20,22,24,31,34,38,42];
num=input('please enter the key:','s');
num=num-48;
nn=length(num);
disp('The number of the key is: ');
disp(nn);
number=zeros(nn,length(t));
for i=1:nn
switch num(i)
case 1
number(i,1:N)=k1;
case 2
number(i,1:N)=k2;
case 3
number(i,1:N)=k3;
case 4
number(i,1:N)=k4;
case 5
number(i,1:N)=k5;
case 6
number(i,1:N)=k6;
case 7
number(i,1:N)=k7;
case 8
number(i,1:N)=k8;
case 9
number(i,1:N)=k9;
case 0
number(i,1:N)=k0;
case 49
number(i,1:N)=ka;
case 50
number(i,1:N)=kb;
case 51
number(i,1:N)=kc;
case 52
number(i,1:N)=kd;
case -6
number(i,1:N)=km;
case -13
number(i,1:N)=kj;
otherwise
error('The key is not right!');
end
end
disp('The key is: ');
for i=1:nn
xgk=goertzel(number(i,1:N),k+1);
figure;
x=[697,770,852,941,1209,1336,1477,1633];
stem(x,abs(xgk),'.');
xlabel('Hz');
title('频谱分析');
zb=find(abs(xgk)>50);
disp(key(zb(1),zb(2)-4));
end

结果: 随机输入数字:12378