Как найти изменения на спутниковых снимках и опубликовать результат

Привет! Я решил проверить, можно ли с помощью ИИ-агента собрать рабочий процесс мониторинга территории: найти открытые спутниковые снимки, выделить изменения за год и опубликовать результат на картографическом сервисе, чтобы его могли посмотреть коллеги.

Для эксперимента выбрал участок Новой Москвы. Задача — найти места, где могли появиться новые строительные площадки и искусственные поверхности. Вместо просмотра всей территории вручную мне хотелось получить предварительный список участков для проверки и сохранить его на карте вместе с исходными данными.

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

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

Территория и условия эксперимента

Для проверки взял прямоугольный участок около Троицка, Десны, Красной Пахры и Щапова. Исходные границы заданы bounding box [37.18, 55.39, 37.46, 55.61] в WGS84. После построения охватывающей сетки в UTM площадь составила 455,60 км². Прямоугольник охватывает часть Новой Москвы и прилегающую территорию.

Для эксперимента взял открытые снимки Sentinel-2 и модели с готовыми весами. Результат планировал передать коллегам по ссылке: с возможностью сравнить снимки, включить нужный слой и посмотреть свойства выделенного участка.

1. Поиск и отбор снимков

Снимки искал в Earth Search от Element 84. Каталог поддерживает STAC API — стандартный интерфейс поиска геопространственных данных: в запросе можно задать территорию, даты и облачность, а в ответе получить метаданные сцен и ссылки на каналы.

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

Пример запроса поиска для первого сезона:

import requests

response = requests.post(
    "https://earth-search.aws.element84.com/v1/search",
    json={
        "collections": ["sentinel-2-l2a"],
        "bbox": [37.18, 55.39, 37.46, 55.61],
        "datetime": "2025-06-01T00:00:00Z/2025-08-31T23:59:59Z",
        "query": {"eo:cloud_cover": {"lte": 15}},
        "limit": 100,
    },
    timeout=60,
)
response.raise_for_status()
items = response.json()["features"]
# Если выдача разбита на страницы, нужно обработать ссылку rel="next".

STAC-описания выбранных сцен сохранил вместе с данными: в них остаются ссылки на исходные файлы и параметры каналов.

Низкая облачность всей сцены ещё не гарантирует, что нужный участок открыт. Для локальной проверки использовал SCL — Scene Classification Layer, слой классификации сцены. Он входит в продукт Sentinel-2 L2A и содержит, в частности, метки облаков, теней, снега и отсутствующих данных. По этим меткам можно исключить непригодные пиксели. Описание SCL и уровней обработки Sentinel-2.

После проверки подошли три сцены тайла 37UCB:

Проценты в таблице рассчитаны по меткам SCL. Основную пару разделяют 355 дней, два снимка 2026 года — 20 дней.

2. Подготовка каналов и общей сетки

Чтобы модели сравнивали одну и ту же поверхность, снимки нужно привести к общей сетке и подготовить каналы в том виде, который ожидают их входы.

У Sentinel-2 есть несколько уровней обработки. L1C содержит значения отражения на верхней границе атмосферы, L2A — после атмосферной коррекции, на уровне поверхности. Для просмотра взял RGB из L2A, оттуда же — SCL для маски качества. Спектральные входы выбранных нейросетей подготовил из соответствующих сцен L1C.

Для всех дат задал общую сетку:

• система координат EPSG:32637, UTM 37N;

• размер 1829 × 2491 пиксель;

• шаг 10 м;

• одинаковые границы и положение пикселей для всех дат.

Спектральные каналы пересчитал билинейной интерполяцией. Для категорий SCL использовал метод ближайшего соседа, чтобы не смешивать коды классов. Часть каналов имеет исходное разрешение 20 и 60 м. На общей сетке их значения распределены по более мелким пикселям, но детализация остаётся прежней.

Значения каналов преобразовал с учётом scale и offset из метаданных. Входы L1C привёл к гармонизированным DN — числовым значениям в масштабе, который использует подготовка данных модели, — затем нормализовал для инференса.

Облачные пиксели и тени исключил вместе с буфером 20 м вокруг них. Из анализа также убрал дефекты, пропуски и, для сравнения дат, SCL-класс снега. Общая пригодная для сравнения площадь трёх дат составила 455,46 км².

