three.jsと標高タイルで3Dマップを生成する

今日はクリスマスイブですね．皆様いかがお過ごしでしょうか．かく言う私はお布団にくるまり http://qiita.com/advent-calendar/2016/no-girlfriend を眺めています．お察しください．

ブラウザ上で現在位置を中心とした 3D マップを 国土地理院の標高タイル と three.js で作ることができたのでメモ．今月の 17，18 日に開催された 林業応援ハッカソン で作った．
github.com
流れとしてはまず現在位置を取得し，その取得した現在位置からタイル座標を求める．ちなみに私は初めてタイル座標というものを知った．このタイル座標や後々利用するピクセル座標については TrailNote : 座標の変換（世界座標、ピクセル座標、タイル座標、緯度・経度） が参考になる．次に求めたタイル座標を元に 国土地理院の標高タイル から標高データを取得し three.js で描画する．

タイル座標を求める

緯度と経度からタイル座標を求める．Showing pixel and tile coordinates  |  Google Maps JavaScript API  |  Google Developers を参考にした．

const createCoordinate = ( lat, lng ) => {
    const TILE_SIZE  = 256,
          ZOOM_LEVEL = 14
    return createInfoWindowContent( ( new google.maps.LatLng( lat, lng ) ), ZOOM_LEVEL, TILE_SIZE )
}

const createInfoWindowContent = ( latLng, zoom, tileSize ) => {
    const scale = 1 << zoom
    
    const worldCoordinate = ( ( latLng, tileSize ) => {
        const siny = Math.min( Math.max( ( Math.sin(latLng.lat() * Math.PI / 180) ), -0.9999 ), 0.9999 )
        return new google.maps.Point( tileSize * ( 0.5 + latLng.lng() / 360 ), tileSize * ( 0.5 - Math.log( ( 1 + siny ) / ( 1 - siny ) ) / ( 4 * Math.PI ) ) )
    } )( latLng, tileSize )

    const pixelCoordinate = new google.maps.Point( Math.floor( worldCoordinate.x * scale ), Math.floor( worldCoordinate.y * scale ) ),
          tileCoordinate  = new google.maps.Point( Math.floor( worldCoordinate.x * scale / tileSize ), Math.floor( worldCoordinate.y * scale / tileSize ) )

    return {
        lat            : latLng.lat(),
        lng            : latLng.lng(),
        zoomLevel      : zoom,
        pixelCoordinate: { x: pixelCoordinate.x, y: pixelCoordinate.y },
        worldCoordinate: { x: worldCoordinate.x, y: worldCoordinate.y },
        tileCoordinate : { x: tileCoordinate.x, y: tileCoordinate.y }
    }
}

標高データを取得する

タイル座標を使い，標高データを CSV で取得する．また，タイル座標に対応したテクスチャも取得することができる．ただし，国土地理院の標高タイルではズームレベルが 14 で固定されており，他のズームレベルでは標高データを取得することができないことに注意する．

標高データ
http://cyberjapandata.gsi.go.jp/xyz/dem/ZOOM_LEVEL/TILE_COORDINATE_X/TILE_COORDINATE_Y.txt
http://cyberjapandata.gsi.go.jp/xyz/dem/14/14547/6463.txt

テクスチャ
http://cyberjapandata.gsi.go.jp/xyz/relief/ZOOM_LEVEL/TILE_COORDINATE_X/TILE_COORDINATE_Y.png
http://cyberjapandata.gsi.go.jp/xyz/relief/14/14547/6463.png

取得した標高データは適当に処理して配列に変換する．

mapData = mapData.substring( 0, mapData.length-1 ).split( '\n' ).map( ( row ) => {
    return row.split( ',' ).map( ( height ) => {
        return parseFloat( height ) || -1
    } )
} )

あとは three.js の PlaneGeometry で元となる Plane を作成し，取得した標高データを元に各点の高さを変更する．

let geometry = new THREE.PlaneGeometry( 0, 0, xLength-1, yLength-1 )

for( let y=0; y<yLength; y++ )
    for( let x=0; x<xLength; x++ )
        geometry.vertices[x + ( y * xLength )] = { x: x, y: y, z: mapData[y][x] }

自分を中心とした 3D マップを生成する

