Отслеживание основного тона с использованием нескольких оценок основного тона и HMM

Открыть сценарий в реальном времени

В этом примере показано, как выполнять отслеживание основного тона с использованием нескольких оценок основного тона, октавного и медианного сглаживания и скрытой модели Маркова (HMM).

Введение

Обнаружение основного тона является основным строительным блоком в обработке речи, кодировании речи и поиске музыкальной информации (MIR). При распознавании речи и говорящего шаг используется как признак в системе машинного обучения. Для колл-центров шаг используется для обозначения эмоционального состояния и пола клиентов. В логопедической терапии с помощью шага выявляют и анализируют патологии и диагностируют физические дефекты. В MIR, pitch используется для классификации музыки, для систем запроса по гудению, и в качестве основной функции в системах идентификации песни.

Обнаружение основного тона для чистой речи в основном считается решенной проблемой. Обнаружение основного тона с шумом и многолучевым отслеживанием остаются трудными проблемами. Существует много алгоритмов, которые широко описаны в литературе с известными компромиссами между вычислительной стоимостью и устойчивостью к различным типам шума.

Обычно алгоритм обнаружения основного тона (PDA) оценивает основной тон для данного момента времени. Затем оценка основного тона проверяется или корректируется в системе отслеживания основного тона. Системы отслеживания основного тона обеспечивают непрерывность оценок основного тона во времени.

В этом примере приводится пример функции: HelperPitchTracker, которая реализует систему отслеживания основного тона. В примере рассматривается алгоритм, реализованный HelperPitchTracker функция.

Сводка проблем

Загрузите аудиофайл и соответствующий опорный шаг для аудиофайла. Опорный тон сообщается каждые 10 мс и определяется как среднее из нескольких сторонних алгоритмов в файле чистой речи. Регионы без озвученной речи представлены как nan.

[x,fs] = audioread('Counting-16-44p1-mono-15secs.wav');
load TruePitch.mat truePitch

Используйте pitch функция для оценки основного тона звука во времени.

[f0,locs] = pitch(x,fs);

При определении ошибки основного тона обычно сообщаются две метрики: грубая ошибка основного тона (GPE) и ошибка принятия речевого решения (VDE). Так как алгоритмы основного тона в этом примере не предоставляют решение озвучивания, сообщается только GPE. В этом примере ГПЭ рассчитывается как процент оценок основного тона за пределами $\pm 10$ % опорного основного тона на протяжении диапазона вокализированных сегментов.

Рассчитайте ГПЭ для областей речи и постройте график результатов. Прослушать чистый звуковой сигнал.

isVoiced = ~isnan(truePitch);
f0(~isVoiced) = nan;

p = 0.1;
GPE = mean(abs(f0(isVoiced)-truePitch(isVoiced)) > truePitch(isVoiced).*p).*100;

t = (0:length(x)-1)/fs;
t0 = (locs-1)/fs;
sound(x,fs)

figure(1)
tiledlayout(2,1)
nexttile
plot(t,x)
ylabel('Amplitude')
title('Pitch Estimation of Clean Signal')

nexttile
plot(t0,[truePitch,f0])
legend('Reference','Estimate','Location','NW')
ylabel('F0 (Hz)')
xlabel('Time (s)')
title(sprintf('GPE = %0.1f%%',GPE))

Смешайте речевой сигнал с шумом при $SNR -5$ дБ.

Используйте pitch функция шумного звука для оценки основного тона во времени. Вычислите ГПЭ для областей озвученной речи и постройте график результатов. Слушайте шумный звуковой сигнал.

desiredSNR = -5;
x = mixSNR(x,rand(size(x)),desiredSNR);

[f0,locs] = pitch(x,fs);
f0(~isVoiced) = nan;
GPE = mean(abs(f0(isVoiced) - truePitch(isVoiced)) > truePitch(isVoiced).*p).*100;

sound(x,fs)

figure(2)
tiledlayout(2,1)

nexttile
plot(t,x)
ylabel('Amplitude')
title('Pitch Estimation of Noisy Signal')

