Дискретные вейвлет-преобразования

Непрерывное вейвлет-преобразование

Свойства вейвлет преобразования

Требования к вейвлетам

Для осуществления вейвлет-преобразования вейвлет-функции должны удовлетворять следующим критериям:

1. Вейвлет должен обладать конечной энергией:

2. Если фурье-преобразование для, то есть

тогда должно выполняться следующее условие:

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

3. Дополнительный критерий предъявляется для комплексных вейвлетов, а именно, что для них Фурье-преобразование должно быть одновременно вещественным и должно убывать для отрицательных частот.

4. Локализация: вейвлет должен быть непрерывным, интегрируемым, иметь компактный носитель и быть локализованным как во времени (в пространстве), так и по частоте. Если вейвлет в пространстве сужается, то его средняя частота повышается, спектр вейвлета перемещается в область более высоких частот и расширяется. Этот процесс должен быть линейным – сужение вейвлета вдвое должно повышать его среднюю частоту и ширину спектра также вдвое.

1. Линейность

2. Инвариантность относительно сдвига

Сдвиг сигнала во времени на t0 приводит к сдвигу вейвлет-спектра также на t0.

3. Инвариантность относительно масштабирования

Растяжение (сжатие) сигнала приводит к сжатию (растяжению) вейвлет-спектра сигнала.

4. Дифференцирование

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

Вейвлет преобразование для непрерывного сигнала относительно вейвлет функции определяется следующим образом:

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

Весовая функция.

Мы можем определить нормированную функцию следующим образом

что означает временной сдвиг на b и масштабирование по времени на a. Тогда формула вейлет-преобразования изменится на

Исходный сигнал может быть восстановлен по формуле обратного преобразования

В дискретном случае, параметры масштабирования a и сдвига b представлены дискретными величинами:

Тогда анализирующий вейвлет имеет следующий вид:

где m и n - целые числа.

В таком случае для непрерывного сигнала дискретное вейвлет-преобразование и его обратное преобразование запишутся следующими формулами:

Величины также известны как вейвлет-коэффициенты.

есть постоянная нормировки.

Некоторые идеи теории вейвлетов появились очень давно. Например, уже в 1910 году А.Хаар опубликовал полную ортонормальную систему базисных функций с локальной областью определения (теперь они называются вейвлетами Хаара). Первое упоминание о вейвлетах появилось в литературе по цифровой обработке и анализу сейсмических сигналов (работы А.Гроссмана и Ж.Морле).

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

Несмотря на то, что теория вейвлет-преобразования уже в основном разработана, точного определения, что же такое "вейвлет", какие функции можно назвать вейвлетами, насколько мне известно, не существует. Вейвлеты могут быть ортогональными, полуортогональными, биортогональными. Эти функции могут быть симметричными, асимметричными и несимметричными.

Различают вейвлеты с компактной областью определения и не имеющие таковой. Некоторые функции имеют аналитическое выражение, другие – быстрый алгоритм вычисления связанного с ними вейвлет-преобразования. Попробуем дать вначале неформальное определение вейвлет-преобразования, а затем – его точное математическое обоснование.

Вейвлеты и многомасштабный анализ

Рассмотрим задачу, которая очень часто встречается на практике: у нас есть сигнал (а сигналом может быть все, что угодно, начиная от записи показаний датчика и кончая оцифрованной речью или изображением). Идея многомасштабного анализа (multiscale analysis, multiresolutional analysis) заключается в том, чтобы взглянуть на сигнал сначала вплотную – под микроскопом, затем через лупу, потом отойти на пару шагов, потом посмотреть издалека (рис.1).

Что это нам дает? Во-первых, мы можем, путем последовательного огрубления (или уточнения) сигнала выявлять его локальные особенности (ударение в речи или характерные детали изображения) и подразделять их по интенсивности. Во-вторых, таким образом обнаруживается динамика изменения сигнала в зависимости от масштаба.

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

