RTL-SDRでエアレース機のADS-Bメッセージを受信し、飛行軌跡をJupyterとpandasでGeoJSONに変換してプロットしてみる
機材のセットアップですが、受信機にはいつものRTL2832UドングルをノートPC(Mac Book Pro)に接続しています。受信用ソフトウェアはdump1090をコマンドラインから使います。dump1090は受信したADS-Bメッセージをデコードした結果をTCPのポート30003で出力する機能があります。これを活用してロギングしてみることにします。 UNIX系のOSでは、nc(NetCat)というコマンドで、任意のホスト、ポートにTCP接続をし、受け取ったデータをそのまま標準出力に出すことができます。これをリダイレクトすることでファイルに記録できます。この出力フォーマットはCSVのテキスト形式ですので、データの処理は難しくありません。
MacOSの場合、dump1090コマンドはソースからビルドする必要があります。準備として、RTL2832U用のライブラリを用意しておく必要があります。ライブラリはbrewでインストールしておくのが簡単です。あとはdump1090のソースコードをgithubから落としてきて、makeするだけでOKです。Linux(Raspberry Piなど)の場合もbrewの代わりにapt-get等を使うことで同様にビルドできます。
$ brew install rtl-sdr $ git clone https://github.com/MalcolmRobb/dump1090 $ cd dump1090 $ make
RTL-SDRドングルをPCに接続し、dump1090コマンドを起動します。アンテナは付属のものをつかいます。受信に成功したらログメッセージが流れていきます。次にやるTCPでデータを取れるように--netオプションを付けて起動しておきます。
$ ./dump1090 --net
この状態で別のターミナルを開いて、ローカルホスト(127.0.0.1)、ポート30003からデータを読むため、下記のコマンドを入力します。
$ nc 127.0.0.1 30003
そうすると、テキストでCSV形式のデータが表示されるはずです。デコードできたメッセージが、1メッセージ=1行で出力されます。あとは、これをリダイレクトでファイルに落とせば記録できることになります。
$ nc 127.0.0.1 30003 > adsb.csv
ところでdump1090には、いくかバリエーションがあります。検索で最初に引っかかる https://github.com/antirez/dump1090 で公開されているバージョンでは、30003ポートから下記のような出力が得られます。
MSG,8,,,84BB74,,,,,,,,,,,,,,,,, MSG,8,,,84B79C,,,,,,,,,,,,,,,,, MSG,3,,,D9E098,,,,,,,40325,,,,,,,0,0,0,0 MSG,5,,,84B79C,,,,,,,2400,,,,,,,,,, MSG,8,,,84BB74,,,,,,,,,,,,,,,,, MSG,5,,,84B79C,,,,,,,2450,,,,,,,,,, MSG,5,,,84B79C,,,,,,,2475,,,,,,,,,, MSG,5,,,84B79C,,,,,,,2500,,,,,,,,,, MSG,6,,,84B79C,,,,,,,,,,,,,7301,0,0,0,0 MSG,5,,,84B79C,,,,,,,2525,,,,,,,,,, MSG,8,,,84B79C,,,,,,,,,,,,,,,,,
カンマの連続が目立ちます。空のフィールド多いようです。
このバージョンの他に、https://github.com/MalcolmRobb/dump1090 にも公開されています。こちらは種々の改良が加えられており、CSV出力の形式形式が違っています。こんな出力が得られます。
MSG,7,111,11111,8675BC,111111,2017/06/04,10:18:14.261,2017/06/04,10:18:14.252,,7375,,,,,,,,,,0 MSG,4,111,11111,8675BC,111111,2017/06/04,10:18:14.316,2017/06/04,10:18:14.254,,,271,217,,,3968,,,,,0 MSG,7,111,11111,86752C,111111,2017/06/04,10:18:14.356,2017/06/04,10:18:14.318,,11025,,,,,,,,,,0 MSG,7,111,11111,862248,111111,2017/06/04,10:18:14.530,2017/06/04,10:18:14.514,,22525,,,,,,,,,,0 MSG,5,111,11111,86752C,111111,2017/06/04,10:18:14.538,2017/06/04,10:18:14.515,,11025,,,,,,,0,,0,0 MSG,8,111,11111,8675BC,111111,2017/06/04,10:18:15.106,2017/06/04,10:18:15.104,,,,,,,,,,,,0 MSG,8,111,11111,8675BC,111111,2017/06/04,10:18:15.115,2017/06/04,10:18:15.104,,,,,,,,,,,,0 MSG,8,111,11111,8675BC,111111,2017/06/04,10:18:15.124,2017/06/04,10:18:15.105,,,,,,,,,,,,0 MSG,7,111,11111,85D7C0,111111,2017/06/04,10:18:15.195,2017/06/04,10:18:15.170,,2550,,,,,,,,,,0 MSG,8,111,11111,8675BC,111111,2017/06/04,10:18:15.195,2017/06/04,10:18:15.170,,,,,,,,,,,,0 MSG,3,111,11111,8675BC,111111,2017/06/04,10:18:15.226,2017/06/04,10:18:15.172,,7425,,,35.44315,139.86805,,,,,,0
こちらは埋められたフィールドが増えています(一部のフィールドはダミーのようです)。重要なのは日時についての情報が入っていることです。航跡を時系列データとして記録するためには、MalcolmRobb版を使った方が良さそうです。
以上のような準備をしたえうえで、当日現地で、屋外でノートPC(MBP)を広げて、データを記録しました。日差しが強く、ノートPCの画面が見えにくく少し難儀しました。
さて、得られたデータを観察してみると、行によって内容が違っています。データフォーマットはこちらに解説がありました。2 番目のフィールド2=Transmission Typeでメッセージの種別が区別されているようです。フィールド5にAircraftIDがあり、機体(メッセージの発信者)を識別するようです。
緯度経度はTransmission Typeすなわちタイプ3メッセージのカラム15,16に記録されています。また、 高度はタイプ2,3,5,6,7のいずれか、カラム12のようです。これらを取り出せば軌跡が得られるはずです。
さてこの手順でデータを処理して、プロットすることが目標です。何かツールがあるかもしれませんが、テキストファイルですので、自力で処理してみます。時系列データですので、PythonとPandasを使ってみます。最終的にGeoJSONに変換することで、Google Mapをはじめ、いろんなツールで取り扱えるはずです。
まずはpandasをimportし、CSVファイルを読み込みます。
import pandas as pd dat = pd.read_csv('20170604-1252.csv', header=None, index_col=9)
フィールド9にある時刻をインデックスに指定しています。こうすることで取り扱い易く易くなります。
読み込んだデータはJupyter notebookのうえではこんなふうに表示されます。
NaNが多く、意外とスカスカであることがわかります。
メッセージタイプ=3のレコードから、AirCraft IDと、時間、高度、緯度、経度を取り出してみます。
メッセージタイプ5を取り出します。Aircraft ID,時刻,コールサイン,高度(FT)のフィールドだけを表示しています。これでコールサイン=便名とAircraft IDの対応がわかります。
高い高度で飛んでいるのは商用機です。羽田の近くなので飛行機はたくさん飛んでいます。
エアレース機は低空を飛んでいますので、1000FT以下のメッセージだけを取り出してみます。NaNになっているレコードが多いのでそれも取り除きます。
エアレース機らしいコールサインが見つかります。このうち、N31YMというのがFalken室屋機です。Aircraft IDが A346AA ということがわかります。
Aircraft ID(フィールド4)が’A346AA’のレコードを取り出します。
結構たくさん出てきます。930行ありました。これならプロットできそうです。
matplotlibでグラフ化してみます。jupyter notebookの中で%matplotlibディレクティブでinlineを指定します。
見事グラフになりました。上下に大きく変化しているのがわかります。
プロットがブツ切れになっているのは、NaNが混じっているせいです。dropna関数を使って取り除きます。
いい感じになりました。途中怪しげなジャンプがありますが、何かあるのかもしれません。
さて次は緯度と経度です。メッセージタイプ3の、フィールド14,15に入っています。
こちらは186個取れてました。軌跡は描けそうです。
さて、このデータをgeojsonに変化してみます。まずはpandasから緯度経度を普通のデータに直します。イテレータとリスト包括表記で、経度緯度のタプルの配列に変換します。
うまくいったら、今度はgeojsonというモジュールを使ってgeojson形式に変換します。緯度経度の配列からLineStringを生成、それをFeatureとしてgeojsonに追加します。geojsonのデータができたら、こんどはjsonモジュールを使って、jsonファイルを生成します。
ファイルの冒頭を見ると、それらしいデータが入っています。
このgeojsonファイルを地図の上にプロットしてみます。今回は簡単のためgeojson.io というWebサイトを利用します。ブラウザで http://geojson.io/ を開き、geojsonファイルをアップロードするだけで表示させることができます。
稲毛海岸で待機旋回している様子や、会場前で競技の折り重なった飛跡、そして終わった後、浦安方面を戻る様子が記録できていました。
手順は省略しますが、他機や商用機もすべて含めたgeojsonを作ってみました。geojson.ioは、GITHub/Gistに置いたデータを開くことができます。置いておきましたので興味があればご覧ください。こちらのgeojsonには高度情報も含めています。
本来であれば、Jupyterの上で表示させると良かったのですが、ごにょごにょと面倒なので、ここではGeoJsonファイルを作ってWebサイトで表示させてみました。コードの断片をスクリプトファイルにまとめることも可能です。是非お試しください。
ちなみにGeoPandasという地理情報用のPandasパッケージもあります。geojsonファイルを直接扱うことができ、簡単なプロットも可能です。shapelyというモジュールで、いきなりセルの上でプロットされてしまいます。こちらも便利そうです。
補足
Jupyterを使える環境にある場合は問題ないと思いますが、pandasやgeojsonはモジュールとして追加が必要です。pipでインストールすることが可能です。
pip install pandas pip install geojson
参考
使用したdump1090(ADS-Bデコーダ) https://github.com/MalcolmRobb/dump1090
ADS-Bメッセージの解説 http://woodair.net/sbs/Article/Barebones42_Socket_Data.htm
Jupyter http://jupyter.org/
Pandas http://pandas.pydata.org/
GeoPandas http://geopandas.org/
geojson.io http://geojson.io/
受信したデータをプロットしたものhttp://geojson.io/#id=gist:/91a831b32253595ff02777cee42cbb1e
