Основные геоморфометрические параметры: теория

Эта страница является черновиком статьи.

Расчет (Zevenbergen-Thorne, 1987) и интерпретация уклона, экспозиции, кривизн рельефа земной поверхности

Геморфометрический анализ растровых ЦМР базируется на двух исходных положениях. Первое основывается на математической формализации земной поверхности, а второе предусматривает расчет показателя в точке (пикселе) с учетом окружения.

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

$z=f(x,y)$

(1)

где $z$ – значение высоты в точке с географическими координатами $(x,y)$, которое для лучшей аппроксимации рельефа может быть выражено более сложными функциями, например – полиномиальными (или многочленами). В таком случае многочлен 2-го порядка, используемый для аппроксимации земной поверхности, может иметь следующий вид:

$z=Ax^{2}+By^{2}+Cxy+Dx+Ey+F$

(2)

где $x$ и $y$ географические координаты точки, высоту $z$ которой необходимо определить, $A...F$ – коэффициенты уравнения аппроксимирующей поверхности 2-го порядка. Многочлены являются одними из наиболее простых и хорошо изученных функций в математике. Они характеризуются такими свойствами как непрерывность и сглаженность, благодаря чему их легко можно интегрировать и дифференцировать. Это открывает возможности использования математического анализа не только для более совершенного представления земной поверхности, но и для изучения ее свойств, например, на основе производных разных порядков.

Согласно второму положению, основной аналитической операцией в ГИС, которая используется для расчета большинства параметров на основе растровых ЦМР является анализ окружения. Он позволяет количественно описать связь между точкой (пикселем) и его ближайшим окружением, применяя для расчета локальное (чаще всего, размером 3×3 пиксела) скользящее окно (рис. 1).

рисунок 1 Рис. 1 Расчет большинства геоморфометрических параметров как правило производится на основе скользящего окна размером 3×3 пиксела (Geomorphometry…, 2008)

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

В данной статье мы будем рассматривать особенности расчета основных геморфометрических параметров на примере алгоритма Zevenbergen-Thorne (Zevenbergen, Thorne 1987), который характеризуется расчетной эффективностью и высокой достоверностью результатов (Skidmore 1989, Jones 1998, Zhou, Liu 2004, Rodríguez, Suarez 2010). Кроме того, он реализован как в Открытых (SAGA), так и проприетарных ГИС (кривизны в ArcGIS, расширение для ArcGIS DEM Surface Tools от Jenness Enterprises).

Алгоритм Zevenbergen-Thorne использует модификацию (2) следующего вида:

$z=Ax^{2}y^{2}+Bx^{2}y+Cxy^{2}+Dx^{2}+Ey^{2}+Fxy+Gx+Hy+I$

(3)

где $A...I$ – коэффициенты аппроксимации, рассчитанные с помощью полиномов Лагранжа на основе 9 значений $z$ в ячейках окна 3×3. Геоморфометрические параметры получаются в результате дифференциациии (3) и решения соответствующих уравнений для центральной ячейки квадратной матрицы 3×3.

Основные геоморфометрические параметры, рассчитываемые на основе производных первого порядка

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

Уклон поверхности (Slope)

Понятие

Уклон поверхности – угол наклона в точке пересечения между горизонтальной плоскостью и плоскостью касательной к земной поверхности; фиксирует интенсивность перепада высот (градиент) между двумя заданными точками (рис. 2)

рисунок 2 Рис. 2 XXX

Если земная поверхность представлена функцией $z=f(x,y)$, то уклон рассчитывается с учетом изменений значений $z$ в двух направлениях как ${\frac {\partial z}{\partial xy}}$:

$\alpha =\tan ^{-1}\left({\sqrt {\left({\frac {\partial z}{\partial x}}\right)^{2}+\left({\frac {\partial z}{\partial y}}\right)^{2}}}\right)$

(4)

где ${\frac {\partial z}{\partial x}}$ и ${\frac {\partial z}{\partial y}}$ - производные первого порядка, представляющие изменение значений абсолютной высоты $z$ с запада на восток ($x$) и с севера на юг ($y$).

Расчет по Zevenbergen-Thorne

Процедура определения уклона поверхности по алгоритму Zevenbergen-Thorne сводится к следующим шагам (рис. 3):

рисунок 3 Рис. 3 Определение уклона поверхности (Zevenbergen, Thorne 1987) XXX

1. определение уклона поверхности в направлении с востока на запад:

$G={\frac {\left(z_{6}+z_{4}\right)}{2l}}$

(5)

2. определение уклона поверхности в направлении с севера на юг:

$G={\frac {\left(z_{2}+z_{8}\right)}{2l}}$

(6)

где $z_{2...8}$ - высотные отметки в соответствующих ячейках растра, а $l$ - расстояние между индивидуальными элементами матрицы высот, другими словами – пространственное разрешение растра. При этом предусматривается, что единицы измерения абсолютной высоты и пространственного разрешения идентичны (как правило, метры);

3. определение интегрального значения уклона поверхности для центральной ячейки скользящего окна:

$\alpha =\tan ^{-1}\left({\sqrt {G^{2}+H^{2}}}\right)$

(7)

Рассмотрим алгоритм на примере Рис. 4:

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

Пример расчета уклона поверхности

     
 Рассчитаем уклон по 5: 
$G={\frac {3-10}{2\times 10}}=-0.35$
Рассчитаем уклон по 6: 
$H={\frac {6-7}{2\times 10}}=-0.05$
Тогда уклон поверхности в центральной ячейке, рассчитанный по 7 составит:
$\alpha =\tan ^{-1}\left({\sqrt {\left(-0.35\right)^{2}+\left(-0.05\right)^{2}}}\right)=\tan ^{-1}\left({\sqrt {0.125}}\right)=\tan ^{-1}\left(0.3536\right)=19.47^{\circ }$

Интерпретация

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

• поверхностный сток и дренирование – чем более крутой склон, тем интенсивнее поверхностный сток и меньше инфильтрация влаги в почвенную толщу. Таким образом, уклон имеет принципиальное значение для режима увлажнения почв, особенно – верхних слоев;
• эрозия – интенсивность эрозии растет экспоненциально с увеличением уклона. Это объясняется тем, что с увеличением градиента кинетическая энергия осадков остается постоянной, но транспорт ускоряется в направлении подножья. В результате, кинетическая энергия стока превышает кинетическую энергию осадков, когда склон переходит отметку 8,5°, что и способствует проявлению эрозионных процессов;
• мощность почвенного профиля на склоне закономерно изменяется в соответствии с уклоном и относительной высотой. Как правило, почвенная толща меньше на возвышенных наклонных участках вследствие эрозионных процессов и гравитационного перемещения материала, и постепенно увеличивается в направлении пониженных участков с небольшим уклоном;
• количество солнечной энергии также зависит от уклона, поскольку он определяет угол падения солнечных лучей на земную поверхность. Увеличение уклона поверхности в направлении поступления солнечных лучей увеличивает угол их падения, а значит – количество энергии, которое получает поверхность. Это определяет микроклиматические особенности участка, в частности температуру, эвапотранспирацию и влажность верхних слоев почвы;
• особенности растительного покрова совокупно отражают все вышеперечисленные характеристики, поскольку они прямо или косвенно влияют на такие эдафические факторы как водный и температурный режим почвы, механический состав корнесодержащего слоя, содержание питательных элементов и т.д.

Простота расчета и информативность делают уклон поверхности наиболее употребимым показателем в моделировании процессов перераспределения поверхностного и внутрипочвенного стока, эрозии, определении эдафических условий, индикационном картографировании в физической географии и близких отраслях. Как правило, значения показателя измеряются в градусах (также это могут быть проценты или радианы) и колеблются в диапазоне от 0° (горизонтальная плоскость) до 90° (вертикальная плоскость). Часто используются следующие градации уклона поверхности:

Для равнинных территорий (Жучкова, Раковская, 2004)		Для горных территорий
		(Жучкова, Раковская, 2004)		(Міллер, 1996)
меньше 1°	плоские (субгоризонтальные) равнины	меньше 4°	плоские и почти плоские поверхности	меньше 3°	очень пологие склоны
1-3°	слабонаклонные равнины (очень пологие склоны)	4-10°	пологие склоны	3-6°	пологие склоны
3-5°	пологие склоны (наклонные равнины)	10-20°	покатые склоны	6-9°	слабопокатые склоны
5-7°	слабопокатые склоны	20-30°	склоны средней крутизны	9-12°	покатые склоны
7-10°	покатые склоны	30-45°	крутые склоны	12-15°	сильнопокатые склоны
10-15°	сильнопокатые склоны	45-60°	очень крутые склоны	15-30°	крутые склоны
15-20°	крутые склоны	больше 60°	скалистые (обрывистые) склоны	30-45°	очень крутые склоны
20-40°	очень крутые склоны			больше 45°	обрывистые склоны
больше 40°	обрывистые склоны

