傾斜角 (Slope) の仕組み—ArcGIS Pro

Spatial Analyst のライセンスで利用可能。

3D Analyst のライセンスで利用可能。

[傾斜角 (Slope)] ツールはラスターサーフェスの各セルの傾斜を特定します。傾斜角の値が小さくなるほど地表は平らになり、傾斜角の値が大きくなるほど地表が急勾配になります。

出力傾斜角ラスターは、「度」と「パーセント (傾斜率) 」の 2 種類の単位で計算できます。傾斜率は、「高低差÷距離 × 100」で求めることができます。たとえば、下図の三角形 B の場合、角度 45 度なので「高低差」と「距離」が等しく、傾斜率は 100% になります。三角形 C のように、傾斜角が垂直 (90 度) に近づくほど傾斜率が大きくなります。

多くの場合、[傾斜角 (Slope)] ツールは、次の図に示すように標高データセットに対して実行されます。この出力傾斜角ラスターでは、傾斜の急な斜面が暗褐色で表示されています。

人口データなど、他の種類の連続データでこのツールを使用して、値 (人口) の急激な変化を調べることもできます。

計算方法とエッジ効果

傾斜角を計算する方法は 2 つあります。[方法] パラメーターを使用して [平面] 計算を実行するか、[測地線] 計算を実行するかを選択できます。

平面方法の場合、傾斜角は、特定のセルの値と、その隣接セルの値とを比較したときの最大変化率として計測されます。二次元の直交座標系を使用し、投影された平面上で計算が実行されます。傾斜角の値は平均最大法 (Burrough、1998) を用いて計算されます。

測地線方法では、地球の形を楕円体と見なすことで三次元の直交座標系で計算が実行されます。傾斜角の値は、地形サーフェスと参照先の測地基準系との間の角度を計測することで計算されます。

平面計算と測地線計算のいずれもが、3 セル X 3 セルの近傍 (移動する枠) の使用によって実行されます。各近傍で、処理 (中央の) セルが NoData の場合、出力は NoData になります。また、計算を実行するには、処理セルに隣接しているセルのうち、少なくとも 7 つのセルが有効な値を保持している必要があります。有効なセルの数が 7 未満である場合、計算は実行されず、その処理セルの出力は NoData になります。

出力ラスターの最も外側のロウとカラムにあるセルは NoData になります。これは、これらのセルが、入力データセットの境界沿いにあり、必要なだけの有効な隣接セルに囲まれていないためです。

平面方法

このツールは、対象セルとその近隣セルとを比較し、それらの値の最大変化率を計算します。原則として、対象セルとその 8 つの近隣セルの高度を比較したとき、距離に対するそれらの最大高度変化が、対象セルからの最も急な降下傾斜となります。

平面上での傾斜角のアルゴリズム

中央セルから水平方向 (dz/dx) および垂直方向 (dz/dy) へのサーフェスの変化率 (差分) が傾斜角を決定します。傾斜角 (slope_radians) を計算するための基本的なアルゴリズムは次のとおりです。

 slope_radians = ATAN ( √ ([dz/dx]² + [dz/dy]²) )

通常、傾斜角は「度」で計測され、次のアルゴリズムでもこの単位を使用しています。

 slope_degrees = ATAN ( √ ([dz/dx]² + [dz/dy]²) ) * 57.29578

メモ:

ここに示されている値 57.29578 は、180/pi の算出値を切り捨てたものです。

傾斜角アルゴリズムは次のように表すこともできます。

 slope_degrees = ATAN (rise_run) * 57.29578

ここで、
```
 rise_run = √ ([dz/dx]² + [dz/dy]²]
```

中央セルの値、およびそれを囲む 8 つのセルの値は、横方向の差分と縦方向の差分を決定します。隣接セルは a ～ i の文字で識別され、傾斜方向の計算対象となるセルには e が割り当てられます。

セル e の X 方向の変化率は、次のアルゴリズムで計算されます。

  [dz/dx] = ((c + 2f + i)*4/wght1 - (a + 2d + g)*4/wght2) / (8 * x_cellsize)

ここで、
wght1 と wght2 は有効なセルの水平加重計数値です。

次に例を示します。
- c、f、i の値がすべて有効である場合、wght1 = (1+2*1+1) = 4
- i が NoData の場合、wght1 = (1+2*1+0) = 3
- f が NoData の場合、wght1 = (1+2*0+1) = 2
隣接位置が a、d、g である場合を除き、同じようなロジックが wght2 にも適用されます。

セル e の Y 方向の変化率は、次のアルゴリズムで計算されます。

  [dz/dy] = ((g + 2h + i)*4/wght3 - (a + 2b + c)*4/wght4) / (8 * y_cellsize)

ここで、
wght3 と wght4 は [dz/dx] の計算で見られたものと同じ概念を表しています。

平面上での傾斜角の計算例

例として、次のような移動ウィンドウの中心セルの傾斜角値を計算します。

中央セル e の X 方向の変化率は、次のように計算されます。

  [dz/dx] = ((c + 2f + i)*4/wght1 - (a + 2d + g)*4/wght2) / (8 * x_cellsize)          = ((50 + 60 + 10)*4/(1+2+1) - (50 + 60 + 8)*4/(1+2+1)) / (8 * 5)          = (120 - 118) / 40          = 0.05

