MADOCALIBによる高精度GPS測位

category: gnss

はじめに

米国のGPSをはじめとした測位衛星からの電波を受信することで、現在位置を知ることができます。このような測位衛星システム(GNSS: global navigation satellite system)には、GPSのほかに、ロシアのGLONASS、欧州のGalileo、日本の準天頂衛星みちびき(QZSS: quasi-zenith satellite system)、中国のBeiDou、インドのNavIC(navigation with Indian constellation)があります。

この測位電波には、衛星から受信機までの擬似距離を測定するための測距信号と、衛星位置や時刻や電波経路遅延量を伝送するための航法メッセージ信号とが含まれます。受信機は、この衛星からのこの2つの信号を参考に、自らの座標を推定します。

しかしながら、航法メッセージ中にある衛星座標、衛星時刻、電波経路上の電離層遅延量が、モデル化の過程などにより、現実と異なる可能性があります。このような現実との相違は、測位誤差の発生につながります。地上での衛星電波観測網によりこれらの補正量を推定して受信機に伝達することにより、このような誤差を軽減して測位精度を高めることは補強(オーギュメンテーション)と呼ばれ、このようにして実現される高精度測位はPPP(精密単独測位、precise point positioning)と呼ばれます。

日本の準天頂衛星みちびきは、3種類の補強信号を放送しています。この一つがMADOCA-PPP(マドカ・ピーピーピー、multi-GNSS advanced orbit and clock augmentation - precise point positioning)です。そして、MADOCA-PPPを利用するための参考プログラムがMADOCALIB(マドカリブ、MADOCA-PPP test library)です。

2024年7月1日に、このMADOCALIBバージョンアップ版が公開され、将来のみちびき衛星から放送される補強信号を処理できるようになりました。この新しいMADOCALIBを用いて、将来の高精度測位を体験してみます。

Start of the Internet distribution of MADOCA-PPP and release of MADOCALIB v1.2

測位衛星から放送される補強情報には、みちびきのCLAS(シーラス、centimeter level augmentation service)とMADOCA-PPPのほかに、GalileoのHAS(ハス、high accuracy service)、BeiDouのPPP-B2bがあります。CLASやMADOCA-PPPでは、補強情報に搬送波位相バイアスを含むため、AR(amibiguity resolution)という方法で、PPP解よりも高精度なFix解を得られる見込みがあります。

MADOCALIBのダウンロードと準備

ここでは、Windows PCを用いて、MADOCALIBにあるサンプルデータを処理します。

内閣府は、2024年6月5日より、MADOCALIBをオープンソースソフトウェアとしてGitHub上に公開しています。2024年7月1日に、MADOCALIBはバージョン1.2に更新されました。下記リンクをクリックして、CodeDownload ZIPの順にクリックすることで(またはgitクライアントソフトウェアを用いることで)、MADOCALIBをダウンロードできます。

MADOCALIB: MADOCA-PPP Test Library

たとえばデスクトップにmadocalibというフォルダを作成し、必要ファイルをここに集めます。ダウンロードしたMADOCALIBから、

  • binフォルダにあるrnx2rtkp.exe
  • app\consapp\rnx2rtkp\gcc_mingwフォルダにあるsample.conf (のちに内容を編集します)
  • sample_dataフォルダにある2024162all.200.l6など全部(batのつくファイルはなくてもかまいません)

をこのフォルダにコピーします。

preparation for MADOCALIB execution

docフォルダにあるMADOCALIB_Detailed manual_ver001.pdfはMADOCALIBの説明書で、MADOCALIB_manual_ver002.pdfはv.1.2で新しく追加された電離層遅延量推定機能の説明書です。

そこで、後者の説明書に書かれた手順にて、このsample.confを修正します

また、データ変換や表示のために、ここではさらにRTKLIB v.2.4.3b34のWindowsバイナリーファイルにあるrtkplot.exertkconv.execonvbin.exeも利用します。これらも、このフォルダにコピーしておきます。

同梱サンプルの実行

その後、コマンドプロンプトを開き、このフォルダに移動し、測位を行います。rnx2rtkp.exeを用いて、測距データTSK200JPN_S_20241620000_01D_30S_MO.rnxと航法データTSK200JPN_S_20241620000_01D_MN.rnxのみで測位できることを確認します。このような方法を単独測位と呼びます。

rnx2rtkp.exe -p 0 -sys G,R,E,J -f 2 TSK200JPN_S_20241620000_01D_30S_MO.rnx TSK200JPN_S_20241620000_01D_MN.rnx

