tomopinのブログ

ジャンク好き。おもちゃ好き。二人の孫はもっと好き。

QGisを用いて鉄道データの座標値を抽出

前回まで、国土数値情報の鉄道データ(路線、駅)を少しいじってみました。今回はこれらのデータから、普段は見えない座標値を抽出してみたいと思います。何のためかというと、駅の座標値がわかるとGPSへ送ったり、路線の連続直線の交点の座標間距離の加算から駅間距離を算出したりすることができます。あまり意味がないかもしれませんが、こんなこともできるというGisソフトのひとつの機能紹介です。

今回の話題では、全国の鉄道データを使うとデータ量が膨大になるので、前回及び前々回で使った「わたらせ渓谷鉄道」の駅および路線のデータを使用したいと思います。

路線別のデータおよびQGisそのものがない場合には、前々回の「Gisを使って特定の路線表示」をご参照いただき、路線別の路線および駅データを作成して下さい。

まずは、QGisを新規立ち上げ、もしくは前々回作成した「わたらせ渓谷鉄道」のプロジェクトファイルを保存してあればそれを起動させてもOKです。

新規立ち上げの場合には基盤地図はどれでもよいのですが、「地理院地図」などを表示させておいてください。

次に路線別に分離した路線データおよび駅データのフォルダから、それぞれ「わたらせ渓谷鉄道」の路線データと駅データを読み込みます。

わたらせ渓谷鉄道
見やすいように路線データの線を少し太めの赤線にしています。

本日のご説明はここからです。

まずは簡単な駅データから座標値を作っていきます。

 

1. 駅データの座標値の抽出

 わたらせ渓谷鉄道の駅データ「N05_002_わたらせ渓谷線」のアイコンが「◯」になっている方の名前を右クリックして、「属性テーブルを開く(A)」をクリックします。

「N05_002_わたらせ渓谷線」の駅データの「属性テーブルを開く(A)」
データの内容を表示します。

駅データの属性テーブルは初期状態では以下のようになっています。

「駅データ」の属性テーブル

駅データの属性テーブルの項目列は1から11まであり、11列目に「駅名」が入っています。路線データは10までです。
ここの右端に「X」と「Y」の列を追加し、追加した列に隠れている座標値を表示させます。

このデータ表示画面で、この例ではアイコンの並んだところの右から4番め、そろばんを縦にしたアイコンで「フィールド計算機を開く(Ctrl+l)」をクリックします。

「フィールド計算機を開く(Ctrl+l)」

すると、「フィールド計算機」の設定画面が開きます。

「フィールド計算機」設定画面
左が式の入力画面、真ん中が関数テーブル、右は備考欄です。
式やデータ値などはほとんど、関数テーブルの中から選択します。

まずは、列の追加から、最初はX列を追加します。

X列の追加、定義

上の赤枠内をご参照ください。

まずは「▢ 新規フィールドを作成」の▢にチェック、

次の次の行、

フィールド名;「X」と入力(Eでも構いません)、

フィールド型;「小数点付き数値(real)」を選択、

フィールド長;「10」を選択、

精度;「6」を選択

ここまで設定が終ったら、次は真ん中の関数テーブルの中から、「▼ ジオメトリ」の▼をクリックするとたくさん関数が出てきますので、

「フィールド計算機」「ジオメトリ」一覧を開く

アルファベット順に並んでいるので、ずーっと下の方、「$x」を見つけたらその行をダブルクリックしますと、式の欄に「$x」と入力されます。

「フィールド計算機」「ジオメトリ」「$x」を選択
ダブルクリックすると左の空欄に「$x」と入力されます。

式が入力されている状態で、式の欄の下、「プレビュー:」のところに数値が表示されていればOKです。数式が成立していない場合にはここに「不正」とかなんか、忌まわしい単語が並びますので、余計な記号や数値がないか、式を見直します。

「プレビュー:」右が数字であれば、ここで右下の「OK」をクリックします。すると先程の属性テーブルに戻ります。

属性テーブルに「X」列が追加されました。
駅の経度が度単位で表示されています。
上で小数点以下を「6」にしていないと値に誤差がでます。

同じようにして、Y列を追加します。再び、右上の縦そろばんアイコン「フィールド計算機」をクリックします。

出てきたフィールド計算機設定画面で、

「▢ 新規フィールドを作成」の▢にチェック、

次の次の行、

フィールド名;「Y」と入力(Nでも構いません)、

フィールド型;「小数点付き数値(real)」を選択、

フィールド長;「10」を選択、

精度;「6」を選択

ここまで設定が終ったら、次は真ん中の関数テーブルの中から、「▼ ジオメトリ」の▼をクリックするとたくさん関数が出てきますので、「$y」をダブルクリックします。

下の画面のように、式の欄に「$y 」が表示されます。同時に一番下の「プレビュー:」に数値が表示されていれば「OK」です。数値でなく不吉な文字が現れていれば、数式欄に余分な記号等が混じってないか確認します。

Y列を追加

「OK」をクリックすると属性テーブルに戻ります。

