チュートリアル: DTMから地形湿潤指数（Topographic Wetness Index）を計算する

1. 概要

2. DTMのボイド（Voids）を補間する

提供されているDTMは、オランダのLIDARデータ（AHN3）に由来しています。オランダのAmersfoort投影（EPSG: 28992）のデータです。しかし、DTMにはボイドがあります。ボイドとは、データのないピクセルのことです。これらのボイドは、次に進む前に補間する必要があります。

1. QGISを起動します。

2. 提供されたDTM (dtm.tif) をブランクプロジェクトに追加します。

ボイドがあることがはっきりと見てとれます。

3. プロセッシングツールボックスにいきます: メインメニューのプロセッシング | ツールボックスを選びます。

4. プロセッシングツールボックスで、GDAL | nodata値を内挿値で埋めるを選択します。

5. nodata値を内挿値で埋めるダイアログでdtmを入力レイヤとして選びます。その他の設定はデフォルトのままにしておきます（もしボイドがより大きかったらそれらの設定を変えても良いでしょう）。dtm_voidfilled.tifとして結果を保存します。

6. 実行をクリックします。プロセッシングの後ダイアログを閉じます。

7. レイヤパネルからdtmレイヤを削除します。

結果を見ると、ボイドが補間されていることがはっきりとわかります。

8. DTMを単バンド疑似カラーでスタイル設定します。 カタログのcpt-cityからTopographyランプを使用できます。DTMを陰影図とブレンドします（この方法については他のチュートリアルを参照してください）。

これでDTMの準備ができたので、傾きの計算を進めることができます。

注：QGIS 3.14 では、各計算後にラスター レイヤーの投影を設定する必要があることに注意してください。これは 3.10 LTR バージョンでは異なります。

3. 傾きを度数とラジアンで計算する

TWIを計算するには、ラジアン単位の傾きが必要です。ラスタ解析ツールを使用すると、度単位とパーセンテージ単位の傾きを計算することができます。そこで、まず度単位の傾きを計算し、その値をラジアンに変換します。

1. メインメニューのラスター  | 解析 | 傾斜（slope）...と進みます。

2.  傾きダイアログで、入力レイヤーとして dtm_voidfilled を選択します。Hornの公式ではなく、Zevenbergen Thorneの公式を使用するにチェックを入れます。Zevenbergen Thorneの公式はレリーフが少ない場合に有効ですが、急峻な部分にはHornの式を使用します。出力する傾斜（slope）レイヤとして slopedegrees.tif を選択します。

3. 実行をクリックします。処理が終わったらダイアログを閉じます。

後に自然対数 (ln) を計算することになるので、傾きが 0 度の画素があるかどうかをチェックする必要があります。

4. slopedegreesを右クリックし、プロパティを選びます。

5.レイヤプロパティダイアログで情報タブに進んで、最小値(STATISTICS_MINIMUM)をチェックします。

この場合、最小勾配は0よりも大きいので、次に進むことができます。しかし、最小勾配が0の場合は、以下の手順を実行する必要がありますが、それは何の影響もないので、このデータにも適用できます。

6. キャンセルをクリックして、レイヤープロパティダイアログを閉じます。

7. メインメニューからラスタ | ラスタ計算機と進みます。

8. ラスタ計算機では次の公式を作成します:

( "slopedegrees@1" <= 0 ) * 1 + ( "slopedegrees@1" > 0 ) * "slopedegrees@1"

これは次のことを意味します: セルの傾きが0以下の場合は、値1 (True)を与えて1を掛けます。 0より大きいセルについては、値1 (True)を与えて、その傾きの値を掛けます。結果は、セル <= 0 の場合は 1 になり、セル > 0 の場合は元の傾きになります。

9. 結果をslopedegrees_modified.tifに保存し、OKをクリックします。

度からラジアンに変換するには、ピクセル値に 0.01745 を乗算する必要があります。

10. ラスター計算機に戻り、次の式を計算します。:

"slopedegrees_modified@1" * 0.01745

11. 結果をsloperadians.tifに保存し、OKをクリックします。

12. 結果を黄色から赤にスタイリングします。

13. 他の傾きのレイヤを削除します。

次のステップでは、DTMの各ピクセルの上りの傾斜の集水域を計算します。

4. 上りの斜面の集水域の面積を計算する

ここでは、DTMの各ピクセルについて、上り斜面の集水域の面積を計算していきます。

1. プロセッシングツールボックスに進んで、SAGA | Simulation | Flow accumulation (qm of esp)を選びます。

2. ダイアログでDEMとしてdtm_voidfilledを選択し、デフォルトの前処理の fill sinks temporarilyのままにして、集水域をcontribution_upslope_area.sdatに保存します。

3. 実行をクリックします。処理が終わったらダイアログを閉じます。

結果は、各画素について、上流側の画素の量を表しています。

これでTWIを計算する準備ができました。

5. TWIを計算する

これで、ラジアン単位での傾きと斜面の集水域面積（ピクセル単位）が得られたので、ラスター計算機で地形湿潤指数（TWI)を計算することができます。

1. ラスタ計算機に進むます。

2. 次の公式を書きます:

ln ( ( "contributing_upslope_area@1" * 5 * 5 ) / tan ( "sloperadians@1" ) )

なお、上り斜面の集水域面積（画素量）に、画素の面積（５ｍ×５ｍ）を乗じたものを使用しています。

3. 結果を twi.tif に保存します。

出力CRSが不明な場合は、入力レイヤーに設定する必要があるので注意してください。

4. OKをクリックします。

5. 赤から黄色、青へのカラーランプでスタイルを作成し、陰影レイヤーとブレンドさせるようにします。

結果を3Dで可視化することも可能です。

6. 3DビューでTWIを可視化する