山体阴影工具的工作原理

山体阴影工具通过为栅格中的每个像元确定照明度,来获取表面的假定照明度。 通过设置假定光源的位置并计算每个像元相对于相邻像元的照明度值来实现此目的。 它可以显著增强用于分析或图形显示的表面的可视化效果,尤其是在使用透明度时。

默认情况下,阴影和光是与从 0 到 255(从黑色到白色递增)的整数关联的灰色阴影。

山体阴影参数

为任何特定位置创建山体阴影地图时,主要因素是太阳在天空中的位置。

方位角

方向角指的是太阳的角度方向,是以北为基准方向在 0 到 360 度范围内按顺时针进行测量的。 90 度的方位角指向为东。 默认方位角为 315 度 (NW)。

太阳方位角为 315º
山体阴影的默认太阳方位角(方向)为 315º。

高度

高度指的是照明源高出地平线的角度或坡度。 高度的单位为度,范围为 0(位于地平线上)到 90(位于头上)之间。 默认值为 45 度。

太阳高度角为 45º
山体阴影的默认太阳高度角为 45º。

山体阴影示例

下面的山体阴影示例的方位角为 315 度,高度角为 45 度。

山体阴影输出示例
山体阴影输出示例

使用山体阴影进行显示

通过将高程栅格放置在山体阴影栅格的上方,然后对高程栅格的透明度进行调整,可创建外观精美的地表地貌图。 要了解如何调整显示和外观,请参阅影像外观

使用透明度将高程栅格和山体阴影相结合
通过使用透明度将高程和山体阴影栅格相结合来创建晕渲地貌栅格。

您可以添加其他图层,例如土地利用类型、植被、道路或河流,以进一步增加显示中的信息内容。

在分析中使用山体阴影

通过对阴影建模(默认选项),您可以计算局部照明以及像元是否在阴影中。

通过对阴影进行建模,可识别一天的特定时间内将要落入其他像元阴影的各个像元。 位于其他像元阴影中的像元编码为 0;所有其他像元的编码为介于 1 和 255 之间的整数。 可将所有大于 1 的值重分类为 1,从而生成二元输出栅格。 在以下示例中,黑色区域位于阴影中。 这两张图像中的方位角相同,只是太阳角度(高度角)已修改。

太阳角度较低的阴影
太阳角度较低会产生更多的阴影区域。
太阳角度较高的阴影
太阳角度较高会产生更少的阴影区域。

山体阴影计算

要计算阴影值,首先需要提供光照源的高度角和方向角。 在计算坡度和坡向的同时会处理这些值,从而确定输出栅格中每个像元的最终山体阴影值。

山体阴影算法

用于计算山体阴影值的算法如下:

(1) Hillshade = 255.0 * ((cos(Zenith_rad) * cos(Slope_rad)) + (sin(Zenith_rad) * sin(Slope_rad) * cos(Azimuth_rad - Aspect_rad)))

请注意,如果山体阴影值的计算结果小于 0,则输出像元值将等于 0。

计算入射角度

以高于地平线的角度指定照明源的高度角。 但是,用于计算山体阴影值的公式要求该角度以弧度表示,并且是相对于垂直线的偏转。 将从表面垂直向上的方向(正上方)标注为天顶。 天顶角为天顶点与光源方向的夹角,是高度角的 90 度余角。 要计算照明角度,首先要将高度角转换为天顶角。 接下来,将角度转换为弧度。

将高度角转换为天顶角:

(2) Zenith_deg = 90.0 - Altitude

转换为弧度:

(3) Zenith_rad = Zenith_deg * pi / 180.0

计算入射方向

以度为单位指定照明源的方向和方位角。 山体阴影公式要求方位角采用弧度作为单位。 首先,将方位角从地理单位(罗盘方向)转换为数学单位(直角)。 接下来,将方位角转换为弧度。

转换方位角的方法:

(4) Azimuth_math = 360.0 - Azimuth + 90.0

请注意,如果 Azimuth_math >= 360.0,则将出现以下情况:

(5) Azimuth_math = Azimuth_math - 360.0

转换为弧度:

(6) Azimuth_rad = Azimuth_math * pi / 180.0

计算坡度和坡向

移动 3 x 3 窗口将移动至输入栅格中的每个像元,对于位于窗口中心的每个像元,其坡向和坡度值将通过一种将纳入 8 个相邻像元值的算法进行计算。 这些像元使用字母 ai 进行标识,其中 e 表示当前正在计算坡向的像元。

