Claude Media
Claude Codeで標高APIを引いて縦断図を描く — 国土地理院のデータなし点の扱い

Claude Codeで標高APIを引いて縦断図を描く — 国土地理院のデータなし点の扱い

国土地理院の標高APIを一定間隔で呼び、路線沿いの地盤高から概略の縦断図を描きます。「-----」が返るデータなしの点を欠測として残す、補間・取得・描画の3本のスクリプトを載せます。

道路や水路の検討の初期段階では、測量成果が手元にない状態で「この路線は登りか下りか、どこが谷か」を知りたい場面があります。国土地理院の標高APIは、緯度経度を渡すと地盤の標高を返すWebエンドポイントです。中心線の座標列を一定間隔に刻んで順に問い合わせれば、概略の縦断図になります。

作るのは3つです。座標列を距離で刻むスクリプト、標高を取って表にするスクリプト、グラフにするスクリプト。ここで一番のつまずきは、海上や範囲外の点で標高の代わりに文字列「-----」が返ることです。これを0として線につなぐと、東京湾のど真ん中に深さ0mの谷ができます。データなしを欠測として残す書き方を中心に、Claude Codeに回させる設定まで通して扱います。

国土地理院の標高APIは何を返すのか

APIは「サーバサイドで経緯度から標高を求めるプログラム」として、地理院地図の技術情報のページに載っています。エンドポイントは getelevation.php で、必須パラメータは経度の lon と緯度の lat(いずれも度の10進法)です。結果の形式は、JSONPなら callback、JSONなら outtype=JSON で選び、2つは同時に指定しません。

返る項目は2つです。

項目内容エラー時
elevation内容標高値エラー時文字列「-----」
hsrc内容標高データのデータソースエラー時文字列「-----」

hsrc は「1m(レーザ)」「5m(レーザ)」「5m(写真測量)」「5m(写真測量5C)」「10m」のどれかです。その地点で最も精度のよいデータソースの値が返ります。標高値の桁数もソースで変わり、1mの元データなら0.1mの位まで、10mの元データなら1mの位までです。

緯度と経度の取り違えには注意が要ります。URLのパラメータには名前が付いていますが、手元の座標列が「緯度, 経度」の並びのとき、列を入れ替えずに渡すと別の場所を引きます。座標列を読む段階で列名を lon と lat に固定しておけば、この取り違えは起きません。

公式ページが示す利用上の約束

使用上の注意には2点あります。ひとつは、サーバに過度の負担を与えるアクセスは予告なく遮断される場合があること。もうひとつは、プログラムの内容変更や提供停止が予告なく行われる可能性があることです。

どちらも、縦断図のスクリプトの作り方に直結します。遮断されないよう1秒に1回に絞り、返り値の形が変わったときは黙って通さず例外で止めます。何点まで許されるかは書かれていないので、数百点規模の一括取得を日常的に回す用途には向きません。

実際に叩いた結果と、空の応答への備え

東京駅付近(経度139.7671、緯度35.6812)に問い合わせると、次の応答が返りました。

curl -s "https://cyberjapandata2.gsi.go.jp/general/dem/scripts/getelevation.php\
?lon=139.7671&lat=35.6812&outtype=JSON"
{"elevation":3.6,"hsrc":"1m(レーザ)"}

hsrc の日本語は \uXXXX のエスケープで返るため、JSONとして読めば「1m(レーザ)」に戻ります。文字列のまま比較せず、必ずJSONパーサーを通してください。

東京湾上(経度139.85、緯度35.50)では、次の応答でした。

{"elevation":"-----","hsrc":"-----"}

HTTPステータスは200のままで、エラーは本文の文字列でしか分かりません。ステータスコードだけを見るコードは、これを見逃します。数値であるはずの elevation が文字列で返る、という型の揺れがあるので、float() に直接渡すと例外になるか、「-----」を握りつぶす書き方をすると0に化けます。

もうひとつ、手元で見た挙動があります。間隔を置かずに curl を続けて叩いたとき、HTTP 200のまま本文が空の応答が返ることがありました。数秒空けて同じ座標を叩き直すと、正常な応答が返りました。仕様に書かれた挙動ではなく、手元で数回遭遇した観察です。原因は分かっていませんが、空の本文をJSONとして読もうとして落ちる実装は避けたほうが安全なので、後述のスクリプトでは空なら待って再試行します。

標高タイル方式との違い

同じ技術情報のページには、クライアントサイドで標高を求める方式も載っています。標高タイル(PNG)をダウンロードし、ピクセルの画素値(RGB)から標高を計算するもので、タイルの種別はDEM1A、DEM5A、DEM5B、DEM5C、DEM10B、DEMGMの6種類です。海部などデータがない場所のピクセル値は(128,0,0)になる、という説明も同じページにあります。