Идея применения вейвлетов для многомасштабного анализа заключается в том, что разложение сигнала производится по базису, образованному сдвигами и разномасштабными копиями функции-прототипа (то есть вейвлет-преобразование по своей сути является фрактальным). Такие базисные функции называются вейвлетами (wavelet ), если они определены на пространстве L 2 (R) (пространство комплекснозначных функций f(t) на прямой с ограниченной энергией), колеблются вокруг оси абсцисс и быстро сходятся к нулю по мере увеличения абсолютного значения аргумента (рис.2).

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

Согласно принципу неопределенности, чем лучше функция сконцентрирована во времени, тем больше она размазана в частотной области. При перемасштабировании функции произведение временного и частотного диапазонов остается постоянным и представляет собой площадь ячейки в частотно-временной (фазовой) плоскости.

Преимущество вейвлет-преобразования перед, например, преобразованием Габора заключается в том, что оно покрывает фазовую плоскость ячейками одинаковой площади, но разной формы (рис.3). Это позволяет хорошо локализовать низкочастотные детали сигнала в частотной области (преобладающие гармоники), а высокочастотные – во временной (резкие скачки, пики и т.п.).

Более того, вейвлет-анализ позволяет исследовать поведение фрактальных функций – то есть не имеющих производных ни в одной своей точке!

Ортогональное вейвлет-преобразование

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

Вообще говоря, для точного восстановления сигнала достаточно знать его вейвлет-преобразование на некоторой довольно редкой решетке в фазовой плоскости (например, только в центре каждой ячейки на рис.3). Следовательно, и вся информация о сигнале содержится в этом довольно небольшом наборе значений.

Идея здесь заключается в том, чтобы масштабировать вейвлет в некоторое постоянное (например, 2) число раз, и смещать его во времени на фиксированное расстояние, зависящее от масштаба. При этом все сдвиги одного масштаба должны быть попарно ортогональны – такие вейвлеты называются ортогональными.

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

Последовательно применяя такое преобразование, мы можем получить результат нужной нам степени детальности (гладкости) и набор деталей на разных масштабах – то, о чем говорили в начале статьи. Более того, применив вейвлет-преобразование к заинтересовавшей нас детали сигнала, мы можем получить ее "увеличенное изображение". И наоборот, отбросив несущественные детали и выполнив обратное преобразование, мы получим сигнал, очищенный от шумов и случайных выбросов (например, "убрать" случайно попавшую в кадр птицу на фотографии здания).

Дискретное вейвлет-преобразование и другие направления вейвлет-анализа

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

К счастью, И.Добеши удалось найти метод, позволяющий построить (бесконечную) серию ортогональных вейвлетов, каждый из которых определяется конечным числом коэффициентов. Стало возможным построить алгоритм, реализующий быстрое вейвлет-преобразование на дискретных данных (алгоритм Малла). Достоинство этого алгоритма, помимо всего вышесказанного, заключается в его простоте и высокой скорости: и на разложение, и на восстановление требуется порядка cN операций, где с – число коэффициентов, а N – длина выборки.

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

Применение вейвлет-преобразования

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

