今回は最近諸事情により関わることになった「地域メッシュコード」の算出方法を紹介したいと思います。
詳細については総務省統計局によるこちらの PDF をご参照ください。 → 第1章 地域メッシュ統計の特質・沿革 – gaiyo1.pdf
地域メッシュとはその名の通り網の目のように地域を分割する規格で、その区域ごとに一意のコードが割り当てられています。以下のように網の目が細かくなるごとに桁数が増えていくのが特徴です :
この記事では、これらのメッシュコードがどのように算出されているのかを具体的なソースコードを見ながら追っていきます。
まずはメッシュコードが割り当てられる個々の区画の定義です :
case class GridCell(
meshCode: BigInt,
// セル西端の経度
longitude: Double,
// セル南端の緯度
latitude: Double,
)
(メッシュという言葉の指しているものが網なのか網の目なのかがイマイチ不明瞭だったので、ここでは意図的にセルという別名を与えています。)
次にこの GridCell を提供するためのクラスも定義しましょう :
case class CoordinateCells(longitude: Double, latitude: Double) {
// 一次メッシュ
val first: GridCell = ???
// 二次メッシュ
lazy val second: GridCell = ???
// 三次メッシュ
lazy val third: GridCell = ???
// 四次メッシュ
lazy val fourth: GridCell = ???
}
この ??? を埋めるのがこの記事の趣旨となります。
まずは全てのメッシュコードの基準となる一次メッシュからです。冒頭の PDF によると正式には第1次地域区画と呼ぶそうですが、この記事では、慣用的に広く使われている一次メッシュという呼称を用います。(二次以降も同様です。)
このサンプルが少々厄介で、以下のようにセルの間隔については明示されているものの… :
経度の間隔 | 緯度の間隔 |
---|---|
1度 | 40分 |
PDF には (東経 138, 北緯 36) の座標から 5438 というメッシュコードを求めるための謎の計算式?しか載っていません :
地域メッシュ・コードの付け方南端緯度×1.5 [36×1.5 = 54]西端経度の下2けた [138 → 38]
PDF を少し遡って以下の文章に辿り着くと :
標準地域メッシュが設定される範囲は,北緯 20 度から 46 度まで,東経 122 度から 154 度までの地域で,東端が南鳥島(当該第1次地域区画の地域メッシュ・コードは 3653),西端が与那国島(同 3622),南端が沖ノ鳥島(同 3036),北端が択捉島(同 6848)になります。
ここで一次メッシュが (東経 100 度, 北緯 0 度) という座標を基点としていることが読み取れます。 そう明示してほしいものですね…
つまりこの不思議な計算式の正体は、マジックナンバーを展開するとこんな形になっています :
(36 - 0) / (40 / 60) → 54
(138 - 100) / 1 → 38
あとは上記の計算式を書き下すだけです :
val first: GridCell = {
// 一次メッシュの緯度差は 40 分、基点は北緯 0 度
val degreeY = 40D / 60
val y = (latitude / degreeY).floor.toInt
// 一次メッシュの経度差は 1 度、基点は東経 100 度
val east = longitude.floor.toInt
val x = east - 100
GridCell(
meshCode = x + 100 * y,
longitude = east,
latitude = y * degreeY,
)
}
計算方法は一次メッシュとほぼ同じですが、基点が (東経 100 度, 北緯 0 度) から (一次メッシュの西端, 南端) に変わり (一次メッシュを縦横に 8 等分したものが二次メッシュなので) 緯度・経度の間隔も以下のように変わります。
経度の間隔 | 緯度の間隔 |
---|---|
7分30秒 | 5分 |
二次メッシュが必要になった場合にのみ計算が走るように lazy 修飾子が付いています :
lazy val second: GridCell = {
// 二次メッシュの緯度差は 5 分、基点は一次メッシュの南端
val degreeY = 5D / 60
val y = ((latitude - first.latitude) / degreeY).floor.toInt
// 二次メッシュの経度差は 7 分 30 秒、基点は一次メッシュの西端
val degreeX = 7.5D / 60
val x = ((longitude - first.longitude) / degreeX).floor.toInt
GridCell(
meshCode = first.meshCode * 100 + y * 10 + x,
longitude = first.longitude + x * degreeX,
latitude = first.latitude + y * degreeY,
)
}
三次メッシュは二次メッシュを更に縦横に 10 等分したものです。要領は全く同じです。
経度の間隔 | 緯度の間隔 |
---|---|
45 秒 | 30 秒 |
lazy val third: GridCell = {
// 三次メッシュの緯度差は 30 秒、基点は二次メッシュの南端
val degreeY = 30D / 60 / 60
val y = ((latitude - second.latitude) / degreeY).floor.toInt
// 三次メッシュの経度差は 45 秒、基点は二次メッシュの西端
val degreeX = 45D / 60 / 60
val x = ((longitude - second.longitude) / degreeX).floor.toInt
GridCell(
meshCode = second.meshCode * 100 + y * 10 + x,
longitude = second.longitude + x * degreeX,
latitude = second.latitude + y * degreeY,
)
}
四次メッシュ以降は少しだけ特殊です。三次メッシュを 4 等分して (南西, 南東, 北西, 北東) の順に割り当てられた 1,2,3,4 の番号がメッシュコードの末尾になります。
経度の間隔 | 緯度の間隔 |
---|---|
22.5 秒 | 15 秒 |
ちなみに PDF によると、これは2分の1地域メッシュというのが正式な呼称のようです。
lazy val fourth: GridCell = {
// 四次メッシュの緯度差は 15 秒、基点は三次メッシュの南端
val degreeY = 15.0D / 60 / 60
val y = ((latitude - third.latitude) / degreeY).floor.toInt
// 四次メッシュの経度差は 22.5 秒、基点は三次メッシュの西端
val degreeX = 22.5D / 60 / 60
val x = ((longitude - third.longitude) / degreeX).floor.toInt
// 四次メッシュの末尾は (南西: 1, 南東: 2, 北西: 3, 北東: 4)
val code = (x + 1) + y * 2
GridCell(
meshCode = third.meshCode * 10 + code,
longitude = third.longitude + x * degreeX,
latitude = third.latitude + y * degreeY,
)
}
さらに粒度の細かい4分の1地域メッシュや8分の1地域メッシュについては全く同様の計算になるので割愛します。
実行してみましょう。
val cells = CoordinateCells(
longitude = 139.733231,
latitude = 35.680916
)
println(
cells.first.meshCode,
cells.second.meshCode,
cells.third.meshCode,
cells.fourth.meshCode)// 5339, 533945, 53394518, 533945184
無事に計算できているようです。
これらのメッシュコードについてはマスタデータのようなものは見つからなかったため、下記の複数個のサービスで同一の結果が得られることを確認させて頂きました 🙏
地域メッシュは日本の座標しか考慮されていません。また周辺ライブラリも少なく、今回のように簡単な計算であっても自前で用意する必要が多々あります。もしも特別な要件がなければ Geohash や S2 Geometry といった規格の採用を検討しましょう………
Mobility Technologies では共に日本のモビリティを進化させていくエンジニアを募集しています。話を聞いてみたいという方は、是非 募集ページ からご相談ください!