На что влияет размер (глубина) семпла?

В этом посте я разберу ещё одну характеристику представления аудио на компьютере — размер или глубину семпла. В теореме Шеннона, с которой мы ознакомились в первом посте, сказано, что семплы нужно брать с определенной частотой, чтобы восстановить исходный сигнал. Однако, в ней ничего не сказано о том, сколько памяти компьютера нужно выделить на один семпл. По сути  ней семплы — это обычные вещественные числа. В компьютере же семплы чаще всего представляются целыми числами в диапазоне от $-2^{n-1}$ до $2^{n-1}-1$ (целые числа со знаком), где $n$ — количество бит на один семпл (bps, bits per sample). С количеством бит на семпл связана другая величина — динамический диапазон (или отношение сигнал-шум, что в данном случае то же самое). Он измеряется в децибелах определяется так:
$DR = 20 \lg (2^{n}) \approx 6 n$

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

Теперь, что будет, если мы выберем $n$ маленьким? Если взять синусоидальную волну с амплитудой 8 и представить её в виде целых чисел, округлив семплы в меньшую сторону, и взять разность оригинального сигнала и округленного (найти ошибку округления), увидим примерно следующую картину:

Видно, что ошибка округления сама является источником колебаний с несколькими частотами, большими, чем наша синусоида. Чем меньше число бит на семпл, тем больше ошибка. В конце концов она становится слышимой как неприятный металлический звон. Эту ошибку ещё называют ошибкой квантования (под квантованием подразумевают представление вещественного числа в виде "квантов", представимых в виде конечного числа бит, будь то целые числа, числа с плавающей запятой итд).

В том случае, если мы проводим операцию понижения числа бит на семпл (например, с 24 до 16) или проводим обработку аудио в несколько стадий (например, последовательно применяя несколько фильтров) с высоким значением bps (24-32) с последующим понижением до "обычного" (16), мы можем повлиять на ситуацию, превратив неприятный шум в приятный (сия процедура называется dithering). О том, как это сделать, и будет данный пост.

Как понизить bps?

Прежде чем разобраться, что делать с этим неприятным звоном, надо разобраться, как лучше всего понижать bps. В случае семплов, представленных целыми числами, можно делать арифметический сдвиг вправо на величину $old bps - newbps$. Также можно делать деление семлпа $s / 2^{oldbps - newbps}$ с последующим округлением результата до целого.

Есть четыре вида округления — в сторону минус бесконечности (floor, представлен на картинке выше), в сторону плюс бесконечности (ceil), в сторону нуля (truncate) и в сторону ближайшего  целого (round). Варианты floor и ceil мало отличаются между собой. На картинке выше видно, что ошибка сосредоточена вокруг значения $1/2$, то есть в сигнале среди прочих колебаний появляется константное значение. Этот факт будет не желателен в дальнейшем, так что эти способы округления я отмету.

У truncate ошибка сосредоточена вокруг нуля, но её колбасит с вдвое большей амплитудой и как-то вообще слишком страшно, поэтому я отмету и этот вариант:

Остается round. У него и ошибка вокруг нуля, и диапазон значений в отрезке $\left[-1/2, 1/2 \right]$

Теперь узнаем, как убрать неприятный шум.

Dithering.

Всё просто — чтобы убрать неприятный шум, нужно ввести приятный, маскирующий неприятный. Для этого перед понижением количества бит на семпл к сигналу добавляют некую случайную величину. Тогда ошибка квантования как бы декоррелируется от сигнала. Мы будем брать такую случайную величину $\xi$, для которой $M\xi = 0$ по той же причине, по какой выбрали метод округления — чтобы (систематически) не смещать исходный сигнал. Самый простой вид такой величины — случайная величина, равномерно распределенная в диапазоне $\left[ -\Delta, \Delta \right]$, где $\Delta$ характеризует наименее значащий бит. Например, при понижении bps с 24 до 16 $\Delta = 2^{24-16}/2 = 128$. Для оценки того, как этот шум повлияет на исходный сигнал, ему можно приписать "мощность". Вспомним пружинный маятник из курса физики. Его энергия зависит от отклонения маятника от положения равновесия: $E = \frac{k x^{2}}{2}$. Тогда разумно будет приписать ему "мощность" $W = \frac{\int_{0}^{T} E(t) dt}{T}$, где $T$ — период колебания. Для несмещенной случайной величины разумно оценить "мощность" как второй центральный момент $W = D = M(\xi - M\xi)^{2} = M\xi^{2} - (M\xi)^{2} = M\xi^{2}$. Для равномерно распределенной случайной величины $D = \int_{-\Delta}^{\Delta} x^{2} p(x) dx = \frac{1}{2 \Delta} \int_{-\Delta}^{\Delta} x^{2} dx = \Delta^{2}/3$. Интуитивно можно ожидать, что чем меньше мощность шума, тем тише он будет звучать на фоне исходного сигнала.