1. Обработка экспериментальных данных. Поскольку вейвлеты появились именно как механизм обработки экспериментальных данных, их применение для решения подобных задач представляется весьма привлекательным до сих пор. Вейвлет-преобразование дает наиболее наглядную и информативную картину результатов эксперимента, позволяет очистить исходные данные от шумов и случайных искажений, и даже "на глаз" подметить некоторые особенности данных и направление их дальнейшей обработки и анализа. Кроме того, вейвлеты хорошо подходят для анализа нестационарных сигналов, возникающих в медицине, анализе фондовых рынков и других областях.
2. Обработка изображений. Наше зрение устроено так, что мы сосредотачиваем свое внимание на существенных деталях изображения, отсекая ненужное. Используя вейвлет-преобразование, мы можем сгладить или выделить некоторые детали изображения, увеличить или уменьшить его, выделить важные детали и даже повысить его качество!
3. Сжатие данных. Особенностью ортогонального многомасштабного анализа является то, что для достаточно гладких данных полученные в результате преобразования детали в основном близки по величине к нулю и, следовательно, очень хорошо сжимаются обычными статистическими методами. Огромным достоинством вейвлет-преобразования является то, что оно не вносит дополнительной избыточности в исходные данные, и сигнал может быть полностью восстановлен с использованием тех же самых фильтров. Кроме того, отделение в результате преобразования деталей от основного сигнала позволяет очень просто реализовать сжатие с потерями – достаточно просто отбросить детали на тех масштабах, где они несущественны! Достаточно сказать, что изображение, обработанное вейвлетами, можно сжать в 3-10 раз без существенных потерь информации (а с допустимыми потерями – до 300 раз!). В качестве примера отметим, что вейвлет-преобразование положено в основу стандарта сжатия данных MPEG4.
4. Нейросети и другие механизмы анализа данных. Большие трудности при обучении нейросетей (или настройке других механизмов анализа данных) создает сильная зашумленность данных или наличие большого числа "особых случаев" (случайные выбросы, пропуски, нелинейные искажения и т.п.). Такие помехи способны скрывать характерные особенности данных или выдавать себя за них и могут сильно ухудшить результаты обучения. Поэтому рекомендуется очистить данные, прежде чем анализировать их. По уже приведенным выше соображениям, а также благодаря наличию быстрых и эффективных алгоритмов реализации, вейвлеты представляются весьма удобным и перспективным механизмом очистки и предварительной обработки данных для использования их в статистических и бизнес-приложениях, системах искусственного интеллекта и т.п.
5. Системы передачи данных и цифровой обработки сигналов. Благодаря высокой эффективности алгоритмов и устойчивости к воздействию помех, вейвлет-преобразование является мощным инструментом в тех областях, где традиционно использовались другие методы анализа данных, например, преобразование Фурье. Возможность применения уже существующих методов обработки результатов преобразования, а также характерные особенности поведения вейвлет-преобразования в частотно-временной области позволяют существенно расширить и дополнить возможности подобных систем.

И это еще далеко не все!

Заключение

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

Литература

1. Добеши И. Десять лекций по вейвлетам. Москва, "РХД", 2001 г.
2. Воробьев В.И., Грибунин В.Г. Теория и практика вейвлет-преобразования. С.-Петербург, ВУС, 1999 г.
3. Mallat S. A theory for multiresolutional signal decomposition: the wavelet representation. IEEE Trans. Pattern Analysis and Machine Intelligence, 1989, N7, p.674-693.

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

где * - символ комплексной сопряженности и функция ψ - некоторая функция. Функция может быть выбрана произвольным образом, но она должна удовлетворять определённым правилам.

Как видно, вейвлет-преобразование на самом деле является бесконечным множеством различных преобразований в зависимости от оценочной функции, использованной для его расчёта. Это является основной причиной, почему термин «вейвлет-преобразование » используется в весьма различных ситуациях и для различных применений. Также существует множество типов классификации вариантов вейвлет-преобразования. Здесь мы покажем только деление, основанное на ортогональности вейвлетов. Можно использовать ортогональные вейвлеты для дискретного вейвлет-преобразования и неортогональные вейвлеты для непрерывного. Эти два вида преобразования обладают следующими свойствами:

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

Дополнительные подробности о вейвлет-преобразовании доступны на тысячах интернет-ресурсов о вейвлетах в сети, или, например, здесь .

В библиотеке обработки данных Gwyddion реализованы оба этих преобразования и использующие вейвлет-преобразование модули доступны в меню Обработка данных → Интегральные преобразования .

Дискретное вейвлет-преобразование

Дискретное вейвлет-преобразование (DWT) - реализация вейвлет-преобразования с использованием дискретного набора масштабов и переносов вейвлета, подчиняющихся некоторым определённым правилам. Другими словами, это преобразование раскладывает сигнал на взаимно ортогональный набор вейвлетов, что является основным отличием от непрерывного вейвлет-преобразования (CWT), или его реализации для дискретных временных рядов, иногда называемой непрерывным вейвлет-преобразованием дискретного времени (DT-CWT).

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

