世界の測量

Sibling of "Relevant, Timely, and Accurate, " but much lighter and shorter ※自らの所属する組織の見解を示すものでない

Tangram で国土地理院ベクトルタイル提供実験のデータを見る

Tangram で国土地理院ベクトルタイル提供実験のデータ(基盤地図情報)を見るコードを作成してみました。

f:id:hfu:20150520053831p:plain

下記 gsi.html と gsi.yaml からなっています。ファイルは、GitHub レポジトリ https://github.com/handygeospatial/simple-demo においています。

なお、ローカルのファイルシステムにあるものをブラウザでそのまま開いても動作しないかもしれません。python -m SimpleHTTPServer するなどして、HTTPでサーブししたものを開く必要があります。

なお、完成品を http://handygeospatial.github.io/simple-demo/gsi.html にも置いています。

gsi.html

<!doctype html>
<html>
<head>
<meta charset="utf-8">
<meta name="viewport" content="width=device-width, initial-scale=1.0, maximum-scale=1.0, user-scalable=no">
<title>Minimalist Tangram</title>
<link rel="stylesheet" href="lib/leaflet/leaflet.css" />
<style>
body {margin: 0px; border: 0px; padding: 0px;}
#map {height: 100%; width: 100%; position: absolute;}
</style>
</head>
<body>
<div id="map"></div>
<script src="lib/leaflet/leaflet.js"></script>
<script src="lib/leaflet-hash.js"></script>
<script src="https://mapzen.com/tangram/0.1/tangram.min.js"></script>
<script>
var map = L.map('map', {maxZoom: 18, minZoom: 18});
var layer = Tangram.leafletLayer({
    scene: 'gsi.yaml',
    attribution: '国土地理院ベクトルタイル提供実験'
});
layer.addTo(map);
map.setView([36.08600, 140.11373], 18);
var hash = new L.Hash(map);
</script>
</body>
</html>

gsi.yaml

こちらでデータの処理(スタイルなど)を指定します。文法は https://github.com/tangrams/tangram/wiki が参考になると思います。

cameras:
    camera1:
        type: isometric

lights:
    light1:
        type: point
        direction: [0, 0, -10800]
        origin: camera
        diffuse: 1.
        ambient: .3
        specular: .2

sources:
    gsi:
        type: GeoJSONTiles
        url:  //cyberjapandata.gsi.go.jp/xyz/experimental_fgd/{z}/{x}/{y}.geojson
        max-zoom: 18

layers:
    blda:
        data: {source: gsi}
        filter:
            class: BldA
        draw:
            polygons:
                order: 2
                color: '#ffe6be'
        3d-buildings:
            draw:
                polygons:
                    extrude: function() {return true;}
              
    xray:
        data: {source: gsi}
        draw:
            lines:
                order: 0
                color: '#040'
                width: .5

今後

yaml ファイルをさらにカスタマイズしていくことで、地図表現は改善していくことができます。

地図は「まちあるきファシリテーションツール」であり「意見交換プラットフォーム」(「国土交通」都市防災)

「国土交通」No.131(2015.4-5)「特集 都市防災」に「まちを歩いて防災マップをつくろう!!」という1ページ記事がある。

「国土交通」No.131(2015.4-5) 特集 都市防災 - 国土交通省

これを読んで思ったのだが、「マップづくり」という活動の主目的は、地図づくりとは限らないし、地図づくりが主目的ではない「マップづくり」こそが重要なのかもしれない。つまり「まちあるきをファシリテートすること」だとか「意見交換のプラットフォームとなること」をこそ主目的とする「マップづくり」活動について、よく考察しておく必要があるのではないかと思う。

全国児童生徒地図優秀作品展の構造

第18回全国児童生徒地図優秀作品展の構成作品展の一覧を整理することを通じて、当該作品展の構造を考察する。

