(57)【要約】【課題】航法衛星および衛星信号受信機のそれぞれに起因する周波数間バイアスを推定する。【解決手段】複数周波数の衛星信号を送信する航法衛星から放送される群遅延情報から得られる放送衛星バイアスにより、各航法衛星までの位相距離の周波数間の差である周波数間差分を補正し、各航法衛星における補正した周波数間差分の最小値である最小周波数間差分のうち最頻出の値を基準衛星バイアスとし、各航法衛星における最小周波数間差分と基準衛星バイアスとの差により衛星放送バイアスを補正して衛星バイアス推定値を算出し、この衛星バイアス推定値、および電離層電子密度分布モデルから推定されるＴＥＣ値から受信機バイアス推定値を算出する。【選択図】図１
2010025679.tif 000002
【特許請求の範囲】
【請求項１】
複数の航法衛星からそれぞれ送信される第１および第２の周波数の衛星信号に含まれる群遅延情報に基づいて、前記航法衛星内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差の初期値である衛星放送バイアスを前記航法衛星ごとに算出し、複数の観測時刻において地上で受信した前記第１および第２の周波数の衛星信号の伝搬経路における伝搬遅延量と位相量とに基づいて算出される前記各航法衛星までの位相距離の周波数間の差である周波数間差分を、前記衛星放送バイアスを用いて補正するステップと、
前記複数の観測時刻における補正した前記周波数間差分の最小値である最小周波数間差分を前記航法衛星ごとに算出し、前記各航法衛星における前記最小周波数間差分のうち最頻出の値である基準衛星バイアスを算出するステップと、
前記各航法衛星における前記最小周波数間差分と前記基準衛星バイアスとの差を用いて前記衛星放送バイアスを補正し、この補正した前記衛星放送バイアスを、前記航法衛星内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差である衛星バイアスの推定値である第１の衛星バイアス推定値として前記航法衛星ごとに算出するステップと、
電離層電子密度分布モデルを用いて、前記最小周波数間差分が前記基準衛星バイアスに最も近い前記航法衛星において、前記最小周波数間差分をとる観測時刻における前記衛星信号の受信位置から当該航法衛星までの総電子数を算出するステップと、
前記基準衛星バイアスと前記総電子数との差を、前記衛星信号を受信した衛星信号受信機内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差である受信機バイアスの推定値である受信機バイアス推定値として算出するステップと
を含むことを特徴とする衛星バイアスおよび受信機バイアスの推定方法。
【請求項２】
複数日を含む第１の期間内の複数の第２の期間ごとの前記受信機バイアス推定値の前記第１の期間内における平均値である補正受信機バイアス推定値を算出するステップと、
前記第１の期間内の複数の観測時刻における前記各航法衛星までの位相距離の周波数間の差、前記第１の期間内の観測時刻ごとに電離層電子密度分布モデルを用いて算出される前記衛星信号の受信位置から前記各航法衛星までの総電子数、および前記補正受信機バイアス推定値を用いて、前記第１の期間内の観測時刻ごとの前記衛星バイアスの推定値である第２の衛星バイアス推定値を前記航法衛星ごとに算出し、前記第２の衛星バイアス推定値の前記第１の期間内における平均値である補正衛星バイアス推定値を前記航法衛星ごとに算出するステップと
をさらに含むことを特徴とする請求項１に記載の衛星バイアスおよび受信機バイアスの推定方法。
【請求項３】
複数の航法衛星からそれぞれ送信される第１および第２の周波数の衛星信号に含まれる群遅延情報に基づいて、前記航法衛星内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差の初期値である衛星放送バイアスを前記航法衛星ごとに算出し、複数の観測時刻において地上で受信した前記第１および第２の周波数の衛星信号の伝搬経路における伝搬遅延量と位相量とに基づいて算出される前記各航法衛星までの位相距離の周波数間の差である周波数間差分を、前記衛星放送バイアスを用いて補正するステップと、
前記複数の観測時刻における補正した前記周波数間差分の最小値である最小周波数間差分を前記航法衛星ごとに算出し、前記各航法衛星における前記最小周波数間差分のうち最頻出の値である基準衛星バイアスを算出するステップと、
前記各航法衛星における前記最小周波数間差分と前記基準衛星バイアスとの差を用いて前記衛星放送バイアスを補正し、この補正した前記衛星放送バイアスを、前記航法衛星内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差である衛星バイアスの推定値である第１の衛星バイアス推定値として前記航法衛星ごとに算出するステップと、
電離層電子密度分布モデルを用いて、前記最小周波数間差分が前記基準衛星バイアスに最も近い前記航法衛星において、前記最小周波数間差分をとる観測時刻における前記衛星信号の受信位置から当該航法衛星までの総電子数を算出するステップと、
前記基準衛星バイアスと前記総電子数との差を、前記衛星信号を受信した衛星信号受信機内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差である受信機バイアスの推定値である受信機バイアス推定値として算出するステップと、
複数日を含む第１の期間内の複数の第２の期間ごとの前記受信機バイアス推定値の前記第１の期間内における平均値である補正受信機バイアス推定値を算出するステップと、
前記第１の期間内の複数の観測時刻における前記各航法衛星までの位相距離の周波数間の差、前記第１の期間内の観測時刻ごとに電離層電子密度分布モデルを用いて算出される前記衛星信号の受信位置から前記各航法衛星までの総電子数、および前記補正受信機バイアス推定値を用いて、前記第１の期間内の観測時刻ごとの前記衛星バイアスの推定値である第２の衛星バイアス推定値を前記航法衛星ごとに算出し、前記第２の衛星バイアス推定値の前記第１の期間内における平均値である補正衛星バイアス推定値を前記航法衛星ごとに算出するステップ
との組み合わせを単位工程とし、２回目以降の前記単位工程では前記衛星放送バイアスを前回の前記単位工程で算出した前記補正衛星バイアス推定値に置き換え、前記単位工程で算出される前記補正受信機バイアス推定値、および少なくとも１つの前記航法衛星に対応する前記補正衛星バイアス推定値が、前回の前記単位工程で算出した値との差がそれぞれ予め設定した閾値よりも小さくなるまで前記単位工程を繰り返すことを特徴とする衛星バイアスおよび受信機バイアスの推定方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、航法衛星および衛星信号受信機に起因する複数周波数の衛星信号の周波数間バイアスを推定する衛星バイアスおよび受信機バイアスの推定方法に関する。
【背景技術】
【０００２】
ＧＰＳ（Global Positioning System）衛星、Galileo衛星、準天頂衛星などの航法衛星から送信される複数周波数の衛星信号における各周波数の伝搬遅延量を観測することにより、衛星信号の通過経路に存在する総電子数（ＴＥＣ：Total Electron Content）の推定値を算出することができる（例えば、非特許文献１参照）。
【０００３】
このように算出されるＴＥＣは、例えば、ＨＦ帯などの電波通信における電波の伝搬経路を算出する際に用いるＩＲＩ（International Reference Ionosphere）モデル等の電離層電子密度分布モデルを修正し、良好な通信環境を構築するために用いられる。
【非特許文献１】Pratap Misra,PerEnge,“GLOBAL POSITIONING SYSTEM Signals,Measurements,and Performance．”Ganga-Jamuna Press,2001
【発明の開示】
【発明が解決しようとする課題】
【０００４】
しかしながら、上述のＴＥＣの推定値には、航法衛星内部および衛星信号受信機内部でそれぞれ生じる複数周波数の衛星信号の周波数間の遅延時間差（周波数間バイアス）である衛星バイアスおよび受信機バイアスの影響が含まれてしまう。
【０００５】
衛星バイアスおよび受信機バイアスを求めることができれば、その影響を除いたＴＥＣの推定値を算出することができるため、衛星バイアスおよび受信機バイアスを推定する方法が要望されていた。
【０００６】
本発明は上記に鑑みてなされたもので、航法衛星および衛星信号受信機のそれぞれに起因する周波数間バイアスを推定する衛星バイアスおよび受信機バイアスの推定方法を提供することを目的とする。
【課題を解決するための手段】
【０００７】
上記目的を達成するため、本発明の衛星バイアスおよび受信機バイアスの推定方法は、複数の航法衛星からそれぞれ送信される第１および第２の周波数の衛星信号に含まれる群遅延情報に基づいて、前記航法衛星内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差の初期値である衛星放送バイアスを前記航法衛星ごとに算出し、複数の観測時刻において地上で受信した前記第１および第２の周波数の衛星信号の伝搬経路における伝搬遅延量と位相量とに基づいて算出される前記各航法衛星までの位相距離の周波数間の差である周波数間差分を、前記衛星放送バイアスを用いて補正するステップと、前記複数の観測時刻における補正した前記周波数間差分の最小値である最小周波数間差分を前記航法衛星ごとに算出し、前記各航法衛星における前記最小周波数間差分のうち最頻出の値である基準衛星バイアスを算出するステップと、前記各航法衛星における前記最小周波数間差分と前記基準衛星バイアスとの差を用いて前記衛星放送バイアスを補正し、この補正した前記衛星放送バイアスを、前記航法衛星内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差である衛星バイアスの推定値である第１の衛星バイアス推定値として前記航法衛星ごとに算出するステップと、電離層電子密度分布モデルを用いて、前記最小周波数間差分が前記基準衛星バイアスに最も近い前記航法衛星において、前記最小周波数間差分をとる観測時刻における前記衛星信号の受信位置から当該航法衛星までの総電子数を算出するステップと、前記基準衛星バイアスと前記総電子数との差を、前記衛星信号を受信した衛星信号受信機内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差である受信機バイアスの推定値である受信機バイアス推定値として算出するステップとを含むことを特徴とする。
【０００８】
また、本発明の衛星バイアスおよび受信機バイアスの推定方法は、複数日を含む第１の期間内の複数の第２の期間ごとの前記受信機バイアス推定値の前記第１の期間内における平均値である補正受信機バイアス推定値を算出するステップと、前記第１の期間内の複数の観測時刻における前記各航法衛星までの位相距離の周波数間の差、前記第１の期間内の観測時刻ごとに電離層電子密度分布モデルを用いて算出される前記衛星信号の受信位置から前記各航法衛星までの総電子数、および前記補正受信機バイアス推定値を用いて、前記第１の期間内の観測時刻ごとの前記衛星バイアスの推定値である第２の衛星バイアス推定値を前記航法衛星ごとに算出し、前記第２の衛星バイアス推定値の前記第１の期間内における平均値である補正衛星バイアス推定値を前記航法衛星ごとに算出するステップとをさらに含むことを特徴とする。
【０００９】
また、本発明の衛星バイアスおよび受信機バイアスの推定方法は、複数の航法衛星からそれぞれ送信される第１および第２の周波数の衛星信号に含まれる群遅延情報に基づいて、前記航法衛星内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差の初期値である衛星放送バイアスを前記航法衛星ごとに算出し、複数の観測時刻において地上で受信した前記第１および第２の周波数の衛星信号の伝搬経路における伝搬遅延量と位相量とに基づいて算出される前記各航法衛星までの位相距離の周波数間の差である周波数間差分を、前記衛星放送バイアスを用いて補正するステップと、前記複数の観測時刻における補正した前記周波数間差分の最小値である最小周波数間差分を前記航法衛星ごとに算出し、前記各航法衛星における前記最小周波数間差分のうち最頻出の値である基準衛星バイアスを算出するステップと、前記各航法衛星における前記最小周波数間差分と前記基準衛星バイアスとの差を用いて前記衛星放送バイアスを補正し、この補正した前記衛星放送バイアスを、前記航法衛星内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差である衛星バイアスの推定値である第１の衛星バイアス推定値として前記航法衛星ごとに算出するステップと、電離層電子密度分布モデルを用いて、前記最小周波数間差分が前記基準衛星バイアスに最も近い前記航法衛星において、前記最小周波数間差分をとる観測時刻における前記衛星信号の受信位置から当該航法衛星までの総電子数を算出するステップと、前記基準衛星バイアスと前記総電子数との差を、前記衛星信号を受信した衛星信号受信機内部で生じる前記第１および第２の周波数の衛星信号の周波数間の遅延時間差である受信機バイアスの推定値である受信機バイアス推定値として算出するステップと、複数日を含む第１の期間内の複数の第２の期間ごとの前記受信機バイアス推定値の前記第１の期間内における平均値である補正受信機バイアス推定値を算出するステップと、前記第１の期間内の複数の観測時刻における前記各航法衛星までの位相距離の周波数間の差、前記第１の期間内の観測時刻ごとに電離層電子密度分布モデルを用いて算出される前記衛星信号の受信位置から前記各航法衛星までの総電子数、および前記補正受信機バイアス推定値を用いて、前記第１の期間内の観測時刻ごとの前記衛星バイアスの推定値である第２の衛星バイアス推定値を前記航法衛星ごとに算出し、前記第２の衛星バイアス推定値の前記第１の期間内における平均値である補正衛星バイアス推定値を前記航法衛星ごとに算出するステップとの組み合わせを単位工程とし、２回目以降の前記単位工程では前記衛星放送バイアスを前回の前記単位工程で算出した前記補正衛星バイアス推定値に置き換え、前記単位工程で算出される前記補正受信機バイアス推定値、および少なくとも１つの前記航法衛星に対応する前記補正衛星バイアス推定値が、前回の前記単位工程で算出した値との差がそれぞれ予め設定した閾値よりも小さくなるまで前記単位工程を繰り返すことを特徴とする。
【発明の効果】
【００１０】
本発明の衛星バイアスおよび受信機バイアスの推定方法によれば、航法衛星および衛星信号受信機のそれぞれに起因する周波数間バイアスを推定することができる。
【発明を実施するための最良の形態】
【００１１】
以下、本発明を実施するための最良の形態について、図面を参照して説明する。
【００１２】
図１は、本発明の実施の形態に係る衛星バイアスおよび受信機バイアスの推定方法を実施する地上システムの構成を示すブロック図である。図１に示す地上システムは、衛星信号受信アンテナ１と、衛星信号受信機２と、インターネットデータ処理装置３と、バイアス推定処理装置４と、データ記録装置５とを備える。上記各装置は、ＬＡＮ（Local Area Network）等からなるネットワーク６を介して互いに接続されている。
【００１３】
衛星信号受信機２は、ＧＰＳ衛星、Galileo衛星、準天頂衛星などの航法衛星から送信される複数周波数の衛星信号を、衛星信号受信アンテナ１を介して受信し、この衛星信号をネットワーク６を介してバイアス推定処理装置４へ供給する。
【００１４】
インターネットデータ処理装置３は、ＧＥＯＮＥＴ収集データ処理部３１と、中継部３２とを備える。
【００１５】
ＧＥＯＮＥＴ収集データ処理部３１は、一般に公開されている電離層関連の情報や、国土地理院が管理するＧＰＳ受信観測網（ＧＥＯＮＥＴ）のＧＰＳ観測データや、国際的にＧＰＳ観測結果を公開しているＩＧＳ（International GPS Service for Geodynamics）のＧＰＳ観測データなどをインターネット７を介して取得する。
【００１６】
中継部３２は、スイッチングハブあるいはルータにより構成され、ＧＥＯＮＥＴ収集データ処理部３１で取得した各種データを、ネットワーク６を介してバイアス推定処理装置４およびデータ記録装置５へ出力する。
【００１７】
なお、インターネットデータ処理装置３は、外部とのつながりあるため、ＧＥＯＮＥＴ収集データ処理部３１および中継部３２は、ファイアウォール機能を有するものとする。
【００１８】
バイアス推定処理装置４は、衛星信号受信機２で受信した衛星信号、ＧＥＯＮＥＴ収集データ処理部３１で取得したＧＰＳ観測データ等を用いて、航法衛星内部で生じる複数周波数の衛星信号の周波数間の遅延時間差である衛星バイアスの推定値、および衛星信号受信機２内部で生じる複数周波数の衛星信号の周波数間の遅延時間差である受信機バイアスの推定値等を算出する処理を行う。
【００１９】
データ記録装置５は、上記各装置で得られた各種演算結果や各種データを、ネットワーク６を介して受け取り、これらを記憶する。
【００２０】
上記構成の地上システムにおいて、航法衛星から送信される複数周波数の衛星信号は、アンテナ１を介して衛星信号受信機２により受信され、ネットワーク６を介してバイアス推定処理装置４へ供給される。バイアス推定処理装置４は、受信した衛星信号、ＧＥＯＮＥＴ収集データ処理部３１で取得したＧＰＳ観測データ等を用いて、以下に説明する手順により、衛星バイアスの推定値、受信機バイアスの推定値等を算出する処理を行う。
【００２１】
ここでは、衛星信号としてＧＰＳのＬ１周波数（１５７５．４２ＭＨｚ）、Ｌ２周波数の２周波（１２２７．６０ＭＨｚ）の信号を用いる場合について説明する。なお、電離層がダイナミックに変化するような時間帯における衛星の観測結果を用いると推定値の算出が複雑になるため、ローカル時刻（ＳＬＴ：Sun Local Time）の深夜２時頃の観測結果を主に使用する。
【００２２】
以下の（数式１）により衛星測位における擬似距離（コード距離、シュードレンジ）が算出され、（数式２）により位相が算出される。なお、（数式１），（数式２）における添え字中の“Ｌ１ｏｒＬ２”はＬ１，Ｌ２のいずれかを示す。
【数１】
2010025679.tif 000003
【００２３】
２つの周波数間の擬似距離の差分Δρ、位相距離の差分（周波数間差分）ΔＦは、それぞれ（数式３），（数式４）のように示される。
【数２】
2010025679.tif 000004
【００２４】
ここで、ρは擬似距離、φは位相、ｒは真の距離、ｃは光速、ｆは衛星信号の周波数、λは衛星信号の波長、δｔuは受信機時刻誤差、δｔsは衛星時刻誤差、δｔu,L1orL2biasは受信機周波数依存ハードウェアバイアス、δｔs,L1orL2biasは衛星ハードウェア依存バイアス、Ｉは電離層伝搬遅延量、Ｔは対流圏伝搬遅延量、Ｎambは整数不確定値、εは観測誤差である。
【００２５】
上記（数式４）には、観測値からは分からない不確定値（λL1ΔNL1,amb−λL2ΔNL2,amb）が含まれている。そこで、（数式３）と（数式４）とを組み合わせ、（数式５）により上記不確定値を消去する。
【数３】
2010025679.tif 000005
【００２６】
ここで、添え字ｕはデータの番号、Ｍは連続的に複数の観測時刻で収集できたサンプル数の合計、ＴＥＣk,uはＧＰＳ衛星ｋの伝搬経路上の総電子数（ 個／ｍ2 ）、ΔRx,biasは受信機バイアス、Δksat,biasはＧＰＳ衛星ｋの衛星バイアスである。衛星信号受信機２とＧＰＳ衛星によるバイアス項Δkbias（＝ΔRx,bias＋Δksat,bias）は、数ヶ月程度は一定であると仮定する。
【００２７】
次に、衛星バイアスの推定値を算出する手順を説明する。
【００２８】
ＧＰＳ衛星の場合、衛星信号の航法メッセージに群遅延（Tgd）情報が放送されている。群遅延情報には、衛星の製造時における衛星内部で生じるＬ１周波数とＬ２周波数との間の遅延時間差を示す情報が含まれている。
【００２９】
この群遅延情報を用いて、各衛星内部で生じるＬ１周波数とＬ２周波数の衛星信号の周波数間の遅延時間差の初期値である衛星放送バイアスΔL1-L2,Broadを、以下の（数式６）により求めることができる。
【数４】
2010025679.tif 000006
【００３０】
各衛星における衛星放送バイアスΔL1-L2,Broad,kを用いて、（数式３）、（数式５）で求められた値を、以下の（数式７）、（数式８）により補正する。
【数５】
2010025679.tif 000007
【００３１】
上記（数式８）で算出した補正後の周波数間差分ΔＦM1k,uについて、衛星ごとに最小値を算出し、各衛星の最小値（最小周波数間差分）をΔＦM1k,minとする。
【００３２】
そして、例えば図２のように分布する各衛星の最小周波数間差分ΔＦM1k,minの中で、最頻出の値を基準衛星バイアスΔs0とする。
【００３３】
さらに、（数式９）により、各衛星における最小周波数間差分ΔＦM1k,minの、基準衛星バイアスΔs0からの偏差を、衛星放送バイアスΔL1-L2,Broad,kを補正するための補正値ΔmodL1-L2,kとして算出する。
【数６】
2010025679.tif 000008
【００３４】
この補正値ΔmodL1-L2,kを用いて衛星放送バイアスΔL1-L2,Broad,kを補正した値（ΔL1-L2,Broad,k＋ΔmodL1-L2,k）を、衛星バイアス推定値とする。
【００３５】
なお、衛星放送バイアスがない場合は、（数式３）と（数式７）、（数式５）と（数式８）はそれぞれ同じになる。この場合、（数式５）の周波数間差分ΔＦ´k,uの各衛星における最小値の中で最頻出の値を衛星バイアス推定値とする。
【００３６】
次に、受信機バイアスの推定値を算出する手順を説明する。本実施の形態では、電離層電子密度分布モデルを用いて垂直方向のＴＥＣを算出し、これを用いて受信機バイアスの推定値を算出する。
【００３７】
まず、最小周波数間差分ΔFM1k,minが基準衛星バイアスΔs0値に最も近い衛星を見つける。この衛星をｋ＝ｖとする。
【００３８】
次に、この衛星ｖが最小周波数間差分ΔFM1k,minをとる観測時刻を求め、この観測時刻について修正した電離層電子密度分布モデルにより、衛星ｖまでの垂直方向のＴＥＣを計算する。なお、時刻としては、国際標準時（ＵＴ）、日本標準時（ＪＳＴ）等を用いることができる。また、電離層電子密度分布モデルとしては、ＩＲＩモデルやGallagherのモデルを用いることができる。他のモデルでも以降の処理は同様である。
【００３９】
電離層電子密度分布モデルを用いた垂直方向のＴＥＣ計算では、高さＨ＝４００ｋｍ以下とそれ以上の区間に分けてＴＥＣを計算する。
【００４０】
高さ４００ｋｍ以下における垂直方向のＴＥＣ値であるＴＥＣ（Ｈ≦４００ｋｍ）modelは、次の（数式１０）により求められる。
【数７】
2010025679.tif 000009
【００４１】
ここで、ＩＯＮmodel,uは、高さΔＨごとに電離層電子密度分布モデルにより推定した電離層電子密度である。
【００４２】
次に、高さ４００ｋｍから衛星高度Ｈsat（約２００００ｋｍ）までのＴＥＣ値を算出する手順について説明する。
【００４３】
次の（数式１１）により、電離層電子密度分布モデルにより推定した高さ４００ｋｍから１０００ｋｍまでの電離層電子密度ＩＯＮmodel,u（１０００≧Ｈ＞４００ｋｍ）を指数関数近似する。
【数８】
2010025679.tif 000010
【００４４】
ここで、ＩＯＮmodel（Ｈ＝４００ｋｍ）は、電離層電子密度分布モデルにより推定した高さ４００ｋｍにおける電離層電子密度である。
【００４５】
上記（数式１１）を直線関数で近似する。近似関数を求める方法としては、最小２乗法等がある。これにより求められた係数をａ，ｂとすると、修正モデル推定値ＩＯＮ´model,uは、次の（数式１２）で表される。
【数９】
2010025679.tif 000011
【００４６】
この指数関数（数式１２）を、（数式１３）に示すように、高さ４００ｋｍから衛星高度Ｈsatまで解析的に積分する。
【数１０】
2010025679.tif 000012
【００４７】
この（数式１３）により求められた値ＴＥＣ（Ｈsat≧Ｈ＞４００ｋｍ）modelを、高さ４００ｋｍから衛星高度Ｈsatまでの垂直方向モデルＴＥＣ値とする。
【００４８】
総合したＴＥＣ値を（数式１４）に示す。
【数１１】
2010025679.tif 000013
【００４９】
ここで、上記（数式１４）のＴＥＣ（Ｈsat≧Ｈ＞０ｋｍ）modelは、純粋に垂直方向のＴＥＣ値を示している。しかしながら、衛星が必ずしも垂直方向に位置しているとは限らないため、次の（数式１５）に示すように、衛星仰角分の補正係数αを（数式１４）に作用して調整する。
【数１２】
2010025679.tif 000014
【００５０】
補正係数αの一例を（数式１６）に示す。
【数１３】
2010025679.tif 000015
【００５１】
なお、補正係数α＝１としてもよい。
【００５２】
次いで、（数式１７）に示すように、前述の基準衛星バイアスΔs0と、上記（数式１５）で算出されるＴＥＣ（Ｈsat≧Ｈ＞０ｋｍ）´modelとの差を受信機バイアス推定値Δ´Rx,biasとする。すなわち、基準衛星バイアスΔs0が衛星信号受信機によるバイアスによる効果とみなす。
【数１４】
2010025679.tif 000016
【００５３】
そして、（数式１８）に示すように、（数式５）に放送衛星バイアスΔL1-L2,Broad,k、放送衛星バイアスの補正値ΔmodL1-L2,k、および受信機バイアス推定値Δ´Rx,biasを代入してバイアス項Δkbiasを補正する。
【数１５】
2010025679.tif 000017
【００５４】
この（数式１８）より、最終的なＴＥＣの推定値ＴＥＣ´k,uは、次の（数式１９）のようになる。
【数１６】
2010025679.tif 000018
【００５５】
次に、受信機バイアス推定値を補正する手順を説明する。
【００５６】
電離層電子密度分布モデルは、月平均レベルの精度のデータを示すものであり、電離層電子密度分布モデルを用いて算出した上述の垂直方向のＴＥＣ値に現実の値との誤差が生じる可能性がある。
【００５７】
そこで、数日間〜数ヶ月間の期間（第１の期間）において、第１の期間内の複数の期間（第２の期間）ごとに、上述の手順により受信機バイアス推定値を算出し、以下の（数式２０）により、第１の期間において算出した複数の受信機バイアス推定値の平均値をとることにより、受信機バイアス推定値の補正を行う。
【数１７】
2010025679.tif 000019
【００５８】
ここで、ｗは複数ある衛星信号受信機の識別子であり、Ｍは衛星信号受信機ｗで第１の期間内に算出した受信機バイアス推定値の数である。この平均操作により補正受信機バイアス推定値ΔwRx,biasが算出される。
【００５９】
なお、第１の期間内の複数の第２の期間ごとに算出した複数の受信機バイアス推定値の頻度分布を求め、最大頻度となる付近の相当する値を衛星信号受信機ｗの補正受信機バイアス推定値としてもよい。
【００６０】
次に、上述のように算出した補正受信機バイアス推定値ΔwRx,biasを用いて、補正衛星バイアス推定値を算出する。
【００６１】
衛星信号受信機ｗで数日間〜数ヶ月間の期間（第１の期間）内の複数の観測時刻で観測した衛星信号に基づいて、次の（数式２１）の値を衛星ごとに計算する。
【数１８】
2010025679.tif 000020
【００６２】
ここで、ｋは衛星の識別子、ｗは衛星信号受信機の識別子（受信した場所）、Ｍ´は衛星信号受信機の数、ｕは衛星ごとのデータの識別子（１日単位程度）、Ｎは衛星ごとのデータの数である。（数式２１）に示された衛星バイアス値Δw,ksat,bias,uは、受信した場所ｗやデータの識別子ｕに依存しないと仮定する。数日〜数ヶ月で衛星バイアスは大きく変化しないからである。そこで、（数式２１）の値を、以下の（数式２２）により、衛星ごとに衛星信号受信機およびデータで平均をとる。
【数１９】
2010025679.tif 000021
【００６３】
この平均操作により補正衛星バイアス推定値Δksat,bias,estが算出される。
【００６４】
なお、衛星ごとに（数式２１）による衛星バイアス値Δw,ksat,bias,uの頻度分布を求め、最頻度値をその衛星の補正衛星バイアス推定値Δksat,bias,estとしてもよい。
【００６５】
このようにして算出した補正衛星バイアス推定値Δksat,bias,estを、前述の衛星放送バイアスΔL1-L2,Broad,kと置き換えて（数式８）に代入し、以降の処理を上述と同様に実施し、
補正受信機バイアス推定値ΔwRx,bias、および補正衛星バイアス推定値Δksat,bias,estを算出する。
【００６６】
この一連の処理（単位工程）を複数回繰り返し、補正受信機バイアス推定値ΔwRx,bias、および少なくとも１つの衛星の補正衛星バイアス推定値Δksat,bias,estが、前回の単位工程で算出した値との差がそれぞれ予め設定した閾値よりも小さくなった時点で処理を止め、最終的な受信機バイアス推定値および衛星バイアス推定値とする。
【００６７】
以上説明したように本実施の形態では、衛星から放送される群遅延情報から得られる放送衛星バイアスにより、各衛星までの位相距離の周波数間の差である周波数間差分を補正し、各衛星における補正した周波数間差分の最小値である最小周波数間差分のうち最頻出の値を基準衛星バイアスとし、各衛星における最小周波数間差分と基準衛星バイアスとの差により衛星放送バイアスを補正して衛星バイアス推定値を算出し、この衛星バイアス推定値、および電離層電子密度分布モデルから推定されるＴＥＣ値から受信機バイアス推定値を算出することで、正確な衛星バイアス推定値および受信機バイアス推定値が得られる。
【００６８】
また、数日間〜数ヶ月間の期間にわたって観測した衛星信号に基づいて、それぞれ複数の受信機バイアス推定値および衛星バイアス推定値を算出し、これらの平均値をとることにより、より正確な値に補正した補正衛星バイアス推定値および補正受信機バイアス推定値を得ることができる。また、算出した補正衛星バイアス推定値を衛星放送バイアスと置き換えて一連の処理を繰り返すことで、さらに精度を高めることができる。
【００６９】
その結果、電離層電子密度分布モデルを用いたＴＥＣ推定精度を向上させ、ＨＦ帯などの電波通信における電波の伝搬経路を正確に算出することができ、良好な通信が可能となる。また、航法衛星を用いた測位精度も向上する。
【００７０】
なお、本発明は上記実施の形態そのままに限定されるものではなく、実施段階ではその要旨を逸脱しない範囲で構成要素を変形して具体化できる。また、上記実施の形態に開示されている複数の構成要素の適宜な組み合せにより、種々の発明を形成できる。例えば、実施の形態に示される全構成要素から幾つかの構成要素を削除してもよい。
【図面の簡単な説明】
【００７１】
【図１】本発明の実施の形態に係る衛星バイアスおよび受信機バイアスの推定方法を実施する地上システムの構成を示すブロック図である。
【図２】最小周波数間差分の分布を示す図である。
【符号の説明】
【００７２】
１  衛星信号受信アンテナ
２  衛星信号受信機
３  インターネットデータ処理装置
４  バイアス推定処理装置
５  データ記録装置
６  ネットワーク
３１  ＧＥＯＮＥＴ収集データ処理部
３２  中継部
【図１】
2010025679.tif 000022
【図２】
2010025679.tif 000023