где S - фактор масштаба (обычно выбирается как 2). Более того, площадь под функцией должна быть нормализована и функция масштабирования должна быть ортогональна к своим численным переносам, т.е.

После введения некоторых дополнительных условий (поскольку вышеупомянутые ограничения не приводят к единственному решению) мы можем получить результат всех этих уравнений, т.е. конечный набор коэффициентов a k которые определяют функцию масштабирования, а также вейвлет. Вейвлет получается из масштабирующей функции как N где N - чётное целое. Набор вейвлетов затем формирует ортонормированный базис, который мы используем для разложения сигнала. Следует отметить, что обычно только несколько коэффициентов a k будут ненулевыми, что упрощает расчёты.

На следующем рисунке показаны некоторые масштабирующие функции и вейвлеты. Наиболее известным семейством ортонормированных вейвлетов явлется семейство Добеши. Её вейвлеты обычно обозначаются числом ненулевых коэффициентов a k , таким образом, мы обычно говорим о вейвлетах Добеши 4, Добеши 6, и т.п. Грубо говоря, с увеличением числа коэффициентов вейвлета функции становятся более гладкими. Это явно видно при сравнении вейвлетов Добеши 4 и 20, представленных ниже. Другой из упомянутых вейвлетов - простейший вейвлет Хаара, который использует прямоугольный импульс как масштабирующую функцию.

Функция масштабирования Хаара и вейвлет (слева) и их частотные составляющие (справа).

Функция масштабирования Добеши 4 и вейвлет (слева) и их частотные составляющие (справа).

Функция масштабирования Добеши 20 и вейвлет (слева) и их частотные составляющие (справа).

Существует несколько видов реализации алгоритма дискретного вейвлет-преобразования. Самый старый и наиболее известный – алгоритм Малла (пирамидальный). В этом алгоритме два фильтра – сглаживающий и несглаживающий составляются из коэффициентов вейвлета и эти фильтры рекуррентно применяются для получения данных для всех доступных масштабов. Если используется полный набор данных D = 2 N и длина сигнала равна L , сначала рассчитываются данные D /2 для масштаба L /2 N - 1 , затем данные (D /2)/2 для масштаба L /2 N - 2 , … пока в конце не получится 2 элемента данных для масштаба L /2 . Результатом работы этого алгоритма будет массив той же длины, что и входной, где данные обычно сортируются от наиболее крупных масштабов к наиболее мелким.

В Gwyddion для расчёта дискретного вейвлет-преобразования используется пирамидальный алгоритм. Дискретное вейвлет-преобразование в двумерном пространстве доступно в модуле DWT.

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

где Y ij соответствует всем коэффициентам наиболее высокого поддиапазона масштаба разложения (где, как предполагается, должна присутствовать большая часть шума). Или же дисперсия шума может быть получена независимым путём, например, как дисперсия сигнала АСМ, когда сканирование не идёт. Для наиболее высокого поддиапазона частот (универсальный порог) или для каждого поддиапазона (для адаптивного по масштабу порога) или для окружения каждого пикселя в поддиапазоне (для адаптивного по масштабу и пространству порога) дисперсия рассчитывается как

Значение порога считается в конечном виде как

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

Удаление шума DWT доступно в меню Обработка данных → Интегральные преобразования → Удаление шума DWT .

Непрерывное вейвлет-преобразование

Непрерывное вейвлет-преобразование (CWT) - реализация вейвлет-преобразования с использованием произвольных масштабов и практически произвольных вейвлетов. Используемые вейвлеты не ортогональны и данные, полученные в ходе этого преобразования высоко коррелированы. Для дискретных временных последовательностей также можно использовать это преобразование, с ограничением что наименьшие переносы вейвлета должны быть равны дискретизации данных. Это иногда называется непрерывным вейвлет-преобразованием дискретного времени (DT-CWT) и это наиболее часто используемый метод расчёта CWT в реальных применениях.

