Метод сжатия гиперспектральных изображений на основе блочного кодирования с преобразованием тема диссертации и автореферата по ВАК РФ 05.13.17, кандидат наук Юзькив Руслан Романович

• Юзькив Руслан Романович
• кандидат наук
• 2019

Оглавление диссертации кандидат наук Юзькив Руслан Романович

ВВЕДЕНИЕ

РАЗДЕЛ 1 Задача сжатия гиперспектральных изображений

1.1 Основные характеристики методов сжатия

1.2 Общая схема трёхмерного метода кодирования с преобразованием

1.3 Задачи, требующие решения при использовании трёхмерного метода кодирования с преобразованием

1.3.1 Выбор преобразования А

1.3.2 Выбор коэффициентов квантования Q

1.3.3 Выбор способа развёртки АС-коэффициентов

1.4 Модель автокорреляционной функции гиперспектральных изображений

1.4.1 Построение модели автокорреляционной функции

1.4.2 Оценка параметров модели

1.4.3 Экспериментальное исследование модели

1.5 Выводы и результаты

РАЗДЕЛ 2 Метод кодирования с преобразованием для сжатия гиперспектральных изображений

2.1 Вычисление коэффициентов квантования Q

2.1.1 Вычисление начальных значений коэффициентов квантования

2.1.2 Итеративное уменьшение значений коэффициентов квантования

2.2 Выбор способа развёртки АС-коэффициентов в одномерный массив

2.3 Стационаризация блоков по пространственным координатам

2.3.1 Стационаризация на основе бинаризации

2.3.2 Стационаризация на основе кластеризации

2.4 Выводы и результаты

РАЗДЕЛ 3 Программный комплекс

3.1 Архитектура разработанного программного комплекса

3.1.1 Вспомогательные модули

3.1.2 Модуль преобразования

3.1.3 Модуль компенсации средних значений

3.1.4 Модуль вычисления коэффициентов квантования

3.1.5 Модуль вычисления схемы развёртки AC-коэффициентов

3.1.6 Модуль статистического кодирования

3.2 Консольный интерфейс пользователя

3.2.1 Кодирование данных

3.2.2 Декодирование данных

3.2.3 Сравнение кодированных и декодированных данных

3.3 Выводы и результаты

РАЗДЕЛ 4 Экспериментальные исследования

4.1 Сравнение алгоритмов в рамках метода без стационаризации

4.1.1 Сравнение начальных приближений Q и способов развёртки AC-

коэффициентов

4.1.2 Сравнение алгоритмов вычисления ошибок квантования

4.1.3 Влияние размера блока на эффективность сжатия

4.2 Сравнение с известными решениями

4.2.1 Сравнение с послойным применением JPEG

4.2.2 Сравнение с алгоритмами на основе работ Haiyan, Qiao

4.2.3 Сравнение с методами на основе вейвлет-преобразований

4.3 Исследование эффективности стационаризации

4.3.1 Стационаризация на основе бинаризации

4.3.2 Стационаризация на основе кластеризации

4.4 Выводы и результаты

ЗАКЛЮЧЕНИЕ

ПЕРЕЧЕНЬ СОКРАЩЕНИЙ

ПЕРЕЧЕНЬ УСЛОВНЫХ ОБОЗНАЧЕНИЙ

СПИСОК ЛИТЕРАТУРЫ

ПРИЛОЖЕНИЕ А. Тестовый набор гиперспектральных изображений

ПРИЛОЖЕНИЕ Б. Числовые параметры построенных моделей автокорреляционных функций

ПРИЛОЖЕНИЕ В. Исходный код ключевых алгоритмов

ПРИЛОЖЕНИЕ Г. Результаты экспериментальных исследований

ПРИЛОЖЕНИЕ Д. Акт об использовании результатов диссертации в АО «Самара-Информспутник»

ВВЕДЕНИЕ

Введение диссертации (часть автореферата) на тему «Метод сжатия гиперспектральных изображений на основе блочного кодирования с преобразованием»

Актуальность темы исследования

Гиперспектральное изображение (ГСИ) представляет собой трёхмерный массив данных, в который входит двумерная пространственная информация и одномерная спектральная. В отличие от обычных изображений, каждая пространственная точка описывается не одним-четырьмя спектральными каналами, а целыми сотнями (например, на гиперспектрометре АУЖИБ их количество равно 224 [1]). Это предоставляет обширные возможности по обнаружению и анализу свойств наблюдаемых объектов и применению ГСИ в самых различных областях, таких как картография, обнаружение и мониторинг чрезвычайных ситуаций, анализ экологической обстановки, поиск полезных ископаемых и так далее [2].

Однако за наличие подробной спектральной информации приходится платить объёмом данных. При этом использование методов сжатия данных без потерь, как правило, позволяет уменьшить объём не более чем в два раза [3]. Это приводит к целесообразности рассмотрения более эффективных методов сжатия с потерями, позволяющих получить высокие коэффициенты сжатия при некоторых допустимых искажениях. Существует множество известных подходов [4, 5], которые можно условно разделить на две группы.

Первая группа использует в своей основе дифференциальное кодирование [6, 7, 8]. Основная идея заключается в том, что соседние отсчёты данных используются для предсказания текущего отсчёта. Ошибка между предсказанными и изначальными значениями квантуется. Далее эта ошибка кодируется вместе с добавочной информацией, необходимой для осуществления процесса декодирования.

Вторая группа основана на использовании кодирования с преобразованием [9, 10, 11]. При этом изначальные данные обычно представляются через коэффициенты в новом базисном пространстве и квантуются. Как правило, этот базис подбирается таким образом, чтобы основной объём информации в новом представлении оказывался сконцентрированным в как можно меньшем числе отсчётов. При этом

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

Известным недостатком методов сжатия на основе кодирования с преобразованием является вычислительная сложность. Однако экспоненциальный рост эффективности вычислительных средств (с 1972 по 2007 год, [12]) в какой-то степени компенсирует алгоритмическую сложность данного подхода и позволяет рассматривать методы этой группы как эффективную (с точки зрения получаемых коэффициентов сжатия) альтернативу методам, основанным на дифференциальном кодировании.

Степень разработанности темы исследования

Для сжатия обычных (не гиперспектральных) изображений известными стандартами являются основанный на дискретном косинусном преобразовании (ДКП) JPEG [13] и основанные на вейвлет-преобразованиях JPEG 2000 [14] и ICER [15]. Последний был специально разработан для марсианских миссий Национального управления США по аэронавтике и исследованию космического пространства.

При переходе от двумерных изображений к ГСИ возникает необходимость учёта корреляции вдоль третьего спектрального измерения. Такой учёт может быть реализован различными способами. Например, для спектрального измерения можно использовать дифференциальное кодирование, а двумерные пространственные слои сжимать классическим JPEG [16]. Другим подходом является использование трёхмерных преобразований.

Примером использования трёхмерных вейвлет-преобразований является стандарт ICER 3D [17], созданный на основе своего двумерного предшественника ICER. Помимо этого, теме кодирования ГСИ с помощью трёхмерных вейвлет-преобразований посвящено множество работ (например, [18, 19, 20, 21]).

Что касается методов на основе ДКП, то в работах [22, 23, 24] было показано, что использование трёхмерных блоков в задачах сжатия видео позволяет существенно увеличить эффективность сжатия по сравнению с двумерными блоками. Аналогичный эффект можно ожидать и при сжатии ГСИ, являющихся по сути трёхмерными.

Несмотря на явный уклон интереса исследователей в последнее время к вейвлет-преобразованиям, потенциал использования методов, основанных на ДКП, раскрыт ещё не полностью. Во-первых, независимая обработка блоков позволяет распараллелить вычисление преобразований и дальнейшее кодирование блоков. Во-вторых, несмотря на то, что для двумерных изображений технологии сжатия на основе ДКП уже давно известны и хорошо проработаны, при сжатии трёхмерных данных возникает ряд новых задач. В частности, подбор адекватной модели ГСИ и улучшение различных алгоритмов в рамках метода до сих представляют собой поле для возможных исследований, что показывают работы [25, 26, 27, 28]. В-третьих, практически все предложенные на данные момент алгоритмы вычисления коэффициентов квантования зависят от абстрактных параметров качества, не позволяющих непосредственно контролировать получаемую ошибку восстановления.

