傾斜方向の仕組み

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

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

傾斜方向 ツールは、下りの傾斜角が指している方向を特定します。

出力ラスターの各セルの値は、その位置でサーフェスが向いている方位方向を示しています。 0 度 (真北) ～ 360 度 (真北) で右回りに計測され、完全に円になります。下り方向のない平坦な領域には -1 という値が割り当てられます。

以下の図に、入力標高データセットと出力傾斜方向ラスターを示します。

傾斜方向の出力例

注意:

サーフェスパラメーター ツールは、傾斜方向を新たに実装するので、傾斜方向 ツールの代わりに使用することをおすすめします。 傾斜方向 ツールは平面を 9 個のローカルセルに収めますが、平面は地形の記述子としてはあまり役に立たないことが多く、対象の自然変動がマスクされたり、強調されたりすることがあります。 サーフェスパラメーター ツールは、サーフェスを平面ではなくセルの近傍に収めるため、地形をより自然に適合させることができます。

傾斜方向 ツールは、3 x 3 のセルウィンドウを使用して値を計算するのに対し、サーフェスパラメーター ツールは 3 x 3 から 15 x 15 のセルまでのウィンドウサイズを使用できます。高解像度の標高データには、大きいウィンドウサイズが便利です。これにより、地表面の処理を適切な縮尺で捉えられるようになります。 サーフェスパラメーター ツールにはアダプティブウィンドウオプションもあります。このオプションは地形のローカル変動を評価し、セルごとに最大の適切な近傍サイズを特定します。これは、河川、道路、傾斜の急変によって途切れる、緩やかな同種の地形において便利です。

以前のツールの実行と正確に一致する結果が必要である場合、またはアルゴリズムの正確さよりも高速な実行の方が重要である場合には、傾斜方向 ツールの従来のアプローチを引き続き使用できます。

傾斜方向 ツールを使用すると、次の操作を実行できます:

スキーに最適な斜面を検索する際に、ある山の北側斜面をすべて検出する。
各地の生活の多様性を調べる研究の一環として、ある地域の各場所における直射日光照度を計算する。
雪解け水が最初に流れ込む可能性がある住宅地域を特定する研究の一環として、山岳地域の南に面する傾斜のうち、雪が最初に融ける場所を検出する。
緊急時に航空機が着陸できる平坦地を特定する。

計算方法とエッジ効果

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

平面方法では、2D 直交座標系を使用し、投影された平面上で計算が実行されます。測地線方法では、地球の形を楕円体と見なすことで 3D 直交座標系で計算が実行されます。

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

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

平面方法

平面方法は、傾斜方向を計算する従来の方法です。

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

移動する 3 x 3 の枠が入力ラスターの各セルを巡り、枠の中央にある各セルについて、そのセルの 8 個の近傍セルの値を加味したアルゴリズムを使用して傾斜方向の値が計算されます。これらのセルは \(a\) ～ \(i\) の文字で識別され、傾斜方向の計算対象となるセルには \(e\) が割り当てられます。

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

\([dz/dx] = \Large \frac {((c + 2f + i) * 4 / wght1) - (a + 2d + g) * 4 / wght2)}{8}\)

ここで:

\(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/dx] = \Large \frac {((g + 2 h + i) * 4 / wght3) - (a + 2b + c) * 4 / wght4)}{8}\)

ここで:

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

セル \(e\) の X 方向と Y 方向の両方の変化率を考慮して、傾斜方向は次のように計算されます:

\(aspect = 57.29578 * atan2 ([dz/dy], -[dz/dx])\)

傾斜方向の値は、次の規則に従って、コンパス方向の値 (0 ～ 360 度) に変換されます:

if aspect < 0
cell = 90.0 - aspect
else if aspect > 90.0
cell = 360.0 - aspect + 90.0
else
cell = 90.0 - aspect

平面上での傾斜方向の計算例

例として、移動する枠の中央セルについて、平面上での傾斜方向値を計算します。

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

\([dz/dx] = \Large \frac {((c + 2f + i) * 4 / wght1) - (a + 2d + g) * 4 / wght2)}{8}\)