nexttile
plot(t0,[truePitch,f0])
legend('Reference','Estimate','Location','NW')
ylabel('F0 (Hz)')
xlabel('Time (s)')
title(sprintf('GPE = %0.1f%%',GPE))

Этот пример показывает, как улучшить оценку основного тона шумных речевых сигналов с использованием генерации множественных кандидатов основного тона, октавного сглаживания, медианного сглаживания и HMM.

Алгоритм, описанный в этом примере, реализован в функции примера HelperPitchTracker. Чтобы узнать о HelperPitchTracker функция, введите help HelperPitchTracker в командной строке.

help HelperPitchTracker

 HelperPitchTracker Track the fundamental frequency of audio signal
    f0 = HelperPitchTracker(audioIn,fs) returns an estimate of the
    fundamental frequency contour for the audio input. Columns of the
    input are treated as individual channels. The HelperPitchTracker
    function uses multiple pitch detection algorithms to generate pitch
    candidates, and uses octave smoothing and a Hidden Markov Model to
    return an estimate of the fundamental frequency.
 
    f0 = HelperPitchTracker(...,'HopLength',HOPLENGTH) specifies the number
    of samples in each hop. The pitch estimate is updated every hop.
    Specify HOPLENGTH as a scalar integer. If unspecified, HOPLENGTH
    defaults to round(0.01*fs).
 
    f0 = HelperPitchTracker(...,'OctaveSmoothing',TF) specifies whether or
    not to apply octave smoothing. Specify as true or false. If
    unspecified, TF defaults to true.
 
    f0 = HelperPitchTracker(...,'EmissionMatrix',EMISSIONMATRIX) specifies
    the emission matrix used for the HMM during the forward pass. The
    default emission matrix was trained on the Pitch Tracking Database from
    Graz University of Technology. The database consists of 4720 speech
    segments with corresponding pitch trajectories derived from
    laryngograph signals. The emission matrix corresponds to the
    probability that a speaker leaves one pitch state to another, in the
    range [50, 400] Hz. Specify the emission matrix such that rows
    correspond to the current state, columns correspond to the possible
    future state, and the values of the matrix correspond to the
    probability of moving from the current state to the future state. If
    you specify your own emission matrix, specify its corresponding
    EMISSIONMATRIXRANGE. EMISSIONMATRIX must be a real N-by-N matrix of
    integers.
 
    f0 = HelperPitchTracker(...,'EmissionMatrixRange',EMISSIONMATRIXRANGE)
    specifies how the EMISSIONMATRIX corresponds to Hz. If unspecified,
    EMISSIONMATRIXRANGE defaults to 50:400.
 
    [f0,loc] = HelperPitchTracker(...) returns the locations associated
    with each pitch decision. The locations correspond to the ceiling of
    the center of the analysis frames.
 
    [f0,loc,hr] = HelperPitchTracker(...) returns the harmonic ratio
    associated with each pitch decision.
 
  See also pitch, voiceActivityDetector

Описание системы отслеживания тангажа

На рисунке представлен обзор системы отслеживания тангажа, реализованной в функции примера. Следующий код проходит через внутренние рабочие процессы HelperPitchTracker примерная функция.

1. Создание кандидатов на несколько шагов

На первом этапе системы отслеживания основного тона создается несколько кандидатов основного тона с использованием нескольких алгоритмов обнаружения основного тона. Основные кандидаты основного тона, которые обычно являются более точными, генерируются с использованием алгоритмов, основанных на алгоритме суммирования остаточных гармоник (SRH) [2] и алгоритме фильтра оценки основного тона с амплитудным сжатием (PEFAC) [3].

Буферизация шумного входного сигнала в перекрывающиеся кадры, а затем использование audio.internal.pitch.SRH для генерации 5 кандидатов основного тона для каждого транзитного участка. Также возвращает относительную уверенность каждого кандидата основного тона. Постройте график результатов.

RANGE = [50,400];
HOPLENGTH = round(fs.*0.01);

% Buffer into required sizes
xBuff_SRH = buffer(x,round(0.025*fs),round(0.02*fs),'nodelay');