Совмещение проверил на восьми фрагментах двух пар: сдвиг изображения мог бы дать ложные изменения вдоль границ объектов. На этих фрагментах изображения совпали с точностью до целого пикселя. Возможные смещения на доли пикселя стоит учитывать при просмотре тонких полос вдоль границ объектов.

3. Выбор моделей

Для поиска изменений я рассмотрел открытые реализации нейросетей и доступные предобученные веса. Критерии были практическими: работа с Sentinel-2, возможность запустить всю модель без дообучения и доступ к оценкам на уровне пикселей.

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

Классы поверхности: Dynamic World

Dynamic World — проект классификации земной поверхности по Sentinel-2, разработанный Google совместно с World Resources Institute. Его модель различает девять классов: воду, древесную и травянистую растительность, затопленную растительность, сельхозугодья, кустарники, застройку, открытый грунт, снег и лёд.

Для этой задачи модель удобна тем, что возвращает оценки всех классов на каждом пикселе. Это позволяет отслеживать рост оценки застройки между датами и использовать его как условие отбора участков. Для локального запуска выбрал PyTorch-порт Dynamic World.

На вход подал девять каналов: B02, B03, B04, B05, B06, B07, B08, B11 и B12. Для каждого пикселя сохранил наиболее вероятный класс (argmax), оценки всех девяти классов и максимальную из них. Эти оценки затем использовал при отборе изменений.

Чтобы обрабатывать снимок частями и уменьшить влияние краёв фрагментов, инференс выполнялся окнами 320 × 320 пикселей с дополнительным контекстом 32 пикселя. В итоговый растр переносилась центральная часть предсказания.

На снимке 6 августа у 23,39% пикселей максимальная оценка класса была ниже 0,5. На таких участках модель слабо выделяет один класс среди остальных. Карту этих оценок удобно смотреть вместе с сегментацией: она помогает выбрать спорные места для проверки по исходному снимку.

Изменения между датами: FC-Siam-diff

FC-Siam-diff — полносвёрточная сиамская сеть для обнаружения изменений. Две ветви с общими весами извлекают признаки из снимков разных дат. Разности этих признаков передаются в декодер, который строит попиксельную карту изменения. В отличие от сегментации отдельного снимка, модель сразу обрабатывает пару изображений.

Такой вариант подходит для основной задачи: выделить места, которые стоит сравнить вручную. Я нашёл реализацию семейства FC-Siamese в каталоге TorchGeo, а готовые веса всей модели — в DeepEnsembles4ChangeDetection. Checkpoint обучен на OSCD — наборе пар снимков Sentinel-2 с разметкой изменений.

Для запуска выбранных весов использовал исходный код FC-Siam-diff. Так удалось собрать детектор из готовой реализации и предобученной модели. TorchGeo помог на этапе поиска подходящей архитектуры.

Модель содержит 1 351 586 параметров и получает 13 каналов каждой даты в порядке:

B04, B03, B02, B08, B05, B06, B07, B8A, B11, B12, B01, B09, B10

Checkpoint загрузил со строгой проверкой соответствия всех слоёв. Для нормализации каждой даты использовал одно среднее и стандартное отклонение по всем каналам на пригодных пикселях. Рабочая область фрагмента — 256 × 256 пикселей, контекст — 64 пикселя.

В эксперименте использовал готовые веса OSCD. Размер участка, статистика нормализации и интерполяция каналов отличались от авторского примера. Для контроля подал один и тот же фрагмент на оба входа: во всём контрольном фрагменте оценки остались ниже порога изменения 0,5.

4. Отбор участков с изменениями

Чтобы сузить выборку до возможной новой застройки, задал два условия: FC-Siam-diff должна показать изменение на обеих парах дат, а Dynamic World — рост оценки класса застройки на обоих снимках 2026 года.

Сокращённая логика отбора:

# change_july, change_august — оценки FC-Siam-diff относительно 2025 года
# built_2025, built_july, built_august — оценки класса застройки Dynamic World
persistent_change = (change_july > 0.5) & (change_august > 0.5)

persistent_built_gain = (
    (built_2025 < 0.35)
    & (built_july > 0.55)
    & (built_august > 0.55)
    & ((built_july - built_2025) > 0.25)
    & ((built_august - built_2025) > 0.25)
)

candidates = persistent_change & persistent_built_gain & valid_all_dates

Пороги для этого эксперимента задал вручную. В полученной маске оставил связные области от девяти пикселей — 900 м², чтобы убрать мелкие пятна. Этот порог относится к площади выделенной области; заметность реального объекта на снимке зависит также от его формы и окружения.

