コンセプトに基づく改良

コンセプトに基づく改良

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)

# 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()   # 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))


# 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"]

# ▶ 互換呼び分け
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)


# 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 を作成")

の結果

---------------------------------------------------------------------------
CPLE_AppDefinedError                      Traceback (most recent call last)
File pyogrio/_io.pyx:2344, in pyogrio._io.create_ogr_dataset_layer()

File pyogrio/_err.pyx:218, in pyogrio._err.check_pointer()

CPLE_AppDefinedError: Layer roads_all already exists, CreateLayer failed. Use the layer creation option OVERWRITE=YES to replace it.

During handling of the above exception, another exception occurred:

DataLayerError                            Traceback (most recent call last)
Cell In[2], line 129
    127 shelters.to_file(GPKG, layer="shelters_all", driver="GPKG")
    128 safe_shelters.to_file(GPKG, layer="shelters_safe", driver="GPKG")
--> 129 edges.to_file(GPKG, layer="roads_all", driver="GPKG")
    130 edges_keep.to_file(GPKG, layer="roads_safe", driver="GPKG")
    132 # 10) Folium地図

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/geopandas/geodataframe.py:1637, in GeoDataFrame.to_file(self, filename, driver, schema, index, **kwargs)


   1543 """Write the ``GeoDataFrame`` to a file.
   1544 
   1545 By default, an ESRI shapefile is written, but any OGR data source
   (...)   1633 
   1634 """
   1635 from geopandas.io.file import _to_file
-> 1637 _to_file(self, filename, driver, schema, index, **kwargs)

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/geopandas/io/file.py:731, in _to_file(df, filename, driver, schema, index, mode, crs, engine, metadata, **kwargs)
    728     raise ValueError(f"'mode' should be one of 'w' or 'a', got '{mode}' instead")
    730 if engine == "pyogrio":
--> 731     _to_file_pyogrio(df, filename, driver, schema, crs, mode, metadata, **kwargs)
    732 elif engine == "fiona":
    733     _to_file_fiona(df, filename, driver, schema, crs, mode, metadata, **kwargs)

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/geopandas/io/file.py:793, in _to_file_pyogrio(df, filename, driver, schema, crs, mode, metadata, **kwargs)
    790 if not df.columns.is_unique:
    791     raise ValueError("GeoDataFrame cannot contain duplicated column names.")
--> 793 pyogrio.write_dataframe(df, filename, driver=driver, metadata=metadata, **kwargs)

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pyogrio/geopandas.py:710, in write_dataframe(df, path, layer, driver, encoding, geometry_type, promote_to_multi, nan_as_null, append, use_arrow, dataset_metadata, layer_metadata, metadata, dataset_options, layer_options, **kwargs)
    707         field_data.append(values)
    708         field_mask.append(None)
--> 710 write(
    711     path,
    712     layer=layer,
    713     driver=driver,
    714     geometry=geometry,
    715     field_data=field_data,
    716     field_mask=field_mask,
    717     fields=fields,
    718     crs=crs,
    719     geometry_type=geometry_type,
    720     encoding=encoding,
    721     promote_to_multi=promote_to_multi,
    722     nan_as_null=nan_as_null,
    723     append=append,
    724     dataset_metadata=dataset_metadata,
    725     layer_metadata=layer_metadata,
    726     metadata=metadata,
    727     dataset_options=dataset_options,
    728     layer_options=layer_options,
    729     gdal_tz_offsets=gdal_tz_offsets,
    730     **kwargs,
    731 )

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pyogrio/raw.py:723, in write(path, geometry, field_data, fields, field_mask, layer, driver, geometry_type, crs, encoding, promote_to_multi, nan_as_null, append, dataset_metadata, layer_metadata, metadata, dataset_options, layer_options, gdal_tz_offsets, **kwargs)
    718 # preprocess kwargs and split in dataset and layer creation options
    719 dataset_kwargs, layer_kwargs = _preprocess_options_kwargs(
    720     driver, dataset_options, layer_options, kwargs
    721 )