В принципеЮ непрерывное вейвлет-преобразование работает используя напрямую определение вейвлет-преобразования, т.е. мы рассчитываем свёртку сигнала с масштабированным вейвлетом. Для каждого масштаба мы получаем этим способом набор той же длины N , что и входной сигнал. Используя M произвольно выбранных масштабов мы получаем поле N×M , которое напрямую представляет плоскость время-частота. Алгоритм, используемый для этого расчёта может быть основан на прямой свёртке или на свёртке посредством умножения в Фурье-пространстве (это иногда называется быстрым вейвлет-преобразованием).

Выбор вейвлета для использования в разложении на время-частоту является наиболее важной вещью. Этим выбором мы можем влиять на разрешение результата по времени и по частоте.Нельзя изменить этим путём основные характеристики вейвлет-преобразования (низкие частоты имеют хорошее разрешение по частотам и плохое по времени; высокие имеют плохое разрешение по частотам и хорошее по времени), но можно несколько увеличить общее разрешение по частотам или по времени. Это напрямую пропорционально ширине используемого вейвлета в реальном и Фурье-пространстве. Если, например, использовать вейвлет Морле (реальная часть – затухающая функция косинуса), то можно ожидать высокого разрешения по частотам, поскольку такой вейвлет очень хорошо локализован по частоте. наоборот, используя вейвлет Производная Гауссиана (DOG) мы получим хорошую локализацию по времени, но плохую по частоте.

Непрерывное вейвлет-преобразование реализовано в модуле CWT, который доступен в меню Обработка данных → Интегральные преобразования → CWT .

Источники

A. Bultheel: Bull. Belg. Math. Soc.: (1995) 2

S. G. Chang, B. Yu, M. Vetterli: IEEE Trans. Image Processing, (2000) 9 p. 1532

S. G. Chang, B. Yu, M. Vetterli: IEEE Trans. Image Processing, (2000) 9 p. 1522

• Tutorial

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

Попробуем на простых примерах разобраться, откуда же вообще берутся вейвлеты и как их можно использовать при сжатии. Предполагается, что читатель знаком с основами линейной алгебры, не боится слов вектор и матрица, а также умеет их перемножать. (А во даже попробуем что-то запрограммировать.)

Сжатие изображений

Упрощённо, изображение представляют собой таблицу, в ячейках которой хранятся цвета каждого пикселя. Если мы работаем с чёрно-белым (или, точнее, серым) изображением, то вместо цвета в ячейки помещают значения яркости из отрезка . При этом 0 соответствует чёрному цвету, 1 - белому. Но с дробями работать неудобно, поэтому часто значения яркости берут целыми из диапазона от 0 до 255. Тогда каждое значение будет занимать ровно 1 байт.

Даже небольшие изображения требуют много памяти для хранения. Так, если мы кодируем яркость каждого пикселя одним байтом, то изображение одного кадра формата FullHD (1920×1080) займёт почти два мегабайта. Представьте, сколько памяти потребуется для хранения полуторачасового фильма!

Поэтому изображения стремятся сжать. То есть закодировать таким образом, чтобы памяти для хранения требовалось меньше. А во время просмотра мы декодируем записанные в память данные и получаем исходный кадр. Но это лишь в идеале.

Существует много алгоритмов сжатия данных. О их количестве можно судить по форматам, поддерживаемым современными архиваторами: ZIP, 7Z, RAR, ACE, GZIP, HA, BZ2 и так далее. Неудивительно, что благодаря активной работе учёных и программистов в настоящее время степень сжатия данных вплотную подошла к теоретическому пределу.

Плохая новость в том, что для изображения этот теоретический предел не так уж и велик. Попробуйте сохранить фотографию (особенно с большим количеством мелких деталей) в формате PNG - размер получившегося файла может вас расстроить.

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

Алгоритмы сжатия «любят», когда в данных есть закономерность. Лучше всего сжимаются длинные последовательности нулей (закономерность тут очевидна). В самом деле, вместо того, чтобы записывать в память 100 нулей, можно записать просто число 100 (конечно, с пометкой, что это именно количество нулей). Декодирующая программа «поймёт», что имелись в виду нули и воспроизведёт их.

