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

На сегодня реализовано большое количество алгоритмов сжатия, используемых для растровых данных (изображений). Они делятся на алгоритмы с потерей и без потери качества. В областях, связанных с обработкой данных ДЗЗ, как правило применяются алгоритмы без потери качества, чтобы избежать искажения данных при расчётах. Среди таких алгоритмов наибольшей степенью сжатия отличается алгоритм JPEG2000, но для него характерна невысокая скорость сжатия и во много раз большее время чтения, чем для большинства алгоритмов этого класса. Поэтому данный алгоритм сжатия плохо подходит для систем, обеспечивающих оперативный (особенно интерактивный) доступ к данным ДЗЗ. Большинство остальных алгоритмов, таких как LZW, DEFLATE, LZMA, ZSTD и другие, построено на базе алгоритмов Лемпеля-Зива-Уэлча и кодов Хаффмана. На текущем этапе развития их эффективность и быстродействие более-менее сопоставимы и зависят от природы сжимаемых данных. Практически все подобные алгоритмы имеют средний предел эффективности степени сжатия порядка 2,3 для данных разрядности Int16 и порядка 1,3 для данных разрядности Float32, выше которого можно перейти только в редких случаях (см., например https://kokoalberti.com/articles/geotiff-compression-optimization-guide).

Использование классических алгоритмов сжатия с потерей качества, таких, как JPEG, позволяет существенно увеличить достижимую степень сжатия данных по сравнению с предыдущими алгоритмами, однако область их применения ограничивается данными, которые используются только для визуализации и не предполагают «количественной» обработки и анализа. Это связано с тем, что ошибка, вносимая в данные алгоритмом JPEG и подобными, непредсказуема и может в конечном итоге привести к существенным искажениям как данных, так и полученных результатов их обработки.
 
В 2015 г. стал доступен новый алгоритм сжатия с потерями — LERC (Limited Error Raster Compression, https://esri.github.io/lerc/, Becker, 2016), который позволяет задать и гарантировать максимальную вносимую ошибку в значение каждого пиксела  (MaxError), тем самым сохраняя радиометрические свойства изображения с определённой погрешностью. Алгоритм основан на том, что данные разбиваются на небольшие блоки, в каждом из блоков минимальное значение пиксела берётся за основу, вычисляется разность остальных значимых пикселов с ним, полученная разность делится на 2*MaxError (где MaxError — заданный уровень допустимой ошибки) и округляется до целого, после чего полученный блок сжимается без потерь. Таким образом, достигается уменьшение размерности данных за счет перевода данных в целочисленные значения. Большая часть значений при этом изменяется, но не более чем на величину установленного порога допустимой ошибки. Результатом такого изменения является уменьшение словарей для последующего сжатия без потерь. При этом естественно, что повышение порога допустимой ошибки приводит к более существенному сокращению словарей и, как следствие, к более высоким степеням сжатия. Опыт использования данного алгоритма показывает, что даже при относительно низких допустимых значениях допустимой ошибки достигаются значительные уровни сжатия (https://kokoalberti.com/articles/geotiff-compression-optimization-guide/).

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

В качестве основного набора тестовых данных были использованы восстановленные композитные изображения, очищенные от облачности, помех и шумов, полученные на основе данных прибора MSI, установленного на спутниках серии Sentinel-2. Такой выбор был обусловлен тем, что эти данные активно используются в работах, идущих в ЦКП «ИКИ-Мониторинг», уже занимают достаточно большой объем в архивах, а в будущем планируется многократное увеличение количества таких продуктов. Кроме этого, эти данные менее эффективно сжимаются стандартным алгоритмом DEFLATE. Данные хранятся в формате Int16, средняя степень сжатия для них составляет всего 1,33, в то время как степень сжатия исходных сцен данных этого прибора в среднем превышает 2, что, вероятно, объясняется наличием в исходных данных однородных участков, покрытых облачностью. Следует также отметить, что проведение анализа именно на предлагаемом типе данных позволяет оценить возможности использования алгоритма LERC при сжатии данных, свободных от облачности, которые в основном и используются при изучении различных процессов и явлений, происходящих на поверхности Земли.

Для построения выборки использовались фрагменты композитных изображений по каналам NIR и RED, а также композита вегетационного индекса NDVI. В качестве дополнительного тестового набора использовались аналогичные информационные продукты, полученные по данным прибора MODIS (спутники TERRA, AQUA), отличающегося существенно меньшим пространственным разрешением и намного более широкой полосой обзора.

Для каждого из этих информационных продуктов были подготовлены выборки, состоящие из примерно 500 фрагментов за разные даты и с разным географическим расположением. Затем для каждого из экземпляров данных были получены сжатые алгоритмом LERC изображения с максимальной ошибкой от 0 до 5 процентов от исходного диапазона значений. После этого были получены зависимости среднего объема сжатых данных от используемой ошибки LERC. Пример такой зависимости для восстановленных композитных изображений индекса NDVI (диапазон значений от 0 до 10000), полученных по данным прибора MSI, приведен на рисунке 1. На графике красной пунктирной линией показан объем несжатых данных, принятый за 100%, зеленой — объем данных, сжатых при помощи алгоритма DEFLATE, а синей — процент от исходного объема, получающийся при сжатии алгоритмом LERC при разных значениях максимальной ошибки. По оси абсцисс указана ошибка, как процент от возможного диапазона значений в исходных файлах. Также на графике оранжевым цветом отображена кривая, показывающая объем, получающийся после перепаковки алгоритмом DEFLATE данных, которые были сжаты алгоритмом LERC. Интерес к этой зависимости вызван тем, что алгоритм сжатия LERC поддерживается только в актуальных версиях программных пакетов работы с изображениями (прежде всего это пакет GDAL, см. https://gdal.org/en/latest/), что может привести к проблемам с совместимостью для более старых версий ПО работы с данными. Из графика видно, что при использовании алгоритма LERC c максимальной ошибкой в 1% от диапазона значений достигается экономия дискового пространства более чем в 3 раза относительно данных, сжатых алгоритмом DEFLATE, и даже при использовании компромиссного варианта — более чем в два раза относительно них же.

Рисунок 1 — Эффективность сжатия ежедневных композитов индекса NDVI, полученного по данным прибора MSI

Аналогичная зависимость для данных композитов индекса NDVI, полученных по данным прибора MODIS, приведена на рисунке 2. Степень сжатия алгоритмом DEFLATE у этих данных оказалась существенно лучше, в то время как эффективность сжатия при помощи алгоритма LERC отличается только для малых значений допустимой ошибки LERC. Эти композиты получаются по аналогичной методике, а отличия, судя по всему, связаны с существенно более низким пространственным разрешением этих данных.

Рисунок 2 — Эффективность сжатия ежедневных композитов индекса NDVI, полученного по данным прибора MODIS

На приведённых графиках хорошо видно, что даже относительно небольшие значения вносимой ошибки позволяют получить степени сжатия порядка 4–5. Тем не менее для практического применения алгоритма LERC необходимо оценить порог допустимой максимальной ошибки, при которой вносимые при сжатии искажения не будут критичными при обработке данных. С точки зрения радиометрических характеристик изображения, такую ошибку можно получить на основе оценки радиометрической точности исходных данных, задав её значение, не превышающее половины исходной точности данных. Так, например, для исходных сцен данных спутников Sentinel-2 она составляет 3–5% (см. ежегодный отчёт о функционировании системы Sentinel-2 за 2023 г., https://sentiwiki.copernicus.eu/__attachments/1673423/OMPC.CS.APR.003%20-%20S2%20MSI%20Annual%20Performance%20Report%202023%20-%201.0.pdf?inst-v=5d892f0e-0a4f-4a11-88a9-f6f3caa8d5c0), а для используемых в выборке восстановленных композитных изображений - порядка 5%. Таким образом, для того, максимальная ошибка сжатия для этих данных будет составлять 1,5–2,5 и 2,5% соответственно.

Вторым важным фактором при оценке допустимой ошибки сжатия является сохранение структуры изображения, так как при высоких степенях потери в силу особенностей работы алгоритма с блоками данных она заведомо портится, что влияет на решение таких задач, как детектирование объектов на изображении, расчёт средних по объектам и тому подобное. Для выявления пространственных искажений анализировались изображения, получаемые как разность между исходными и сжатыми данными. Если разность близка к шуму, искажений в пространственной структуре сжатого изображения нет. Наличие в разности явной структуры является признаком того, что пространственная структура изображения была искажена. Существенно, что характеристики таких изображений зависят не только от конкретных типов спутниковых данных, но и от конкретного региона наблюдения и времени наблюдения, однако на качественном уровне все закономерности сохраняются. Для иллюстрации того, как меняется структура изображения при разных ошибках сжатия, ниже анализируется фрагмент изображения индекса NDVI размером 1000x1000 пикселов, полученного по данным спутников Sentinel-2.

Полученная зависимость разности между оригинальным и сжатыми изображениями от допустимой ошибки показана на рисунке 3. В первой колонке таблицы указана заданная ошибка LERC (как процент от диапазона значений), во второй полученная разность, в третьей — гистограмма значений в ней, а в четвертой — соответствующее сжатое изображение. При ошибке в 1% изображение, соответствующее разности, близко к шуму, что подтверждается его гистограммой. Начиная с ошибки в 2%, на разности начинают проявляться однородные области (структура), а гистограмма перестает быть равномерной. При последующем увеличении ошибки неравномерность гистограммы также всё более увеличивается. При этом в разнице изображений появляется текстура, аналогичная текстуре, присутствующей на несжатом изображении, что говорит о том, что часть пространственной информации была потеряна при сжатии. Сами сжатые изображения при этом остаются визуально практически неотличимыми от оригинального, а отличия заметны только при переходе на пиксельный уровень. Также при увеличении допустимой ошибки начинает проявляться эффект занижения значений в сжатом изображении, который усиливается с увеличением заданной ошибки. Такой эффект наблюдается при сжатии алгоритмом LERC изображений, в которых присутствуют сравнительно однородные области, значения в которых варьируются меньше, чем на заданную ошибку.

Таким образом, несмотря на то что с точки зрения сохранения радиометрических свойств изображения допустимая ошибка в 2,5–3% не является существенной, она уже искажает пространственную текстуру оригинального изображения, поэтому именно критерий сохранения пространственной структуры изображений является определяющим при выборе максимальной ошибки при использовании алгоритма сжатия LERC.

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

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

Поэтому были выбраны статистические критерии для оценки искажений, вносимых в пространственную структуру изображений. Как было показано выше, при появлении пространственной структуры в разностях между исходным и сжатым изображениями их гистограммы «перекашиваются» вправо, поэтому в качестве первого критерия было выбрано среднее значение по разнице. Однако такой показатель характеризует только изображение в целом. В качестве более детального способа выявления структуры в разнице между изображениями был использован индекс автокорреляции изображения со сдвигом 1. Так как размеры изображений в пикселах, как правило, достаточно велики, искажения, связанные с использованием одномерного массива значений вместо двумерного, можно считать несущественными. По сути, этот индекс характеризует вероятность того, что соседние пиксели в сжатом изображении были смещены схожим образом, что и означает наличие текстуры в анализируемой разнице.

Рисунок 3 —  Разности между оригинальным и сжатым изображениями при использовании разных допустимых ошибок LERC и соответствующие им гистограммы значений. Для изображений используется палитра, которая пропорционально распределяется от минимальных до максимальных значений, присутствующих в пикселях изображений

Для проведения экспериментов использовалась описанная выше выборка восстановленных композитов индекса NDVI (диапазон значений от нуля до 10000), полученных по данным прибора MSI, установленного на спутниках серии Sentinel-2. Для каждого сжатого изображения были рассчитаны средние значения по разности, индексы автокорреляции, а также степень сжатия. Потом для каждого из значений максимальных ошибок сжатия были получены средние (агрегированные) значения этих показателей. Для полученных агрегированных значений была построена корреляционная матрица, которая показала очень высокую корреляцию (около 0,98) между средним по разности и автокорреляцией со сдвигом 1. Это позволяет предположить, что среднее значение по разности между исходными и сжатыми изображениями так же хорошо характеризует вносимые пространственные искажения, как и индекс автокорреляции со сдвигом 1. Примерно такая же высокая степень корреляции была обнаружена у рассматриваемых показателей и со степенью сжатия изображений, вычисляемой как отношение размера исходного изображения к размеру сжатого. Поскольку все выбранные показатели хорошо коррелируют, для оценки пространственных искажений в дальнейшем достаточно выбрать какой-либо один из них, например – среднее по разности, как наиболее простое в вычислении. На рисунке 4 приведена двумерная гистограмма зависимости между индексом автокорреляции со сдвигом 1 (представлен по оси X) и величиной средних значений по разности между исходными и сжатыми изображениями (представлены по оси Y), построенная для наборов данных с максимальной ошибкой сжатия в 0,5, 1, 1,5, 2 и 3%. На представленной гистограмме видно, что для малых значений автокорреляции и среднего, которые как раз представляют основной интерес, эти два показателя действительно хорошо коррелированы между собой, на больших же значениях параметров корреляция падает.

Рисунок 4 — Гистограмма корреляции индекса автокорреляции и среднего по разности изображений

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

Рисунок 5 — MSI NDVI. Зависимость средних значений по разности изображений от максимальной ошибки LERC

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

Для подтверждения выбранного порогового значения нами были рассмотрено поведение гистограмм и агрегированных средних по разности для значений максимальной ошибки сжатия в диапазоне от 0,25 до 1,5%, представленных на рисунке 6. Ошибка дана в процентах от исходного диапазона, средние — в абсолютных значениях. Нулевые значения были выколоты для большей детальности приведённых гистограмм.

Рисунок 6 — Поведение гистограмм разностей и средних значений для ошибок LERC в диапазоне 0,25–1,5%

Хорошо видно, что при малых ошибках сжатия гистограммы значений в разностях остаются стабильными и отражают распределение значений, близкое к шуму, среднее значение по разности так же остаётся стабильным. После превышения выбранного порога в 1% гистограмма перестаёт быть стабильной (см. также рис. 3), средние значения по разности также резко начинают расти (для ошибки в 1,5% среднее по разности составляет 1,4, для ошибки в 2% — 4,5).  Степень сжатия алгоритмом LERC с использованием такой максимальной ошибки относительно изображений, сжатых при помощи алгоритма DEFLATE, равна примерно 3,3 (т.е. данные сжимаются более, чем в три раза).

Для дополнительного обоснования такого выбора была также проанализирована гистограмма распределения средних ошибок по разницам значений для всех файлов в выборке для разных ошибок LERC, приведенная на рисунке 7. При ошибках больше, чем 1%, резко увеличивается дисперсия распределения. Это означает, что для значительной части изображений в тестовой выборке пространственные искажения существенно выше, чем в среднем, и выбранный нами порог корректен.

Рисунок 7 — Распределение средних ошибок по разницам изображений в выборке для разных допустимых ошибок LERC

Для получения методики также необходимо было оценить корректность подхода и достаточность созданной выборки из 500 элементов для получения распределения значений индекса автокорреляции. Для этого была использована упрощенная формула для оценки минимального размера выборки для достижения требуемой точности при определении среднего по ней значения: MSS=((2*S)/P)², где S — стандартное отклонение, а P — требуемая точность.
    Из этой формулы можно получить выражение для получения точности для известного размера выборки:

Для распределения средних значений, соответствующих выявленному пороговому значению ошибки в 1% (абсолютное значение ошибки 100), дисперсия средних по разницам изображений равна примерно равна 0,4. Согласно приведенной формуле, выборка из 500 элементов позволяет получить среднее по этой величине с точностью лучше, чем 0,036. Для достижения же допустимой точности в 0,05 достаточно и вдвое меньшей выборки. Можно предположить, что и для других типов спутниковых изображений оценки необходимого размера выборки будут похожими.

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

 

1. Подготовка случайной по времени и географическому расположению выборки спутниковых изображений анализируемого типа, состоящую из нескольких сотен экземпляров данных;
2. Получение набора сжатых алгоритмом LERC изображений для ряда допустимых ошибок LERC в диапазоне от 0,5 до меньшего из половины диапазона радиометрической точности данных или 2% от диапазона значений в изображениях для нескольких экземпляров из выборки;
3. Получение зависимости агрегированного среднего по разностям между оригинальными и сжатыми изображениями от максимальной ошибки LERC и выбор значения, при котором эта величина максимально близка к нулю, в качестве порогового;
4. Повтор шагов 2 и 3 для всей выборки для диапазона максимальных ошибок сжатия +/- 0,5% от полученного порогового значения для уточнения полученного порога.

Предложенная методика была применена для получения порогового значения сжатия для восстановленного композита канала RED по данным спутников Sentinel-2. Диапазон используемых значений для этого продукта существенно уже (около 2500), а полученное пороговое значение ошибки LERC равно примерно 1,2% (около 30). При этом достигается относительная степень сжатия около 3,32. На рисунке 8 приведена зависимость агрегированного среднего значения по разнице между оригинальными и сжатыми изображениями для восстановленного композита канала RED по данным спутников Sentinel-2.

Также методика была применена для восстановленных композитов индекса NDVI по данным прибора MODIS, реализованных в том же диапазоне значений, что и аналогичный продукт по данным спутников Sentinel-2. На рисунке 9 приведена аналогичная зависимость для выборки. Полученное пороговое значение ошибки LERC равно примерно 1,2% от диапазона значений, при этом достигается относительная степень сжатия данных около 2,9. График зависимости средних по разности от максимальной ошибки сжатия приведён на рисунке 9.

Рисунок 8 — Зависимость агрегированных средних значений по разнице изображений от максимальной ошибки LERC для данных Sentinel 2 RED

 
Рисунок 9 — Зависимость агрегированных средних значений по разности изображений от максимальной ошибки LERC для данных MODIS NDVI

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

Таким образом, ключевым критерием, больше всего влияющим на выбор порогового значения ошибки для этого алгоритма, является сохранение пространственной структуры спутниковых изображений. Разработанная для определения допустимого порогового значения ошибки методика основана на статистических критериях для оценки вносимых в сжатые изображения искажений. Результаты применения созданной методики для созданных тестовых выборок показали, что использование алгоритма LERC с полученными значениями максимальной ошибки позволяет сократить объем данных в архиве более чем в три раза относительно данных, сжатых при помощи алгоритма DEFLATE, не потеряв при этом радиометрических и структурных свойств данных, а сам алгоритм LERC имеет широкие перспективы применения для реальных задач обработки и анализа данных при условии соблюдения критериев, обеспечивающих сохранность свойств изображения. Так, его применение к накопленным в ЦКП «ИКИ-Мониторинг» архивам исходных данных КА Sentinel-2 уровней обработки L1С и L2A, объём которых превышает 3 Пб, потенциально может уменьшить эту величину примерно на 1.5–2 Пб.