TL;DR
以下のコマンドでCS立体図を作成出来ます。
pip install csmap-py # インストール
csmap dem.tif csmap.tif # 処理実行
あまりよくテスト出来ていないので、フィードバックお待ちしております（Issueまで)。
CS立体図
CS立体図とは、DEMをもとにした地形判読のための地形表現図法です。以下の林野庁の資料が詳しいです。
- 赤い箇所は凸地形、青い箇所は凹地形
- 傾斜率が高いと、うっすら茶色っぽい
ざっくりこのような着色をしておくことで、平面上の図面でも地形判読が容易となります。
CS立体図の作成方法
CS立体図を作成するには、オープンなものだと上記のQGISプラグインに限られていると思われます（弊社代表の作品）。一方、こちらのツールはQGISに内蔵のSAGAに依存しており、そして最近のQGISはSAGAが内蔵されなくなったという問題がありました。なので新しく作りました。
csmap-py
- シンプルな依存関係（rasterioとnumpy）
- 高速
- 並列処理
- 省メモリ
使い方
pip install csmap-py # インストール
csmap dem.tif csmap.tif # 処理実行
csmap dem.vrt csmap.tif # 仮想ラスターや、GDALがサポートしているラスターを読める
csmap dem.tif csmap.tif --max_workers 4 # max_workersオプションで並列処理数を設定可能
csmap dem.tif csmap.tif --gf_size 3 --gf_sigma 1.0 --curvature_size 1 # CS立体図作成にまつわるパラメータを変更できる
仕様
- 入力のDEMに対し、外縁部の一定部分（＝ガウシアンフィルタの影響範囲）はトリミングされます
- No-Dataの部分はあまり気を使っておらず、一定の計算が実行されます
- 着色はCSMapMakerを参考にしていますが、厳密な再現はしていません
- パラメータのデフォルト値はCSMapMakerに準拠しています
ベンチマーク
- 1mDEM, 解像度17602x6794, 10個のTIFFファイルから成る仮想ラスターを処理
- Macbook Air M2
time csmap merged.vrt csmap.tif --max_worker 1
real 0m29.345s
user 0m27.055s
sys 0m3.457s
time csmap merged.vrt csmap.tif --max_worker 16
real 0m10.848s
user 0m43.075s
sys 0m5.517s
謝辞
本ツールは、「SAGA依存しないでNumpyで書けばいいんじゃない、しらんけど（意訳）」というわたしの雑なSlack投稿を拾ってコア部分（SAGAでやってる処理のNumpy実装）を書いてくれた@ciscornさんとのゆるいコールアンドレスポンスにより出来ています、感謝。