Однако если в нашей последовательности в середине вдруг окажется единица, то одним числом 100 ограничится не удастся.

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

Преобразование Хаара

Итак, наша цель - преобразовать изображение так, чтобы оно хорошо сжималось классическими алгоритмами. Подумаем, как нужно изменить его, чтобы получить длинные цепочки нулей.

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

Рассмотрим фрагмент первой строки яркостей из известного изображения «Lenna» (на рисунке).

154, 155, 156, 157, 157, 157, 158, 156

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

Получаем:

154, 1, 1, 1, 0, 0, 1, -2.

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

Попробуем поступить иначе. Не будем пытаться сразу получить хорошую последовательность, попробуем улучшить её хотя бы немного.

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

(154, 155), (156, 157), (157, 157), (158, 156)
(154.5, 0.5), (156.5, 0.5), (157, 0.0), (157, -1.0)

Почему именно полусуммы и полуразности? А всё очень просто! Полусумма - это среднее значение яркости пары пикселей. А полуразность несёт в себе информацию об отличиях между значениями в паре. Очевидно, зная полусумму a и полуразность d можно найти и сами значения:
первое значение в паре = a - d,
второе значение в паре = a + d.

Это преобразование было предложено в 1909 году Альфредом Хааром и носит его имя.

А где же сжатие?

Полученные числа можно перегруппировать по принципу «мухи отдельно, котлеты отдельно», разделив полусуммы и полуразности:

154.5, 156.5, 157, 157; 0.5, 0.5, 0.0, -1.0.

Числа во второй половине последовательности как правило будут небольшими (то, что они не целые, пусть пока не смущает). Почему так?

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

В самом деле, рассмотрим первые 2000 пар соседних пикселей и каждую пару представим на графике точкой.

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

А теперь рассмотрим график, точками в котором будут полусуммы и полуразности.

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

Применим математику!

Попробуем записать математические выражения, описывающие преобразование Хаара.

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

Такое преобразование описывается матрицей .

В самом деле , что нам и требовалось.

Внимательный читатель наверняка заметил, что рисунки из точек на двух последних графиках одинаковы. Разница лишь в повороте на угол в 45°.

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

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

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

Получаем в итоге матрицу

А как декодировать?

Как известно, если у матрицы определитель не равен нулю, то для неё существует обратная матрица, «отменяющая» её действие. Если мы найдём обратную матрицу для H, то декодирование будет заключаться просто в умножении векторов с полусуммами и полуразностями на неё.

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

Рассмотрим поближе нашу матрицу. Она состоит из двух вектор-строк: и . Назовём их v 1 и v 2 .

Они обладают интересными свойствами.

Во-первых, их длины равны 1, то есть . Здесь буква T означает транспонирование. Умножение вектор-строки на транспонированный вектор-строку - это скалярное произведение.

Во-вторых, они ортогональны, то есть .

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

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

Мы любим ортогональные матрицы!

Увеличиваем число точек

Всё сказанное хорошо работает для двух точек. Но что делать, если точек больше?

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

То есть. исходный вектор просто обрабатывается независимо по парам.

Фильтры

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

Полученные «полусуммы» (из-за того, что делим не на 2, приходится использовать кавычки) - это, как мы уже выяснили, средние значения в парах пикселей. То есть, фактически, значения полусумм - это уменьшенная копия исходного изображения! Уменьшенная потому, что полусумм в два раза меньше, чем исходных пикселей.

Но что такое разности?

Полусуммы усредняют значения яркостей, то есть «отфильтровывают» случайные всплески значений. Можно считать, что это некоторый частотный фильтр.

Аналогично, разности «выделяют» среди значений межпиксельные «всплески» и устраняют константную составляющую. То есть, они «отфильтровывают» низкие частоты.

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

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

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

