カテゴリー: FamilySafetyMap

tsunami_union を保存 → 交差チェックをそれに切替
上位N→ネットワーク最短の選定を導入（経路品質UP）
アイソクローンのレイヤをGPKG/地図に追加（“どこまで行けるか”の可視化）
警告フラグをUIに表示
GraphDB CSV出力で家族アプリへの橋渡し

を行う

1) tsunami_union をGPKGにキャッシュ（毎回のunion計算を省略）
これで気になったのが
このtsunami_union をGPKGにキャッシュしたものは現在地、もしくはすでに計測したもののみ有効？
もし新規の場所で避難することがあれば対応したいので
すでにあるルートならキャッシュ、
そうでないのなら新規ルートを使うようにしたい

結論：tsunami_union は“現在地に依存しない”全域キャッシュです。
静岡県の「津波浸水想定（A40）」レイヤ全体を union したものなので、同じデータ（A40）がソースならどの場所でも再利用可能。新規の場所でもそのまま使えます。再生成が必要なのは、A40を更新・差し替えたときだけです。

A) tsunami_union（全域）のキャッシュ方針
* 有効範囲：A40レイヤ（静岡版）全域。現在地・ルートに非依存。
* 再生成条件：A40が変わった可能性があるとき（例：ファイル更新、feature数やBBOXが変化）。
* 高速化：ルート検証時は tsunami_union をその場のAOIでクリップして使うと軽いです。

例：キャッシュの検証＆必要なら再生成

import geopandas as gpd
from fiona import listlayers

GPKG = "outputs/shizuoka_hazard.gpkg"

tsunami = gpd.read_file(GPKG, layer="tsunami_inundation")
feat_count = len(tsunami)
bbox = tuple(tsunami.total_bounds)  # (xmin,ymin,xmax,ymax)

# 既存 union があるか＆メタデータが一致するか
need_rebuild = True
if "tsunami_union" in listlayers(GPKG):
    union = gpd.read_file(GPKG, layer="tsunami_union")
    # 簡易同一性チェック（必要なら厳格化）
    if "meta_feat" in union.columns and "meta_bbox" in union.columns:
        need_rebuild = not (int(union.loc[0,"meta_feat"]) == feat_count and
                            tuple(map(float, union.loc[0,"meta_bbox"].split(","))) == bbox)

if need_rebuild:
    flood_union = gpd.GeoDataFrame(
        {"meta_feat":[feat_count], "meta_bbox":[",".join(map(str,bbox))]},
        geometry=[tsunami.geometry.union_all()], crs=tsunami.crs
    )
    flood_union.to_file(GPKG, layer="tsunami_union", driver="GPKG")
else:
    flood_union = gpd.read_file(GPKG, layer="tsunami_union")

例：ルート検証時は AOI で軽量化

# route_gdf: ルート LineString (EPSG:4326)
AOI_m = 15000  # 15km
aoi = route_gdf.to_crs("EPSG:32654").buffer(AOI_m).to_crs(4326)
flood_union_aoi = gpd.overlay(flood_union, gpd.GeoDataFrame(geometry=aoi, crs=4326), how="intersection")
# 以後は flood_union_aoi で交差・近接判定（軽い）

B) 「既に計測したルート」はルート側をキャッシュ（任意）
tsunami_union は“全域・一回”でOKですが、ルートの結果は現在地や避難所・パラメータに依存するので、キー付きでGPKGに保存して再利用が有効です。
例：ルート結果のキャッシュ（キー = 端点・モード・半径）

import hashlib, json, geopandas as gpd

def route_key(home_wgs, dest_wgs, mode="walk", radius_m=12000):
    payload = {"home":home_wgs, "dest":dest_wgs, "mode":mode, "r":radius_m}
    return hashlib.md5(json.dumps(payload, sort_keys=True).encode()).hexdigest()

key = route_key(home_wgs, (dest_pt.y, dest_pt.x), mode="walk", radius_m=GRAPH_RADIUS_M)

# 取り出し（あれば使う）
try:
    cache = gpd.read_file(GPKG, layer="route_cache")
    hit = cache.loc[cache["rid"]==key]
    if len(hit):
        route_gdf = hit.set_crs(4326)  # そのまま使える
    else:
        raise IOError
except Exception:
    # 計算 → route_gdf を作ったら保存
    route_gdf["rid"] = key
    route_gdf.to_file(GPKG, layer="route_cache", driver="GPKG", mode="a")  # 追記

A40をアップデートしたら、route_cacheは無効化（レイヤを削除、またはキーにA40のメタ情報を混ぜる）するのが安全です。

A40をアップデートしたら、route_cacheは無効化（レイヤを削除、またはキーにA40のメタ情報を混ぜる）するのが安全です。

ということで

ルート再利用性を高めたい場合は、ルート結果をキー付きでGPKGにキャッシュする機能を追加したい

import hashlib, json, geopandas as gpd

def route_key(home_wgs, dest_wgs, mode="walk", radius_m=12000):
    payload = {"home":home_wgs, "dest":dest_wgs, "mode":mode, "r":radius_m}
    return hashlib.md5(json.dumps(payload, sort_keys=True).encode()).hexdigest()

key = route_key(home_wgs, (dest_pt.y, dest_pt.x), mode="walk", radius_m=GRAPH_RADIUS_M)

# 取り出し（あれば使う）
try:
    cache = gpd.read_file(GPKG, layer="route_cache")
    hit = cache.loc[cache["rid"]==key]
    if len(hit):
        route_gdf = hit.set_crs(4326)  # そのまま使える
    else:
        raise IOError
except Exception:
    # 計算 → route_gdf を作ったら保存
    route_gdf["rid"] = key
    route_gdf.to_file(GPKG, layer="route_cache", driver="GPKG", mode="a")  # 追記

これはどこに追記する？
入れる場所は**「目的地（dest_pt）が決まった直後、経路計算（nx.shortest_path）に入る前」です。
ヒットしたらそのまま使い、ミスなら今までどおり計算→検証→保存しつつキャッシュへ追記**します

# ==== ルートキャッシュ（GPKG: route_cache） ====
import hashlib, json, geopandas as gpd
from fiona import listlayers