% Define pitch parameters
params_SRH = struct('Method','SRH', ...
    'Range',RANGE, ...
    'WindowLength',round(fs*0.06), ...
    'OverlapLength',round(fs*0.06-HOPLENGTH), ...
    'SampleRate',fs, ...
    'NumChannels',size(x,2), ...
    'SamplesPerChannel',size(x,1));
multiCandidate_params_SRH = struct('NumCandidates',5,'MinPeakDistance',1);

% Get pitch estimate and confidence
[f0_SRH,conf_SRH] = audio.internal.pitch.SRH(xBuff_SRH,x, ...
                                             params_SRH, ...
                                             multiCandidate_params_SRH);

figure(3)
tiledlayout(2,1)
nexttile
plot(t0,f0_SRH)
ylabel('F0 Candidates (Hz)')
title('Multiple Candidates from SRH Pitch Estimation')
nexttile
plot(t0,conf_SRH)
ylabel('Relative Confidence')
xlabel('Time (s)')

Создайте дополнительный набор кандидатов первичного основного тона и связанную с ним уверенность с помощью алгоритма PEF. Формирование резервных кандидатов и связанных с ними доверительных отношений с использованием алгоритма нормализованной корреляционной функции (NCF) и алгоритма определения основного тона (CEP). Зарегистрируйте только наиболее достоверную оценку кандидатов на резервное копирование.

xBuff_PEF = buffer(x,round(0.06*fs),round(0.05*fs),'nodelay');
params_PEF = struct('Method','PEF', ...
    'Range',RANGE, ...
    'WindowLength',round(fs*0.06), ...
    'OverlapLength',round(fs*0.06-HOPLENGTH), ...
    'SampleRate',fs, ...
    'NumChannels',size(x,2), ...
    'SamplesPerChannel',size(x,1));
multiCandidate_params_PEF = struct('NumCandidates',5,'MinPeakDistance',5);
[f0_PEF,conf_PEF] = audio.internal.pitch.PEF(xBuff_PEF, ...
                                             params_PEF, ...
                                             multiCandidate_params_PEF);

xBuff_NCF = buffer(x,round(0.04*fs),round(0.03*fs),'nodelay');
xBuff_NCF = xBuff_NCF(:,2:end-1);
params_NCF = struct('Method','NCF', ...
    'Range',RANGE, ...
    'WindowLength',round(fs*0.04), ...
    'OverlapLength',round(fs*0.04-HOPLENGTH), ...
    'SampleRate',fs, ...
    'NumChannels',size(x,2), ...
    'SamplesPerChannel',size(x,1));
multiCandidate_params_NCF = struct('NumCandidates',5,'MinPeakDistance',1);
f0_NCF = audio.internal.pitch.NCF(xBuff_NCF, ...
                                  params_NCF, ...
                                  multiCandidate_params_NCF);

xBuff_CEP = buffer(x,round(0.04*fs),round(0.03*fs),'nodelay');
xBuff_CEP = xBuff_CEP(:,2:end-1);
params_CEP = struct('Method','CEP', ...
    'Range',RANGE, ...
    'WindowLength',round(fs*0.04), ...
    'OverlapLength',round(fs*0.04-HOPLENGTH), ...
    'SampleRate',fs, ...
    'NumChannels',size(x,2), ...
    'SamplesPerChannel',size(x,1));
multiCandidate_params_CEP = struct('NumCandidates',5,'MinPeakDistance',1);
f0_CEP = audio.internal.pitch.CEP(xBuff_CEP, ...
                                  params_CEP, ...
                                  multiCandidate_params_CEP);

BackupCandidates = [f0_NCF(:,1),f0_CEP(:,1)];

2. Определение долгосрочной медианы

Долгосрочная медиана кандидатов основного тона используется для уменьшения числа кандидатов основного тона. Чтобы рассчитать долгосрочный медианный шаг, сначала вычислите отношение гармоник. Оценки основного тона действительны только в областях звонкой речи, где отношение гармоник высокое.

