Сегментация спутникового снимка, часть 1

Это первый пример из цикла скриптов,
описывающих типичные задачи цифровой обработки
изображений. В первой части мы познакомимся с
базовыми подходами и алгоритмами обработки изображений
для задачи сегментации спутникового снимка.
Сегментация - это сопоставление пикселей
изображения различным классам.
В самом простом случае можно отделять
объекты на условном "переднем плане"
(foreground) от объектов условного "фона" (background).
В этом случае результатом сегментации
может служить бинарная маска - матрица таких же
размеров, что и исходное изображение, где
элементы 1 соответствуют интересующим нас
объектам, а элементы 0 - фону.
Рассматриваемые темы:

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

Загрузка необходимых библиотек, импорт данных

Следующую строчку кода стоит разкомментировать,
если чего-то из применяемого не хватает:

# Pkg.add(["Images", "ImageShow", "ImageContrastAdjustment", "ImageBinarization", "ImageMorphology"])

Подключение загруженных библиотек:

using Images, ImageShow, ImageContrastAdjustment, ImageBinarization, ImageMorphology

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

I = load("$(@__DIR__)/map_small.jpg")       # загрузка исходного цветного изображения

Предобработка

Как правило, в задачи предобработки
входят анализ изображения для выделения
характерных признаков объектов (отличия их от фона),
изменения контраста и/или яркости,
фильтрация для удаления шумов или сглаживания текстуры,
и проч.

Рассмотрим, какие из алгоритмов предобработки
подойдут в нашем случае. Учитывая, что изображение
цветное, можно попробовать рассмотреть отдельные
каналы (красный, зелёный и синий). Для этого
воспользуемся функцией channelview:

h = size(I, 1);         # высота изображения, пикс
CV = channelview(I);    # разложение изображения на каналы
[ RGB.(CV[1,:,:], 0.0, 0.0) RGB.(0.0, CV[2,:,:], 0.0) RGB.(0.0, 0.0, CV[3,:,:]) ]

Для того, чтобы получить бинарную маску,
нам нужно работать с изображением в
градациях серого (так называемой матрицей интенсивности).
Мы можем перевести цветное изображение в
матрицу интенсивности функцией Gray, или же
рассмотреть один из каналов. Как мы видим,
в синем канале область леса на изображении
сильнее отличается по интенсивности от области
города. И задача сегментации сводится к отделению
более ярких областей пикселей от более тёмных.

[RGB.(CV[3,:,:]) simshow(ones(h,20)) Gray.(I)]

Теперь воспользуемся линейным сглаживающим фильтром,
чтобы выровнять интенсивность близко расположенных
пикселей в синем канале. Сравним два линейных
фильтра - фильтр Гаусса со
скользящим окном размером 9 на 9 пикселей, а также
простой осредняющий фильтр с окном размером 7 на 7.
Применим их к матрице BLUE при помощи функции
imfilter, которой мы передадим разные кернелы.

BLUE = CV[3,:,:];                                       # интенсивность синего канала RGB
gaussIMG = imfilter(BLUE, Kernel.gaussian(2));          # линейный фильтр Гаусса
kernel_average = ones(7,7) .* 1/(7^2);                  # кернел среднего арифметического
avgIMG = imfilter(BLUE, kernel_average);                # линейный фильтр среднего
[simshow(gaussIMG) simshow(ones(h,20)) simshow(avgIMG)]