海の判別には、衛星データを使って水面を検出するのが普通だと思いますが、ここではお手軽に標高データから海を判別します。流域界を求める際，海をマスクする必要がありますし、逆に衛星データの海の部分だけ使いたい時もあると思います。こうした際に必要になる技術です．
QGISはフリーでオープンソースの地理情報システムです。
DEMの準備
例として使用するのは国土地理院からダウンロードした与論島10mDEMです。本稿では他の目的で既にfill処理（窪地を埋める）してあるものを用いますが、海水面より低い土地でもなければ生のDEMを用いても問題ありません。
DEMをポリゴンに変換する
まず、GDALの「ラスタのベクタ化 (polygonize)」ツールを使います。
元のDEMは実数ですが、この処理で生成されたポリゴンのDNフィールドは整数で、DEMの数値を四捨五入したものになっています。下図は、DNが小さなものだけ色分け表示したものです。このDEMではDNが0から１の範囲が海と考えてよさそうです。本稿ではこのまま進みますが、もし厳密に海岸線を出したければ、海と陸を分ける標高を調べて、水際が四捨五入の分かれ目になるようにラスタ計算でDEMの値を調整してやる必要があるかもしれません。
ポリゴンから海を取り除く
islandポリゴンの属性テーブルのツールバーから「式による地物選択」ツールを開き、「”DN” < 2」という式を入力して「地物を選択」ボタンをクリックします。
すると下図のように海のポリゴンが選択されるので、編集をオンにして、選択されたポリゴンを削除し、結果を保存します。
すると、2m以上のポリゴンが残ります。（背景の白黒模様は別の処理の結果なので気にしないでください。）
残ったポリゴンを融合(dissolve)する
GDALの「融合(dissolve)」ツールを使います。
QGISのdissolveはArcGISと比べて遅いのが難点ですが、一服するうちには終了します。（下図の背景の黄黒縞模様も、関係ありませんので気にしないでください。）
脱線：やっぱりdissolveが遅い件
やはりdissolveが我慢ならないくらい時間がかかる場合があるようです。実際私も今、もう少し広く複雑な場所で同じ作業をしていましたが、３時間ぐらいdissolveしても終わらないのでArcGISでdissolveをかけたら10分もせずに終わりました。そのような時にあくまでQGISで立ち向かうなら、むしろフィーチャ数が少ない海のポリゴンを選択してエクスポートし、領域全体を包む単一フィーチャの矩形ポリゴンを新たに作って海のポリゴンでClipした方が早いかもしれません。dissolveの速度についてはずいぶん昔から遅いとか遅くないとか議論が続いていることが英語で検索するとわかります。無料ソフトなので工夫して切り抜けましょう！
本稿は以上です。Have fun!
コメント