Практика 9. Паншарпенинг

Скачивание данных Landsat 8

Откройте плагин: Модули → STAC API Browser Plugin → Open STAC API Browser.

Рис. П9.1 — путь к плагину: Модули → STAC API Browser Plugin → Open STAC API Browser

В источнике Earth Search найдите коллекцию Landsat-8 L1 Collection-1, задайте диапазон дат и охват по территории интереса:

Рис. П9.2 — панель поиска STAC API Browser: источник Earth Search, коллекция Landsat-8 L1 Collection-1, фильтр по дате и охвату

Скачайте по интересующей территории каналы 2 (Blue), 3 (Green), 4 (Red) и 8 (Panchromatic) — панхроматический канал понадобится для повышения пространственного разрешения.

Рис. П9.3 — список доступных каналов Landsat 8 в STAC API Browser — выбраны B2, B3, B4 и B8

Канал	Длины волн, мкм	Разрешение
Диапазоны OLI (Operational Land Imager)
1 — Побережья и аэрозоли (Coastal/Aerosol)	0,433–0,453	30 м
2 — Синий (Blue)	0,450–0,515	30 м
3 — Зелёный (Green)	0,525–0,600	30 м
4 — Красный (Red)	0,630–0,680	30 м
5 — Ближний ИК (NIR)	0,845–0,885	30 м
6 — Коротковолновый ИК (SWIR 2)	1,560–1,660	30 м
7 — Коротковолновый ИК (SWIR 3)	2,100–2,300	30 м
8 — Панхроматический (Panchromatic)	0,500–0,680	15 м
9 — Перистые облака (Cirrus)	1,360–1,390	30 м
Диапазоны TIRS (Thermal Infrared Sensor)
10 — Дальний ИК (TIR1)	10,30–11,30	100 м
11 — Дальний ИК (TIR2)	11,50–12,50	100 м

Панхроматический канал специально снимается с более высоким пространственным разрешением (15 м), чем цветные каналы (30 м) — именно эта разница в разрешении и используется методами паншарпенинга.

Паншарпенинг

Паншарпенинг (pan-sharpening, panchromatic sharpening) — метод увеличения пространственного разрешения растрового изображения. Суть метода — использование канала панхроматического диапазона (снятого с более высоким разрешением) для «резкости» многоканального набора данных: результат сохраняет цвета исходных каналов, но с детализацией панхроматического канала.

Классификация методов

Методы повышения пространственного разрешения делятся на несколько крупных групп в зависимости от того, что именно используется для переноса детализации панхроматического канала на цветное изображение:

Основаны на замене компонент
IHS, PS, GS

Основаны на спектральном вкладе
Brovey, IM, P+XS

Основаны на внедрении деталей
HPF, HPM

Основаны на статистике изображения

Основаны на кратном увеличении разрешения

Нейросетевые методы
SRCNN, SRGAN, EDSR

Рис. П9.1 — основные группы методов паншарпенинга; в этой практике разбираются метод IHS (замена компонент) и метод Brovey (спектральный вклад)

Метод IHS (замена компонент)

Метод IHS (Intensity-Hue-Saturation, интенсивность-оттенок-насыщенность) раскладывает цветное изображение на пространственную и хроматическую составляющие: интенсивность (I) отвечает за детализацию и яркость, а оттенок и насыщенность (H, S) — за собственно цвет.

Рис. П9.4 — цветное изображение (слева) раскладывается на хроматическую компоненту — оттенок и насыщенность (справа сверху) — и пространственную компоненту — интенсивность (справа снизу)

Алгоритм метода:

Разложить цветное изображение на компоненты I, H, S.
Заменить компоненту интенсивности I на панхроматический канал (предварительно растянутый по гистограмме под диапазон значений I).
Выполнить обратное преобразование H, S и нового I обратно в цветовое пространство RGB.

В QGIS:

Объедините цветные каналы (кроме панхроматического) инструментом Растр → Прочее → Слияние растров, отметив каналы 2, 3, 4 и включив опцию «Поместить каждый исходный файл в отдельный канал» — иначе каналы наложатся друг на друга вместо того, чтобы сформировать RGB-композит (см. также Практику 3).

Рис. П9.5 — диалог «Слияние растров»: выбрано 3 исходных слоя, отмечена опция «Поместить каждый исходный файл в отдельный канал»

Откройте панель инструментов (Анализ → Панель инструментов) и найдите алгоритм «Повышение пространственного разрешения».
В качестве спектрального (цветного) набора данных укажите полученный композит, а в качестве панхроматического — канал 8; алгоритм передискретизации можно оставить кубическим (см. Практику 4 про методы передискретизации).