Экспозиция (Aspect)

Понятие

Экспозиция поверхности – угол по часовой стрелке между определенным направлением (как правило, на север) и проекцией уклона на горизонтальную плоскость; фиксирует направление (азимут) максимального уклона (градиента) земной поверхности (рис. 5).

рисунок 5 Рис. 5 ?Экспозиция с математической точки зрения?

Для земной поверхности представленной функцией $z=f(x,y)$ экспозиция рассчитывается как угол между двумя производными по формуле:

$\beta =-\tan ^{-1}\left({\frac {\partial z/\partial x}{\partial z/\partial y}}\right)$

(8)

где ${\frac {\partial z}{\partial x}}$ и ${\frac {\partial z}{\partial y}}$ - производные первого порядка, представляющие изменение значений абсолютной высоты $z$ с запада на восток ($x$) и с севера на юг ($y$).

Расчет по Zevenbergen-Thorne

Интегральное значение экспозиции поверхности в центральной ячейке скользящего окна по алгоритму Zevenbergen-Thorne определяется по формуле:

$\beta =-\tan ^{-1}\left({\frac {\pm H}{\pm G}}\right)$

(9)

В данном случае ± определяет в какой четверти находится β по отношению к выбранному направлению. Как правило, экспозиция отсчитывается по часовой стрелке от северного направления географического меридиана, потому (9) приобретает следующий вид:

$\beta =180^{\circ }-\tan ^{-1}\left({\frac {\pm H}{\pm G}}\right)+90^{\circ }\left({\frac {\pm G}{\left|G\right|}}\right)$

(10)

Пример расчета экспозиции поверхности

     
 Рассмотрим алгоритм на примере рис. 4. Уже известно, что по (5) G=-0.35, а по (6) H=-0.05. Тогда экспозиция в центральной ячейке, рассчитанная по (10) составит: 
$\beta =180^{\circ }-\tan ^{-1}\left({\frac {-0.05}{-0.35}}\right)+90^{\circ }\left({\frac {-0.35}{\left|0.35\right|}}\right)=180^{\circ }-\tan ^{-1}\left(0.1429\right)+90^{\circ }\left(-1\right)=180^{\circ }-8.1325^{\circ }-90^{\circ }=81.87^{\circ }$

Интерпретация

Функциональная интерпретация экспозиции может вестись в нескольких направлениях, поскольку она характеризует:

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

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

Отмеченные выше особенности распределения тепла дополнительно усложняются тем, что склоны восточной и западной ориентации формально имея одинаковую экспозицию (по 90°), также демонстрируют контрастность тепловых условий. Это объясняется тем, что на восточные склоны солнечные лучи попадают в утренние прохладные часы и тратятся на прогревание поверхности, а западне склоны освещаются во второй половине дня, когда поверхность уже прогрета. В результате, западне склоны несколько тепле и суше, чем восточные. М. Гродзинський (2013) отмечает, что при весьма значительной крутизне склонов (особенно в горных ландшафтах) склоны юго-западной экспозиции оказываются тепле и суше, чем южной экспозиции, а склоны северо-восточной экспозиции – холоднее и влажнее, чем северные. Еще раньше эту особенность отметил Р. Уиттекер (1980), который самым холодным считает не северное, а северо-восточное местоположение, а наиболее теплым – юго-восточное.

Типология инсоляционных местоположений проводится по сторонам горизонта (румбам)

северное	N	0-22.5°; 337.5-360°
северо-восточное	NE	22.5-67.5°
восточное	E	67.5-112.5°
юго-восточное	SE	112.5-157.5°
южное	S	157.5-202.5°
юго-западное	SW	202.5-247.5°
западное	W	247.5-292.5°
северо-западное	NW	292.5-337.5°

Определение соотношения по теплообеспеченности между этими типами может трактоваться по-разному с учетом того, как будут рассматриваться склоны восточной и западной экспозиций. Поэтому возможно несколько вариантов упорядочивания местоположений в инсоляционный ряд (табл. 1, рис. 4).

Варианты рядов тепло(влаго)обеспеченности по инсоляционной экспозиции (Гродзинський, 2013)

	Холодно => Тепло Влажно => Сухо
Ориентационный ряд	N → NE=NW → E=W → SE=SW → S
Традиционный "компасный" ряд	N → NE → NW → E → W → SE → SW → S
Ряд Уиттекера	NE → N → NW → E → W → SE → S → SW