GPKG = "outputs/shizuoka_hazard.gpkg"

def route_key(home_wgs, dest_wgs, mode="walk", radius_m=GRAPH_RADIUS_M, hazard_meta=None):
    payload = {"home": home_wgs, "dest": dest_wgs, "mode": mode, "r": radius_m, "hz": hazard_meta}
    return hashlib.md5(json.dumps(payload, sort_keys=True).encode()).hexdigest()

# A40が変わったら無効化できるように簡易メタを混ぜる（件数＋BBOX）
hazard_meta = None
try:
    if "tsunami_union" in listlayers(GPKG):
        # union保存時にメタ列を持たせていればそれを使う／無ければ tsunami_inundation から作る
        tu = gpd.read_file(GPKG, layer="tsunami_union")
        if {"meta_feat","meta_bbox"} <= set(tu.columns):
            hazard_meta = f'{int(tu.loc[0,"meta_feat"])}|{tu.loc[0,"meta_bbox"]}'
    if hazard_meta is None:
        ti = gpd.read_file(GPKG, layer="tsunami_inundation")
        hazard_meta = f'{len(ti)}|{",".join(map(str, ti.total_bounds))}'
except Exception:
    pass

key = route_key(home_wgs, (dest_pt.y, dest_pt.x), mode="walk", radius_m=GRAPH_RADIUS_M, hazard_meta=hazard_meta)

route_gdf_cached = None
if "route_cache" in listlayers(GPKG):
    cache = gpd.read_file(GPKG, layer="route_cache")
    hit = cache.loc[cache["rid"] == key]
    if len(hit):
        route_gdf_cached = hit.set_crs(4326)
        print("↩︎ ルート: キャッシュヒット")

# ヒット時の扱い：後続の検証/表示にそのまま使えるよう最低限の変数を整備
if route_gdf_cached is not None:
    route_gdf = route_gdf_cached[["meters","eta_min_walk","eta_min_walk_safe","eta_min_drive","geometry"]].copy()
    # route_latlon / meters 等は既存列から復元（既に保存してある前提）
    meters = float(route_gdf["meters"].iloc[0])
    eta_walk_min  = float(route_gdf["eta_min_walk"].iloc[0])
    eta_walk_safe = float(route_gdf["eta_min_walk_safe"].iloc[0])
    eta_drive_min = float(route_gdf["eta_min_drive"].iloc[0])
    # folium描画用
    coords = list(route_gdf.geometry.iloc[0].coords)
    route_latlon = [(lat, lon) for lon, lat in coords]
    # 以降の「堅牢チェック→保存」はスキップしてOK（または最新A40で再チェックしたいならそのまま通しても良い）
    use_cached_route = True
else:
    use_cached_route = False

その後のフローはこうしてください：
* use_cached_route == True のとき：
* すでに route_gdf / meters / eta_* / route_latlon が揃っているので、堅牢チェックと保存をスキップして Folium 表示だけ行う（必要なら最新 tsunami_union で再チェックだけ通してもOK）。
* use_cached_route == False のとき：
* これまでどおり nx.shortest_path → 堅牢チェック → ETA計算 → route_latest 保存。
* さらに route_cache に追記します（下の1行を“保存直後”に追加）。

# route_latest を保存した直後にキャッシュへ追記
route_gdf_cache = route_gdf.copy()
route_gdf_cache["rid"] = key
# append（GeoPandas 0.14+）：if_exists="append" が使えます
try:
    route_gdf_cache.to_file(GPKG, layer="route_cache", driver="GPKG", if_exists="append")
except TypeError:
    # 古い版は append がないので自前で結合して上書き
    if "route_cache" in listlayers(GPKG):
        old = gpd.read_file(GPKG, layer="route_cache")
        gpd.GeoDataFrame(pd.concat([old, route_gdf_cache], ignore_index=True), crs=4326)\
           .to_file(GPKG, layer="route_cache", driver="GPKG")
    else:
        route_gdf_cache.to_file(GPKG, layer="route_cache", driver="GPKG")

まとめ（差し込む位置の目安）
1. 目的地 dest_pt 決定直後に「キャッシュ照会ブロック」を挿入
2. use_cached_route が False なら従来の「計算→検証→保存」ブロックへ進む
3. 保存直後に route_cache へ追記
4. True なら保存はスキップし、地図表示へ直行（または再検証のみ実施）
これで「既知の組み合わせは即復帰／未知は新規計算」の運用ができます。

これはモジュールにして読み込んだ方が良い？

結論：はい。モジュール化して読み込むのが吉です。
理由は「再利用・テスト・差し替え（例：県を静岡→他県）」が一気に楽になるから。

貼り付け用の完成ブロック（そのまま置換OK）

# === 検証→ETA→保存（検証に通ったルートだけ保存・表示） ===
from shapely.geometry import LineString
import geopandas as gpd

# 1) ルートのLineString（まだ属性は付けない）
route_lonlat = [(G2.nodes[n]["x"], G2.nodes[n]["y"]) for n in route]  # (lon,lat)
route_geom = LineString(route_lonlat)
route_gdf  = gpd.GeoDataFrame(geometry=[route_geom], crs="EPSG:4326")

# 2) 堅牢チェック（交差ゼロ＆境界近接なし）
#   union_all が使えない環境向けにフォールバックも用意
try:
    flood_union = gpd.GeoDataFrame(geometry=[tsunami.geometry.union_all()], crs="EPSG:4326")
except Exception:
    flood_union = gpd.GeoDataFrame(geometry=[tsunami.unary_union], crs="EPSG:4326")

METRIC = "EPSG:32654"  # 静岡周辺（必要なら地域に合わせて変更）
route_m = route_gdf.to_crs(METRIC)
flood_m = flood_union.to_crs(METRIC)

# 交差長[m]
inter = gpd.overlay(route_m, flood_m, how="intersection")
inside_len_m = float(inter.geometry.length.sum()) if len(inter) else 0.0

# 近接（境界 ±2m 以内をNGにする例）
near = not gpd.overlay(route_m.buffer(2.0), flood_m, how="intersection").empty

