A new algorithm for coregistration of digital elevationmodels (ilem)

Cover Page

Cite item

Full Text

Abstract

This paper proposes a new algorithm that allows performing a high-precision fitting of multi-temporal digital elevation models, which do not have appropriate geographic reference, in order to calculate the difference in elevation over a known time interval. Similar algorithms exist, the proposed algorithm is based on different principles, and therefore it can complement the toolkit for spatial data coregistration. The paper describes the stages of the algorithm operation, which in generalized form includes first the adjustment of the registered model to the reference model in plan, then – in vertical direction. The algorithm was tested on 2 sites and different kinds of data: 1) the 2014 landslide site in the valley of the Geysernaya River in Kamchatka using space imagery and stereo photogrammetry (ArcticDEM), and 2) an erosion monitoring site in the Gitche-Gizhgit catchment in the Greater Caucasus using aerial photography and a structure-from-motion approach (UAV). The proposed algorithm is effectively applicable to data of different origin, detail, spatial coverage. Conditions for its effective application: 1) presence of any significant areas with unchanged relief, 2) presence of a pronounced pattern of topographic dissection (texture of image / digital elevation model). It is shown that the refinement of the geographical reference of the registered elevation model significantly improves estimates of the volumes of denuded and accumulated material, which is especially important in the tasks of dynamic geomorphology. In the given examples, the registration error of digital elevation models decreased from 3–4 to almost 70 times. The volumes of surface changes in the areas of reliably prevailing denudation were corrected both in magnitude (as a rule, downward) and in sign.

Full Text

1. ВВЕДЕНИЕ

Сравнение разновременных цифровых моделей высот (ЦМВ), понимая под ними как цифровые модели рельефа твердой земной поверхности (ЦМР), так и модели высот интегральной поверхности грунта, строений, крон деревьев [цифровая модель поверхности или ЦМП, по (Чибуничев, 2022)] дает информацию об изменениях исследуемой поверхности за фиксированный промежуток времени. Эта информация не позволяет напрямую судить о динамике рельефа этой поверхности, т.е.говорить о механизмах рельефообразования, но, как минимум, позволяет анализировать кинематику рельефа (Девдариани, 1950), т.е. объемы и площади произошедших изменений.

Получение разновременных ЦМВ сейчас сопряжено, чаще всего, с фотограмметрической обработкой стереофотоизображений (Shean, 2016), в том числе с использованием подхода “структура из движения” (Structure from Motion (Westoby, 2012)), радарной интерферометрией (Crosetto, 2002). Географическую привязку получаемых материалов обеспечивает проведение наземных геодезических измерений. Вместе с тем во многих случаях получение наземных опорных точек в результате полевых работ невозможно или нецелесообразно. Помимо этого, даже абсолютная географическая привязка не всегда обеспечивает требуемую точность (субдециметровую/субсантиметровую) регистрации в пространстве всей сцены съемки, которая бы позволила отслеживать произошедшие изменения рельефа.

Существует несколько наиболее популярных алгоритмов корегистрации ЦМВ, в основе которых лежат различные подходы (в частности, глобального или локального совмещения ЦМВ) и метрики оценки качества совмещения (Van Niel, 2008). Глобальная корегистрация предполагает, что к регистрируемой модели высот (ко всей ее площади разом) применяется одно и то же математическое преобразование таким образом, чтобы совместить ее наилучшим образом с референсной моделью высот. Наилучшее совмещение такое, чтобы произвольно взятая точка (X, Y) на обеих моделях соответствовала строго одной точке поверхности на местности. Применение глобального преобразования требуется в случаях, когда такое совмещение удается получить одним для всей площади сдвигом по трем осям, углами поворота (в горизонтальной плоскости или в трехмерном пространстве), коэффициентом масштабирования регистрируемой модели, в общем – аффинными преобразованиями. В случае, если ошибка совмещения двух ЦМВ меняется по площади, при этом не имеет строгого тренда (т. е. нельзя однозначно вычислить ошибку совмещения как функцию от X и Y), то требуется локальное преобразование. Частным случаем локального преобразования является регистрация по наземным опорным точкам, при котором почти всегда происходит не только сдвиг, но и некоторое гладкое искажение регистрируемой ЦМВ.

Опыт автора в вопросе точного совмещения ЦМВ, полученных как по данным космической съемки, так и по результатам использования бюджетных БПЛА, показывает, что ошибка и в плане, и по высоте часто именно плавающая (Харченко, 2023). В результате ставшие уже стандартными алгоритмы [Nuth&Kaab (Nuth, 2011), Ames Stereo Pipeline (Beyer, 2018)] совмещения моделей высот в ряде сложных случаев не могут дать удовлетворительные результаты. Алгоритмы локального совмещения [например, CODEM (NCALM-UH/CODEM, 2023)] могут давать впечатляющие результаты, в то же время в нашей работе было несколько примеров, когда и этот передовой метод не справлялся с задачей. Большая группа способов корегистрации основана на использовании разных модификаций метода ICP (iterative closest point) (Besl, 1992). Его применяют как, например, в CloudCompare (Girardeau-Montaut, 2016), так и на этапе точного совмещения в упомянутом выше CODEM и т.д. Хорошие результаты данный метод дает именно в том случае, если грубое совмещение уже проведено (т.е. для растров – на этапе субпиксельного совмещения), но никак не в противном случае. Отдельную группу методов корегистрации ЦМВ составляют алгоритмы, работающие на основе SIFT (scale-invariant feature transform). Чаще эти алгоритмы используются для поиска общих точек на фотоизображениях, в том числе в составе рабочего процесса Structure from Motion для построения разреженного облака точек. Однако применяются они и для корегистрации ЦМВ (Aguilar, 2012; Sedaghat, 2018) в случае, если текстура этих одноканальных изображений дает достаточное количество пар общих точек.

