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つ問題がある。
一部の描画がおかしくなった。ポリゴンの外側が塗りつぶされているようだ。
第二に、ポリゴンを適当に分割しないと大陸ポリゴンが大きすぎる。 世界地図は間引きをした低ズーム用のみとして、日本地図は間引き無しとする案があるが、 パフォーマンス上は分割が望ましい。
なんとか分割方法を考えたい。
時計回りでないものがあった。これを反転するのを忘れていた。この追加修正で描画エラーはなくなった。