Цели и задачи

Целью данной работы является повышение эффективности сжатия ГСИ методом кодирования на основе преобразования. Для достижения этой цели в работе решаются следующие задачи:

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

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

3. Экспериментальное исследование эффективности разработанного метода сжатия ГСИ, выработка рекомендаций по рациональному выбору его параметров.

Методология и методы исследования

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

Научная новизна

1. Предложен новый метод сжатия ГСИ на основе блочного кодирования с преобразованием, отличающийся от известных методов введением этапа предварительной обработки данных, способами квантования и развёртки спектральных компонент трёхмерных блоков ГСИ.

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

3. Предложена и экспериментально обоснована статистическая модель ГСИ в виде трёхмерного стационарного (пространственно-однородного) случайного поля с автокорреляционной функцией, двухкомпонентной по пространственным координатам и однокомпонентной - по оптическому спектру.

4. Разработан алгоритм эвристической оценки параметров автокорреляционной функции случайного поля в рамках предложенной модели ГСИ.

5. Предложен алгоритм предварительной обработки (стационаризации) ГСИ на основе послойной бинаризации блоков методом Оцу.

6. Предложен алгоритм предварительной обработки (стационаризации) ГСИ на основе кластеризации всего ГСИ методом к-средних.

7. Предложен новый способ развёртки компонент трёхмерного спектрального преобразования блоков ГСИ в одномерный массив, необходимый для дальнейшего статистического кодирования данных.

8. Проведено экспериментальное исследование предложенного метода сжатия ГСИ применительно к обработке гиперспектральных данных аэрокосмического дистанционного зондирования Земли, подтвердившее работоспособность нового метода и демонстрирующее его преимущества перед известными.

Теоретическая и практическая значимость работы

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

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

позволяет использовать его для дальнейшей разработки и экспериментальных исследований различных алгоритмов сжатия ГСИ в рамках методологии блочного кодирования с преобразованием.

Положения, выносимые на защиту

1. Новый метод сжатия ГСИ на основе блочного кодирования с преобразованием.

2. Итеративный алгоритм вычисления коэффициентов квантования компонент трёхмерного спектрального преобразования блоков ГСИ.

3. Статистическая модель ГСИ и алгоритм оценки параметров её автокорреляционной функции.

4. Алгоритмы предварительной обработки данных на основе бинаризации и кластеризации для учёта особенностей статистической модели ГСИ при сжатии.

5. Способ развертки компонент трёхмерного спектрального преобразования блоков ГСИ в одномерный массив.

6. Экспериментальное подтверждение эффективности предложенного метода сжатия ГСИ применительно к обработке гиперспектральных данных ДЗЗ, рекомендации по рациональному выбору параметров метода.

Степенъ достоверности и апробация результатов

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

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

- на 10-ой международной научной конференции «Интеллектуализация обработки информации» (ИОИ 2014), Греческая Республика, Крит, 4-11.10.2014;

- на 17-ой международной конференции «Цифровая обработка сигналов и её применение» ДОРА 2015), Россия, Москва, 24-27.03.2015;

- на 4-ой международной конференции «Анализ изображений, сетей и текстов» (АИСТ 2015), Россия, Екатеринбург, 9-11.04.2015;

- на 3-ей международной конференции «Информационные технологии и нанотехнологии» (ИТНТ 2017), Россия, Самара, 25-27.04.2017;

- на международной конференции «Перспективные информационные технологии» (ПИТ 2019), Россия, Самара, 24-26.06.2019.

Реализация результатов работы

Результаты диссертации были использованы при выполнении следующих работ:

- в проекте Российского научного фонда № 14-31-00014 «Создание лаборатории прорывных технологий дистанционного зондирования Земли»;

- в гранте Российского фонда фундаментальных исследований № 16-29-09494-офи-м «Методы компьютерной обработки мультиспектральных данных дистанционного зондирования Земли для определения ареалов растений в специальных криминалистических экспертизах»;

- в государственном задании Института систем обработки изображений РАН

- филиала ФНИЦ «Кристаллография и фотоника» РАН, проект № 0026-2014-0106 «Методы кодирования данных дистанционного зондирования Земли»;

- в АО «Самара-Информспутник» (договор № 48/17 от 15.07.2017, заказчик

- ФГУП «ГосНИИПП»; акт об использовании результатов приведён в приложении Д).

Публикации

Основные результаты исследования опубликованы в девяти работах, из них две научные статьи [29, 30] опубликованы в рецензируемом научном издании, рекомендуемом ВАК; две работы [31, 32] опубликованы в англоязычных сборниках трудов конференций, индексируемых в Scopus и Web of Science; четыре работы [33,

34, 35, 36] опубликованы в сборниках трудов международных конференций. Получено одно свидетельство об официальной регистрации программы для ЭВМ [37].

РАЗДЕЛ 1

Задача сжатия гиперспектральных изображений

Введём формальное определение ГСИ следующим образом:

f - /(4-1,42,43),

0<n1<N1, 0 <п2 < N2, 0 <п3 < N3, где f - трёхмерный массив значений яркости;

п1, п2 и п3 - целочисленные индексы соответственно по вертикальной, горизонтальной и спектральной координате;

N1, N2, N2 - соответственно высота, ширина и глубина (количество спектральных каналов) изображения.

Для краткости и удобства записи также будем использовать векторное обозначение координат:

f(n) = ffai, П2, п3), п - (щ, щ, п3), neN = {0,...,N1-l}x{0,...,N2-l}x{0,...,N3-l},

где п - трёхмерный вектор целочисленных индексов;

N - множество допустимых значений п с INI — N1 • N2 • N3 элементами.

1.1 Основные характеристики методов сжатия

Сжатие (компрессия) данных - алгоритмическое преобразование данных, производимое с целью уменьшения занимаемого ими объёма. Самым важным показателем эффективности является коэффициент сжатия [38]:

К -'о

где 10 - объём данных до сжатия; I - объём данных после сжатия. Все методы сжатия данных можно условно разделить на два класса: сжатие без потерь и сжатие с потерями.

Сжатие без потерь подразумевает, что закодированные данные могут быть однозначно восстановлены с точностью до бита. Яркими примерами являются кодирование Хаффмана [39] и арифметическое кодирование [40].

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

Существует множество известных критериев - например, максимальная ошибка, вероятностно-зональный критерий [41], среднеквадратическая ошибка, отношение сигнала к шуму, отношение пикового уровня сигнала к шуму (см. с. 185 в [42], с. 613 в [43]). Каждый из перечисленных критериев обладает своими преимуществами и недостатками.

Самым строгим критерием является абсолютное (максимальное) отклонение:

В данной работе в качестве основного показателя потери качества будет использоваться среднеквадратическая ошибка:

При необходимости использования критерия (1) существует известный подход [44], связанный с «достраиванием» любых методов сжатия за счёт дополнительной обработки поля ошибок и, таким образом, обеспечения необходимой максимальной ошибки. Поэтому в рамках данной работы ограничимся рассмотрением только среднеквадратического критерия (2), полагая, что в случае надобности критерий (1) может быть достигнут путём дополнительной обработки поля ошибок.

£абс = тах1Т(п)-Т(п где f(n) - отсчёты изначального изображения; f(n) - отсчёты восстановленного изображения.

|/(п)-/(п)|.

(1)

1.2 Общая схема трёхмерного метода кодирования с преобразованием

На рисунке 1 приведена общая схема сжатия трёхмерного массива f методом кодирования с преобразованием.

С1 Ск-1

Разбиение на блоки

Этап 1

Квантование

С1 =

Этап 3

йс

Преобразование В1 = АЪ,

Этап 2

Развёртка АС ас^п)) = С1(п)

