Шумодав для рацайки. Часть 4. Шумодав.

Предварительные замечания.


В предыдущей части я описал внешние и сопутствующие процедуры для вызова, собственно, шумоизнечтожителя.
В этой части я постараюсь описать последовательность процедур и методов, которые могут привести к желаемому результату. Однако сразу замечу, что я не претендую на какую-либо научность или точность результатов. Я просмотрел массу материалов по шумодавлению в сети и, как я понял, все методы, за исключением случая когда сигнал известен, не могут дать заранее известного результата. А поскольку результат заранее неизвестен, то будем опираться на собственное восприятие того, как должно или как хорошо могло бы быть.


Вообще говоря я даже и не знаю с чего и начать. Наверное с изучения того что и каким методом мы будем давить. В любом случае наш сигнал можно представить аддитивной смесью шума и голоса — x(k) = s(k) + n(k)

Для того, чтобы выделить сигнал, нужно знать или сигнал или шум. Но мы не знаем ни того, ни другого. Ну, шум процесс случайный. Его можно характеризовать плотностью вероятности или спектральной плотностью шума. Но, несколько сделанных записей показывают, что характеристики шума могут значительно меняться во времени и, кроме того, он может меняться достаточно быстро. Более того, сам сигнал тоже крайне непредсказуем. Далее — не совсем понятно каким способом отслеживать узкополосные помехи (атмосферные, например) и что с ними вообще делать.
Кроме того, некоторые записи показывают, что шум не аддитивный! В некоторых случаях прослеживается некоторая зависимость шума от сигнала. Вообще шляпа!

Общие положения.
За исключением способов шумоочистки, когда характеристики сигналов известны, для шумоизничтожения применяются или же винеровская фильтрация или же гребенчатые фильтры (его называют еще оптимальный фильтр, метод спектрального вычитания).
Что такое первое, я не смог до конца разобраться.
Второй способ базируется на предположении нестационарности голосового сигнала. Т.е. амплитуда сигнала будет довольно-таки быстро меняться. Это несколько контрастирует с тем, что я говорил раньше — шум тоже может быть крайне нестационарным. Но… потом посмотрим…
Метод заключается в следующем:
  • Сигнал разбивается на множество узкополосных сигналов таким образом, что потом исходный сигнал из узкополосных можно восстановить без потерь.
  • В каждой полосе каким-либо образом определяется уровень или мощность сигнала.
  • Если этот уровень сильно и быстро меняется, то на основании чего-либо делается вывод о наличии полезного сигнала в этой полосе.
  • Далее зануляются сигналы в тех полосах, где не детектируется активность.
  • Из ненулевых сигналов собирается очищенный сигнал.


Таким образом, если мы правильно распилим входной сигнал на много узкополосных, правильно определим те области, где есть полезный сигнал и где есть шум, правильно соберем потом этот сигнал, то получим некоторый аналог следящего гребенчатого фильтра. Ключевое слово здесь "правильно".

На самом деле здесь не совсем корректно говорить о шумоочистке — сам шум из полезного сигнала не вычитается. Увеличивается лишь соотношение сигнал-шум широкополосного сигнала путем уменьшения коэффициента усиления на тех частотах, где полезного сигнала нет.

Недостатки этого метода следующие:
  • Появление музыкального шума — эффекта льющейся воды или сливного бачка :) ввиду того, что в некоторых случаях колебания амплитуды сигнала в полосах могут ошибочно приниматься за полезный сигнал.
  • Ухудшение разборчивости ввиду того, что в некоторых случаях полезный сигнал может быть «недоопределен» т.е. приниматься за шум.
  • Резкие колебания шума исходного сигнала может давать эффект выстрела — ба-бах!
  • Ну и еще несколько других...

Тем не менее, мы постараемся сделать так, чтобы было хорошо.

В основном я ссылался на статью Acoustic Noise Suppression for Speech Signals using Auditory Masking Eects, но есть еще несколько в сети. Даже на русском языке.

Замечу, что теория кончается очень быстро, после располосовки… Далее я буду опираться исключительно на здравый смысл и моё представление сигналов.

Схема все обработки выглядит так:


Т.е с помощью оконного Фурье располосовываем сигнал. Вычисляем есть или нет голоса. Вычисляем усиления в каждой из полос. Применяем их. Сворачиваем взад.

Параметры нашей обработки будут такими:
  • Частота дискретизации 8кГц
  • Длина кадра 32мс (256) отсчетов
  • Период кадров — 16мс
  • Хотим иметь ширину узкой полосы при обработке герц эдак 30-40