現在位置からタイル座標を求める場合，自分の現在位置に対応したタイル座標を得ることができる．だがタイル上のどの位置にいても同じタイル座標を得てしまう．まあタイル座標というくらいだし当たり前なのだろうけど．またまた TrailNote : 座標の変換（世界座標、ピクセル座標、タイル座標、緯度・経度） を参考にした．自分の現在位置を元に 3D マップを作成しても自分を中心としたマップは作成されないのでいい感じにゴニョゴニョする．

仕方がないので差分をとって取得した標高データを合成する．タイル座標からピクセル座標を求めて現在位置のピクセル座標と比較する．自分がタイルのどの辺りにいるのかを調べて足りない標高データを取得するようにした．

もっとスマートにやりたかった．恥ずかしい．mori_hack/main.js at master · calmery/mori_hack · GitHub に頑張ったあとがある．

const targetTile = createCoordinate( position.lat, position.lng )

const difference = ( ( targetTile ) => {
    let t = targetTile.tileCoordinate,
        p = targetTile.pixelCoordinate

    let tileToLatLng = ( tx, ty, zoom ) => {
        let numOfTiles = Math.pow( 2, zoom ),
            x = tx / numOfTiles,
            y = ty / numOfTiles,
            latRadians = ( y - ( 1 / 2 ) ) / - ( 1 / ( 2 * Math.PI ) )
        return { lat: ( 2 * Math.atan( Math.exp( latRadians ) ) - Math.PI / 2 ) / Math.PI * 180,  lng: ( x - ( 1 / 2 ) ) / ( 1 / 360 ) }
    }

    let latLng = tileToLatLng( t.x, t.y, t.z ),
        tile = createCoordinate( latLng.lat, latLng.lng )

    return { x: p.x - tile.pixelCoordinate.x, y: p.y - tile.pixelCoordinate.y }
} )( targetTile )

const essentialTile = ( ( targetTile, difference ) => {
    let tx = targetTile.tileCoordinate.x,
        ty = targetTile.tileCoordinate.y

    let tiles = [ [], [] ],
        differenceFromBase = {}

    if( difference.x <= 128 && difference.y <= 128 ){
        tiles[0][0] = { x: -1, y: -1 }
        tiles[0][1] = { x:  0, y: -1 }
        tiles[1][0] = { x: -1, y:  0 }
        tiles[1][1] = { x:  0, y:  0 }

        differenceFromBase.x = 256 + difference.x
        differenceFromBase.y = 256 + difference.y
    } else if( difference.x > 128 && difference.y <= 128 ){
        ...
    } else if( difference.x <= 128 && difference.y > 128 ){
        ...
    } else if( difference.x > 128 && difference.y > 128 ){
        ...
    }

    return { tiles: tiles, difference: differenceFromBase }
} )( targetTile, difference )

足りない標高データを補うため合計で 4 回，標高データを取得する必要がある．配列に現在位置より求めたタイル座標からどの程度ずれているかの情報を格納している．上の図の場合，こんな配列になる．

[
    [{ x: -1, y: -1 }, { x:  0, y: -1 }],
    [{ x: -1, y:  0 }, { x:  0, y:  0 }]
]

この配列を元に標高データを取得し，あとはタイル座標から求めたピクセル座標と，現在位置のピクセル座標の差分を元に，足りない標高データを補った配列を作る．テクスチャも同様に 4 枚取得し canvas に貼り付けて切り取る．人吉駅を中心とするとこんな感じの 3D マップができる．
ちょっとわかりにくいが中心に人吉駅がある．Google Map で見るとちゃんと中心にきていることがわかるかも．

まとめ

Unity を使って OpenStreetMap Japan | 自由な地図をみんなの手に/The Free Wiki World Map とかと組み合わせると面白そうだと思った．mori_hack/sample at master · calmery/mori_hack · GitHub にすぐに動かせるものがあるのでお試しあれ．

参考

TrailNote : 座標の変換（世界座標、ピクセル座標、タイル座標、緯度・経度）
Showing pixel and tile coordinates  |  Google Maps JavaScript API  |  Google Developers
three.jsで3Dを作りたい！ - Qiita
地理院タイルとOpenStreetMapのデータを使ってthree.jsで立体地図を描画する - Qiita
国土地理院タイル地図をthree.jsで3D表示する