Этап 4

| | 11 I-• -п ас0 гггп-п ас{ 11 | 11-• -п аск_1

Разностное кодирование

Статистическое кодирование

Этап 5

Сжатые данные

Рисунок 1 - Схема сжатия трёхмерных данных методом кодирования с преобразованием

Первый этап

Пусть дано ГСИ f размера Ы1 X Ы2 X Ы3. На первом этапе происходит разбиение изображения на непересекающиеся трёхмерные блоки размера М1х М2х М3. Обозначим эти блоки как

Ь1, 0<КК, где К = К1 • К2 • К3 - итоговое количество блоков;

К± =

N1

мл

К? =

Я?

Мп

К* =

N1

мл

- количество блоков вдоль вертикальной, гори-

зонтальной и спектральной оси соответственно;

¡•] - оператор округления вверх до ближайшего целого числа.

Аналогично множеству N допустимых значений вектора координат для ГСИ, введём множество допустимых значений вектора координат для отдельно рассматриваемого блока (т. е. множество допустимых значений для относительных координат в пределах одного блока):

М = [0,...,М1-1}х[0,...,М2-1}х[0,...,М3-11 \М\ = М1^М2^М3.

Второй этап

На втором этапе каждый из блоков Ь^ подвергается преобразованию:

Вь(т) = ^ Ьь(п) • А(т,п), теМ, 0<КК,

пем

где А - ядро некоторого преобразования. По аналогии с одномерным случаем [38] будем называть получаемые значения В^ значениями в обобщённых координатах, или трансформантами. В рамках рассматриваемой задачи сжатия данных на преобразование А естественным образом накладываются следующие ограничения. Во-первых, оно должно быть обратимым, то есть должно существовать ядро обратного преобразования А-1. Во-вторых, в результате преобразования основной объём информации должен концентрироваться по возможности в меньшем числе обобщённых координат для обеспечения эффекта сжатия.

Третий этап

На третьем этапе происходит квантование трансформант Bt:

\Bt(n)]

С tin) =

Q(n)

пем, 0<i< к,

где Q - массив размера М1Х М2 х М3 заданных коэффициентов квантования (этот массив один и тот же для всех блоков ГСИ);

[•] - оператор округления до ближайшего целого числа.

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

Четвёртый этап

На четвёртом этапе полученные целочисленные значения трёхмерных блоков С t развёртываются в одномерные массивы. Для большинства видов дискретных преобразований (например, Фурье и косинусного) первый коэффициент в обобщённых координатах Сt ( 0), где 0 = (0,0,0), соответствует среднему значению яркости и, как правило, его значение сильно отличается от значений других коэффициентов Сt (п), п ф 0, соответствующих более высоким частотным составляющим кодируемого сигнала. Вследствие этого упомянутые значения кодируются раздельно. Для унификации дальнейшего процесса изложения положим, что используемое преобразование А также имеет указанную особенность, а именно: самый первый коэффициент в обобщённых координатах С ( 0) соответствует среднему значению яркости изначального сигнала.

Значения Сt ( 0) получили название DC-коэффициентов (от англ. «Direct Current»). На четвёртом этапе DC-коэффициенты из всех блоков извлекаются в отдельный массив d с:

dc(i) = С(0), 0< К к.

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

Остальные значения Ct(n), п ф 0, получили названия AC-коэффициентов (от англ. «Alternating Current»). На четвёртом этапе AC-коэффициенты каждого блока С t извлекаются в соответствующие одномерные массивы act, каждый из которых имеет размер М± • М2 • М3 — 1:

act(s(n)) = С (п), п е М\{0}, 0<i<K. где s: М\{0} ^ {0,..., 1М1 — 2} - некоторое правило развёртки трёхмерного блока в одномерный массив. На следующем этапе AC-коэффициенты будут кодироваться длинами серий нулей, поэтому правило развёртки должно быть выбрано таким образом, чтобы как можно больше ненулевых значений оказалось сконцентрированным в начале массива и как можно больше нулевых значений - в конце массива.

Пятый этап

На пятом этапе полученные массивы d и a t подвергаются статистическому кодированию. Отметим, что на данном этапе размерность кодируемых данных уже не имеет значения, так как изначальная задача сведена к кодированию одномерных массивов. Поэтому на данном этапе могут использоваться те же самые алгоритмы, что и для двумерных данных в стандарте JPEG [13]. Как следствие, рассмотрение пятого этапа как такового выходит за рамки данной работы. Однако для полноты картины кратко опишем особенности статистического кодирования полученных одномерных массивов в упомянутом стандарте. Эвристические соображения по особенностям используемых алгоритмов с подробными примерами могут быть найдены в [45].

ОС-коэффициенты, представляющие собой средние значения блоков, кодируются дифференциальным методом. То есть вместо самих значений йс рассматриваются их разности dif^.

= йс(0), = йс(1) - йф -1), 1<КК.

Каждое значение dif(i) дальше представляется парой «длина бинарного кода; бинарный код» в соответствии с кодовой таблицей 1.

Таблица 1 - Коды, используемые для кодирования коэффициентов DC и AC

Значения Длина бинарного кода Бинарные коды

0 0 Нет

-1 1 0

1 1

-3, -2, 2 00, 01

2, 3 10, 11

-7, -6, -5, -4 3 000, 001, 010, 011

4, 5, 6, 7 100, 101, 110, 111

Длина бинарного кода записывается с помощью кода Хаффмана [39]. Таким образом, каждый DC-коэффициент dc(i) кодируется парой «код Хаффмана для длины бинарного кода dif(i); бинарный код dif(i)».

Массивы AC-коэффициентов кодируются длинами серий нулей. Рассмотрим процедуру кодирования на примере одного массива act. Каждый ненулевой AC-коэффициент aCi(j) представляется тройкой «количество нулевых AC-коэффициентов перед текущим AC-коэффициентом; длина бинарного кода aci(j); бинарный код aCi(j)». Первые два значения объединяются в комбинированную пару. В итоге ненулевой AC-коэффициент aci(j) кодируется парой: «код Хаффмана комбинированной пары, состоящей из количества предшествующих нулевых AC-коэффициентов и длины бинарного кода ас^^); бинарный код aCi(j)». После кодирования последнего ненулевого AC-коэффициента в массиве act используется специальная комбинированная пара «0; 0», означающая, что все оставшиеся значения в массиве равны нулю.

В JPEG комбинированная пара занимает один байт: четыре бита на количество нулевых значений и четыре бита на длину бинарного кода AC. Количество предшествующих нулей может оказаться большим, чем можно представить с помощью четырёх бит (т. е. больше 15). Поэтому максимальное значение 15 для количества нулей трактуется как переполнение. В случае переполнения для представления количества нулей используется несколько комбинированных пар. Например, для 18 нулей будет использовано две пары: «15; 0», «3; длина бинарного кода».

При кодировании как AC-коэффициентов, так и DC-коэффициентов, вместо кода Хаффмана может использоваться арифметическое кодирование [40].

На рисунке 2 приведена общая схема восстановления трёхмерного массива данных методом кодирования с преобразованием, которая в целом представляет собой схему, обратную схеме кодирования, приведённой на рисунке 1.

На первом этапе декодируются массивы йс и аС1.

На втором этапе по массивам йс и аС1 восстанавливаются блоки С^.

На третьем этапе по квантованным значениям С восстанавливаются значения обобщённых координат В^:

На четвёртом этапе для каждого блока В^ выполняется обратное преобразо-

Схема декодирования

Èi(ri) = Ci(n) • Q(n), ne M, 0<i<K.

вание для получения восстановленных блоков изображения bi :

b)i(n) = Êi(m) • A-1(m,n), ne M, 0<i<K.

тем

На пятом этапе на основе набора блоков Ь осуществляется сборка восстановленного куба значений f.

Сжатые данные

Разностное декодирование

□ •••□•••□ йс

Реконструкция БС С1(0) = йс(1)