Что мы можем себе позволить:
  • Ресурсов проца у нас вроде достаточно. Т.е. не на скорости не на памяти пока что можно не экономить.
  • Можем все делать в плавающей арифметике и не париться с дробной.
  • Поскольку у нас связь симплексная, то мы можем не особо беспокоиться о задержке обработки. Нет, конечно же не секунды, но сотни миллисекунд задержки никто не заметит.
  • Можем не следовать стандартам, а делать как удобно. Вообще не смотрим в стандарты.
  • Сконцентрируемся пока что на адекватности различения шума, помех и голоса.


Теперь по-порядку:

Располосовка и синтез.

Вообще говоря для этих целей применяют банк фильтров с идеальной реконструкцией. Таковым является ФФТ с перекрывающимися окнами. Т.е. в нашем случае:
  • Каждые 16мс берутся последние 256 отсчетов.
  • Каждый отсчет умножается на значение коэффициента окна.
  • Выполняется преобразование Фурье над полученными данными.

Схема такая:


На выходе каждые 16мс будем получать значения гармоник (бины) входных отсчетов. Строго говоря это не совсем так и вот почему:
Представим что нас интересует тоько нулевая полоса (нулевой бин). Значение на выходе нулевого бина есть свертка исходного сигнала с фильтром, импульсная характеристика которого равна нашему окну. Поскольку наше окно подобрано так, что ширина полосы пропускания образованного им фильтра есть Fs/N, где Fs — частота дискретизации, N — количество отсчетов (длина Фурье), то полоса сигнала на выходе этого фильтра будет равна Fs/N. Тогда для этого «выходного» сигнала мы можем понизить частоту дискретизации до 2*Fs/N без потери информации.

Для бина n выходной сигнал есть свертка исходного и фильтром, импульсная х-ка которого есть наше окно, умноженное на комплексную экспоненту exp(j*2*Pi/N*k*n), где к — номер отсчета или, переходя в «аналоговую область» — время… Согласно свойствам пр.фурье это эквивалентно сдвигу АЧХ фильтра в частотной области на exp(j*2*Pi/N*n). Ширина полосы получаемого сигнала есть Fs/N как и для нулевого бина, а вот центральная частота есть Fs/N*n. Тогда, чтобы понизить частоту дискретизации до 2*Fs/N, нам надо умножить выходной отсчет на exp(-j*2*Pi/N*k*n), т.е. скинуть частоту в ноль. При этом k=N/2*m, где m — номер кадра. Тода эта экспонента превращается в exp(-j*2*Pi/N*N/2*n) == exp(-j*Pi*n*m) == (-1)^[n*m]. Или же это означает, что каждый нечетный бин надо умножать на -1 каждый нечетный кадр (фрейм).

Поскольку нам нужно знать только амплитуду сигнала в полосе, то мы не будем заморачиваться с умножением на -1.

Замечу так же, что не важно какой знак комплексной экспоненты — то ли это +j то ли -j. Главное, чтобы у обратного Фурья был обратный знак.

Поскольку входные данные вещественны, то получившиеся значения амплитуд бинов после фурье будут симметричны относительно нуля, а фаза их будет функцией нечетной.
В этой связи мы будем работать только с первыми N/2 == 128 бинами, которые нам являют сигналы на центральных частотах, определяемых из соотношения Fn = Fs/N*n == 31.25*n [Гц]. При синтезе-сборке мы восстановим вторую половину спектра по первой с учетом вышесказанного.

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

Схема:



Т.е. первую половину выхода мы складываем со второй половиной предыдущего.

В качестве окна выбираем окно Ханна длины 256 (по числу отсчетов) и с центром в точке 127.5.
Его АЧХ такая:


Более детально в области низких частот:


Видим, что у нас довольно-таки существенное перекрытие полос — алиасинг, а первый ноль АЧХ приходится на полосу через одну от текущей. Уровень первого бокового лепестка есть -32дБ. Это значит, что сигнал из соседних полос будет влиять на сигнал текущей. Это, наверное, не очень хорошо. Тем более это идет несколько вразрез с мыслями выше…
Можно, конечно же с этим бороться путем удлинения окна или подбором другого окна, но, как показывает практика, в некотором смысле даже хорошо что полосы перекрываются. Об этом позже.

Если есть алиасинг, то возникает вопрос — почему таким образом располосованный сигнал будет идеально восстанавливаться на выходе? Ответ проще искать не в частотной, а во временной области. Предположим, что мы ничего не делаем с полосовыми сигналами. Тогда мы как-будто бы неискаженный сигнал умножаем дважды на sin(Pi / N * n) — получается sin^2(Pi / N * n) и потом складываем это с задержанным на четверть периода синуса — т.е. с cos^2. sin^2 + cos^2 == 1. Значит ничего как будто бы и не произошло!

А вот если делаем с узкополосными сигналами, то сложнее… Здесь пока что не буду распыляться на эту тему… Отмечу лишь что, например, если мы занулим все гармоники за исключением одной (или вообще руками заменим значение соответствующего бина константой), то на выходе получим слегка амплитудно модулированную синусоиду с периодом модуляции близким к Fs/N.