В этой статье предлагается алгоритм локальной корегистрации, реализованный в качестве скрипта рабочего процесса и набора отдельных функций на языке R и названный ILEM (Iterative Local fitting of Elevation Models). Как и перечисленные выше алгоритмы, ILEM не совершенен, его ограничения будут описаны ниже. Однако он может быть альтернативой в случаях, если иные использованные методы не привели к получению желаемых результатов. Его разработка стала реакцией именно на невозможность получить удовлетворительные результаты имеющимся инструментарием.

2. МЕТОДИКА

Метод основан на итеративном локальном совмещении (независимо друг от друга) разных частей регистрируемой и референсной ЦМВ при условии минимизации метрики ошибки. Вся площадь сцены разбивается гексагональной решеткой с выбранным шагом. Шаг подбирается таким образом, чтобы, по возможности, каждый полигон охватывал предположительно стабильные площадки, наряду с динамичными. Например, если регистрируется ЦМВ участка многорукавного русла реки с островами и осередками, то размер полигонов должен позволять охватывать как участки отдельных проток, так и соседние поверхности островов, а не только лишь протоки, подмывающие берега и меняющие рельеф. Строго говоря, шаг сетки можно сделать сколь угодно малым (и в ряде случае это позволяет иметь большую по размеру выборку точек привязки), но если задать коэффициент размера полигонов > 1 (значение по умолчанию), тогда гексагоны будут в поперечнике больше шага сетки (рис. 1).

 

Рис. 1. Варианты генерации гексагональной сетки для поиска устойчивых площадок на участке многорукавного русла р. Гейзерной. Во всех случаях шаг сетки D = 15 м. Размеры полигонов: (а) – 1D, (б) – 0.7D, (в) – 1.5D, (г) – 2D.

Fig. 1. Variants of hexagonal grid generation to search for stable sites in the braided section of the Geysernaya River. In all cases the grid step D = 15 m. Polygon size: (а) – 1D, (б) – 0.7D, (в) – 1.5D, (г) – size 2D.

 

В качестве метрики ошибки используется среднеквадратическое отклонение разности высот по каждому из гексагонов при различных сдвигах по X и Y. При том сдвиге, при котором метрика достигает минимума, как правило, достигается и оптимальное совмещение. При этом сама разность высот не важна: если рельеф на двух ЦМВ оказывается морфологически подобен (относительные высоты и их пространственный рисунок эквивалентны), то при любой абсолютной разности высот среднеквадратическое отклонение окажется равным нулю. Иными словами, это позволяет находить позиции с одинаковой морфологией рельефа при любой ошибке абсолютной высоты. Использование подобной метрики дает возможность понять, есть ли в границах конкретной площадки (гексагона) сколько-нибудь значительная по площади стабильная поверхность, т.е. метод локально разыскивает не наземные опорные точки, а условно неизменные площадки, в чем состоит одно из его преимуществ.

Сначала для каждого гексагона получаем оптимальные сдвиги с грубым шагом, а вместе с ними – саму достигнутую (минимальную из полученных при разных сдвигах) величину метрики ошибки. На участках с активными изменениями рельефа, с растительностью, с возникшими между съемками антропогенными объектами эта метрика даже при оптимальном сдвиге окажется велика и позволит отсепарировать данные участки, которые невозможно корректно привязать по текстуре поверхности, так как внутри соответствующих гексагонов не найдено условно стабильных площадок. Затем итеративно для оставшихся после сепарации гексагонов оптимальный сдвиг уточняется до субпиксельного уровня.

Если поверхность имеет ярко выраженную текстуру, хорошо различимые осложняющие отрицательные и положительные формы рельефа, то в результате наложения исходных центров гексагонов и их положения после оптимального сдвига проявится тенденция, в чем именно состоит искажение регистрируемой ЦМВ относительной референсной. Проявится, есть ли сдвиг, одинаков по величине и однонаправлен ли он по площади, наблюдается ли вращение. Одновременно с этим проявятся те гексагоны, которые явно выбиваются из общих закономерностей по перечисленным выше параметрам – их необходимо удалить, а поворот, в случае его наличия, устранить. Затем по точкам центров оставшихся гексагонов регистрируемая модель совмещается с референсной в плане. Эту же процедуру по наборам точек, полученным при анализе ЦМВ, можно применить к ортофотоплану и совместить его с референсным.

Так как после сепарации в сетке гексагонов остаются только площадки с низкой величиной метрики ошибки, эти участки можно принять за условно стабильные. Следовательно, если внутри каждого гексагона ни денудация, ни аккумуляция не занимает больше половины площади, то медианная разность высот должна стремиться к нулю. Среднеквадратическое отклонение медианных разностей высот двух ЦМВ – есть ошибка высот после совмещения этих моделей в плане. Чаще всего эта ошибка имеет выраженный пространственный тренд: либо нарастает/убывает в каком-либо направлении на плоскости, либо меняется по высоте, например меньше на междуречьях и больше в днищах долин. Интерполяция тренда этой ошибки и последующее вычитание тренда из регистрируемой ЦМВ после ее совмещения с референсной ЦМВ в плане, дает совмещенную в плане и по высоте зарегистрированную ЦМВ. Повторный расчет среднеквадратического отклонения медианных разностей высот для условно стабильных площадок дает остаточную ошибку, которая, разумеется, никогда не оказывается равной нулю.