「属性テーブル」の右端に「Y」列が追加
駅の緯度が度単位で表示されました。
最後の3列だけでも使えそうです。

ここで忘れずに、修正した属性テーブルを保存します。

修正した「属性テーブル」の保存
左上の「フロッピー」アイコンをクリック、
修正がない場合には押せないように薄くなっています。

これで「駅」データのshp形式ファイルまで修正保存されます。

テーブル右上端の「✕」をクリックして閉じます。

これで、駅データをcsv形式(shift-JIS)でエクスポートするとエクセルやテキストエディタで編集できます。エクスポートの手順は次の項目でご説明します。

 

2. 路線データの座標値の抽出

で、次は路線データです。路線データは細かい直線の連続なので、その交点が座標値として保存されています。

ただ、駅データのようにすぐには編集できません。

まずは路線データの属性テーブルを開いてみます。

路線データの属性テーブルを開く

開いてみると、

路線データの属性テーブル
1路線で1行のデータです。

1路線では1行しかデータがありません。1本の連続直線ということです。

まずは、この連続直線を構成する各直線の交点を抽出します。

一旦、属性テーブルを閉じます。

路線データのデータ名「N05_002_わたらせ渓谷線」を選択した状態で、上のメニューバーから「ベクタ(O)」→「ジオメトリツール(E)」→「頂点を抽出...」とクリックしていきます。

既出ですがメニューバーに「ベクタ」がない場合には、メニューバーの「プラグイン」から「Processing」プラグインをインストールしてください。「未インストールのプラグイン」一覧から探してインストールするだけです。

「ベクタ(O)」→「ジオメトリツール(E)」→「頂点を抽出...」

すると、頂点を抽出する設定画面が開きます。

「頂点を抽出」設定画面

「入力レイヤ」には路線データ名「N05_002_わたらせ渓谷線」が選択された状態で出てきます。

「出力」欄は「[一時レイヤを作成]」のままです。

座標が必要ない場合等はここで、CSV形式などで保存したりもできますが、ここではそのまま保存してもほとんど意味がないので、とりあえず「一時レイヤ」です。

駅データ同様に「X」列「Y」列を追加してから保存します。

「▢ アルゴリズムの終了後、出力ファイルを開く」の▢にチェックを入れます。

で、右下の「OK」をクリックすると、路線データの頂点が抽出されます。

「頂点の座標抽出」作業終了

作業が完了すると下の経過を示す棒グラフのところに「完了」と表示され、

後の画面のレイヤ欄に「出力レイヤ」として、点のデータが新たに作成されています。

この時点での出力レイヤ名の右端に「ICチップ」みたいなマークが出ていますが、このマークはまだどこにも保存されていないデータという意味で、この状態で作業を終了すると、このデータは消えてしまいます。いわゆる「保存してケロ」要求マークです。

ここではこのデータを後で保存するので、そのまま続行です。

作業画面を終了させて、「出力レイヤ」の属性テーブルを表示させてみます。

「出力レイヤ」選択

「出力レイヤ」右クリック→「属性テーブルを開く(A)」

そうすると、始点から始まって、直線の交点が表示されています。この路線の場合には1,037個の交点からなる、すなわち、1,036本の直線から構成されています。

「出力レイヤ」の「属性テーブル」
右の方に「distance」の欄がありますが、これは
角度単位なので、実際の距離ではありません。
地球の楕円体により、角度は方位によって
長さが異なりますので、
簡単には距離単位には変換できません。

ここから、駅のデータと全く同じ手順で、座標値のX列とY列を追加していきます。

ツールアイコンの縦ソロバンアイコン「フィールド計算機を開く(Ctrl+l)」をクリックします。

「属性テーブル」の「フィールド計算機」

 

「属性テーブル」の「フィールド計算機」「X」列を追加

 

駅データのときと全く同様に、「フィールド計算機」「X」列を追加します。

まずは「▢ 新規フィールドを作成」の▢にチェック、

次の次の行、

フィールド名;「X」と入力(Eでも構いません)、

フィールド型;「小数点付き数値(real)」を選択、

フィールド長;「10」を選択、

精度;「6」を選択

ここまで設定が終ったら、次は真ん中の関数テーブルの中から、「▼ ジオメトリ」の▼をクリックするとたくさん関数が出てきますので、「$x」をダブルクリックすると、式の欄に「$x」と入力され、そのずーっと下の「プレビュー:」欄右に数値が表示されていればOKです。

「属性テーブル」の「フィールド計算機」→「X」列の追加

ここで右下のOKをクリックすると属性テーブルに戻ります。

「出力レイヤ」の「属性テーブル」右端に「X」列が追加された

同様に「Y」列もフィールド計算機から追加していきますが、手順は省略して、

「出力レイヤ」の「属性テーブル」の右端に「X」列と「Y」列が追加された

これで路線データの連続直線の交点が抽出されましたので、左上の「フロッピー」マークを押してから「属性テーブル」を閉じます。