Ну, други мои… Код здесь публиковать не буду… Я уверен, что каждый может соорудить линию задержки и сложить и перемножить 2 числа в цикле. Фурье быстрое я публиковал в нескольких местах, например на электрониксе. Правда там он целочисленный. Его переделать во флоты — 5 минут.

Еще один момент. Поскольку у нас входной и выходной сигнал вещественен, то нам достаточно только половины его спектра (от 0 до N/2-1). За сим Фурье длины 256 можно выполнять не над вещественной последовательностью длиной 256, а над комплексной длины 128, получив ее их исходной. При этом длина Фурье будет так же 128, что, скорее всего, быстрее. Как это сделать написано в Numerical Recipes in C, Second Edition (1992), главы 12 и 13.

В любом случае. пишем программу и смотрим, сколько времени занимает одно Фурье + умножение — 950мкс. Терпимо. Во всяком случае, пока что влезает.

Вот… после располосовки нам надо вычислить квадраты амплитуд бинов — re^2 + im^2.

Еще момент… На входе прямое фурье, на выходе обратное. На входе-выходе сигнал вещественный. Обратное фурье от вещественного сигнала можно получить из прямого путем инвертирования знака мнимой части результата фурье преобразования. Т.е. _везде_ мы будем использовать одно Фурье, но после располосовки поменяем знак мнимой части каждого бина.

VAD адекватный.

VAD — Voice Activity Detector или детектор голосовой активности (ДГА). Как его сделать…

Коррелятор.

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

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

Посмотрим на типичную фонему:


Видим, что форма сигнала повторяется каждые 5-6мс. При условии, что средняя длина фонемы есть примерно 80мс, а длинна обрабатываемых нами данных на фрейме есть 256 отсчетов или 32мс, то, скорее всего, мы не пропустим именно голос, хотя можем и пропустить всякие посторонние звуки. Ну и ладно.
Так вот — сигнал повторяется 5-6мс (это женский голос. мужской — чаще). Значит у нас будут локальные максимумы у Кф при задержки сигналов каждые 40-50 отсчетов. Но! Нас интересует не столько «когда» появился максимум Кф, сколько его наличие при ненулевых задержках. Поскольку у нас на входе 32мс данных, то мы можем себе позволить задержки от почти нуля до 32мс! Нас, конечно же, не интересуют большие задержки. Предполагая, что с учетом АЧХ детекторного тракта, максимальная гармоника речевого сигнала лежит в диапазоне 250-1000Гц, нас будут интересовать задержки от 7-8 до 30-40. Если мы нашли локальный максимум Кф в этом диапазоне, то у нас есть тоновый сигнал.

Например:


Вверху — сигнал, внизу Кф*10000. Как видим, значение Кф достаточно адекватно повторяет фонемы. Значит, его, коррелятор, на первый взгляд можно использовать в качестве VADa. Но это только на первый взгляд…

Рассмотрим такой сигнал (атмосферная помеха):


Вверху — спектр сигнала, под ним — его вид и в самом низу — Кф.

Видно что помеха дает хоть и не большой, но все-же пик Кф. Значит, если порог VADa по каким-либо причинам оказался заниженным, то мы ошибочно примем атмосферку за голос и выдадим неприятный треск в уши.

Путей решения проблемы несколько — 1. Завысить порог, но тогда мы можем пропускать наш полезный сигнал. 2. Сузить диапазон частот (отфильтровать сигнал) в котором мы рассматриваем Кф, но тогда, согласно теореме Винера-Хинчина мы за голос можем ошибочно принимать колебания шума. 3. Скомбинировать коррелятор с еще каким-либо независимым методом и на основании уже нескольких методов принимать решение о наличии сигнала.

Что ж. Будем считать, что коррелятор нам пригодится, поэтому будем его реализовывать.

Нам надо взять входной сигнал, родить его смещенную копию, перемножить и сложить. Сложность такой операции будет составлять O(N*m), где m — кол-во смещенной копии сигнала. Оч тяжелые вычисления будут. Однако, посмотрим опять на Numerical Recipes in C, Second Edition (1992) раздел 13.2.
И от себя — поскольку Кф есть свертка, то этой свертке в частотной области соответствует произведение спектра одного сигнала на комплексно-сопряженный спектр другого. Обратное Фурье от этого произведения нам даст Кф.
Но! у нас только один сигнал. Тогда произведение спектра одного сигнала на комплексно-сопряженный спектр другого есть квадрат модуля спектра исходного сигнала, а он у нас уже вычислен раньше, сразу за располосовкой. Осталось только сделать обратное фурье, посчитать мощность, найти максимум и поделить значение его на мощность сигнала.
Поскольку на входе обратного фурье гармоники чисто вещественные и симметричные относительно нуля, то и на выходе получим тоже чисто вещественный сигнал — значение Кф. А тогда нам знак комплексной экспоненты не важен. Т.е. мы можем использовать то же самое фурье что и раньше.

