如何使用matlab?MFCC+GMM的安全事件聲學(xué)檢測(cè)系統(tǒng)

發(fā)布時(shí)間：2022-02-08 09:30:10 來(lái)源：億速云閱讀：135 作者：小新欄目：開發(fā)技術(shù)

小編給大家分享一下如何使用matlab MFCC+GMM的安全事件聲學(xué)檢測(cè)系統(tǒng)，希望大家閱讀完這篇文章之后都有所收獲，下面讓我們一起去探討吧！

一、安全事件聲學(xué)檢測(cè)簡(jiǎn)介（附lunwen）

1 選題背景

公共安全問(wèn)題是社會(huì)安全穩(wěn)定所聚焦的話題之一。近年來(lái)，檢測(cè)技術(shù)與監(jiān)控自動(dòng)化正深刻地改變著人們的生活。尤其在安防領(lǐng)域，閉路電視CCTV（Closed Circuit Television）、視頻流分析、智能監(jiān)控等新技術(shù)得到了廣泛應(yīng)用，大大提高了安防監(jiān)控的管理效率。然而值得注意的是，基于視頻流的監(jiān)控手段不可避免地也具有一定的先天性缺漏，例如存在視野盲區(qū)、易受光照影響等問(wèn)題，對(duì)于事件檢測(cè)，還可能存在語(yǔ)義不明的問(wèn)題，監(jiān)控手段不夠全面。純視頻手段在槍擊、爆炸、暴恐襲擊、人群恐慌等具有較強(qiáng)語(yǔ)義性的突發(fā)公共安全事件中，往往不如聲學(xué)檢測(cè)分析手段敏感和有效。聲學(xué)事件檢測(cè)主要是使用一些聲學(xué)處理方法，刻畫現(xiàn)場(chǎng)音頻流的聲學(xué)特征，再結(jié)合適當(dāng)?shù)姆诸惼鬟M(jìn)行檢測(cè)分類，從而實(shí)現(xiàn)對(duì)音頻流中出現(xiàn)的聲學(xué)事件進(jìn)行檢測(cè)分析?；诼晫W(xué)的公共安全事件檢測(cè)在反恐、維穩(wěn)、社會(huì)治安等多個(gè)領(lǐng)域具有廣泛的使用價(jià)值和應(yīng)用前景。本課題重點(diǎn)針對(duì)槍擊與爆炸兩類突發(fā)公共安全事件，對(duì)相應(yīng)的聲學(xué)檢測(cè)方法進(jìn)行了研究。

2 研究現(xiàn)狀

