2

我在信号处理任务中遇到了一些问题。这个想法是分析声音信号,其中包含有人在老式模拟电话线上拨打电话号码的声音。函数dtmf_84125P将信号的文件名作为输入,读取它并返回一个包含所拨号码的字符串,例如

number = dtmf_84125P('dtmf_84125P.wav')

然后例如,数字 = '0504063452'。

我试图做的如下所示:

function [output_seq] = dtmf_84125P(input_file)
%dtmf_84125P Decodes DTMF-signal to numbers.
%   input_file: string with filename of DTMF-signal
%   output_seq: string with decoded numbers from the DTMF-signal
%
%   Mikael Eriksson
%   84125P
%   mikael.eriksson@aalto.fi

[x,fs] = wavread(input_file);
N = 205;    % N-point Fourier transform
L = 5; %Frame length for energy calculations
MINPEAK = 20;
f = linspace(0,fs/2,N);
f2 = linspace(0,fs/2,ceil(N/L));
soundlength = 0.07;     
breaklength = 0.04; 
%soundsc(x)
output_seq = [];
harm1 = [697,770,852,941];
harm2 = [1209,1336,1477];
%           1209 Hz     1336 Hz     1477 Hz
%   697 Hz  1           2           3       
%   770 Hz  4           5           6 
%   852 Hz  7           8           9
%   941 Hz              0    


%Test
% t = linspace(0,1,8000);
% fr1 = 500;
% fr2 = 900;
% x = 5*sin(2*pi*fr1*t)+2*sin(2*pi*fr2*t);
% soundsc(x) 

for i = 1:round((soundlength+breaklength)*fs):length(x)-round((soundlength+breaklength)*fs)

    %Choose subset of signal i time-domain, moving the window forward in
    %every loop. Perform Fourier-transform on this subsignal...
    num = '';
    s = x(i:i+round((soundlength+breaklength)*fs)-1);
    S = fft(s,N);
    S_abs = abs(S);        
    E = zeros(ceil(N/L),1);
    ind  = 0;

    %...and calculate its energy.
    for k = 1 : L : N-L
        ind = ind+1;
        E(ind) = sum(S_abs(k:k+L-1).^2)/L;
    end;

   %Plot the Fourier-transform and the energy of the signal in the current
   %window.
   plot(f,S_abs,f2,E,'ro')
   soundsc(s)
   %pause(0.2)

    % See if there are energy peaks at certain frequencies. First find all
    % peaks, then find the two different frequerencies, and finally the
    % number of this frequency-pair.
    [pks,locs] = findpeaks(E(2:end),'MINPEAKHEIGHT', min([MINPEAK, 0.99*max(E(2:end))]));
    freq = f2(locs);

    if length(pks) >= 2
        [amplitude_sorted, index] = sort(pks,'descend');
        freq = freq(index);

        if freq(1) < freq(2) 
            [tmp,ind] = min(abs(harm1-freq(1)));
            freq1 = harm1(ind);
            [tmp,ind] = min(abs(harm2-freq(2)));
            freq2 = harm2(ind);
        else
            [tmp,ind] = min(abs(harm1-freq(2)));
            freq1 = harm1(ind);
            [tmp,ind] = min(abs(harm2-freq(1)));
            freq2 = harm2(ind);
        end

        switch freq1
            case 697
                switch freq2
                    case 1209
                        num = '1';                
                    case 1336
                        num = '2';
                    case 1477
                        num = '3';
                end
            case 770
                switch freq2
                    case 1209
                        num = '4';                
                    case 1336
                        num = '5';
                    case 1477
                        num = '6';
                end
            case 852
                switch freq2
                    case 1209
                        num = '7';                
                    case 1336
                        num = '8';
                    case 1477
                        num = '9';
                end
        end
    elseif length(pks) == 1
        num = '0';
    end

    % If a number was dialed in this window of the signal and this was
    % detected, add it to the output sequence.
    if num ~= ''
        output_seq = [output_seq,num];
    end          

end

end

所以想法是对信号进行小部分分析,知道在 0.07+0.04 s 的窗口内拨打的号码永远不会超过一个。所以我循环遍历原始信号x,每次信号的分析部分都是这个长度并命名为s。I 将其傅里叶变换为 S,并分析光谱不同点的能量。如果存在显着能量的峰值,我会尝试根据代码开头作为注释给出的表格将这些峰值分类为一个数字。

但是,此代码生成的输出只是一个空字符串。我对我在每个窗口中绘制的频谱持怀疑态度,因为它通常包含 4 个峰值而不是 2 个。此外,与桌子。我实际上不确定这种信号能量的做事方式是否是最好的方法,所以我愿意接受建议,但是已经实现了这样的代码,如果可以使用它当然会很方便。它至少可以为我节省一些时间。

我现在不会再提出这个问题了,但请不要犹豫,询问您是否需要更多信息。谢谢您的帮助!

-米凯尔

4

0 回答 0