Рис. П9.6 — диалог «Повышение пространственного разрешения»: спектральный набор данных, панхроматический набор данных (канал B8), алгоритм передискретизации — кубическая (ядро 4×4)

Типичная ошибка

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

Рис. П9.7 — результат паншарпенинга методом IHS с заметной голубой дымкой над пустынной территорией (Landsat 8)

Рис. П9.8 — та же территория по данным Sentinel-2A — естественный цвет без постороннего оттенка, для сравнения

Причина — панхроматический канал Landsat 8 (0,500–0,680 мкм) захватывает более широкий и немного иной спектральный диапазон, чем комбинация каналов 2–4, поэтому подстановка его как новой интенсивности слегка искажает цветовой баланс. Исправить это можно вручную через тонирование слоя: Свойства слоя → Стиль → Тонирование, выбрав цвет, противоположный паразитному оттенку (для голубой дымки — тёплый, например оранжево-жёлтый), с интенсивностью порядка 20–25%.

Рис. П9.9 — диалог свойств слоя: включена опция «Тонировать», выбран корректирующий цвет, интенсивность 23%

После коррекции результат заметно чище:

Рис. П9.10 — результат паншарпенинга после тонирования: чёткая застройка без постороннего цветового оттенка

Метод Brovey (спектральный вклад)

Метод Brovey (Brovey Transform, BT) основан на нормализации цветовых характеристик пикселов — пропорциональном распределении интенсивности между каналами так, чтобы их сумма соответствовала панхроматическому значению. Сначала для каждого пиксела вычисляется общая интенсивность \( I \) как взвешенная сумма цветовых каналов:

\[ I = \frac{0{,}42 \cdot Blue + 0{,}98 \cdot Green + 0{,}6 \cdot Red}{2}, \]

а затем каждый цветовой канал пересчитывается пропорционально отношению панхроматического значения к этой интенсивности:

\[ MS_{pan} = MS \cdot \frac{PAN}{I}, \]

где \( MS \) — исходное значение пиксела в одном из цветных (многоспектральных) каналов, а \( PAN \) — значение того же пиксела в панхроматическом канале. Веса при Blue, Green и Red в формуле для \( I \) подобраны так, чтобы учесть разный вклад каналов в яркость, воспринимаемую панхроматическим сенсором; при умножении на \( PAN/I \) цветовые пропорции между каналами сохраняются, но яркость и детализация приходят от панхроматического канала.

Реализация в QGIS

В отличие от IHS, для метода Brovey в QGIS нет готового алгоритма — его собирают вручную в Конструкторе моделей (см. также Практику 5) из калькуляторов растров и инструментов передискретизации:

Рис. П9.11 — схема модели Brovey в конструкторе моделей QGIS: входы Red, Green, Blue, Pan → передискретизация → калькулятор интенсивности I → калькуляторы MS·PAN/I по каждому каналу → объединение в виртуальный растр Res

Каждый цветной канал и панхроматический канал сначала приводятся к единому охвату и разрешению узлом «Деформация (перепроецирование)».
Отдельный калькулятор растров считает интенсивность \( I \) по приведённой выше формуле.
Для каждого цветного канала отдельный калькулятор растров считает \( MS \cdot PAN / I \).
Три результата объединяются инструментом «Создать виртуальный растр» в конечное изображение.

Типичная ошибка

Во всех узлах-калькуляторах модели нужно вручную задать одинаковый «Целевой охват» (Target extent) — если оставить это поле пустым, разные узлы могут рассчитать результат по разным охватам, и при объединении в конце цепочки каналы не совпадут по геометрии.

Рис. П9.12 — диалог калькулятора растров с выделенным полем «Целевой охват», которое нужно задать одинаковым во всех узлах модели

Рис. П9.13 — исходный композит с разрешением 30 м (слева) и результат паншарпенинга методом Brovey с разрешением 15 м (справа)

Типичная ошибка

Метод Brovey хорошо сохраняет цветовой баланс, но не гарантирован от искажений там, где интенсивность \( I \) близка к нулю (очень тёмные пикселы) — деление на малое число может дать выброс яркости. На практике такие пикселы стоит проверять отдельно, а не считать результат метода одинаково надёжным по всему снимку.

Домашнее задание

Подготовленные в этой практике паншарпенированные снимки Landsat (лето 2016 и лето 2021 года) используются в Лабораторной работе №2: на их основе собирается обучающая выборка в SCP, обучается классификатор SVM (см. Практики 6–8) для нескольких классов местности, и по итогам классификации обоих лет делается вывод о динамике развития выбранной территории — какие классы приросли, а какие сократились за пять лет.