行事名 地区 主催
札幌市児童生徒社会研究作品展 札幌地区 札幌市教育地図研究会
私たちの身の回りの環境地図作品展 全国地区 環境地図教育研究会北海道教育大学・北海道地図ほか旭川関係)
仙台市中学校生徒地図作品展 仙台地区 仙台市中学校社会科研究会
いばらき児童生徒地図作品展 茨城地区 いばらき児童生徒地図研究会
多摩市身のまわりの環境地図作品展 多摩地区 多摩市身のまわりの環境地図作品展運営委員会
とやまみんなの地図作品展 富山地区 とやま地図作品研究会
児童生徒地図作品展 岐阜地区 岐阜県図書館(世界分布図センターを併設)
地域の安全安心マップコンテスト 全国地区 立命館大学歴史都市防災研究所
神戸市小学校社会科作品展 神戸地区 神戸市小学校教育研究会社会科部会
あこう絵マップコンクール 赤穂地区 あこう絵マップコンクール実行委員会
鳥取県児童生徒地域地図発表作品展 鳥取地区 鳥取県地域社会研究会
しまね児童・生徒地図作品展 島根地区 島根地理学会島根大学
地図ならびに地理作品展 広島地区 広島県地理作品展運営委員会
児童・生徒の地図作品展 徳島地区 徳島地理学会鳴門教育大学
大分市児童地図作品展 大分地区 大分市小学校教育研究会社会科部会
社会科児童作品展 滋賀地区 滋賀県小学校教育研究会社会科部会

Geo+Webの新局面を象徴する「Turf for local gov」について


Turf for local gov: potholes and parking meters | Mapbox

ウェブ地図が単なる「プレースマークに吹き出し」だけの段階を脱する、それはウェブ地図の未来像の中に確実に描かれていることなのですが、その方向に力が入り始めてきました。呼応する分野の選択も、さすがです。

属性表示を吹き出し表示に持たせず、地図の外のしっかりデザインされた場所に持たせる、それは技術的には大したことではなく、また、以前からも行われたことですが、ウェブ地図が趣味ではなくて生活に使われるためには、そういった方向のデザインが死活的に重要です。そういうところをしっかり行っていることが重要です。

もうひとつは、「ウェブ地図は簡易GISである」という古い偏見を脱して、ウェブにまともなジオ技術を入れることを怠っていないこと。ウェブブラウザは制限された計算機環境というよりはむしろ、ITのメジャープレイヤーがしのぎを削って新しい計算機の使い方を切り開いていっている、最新の計算機環境です。ここでGISが数十年にわたって育ててきた技術を花開かせること、その仕事をしっかりと行っています。ちなみにデータのフローがデスクトップとウェブでは大きく違う、というのが致命的な違いで、データフローについてデスクトップと同じ発想をしていると、根本的に誤ります。そこを誤っていない、というのがこの分野の実力者共通の特徴だと思っています。

「UIないし可視化の問題」と、「ブラウザサイドのGIS処理の問題」が、方向を開くために外さなければならない「かんぬき」なのですが、かれらの動きは、このかんぬきを2つ同時に外す動きに見えますね。
2つ同時に外すことで、革新の外見を得られますし、他が追いつく可能性が減少する。オープンソースの手法を取りながら他の追従を許さないには、相応の実力が必要ですが、今の彼らは実力者の集中に成功しているように見えます。

自分がこの技術潮流の中でどう動くか考える以前に、ジオ界がこれだけまともに技術の発展を迎えられたことを見るのが何より嬉しいです。
ジオがこの方向に進むこと、これがここ1年半くらいで大切な、ジオのブレークスルーになり得ると思います。
こういった技術を使えるプレイヤーを、早期に一者でも増やせなければ、我々はジオの発展の傍観者になってしまいます。

もうひとつは、仕事のかたち、の問題ですね。この先4年くらいは、技術が事業の形、ひいては制度をどんどん変えてくる状況になるはず。2000年代に制度が技術を誘導できると考えた(XML, UML, CALS etc.)のと真逆の動きですね。それに対応できるためには、技術と事業を同時に考える形、動きの早い技術潮流の中に継続的に身を浸しながら事業の軌道修正をかけられる形を取れるひとが強いです。(勿論、形以上に中身が大切です。)クラシック(オーソドックス)なコンサルティングや受託開発の形では、この形を取るのは非常に困難です。ソーシャルコーディングを使える事業体、というのは、既存の慣行のなかでとり得る一つの良好な形だと思います。

FOSS4G Advent Calendar 2014 に投稿

FOSS4G Advent Calendar 2014 http://qiita.com/advent-calendar/2014/foss4goss4g の 2014-12-15 のエントリとして、http://handygeospatial.github.io/foss4g-advent-calendar-2014/ を作成した。