if inside_len_m > 0 or near:
    raise RuntimeError(f"ルート不適合: 交差長 {inside_len_m:.2f} m / 近接={near}. "
                       "別避難所 or 半径拡大で再探索してください。")

# 3) 距離[m]→ETA（バージョン互換）
try:
    route_edges = ox.routing.route_to_gdf(G2, route, fill_edge_geometry=True)
    meters = float(route_edges["length"].sum())
except Exception:
    meters = 0.0
    for u, v in zip(route[:-1], route[1:]):
        ed = G2.get_edge_data(u, v)
        meters += min(d.get("length", 0.0) for d in ed.values())

walk_kmh = 4.8
drive_kmh = 20.0
eta_walk_min  = meters / (walk_kmh  * 1000 / 60)
eta_drive_min = meters / (drive_kmh * 1000 / 60)
eta_walk_safe = eta_walk_min * 1.3  # 安全係数（例：+30%）

# 4) GPKG に“合格ルート”だけ属性付きで保存
out_gpkg = "outputs/shizuoka_hazard.gpkg"
route_gdf = gpd.GeoDataFrame(
    {
        "meters":[meters],
        "eta_min_walk":[eta_walk_min],
        "eta_min_walk_safe":[eta_walk_safe],
        "eta_min_drive":[eta_drive_min],
    },
    geometry=[route_geom], crs="EPSG:4326",
)
route_gdf.to_file(out_gpkg, layer="route_latest", driver="GPKG")
print("✅ 交差なし・保存完了")

# 5) （任意）Foliumに表示するならここで
route_latlon = [(G2.nodes[n]["y"], G2.nodes[n]["x"]) for n in route]
folium.PolyLine(
    route_latlon, weight=6, opacity=0.9,
    tooltip=f"距離 {meters:.0f} m / 徒歩 {eta_walk_min:.1f} 分 (安全側 {eta_walk_safe:.1f} 分)"
).add_to(m)

を実行すると

---------------------------------------------------------------------------
NotImplementedError                       Traceback (most recent call last)
Cell In[12], line 26
     23 inside_len_m = float(inter.geometry.length.sum()) if len(inter) else 0.0
     25 # 近接（境界 ±2m 以内をNGにする例）
---> 26 near = not gpd.overlay(route_m.buffer(2.0), flood_m, how="intersection").empty
     28 if inside_len_m > 0 or near:
     29     raise RuntimeError(f"ルート不適合: 交差長 {inside_len_m:.2f} m / 近接={near}. "
     30                        "別避難所 or 半径拡大で再探索してください。")

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/geopandas/tools/overlay.py:273, in overlay(df1, df2, how, keep_geom_type, make_valid)
    270     raise ValueError(f"`how` was '{how}' but is expected to be in {allowed_hows}")
    272 if isinstance(df1, GeoSeries) or isinstance(df2, GeoSeries):
--> 273     raise NotImplementedError(
    274         "overlay currently only implemented for GeoDataFrames"
    275     )
    277 if not _check_crs(df1, df2):
    278     _crs_mismatch_warn(df1, df2, stacklevel=3)

NotImplementedError: overlay currently only implemented for GeoDataFrames

となる

原因を調べたら
route_m.buffer(2.0) は GeoSeries なので、
gpd.overlay()（＝GeoDataFrame同士専用）に渡すとエラーになる

バッファをGeoDataFrame化してから overlay するか、sjoinで近接判定に変える

今後も考慮してGPTに質問する

上位N候補の避難所からネットワーク距離最短を自動選定
10/20分アイソクローンの可視化
現在地が浸水想定内かの警告フラグ
GraphDB用にノード/リレーション出力（Home, Shelter, RoadNode/Edge, TsunamiPolygon, ROUTE_TO など）
も今後行うなら

今回変更するのは
overlay するか、sjoin
どちらが適切？

短く結論：
* 基本は sjoin 推奨（速い・メモリ軽い・属性付与しやすい）。
* 幾何そのものが必要なときだけ overlay（交差部分の長さ/面積を測る、クリップして可視化したい等）。
1. 道路を浸水域で“遮断”（真偽判定・フラグ付け）
* sjoin で十分（predicate=”intersects”）。
* 例：edges[“blocked”] = gpd.sjoin(edges, tsunami[[“geometry”]], predicate=”intersects”, how=”left”)[“index_right”].notna()
2. **近接NG（例：2m以内）**のチェック
* ルートを buffer(2m) → sjoin（軽い）。
* 例：near = not gpd.sjoin(route_m.buffer(2), flood_m, predicate=”intersects”).empty
3. 交差“長さ/面積”が要る（ルートが浸水域内に何 m 入ってる？）
* overlay(…, how=”intersection”) で交差ジオメトリを作り、length/area を合計。
* 今回の「inside_len_m」は overlay が最適。
4. AOIクリップや可視化用の切り出し
* 結果ジオメトリが要るので overlay。
5. 上位N避難所のネットワーク最短選定
* 候補抽出は直線距離（to_crs(UTM)）→ N件
* その後は NetworkX で経路長を比較。
* 候補フィルタ（浸水域内除外）は sjoin でOK。
6. 10/20分アイソクローン
* 等時間ポリゴン生成自体はネットワーク処理（eg. 到達可能ノードの凸包/αシェイプ等）。
* 仕上げの「浸水域で塗り分け」や表示用クリップは overlay。
7. GraphDB 連携（ノード/エッジに属性フラグ）
* IN_FLOOD, NEAR_FLOOD, BLOCKED 等のブール属性付与は sjoin が最適。
* 「交差長」など数値を持たせたい関係は overlay で測って格納。

いまのコードへの具体反映
* 継続：inside_len_m の算出 → overlay（長さが欲しいから）
* 修正済：近接チェック → sjoin（軽い）





buf_gdf = gpd.GeoDataFrame(geometry=route_m.buffer(2.0), crs=METRIC)
near = not gpd.sjoin(buf_gdf, flood_m, how="inner", predicate="intersects").empty

この方針なら、今後の「上位N・アイソクローン・GraphDB」までスムーズに拡張できます。