Статистическое декодирование

Этап 1

гттп - р ас0 гггп-п ась гггп - п аск_г

Реконструкция АС С^п) = ас^п))

Этап 2

Восстановление В1(п) = С1(п) • Q(n)

Этап 3

Преобразование Ь1 = А-1В,

Этап 4

Сборка куба

Этап 5

Рисунок 2 - Схема восстановления трёхмерных данных методом блочного кодирования с преобразованием

1.3 Задачи, требующие решения при использовании трёхмерного метода кодирования с преобразованием

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

1) выбор трёхмерного преобразования А, которое определяет, насколько эффективно происходит концентрация информации при переходе в спектральную область;

2) выбор коэффициентов квантования Q которые определяют уровень допустимых потерь;

3) выбор способа развёртки 5 трёхмерного блока АС-коэффициентов С I в одномерный массив ас^, который определяет эффективность статистического кодирования на последнем этапе.

Рассмотрим каждую из этих задач подробнее.

1.3.1 Выбор преобразования А

В подавляющем числе работ по сжатию двумерных данных методом кодирования используется ДКП. Популярность косинусного преобразования обуславливается его близостью к идеально возможному преобразованию Хотеллинга (разложению по собственным векторам АКФ) для сигналов с экспоненциальной АКФ [46] и наличием быстрых алгоритмов его вычисления.

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

Таким образом, первой возникающей задачей при использовании трёхмерного метода кодирования с преобразованием является построение АКФ ГСИ.

1.3.2 Выбор коэффициентов квантования Q

В двумерном алгоритме JPEG возможные значения коэффициентов квантования определены самим стандартом и определяются допустимым уровнем потерь качества (значения коэффициентов были подобраны экспериментальным путём на основе свойств восприятия зрительной системы человека [13]). При этом алгоритм позволяет использовать произвольные значения коэффициентов квантования и существуют работы, позволяющие вычислять эти коэффициенты по некоторому заданному критерию. Например, в [47] используется вычисление коэффициентов квантования, адаптированное к кодируемым данным, с целью минимизации ошибки восприятия, которая учитывает особенности контраста и освещённости изображения.

Список литературы диссертационного исследования кандидат наук Юзькив Руслан Романович, 2019 год

СПИСОК ЛИТЕРАТУРЫ

1. NASA. AVIRIS Data Portal. URL: https://aviris.jpl.nasa.gov/alt_locator/ (дата обращения: 01.09.2016).

2. Перспективные информационные технологии дистанционного зондирования Земли / М. В. Гашников [и др.]; под ред. В. А. Сойфера. - Самара: Новая техника, 2015. - 256 с.

3. Ефремов, В. Ю. Технологии построения автоматизированных систем хранения спутниковых данных / В. Ю. Ефремов, Е. А. Лупян, А. А. Мазуров, А. А. Прошин, Е. В. Флитман // Современные проблемы дистанционного зондирования Земли из космоса. - 2004. - Т. 1. - С. 437-443.

4. Yu, G. Image compression systems on board satellites / G. Yu, T. Vladimirova, M. N. Sweeting // Acta Astronautica. - 2009. - Vol. 64. - P. 988-1005.

5. Babu, S. Hyperspectral image compression algorithms / S. Babu, K. K. Thyaghara-jan // Advances in Intelligent Systems and Computing. - 2014. - Vol. 325. -P. 127-138.

6. Mielikainen, J. S. Lossless hyperspectral image compression via linear prediction / J. S. Mielikainen, A. Kaarna, P. J. Toivanen // Algorithms and Technologies for Multispectral, Hyperspectral, and Ultraspectral Imagery VIII. - 2002. - Vol. 4725. - P. 600-608.

7. Wang, H. Lossless hyperspectral image compression using context-based conditional averages / H. Wang, S. D. Babacan, K. Sayood // IEEE Geoscience and Remote Sensing Letters. - 2007. - Vol. 45(12). - P. 4187-4193.

8. Гашников, М. В. Иерархичная сеточная интерполяция при сжатии гиперспектральных изображений / М. В. Гашников, Н. И. Глумов // Компьютерная оптика. - 2014. - Т. 38, № 1. - С. 87-93.

9. Прэтт, У. К. Цифровая обработка изображений. Книга 2 / У. К. Прэтт. - М.: Мир, 1982. - 480 с.

10. Методы компьютерной обработки изображений. 2-е испр. изд. / М. В. Гашников [и др.]; под ред. В. А. Сойфера. - М.: Физматлит, 2003. - 784 с.

11. Гашников, М. В. Методы сжатия цифровых сигналов и изображений / М. В. Гашников, Н. И. Глумов, С. Б. Попов, В. В. Сергеев. - Самара: Самарский государственный аэрокосмический университет, 2006. - 88 с.

12. Хокинг, С. О Вселенной в двух словах / С. Хокинг. - М.: АСТ, 2017. - 224 с.

13. ISO/IEC 10918-1. Information technology. Digital compression and coding of continuous-tone still images. Requirements and guidelines / ISO/IEC JTC 1/SC 29 Coding of audio, picture, multimedia and hypermedia information. - 1994. - 182 p.

14. ISO/IEC 15444-1. Information technology. JPEG 2000 image coding system. Core coding system / ISO/IEC JTC 1/SC 29 Coding of audio, picture, multi-media and hypermedia information. - 2016. - 221 p.

15. Kiely, A. The ICER Progressive Wavelet Image Compressor / A. Kiely, M. Klimesh // The Interplanetary Network Progress Report. - 2003. - Vol. 42155 (J). - 46 p.

16. Abousleman, G. P. Compression of hyperspectral imagery using hybrid DPCM/DCT and entropy-constrained trellis coded quantization / G. P. Abousleman // Proceedings of Data Compression Conference. - 1995. -P. 322-331.

17. Kiely, A. ICER-3D: A Progressive Wavelet-Based Compressor / A. Kiely, M. Klimesh, H. Xie, N. Aranki // The Interplanetary Network Progress Report. -2006. - Vol. 42-164 (A). - 21 p.

18. Penna, B. Progressive 3-D Coding of hyperspectral images based on JPEG 2000 / B. Penna, T. Tillo, E. Magli, G. Olmo // IEEE Geoscience and Remote Sensing Letters. - 2006. - Vol. 3(1). - P. 125-129.

19. Penna, B. Transform coding techniques for lossy hyperspectral data compression / B. Penna, T. Tillo, E. Magli, G. Olmo // IEEE Geoscience and Remote Sensing Letters. - 2007. - Vol. 45(5). - P. 1408-1421.

20. Tang, X. Hyperspectral image compression using three-dimensional wavelet coding / X. Tang, W. A. Pearlman, J. W. Modestino // Image and Video Communications and Processing. - 2003. - Vol. 5022. - P. 1037-1047.

21. Du, Q. Hyperspectral image compression using JPEG2000 and principal component analysis / Q. Du, J. E. Fowler // IEEE Geoscience and Remote Sensing Letters. - 2007. - Vol. 4. - P. 201-205.

22. Westwater, R. Three-dimensional DCT video compression technique based on adaptive quantizers / R. Westwater, B. Furht // Engineering of Complex Computer Systems Proceedings. - 1996. - P. 189-198.

23. Lee, M. C. Quantization of 3D-DCT coefficients and scan order for video compression / M. C. Lee, R. K. W. Chan, D. A. Adjeroh // Journal of Visual Communication and Image Representation. - 1997. - Vol. 8(4). - P. 405-422.

24. Bozinovic, N. Scan order and quantization for 3D-DCT coding / N. Bozinovic, J. Konrad // Proceedings of Visual Communications and Image Processing. - 2003.

- Vol. 5150. - P. 1204-1215.

25. Karami, A. Hyperspectral image compression using 3D discrete cosine transform and support vector machine learning / A. Karami, S. Beheshti, M. Yazdi // Proceedings of Information Science, Signal Processing and their Applications (IS-SPA). - 2012. - P. 826-829.