hr = harmonicRatio(xBuff_PEF,fs, ...
                  'Window',hamming(size(xBuff_NCF,1),'periodic'), ...
                  'OverlapLength',0);

figure(4)
tiledlayout(2,1)
nexttile
plot(t,x)
ylabel('Amplitude')
nexttile
plot(t0,hr)
ylabel('Harmonic Ratio')
xlabel('Time (s)')

Используйте отношение гармоник к пороговым областям, которые не включают озвученную речь в долгосрочный медианный расчет. После определения долгосрочной медианы рассчитайте нижнюю и верхнюю границы для кандидатов основного тона. В этом примере нижняя и верхняя границы были определены эмпирически как 2/3 и 4/3 медианного шага. Кандидаты за пределами этих границ наказываются на следующей стадии.

idxToKeep = logical(movmedian(hr>((3/4)*max(hr)),3));
longTermMedian = median([f0_PEF(idxToKeep,1);f0_SRH(idxToKeep,1)]);
lower = max((2/3)*longTermMedian,RANGE(1));
upper = min((4/3)*longTermMedian,RANGE(2));

figure(5)
tiledlayout(1,1)
nexttile
plot(t0,[f0_PEF,f0_SRH])
hold on
plot(t0,longTermMedian.*ones(size(f0_PEF,1)),'r:','LineWidth',3)
plot(t0,upper.*ones(size(f0_PEF,1)),'r','LineWidth',2)
plot(t0,lower.*ones(size(f0_PEF,1)),'r','LineWidth',2)
hold off
xlabel('Time (s)')
ylabel('Frequency (Hz)')
title('Long Term Median')

3. Сокращение числа кандидатов

По умолчанию кандидаты, возвращаемые алгоритмом обнаружения основного тона, сортируются в порядке убывания достоверности. Уменьшите доверие любого основного кандидата за пределами нижней и верхней границ. Уменьшите доверие в 10 раз. Выполните повторную сортировку кандидатов для алгоритмов PEF и SRH, чтобы они находились в порядке убывания достоверности. Объедините кандидатов, сохранив только двух наиболее уверенных кандидатов из каждого алгоритма.

Постройте график сокращенных кандидатов.

invalid = f0_PEF>lower | f0_PEF>upper;
conf_PEF(invalid) = conf_PEF(invalid)/10;
[conf_PEF,order] = sort(conf_PEF,2,'descend');
for i = 1:size(f0_PEF,1)
    f0_PEF(i,:) = f0_PEF(i,order(i,:));
end

invalid = f0_SRH>lower | f0_SRH>upper;
conf_SRH(invalid) = conf_SRH(invalid)/10;
[conf_SRH,order] = sort(conf_SRH,2,'descend');
for i = 1:size(f0_SRH,1)
    f0_SRH(i,:) = f0_SRH(i,order(i,:));
end

candidates = [f0_PEF(:,1:2),f0_SRH(:,1:2)];
confidence = [conf_PEF(:,1:2),conf_SRH(:,1:2)];

figure(6)
plot(t0,candidates)
xlabel('Time (s)')
ylabel('Frequency (Hz)')
title('Reduced Candidates')

4. Сделать отличительным

Если два или более кандидатов находятся в пределах заданного диапазона 5 Гц, установите их среднее значение и суммируйте их уверенность.

span = 5;
confidenceFactor = 1;
for r = 1:size(candidates,1)
    for c = 1:size(candidates,2)
        tf = abs(candidates(r,c)-candidates(r,:)) < span;
        candidates(r,c) = mean(candidates(r,tf));
        confidence(r,c) = sum(confidence(r,tf))*confidenceFactor;
    end
end
candidates = max(min(candidates,400),50);

5. Прямая итерация HMM со сглаживанием октавы

Теперь, когда кандидаты были сокращены, вы можете передать их в HMM для обеспечения ограничений непрерывности. Контуры основного тона, как правило, непрерывны для речевых сигналов при анализе в 10 мс прыжков. Вероятность перемещения основного тона из одного состояния в другое во времени называется вероятностью излучения. Вероятности излучения могут быть инкапсулированы в матрицу, которая описывает вероятность перехода от любого значения основного тона в заданном диапазоне к любому другому в заданном диапазоне. Матрица выбросов, использованная в этом примере, была создана с использованием базы данных Граца. [1]

