WavUtils - библиотека инструментов вейвлет-преобразования
WavUtils - библиотека инструментов вейвлет-преобразования
Вейвлет-преобразование является инструментом кратномасштабного анализа: с его помощью сигнал представляется в виде некоторой аппроксимации и набора деталей, отличающих данную аппроксимацию от исходного сигнала. Такое представление позволяет "отсеять" незначительные изменения сигнала, сосредоточившись на изучении его глобальных особенностей, либо, напротив, рассмотреть его "мелкие детали", неразличимые на фоне "крупномасштабных" изменений.
Кроме того, поскольку шум представляет собой, как правило, незначительные высокочастотные ("мелкомасштабные") изменения сигнала, вейвлет-преобразование можно использовать для очистки сигнала от шума путем удаления из него таких "незначительных" деталей.
В предыдущих статьях на нашем сайте был дан краткий обзор вейвлет-анализа и описан алгоритм выполнения дискретного вейвлет-преобразования. Теперь настало время представить инструментарий, реализующий на практике некоторые возможности использования вейвлетов.
Модуль WavUtils предоставляет собой библиотеку инструментов для разложения сигнала при помощи вейвлет-преобразования, очистки сигнала от шума путем удаления незначащих деталей и восстановления сигнала по его разложению.
В качестве инструмента вейвлет-преобразования были выбраны вейвлеты Добеши – как наиболее хорошо изученные и удобные в вычислительном плане. Использование вейвлетов до 10 порядка включительно представляется достаточным для решения подавляющего большинства прикладных задач, связанных с цифровой обработкой сигналов.
Модуль WavUtils экспортирует функции разложения (Decompose) и восстановления (Reconstruct) сигнала, а также очистки сигнала от шума путем удаления несущественных деталей (Denoise).
Обработка ошибок может быть реализована при помощи вызова исключений или путем возврата кода завершения функции – в зависимости от определения символа UseExceptions. TVector – базовый тип для представления данных
Поскольку библиотеке WavUtils ориентирована на обработку одномерных дискретных сигналов (векторов данных), базовым типом в ней является TVector – одномерный динамический массив чисел с плавающей точкой. В библиотеке тип Tvector описан следующим образом:
type
TVector = array of Double;
Однако, это определение может быть изменено, например, для адаптации к требуемой точности и скорости вычислений. Поскольку все последующие определения базируются на типе TVector, для изменения точности вычислений достаточно внести следующие изменения:
TVector = array of Single (ухудшает точность, но повышает скорость вычислений)
или
TVector = array of Extended (наивысшая точность, но крайне неоптимальный код).
Кроме того, тип TVector может быть переопределен как псевдоним существующего типа для обеспечения совместимости с уже имеющимся кодом, например:
uses TypesBase;
type TVector = TypesBase.TVectorFloat;
В целом, общее требование к типу TVector – он должен представлять собой одномерный динамический массив действительных чисел. TDecomposition – информация о разложении сигнала
Как уже было отмечено, сущность вейвлет-преобразования заключается в разложении сигнала на некоторую его аппроксимацию и набор деталей, отличающих эту аппроксимацию от исходного сигнала. Как правило, считают, что аппроксимация несет достаточно полезной информации об исходном сигнале, а детали представляют собой шум и могут быть отброшены. Однако, даже аппроксимация сигнала может еще содержать значительное количество шума. В этом случае процесс вейвлет-разложения применяют уже к аппроксимации сигнала – и повторяют этот процесс до тех пор, пока вся незначащая информация не будет выделена в наборы деталей. После этого коэффициенты детализации могут быть обработаны известными статистическими методами, или даже попросту отброшены. Затем должно быть произведено восстановление сигнала – однако для этого необходимо иметь сведения о том, каким образом выполнялось разложение. Эта информация хранится в структуре TDecomposition, определенной следующим образом:
type
TDecomposition = record WaveletOrder: Integer; SignalLength: Integer; Approx: TVector; Details: array of TVector; end;
Рассмотрим эту структуру подробнее.
Для восстановления сигнала должен применяться вейвлет, соответствующий использованному для разложения. Поскольку в библиотеке WavUtils используются вейвлеты Добеши, для восстановления сигнала достаточно знать только порядок вейвлета (т.е. количество его коэффициентов). Это значение хранится в поле WaveletOrder. Заметьте, что, поскольку вейвлеты Добеши всегда имеют четное число коэффициентов, их количество в два раза больше, чем порядок вейвлета. То есть, вейвлет Добеши 1-го порядка имеет 2 коэффициента, 2-го порядка – 4 коэффициента и т.д.
Коэффициенты аппроксимации и детализации хранятся в полях Approx и Details соответственно. При этом, поскольку на каждом шаге преобразования мы выполняем разложение очередной аппроксимации, достаточно хранить только коэффициенты последней аппроксимации, хотя детальные коэффициенты должны сохраняться для каждого уровня разложения. Тогда восстановление сигнала в обратном порядке позволит шаг за шагом получать все более точные аппроксимации реконструируемого сигнала – и в конце концов сам этот сигнал. Отметим, что глубину разложения (т.е. количество уровней, или итераций вейвлет-преобразования) можно отдельно не сохранять – эта величина соответствует размеру массива Details – ведь коэффициенты детализации сохраняются для каждого уровня разложения.
Дискретное вейвлет-преобразование обладает свойством децимации в 2 раза – т.е. получаемые вектора коэффициентов аппроксимации и детализации короче исходного вектора в 2 раза. При этом, если исходный вектор имеет нечетную длину, он должен быть дополнен до четного числа отсчетов. Кроме того, для полного вычисления свертки сигнала и вейвлет-фильтра, исходный вектор должен быть продолжен соответствующим образом (по крайней мере, на длину фильтра). В библиотеке используется метод периодического дополнения – т.е. считается, что за последним отсчетом сигнала вновь следует первый, второй и т.д. отсчеты сигнала. Если число отсчетов нечетное, то между последним значением сигнала и первым значение дополнения вставляется отсчет, равный их среднему арифметическому – это позволяет снизить влияние краевых искажений для непериодических сигналов. Однако, это приводит к неопределенности при восстановлении сигнала. Конкретно, неизвестно, какова должны быть его длина – четная (в два раза больше длины векторов соответствующих коэффициентов) или нечетная (на единицу меньше, т.е. последний отсчет был вставлен искусственно и может быть отброшен). На промежуточных уровнях разложения длина очередного вектора принимается равной длине вектора детализации следующего уровня – поскольку длина векторов аппроксимации (синтезируемого) и детализации (сохраненного) должна быть равной. Однако на последнем шаге восстановления это значение неизвестно, и длина исходного сигнала также хранится в структуре разложения (поле SignalLength).
Таким образом, получив структуру разложения сигнала, при его последующей обработке допускается изменять значения коэффициентов аппроксимации и детализации, однако категорически недопустимо изменять значения полей WaveletOrder и SignalLength или размерности массивов Approx и Details – в противном случае правильная реконструкция сигнала будет невозможна. Разложение и восстановление сигнала – функции Decompose и Reconstruct
Функции Decompose и Reconstruct используются для разложения и восстановления сигнала – и работают "в паре".
Функция Decompose принимает на входе вектор значений исходного сигнала (параметр Signal), порядок вейвлета (параметр Order) и глубину разложения (параметр Level) и возвращает структуру разложения данного сигнала с указанными параметрами (как значение функции либо через переменную Decomp – в зависимости от способа обработки ошибок). При этом могут возникать два вида ошибочных ситуаций: Указан неправильный порядок вейвлета. Поскольку в библиотеке реализованы вейвлеты Добеши порядка от 1 до 10, параметр Order может принимать значения только из этого диапазона. Длина сигнала меньше длины фильтра. Поскольку на каждом шаге разложения количество коэффициентов уменьшается вдвое, при большой глубине разложения количество коэффициентов аппроксимации может стать меньше, чем количество коэффициентов фильтра – а это недопустимо. При возникновении такой ситуации надо или уменьшить глубину разложения, или использовать более короткие фильтры (вейвлеты меньшего порядка). Вообще, не рекомендуется использовать глубину разложения (значение параметра Level) большие 9. Разумеется, отрицательные значения этого параметра также недопустимы.
Функция Reconstruct, в свою очередь, принимает на входе структуру разложения Decomp и выполняет восстановление сигнала по хранящейся в ней информации. Реконструированный сигнал возвращается либо как значение функции, либо в переменной SynthSig – опять же в зависимости от способа обработки ошибок. При нарушении указанных выше (в описании типа TDecomposition) ограничений при выполнении функции могут возникать ошибочные ситуации, связанные с неправильным значением порядка вейвлета или несоответствием длины синтезированного сигнала длине исходного. Функция Denoise – очистка сигнала от шума
Функция Denoise выполняет очистку сигнала от шума путем удаления незначащих деталей. Входные параметры функции – такие же, как и у Decompose (исходный сигнал, порядок вейвлета и глубина разложения), выходное значение – синтезированный сигнал (как у Reconstruct). Алгоритм работы функции достаточно прост: выполняется разложение исходного сигнала с указанными параметрами, затем все коэффициенты детализации обнуляются и выполняется реконструкция сигнала. Эксперименты показали, что такой метод хорошо работает со многими реальными данными, в том числе и сигналами, имеющими нестационарные характеристики шума, однако оптимальные значения параметров разложения зависят от характеристик сигнала и должны подбираться экспериментально. Дадим следующие рекомендации: Глубина разложения (параметр Level) определяет "масштаб" отсеиваемых деталей: чем больше эта величина, тем более "крупные" изменения сигнала будут отброшены. При достаточно больших значениях параметра Level (порядка 7-9) выполняется не только очистка сигнала от шума, но и его сглаживание ("обрезаются" пики сигнала). И не забывайте, что использование слишком больших значений глубины разложения может привести к возникновению ошибки "Длина сигнала меньше длины фильтра". Порядок вейвлета определяет гладкость восстановленного сигнала: чем меньше значение параметра Order, тем ярче будут выражены "пики" сигнала, и наоборот – при больших значения параметра Order пики сигнала будут сглажены. Для первого эксперимента рекомендуется выбрать значение параметров Order и Level равные 3, а затем изменять их, в соответствии с приведенными рекомендациями, до достижения наилучшего результата, т.е. пока детали разложения будут содержать шумоподобную составляющую, а аппроксимация – достаточно хорошо описывать исходный сигнал. Обработка ошибок: символ UseExceptions и тип TWaveletResult
Рекомендуемым способом обработки ошибочных ситуации при программировании на Delphi является вызов исключений (класс Exception и его потомки). Однако, в некоторых случаях, желательно, чтобы функция не вызывала исключение, а просто возвращала некоторое значение, сигнализирующее об ошибке. Библиотека WavUtils позволяет использовать оба метода обработки ошибок – в зависимости от определения символа условной компиляции UseExceptions. Если символ определен (непосредственно после заголовка модуля), то функции библиотеки будут возвращать своих выходные значения в качестве результата и генерировать исключения при возникновении ошибок. Однако, если символ UseExceptions не определен, то выходное значение будет передаваться через переменную, указанную при вызове функции, а ее результатом будет значение типа TWaveletResult:
type
TWaveletResult = ( wavSuccess, // Удачное завершение операции wavErrorWrongOrder, // Указан неправильный порядок вейвлета wavErrorSignalTooShort, // Длина сигнала меньше длины фильтра wavErrorWrongLength // Длина восстановленного сигнала не соответствует длине исходного
);
Справка и демонстрационный пример
Библиотека WavUtils доступна на сайте BaseGroup.ru в виде Delphi-модуля. Там же находится справка по библиотеке в формате Windows Help и демонстрационный пример, позволяющий загружать данные из текстового файла, выполнять их вейвлет-преобразование с указанными параметрами и синтезировать сигнал, очищенный от шума описанным методом.
