一误澳、什么是嘯叫（Howling）

在本地擴聲系統(tǒng)中兼蕊，麥克風拾取揚聲器播放的音頻后又傳回給揚聲器播放，當擴音的增益足夠大贸毕，在某些頻率就會產(chǎn)生自激振蕩郑叠，從而產(chǎn)生嘯叫。

圖1 聲音反饋回路崖咨，x(n)為麥克風采集的音頻锻拘，y(n)為揚聲器播放的音頻

二心包、相移器

相移器的原理圖如圖2所示，其中粗線代表復數(shù)信號馒铃。輸入信號是蟹腾，希爾伯特變換（Hilbert transformer）會將輸入信號轉(zhuǎn)換為。將該信號與（θ是期望的相位偏移）相乘区宇，得到復數(shù)信號娃殖。取該復數(shù)信號的實部就得到我們想要的輸出：。

圖2 相移器原理

圖3 相移器的實現(xiàn)

    % 31-tap Hilbert transformer
    b= 2/pi * [-1/15 0 -1/13 0 -1/11 0 -1/9 0 -1/7 0 -1/5 0 -1/3 0 -1 0 1 0 … 
                1/3 0 1/5 0 1/7 0 1/9 0 1/11 0 1/13 0 1/15];
    b1= b.*hamming(31)';                          % window the coefficients
    b2= [0 0 0 0  0 0 0 0  0 0 0 0  0 0 0 1];     % delay of 15 (center tap of HT)

接下來創(chuàng)建一個6Hz的正弦輸入信號，其采樣率為100Hz：

    fs= 100;                    % Hz sample frequency
    Ts= 1/fs;
    N= 128;
    n= 0:N-1;
    f0= 6;                       % Hz freq of input signal
    x= cos(2*pi*f0*n*Ts);        % input signal

現(xiàn)在對執(zhí)行希爾伯特變換：

    I= filter(b2,1,x);        % I= center tap of HT
    Q= filter(b1,1,x);        % Q= output of HT

最后根據(jù)公式1得到相移的輸出：

    theta= -pi/3;                          % rad phase shift with respect to I
    y= I*cos(theta) - Q*sin(theta);        % phase shifter output

在圖3中我們假設希爾伯特變換不會引起相移胜宇，但事實上在實現(xiàn)中會引入15個采樣點的延遲耀怜，換句話說我們是在希爾伯特轉(zhuǎn)換輸出的的基礎上增加的相移得到輸出。圖4打印了和的32個樣點掸屡，可以看到和相位相差封寞。

圖4 相移器的輸入和輸出，f0=6Hz仅财，fs=100Hz狈究，?= - π/3，藍線=I盏求，綠線=y

三抖锥、頻移器

頻移器的實現(xiàn)如圖5所示，它和相移器類似碎罚，只是我們用替換公式1中的磅废，其中為期望的頻率偏移。

圖5 頻移器的實現(xiàn)

    % 31-tap Hilbert transformer
    b= 2/pi * [-1/15 0 -1/13 0 -1/11 0 -1/9 0 -1/7 0 -1/5 0 -1/3 0 -1 0 1 0 … 
               1/3 0 1/5 0 1/7 0 1/9 0 1/11 0 1/13 0 1/15];
    b1= b.*blackman(31)';                         % window the coefficients
    b1= round(b1*2^12)/2^12;                      % quantize coefficients
    b2= [0 0 0 0  0 0 0 0  0 0 0 0  0 0 0 1];     % delay of 15 (center tap of HT)

生成輸入信號，并對它執(zhí)行希爾伯特變換导绷。這次的輸入信號不用正弦信號犀勒，而是選擇有矩形頻譜的脈沖調(diào)制信號。

    fs= 100;                       % Hz sample frequency
    Ts= 1/fs;
    N= 2048;
    n= 0:N-1;
    fc= 12;                        % Hz carrier frequency
    bw= 2;                         % Hz -3 dB bandwidth of modulated pulse
    x= .5*modpulse(fc,bw,N,fs);    % modulated pulse with approx rect spectrum
    % Apply modulated pulse to Hilbert transformer
    I= filter(b2,1,x);             % I= center tap of HT filter
    Q= filter(b1,1,x);             % Q= output of HT filter

接著通過NCO（數(shù)字振蕩器）生成正弦信號與余弦信號妥曲，并利用它們實現(xiàn)頻移贾费。

    df= 3;                              % Hz desired frequency shift
    u= mod(df*n*Ts,1);                  % phase accumulator
    theta= 2*pi*u;                      % phase
    y= I.*cos(theta) –Q.*sin(theta);    % final output

現(xiàn)在我們已經(jīng)可以得到實數(shù)輸出的頻譜了，但如果我們將圖5中所有的信號放在一起展示會看到更多有趣的內(nèi)容檐盟。圖5有實數(shù)信號褂萧、和三個復數(shù)信號：希爾伯特變換的輸出、NCO的輸出遵堵、這兩者相乘的輸出箱玷。這些復數(shù)信號可以通過以下代碼得到：

    xc= I + j*Q;            % complex Hilbert transformer output
    nco= exp(j*theta);      % complex nco output
    yc= xc.*nco;            % complex product