26. Qiao, T. Effective compression of hyperspectral imagery using an improved 3D DCT approach for land-cover analysis in remote-sensing applications / T. Qiao, J. Ren, M. Sun, J. Zheng, S. Marshall // International Journal of Remote Sensing. -2014. - Vol. 35(20). - P. 7316-7337.

27. Mankar, P. R. Image compression based on 3D-DCT / P. R. Mankar, S. S. Rane, A. E. Patil // International Journal of Research in Science & Engineering. - 2017.

- Vol. 3(2). - P. 16-22.

28. Haiyan, T. Research on quantization and scanning order for 3-D DCT video coding / T. Haiyan, S. Wenbang, G. Bingzhe, Z. Fengjing // International Conference on Computer Science and Electronics Engineering (ICCSEE) Proceedings. - 2012. -Vol. 1. - P. 200-204.

29. Чичёва, М. А. Сжатие гиперспектральных данных на основе метода кодирования с преобразованием / М. А. Чичёва, Р. Р. Юзькив // Компьютерная оптика. - 2014. - Т. 38, № 4. - С. 798-803.

30. Сергеев, В. В. Параметрическая модель автокорреляционной функции космических гиперспектральных изображений / В. В. Сергеев, Р. Р. Юзькив // Компьютерная оптика. - 2016. - Т. 40, № 3. - С. 416-421.

31. Chicheva, M. A. Transform coding method for hyperspectral data: influence of block characteristics to compression quality / M. A. Chicheva, R. R. Yuzkiv // Communications in Computer and Information Science. - 2015. - Vol. 542. -P. 51-56.

32. Yuzkiv, R. R. Transform-based coding method for remote sensing hyperspectral data compression / R. R. Yuzkiv, V. V. Sergeev // Procedia Engineering. - 2017. -Vol. 201. - P. 249-257.

33. Чичёва, М. А., Юзькив, Р. Р. Сжатие гиперспектральных данных на основе метода кодирования с преобразованием // Тез. докл. междунар. науч. конф. «Интеллектуализация обработки информации» (И0И-10), 4-11 октября 2014 г., Крит. - М.: Торус Пресс, 2014. - С. 226-227.

34. Чичёва, М. А., Юзькив, Р. Р. Влияние характеристик блоков на качество сжатия гиперспектральных данных методом кодирования с преобразованием // Труды Российского научно-технического общества радио-техники, электроники и связи имени А. С. Попова, серия «Цифровая обработка сигналов и её применение» (DSPA-2015), 24-27 марта 2015 г., Москва. - М.: РНТОРЭС, 2015. - Выпуск XVII, Т. 2. - С. 739-742.

35. Юзькив, Р. Р., Сергеев, В. В. Метод кодирования с преобразованием для сжатия космических гиперспектральных изображений // Сб. тр. междунар. конф. «Информационные технологии и нанотехнологии» (ИТНТ-2017), 25-27 апреля 2017 г., Самара. - Самара: Предприятие «Новая техника», 2017. - С. 522528.

36. Юзькив, Р. Р. Кластеризация гиперспектральных изображений для повышения эффективности сжатия методом кодирования с преобразованием // Сб. науч. тр. междунар. науч. -техн. конф. «Перспективные информационные технологии» (ПИТ-2019), 24-26 июня 2019 г., Самара. - Самара: Издательство Самарского научного центра РАН, 2019. - С. 313-317.

37. Юзькив, Р. Р. Модуль сжатия гиперспектральных изображений на основе кодирования с преобразованием. Свидетельство об официальной регистрации программы для ЭВМ № 2018661082 от 11 сентября 2018 г.

38. Ольховский, Ю. Б. Сжатие данных при телеизмерениях / Ю. В. Ольховский, О. Н. Новоселов, А. П. Мановцев. - М: Советское радио, 1971. - 304 с.

39. Huffman, D. A. A method for the construction of minimum-redundancy codes / D. A. Huffman // Proceedings of the IRE. - 1952 - Vol. 40(9). - P. 1098-1101.

40. Rissanen, J. J. Arithmetic coding / J. J. Rissanen, G. G. Langdon // IBM Journal of Research and Development. - 1979. - Vol. 23(2). - P. 149-162.

41. Сойфер, В. А. Введение в цифровую обработку сигналов и изображений: критерии качества изображений и погрешности их дискретного представления : учеб. пособие / В. А. Сойфер, В. В. Сергеев, С. Б. Попов [и др.]. - Самара: Изд-во Самар. гос. аэрокосм. ун-та, 2006. - 36 с.

42. Прэтт, У. К. Цифровая обработка изображений. Книга 1 / У. К. Прэтт. - М.: Мир, 1982. - 312 с.

43. Гонсалес, Р. Цифровая обработка изображений / Р. Гонсалес, Р. Вудс. - М.: Техносфера, 2005. - 1072 с.

44. Тимбай, Е. И. Повышение эффективности компрессии изображений на основе анализа поля ошибок : дис. ... канд. техн. наук : 05.13.17 / Тимбай Елена Ивановна; [Место защиты: Сам. гос. аэрокосм. ун-т им. С. П. Королева]. -Самара, 2010. - 159 с.: РГБ ОД, 61 11-5/328.

45. Изобретаем JPEG. URL: https://habr.com/ru/post/206264/ (дата обращения: 01.10.2016).

46. Ахмед, Н. Ортогональные преобразования при обработке цифровых сигналов / Н. Ахмед, К. Р. Рао. - М.: Связь, 1980. - 248 с.

47. Rosenholtz, R. Perceptual adaptive JPEG coding / R. Rosenholtz, A. B. Watson // Image Processing Proceedings. - 1996. - Vol. 1. - P. 901-904.

48. Yeo, B.-L. Volume rendering of DCT-based compressed 3D scalar data / B.-L. Yeo, B. Liu // IEEE Transactions on Visualization and Computer Graphics. -1995. - P. 29-43.

49. Mastriani, M. 3D zigzag for multislicing, multiband and video processing / M. Mastriani // The Computing Research Repository (CoRR). - 2016. - Vol. 6. -13 p.

50. Сергеев, В. В. Методы цифрового моделирования оптико-электронных систем дистанционного формирования и обработки изображений : дис. ... докт. техн. наук : 05.13.16 / Сергеев Владислав Викторович; [Место защиты: Сам. гос. аэрокосм. ун-т им. С. П. Королева]. - Самара, 1993. - 432 с.

51. Порфирьев, Л. Ф. Основы теории преобразования сигналов в оптико-электронных системах / Л. Ф. Порфирьев. - Санкт-Петербург: Лань, 2013. - 400 с.

52. Сергеев, В. В. Теория цифровой обработки сигналов и изображений / В. В. Сергеев, М. А. Чичёва. - Самара: Изд-во Самарского государственного аэрокосмического университета, 2013. - 206 с.

53. Джайн, А. К. Успехи в области математических моделей для обработки изображений / А. К. Джайн // ТИИЭР. - 1981. - Т. 69, № 5. - С. 9-39.

54. Белокуров, А. А. Стохастические модели в задачах анализа и обработки изображений / А. А. Белокуров, В. В. Сечко // Зарубежная радиоэлектроника. -1989. - Т. 5. - С. 3-18.

55. Сергеев, Г. А. Статистические методы исследования природных объектов / Г. А. Сергеев, Д. А. Янутш. - Ленинград: Гидрометеоиздат, 1973. - 300 с.

56. Мельканович, А. Ф. Фотографические средства и их эксплуатация / А. Ф. Мельканович. - М.: Изд-во Министерства обороны СССР, 1984. - 576 с.

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

58. Вентцель, Е. С. Теория вероятностей: Учеб. для вузов / Е. С. Вентцель. — 6-е изд. стер. — М.: Высш. шк., 1999. — 576 с.

59. Харалик, Р. М. Статистический и структурный подходы к описанию текстур / Р. М. Харалик // ТИИЭР. - 1979. - Т. 67, № 5. - С. 98-120.