Можно взять ещё менее мощный шум, лежащий в том же  диапазоне $\left[ -\Delta, \Delta \right]$. Например шум, распределенный "треугольно":

Гистограмма выборки треугольного шума и его функция плотности вероятности для $\Delta=128$.

Функция плотности вероятности для этого распределения для $-\Delta \le x \le \Delta$ равна $p(x) = ( \Delta - \left|x\right|) / \Delta^{2}$, а дисперсия (в нашем случае равная "мощности") $W = D = 2\int_{0}^{\Delta} x^{2} p(x) dx = \Delta^{2}/6$. Этот шум так же хорошо "глушит" ошибку квантизации, но звучит чуть тише на фоне "полезного" сигнала.

На компьютере случайную величину с таким распределением можно взять, сложив 2 числа, полученных с помощью обычной функции random. Проверим это. Пусть $\xi_{1}, \xi_{2}$ — случайные независимые величины с равномерным распределением в отрезке $[0, \Delta]$. Найдем функцию вероятности $F(y)$ того, что сумма этих величин $\xi = \xi_{1} + \xi_{2}$ меньше (или меньше или равна) $y$. Очевидно, что ни одна из случайных величин не может быть больше, чем $y$. Используя этот факт, найдём:

Для $0 \le y \le \Delta$: $F(y) = \int_{0}^{y-\xi_{1}} d\xi_{2} \int_{0}^{y} d\xi_{1} p(\xi_{1}, \xi_{2})$
Для $\Delta \le y \le 2\Delta$: $F(y) = \int_{0}^{y-\xi_{1}}d\xi_{2} \int_{y-\Delta}^{\Delta}d\xi_{1} p(\xi_{1}, \xi_{2}) + \int_{0}^{\Delta}d\xi_{2} \int_{0}^{y-\Delta}d\xi_{1} p(\xi_{1}, \xi_{2})$

В силу независимости и равномерного распределения $\xi_{1}$ и $\xi_{2}$, $p(\xi_{1}, \xi_{2}) = p(\xi_{1})p(\xi_{2}) = 1/\Delta^{2}$

Итого для $0 \le y \le \Delta$: $F(y) = \frac{1}{2\Delta^{2}} y^{2}$
Для $\Delta \le y \le 2\Delta$: $F(y) = \frac{2y}{\Delta} - \frac{y^{2}}{2\Delta^2} -1$

Продифференцировав $F(y)$, убедимся, что случайная величина $\xi = \xi_{1} + \xi_{2}$ распределена "треугольно".

Вот как я получаю случайные величины нужного распределения в диапазоне $[-2^{magnitude-1}, 2^{magnitude-1}-1]$ на компьютере (язык Common Lisp):

(defun uniform-noise (magnitude)
  (- (random (ash 1 magnitude))
     (ash 1 (1- magnitude))))

(defun triangle-noise (magnitude)
  (+ (uniform-noise (1- magnitude))
     (uniform-noise (1- magnitude))))

Так я применяю dithering к входному сигналу (в качестве noise-fn могут выступать uniform-noise или triangle-noise.

(defun apply-simple-noise (noise-fn)
  (lambda (magnitude x)
    (+ x (funcall noise-fn magnitude))))

(defun apply-dither (array magnitude noise-fn)
  (let ((new-array (make-array (length array))))
    (map-into new-array
              (lambda (x)  (funcall noise-fn magnitude x))
              array)))

После этого можно понижать bps (здесь я захардкодил понижение с 24 бит до 16, семплы представлены как целые числа со знаком):

(defun reduce-bit-depth (array)
  (let ((new-array (make-array (length array) :element-type '(signed-byte 32))))
    (map-into new-array (lambda (x) (min 32767 (max -32768  (round x 256)))) array)))

Вот что получается на деле. Возьмем сигнал в виде гармонического осциллятора с переменной амплитудой и частотой колебания 440 Гц:

Sine wave 440 Hz (44100 Hz, 24 bps)

Понизим количество бит на семпл до 7 (семплы представлены как 16-битные целые числа со сброшенными в 0 младшими битами):

Sine wave (no dither) 440 Hz (44100 Hz, 7 bps)

Слышно, что второй файл звучит более "грязно". Сравним спектр исходного сигнала с "урезанным":

vs

Теперь возьмем вариант, где до понижения bps был применен dithering с равномерно распределенным шумом: Sine wave (uniform dither) 440Hz (44100 Hz, 7 bps). Звучит лучше? Вот его спектр:

Как мы видим, dithering сгладил все прочие шумы, так что фоновый шум стал почти белым.

Триангулярно распределенный шум почти не отличается от равномерно распределенного (хотя если чередовать их в audacity, триангулярный вроде действительно звучит тише):
Sine wave (triangular dithering) 440 Hz (44100 Hz, 7 bps)

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