рисунок 6 Рис. 6 Варианты типологии инсоляционных местоположений (Гродзинський, 2013): а – по ориентационному, б – по «компасному», в – по экологическому принципу (топографический ряд Уиттекера); 1-8 – типы местоположений (увеличение порядкового номера указывает на увеличение сухости местоположения). Стрелкой показано направление движения солнца в светлую пору дня.

Основные геоморфометрические параметры, рассчитываемые на основе производных второго порядка

Для рассмотренных выше производных функции $z=f(x,y)$ в свою очередь могут быть рассчитаны производные. И если первая производная описывала градиент поверхности (его величину и направление), то вторая фиксирует меру изменений этого градиента, т.е. является градиентом первой производной в заданном направлении. На производных второго порядка основывается система морфометрических кривизн, описывающих форму поверхности.

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

Горизонтальная (плановая) кривизна (Plane Curvature)

Понятие

Горизонтальная (плановая) кривизна – кривизна линии, образованной пересечением земной поверхности с плоскостью, перпендикулярной к направлению ориентации максимального градиента (экспозиции). Как производная второго порядка, горизонтальная кривизна описывает градиент экспозиции вдоль заданного контура (рис. 7).

рисунок 7 Рис. 7 Горизонтальная кривизна

Если земная поверхность представлена функцией $z=f(x,y)$, то плановая кривизна является функцией ее частичных производных:

$PLANC={\frac {\left({\frac {\partial z}{\partial y}}\right)^{2}{\frac {\partial ^{2}z}{\partial x^{2}}}-2{\frac {\partial z}{\partial x}}{\frac {\partial z}{\partial y}}{\frac {\partial ^{2}z}{\partial x\partial y}}+\left({\frac {\partial z}{\partial x}}\right)^{2}{\frac {\partial ^{2}z}{\partial y^{2}}}}{\sqrt {\left(\left({\frac {\partial z}{\partial x}}\right)^{2}+\left({\frac {\partial z}{\partial y}}\right)^{2}\right)^{3}}}}$

(11)

Расчет по Zevenbergen-Thorne (1987)

Для получения значения горизонтальной кривизны на основе алгоритма Zevenbergen-Thorne используют следующее выражение:

$PLANC={\frac {-2\left(DH^{2}+EG^{2}-FGH\right)}{\left(G^{2}+H^{2}\right)}}\times 100^{optionally}$

(12)

В данном случае G и Н рассчитываются по выражениям (5) и (6) соответственно, а другие частные производные как:

$D={\frac {\left[{\frac {\left(z_{4}+z_{6}\right)}{2}}-z_{5}\right]}{l^{2}}}$

(13)

$E={\frac {\left[{\frac {\left(z_{2}+z_{8}\right)}{2}}-z_{5}\right]}{l^{2}}}$

(14)

$F={\frac {\left(z_{3}-z_{1}+z_{7}-z_{9}\right)}{4l^{2}}}$

(15)

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

     
 Рассмотрим алгоритм на примере рис. 4. Уже известно, что по (5) G=-0.35, а по (6) H=-0.05. Рассчитаем другие частные производные:
по (13) $D={\frac {\left[{\frac {\left(10+3\right)}{2}}-4\right]}{10^{2}}}={\frac {2.5}{100}}=0.025$
по (14) $E={\frac {\left[{\frac {\left(6+7\right)}{2}}-4\right]}{10^{2}}}={\frac {2.5}{100}}=0.025$
по (15) $F={\frac {\left(9-4+8-1\right)}{4\times 10^{2}}}={\frac {12}{400}}=0.03$
Тогда горизонтальная кривизна для центральной ячейки, полученная по (12) составит:
Невозможно разобрать выражение (синтаксическая ошибка): {\displaystyle PLANC = \frac{-2 \left (0.025 \times \left (-0.05 \right )^2+0.025 \times \left (-0.35 \right )^2-0.03 \times left (-0.35 \right ) \times left (-0.05 \right ) \right )}{left (-0.35 \right )^2+left (-0.05 \right )^2}=}

Невозможно разобрать выражение (SVG (MathML можно включить с помощью плагина для браузера): Недопустимый ответ («Math extension cannot connect to Restbase.») от сервера «https://wikimedia.org/api/rest_v1/»:): {\displaystyle PLANC = \frac{-2 \left ( \right )}{left (-0.35 \right )^2+left (-0.05 \right )^2}=}