Полученные 28 полигонов — участки для проверки. Внутри одного могут оказаться несколько зданий, дорога или часть строительной площадки. Для подсчёта новых зданий и определения их состояния потребуется более детальная съёмка и сведения о самих объектах.

Для визуальной проверки выбрал шесть крупнейших кандидатов и сопоставил их на трёх снимках Sentinel-2. Во всех шести случаях заметны изменения искусственных поверхностей и площадок. Остальные 22 участка остаются в очереди на проверку.

Почему на снимках видно больше изменений, чем на карте

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

У такого правила есть три ограничения:

• Застройка на уже освоенной территории. Если на площадке появился ещё один объект, а в 2025 году модель уже относила её к застройке, условие built_2025 < 0.35 отбросит эти пиксели. Отбор ориентирован на переход к застройке с другого типа поверхности.

• Подтверждение на двух датах. Для каждого пикселя нужны и оценка изменения FC-Siam-diff, и рост оценки застройки Dynamic World — причём на обоих снимках 2026 года. Достаточно одной оценки ниже порога, чтобы пиксель выпал из результата. Поэтому контур может покрывать только часть площадки.

• Мелкие и разрозненные изменения. После пересечения условий остаются отдельные пятна. Все связные области меньше 900 м² удаляются, даже если вместе они относятся к одной большой стройке.

По отдельности методы дали 144 области изменений от FC-Siam-diff и 339 областей роста оценки застройки от Dynamic World. Эти промежуточные результаты не опубликованы. На карте показано пересечение двух условий — 28 участков для первоочередной проверки.

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

5. Векторизация

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

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

Для этой операции использовал rasterio.features.sieve с порогом девять пикселей и четырёхсвязной окрестностью. Соседними считаются пиксели с общей стороной. Фильтрация работает с областями растра, ещё до построения векторных границ. В бинарной маске изменений маленькие области просто исключались.

Полигоны построил через rasterio.features.shapes:

clean = sieve(labels, size=9, mask=valid, connectivity=4)

for geometry, class_id in shapes(
    clean, mask=valid, transform=transform, connectivity=4
):
    polygon = shape(geometry)          # геометрия в UTM
    area_ha = polygon.area / 10_000    # площадь в метрической системе
    # Затем геометрия преобразуется в WGS84 для GeoJSON.

Площади вычислил в UTM, где координаты заданы в метрах. Геометрию для GeoJSON преобразовал в WGS84, а рассчитанную площадь сохранил в атрибутах.

Для снимка 6 августа фильтрация изменила 0,1691% пикселей. После векторизации получилось 2 709 полигонов однородных классов. Обе версии растра — до фильтрации и после — сохранены.

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

6. Загрузка данных и публикация карты

После расчёта у меня были GeoTIFF со снимками и GeoJSON с полигонами. Для обсуждения с коллегами я хотел собрать их в одной карте: открыть исходное изображение, наложить результат модели и проверить конкретный участок. Для этого подойдёт любой картографический сервис, который поддерживает загрузку GeoTIFF и GeoJSON и публикацию слоёв.

Я использовал ORBISmap: здесь можно бесплатно опубликовать свои данные и продолжить работу с ними — добавлять слои, обновлять объекты и менять оформление карты. Далее покажу загрузку через MCP и публикацию через API этого сервиса.

В карту вошли три RGB-снимка, 28 кандидатов изменений, сегментация небольшого участка на три даты и полная сегментация от 6 августа. Получилось восемь слоёв без учёта подложки. Карта использует EPSG:3857, растры сохраняют исходную геопривязку UTM. Для двух RGB-файлов подготовил копии с JPEG-сжатием качества 95, чтобы сократить объём загрузки. Эти копии нужны для просмотра; исходные данные для расчётов хранятся отдельно.

Как агент работает с ГИС через MCP

MCP — Model Context Protocol — позволяет ИИ-агенту вызывать предоставленные сервисом инструменты с заданными параметрами. Через MCP-сервер доступны проекты, карты, таблицы, растры и другие сущности ГИС. Для подключения используется адрес https://my.orbismap.ru/mcp и авторизация заголовком; ключ хранится в конфигурации клиента.

Сначала агент получает список проектов через projects-list и создаёт карту инструментом create-map. Затем для каждого файла вызывает upload-init, передаёт содержимое по полученному адресу и подтверждает завершение через upload-confirm.

