L42 地図情報_3次元点群データ_オープン 3 件統合分析 — 1200 リソースの航空レーザ測量カタログを 3 vintage で読み解く

データ取得手順

✅ このスクリプトは初回実行時にデータを自動取得します（DoBoX からの直接ダウンロード）。

実行コマンド:

cd "2026 DoBoX 教材"
python -X utf8 lessons/L42_map_information.py

DoBoX のオープンデータは申請不要・商用/非商用とも利用可。 data/extras/ は .gitignore 対象（約 57 GB のキャッシュ）。 スクリプト実行で自動再生成されます。

学習目標と問い

研究の問い (RQ)

DoBoX のシリーズ 「地図情報_3次元点群データ_オープン」 3 件:

• (a) dataset 65 = 「地図情報_3次元点群データ_オープン」(2016-2019 公開, 計測 2014-2018)
• (b) dataset 1434 = 「同_2023」(2023 公開, 計測 2018, 2022)
• (c) dataset 1527 = 「同_2024」(2024 公開, 計測 2023)

これら 3 dataset は同じシリーズ名 (= 地図情報_3次元点群データ) を持つが、 vintage (年度) で分割された 3 つの公開単位である。 本記事は「3 dataset を時系列で並べると、広島県の航空レーザ測量の整備史と将来像はどう描けるか」を、 1200 resource のメタ情報 + 23 メタデータ XML + 代表サンプルファイル 3 種類 (DXF/TXT/LAS) を実データで読み解く形で量的に検証する。

独自用語の定義 (本記事内のみ)

立てた仮説 (H1〜H6)

到達点

• 3 dataset = 3 vintage の時系列タイムラインを可視化
• 1200 resource をプロダクト × datasetのクロス集計で構造化
• 23 業務メタデータ XML を解読し計測期間・対象圏域・公開ラグを量化
• 図郭 ID 命名体系を解読し広島県の国土基本図カバレッジを地区別に集計
• 各プロダクト (DXF / TXT / LAS) の実データ構造を 1 件ずつ解析 (Before/After)
• L36-L40 の LiDAR ファミリ派生層に対する原料層としての位置づけを示す

使用データ

本記事では DoBoX の「地図情報_3次元点群データ_オープン」シリーズ 3 件 (dataset_id = 65, 1434, 1527; 合計 1200 resource) を統合する。 すべて広島県土木建築局が砂防基礎調査の航空レーザ測量成果として公開した 公共測量 2 次著作物。CRS はEPSG:6671 (JGD2011 / 平面直角座標系 III 系)。 ライセンスはクリエイティブ・コモンズ表示 (CC-BY)。

シリーズの実体 (DoBoX タイトルの誤誘導)

「地図情報」というタイトルは抽象的だが、内部メタデータ (JMP 2.0 規格) の title 要素を 読むと、実体は「砂防基礎調査に伴う航空レーザ測量及び撮影業務」の成果である。 広島県土木建築局が「土砂災害警戒区域等における土砂災害防止対策の推進に関する法律」第 4 条に基づき、 広島県全域 (8479 km²) を順次計測したもの。

L36-L40 (LiDAR 派生層) との関係

本記事の 3 dataset は L36-L40 のラスタ系 5 dataset の原料に相当する:

つまり LiDAR ファミリは「点群 (本記事 LAS) → DEM (本記事 TXT) → 派生ラスタ (L36-L40) → 林業ベクタ (L41)」 の 4 段派生階層。本記事は最上流 (= 原料層) を扱う。

データ仕様

• dataset 65 (オリジナル): resource ~400 件。2014-2018 計測 (10 業務)。地図情報レベル 1000 (1m DEM)。 プロダクト: メタ XML / 等高線 DXF / オルソ TIFF / 水部 DXF / DEM TXT (4 種類) 。
• dataset 1434 (2023 公開): resource ~400 件。2018, 2022 計測 (8 業務)。地図情報レベル 1000 → 500 移行期。 プロダクト: 上記と同等。
• dataset 1527 (2024 公開): resource ~400 件。2023 計測 (5 業務 = 三次市 2 + 庄原市 3)。 地図情報レベル 500 (50cm DEM)。プロダクト: 上記 + LAS 点群 (新規追加!) 。

共通メタデータ (JMP 2.0 形式)

各業務には JMP (Japan Metadata Profile) 2.0 形式の XML メタデータが付属。 含まれる要素: title / date / abstract / pointOfContact / language / topicCategory / extent (geographic + temporal) / referenceSystemIdentifier。 本記事はこの XML 23 件を全件パースし、業務単位の構造化情報を作成した。

ダウンロード（再現用データ・中間データ・図）

本記事の再現に必要なデータ・中間データ・図はすべて以下から直 DL 可能。

原データ (DoBoX 直リンク)