Остаточная ошибка случайна и чаще всего распределена по нормальному статистическому закону со средним значением около 0 м. Вычитанием референсной и зарегистрированной ЦМВ получаем растр разности высот, так называемый DEM of Differences, или DoD. Пиксели DoD со значениями в интервале ±2 остаточные ошибки – заполняются отсутствующими данными (No Data), таким образом устраняются данные с вероятностью 95% (т.е. в интервале ±2σ), являющиеся результатом ошибки привязки. Все остальные значения разности высот с той же вероятностью (или, напротив, вероятностью ошибки 5%) являются содержательным, полезным сигналом. Устраняемый из DoD слой в литературе называется слоем неопределенности (uncertainty layer), а половина его ширины – порогом обнаружения (level of detection) произошедших изменений высот (Bishop, 2006).

Помимо описанного, алгоритм включает опциональные функции, одной из важнейших из которых является быстрая сепарация территории по визуальным признакам на покрытую растительностью (под ней часто протекают геоморфологические процессы, но они, как правило, медленны, а их следы не регистрируются на ЦМВ, полученных с воздуха) и оголенную, создание маски растительности и вырезание из итоговой DoD или из исходных ЦМВ участков, заведомо непригодных. Общая схема алгоритма показана на рис. 2.

 

Рис. 2. Процедура корегистрации ЦМВ с использованием алгоритма ILEM. Пронумерованы обязательные пункты. Сокращения: ОФП – ортофотоплан, реф. – референсный (ая), рег. – регистрируемый (ая), стат. – статистическое.

Fig. 2. DEM co-registration workflow using the ILEM algorithm. Mandatory items are numbered. Abbreviations: ОФПorthophoto, реф. – reference, рег. – registered, стат. – statistical.

 

Алгоритм ILEM – полуавтоматический, требует от пользователя наблюдения за промежуточными результатами (в первую очередь, статистическим распределением полученных значений оптимального сдвига) и принятия решений по поводу дальнейших этапов обработки. Пользовательские параметры алгоритма (в порядке назначения):

1) имена файлов с ЦМВ;
2) шаг гексагональной сетки в метрах;
3) коэффициент для установления размера гексагонов относительно шага сетки (см. рис. 1);
4) ряд значений потенциального сдвига по X и Y в метрах (от-до-шаг для итераций от 1 до N);
5) после первой итерации и ознакомления с гистограммами значений оптимального сдвига и величины остаточной ошибки – пороговые значения сдвига и предельно допустимой ошибки, определяемые визуально по графикам, в метрах.

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

В остальном скрипт с алгоритмом может использоваться без редактирования.

Несмотря на то, что в качестве вариантов движения используется лишь сдвиг по двум осям, ILEM позволяет корректировать и угол поворота в горизонтальной плоскости (а по двум другим осям наклон корректируется на этапе оценки тренда ошибки высоты). При наличии небольшого угла поворота (до первых десятков градусов) регистрируемой ЦМВ относительно референсной, данный факт проявится при визуализации двух наборов точек – центроидов исходных гексагонов и сдвинутых при поиске оптимального совмещения моделей высот. Затем при необходимости можно вычислить и центр вращения, и медианное значение угла поворота, внести коррекцию, сохранить файл с устраненной ошибкой поворота и провести процедуру корегистрации заново.

Более подробно этапы работы алгоритма и ряд важных технических деталей, которые нужно понимать для его эффективного использования, объясняются в видеоинструкции (https://github.com/sergeikharchenko/ILEM).

3. АПРОБАЦИЯ НА ТЕСТОВЫХ УЧАСТКАХ

Методика в равной мере применима как к ЦМР, так и к ЦМП, в любых масштабах рассмотрения. В настоящей работе приведем два примера корегистрации разновременных моделей высот – 1) ArcticDEM для Камчатки; 2) собственные сверхдетальные ЦМВ для Северного Кавказа.

3.1. Оползень 2014 г. на левом борту долины р. Гейзерной (Камчатка)

В начале января 2014 г. в прибровочной части левого борта долины р. Гейзерной, на переходе к лавовому плато у подножия сопки Желтой активизировался блок, отделенный от остальной поверхности ранее сформировавшимися трещинами отседания (регистрировались уже в 70-х гг. XX в. (Леонов, 2014)). Возник обвал, под воздействием которого нижележащий склон также пришел в движение, и процесс из обвального трансформировался в оползневой. Примечательно, что на этом же борту, буквально в нескольких сотнях метров выше по течению реки, в октябре 1981 г. также сходил крупный оползень (Лебедева, 2020). По данным сотрудников Кроноцкого заповедника, в ведении которого находится и территория Долины Гейзеров, и сотрудников центра “СканЭкс”, общий объем мобилизованного вещества составил 2.5 млн м3, из которых большая часть сформировала плотину, приведшую к возникновению подпрудного озера (рис. 3, (а)), а 0.75 млн м3 прошло вниз по долине в виде селя; впоследствии часть материала была вынесена речной эрозией.

 

Рис. 3. Участок долины р. Гейзерной (Камчатка). (а) – фотоплан (ESRI Imagery); разности высот за 2012–2020 гг.: (б) – до корегистрации и (в) – после нее, (г) и (д) – то же за 2012–2022 гг.