中央セル e の Y 方向の変化率は、次のように計算されます。

  [dz/dy] = ((g + 2h + i)*4/wght3 - (a + 2b + c)*4/wght4) / (8 * y_cellsize)          = ((8 + 20 + 10)*4/(1+2+1) - (50 + 90 + 50)*4/(1+2+1)) / (8 * 5)          = (38 - 190 ) / 40          = -3.8

X 方向と Y 方向の変化率を基に、中央セル e の傾斜角を次のように計算します。

  rise_run = √ ([dz/dx]² + [dz/dy]²)           = √ ((0.05)² + (-3.8)²)           = √ (0.0025 + 14.44)           = 3.80032

  slope_degrees = ATAN (rise_run) * 57.29578                = ATAN (3.80032) * 57.29578                = 1.31349 * 57.29578                = 75.25762

セル e の傾斜角値 (整数) は 75 度になります。

測地線方法

測地線方法では、地球の形を楕円体と見なすことで、地心三次元座標系 (Earth Centered、Earth Fixed (ECEF) 座標系とも呼ばれる) での傾斜角を計測します。データセットの投影方法は、この計算結果に影響しません。入力ラスターの Z 単位が空間参照で定義されている場合は、この単位が使用されます。入力の空間参照で Z 単位が定義されていない場合は、Z 単位パラメーターを使用して定義する必要があります。測地線方法では、平面方法と比べて、より精度の高い傾斜角が得られます。

測地座標変換

ECEF 座標系は、地球の中心を原点とする三次元の右手直交座標系です。この座標系では、すべての位置が X、Y、Z 座標で表現されます。次の図は、地心座標により表された目標位置 T の例です。

測地線計算では、測地座標 (緯度 φ、経度 λ、高度 h) に基づいて計算される X, Y, Z 座標を使用します。入力サーフェスラスターの座標系が投影座標系 (PCS) である場合は、最初に、各位置に測地座標が含まれているラスターが地理座標系 (GCS) に再投影され、その次に ECEF 座標系に変換されます。高度 h (Z 値) は楕円形表面を基準とする楕円体の高さです。次の図をご参照ください。

ECEF 座標を測地座標 (緯度 φ、経度 λ、高度 h) に変換するには、次の式を使用します。

X = (N(φ)+h)cosφcosλ

Y = (N(φ)+h)cosφsinλ

Z = (b²/a²*N(φ)+h)sinφ

ここで、
- N( φ ) = a²/ √(a²cosφ²+b²sinφ²)
- φ = 緯度
- λ = 経度
- h = 楕円体の高さ
- a = 楕円体の長軸
- b = 楕円体の短軸

上記の式で、楕円体の高さはメートル単位です。ご使用の入力ラスターの Z 単位がこれ以外の単位で指定されている場合は、内部的にメートルに変換されます。

傾斜角の計算

測地線傾斜角は、地形サーフェスと楕円体表面との間に形成される角度を表しています。楕円体表面に平行であるサーフェスの傾斜角は 0 になります。各位置の傾斜角を計算するために、3 セル X 3 セルの近傍平面が、最小二乗法 (LSM) を使用して、各処理セルの周りでフィッティングされます。LSM による最適なフィットの位置では、実際の Z 値とフィッティング後の Z 値の差の二乗和 (SSD: Sum of Squared Differece) (dz_i) が最小化されます。例については、下の図をご参照ください。

平面は z = Ax + By + C で表されます。各セルの中心で、dz_i は実際の Z 値とフィッティング後の Z 値の差を示します。

∑⁹_i=1dz_i² が最小になったときに平面が最適にフィッティングされます。

平面のフィッティング後、そのセルの位置でサーフェスの法線が計算されます。また、同じ位置で、楕円表面の接平面に垂直である楕円体の法線も計算されます。

度単位の傾斜角は、楕円体の法線と地形サーフェスの法線との間の角度 (ここでは β で表す) から計算されます。上図の角度 α は測地線傾斜角であり、ジオメトリ一致の法則に従って角度 β と同じになります。

勾配 (%) として傾斜角を計算するには、次の式を使用します。

Slope_PercentRise = ATAN(β) * 100%

GPU の使用

測地線方法の場合、このツールは、特定の GPU ハードウェアをシステムにインストールしている場合に、パフォーマンスの向上が可能です。このツールをサポートする方法、構成する方法、および有効化する方法の詳細については「Spatial Analyst での GPU 処理」をご参照ください。

参照

Burrough, P. A., and McDonell, R. A., 1998. Principles of Geographical Information Systems (Oxford University Press, New York), 190 pp.

Marcin Ligas, and Piotr Banasik, 2011. Conversion between Cartesian and geodetic coordinates on a rotational ellipsoid by solving a system of nonlinear equations (GEODESY AND CARTOGRAPHY), Vol. 60, No 2, 2011, pp. 145-159

B. Hofmann-Wellenhof, H. Lichtenegger and J. Collins, 2001. GPS - theory and practice. Section 10.2.1. p. 282.

David Eberly 1999. Least Squares Fitting of Data (Geometric Tools, LLC), pp. 3.