本記事のAPI方式は、この計算をサーバ側が済ませて数値で返してくれる点が違います。数点から数十点を引く用途ならAPIが手軽で、数千点を超えるなら、標高タイルを取得して手元で計算する方式も選択肢になります。

縦断図に必要なのは「一定間隔の点列」と「欠測の区別」

縦断図は、横軸が路線に沿った距離、縦軸が標高のグラフです。標高APIは点ごとにしか答えないので、中心線の折れ線を一定の間隔で刻んだ点列をこちらで作ります。ここで3つの設計判断があります。

  • 距離は緯度経度の差ではなく、地球上の2点間距離(球面上の大円距離)で測る。緯度経度の差を足すと、東西と南北で長さの尺度が違うため、距離軸が歪む
  • 終点は必ず点列に含める。間隔で割り切れない端数の区間でも、終点を落とさない
  • データなしの点は空欄にして、hsrc に「nodata」と書く。0や直前の値で埋めない

3つ目が縦断図の信頼性を決めます。グラフでは欠測の前後で線を切り、表では空欄のまま残すことで、読む人は「そこは測れていない」と分かります。

取得元が混在する路線の読み方

1本の路線でも、地点ごとに hsrc が変わります。市街地は1m(レーザ)でも、郊外や臨海部は10mに落ちることがあるためです。実際に東京駅から南東へ引いた試行では、1m(レーザ)の点の並びの途中に10mの点が1つ混ざりました。

10mのデータは等高線から作られた元データなので、近くの1m(レーザ)の点と並べると、細かい起伏の出方が変わります。縦断図の1本の線に精度の違う点が混在することは避けられないため、表には hsrc の列を残し、後から「この区間は10m」と分かるようにしておきます。

標高を取って表にするスクリプト

次の profile.py は、中心線のCSV(lon,lat の2列)を読み、指定間隔で刻み、標高を取ってCSVに書き出します。標準ライブラリだけで動きます。

import csv, json, math, sys, time, urllib.parse, urllib.request
 
API = "https://cyberjapandata2.gsi.go.jp/general/dem/scripts/getelevation.php"
 
def haversine(a, b):
    r = 6371008.8
    p1, p2 = math.radians(a[1]), math.radians(b[1])
    dp, dl = p2 - p1, math.radians(b[0] - a[0])
    h = math.sin(dp / 2) ** 2 + math.cos(p1) * math.cos(p2) * math.sin(dl / 2) ** 2
    return 2 * r * math.asin(math.sqrt(h))
 
def resample(points, step):
    """(lon, lat) の折れ線を step[m] ごとの点列にする。端点は必ず含める。"""
    out, total, next_at = [(points[0], 0.0)], 0.0, step
    for a, b in zip(points, points[1:]):
        seg = haversine(a, b)
        while next_at <= total + seg:
            t = (next_at - total) / seg
            out.append(((a[0] + (b[0] - a[0]) * t,
                         a[1] + (b[1] - a[1]) * t), next_at))
            next_at += step
        total += seg
    if out[-1][1] < total - 1e-6:
        out.append((points[-1], total))
    return out
 
def elevation(lon, lat):
    q = urllib.parse.urlencode(
        {"lon": f"{lon:.6f}", "lat": f"{lat:.6f}", "outtype": "JSON"})
    for _ in range(3):
        with urllib.request.urlopen(f"{API}?{q}", timeout=20) as res:
            text = res.read().decode("utf-8")
        if text.strip():
            break
        time.sleep(5)              # 空の応答は待って再試行
    else:
        raise RuntimeError("3 回続けて空の応答")
    body = json.loads(text)
    if set(body) != {"elevation", "hsrc"}:
        raise RuntimeError(f"想定外の応答: {body}")
    if body["elevation"] == "-----":
        return None, None          # データなし(海上・範囲外など)
    return float(body["elevation"]), body["hsrc"]
 
def main(route_csv, step, out_csv):
    with open(route_csv, newline="") as f:
        pts = [(float(r["lon"]), float(r["lat"])) for r in csv.DictReader(f)]
    rows = []
    for (lon, lat), d in resample(pts, step):
        e, src = elevation(lon, lat)
        rows.append({"distance_m": round(d), "lon": f"{lon:.6f}",
                     "lat": f"{lat:.6f}",
                     "elevation_m": "" if e is None else e,
                     "hsrc": src or "nodata"})
        time.sleep(1)              # 1 秒 1 リクエスト
    with open(out_csv, "w", newline="") as f:
        w = csv.DictWriter(f, fieldnames=rows[0].keys())
        w.writeheader()
        w.writerows(rows)
 