60. Виттих, В. А. Обработка изображений в автоматизированных системах научных исследований / В. А. Виттих, В. В. Сергеев, В. А. Сойфер. - М.: Наука, 1982. - 214 с.

61. Чочиа, П. А. Двухмасштабная модель изображения / П. А. Чочиа // Кодирование и обработка изображений. - М.: Наука, 1988. - С. 69-87.

62. Sergeev, V. V. Spectral energy identification method of the linear observation model for remote sensing of the Earth / V. V. Sergeev, A. Yu. Denisova // Pattern Recognition Image Analysis. - 2011. - Vol. 21(2). - P. 321-323.

63. Chakrabarti, A. Statistics of real-world hyperspectral images / A. Chakrabarti, T. Zickler // IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR). - 2011. - P. 193-200.

64. Деза, Е. И., Деза, М. М. Энциклопедический словарь расстояний / Е. И. Деза, М. М. Деза. - М: Наука, 2008. - 444 c.

65. Otsu, N. A treshold selection method from gray-level histogram / N. Otsu // IEEE Transactions on Systems, Man, and Cybernetics. - 1979. - Vol. 9(1). - P. 62-66.

66. MacQueen, J. Some methods for classification and analysis of multivariate observations / J. MacQueen // Proceedings of 5th Berkeley Symposium on Mathematical Statistics and Probability. - 1967. - P. 281-297.

67. Arthur, D. k-means++: the advantages of careful seeding / D. Arthur, S. Vassil-vitskii // Proceedings of the 18th annual ACM-SIAM symposium on Discrete algorithms. - 2007. - P. 1027-1035.

68. International Standard ISO/IEC 14882:2014(E) Programming Language C++.

69. TIFF 6.0 Specification // Adobe Developers Association. - 1992. - 121 p.

70. ENVI Header Files. URL: http://www.harrisgeospatial.com/docs/ENVIHeader-Files.html (access date: 10.03.2017).

71. Чернов, В. М. Быстрые алгоритмы дискретных косинусных преобразований коротких длин с минимальной вычислительной сложностью / В. М. Чернов, М. А. Чичёва // Известия Самарского научного центра РАН. - 1999. - Т. 2. -С. 241-248.

72. Narasimha, M. J. On the computation of the discrete cosine transform / M. J. Nara-simha, A. M. Peterson // IEEE Transactions on Communications. - 1978. -Vol. COM-26(6). - P. 934-936.

73. Fenwick, P. M. A new data structure for cumulative frequency tables / P. M. Fen-wick // Software: Practice and Experience. - 1994 - Vol. 24(3). - P. 327-336.

74. CCSDS 122.0-B-2. Image data compression. Issue 2 // The Consultative Committee for Space Data Systems. - 2017. - 65 p.

75. CCSDS 122.1-B-1. Spectral preprocessing transform for multispectral and hyperspectral image compression. Issue 1 // The Consultative Committee for Space Data Systems. - 2017. - 88 p.

76. QccPack - Quantization, compression and coding Library. URL: http://qccpack.sourceforge.net/ (access date: 10.08.2017).

77. Pearlman, W. A. Efficient, low-complexity image coding with a set-partitioning embedded block coder / W. A. Pearlman, A. Islam, N. Nagaraj, A. Said // IEEE Transactions on Circuits and Systems for Video Technology. - 2004. -Vol. 14(11). - P. 1219-1235.

78. Islam, A. An embedded and efficient low-complexity hierarchical image coder / A. Islam, W. A. Pearlman // Visual Communications and Image Processing, K. Ai-zawa, R. L. Stevenson, Y.-Q. Zhang, Eds., San Jose, CA. - January 1999. -Vol. SPIE 3653. - P. 294-305.

79. Lu, Z. Wavelet coding of video object by object-based SPECK algorithm / Z. Lu, W. Pearlman // Proceedings of the 22nd Picture Coding Symposium (PCS '01), Seoul, Korea. - April 2001. - P. 413-416.

80. Said, A. A new, fast, and efficient image codec based on set partitioning in hierarchical trees / A. Said, W. A. Pearlman // IEEE Transactions on Circuits and Systems for Video Technology. - 1996. - Vol. 6(3). - P. 243-250.

81. Kim, B.-J. An embedded wavelet video coder using three-dimensional set partitioning in hierarchical trees (SPIHT) / B.-J. Kim, W. A. Pearlman // DCC '97 Proceedings of the Conference on Data Compression. - 1997. - P. 251-260.

82. Kim, B.-J. Low bit-rate scalable video coding with 3-D set partitioning in hierarchical trees (3-D SPIHT) / B.-J. Kim, Z. Xiong, W. A. Pearlman // IEEE Transactions on Circuits and Systems for Video Technology. - 2000. - Vol. 10(8). -P. 1374-1387.

ПРИЛОЖЕНИЕ А Тестовый набор гиперспектральных изображений

В качестве материала для экспериментальных исследований были использованы участки (количество строк - 512, количество столбцов - 512, количество каналов - 224) ГСИ, полученных сканером «Aviris» и находящихся в открытом доступе в интернете [1]. Уникальные идентификаторы используемых массивов данных и координаты взятых участков приведены в таблице А.1. Все значения были приведены к диапазону [0; 255]. Сечения изображений, взятые по спектральному каналу № 64, приведены на рисунке 3 на странице 26 данной работы.

Таблица А.1 - Набор тестовых ГСИ

Изображение (участок данных) Уникальный идентификатор массива данных Координаты участка (х, у) в массиве данных

К1 (КаИики-1) ГО70203Ю1р00г02 (200, 6 300)

К2 (КаИики-2) (200, 7 000)

КЗ (КаИики-3) (140, 11 000)

8Б1 (8апБ1едо-1) f111116t01p00r19 (325, 300)

8Б2 (8апБ1едо-2) (300, 1 100)

8Б3 (8апБ1едо-3) (330, 1 750)

81 (8^Иопе-1) f151026t01p00r20 (60, 700)

82 (8^Иопе-2) (530, 6 200)

83 (8^Иопе-3) (600, 6 800)

У1 (Уе11о^^опе-1) fD60925t01p00r10 (100, 100)

У2 (Ye11owstoпe-2) (100, 2 000)

У3 (Ye11owstoпe-3) (50, 4 700)

ПРИЛОЖЕНИЕ Б

Числовые параметры построенных моделей автокорреляционных функций

В таблице Б.1 приведены значения параметров построенных моделей АКФ вида (6). Количество рассматриваемых отсчётов по каждой из координат ограничено значением 50. В качестве точек (х0,у0), (х1,у1), (х2,у2) для оценки параметров пространственной части взяты значения (0,0), (20,20), (21,21). Последний параметр гс является радиусом доверительного интервала, построенного при уровне доверия 0,95 для первых 50 отсчётов по любому из направлений (т. е. для 503 значений АКФ).

Таблица Б.1 - Числовые параметры построенных моделей АКФ

Изображение Тараметры построенной модели

а 0 У Р гг

КаИики-1 0,0032 0,4256 0,0235 0,7733 0,05

КаИики-2 0,0032 0,4282 0,0242 0,8181 0,06

КаИики-Э 0,0038 0,2893 0,0249 0,8165 0,06

8апБ1едо-1 0,0042 0,3668 0,0250 0,8032 0,03

8апБ1едо-2 0,0035 0,3810 0,0264 0,7876 0,04

8апБ1едо-3 0,0021 0,2415 0,0238 0,7938 0,07

БИовИопе^ 0,0065 0,1168 0,0247 0,9148 0,08

БИовИопе^ 0,0037 0,4899 0,0292 0,9943 0,05

БИовИопе-З 0,0047 0,1011 0,0251 0,8921 0,04

Yellowstone-1 0,0056 0,1770 0,0195 0,7695 0,08

Yellowstone-2 0,0045 0,2011 0,0351 0,9030 0,16

Yellowstone-3 0,0066 0,2236 0,0201 0,8835 0,07