geojson-vt のサンプルサイトを作成したという話であるが、その話にソーシャルコーディングツールを使った作業のチュートリアルを絡めている。geojson-vt の目的の分析、といったような話は割愛している。ソーシャルコーディングツールを使うプロセスの普及を、分析結果の共有に優先したことになる。

Hadoop を使った基盤地図情報の GeoJSON タイル変換 Pt. 1

本物の MapReduce を使った基盤地図情報の GeoJSON タイル変換の流れができた。ソースデータごとに GeoJSONL に変換するステップの後、一度に MapReduce を使って変換する方法である。ソースデータの一部が変更された時に、変更を受けていないソースデータに対応する GeoJSONL 変換をしなくて済むようになっている。

Hadoop のインストール(OS X

MapReduce といっても、特段分散計算機環境を設定しない、HomeBrew のデフォルト状態の Hadoop を使っている。それでも、特にソースデータが大規模の場合に、Unix sort よりは高速になる。Hadoop をきちんと設定して高速化する作業は、あとの楽しみにとっておく。

$ brew install mapreduce

現在、手元では hadoop の 2.5.1 がインストールされている。
OS XJava はターミナルで日本語を CP932 で出力して文字化けさせてしまう癖がある。これを回避して日本語を綺麗に表示するためには、次のようにする。Haddop のステータス出力には日本語が混じるので、設定しておいたほうがよい。

$ export _JAVA_OPTIONS="-Dfile.encoding=UTF-8"

また、大規模なデータセットを変換する場合には、heap space 不足でエラーを出す場合があるので、これを予防するために、次のように環境変数を設定しておくと良い。

$ export HADOOP_CLIENT_OPTS="-Xmx2048m $HADOOP_CLIENT_OPTS"

これらは、.bash_profile あたりに書き込んでおくとよい。

ディレクトリの構成

ディレクトリの構成は、次のようにしている。

. -+- source/ : 基盤地図情報XMLファイルを格納
   +- input/ : convert.rb を使って基盤地図情報をGeoJSONL.gzにしたものを格納
   +- output/ : hadoop を使って、GeoJSONLをタイルごとにreduceしたファイルを格納
   +- dst/ : deploy.rb を使って、output をタイルに展開したものを格納
   +- Rakefile : タスクを記述
   +- convert.rb
   +- reduce.rb
   +- deploy.rb

Rakefile

task :default do
  sh "rm -r output" if File.directory?('output')
  Dir.glob('source/*.xml') {|path|
    dst_path = "#{path.sub('source', 'input').sub('xml', 'geojsonl.gz')}"
    next if File.exist?(dst_path)
    sh "ruby convert.rb #{path} | gzip -c > #{dst_path}"
  }
  sh "hadoop jar /usr/local/Cellar/hadoop/2.5.1/libexec/share/hadoop/tools/lib/hadoop-streaming-2.5.1.jar -input input -output output -mapper cat -reducer reduce.rb"
  sh "ruby deploy.rb"
end

convert.rb

現在のところ、数値標高データ(DEM5A及びDEM10B)についてXMLデータをGeoJSONLに変換する部分ができている。

# coding: utf-8
require 'json'

# Meshcode is probably Japanese English.
module Meshcode
  def self.width(code)
    case code.size
    when 8
      45.0 / 60 / 60
    else
      raise 'not implemented.'
    end
  end

  def self.height(code)
    case code.size
    when 8
      30.0 / 60 / 60
    else
      raise 'not implemented.'
    end
  end

  def self.lefttop(code)
    case code.size
    when 6
      [(code[2..3].to_f + code[5].to_f / 8) + 100, 
       (code[0..1].to_f + (code[4].to_f + 1) / 8) * 2 / 3]
    when 8
      [(code[2..3].to_f + code[5].to_f / 8 + code[7].to_f/ 80) + 100, 
       (code[0..1].to_f + code[4].to_f / 8 + (code[6].to_f + 1) / 80) * 2 / 3]
    else
      raise 'not implemented.'
    end
  end
end

module Math
  def self.sec(x)
    1.0 / cos(x)
  end
end

class DEM
  def parse(params)
    (left, top) = Meshcode::lefttop(params[:meshcode])
    skip = true
    count = 0
    File.foreach(params[:path], encoding: 'cp932') {|l|
      if l.include?('<gml:tupleList>')
        skip = false
        next
      elsif l.include?('</gml:tupleList>')
        skip = true
        next
      elsif !skip
        (i, j) = [count % @n_lng, count / @n_lng]
        lng = left + @d_lng * (i + 0.5)
        lat = top - @d_lat * (j + 0.5)
        lat_rad = lat * 2 * Math::PI / 360
        n = 2 ** params[:z]
        x = n * ((lng + 180) / 360)
        y = n *
          (1 - (Math::log(Math::tan(lat_rad) + Math::sec(lat_rad)) /
          Math::PI)) / 2
        x = x.to_i
        y = y.to_i
        (type, height) = l.encode('UTF-8').strip.split(',')
        next unless type
        f = {:type => 'Feature', 
          :geometry => {:type => 'Point', :coordinates => [lng, lat]},
          :properties => {
            :type => type, :height => height.to_f,
            :datePublished => params[:datePublished]
        }}
        print "#{params[:z]}/#{x}/#{y}.geojson\t#{JSON::dump(f)}\n"
        count += 1
      end
    }
  end
end

class DEM5A < DEM
  def initialize
    @n_lng = 225
    @n_lat = 150
    @d_lng = 1.0 / 80 / @n_lng
    @d_lat = 2.0 / 3 / 80 / @n_lat
  end
end

class DEM10B < DEM
  def initialize
    @n_lng = 1125
    @n_lat = 750
    @d_lng = 1.0 / 8 / @n_lng
    @d_lat = 2.0 / 3 / 8 / @n_lat
  end
end

ARGV.each {|path|
  next unless /xml$/.match path
  r = File.basename(path, '.xml').split('-')
  r.pop if r[-1].size == 4
  next unless r.shift == 'FG'
  next unless r.shift == 'GML'
  datePublished = r.pop.insert(4, '-').insert(7, '-')
  type = r.pop
  meshcode = r.join
  params = {:path => path, :type => type, :meshcode => meshcode, :z => 18,
    :datePublished => datePublished}
  case type
  when 'DEM5A'
    Kernel.const_get(type).new.parse(params)
  when 'dem10b'
    DEM10B.new.parse(params)
  else
    # print "converter for #{type} not implemented.\n"
  end
}

reduce.rb

hadoop の reducer に当てて、タイルごとに GeoJSONL のカンマ区切りリストを作成するスクリプト。reduce.rb で GeoJSON に組み立てることも可能であるが、reduce を多段に実施する場合を考えて、FeatureCollection を当てて GeoJSON に組み立てる作業は deploy.rb に移譲している。

#!/usr/bin/env ruby
require 'fileutils'
require 'json'

last = nil
geojsonl = nil

def write(geojsonl, path)
  print "#{path}\t#{JSON.dump(geojsonl)[1..-2]}\n"
end

while gets
  r = $_.strip.split("\t")
  current = r[0]
  if current != last
    write(geojsonl, last) unless last.nil?
    geojsonl = []
  end
  geojsonl << JSON.parse(r[1])
  last = current
end
write(geojsonl, last)

deploy.rb

上記 reduce の結果をファイルシステムにタイルセットとして配備するスクリプト。ここで Feature のカンマ区切りリストを FeatureCollection に組み上げ、valid な GeoJSON ファイルにしている。

require 'fileutils'
Dir.glob('output/part*') {|path|
  File.foreach(path) {|l|
    (path, geojsonl) = l.split("\t")
    path = "dst/#{path}"
    FileUtils.mkdir_p(File.dirname(path))
    File.open(path, 'w') {|w|
      geojson = <<-EOS
{"type": "FeatureCollection", "features": [#{geojsonl}]}
      EOS
      print "#{geojson.size} characters to #{path}\n"
      w.print geojson
    }
  }
}

現在までのテスト状況

DEM5A と DEM10B データについて、一次メッシュ数個程度をまとめて変換する作業が無事終了することが分かっている。UNIX sort を使う場合に比べて高速であるようだが、ディスクを大量に使用することに留意する必要がある。

今後の予定

今回は、数値標高データのポイントデータタイル変換を行ったが、今後、ベクトルデータの変換も行いたい。

8月16日からの大雨(広島市内)正射画像をCesium 1.0で見る

地理院タイル「8月16日からの大雨(広島市内)正射画像(2014年8月28日撮影)」をCesium 1.0で見るサイトを作成しました。
サイト:8月16日からの大雨(広島市内)正射画像(2014年8月28日撮影)

f:id:hfu:20140831052612p:plain