Fig. 3. The Geysernaya river valley (Kamchatka). (а) – photoplan of the site (ESRI Imagery), (б)–(в) – elevation differences for 2012–2020 before and after co-registration, (г)–(д) – the same for 2012–2022.

 

Сотрудники Кроноцкого заповедника получили эту оценку фотограмметрическим методом путем обработки космических стереоизображений для состояния территории “до оползня” и “после оползня”, осуществив соответствующую привязку изображений. Доступные сейчас разновременные (с декабря 2007 г. по октябрь 2022 г.) данные ArcticDEM на всю территорию Арктики, также созданные из стереофотоизображений WorldView – 1, 2, 3, позволяют за последние 10–15 лет отслеживать значительные преобразования рельефа поверхности на огромной площади. Однако отдельные полосы (strips) ArcticDEM не имеют точной абсолютной и относительной привязки. Ошибка в плане часто достигает нескольких метров (при разрешении моделей 2 м). На рассматриваемую территорию оползня 2014 г. и ближайших окрестностей имеются качественные (в частности, без облаков на исходных снимках, из которых получены ЦМВ, а также за бесснежный период) изображения от сентября 2012 г., сентября 2020 г. и июня 2022 г. При сравнении ЦМВ до оползня (2012 г.) и после оползня (2020 г.), установлено, что характерная величина планового сдвига двух моделей 4–5 м, что в условиях значительной крутизны поверхности на рассматриваемой территории при использовании исходной привязки приводит к получению сильно искаженных результатов оценки смещений поверхности. Так, на склоне крутизной 30° при условии неизменности его морфологии на два временные среза, но относительном плановом сдвиге двух ЦМВ на 4 м, будет получена величина слоя денудации или аккумуляции (в зависимости от направления планового сдвига) по вертикали, равная 2.3 м. Таким образом, в горных регионах точная корегистрация ЦМВ для оценки изменений рельефа особенно важна. Для второй сравниваемой пары ЦМВ (2012 и 2022 гг.) ошибка планового сдвига еще больше – порядка 6–8 м.

При корегистрации обеих пар ЦМВ среднеквадратическое отклонение высот условно стабильных площадок составило около 4 см, что является весьма успешным результатом при совмещении данных с разрешением 2 м, полученных по космическим снимкам. Однако столь успешный результат обусловлен прежде всего обилием условно стабильных площадок на междуречных пространствах, и не отвечает реальной ошибке привязки динамичного днища и склонов долины. Поэтому при рассмотрении небольшой площади участка оползня 2014 г. ошибка корегистрации была оценена нами повторно. Она составила для исходных данных 0.79–1.87 м, а для итоговых данных (ЦМВ 2020 г. была совмещена с ЦМВ 2012 г.) 0.22–0.32 м (табл. 1). Таким образом, ширина слоя неопределенности – 0.44–0.64 м в двух случаях, корегистрация снизила ее в 3.6–5.9 раза.

До регистрации моделей медианное значение разности высот составляло 2.6–3.3 м, что означает, что больше чем на половине площади произошло увеличение относительных высот как минимум на указанную величину. Так как осуществляется относительная привязка одной модели к другой по условно стабильным площадкам, то этот результат не может быть следствием, например, вертикальных тектонических движений (даже если принять, что подъем земной коры на такие величины за несколько лет без разломной тектоники возможен). При медленных вертикальных движениях на всей площади участка условно стабильные площадки движутся примерно на ту же величину и в том же направлении, что и остальная часть площади. Подобные движения не могут быть оценены при взаимной привязке ЦМВ. Таким образом, при относительной привязке (корегистрации) ЦМВ такие сдвиги могут быть следствием только экзогенных процессов (денудации поверхности и сопряженной аккумуляции наносов на более низких уровнях). При этом рассматриваемая территория с позиций морфодинамики – скорее денудационная зона, с вышележащего плато (из-за границ сцены) на склоны вряд ли поступает много материала, с учетом незначительных уклонов и бронирования поверхности лавами. Общее повышение поверхности может происходить только за счет заполнения днища долины, но оно не занимает больше половины площади участка, а следовательно, на медианное значение изменения высот никак не могут приходиться величины выше 0 м (вернее, выше уровня обнаружения).

 

Таблица 1. Изменения высот на участке оползня 2014 г. в долине р. Гейзерной за два временных среза

Table 1. Elevation changes at the site of the 2014 landslide in the valley of the Geysernaya River for two time slices

Параметры

09.29.2012–09.17.2020*

09.29.2012–07.14.2022*

Ошибка корегистрации, м

1.87 / 0.32

0.79 / 0.22

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

Медиана, м

2.63 / 0.37

3.34 / –0.12

Максимум плотности распределения, м

1.75 / 0.42

2.95 / –0.26

Нижний квартиль (25%), м

0.73 / –1.08

0.71 / –1.17

Верхний квартиль (75%), м

5.67 / 2.35

5.98 / 1.89

Минимум, м

–102.8 / –105.1

–103.4 / –106.9

Максимум, м

38.7 / 36.1

38.6 / 36.9

Объемы отрицательных и положительных изменений рельефа

Объем вынесенного материала, млн м3

–2.52 / –3.24

–2.76 / –3.30

Объем отложенного материала, млн м3

3.50 / 2.69

4.51 / 2.66