--> 723 ogr_write(
    724     path,
    725     layer=layer,
    726     driver=driver,
    727     geometry=geometry,
    728     geometry_type=geometry_type,
    729     field_data=field_data,
    730     field_mask=field_mask,
    731     fields=fields,
    732     crs=crs,
    733     encoding=encoding,
    734     promote_to_multi=promote_to_multi,
    735     nan_as_null=nan_as_null,
    736     append=append,
    737     dataset_metadata=dataset_metadata,
    738     layer_metadata=layer_metadata,
    739     dataset_kwargs=dataset_kwargs,
    740     layer_kwargs=layer_kwargs,
    741     gdal_tz_offsets=gdal_tz_offsets,
    742 )

File pyogrio/_io.pyx:2463, in pyogrio._io.ogr_write()

File pyogrio/_io.pyx:2360, in pyogrio._io.create_ogr_dataset_layer()

DataLayerError: Layer roads_all already exists, CreateLayer failed. Use the layer creation option OVERWRITE=YES to replace it.

となった

原因は
エラーの原因は、そのレイヤー名(roads_all)が GPKG 内にすでに存在しているのに、上書きオプションを付けずに to_file() しているためです。対処は「上書き」するか「追記」するかの2択です。再実行に強いのは上書きです。

いくつか方法を提示されたが、コンセプトに基づき改造しないと後で全部修正になる

使用する用途のコンセプト また震源地から津波が到達するまでの時間と津波の高さを計算し GIS情報を使い避難するエリアを算出しルート案内できるようにしたい 家族の避難経路・地域防災マップをGraphDBで整理 避難所や重要連絡先も紐付け 使い方例: 「地震が来たら最寄りの避難所はどこ?」 「非常食の賞味期限が近いものは?」 など この用途で使うなら 今回の場合は (上書き) 既存レイヤーを削除してから書く 同じレイヤーにレコードを追加 どれを推奨するか?

結論から:
* 今回の用途(津波到達時間・高さを計算 → 安全エリア算出 → その時点で最良の避難経路を提示)では、
“計算で毎回作り直す派生レイヤー”は原則《上書き》を推奨です。
* 履歴を残したいテーブル(実行ログ・過去の経路・賞味期限チェック結果など)は《追記》を併用します。
* 《既存レイヤー削除→書き込み》は、スキーマを頻繁に変えるときやエンジン差の不具合回避用のセーフティとして使う位置づけ。

なぜ《上書き》が最適か
あなたのワークフローは「同じ入力(ベースの津波想定・OSM・避難所)→同じ処理→決定論的な派生成果物(roads_safe, shelters_safe, route_latest など)」という再現可能・冪等パイプラインです。
このタイプは毎回クリーンに作り直せることが最重要なので、レイヤー名固定+上書きが最もシンプルで壊れにくいです。
* 冪等:再実行で結果が混ざらない
* 運用:下流(Folium表示/GraphDB同期)が常に同じレイヤー名を参照
* 事故防止:前回分の残骸や列増減の不一致で壊れない

まとめ(今回どれを勧める?)
* メインは《上書き》一択(派生レイヤーは常に最新をクリーンに再生成)
* 履歴が必要なものだけ《追記》を別レイヤーで併用
* スキーマ変更直後や不具合回避には《削除→書き込み》をスポット的に
この方針にすると、
「地震が来たら最寄りの避難所は?」→ shelters_safe+route_latest を即参照、
「非常食の賞味期限が近いものは?」→ food_expiry_alerts(追記ログ)や在庫DB を参照、
という形で運用がシンプル&堅牢になります。

# 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")

# 9) GPKGへ書き出し(毎回クリーンに上書き)
common_kwargs = dict(
    driver="GPKG",
    engine="pyogrio",
    layer_options={"OVERWRITE": "YES"},  # 既存レイヤーを置き換える
)

shelters.to_file(GPKG, layer="shelters_all",   **common_kwargs)
safe_shelters.to_file(GPKG, layer="shelters_safe", **common_kwargs)
edges.to_file(GPKG,    layer="roads_all",     **common_kwargs)
edges_keep.to_file(GPKG, layer="roads_safe",  **common_kwargs)

というようにコードを変更

