УДК 681.3
РАЗВИТИЕ МЕТОДОВ БЫСТРОГО ВЕЙВЛЕТ-ПРЕОБРАЗОВАНИЯ С ПОМОЩЬЮ ФИЛЬТРОВ ДОБЕШИ
Н. И. ЧЕРВЯКОВ Ю. В. КОНДРАШОВ
Ставропольский военный институт связи ракетных войск
e-mail: yrii-fifa83@mail.ru
В работе показана возможность построения высокоэффективных вейвлет — преобразований ортогональными вейвлет — фильтрами Добеши по алгоритму Малла на базе свертки, вычисленной в системе остаточных классов.
Ключевые слова: вейвлет — преобразование, вейвлет — фильтр, вейвлет — функция, алгоритм Малла, свертка, цифровая обработка сигналов, система остаточных классов.
В настоящее время вейвлет — преобразование широко применяется в задачах обработки и кодирования сигналов и изображений самой различной природы (речь, спутниковые изображения, рентгенограммы внутренних органов и др.), распознавания образов, при изучении свойств поверхностей кристаллов и нанообъектов и во многих других случаях.
Появление в 1988 году ортогональных вейвлет-фильтров Добеши или так называемых фильтров с компактным носителем в значительной мере усилило интерес к вейвлет-анализу, т.к. открылись новые возможности не только для теоретического, но и для практического применения вейвлет-преобразования.
Важно отметить то, что вейвлет-фильтры Добеши строятся, исходя из критерия длины фильтров и, следовательно, являются фильтрами с конечным числом коэффициентов [1]. Вейвлет-функции /(t ) фильтров Добеши принято обозначать литерой D с добавлением цифры, соответствующей длине вейвлет-фильтра Добеши, т.е. D2, D4, D6 и т.д.
Пусть даны два фильтра h и g с ненулевыми элементами [2]:
, 1 + >/з , 3 + >/з , 3-л/3 , 1-V3
h0 =----------------, h1 =-, h2 =---, h3 =----------------------------------; (1)
0 8 1 8 2 8 3 8 (l)
1 -л/3 3-л/3 3 + л/3 1 + л/3
g 0 = 8 , g 1 = 8 ’ g 2 = 8 ’ g 3 = 8 '
Отметим соотношения между коэффициентами этих фильтров [2]:
ho + h + h2 + h3 = 1;
g 0 + g1 + g 2 + g 3 =0; (2)
g 0 = h3, g1 = — h2 , g 2 = h1, g 3 = — h0; h0 + h2 = h + h3 = 1/2, 2 h2 = h + 3h3.
Найдём передаточные функции H(z) и G(z) в z -представлении:
H ( z) = h0 + h^z + h2 z2 + h3 z3; (3)
G(z) = g0 + gz + g2z2 + gz3 = h3 - h2z + V2 - h0z3 = (4)
= - z3 (-h3 z ~3 + h2 z-2 - hz- + h0 ) = - z 3H (-z _1).
Таким образом, мы получили
G ( z ) = - z3 H (-z-1). (5)
Для восстановления сигнала требуются дополнительные фильтры Н (2) и G(г) . Определим их как сопряжённые квадратурные фильтры по формулам:
Н (г) — Н (г“1), О (г) — О( г~1) = -г“3 Н (-г) (6)
Тогда второе соотношение Н (г)Н (-г) + О(г)О(-г) = 0 из (3-4) выполняется. Первое соотношение Н ( г )Н ( г) + О ( г )О( г) = 1 принимает вид:
Н (г~1) Н (г) + О( г-1 )О( г) = 1. (7)
Вернёмся к частотной переменной г = е 1С0. Поскольку коэффициенты {Пи} —
вещественные, то Н(г~') = Н(а). Поэтому последнее соотношение принимает вид
I I? I I?
\Н (а)|2 + |О (а)|2 — 1.
Найдём коэффициенты фильтров восстановления Н (г) и (О (г) из их определения Н (г) = Н (г_1), О(г) = О(г_1) :
г 1 — л/3 г 3 - л/3 Г 3 + л[3 1 + у/3
¡1 3 — , ¡1 2 — 5 П 1 — 5 По — ; / О \
-3 & -2 & -1 & 0 8 (8)
_ 1 + >/3 _ _3^/3 _ 3 - л/э _ 1-л/э
<^?-з 8 5 <^?-2 & 5 <?-1 & 5 <?° & •
Таким образом, метод одномерного дискретного вейвлет — преобразования (ДВП) N -го порядка последовательности Хп определяется следующими рекуррентными соотношениями:
N— 1
аП) = 1 2Л i = 1,2,..., J;
k = 0 N-1
(9)
/(i) = У h,a2i_ 1) a(0) = х ,
n / у k 2n-k n n '
k = 0
где а^} и ^^(i) являются аппроксимирующими и детализирующими коэффициентами i -го уровня, а gk и hk (k = 0,1,.,N — 1) — коэффициенты низкочастотного и высокочастотного анализирующих фильтров, соответственно.
С другой стороны, сигнал Хп может быть восстановлен по коэффициентам
janJ), d(n J), d(n J 1),., | путём последовательной итерации по формулам:
a{i—1) =
m
Zg^.a(i) + У h2kd(i) , m чётное
S2k m-k^ ¿k m—k (10)
2k m
k=0 2 ^ k=0 2
N /2—1__ N/2—1_ _
У g 2k+1am—1 —, + У h2k+1dm—1 , m нечётное
k=0 “2 k k=0 “2 k
где gk и являются коэффициентами низкочастотного и высокочастотного синтезирующих фильтров, соответственно.
Для того, чтобы восстановленный сигнал соответствовал исходному, должны быть соответствующим образом подобраны анализирующий (раскладывающий) и синтезирующий (собирающий) фильтры.
Для вейвлет-преобразования функции /(х) необходимо вычислить серию
коэффициентов |ап5йп5йп^5...5^1|, где ап - аппроксимация функции, di - детази-
зирующие коэффициенты функции, i — 15.5п . Каждый коэффициент находится интегрированием (11, 12):
N /2 1
N /2 1
114
НАУЧНЫЕ ВЕДОМОСТИ
№ 15(70) 2009
а
а
J - т ,к
1 (I, VJ - N ,к ) = | I( Х - N ,к ( Х )аХ ;
Я
( 1 У J - т ,к )= | 1 ( Х )/ J - т ,к ( Х )ах , т = 1,2,-‘- , N ■
(11)
(12)
Возникает проблема вычисления большого количества интегралов с необходимой точностью. Следует также учитывать, что при высоком уровне разрешения J носители функций (р^1к (х) и у к (х) становятся малыми порядка } .
Быстрое вейвлет-преобразование, предложенное Мала позволяет решить эту проблему. Алгоритм Малла даёт возможность вычислять коэффициенты вейвлет-разложения без интегрирования, используя алгебраические операции на основе свёртки:
N-1
3П‘) = Е ёк^П-к
(13)
,(¡-1)
где а() и
ми
) являются аппроксимирующими и детализирующими коэффициента-
I -го уровня, а и Ик (к = 0,1,...,N - 1) - коэффициенты низкочастотного и
высокочастотного анализирующих фильтров, соответственно; Хп — исходный сигнал;
N — порядок фильтра.
Эти равенства обеспечивают быстрые алгоритмы вычисления вейвлет-коэффициентов (каскадные алгоритмы, алгоритмы Малла). Термин «быстрые» означает не только, что в (13) используются более быстрые алгебраические процедуры, но и то, что при каждом преобразовании общее число новых коэффициентов не увеличивается в два раза, а остаётся прежним.
Схема разложения сигнала по алгоритму Малла приведена на рис. 1.
Рис. 1. Последовательность получения вейвлет-коэффициентов третьей октавы; Н (2) и 0( г ) , соответственно, высокочастотные и низкочастотные анализирующие фильтры
в г -представлении
Единственное отличие фильтрации в алгоритме Малла от классического КИХ-фильтра, задаваемого уравнением ^) = ь х(к - ¡) [3], заключается в том, что значе-
¡=0
ния фильтруемого ряда выбираются через один — индекс 2п - к в а2П-). Это и есть
децимация 21 — исключение из обработки каждого второго элемента.
Для двумерных сигналов — изображений — алгоритм разложения аналогичен тому, что применяется в одномерном случае (13). Пусть (р(х) — масштабирующая
вейвлет-функция и / (х) — материнский вейвлет. Как известно, они порождают ба-
к=0
к=0
зисные функции п (х) и к (х). Двумерный сигнал а(п1, п2) раскладывается по
базисным в Ь (К ) функциям PJ,п (х)PJ, т (У ) ' PJ,п (х)/J, т (У ) ' /J,п (x)PJ, т (У ) и /J п (х)/_; т (У) . Соответствующие коэффициенты принято называть следующим образом.
Аппроксимирующие коэффициенты а(3 )(п1, п2) получаются как коэффициенты разложения по вейвлет-базису рJ п (х)рJ т (у). На рис. 2 (а) показано распределение пикселов после пошаговой обработки исходного изображения банком фильтров.
Горизонтальные детализирующие коэффициенты ¿23)(п1, п2) получаются как коэффициенты разложения по вейвлет-базису рJ п (х)/J т (у) .
Вертикальные детализирующие коэффициенты ¿13)(п1, п2) получаются как коэффициенты разложения по вейвлет-базису /J п (х)р_1 т (у) .
Диагональные детализирующие коэффициенты ¿33)(п1, п2) получаются как коэффициенты разложения по вейвлет-базису /J п (х)у_/ т (у) .
Схема разложения сигнала а°(п1, п2) изображена на рис. 2 (б).
а(1)(п1, и2) и2)
а<0)(п1, и2) —
d¡1)(n¡, и2) йз(1)(п, пг)
а(2) (п1, п. й2(2)(п, п) п2)
^Ц, п2) íз<2>(n, п2)
^Ч^ п2) йз(1)(п1, ^2 )
(а)
роды
колонки
ряды ,
“I
<'(£} |—— -@ъ
[МЬ@^СЦ
т=) Н@)-| *■
I--► (Л,, Яг)
(б)
Рис. 2. Последовательность получения вейвлет-коэффициентов третьей октавы для двумерного сигнала: (а) —распределение пикселов после пошаговой обработки исходного изображения банком фильтров; (б) — в виде последовательности фильтров;
Н (г) и С(г) , соответственно, высокочастотные и низкочастотные
анализирующие фильтры в г -представлении
В аналитическом виде разложение двумерного сигнала фильтрами можно записать следующим образом:
116
a(1+1)(n,,n
НАУЧНЫЕ ВЕДОМОСТИ
N-1 N-1
(ni,n2) = Z Z g(ki)g(k2)a(')(2П1 -ki,2n2 - к2);
№ 15(70) 2009
к! =0 к2 =0 N 1N 1
d(‘+!)(ni,n2) = Z Z g(ki)h(к2)a(')(2n -ki,2n2 -к2);
к1 =0 к2 =0 N-1 N-1
d2,+1)(ni,n2) = Z Z h(к\)g(к2)a(’)(2ni -^,2n2 -к2);
к1 =0 к2 =0
N-1 N-1
J3(I+1)(n1,n2) = Z Z h(к1)h(к2)a(i)(2n1 - к1,2n2 - к2).
(14)
В качестве собственно фильтров могут использоваться фильтры Добеши D4 четвёртого порядка. Вейвлеты Добеши являются вейвлетами с компактным носителем, что обеспечивает хорошие свойства приближения вейвлет-разложений. Они не имеют эксплицитного (явного) выражения, а задаются коэффициентами фильтрации. Анализирующие (разлагающие) высокочастотные (^ и низкочастотные ^) коэффициенты фильтра Добеши D4 задаются следующими коэффициентами [2]:
g 0 =
1 - л/3
4V2 :
hi =
g 1 =
3 + Уз 4V2 , з - л/з 4V2 :
h 2 =
g 2 =
з-Уз 4V2 , з + л/з 4V2 ,
, 1 - л/з"
hз = “WT;
(15)
g з =
1 + л/з
4V2 .
Графики вейвлетов Добеши D4 (db4) в среде MATLAB можно увидеть следующим образом (рис. 3):
[phi,psi,x]=wavefun('db4',10); subplot(121); plot(x,phi); title('y=\phi(x)'); axis square; grid on;
subplot(122); plot(x,psi); title('y=\psi(x)'); axis square; grid on;
Рис. 3. Масштабирующая вейвлет-функция и материнский вейвлет Добеши Э4
С целью повышения скорости вычисления свертки (13) предлагается её вычислять в системе остаточных классов, тогда выбирая модуль р свертка может быть выражена как:
a(‘) = n N-1, Z g a (i-1) Аи^к^12n-k|
к=0 p
N - 1
d(1) = n Z\KaT | к 2n- к к=0
i=1,2,..., J,
(16)
a
(0)
X .
n
к =0 к =0
12
p
p
Pi
Система остаточных классов и модулярные вычисления являются практически идеальным инструментом реализации линейной свертки, поскольку операции сложения, вычитания и умножения выполняются очень просто, а именно, если даны два числа A и B, представленные в системе остаточных классов (с набором взаимно простых оснований т, т2,...,тС) следующим образом:
A = (a1,a2,..., aL) : A = a1(mod m1), A = a2(mod m2),..., A = aL (mod mL),
B = (b1,b2,..., bL): B = b1(mod m 1), B = b2(mod m2),..., B = bL (mod mL) (17)
то
A±B = ^a2,...0L)±(bi,b2,...bL) = (a ±bim),(a ±¿21 m),...(aL Щщ) (18)
и
A-B=a ^...aO-(bl,b2,... bL)=( a-bi| m),(a-¿J щ),...(a \\щ) (19)
Математические модели (17-19) вычисляются на основе использования нейронных сетей конечного кольца [3], число которых определяется рядом каналов по числу оснований, работающих независимо друг от друга и параллельно во времени. Если каждую нейронную сеть конечного кольца отожествить с отдельным основанием системы остаточных классов, то образованная совокупность каналов будет представлять собой арифметическое устройство выполняющее с большой эффективностью вейвлет — преобразование сигналов.
Итак, система остаточных классов является наиболее подходящей технологией для реализации высокоэффективного вейвлет - преобразования для задач цифровой обработки сигналов.
Для эффективной реализации операций вейвлет - преобразования по алгоритму Малла на основе свертки предлагается использовать математическую модель вычислительного объекта, оперирующую числами, представленными в системе остаточных классов.
Литература
1. Добеши И. Десять лекций по вейвлетам. - М.: Ижевск: РХД, 2001.
2. Daubechies I. The Wavelet Transform, Time-Frequency Localization and Signal Analysis // IEEE Trans. Inform. Theory, 1990, № 5. P. 961-1005.
3. Червяков Н. И., Сахнюк П. А., Шапошников А. В., Макоха А. Н. Под редакцией А. И. Галушкина. Учебное пособие для ВВУЗов. - М.: Радиотехника, 2003. - 272 с.
4. Сергиенко А.Б. Цифровая обработка сигналов. - СПб.: Питер, 2003. - 604 с.
5. Астафьева Н.М. Вейвлет-анализ: основы теории и некоторые приложения // Успехи физических наук, 1996, № 11. С. 1145-1170.
6. Goswami J.C., Chan A.K. Fundamentals of Wavelets. Theory, Algorithms, and Applications. Wiley, 2000. - 306 p.
7. Прокис Дж. Цифровая связь. - М.: Радио и связь, 2000. - 800 с.
8. Желудев В.А. О вейвлетах на базе периодических сплайнов / / Докл. РАН, 1994, № 1. С. 9-13.
9. Новиков Л.В. Основы вейвлет-анализа сигналов. Учебное пособие. - СПб., 1999. - 152 с.
DEVELOPMENT OF METHODS FAST WAVELET-TRANSFORMATION BY MEANS OF FILTERS DOBESHI
N. I. CHERVYKOV Y. V. KONDRASHOV
Stavropol military institute of communication of rockets armies
e-mail: yrii-fifa83@mail.ru
In work possibility of construction highly effective wavelet - transformations orthogonal wavelet - filters Dobeshi on algorithm of Mull on the basis of the convolution calculated in system of residual classes is shown.
Key words: wavelet - transformation, wavelet - the filter, wavelet -functions, algorithm of Mull, convolution, system of residual classes, digital processing of signals.