【BigQuery 】緯度経度から50mメッシュコードを算出する方法

はじめに

緯度経度から50mメッシュ(1/20メッシュ)を算出するBigQueryのSQLです。

ぐぐっても見つからなかったので自分で書きました。

ただ、未検証なので間違っている可能性があります。

もし使う場合はご了承のほどよろしくお願いいたします。

※別でやっているブログ(データ分析備忘録)でも記載しています。

BigQuery

WITH
  t0 AS (
  SELECT
    139.6917337 AS longitude,
    35.6895014 AS latitude ),
  t1 AS (
  SELECT
    latitude,
    longitude,
    /* FLOOR X 以下で最大の整数値を返します。*/ 
    /*round X のみが存在する場合、X を最も近い整数に丸めます。N が存在する場合、X を小数点以下の N 桁に丸めます。*/ 
    /*MOD(X, Y) Y による X の除算の剰余を返します。*/ 
    /* 1次メッシュ 辺の長さ:約80km 経度差:1度 緯度差:40分*/ 
    CONCAT( SAFE_CAST(SAFE_CAST(FLOOR(latitude*60/40) AS INT64) AS STRING), SAFE_CAST(SAFE_CAST(FLOOR(longitude)-100 AS INT64) AS STRING) ) AS mesh1d,
    /* 2次メッシュ 辺の長さ:約10km 経度差:7分30秒 緯度差:5分*/ 
    CONCAT( SAFE_CAST(SAFE_CAST(FLOOR(MOD(SAFE_CAST(ROUND(latitude*60, 7) AS NUMERIC), 40.0)/5) AS INT64) AS STRING), SAFE_CAST(SAFE_CAST(FLOOR(MOD(SAFE_CAST(longitude-100 AS NUMERIC), 1.0)*60/7.5) AS INT64) AS STRING) ) AS mesh2d,
    /*3次メッシュ 辺の長さ:約1km 経度差:45秒 緯度差:30秒 */ 
    CONCAT( SAFE_CAST(SAFE_CAST(FLOOR(MOD(MOD(SAFE_CAST(ROUND(latitude*60, 7) AS NUMERIC),40.0),5.0)* 60/30) AS INT64) AS STRING), SAFE_CAST(SAFE_CAST(FLOOR(MOD(MOD(SAFE_CAST(longitude AS NUMERIC),1.0)*60,7.5) * 60/45) AS INT64) AS STRING) ) AS mesh3d,
    /* 2分の1地域メッシュ 辺の長さ:約500m 経度差:22.5秒 緯度差:15秒 三次メッシュを緯線方向、経線方向に2等分してできる区域*/ 
    CONCAT( SAFE_CAST(SAFE_CAST(FLOOR(MOD(MOD(MOD(SAFE_CAST(ROUND(latitude*60, 7) AS NUMERIC),40.0),5.0),0.5)*4)*2+ FLOOR(MOD(MOD(MOD(SAFE_CAST((longitude-100) AS NUMERIC),1.0)*60,7.5),0.75)* 60/22.5)+1 AS INT64) AS STRING) ) AS mesh4d,
    /* 4分の1地域メッシュ 辺の長さ:約250m 経度差:11.25秒 緯度差:7.5秒 */ 
    CONCAT( SAFE_CAST(SAFE_CAST(FLOOR(MOD(MOD(MOD(MOD(SAFE_CAST(ROUND(latitude*60, 7) AS NUMERIC),40.0),5),0.5),0.25)*8)*2+ FLOOR(MOD(MOD(MOD(MOD(SAFE_CAST((longitude-100) AS NUMERIC),1)*60,7.5),0.75),0.375)*60/11.25)+1 AS INT64) AS STRING) ) AS mesh5d,
    /* 8分の1地域メッシュ 辺の長さ:約125m 経度差:5.625秒 緯度差:3.75秒 */
     CONCAT( SAFE_CAST(SAFE_CAST(FLOOR(MOD(MOD(MOD(MOD(MOD(SAFE_CAST(ROUND(latitude*60, 7) AS NUMERIC),40.0),5.0),0.5),0.25),0.125)*16)*2+ FLOOR(MOD(MOD(MOD(MOD(MOD(SAFE_CAST((longitude-100) AS NUMERIC),1)*60,7.5),0.75),0.375),0.1875)* 60/5.625)+1 AS INT64) AS STRING) ) AS mesh6d,
    /* 10分の1地域メッシュ 辺の長さ:約100m 経度差:4.5秒 緯度差:3秒 三次メッシュを緯線方向、経線方向に10等分してできる区域*/ 
    CONCAT( SAFE_CAST(SAFE_CAST(FLOOR(MOD(MOD(MOD(SAFE_CAST(ROUND(latitude*60, 7) AS NUMERIC),40.0),5.0),0.1)*20)*2+ FLOOR(MOD(MOD(MOD(SAFE_CAST((longitude-100) AS NUMERIC),1.0)*60,7.5),0.75)* 60/4.5)+1 AS INT64) AS STRING) ) AS mesh7d,
    /* 20分の1地域メッシュ 辺の長さ:約50m 経度差:2.25秒 緯度差:1.5秒 1/10メッシュを緯線方向、経線方向に2等分してできる区域*/ 
    CONCAT( SAFE_CAST(SAFE_CAST(FLOOR(MOD(MOD(MOD(MOD(SAFE_CAST(ROUND(latitude*60, 7) AS NUMERIC),40.0),5.0),0.1),0.05)*40)*2+ FLOOR(MOD(MOD(MOD(MOD(SAFE_CAST((longitude-100) AS NUMERIC),1.0)*60,7.5),0.75),0.075)* 60/2.25)+1 AS INT64) AS STRING) ) AS mesh8d
  FROM
    t0 )
SELECT
  latitude,
  longitude,
  CONCAT(mesh1d,mesh2d,mesh3d,mesh4d,mesh5d,mesh6d,mesh7d,mesh8d) AS mesh50m
FROM
  t1

参考サイト

タイトルとURLをコピーしました