Применение теории вейвлетов при сжатии и фильтрации геоинформации

УДК 528.2, 004:528

ПРИМЕНЕНИЕ ТЕОРИИ ВЕЙВЛЕТОВ ПРИ СЖАТИИ И ФИЛЬТРАЦИИ ГЕОИНФОРМАЦИИ

А.С.ЯРМОЛЕНКО1, О.В.СКОБЕНКО2

1 Новгородский государственный университет им. Ярослава Мудрого, Великий Новгород, Россия

2 Белорусская государственная сельскохозяйственная академия, г. Горки, Могилевская область, Республика Беларусь

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

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

Ключевые слова: базис вейвлет-преобразования; свертка; КМА-анализ; фильтрация; сжатие; точность; программа; изображение; поле данных

Как цитировать эту статью: Ярмоленко А.С. Применение теории вейвлетов при сжатии и фильтрации геоинформации / А.С.Ярмоленко, О.В.Скобенко // Записки Горного института. 2018. Т.234. С. 612-623. DOI: 10.31897/PMI.2018.6.612

Введение и постановка задачи. Сжатие геоинформации и ее фильтрация являются актуальной проблемой теории математической обработки как геодезической информации (измерений), так и изображений. При этом необходимо максимально использовать полученную информацию и получать конечные результаты с достаточной точностью и минимальными затратами при ее хранении, что связано со сжатием информации. Это важно при работе с геоинформацией в геодезии, землеустройстве, природообустройстве, мониторинге земель, при ведении точного сельского хозяйства [20, с.200]. При этом алгоритмы должны обладать простотой и точностью вычислений.

В настоящее время для моделирования объектов сжатия и фильтрации информации широко применяются преобразования Фурье [6, 9-11, 16-19, 23]. Однако даже в своем быстром варианте (быстрое преобразование Фурье - БПФ) оно сопряжено с большим количеством вычислений. В свою очередь в работах [11, 18] отмечается, что в отличие от преобразований Фурье возможны преобразования с другими базисами, восстанавливающие дискретные и непрерывные функции, но значительно сокращающие вычисления. Одним из таких базисов является вейвлетный. В настоящее время он постепенно находит практическое применение [7-15]. Следует при этом отметить, что работы [5-7, 21, 22] носят лишь ознакомительный характер, в основательной работе [21] о применении вейвлетов ничего не сказано, но детально описывается большинство фильтров подавления шумов, которые могут использоваться при моделировании процессов на основе вейвлетного базиса. В трудах [3, 17] приводятся основные теоретические положения по теории вейвлетов, базирующиеся на основе зарубежных исследований. Для создания технологий обработки геоинформации на основе этих работ необходимы дополнительные исследования. В работах [1, 11] выполнены значительные исследования по описанию вейвлетами гравитационного поля Земли. Показано, без детального описания технологии, что применение сжатия информации не всегда эффективно при работе с вейвлетами, хотя при использовании всей информации осуществляется быстрое и полное восстановление поля. В [24, 25] детально

Ж АС.Ярмоленко, О.В.Скобенко 001: 10.31897/РМ1.2018.6.612

:,и  Применение теории вейвлетов...

изложена теория лишь первого этапа вейвлет-преобразования, при этом она не доведена до универсального алгоритма. В [14, 15] выдвигается оригинальная идея комбинирования базиса Фурье и вейвлет-базиса при моделировании гравитационного поля Земли (ГПЗ), приводятся существующие вейвлеты, представлена лишь общая схема вейвлет-анализа-синтеза. В [4] в вейвлет-преобразовании используются зарубежные пакеты при неизвестном алгоритме.

Перечисленные труды имеют практическое и теоретическое значение при обработке геоинформации, но в них теория вейвлетов не доведена до инженерного уровня их применения. На основе отмеченного в данной статье с целью разработки детальной и доступной технологии применения вейвлетов в обработке геоинформации ставятся и решаются следующие задачи: 1) построение ортонормированных базисов вейвлет-преобразований и фильтрация; 2) создание вейвлет-фильтров; 3) построение фильтров последовательностью сверток; 4) применение КМА-анализа для построения ортонормированного базиса вейвлет-преобразования; 5) вейвлет-фильтрация и сжатие данных на основе составленных конкретных программ моделирования полей данных геонформации.

Базовые положения. В соответствии с [18] приведем следующие определения. Оператором циклического сдвига последовательности Z [18, с.125] на к-позиций вправо является оператор Як, осуществляющий создание новой последовательности Як2 по формуле

где п - номер элемента в создаваемой последовательности. Пусть имеется ряд