Загрузите матрицу выбросов и соответствующий диапазон. Постройте график функции плотности вероятности (PDF) шага в состоянии 150 Гц.

load EmissionMatrix.mat emissionMatrix emissionMatrixRange

currentState = 150;

figure(7)
plot(emissionMatrixRange(1):emissionMatrixRange(2),emissionMatrix(currentState - emissionMatrixRange(1)-1,:))
title(sprintf('Emission PDF for %d Hz',currentState))
xlabel('Next Pitch Value (Hz)')
ylabel('Probability')

HMM, используемый в этом примере, объединяет вероятности излучения, которые обеспечивают непрерывность, и относительную уверенность основного тона. На каждом скачке вероятности излучения объединяются с относительной уверенностью для создания доверительной матрицы. Наилучший выбор для каждого пути определяется как максимум доверительной матрицы. HMM, используемый в этом примере, также предполагает, что только один путь может быть назначен заданному состоянию (предположение алгоритма Витерби).

В дополнение к HMM в этом примере отслеживаются октавные скачки относительно краткосрочной медианы траекторий тангажа. Если обнаружен переход на октаву, то в качестве опций для HMM добавляются кандидаты на резервный тон.

% Preallocation
numPaths = 4;
numHops = size(candidates,1);
logbook = zeros(numHops,3,numPaths);
suspectHops = zeros(numHops,1);

% Forward iteration with octave-smoothing
for hopNumber = 1:numHops
    nowCandidates = candidates(hopNumber,:);
    nowConfidence = confidence(hopNumber,:);
    
    % Remove octave jumps
    if hopNumber > 100
        numCandidates = numel(nowCandidates);
        
        % Weighted short term median
        medianWindowLength = 50;
        aTemp = logbook(max(hopNumber-min(hopNumber,medianWindowLength),1):hopNumber-1,1,:);
        shortTermMedian = median(aTemp(:));
        previousM = mean([longTermMedian,shortTermMedian]);
        lowerTight = max((4/3)*previousM,emissionMatrixRange(1));
        upperTight = min((2/3)*previousM,emissionMatrixRange(2));
        numCandidateOutside = sum([nowCandidates < lowerTight, nowCandidates > upperTight]);
        
        % If at least 1 candidate is outside the octave centered on the
        % short-term median, add the backup pitch candidates that were
        % generated by the normalized correlation function and cepstrum
        % pitch determination algorithms as potential candidates.
        if numCandidateOutside > 0
            % Apply the backup pitch estimators
            nowCandidates = [nowCandidates,BackupCandidates(hopNumber,:)];%#ok<AGROW>
            nowConfidence = [nowConfidence,repmat(mean(nowConfidence),1,2)];%#ok<AGROW>
            
            % Make distinctive
            span = 10;
            confidenceFactor = 1.2;
            for r = 1:size(nowCandidates,1)
                for c = 1:size(nowCandidates,2)
                    tf = abs(nowCandidates(r,c)-nowCandidates(r,:)) < span;
                    nowCandidates(r,c) = mean(nowCandidates(r,tf));
                    nowConfidence(r,c) = sum(nowConfidence(r,tf))*confidenceFactor;
                end
            end
        end
    end
    
    % Create confidence matrix
    confidenceMatrix = zeros(numel(nowCandidates),size(logbook,3));
    for pageIdx = 1:size(logbook,3)
        if hopNumber ~= 1
            pastPitch = floor(logbook(hopNumber-1,1,pageIdx)) - emissionMatrixRange(1) + 1;
        else
            pastPitch = nan;
        end
        for candidateNumber = 1:numel(nowCandidates)
            if hopNumber ~= 1
                % Get the current pitch and convert to an index in the range
                currentPitch = floor(nowCandidates(candidateNumber)) - emissionMatrixRange(1) + 1;
                confidenceMatrix(candidateNumber,pageIdx) = ...
                    emissionMatrix(currentPitch,pastPitch)*logbook(hopNumber-1,2,pageIdx)*nowConfidence(candidateNumber);
            else
                confidenceMatrix(candidateNumber,pageIdx) = 1;
            end
        end
    end
    
    % Assign an estimate for each path
    for pageIdx = 1:size(logbook,3)
        % Determine most confident transition from past to current pitch
        [~,idx] = max(confidenceMatrix(:));
        
        % Convert confidence matrix index to pitch and logbook page
        [chosenPitch,pastPitchIdx] = ind2sub([numel(nowCandidates),size(logbook,3)],idx);
        
        logbook(hopNumber,:,pageIdx) = ...
            [nowCandidates(chosenPitch), ...
            confidenceMatrix(chosenPitch,pastPitchIdx), ...
            pastPitchIdx];
        
        % Remove the chosen current pitch from the confidence matrix
        confidenceMatrix(chosenPitch,:) = NaN;
    end
    % Normalize confidence
    logbook(hopNumber,2,:) = logbook(hopNumber,2,:)/sum(logbook(hopNumber,2,:));