これでエラーは解決

しかし

# # === 検証→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)

# === ルート検証 → ETA → 保存(キャッシュ対応) ===
from shapely.geometry import LineString, Point
import geopandas as gpd
import hashlib, json
from datetime import datetime
from fiona import listlayers
import osmnx as ox
import networkx as nx

GPKG = "outputs/shizuoka_hazard.gpkg"
METRIC = "EPSG:32654"  # 静岡周辺

# 0) 津波 union はキャッシュを優先(なければ作る)
if "tsunami_union" in listlayers(GPKG):
    flood_union = gpd.read_file(GPKG, layer="tsunami_union")
else:
    flood_union = gpd.GeoDataFrame(geometry=[tsunami.geometry.union_all()], crs=tsunami.crs)
    flood_union.to_file(GPKG, layer="tsunami_union", driver="GPKG")

# 0') ハザードの簡易メタ(A40が変わったらルートキャッシュ無効化に使う)
ti = gpd.read_file(GPKG, layer="tsunami_inundation")
hazard_meta = f'{len(ti)}|{",".join(map(str, ti.total_bounds))}'

# 1) キャッシュキー生成(端点・モード・半径・ハザードメタ)
def route_key(home_wgs, dest_wgs, mode="walk", radius_m=12000, 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()

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

# 2) ルートキャッシュにヒットすれば使う
route_gdf = None
use_cached_route = False
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 = hit.set_crs(4326)
        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])
        coords = list(route_gdf.geometry.iloc[0].coords)
        route_latlon = [(lat, lon) for lon, lat in coords]  # folium用 (lat,lon)
        use_cached_route = True
        print("↩︎ ルート: キャッシュヒット")

if not use_cached_route:
    # 3) これまで通りの経路('route')から 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")

    # 4) 堅牢チェック(交差ゼロ&境界近接なし)
    route_m = route_gdf.to_crs(METRIC)
    flood_m = flood_union.to_crs(METRIC)

    inter = gpd.overlay(route_m, flood_m, how="intersection")
    inside_len_m = float(inter.geometry.length.sum()) if len(inter) else 0.0

    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 半径拡大で再探索してください。")

    # 5) 距離[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

    # 6) route_latest へ保存(検証済みのみ)
    out_gpkg = 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],
            "created_at":[datetime.utcnow().isoformat(timespec="seconds")+"Z"],
        },
        geometry=[route_geom], crs="EPSG:4326",
    )
    route_gdf.to_file(out_gpkg, layer="route_latest", driver="GPKG")
    print("✅ 交差なし・保存完了")

    # 7) ルートキャッシュに追記
    route_cache_rec = route_gdf.copy()
    route_cache_rec["rid"] = key
    try:
        route_cache_rec.to_file(out_gpkg, layer="route_cache", driver="GPKG", if_exists="append")
    except TypeError:
        # 古いGeoPandasのためのフォールバック
        if "route_cache" in listlayers(out_gpkg):
            old = gpd.read_file(out_gpkg, layer="route_cache")
            gpd.GeoDataFrame(pd.concat([old, route_cache_rec], ignore_index=True), crs=4326)\
               .to_file(out_gpkg, layer="route_cache", driver="GPKG")
        else:
            route_cache_rec.to_file(out_gpkg, layer="route_cache", driver="GPKG")

    # 8) folium用座標
    route_latlon = [(G2.nodes[n]["y"], G2.nodes[n]["x"]) for n in route]

# 9) Folium表示(キャッシュ・新規どちらでも)
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)

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[8], line 175
    173 route_cache_rec["rid"] = key
    174 try:
--> 175     route_cache_rec.to_file(out_gpkg, layer="route_cache", driver="GPKG", if_exists="append")
    176 except TypeError:
    177     # 古いGeoPandasのためのフォールバック
    178     if "route_cache" in listlayers(out_gpkg):

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/geopandas/geodataframe.py:1637, in GeoDataFrame.to_file(self, filename, driver, schema, index, **kwargs)
   1543 """Write the ``GeoDataFrame`` to a file.
   1544 
   1545 By default, an ESRI shapefile is written, but any OGR data source
   (...)   1633 
   1634 """
   1635 from geopandas.io.file import _to_file