2 =

при N = 4. Тогда при к = 1 (Я2)(0) = 2(0 - 1) = 2(-1) = - 1) = Z(4 - 1) = 2(3) = 1; (Я2)(1) = = 2(1 - 1) = 2(0) = 1; (Я2)(2) = 2(2 - 1) = 2(1) = 0; (Я2)(3) = 2(3 - 1) = 2(2) = 0. Таким образом получим новую последовательность-вектор

Я2 =

При к = 4 получаем снова исходный вектор

Я2 =

Согласно определению 3.7 и теоремы 3.8 [18, с.176] при некоторых заданных векторах и и V, принадлежащих тому же пространству элементов, что и вектор 2 , например (2), возможно построение ортонормированного базиса вида

В = {VМ1 и {ВДМ-1 = {V, Я?, Я^,..., RN-V, и, Я2и, Я4и,..., RN-2и},

где и, согласно [18], можно назвать отцовским вейвлетом, а V- материнским; символ и - логическое объединение множеств, М = N/2.

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

Ортонормированность (5) возможна тогда лишь [13, с.250-258, лемма 7.1; 18, теорема 3.8], когда система матриц

( Т-7/„Ч т?/„ч

А(п)=т?

и

для п = 0, 1, 2, ..., М- - унитарна.

Под унитарной понимается такая матрица [18, с.100], для которой

А-1 = А*,

где А-1 - обратная к А матрица, а А* - сопряженная к А матрица, получаемая комплексно-сопряженной значением от всех элементов матрицы Аг, транспонированной к А. Комплексно-сопряженным к числу г = х + ¡у есть число г = х + ¡у. Для матриц вещественных чисел А = Аг,

и(п), V (п) _ соответственно преобразования Фурье заданных векторов и и V. Считается, что пара и и Vпорождает ортонормированный базис. В [18] ^мерное пространство элементов вектора 2 обозначено как ¡2(2м) пространство значений, в том числе и комплексных, на котором определена /2-норма [18, с.112]:

А = [£ I2 (к )|:

Множество индексов последовательностей обозначается как = {0,1, ..., N-1} [18, с.111]. Обобщением для двумерного случая является МД^мерное пространство ^(АдуАд^). Если принять, что {В0, В1, ..., ВД1_1} есть ортонормированный базис [18, с.134] ¡2(2Д]) и {С0, С\, ..., Сд-1} есть ортонормированный базис /2(2дт2), то ортонормированный базис пространства I2(2Д1-2Д2) будет Dm1,m2(nl, п2) = Вщ(п1)Ст2(п2), где 0 < Ш1 < N1 _ 1; 0 < щ < N1 _ 1; 0< Ш2 < N2 _ 1, 0 < П2 < N1 _ 1. Например, при

1

1 1 ^ V1 -

1 Г1 1 ^

V1 - Ъ

получим

Dm1,m2(nl, П2) --

Г1 1 1 1 ^

1 -11 -1

1 1 -1 -1

V1 -1 -1 1 У

Скалярным произведением (комплексным) векторов 7 и w называется выражение [18, с.91]

п

1-1

где есть число, комплексно-сопряженное Wj. Сверткой z*w называется вектор с компонентами

N-1

г * w(m) - ^ г(т - п^(п).

п-0

Для использования аппарата свертки в [18, с.172] вводится сопряженное отражение ю вектора ю:

ю(к - п) - ю(п - к).

тогда

N-1

N-1

г *ю(к) - £ 2 (п) *5(к - п)2 (п) * ю(п - к).

п-0

п-0

Построение ортонормированных базисов вейвлет-преобразований и фильтрация. В качестве примера построения ортонормированного базиса примем векторы Хаара [13, с.270; 18, с.190]:

и-| * , * ,0,...,0 | ;

л/2 л/2

V-|-^,0,...,0 | .

л/2 л/2

Поскольку в [18] установлено, что матрицы А(п) (6) для векторов и и V унитарны, то пользуясь правилом (5), построим ортонормированный базис для вейвлет-преобразования вектора 2, например (2). Здесь N = 4 и М = N/2 = 2.

2

А.С.Ярмоленко, О.В.Скобенко

Применение теории вейвлетов...

Тогда в соответствии с (1), (5) на основании (8), (9) можно составить ортонормированный вейвлет-базис:

В -

( 1 0 10 ^ -10 10 0 10 1

ч 0 -10 1 у

Тогда вейвлет-преобразование на данном этапе будет

В нашем примере

= В 2 .

( 1 0 10 у 1 ^

2 - -Г

л/2

-10 10 0 10 1

V0 -10 1А1 у

1

72

(1 ^ -1 1

V 1 у

а обратное преобразование, следуя (11), имеет вид

В примере это

2 - В2 - В*2 - ВТ2.

(1 -10 0 у 1 ^

2 - -Г

л/2

0 0 1 -1 110 0

V0 0 1 1 Л 1 у

-1 1

1

И

(1 ^ 0 0

V1 у

Теперь осуществим фильтрацию сигнала, т.е. его разложение по отцовским и материнским составляющим вейвлет-базиса. Согласно [13, с.26; 18, с.99] в пространстве с комплексным скалярным произведением <,> и ортонормированным базисом Я={ы\, и2, ..., ип} для любого V этого пространства справедливо

V-2

1 -1

V,и]и].

Здесь V соответствует вектору данных 2, а базисные векторы иу - всем столбцам матрицы (10). После подстановки их в (14) получим

где

2 - Q(г) + Р(2),

^2) - 2

(1 У 1 ^

-1

V 0 у

-1

V 0 у

( 0 у 0 У ( ( 1 ^ ( 0 У ( 1 ^

+

0 1

v- 1у

0 1

V- 1уу

-1

V 0 у

+

0 1

V- 1уу

-1 -1

V 1 у

р( 2) - ^

(1У1 ^

(0 У 0 у

+

0 1

V1 уууу

( (1 ^ 1 0

V V °у

+1"

( 0 У 0 1

V1 уу

ч1У

Из примера видно, что Q(Z) и Р(2 составляющие сигнала высоких и низких частот соответственно. Следуя [18, с.180], составляющую высоких частот сигнала представим в виде формулы

*

1

0

0

0

А.С.Ярмоленко, О.В.Скобенко

Применение теории вейвлетов...

( N/2)—1

Q( 2) = £ (2, я2у)я2у = ((2, яу) )я0У

+

к=0

+

((2, яу) )Я2У +... +

( ,

\

2, V

\\ 2 /у

Ям

2

а составляющую низких частот в виде

(N/2)—1

Р(2) = £(2, Я^У = ((2,Щ)яи+

к=0

-((2, Я2и) +... +

Л

2, и

и.

2

Примеры фильтров на основе базисов Шеннона с вещественными и комплексными числами. Если воспользоваться составляющими векторов и и V (преобразований Фурье) в [18, с.179, 181], то при N = 4 с использованием обратных преобразований Фурье можно по аналогии с (14)-(17) определить высоко- и низкочастотные составляющие для данных примеров и в этих базисах (табл.1).

Таблица 1

Разложение вектора ^ по отцовскому и материнскому вейвлетам в различных базисах

Исходный Базисы

вектор 2 Хаара Шеннона 1 Шеннона 2

Р

1 1 3 1 . 1

2 1 4 4 2

0 1 2 11 — - — г 4 4 0

0 1 11 — + — г 44 1

2 2

1 1 2 3 1 — + — г 44 1

6

1 1 11 — + — г 44 1

2 2

0 1 2 11 - — + — г 44 0

0 1 1 1 1

2 — — г 44 2

1 1 2 11 — — — г 44 0

Из сравнения базисов, приведенных в табл.1, а также базисов Мейера, Баттла - Лемарье, Добеши [13, с.207] преимущество следует отдавать базису Хаара по следующим причинам:

• базис Хаара прост в вычислении;

• понятна на его основании и фильтрация сигнала; так, низкочастотная часть на первом этапе равна нулевому коэффициенту преобразования Фурье, высокочастотная - соответствует отклонениям сигнала от его середины;

• в базисах, отличных от Хаара, задаются дополнительные требования к числу N. Например, в базисах Шеннона оно должно быть кратным 4, а в базисе Хаара кратности 2р при р = числом N определяется и вид базиса Добеши [18].

Ж АС.Ярмоленко, О.В.Скобенко DOI: 10.31897/PMI.2018.6.612

:,и  Применение теории вейвлетов...

Построение фильтров последовательностью сверток. В вейвлет-разложениях (16), (17) замечено, что коэффициенты при R2kV, R2kU являются соответствующими скалярными произведениями, которые в виде сверток можно записать так [13, 18, 19]:

Z, R2kv) = ZV~(2k), Z, RiU) = ZU~(2k).

Настоящие выражения являются теоретической основой быстрого вейвлет-преобразования. Тогда фильтр на основе такой свертки строится следующим образом [13, с.279; 18, с.187]: 1) образуются векторы U, V (описаны в конце базовых положений); 2) осуществляются свертки Z*U, Z*V; 3) вводится оператор децимации - удаления составляющих с нечетными номерами D(ZU), D(ZV); 4) реализуется оператор разрежающей выборки - удваивания размера вектора вставкой нуля между двумя смежными значениями U(D(Z*U)), U(D(Z*V)); 5) осуществляется фильтрация -строится вектор низкочастотной составляющей в виде свертки P(Z) = U*U(D(Z*U)) и высокочастотной Q(Z) = V*U(D(Z*V)); 6) восстанавливается сигнал Z = P(Z) + Q(Z).

На этом заканчивается первый этап разложения сигнала на высоко- и низкочастотные составляющие. Число всех этапов определяется по формуле р = log2N. Каждый последующий этап состоит из фазы анализа и фазы синтеза. В фазе анализа на этапе n осуществляется: 1) ввод вектора ZTn_x *Un_j и его децимация Zn = D(ZTn_x *Un_x); 2) ввод векторов Un, Vn нормированного базиса размерности N\ = N/2"-1; 3) осуществляются децимация и разреживание сверток ZTn *Un, ZTn *Vn; UDUn = UD(ZTnUn), UDVn = UD(ZTnVn); 4) вычисляются высокочастотная и низкочастотная составляющие вектора Z„: Q(Z„) = Vn*U(D(Zn* Vn)); P(Z„) = U*U(D(Zn*U„)).

В фазе синтеза осуществляется: 1) разреживание векторов Q(Zn), P(Zn) (при этом разреживание выполняется числом n - 1-раз до достижения размерности исходного сигнала): U(Q(Zn)), U(P(Zn)); 2) операциями свертки получают высоко- и низкочастотную составляющие сигнала на этапе n: Qn(Z) = U(Q(Z„))*Ub P„(Z) = U(P(Zn))*Ul.

Окончательный результат анализа синтеза:

Z = P( Z) + 2 Qp_x( Z).

i=0

Настоящий алгоритм реализован в специально составленной авторами программе Sub Вейв-лет Анализ Синтез ().

Сигнал Z можно представить в виде разложения

Z = c0b0 +C1b1 + ... + CN _1bN

по базису B, который представляет собой совокупность ортонормированных базисов [18, с.176, с.185-186] B = {bi }N=01 .Очевидно, что в этом случае ci = (ZbT). Задача заключается в том, чтобы определить ортонормированные векторы в разложении (18). В теории вейвлетов [18, с.209-225, определение 3.28] вместо (18) принимается следующая запись:

+ 2 C2,k p,k + 2 С3

k=0 k=0 k=0 k=0

Z = 2 c1,k Ф_ p,k + 2 c2,k p,k + 2 c3,k V_

Если принять р = то (19) перепишется в виде

Z = С1,0ф_р,0 + С2,0ф-р,0 + С2,оф-р,1 + • • • + + С#/2,0ф-1,0 + сМ2,1ф-1,1 + — + СМ2,М2-1ф-1,М2-1.

Здесь все векторы базиса В ф_д записаны слева направо по степени детализации вектора 2. В [18] они записаны наоборот справа налево. Построение базисных векторов ф_д осуществляется в следующем порядке [18, с.225]: 1) применяется последовательность вейвлет-фильтров и1, У\, и2, У2; ...; ир, Ур; при этом и1, У1 е 12( например, исходя из (8), (9) при 1 = 2 будет

U2 1 l V2 j

_ 1

V2

Г1 1 0 0^ 1 _ 10 0

; 2) каждый из базисов строится по формуле [18, с.225]

Ж АС.Ярмоленко, О.В.Скобенко 001: 10.31897/РМ1.2018.6.612

:,и  Применение теории вейвлетов...

V— ;,к = Я2 >к/1 ,

где/ = gl-lUl-1(Vl), а gl = g/-1U/-1(Ul) при начальных значениях, равных/1 = VI, gl = Ц\. В приведенных выражениях оператор Ul-1(Vl) означает /-1-кратное разреживание вектора Vl или и1. Элементы сверток/ (п), g/ (п) для п = 0, 1, ..., N-1 вычисляются по формулам

N-1 N-1

fx(«) = Z(gi-i(m)Ul-1(Vl)(n - m)); gl(n) = £(gl-l(m)Ul-1(Ul)(n - m)) .

m=0 m=0

КМА-анализ для построения ортонормированного базиса вейвлет-преобразования. В фундаментальных работах [13, 18, 19] вейвлет-разложение предполагает наличие масштабирующей и уточняющей функций. В [18, определение 5.30] масштабирующая функция еще называется отцовским вейвлетом. В [13, с.103; 19, с.44] термин «отцовский вейвлет» в определении масштабирующей функции не применяется. В тех же работах [18, с.198; 19, с.39-42, гл.5] уточняющая функция (материнский вейвлет) просто названа вейвлетом.

В системе Хаара в соответствии с [18, (5.49)] отцовский вейвлет записан в виде

Г1, 0 < х < 1;

Ф

[0, иначе.

Материнский вейвлет по формуле [18, с.380]

у

и, что легко показать, имеет вид

у( х) =

-1, 0 < х <-;

2

1, 1 < х < 1;

2

0, иначе.

Отцовский и материнский вейвлеты позволяют построить ортонормированный базис вейв-лет-преобразования для представления дискретно-заданной геоинформации с целью ее обработки. Для этого в вейвлет-теории разработан так называемый кратно-масштабный анализ (КМА-анализ). На его основе создается удобный в использовании алгоритм построения ортонормированного базиса вейвлет-преобразования. КМА-анализ базируется на функциях

Фм = 2-2 ф

у;,k = 2-* у

Эти формулы приведены в [5, с.113, с.128; 13, с.241; 18, с.368; 19, с.193]. При этом в приведенных работах показателям степени в (24), (25) приписывается как положительный, так и отрицательный знак, как в нашем случае. В случае отрицательной степени график функции растягивается по оси х, а при положительной - сжимается. В работе нас интересует растяжение по оси х, поэтому принята запись степени с отрицательным знаком. При такой записи уточняются значения функций в зависимости от числа ортонормированных векторов базиса вейвлет-разложения, т.е. осуществляется увеличение деталей анализируемой информации или увеличение разрешения. Поэтому в зарубежной литературе [5, с.113] КМА-анализ справедливо называется многоразрешающим анализом (multiresolution analysis). Построение вейвлет-базисов в системе Хаара будем вести на основе (21)-(25) в следующем порядке:

1) построение базисного вектора нулевого приближения ф;

2) построение последующих уточняющих базисных векторов.

1. Построение базисного вектора нулевого приближения. В основу построения всех базисных векторов положим вейвлет-базис вида [18, (3.67), (3.84)], применяемый нами в разложении (20). Хотя в теории вейвлетов [3, 5, 13, 17-19] допускается возможность построения нескольких

А.С.Ярмоленко, ОЯ.Скобенко

Применение теории вейвлетов...

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

В соответствии с разложением (20), принятым здесь за основу, базисный вектор ф-р, с единствен, так как С «пробегает» значения от 0 до ((N2) - 1) [18, (3.67)]. (В КМА-анализе обозначениям ф_р, с, у-р, с разложения (20) соответствуют ф-р, с, у-р, с.) Поскольку N = 2Р, то С = 0. Следовательно, в (24) > = р, С = 0. В соответствии с (21) 0 < 2-рх - С < 1, или С < х < (1 + С)/2-р. При С = 0 будет 0 < х < 1/2-р и составляющие этого вектора определятся по формуле

Фх)=2—р я,0 <х <2 р;

.0,

Таким образом, на основе (26) получим вектор

иначе.

—р_

фр,0 = 2 2

гр,0 -

с количеством одинаковых членов 2р.

2. Построение последующих уточняющих базисных векторов. Последующие уточняющие базисные векторы вычисляются по (25) с учетом (23). Порядок множества непересекающихся базис. векторов [18, (3.67)] и{у—> С 2 ^ определяется значением>. Это значение изменяется отр - 1

ных

до 0 с шагом - 1. Тогда число базисных уточняющих векторов порядка > составит величину N/21. Значение С будет изменяться от 0 до N/21 — 1. Например, для > = р будем иметь следующий уточняющий вектор (единственный): у-р, 0. Для> = / уточняющими векторами порядка / будут: у^, 0; у^, 1 ; у—/ ^21 )-1. Так, при N = 8 р = 3 уточняющими векторами порядка / = 1 при верхнем пределе С, равном

N/21 - 1 = 2р/2-1= 3, будут: у-1,0; у-1,1; у-1,2; Уу. Каждый из этих векторов при определенном С также определяется по (25) с учетом материнского вейвлета (23). Тогда можно записать

У/с = 2 2

— 1, 0 < 2—/х — С < 1 2

1, 1 < 2—/х — С < 1; 2

0, иначе

или

У/с = 2 2

— 1, С21 < х <|1 + С |21;

1, ^1 + С|21 < х < (1 + С)21;

2

0, иначе.

В качестве примера возьмем р = 3, / = 1, С = 2. Тогда

— 1, 4 < х < 5;

1

= 21

1, 5 < х < 6; 0, иначе.

В табл.2 дан пример вейвлет-разложения как методом сверток, так и КМА-методом.

Таблица 2

Разложение вектора ^ по составляющим ф и у

2

Исходный Составляющие

вектор с1,0-ф-3,0 с1,0-У-3,0 с1,0-У-3,0 с3,ГУ-2,1 04ЖУ-1,0 с4,ГУ-1,1 с4,2-У-1,2 с4,3-У-1,3

4 7,25 -3,25 -1 0 1 0 0 0

2 7,25 -3,25 -1 0 -1 0 0 0

3 7,25 -3,25 1 0 0 -2 0 0

7 7,25 -3,25 1 0 0 2 0 0

10 7,25 3,25 0 -1,5 0 0 1 0

8 7,25 3,25 0 -1,5 0 0 -1 0

10 7,25 3,25 0 1,5 0 0 0 -2

14 7,25 3,25 0 1,5 0 0 0 2

Ж АС.Ярмоленко, О.В.Скобенко DOI: 10.31897/PMI.2018.6.612

:,и  Применение теории вейвлетов...

Каждая составляющая вычислялась по формуле:

для низкочастотной части

N/2 p - j

Pu = 1/2p-j ZP,j-1 i=1

и высокочастотной

Qij = Py-1 - Pij

Здесь j - номер этапа; i - номер составляющей сигнала в группе (группа состоит соответственно из 2, 4,... 2р-элементов в зависимости от порядкового номера этапа): p - число этапов. Окончатель-

p-i

ный результат: Z = P1 + ZJQp-i , где P1 - составляющая нулевой частоты (вектор у-3,0 в примере

i=0

табл.1); Qp-,i - уточняющие составляющие. В приведенных преобразованиях на каждом этапе при нецелом p возможен остаток элементов, число которых меньше числа 2ip при ip = 1, 2, ..., 2p. По этому остатку находится также среднее значение TZ, которое записывается в старшей строке pip на место этого остатка элементов. В последней строке, следующей сразу после строки с номером целой части p, обозначаемой через pf, находится среднее по всем элементам предыдущей строки Ppf. Это среднее является одним и тем же для строки Ppf + 1. Подстрока qpf + 1 вычисляется в общем порядке:

qpf+1 = Ppf- Ppf + 1

Исследование эффективности сжатия геоинформации и фильтрации шумов с помощью вейвлетов. В основу исследований положены специально составленные авторами программы на языках VISUAL BASIC Excel (VBE) и IDL системы ENVI по сжатию и фильтрации геоинформации. В качестве первого объекта исследований принята модель рельефа, приведенная в [23]. В VBE-программе Sub МакросВЕЙВЛсИстОш() истинные высоты точек представлены массивом сс1(), а высоты, отягощенные случайными ошибками, - массивом cm(). Так же, как и в [23], связь настоящих массивов определяется формулой cm(i) = cc1(i) + delta, где i изменяется от 0 до N - 1, а delta = Randbetwen(-t, +t)Std - функция языка VBE-генерирования случайного числа в интервале значений квантиля от -t до +t; Std - задаваемый стандарт случайных ошибок (шумов). Именно высоты массива cm() подвержены вейвлет-разложению по частотам. При вейвлет-сжатии и соответственно фильтрации оставались самые большие по амплитуде члены разложения на всех частотах. Фильтрационный вектор представлен массивом Filt = Array(1; 0,8; 0,6; 0,5; 0,25; 0,15; 0,1; 0). Например, при значении фильтра Filter = Filt(1) = 0,8 на всех частотах разложения оставались значения более 0,8, остальные значения обнулялись. На выходе формировался отфильтрованный массив высот Tw(). Среднее квадратическое отклонение (СКО) отфильтрованных высот определялось по формуле

CKO =

N

Z (Tw(i) - cc(i))2

i =1_

N

Таблица 3

Средние квадратические отклонения (СКО) по каждому фильтру (в метрах) в зависимости от стандартов (Std) распределения случайных ошибок высот

Filter Оставшаяся информация, % Значения Std

0,05 0,1 0,3 0,4 0,5

1 25 0,66 0,88 0,71 0,93 0,94

0,8 33 0,55 0,35 0,63 0,86 0,88

0,6 40 0,23 0,26 0,63 0,86 0,94

0,5 48 0,23 0,26 0,64 0,86 0,94

0,25 51 0,19 0,24 0,63 0,86 0,99

0,15 66 0,15 0,20 0,64 0,86 0,99

0,1 70 0,11 0,17 0,63 0,86 0,99

0 100 0,09 0,18 0,63 0,86 0,99

В табл.3 приведены значения СКО в зависимости от значения фильтра - Filter и точности отмоделированных высот - величины Std. Оставшаяся информация (в процентах) вычислена как процентное отношение числа наибольших оставшихся коэффициентов к их общему числу до сжатия, равному 15.

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

СКО по каждому фильтру при высотах, не отягощенных случайными ошибками:

Filter 0,55 0,34 0,25 0,23 0,15 0,10 0,06 0 СКО 1 0,8 0,6 0,5 0,25 0,15 0,1 0

Из экспериментальных исследований следует:

• для сохранения максимально высокой точности исходных высот при их вейвлет-разложении сжатие информации недопустимо; даже при отсутствии сжатия (Filter = 0) и наличии случайных ошибок среднее квадратическое отклонение результирующих высот больше их стандартного отклонения на входе;

• сжатие информации можно допускать, но при этом следует учитывать порог понижения точности изображения рельефа; например, при трехкратном сжатии (Filter = 0,8) и Std = 0,1 среднее квадратическое отклонение высот на выходе СКО = 0,35 больше Std в 3,5 раза, в остальных случаях оно больше в 2 раза. На рис.1 показан рельеф, построенный по истинным высотам, а на рис.2 - рельеф, полученный при трехкратном вейвлет-сжатии (Filter = 0,8) и стандарте ошибок высот 0,3 м (Std = 0,3 м);

• закономерности влияния вейвлет-сжатия на точность получаемой на выходе информации такие же, как и при фильтрации информации в рядах Фурье [23], но здесь объем вычислений ничтожно мал по сравнению с преобразованиями Фурье.

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

1. По известной методике определяется оценка стандарта высот рельефа (Std) в виде средней квадратической ошибки съемки рельефа.

2. Для данного объекта определяются средние квадратические отклонения по каждому фильтру.

3. Если значение средней ошибки рельефа, полученное как 0,8Std [2], меньше трети высоты сечения рельефа, то по полученной экспериментально можно выбрать порог (величину СКО) понижения точности сжатой информации о рельефе при определенном значении фильтра (Filter). Если средняя ошибка рельефа больше трети высоты сечения рельефа и приближается к величине СКО, получаемой при 2530 %-ном сжатии информации (т.е. трех- и четырехкратном сжатии), то в таком случае значение порога (величины СКО) понижения точности сжатой информации о рельефе не суть важно.

Одновременно авторами составлена программа pro oroi_data_ corr24bitWAVE на алгоритмическом языке IDL системы ENVI вейвлет-

Рис. 1.Рельеф в горизонталях, построенный по высотам, принятым в качестве исходных [23]

Рис.2.Отфильтрованный рельеф при трехкратном вейвлет-сжатии информации

Рис.3. Изображения - исходное (восстановленное по всем вейвлет-элементам) (а) и сжатое в 2,5 раза (б) с отфильтрованными значениями (менее 4 единиц) составляющих вейвлет-разложения по модулю

фильтрации (сжатия) изображений. На рис.3 приведены изображения - исходное (и восстановленное по всем вейвлет-элементам) и сжатое в 2,5 раза с отфильтрованными значениями (менее 4 единиц) составляющих вейвлет-разложения по модулю.

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

Выводы

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

2. При фильтрации и сжатии изображений рекомендуется подобрать вначале фильтр в соответствии с пунктами 1-3 в случае рельефа. При этом необходимо задавать СКО значения пиксела, при котором допустимо сжатие и отсутствуют видимые несовпадения.

ЛИТЕРАТУРА

1. БагровА.А. Разложение сферических функций по вейвлетам Хаара / А.А.Багров, А.С.Багрова // Известия вузов. Геодезия и аэрофотосъемка. 2008. № 4. С. 6-8.

2. Большаков В.Д. Теория математической обработки геодезических измерений / В.Д.Большаков, П.А.Гайдаев. М.: Недра, 1977. 367 с.

3. Воробьев В.И. Теория и практика вейвлет-преобразования / В.И.Воробьев, В.Г.Грибунин. СПб: ВУС. 1999. 204 с.

4. Гонжа Е.А. О вейвлет-фильтрации цифрового изображения земной поверхности // Известия вузов. Геодезия и аэрофотосъемка. 2017. № 3. С. 105-110.

5. ДьяконовВ.П. Вейвлеты. От теории к практике. М.: Солон-Р, 2002. 448 с.

6. Журкин И.Г. Автоматизированная обработка данных дистанционного зондирования / И.Г.Журкин, Н.К.Шавенько. М.: Диона, 2013. 456 с.

7. Красильников Н.П. Цифровая обработка 2D- и 3D-изображений. СПб: БХВ-Петербург, 2011. 608 с.

8. Лапшин А.Ю. Разработка и исследование методики вычисления гравиметрической высоты квазигеоида и составляющих уклонения отвеса на основе вейвлет-преобразования: Автореф. дис. ... канд. техн. наук / Московский университет геодезии и картографии. М., 2011. 207 с.

9. Мазурова Е.М. Алгоритмы быстрого преобразования Фурье // Известия вузов. Геодезия и аэрофотосъемка. 2004. № 3. С. 18-35.

10. Мазурова Е.М. Двумерное и матричное представление быстрого преобразования Фурье // Известия вузов. Геодезия и аэрофотосъемка. 2004. № 4. С. 3-12.

11. Мазурова Е.М. К вопросу о вычислении аномалии высоты на основе вейвлет-преобразования и быстрого преобразования Фурье в плоской аппроксимации / Е.М.Мазурова, А.С.Багрова // Известия вузов. Геодезия и аэрофотосъемка.

2008. № 4. С. 6-8.

12. Малинников В.А. Анализ методов формирования мультифрактальной меры, основанных на вейвлет-обработке экспериментальных данных / В.А.Малинников, Д.В.Учаев // Известия вузов. Геодезия и аэрофотосъемка. 2007. № 6. С. 57-61.

13.Мала С. Вейвлеты в обработке сигналов: Пер с англ. М.: Мир, 2005. 671 с.

14. Нейман Ю.М. Об адаптации глобальной модели геопотенциала к региональным особенностям. Ч.1 / Ю.М.Нейман, Л.С.Сугаипова // Известия вузов. Геодезия и аэрофотосъемка. 2014. № 3. С. 3-12.

15. Нейман Ю.М. Основы разномасштабной аппроксимации геопотенциала / Ю.М.Нейман, Л.С.Сугаипова. М.: Изд-во МИИГАиК, 2016. 218 с.

16. ПискуновН.С. Дифференциальные и интегральные исчисления: В 2-х т. М.: Наука. 1978. Т.2. 575 с.

17. Уэлстид С. Фракталы и вейвлеты для сжатия изображений в действии. М.: Изд-во Триумф, 2003. 320 с.

18. Фрейзер М. Введение в вейвлеты в свете линейной алгебры: Пер. с англ. М.: Бином. Лаборатория знаний. 2008. 487 с.

19. Чуи К Введение в вейвлеты: Пер. с англ. М.: Мир, 2001. 412 с.

20. Шпарр Д. Точное сельское хозяйство (Precision agriculture) / Д.Шпаар, А.В.Захаренко, В.П.Якушев /Федеральное министерство продовольствия, сельского хозяйства и защиты прав потребителя Германии. Санкт- Петербург - Пушкин.

2009. 398 с.

21. Шовенгердт Р.А. Дистанционное зондирование, модели и методы обработки изображений. М.: Техносфера, 2010.

560 с.

22. ЯковлевА.Н. Введение в вейвлет-преобразование. Новосибирск: Изд-во НГТУ, 2003. 104 с.

23. ЯрмоленкоА.С. Фильтрация геоинформации в рядах Фурье / А.С.Ярмоленко, О.В.Скобенко // Известия вузов. Геодезия и аэрофотосъемка. 2016. № 1. С. 107-113.

24. ЯрмоленкоА.С. Использование вейвлетов в аналитическом представлении дискретных функций графической информации // Известия вузов. Геодезия и аэрофотосъемка. 2008. № 3. С. 20-30.

25. Ярмоленко А.С. Вейвлет-преобразования в кодировании графической информации // Известия вузов. Геодезия и аэрофотосъемка. 2010. № 4. С. 18-25

Авторы: АС.Ярмоленко, д-р техн. наук, профессор, yarmolenko_alex@mail.ru (Новгородский государственный университет им. Ярослава Мудрого, Великий Новгород, Россия), О.В.Скобенко, ассистент, olga-skobenko@mail.ru (Белорусская государственная сельскохозяйственная академия, г. Горки, Могилевская область, Республика Беларусь). Статья поступила в редакцию 19.04.2018. Статья принята к публикации 21.09.2018.