И еще — в нашем методе входной сигнал умножен на окно. Это проявится как размытие максимума Кф. И черт с ним. Он никуда не денется!

Стандартное отклонение и колебания спектра.

Посмотрим внимательно на спектр шумового сигнала:


И сравним его со спектром сигнала, где присутствует голос:


Видим, что если мы аппроксимируем шум какой-либо прямой с учетом АЧХ детекторного тракта, то заметим что колебания спектра отнесенные к мощности сигнала будут малы. В то же время, проделав ту же процедуру с голосовым сигналом получим большие колебания спектра. Тогда, если мы знаем, что основная энергия голосового сигнала сосредоточена в диапазоне 300-1000 Гц, то найдя эти колебания мы получим величину пропорциональную… эээ… типа энергии голоса в полосе частот.

Но тут есть хитрость. У меня не получилось адекватного результата… Я читал это, читал вот это, потом еще… Ничего хорошего не получилось. Не, получилось, но мне показалось, что так лучше:

1. В полосе частот от 150 до 1000 Гц ищем сумму квадратов разностей амплитуд текущей гармоники и предыдущей. Относим эту сумму к мощности сигнала всего. Смотрим на картинку:



Видим, что в большинстве случаев это работает. Правда я из результата вычел волшебную константу 0.4 — это необходимо для компенсации наклона АЧХ тракта. Если наклон будет другой, то и константа будет другая.

Однако, поскольку мы ограничились лишь 1000 Гц, то мы таким образом можем пропустить звуки шипящие, основная энергия которых сосредоточена выше 1500Гц. Чтобы этого не случилось, мы проделаем похожую процедуру в области высоких частот. НО! Здесь мы начнем считать «среднеквадратичное отклонение» (СОС) спектра от 1000 до 2500Гц.

Вот:


Здесь волшебная константа равна 0.76.

Видим, что эта картинка более адекватна. Однако, при наличии нестационарного шума именно СОС лажает больше всего.
Поэтому, блин, нужно учитывать и первый, эфимерный метод.

Код публиковать этого не буду — слишком элементарные операции.

Четвертый метод вычисления VAD

Три метода выше пытаются вычислить VAD рассматривая только входной сигнал без какой-либо обработки. Т.е. учитываются и шумы входа, их колебания, сигнал и т.д. Однако мы ж хотим этот шум задавить — т.е. у нас есть еще потенциальная возможность оценить VAD на задавленном сигнале. Только как? Позже. когда отфильтруем шумы.

VAD окончательный.

Вопрос — как использовать все 4 описанные выше оценки. Первое, что приходит в голову — сложить. А можно ли просто сложить? Очевидного ответа нет! Посмотрим на все 4 графика для хорошего сигнала (картинки выше):


Svar — девиация спектра в области НЧ
STD — девиация спектра в области ВЧ
Kf — Кф
VAD — четвертый метод :)

Как видим, все 4 зависимости примерно одинаковы и одинаково себя ведут при наличии сигнала. В отсутствии сигнала эти графики представляют собой некоторую шумовую дорожку. С одной стороны, поскольку любой метод работает с одним и тем же сигналом, то результат его работы должен соответствовать другим. С другой стороны, поскольку все методы независимы, то и результаты должны быть независимыми. А это значит, что при наличии сигнала, окончательное значение VADа может быть суммой всех 4х оценок, а при отсутствии — они будут складываться как шум, т.е. пропорционально корню размаха. Складываем, смотрим:


В натуре получилось! В вокальных местах уровень на 4 умножился, а в других только в 2 раза.

Конечно же, одна картинка ничего не доказывает. Но я пересмотрел массу записей и поведение моего VADa примерно везде одинаковое. На момент написания сего опуса все еще остаются несколько типа 150мс фрагментов, на которых VAD не срабатывает.

Еще один момент — на вокализованных участках все 4 оценки должны быть примерно одинаковыми. Тогда, если какая-либо из оценок сильно отличается от суммы других, то ее исключаем из рассмотрения. Зачем? Затем: Наличие сильного узкополосного тона приведет к скачку Кф, в то время как остальные оценки подрастут не сильно. Наличие атмосферки широкополосной увеличит девиацию спектра и оставит Кф прежней. Резкий скачок шума может дать колебание Кф, девиацию не тронет. И т.д.