if __name__ == "__main__":
    main(sys.argv[1], float(sys.argv[2]), sys.argv[3])

読みどころは elevation 関数の3か所です。応答のキーが elevation と hsrc の2つ以外なら例外で止める、データなしは None で返す、空の応答は待って再試行する。time.sleep(1) は、前節の利用上の注意に合わせた間隔です。

東京駅から上野へ、500mごとに引いた結果

東京駅(経度139.7671、緯度35.6812)から途中の点を経て上野(経度139.7774、緯度35.7141)に至る3点のCSVを作り、500mおきに引いた実行結果です。

python3 -I profile.py route.csv 500 out.csv
distance_melevation_mhsrc
0elevation_m3.6hsrc1m(レーザ)
500elevation_m3.0hsrc1m(レーザ)
1000elevation_m4.3hsrc1m(レーザ)
1500elevation_m4.4hsrc1m(レーザ)
2000elevation_m3.6hsrc1m(レーザ)
2500elevation_m2.5hsrc1m(レーザ)
3000elevation_m2.3hsrc1m(レーザ)
3500elevation_m3.2hsrc1m(レーザ)
3791elevation_m11.1hsrc1m(レーザ)

最後の行が終点の上野で、間隔で割り切れない端数(3500mから291m)の区間を、終点として残しています。ここで標高が3.2mから11.1mへ約8m上がるため、500m間隔の点だけでは、この台地への立ち上がりの形は分かりません。

海を横切る路線でのデータなし点

東京駅から東京湾上の座標(経度139.85、緯度35.50)へ向けて4000mおきに引くと、次の結果でした。

distance_melevation_mhsrc
0elevation_m3.6hsrc1m(レーザ)
4000elevation_m(空欄)hsrcnodata
8000elevation_m(空欄)hsrcnodata
12000elevation_m(空欄)hsrcnodata
16000elevation_m(空欄)hsrcnodata
20000elevation_m(空欄)hsrcnodata
21498elevation_m(空欄)hsrcnodata

起点以外の6点がすべてデータなしです。この表を縦断図にすると、グラフに残るのは起点の1点だけで、線は引けません。それが正しい出力です。陸から海に出る路線の縦断図は、海岸線から先は描かない、が正解になります。

落とし穴

データなしを0として描くと、海面に谷ができる

「-----」を float() で落とさないために、try ... except で0を入れる書き方を見かけます。これは縦断図では誤りです。標高0mは実在し得る値で、データなしとは別物です。欠測は空欄にし、グラフでは線を切ってください。

縦断図に描くスクリプト

plot.py は、profile.py が書き出したCSVからグラフを描きます。欠測は nan にして、その前後で線が切れるようにします。matplotlibが必要で、未導入なら pip install matplotlib で入れます。

import csv, sys
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
 
def main(src, dst):
    with open(src, newline="") as f:
        rows = list(csv.DictReader(f))
    dist = [float(r["distance_m"]) for r in rows]
    elev = [float(r["elevation_m"]) if r["elevation_m"] else float("nan")
            for r in rows]
    fig, ax = plt.subplots(figsize=(9, 3.5))
    ax.plot(dist, elev, marker="o", markersize=3, linewidth=1)
    for r in rows:
        if r["hsrc"] == "nodata":
            ax.axvline(float(r["distance_m"]), color="gray",
                       linewidth=0.5, linestyle=":")
    ax.set_xlabel("distance [m]")
    ax.set_ylabel("elevation [m]")
    fig.tight_layout()
    fig.savefig(dst, dpi=150)
 
if __name__ == "__main__":
    main(sys.argv[1], sys.argv[2])

データなしの位置には灰色の点線を縦に引き、線が途切れた理由をグラフ上で示しています。日本語のラベルを使う場合は、環境に入っている日本語フォントを font.family に指定する必要があるので、スクリプトの軸ラベルは英字にしてあります。

Claude Codeに回させる設定

ここまでのスクリプトは、人の手でも動かせます。Claude Codeに任せる利点は、中心線のCSVの整形、間隔の調整、実行、結果の読み取りまでを一続きの会話で回せることです。そのために、守らせたいことをCLAUDE.mdに書きます。

CLAUDE.mdは毎セッションの開始時に読み込まれる指示ですが、強制される設定ではなく文脈として扱われます。公式のメモリのページには、指示は具体的で短いほど守られやすいと書かれています。次のように、検証できる形で書きます。

# 標高の縦断図
 
