ジオイドモデルの内挿計算 — 格子からNを求める

⚠ 本記事は、「ジオイド2024のISG形式とGML形式 — 楕円体高→標高で押さえるデータの違い」で紹介したジオイドデータの計算版です。前回はISG・GML形式の中身（ファイル構造の違い）を扱いましたが、今回は「そのファイルから、任意地点のジオイド高 $N$ をどう計算するか」という内挿処理そのものを掘り下げます。

考え方自体は、グリッド補正で使うバイリニア補間と同じです。TKY2JGD・PatchJGDが「補正量」を格子に持つのに対し、ジオイドモデルは「ジオイド高 $N$」を格子に持つ、という違いだけで、内挿の数式は共通しています。

おさらい：H = h − N の N はどこから来るか

GNSS測位が直接出せるのは楕円体高 $h$（GRS80楕円体面からの高さ）です。私たちが日常使う標高 $H$（東京湾平均海面からの高さ）を得るには、その地点のジオイド高 $N$（ジオイド面と楕円体面のズレ）を引く必要があります。

$$
H = h – N
$$

楕円体高から標高へ（ジオイド2024）で紹介したこの関係式の、$N$ を求める部分が今回のテーマです。

楕円体高・ジオイド高・標高の関係

ジオイド格子の構造

ジオイド2024（JPGEO2024）は、日本全国を覆う緯度経度の格子上に、あらかじめ計算されたジオイド高 $N$ の値を格納したモデルです。緯度方向はおおむね北緯15度〜50度、経度方向はおおむね東経120度〜160度の範囲をカバーし、格子の分解能は約1〜2km相当（分単位のメッシュ）で構成されています（具体的な分解能・範囲はISG/GMLファイルのヘッダに記載されており、バージョンによって多少の差があります）。

ジオイド高2024の分布図

この格子1つ1つの交点に「その地点でのジオイド高」が数値として格納されており、対象地点がちょうど交点に乗ることはまずないため、周囲の格子点から補間して求める必要があります。

内挿計算 — バイリニア＋端部処理

任意地点のジオイド高 $N$ を求める基本手順は、グリッド補間の実装で扱ったバイリニア補間と同じ流れです。

対象地点の緯度経度から、それを取り囲む格子の4隅（南西・南東・北西・北東）を特定する
4隅のジオイド高 $N_{00}, N_{10}, N_{01}, N_{11}$ を格子データから取得する
対象地点が格子内のどの位置にあるかを、比率 $u, v$（0〜1）で正規化する
次の式でジオイド高を補間する

$$
N(u, v) = N_{00}(1-u)(1-v) + N_{10}\,u(1-v) + N_{01}(1-u)v + N_{11}\,uv
$$

グリッド補正のバイリニア補間との違いは、格子に格納されている値の意味（補正量 vs ジオイド高）だけで、計算そのものは同一です。

端部処理の注意点

ジオイドモデルは日本全国をほぼ覆っていますが、沿岸部や離島の一部では格子の4隅すべてが有効な値を持たないケースがあります。この場合、有効な点だけで代用するか、内挿不可としてエラーを返すかは実装判断になります。全テストポイントを国土地理院の変換サイトと比較するで扱った検証でも、境界域は特に注意が必要な領域として扱っています。

ISG/GML形式それぞれの読み方（要約）

ジオイドデータの配布形式には主にISG形式（International Society for the Geoid標準）とGML形式（国土地理院配布）があり、ヘッダ情報の持ち方やバイナリ/テキストの違いなど、ファイル構造そのものの解説は前回記事に譲ります。今回の内挿計算は、どちらの形式でパースしたとしても「緯度経度→格子インデックス→4隅の値→バイリニア補間」という流れ自体は変わりません。

精度検証について

内挿誤差は、格子の分解能が細かいほど、またジオイド面の変化が滑らかであるほど小さくなります。ジオイド面は地殻変動補正の格子に比べて空間的な変化が緩やかな性質を持つため、標準的な分解能であれば内挿による誤差は実務上小さい水準に収まるよう設計されています。GeoCoreJPでの実装値が国土地理院の公式計算結果とどの程度一致するかについては、個別地点ごとの厳密な数値比較は継続的な検証項目としており、断定的な誤差量はここでは示しません。気になる方は国土地理院のジオイド高計算サービスとの突き合わせを推奨します。