Примечание: * – в числителе – значения до корегистрации ЦМВ, в знаменателе – после нее.

 

При условии, если на участке есть хоть какой-то ареал с неизменными отметками высот, а остальная территория (пусть даже больше половины площади) – динамична, модальное значение разности высот будет приходиться на 0 м и соответствовать данному стабильному участку. При неточной корегистрации это значение может быть любым, а при точной оно должно попадать в слой неопределенности (т.е. находиться в пределах погрешности совмещения ЦМВ). В табл. 1 показано, что в обоих случаях проведенная корегистрация позволяет достичь требуемого результата. Аналогичные выводы позволяют сделать значения нижнего и верхнего квартилей распределения изменений высот за 2012–2020 (2022) гг. В обоих случаях (до корегистрации) они говорят о том, что здесь происходила крайне мощная аккумуляция, что неверно. Минимальные и максимальные разности высот, на фоне их больших величин, почти не отличаются во всех четырех случаях (по 2 сравнения для каждого из двух временных срезов). В центральной части обвалившегося блока на склоне понижение высот составило более 100 м, а максимальные значения мощности оползневой плотины находятся в интервале 36–39 м, хотя камчатскими коллегами были оценены лишь в 22 м (Леонов, 2014).

 

Таблица 2. Изменения высот на участке развития плоскостной и ручейковой эрозии на водосборе пруда Гитче-Гижгит за 2020–2022 гг.

Table 2. Elevation changes at the site of shield and ephemeral gullies erosion in the Gitche-Gijgit Pond watershed for 2020–2022

Параметры

08.18.2020–07.15.2022*

Ошибка корегистрации, мм

367 / 5.3

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

Медиана, мм

288 / –18

Максимум плотности

распределения, мм

48 / 2.2

Нижний квартиль (25%), мм

5.4 / –45

Верхний квартиль (75%), мм

740 / 31

Минимум, м

–2.0 / –1.59

Максимум, м

3.08 / 1.26

Объемы отрицательных и положительных изменений рельефа

Объем вынесенного материала, м3

–140.4 / –53.2

Объем отложенного материала, млн м3

970.5 / 27.8

Примечание. * – в числителе – значения до корегистрации ЦМВ,
в знаменателе – после нее.

 

И по экстремальным значениям изменений высот, и по рисунку с пространственными распределениями величины изменений высот до корегистрации и после нее может возникнуть впечатление, что более точная привязка не меняет результат принципиально. Действительно, максимальная мощность обвала, предшествовавшего оползню, во всех случаях отличается на первые проценты. Однако при проведении геоморфологической интерпретации результатов очень важны могут быть участки с невысокими значениями слоя аккумуляции или денудации (дециметры и первые метры), занимающие, как правило, гораздо большие площади, нежели участки катастрофических изменений. В подтверждение этого говорят объемы вынесенного и переотложенного материала на участке оползня 2014 г. – до корегистрации на 1–1.75 млн м3 преобладали положительные изменения высот, а после точной привязки лидирует денудация, сальдо ее объемов и объемов аккумуляции – 0.64–0.72 млн м3. Во-первых, очевидно уменьшение разброса значений за 2 близких периода времени (а значит, и близких по результату совокупности рельефообразующих процессов). Во-вторых, полученные величины безвозвратного удаления наносов с участка весьма близки к оценкам, произведенным местными специалистами, в том числе по результатам полевых обследований – 0.75 млн м3 (Леонов, 2014). В-третьих, это соответствует общей морфодинамической схеме участка, когда крутые и весьма крутые склоны вероятнее являются денудационными позициями, и лишь днище долины (а также узкие днища ложбин на склоне по левому борту долины Гейзерной) – аккумулятивными. Снижение среднеквадратической ошибки взаимной регистрации двух ЦМВ в нашем случае привело к смене результата на прямо противоположный, а именно что весь участок в целом скорее денудационный. Средний слой денудации (то есть среднее направленное понижение поверхности за взятый период) на всей рассмотренной площади участка (рис. 1, (а)) в 2.53 км², включая условно стабильные задернованные и залесенные склоны, составил 25.3–28.4 см за 8–10 лет.

3.2. Участок развития линейной и плоскостной эрозии Гитче-Гижгит (Кабардино-Балкария)

Участок расположен на водосборе пруда Гитче-Гижгит, являющегося обособившейся частью пруда-отстойника Тырныаузского горно-металлургического комбината, на южном, аструктурном склоне Скалистого хребта Большого Кавказа. В целом склон водосбора ступенчатый, уступы ступеней представляют собой обвально-осыпные склоны, переходящие в коллювиально-пролювиальные шлейфы, подверженные медленным массовым смещениям наносов и эрозионным процессам. На одном из таких склонов в 2020 г. организован ежегодный (1–2 раза в год) мониторинг поверхности путем ее съемки с БПЛА и восстановления из данных съемки ЦМВ и ортофотоплана. При этом по ряду причин (включая, например, невозможность работать GNSS-оборудованием в RTK-режиме вне зоны устойчивого приема радиосвязи, а также удаленность участка от ближайшей базовой станции в г. Нальчик для применения поправок к координатам в постобработке) было решено использовать относительную, а не абсолютную привязку разновременных данных.