ということなので
Sjoin　に変更する

# === 検証→ETA→保存（検証に通ったルートだけ保存・表示） ===
from shapely.geometry import LineString
import geopandas as gpd

# 1) ルートのLineString（まだ属性は付けない）
route_lonlat = [(G2.nodes[n]["x"], G2.nodes[n]["y"]) for n in route]  # (lon,lat)
route_geom = LineString(route_lonlat)
route_gdf  = gpd.GeoDataFrame(geometry=[route_geom], crs="EPSG:4326")

# 2) 堅牢チェック（交差ゼロ＆境界近接なし）
#   union_all が使えない環境向けにフォールバックも用意
try:
    flood_union = gpd.GeoDataFrame(geometry=[tsunami.geometry.union_all()], crs="EPSG:4326")
except Exception:
    flood_union = gpd.GeoDataFrame(geometry=[tsunami.unary_union], crs="EPSG:4326")

METRIC = "EPSG:32654"  # 静岡周辺（必要なら地域に合わせて変更）
route_m = route_gdf.to_crs(METRIC)
flood_m = flood_union.to_crs(METRIC)

# 交差長[m]
inter = gpd.overlay(route_m, flood_m, how="intersection")
inside_len_m = float(inter.geometry.length.sum()) if len(inter) else 0.0

# 近接（境界 ±2m 以内をNGにする例）
# near = not gpd.overlay(route_m.buffer(2.0), flood_m, how="intersection").empty
buf_gdf = gpd.GeoDataFrame(geometry=route_m.buffer(2.0), crs=METRIC)
near = not gpd.sjoin(buf_gdf, flood_m, how="inner", predicate="intersects").empty


if inside_len_m > 0 or near:
    raise RuntimeError(f"ルート不適合: 交差長 {inside_len_m:.2f} m / 近接={near}. "
                       "別避難所 or 半径拡大で再探索してください。")

# 3) 距離[m]→ETA（バージョン互換）
try:
    route_edges = ox.routing.route_to_gdf(G2, route, fill_edge_geometry=True)
    meters = float(route_edges["length"].sum())
except Exception:
    meters = 0.0
    for u, v in zip(route[:-1], route[1:]):
        ed = G2.get_edge_data(u, v)
        meters += min(d.get("length", 0.0) for d in ed.values())

walk_kmh = 4.8
drive_kmh = 20.0
eta_walk_min  = meters / (walk_kmh  * 1000 / 60)
eta_drive_min = meters / (drive_kmh * 1000 / 60)
eta_walk_safe = eta_walk_min * 1.3  # 安全係数（例：+30%）

# 4) GPKG に“合格ルート”だけ属性付きで保存
out_gpkg = "outputs/shizuoka_hazard.gpkg"
route_gdf = gpd.GeoDataFrame(
    {
        "meters":[meters],
        "eta_min_walk":[eta_walk_min],
        "eta_min_walk_safe":[eta_walk_safe],
        "eta_min_drive":[eta_drive_min],
    },
    geometry=[route_geom], crs="EPSG:4326",
)
route_gdf.to_file(out_gpkg, layer="route_latest", driver="GPKG")
print("✅ 交差なし・保存完了")

# 5) （任意）Foliumに表示するならここで
route_latlon = [(G2.nodes[n]["y"], G2.nodes[n]["x"]) for n in route]
folium.PolyLine(
    route_latlon, weight=6, opacity=0.9,
    tooltip=f"距離 {meters:.0f} m / 徒歩 {eta_walk_min:.1f} 分 (安全側 {eta_walk_safe:.1f} 分)"
).add_to(m)

というように

# 近接（境界 ±2m 以内をNGにする例）
# near = not gpd.overlay(route_m.buffer(2.0), flood_m, how="intersection").empty
を
buf_gdf = gpd.GeoDataFrame(geometry=route_m.buffer(2.0), crs=METRIC)
near = not gpd.sjoin(buf_gdf, flood_m, how="inner", predicate="intersects").empty

にすることで

✅ 交差なし・保存完了
<folium.vector_layers.PolyLine at 0x31c488b10>

となった

ここまでで「安全確認→距離/ETA→保存→表示」の最小ループは完成

tsunami_union を保存 → 交差チェックをそれに切替
上位N→ネットワーク最短の選定を導入（経路品質UP）
アイソクローンのレイヤをGPKG/地図に追加（“どこまで行けるか”の可視化）
警告フラグをUIに表示
GraphDB CSV出力で家族アプリへの橋渡し

を行っていく

ルート距離・ETA計算 → GPKG保存 → 追加の堅牢化

1) ルート長とETA（徒歩/車）を出す

# ルートの延長[m]は G2 の edge 属性 'length' を合計（OSMnxが付けてくれてます）
meters = sum(ox.utils_graph.get_route_edge_attributes(G2, route, 'length'))
walk_kmh = 4.8   # 歩行速度（必要に応じて変更）
drive_kmh = 20.0 # 災害時の車想定の控えめ値（使うなら）

eta_walk_min  = meters / (walk_kmh  * 1000 / 60)
eta_drive_min = meters / (drive_kmh * 1000 / 60)
print(f"距離: {meters:.0f} m, 徒歩ETA: {eta_walk_min:.1f} 分, 車ETA(仮): {eta_drive_min:.1f} 分")

を実行すると

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In[7], line 2
      1 # ルートの延長[m]は G2 の edge 属性 'length' を合計（OSMnxが付けてくれてます）
----> 2 meters = sum(ox.utils_graph.get_route_edge_attributes(G2, route, 'length'))
      3 walk_kmh = 4.8   # 歩行速度（必要に応じて変更）
      4 drive_kmh = 20.0 # 災害時の車想定の控えめ値（使うなら）

AttributeError: module 'osmnx' has no attribute 'utils_graph'

これは OSMnx のAPI変更が原因です。utils_graph.get_route_edge_attributes は環境によって存在しない

# ▼ 距離[m]の取り方をバージョン互換で
try:
    # 新しめのOSMnx（推奨）: ルートをGDF化して length を合計
    route_edges = ox.routing.route_to_gdf(G2, route, fill_edge_geometry=True)
    meters = float(route_edges["length"].sum())