Степень сжатия можно увеличить, применяя преобразование Хаара многократно. В самом деле, высокочастотная составляющая - это всего лишь половина от всего набора чисел. Но что мешает применить нашу процедуру ещё раз к низкочастотным данным? После повторного применения, высокачастотная информация будет занимать уже 75%.

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

После многократного применения к, например, фотографии замка Лихтенштейн, получим следующий рисунок.

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

Этот процесс называется квантованием. И именно на этом этапе происходит потеря части информации. (К слову, такой же подход используется в JPEG, только там вместо преобразования Хаара используется дискретное косинус-преобразование.) Меняя число обнуляемых коэффициентов, можно регулировать степень сжатия!

Конечно, если обнулить слишком много, то искажения станут видны на глаз. Во всём нужна мера!

После всех этих действий у нас останется матрица, содержащая много нулей. Её можно записать построчно в файл и сжать каким-то архиватором. Например, тем же 7Z. Результат будет неплох.

Декодирование производится в обратном порядке: распаковывем архив, применяем обратное преобразование Хаара и записываем декодированную картинку в файл. Вуаля!

Где эффективно преобразование Хаара?

Когда преобразование Хаара будет давать наилучший результат? Очевидно, когда мы получим много нулей, то есть, когда изображение содержит длинные участки одинаковых значений яркости. Тогда все разности обнулятся. Это может быть, например, рентгеновский снимок, отсканированный документ.

Говорят, что преобразование Хаара устраняет константную составляющую (она же - момент нулевого порядка), то есть переводит константы в нули.

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

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

Преобразование Добеши

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

То есть, если исходная последовательность - 1, 2, 3, 4, 5, 6,…, N-1, N, то будем брать четвёрки (1, 2, 3, 4), (3, 4, 5, 6) и т. д. Последняя четвёрка «кусает последовательность за хвост»: (N-1, N, 1, 2).

Точно так же попробуем построить два фильтра: высокочастотный и низкочастотный. Каждую четвёрку будем заменять на два числа. Так как четвёрки перекрываются, то количество значений после преобразования не изменится.

Пусть значения яркостей в четвёрке равны x, y, z, t. Тогда первый фильтр запишем в виде

Четыре коэффициента, образующих вектор-строку матрицы преобразования, пока нам неизвестны.

Чтобы вектор-строка коэффициентов второго фильтра был ортогонален первому, возьмём те же коэффициенты но переставим их и поменяем знаки:

Матрица преобразования будет иметь вид.

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

Векторы должны иметь единичную длину (иначе определитель будет не единичным):

Преобразование должно обнулять цепочку одинаковых значений (например, (1, 1, 1, 1)):

Преобразование должно обнулять цепочку линейно растущих значений (например, (1, 2, 3, 4)):

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

Получили 4 уравнения, связывающие коэффициенты. Решая их, получаем:

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

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

Это преобразование получило название вейвлета D4 (читателю предлагается самостоятельно разгадать тайну этого буквенно-цифрового названия).

Другие вейвлеты

Мы, конечно, можем не остановиться на этом, и потребовать устранения параболической составляющей (момент 2-го порядка) и так далее. В результате получим вейвлеты D6, D8 и другие.

Чтобы не считать всё вручную, коэффициенты можно посмотреть в википедии .

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

Домашнее задание

Чтобы окончательно разобраться с основами, предлагаю написать на вашем любимом языке программу, которая открывает изображение, выполняет преобразование Хаара (или даже D4), квантует результат, а потом сохраняет результат в файл. Попробуйте сжать этот файл своим любимым архиватором. Хорошо сжимается?

Попробуйте выполнить обратное преобразование. Как вы объясните характер артефактов на изображении?

Заключение

Конечно, в этой статье не были рассмотрены очень многие интересные математические детали и практические применения вейвлет-преобразований. Но нельзя объять необъятное. Да и многое сложно объяснить не повышая градус матана. Надеюсь, что и написанное оказалось кому-то полезным. Добавить метки

Известно, что произвольный сигнал , для которого выполняется условие может быть представлен ортогональной системой функций :