При проведении аэрофотосъемки бюджетными малыми БПЛА, широко распространившимися в последние годы, географическая привязка фотографий и, в конечном счете, ЦМП/ЦМР и фотоплана осуществляется с помощью встроенного навигационного модуля. Однако его погрешность в плане составляет несколько метров, а по высоте – до нескольких десятков метров. Это приводит к тому, что модели, построенные без наземного геодезического обоснования (а иногда – и с ним), не позволяют рассчитывать разность высот на интересующих участках простым вычитанием растровых изображений.

Общий характер рельефа: вогнутый в плане склон на высотах 1460–1498 м со средней крутизной 36±13°, осложнен одной крупной (глубиной до 1.5 м) и серией мелких промоин, сформированных ручейковой эрозией. Остальная часть площади оголенная, подвержена плоскостному смыву во время редких сильных ливней и снеготаяния. Участок пересекает скотопрогонная тропа с крутизной вплоть до горизонтальной, имеются и отвесные (до 89.2°) участки. Общая площадь – 2170 м2 (рис. 4). В силу того, что участок практически полностью свободен от растительности, имеем дело с ЦМР, а не с ЦМП.

Для апробации метода взяты результаты съемок участка от августа 2020 г. и июля 2022 г.; обе ЦМР с пространственным разрешением 0.6 см. Исходная регистрируемая модель рельефа (на 2022 г.) имеет плановый сдвиг ~1.16 м по азимуту 220°. В силу этого, склон к югу от наиболее крупной промоины оказывается мнимо приподнят за 2 года на величину до 3 м (зеленые оттенки на рис. 4, (в)), а склон непосредственно к северу мнимо опустился до 2 м. При отсутствии точной привязки оказывается, что аккумуляция на данном участке почти в 7 раз по объемам превышает денудацию, причем средний слой превышения аккумуляции над денудацией – 38 см за 2 года. Однако полученные значения лишь следствие взаимного планового сдвига двух ЦМВ.

 

Рис. 4. Участок Гитче-Гижгит (Кавказ). (а) – ортофотоплан участка на 2020 г., (б) – цифровая модель рельефа на 2020 г., (в) – разность высот за 2020–2022 гг. без корегистрации, (г) – разность высот за 2020–2023 гг. после корегистрации.

Fig. 4. The Gitche-Gizhgit site (Greater Caucasus). (а) – orthophotoplan of the site for 2020, (б) – digital elevation model for 2020, (в) – elevation difference for 2020–2022 without co-registration, (г) – elevation difference for 2020–2023 after co-registration.

 

Корегистрация для этого участка позволяет снизить ошибку привязки с 367 до 5.3 мм (табл. 2), т.е. почти в 70 раз, чему, безусловно, способствует четкий рисунок ручейкового расчленения, используемый ILEM для точного взаимного совмещения ЦМВ. Медианное значение разности высот составляет 18 мм, что выходит за рамки слоя неопределенности, но в данном случае легко объяснимо: почти весь этот участок денудационный, соответственно, центральное значение в ряду разностей высот будет отрицательным, так как очаги денудации покрывают больше половины площади участка. Ясно проявилась красным цветом (рис. 4, (г)) ручейковая сеть, максимальные величины вреза в крупнейшей из промоин составили 1.56 м за 2 года. Максимальные величины аккумуляции (1.23 м) характерны для ее же днища, они связаны со смещением одиночной крупной глыбы и последовавшей аккумуляцией на ней мелкозема. Прочие очаги аккумуляции опираются на дорогу и сконцентрированы у тылового шва площадки дороги, а ближе к бровке значения разности высот находятся в пределах погрешности, следовательно, эту часть дороги считаем стабильной (белый цвет вдоль дороги на рис. 4, (г)).

Как результат, объем денудированного материала на участке составляет 53.2 м3, а аккумулированного (как принесенного с вышележащего склона, так и переотложенного в пределах участка) – 27.8 м3. Для участка в целом это соответствует понижению поверхности на 1.2 см за 2 года сверх погрешности привязки.

4. ОГРАНИЧЕНИЯ МЕТОДА

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

1) в первую очередь важно подчеркнуть, что качество корегистрации ЦМВ напрямую зависит от качества исходных данных: наличие любых артефактов в силу ограничений родительских по отношению к ILEM методов получения ЦМВ непременно лимитирует максимально достижимое качество результата;
2) в отличие от ряда прямых аналогов, данный метод полуавтоматический, а не автоматический. На ключевых этапах процедура регулируется пользователем, что, в общем, можно воспринимать и как достоинство, так как позволяет контролировать качество на всем протяжении процедуры корегистрации;
3) для участков без выраженного рисунка расчленения метод может давать неудовлетворительные результаты, даже если имеются одиночные удобные для корегистрации точки, например, выступающие над плоской поверхностью вершины валунов, занесенных мелкодисперсным грунтом;
4) без внимательного визуального контроля качества итогового набора гексагонов (на предмет того, нет ли выбивающихся из общего ряда участков с принципиально иным характером искажений, сдвигом и т.п. относительно окружающих участков) возможны очень сильные ложные деформации регистрируемой ЦМВ;
5) итоговая область корегистрации ограничивается выпуклой оболочкой для центров гексагонов, в которых найдены условно стабильные площадки. Прочая площадь отсекается.

Тем не менее внимание к деталям при использовании алгоритма часто позволяет достичь впечатляющих результатов, превышающих по качеству полностью ручную привязку.

5. ЗАКЛЮЧЕНИЕ

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

