QGISラスタ計算機で斜面方位をカテゴリー化

考古学でよく使う斜面方位ラスタを4方位のベクタ地図に変換します。

こんな斜面方位ラスタ

4方位に区分したベクタに変換したい

ラスタ計算機で斜面方位ラスタをカテゴリー化

ラスタ演算式は次のとおりです。なお、東が0で半時計回りに増加するラスタ地図です。GRASS GISで作成した傾斜方位地図などがデフォルトで東が0の半時計回りとなっています。

完成したラスタ地図は

• 東10
• 北20
• 西30
• 南40

が代入されます。

数式説明

("Aspect@1">0)*("Aspect@1"<=45)*10+
("Aspect@1">45)*("Aspect@1"<=135)*20+
("Aspect@1">135)*("Aspect@1"<=225)*30+
("Aspect@1">225)*("Aspect@1"<=315)*40+
("Aspect@1">315)*10

• 「"Aspect@1"」はAspectレイヤのバンド1を意味します。
• 「"Aspect@1">0」が真なら計算機は「1」を返し、偽なら「0」を返します。
• 続く「"Aspect@1"<=45」が真なら同様に「1」を、偽なら「0」を返します。
• 「("Aspect@1">0)*("Aspect@1"<=45)」では0より大きく、45以下でない値はすべて0が返されます。
• 「("Aspect@1">0)*("Aspect@1"<=45)*10」では、0より大きく45以下の値が真の場合に10をかけるので、0〜45の値をとるピクセルには10が代入されます。
• 同様に45〜135（北）では20が代入され、135〜225では30が代入され、225〜315では40が代入されます。
• 真となるのはどれか一つの項なので全部の項を足し合わせると真となる項の数字だけが残ります。

上記の計算式は「カッパ出没マップを作成する」を参照して作成しました。合わせてどうぞ。

ラスタのポリゴン化

gdal_polygonize.py コマンドを使用します。

# 引数は入力ファイル（Aspect_reclass.tif）、出力ファイル（Aspect.shp）の順に指定
# -b オプションはバンドの指定「1」を指定します。
gdal_polygonize.py Aspect_reclass.tif Aspect.shp -b 1

ラスタ計算機で対応する数字を文字列に変換

方位を表す整数値が「DN」フィールドに入力されているので、文字列（東西南北）に変換します。

計算式は以下のとおりです。特に説明することはないと思います。

CASE 
WHEN "DN"=10 THEN  '東'
WHEN "DN"=20 THEN  '北'
WHEN "DN"=30 THEN  '西'
WHEN "DN"=40 THEN  '南'
END

aspectフィールドに文字列の方位が入力されていれば完成。