, (18)

коэффициенты определяются из соотношения

,

где - квадрат нормы или энергия базисной функции . Ряд (18) называется обобщенным рядом Фурье. При этом произведения вида , входящие в ряд (18), представляют собой спектральную плотность сигнала , а коэффициенты - спектр сигнала. Суть спектрального анализа сигнала состоит в определении коэффициентов . Зная эти коэффициенты возможен синтез (аппроксимация) сигналов при фиксированном числе ряда:

.

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

.

Известные преобразования (Адамара, Карунена-Лоэва, Фурье) «плохо» представляют нестационарный сигнал в коэффициентах разложения. Покажем это на следующем примере. Пусть дана нестационарная функция

и ее преобразование Фурье (рис. 9).

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

Сложно провести анализ временного сигнала по его Фурье образу;

Приемлемая аппроксимация временного сигнала возможна при учете большого числа высокочастотных коэффициентов;

Плохое визуальное качество реальных изображений восстановленных по низкочастотным коэффициентам; и т.п.

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

Рис. 9. Преобразование Фурье синусоидального сигнала с небольшими ступеньками при переходе через нуль

ВП одномерного сигнала – это его представление в виде обобщенного ряда Фурье или интеграла Фурье по системе базисных функций локализованных как в пространственной, так и в частотной областях. Примером такой базисной функции может служить вейвлет Хаара, который определяется выражением

(20)

Графически вейвлет Хаара представляется следующим образом:

Рис. 10. Базисная функция вейвлета Хаара

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

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

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

получим следующую аппроксимацию:

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

Рис. 11. Преобразование пар случайных величин

В результате исходный сигнал точно описывается коэффициентами вейвлет-преобразования Хаара. Вейвлет-коэффициенты сигнала (19) показаны на рис. 10.

Из приведенного рисунка видно, что нестационарности сигнала (резкие перепады) локализуются в малом числе вейвлет-коэффициентов. Это приводит к возможности лучшего восстановления нестационарного сигнала по неполным данным.

Рис. 12. Вейвлет-коэффициенты одного периода функции (19)

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

,

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

.

Прямое и обратное дискретные ВП вычисляются по формулам

,

.

Следует отметить, что если число отсчетов , то максимальное значение равно . Наибольшее значение для текущего равно .

Для непрерывных сигналов будут справедливы следующие интегральные выражения:

,

.

Таким образом, задавая вейвлет-функции, можно выполнять разложение сигнала по вейвлет-базису непрерывных или дискретных сигналов.

Рис. 13. Распределение базисных функций Хаара при анализе сигнала

Функция может образовывать вейвлет-базис, если она удовлетворяет следующим условиям:

1. Ограниченность нормы:

.

2. Вейвлет-функция должна быть ограничена и по времени и по частоте:

и , при .

Контрпример: дельта-функция и гармоническая функция не удовлетворяют данному условию.

3. Нулевое среднее:

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

В качестве примера приведем следующие известные вейвлет-функции:

, .

Для ВП, также как и для ДПФ существует алгоритм быстрого преобразования. Рассмотрим снова ВП Хаара. Из рис. 13 видно, что функции с малым масштабным коэффициентом используют те же отсчеты сигнала для вычисления коэффициентов, что и функции с большим масштабным коэффициентом. При этом операция суммирования одних и тех же отсчетов повторяется неоднократно. Следовательно, для уменьшения объема вычислений целесообразно вычислять ВП с самого малого масштабного коэффициента. В результате получаем вейвлет-коэффициенты, представляющие собой средние значения и разности . Для коэффициентов повторяем данную процедуру. При этом усреднение коэффициентов будет соответствовать усреднению четырех отсчетов сигнала, но при этом расходуется одна операция умножения и одна операция сложения. Процесс разложения повторяется до тех пор, пока не будут вычислены все коэффициенты спектра .

Запишем алгоритм быстрого вейвлет-преобразования Хаара в матричном виде. Пусть дан вектор размером 8 элементов. Матрица преобразования Хаара запишется в виде