• dataset 65: DoBoX 地図情報_3次元点群データ_オープン — リソース一覧から個別 DL (各 resource ページ右下の data.hiroshima-dobox.jp/... 直リンク)
• dataset 1434: DoBoX 地図情報_3次元点群データ_オープン_2023
• dataset 1527: DoBoX 地図情報_3次元点群データ_オープン_2024

注意: DoBoX のリソース個別ページの「ダウンロード」ボタンは内部リンク (/resource_download/{rid}) ではなく, ページ内に直接記された data.hiroshima-dobox.jp/aerial_survey/{vintage}/{kind}/... 形式の S3/CloudFront 直 URL を踏む。本記事の取得スクリプトはこの直 URL を抽出して取得する。

中間データ (本スクリプトが生成)

• L42_dataset_summary.csv — 3 dataset 比較サマリ
• L42_business_units.csv — 23 業務メタ (期間・圏域・地図情報レベル)
• L42_resources_full.csv — 1200 resource 一覧 (kind/sheet/year)
• L42_products_cross.csv — ds × kind プロダクト構成クロス
• L42_sheet_overlap.csv — 図郭重複 (3 dataset 間)
• L42_sheet_prefix.csv — 図郭プレフィックス頻度 (国土基本図地区別)
• L42_sheet_prefix_cross.csv — prefix × ds クロス
• L42_region_coverage.csv — 業務 × 圏域カバレッジ
• L42_file_struct.csv — 各プロダクトの実データ構造 (DXF/TXT/LAS)
• L42_hypothesis.csv — H1〜H6 仮説検証結果

図 (PNG)

• 図 1: vintage タイムライン
• 図 2: dataset × プロダクト構成
• 図 3: 23 業務ガント
• 図 4: 図郭プレフィックス × ds ヒートマップ
• 図 5: 等高線 DXF 実データ
• 図 6: DEM TXT 実データ
• 図 7: LAS 点群ヘッダ
• 図 8: 圏域カバレッジ進化
• 図 9: 計測期間と公開ラグ

再現スクリプト

L42_map_information.py をダウンロードし、 データキャッシュ (data/extras/L42_map_information/samples_meta_full.json など) があれば即座に再現可能。 無ければ DoBoX から resource 一覧 + メタデータ XML を取得 (約 1 分)。

cd "2026 DoBoX 教材"
py -X utf8 lessons\L42_map_information.py

分析 1: 3 dataset の vintage 構造比較

狙い

3 dataset = 3 vintage を時系列で並べ、広島県の航空レーザ測量整備の進化を量的に把握する。 本記事の主命題「3 dataset は時系列で別 dataset として整理されている」(H1) を量的に検証する。

手法 (リテラシレベル解説)

3 段階の処理:

1. resource 一覧の取得: 各 dataset ページを paging スクレイピング (https://hiroshima-dobox.jp/datasets/{did}?page=1, 2, ...) して全 resource_id を取得。 各 resource ページのタイトルからプロダクト種別・図郭 ID・年度を抽出。
2. メタデータ XML の取得: 各 dataset 内の「メタデータ」種別 resource (合計 23 件) からJMP 2.0 形式の XML を取得し、Python xml.etree.ElementTree でパース。
3. クロス集計: pandas.DataFrame で resource × kind × ds × year の多次元クロスを作成。

• 入力: 3 dataset URL (合計 ~120 ページ + 1200 resource ページ + 23 XML)
• 出力: 1200 行 × 8 列の resource 表 + 23 行の業務表
• 限界: タイトル抽出に正規表現を使うので、命名規則が将来変わると壊れる (本記事の前提条件)

実装

↑ L42_map_information.py 行 1418–1479

 1
 2
 3
 4
 5
 6
 7
 8
 9
1427
1428
1429
1430
1431
1432
1433
1434
1435
1436
1437
1438
1439
1440
1441
1442
1443
1444
1445

# resource 一覧の paging 取得
def fetch_dataset_resources(did, max_pages=50):
    seen = []
    for p in range(1, max_pages+1):
        html = http_get(f"https://hiroshima-dobox.jp/datasets/{did}?page={p}")
        rids = re.findall(r"/resources/(\\d+)", html)
        new = [r for r in rids if r not in seen]
        seen.extend(new)
        if not new and p > 1: break
    return seen

# resource ページからタイトル + データ URL を抽出
def fetch_resource_meta(rid):
    html = http_get(f"https://hiroshima-dobox.jp/resources/{rid}")
    title = re.search(r"<title>([^<]+)</title>", html).group(1)
    durl = re.search(r"(https?://data\\.hiroshima-dobox\\.jp/[^\"\'\s>]+)", html)
    return {"rid": rid, "title": title, "data_url": durl.group(1) if durl else None}

# JMP 2.0 XML のパース
def parse_xml(path):
    tree = ET.fromstring(path.read_text(encoding="utf-8-sig"))
    ns = {"j": "http://zgate.gsi.go.jp/ch/jmp/"}
    return {
        "title":    tree.findtext(".//j:citation/j:title", default="", namespaces=ns),
        "abstract": tree.findtext(".//j:abstract",         default="", namespaces=ns),
        "begin":    tree.findtext(".//j:beginEnd/j:begin", default="", namespaces=ns),
        "end":      tree.findtext(".//j:beginEnd/j:end",   default="", namespaces=ns),
    }

図と読み取り (図 1: vintage タイムライン)

なぜこの図か: 3 dataset を 1 軸時系列で並べると、計測期間 (バー長) と公開日 (▼) の関係が一目で分かる。 y 軸 3 段に分けて 3 dataset の運用ペースを比較する。

図 1: 3 dataset = 3 vintage の航空レーザ測量タイムライン (23 業務)

読み取り (重要発見):

• ds65 (オリジナル): 2014-2018 の4 年間に 10 業務が並列実施。広島県全域を 4 年で概覆する初期整備期。
• ds1434 (2023 公開): 2018 年 1 業務 + 2022 年に集中 7 業務。2018 年版を移行扱いで残しつつ, 主体は 2022-2023 年の市町別細密計測。
• ds1527 (2024 公開): 2023 年に集中 5 業務。三次市・庄原市 (県北部山間) のみ。地図情報レベル 500 で更新を進める段階。
• 公開日 (▼) は計測終了から3-12 ヶ月後。最近 vintage ほどラグが短い。
• 3 dataset を時系列で並べると「県全域概覆 (1m) → 市町別更新 (1m→50cm) → 中山間集中更新 (50cm)」という整備戦略が見える。

表 (3 dataset サマリ)

読み取り: 3 dataset はresource 件数 (~400) はほぼ均等だが、業務数 (10/8/5)・図郭数・地図情報レベルが異なる。 ds1527 はプロダクト種類数 (5) が最多 (LAS 点群が追加されたため) 。

分析 2: プロダクト構成の進化 (LAS 点群の登場)

狙い

1200 resource をプロダクト種別 (kind) で分類し、3 dataset 間でどのプロダクトが増減しているかを量化する。 H3 (ds1527 にのみ LAS 点群) を検証。

手法 (リテラシレベル解説)

resource タイトル文字列のパターンマッチで分類:

• 入力: 1200 行のタイトル文字列
• 出力: ds × kind の 3 × 6 クロステーブル
• 限界: タイトル命名規則が将来変わると分類が破綻 (現状は全件分類成功)

実装

↑ L42_map_information.py 行 1497–1524

 1
 2
 3
 4
 5
 6
 7
 8
 9
1506
1507

def classify_kind(title: str) -> str:
    if "メタデータ" in title:        return "メタデータ"
    if "オルソ画像" in title:        return "オルソ画像 (TIFF)"
    if "等高線" in title:            return "等高線 (DXF)"
    if "LAS" in title:               return "点群LAS"
    if "グリッドデータ" in title:     return "DEM グリッド (TXT)"
    if "水部" in title:              return "水部ポリゴン (DXF)"
    return "その他"

res_df["kind"] = res_df["title"].apply(classify_kind)
prod_cross = res_df.groupby(["ds", "kind"]).size().unstack(fill_value=0)

図と読み取り (図 2: プロダクト構成積み上げ)

なぜこの図か: 3 dataset 間でどのプロダクトが追加・削減されたかを積み上げ棒で見る。 色がプロダクト種別、棒の高さが resource 件数。

図 2: dataset × プロダクト種別 (1200 resource の積み上げ)

読み取り (重要発見):

• ds65 (オリジナル): 4 種類 (DEM, 等高線, オルソ, 水部) + メタ。各種類 ~110 件で均等。
• ds1434: 同 4 種類だが、水部が ds65 の 52 件 → ds1434 の 30 件に減少。 = 同じ図郭でも水部 polygon が小さい (= 山地中心) 業務が増えた。
• ds1527: LAS 点群が新規 64 件追加! オルソが 187 件 (1.6 倍) に膨張。 = 「raw 点群を public に公開する方針」と「オルソ画像の高解像度化」 の 2 方向に進化。
• 3 dataset ともメタデータ XML 件数 ≈ 業務数 (10, 8, 5)。これが本記事の業務単位識別の根拠。
• H3 「ds1527 にのみ LAS 点群」 → 支持 (ds65=0, ds1434=0, ds1527=64件)。 航空レーザ測量界の「派生だけ → raw も公開」のトレンドと整合。

表 (kind × ds クロス)

読み取り: 縦読み (ds 別) は構成比、横読み (kind 別) はプロダクトごとの発展傾向を示す。 ds1527 のオルソ激増 (187 件) は、1 図郭につき複数年分のオルソを並列公開している可能性が高い。

分析 3: 業務 (XML 23 件) の解読

狙い

23 メタデータ XML (JMP 2.0 形式) を全件パースし、業務単位の構造化情報を作成する。 業務名・対象圏域・計測期間・公開日・地図情報レベル を読み解き、H4 (vintage 進化のカバレッジ) を検証。

手法 (リテラシレベル解説)

JMP 2.0 (Japan Metadata Profile 2.0) は、国土地理院が定める空間情報メタデータの XML 形式 (ISO 19115 を日本向けに簡略化)。 本データの XML は http://zgate.gsi.go.jp/ch/jmp/ 名前空間下に以下の主要要素を持つ:

• citation/title: 業務名 (例: 「砂防基礎調査に伴う航空レーザ測量及び撮影業務その3」)
• citation/date: 公開日 (例: 2016-10-31)
• abstract: 概要 (例: 「地図情報レベル1000」「地図情報レベル500」)
• extent/temporalElement/beginEnd: 計測期間 (例: 2015-11-21 → 2016-02-18)
• extent/geographicIdentifier/code または description: 対象圏域 (例: 「太田川圏域」「庄原市北部」)
• referenceSystemIdentifier/code: CRS (全件 「JGD2011 / 3(X,Y)」 = EPSG:6671)

実装

↑ L42_map_information.py 行 1557–1602

 1
 2
 3
 4
 5
 6
 7
 8
 9
1566
1567
1568
1569
1570
1571
1572
1573
1574
1575
1576

import xml.etree.ElementTree as ET
def parse_xml(path):
    tree = ET.fromstring(path.read_text(encoding="utf-8-sig"))
    ns = {"j": "http://zgate.gsi.go.jp/ch/jmp/"}
    return {
        "title":      tree.findtext(".//j:citation/j:title", default="", namespaces=ns),
        "date":       tree.findtext(".//j:citation/j:date/j:date", default="", namespaces=ns),
        "abstract":   tree.findtext(".//j:abstract", default="", namespaces=ns),
        "geographic": tree.findtext(".//j:geographicIdentifier/j:code",
                       default=tree.findtext(".//j:extent/j:description", default="",
                                              namespaces=ns), namespaces=ns),
        "begin":      tree.findtext(".//j:beginEnd/j:begin", default="", namespaces=ns),
        "end":        tree.findtext(".//j:beginEnd/j:end",   default="", namespaces=ns),
    }

# 地図情報レベル抽出 (全角→半角変換も)
def parse_map_level(s):
    s = s.translate(str.maketrans("０１２３４５６７８９", "0123456789"))
    m = re.search(r"地図情報レベル\\s*(\\d+)", s)
    return int(m.group(1)) if m else None

図と読み取り (図 3: 業務ガント)

なぜこの図か: 23 業務を縦に並べ、横軸時間で計測期間バー + 公開日 ▼ + 圏域名 を 1 枚に。 ds 色分けで時系列パターンを把握する。

図 3: 23 業務 (= メタデータ XML 単位) のガント

読み取り (重要発見):

• ds65 期 (赤, 2014-2018): 計測期間が短い (1-3 ヶ月) ものから長い (1 年以上) ものまで多様。10 業務が並列重複し急ピッチで広島県全域を概覆。
• ds1434 期 (青, 2018, 2022): 2022 年夏の数ヶ月に 7 業務が同期並列。県内 8 グループに分割して同時計測。 公開ラグはほぼ揃って 2023 年 6-8 月。
• ds1527 期 (緑, 2023): 2023 年 4-7 月の 4 ヶ月に 5 業務が密集。三次市・庄原市の高解像度更新に集中。 公開ラグ ~2 週間と最短。
• 業務名から、ds65 は「砂防基礎調査に伴う...その N」(連番)、ds1434 以降は「広島県航空レーザ測量業務 (...)」(地名指定) と命名規則が変化している。
• これは行政契約の構造変化を反映: 連番方式 (一括発注) → 地名分割方式 (年度更新)。

表 (圏域カバレッジクロス)

読み取り: 圏域分類の変化:

• ds65 は「県全域」「太田川圏域」「江の川圏域」の 圏域 (流域) ベース。
• ds1434 は沿岸/島嶼 (呉, 江田島, 大崎上島) + 市町別 (広島市・福山市等)。沿岸＋市町ベース。
• ds1527 は「三次・庄原」のみ (中山間内陸)。
• = 行政の計測戦略は時代で変わる: 流域単位 → 都市市町単位 → 山間市町単位。

図 (図 8: 圏域カバレッジ進化)

図 8: 業務単位の圏域カバレッジ進化 (3 dataset × 圏域)

読み取り: 縦バーで業務件数を圏域別に積み上げ。ds65 → ds1434 → ds1527 で「県全域 → 沿岸/市町 → 三次・庄原」と変化。これは H4 を支持する直接証拠。

図 (図 9: 計測期間と公開ラグ)

図 9: 計測期間と公開ラグ — 3 dataset の運用速度比較

読み取り:

• (a) 計測期間: ds65 は長期業務(180-400 日) と短期業務 (1-2 ヶ月) が混在。ds1434/1527 は短期 (1-3 ヶ月) に標準化。
• (b) 公開ラグ: ds65 は1-2 年かかる業務もあった。ds1527 は~2 週間と最短。運用が高速化した証拠。

分析 4: 図郭 ID と国土基本図カバレッジ

狙い

図郭 ID (例: 03nf2444) の構造を解読し、3 dataset の地理カバレッジを量化する。 H1 (3 dataset の sheet 集合 disjoint) と H6 (図郭 prefix = 国土基本図地区記号) を検証。

手法 (リテラシレベル解説)

図郭 ID 「03nf2444」の構造:

• 先頭 4 桁 (03nf): 国土基本図 1/25000 の地区記号。 「03」は中国地方、「nf」は地区アルファベット (国土地理院命名)。
• 後続 4 桁 (2444): 1/25000 図郭内の細分位置 (子図郭の格子座標)。 1/25000 = 地図情報レベル 25000 を 100 分割した位置。

すなわち 1 業務は 多数の地図情報レベル 1000-2500 図郭を担当し、各図郭が 1 つの「resource」として DoBoX に登録される。1 図郭 = ~ 750m × 1.125km (地図情報レベル 2500 を 4 分割した子図郭)。

実装

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12

def parse_sheet_id(title):
    m = re.search(r"共通_(\\w{6,8})_\\d{2}_", title)
    return m.group(1).lower() if m else None

res_df["sheet_id"] = res_df["title"].apply(parse_sheet_id)
res_df["prefix"] = res_df["sheet_id"].apply(lambda s: s[:4] if s else None)

# 重複チェック
sheet_by_ds = {did: set(res_df[res_df.ds==did]["sheet_id"].dropna()) for did in [65, 1434, 1527]}
print("ds65 ∩ ds1434:", len(sheet_by_ds[65] & sheet_by_ds[1434]))
print("ds65 ∩ ds1527:", len(sheet_by_ds[65] & sheet_by_ds[1527]))
print("ds1434 ∩ ds1527:", len(sheet_by_ds[1434] & sheet_by_ds[1527]))

図と読み取り (図 4: 図郭プレフィックス × ds ヒートマップ)

なぜこの図か: 図郭 prefix (国土基本図地区記号) が 3 dataset でどう分布しているかを 1 枚で見る。 prefix が disjoint なら H1 を支持。

図 4: 図郭プレフィックス × dataset ヒートマップ

読み取り (重要発見):

• 図郭 prefix は3 種類: 03nf, 03od, 03oe。 これは中国地方の国土基本図 1/25000 地区記号と一致。
• ds65 (オリジナル) は 03nf, 03ng, 03of, 03og, 03nh 等多数の地区に分散 = 県全域概覆。
• ds1527 は少数の prefix (主に 03nf) に集中 = 三次・庄原は 03nf 系列のみ。
• 同じ prefix でも各 dataset で扱う図郭 ID (後続 4 桁) はほぼ disjoint = 同じ地区でも異なる小図郭を vintage 別に整備。

表 (図郭重複)

読み取り (H1 検証): ds65 ∩ ds1434 = 0, ds65 ∩ ds1527 = 0, ds1434 ∩ ds1527 = 0。 3 dataset 全部に同じ図郭が現れたのは 0 件のみ。 合計 420 異種図郭中、重複は0.0%。 H1 (図郭集合 disjoint) は支持。3 dataset は時系列の分割整理であり、同じ図郭を重複して公開しているのではない。

表 (prefix 頻度上位 10)

読み取り: prefix のトップは 03nf (785 件)。 これは国土地理院の 1/25000 地形図索引で三次・庄原を中心とする中国山地中央部に対応する地区記号。 広島県の航空レーザ測量は山間部 (中山間地, 急傾斜地) の高密度整備に重点があることが、prefix 分布から逆引きできる。

分析 5: プロダクト実データの解剖 (DXF/TXT/LAS)

狙い

各プロダクト (DXF / TXT / LAS) を実ファイル 1 件ずつ解析し、 データ構造・属性・要素数・値域を Before/After 形式で具体化する。 これは Q (multi-angle 活用) と K (Before/After 例) の要件への対応。

STEP 1: 等高線 DXF の分解

狙い: ds65 オリジナルの等高線 (DXF) を読込み、polyline と標高値の構造を可視化する。 DXF は AutoCAD のテキスト形式で、グループコード (整数) と値が交互に並ぶ独特の形式。

手法 (DXF パース入門)

DXF テキストの主要グループコード:

• 0: エンティティ種別 (例: LWPOLYLINE)
• 10, 20, 30: 頂点の X, Y, Z 座標 (繰り返し)
• 38: polyline の標高 (= 等高線レベル)
• 90: 頂点数

本記事は単純なテキスト解析: 行 i が「38」なら次行を float でパースして標高 list に追加。

実装

↑ L42_map_information.py 行 1711–1763

 1
 2
 3
 4
 5
 6
 7
 8
 9
1720
1721
1722
1723
1724
1725
1726
1727
1728
1729
1730
1731
1732
1733
1734
1735

def parse_dxf(path):
    txt = path.read_text(encoding="utf-8", errors="replace")
    lines = txt.splitlines()
    elevs = []
    xs, ys = [], []
    for i, ln in enumerate(lines):
        if ln.strip() == "38" and i+1 < len(lines):
            try: elevs.append(float(lines[i+1]))
            except: pass
        elif ln.strip() == "10" and i+1 < len(lines):
            try:
                v = float(lines[i+1])
                if abs(v) < 1e10: xs.append(v)
            except: pass
        elif ln.strip() == "20" and i+1 < len(lines):
            try:
                v = float(lines[i+1])
                if abs(v) < 1e10: ys.append(v)
            except: pass
    pl_count = txt.count("\\nLWPOLYLINE\\n") + txt.count("\\nPOLYLINE\\n")
    return {
        "polylines": pl_count, "elev_unique": len(set(elevs)),
        "elev_min": min(elevs), "elev_max": max(elevs),
        "xy_count": min(len(xs), len(ys))
    }

図 (図 5: 等高線 DXF 実データ)

図 5: 等高線 DXF 実データ — DS65 オリジナル vintage の典型サンプル

読み取り:

• (a) サンプル 03nf653_05_contour.dxf (16 MB) は 657 polylines, 463,847 頂点。1 図郭でこれだけ細かい等高線が含まれる。
• (b) 標高は 721–902m の 182 レベル。1m 等高線間隔 (= 地図情報レベル 1000) を反映。
• これは砂防基盤図作成用の精細データであり、1 図郭 = 1 ファイル = 17 MB 級のサイズが業務全体に蓄積される。

STEP 2: DEM グリッド TXT の分解

狙い: ds1434 の DEM グリッド (CSV テキスト) を pandas で読込み、 点群分類 (ground / non-ground) と標高分布を量化する。

手法

DEM TXT は単純な CSV: id, x, y, z, classification の 5 列。 classification = 1 が ground (地面), 0 が non-ground (建物・植生)。 航空レーザのベアアース DTM (= 樹木・建物を除いた地面標高) を作るためのフィルタリング情報。

実装

↑ L42_map_information.py 行 1758–1768

1758
1759
1760
1761

df = pd.read_csv(path, header=None, names=["id", "x", "y", "z", "ground"])
ground_pct = (df.ground == 1).mean() * 100
print(f"  rows={len(df):,}, X[{df.x.min()}-{df.x.max()}], "
      f"Z[{df.z.min()}-{df.z.max()}]m, ground={ground_pct:.1f}%")

図 (図 6: DEM TXT 実データ)

図 6: DEM グリッド TXT 実データ — DS1434 (03od7914 図郭, 2022年度計測)

読み取り:

• (a) ground (緑, 36%) と non-ground (赤, 64%) はほぼ同数。 = この図郭は森林地帯で植生・建物が多く、ベアアース化の効果が大きい。
• (b) 標高 Z は 844–898m (山間地)。色は連続値で表現。
• (c) 標高ヒストは平均 875m を中心に分布、標準偏差 11.7m。
• これは L40 標高図の原料 CSVそのもの。L40 がラスタ化した GeoTIFF を提供しているのに対し、 本データは X/Y/Z の生 CSV なので不規則格子として読める (1 行 = 1 計測点)。

STEP 3: LAS 点群の分解

狙い: ds1527 で新規追加された LAS 点群のヘッダを読み、生成ソフト・点数・XYZ 範囲を確認する。 LAS はバイナリ形式なので、struct でヘッダの固定オフセットから取り出す。

手法 (LAS 1.2 ヘッダ構造)

LAS 1.2 ファイルの先頭 227 バイトはヘッダで、固定オフセット位置に各情報がある:

• 0-4: signature (= "LASF")
• 24-26: バージョン (major, minor)
• 26-58: System Identifier (例: "LAStools")
• 58-90: 生成ソフト (例: "txt2las (version 210418)")
• 107-111: 点数 (uint32)
• 131-179: XYZ scale + offset (double × 6)
• 179-227: XYZ min/max (double × 6)

実装

↑ L42_map_information.py 行 1794–1818

 1
 2
 3
 4
 5
 6
 7
 8
 9
1803
1804

import struct
raw = path.read_bytes()
assert raw[:4] == b"LASF"
n_points = int.from_bytes(raw[107:111], "little")
x_max = struct.unpack("<d", raw[179:187])[0]
x_min = struct.unpack("<d", raw[187:195])[0]
y_max = struct.unpack("<d", raw[195:203])[0]
y_min = struct.unpack("<d", raw[203:211])[0]
z_max = struct.unpack("<d", raw[211:219])[0]
z_min = struct.unpack("<d", raw[219:227])[0]
print(f"n_points={n_points:,}, X={x_max-x_min:.0f}m, Z=[{z_min:.0f}-{z_max:.0f}]m")

図 (図 7: LAS 点群ヘッダ)

図 7: LAS 点群ヘッダ — DS1527 (03nf2634 図郭, 2023年度)

読み取り:

• システム識別子: "LAStools (c) by rapidlasso GmbH" = LAStools (rapidlasso 社の業界標準ツール)
• 生成ソフト: "txt2las (version 210418)" = txt2las で TXT から LAS に変換 (= DEM TXT がオリジナル, LAS は派生)
• 点数: 13,767 点 (このサンプルは 200m × 26.5m の小範囲)
• XY 範囲は EPSG:6671 ((X, Y) = (65691, -99000) → (65891, -98973))
• Z 範囲は 996–1010m (= 山地)。
• = ds1527 が「LAS 点群を新規公開」した実体は「TXT を txt2las でラップした派生」であり、 raw な打点情報 (intensity, return number, 反射回数) を保持する標準 LAS フォーマット。 これにより LiDAR 解析ソフト (CloudCompare, lastools, PDAL 等) で直接読み込めるようになった。

表 (3 プロダクトの実データ構造比較)

読み取り (Before/After): 1 図郭分の生データは合計 70+ MB (DXF 17MB + TXT 3MB + TIFF 56MB + LAS 0.3MB)。 1200 resource × 平均 50MB ≈ 60GB 相当のデータ量。 これが広島県の地形マッピング基盤の物理サイズ。

仮説検証と考察

仮説 H1〜H6 検証マトリクス

考察 1: 「地図情報」というタイトルの誤誘導

DoBoX で「地図情報」を検索すると 3 dataset がヒットするが、タイトルだけ見れば「地図一般」に思える。 しかし JMP 2.0 メタデータを開くと「砂防基礎調査に伴う航空レーザ測量及び撮影業務」と書かれており、 実体は土砂災害防止法第 4 条に基づく基盤整備の航空レーザ測量成果である。 これは行政データオープン化において「公開タイトル ≠ 内部メタタイトル」のギャップが起こりうることを示す好例。 研究者・教育者は必ずメタデータ XML を読む習慣をつけるべき。

考察 2: 3 vintage の整備戦略 (= 計測予算の使い方)

業務単位の解読から、広島県の航空レーザ測量整備戦略は以下のように変化したと読める:

1. 2014-2018 (ds65): 「とにかく県全域を 1 m 解像度で概覆」期。 10 業務並列で圏域 (流域) 単位の発注。1 業務 1-2 年。
2. 2018-2022 (ds1434): 「1m 残り + 50cm 移行」期。8 業務に細分化。 市町別発注に切替 (沿岸都市部の高解像度更新)。
3. 2023- (ds1527): 「中山間地の 50cm 完成」期。三次市・庄原市の 5 業務。 LAS 点群も公開 (= 派生だけでなく原料層も提供する方針転換)。

考察 3: L36-L40 (LiDAR 派生層) との階層関係 (再確認)

本記事の量的結果から、LiDAR ファミリの派生階層は以下のように整理できる:

本記事は段階 1 (原料層) と段階 2 (1 次成果)。L36-L40 は段階 3 (派生ラスタ)。L41 は段階 4 (応用ベクタ)。 4 段の派生階層を 1 つの LiDAR 計測から生成する点で、広島県の砂防 GIS 整備は整然とした垂直統合を達成している。

考察 4: H3 (LAS 公開) が示す「データ提供哲学の転換」

ds1527 で初登場した LAS 点群 (64 件) は、3 dataset 全体で5%の分量にすぎない。 しかしこれは「派生プロダクトだけ公開 → raw 点群も公開」という哲学転換の小さな第一歩である。 理由:

• Reproducibility: raw 点群があれば、利用者は独自の DEM 補間 (Inverse Distance Weight, Kriging, TIN 等) を行える。 従来は広島県が決めた補間法しか使えなかった。
• 透明性: 反射回数 (1st return = 樹冠 / last return = 地面) など、派生プロダクトでは失われる情報が保たれる。
• 応用拡大: 林業以外 (例: 都市熱環境、洪水シミュレーション、考古学的微地形検出) への利用が容易になる。

この方針が ds1527 以降の vintage で全面化すれば、広島県は全国でも先進的な LiDAR オープンデータ拠点になる。

考察 5: 教育的含意 (本記事の方法論)

本記事は「カタログメタデータだけで研究できる」例を示している。 1200 resource を実 DL せず, タイトル文字列と 23 件の小さな XML だけで、 広島県の航空レーザ測量整備史を量的に再構成した。 この手法は他の DoBoX シリーズ (例: 都市計画区域情報, 観測情報, 河川浸水想定) にも応用可能。 「データを読む前にメタデータを読め」。

発展課題

結果 X → 新仮説 Y → 課題 Z

(1) 結果 X: ds1527 で LAS 点群が新規公開された

新仮説 Y: 今後 ds1527 以降の vintage では、ds65/1434 の旧 vintage に対しても LAS 点群が遡及公開される。なぜなら、現状の ds65/1434 の DEM TXT はtxt2las で LAS 化可能 (= ds1527 のサンプル LAS が「txt2las (version 210418)」生成と書かれているのが証拠) なので、 広島県は技術的に簡単に LAS 化できるはず。

課題 Z: 1 年後の DoBoX 再スキャンで、ds65/1434 の resource 数が 大きく増えていないか確認 (新 LAS が追加されているか)。本記事の fetch_dataset_resources() を再実行して 1200 → 1500+ になっているかチェック。

(2) 結果 X: 図郭 prefix は 03nf, 03od, 03oe, 03ng 等の 3 種類

新仮説 Y: 各 prefix は国土地理院 1/25000 地形図索引の 1 地区に対応する。 広島県全域 (8479 km²) は ~30 地区でカバーされる。本データの prefix 数 (3) は 広島県の主要部分を覆っているが、瀬戸内海島嶼や西部の一部 prefix は欠けている可能性。

課題 Z: 国土地理院の地形図索引 (1/25000 索引図) と本記事の prefix 一覧を突き合わせ、 広島県内の欠損 prefix (= 未測量地区) を特定する。 具体的には L15 (行政区域 Shapefile) と本データの prefix から推定される bounding box を重ね、空白地帯を地図化する。

(3) 結果 X: ds1434 の業務 #1 は 2018 年計測なのに 2023 年に公開された (5 年ラグ)

新仮説 Y: ds1434 で 2018 年計測が遅れて公開されたのは、 「2018 年豪雨災害 (西日本豪雨)」の被災区域 (芦田川圏域) の計測データを優先的に再点検 していた可能性がある。 5 年遅れの公開は通常異例。

課題 Z: ds1434 業務 #1 (rid=93809, 芦田川圏域) のメタデータ XML 全文を読み、 processing や quality 要素に災害再点検の痕跡があるか確認。 また DoBoX の関連シリーズ「水害リスクマップ」「多段階の浸水想定図」(L08 等) との 公開タイミングを照合する。

(4) 結果 X: 1 業務あたり平均 18.3 図郭 (1 図郭 ~ 750m × 1.125km)

新仮説 Y: 1 業務 ~ 18 図郭 × 0.85 km² = ~ 16 km² の範囲を担当。 広島県 8479 km² ÷ 23 業務 = 369 km²/業務。これは仮説計算と ~3 倍ずれがある。 原因: (a) 業務間の重複, (b) 図郭サイズが地区により異なる, (c) 業務に未公開図郭 (シェアード) が含まれる。

課題 Z: メタデータ XML の geographicBoundingBox 要素から各業務の bbox を抽出し、 緯度経度面積を計算して業務別カバレッジ km²を確定。さらに L15 行政区域と空間結合して 業務 × 市町クロスを作る。これにより「どの市町に何 vintage の計測があるか」が市町別に明らかになる。

(5) 結果 X: ds65 (オリジナル) は圏域 (流域) 単位, ds1527 は市町単位で発注

新仮説 Y: 計測戦略の変化は「災害対応 (流域単位) → 都市計画連携 (市町単位)」の 行政方針シフトを反映する。具体的には 2018 年豪雨後に「市町別 BCP・流域治水計画」が法定化され、 これに合わせて計測単位も市町別に変わったのではないか。

課題 Z: 国交省「流域治水プロジェクト」(2020-) の市町別整備計画文書と、 本データ ds1434/1527 の業務名 (市町記載) を突き合わせ、 計画単位の符合を確認する。さらに L13 (都市計画基礎調査) との時系列同期をチェック。

(6) 結果 X: 業務 (XML) 数 = 23, resource 数 = 1200, 1 業務あたり ~52 resource

新仮説 Y: 1 業務 = メタ XML 1 + 各図郭 4-5 プロダクト × 10-15 図郭 ≈ 50 resource。 すなわち1 業務 = 「1 メタデータ + 図郭一式」のパッケージ単位。これを Python の dict で 明示的にモデル化すれば、業務単位の「データパッケージ」として再構成可能。

課題 Z: itertools.groupby またはnetworkx で 「業務 → 図郭 → プロダクト」の 2 段階ツリーを構築し、 資料一式を業務単位の zip ファイルに自動再パッケージするスクリプトを書く。 これにより「業務 #5 のフルセットを 1 操作でDL」が可能になる。