В первом случае показано, что исходная точность привязки ArcticDEM явно недостаточна для анализа временных рядов ЦМВ, вплоть до того, что результатом такого анализа на явно денудационных площадках может стать мнимое преобладание процессов аккумуляции. Проанализирован участок левого борта долины Гейзерной, на котором в 2014 г. возник обвал скальной стенки, перешедший в оползень ниже по склону, который, в свою очередь, сформировал в днище долины плотину и подпрудное озеро. Показано, как снижение ошибки корегистрации ЦМВ с 0.8 до 0.22 м меняет итоговую величину аккумуляции материала (в основном в теле этой плотины) с 4.5 до 2.6 млн м3, т.е. почти в 2 раза. Результаты оценки объемов аккумуляции и денудации после корегистрации за 2 временных интервала, каждый из которых покрывает 2014 г., близки друг с другом и с данными других исследователей. В то же время результаты оценок, основанные на исходных полосах ArcticDEM, визуально вполне неплохо совмещенных друг с другом, не соответствуют ни географической логике, ни контрольным данным, ни даже друг другу для двух интервалов времени. Процедура корегистрации позволила уменьшить ошибку взаимной привязки на первом тестовом участке в 3.6–5.9 раза.

Для участка Гитче-Гижгит в горах Большого Кавказа корегистрации подверглись разновременные ЦМР, полученные по данным низковысотной (20–35 м) съемки с БПЛА. Цель подобных съемок – мониторинг медленных геоморфологических процессов, в частности, плоскостного смыва и ручейковой эрозии со скоростями по вертикали от первых сантиметров в год. При попытке сопоставить разновременные ЦМР без предварительного совмещения отлично проявляется проблема влияния их планового совмещения на итоговую разность высот. Одна из фасеток склона оказывается мнимо приподнята на 3 м, другая опущена на 2 м, что явно не отвечает реальным изменениям рельефа. Корегистрация позволила понизить среднеквадратическую ошибку привязки с 367 до 5.3 мм. На итоговой модели DoD хорошо проявляется пространственный рисунок эрозионно-аккумулятивного процесса, становятся четко различимы все промоины, сформированные ручейковой эрозией, а также очаги временной аккумуляции материала. Объемы денудации превышают объемы аккумуляции примерно в 2 раза.

БЛАГОДАРНОСТИ

Работа выполнена за счет гранта РНФ, проект № 23-77-01027. Автор выражает благодарность за данные, используемые для тестирования алгоритма, В.Н. Голосову (Большой Кавказ) и Е.В. Лебедевой (Камчатка).

ACKNOWLEDGMENTS

The study was carried out under support of The Russian Science Foundation (RSF), project № 23-77-01027. The author would like to thank for the data used to test the algorithm V.N. Golosov (Greater Caucasus) and E.V. Lebedeva (Kamchatka) for the data used to test the algorithm.

×

About the authors

S. V. Kharchenko

Lomonosov Moscow State University; Institute of Geography, RAS

Author for correspondence.
Email: xar4enkkoff@yandex.ru

Faculty of Geography

Russian Federation, Moscow; Moscow

References

  1. Aguilar F.J., Aguilar M.A., Fernandez I., et al. (2012). A New Two-Step Robust Surface Matching Approach for Three-Dimensional Georeferencing of Historical Digital Elevation Models. IEEE Trans. Geosci. Electron. № 9. P. 589–593. https://doi.org/10.1109/LGRS.2011.2175899.2012
  2. Besl P.J., McKay N.D. (1992). Method for registration of 3-D shape. Sensor fusion IV: control paradigms and data structures. V. 1611. P. 586–606.https://doi.org/10.1109/34.121791
  3. Beyer R.A., Alexandrov O., McMichael S. (2018). The Ames Stereo Pipeline: NASA’s open source software for deriving and processing terrain data. Earth and Space Sci. V. 5. Iss. 9. P. 537–548.https://doi.org/10.1029/2018EA000409
  4. Bishop T.F., Minasny B., McBratney A.B. (2006). Uncertainty analysis for soil‐terrain models. Int. J. of Geographical Information Sci. V. 20. Iss. 2. P. 117–134.https://doi.org/10.1080/13658810500287073
  5. Chibunichev A.G. (2022). Fotogrammetrija (Photogrammetry). M.: MIIGAiK (Publ.). 328 p.
  6. Crosetto M. (2002). Calibration and validation of SAR intererometry for DEM generation. ISPRS J. of Photogrammetry and Remote Sensing. V. 57. Iss. 3. P. 213–227. https://doi.org/10.1016/S0924-2716(02)00107-7
  7. Devdariani A.S. (1950). Kinematics of terrain. In: Voprosy geografii. Sbornik 21. Geomorfologiya. Moscow: Geografgiz (Publ.). P. 55–85. (in Russ.)
  8. Girardeau-Montaut D. (2016). CloudCompare. France: EDF R&D Telecom ParisTech. V. 11. P. 5.
  9. Kharchenko S.V. (2023). The method for co-registration of digital terrain data to obtain hydrologically correct model of the earth’s surface. Geomorfologiya i Paleogeografiya. V. 54. № 3. P. 150–164. (in Russ.) https://doi.org/10.31857/S2949178923030039
  10. Lebedeva E.V., Sugrobov V.M., Chizhova V.P. et al. (2020). The Valley of the River Geyzernaya (Kamchatka): Hydrothermal Activity and Features of Relief Forming. Geomorfologiya. № 2. P. 60–73. (in Russ.) https://doi.org/10.31857/S0435428120020066
  11. Leonov V.L. (2014). Landslide and landslide that occurred on January 4, 2014 in The Valley of Geysers, Kamchatka, and their consequences. Vestnik KRAUNTS. Nauki o Zemle. V. 23. № 1. P. 7–20. (in Russ.)
  12. NCALM-UH/CODEM [Electronic data]. Access way: https://github.com/NCALM-UH/CODEM/tree/main (access date: 09.09.2023).
  13. Nuth C., Kääb A. (2011). Co-registration and bias corrections of satellite elevation data sets for quantifying glacier thickness change. The Cryosphere. V. 5. Iss. 1. P. 271–290. https://doi.org/10.5194/tc-5-271-2011
  14. Sedaghat A., Naeini A.A. (2018). DEM orientation based on local feature correspondence with global DEMs. GIS cience & Remote Sensing. № 55. P. 110–129. https://doi.org/10.1080/15481603.2017.1364879
  15. Shean D.E., Alexandrov O., Moratto Z.M. et al. (2016). An automated, open-source pipeline for mass production of digital elevation models (DEMs) from very-high-resolution commercial stereo satellite imagery. ISPRS J of Photogrammetry and Remote Sensing. № 116. P. 101–117. https://doi.org/10.1016/j.isprsjprs.2016.03.012
  16. Van Niel T.G., McVicar T.R., Li L. et al. (2008). The impact of misregistration on SRTM and DEM image differences. Remote Sensing of Environment. V. 112. Iss. 5. P. 2430–2442. https://doi.org/10.1016/j.rse.2007.11.003
  17. Westoby M.J., Brasington J., Glasser N.F. et al. (2012). ‘Structure-from-Motion’ photogrammetry: A low-cost, effective tool for geoscience application. Geomorphology. № 179. P. 300–314. https://doi.org/10.1016/j.geomorph.2012.08.021

