(57)【特許請求の範囲】
【請求項１】
複数の測位衛星から送信される測位用信号をそれぞれに受信し、互いに予め設定された間隔で配置された複数のアンテナと、
複数のアンテナ毎に設置され、入力された測位用信号をそれぞれに復調する複数の受信機と、
各受信機から得られる衛星情報から各受信機と各測位衛星との間のコード擬似距離およびキャリア位相積算値を観測し、各受信機に対する整数値バイアス、各受信機の受信機時計誤差、および各受信機位置を未知数として含み、該各受信機位置を過去の推定演算値で、衛星軌道位置を放送軌道情報による推定軌道位置で１次テイラー級数展開することで線形近似し、前記整数値バイアス、前記受信機時計誤差、前記受信機位置、衛星時計誤差、衛星軌道位置を説明変数に含み、前記キャリア位相積算値およびコード擬似距離を目的変数とし、前記線形近似により得られる近似行列を演算要素とする１つの線形回帰方程式を構成し、該線形回帰方程式から前記整数値バイアスと、前記各受信機の位置と、を推定演算して、該各受信機の位置と各受信機間の間隔とから測位を行う測位演算手段と、
を備えたことを特徴とする測位装置。
【請求項２】
前記航法メッセージまたはオフライン処理により衛星軌道情報を取得する衛星軌道情報取得手段を備え、
前記測位演算手段は、前記衛星軌道情報を前記線形回帰方程式に加えて、測位衛星位置をも推定演算する請求項１に記載の測位装置。
【請求項３】
前記測位演算手段は、前記線形回帰方程式にカルマンフィルタを適用する請求項１または請求項２に記載の測位装置。
【請求項４】
前記測位演算手段は、前記複数のアンテナの幾何学的中心位置を推定演算する請求項１〜請求項３のいずれかに記載の測位装置。
【請求項５】
前記アンテナおよび受信機が２つである請求項１〜請求項４のいずれかに記載の測位装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、測位用信号を用いて受信機の単独測位を行う測位装置に関するものである。
【背景技術】
【０００２】
従来、測位衛星から送信される測位用信号を用いて測位を行う装置が各種開示されており、測位方法から単独測位と相対測位とが存在する。単独測位の測位演算は、従来、コード擬似距離を用いて受信機の３次元の位置（誤差）と受信機の時計誤差とを未知数とした非線形連立方程式を構成することにより算出される。この際、非線形連立方程式にニュートン法や拡張カルマンフィルタを適用することで測位演算は行われる。また、電離層遅延による誤差や対流圏遅延による誤差の影響を除去するため、これらの遅延による影響を初期状態から「０」に設定したり、１重位相差等の演算処理を行って、これらの遅延量を推定していた。
【非特許文献１】日本測地学会編著，「新訂版ＧＰＳ−人工衛星による精密測位システム−」，社団法人日本測量協会，１９８９年１１月１５日，ｐ．１２１−１４０
【発明の開示】
【発明が解決しようとする課題】
【０００３】
ところが、従来のコード擬似距離のみを用いた測位装置および測位方法では、現実に存在する電離層遅延、対流圏遅延、衛星軌道誤差、および衛星時計誤差を無視したり、あるいは十分とはいえない精度で推定演算を行うので、測位の推定値のばらつきが大きかった。このため、高精度な測位を行うことができなかった。また、コード擬似距離を用いるとともに、電離層遅延情報や対流圏遅延情報を用いて１重位相差処理等を行う場合、前述のコード擬似距離のみを用いる場合よりも推定精度は向上する。しかしながら、推定演算内で１重位相差処理を行う等、推定演算処理が複雑になる。また、１重位相差処理を行った後のノイズの設定により大きく推定演算結果が異なる。このため、推定演算の処理内容が複雑な割には、推定演算結果の精度の向上性が悪く、さらに推定演算処理速度が遅くなる可能性があった。【０００４】
したがって、この発明の目的は、推定演算処理を複雑にすることなく確実に高精度測位を実現することができる測位装置を提供することにある。
【課題を解決するための手段】
【０００５】
この発明の測位装置は、複数の測位衛星から送信される測位用信号をそれぞれに受信し、互いに予め設定された間隔で配置された複数のアンテナを備える。この発明の測位装置は、これら複数のアンテナ毎に設置され、入力された測位用信号をそれぞれに復調する複数の受信機を備える。そして、この発明の測位装置は、各受信機から得られる衛星情報から各受信機と各測位衛星との間のコード擬似距離およびキャリア位相積算値を観測し、各受信機に対する整数値バイアス、各受信機の受信機時計誤差、および各受信機位置を未知数として含み、該各受信機位置を過去の推定演算値で、衛星軌道位置を放送軌道情報による推定軌道位置で１次テイラー級数展開することで線形近似し、整数値バイアス、受信機時計誤差、受信機位置、衛星時計誤差、衛星軌道位置を説明変数に含み、キャリア位相積算値およびコード擬似距離を目的変数とし、線形近似により得られる近似行列を演算要素とする１つの線形回帰方程式を構成し、該線形回帰方程式から整数値バイアスと各受信機の位置とを推定演算して、該各受信機の位置と各受信機間の間隔とから測位を行う測位演算手段と、を備えたことを特徴としている。【０００６】
この構成では、各測位衛星と各受信機とのキャリア位相積算値およびコード擬似距離を観測値とし、受信機の３次元位置、受信機の時計誤差、測位衛星毎の衛星時計誤差、衛星軌道位置、整数値バイアスを未知数として、前記観測値を目的変数とし、前記未知数を説明変数とする線形回帰方程式が構成される。この際、線形回帰方程式の演算要素に、各受信機位置を過去の推定演算値で、衛星軌道位置を放送軌道情報による推定軌道位置で１次テイラー級数展開することで線形近似した近似行列を用いる。そして、この線形回帰方程式をエポック単位で順次演算することで整数値バイアス、受信機時計誤差や受信機位置等の未知数が推定演算される。この際、受信機数が複数であっても用いる線形回帰方程式は１つである。【０００７】
また、この発明は、航法メッセージまたはオフライン処理により衛星軌道情報を取得する衛星軌道情報取得手段を備え、測位演算手段で、衛星軌道情報を線形回帰方程式に加えて、測位衛星位置をも推定演算することを特徴としている。【０００８】
この構成では、衛星軌道情報を観測値として目的変数に加えるとともに、測位衛星位置を未知数として説明変数に加えることで、線形回帰方程式から測位衛星位置が推定演算される。【０００９】
この発明は、測位演算手段で、線形回帰方程式にカルマンフィルタを適用することを特徴としている。【００１０】
この構成では、カルマンフィルタを適用することで、移動する受信機位置の推定も可能となる。【００１１】
また、この発明は、測位演算手段で、複数のアンテナの幾何学的中心位置を推定演算することを特徴としている。【００１２】
この構成では、複数のアンテナの幾何学的中心位置を用いることで、推定される位置の分散および標準偏差が低くなる。【００１３】
また、この発明は、アンテナおよび受信機を２つにすることを特徴としている。【００１４】
この構成では、アンテナおよび受信機の数が２つであるので、１つずつしか備えない場合よりも高い精度を得られ、３つ以上の備えた場合よりも装置が小型化される。
【発明の効果】
【００１５】
この発明によれば、１つの回帰方程式のみを用いて、それぞれの位置間隔が予め決まっている複数の受信機の位置を高精度に推定演算することができるので、これらの測位結果と位置関係とを用いて、高精度な測位を行うことができる。この際、１つの回帰方程式だけしか用いないため、演算処理が簡素化され、素早くかつ少ない演算で測位を行うことができる。【００１６】
また、この発明によれば、衛星軌道情報を用いることにより、測位衛星位置の推定演算が可能となるので、受信機と測位衛星との間の距離をより高精度に測位することができる。これにより、受信機位置をより高精度に測位することができる。【００１７】
また、この発明によれば、カルマンフィルタを用いることで、固定された受信機のみでなく移動する受信機の位置を推定演算することができる。【００１８】
また、この発明によれば、複数のアンテナの幾何学的中心を用いることで、より高精度に測位することができる。【００１９】
また、この発明によれば、アンテナと受信機との数を２つずつにすることで、小型でありながら高精度な測位演算結果が得られる測位装置を構成することができる。
【図面の簡単な説明】
【００２０】
【図１】本発明の実施形態の測位装置の概略構成を示すブロック図である。
【図２】本発明の実施形態の測位装置の処理フローを示すフローチャートである。
【符号の説明】
【００２１】
１０−測位装置
１Ａ，１Ｂ−ＧＰＳ受信機
２−航法メッセージ解析部
３−衛星情報処理部
４−測位演算部
【発明を実施するための最良の形態】
【００２２】
本発明の実施形態に係る測位装置について図を参照して説明する。なお、以下の説明では、ＧＰＳを用いた単独測位装置について説明するが、他の全てのＧＮＳＳ（全地球的航法衛星システム）に適用することができる。
図１は本実施形態の単独測位装置の概略構成を示すブロック図である。
また、図２は本実施形態の単独測位装置の処理フローを示すフローチャートである。【００２３】
図１に示すように、単独測位装置１０は、２個のＧＰＳ受信機１Ａ，１Ｂ、航法メッセージ解析部２、衛星情報処理部３、および測位演算部４を備える。２個のＧＰＳ受信機１Ａ，１ＢはそれぞれＧＰＳアンテナ２１，２２に接続している。ＧＰＳアンテナ２１，２２との間隔ｄは予め所定値に設定され、図示しない記憶部に記憶されている。【００２４】
ＧＰＳアンテナ２１，２２は、ＧＰＳ衛星（測位衛星）からのＧＰＳ信号（測位用信号）を受信して、それぞれＧＰＳ受信機１Ａ，１Ｂに出力する。各ＧＰＳ受信機１Ａ，１Ｂは、入力されたＧＰＳ信号から既知の方法で搬送波のキャリア位相積算値、コード擬似距離を観測する。また、各ＧＰＳ受信機１Ａ，１Ｂは搬送波に重畳された航法メッセージを取得する（Ｓ１→Ｓ２）。各ＧＰＳ受信機１Ａ，１Ｂは、それぞれが取得した航法メッセージを航法メッセージ解析部２に出力し、キャリア位相積算値とコード擬似距離とを測位演算部４に出力する。【００２５】
航法メッセージ解析部２は、入力された航法メッセージを解析して、電離層遅延情報、各ＧＰＳ衛星の時計誤差、軌道情報を取得する。航法メッセージ解析部２は、所定の推定演算モデルを用いて対流圏遅延情報を取得する。航法メッセージ解析部２は、電離層遅延情報、各ＧＰＳ衛星の時計誤差、軌道情報および対流圏遅延情報を衛星情報処理部３に出力する（Ｓ２）。【００２６】
衛星情報処理部３は、ＧＰＳ衛星のエフェメリス情報を用いて測位に用いるＧＰＳ衛星を選定する。この際、衛星情報処理部３は、それぞれのＧＰＳ受信機１Ａ，１Ｂ（アンテナ２１，２２）から出力された各ＧＰＳ衛星の情報が整合するように、ＧＰＳ衛星を選定する。すなわち、受信機１Ａ，１Ｂから出力された衛星情報のうち同じＧＰＳ衛星からのものだけを選択する（Ｓ３）。また、衛星情報処理部３は、選定したＧＰＳ衛星に関する時計遅延情報、軌道情報、電離層遅延情報、および対流圏遅延情報を測位演算部４に出力する。【００２７】
測位演算部４は、入力された各観測値、すなわち、選定された各ＧＰＳ衛星に対するキャリア位相積算値、コード擬似距離、時計誤差情報、電離層遅延情報、対流圏遅延情報を用いて、それぞれの受信機に対して後述する線形回帰方程式を用い、最小二乗法を適用させて、それぞれのＧＰＳ受信機とＧＰＳ衛星に対する整数値バイアスＮ１，Ｎ２と、各ＧＰＳ受信機位置ｕ１，ｕ２とを推定演算する。この際、この線形回帰方程式に対しては、カルマンフィルタを適用することにより、ＧＰＳ衛星、受信機が移動しても位置推定が可能になる。この場合、線形回帰方程式の説明変数に相当する情報に対しては、既知のモデルにより動的方程式を適応する。この推定演算は、推定値の変化が予め設定された所定閾値以下になるまで繰り返し行われる。そして、推定値の変化が所定閾値以下になった時点で、推定値が収束したとして、この推定値を測位演算結果として取得する（Ｓ５）。測位演算部４は、このようにして得られた２つのＧＰＳ受信機の位置と、予め記憶された間隔ｄとを用いて、例えば、アンテナ２１，２２の中間点の位置等を、これらが備えられた測位装置の位置として算出し測位結果として出力する（Ｓ６）。【００２８】
なお、前述の説明では、各ＧＰＳ受信機位置（アンテナ位置）ｕ１，ｕ２をそれぞれ推定演算する方法を示したが、受信機位置ｕ１，ｕ２間と幾何学距離ｄとの関係から拘束条件を設定することで、受信機位置ｕ１，ｕ２の幾何学的中心位置を推定演算することができる。そして、このような２つのＧＰＳ受信機位置の中心位置を推定演算することにより、推定される中心位置の分散は約半分になるので、より高精度に測位演算することができる。【００２９】
ここで、航法メッセージ解析部２、衛星情報処理部３、測位演算部４は、以下のアルゴリズムを実現するマイクロプロセッサ等の数値演算処理器からなる。これらの機能部は複数の数値演算処理器により形成してもよく、１つの数値演算処理器により形成してもよい。【００３０】
このような構成とすることで、測位装置に備えられた２つの受信機（アンテナ）で高精度な単独測位を行い、これらの２つの高精度な単独測位結果を用いてこれらの中間点に対して、さらに高精度の単独測位結果を得ることができる。【００３１】
次に、前述の整数バイアスＮおよび各ＧＰＳ受信機位置ｕ１，ｕ２を推定演算するアルゴリズムについて詳述する。なお、以下の説明では、キャリア位相積算値を単に「キャリア位相」と称す。【００３２】
一般に、受信機ｕ、ＧＰＳ衛星ｐに対するキャリア位相φpL,uの観測方程式は式（１）で表すことができ、コード擬似距離（擬似距離）ρpc,uの観測方程式は式（２）で表される。ここで、マルチパス誤差は微少として無視する。【００３３】
【数１】
0004108738.tif 000002
【００３４】
ここで、λＬはＬ波の波長を示し、ｒｐｕ（ｔ，ｔ−τｐｕ）は時刻ｔでの受信機ｕと時刻（ｔ−τｐｕ）でのＧＰＳ衛星ｐとの距離を示し、δＩｐｕ（ｔ）はＬ１波の電離層遅延を示し、δＴｐｕ（ｔ）はＬ１波、Ｌ２波の対流圏遅延を示し、δｔｕ（ｔ）は真の時刻ｔでの受信機ｕの時計誤差を示し、δｔｐ（ｔ−τｐｕ）は時刻（ｔ−τｐｕ）でのＧＰＳ衛星ｐの時計誤差を示し、ＮｐＬ，ｕはＬ波における受信機ｕと衛星ｐとの間の整数値バイアスを示し、εｐＬ，ｕ（ｔ），ｅｐｃ，ｕ（ｔ）はそれぞれ観測雑音を示す。【００３５】
ここで、ΦpL,u＝λLφpL,uとし、Ｌ１波の周波数をｆL1とし、Ｌ２波の周波数をｆL2とする。これにより、式（１）、式（２）から、ＣＡコードによるコード擬似距離ρpCA,uの観測方程式は式（３）となり、Ｐ，Ｙコードによるコード擬似距離ρpPY,uの観測方程式は式（４）となり、Ｌ１波におけるキャリア位相ΦpL1,uの観測方程式は式（５）となり、Ｌ２波におけるキャリア位相ΦpL2,uの観測方程式は式（６）となる。【００３６】
【数２】
0004108738.tif 000003
【００３７】
ところで、受信機ｕとＧＰＳ衛星ｐとの距離ｒｐｕ（ｔ，ｔ−τｐｕ）は、式（７）で表すことができる。【００３８】
【数３】
0004108738.tif 000004
【００３９】
ここで、未知数である受信機位置ｕ（ｔ）≡［ｘｕ（ｔ），ｙｕ（ｔ），ｚｕ（ｔ）］Ｔを先験的な推定受信機位置ｕ（ｊ）（ｔ）≡［ｘｕ（ｊ）（ｔ），ｙｕ（ｊ）（ｔ），ｚｕ（ｊ）（ｔ）］Ｔのまわりで、衛星軌道位置ｓpを航法メッセージに含まれる放送軌道情報による推定軌道位置ｓｅpのまわりで１次のテイラー級数展開を行い、ｒｐｕ（ｔ）を線形近似すると、式（８）、（９）となる。【００４０】
【数４】
0004108738.tif 000005
【００４１】
この線形近似式を式（３）〜（６）に代入して、さらに、式（１０）、（１１）の関係を適用する。【００４２】
【数５】
0004108738.tif 000006
【００４３】
この結果として次式が得られる。【００４４】
【数６】
0004108738.tif 000007
【００４５】
このように前述の式（１２）〜（１５）は、キャリア位相、コード擬似距離を目的変数とし、受信機位置、電離層遅延、対流圏遅延、整数値バイアスを説明変数とする近似的な線形回帰方程式となる。【００４６】
ここで、線形近似要素を示す行列を、式（１６）、（１７）で表す。【００４７】
【数７】
0004108738.tif 000008
【００４８】
これにより、式（１２）〜（１５）は次に示す行列式で表される。【００４９】
【数８】
0004108738.tif 000009
【００５０】
この線形回帰方程式に対して、コード擬似距離、キャリア位相積算量の観測値を用いて最小二乗法から受信機位置ｕ、整数バイアスＮを推定演算する。この際、測位衛星位置ｓ、衛星時計誤差ｃδｔsについても同時に推定演算することができる。【００５１】
ここで、この線形回帰方程式を２つのＧＰＳ受信機ｕ１，ｕ２に適用すると、式（１８）から次式が得られる。この際、ＧＰＳ受信機ｕ１，ｕ２は単一の測位装置に設置されているので、電離層遅延や対流圏遅延は共通なものと見なすことができる。したがって、【００５２】
【数９】
0004108738.tif 000010
【００５３】
このような線形回帰方程式であっても、各受信機ｕ１，ｕ２におけるコード擬似距離、キャリア位相積算量を用いることで、前述の式（１８）と同様に受信機位置ｕ１，ｕ２およびそれぞれの整数値バイアスＮを推定演算することができる。すなわち、１つの線形回帰方程式で、２つの受信機位置を推定演算することができる。【００５４】
ここで、式（１９）においてＬ１波のみ関係する項のみを残す。すなわち、Ｌ１波のみが捕捉、追尾できて、ＣＡコードのコード擬似距離とＬ１波のキャリア位相積算量とが得られる場合を考えると、式（１９）より、式（２０）が得られる。【００５５】
【数１０】
0004108738.tif 000011
【００５６】
この際、受信機ｕ１，ｕ２で同一の時計を用いることにより、受信機ｕ２の時計誤差は受信機ｕ１の誤差による得られるものとして線形回帰方程式には使用しない。【００５７】
ここで、受信機ｕ１，ｕ２間の距離をｄと設定し、設置誤差をεdとすると、次式の拘束条件を加える。【００５８】
【数１１】
0004108738.tif 000012
【００５９】
そして、適当な推定値Ｘｅ１，Ｙｅ１，Ｚｅ１，Ｘｅ２，Ｙｅ２，Ｚｅ２の周りでテイラー展開することにより、線形回帰方程式の観測値を１つ増加させることができる。【００６０】
このように、式（２０）に示す線形回帰方程式を用いることで、Ｌ１波の捕捉、追尾のみで受信機ｕ１，ｕ２の位置、測位衛星位置ｓ、衛星時計誤差ｃδｔsおよび整数値バイアスＮを推定演算することができる。すなわち、式（１９）に示す線形回帰方程式を用いる場合のように、Ｌ１波、Ｌ２波ともに捕捉、追尾を行わなければならないのではなく、Ｌ１波のみでよいので、現状で一般に入手可能なＧＰＳ受信機を利用して、高精度な測位を行うことができる。さらに、線形回帰方程式の演算数が減少するので、より高速に推定演算を行うことができる。この際、測位衛星数ｎｓとすると、未知数である説明変数が８ｎｓ＋７個であり、観測値である目的変数が１０ｎｓ個、拘束条件が１個であるので、８ｎｓ＋６≦１０ｎｓで、ｎｓ≧３ならば、すなわち、連続して観測可能な衛星数が３個以上あれば推定演算することができる。【００６１】
さらに、式（２０）より、第１行要素と第３行要素の差分、第２行要素と第４行要素との差分を用いて変形すると、式（２１）になる。【００６２】
【数１２】
0004108738.tif 000013
【００６３】
さらに、整数値バイアス項を変形して、式（２２）が得られる。【００６４】
【数１３】
0004108738.tif 000014
【００６５】
このように、線形回帰方程式をより簡素化することができる。これにより、さらに高速に推定演算を行うことができ、測位スピードを向上することができる。【００６６】
以上のように、本実施形態の構成を用いることで、１つの線形回帰方程式だけで、２つの受信機位置、測位衛星位置、整数値バイアスを推定演算することができ、簡素な構造および処理演算で高精度な測位を高速に行うことができる。【００６７】
なお、前述の説明では、ＧＰＳ受信機数を２つとしたが、３つ以上のＧＰＳ受信機を用いてもよく、使用するＧＰＳ受信機数が増加するほど、処理演算量は増加するものの、より高精度な測位を行うことができる。また、前述の説明では、Ｌ１波およびＬ２波を用いて推定演算する例を示したが、今後増加する、Ｌ５波等の搬送波を演算に加えてもよい。
【図１】
0004108738.tif 000015
【図２】
0004108738.tif 000016