-> 1637 _to_file(self, filename, driver, schema, index, **kwargs)

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/geopandas/io/file.py:731, in _to_file(df, filename, driver, schema, index, mode, crs, engine, metadata, **kwargs)
    728     raise ValueError(f"'mode' should be one of 'w' or 'a', got '{mode}' instead")
    730 if engine == "pyogrio":
--> 731     _to_file_pyogrio(df, filename, driver, schema, crs, mode, metadata, **kwargs)
    732 elif engine == "fiona":
    733     _to_file_fiona(df, filename, driver, schema, crs, mode, metadata, **kwargs)

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/geopandas/io/file.py:793, in _to_file_pyogrio(df, filename, driver, schema, crs, mode, metadata, **kwargs)
    790 if not df.columns.is_unique:
    791     raise ValueError("GeoDataFrame cannot contain duplicated column names.")
--> 793 pyogrio.write_dataframe(df, filename, driver=driver, metadata=metadata, **kwargs)

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pyogrio/geopandas.py:710, in write_dataframe(df, path, layer, driver, encoding, geometry_type, promote_to_multi, nan_as_null, append, use_arrow, dataset_metadata, layer_metadata, metadata, dataset_options, layer_options, **kwargs)
    707         field_data.append(values)
    708         field_mask.append(None)
--> 710 write(
    711     path,
    712     layer=layer,
    713     driver=driver,
    714     geometry=geometry,
    715     field_data=field_data,
    716     field_mask=field_mask,
    717     fields=fields,
    718     crs=crs,
    719     geometry_type=geometry_type,
    720     encoding=encoding,
    721     promote_to_multi=promote_to_multi,
    722     nan_as_null=nan_as_null,
    723     append=append,
    724     dataset_metadata=dataset_metadata,
    725     layer_metadata=layer_metadata,
    726     metadata=metadata,
    727     dataset_options=dataset_options,
    728     layer_options=layer_options,
    729     gdal_tz_offsets=gdal_tz_offsets,
    730     **kwargs,
    731 )

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pyogrio/raw.py:719, in write(path, geometry, field_data, fields, field_mask, layer, driver, geometry_type, crs, encoding, promote_to_multi, nan_as_null, append, dataset_metadata, layer_metadata, metadata, dataset_options, layer_options, gdal_tz_offsets, **kwargs)
    711     warnings.warn(
    712         "'crs' was not provided.  The output dataset will not have "
    713         "projection information defined and may not be usable in other "
    714         "systems.",
    715         stacklevel=2,
    716     )
    718 # preprocess kwargs and split in dataset and layer creation options
--> 719 dataset_kwargs, layer_kwargs = _preprocess_options_kwargs(
    720     driver, dataset_options, layer_options, kwargs
    721 )
    723 ogr_write(
    724     path,
    725     layer=layer,
   (...)    741     gdal_tz_offsets=gdal_tz_offsets,
    742 )

File ~/.pyenv/versions/3.11.0/lib/python3.11/site-packages/pyogrio/raw.py:514, in _preprocess_options_kwargs(driver, dataset_options, layer_options, kwargs)
    512             layer_kwargs[k] = v
    513         else:
--> 514             raise ValueError(f"unrecognized option '{k}' for driver '{driver}'")
    516 return dataset_kwargs, layer_kwargs

ValueError: unrecognized option 'IF_EXISTS' for driver 'GPKG'

原因は
if_exists=”append” は pyogrio では非対応 なのでエラーになっています。
to_file() は engine=”pyogrio” + append=True を使ってください。ついでに route_latest など「毎回作り直す」レイヤーは OVERWRITE で上書きするのが安全です

要点
* pyogrio では if_exists=”append” は使えません → append=True を使う。
* 毎回再計算するレイヤーは layer_options={“OVERWRITE”:”YES”} で確実に上書き。
* 追記でスキーマ差が出たら、読み出し→列合わせ→連結→上書き で堅牢化。