except Exception:
    # フォールバック: NetworkXからエッジ属性 length を合計
    meters = 0.0
    for u, v in zip(route[:-1], route[1:]):
        ed = G2.get_edge_data(u, v)  # MultiDiGraph: {key: attr_dict}
        if isinstance(ed, dict) and len(ed):
            meters += min(d.get("length", 0.0) for d in ed.values())

# ETA（必要なら速度を調整）
walk_kmh = 4.8
drive_kmh = 20.0
eta_walk_min  = meters / (walk_kmh  * 1000 / 60)
eta_drive_min = meters / (drive_kmh * 1000 / 60)
print(f"距離: {meters:.0f} m, 徒歩ETA: {eta_walk_min:.1f} 分, 車ETA(仮): {eta_drive_min:.1f} 分")

# Folium のツールチップにも表示
route_latlon = [(G2.nodes[n]["y"], G2.nodes[n]["x"]) for n in route]
folium.PolyLine(
    route_latlon, weight=6, opacity=0.9,
    tooltip=f"距離 {meters:.0f} m / 徒歩 {eta_walk_min:.1f} 分"
).add_to(m)

に変更

これで

距離: 3017 m, 徒歩ETA: 37.7 分, 車ETA(仮): 9.1 分
<folium.vector_layers.PolyLine at 0x33d5e6650>

となった

約3.0 km / 徒歩 ~38分 / 車 ~9分（仮） は現実的な値

とりあえずやることは
1) ルート長とETA（徒歩/車）を出す
2) ルートを GPKG に属性付きで保存
3) “ネットワーク的に最も近い”安全避難所を選ぶ（推奨）
いまは直線距離で1件選んでいます。候補を直線距離Top Nに絞ってから実際の経路長で最短を選ぶと堅牢
4) 追加の安全策（おすすめ）
* 自宅が浸水域内かを先に判定：内なら「即時・垂直避難」をUIに強調し、経路探索は参考扱いに。
5) 動作チェックのミニ・リスト
* evac.html を開いて、赤いルートが青の浸水域と交差していないか。
* shelters_safe に目的地が含まれているか（含まれなければCSV側の位置精度を見直し）。
* QGISで roads_safe と tsunami_inundation を重ね、視覚的にもOKか。

まずは
1) ルートをGPKGに属性付きで保存

from shapely.geometry import LineString
import geopandas as gpd

route_lonlat = [(G2.nodes[n]["x"], G2.nodes[n]["y"]) for n in route]
route_gdf = gpd.GeoDataFrame(
    {"meters":[meters], "eta_min_walk":[eta_walk_min], "eta_min_drive":[eta_drive_min]},
    geometry=[LineString(route_lonlat)], crs="EPSG:4326"
)
route_gdf.to_file("outputs/shizuoka_hazard.gpkg", layer="route_latest", driver="GPKG")

2) ルートが浸水域を跨いでないか最終チェック

import geopandas as gpd
ok = not route_gdf.intersects(tsunami.unary_union).iloc[0]
print("交差なし✅" if ok else "⚠️ 交差あり（パラメータ見直し）")

実行結果は

/var/folders/db/tkl5w9dd3kn3h53ctyl4s6180000gn/T/ipykernel_49497/3699053437.py:2: DeprecationWarning: The 'unary_union' attribute is deprecated, use the 'union_all()' method instead.
  ok = not route_gdf.intersects(tsunami.unary_union).iloc[0]

交差なし✅

となった

DeprecationWarning は unary_union が将来的に非推奨という通知だけなので、動作上は問題ありません。念のため 推奨APIに差し替える

import geopandas as gpd
# ok = not route_gdf.intersects(tsunami.unary_union).iloc[0]
# print("交差なし✅" if ok else "⚠️ 交差あり（パラメータ見直し）")

# 旧: tsunami.unary_union  ← 警告が出る
# 新:
flood_union = tsunami.geometry.union_all()
ok = not route_gdf.intersects(flood_union).iloc[0]
print("交差なし✅" if ok else "⚠️ 交差あり（パラメータ見直し）")

これで警告は消える

これを修正したらGPTから
もう一段だけ堅牢に（任意）の項目が提案された
この部分は
3) 安全バッファ付きETA（渋滞・混雑を考慮）
4) Foliumツールチップに表示
の前に行う処理？
と聞いてみる

結論：はい、“堅牢チェック”は 3)安全バッファ付きETA・4)ツールチップ表示の前 にやるのがベストです。
理由：交差や境界ニアミスが見つかったらその経路は捨てて再探索するので、ETAや表示は検証済みルートに対してだけ出すのが安全

おすすめの処理順
1. ルート算出（orig→dest）
2. 堅牢チェック
* union_all() で津波ポリゴン結合
* 交差の有無 → あれば再探索（別避難所 or 半径拡大）
* 交差長[m]（=0が理想）
* 近接判定（例：ルート1–5mバッファが浸水域に当たるか）
3. （ここで初めて）距離[m]→ETA計算（徒歩/車＋安全係数）
4. Folium描画（ツールチップに距離/ETAを表示）
5. GPKGへ route_latest などを保存

チェック（交差・近接の堅牢判定）→ OKなら保存の順がベスト
不適切なルートをGPKGに残さないため。保存は“合格”ルートだけにしましょう。
とのこと

道路→浸水域回避→最寄り安全避難所への経路の作成

import pandas as pd
import geopandas as gpd
from shapely.geometry import Point, LineString
import osmnx as ox
import networkx as nx
import folium
import os

GPKG = "outputs/shizuoka_hazard.gpkg"

# 1) 既存レイヤ読み込み（両方とも 4326）
coast   = gpd.read_file(GPKG, layer="coastline")
tsunami = gpd.read_file(GPKG, layer="tsunami_inundation")

# ----- パラメータ -----
# 自宅座標（lat, lon）←必ず置き換え
home_wgs = (34.728, 137.978)
# OSM道路の取得半径[m]（到達不能なら大きく：15000, 20000…）
GRAPH_RADIUS_M = 12000
# 避難所CSV（cp932想定）
SHELTER_CSV = "data/20250801_指定緊急避難場所一覧.CSV"
# ----------------------