このファイル名はとても長いのですが、コマンドプロンプトではファイル名の一部分を入力した上でTabキーを押すと、ファイル名の残りが補完されます。異なるファイル名が補完されたら、さらにTabキーを押すことを正しいファイル名が補完されるまでくり返します。

ここで指定したオプションは次のとおりです。sample.confを使わずに単独測位を行うためのオプション指定です。

  • -p 0:単独測位を指定します。
  • -sys G,R,E,J:GPS(GPSのG)、GLONASS(RussiaのR)、Galileo(EuropeのE)、QZSS(JapanのJ)を用いて測位します。sample.confで指定したファイルと同じ衛星システムを利用するために、このオプションを指定します。
  • -f 2:周波数帯としてL1とL2を使います。このオプションはデフォルト動作なので、指定しなくてもかまいません。

この測距データは、2024年6月10日(1月1日からの日数が162日)に30秒ごとに衛星電波を観測したもので、24時間分のデータです。このファイル名から、国土地理院の茨城県つくば市のGNSS観測点(座標36.1056N、140.0871E)のデータであると推定できます。Miraiにアカウントを作ると、他地点や他日時のリアルタイムデータやアーカイブデータもダウンロードできます。Miraiのアカウントは無料で作成できます

Mirai system screenshot

これらのファイルはテキスト形式であるRINEX(ライネックス、receiver independent exchange format)形式で書かれているので、テキストエディタで内容を見ることができます。

たくさんの文字と共に、36.105576920や140.087107164などの座標らしきものが表示されたら、正しく動作していると言えます。-o single.posオプションを追加して、この結果をファイルsingle.posに保存します。

次に、衛星軌道と衛星クロックの補強情報2024162all.204.l6と、v.1.2から新しく利用できるようになった電離層遅延の補強情報2024162all.200.l6を追加して、測位します。

rnx2rtkp.exe -k sample.conf TSK200JPN_S_20241620000_01D_30S_MO.rnx TSK200JPN_S_20241620000_01D_MN.rnx 2024162all.204.l6 -mdciono 2024162all.200.l6 -o mdcar200.pos

ここで指定したオプションの意味は次のとおりです。

  • -k sample.conf:既定値の代わりに、sample.confに記述されたオプションを読み込みます。
  • -mdciono 2024162all.200.l6:電離層(ionosphere)遅延の補強情報として2024162all.200.l6を利用します。
  • -o mdcar200.pos:結果を、画面の代わりに、ファイルmdcar200.posに出力します。

このサンプルでは、電離層遅延の補強情報として2024162all.200.l62024162all.201.l6とが用意されています。これらの2つの情報を同時に利用することはできないようですので、はじめに前者を用います。

これらの結果をプロットして比較します。rtkplot.exeを起動し、FileOpen solution-1...メニューにてsingle.posを、Open solution-2...メニューにてmdcar200.posを、指定します。

comparison of single positioning and MADOCA positioning

ここで、罫線は50センチメートル間隔、赤いプロットは単独測位結果、薄青色のプロットは新MADOCA-PPPによる測位結果です。 ここでは、右上の統計情報表示を有効にするために、EditOptions...メニューの中のShow Statisticsをオンにしています。

RTKPLOT option of enabling statistics display

単独測位プロットは24時間の間に大きくばらつくのに対して、MADOCA-PPPの結果プロットは集まっています。統計情報によると、座標は36.105591478N、140.087128046Eとなっています。観測点座標の正解はわかりませんが、MADOCA-PPP補強情報によりばらつきは小さくなりました。

次に時間変化を見るために、RTKPLOTの左上のメニューをGnd Trk(地面軌跡: ground track)からPosition(位置)に変更します。

position display in RTKPLOT

このモード変更により、時間に対する3本のグラフが表示されます。上から、東西方向(easting)、南北方向(northing)、および、上下方向(up-down)で、単位はメートルです。横軸の罫線は2時間単位です。単独測位では位置の時間変化が大きく、MADOCA-PPP測位のそれは初期状態を除いてほとんどありません。

position comparison between single point positioning and MADOCA-PPP positioning

ここで、MADOCA-PPPの結果のみを表示するために、1のボタンをクリックします。測位開始から15分程度は薄い青色(PPP解)ですが、それ以降は濃い青色(Fix解)になっています。Fix解は、ARにより測距補強の矛盾を解決し続けて、補強情報による修正を行うPPPよりもさらに高精度測位が達成されたと自己評価されている状態です。

position of MADOCA-PPP results

