紙の図面から図形の頂点座標を読み取ってデータ化したい時はありませんか？ここでは、できるだけ簡単に河川横断図をデータ化するためにQGISを使って行う方法を解説します。ArcGISなど他のGISソフトでも同じようにできますよ。
QGISはフリーでオープンソースの地理情報システムです。
図面をQGISに読み込む
PDF図面なら画面キャプチャ、紙図面ならスキャナなどを用いて、あらかじめ図面を画像ファイルにしておきます。
QGISを起動し、画像ファイルをラスタとしてレイヤに登録します。画像ファイルには座標系が設定されていませんが、他にレイヤがなければ表示されるはずです。
表示が確認できれば、プロジェクトの座標系とレイヤの座標系を、適当な投影座標系（例えばEPSG6671）に設定しておきます。
設定する座標系は、平面直角座標系やUTM座標系ならなんでもいいよ。m単位のX-Y座標が欲しいだけなんだ。
トレースするためのシェープファイルを２つ用意する
読み込んだラスタレイヤをポイントシェープファイルを使ってトレースします。
QGISの座標を河川横断図の座標に変換する式を作成する必要があるので、横断図中の座標が明確に判読できる点をトレースするためのreference.shpと、河床の変曲点をトレースするためのxsection.shpの２つのシェープファイルを作成します。
「レイヤ」ー「レイヤを作成」ー「新規シェープファイルレイヤ」でダイアログを開き、以下のように設定してOKしましょう。
図面をトレースする
まず、referenceレイヤを編集しましょう。
横断図のX軸の値は河川中央部からの距離として与えられているので、河川中央をX軸の原点とし、右側を正、左側を負と考えます。Y軸の値は緑色の「計画河床高」と紫色の「計画H.W.L.」で与えられています。
以上から、本稿では下図の赤丸の位置の４点をreferenceレイヤに作成しました。
次に、xsectionレイヤとしては、横断図の変曲点に下図の青丸のようにポイントを作成します。
トレースした点の座標をフィールド計算機で読み取る
トレースした点のQGIS上の座標を、フィールド計算機を使って与えます。
下図では、フィールド計算機でX座標を与える方法を示しています。X座標は$x関数で調べることができます。$x関数はフィールド計算機下部中央の関数リストのジオメトリグループ内にあります。同じようにして、Y座標も$y関数を用いて与えます。
以上の作業を、referenceレイヤ、xsectionレイヤ両方に対して行い、作業結果を保存しましょう。
座標をExcelに読み込む
シェープファイルは複数のファイルで構成されており、フィールドデータは拡張子が.dbfのファイルに保存されています。.dbfファイルはExcelを空白のブックで起動した後、「ファイル」ー「開く」でdBASEファイルを指定することで開くことができますので、別の空白のブックにコピぺしましょう。
下図では、referenceシートとxsectionシートを作成してコピペしてあります。
.dbfファイルをExcelで開いたら、そのまま作業しないでかならず別のExcelファイルにコピーしてね。.dbfファイルを編集して保存してしまうと、シェープファイルが壊れてしまうわよ！
アフィン変換で正しい座標系に投影する
前項で読み取ったX, Y座標はあくまでもQGIS上の座標なので、referenceシートを使ってアフィン変換の係数を求めます。
まず、X, Yの右にtargetX, targetYカラムを作成し、横断図から読み取った座標を入力していきます。下図の下段に、座標値の根拠を記載してありますので、じっくり見てください。
次に、左上に6×2のスペースができるよう行・列を挿入し、下図左側の紫枠のように値を設定します。これらの値は、その右に記載したアフィン変換式の係数です。未知数が６個ですので、referenceは最低3点あれば係数値を求めることができます。今回は４点あるので大丈夫ですね。
アフィン変換を使ったのは、縦横比の違う断面図やスキャンした時に傾いてしまった画像データも扱うためなんだ。断面図が実際の地形と相似であることが確実なら、ヘルマート変換（e=a, d=-b)などもっと係数が少ない、つまりreferenceの点数が少なくても良い変換方法を使うこともできるよ。ヘルマート変換の未知数は４個だから、referenceは最低２点でいいんだ。
次に、projX, projY, residual^2の数式をシートに入力し、residual^2のカラムの合計値をSUM residualとする計算式を入力します。
最後に、ソルバーを使ってSUM residualを最小にするa, b, c, d, e, fの値を求めます。下図右のようにソルバーのパラメータを設定して「解決」ボタンをクリックします。
ソルバーがメニューに見つからない人もいるかもね。ソルバーはアドインなので、見つからなければ自分で読み込まないと使えないのよ。すごく便利な機能なので、是非活用してね！
以上の操作で、下図左側のようにSUM residualが0.001程度にまで小さくなり、アフィン変換の係数が赤枠内のように求まりました。targetX,YとprojX,Yがほぼ同じ値になっていることを確認してください。
今度は、xsectionシートのX, Yに対して同じ係数値を用いてprojX, projYを計算します。projX, projYが、横断図の座標系で見た青丸の点の座標になります。
念のために、セルに入力した数式を下図に示しておきますね。下図ではreferenceシートのtargetX, Yのセルにもいくつか数式が入っていますが、Excelを電卓替わりにして計算したからです。G12やF10は断面図に値がそのまま記載されているので式にする必要なかったですね。
工程は多く見えますが、昔はモノサシで図面の寸法を測り、電卓で投影計算していたことを思えば、圧倒的にラクな方法ではないでしょうか。いろいろ応用する場面があると思いますので、皆さんで工夫してください。Have fun!
コメント