對(duì)于槍聲的研究聚焦在對(duì)于膛口波和彈道波的研究分析上。吳松林等深基于彈丸的空氣動(dòng)力學(xué)模型，深入分析了彈道激波的成因和理論波形；蔣灝等分析了小口徑武器發(fā)射的膛口波和彈道激波，并設(shè)計(jì)了基于膛口激波的DOA模型對(duì)彈丸彈道軌跡進(jìn)行估計(jì)；盧慧洋分析了彈道波和膛口波在槍聲檢測(cè)與定位中的作用，并設(shè)計(jì)了一套基于正三角形麥克風(fēng)陣列的槍聲定位與測(cè)距軟硬件系統(tǒng)。
對(duì)于聲信號(hào)處理和聲學(xué)特征的研究，趙力等給出了常用的信號(hào)加窗成幀、端點(diǎn)檢測(cè)以及常用聲學(xué)特征的計(jì)算方法；韓紀(jì)慶等對(duì)聲學(xué)事件檢測(cè)技術(shù)與常用模型做了綜述性介紹。徐大為等對(duì)比了基于不同聲信號(hào)特征的端點(diǎn)檢測(cè)方法，并分析了他們對(duì)噪聲的抵制能力和運(yùn)算實(shí)時(shí)性。
針對(duì)槍聲的信號(hào)處理與聲學(xué)事件檢測(cè)研究中，蔣小為和張文等[7]通過(guò)低通濾波和譜減法針對(duì)膛口波進(jìn)行去噪處理，在實(shí)驗(yàn)中得到了與理論波形高度相近的膛口波信號(hào)波形，如圖1.22所示，并提出可以使用相關(guān)分析進(jìn)行槍聲檢測(cè)。張克剛等人研究了基于短時(shí)能量分析對(duì)槍聲信號(hào)進(jìn)行端點(diǎn)檢測(cè)的方法，并提出使用持續(xù)時(shí)間處理來(lái)剔除瞬時(shí)大能量噪聲。張濤、張文、朱強(qiáng)強(qiáng)等人的研究中指出，可以采用MFCC作為目標(biāo)片段的特征，用于進(jìn)一步給分類器進(jìn)行分類檢測(cè)。
對(duì)于聲學(xué)事件檢測(cè)的分類器，Clavel等討論了監(jiān)控環(huán)境中的槍聲檢測(cè)，并通過(guò)PCA選擇13維特征作為GMM模型的輸入特征；劉力維等提出使用10階中值濾波處理端點(diǎn)檢測(cè)中的能量序列，并用GMM對(duì)目標(biāo)片段的按MFCC特征進(jìn)行分類；朱強(qiáng)強(qiáng)分析了Logo、FFS、Adaboost三種特征選擇算法，用特征選擇算法對(duì)時(shí)域特征、頻域特征、感知域特征、基于自相關(guān)函數(shù)的特征等共計(jì)9個(gè)特征組成的特征全集進(jìn)行特征選擇，并最后輸入到GMM中進(jìn)行分類；Pimentel等人提出了通過(guò)分析聚類過(guò)程中的WSS指標(biāo)來(lái)確定聚類算法中聚類中心數(shù)目的方法。
關(guān)于聲學(xué)事件數(shù)據(jù)庫(kù)，F(xiàn)onseca等人所在的龐培法布拉大學(xué)（Universitat Pompeu Fabra, Barcelona）音樂(lè)技術(shù)研究小組為了解決目前數(shù)據(jù)驅(qū)動(dòng)型（data-driven）聲學(xué)計(jì)算研究所遇到的瓶頸和困難，發(fā)起了Freesound Datasets項(xiàng)目，并建立了一個(gè)基于眾包（crowdsourcing）、規(guī)模宏大、音頻種類較齊全的大型公開數(shù)據(jù)庫(kù)Freesound；坦佩雷理工大學(xué)(Tampere University of Technology，TUT)信號(hào)處理學(xué)系Mesaro等人發(fā)起了事件檢測(cè)挑戰(zhàn)TUT Sound Events Challenge與聲學(xué)場(chǎng)景檢測(cè)挑戰(zhàn)Acoustic Scene Classification Challenge，加速了基于聲學(xué)的事件檢測(cè)和場(chǎng)景分析的相關(guān)研究。

3 算法流程設(shè)計(jì)

如何使用matlab?MFCC+GMM的安全事件聲學(xué)檢測(cè)系統(tǒng)

二、部分源代碼