# 2) 避難所CSV→点化（列名は自動検出）
df = pd.read_csv(SHELTER_CSV, encoding="cp932")
lat_candidates = 
lon_candidates = 
if not lat_candidates or not lon_candidates:
    raise ValueError("CSVに緯度/経度列が見つかりません。列名を教えてください。")
lat_col, lon_col = lat_candidates[0], lon_candidates[0]
df[lat_col] = pd.to_numeric(df[lat_col], errors="coerce")
df[lon_col] = pd.to_numeric(df[lon_col], errors="coerce")
df = df.dropna(subset=[lat_col, lon_col])
shelters = gpd.GeoDataFrame(df.copy(),
    geometry=gpd.points_from_xy(df[lon_col], df[lat_col], crs="EPSG:4326")
)

# 3) 安全な避難所 = 浸水域と非交差
#   （unionはメモリ負荷になるので、空間結合で判定）
safe_join = gpd.sjoin(shelters, tsunami[["geometry"]], how="left", predicate="intersects")
shelters["safe"] = safe_join["index_right"].isna()
safe_shelters = shelters[shelters["safe"]].copy()

# 4) OSM道路取得（徒歩網）
G = ox.graph_from_point(home_wgs, dist=GRAPH_RADIUS_M, network_type="walk", simplify=True, retain_all=False)
nodes, edges = ox.graph_to_gdfs(G, nodes=True, edges=True)
if edges.crs is None: edges.set_crs(4326, inplace=True)
elif edges.crs.to_epsg() != 4326: edges = edges.to_crs(4326)

# 5) 浸水域と交差する道路を除去（空間結合で高速化）
tag = gpd.sjoin(edges[["u","v","key","geometry"]], tsunami[["geometry"]],
                how="left", predicate="intersects")
tag["blocked"] = tag["index_right"].notna()
edges_keep = edges.loc[~tag["blocked"].values].copy()

# 6) サブグラフ再構築（通行可のみ）
used_nodes = pd.unique(edges_keep[["u","v"]].values.ravel())
nodes_keep = nodes.loc[nodes.index.isin(used_nodes)].copy()
G2 = ox.utils_graph.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=G.graph)

# 7) 出発・目的地ノード（最寄り“安全”避難所）
home_pt = Point(home_wgs[1], home_wgs[0])
# 最寄り安全避難所（直線距離）…なければ全体から選ぶ
TARGET = safe_shelters if len(safe_shelters) else shelters
METRIC = "EPSG:32654"  # 静岡周辺 UTM
hp = gpd.GeoSeries([home_pt], crs=4326).to_crs(METRIC).iloc[0]
tp = TARGET.to_crs(METRIC).copy()
tp["__dist_m__"] = tp.geometry.distance(hp)
dest_pt = TARGET.loc[tp["__dist_m__"].idxmin(), "geometry"]

orig = ox.distance.nearest_nodes(G2, X=home_pt.x, Y=home_pt.y)
dest = ox.distance.nearest_nodes(G2, X=dest_pt.x, Y=dest_pt.y)

# 8) 最短経路（通行可能グラフ）
try:
    route = nx.shortest_path(G2, orig, dest, weight="length")
except nx.NetworkXNoPath:
    raise RuntimeError("到達経路なし。GRAPH_RADIUS_Mを増やすか、別の避難所で再試行を。")

# 9) GPKGへ追記
#    既にある場合は上書きしたいなら driverオプションを切り替える（ここは追記）
shelters.to_file(GPKG, layer="shelters_all", driver="GPKG")
safe_shelters.to_file(GPKG, layer="shelters_safe", driver="GPKG")
edges.to_file(GPKG, layer="roads_all", driver="GPKG")
edges_keep.to_file(GPKG, layer="roads_safe", driver="GPKG")

# 10) Folium地図
m = folium.Map(location=[home_wgs[0], home_wgs[1]], zoom_start=12, control_scale=True)
folium.GeoJson(tsunami.__geo_interface__, name="津波浸水想定(A40)").add_to(m)
folium.GeoJson(edges_keep.__geo_interface__, name="道路(浸水回避)").add_to(m)

# 経路を線で描画
route_latlon = [(G2.nodes[n]["y"], G2.nodes[n]["x"]) for n in route]
folium.PolyLine(route_latlon, weight=6, opacity=0.9, tooltip="避難経路").add_to(m)

# 出発・目的地
folium.Marker([home_wgs[0], home_wgs[1]], popup="現在地", icon=folium.Icon(color="green")).add_to(m)
folium.Marker([dest_pt.y, dest_pt.x], popup="避難所", icon=folium.Icon(color="blue")).add_to(m)

# 避難所レイヤ（安全のみ強調）
folium.GeoJson(safe_shelters.__geo_interface__, name="安全な避難所",
               marker=folium.CircleMarker(radius=4)).add_to(m)
folium.LayerControl().add_to(m)

# 保存（任意）
m.save("outputs/evac.html")
print("✓ 完了: GPKG追記 & outputs/evac.html を作成")

実行すると

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
Cell In[1], line 51
     48 elif edges.crs.to_epsg() != 4326: edges = edges.to_crs(4326)
     50 # 5) 浸水域と交差する道路を除去（空間結合で高速化）
---> 51 tag = gpd.sjoin(edges[["u","v","key","geometry"]], tsunami[["geometry"]],
     52                 how="left", predicate="intersects")
     53 tag["blocked"] = tag["index_right"].notna()
     54 edges_keep = edges.loc[~tag["blocked"].values].copy()

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/geopandas/geodataframe.py:1896, in GeoDataFrame.__getitem__(self, key)
   1890 def __getitem__(self, key):
   1891     """
   1892     If the result is a column containing only 'geometry', return a
   1893     GeoSeries. If it's a DataFrame with any columns of GeometryDtype,
   1894     return a GeoDataFrame.
   1895     """