GeoJSON импортируется в таблицу вызовом import-create. После окончания импорта можно проверить число объектов через objects-count. Для GeoTIFF создаётся растровый слой (create-raster), к которому добавляется файл (add-raster-object) с датой съёмки и порядком RGB-каналов. Названия, атрибуты и стили слоёв настраиваются через edit-table и edit-column.

Так агент переносит результаты обработки в ГИС вместе с датами, классами и площадями. Эти данные остаются доступны для последующего просмотра и обновления. Для новой даты можно повторить загрузку, сохранив принятую структуру слоёв и цвета классов.

Публикация через API

Для публикации использовал HTTP API сервиса. Запрос POST /api/v1/publications/{project}/publications создаёт публикацию: в нём указываются карта, доступные слои, начальный масштаб и настройки интерфейса.

Например, сокращённое тело запроса для одной карты и слоя сегментации выглядит так. Значения идентификаторов в примере относятся к этому эксперименту:

{
  "code": "new_moscow_changes_2025_2026",
  "name": "Новая Москва — изменения 2025–2026",
  "enabled": true,
  "maps": [{
    "map_id": 286173,
    "selected": true,
    "expanded": true,
    "containers": [{
      "container_id": 8,
      "entity": "table",
      "code": "landcover_summer_20260806",
      "selected": true,
      "params": {
        "col_name": "class_ru",
        "popup": true,
        "objects_list": true
      }
    }]
  }]
}

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

7. Цвета классов и информация об участке

Чтобы сегментацию можно было читать поверх снимка, каждому классу назначил свой цвет: воде — синий, древесной растительности — зелёный, застройке — красный. Полупрозрачная заливка сохраняет видимость исходного RGB.

Цвет выбирается по class_id. Этот атрибут должен передаваться в векторный тайл, поэтому в стиле явно указан список columns. Поле class_ru используется для подсказки:

{
  "columns": ["id", "class_id", "class_ru"],
  "tooltip": "class_ru",
  "layers": [{
    "type": "fill",
    "paint": {
      "fill-opacity": 0.65,
      "fill-color": [
        "match", ["get", "class_id"],
        0, "#419bdf", 1, "#397d49", 2, "#88b053",
        3, "#7a87c6", 4, "#e49635", 5, "#dfc35a",
        6, "#c4281b", 7, "#a59b8f", 8, "#b39fe1",
        "#999999"
      ]
    }
  }]
}

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

Посмотреть снимки, изменения и сегментацию на карте.

8. Аналитика по результатам сегментации

По векторам сегментации можно получить сводку по всей территории. Для снимка 6 августа я локально посчитал площадь и число полигонов каждого класса после фильтрации:

Это статистика классификации. Её ошибки переходят и в расчёт площадей: например, модель выделила 0,19 га снега и льда на летнем снимке. Эти два полигона требуют проверки.

После загрузки в ГИС векторы хранятся в PostgreSQL с расширением PostGIS. Сводки, выборки и пространственные расчёты можно выполнять прямо внутри сервиса, создавая на их основе новые таблицы и виртуальные слои. Например, площадь и число полигонов по классам рассчитываются следующим запросом; вместо <таблица_сегментации> нужно указать имя загруженной таблицы:

SELECT class_ru,
       COUNT(*) AS polygon_count,
       ROUND(SUM(area_ha)::numeric, 2) AS area_ha
FROM <таблица_сегментации>
GROUP BY class_ru
ORDER BY area_ha DESC;

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

Через MCP для этого используется create-table: параметр type=view создаёт виртуальный слой на основе SQL-запроса, а type=table — таблицу с результатами выборки. Виртуальный слой получает данные из исходных таблиц при обращении к нему; отдельная таблица сохраняет результат на момент расчёта. Новому слою можно назначить стиль и добавить его в публикацию.

Для продолжения эксперимента полезны три расчёта:

• Распределение классов по районам. Пересечь сегментацию с границами и пересчитать площади частей полигонов в метрической системе координат. Сумма исходных area_ha даст завышение для объектов, пересекающих границу.

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

• Очередь проверки изменений. Ранжировать кандидатов по площади и расположению, добавить статусы ручной проверки.

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

Результат

В результате я получил карту с тремя снимками, 28 участками для проверки и летней сегментацией из 2 709 полигонов. Снимки и результаты моделей собраны на одной карте: можно сравнить даты и посмотреть атрибуты изменений. Если появится новая сцена или уточнённая классификация, их можно добавить в ту же карту.

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