Еще отметим следующее — если принимаемый голос тонет в шуме, то такой VAD его определяет хорошо. Однако при этом невозможно разобрать что там говорится. Поэтому, в большинстве случаев нет смысла устанавливать порог VADa чуть выше уровня шума — будет слышно много лажи. Так же нельзя устанавливать порог VADa слишком высоким — начало фразы будет теряться и будет слышно только булькание. Поэтому, мы будем не просто вычислять VAD, а вычислять его заранее! Т.е. вычислили VAD для текущего фрейма и задержали его, а вычисленный VAD применили уже к ранее задержанному фрейму.

В общем, я не претендую на научность и строгую обоснованность выбранных методов… Скажем так — меня устраивает как оно работает.

SQUELSH детектор

В общем-то это просто пиковый детектор порога VAD. Т.е. если Vad больше какого-либо установленного уровня, то Sql открываем, если меньше, то постепенно уровень SQL затухает до определенного порога, после чего закрывается.
Сюда же привинтим и уровень входного сигнала — если уровень большой, то Sql затухает медленно, меленький — быстро. Если обнаруживаем резкие скачки уровня «вниз», то затыкаем SQL, если вверх — открываем. Все просто. Вот схема:


Черт, перепутал. Нижний график — выход ноги 2!

Итак… Чем больше уровень Vad и RSSI, тем больше сопротивление резистора. Соответственно кондер разряжается медленно. Чем меньше, тем сопротивление меньше — кондер разряжается быстрее. Побором емкости и сопротивления можно регулировать время затухания.

Однако, на сегодняшний день у меня не реализовано RSSI, ибо Алан жалко, а Мегажуть пока что в машине. Позже…

Детектор уровня шума.


Собственно говоря вот… Это главное. Надо адекватно определить т.н. «шумовой профиль». Т.е. некий уровень, ниже которого сигнал считается шумом, а выше — сигналом. Это сделать довольно-таки просто если шум стационарен. Практически невозможно, если его колбасит. А надо.

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

Итак, для начала считаем шумы в полосах величинами случайными и независимыми. Тогда, если в полосе присутствует только шум, то усредненный уровень сигнала в полосе (амплитуды бина) даст нам уровень шума в полосе.

Если шум меняется, скажем его уровень уменьшается, то тогда скорость затухания оценки уровня шума будет зависеть от постоянной времени интегрирования. Соответственно, там где голоса нет, постоянна времени может быть одной, там где есть — другой. Какой? х3. (х3 — Хочу Знать).

Жажда знаний заставила меня набрать в гугле speech noise floor estimation и почитать что мне предложат. Из этого я натаскал несколько статей:



Во! отсюда начинаем вкуривать…

Вообще говоря, столько одновременно интересных, научно и не очень фантастических, откровенно смешных и грустных статей я не читал. Я честно говоря не понимаю, как можно делать вывод a priori если априори же неизвестно про сигнал ВООБЩЕ ничего!!!

Тем не менее. Включаем голову и на основании прочитанного делаем следующие выводы относительно общей схемы оценки уровня шума:
  • Если в полосе (отсюда и далее «в полосе» — по умолчанию) отсутствует полезный сигнал, то при условии стационарности шума, его уровень есть средний уровень в полосе.
  • Ср. уровень шума (СУШ) можно получить путем интегрирования сигнала за довольно длительный промежуток времени. При этом предполагается, что сигнал не меняется никогда, т.е. шум стационарен.
  • Если шум меняется, то это должно отслеживаться внешними средствами.
  • В зависимости от детектированного состояния сигнала (есть голос/нет голоса/изменился шум/и т.д.) можно менять параметры интегрирования.
  • Возможно так же менять как-либо уровень шума в зависимости от текущего состояния сигнала.
  • Шум в полосах можно считать независимым.
  • Полезный сигнал тоже необходимо считать независимым ибо он случаен.
  • Если полезный сигнал больше уровня шума на сколько-то-там дБ, то пропускаем, если меньше, то зануляем.


Вот. Т.е. если удается детектировать паузы между буквами :), то мгновенный уровень сигнала в паузе будет именно уровень шума в полосе как-либо усредненный ранее. Пока что мы не знаем как детектировать паузы, хоть у нас и есть VAD. Так же мы пока что не знаем, что делать, если пауза не продетектировалась. Ведь VAD то работает в широкой полосе, а в узкой полосе может быть шум, который надо давить. Так же не понятно что делать, если шум в полосе меняется во времени… Но пока что оставим этот вопрос…

Все работы по оценкам уровня шума что я прочитал можно разделить на 2 (или даже 3) категории:
  1. Многа букаф. Ничего не понятно, но чувствуется что авторы круты. Что они делали тоже не понятно, но их результат сильно лучше всех остальных. Как сделать то, что они придумали они и сами не знают :). От себя — я, честно говоря, думал что сие свойственно только нашим академическим авторам. Ан нет, ошибся. Им тоже...
  2. Эмпирика интегрально-счетная основанная на определении уровня сигнала в паузах между буквами. Прародитель — Rainer Martin, Senior Member, IEEE
  3. Эмпирика, основанная на предыдущей + вариации. Например… Прикольно и, главное, более понятно.