for ii = 1:24 % 7 8 wav分別為背景聲和背景+槍聲，21-24是爆炸聲
    if 8 < ii &&  ii < 21
        continue;
    end
    % 取信號(hào)
    file_name = strcat('gun',num2str(ii));
    file_name = strcat(file_name,'.wav');
    fprintf('reading %s...\n',file_name);
    [y,fs] = audioread(file_name);   
    sz = size(y);
    gun = (y(:,1))'; % 單聲道
    % 原信號(hào)
    figure(ii);
    p2 = abs(fft(gun)/length(gun));
    % size(gun)
    % size(1:length(gun)/2+1)
    gun_fft = p2(1:length(gun)/2+1);
    gun_fft(2:end-1) = 2*gun_fft(2:end-1);
    f = fs*(0:(length(gun)/2))/length(gun);
    subplot(3,2,1);plot(gun);xlabel('t / s');title('signal');
    subplot(3,2,2);plot(f,gun_fft);title('spectrum');xlabel('frequency / Hz');
    % 短時(shí)能量分析
    N = 300; % 窗寬（張克剛）
    inc = 100; % 幀移（張克剛）
    win = hamming(N);
    % frameout: num x N 
    % t: num x 1, centers of frames
    % energy: 1 x num
    [frameout,t,energy]=enframe(y,win,inc);
    t = t';
    % 自適應(yīng)短時(shí)能量閾值分割
    %size(energy)
    threshold = min(energy)+0.2*(max(energy)-min(energy));
    processed_energy = energy;
    for i = 1:length(energy)
        processed_energy(i) = 0;
        if energy(i) >= threshold 
            processed_energy(i) = 1;
        end
        %fprintf('%d: %f > %f = %d\n',ii,energy(i),threshold,processed_energy(i));
    end
    subplot(3,2,3); plot(energy,'b');title('energy');
    hold on;plot(threshold*ones(size(energy)),'g');
    subplot(3,2,5); plot(processed_energy);title('binarized energy')
    % 持續(xù)時(shí)間分析
    thr = 30; % 持續(xù)采樣點(diǎn)
    cnt = 0;
    for i = 1:length(processed_energy)
        if processed_energy(i) == 1
            if cnt > 0
                cnt = cnt+1; %計(jì)數(shù)器累加
            elseif cnt == 0
                cnt = 1; %初始化計(jì)數(shù)器
            end
            if i == length(processed_energy) && cnt < thr
                processed_energy((i-cnt):i) = 0;
            end
        elseif processed_energy(i) == 0
            if cnt > 0
                if cnt < thr
                    processed_energy((i-cnt):i) = 0;
                end
            end
            cnt = 0;
        end
        %fprintf('%f, %f\n',i,processed_energy(i));
    end
    
    subplot(3,2,3);hold on;plot(processed_energy*max(energy),'r');hold off;
    subplot(3,2,5);hold on;plot(processed_energy,'r');hold off;
end
function [f,t,eng,zcr]=enframe(x,win,inc)
%ENFRAME split signal up into (overlapping) frames: one per row. [F,T]=(X,WIN,INC)
%
%	F = ENFRAME(X,LEN) splits the vector X(:) up into
%	frames. Each frame is of length LEN and occupies
%	one row of the output matrix. The last few frames of X
%	will be ignored if its length is not divisible by LEN.
%	It is an error if X is shorter than LEN.
%
%	F = ENFRAME(X,LEN,INC) has frames beginning at increments of INC
%	The centre of frame I is X((I-1)*INC+(LEN+1)/2) for I=1,2,...
%	The number of frames is fix((length(X)-LEN+INC)/INC)
%
%	F = ENFRAME(X,WINDOW) or ENFRAME(X,WINDOW,INC) multiplies
%	each frame by WINDOW(:)
%
%   The second output argument, T, gives the time in samples at the centre
%   of each frame. T=i corresponds to the time of sample X(i). 
%
nx=length(x);
nwin=length(win);
if (nwin == 1)
   len = win;
else
   len = nwin;
end
if (nargin < 3)
   inc = len;
end
len = nwin;
nf = fix((nx-len+inc)/inc);
f=zeros(nf,len);
indf= inc*(0:(nf-1)).';
inds = (1:len);
f(:) = x(indf(:,ones(1,len))+inds(ones(nf,1),:));
if (nwin > 1)
    w = win(:)';
    f = f .* w(ones(nf,1),:);
end
t = floor((1+len)/2)+indf;
%fprintf('size of f\n');
szf = size(f);
% ff = f(:).*f(:);
for i = 1:szf(1)
    %ff = f(i,:).*f(i,:)
%     ff = abs(f(i,:));
%     eng(i) = sum(ff);
    eng(i) = 0;
    zcr(i) = 0;
    for j = 1:szf(2)
        eng(i) = eng(i)+abs(f(i,j));
        if j+1 <= szf(2)
            zcr(i) = zcr(i)+abs(sign(f(i,j+1))-sign(f(i,j)));
        end
    end
    zcr(i) = 0.5*zcr(i);
end