- 標高は国土地理院のAPI(getelevation.php)だけで取る。他の標高ソースに切り替えない
- リクエストは1秒に1回まで。並列に呼ばない
- 応答の elevation が文字列「-----」の点は欠測。0や前後の値で補わず、CSVでは空欄、hsrc は nodata にする
- 表には必ず hsrc の列を残す。精度の違う点が混ざっていることを後から追えるようにする
- 中心線CSVは lon,lat の2列。緯度経度の並びを入れ替えない
- 出力する縦断図には「概略。測量成果ではない」と注記する

CLAUDE.mdが「守られるべきこと」を伝える層なら、権限は「実行できること」を決める層です。非対話モードで回すなら、--allowedTools で許可を絞ります。許可の考え方はClaude Codeの/permissionsコマンドで権限ルールとauto mode拒否を管理するにまとめています。

claude -p "route.csv を500m間隔で刻んで標高を取り、out.csv と profile.png を作って。\
データなしの点の数も報告して" \
  --allowedTools "Read,Bash(python3 profile.py *),Bash(python3 plot.py *)"

Bash に何も付けずに許可すると、あらゆるシェルコマンドが通ります。ここでは、2つのスクリプトの呼び出しだけに絞っています。claude -p の使い方と、--bare で設定の読み込みを止める判断はclaude -pオプションとはで扱っています。--bare を付けるとCLAUDE.mdも読まれなくなる点は、この用途では「欠測の扱いを書いた指示が効かなくなる」ことを意味します。付けずに回すか、指示をプロンプトに含めてください。

書き方そのものの型はCLAUDE.mdの書き方に詳しいので、ここでは標高に固有の項目だけを挙げました。

Claudeに返させる報告の形

実行後の報告は、グラフの出来を褒めさせるより、機械的に数えられる事実を返させるほうが確実です。

  • 点の総数、取得できた点の数、データなしの点の数
  • hsrc ごとの点数(1m(レーザ)が何点、10mが何点か)
  • 標高の最小値・最大値と、その位置の距離

これはスクリプトの出力から数えられるので、数字が食い違えばCSVを開いて確かめられます。住所のような文字列と違い、標高は桁が合っていても意味が違うことがあるため、点数の内訳を毎回出させる運用が検証の足がかりになります。座標の検証という点では、国土地理院の住所検索APIでClaude抽出の住所を検証するの手順で、中心線を作る前の地名の取り違えを先に潰せます。

概略として使うための限界

この方法で得られるのは、地盤の標高の概略です。使う前に、次の限界を頭に置いてください。

  • 点と点の間は分からない。500m間隔で引くと、数十mの幅の谷や盛土は点列に出ない。前掲の上野では、端数の291m区間だけで約8m上がっていた
  • 公式のページには、橋梁やトンネル、盛土・切土の区間で何を返すかの説明がない。道路や鉄道の「計画高」と、APIが返す地盤の標高は別のものとして扱う
  • 返り値の桁数はデータソースで違う。1mの元データなら0.1m、10mの元データなら1mの位までなので、2桁目の小数を比べて差を論じない
  • 提供が予告なく変わる可能性が公式に書かれている。業務の成果品に直接組み込まず、検討用の下図に留める

間隔を狭めたいときは、リクエスト数との兼ね合いです。10kmの路線を100m間隔で刻むと約100点で、1秒1回なら2分弱で終わります。10m間隔にすると約1000点、16分強になり、サーバへの負担も増えます。点を増やす前に、本当にその細かさが検討に必要かを確かめてください。

要点

標高APIで縦断図を作るときの分かれ目

  • 標高APIは、海上などで「-----」を返す。ステータスは200のままなので、本文で判定する
  • データなしは空欄にして、グラフの線を切る。0で埋めると、存在しない谷ができる
  • 距離は球面上の2点間距離で測り、終点を必ず点列に残す
  • hsrc を表に残し、10mなど精度の違う点を後から追えるようにする
  • 1秒1回に絞り、応答の形が変わったら例外で止める

縦断図は、グラフの見栄えよりも、どの点が測れていてどの点が測れていないかを読み手が追えるかで価値が決まります。データなしを0にしない、取得元を残す、概略と明記する。この3つが守られていれば、スクリプトを数十行で作っても、検討の下図として使える形になります。

参照した一次ソース

この記事の記述は、次の一次ソースと照合しています。更新時に参照先が変わった場合は一覧も更新します。

  1. maps.gsi.go.jp/development/elevation_s.html
  2. maps.gsi.go.jp/development/elevation.html
  3. code.claude.com/docs/en/headless
  4. code.claude.com/docs/en/memory
この記事を共有:XはてブLinkedIn