end

6. Отслеживание HMM

Как только прямая итерация HMM завершена, окончательный контур основного тона выбирается как наиболее уверенный путь. Для определения контура тангажа, выводимого с помощью HMM, перейдите в журнал назад. Рассчитайте ГПЭ и постройте график нового контура шага и известного контура.

numHops = size(logbook,1);
keepLooking = true;
index = numHops + 1;

while keepLooking
    index = index - 1;
    if abs(max(logbook(index,2,:))-min(logbook(index,2,:)))~=0
        keepLooking = false;
    end
end

[~,bestPathIdx] = max(logbook(index,2,:));
bestIndices = zeros(numHops,1);
bestIndices(index) = bestPathIdx;

for ii = index:-1:2
    bestIndices(ii-1) = logbook(ii,3,bestIndices(ii));
end

bestIndices(bestIndices==0) = 1;
f0 = zeros(numHops,1);
for ii = (numHops):-1:2
    f0(ii) = logbook(ii,1,bestIndices(ii));
end

f0toPlot = f0;
f0toPlot(~isVoiced) = NaN;
GPE = mean( abs(f0toPlot(isVoiced) - truePitch(isVoiced)) > truePitch(isVoiced).*p).*100;
figure(8)
plot(t0,[truePitch,f0toPlot])
legend('Reference','Estimate')
ylabel('F0 (Hz)')
xlabel('Time (s)')
title(sprintf('GPE = %0.1f%%',GPE))

7. Перемещение фильтра медианы

В качестве последнего шага постобработки примените подвижный медианный фильтр с длиной окна три прыжка. Рассчитайте окончательный ГПЭ и постройте график конечного контура шага и известного контура.

f0 = movmedian(f0,3);
f0(~isVoiced) = NaN;

GPE = mean(abs(f0(isVoiced) - truePitch(isVoiced)) > truePitch(isVoiced).*p).*100;
figure(9)
plot(t0,[truePitch,f0])
legend('Reference','Estimate')
ylabel('F0 (Hz)')
xlabel('Time (s)')
title(sprintf('GPE = %0.1f%%',GPE))

Оценка эффективности

HelperPitchTracker функция использует HMM для применения ограничений непрерывности к контурам основного тона. Матрица излучения HMM может быть установлена непосредственно. Лучше всего тренировать матрицу излучения на источниках звука, подобных тем, которые вы хотите отслеживать.

В этом примере используется база данных отслеживания основного тона Грацского технологического университета (PTDB-TUG) [4]. Набор данных состоит из 20 англоязычных носителей, читающих 2342 фонетически насыщенных предложения из корпуса TIMIT. Загрузите и извлеките набор данных. В зависимости от системы загрузка и извлечение набора данных может занять около 1,5 часов.

url = 'https://www2.spsc.tugraz.at/databases/PTDB-TUG/SPEECH_DATA_ZIPPED.zip';
downloadFolder = tempdir;
datasetFolder = fullfile(downloadFolder,'PTDB-TUG');

