基于matlab的DTMF信号的产生和检测（1）

DTMF是什么？Double Tone Multiple Frequency就是这个，双音多频信号，不仅在我们现在的音频电话中使用，还在数据通信中处理各种数据流等信息。

一个DTMF信号由两个单音频正弦信号组成，有一个低频组（697hz，770hz，852hz，941hz），还有一个高频群（1209hz，1366hz，1477hz，1633hz），特点：

1.所有音频都位于人的可听范围内，因此按键下去时人可以听到。

2.8个频率中没有一个频率是其他任意一个频率的倍数。

3.任意两个频率的组合，相加或相减都不等于其他任意一个频率。

这些特性不仅简化了DTMF信号的解码同时也降低了DTMF误检的概率。

记得在《名侦探柯南》剧场版12里面柯南和怜子被困在水库里，柯南用足球把固定电话的听筒给打下来后和怜子分别发出697hz，1209hz以及941hz，1366hz的声音，最终报警成功。我也是醉了，柯南不是音痴吗？怎么会有绝对音感？且不说这个神奇的现象。拨号电话能够这样通过两个频率的声音组合作为DTMF信号拨出去，还是让我对这个提起了一点兴趣。

但是规定DTMF生成模块中信号的采样频率为1M多，我不知道什么音频信号会达到那么高，但是现在用matlab仿真的话的确还是有可能的。

产生的信号频率精度优于1%，这个之后再实际硬件上再实现就好了。

DTMF用下式表示：

x(t) = sin(2*pi*f1*t) + sin(2*pi*f2*t)

x = sin(2*pi*n*f1(p)/fs) + sin(2*pi*n*f2(q)/fs)

这个就是在matlab上面写的，但是采样频率这里还有待考虑。

输入的时候还需要加上高斯白噪声，

SNR = input(输入信道信噪比)；w = awgn(x,SNR);

这个函数就有用。

tm=[1,2,3,65;4,5,6,66;7,8,9,67;42,0,35,68];   % DTMF信号代表的16个数

N=205;K=[18,20,22,24,31,34,38,42];

f1=[697,770,852,941];                   % 行频率向量

f2=[1209,1336,1477,1633];               % 列频率向量

TN=input('键入11位电话号码= ');          % 输入11位数字

SNR=input('键入信道信噪比= ');           % 输入信噪比

TNr=0;                                  %接收端电话号码初值为零

for l=1:11;

d=fix(TN/10^(11-l));

TN=TN-d*10^(11-l);

for p=1:4;

for q=1:4;

if tm(p,q)==abs(d); break,end       % 检测码相符的列号q

end

if tm(p,q)==abs(d); break,end      % 检测码相符的行号p

end

n=0:1023;                               % 为了发声，加长序列

x = sin(2*pi*n*f1(p)/fs) + sin(2*pi*n*f2(q)/fs);% 构成双频信号

sound(x,8000);                               % 发出声音（真的就是发出话音信号的声音了）

pause(0.15)

很多代码需要解释，你产生了信号，当然需要检测出来，因为matlab仿真，需要显示出你所按的键的图样。这里用到了goertzel算法。

这个算法实质上是以一个二阶IIR带通滤波器为基础，可由两个差分方程表示，其实是由DFT实现的，其运算次数比FFT在点数多的时候方便，因为可以有选择地计算个别点出的DFT，而对不需要的点不进行运算。

X=goertzel(w(1:N),K+1);              % 用Goertzel算法计算八点DFT样本

val = abs(X);                           % 列出八点DFT向量

subplot(4,3,l);

stem(K,val,'.');grid;xlabel('k');ylabel('|X(k)|') % 画出DFT(k)幅度

axis([10 50 0 120])

limit = 80;                 %

for s=5:8;

if val(s) > limit, break, end       % 查找列号

end

for r=1:4;

if val(r) > limit, break, end       % 查找行号

end

TNr=TNr+tm(r,s-4)*10^(11-l);

end

disp('接收端检测到的号码为：')         % 显示接收到的字符

%disp(char(65))

disp(TNr)

这里需要DFT，所以要涉及到N，和k，N就是DFT所需要的采样的，所变换的点数。k就是要求计算的DFT[x]的频点序号向量。

检测DTMF信号的DFT参数：

确定三个：采样频率fs，DFT的变换点数，需要对信号的观察时间的长度Tp。

对信号频谱分析也有要求：频率分辨率，谱分析的频谱范围，检测频率的准确性。

1.相邻间隔频率最小的是73hz，Fmin =73hz，Tpmin = 13.7ms，要求按键时间>40ms

2.697~1633hz，这里不考虑二次谐波的干扰。

3.序列的N点DFT是对序列频谱函数0~2pi区间的N点等间隔采样，若是一个周期序列，截取周期序列的整数倍周期，进行DFT，其采样点刚好在周期信号的频率上，DFT的幅度最大处就是信号的准确频率。

对应模拟频率的采样点频率为fk = fs*k/N（又叫频率分辨率）

希望选一个合适的N，使该公式算出的fk能接近要检测的频率。或者用8个频率中任意一个fk’带入进去，得到k值，然后k取整数。

ITUT规定，每秒最多按10个键，即每个键时隙最短为100ms，其中音频实际持续时间至少为45ms，不大于55ms。因为我们生成DTMF信号的采样频率巨大，意思就是1s采样100多万个点哦，100ms采样15万多个个点，即我们自己设置78125个点产生DTMF信号，剩下的78125个点就是静音的点。都是50ms。

所以呢，要把N点求出来，从而求出k值。

fk/fs = k/N;

假设fs巨大，fk最多也就是1633hz吧，那么k/N就得很小很小，那么N就要很大，因为k是要从图上面看的检测DFT幅度显示所对应的频率的系数的。所以N就很大，那么就加大N，捕捉N个样值所需要的时间就增加了，而这个时间的增加会对检测的效果产生影响

弄不懂的是，DTMF产生时的采样信号和检测时的采样信号不同，而且频率偏移度根据ITU规定1.5%以内正常，超出3.5%就拒绝，这样的话697hz的拒绝误差就比1633hz的允许误差要小，因为低频群和高频群的绝对误差是不同的，应该要区别处理，到底要不要呢？要就应该把N分开用，高频群一个N，低频群一个N。

今天暂时分析到这里.。。