OpenStreetMapの陸地のレンダリングは一般にはOSMデータから作り出されたシェイプファイル[3]を使う。
自作OSM地図では、現在は Mapnik を使わないため、シェイプファイルのままでは使えない。 これまでは、このシェイプファイルをPostgreSQL/PostGISデータベースに一旦インポートして、 OpenGIS の WKT(Well Known Text)形式でエクスポートしたものを使ってきた。
手間はかかるが、プログラム的には楽である。
今回は、より深く OpenStreetMap を理解するために、 OSMデータから直接陸地ポリゴンデータを取り出すことを試みた。
当面、世界地図は陸地の描画だけであるが、いずれ、国境、国名、大都市名なども描画したい。 Overpass APIには制約があるため、OSMデータから抽出する方がよい。 このためにも、OSMデータの扱いに慣れておきたい。
最初に、nodeセクションをスキップして wayセクションから 海岸線データを得る。europe と asia にまたがる海岸線データは同じものが両者にある。
europe の海岸線データだけでは、asia にまたがる場合は単独で繋ぎ合わせても閉ループとはならない。
繋ぐ合わせは行わず、wayオブジェクト単位の海岸線データをファイル出力しておく。 xml形式であっても、世界全体で 2GB弱の大きさであるから、CSVファイルでいいだろう。 行先頭を way id として、その後を node id 列とすればよい。
osmconvert で解凍した xml(テキスト)をパイプで受ける。 全体のパフォーマンスは osmconvert の解凍時間に支配される。日本地図データでは 10分弱である。
nodeセクションは空読みでいいが、解凍に時間がかかる。
世界全体では合計3~4時間程度かかった。同時に2、3の並列実行とした。
osmconvert d:/downloads/asia-latest.osm.pbf | java -Dfile.encoding=UTF-8 -Xmx5000m LandPolygon asia
// xmlファイルは標準入力(パイプ)から static void parseWay(String name) { char quote = '"'; String pathout = "c:/osm/" + name + "_coastline.csv"; List<Long> listNodes = new ArrayList<>(); long idWay = 0; try (BufferedReader br = new BufferedReader(new InputStreamReader(System.in)); FileWriter fw = new FileWriter(new File(pathout)) ) { String line; boolean fCoastline = false; boolean fWaySection = false; while ((line = br.readLine()) != null) { int ix = line.indexOf('<'); if (ix < 0) continue; String s = line.substring(ix); // 行先頭のスペースを削除する if (fWaySection && s.startsWith("<nd ")) { int bgn = 9; int end = s.indexOf("\"", bgn); String ref = s.substring(bgn, end); long idRef = Long.parseLong(ref); listNodes.add(idRef); } else if (fWaySection && s.startsWith("<tag ")) { int kbgn = 8; int kend = s.indexOf(quote, kbgn); String key = s.substring(kbgn, kend); int vbgn = kend + 5; int vend = s.indexOf(quote, vbgn); String val = s.substring(vbgn, vend); if (key.equals("natural") && val.equals("coastline")) { fCoastline = true; } } else if (s.startsWith("<way ")) { // wayオブジェクト始まり if (!fWaySection) fWaySection = true; listNodes.clear(); fCoastline = false; int bgnId = 9; int endId = s.indexOf(quote, bgnId); String strId = s.substring(bgnId, endId); idWay = Long.parseLong(strId); } else if (fWaySection && s.startsWith("</way")) { // wayオブジェクト終わり if (fCoastline) { fw.write("" + idWay); for (int i = 0; i < listNodes.size(); i++) { fw.write("," + listNodes.get(i)); } fw.write("\n"); } } } } catch (Exception e) { e.printStackTrace(); } }
node の id をキー、lon, lat を値として HashMap で管理するのはノートパソコンでは メモリが足りない。記事[7] に示したハッシュ法を採用した。
LandPolygon は nodeセクションが終わると、実行を終了するため、 前段の osmconvert は write error で終了するが、問題はない。
c:\gisa>osmconvert d:/downloads/antarctica-latest.osm.pbf | java -Dfile.encoding=UTF-8 -Xmx5000m LandPolygon -node antarctica 実行時間=0.6分 osmconvert Error: write error.
static void parseNode(String name) { lons = new int[MAXSIZE]; lats = new int[MAXSIZE]; String file = "c:/osm/" + name + "_coastline.csv"; try (BufferedReader br = new BufferedReader(new FileReader(file))) { String line; while ((line = br.readLine()) != null) { String[] nids = line.split(","); for (int n = 1; n < nids.length; n++) { put(Long.parseLong(nids[n]), 0, 0); } } } catch (IOException e) { e.printStackTrace(); return; } char quote = '"'; int cntFound = 0; String pathOut = "c:/osm/" + name + "_idlonlat.dat"; // 標準入力 try (BufferedReader br = new BufferedReader(new InputStreamReader(System.in)); DataOutputStream dos = new DataOutputStream(new FileOutputStream(pathOut))) { String line; while ((line = br.readLine()) != null && ixNode < num_nodes) { int ix = line.indexOf('<'); if (ix < 0) continue; String s = line.substring(ix); // 行先頭のスペースを削除する if (s.startsWith("<node ")) { // nodeオブジェクト int bgnId = 10; // <node id="31252942" lat="35.6528795" lon="139.7573991" ... > int endId = s.indexOf(quote, bgnId); String strId = s.substring(bgnId, endId); long idNode = Long.parseLong(strId); if (get(idNode) > 0) { // 海岸線ノード cntFound++; int bgnLon = s.indexOf("lon=", endId); int endLon = s.indexOf(quote, bgnLon+5); String strLon = s.substring(bgnLon+5, endLon); int bgnLat = s.indexOf("lat=", endId); int endLat = s.indexOf(quote, bgnLat+5); String strLat = s.substring(bgnLat+5, endLat); int lon = (int)(Double.parseDouble(strLon)*E7); int lat = (int)(Double.parseDouble(strLat)*E7); dos.writeLong(idNode); dos.writeInt(lon); dos.writeInt(lat); } } else if (s.startsWith("<way ")) { // wayオブジェクト break; } } } catch (Exception e) { e.printStackTrace(); } }
同時に2、3の並列実行で約2時間かかった。
海岸線データ(wayデータ)を繋ぎ合わせ、 ノード番号を座標値に変えて、ポリゴンデータとする。
経度 -180度は 180度と同じであることを考慮して繋ぎ合わせている。
プログラムソースコード: LandPolygon.java
ノード座標値が得られないものが1点あった。除外しても問題はないことを JOSM で確認した。
閉ループにならない海岸線がひとつだけあったが、 ノードは3点だけで、-180度の位置であるから、無視して問題はない。
全体で約73万ポリゴンにもなるのにエラーがわずかこの二つに過ぎないのが驚きである。
procCoastLine antarctica mapCoastLines エントリ数=0 procCoastLine africa mapCoastLines エントリ数=66 procCoastLine europe mapCoastLines エントリ数=74 procCoastLine asia mapCoastLines エントリ数=4 procCoastLine north-america Error 座標値が得られない Id=7670000000 mapCoastLines エントリ数=8 procCoastLine central-america mapCoastLines エントリ数=6 procCoastLine south-america mapCoastLines エントリ数=4 procCoastLine australia-oceania mapCoastLines エントリ数=2 3 (-1800000000 675522267) (-1800000000 669296325) ポリゴン数=729300 実行時間=1.53分
結果には2つ問題がある。
一部の描画がおかしくなった。ポリゴンの外側が塗りつぶされているようだ。
第二に、ポリゴンを適当に分割しないと大陸ポリゴンが大きすぎる。 世界地図は間引きをした低ズーム用のみとして、日本地図は間引き無しとする案があるが、 パフォーマンス上は分割が望ましい。
なんとか分割方法を考えたい。
時計回りでないものがあった。これを反転するのを忘れていた。この追加修正で描画エラーはなくなった。