電離層遅延の補強データとして2つのファイルがあり、今まではPRN(pseudo random noise)番号200のものを利用してきました。同様にPRN201のものでも測位計算を行い、これら2本の位置変化を同一グラフにプロットします。時間軸のところでマウスホイールを回転させて、時間軸を拡大します。横軸の罫線間隔を5分に設定しました。測位開始から10分間のPPP測位状態で、紫色のプロットがPRN200のもの、緑色のプロットがPRN201のものです。

position of MADOCA-PPP results

測位結果座標の時間変化は、PRN200を用いたときとPRN201を用いたときとで異なりました。一方、10分程度の後、両者はほぼ同一の座標を示しました。

自らの観測データによるMADOCALIB高精度測位

次に、私が職場に設置しているGNSS観測局データと、内閣府のみちびきホームページで公開されているL6信号アーカイブを用いて、MADOCA-PPP測位を行います。

ここでは、2台の受信機のうちのNovAtel(ノバテル)OEM729受信機にて、2023-07-03 00:00:00 UTC(UTCは世界協定時coordinated universal timeで、日本標準時よりも9時間遅れています)から1時間の衛星信号を収録した生データ20240703a.gpsを用います。この生データには、航法データと測距データの両方が含まれます。u-blox ZED-F9P受信機の生データでも試しましたが、8衛星程度しか利用できなくて、こちらは断念しました。

また、MADOCA-PPPの同時間帯の衛星軌道・クロックデータと電離層遅延データが必要です。これらは、みちびきホームページのアーカイブからダウンロードします。データのダウンロードには利用条件の同意が必要です。日付検索を行い、必要ファイルをダウンロードします。

Multi-GNSS Advanced Orbit and Clock Augmentation - Precise Point Positioning (MADOCA-PPP) Service

このアーカイブデータのうち、衛星軌道・クロックに関するものは現在のみちびきからも放送されていますが、電離層遅延に関するものはまだみちびきから放送されていません。信号仕様書IS-QZSS-MDC-002の6ページにあるPRN番号と内容対応表からCLASの情報を除き、みちびき衛星番号を加えたものを次の表にまとめます。

QZSS satellitePRN and MADOCA-PPP contents
QZS-2orbit&clock (PRN204)
QZS-4orbit&clock (PRN205)
QZS-1Rorbit&clock (PRN206)
(QZS-5?)orbit&clock (PRN207)
QZS-3orbit&clock (PRN209)
(QZS-6?)orbit&clock (PRN210), ionosphere (PRN200)
(QZS-7?)orbit&clock (PRN211), ionosphere (PRN201)

みちびき5号機、6号機、7号機はまだ打ち上げられていませんので、考えられる衛星番号を括弧書きで記しました。衛星軌道・衛星クロックの補強情報内容は、現在のところ、すべての衛星間で同一ですので、PRN204からPRN211までのいずれか一つのものをダウンロードします。ここではPRN209のものである2024185A.209.l6をダウンロードしました。ここで、2024は2024年を、185は1月1日から185日目の2023-07-03を、また、Aは1時間分データの開始時刻00:00:00 UTCからのものであることを、209はPRN209を表しています。

一方、電離層遅延の補強情報については、PRNにより異なる内容が配信されるようです。ここでは、この両方の2024185A.200.l62024185A.201.l6をダウンロードします。電離層遅延情報を日付検索して、電離層遅延データの提供開始が2024年6月28日 00:00 UTCであることがわかりました。

ダウンロードしたこれらのファイルを先ほどのフォルダmadocalibにまとめておきます。

次に、rtkconv.exeを用いて、20240703a.gpsから航法メッセージファイルと測距ファイルを抽出します。Convertボタンを押して作成される20240703a.navが航法データファイル、20240703a.obsが測距データファイルです。

また、convbin.exeを用い、コマンドプロンプト上でこの処理を行うことも可能です。ここでは、-yオプションにて、GLONASS、SBAS、BeiDou、NavICのデータを除外しました。

convbin.exe -y R,S,C,I 20240703a.gps

この2つのファイルを用い、単独測位します。マスク角は既定値の10度、周波数帯も規定値のL1帯とL2帯です。衛星システムとして、GPS、Galileo、みちびきの3システムに限定しました。

rnx2rtkp -p 0 -t -sys G,E,J 20240703a.nav 20240703a.obs -o 20240703a-single.pos

この出力ファイル20240703a-single.posをRTKPLOTにてプロットしました。

single point positioning with own data

1時間の間に複数の座標の不連続があります。よい観測データとはいえませんでした。次に、座標の時間変化を表示します。

single point positioning with own data

座標の不連続がありますが、このデータに対してMADOCA-PPPの補強情報を利用します。はじめに、PRN209の衛星軌道・クロックの補強情報のみを用いて測位計算し、プロットします。