Supplementary files

Supplementary Files
Action
1. JATS XML
2. Fig. 1. Variants of hexagonal grid generation to search for stable sites in the braided section of the Geysernaya River. In all cases the grid step D = 15 m. Polygon size: (а) – 1D, (б) – 0.7D, (в) – 1.5D, (г) – size 2D.

Download (1MB)
3. Fig. 2. DEM co-registration workflow using the ILEM algorithm. Mandatory items are numbered. Abbreviations: ОФП – orthophoto, реф. – reference, рег. – registered, стат. – statistical.

Download (1MB)
4. Fig. 3. The Geysernaya river valley (Kamchatka). (а) – photoplan of the site (ESRI Imagery), (б)–(в) – elevation differences for 2012–2020 before and after co-registration, (г)–(д) – the same for 2012–2022.

Download (1MB)
5. Fig. 4. The Gitche-Gizhgit site (Greater Caucasus). (а) – orthophotoplan of the site for 2020, (б) – digital elevation model for 2020, (в) – elevation difference for 2020–2022 without co-registration, (г) – elevation difference for 2020–2023 after co-registration.

Download (882KB)

Copyright (c) 2024 Russian Academy of Sciences

Согласие на обработку персональных данных с помощью сервиса «Яндекс.Метрика»

1. Я (далее – «Пользователь» или «Субъект персональных данных»), осуществляя использование сайта https://journals.rcsi.science/ (далее – «Сайт»), подтверждая свою полную дееспособность даю согласие на обработку персональных данных с использованием средств автоматизации Оператору - федеральному государственному бюджетному учреждению «Российский центр научной информации» (РЦНИ), далее – «Оператор», расположенному по адресу: 119991, г. Москва, Ленинский просп., д.32А, со следующими условиями.

2. Категории обрабатываемых данных: файлы «cookies» (куки-файлы). Файлы «cookie» – это небольшой текстовый файл, который веб-сервер может хранить в браузере Пользователя. Данные файлы веб-сервер загружает на устройство Пользователя при посещении им Сайта. При каждом следующем посещении Пользователем Сайта «cookie» файлы отправляются на Сайт Оператора. Данные файлы позволяют Сайту распознавать устройство Пользователя. Содержимое такого файла может как относиться, так и не относиться к персональным данным, в зависимости от того, содержит ли такой файл персональные данные или содержит обезличенные технические данные.

3. Цель обработки персональных данных: анализ пользовательской активности с помощью сервиса «Яндекс.Метрика».

4. Категории субъектов персональных данных: все Пользователи Сайта, которые дали согласие на обработку файлов «cookie».

5. Способы обработки: сбор, запись, систематизация, накопление, хранение, уточнение (обновление, изменение), извлечение, использование, передача (доступ, предоставление), блокирование, удаление, уничтожение персональных данных.

6. Срок обработки и хранения: до получения от Субъекта персональных данных требования о прекращении обработки/отзыва согласия.

7. Способ отзыва: заявление об отзыве в письменном виде путём его направления на адрес электронной почты Оператора: info@rcsi.science или путем письменного обращения по юридическому адресу: 119991, г. Москва, Ленинский просп., д.32А

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

9. Порядок уничтожения персональных данных при достижении цели их обработки или при наступлении иных законных оснований определяется Оператором в соответствии с законодательством Российской Федерации.

10. Я согласен/согласна квалифицировать в качестве своей простой электронной подписи под настоящим Согласием и под Политикой обработки персональных данных выполнение мною следующего действия на сайте: https://journals.rcsi.science/ нажатие мною на интерфейсе с текстом: «Сайт использует сервис «Яндекс.Метрика» (который использует файлы «cookie») на элемент с текстом «Принять и продолжить».