このデータは「出力レイヤ」という一時データなので、このまま終了するとデータは消えてしまいます。なので、ここでデータを「エクスポート」します。

「出力レイヤ」を右クリックして、「エクスポート(x)」→「新規ファイルに地物を保存(A)」をクリックしていきます。

「出力レイヤ」のエクスポート

すると「エクスポート」の設定画面が開きます。

「エクスポート」の設定画面

ここで、

形式;「カンマで区切られた値[csv]」(なんと一番下!)を選択(右端の下向き▼)、

「カンマ区切られた値[CSV]」を選択

CRS;「EPSG:6668 - JGD2011」を選択、

文字コード;「Shift-JIS」を選択、日本でエクセルやテキストエディタで開くならこれです。

を設定してから、再び上に戻って、

「ファイル名」右の「...」をクリックするとファイル保存先を指定するためのブラウザが開きます。

保存先指定のブラウザ

ここでは保存先フォルダを決めて、路線座標がわかるような名前をつけておきます。

保存形式はカンマ区切り形式です。

通常のエクスポートはほとんど「shp」形式で保存しますが、ここでは画面に何も表示されないデータでかつ、出力後の数値の取り扱いができないので、「csv」形式です。

「エクスポート」設定終了画面
書き忘れました、一番下の
「保存されたファイルを地図に追加する」のチェックは外します。

赤枠の箇所を全て設定し終わったら「OK」です。

保存先フォルダを見てみます。

出力ファイルの確認
拡張子「csv」が出力データです。
拡張子「qmd」の方は削除して構いません。

2つファイルが作成されていますが、重要なのは拡張子が「csv」のデータです。拡張子「qmd」のデータはエクスポートの度に出力されますが、使ったことは一度もありません。削除してOKです。ログファイルみたいなものでしょう。

この段階で、忘れずに駅座標の方も同様に「csv」形式でエクスポートしておきましょう。あると便利です。

で、早速出力された路線データの方の「csv」形式ファイルをエクセルで開いてみます。

もちろん、テキストエディタで開いても構いませんが、グラフ化できないのでエクセルが吉です。

エクセルで「路線データに座標データを追加したデータ」を開いてみる
先程の属性テーブルの最終形の内容そのままです。

内容は先程の最終段階の「属性テーブル」と同じです。右端に「X」列と「Y」列として、座標値が追加されています。

この「X」列と「Y」列で試しに散布図(線あり)を描いてみます。

エクセルで「線あり散布図」を選択

シートの上にそのまま描画してみます。

エクセルで「わたらせ渓谷鉄道路線図」
最初の画面と比較して特に曲線に矛盾はありません。

ちゃんと「わたらせ渓谷線」の路線の形になっています。

グラフは試しにかかせてみただけなので、このまま保存せずに終って構いません。

ここから先は目的に応じていろいろと操作してもらうことになります。

駅間距離を出すためには、まずは、度単位の緯度経度データの各直線部の端点2点間距離を出して、これを先に算出している駅座標を参考に、各駅間で累積していくことで駅間距離が出せます。道路のデータも同様の手順で座標抽出が可能と思われますが、莫大なデータ量となるので、そのまま作業するのはおすすめしません。国道番号とか、なんらかの範疇で細切れにしてから作業しましょう。QGisは膨大なデータもとにかく読み込んで計算しようとしてくれます。なかなかエラーを出しません。

駅データはもとの[shp]形式ファイル毎修正されますが、路線データの方は出力レイヤをいじってるので、もとの[shp]形式には影響を及ぼしません。及ぼしては描画のためにかえって不都合が生じますのでこれでOKです。

 

3. 座標データの活用の一例

区間距離の算出について

路線を構成する連続直線の各直線間の距離を累積加算することで、特定の区間の距離を算出することができます。ただし、緯度経度データは地球楕円体上の度単位のデータなので、そのままでは三平方の定理を用いた距離計算はできません。

度単位の2点間距離の計算法はいろいろとあるようですが、ヒュベニの公式と呼ばれる計算方法がギリ、エクセルでの計算が容易かつ精度が高い方法なので、簡単にご紹介します。

・2点間(度単位の緯度経度形式「(X1,Y1)-(X2,Y2)」)の距離

テキストで記載すると何故かスマホでみると化けて出るので、画像ファイルにして添付します。かける意味の「*」が脚注になってしまうようです。なんと。

地球の長半径(赤道半径)と短半径(極半径)は、現在GPSなどで使われている1984年のWGS84という測地系が現時点では最も実測値との差が少なく(場所にもよるらしい)、上のように算出されています。

式に入力ミスがあるかもしれません。検索されて再確認を推奨します。

駅間距離を算出してみたい方は是非とも上の式などをエクセルに突っ込んでやってみてください。

 

一応、QGisを用いたヘビーでソフトな話題はこの辺で終わりにしたいと思います。

最後まで読んでいただいた方、誠にありがとうございます。お役に立てると幸いなのですが。

ライトな話題は今後もちょくちょく出るかもしれません。その場合でもくどくない、ライトな話題にしたいと思います。