ПРИЛОЖЕНИЕ В Исходный код ключевых алгоритмов

Вычисление коэффициентов квантования

// =========

// Quantizer

// =========

void Quantizer::save(DataStream& dataStream, const Array& qm) const { const auto minMax = Statistics::minMax(qm);

const auto bits = System::getHighestBitNumber(static_cast<size_t>(minMax.second)); dataStream.write(bits, CHAR_BIT); for (auto& v : qm) {

dataStream.write(static_cast<size_t>(v), bits);

}

}

Array Quantizer::load(DataStream& dataStream, Size blockSize) const { uint8_t bits = static_cast<uint8_t>(dataStream.read(CHAR_BIT)); Array qm(blockSize); for (auto& v : qm) {

v = static_cast<float>(dataStream.read(bits));

}

return qm;

}

// =================

// AdaptiveQuantizer

// =================

AdaptiveQuantizer::AdaptiveQuantizer(InitAlg initAlg, ErrorAlg errorAlg) : _initAlg(initAlg), _errorAlg(errorAlg) {}

InitAlg AdaptiveQuantizer::getInitAlg() const { return _initAlg;

}

ErrorAlg AdaptiveQuantizer::getErrorAlg() const { return _errorAlg;

}

Array AdaptiveQuantizer::calculate(const Array& spectrum, const Transform& t, const double targetMse, double* outputMse, milliseconds& initTime, milliseconds& loopTime) const {

// Подготовка Array qm(t.size());

GenericArray<double> eCurr(qm.size());

GenericArray<double> eNext(qm.size());

// Вычисление первоначальных коэффициентов квантования

auto start = timer::now();

auto w = calculateWeightArray(t);

double error = 0;

init(spectrum, w, _initAlg, qm); {

//calculateError(spectrum, qm, eCurr, eNext); qm(Point()) = 1;

for (auto it = qm.begin(); it != qm.end(); ++it) { const auto q = static_cast<size_t>(*it); const auto qNext = q / 2 > 1 ? q / 2 : 1;

eCurr(it.abs()) = calculateError(spectrum, qm.size(), it.abs(), q); eNext(it.abs()) = calculateError(spectrum, qm.size(), it.abs(), qNext); error += static_cast<double>(w(it.abs())) * eCurr(it.abs());

}

}

auto end = timer::now();

initTime = duration_cast<milliseconds>(end - start); start = timer::now();

// Итеративное уменьшение коэффициентов квантования while (error > targetMse) {

// Находим коэффициент, который даёт наибольшую ошибку

Point max;

double maxE = 0;

auto noChanges = true;

for (auto it = qm.begin(); it != qm.end(); ++it) { if (*it > 1) {

noChanges = false;

const auto e = w(it.abs()) * (eCurr(it.abs()) - eNext(it.abs())); if (e > maxE) {

max = it.abs(); maxE = e;

}

}

}

if (noChanges || maxE == 0) { break;

}

// Уменьшаем найденный коэффициент квантования

qm(max) = static_cast<float>(static_cast<size_t>(qm(max) / 2));

eCurr(max) = eNext(max);

error -= maxE;

switch (_errorAlg) {

case ErrorAlg::PRECISE: {

const auto qNext = qm(max) / 2 > 1 ? static_cast<size_t>(qm(max) / 2) : static_cast<size_t>(1);

eNext(max) = calculateError(spectrum, qm.size(), max, qNext);

break;

}

case ErrorAlg::COARSE: {

eNext(max) = eNext(max) / 4; break;

}

default:

throw runtime_error("[AdaptiveQuantizer] Неизвестный режим расчёта ошибки.");

}

}

if (outputMse) {

*outputMse = error;

}

end = timer::now();

loopTime = duration_cast<milliseconds>(end - start); return qm;

}

void AdaptiveQuantizer::save(DataStream& dataStream, const Array& qm) const { switch (_initAlg) { case InitAlg::MAX_POWER_OF_TWO: {

const auto minMax = Statistics::minMax(qm);

const auto bits = 1 + static_cast<size_t>(round(log2(minMax.second))); dataStream.write(bits, CHAR_BIT); for (auto& v : qm) {

dataStream.write(static_cast<size_t>(log2(v)), bits);

}

break;

}

case InitAlg::THREE_SIGMA: case InitAlg::MAX: case InitAlg::MAX43: {

Quantizer::save(dataStream, qm); break;

}

default:

throw runtime_error(

"[AdaptiveQuantizer] Неизвестный режим первоначального приближения массива квантования.");

}

}

Array AdaptiveQuantizer::load(DataStream& dataStream, Size blockSize) const { switch (_initAlg) { case InitAlg::MAX_POWER_OF_TWO: {

uint8_t bits = static_cast<uint8_t>(dataStream.read(CHAR_BIT)); Array qm(blockSize); for (auto& v : qm) {

v = static_cast<float>(pow(2, dataStream.read(bits)));

}

return qm;

}

case InitAlg::THREE_SIGMA:

case InitAlg::MAX: case InitAlg::MAX43: {

return Quantizer::load(dataStream, blockSize);

}

default:

throw runtime_error(

"[AdaptiveQuantizer] Неизвестный режим первоначального приближения массива квантования.");

}

}

Array AdaptiveQuantizer::calculateWeight(const VectorTransform& vt) { const auto n = vt.size(); Array buffer(Size(n, 1, 1)); Array w(Size(n, 1, 1));

for (size_t index = 0; index < n; ++index) { for (size_t t = 0; t < n; ++t) {

buffer(Point(t, 0, 0)) = index == t ? 1.0f : 0.0f;

}

auto vector = buffer.vector(Point(), Orientation::HEIGHT, vt.size());

vt.transform(vector, false);

w(Point(index, 0, 0)) = 0;

for (size_t t = 0; t < n; ++t) {

const auto value = buffer(Point(t, 0, 0)); w(Point(index, 0, 0)) += value * value;

}

}

return w;

}

Array AdaptiveQuantizer::calculateWeightArray(const Transform& t) { Array w(t.size());

const auto heightWeights = calculateWeight(t.heightTransform()); const auto widthWeights = calculateWeight(t.widthTransform()); const auto depthWeights = calculateWeight(t.depthTransform()); for (auto it = w.begin(); it != w.end(); ++it) {

const auto wi = heightWeights(Point(it.abs().i, 0, 0)); const auto wj = widthWeights(Point(it.abs().j, 0, 0)); const auto wk = depthWeights(Point(it.abs().k, 0, 0));

*it = static_cast<float>(static_cast<double>(wi) * wj * wk / t.size().total());

}

return w;

}

void AdaptiveQuantizer::calculateError(const Array& spectrum, const Array& qm, GenericArray<double>& eCurr, GenericArray<double>& eNext) {

for (auto it = spectrum.begin(qm.size()); it != spectrum.end(); ++it) { const auto qCurr = static_cast<size_t>(qm(it.rel())); const auto qNext = qCurr / 2 > 1 ? qCurr / 2 : static_cast<size_t>(1); const auto src = static_cast<double>(*it); const auto dstCurr = round(src / qCurr) * qCurr; const auto dstNext = round(src / qNext) * qNext; eCurr(it.rel()) += (dstCurr - src) * (dstCurr - src); eNext(it.rel()) += (dstNext - src) * (dstNext - src);

}

const auto blockN = static_cast<double>((spectrum.size() / qm.size()).total()); eCurr /= blockN; eNext /= blockN;

}

float AdaptiveQuantizer::calculateError(const Array& spectrum, const Size blockSize, const Point p, const size_t q) {

double error = 0;

for (size_t i = p.i; i < spectrum.size().h; i += blockSize.h) {

for (size_t j = p.j; j < spectrum.size().w; j += blockSize.w) {

for (size_t k = p.k; k < spectrum.size().d; k += blockSize.d) {

const auto src = static_cast<double>(spectrum(Point(i, j, k)));

const auto dst = round(src / static_cast<double>(q)) * static_cast<double>(q);

error += (dst - src) * (dst - src);

}

}

}

const auto blockN = (spectrum.size() / blockSize).total(); return static_cast<float>(error / static_cast<double>(blockN));

}