獲取這些信號的頻譜：

    X= fft(x,N);          XdB= 20*log10(abs(x));
    XC= fft(xc,N);        XCdB= 20*log10(abs(XC));
    YC= fft(yc,N);        YCdB= 20*log10(abs(YC));
    NCO= psd(nco,N,fs);   NCOdB= 10*log10(abs(NCO));
    Y= fft(y,N);          YdB= 20*log10(abs(Y));

這些頻譜在圖6中展示。其中圖6a是實數(shù)輸入的頻譜陌宿。圖6b是希爾伯特變換的輸出，可看作解析信號（analytic signal波丰，沒有負頻率分量的復數(shù)信號）壳坪，它的負頻率分量幾乎可以忽略。圖6c是NCO信號掰烟，同樣是個解析信號爽蝴，只在+3Hz處有值。圖6d是和相乘的結(jié)果纫骑，可以看到兩者相乘的結(jié)果是的中心頻率從12Hz移到15Hz蝎亚。最后圖6e是實數(shù)輸出，中心頻率為15Hz先馆。
假如我們只是對兩個實數(shù)信號相乘发框，而不是對兩個復數(shù)信號相乘，那得到的會是中心頻率為9Hz和15Hz的兩個信號疊加的頻譜煤墙，這并非我們期望的結(jié)果梅惯。

圖6 a.實數(shù)輸入X b.希爾伯特輸出XC c.NCO輸出 d.相乘的結(jié)果YC e.最終輸出Y

四、希爾伯特變換

下面介紹希爾伯特變換的基本原理和實現(xiàn)仿野。
一個理想的離散希爾伯特變換的頻率響應為：
該響應可以通過一個濾波器近似得到铣减，該濾波器的脈沖響應為：

該脈沖響應是非因果的，但可以通過截斷來近似脚作，令葫哗，其中為偶數(shù)，且為抽頭數(shù)。抽頭系數(shù)可以乘上一個窗函數(shù)劣针。圖7校镐、圖8顯示7抽頭希爾伯特變換器，其中圖7抽頭系數(shù)序號與公式3保持一致酿秸，而圖8重排了序號灭翔，讓抽頭系數(shù)從開始。通道在濾波器延遲網(wǎng)絡的正中間辣苏，對應著理想希爾伯特濾波器的時刻肝箱。
Q通道的響應對應的相位滯后。因此對于稀蟋，相位滯后煌张，得到為。復數(shù)信號退客。與實數(shù)輸入信號不同骏融，該信號沒有負頻率分量，是一個解析信號萌狂。

圖7 抽頭系數(shù)序號與公式3保持一致

圖8 抽頭系數(shù)從b0開始

在前面的章節(jié)我們實現(xiàn)了一個31抽頭的希爾伯特變換器箫踩，代碼如下：

     % 31-tap Hilbert transformer
    b= 2/pi * [-1/15 0 -1/13 0 -1/11 0 -1/9 0 -1/7 0 -1/5 0 -1/3 0 -1 0 1 0 … 
               1/3 0 1/5 0 1/7 0 1/9 0 1/11 0 1/13 0 1/15];
    b1= b.*blackman(31)';                         % window the coefficients
    b1= round(b1*2^12)/2^12;                      % quantize coefficients
    b2= [0 0 0 0  0 0 0 0  0 0 0 0  0 0 0 1];     % delay of 15 (center tap of HT)

這些濾波器系數(shù)如圖9所示察迟。下面計算頻率響應。為了得到一個可實現(xiàn)的濾波器，我們對輸入信號延遲了15個采樣點得到通道磁椒，因此我們要計算的頻率響應也是延遲了15個采樣點的輸入信號的響應点额。下面Matlab代碼計算得到的頻率響應如圖10所示狈癞，可以看到在低頻和頻率接近的地方夯到，頻率響應與公式2的理想希爾伯特變換的頻率響應有較大偏差，因此我們在執(zhí)行希爾伯特變換的時候售葡，要確保輸入信號的頻率落在圖10中頻率響應的平坦部分看杭。

    N= 1024;
    k= -N/2:N/2-1;
    f= k*fs/N;
    h_Q= fft(b1,N);            % h_Q = response of b1
    h_I= fft(b2,N);            % h_I = response at center tap (delay of 15 samples)
    h= h_Q./h_I;
    h= fftshift(h);            %shift dc to center (swap left and right halves of h)
    plot(f,imag(h))

圖9 31抽頭的希爾伯特變換器系數(shù)

圖10 31抽頭的希爾伯特變換器的頻率響應，fs=100Hz

參考文章
Phase or Frequency Shifter Using a Hilbert Transformer
嘯叫抑制之陷波法