if ~exist(datasetFolder,'dir')
    disp('Downloading PTDB-TUG (3.9 G) ...')
    unzip(url,datasetFolder)
end

Создайте хранилище аудиоданных, указывающее на записи микрофона в базе данных. Установите метку, связанную с каждым файлом, в расположение связанного известного файла основного тона. Набор данных содержит записи 10 самок и 10 самцов говорящих. Использовать subset для выделения 10-го женского и мужского носителей. Обучить матрицу излучения на основе контуров опорного шага как для мужских, так и для женских динамиков с 1 по 9.

ads = audioDatastore([fullfile(datasetFolder,"SPEECH DATA","FEMALE","MIC"),fullfile(datasetFolder,"SPEECH DATA","MALE","MIC")], ...
                     'IncludeSubfolders',true, ...
                     'FileExtensions','.wav');
wavFileNames = ads.Files;
ads.Labels = replace(wavFileNames,{'MIC','mic','wav'},{'REF','ref','f0'});

idxToRemove = contains(ads.Files,{'F10','M10'});
ads1 = subset(ads,idxToRemove);
ads9 = subset(ads,~idxToRemove);

Перетасовка хранилищ аудиоданных.

ads1 = shuffle(ads1);
ads9 = shuffle(ads9);

Матрица излучения описывает вероятность перехода из одного состояния основного тона в другое. На следующем шаге создается матрица эмиссии на основе динамиков с 1 по 9 как для мужчин, так и для женщин. База данных хранит ссылочные значения шага, кратковременную энергию и дополнительную информацию в текстовых файлах с расширением f0. getReferencePitch функция считывает значения основного тона, если кратковременная энергия превышает пороговое значение. Порог был определен эмпирически в тестах на прослушивание. HelperUpdateEmissionMatrix создает двумерную гистограмму на основе текущего состояния основного тона и следующего состояния основного тона. После создания гистограммы она нормализуется для создания матрицы излучения.

emissionMatrixRange = [50,400];
emissionMatrix = [];

for i = 1:numel(ads9.Files)
    x = getReferencePitch(ads9.Labels{i});
    emissionMatrix = HelperUpdateEmissionMatrix(x,emissionMatrixRange,emissionMatrix);
end
emissionMatrix = emissionMatrix + sqrt(eps);
emissionMatrix = emissionMatrix./norm(emissionMatrix);

Определите различные типы фонового шума: белый, атмосферный, двигатель, реактивный и уличный. Выполните повторную выборку до 16 кГц, чтобы ускорить тестирование базы данных.

Определите SNR для тестирования в дБ как 10, 5, 0, -5 и -10.

noiseType = {'white','ambiance','engine','jet','street'};
numNoiseToTest = numel(noiseType);

desiredFs = 16e3;

whiteNoiseMaker = dsp.ColoredNoise('Color','white','SamplesPerFrame',40000,'RandomStream','mt19937ar with seed','BoundedOutput',true);
noise{1} = whiteNoiseMaker();
[ambiance,ambianceFs] = audioread('Ambiance-16-44p1-mono-12secs.wav');
noise{2} = resample(ambiance,desiredFs,ambianceFs);
[engine,engineFs] = audioread('Engine-16-44p1-stereo-20sec.wav');
noise{3} = resample(engine,desiredFs,engineFs);
[jet,jetFs] = audioread('JetAirplane-16-11p025-mono-16secs.wav');
noise{4} = resample(jet,desiredFs,jetFs);
[street,streetFs] = audioread('MainStreetOne-16-16-mono-12secs.wav');
noise{5} = resample(street,desiredFs,streetFs);

snrToTest = [10,5,0,-5,-10];
numSNRtoTest = numel(snrToTest);