rnx2rtkp.exe -k sample.conf 20240703a.nav 20240703a.obs 2024185A.209.l6 -o 20240703a-mdc.pos

MADOCA-PPP positioning with only orbit and clock augmentation

ここで、罫線の間隔は20センチメートルです。測位開始から15分程度は大きく位置変化しますが、その後は一定座標に収束します。次に、PRN200の電離層遅延情報を追加して、測位計算しました。ここで、-x 2はデバッグ用のトレースオプションで、20240703a-mdc200.pos.traceというデバッグファイルが作成されます。

rnx2rtkp.exe -k sample.conf 20240703a.nav 20240703a.obs 2024185A.209.l6 -mdciono 2024185A.200.l6 -o 20240703a-mdc200.pos -x 2

MADOCA-PPP positioning with PRN200 ionosphere corrections

この結果から、電離層遅延補強情報の有無によっても、プロットにほとんどの変化がありませんでした。上述の説明書MADOCALIB_manual_ver002.pdfにデバッグファイルの読み方が書かれています。このデバッグファイル中のmiono_sel_areaを含む行を検索すると、2 miono_sel_area: lat=34.44,lon=132.41と書かれた行しかありませんでした。これは電離層遅延情報が利用されなかったことを表します。

電離層遅延情報として、PRN201のものを同様に試してみます。PRN201のデバッグファイル中には、2 miono_sel_area: lat=34.44,lon=132.41のほかに2 miono_sel_area: closest RegionID=5,AreaNo=5,dist=267.559と書かれた行もありました。電離層遅延情報が利用されたようです。

そこで、PRN200のものと、PRN201のもののGnd Trkを比較します。PRN201の測位結果の方が、より座標変化が小さかったです。

Ground track comparison between PRN200 and PRN201 results of MADOCA-PPP

次に、両者の座標の時間変化の比較します。両者ともに15分程度でほぼ同一座標に収束しますが、この場合はPRN201を用いた方がより座標変化が小さい結果になりました。

Position comparison between PRN200 and PRN201 results of MADOCA-PPP

MADOCA-PPPに新しく追加された電離層遅延情報は、最初から座標の時間変化の少ない測位結果を与える効果があること、与える電離層補強情報により結果が異なる可能性のあることがわかりました。

実は、この実験の過程で、職場のGNSS観測所のデータ収集サーバの時刻同期(NTP)の設定が不完全で、1分程度の時間遅れがあることに気づきました。すでに修正しました。

MADOCA-PPP補強情報のインターネット配信

内閣府は、MADOCALIB v.1.2の提供と同時に、MADOCA-PPPのインターネット配信を始めました。インターネット配信を受信するためには、申請が必要です。配信はMiraiシステムを用いたNTRIP形式にて行われます。利用が許可されると、NTRIPキャスターアドレス、ID、パスワード、マウントポイント名が通知されます。PRNに対する全9箇所のマウントポイントが用意され、同時利用許可セッション数は9です。すなわち、私たちは全PRNの補強情報を同時に利用できます。利用の有効期限は、毎年の3月末日までで、継続利用には申請が必要です。

インターネット配信では、衛星中継機(トランスポンダ)などによる伝送遅延がないため、より早く情報を得ることができます。下記の左画面はインターネット配信のPRN205の、また、右画面はAllystar受信機にてみちびきL6E信号を受信したPRN205の補強メッセージです。インターネット配信の方が9秒程度、早く情報を得ることができました。

たとえば、左画面7行目にあるQZNMA G06(CNAV)は、右画面では16行目にあります。MADOCA-PPPでは1秒間に1データパートが伝送されますので、これは9秒間のL6E信号の遅延に相当します。

おわりに

MADOCALIBがv.1.2にバージョンアップされ、これまでのARによる高精度測位に加えて、短時間で座標を決定できる電離層遅延補強情報を扱えるようになりました。この電離層遅延情報は、打ち上げ予定のみちびき5号機、6号機、7号機から放送される予定ですが、公式アーカイブから取得することができます。

精密単独測位(PPP)では、座標の初期収束のために15分程度の時間が必要で、この間の座標出力は不安定です。MADOCA-PPPでは、新しく利用できるようになった電離層遅延情報により、この座標変化を軽減できるので、より早く座標を得ることができるようになったといえます。今後、2つの電離層遅延量補強データの内容差を解析してみます。

JAXAによるH3ロケット3号機(H3F3)の打ち上げ成功により、みちびき5号機、6号機、7号機の早期利用開始が現実的になりました。みちびきの将来をとても楽しみにしています。