В общем то я верю, что методом спектрального вычитания можно что-то сделать и в блютуснутых гарнитурах именно так и сделано. Не могу сказать, что звук в блютусах мне очень нравится, но все-же…

Итак — пробуем. Я буду пробовать то, что мне пришло в голову. Это хорошо еще и потому, что по дороге найдем чтонидь интересное.

Первый способ. Оценка уровня шума методом скользящего среднего. Наверное, усреднение — это единственно правильный способ оценки уровня шума. Вообще, когда речь идет о случайном процессе, то первое, что приходит мне в голову это усреднение.

Для полосы со средней частотой 500Гц и ранее рисованного голоса это выглядит так —


Зеленая — усреднение. Ну в общем ничего себе так выглядит. Это скользящее среднее по 15 фреймам — 7 вперед и 7 назад или -112ms… +112ms. Очевидно, что мы так не можем делать, ибо нам тогда надо хранить аж 15 фреймов амплитуд + 15 фреймов комплексных бинов и общая задержка будет не менее 112мс. С задержкой мы справимся (нам плевать на нее), а вот с памятью… Считаем: sizeof(float) == 4: 15*4*128 + (4 + 4)*15*128 == 23040 байт. В общем-то для нас не страшно, но что-то тут не то… А если длина фонемы больше 112 миллисекунд? Это она средняя 80мс, а реально может быть и 500мс (например гласная после согласной — «па», «ма» и т.д.). В этом случае наш средний уровень будет ползти вверх на фонемах и шляпа — будем вытирать звук.



Вот. в районе 50ого фрейма и после 70ого потрутся хорошие звуки (файл air16-1.zip).
В принципе, однако, у нас есть VAD, который сработает на 40м фрейме, но становится вопросом тогда за сколько предыдущих фреймов мы должны усреднить уровень, чтобы это был уровень шума? Ну, скажем это может быть треть от максимума на последующих 7 фреймах (после включения VADa) скомбинированная с каким-то количеством предыдущих… Но тогда в полосе, где есть голос это сработает, а в той, где только шум может получиться так, что тогда этот последующий максимум нам даст выброс и мы ничего не почистим!

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

Несмотря на то, что усреднение это, наверное, самый правильный способ оценки х-к случайного процесса, нам скользящее среднее не подходит…

Второй способ. Пиковый детектор. Ну… вообще не вариант:



Нет, конечно же можно как-то там регулировать скорость спада максимума в зависимости от уровня VADa, но маловероятно, что что-либо получится. Да и не совсем понятно, что в этом случае мы получим — уровень чего? Локального пика шума? Или?

Способ третий. Пробуем реализовывать идеи Мартина (это первая статья в списке выше под пунктом 2).

Смотрите, что нам предлагает Мартин сделать.
  • Берем амплитуды гармоник в полосах.
  • Усредняем их каким-либо простым способом. Он предлагает так: A = Aprev * alpha + (1 — alpha) * A, где alpha — параметр усреднения. У него 0.5-0.9. Далее новую А называем сигналом.
  • Далее определяем уровень шума F как минимум сигнала в полосе и предыдущим уровнем шума. Если сигнал меньше, то новый уровень шума есть уровень сигнала и, одновременно, взводим счетчик.
  • Задаемся «сохраненным» уровнем шума Н. Если А > F и одновременно A < H, то H = A.
  • Все время крутим счетчик. Если счетчик досчитал до какого-то определенного значения (им определяется задержа уровня шума), то новый F есть Н. При этом новому H присваиваем какое-то бооооольшое число. Взводим счетчик опять.
  • И так далее.

Т.е. уровень шума отслеживает минимум амплитуды в течении какого-то времени и, если время вышло, становится равный минимуму амплитуды, случившемуся в течении этого времени.

В принципе получается не плохо:



Зеленая — сглаженная амплитуда. Красная — уровень шума.
На картинке показан переход от шума к голосу.

Или вот еще:


Код на матлабе, вычисляющей этот уровень шума, такой:

function martin(x)

N=length(x);
asm = x;
floor = zeros(1,N);
top = floor;
nhold = 1;

cntmax = 20;
cnt  = 0;
f = 1;
h = 1;
aprev = 0;
alpha = 0.85;

for k=2:N
    
    a = x(k);
    a = aprev*alpha + (1-alpha)*a;
    aprev = a;
    asm(k) = a;
    
    if a < f
        
        f = a;
        cnt = 0;
        h = 1;
        
    elseif a < h
 
        h = a; 
        
    elseif cnt >= cntmax
        f = h;
        h = 0.2;
        cnt = 0;
        
    end
    
    cnt= cnt + 1;
    top(k) = h;
    floor(k) = f;
    