Выполните алгоритм обнаружения основного тона для каждого SNR и типа шума для каждого файла. Вычислите среднее значение GPE для речевых файлов. Этот пример сравнивает производительность с популярным алгоритмом отслеживания основного тона: Sawtooth Waveform Inspired Pitch Estimator (SWIPE). Реализацию алгоритма MATLAB ® можно найти в [5]. Чтобы запустить этот пример без сравнения с другими алгоритмами, установитеcompare кому false. Следующее сравнение занимает около 15 минут.

compare = true;
numFilesToTest = 20;
p = 0.1;
GPE_pitchTracker = zeros(numSNRtoTest,numNoiseToTest,numFilesToTest);
if compare
    GPE_swipe = GPE_pitchTracker;
end
for i = 1:numFilesToTest
    [cleanSpeech,info] = read(ads1);
    cleanSpeech = resample(cleanSpeech,desiredFs,info.SampleRate);
    
    truePitch = getReferencePitch(info.Label{:});
    isVoiced = truePitch~=0;
    truePitchInVoicedRegions = truePitch(isVoiced);
    
    for j = 1:numSNRtoTest
        for k = 1:numNoiseToTest
            noisySpeech = mixSNR(cleanSpeech,noise{k},snrToTest(j));
            f0 = HelperPitchTracker(noisySpeech,desiredFs,'EmissionMatrix',emissionMatrix,'EmissionMatrixRange',emissionMatrixRange);
            f0 = [0;f0]; % manual alignment with database.
            GPE_pitchTracker(j,k,i) = mean(abs(f0(isVoiced) - truePitchInVoicedRegions) > truePitchInVoicedRegions.*p).*100;
            
            if compare
                f0 = swipep(noisySpeech,desiredFs,[50,400],0.01);
                f0 = f0(3:end); % manual alignment with database.
                GPE_swipe(j,k,i) = mean(abs(f0(isVoiced) - truePitchInVoicedRegions) > truePitchInVoicedRegions.*p).*100;
            end
        end
    end
end
GPE_pitchTracker = mean(GPE_pitchTracker,3);

if compare
    GPE_swipe = mean(GPE_swipe,3);
end

Постройте график полной ошибки основного тона для каждого типа шума.

for ii = 1:numel(noise)
    figure(9+ii)
    plot(snrToTest,GPE_pitchTracker(:,ii),'b')
    hold on
    if compare
        plot(snrToTest,GPE_swipe(:,ii),'g')
    end
    plot(snrToTest,GPE_pitchTracker(:,ii),'bo')
    if compare
        plot(snrToTest,GPE_swipe(:,ii),'gv')
    end
    title(noiseType(ii))
    xlabel('SNR (dB)')
    ylabel(sprintf('Gross Pitch Error (p = %0.2f)',p))
    if compare
        legend('HelperPitchTracker','SWIPE')
    else
        legend('HelperPitchTracker')
    end
    grid on
    hold off
end

Заключение

Вы можете использовать HelperPitchTracker в качестве основы для оценки производительности GPE вашей системы отслеживания тангажа или адаптации этого примера к вашему приложению.

Ссылки

[1] Г. Пиркер, М. Вулмайр, С. Петрик и Ф. Пернкопф, «Корпус слежения за шагом с оценкой сценария многолучевого слежения», Interspeech, стр. 1509-1512, 2011.

[2] Drugman, Томас и Абир Алван. «Совместное надежное обнаружение озвучивания и оценка основного тона на основе остаточных гармоник». Материалы ежегодной конференции Международной ассоциации речевой коммуникации, ИНТЕРСПИЧ. 2011, стр. 1973-1976.

[3] Гонсалес, Сира и Майк Брукс. 19-я Европейская конференция по обработке сигналов. Барселона, 2011, стр. 451-455.

[4] Лаборатория обработки сигналов и речевой связи. Доступ состоялся 26 сентября 2018 года. https://www.spsc.tugraz.at/databases-and-tools/ptdb-tug-pitch-tracking-database-from-graz-university-of-technology.html.

[5] «Артуро Камачо». Доступ состоялся 26 сентября 2018 года. https://www.cise.ufl.edu/~acamacho/english/.

[6] «Fxpefac». Описание Fxpefac. Доступ состоялся 26 сентября 2018 года. http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html.

Документация