route_latest は毎回クリーンに上書き、
route_cache はappendで追記(pyogrio準拠)、
tsunami_union も必要に応じて上書き耐性で作成、
の3点が組み込んでコードを書き換える

# # === 検証→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)

# === ルート検証 → ETA → 保存(キャッシュ対応) ===
from shapely.geometry import LineString, Point
import geopandas as gpd
import hashlib, json
from datetime import datetime
from fiona import listlayers
import osmnx as ox
import networkx as nx

GPKG = "outputs/shizuoka_hazard.gpkg"
METRIC = "EPSG:32654"  # 静岡周辺

# 0) 津波 union はキャッシュを優先(なければ作る)
if "tsunami_union" in listlayers(GPKG):
    flood_union = gpd.read_file(GPKG, layer="tsunami_union")
else:
    flood_union = gpd.GeoDataFrame(geometry=[tsunami.geometry.union_all()], crs=tsunami.crs)
    flood_union.to_file(GPKG, layer="tsunami_union", driver="GPKG")

# 0') ハザードの簡易メタ(A40が変わったらルートキャッシュ無効化に使う)
ti = gpd.read_file(GPKG, layer="tsunami_inundation")
hazard_meta = f'{len(ti)}|{",".join(map(str, ti.total_bounds))}'

# 1) キャッシュキー生成(端点・モード・半径・ハザードメタ)
def route_key(home_wgs, dest_wgs, mode="walk", radius_m=12000, 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()

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

# 2) ルートキャッシュにヒットすれば使う
route_gdf = None
use_cached_route = False
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 = hit.set_crs(4326)
        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])
        coords = list(route_gdf.geometry.iloc[0].coords)
        route_latlon = [(lat, lon) for lon, lat in coords]  # folium用 (lat,lon)
        use_cached_route = True
        print("↩︎ ルート: キャッシュヒット")

if not use_cached_route:
    # 3) これまで通りの経路('route')から 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")

    # 4) 堅牢チェック(交差ゼロ&境界近接なし)
    route_m = route_gdf.to_crs(METRIC)
    flood_m = flood_union.to_crs(METRIC)

    inter = gpd.overlay(route_m, flood_m, how="intersection")
    inside_len_m = float(inter.geometry.length.sum()) if len(inter) else 0.0

    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 半径拡大で再探索してください。")

    # 5) 距離[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

    # 6) route_latest へ保存(検証済みのみ)
    out_gpkg = 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],
            "created_at":[datetime.utcnow().isoformat(timespec="seconds")+"Z"],
        },
        geometry=[route_geom], crs="EPSG:4326",
    )
    route_gdf.to_file(out_gpkg, layer="route_latest", driver="GPKG")
    print("✅ 交差なし・保存完了")

    # 7) ルートキャッシュに追記
    route_cache_rec = route_gdf.copy()
    route_cache_rec["rid"] = key
    try:
        route_cache_rec.to_file(out_gpkg, layer="route_cache", driver="GPKG", if_exists="append")
    except TypeError:
        # 古いGeoPandasのためのフォールバック
        if "route_cache" in listlayers(out_gpkg):
            old = gpd.read_file(out_gpkg, layer="route_cache")
            gpd.GeoDataFrame(pd.concat([old, route_cache_rec], ignore_index=True), crs=4326)\
               .to_file(out_gpkg, layer="route_cache", driver="GPKG")
        else:
            route_cache_rec.to_file(out_gpkg, layer="route_cache", driver="GPKG")

    # 8) folium用座標
    route_latlon = [(G2.nodes[n]["y"], G2.nodes[n]["x"]) for n in route]

# 9) Folium表示(キャッシュ・新規どちらでも)
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)

を
# === ルート検証 → ETA → 保存(キャッシュ対応, 上書き&追記の方針を反映) ===
from shapely.geometry import LineString, Point
import geopandas as gpd
import pandas as pd
import hashlib, json
from datetime import datetime
from fiona import listlayers
import osmnx as ox
import networkx as nx

GPKG = "outputs/shizuoka_hazard.gpkg"
METRIC = "EPSG:32654"  # 静岡周辺