= ((85 + 170 + 84) * 4 / (1 + 2 + 1) - (101 + 202 + 101) * 4 / (1 + 2 + 1)) / 8
= -8.125

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

\([dz/dy] = \Large \frac {((g + 2h + i) * 4 / wght3 - (a + 2b + c) * 4 / wght4)}{8}\)

= ((101 + 182 + 84) * 4 / (1 + 2 + 1) - (101 + 184 + 85) * 4 / (1 + 2 + 1)) / 8
= -0.375

傾斜方向は、次のように計算されます:

\(aspect = 57.29578 * atan2 ([dz/dy], -[dz/dx])\)

= 57.29578 * atan2 (-0.375, 8.125)
= -2.64

計算された値がゼロより小さいため、最終ルールが次のように適用されます:

\(cell = 90.0 - aspect\)

= 90 - (-2.64)
= 90 + 2.64
= 92.64

中央セル \(e\) の 92.64 という値は、その傾斜方向が東方向であることを示しています。

測地線方法

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

測地座標変換

ECEF (earth-centered earth-fixed、地球中心慣性座標系) 座標系は、地球の中心を原点とする 3D 右手直交座標系です。この座標系では、すべての位置が X、Y、Z 座標で表現されます。次の図は、地心座標により表された目標位置 T の例です。

ECEF 座標系 — サーフェスラスターは、入力座標系から地心 3D 座標系に変換されます。

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

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

\(X = (N(\phi) + h) \space cos \phi \space cos λ\)

\(Y = (N(\phi) + h) \space cos \phi \space sin λ\)

\(Z = (b^2 / a^2 * N(\phi) + h) sin \phi\)

ここで:

\(N(\phi) = a^2 / √(a^2 cos \phi^2 + b^2 sin \phi^2)\)

\(\phi\) は緯度

λ は経度

\(h\) は楕円体の高さ

\(a\) は楕円体の長軸

\(b\) は楕円体の短軸

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

傾斜方向の計算

ある位置における測地線傾斜方向は、楕円体表面に平行である平面上での、北を基準とした下り斜面の方向を表しています。

各位置の傾斜方向を計算するために、3 セル X 3 セルの近傍平面が、最小二乗法 (LSM) を使用して、各処理セルの周りでフィッティングされます。 LSM による最適なフィットの位置では、実際の Z 値とフィッティング後の Z 値の差の二乗和 (SSD: Sum of Squared Differece) (\(dz_i\)) が最小化されます。例については、下の図をご参照ください。

平面は \(z = Ax + By + C\) で表されます。

各セルの中心で、\(dz_i\) は実際の Z 値とフィッティング後の Z 値の差を示します。

\(\large \sum \limits_{i=1}^9\) \(\large dz{_i}{^2}\) が最小になったときに平面が最適にフィッティングされます。

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

楕円体表面の接平面が参照平面と見なされるため、サーフェス法線が平面上に垂直に投影されます。最後に、サーフェス法線の垂直投影と北の間の時計回りの角度 α を計測することにより、測地線傾斜方向が計算されます (上図を参照)。

サーフェスパラメーターツールを使用すべきか

入力ラスター パラメーター値が高解像度でセルサイズが数メートル未満の場合や、特にノイズが多い場合は、このツールの直接 3 × 3 近傍ではなく、サーフェスパラメーター ツールとそのユーザー定義の近傍距離オプションを使用することを検討してください。より大きな近傍を使用することで、ノイズの多いサーフェスの影響を最小限に抑えることができます。また、高解像度のサーフェスを使用する場合は、より大きな近隣を使用することで、地形やサーフェスの特性をより上手く表現することができます。

GPU の使用

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

参考情報

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

E.J.Krakiwsky, and D.E.Wells, 1971. Coordinate Systems In Geodesy (GEODESY AND GEOMATICS ENGINEERING, UNB), LECTURE NOTES, No16, 1971, pp. 18-38

Lancaster, P. and Šalkauskas, K. 1986. Curve and Surface Fitting: An Introduction. London: Academic Press

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.