end

floor(1) = 0;
k = 1:N;
plot(k,x,k,asm,k,floor);


В общем-то все просто и понятно. Не понятен момент как задавать параметр сглаживания. Слишком маленький — будут теряться буквы при переходе от шума к речи. Большой — будут проглатываться буквы, ибо шум будет успевать адаптироваться на фонемы.
Но это не самое страшное. Кстати, метод предполагает параметры оценки независимыми в полосах и даже предлагает какой-то там метод их оценки на основании апостериорных данных отношения С/Ш. Беда в другом. точнее не катастрофа, а то, что мне очень не нравится.

Рассмотрим переход от голоса к шуму:


Замечаем, что время, необходимое для адаптации к новому шуму это примерно время задержки шума (cntmax в коде выше) деленное на alpha и умноженное на 2 (из-за Н — задержанного шума).
Т.е. это будет проявляться как выстрел или эффект сливного бачка, когда респондент отпустит тангенту и АРУ сильно увеличит усиление. Потом, конечно же, этот шум уйдет, да и SQL отработает, но что-то здесь не так. В любом случае размышлизмы позже, а сейчас картинка, иллюстрирующая этот эффект:


Тут даже два эффекта — первый это проглатывание первой фонемы, а второй — это переход от голоса к шуму. Замечу, что эта обработка выполнена с учетом пересчета сглаживающих коэффициентов и времени адаптации — они подобраны наилучшем способом как мне кажется. И тем не менее.

Файлы чтобы посмотреть в архиве martin.zip не загружаются :( Зараза…

И еще один момент — поскольку шум всегда недооценен, то на шумовых участках должен быть установлен порог открывания полосы существенно выше оценки. А это опять-таки приводит к проглатыванию фонем. Если порог занизить, то добавляется раздражающее музыкальное оформление голоса и шума. Заметим, сто мгновенная амплитуда может превышать оценку шума децибел эдак на 20!!!

Этот метод нам, в принципе подходит, но с оговорками.

Четвертый способ основан на статье(ях) Sundarrajan Rangachari, утверждающего, что его кунг-фу сильнее предыдущего.

Этот Хари-Кришна-Хари-Рама предлагает не использовать задержанный шум, а неким нелинейным способом обрабатывать сглаженную амплитуду и постоянно сглаживать оценки шума на основании ошибки оценки. Вообще говоря, это мне что-то как-то отдаленно напоминает фильтр Кальмана, но не будем торопиться с выводами, а просто попробуем, что из этого получится. При этом он разбивает весь диапазон звука на 3 — НЧ, СЧ, ВЧ (цветомузыка, блин!) и в каждой полосе подбирает свои константы.
Пробуем. Код на матлабе такой:


gamma = 0.98;
beta = 0.8;
alpha = 0.8;

p = 0;
Pmin = zeros(1,N);
pprev = 0;
pm = 1;
aprev = 0;
Nn = Pmin;

Nnoi = 1;
alpha1 = 0.8;
alpha2 = 1;
epsilon = 0.8;

for k=1:N
    
    a = x(k);
    a = aprev*alpha + (1-alpha)*a;
    aprev = a;
    asm(k) = a;
    p = a;
    
    if pm < p
        pm = gamma*pm + (1-gamma)/(1-beta)*(p-beta*pprev);
    else
        pm = p;
    end
    
    Pmin(k) = pm;
    
    s = p/pm;
    
    if s < 3
        delta = alpha1;
    else
        delta = alpha2;
    end
    
    if p/Nnoi < 1.3
        Nnoi = epsilon * Nnoi + (1-epsilon)*p;
    else
        Nnoi = delta * Nnoi + (1-delta)*p;
    end
    
    Nn(k)= Nnoi;
    pprev = p;
    
end


Прогоняем те же данные и получаем:


Чуда не произошло. Да, адаптация на шум происходит более быстро (раза в 2 быстрее), но вот что происходит на вокальных участках мне не совсем понятно. Похоже я ошибся с фильтром Кальмана. Не похоже. Что я вижу сдесь это интегратор с переменной постоянной времени, который может замирать на резких переходах и влачить свое жалкое существование все остальное время.
За исключение несколько быстрой адаптации на переходе от голоса к шуму этот метод даже хуже предыдущего в смысле глотаня звуков.

В принципе, может быть у меня не совсем правильный шум. Но какой есть. У них шум от автомобиля и от галдежа…

И остальные методы.
Они все примерно одинаковы, являются комбинацией первых двух и их суть заключается в следующем:
  1. Сглаживаем входные амплитуды.
  2. За уровень шума принимаем минимум амплитуды в перерыве между буквами.
  3. Оцениваем шум как сглаженная амплитуда не более чем локальный минимум.
  4. Если голос есть, то сглаживаем оценку сильно, если нет, то не сильно.

В общем разные методы работают примерно одинаково с некоторыми особенностями на различных шумах и т.п. В принципе их всех можно заставить работать одинаково путем подбора параметров для каждого конкретного случая. Но случаи бывают разные.

Лично мне более всего по душе метод Мартина. На нем, пока что и остановимся.
Допустим, что нам удалось вычислить уровень шума.

Шумооочистка.

Тут все просто: Если уровень сигнала в полосе больше вычисленного уровня шума на какое-то кол-во дБ, то бин в полосе не трогаем.Если меньше, то зануляем. Честно говоря, получается полная фигня! Мне файл не загрузить, но картинка такая:


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

Поэтому, вместо того, чтобы занулять шумные бины, мы будем их ослаблять на какое-то кол-во дБ, которые устанавливаются ручкой переменного резистора который второй слева внизу :)