像元 e 在 x 方向上的变化率将通过以下算法进行计算:

(7) [dz/dx] = ((c + 2f + i) - (a + 2d + g)) / (8 * cellsize)

像元 e 在 y 方向上的变化率使用以下算法进行计算:

(8) [dz/dy] = ((g + 2h + i) - (a + 2b + c)) / (8 * cellsize)

表面上每个像元的最陡坡降是坡度。 结合 z 因子的以弧度为单位的坡度计算算法如下:

(9) Slope_rad = ATAN (z_factor * √ ([dz/dx]2 + [dz/dy]2))

最陡坡坡降面向的方向为坡向。 以弧度为单位的坡向的范围从 0 到 2pi,其中 0 表示朝向东。 坡向是根据以下算法中的规则确定的:

(10) If [dz/dx] is non-zero: Aspect_rad = atan2 ([dz/dy], -[dz/dx]) if Aspect_rad < 0 then Aspect_rad = 2 * pi + Aspect_rad If [dz/dx] is zero: if [dz/dy] > 0 then Aspect_rad = pi / 2 else if [dz/dy] < 0 then Aspect_rad = 2 * pi - pi / 2 else Aspect_rad = Aspect_rad

山体阴影计算示例

在本示例中,将计算移动窗口内中心像元的山体阴影值。

输入高程栅格数据

像元大小为 5 个单位。 将使用 45 度的默认高度角和 315 度的默认方位角

  • 入射角度

    根据等式 2 进行天顶角计算如下:

    (2) Zenith_deg = 90.0 - Altitude = 90.0 - 45.0 = 45.0

    根据等式 3 转换为弧度如下:

    (3) Zenith_rad = Zenith_deg * pi / 180.0 = 45.0 * 3.1415926536 / 180.0 = 0.7853981634

  • 入射方向

    根据等式 4 将方位角从地理角度转换为数学角度的计算如下:

    (4) Azimuth_math = 360.0 - Azimuth + 90.0 = 360.0 - 315.0 + 90.0 = 135.0

    根据等式 6 将方位角转换为弧度如下:

    (6) Azimuth_rad = Azimuth_math * pi / 180.0 = 135.0 * 3.1415926536 / 180 = 2.3561944902

  • 坡度和坡向

    中心像元 e 在 x 方向上的变化率的计算结果如下:

    (7) [dz/dx] = ((c + 2f + i) - (a + 2d + g)) / (8 * cellsize) = ((2483 + 4966 + 2477) - (2450 + 4904 + 2447)) / (8 * 5) = (9926 - 9801) / 40 = 3.125

    中心像元 e 在 y 方向上的变化率的计算结果如下:

    (8) [dz/dy] = ((g + 2h + i) - (a + 2b + c)) / (8 * cellsize) = (2447 + 4910 + 2477) - (2450 + 4922 + 2483) / (8 * 5) = (9834 - 9855) / 40 = -0.525

    坡度角的计算结果如下:

    (9) Slope_rad = ATAN ( z_factor * √ ([dz/dx]2 + [dz/dy]2)) = atan(1 * sqrt((3.125 * 3.125) + (-0.525 * -0.525))) = atan(1 * sqrt(10.04125 + 0.275625)) = atan(1 * 3.1687931457) = 1.2651101670

    根据规则 10 计算 Aspect_rad 角(因为本例中 dz/dx 为非零值)所得出的结果如下:

    Aspect_rad = atan2 ([dz/dy], -[dz/dx]) = atan2(-0.525, -3.125) = -2.9751469600

    由于计算得出的值小于零,因此使用此规则部分:

    Aspect_rad = 2 * pi + Aspect_rad = 2 * 3.1415926536 + -2.9751469600 = 3.3080383471

  • 山体阴影

    山体阴影最终计算如下:

    Hillshade = 255.0 * ((cos(Zenith_rad) * cos(Slope_rad)) + (sin(Zenith_rad) * sin(Slope_rad) * cos(Azimuth_rad - Aspect_rad))) = 255.0 * ((cos(0.7857142857) * cos(1.26511)) + (sin(0.7857142857) * sin(1.26511) * cos(2.3571428571 - 3.310567))) = 153.82

因为输出栅格为整型,所以中心像元 e 的阴影值为 154。

参考资料

Burrough, P. A., and R. A. McDonell. 1998. Principles of Geographical Information Systems. New York: Oxford University Press.