200 m 先の勾配で土石流停止点を抽出

QGISで土石流停止点を抽出する方法

動作確認環境: QGIS 3.28 / DEM10B (1mメッシュ)
最終更新: 2025年7月

1準備するもの

• 流下ライン（Shapefile または GeoPackage形式）
• DEM（数値標高モデル）（GeoTIFF形式）
• QGIS 3.24以降

2Pythonスクリプトの実行

実行手順:
QGISメニューから プラグイン → Python コンソール を開き、下記のスクリプトを

"""
make_points_and_slope.py
  ① ライン上に等間隔ポイント生成（distance フィールド付き）
  ② DEM から標高取得 → *_pts.gpkg に points_z 保存
  ③ 200 m 先の勾配を「200ｍ先の勾配」列へ書込
"""

from qgis.PyQt.QtWidgets import QFileDialog, QInputDialog
import os, math, processing
from qgis.core import QgsProject, QgsVectorLayer, QgsField, edit

# ── 入力 ──
line, _ = QFileDialog.getOpenFileName(None,'ライン選択','','Line layer (*.shp *.gpkg)')
if not line: raise Exception('ライン未選択')
dem , _ = QFileDialog.getOpenFileName(None,'DEM 選択','','GeoTIFF (*.tif)')
if not dem: raise Exception('DEM 未選択')

interval, ok = QInputDialog.getInt(None,'点間隔(m)','1 / 20 / 50 など',1,1,500,1)
if not ok: interval=1
lag_n = max(1, round(200/interval))  # 200m ÷ 点間隔

# ── 出力 gpkg ──
folder = os.path.dirname(line)
base   = os.path.splitext(os.path.basename(line))[0]
gpkg   = os.path.join(folder,f'{base}_pts.gpkg')
if os.path.exists(gpkg): os.remove(gpkg)

# 1️⃣ 等間隔ポイント（distance）
pts = processing.run('native:pointsalonglines',{
    'INPUT':line,'DISTANCE':interval,'INCLUDE_DISTANCE':True,'OUTPUT':'memory:'
})['OUTPUT']

# 2️⃣ DEM サンプリング
pts_z_mem = processing.run('qgis:rastersampling',{
    'INPUT':pts,'RASTERCOPY':dem,'COLUMN_PREFIX':'Z_','OUTPUT':'memory:'
})['OUTPUT']

processing.run('native:savefeatures',{
    'INPUT':pts_z_mem,'OUTPUT':gpkg,'LAYER_NAME':'points_z'
})

uri=f'{gpkg}|layername=points_z'
pts_z=QgsVectorLayer(uri,'points_z','ogr')
if not pts_z.isValid(): raise Exception('points_z 読込失敗')
QgsProject.instance().addMapLayer(pts_z)

# 標高列
zfield = next((f.name() for f in pts_z.fields() if f.name().startswith('Z_')),None)
if not zfield: raise Exception('Z_ 列が無い')

# 勾配列追加
if '200ｍ先の勾配' not in [f.name() for f in pts_z.fields()]:
    with edit(pts_z):
        pts_z.dataProvider().addAttributes([QgsField('200ｍ先の勾配',6,len=10,prec=6)])
    pts_z.updateFields()

slope_idx = pts_z.fields().indexOf('200ｍ先の勾配')

# distance 昇順で計算
f_sorted = sorted(pts_z.getFeatures(),key=lambda f:f['distance'])
zvals = [f[zfield] for f in f_sorted]

with edit(pts_z):
    for i,f in enumerate(f_sorted):
        if i >= lag_n and zvals[i] is not None and zvals[i-lag_n] is not None:
            dz = zvals[i-lag_n]-zvals[i]
            f[slope_idx] = math.degrees(math.atan(dz/200.0))
        else:
            f[slope_idx] = None
        pts_z.updateFeature(f)

print('✅ 完了：points_z に「200ｍ先の勾配」を追加 →',gpkg)

💡 スクリプト実行時に、流下ラインファイル、DEMファイル、点間隔を順次選択するダイアログが表示されます。

3勾配角度による停止点の抽出

生成されたポイントから、勾配角度に応じて停止点を抽出できます。属性テーブルの高度なフィルタで以下の式を使用してください。

抽出条件	フィルタ式
20m間隔のポイント	"distance" % 20 = 0
50m間隔のポイント	"distance" % 50 = 0
20m間隔 & 勾配2°未満	"distance" % 20 = 0 AND "200ｍ先の勾配" < 2

💡 抽出した停止点を別レイヤーに保存する場合は、フィルタ適用後に右クリック → 選択フィーチャをレイヤへ保存 を選択してください。

スクリプトの動作概要

処理の流れ

1. 等間隔ポイント生成: 流下ライン上に指定した間隔でポイントを生成し、各点に距離情報を付与
2. 標高値の取得: DEMからラスタサンプリングを行い、各ポイントの標高値を取得
3. 勾配計算: 各ポイントから200m先（上流側）のポイントとの高度差から勾配角度を計算
4. 結果出力: 「200ｍ先の勾配」フィールドに角度（度）を格納したGeoPackageファイルを出力

勾配の計算方法

勾配角度は以下の式で計算されます：

勾配角度 = arctan((上流側標高 - 現在点標高) / 200) × 180 / π

💡 土石流は勾配が緩やかになると堆積する傾向があります。一般的に勾配2°未満で停止する可能性が高いとされています。

このスクリプトは災害予測や防災計画の基礎資料作成に活用できます。

© 2025 / Creative Commons BY-SA