-> 1896     result = super().__getitem__(key)
   1897     # Custom logic to avoid waiting for pandas GH51895
   1898     # result is not geometry dtype for multi-indexes
   1899     if (
   1900         pd.api.types.is_scalar(key)
   1901         and key == ""
   (...)   1904         and not is_geometry_type(result)
   1905     ):

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pandas/core/frame.py:4108, in DataFrame.__getitem__(self, key)
   4106     if is_iterator(key):
   4107         key = list(key)
-> 4108     indexer = self.columns._get_indexer_strict(key, "columns")[1]
   4110 # take() does not accept boolean indexers
   4111 if getattr(indexer, "dtype", None) == bool:

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pandas/core/indexes/base.py:6200, in Index._get_indexer_strict(self, key, axis_name)
   6197 else:
   6198     keyarr, indexer, new_indexer = self._reindex_non_unique(keyarr)
-> 6200 self._raise_if_missing(keyarr, indexer, axis_name)
   6202 keyarr = self.take(indexer)
   6203 if isinstance(key, Index):
   6204     # GH 42790 - Preserve name from an Index

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pandas/core/indexes/base.py:6252, in Index._raise_if_missing(self, key, indexer, axis_name)
   6249     raise KeyError(f"None of [{key}] are in the [{axis_name}]")
   6251 not_found = list(ensure_index(key)[missing_mask.nonzero()[0]].unique())
-> 6252 raise KeyError(f"{not_found} not in index")

KeyError: "['u', 'v', 'key'] not in index"

原因はコレ：graph_to_gdfs の edges は MultiIndex（u,v,key が“列ではなくインデックス”） なので、edges[[“u”,”v”,”key”,”geometry”]] が KeyError になります。
→ reset_index() で列に起こす＋sjoin の重複行を集計してからフィルタしてください

ということで

# --- 4) OSM道路取得（徒歩網） ---
G = ox.graph_from_point(home_wgs, dist=GRAPH_RADIUS_M, network_type="walk", simplify=True, retain_all=False)
nodes, edges = ox.graph_to_gdfs(G, nodes=True, edges=True)

# u,v,key を“列”にする（重要）
edges = edges.reset_index()

# CRSを4326に
if edges.crs is None:
    edges = edges.set_crs(4326)
elif edges.crs.to_epsg() != 4326:
    edges = edges.to_crs(4326)

# --- 5) 浸水域と交差する道路を除去（sjoin→集計で一意化） ---
# 属性を絞って結合
left = edges[["u","v","key","geometry"]].copy()
right = tsunami[["geometry"]]

tag = gpd.sjoin(left, right, how="left", predicate="intersects")

# エッジごとに「1回でも当たったら blocked=True」
blocked = (
    tag.assign(hit=tag["index_right"].notna())
       .groupby(["u","v","key"], as_index=False)["hit"].max()
       .rename(columns={"hit":"blocked"})
)

# 元の edges に blocked を付与（ない行は False）
edges2 = edges.merge(blocked, on=["u","v","key"], how="left")
edges2["blocked"] = edges2["blocked"].fillna(False)

edges_keep = edges2.loc[~edges2["blocked"]].copy()

# --- 6) サブグラフ再構築 ---
used_nodes = pd.unique(edges_keep[["u","v"]].values.ravel())
nodes_keep = nodes.loc[nodes.index.isin(used_nodes)].copy()
G2 = ox.utils_graph.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=G.graph)

へコード変更

しかしまだ

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In[2], line 80
     78 used_nodes = pd.unique(edges_keep[["u","v"]].values.ravel())
     79 nodes_keep = nodes.loc[nodes.index.isin(used_nodes)].copy()
---> 80 G2 = ox.utils_graph.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=G.graph)
     83 # 7) 出発・目的地ノード（最寄り“安全”避難所）
     84 home_pt = Point(home_wgs[1], home_wgs[0])

AttributeError: module 'osmnx' has no attribute 'utils_graph'

となる

これは
OSMnx のバージョン差です。
utils_graph.graph_from_gdfs は古い呼び方で、いまは ox.graph_from_gdfs（トップレベル）を使います
とのことなので

# --- 6) サブグラフ再構築 ---
used_nodes = pd.unique(edges_keep[["u","v"]].values.ravel())
nodes_keep = nodes.loc[nodes.index.isin(used_nodes)].copy()
G2 = ox.utils_graph.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=G.graph)

を

# --- 6) サブグラフ再構築（通行可のみ） ---
used_nodes = pd.unique(edges_keep[["u","v"]].values.ravel())
nodes_keep = nodes.loc[nodes.index.isin(used_nodes)].copy()

# ▶ ここを置き換え
try:
    # 新しめの OSMnx（1.1+ 〜 2.x）
    G2 = ox.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=getattr(G, "graph", None))
except AttributeError:
    # 旧版 OSMnx 向けフォールバック
    from osmnx import utils_graph as ug
    G2 = ug.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=getattr(G, "graph", None))

に変更

最寄りノード取得もバージョン差があるので互換にしておく

home_pt = Point(home_wgs[1], home_wgs[0])

TARGET = safe_shelters if len(safe_shelters) else shelters
METRIC = "EPSG:32654"  
hp = gpd.GeoSeries([home_pt], crs=4326).to_crs(METRIC).iloc[0]
tp = TARGET.to_crs(METRIC).copy()
tp["__dist_m__"] = tp.geometry.distance(hp)
dest_pt = TARGET.loc[tp["__dist_m__"].idxmin(), "geometry"]

orig = ox.distance.nearest_nodes(G2, X=home_pt.x, Y=home_pt.y)
dest = ox.distance.nearest_nodes(G2, X=dest_pt.x, Y=dest_pt.y)

を

# 7) 出発・目的地ノード
home_pt = Point(home_wgs[1], home_wgs[0])
TARGET = safe_shelters if len(safe_shelters) else shelters
METRIC = "EPSG:32654"
hp = gpd.GeoSeries([home_pt], crs=4326).to_crs(METRIC).iloc[0]
tp = TARGET.to_crs(METRIC).copy()
tp["__dist_m__"] = tp.geometry.distance(hp)
dest_pt = TARGET.loc[tp["__dist_m__"].idxmin(), "geometry"]

