Overpass API で得た海岸線データ(wayオブジェクト) と Geofabrik からダウンロードした europe-latest.osm.pbf などから陸地ポリゴンを作り出すことを試みたが、データ取得日時のずれから、 wayオブジェクトで参照している node の一部が *.osm.pbf には存在しない。
Overpass APIでは、世界規模の海岸線データを node情報を含めて一度に取り出すことはできない。 分割した場合には、時間的なずれがあるため、繋ぎ合わせで問題が出てくる可能性が高い。
Geofabrik の場合、全ての *.osm.pbf がある時刻の planet.osm から切り出されたものと期待する。 もし、europe-latest.osm.pbf と asia-latest.osm.pbf の元となる OSMデータの時刻が異なるならば europe と asia にまたがる海岸線データがぴったり繋がらないケースが出てくるだろう。
たとえ時刻に食い違いはなくても、海岸線データにエラーがあった場合は更新しないといった 対策は必要であるが、まずは、 *.osm.pbf から海岸線wayデータを抽出してみよう。 [2023.1.27]
OpenStreetMapの陸地のレンダリングは一般にはOSMデータから作り出されたシェイプファイル[3]を使う。
自作OSM地図では、現在は Mapnik を使わないため、シェイプファイルのままでは使えない。 これまでは、このシェイプファイルをPostgreSQL/PostGISデータベースに一旦インポートして、 OpenGIS の WKT(Well Known Text)形式でエクスポートしたものを使ってきた。
手間はかかるが、プログラム的には楽である。
今回は、より深く OpenStreetMap を理解するために、 OSMデータから直接陸地ポリゴンデータを取り出すことを試みた。
世界全体の海岸線データを Overpass APIで得ることができる[5]。
wayオブジェクト natural=coastline については以下のようにして、Overpass API で入手できる。 ２GB近いが、９分弱で得られた。
日本全体の osmファイルは 40GB近いので、それに比べれば、5％程度に過ぎないので、 パソコンで解析できる。
c:\gisa>curl --globoff -o c:/osm/world_coastline.xml -d @c:/osm/overpass_coastline.xml http://overpass-api.de/api/interpreter
% Total % Received % Xferd Average Speed Time Time Time Current
Dload Upload Total Spent Left Speed
100 1953M 0 1953M 0 176 3764k 0 --:--:-- 0:08:51 --:--:-- 4379k
<osm-script timeout="900" element-limit="2150000000"><!-- 60 minutes、約 2GB -->
<query type="way">
<has-kv k="natural" v="coastline"/>
</query>
<print/>
</osm-script>
このクエリで得られるのは次のような wayオブジェクトのみである。 狭い範囲であれば、参照されている node オブジェクトも得ることができるが、 世界全体に対しては不可能である。
<way id="131026052" [後略] > <nd ref="1442787079"/> <nd ref="1442787122"/> <nd ref="1442787154"/> <nd ref="1442787072"/> <nd ref="1442787061"/> <nd ref="1442787096"/> <nd ref="1442787105"/> <nd ref="1442787157"/> <nd ref="1442787160"/> <nd ref="1442787150"/> <nd ref="1442787034"/> <tag k="natural" v="coastline"/>
<nd ref="1442787150"> から数値 1442787150 を取り出すのは簡単であり、SAX を使うまでもない。
static void loadCoastlineNodes(String path) {
num_nodes = 0;
try (BufferedReader br = new BufferedReader(
new InputStreamReader(new FileInputStream(path))) ) {
String line;
while ((line = br.readLine()) != null) {
int ix = line.indexOf('<');
if (ix < 0) continue;
String s = line.substring(ix); // 行先頭のスペースを削除する
if (s.startsWith("
取り出された node ID は約 71M であり、時間は１分とかからない。
Overpass APIといえども、一挙に 71M ノードの情報は得られない。 分割すれば Overpass API で入手できるが、Overpass API の場合、 使用制限があるため、あれこれ試行錯誤できない。
今回は、Geofabrikサイトから latest-europe.osm.pbf、latest-asia.osm.pbf など合わせて８つの pbf ファイルから ノード情報を取り出した。
pbfファイルでも合計約 50GB となり、解凍すると１TB 近くにもなるため、 下のようにパイプ接続で node の座標値を取り出す。
osmosis でも解凍できるが、時間がかかりすぎる。
osmconvert d:/downloads/europe-latest.osm.pbf | java LandPolygon europe
プログラム自体は比較的簡単である。LandPolygon は wayセクションのデータしか使用せず、 先に終了するため、osmconvert は最後は write error となる。 LandPolygonの処理が終わった時点で強制終了させてよい。
/*
javac LandPolygon.java
osmosis --rbf d:\downloads\europe-latest.osm.pbf --wx file=- | java LandPolygon europe
osmconvert d:/downloads/europe-latest.osm.pbf | java LandPolygon europe
*/
import java.io.*;
import java.util.*;
import java.time.*;
public class LandPolygon {
static int E7 = 10000000;
static long[] nodes = new long[80000000];
static int num_nodes;
static int ixNode;
static void loadCoastlineNodes(String path) {
num_nodes = 0;
try (BufferedReader br = new BufferedReader(
new InputStreamReader(new FileInputStream(path))) ) {
String line;
while ((line = br.readLine()) != null) {
int ix = line.indexOf('<');
if (ix < 0) continue;
String s = line.substring(ix); // 行先頭のスペースを削除する
if (s.startsWith("<nd ")) { // <nd ref="7943616013"/>
int bgn = 9;
int end = s.indexOf("\"", bgn);
String ref = s.substring(bgn, end);
long idRef = Long.parseLong(ref);
nodes[num_nodes++] = idRef;
}
}
} catch (Exception e) {
e.printStackTrace();
}
}
static void parseXML(String name) {
char quote = '"';
int cntFound = 0;
String pathOut = "c:/osm/" + name + "_idlonlat.dat";
long idCurr = nodes[0]; // 最初に探すノード
ixNode = 0;
// 標準入力
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);
while (idNode > idCurr) {
if (ixNode >= num_nodes) return;
idCurr = nodes[ixNode++];
}
if (idNode < idCurr) {
continue;
} else if (idNode == idCurr) { // 海岸線ノード
cntFound++;
idCurr = nodes[ixNode++]; // 次に探すノード
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);
//System.out.printf("%d: %d %d\n", idNode, lon, lat);
}
}
}
} catch (Exception e) {
e.printStackTrace();
}
}
public static void main(String[] args) {
long start = System.currentTimeMillis();
Instant instant = Instant.ofEpochSecond(start/1000);
ZonedDateTime zonedDateTime = ZonedDateTime.ofInstant(instant, ZoneId.systemDefault());
System.out.println(zonedDateTime);
loadCoastlineNodes("c:/osm/world_coastline00.xml");
System.out.println("海岸線ノード数＝" + num_nodes);
Arrays.sort(nodes, 0, num_nodes-1);
System.out.println("並び替え(昇順)終了");
parseXML(args[0]);
System.out.printf("実行時間＝%.1f分\n",
(System.currentTimeMillis() - start)/60000.0 );
}
}
出力データは node id（８バイト）、lon（４バイト）、lat（４バイト）のバイナリレコードファイルである。
CPU負荷やディスク(SSD)負荷は低いので、同時に２，３ジョブを実行することができる。 全体では 実行時間は 2.5～3 時間となった。
高速化の余地はあるが、陸地ポリゴンデータの更新は精々年に１、２回に過ぎないので、 プログラムを複雑にしない方がよいと思われる。
繋ぎ合わせは、node id のままで実行できる。出力するときに、id を座標値に変える。 最終的にはバイナリレコード形式とするが、まずは、PostGISのエクスポートと同じ WKT形式とする。
元々OSMデータ構造はシンプルであるが、今回は way データのみであるから、SAXも使わず、 自前でパースする。