void AdaptiveQuantizer::init(const Array& spectrum, const Array& w, InitAlg initAlg, Array& qm) { switch (initAlg) {

case InitAlg::THREE_SIGMA: {

//GenericArray<double> means(qm. size()); GenenicAmay<double> vans(qm. size());

//for (auto it = spectrum.begin(qm.size()); it != spectrum.end(); ++it) {

// means(it.rel()) += *it; //}

const auto blockN = spectrum.size().h / qm.size().h * spectrum.size().w / qm.size().w * spectrum.size().d / qm.size().d;

for (auto it = spectrum.begin(qm.sizeO); it != spectrum.end(); ++it) {

//const auto value = static_cast<double>(*it) - means(it.rel()) / static_cast<double>(blockN); const auto value = static_cast<double>(*it); vars(it.rel()) += value * value;

}

for (auto it = qm.begin(); it != qm.end(); ++it) {

const auto threeSigma = 3 * sqrt(vars(it.abs()) / static_cast<double>(blockN - 1)); *it = max(1.0f, static_cast<float>(ceil(threeSigma)));

}

break;

}

case InitAlg::MAX_POWER_OF_TWO: {

for (auto it = spectrum.begin(qm.sizeO); it != spectrum.end(); ++it) { const auto q = abs(*it); if (q > qm(it.relO)) { qm(it.rel()) = q;

}

}

for (auto& q : qm) {

q = max(1.0f, pow(2.0f, floor(log2(q))));

}

break;

}

case InitAlg::MAX: case InitAlg::MAX43: {

for (auto it = spectrum.begin(qm.size()); it != spectrum.end(); ++it) { const auto q = abs(*it); if (q > qm(it.rel())) { qm(it.rel()) = q;

}

}

if (initAlg == InitAlg::MAX43) {

for (auto it = qm.begin(); it != qm.end(); ++it) { *it *= 4.0f / 3.0f;

}

}

for (auto& q : qm) { q = ceil(q);

}

break;

}

default:

throw runtime_error(

"[AdaptiveQuantizer] Неизвестный режим первоначального приближения массива квантования.");

}

}

Стационаризация на основе бинаризации

// =======================

// BinarizationCompensator

// =======================

BinarizationCompensator::BinarizationCompensator() : _ratio(nanf("")) { } BinarizationCompensator::BinarizationCompensator(float ratio) : _ratio(ratio) { }

void BinarizationCompensator::compensate(DataStream & dataStream, EntropyCodec & codec, Array & image, const Transform & tr) const { // Инициализация const auto size = image.size(); const auto h = size.h; const auto w = size.w; const auto d = size.d; const auto block = tr.size(); const auto bh = block.h; const auto bw = block.w; const auto bd = block.d;

// Вычисление srcVars (дисперсий неразделённых блоков M1 x M2 x N3) GenericArray<double> srcVars(Size(h / bh, w / bw, 1)); for (auto it = srcVars.begin(); it != srcVars.end(); ++it) { auto var = 0.0;

for (size_t bk = 0; bk < d; bk += bd) {

const auto offset = Size(it.abs().i * bh, it.abs().j * bw, bk); var += Statistics::var(image, offset, block);

}

*it = var;

}

// Определение слоя с наибольшей дисперсией

size_t layer = 0; {

double maxVar = 0;

for (size_t k = 0; k < d; ++k) {

double var = Statistics::var(image, Point(0, 0, k), Size(h, w, 1)); if (var > maxVar) { layer = k; maxVar = var;

}

}

}

// Построение масок

GenericArray<unsigned char> allMasks(Size(h, w, d)); {

const auto minMax = Statistics::minMax(image);

vector<size_t> buffer;

for (size_t i = 0; i < h; i += bh) {

for (size_t j = 0; j < w; j += bw) { for (size_t k = 0; k < d; ++k) {

otsu(image, allMasks, Point(i, j, k), Size(bh, bw, 1), buffer, minMax);

}

}

}

}

// Построение кумулятивных масок GenericArray<int> cumulativeMasks(Size(h, w, 1)); for (size_t bi = 0; bi < h; bi += bh) {

for (size_t bj = 0; bj < w; bj += bw) { for (size_t k = 0; k < d; ++k) { size_t matches = 0; for (size_t i = 0; i < bh; ++i) {

for (size_t j = 0; j < bw; ++j) {

if (allMasks(Point(bi + i, bj + j, k)) == allMasks(Point(bi + i, bj + j, layer))) { ++matches;

}

}

}

bool inverse = matches < bh * bw / 2; for (size_t i = 0; i < bh; ++i) {

for (size_t j = 0; j < bw; ++j) {

const auto value = allMasks(Point(bi + i, bj + j, k)); if (value == 0 || inverse)

--cumulativeMasks(Point(bi + i, bj + j, 0));

else

++cumulativeMasks(Point(bi + i, bj + j, 0));

// Построение итоговой маски

GenericArray<unsigned char> masks(Size(h, w, 1)); for (auto it = masks.begin(); it != masks.end(); ++it) { *it = cumulativeMasks(it.abs()) >= 0;

}

// Вычисление матожиданий Array means(Size(h / bh, w / bw, d / bd)); Array means1(Size(h / bh, w / bw, d / bd)); Array means2(Size(h / bh, w / bw, d / bd)); for (size_t bi = 0; bi < h; bi += bh) {

for (size_t bj = 0; bj < w; bj += bw) {

// Приведение маски к виду count(0) >= count(1)

size_t falseCount = 0;

for (size_t i = 0; i < bh; ++i) {

for (size_t j = 0; j < bw; ++j) {

if (masks(Point(bi + i, bj + j, 0)) == 0) { ++falseCount;

}

}

}

auto trueCount = bh * bw - falseCount; if (falseCount < trueCount) {

for (size_t i = 0; i < bh; ++i) {

for (size_t j = 0; j < bw; ++j) {

masks(Point(bi + i, bj + j, 0)) = !masks(Point(bi + i, bj + j, 0));

}

}

swap(falseCount, trueCount);

}

// Вычисление матожиданий for (size_t bk = 0; bk < d; bk += bd) { auto suml = 0.0; auto sum2 = 0.0;

for (size_t i = 0; i < bh; ++i) {

for (size_t j = 0; j < bw; ++j) {

for (size_t k = 0; k < bd; ++k) {

const auto value = image(Point(bi + i, bj + j, if (masks(Point(bi + i, bj + j, 0)) == 0) { suml += value;

}

else {

sum2 += value;

}

bk + k));

}

}

const auto p = Point(bi / bh, bj / bw, bk / bd);

means(p) = static_cast<float>((sum1 + sum2) / static_cast<double>(bh * meansl(p) = static_cast<float>(sum1 / static_cast<double>(falseCount means2(p) = (trueCount == 0) ? meansl(p) : static_cast<float>(sum2 /

static_cast<double>(trueCount * bd)); }

}

}

// Вычисление dstVars (дисперсий разделённых блоков Ml x M2 x N3)

GenericArray<double> dstVars(srcVars.size());

for (auto it = dstVars.begin(); it != dstVars.end(); ++it) {

bw * bd)); * bd));

const auto bi = it.abs().i * bh; const auto bj = it.abs().j * bw; auto van = 0.0;

for (size_t bk = 0; bk < d; bk += bd) { auto blockVan = 0.0; for (size_t i = 0; i < bh; ++i) {

for (size_t j = 0; j < bw; ++j) {

for (size_t k = 0; k < bd; ++k) {

const auto value = image(Point(bi + i, bj + j, const auto meanPoint = Point(bi / bh, bj / bw, const auto dif = value - ((masks(Point(bi + i,

bk + k)); bk / bd); bj + j, 0)) == 0) ?

meansl(meanPoint)

means2(meanPoint)) ;

blockVar += static_cast<double>(dif) * dif;

}

}