И еще. Поскольку АЧХ завалена, то было бы прикольно из розового шума сделать белый. Это делается так (только для тех бинов, где определился шум):
  1. Вычисляем мощность шума как сумму оцененных мощностей шумов во всех бинах.
  2. Делим на кол-во бинов.
  3. Сравниваем мощу шума в бине со средним.
  4. Если меньше, то не трогаем этот бин и применяем к нему коэффициент усиления, установленный резистором
  5. Если больше, то находим отношение оценки уровня шума к среднему значению и применяем к этому бину этот коэффициент

Получается:


вот так.

Секретная добавка VAD. Четвертый метод.

Теперь, когда мы получили очищенную речь (если получили), то можем посчитать отношение С/Ш. Т.е. суммируем мощности сигнала в открытых бинах и относим их к сумме в закрытых. Это и будет оно. Оно всегда будет больше или равно 1. Осталось только пронормировать на что-либо. На что? Ну, понятно, что макс. значение будет бесконечность — это когда все бины открыты. Такая ситуация случается, когда шум недоопределен — в основном на переходах от голоса к шуму. Во всех остальных случаях, как правило соотношение открытых полос к закрытым составляет не более 1/2. Значит, в нашей конфигурации максимум будет равен 64. Что-то не то… Слишком эта оценка будет выбиваться из предыдущих.
Сделаем не так. Мы знаем, что основная энергия голоса сосредоточена в полосе 300 — 1500 Гц. В этой полосе мы посчитаем отношение мощности сигнала в открытых полосах к количеству вообще полос. Во! в этом случае получилось хорошо. Картинка есть выше.

Другое дело, что эта оценка получилось на фрейм задержанной, но… и черт с ней.

Заключение

Ну, как-то так… вот… мне файлов сюда почему-то не загрузить. только один удалось. Если кому интересно, могу выслать прошивку для discovery по почте.

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

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

Реализуем идеи выше. Смотрим. Обработка фрейма заняла 4.5мс.

Остался не привинченным RSSI.

… это не конец…
  • +7
  • 17 января 2013, 14:30
  • diwil
  • 1
Файлы в топике: air16-1.zip

Комментарии (7)

RSS свернуть / развернуть
А вот и мясо!
Спасибо, вечером приступлю к трапезе :)
0
Ого! А вообще, поговорю-ка я со знакомыми с кафедры фонетики в Питерском Универе, они точно таким занимались, может подскажут чего ещё почитать, тогда сообщу.

А ещё есть фирма ЦРТ, она такое точно умеет но наверняка ничего не посоветует по понятным причинам :)
0
Интересно, мяса только правда нет… исходников нету. По предыдущим статьям я уж наивно подумал что неужто автор решит поделиться сокровенным, ан нет… :( Зачем тогда в предыдущих статьях было выкладывать исходники которые по большому счету нафиг не сдались, а вот когда перешли непосредственно к теме, увы… печально. В прочем за статью спасибо, тут редкость статьи серьезного уровня.
0
Ну… тут пока что в исходниках то давать нечего. Где взять что-то осязаемое я сказал, как найти максимум последовательности… я не сказал :) Как вычислить сумму квадратов разности двух соседних членов последовательности я тоже не написал :)
Самые интересные вещи это вычисление уровня шума. В принципе я в этой части не определился как его вычислять, за сим, пока что исходников не выкладываю.
0
Ой, блин… ошибся!!!
Вся обработка занимает 1.5мс!!! Это примерно 12МГц!!!
0
Утащил все части в блог «Радио»
0
Красавчик. Цикл статей тянет на хорошую дипломную работу.
+1
Только зарегистрированные и авторизованные пользователи могут оставлять комментарии.