# 0) 津波 union はキャッシュを優先(なければ作る)→ 上書き耐性ありで作成可
if "tsunami_union" in listlayers(GPKG):
    flood_union = gpd.read_file(GPKG, layer="tsunami_union")
else:
    flood_union = gpd.GeoDataFrame(geometry=[tsunami.geometry.union_all()], crs=tsunami.crs)
    flood_union.to_file(
        GPKG, layer="tsunami_union",
        driver="GPKG", engine="pyogrio",
        layer_options={"OVERWRITE": "YES"}  # (3) 上書き耐性
    )

# 0') ハザードの簡易メタ(A40が変わったらルートキャッシュ無効化に使う)
ti = gpd.read_file(GPKG, layer="tsunami_inundation")
hazard_meta = f'{len(ti)}|{",".join(map(str, ti.total_bounds))}'

# 1) キャッシュキー生成(端点・モード・半径・ハザードメタ)
def route_key(home_wgs, dest_wgs, mode="walk", radius_m=12000, 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()

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

# 2) ルートキャッシュにヒットすれば使う
route_gdf = None
use_cached_route = False
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 = hit.set_crs(4326)
        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])
        coords = list(route_gdf.geometry.iloc[0].coords)
        route_latlon = [(lat, lon) for lon, lat in coords]  # folium用 (lat,lon)
        use_cached_route = True
        print("↩︎ ルート: キャッシュヒット")

if not use_cached_route:
    # 3) これまで通りの経路('route')から 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")

    # 4) 堅牢チェック(交差ゼロ&境界近接なし)
    route_m = route_gdf.to_crs(METRIC)
    flood_m = flood_union.to_crs(METRIC)

    inter = gpd.overlay(route_m, flood_m, how="intersection")
    inside_len_m = float(inter.geometry.length.sum()) if len(inter) else 0.0

    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 半径拡大で再探索してください。")

    # 5) 距離[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

    # 6) route_latest へ保存(検証済みのみ)→ 上書き
    out_gpkg = 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],
            "created_at":[datetime.utcnow().isoformat(timespec="seconds")+"Z"],
        },
        geometry=[route_geom], crs="EPSG:4326",
    )
    route_gdf.to_file(
        out_gpkg, layer="route_latest",
        driver="GPKG", engine="pyogrio",
        layer_options={"OVERWRITE": "YES"}  # (1) 上書き
    )
    print("✅ 交差なし・保存完了")

    # 7) ルートキャッシュに追記(pyogrio 流儀)
    route_cache_rec = route_gdf.copy()
    route_cache_rec["rid"] = key
    try:
        route_cache_rec.to_file(
            out_gpkg, layer="route_cache",
            driver="GPKG", engine="pyogrio",
            append=True                     # (2) 追記
        )
    except Exception:
        # フォールバック:スキーマ差や古い環境でも確実に追記する
        if "route_cache" in listlayers(out_gpkg):
            old = gpd.read_file(out_gpkg, layer="route_cache")
            cols = sorted(set(old.columns) | set(route_cache_rec.columns))
            old = old.reindex(columns=cols)
            new = route_cache_rec.reindex(columns=cols)
            merged = gpd.GeoDataFrame(pd.concat([old, new], ignore_index=True), crs=4326)
            merged.to_file(
                out_gpkg, layer="route_cache",
                driver="GPKG", engine="pyogrio",
                layer_options={"OVERWRITE": "YES"}  # マージして上書き
            )
        else:
            route_cache_rec.to_file(
                out_gpkg, layer="route_cache",
                driver="GPKG", engine="pyogrio"
            )

    # 8) folium用座標
    route_latlon = [(G2.nodes[n]["y"], G2.nodes[n]["x"]) for n in route]

# 9) Folium表示(キャッシュ・新規どちらでも)
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)

へ変更

これで動くようになった

!ogrinfo -so outputs/shizuoka_hazard.gpkg route_latest

で動作も問題なし

import geopandas as gpd
gpd.read_file("outputs/shizuoka_hazard.gpkg", layer="route_latest")[["meters","eta_min_walk","eta_min_walk_safe","eta_min_drive"]]

も問題なし

コメントを残す

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です