カシミール3D: アメリカの州境データの軽量化 - 思考と試行の物置部屋(Attic of Thoughts and Experiments)

州境データとして、以前の記事に書いた通り、http://www.npms.rspa.dot.gov/data/data_base.htmを使っていたが、あまりにも詳細なため、カシミール3Dの起動に毎回10秒以上かかっていた。当初これでもよいかと思っていたが、そのうち遅さに耐えられなくなったので、元のデータを間引いて軽量化することにした。

e00データのフォーマットについては、http://avce00.maptools.org/docs/v7_e00_cover.htmlを参照した。

ARCの説明のうち、from_node, to_nodeというのが分かりにくいので補足。e00 データ上でまとまって記録はされていないものの、IDがついたnode(座標値)が多数存在する。それぞれの線分は、これらのnodeを始点・終点として間をIDなしのnodeで結んだ形になっている。この始点・終点のIDを保持しているのが、from_node, to_node という項目である。

カシミール3DのGISツールプラグインは、線分データのみしか扱わないと記載されているので、ARCの部分だけ加工することにした。複数の線分を組み合わせて処理するのは大変そうなので、線分単位で以下の処理を行う。

線分の始点と終点は必ず残す(これにより、from_node, to_nodeの情報を加工しなくて済む)
線分の間に入っている点は、前の点から1マイル以上離れているときに限り残す。

なお、1マイルというのは、普段使っている地形データであるSRTM-30と同程度の粗さになる距離を適当に選んで決めた。

Perlでの実装

Perlで適当に実装する。なお、Cygwin上でperl5.8.7を使用した。

探してみたらGeo::Distanceというモジュールがあったので、インストール。このモジュールがDBIというモジュールを必要とするらしいので、それもインストール。

作成したコードは末尾の通り。自分が使っているe00データを軽量化するためだけに作ったものなので、他のデータに使えるかどうかはほとんど考えていない。また、1度しか実行しないものなので、速度向上の努力もしていない。

これを使って、130MB程度の元データが8MBぐらいに軽量化できた。カシミール3Dの起動時間もまったく気にならないレベルに短縮された。精度もSRTM-30と併用するには十分である(4倍ぐらいまでの拡大ならまったく違いが分からないレベル)。

#!/bin/perl

# e00-simplify: simplify e00 ARC, considering distance between points
# Usage: ./e00-simplify < state.e00 > state_sm.e00

use Geo::Distance;
$geo = new Geo::Distance;
$geo->formula('hsin');

$minimum_dist = 1.0; # in mile

sub distance_between_lines {
    ($line_from, $line_to) = @_;
    $line_from =~ s/^\s*//;
    ($lon1, $lat1) = split(/\s+/, $line_from);
    $line_to =~ s/^\s*//;
    ($lon2, $lat2) = split(/\s+/, $line_to);
    $geo->distance('mile', $lon1, $lat1 => $lon2, $lat2);
}

# process input upto "ARC"
while (<>) {
    last if /ARC/;
    print;
}

print;

# Process ARC part
while (<>) {
    # Search for arc header
    if (/(-*\d+)\s+(\d+)\s+(\d+)\s+(\d+)\s+(\d+)\s+(\d+)\s+(\d+)\s+/) {
        $coverage_num = $1;
        $coverage_id  = $2;
        $from_node    = $3;
        $to_node      = $4;
        $left_polygon = $5;
        $right_poligon= $6;
        $num_coordinates = $7;
        last if ($coverage_num == -1 && $coverage_id == 0);
        # process one (1) arc (set of coordinates preceded by header
        @output_lines = ();
        for ($i = 0; $i < $num_coordinates; $i++) {
            $line[$i] = <>;
        }
        push(@output_lines, $line[0]);
        for ($i = 0; $i < $num_coordinates -1; $i++) {
            $line_from = pop(@output_lines);
            push(@output_lines, $line_from);
            if (distance_between_lines($line_from, $line[$i])
                > $minimum_dist) {
                push(@output_lines, $line[$i]);
            }
        }
        push(@output_lines, $line[$num_coordinates-1]);
        # output one (1) arc
        printf "%8d%8d%8d%8d%8d%8d%8d\n", $coverage_num, $coverage_id,
        $from_node, $to_node, $left_polygon, $right_poligon, $#output_lines+1;
        foreach (@output_lines) {
            print;
        }
    }
}

printf "%8d%8d%8d%8d%8d%8d%8d\n", $coverage_num, $coverage_id,
    $from_node, $to_node, $left_polygon, $right_poligon, $num_coordinates;