# ▶ 互換呼び分け
try:
    orig = ox.distance.nearest_nodes(G2, X=home_pt.x, Y=home_pt.y)
    dest = ox.distance.nearest_nodes(G2, X=dest_pt.x, Y=dest_pt.y)
except AttributeError:
    # 新API名（環境によってはこちらが有効）
    orig = ox.nearest_nodes(G2, X=home_pt.x, Y=home_pt.y)
    dest = ox.nearest_nodes(G2, X=dest_pt.x, Y=dest_pt.y)

に変更

OSMnx のバージョンは import osmnx as ox; print(ox.__version__) で確認できます。
もし整えたいなら

 pip install -U osmnx 

で最新版に寄せるのもアリ
とのこと

--------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
Cell In[3], line 87
     84 # ▶ ここを置き換え
     85 try:
     86     # 新しめの OSMnx（1.1+ 〜 2.x）
---> 87     G2 = ox.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=getattr(G, "graph", None))
     88 except AttributeError:
     89     # 旧版 OSMnx 向けフォールバック
     90     from osmnx import utils_graph as ug

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/osmnx/convert.py:302, in graph_from_gdfs(gdf_nodes, gdf_edges, graph_attrs)
    265 def graph_from_gdfs(
    266     gdf_nodes: gpd.GeoDataFrame,
    267     gdf_edges: gpd.GeoDataFrame,
    268     *,
    269     graph_attrs: dict[str, Any] | None = None,
    270 ) -> nx.MultiDiGraph:
    271     """
    272     Convert node and edge GeoDataFrames to a MultiDiGraph.
    273 
   (...)    300         The converted MultiDiGraph.
    301     """
--> 302     _validate_node_edge_gdfs(gdf_nodes, gdf_edges)
    304     # drop geometry column from gdf_nodes (since we use x and y for geometry
    305     # information), but warn the user if the geometry values differ from the
    306     # coordinates in the x and y columns. this results in a df instead of gdf.
    307     if gdf_nodes.active_geometry_name is None:  # pragma: no cover

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/osmnx/convert.py:241, in _validate_node_edge_gdfs(gdf_nodes, gdf_edges)
    239 edges_index_levels = 3
    240 check1 = gdf_edges.index.nlevels == edges_index_levels
--> 241 uv = set(gdf_edges.index.get_level_values(0)) | set(gdf_edges.index.get_level_values(1))
    242 check2 = uv.issubset(set(gdf_nodes.index))
    243 if not (check1 and check2):  # pragma: no cover

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pandas/core/indexes/base.py:2102, in Index._get_level_values(self, level)
   2066 def _get_level_values(self, level) -> Index:
   2067     """
   2068     Return an Index of values for requested level.
   2069 
   (...)   2100     Index(['a', 'b', 'c'], dtype='object')
   2101     """
-> 2102     self._validate_index_level(level)
   2103     return self

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pandas/core/indexes/base.py:2008, in Index._validate_index_level(self, level)
   2003         raise IndexError(
   2004             "Too many levels: Index has only 1 level, "
   2005             f"{level} is not a valid level number"
   2006         )
   2007     if level > 0:
-> 2008         raise IndexError(
   2009             f"Too many levels: Index has only 1 level, not {level + 1}"
   2010         )
   2011 elif level != self.name:
   2012     raise KeyError(
   2013         f"Requested level ({level}) does not match index name ({self.name})"
   2014     )

IndexError: Too many levels: Index has only 1 level, not 2


原因は、ox.graph_from_gdfs が受け取る edges のインデックス形式です。
この関数は edges のインデックスが MultiIndex（u, v, key） であることを前提にしています。先で reset_index() したので RangeIndex になっており、検証でコケています

# --- 5) 浸水域と交差する道路を除去（sjoin→集計で一意化） ---
left = edges[["u","v","key","geometry"]].copy()   # edges は reset_index() 済み想定
right = tsunami[["geometry"]]

tag = gpd.sjoin(left, right, how="left", predicate="intersects")

blocked = (
    tag.assign(hit=tag["index_right"].notna())
       .groupby(["u","v","key"], as_index=False)["hit"].max()
       .rename(columns={"hit":"blocked"})
)

edges2 = edges.merge(blocked, on=["u","v","key"], how="left")
edges2["blocked"] = edges2["blocked"].fillna(False)

edges_keep = edges2.loc[~edges2["blocked"]].copy()

# ★ ここが重要：graph_from_gdfs に渡す前に MultiIndex に戻す
#   （key が float になっていたら int にしておくと安心）
if edges_keep["key"].isna().any():
    edges_keep["key"] = edges_keep["key"].fillna(0)

edges_keep["key"] = edges_keep["key"].astype(int)
edges_keep = edges_keep.set_index(["u","v","key"]).sort_index()

# --- 6) サブグラフ再構築（通行可のみ） ---
used_nodes = pd.unique(edges_keep.reset_index()[["u","v"]].values.ravel())
nodes_keep = nodes.loc[nodes.index.isin(used_nodes)].copy()

# nodes_keep の geometry 列が欠けていたら補完
if "geometry" not in nodes_keep or nodes_keep["geometry"].isna().any():
    nodes_keep = nodes_keep.set_geometry(
        gpd.points_from_xy(nodes_keep["x"], nodes_keep["y"]), crs=edges_keep.crs
    )

# OSMnx のバージョン差に対応
try:
    G2 = ox.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=getattr(G, "graph", None))
except AttributeError:
    from osmnx import utils_graph as ug
    G2 = ug.graph_from_gdfs(nodes_keep, edges_keep, graph_attrs=getattr(G, "graph", None))

これで

✓ 完了: GPKG追記 & outputs/evac.html を作成

となって

print("edges_keep index levels:", getattr(edges_keep.index, "nlevels", 1))  # → 3 が正解
print(edges_keep.index.names)  # → ['u','v','key'] が正解

の結果は

edges_keep index levels: 3
['u', 'v', 'key']

となる

次に
“差分だけ”足せば、ルート距離・ETA計算 → GPKG保存 → 追加の堅牢化まで整う