ＧＮＳＳ暦の自律的な受信機内予測方法及び装置
全地球航法衛星システム（ＧＮＳＳ）の軌道状態及びクロック状態の自律的な受信機内予測の方法及び装置を記載する。ＧＮＳＳ放送メッセージのみを使用し、周期的な外部通信情報を必要としない。地球方位情報をＧＮＳＳ放送暦から抽出する。地球方位パラメータの正確な推定により、慣性基準フレーム内で、最適に適合するＧＮＳＳ軌道を、時間内に前向きに伝播させることができる。次に、推定した地球方位パラメータを用いて、予測した軌道を地球中心−地球固定（ＥＣＥＦ）座標に変換して、信号の取得に当たりＧＮＳＳ受信機を支援するために使用する。ＧＮＳＳ衛星のクロック状態も放送暦から抽出して、このデータに、パラメータ化したクロック動作のモデルを適合させる。次に、モデル化したクロックを推定し、時間的に前向きに伝播させて、予測した軌道と共に、より迅速なＧＮＳＳ信号取得を可能にする。
【発明の詳細な説明】
【技術分野】
【０００１】
（政府認可書）
本明細書に記載する発明は、ＮＡＳＡ契約に基づく研究の遂行中に行われ、契約者が権利を保有するために選択した公法９６−５１７（３５ＵＳ２０２）の規定に支配される。
【０００２】
本発明は、衛星に基づく測位システム及び方法に関するものである。特に、本発明は、外部データ通信なしに、全地球航法衛星から直接取得した放送暦を用いて、（軌道状態及びクロック状態から成る）衛星暦を予測する方法及び装置に関するものである。
【背景技術】
【０００３】
全地球測位システム（ＧＰＳ：Global Positioning System）及びＧＬＯＮＡＳＳ（Global Navigation Satellite System）のような全地球航法衛星システム（ＧＮＳＳ：Global Navigation Satellite System）からのナビゲーション信号は、ＧＮＳＳ受信機によって用いられて、位置、ナビゲーション（航行）、及びタイミング（ＰＮＴ：Position, Navigation and Timing）サービス、及び多数の関連応用を提供する。
【０００４】
図1に示すように、ＧＮＳＳ受信機（１１０）は一般に、無線周波数（ＲＦ）部（１０２）及びデジタル部（１０５）を備えている。デジタル部（１０５）は更に、中央処理装置（ＣＰＵ：Central Processing unit）（１０７）及びメモリ（１０４）を備えている。一部のＧＮＳＳ受信機は、デジタル部を、携帯電話または個人用ナビゲーション装置のようなホスト装置と共有する。
【０００５】
ほとんど全てのＰＮＴ応用の基になる基本的な測位機能を実行するために、ＧＮＳＳ受信機はまず、十分な数のＧＮＳＳ衛星（１００）からＧＮＳＳナビゲーション信号（１０１）を取得し、次に、この信号（１０１）から必要なナビゲーション情報を抽出する。ＧＮＳＳ信号の取得はＲＦ部（１０２）内で実行され、ナビゲーション情報の抽出は、デジタル部（１０５）内で実行される。ＧＮＳＳナビゲーション信号中に組み込まれ、測位に必要なナビゲーション情報は、i) 受信機から送信衛星までの範囲（疑似範囲としても知られている、というのは、この範囲は衛星及び受信機クロックの不正確な知識に依存するからであり、この知識は、放送ナビゲーション・メッセージからのＧＮＳＳ衛星クロック状態、及び受信機クロック状態の最小二乗推定値を用いて後で調整する必要がある）、及びii) 送信衛星の位置、速度、及びクロック状態から成る。データi) 及びデータii) を用いて、受信機は最小二乗法または同様の推定方法を用いて、自機の位置を三角測量する。
【０００６】
受信機性能の重要な測定基準は初回測位時間（ＴＴＦＦ：Time to first fix）である。この測定基準は、種々の条件下で、受信機がＧＮＳＳ信号を取得して、自機の位置を測定するために要した時間を測定し、これらの条件は、コールドスタート、すなわち自機の位置または時刻、あるいはＧＮＳＳ衛星の現在状態の事前の正確な知識のない状態を含む。最良のユーザ体験のためには、ＴＴＦＦを最小にすることが望ましい。コールドスタート（あるいは、部分的に関係する先験的情報が利用可能な「ウォームスタート」）からのＴＴＦＦのかなりの部分が、受信機のＲＦ部（１０２）によってＧＮＳＳ信号の取得に費やされ、追加的な時間が、デジタル部（１０５）によって、ＧＮＳＳナビゲーション信号のナビゲーション・メッセージ部分から衛星位置及びクロック状態の情報を抽出することに費やされる。
【０００７】
信号取得の時間を要する側面は、所望信号について、範囲及びドップラー値のとり得る空間全体を走査する必要性に起因する。特定のＧＮＳＳ受信機がＧＮＳＳナビゲーション信号を取得する実際の速度は、その設計、位置、及び環境に強く依存する。
【０００８】
ＲＦ部（１０２）が信号を取得した後に、他の時間を要するステップにおいてナビゲーション情報を抽出する必要がある。例えば、ＧＰＳナビゲーション・メッセージは５０ビット／秒のバイナリデータストリームとして符号化され、衛星の軌道状態及びクロック状態についての必要な情報を含む正確な放送暦データは１５００ビットフレーム（ビット列）上に拡散されており、受信機（１１０）が全ての必要な放送暦情報を抽出するために最長３０秒を費やす必要がある。
【０００９】
ＧＮＳＳ受信機のＴＴＦＦを低減する一般的方法は、外部情報で受信機を「支援」することを中心に展開する。こうした情報は一般に、ＧＮＳＳ衛星の予期される軌道及びクロック状態の形をとる。大まかであるにせよ、衛星の位置及び／またはクロック状態の追加的な先験的情報により、受信機はＧＮＳＳ信号の探索範囲及びドップラー値の探索空間を狭めることができ、これらをより迅速に取得することができる。このように受信機に提供される、支援されるＧＮＳＳの軌道状態及びクロック状態（まとめて衛星暦と称する）が十分正確であれば、これらの状態は、等価な情報をナビゲーション信号から抽出する必要性を取り除き、ＴＴＦＦを更に低減する。
【００１０】
この種の支援を受信機（１１０）に提供するために、２つの基本的方法が採用されている。第１の方法は、隣接する受信機によって、あるいは受信機のネットワークによって同時に観測されているＧＮＳＳの軌道状態及びクロック状態を、ほぼリアルタイムで受信機（１１０）に提供することである。この方法は、監視している受信機間の、相当の高帯域を有する連続的な低遅延の通信チャネル、有用な情報を収集し、処理し、広めるための一連のサーバー、及びターゲット（目標）受信機または低軌道衛星を必要とし、ターゲット受信機は、例えば世界中のどこにでもある携帯電話またはカーナビゲーション装置とすることができる。採用されている第２の方法は、将来に向けて延長された期間にわたって有効なＧＮＳＳ衛星についての一組の軌道状態及びクロック状態をターゲット受信機に周期的に転送し、これにより必要な外部通信チャネルの必要性を低減することである。これらの期間延長された軌道状態及びクロック状態は、外部サーバーによって、放送暦及び／または、ＧＮＳＳ追跡受信機のネットワークから長期にわたり収集した範囲測定値を用いて計算される。後者の方法は、Ｔｏｐｅｘ衛星上のＧＰＳ受信機を支援するに当たり、ＪＰＬによって早くも１９９２年に採用され、この受信機ＴＴＦＦを低減している。この方法は、過去または現在の情報から、ＧＮＳＳの軌道状態及びクロック状態を将来に向けて予測または延長する能力を予見されている。（軌道状態及びクロック状態を含む）これらの暦は、「予測」暦または「延長」暦として知られている。後者の方法の変形例では、サーバーが初期状態、あるいは将来に及ぶ軌道状態及びクロック状態の「種」を周期的に計算し、これらの「種」のみを受信機に送信して、必要な通信帯域を低減する。そして受信機は、記憶している衛星の軌道及びクロックの動特性のモデルを用いて、「種」から延長暦を生成する。
【００１１】
支援方法が、ほぼリアルタイムの暦または関連情報（例えばＧＮＳＳナビゲーション・メッセージ全体）によるか、延長暦または関連情報（例えば延長暦モデルのための「種」）によるかに関わらず、これらの方法は全て、追跡受信機、サーバー、及びターゲット受信機と支援情報のサーバーとの間の何らかの通信チャネルから成るネットワークを含む相当量の外部インフラストラクチャ（基盤）を必要とする。これらの全てが、連続的かつ高信頼性で動作しなければならず、設備、用地、通信費、及び人間による監視及び保守のための相当な費用を生じさせる。
【００１２】
ＧＰＳ及びＧＬＯＮＡＳＳのようなＧＮＳＳ衛星の軌跡は、地球の重力及び太陽の放射束のような、完全には知られていない、あるいは完全に予測可能でない多数の物理的要因に影響される。従って、これらの衛星の軌跡を将来に向けて予測することは困難であり、そして常に、こうした予測の精度は時間と共に低下する。これらの衛星上の原子時計の状態を予測することは、更に困難である、というのは、これらの原子時計は、ランダムプロセス及び複雑な環境影響（例えば温度）を共に受けるからである。これに加えて、衛星の軌跡及びクロック状態は時として、これらの衛星のオペレータによって変更されることがあり、軌道状態及びクロック状態の予測を無用にする。これらの理由により、予測暦は相当規則的な間隔で周期的に更新しなければならない。これらの間隔は一般に、衛星クロック予測（予測するのが最も困難な量）の質によって支配され、受信機の位置精度の要求、及び軌道及びクロックの予測アルゴリズムの質に応じて、通常は、ほぼ毎日から週１回の間で変動する。
【００１３】
軌道状態及びクロック状態を予測するための２つの基本的方法が存在する。最も低精度の方法は、衛星状態（軌道及びクロック）の表現を所定時刻に取得して、この状態を先験的モデルと共に転送することである。第２の、より高精度の方法は、衛星状態またはクロック状態の時系列を使用し、一組のモデルパラメータをこの時系列に適合させることである。軌道については、これらのモデルパラメータは、３つの初期位置座標、３つの初期速度座標、太陽圧スケールのようないくつか物理的及び経験的パラメータ、及び定加速度を含む。一旦、これらのモデルパラメータを推定すると、モデルを数値的に将来に向けて伝播させて、より高忠実性で予測した軌道を生じさせることができる。同様のプロセスを続けて、クロック予測を最適化する。モデルパラメータを適合させる軌道状態及びクロック状態の時系列は、例えばネットワークに基づくＧＮＳＳ軌道測定から得られる正確なデータ、あるいは例えば放送暦から得られるより不正確なデータに基づくものとすることができる。それにもかかわらず、ＧＮＳＳの軌道及びクロックの予測における支配的な誤差原因は、伝播モデルにおける不確実性によるものである。データを適合させるソースの正確性は、放送暦でも正確な暦でも、相対的に重要でない。
【００１４】
衛星の動きを支配する物理的モデルは、慣性空間内に固定された座標系（地球中心慣性：ＥＣＩ：Earth-Centered-Inertial）により最良に記述される。従って、ＧＮＳＳ衛星の軌道伝播は、ＥＣＩ座標系内で行われる。しかし、地上ＧＮＳＳ受信機は、実際にはＧＮＳＳ放送暦によって提供される、地球と共に回転する座標系（地球中心−地球固定：ＥＣＥＦ：Earth-Centered-Earth-Fixed）内の位置情報を必要とする。従って、予測軌道は、受信機支援用に利用可能になる前に、ＥＣＩ座標からＥＣＥＦ座標に変換すべきであり、このことは、任意時点における慣性空間内の地球方位の知識を必要とする。しかし、地球方位は、非常に複雑で予測不可能なパターンに従う。地球方位は、観測され、超長基線電波干渉法（ＶＬＢＩ：Very Long Baseline Interferometry）の測地技術、及びＧＰＳを用いた精巧なデータ処理操作の後に報告される。３つの主要な地球方位モデルパラメータ（ＥＯＰ：Earth Orientation Model Parameter）、即ちＸ極運動、Ｙ極運動、及び日長（ＬＯＤ：Length of Day）は、ＮＡＳＡ−ＪＰＬ、国際地球回転観測及び基準座標系事業（ＩＥＲＳ：International Earth Rotation and Reference Systems Service）またはアメリカ国家地球空間情報局（ＮＧＡ：National Geo-Spatial Agency）のようないくつかの機関によって周期的に報告され、地球方位を記述するために広く利用され、ＥＣＩ座標系とＥＣＥＦ座標系との間の変換を可能にする。地球方位を、これらの周期的な外部ソース（情報源）に依存することは、現在、正確で長期的な受信機内の自律的な暦予測及び支援の、乗り越えられない障害の１つである。
【先行技術文献】
【特許文献】
【００１５】
【特許文献１】米国特許第５９６３１６７号明細書
【特許文献２】米国特許第５８２８３３６号明細書
【特許文献３】米国特許第６２９５０２１号明細書
【非特許文献】
【００１６】
【非特許文献１】GPS ICD 200
【非特許文献２】IERS Conventions (2003). Dennis D. McCarthy and Gerard Petit. Verlag des Bundesamts fur Katograpgic und Geodasie, 2004
【非特許文献３】Bate et al, fundamental of Astrodynamics, Dover, 1971, Chapter 5.2
【非特許文献４】Numerical Recipes, Press et al, Cambridge University Press, 1989, Chapter 15.1-15.2
【非特許文献５】IPN Progress Report 42-160 JPL, 2004
【非特許文献６】Bierman, Factorization Methods for Discrete Sequential Estimation, Dover, 2006
【発明の概要】
【発明が解決しようとする課題】
【００１７】
本発明は、信号取得を支援しＴＴＦＦを低減するために受信機が、期間延長された軌道状態及びクロック状態を自律的に生成することを可能にすることによって、高コストを伴うあらゆる外部インフラストラクチャ及び外部操作の必要性を取り除く。
【課題を解決するための手段】
【００１８】
第１の態様によれば、衛星の軌道情報を地球中心−地球固定（ＥＣＥＦ）座標で自律的に予測して、衛星受信機を支援する方法が提供され、この方法は、i) 衛星の放送軌道暦を取得するステップと、ii）収束が得られるまで、衛星の軌道位置、速度、及び動的パラメータを、衛星の動特性モデルに基づいて反復的に推定し、地球方位パラメータ（ＥＯＰ）を、地球方位モデルに基づいて反復的に推定するステップと、iii）推定した軌道位置、速度及び動的パラメータを、衛星の動特性モデルを用いて時間的に前向きに伝播させ、これにより予測軌道を得るステップと、iv）地球方位モデルを用いて、予測軌道を地球中心−地球固定（ＥＣＥＦ）座標に変換するステップと、v）変換した予測軌道を記憶するステップとを含む。
【００１９】
第２の態様によれば、衛星のクロック状態を自律的に予測する方法が提供され、この方法は、衛星についての放送クロック記録を取得するステップと、この衛星の将来のクロック状態を予測モデルに基づいて、収束が得られるまで反復的に推定するステップと、推定したクロック状態を、衛星のクロック状態の自律的予測用に記憶するステップとを含む。
【００２０】
第３の態様によれば、全地球航法衛星システム（ＧＮＳＳ）受信機が、処理装置と、メモリ装置とを備え、このメモリ装置は、衛星の軌道予測用及び衛星のクロック状態予測用のソフトウエアコードを記憶する第１メモリ部分と、記録した衛星の放送軌道及び衛星のクロック状態を記憶する第２メモリ部分と、地球方位パラメータ（ＥＯＰ）値及び軌道の動的パラメータを記憶する第３メモリ部分と、第１メモリ部に記憶されているソフトウエアコードに従って計算して予測した衛星の軌道状態及び衛星のクロック状態を記憶する第４メモリ部分とを備えるように構成されている。
【００２１】
他の好適例は、本願の明細書、図面、及び特許請求の範囲中に示す。
【００２２】
従って、本発明によれば、受信機固有のＧＮＳＳ ＲＦ受信機能以外のあらゆる通信チャネル、外部インフラストラクチャ、及び操作の必要性が取り除かれる。
【００２３】
本発明による方法は、携帯電話をホストとするＧＮＳＳ受信機、個人用ナビゲーション装置、車両用ナビゲーション装置、及び地球軌道内のＧＮＳＳ受信機のような装置上の応用をもたらす。
【００２４】
本発明の第１の利点は、初回測位時間（ＴＴＦＦ）が低減されることにある。
【００２５】
本発明の他の利点は、受信機が放送暦のサンプルを収集する周波数及びタイミングによって表され、必要な電力及びメモリリソースを最小にし、予測精度を最大にする。
【００２６】
本発明の他の利点は、ＧＮＳＳ衛星の軌道が、地球中心−地球固定（ＥＣＥＦ）の基準フレーム内で、限られた期間の衛星のＥＣＥＦ軌道状態のアーカイブ（保存）記録に基づいて、将来に向けて何日間にもわたって正確かつ効率的に予測されることにあり、このアーカイブ記録は、ＧＮＳＳ放送暦または精密な軌道決定解法から得ることができる。
【００２７】
更に、Ｘ極運動、Ｙ極運動、及び日長（ＬＯＤ）という地球方位パラメータ（ＥＯＰ）が、限られた期間の衛星のＥＣＥＦ軌道状態のアーカイブ記録のみに基づいて、正確かつ効率的に予測される。
【００２８】
これに加えて、Ｘ極運動、Ｙ極運動、日長（ＬＯＤ）、及びＵＴ１−ＵＴＣ（ＵＴ１マイナスＵＴＣ、原子標準によって規定されるＵＴＣタイムスケール（時間尺度）と、地球の回転によって規定されるＵＴ１タイムスケールとの差）という地球方位パラメータの値が、これらのパラメータの履歴的に観測された記録から、将来に向けて数年間にわたって予測される。
【００２９】
更に、ＧＮＳＳ衛星のクロック状態を、限られた期間の衛星のクロック状態のアーカイブ記録に基づいて、将来に向けて何日間にもわたって正確に予測することができ、このアーカイブ記録は、ＧＮＳＳ放送暦または正確なクロック決定解法から得ることができる。
【図面の簡単な説明】
【００３０】
【図１】一般的なＧＮＳＳ受信機のブロック図である。
【図２】本発明による、自律的に予測した軌道及びクロックの支援用の受信機要素のブロック図である。
【図３】本発明による軌道伝播プロセスのフローチャートである。
【図４】本発明によるクロック伝播プロセスのフローチャートである。
【発明を実施するための形態】
【００３１】
この開示によれば、工場でインストールされる情報以外の外部提供情報の必要性なしに、受信機固有のリソースのみを用いて、予測したＧＮＳＳの軌道状態及びクロック状態によって、ＧＮＳＳ受信機に信号取得支援を自律的に提供する装置及び方法を説明する。
【００３２】
図２に示すように、本発明による自立的に支援されるＧＮＳＳ受信機（２１０）は、ＧＮＳＳ信号を取得する無線周波数（ＲＦ）部（２０２）と、ＣＰＵ（２０７）及び不揮発性メモリ（２０４）を備えたデジタル部（２０５）とを備えている。デジタル部（２０５）は、個人用ナビゲーション装置（例えば、携帯電話またはカーナビゲーション装置）のようなホスト装置と共有することができる。ＣＰＵ（２０７）は、本発明の適用性を制限することなしに、多数の方法で実現することができる。これらは、マイクロプロセッサ、フールドプログラマブル・ゲートアレイ（ＦＰＧＡ：Field Programmable Gate Array）、または特定用途向け集積回路（ＡＳＩＣ：Application Specific Integrated Circuit）を含むが、これらに限定されない。
【００３３】
不揮発性メモリ（２０４）は、次の部分を備えるように構成されている：
１）図３及び図4を参照して後に説明する方法及びルーチンによる軌道及びクロック予測の部分（２０３）。
２）直近に記録した放送暦及びクロック状態を含む部分（２０６）。
３）地球方位値及び軌道の動的パラメータのデータベース（２０８）。
４）図３及び図4を参照して後に説明する方法を適用した結果として得られる予測した軌道状態及びクロック状態を含む部分（２０９）。
【００３４】
図３は、衛星位置の自律的予測の実施例のステップを示すフローチャートである。
【００３５】
ステップ（Ｓ１）に示すように、衛星毎に、図２に示すＧＮＳＳ放送ナビゲーション・メッセージ（１０１）から放送暦を得る。
【００３６】
電力制限または他の動作上の制限により、ＧＮＳＳ受信機が放送暦データを連続的に観測及び記録できないことがある。こうした問題を克服するために、ＧＮＳＳ受信機の動作制限に整合し、放送データの十分な量及び配分を可能にして軌道とクロックの正確な予測を促進するサンプリング方式を考案することができる。
【００３７】
放送暦記録の有用性を最大にするために、これらの記録の有効な２時間間隔（ＧＰＳの場合、他のＧＮＳＳシステムについては、有効な間隔は異なり得る）を利用して、正確性の損失なしに、各放送暦記録を１５分間ずつの、８つの個別記録に拡張することができる。
【００３８】
また、地上波受信機は常時、ＧＮＳＳ衛星群全体のごく一部（任意時点において、その地理的位置から常時見えるＧＰＳ衛星のおよそ３０％）しか見ておらず、ＧＰＳ衛星の地上航跡は毎日繰り返される傾向がある。本発明の実施例によれば、放送暦を１日当たり少なくとも３回サンプリングして、最適なＧＰＳ衛星群のサンプリングを保証することができる。その代わりに、衛星群全体をサンプリングするまで、２４時間の周期性を回避するサンプリング周波数を用いることができる。こうした代案実施例は、２４時間のサンプリング周期の場合に生じるＥＯＰの可観測性の弱まりを克服する。こうした実施例によれば、可能なサンプリング間隔は、７時間間隔または１０時間間隔とすることができる。７時間間隔は、ＧＰＳ衛星群全体をおよそ１日でサンプリングし、７日間は２４時間の周期性がないが、３日後にほぼ２４時間の周期性を有する。１０時間間隔は、ＧＰＳ衛星群全体をおよそ１．５日でサンプリングし、５日間は２４時間の周期性がない。こうした方法はソフトウエアによって実現することができ、不揮発性メモリに記憶することができる。
【００３９】
予想外の受信機操作を機に得られるあらゆる追加的放送記録を用いて、解決法を強化することができることも、当業者にとって明らかである。
【００４０】
図２に示すように、ステップ（Ｓ１）において観測した全ての放送暦及びクロック記録は、受信機（２１０）によって部分（２０６）に保存される。特定日数（例えば１４日間）よりも古い記録は廃棄することができる。
【００４１】
図３のステップ（Ｓ２）に示すように、次に、取得した放送暦を、現在技術において周知でありGPS ICD 200（非特許文献１）に記載された方法により、ケプラー座標からＥＣＥＦ（地球中心−地球固定）デカルト座標に変換し、不揮発性メモリの部分（２０６）に保存する。GPS ICD 200は、その全文を参考文献として本明細書に含める。
【００４２】
ステップ（Ｓ３）は、放送暦の期間の中頃の一時期についての、履歴的記録からの先験的ＥＯＰ（地球方位モデルパラメータ）値の生成を示す。ＥＯＰ値は、ＵＴ１−ＵＴＣ、ＬＯＤ、Ｘ極運動、Ｘ極運動速度、Ｙ極運動、及びＹ極運動速度を含む。
【００４３】
限定的でない一例として、ＥＯＰ値を予測するためのモデルは、国際地球回転観測及び基準座標系事業（ＩＥＲＳ）によって提供されるもののような、測定したＥＯＰ値の履歴的時系列へのパラメータ適合に基づくことができる。パラメータ毎の予測モデルは、履歴的時間記録（ほぼ５０年にわたる）のパワースペクトル分析を実行し、支配的な周波数を長期的な二次多項式と共に選択することによって生成することができる。次に、こうした予測モデルを用いて、将来に向けた履歴的時系列を効率的に推定することができる。このモデルは、図２の受信機の部分（２０３）に記憶されているソフトウエアとして実現することができる。
【００４４】
ステップ（Ｓ３）の予測ＥＯＰ値は後に、これら及び他のパラメータの最小二乗推定法において先験的値として用いることができる。高精度のＥＯＰ予測は、最小二乗推定法が正しい値に収束し、データ処理時間を最小にして電力を節約することを保証する。完全な自律性を保証するために、ＥＯＰ予測モデルは、一般に数年間である受信機の寿命期間中に有効であるべきである。しかし、このモデルは定期的更新により改良することができる。
【００４５】
図３に示す次のステップ（Ｓ４）は、衛星毎に、a) 初期軌道位置及び速度状態x(ｔ0)、x’(ｔ0)、及びb) 初期軌道の動的パラメータｐdを生成するステップである。
【００４６】
この開示で後に示すように、こうしたＥＯＰ（Ｓ３）、軌道状態（Ｓ４）、及び動的パラメータ（Ｓ４）は反復的に推定する。次の段落は、本発明の実施例による、こうした反復推定を実行する方法を説明する。軌道の動的パラメータｐdは、受信機（２１０）の不揮発性メモリ（２０４）の部分（２０８）に記憶されているデータベースから初期化する。
【００４７】
特に、こうした反復推定は、ＥＣＥＦ衛星状態への最小二乗適合を用いて実行することができる。こうした推定法は、受信機（２１０）の部分（２０３）に記憶されるソフトウエアによる方法とすることができる。簡単のため、以下の推定の説明は単一の衛星を参照し、この衛星は、例えば、その一意的な疑似ランダムノイズ（不規則雑音）番号（ＰＲＮ：Pseudorandom Noise Number）によって識別することができる。マトリクス（行列）代数及び記数法を用いる。この推定法は、共にＥＣＩ座標で表した、衛星のモデル化した軌道と、この衛星の放送暦との間の距離の二乗を、衛星の初期状態、動特性パラメータ、及び３つのＥＯＰを調整することによって最小化することを目標とする。
【００４８】
ｔjを、特定の衛星についての、放送軌道暦を取得した時点とする。索引ｊは、暦の時点を数え、ｊ＝1, 2,...Ｎである。ｒejは、時刻ｔjに放送暦から得られた衛星のＥＣＥＦ位置座標を保持する。ｒ’ejは、放送暦から得られた衛星のＥＣＥＦ速度座標を保持する。ｓjは、各時刻ｔjにおける放送暦位置の推定誤差とする。本発明の実施例によれば、全てのｊについてｓj＝１メートルである。ｓ’jは、各時刻ｔjにおける放送暦速度の推定誤差とする。例えば、全てのｊについてｓ’j＝0.001メートル／秒である。Ｔ＝Ｔ(ｔj,ＥＯＰ)は、ＥＣＥＦからＥＣＩへの回転行列とする。Ｔは、時刻の関数、かつ地球方位パラメータ、ＥＯＰの関数である。換言すれば、Ｔは放送暦の各時期における地球方位を表す。
【００４９】
Ｔを計算する好適な方法は、a) 地球方位を、時刻の関数、かつＥＯＰ（ＵＴ１−ＵＴＣ、ＬＯＤ、Ｘ極運動、Ｙ極運動、及びこれらの速度）の関数として計算し、b) 地球方位のＥＯＰに対する感度を計算する手法を用いる方法である。例えば、手法a)及びb)は、地球の歳差運動、章動、極運動、及び地球の回転のモデル化についての２００３年ＩＥＲＳ議定書に従うことができる（IERS Conventions (2003). Dennis D. McCarthy and Gerard Petit. Verlag des Bundesamts fur Katograpgic und Geodasie, 2004（非特許文献２）参照、その全文を参考文献として本明細書に含める）。これらの方法は、受信機（２１０）の部分（２０３）に記憶されているソフトウエェアルーチンとすることができる。
【００５０】
従って、図３のステップ（Ｓ５）に示すように、現在のＥＯＰ値を用い、現在のＥＯＰ値に対して、ＥＣＩ変換した放送軌道暦の偏導関数（感度）を計算することによって、放送軌道暦をＥＣＥＦからＥＣＩに変換することができる。計算した感度はその後、受信機（２１０）内に記憶する。
【００５１】
反復推定法に戻って説明する。ｘ及びｘ’を三次元ベクトルとし、これらのベクトルはそれぞれ、第１時点ｔ1における衛星のＥＣＩ位置及び速度を保持する。最初に、ｘはＴｒe1によって得ることができ、ｘ’はＴｒ’e1＋Ｔ’ｒe1によって得ることができる。本発明による方法の一実施例では、衛星暦からの速度は使用せず、この場合にｘ’は、Ｔｒe1、Ｔｒe2、ｔ1、ｔ2を用いていわゆるガウス問題を解くことによって得ることができる（Bate et al, fundamental of Astrodynamics, Dover, 1971, Chapter 5.2の227〜231ページ参照、 その全文を参考文献として本明細書に含める）。後に、ｘ及びｘ’を推定し、これを反復する。
【００５２】
ｐdを、ｍ個の軌道動特性パラメータの集合とし、これらのパラメータの先験的な値は十分既知ではないが、ＧＮＳＳ衛星の軌道動特性に対して強い影響を有する。本発明の実施例の１つによれば、ｐdは次の７つのパラメータから成る：太陽圧モデルの尺度（一般に「太陽尺度」と称する）、太陽−地球−衛星の平面に直交する方向の定加速度（一般に「Ｙバイアス」称する）、ラジアル（放射）方向及びクロストラック（軌道交差）方向の、軌道周期当たり２回の加速度を表す４つの調和係数、及びトラック（軌道）内の定加速度。
【００５３】
ｐeopを、ｎ個のＥＯＰパラメータの集合とし、これらのパラメータの先験的な値は十分既知ではないが、ＥＣＥＦフレームからＥＣＩフレームへの回転に対して強い影響を有する。本発明の実施例の１つによれば、ｐeopは次の３つのパラメータを備える：ＬＯＤ、Ｘ極運動、及びＹ極運動。
【００５４】
ｒi＝ｒi(ｘ, ｘ’, ｐd)を、時刻ｔjにおける衛星のモデル化したＥＣＩ位置座標とする。ｒiは、初期状態ｘ及びｘ’の関数、かつ動的パラメータｐdの関数であり、時刻ｔ1から時刻ｔjまでの運動方程式の積分によって得られる。ｒiは、ｐeopに弱くしか依存せず、従って、この依存性は明示的に指定しない。同様に、ｒ’i＝ｒ’i(ｘ, ｘ’, ｐd）を、時刻ｔjにおける衛星のモデル化したＥＣＩ速度座標とする。この速度座標が放送暦から抽出されない際は、これらの速度座標はモデル化もされない。
【００５５】
従って、図３のステップ（Ｓ６）に示すように、衛星位置予測モデルは、ｘ、ｘ’、ｐd及びｐeopの現在値から開始して生成される。一旦、このことを行うと、図３のステップ（Ｓ７）に示すように、パラメータｘ、ｘ’、ｐd、及びｐeopは、ＥＣＩ座標で表されるこうした予測モデル及び放送暦に基づいて推定される。
【００５６】
上述した運動／予測モデルの方程式は、例えば、a) 初期状態（所定時刻における位置及び速度）及び特定のモデルパラメータ（例えば地球の重力定数、衛星の質量）の集合からの衛星軌跡の伝播を可能にし、b) 伝播する状態の、こうした初期状態及びモデルパラメータに対する感度を計算するＧＮＳＳ軌道モデルを用いて得ることができる。別の言い方をすれば、ＧＮＳＳ軌道モデルを提供することができ、このモデルは、運動方程式と変分方程式とを十分な忠実性で統合する。このモデルの忠実性は、結果的に予測される軌道の精度に直接影響を与える。本発明の実施例の１つによれば、ＧＮＳＳ軌道モデルは、４〜５次の可変ステップサイズのルンゲクッタ−フェールベルグ（Runge-Kutta-Fehlberg）法の常微分方程式積分器（例えば、Numerical Recipes, Press et al, Cambridge University Press, 1989, Chapter 15.1-15.2の550〜560ページ参照、その全文を参考文献として本明細書に含める）、及び次の軌道動特性のモデルを採用する：地球、太陽、及び月の相対位置についてのＪＰＬ惑星暦を用いた地球、太陽、及び月についての質点引力；次数１２までの地球重力場；Bar-Sever and Kuangの経験的モデル [IPN Progress Report 42-160 JPL, 2004, http://ipnpr.jpl.nasa.gov/progress_report/42-160/title.htm参照、その全文を参考文献として本明細書に含める]のような太陽圧モデル、このモデルはＧＰＳ衛星の軌道予測を改善することが示されている；米国特許第５９６３１６７号明細書（特許文献１）及び第５８２８３３６号明細書（共に、その全文を参考文献として本明細書に含める）に示されているモデルのような地球と月の円錐形の影のモデル；米国特許第６２９５０２１号明細書（特許文献３、その全文を参考文献として本明細書に含める）に示されているＧＰＳ衛星用のモデルのような衛星姿勢モデル；ＲＦスラストモデル；及び相対性モデル。ＧＮＳＳ軌道モデルは、受信機（２１０）の部分（２０３）に記憶されているソフトウェアルーチンとして実現することができる。
【００５７】
反復推定法に戻って説明し。表記目的で、ｄを、ｘ、ｘ’、ｐdから成る複合的な(６＋ｍ)次元ベクトルとする。即ち、ｄt=[ｘt；ｘ’t；ｐd]であり、上付き文字tは移項演算を示す。
は動特性の勾配演算子であり、次式のように定義される：
【数１】
【００５８】
そして、
は、ｘ、ｘ’、及びｐdに対するｒiの感度の３×(６＋ｍ)行列であり、運動（以上で説明したモデル参照）の変分方程式を時刻ｔ1から時刻ｔjまで積分することによって計算される。必要であれば、
を同様に定義し、これは、運動の変分方程式を積分することによって計算される。
をＥＯＰ勾配演算子とし、次式のように定義する。
【数２】
ここに、
は、上述したＴを得る方法によって計算することができる。
【００５９】
ここで、最小二乗問題を次のように定式化する。最小二乗の意味で、ＥＣＩフレーム内で観測される位置Ｔｒejと、ＥＣＩフレーム内でモデル化した位置ｒj（ｊ＝１,...Ｎ）との間の距離を最小にするパラメータｄとｐeopの集合を見出す。この問題を数値的に解くために、この問題を、解となるパラメータｄ及びｐeopの先験的な値付近で線形化し、結果的な線形最小二乗問題を解いて、これらのパラメータについての調整値Δd及びΔｐeopを見出す。次に、この問題を、解となるパラメータの調整値付近で再線形化し、等々によって、解が収束するまで、この最小二乗の解法を反復する。ここで、線形最小二乗問題は次式のようになる：
【数３】
が
【数４】
を近似する。
【００６０】
未知数を左辺に分離することによって、次式のような最小二乗問題の古典的定式化が生じる：
【数５】
が、Ｔ(ｐeop)ｒej−ｒj(ｄ, ｐeop)；ｊ＝１,...Ｎを近似し、関連するデータノイズｓjを伴う。この数値的な最小二乗問題は、平方根情報フィルタ（ＳＲＩＦ：Square Root Information Filter）因数分解法を用いて解くことができる。可能な方法は、Housholderの正規直交三角測量法、上三角行列の逆行列を求める方法、及びHousholderの正規直交変換を用いてＳＲＩＦ配列を結合する方法[例えば、BiermanによるFactorization Methods for Discrete Sequential Estimation, Dover, 2006（非特許文献６）参照、その全文を参考文献として本明細書に含める]を含む。これらのアルゴリズムは不揮発性メモリに記憶することができる。
【００６１】
図３のステップ（Ｓ８）に示すように、ステップＳ５〜Ｓ７を、a) 最終の反復における衛星軌道が、この衛星の放送暦に十分近くなるまで、あるいはb) 推定したパラメータｘ、ｘ’、ｐd及びｐeopの値が、前回の反復から明らかに変化しなくなるまで反復する。換言すれば、解となる「デルタ（Δ）」の値が所定閾値を下回れば、即ち、小さ過ぎて問題にならなければ、あるいは、モデルと観測値との間の距離が所定閾値未満であれば、最小二乗問題が収束したものと見なし、反復プロセスを停止する。これらの条件を満足せず、かつ反復回数が所定閾値を超えなければ、次の新たな仮のパラメータ値：
ｄ＝ｄ＋Δｄ 及び ｐeop＝ｐeop＋Δｐeop
で最小二乗問題を反復する。
【００６２】
本発明の一実施例によれば、この最小二乗パラメータ推定問題は、各軌道面内に１つの衛星について解かれ、ＧＰＳの場合は６つの衛星（６つの軌道面を有する）となる。他の実施例は、衛星の異なる選択についての問題を解くことができ、ＣＰＵ時間と解の精度とのトレードオフになる。最小二乗問題は、衛星毎に大きく異なる動的パラメータｐdの集合を生じさせることが予期される。しかし、ＥＯＰパラメータの推定値ｐeopも異なり得る。ＥＯＰパラメータの推定値ｐeopの精度を向上させるために、図３のステップ（Ｓ９）に示すように、上述したＳＲＩＦ因数分解法を用いて、各衛星からのＥＯＰパラメータの推定値どうしを、それぞれの共分散行列全体と共に組み合わせることができる。このことは、単一の強固なＥＯＰ値の推定を生じさせる。
【００６３】
一旦、不確定のＥＯＰ値ｐeopが強固に推定されると、全ての衛星についての最小二乗問題を形成して解くことができ、ここでは、初期状態パラメータｘ、ｘ’、及び動的パラメータｐdのみを推定し、ＥＯＰは、最終の組み合わせた推定値に固定して保持する。このことも図３のステップ（Ｓ１０）に示す。このプロセスの結果は、全ての衛星についての、初期状態（ｘ, ｘ’）と動的パラメータｐdとの集合となる。図３のステップ（Ｓ１１）〜（Ｓ１３）に示すように、これらの要素は、（以上で説明したモデルと共に）これらの衛星の軌道の将来に向かう伝播を可能にする。これらの値は受信機の不揮発性メモリ（１０４）に記憶する。こうした値は、信号取得及び位置決めを支援するために用いることのできる予測軌道を表す。
【００６４】
ここで図４を参照し、図４は、本発明の実施例による衛星クロック状態を予測する方法のフローチャートを示す。
【００６５】
ステップ（Ｔ１）に示すように、衛星毎の最新の放送クロック記録をメモリからロードする。好適な実施例では、これらの記録は少なくとも1日にわたる。ステップ（Ｔ２）では、衛星毎に、クロックモデルを放送暦値に適合させることによってクロックのモデルパラメータを決定する。
【００６６】
ＧＮＳＳクロックの仮のモデルは、（時間の二次式）＋（１回転に１回の調波）であり、衛星当たり５つのモデルパラメータを生じさせる（うち３つは二次多項式用、２つは１回転に１回の調波のサイン及びコサイン成分用）。これらのモデルパラメータは、放送暦から抽出した数日分の観測クロック値へのモデルの最小二乗適合を用いて推定することができる。図４の反復ステップ（Ｔ３）を参照。図４のステップ（Ｔ４）に示すように、一旦、衛星毎に個別にモデルパラメータを推定すると、衛星毎のモデルを時間的に前向きに、無限に伝播させて、予測クロック値を生成する。衛星毎に推定したモデルパラメータは、受信機（２１０）のメモリ（２０４）に記憶する。これらをいつでも呼び出して、クロック予測を生成することができる。本発明の一実施例によれば、４日分の放送クロック値を、モデルへの最小二乗適合に用いることができる。ＧＮＳＳのオペレータによって、モデルと不整合なクロックジャンプ（クロックの急変）が時として導入され得るので、最小二乗適合誤差（二乗平均平方根：ＲＭＳ：Root-Mean-Square）を所定閾値に対してチェックし、この閾値は仮に１メートルとする。ＲＭＳ適合がこの閾値を超えれば、適合間隔のうち最も適合が見られない１２時間の部分を試しに除去して、モデル不整合のない間隔を生成することができる。このプロセスは、適合のＲＭＳ値が閾値を下回るまで、あるいは適合間隔が１２時間未満である場合に反復することができる。２日間未満の適合間隔については、線形クロックモデルを適用すべきである（衛星当たり２つのパラメータ）。通常のように、こうしたアルゴリズムモデルは受信機（２１０）の部分（２０３）に記憶することができる。
【００６７】
予測プロセスの最後に、伝播させた軌道及び予測したクロックモデルを利用して、将来のＧＮＳＳ信号の取得に当たり受信機を支援することができる。新たな暦データ記録が利用可能になる毎に、受信機の電力管理の考慮及び他の動作上のトレードオフに応じて、この予測プロセスを反復することができる。
【００６８】
上述したように、また図２を参照して述べたように、本発明の実施例の１つによれば、本発明による受信機は不揮発性メモリ部分（２０４）を備え、この不揮発性メモリ部分は次のものを含む：
１）以上で説明した種々のソフトウェアルーチンのオブジェクトコード、例えばi) ＥＯＰ予測モデル；ii) 回転行列を生成する方法；iii) 運動方程式を統合したＧＮＳＳ軌道モデル；iv) ＳＲＩＦアルゴリズム；v) 軌道初期状態、動的パラメータ、及びＥＯＰの反復推定方法；vi) クロック予測方法；及びvii) ＧＮＳＳ放送暦データの最適なサンプリング方法。
２）最新の記録放送暦。
３）初期化目的用の、軌道の動的パラメータｄのデータベース。
４）上記ＥＯＰ予測モデルを用いて生成した先験的ＥＯＰ値。
５）上記１）で述べた方法を適用した結果として得られた、予測した軌道状態及びクロック状態。
【００６９】
要約すれば、本明細書におけるいくつかの実施例によれば、全地球航法衛星システム（ＧＮＳＳ）の軌道状態及びクロック状態の自律的な受信機内予測の方法及び装置を説明している。ＧＮＳＳ放送メッセージのみを使用し、周期的な外部通信情報を必要としない。地球方位情報をＧＮＳＳ放送暦から抽出する。地球方位パラメータの正確な推定により、慣性基準フレーム内で、最良に適合するＧＮＳＳ軌道を、時間的に前向きに伝播させることができる。次に、推定した地球方位パラメータを用いて、予測した軌道を地球中心−地球固定（ＥＣＥＦ）座標に変換して、信号の取得に当たりＧＮＳＳ受信機を支援するために使用する。ＧＮＳＳ衛星のクロック状態も放送暦から抽出して、このデータに、パラメータ化したクロック動作のモデルを適合させる。次に、モデル化したクロックを推定し、時間的に前向きに伝播させて、予測した軌道と共に、より迅速なＧＮＳＳ信号取得を可能にする。
【００７０】
従って、これまで示してきたものは、ＧＮＳＳ暦の自律的な受信機内予測の方法及び装置である。これらの方法及び装置は、特定実施例及びその応用を用いて説明してきたが、本発明の範囲を逸脱することなしに、当業者が多数の変更及び変形を加え得ることは明らかである。従って、特許請求の範囲内で、本明細書で具体的に説明した以外の方法で本発明を実施することができることは明らかである。
【技術分野】
【０００１】
（政府認可書）
本明細書に記載する発明は、ＮＡＳＡ契約に基づく研究の遂行中に行われ、契約者が権利を保有するために選択した公法９６−５１７（３５ＵＳ２０２）の規定に支配される。
【０００２】
本発明は、衛星に基づく測位システム及び方法に関するものである。特に、本発明は、外部データ通信なしに、全地球航法衛星から直接取得した放送暦を用いて、（軌道状態及びクロック状態から成る）衛星暦を予測する方法及び装置に関するものである。
【背景技術】
【０００３】
全地球測位システム（ＧＰＳ：Global Positioning System）及びＧＬＯＮＡＳＳ（Global Navigation Satellite System）のような全地球航法衛星システム（ＧＮＳＳ：Global Navigation Satellite System）からのナビゲーション信号は、ＧＮＳＳ受信機によって用いられて、位置、ナビゲーション（航行）、及びタイミング（ＰＮＴ：Position, Navigation and Timing）サービス、及び多数の関連応用を提供する。
【０００４】
図1に示すように、ＧＮＳＳ受信機（１１０）は一般に、無線周波数（ＲＦ）部（１０２）及びデジタル部（１０５）を備えている。デジタル部（１０５）は更に、中央処理装置（ＣＰＵ：Central Processing unit）（１０７）及びメモリ（１０４）を備えている。一部のＧＮＳＳ受信機は、デジタル部を、携帯電話または個人用ナビゲーション装置のようなホスト装置と共有する。
【０００５】
ほとんど全てのＰＮＴ応用の基になる基本的な測位機能を実行するために、ＧＮＳＳ受信機はまず、十分な数のＧＮＳＳ衛星（１００）からＧＮＳＳナビゲーション信号（１０１）を取得し、次に、この信号（１０１）から必要なナビゲーション情報を抽出する。ＧＮＳＳ信号の取得はＲＦ部（１０２）内で実行され、ナビゲーション情報の抽出は、デジタル部（１０５）内で実行される。ＧＮＳＳナビゲーション信号中に組み込まれ、測位に必要なナビゲーション情報は、i) 受信機から送信衛星までの範囲（疑似範囲としても知られている、というのは、この範囲は衛星及び受信機クロックの不正確な知識に依存するからであり、この知識は、放送ナビゲーション・メッセージからのＧＮＳＳ衛星クロック状態、及び受信機クロック状態の最小二乗推定値を用いて後で調整する必要がある）、及びii) 送信衛星の位置、速度、及びクロック状態から成る。データi) 及びデータii) を用いて、受信機は最小二乗法または同様の推定方法を用いて、自機の位置を三角測量する。
【０００６】
受信機性能の重要な測定基準は初回測位時間（ＴＴＦＦ：Time to first fix）である。この測定基準は、種々の条件下で、受信機がＧＮＳＳ信号を取得して、自機の位置を測定するために要した時間を測定し、これらの条件は、コールドスタート、すなわち自機の位置または時刻、あるいはＧＮＳＳ衛星の現在状態の事前の正確な知識のない状態を含む。最良のユーザ体験のためには、ＴＴＦＦを最小にすることが望ましい。コールドスタート（あるいは、部分的に関係する先験的情報が利用可能な「ウォームスタート」）からのＴＴＦＦのかなりの部分が、受信機のＲＦ部（１０２）によってＧＮＳＳ信号の取得に費やされ、追加的な時間が、デジタル部（１０５）によって、ＧＮＳＳナビゲーション信号のナビゲーション・メッセージ部分から衛星位置及びクロック状態の情報を抽出することに費やされる。
【０００７】
信号取得の時間を要する側面は、所望信号について、範囲及びドップラー値のとり得る空間全体を走査する必要性に起因する。特定のＧＮＳＳ受信機がＧＮＳＳナビゲーション信号を取得する実際の速度は、その設計、位置、及び環境に強く依存する。
【０００８】
ＲＦ部（１０２）が信号を取得した後に、他の時間を要するステップにおいてナビゲーション情報を抽出する必要がある。例えば、ＧＰＳナビゲーション・メッセージは５０ビット／秒のバイナリデータストリームとして符号化され、衛星の軌道状態及びクロック状態についての必要な情報を含む正確な放送暦データは１５００ビットフレーム（ビット列）上に拡散されており、受信機（１１０）が全ての必要な放送暦情報を抽出するために最長３０秒を費やす必要がある。
【０００９】
ＧＮＳＳ受信機のＴＴＦＦを低減する一般的方法は、外部情報で受信機を「支援」することを中心に展開する。こうした情報は一般に、ＧＮＳＳ衛星の予期される軌道及びクロック状態の形をとる。大まかであるにせよ、衛星の位置及び／またはクロック状態の追加的な先験的情報により、受信機はＧＮＳＳ信号の探索範囲及びドップラー値の探索空間を狭めることができ、これらをより迅速に取得することができる。このように受信機に提供される、支援されるＧＮＳＳの軌道状態及びクロック状態（まとめて衛星暦と称する）が十分正確であれば、これらの状態は、等価な情報をナビゲーション信号から抽出する必要性を取り除き、ＴＴＦＦを更に低減する。
【００１０】
この種の支援を受信機（１１０）に提供するために、２つの基本的方法が採用されている。第１の方法は、隣接する受信機によって、あるいは受信機のネットワークによって同時に観測されているＧＮＳＳの軌道状態及びクロック状態を、ほぼリアルタイムで受信機（１１０）に提供することである。この方法は、監視している受信機間の、相当の高帯域を有する連続的な低遅延の通信チャネル、有用な情報を収集し、処理し、広めるための一連のサーバー、及びターゲット（目標）受信機または低軌道衛星を必要とし、ターゲット受信機は、例えば世界中のどこにでもある携帯電話またはカーナビゲーション装置とすることができる。採用されている第２の方法は、将来に向けて延長された期間にわたって有効なＧＮＳＳ衛星についての一組の軌道状態及びクロック状態をターゲット受信機に周期的に転送し、これにより必要な外部通信チャネルの必要性を低減することである。これらの期間延長された軌道状態及びクロック状態は、外部サーバーによって、放送暦及び／または、ＧＮＳＳ追跡受信機のネットワークから長期にわたり収集した範囲測定値を用いて計算される。後者の方法は、Ｔｏｐｅｘ衛星上のＧＰＳ受信機を支援するに当たり、ＪＰＬによって早くも１９９２年に採用され、この受信機ＴＴＦＦを低減している。この方法は、過去または現在の情報から、ＧＮＳＳの軌道状態及びクロック状態を将来に向けて予測または延長する能力を予見されている。（軌道状態及びクロック状態を含む）これらの暦は、「予測」暦または「延長」暦として知られている。後者の方法の変形例では、サーバーが初期状態、あるいは将来に及ぶ軌道状態及びクロック状態の「種」を周期的に計算し、これらの「種」のみを受信機に送信して、必要な通信帯域を低減する。そして受信機は、記憶している衛星の軌道及びクロックの動特性のモデルを用いて、「種」から延長暦を生成する。
【００１１】
支援方法が、ほぼリアルタイムの暦または関連情報（例えばＧＮＳＳナビゲーション・メッセージ全体）によるか、延長暦または関連情報（例えば延長暦モデルのための「種」）によるかに関わらず、これらの方法は全て、追跡受信機、サーバー、及びターゲット受信機と支援情報のサーバーとの間の何らかの通信チャネルから成るネットワークを含む相当量の外部インフラストラクチャ（基盤）を必要とする。これらの全てが、連続的かつ高信頼性で動作しなければならず、設備、用地、通信費、及び人間による監視及び保守のための相当な費用を生じさせる。
【００１２】
ＧＰＳ及びＧＬＯＮＡＳＳのようなＧＮＳＳ衛星の軌跡は、地球の重力及び太陽の放射束のような、完全には知られていない、あるいは完全に予測可能でない多数の物理的要因に影響される。従って、これらの衛星の軌跡を将来に向けて予測することは困難であり、そして常に、こうした予測の精度は時間と共に低下する。これらの衛星上の原子時計の状態を予測することは、更に困難である、というのは、これらの原子時計は、ランダムプロセス及び複雑な環境影響（例えば温度）を共に受けるからである。これに加えて、衛星の軌跡及びクロック状態は時として、これらの衛星のオペレータによって変更されることがあり、軌道状態及びクロック状態の予測を無用にする。これらの理由により、予測暦は相当規則的な間隔で周期的に更新しなければならない。これらの間隔は一般に、衛星クロック予測（予測するのが最も困難な量）の質によって支配され、受信機の位置精度の要求、及び軌道及びクロックの予測アルゴリズムの質に応じて、通常は、ほぼ毎日から週１回の間で変動する。
【００１３】
軌道状態及びクロック状態を予測するための２つの基本的方法が存在する。最も低精度の方法は、衛星状態（軌道及びクロック）の表現を所定時刻に取得して、この状態を先験的モデルと共に転送することである。第２の、より高精度の方法は、衛星状態またはクロック状態の時系列を使用し、一組のモデルパラメータをこの時系列に適合させることである。軌道については、これらのモデルパラメータは、３つの初期位置座標、３つの初期速度座標、太陽圧スケールのようないくつか物理的及び経験的パラメータ、及び定加速度を含む。一旦、これらのモデルパラメータを推定すると、モデルを数値的に将来に向けて伝播させて、より高忠実性で予測した軌道を生じさせることができる。同様のプロセスを続けて、クロック予測を最適化する。モデルパラメータを適合させる軌道状態及びクロック状態の時系列は、例えばネットワークに基づくＧＮＳＳ軌道測定から得られる正確なデータ、あるいは例えば放送暦から得られるより不正確なデータに基づくものとすることができる。それにもかかわらず、ＧＮＳＳの軌道及びクロックの予測における支配的な誤差原因は、伝播モデルにおける不確実性によるものである。データを適合させるソースの正確性は、放送暦でも正確な暦でも、相対的に重要でない。
【００１４】
衛星の動きを支配する物理的モデルは、慣性空間内に固定された座標系（地球中心慣性：ＥＣＩ：Earth-Centered-Inertial）により最良に記述される。従って、ＧＮＳＳ衛星の軌道伝播は、ＥＣＩ座標系内で行われる。しかし、地上ＧＮＳＳ受信機は、実際にはＧＮＳＳ放送暦によって提供される、地球と共に回転する座標系（地球中心−地球固定：ＥＣＥＦ：Earth-Centered-Earth-Fixed）内の位置情報を必要とする。従って、予測軌道は、受信機支援用に利用可能になる前に、ＥＣＩ座標からＥＣＥＦ座標に変換すべきであり、このことは、任意時点における慣性空間内の地球方位の知識を必要とする。しかし、地球方位は、非常に複雑で予測不可能なパターンに従う。地球方位は、観測され、超長基線電波干渉法（ＶＬＢＩ：Very Long Baseline Interferometry）の測地技術、及びＧＰＳを用いた精巧なデータ処理操作の後に報告される。３つの主要な地球方位モデルパラメータ（ＥＯＰ：Earth Orientation Model Parameter）、即ちＸ極運動、Ｙ極運動、及び日長（ＬＯＤ：Length of Day）は、ＮＡＳＡ−ＪＰＬ、国際地球回転観測及び基準座標系事業（ＩＥＲＳ：International Earth Rotation and Reference Systems Service）またはアメリカ国家地球空間情報局（ＮＧＡ：National Geo-Spatial Agency）のようないくつかの機関によって周期的に報告され、地球方位を記述するために広く利用され、ＥＣＩ座標系とＥＣＥＦ座標系との間の変換を可能にする。地球方位を、これらの周期的な外部ソース（情報源）に依存することは、現在、正確で長期的な受信機内の自律的な暦予測及び支援の、乗り越えられない障害の１つである。
【先行技術文献】
【特許文献】
【００１５】
【特許文献１】米国特許第５９６３１６７号明細書
【特許文献２】米国特許第５８２８３３６号明細書
【特許文献３】米国特許第６２９５０２１号明細書
【非特許文献】
【００１６】
【非特許文献１】GPS ICD 200
【非特許文献２】IERS Conventions (2003). Dennis D. McCarthy and Gerard Petit. Verlag des Bundesamts fur Katograpgic und Geodasie, 2004
【非特許文献３】Bate et al, fundamental of Astrodynamics, Dover, 1971, Chapter 5.2
【非特許文献４】Numerical Recipes, Press et al, Cambridge University Press, 1989, Chapter 15.1-15.2
【非特許文献５】IPN Progress Report 42-160 JPL, 2004
【非特許文献６】Bierman, Factorization Methods for Discrete Sequential Estimation, Dover, 2006
【発明の概要】
【発明が解決しようとする課題】
【００１７】
本発明は、信号取得を支援しＴＴＦＦを低減するために受信機が、期間延長された軌道状態及びクロック状態を自律的に生成することを可能にすることによって、高コストを伴うあらゆる外部インフラストラクチャ及び外部操作の必要性を取り除く。
【課題を解決するための手段】
【００１８】
第１の態様によれば、衛星の軌道情報を地球中心−地球固定（ＥＣＥＦ）座標で自律的に予測して、衛星受信機を支援する方法が提供され、この方法は、i) 衛星の放送軌道暦を取得するステップと、ii）収束が得られるまで、衛星の軌道位置、速度、及び動的パラメータを、衛星の動特性モデルに基づいて反復的に推定し、地球方位パラメータ（ＥＯＰ）を、地球方位モデルに基づいて反復的に推定するステップと、iii）推定した軌道位置、速度及び動的パラメータを、衛星の動特性モデルを用いて時間的に前向きに伝播させ、これにより予測軌道を得るステップと、iv）地球方位モデルを用いて、予測軌道を地球中心−地球固定（ＥＣＥＦ）座標に変換するステップと、v）変換した予測軌道を記憶するステップとを含む。
【００１９】
第２の態様によれば、衛星のクロック状態を自律的に予測する方法が提供され、この方法は、衛星についての放送クロック記録を取得するステップと、この衛星の将来のクロック状態を予測モデルに基づいて、収束が得られるまで反復的に推定するステップと、推定したクロック状態を、衛星のクロック状態の自律的予測用に記憶するステップとを含む。
【００２０】
第３の態様によれば、全地球航法衛星システム（ＧＮＳＳ）受信機が、処理装置と、メモリ装置とを備え、このメモリ装置は、衛星の軌道予測用及び衛星のクロック状態予測用のソフトウエアコードを記憶する第１メモリ部分と、記録した衛星の放送軌道及び衛星のクロック状態を記憶する第２メモリ部分と、地球方位パラメータ（ＥＯＰ）値及び軌道の動的パラメータを記憶する第３メモリ部分と、第１メモリ部に記憶されているソフトウエアコードに従って計算して予測した衛星の軌道状態及び衛星のクロック状態を記憶する第４メモリ部分とを備えるように構成されている。
【００２１】
他の好適例は、本願の明細書、図面、及び特許請求の範囲中に示す。
【００２２】
従って、本発明によれば、受信機固有のＧＮＳＳ ＲＦ受信機能以外のあらゆる通信チャネル、外部インフラストラクチャ、及び操作の必要性が取り除かれる。
【００２３】
本発明による方法は、携帯電話をホストとするＧＮＳＳ受信機、個人用ナビゲーション装置、車両用ナビゲーション装置、及び地球軌道内のＧＮＳＳ受信機のような装置上の応用をもたらす。
【００２４】
本発明の第１の利点は、初回測位時間（ＴＴＦＦ）が低減されることにある。
【００２５】
本発明の他の利点は、受信機が放送暦のサンプルを収集する周波数及びタイミングによって表され、必要な電力及びメモリリソースを最小にし、予測精度を最大にする。
【００２６】
本発明の他の利点は、ＧＮＳＳ衛星の軌道が、地球中心−地球固定（ＥＣＥＦ）の基準フレーム内で、限られた期間の衛星のＥＣＥＦ軌道状態のアーカイブ（保存）記録に基づいて、将来に向けて何日間にもわたって正確かつ効率的に予測されることにあり、このアーカイブ記録は、ＧＮＳＳ放送暦または精密な軌道決定解法から得ることができる。
【００２７】
更に、Ｘ極運動、Ｙ極運動、及び日長（ＬＯＤ）という地球方位パラメータ（ＥＯＰ）が、限られた期間の衛星のＥＣＥＦ軌道状態のアーカイブ記録のみに基づいて、正確かつ効率的に予測される。
【００２８】
これに加えて、Ｘ極運動、Ｙ極運動、日長（ＬＯＤ）、及びＵＴ１−ＵＴＣ（ＵＴ１マイナスＵＴＣ、原子標準によって規定されるＵＴＣタイムスケール（時間尺度）と、地球の回転によって規定されるＵＴ１タイムスケールとの差）という地球方位パラメータの値が、これらのパラメータの履歴的に観測された記録から、将来に向けて数年間にわたって予測される。
【００２９】
更に、ＧＮＳＳ衛星のクロック状態を、限られた期間の衛星のクロック状態のアーカイブ記録に基づいて、将来に向けて何日間にもわたって正確に予測することができ、このアーカイブ記録は、ＧＮＳＳ放送暦または正確なクロック決定解法から得ることができる。
【図面の簡単な説明】
【００３０】
【図１】一般的なＧＮＳＳ受信機のブロック図である。
【図２】本発明による、自律的に予測した軌道及びクロックの支援用の受信機要素のブロック図である。
【図３】本発明による軌道伝播プロセスのフローチャートである。
【図４】本発明によるクロック伝播プロセスのフローチャートである。
【発明を実施するための形態】
【００３１】
この開示によれば、工場でインストールされる情報以外の外部提供情報の必要性なしに、受信機固有のリソースのみを用いて、予測したＧＮＳＳの軌道状態及びクロック状態によって、ＧＮＳＳ受信機に信号取得支援を自律的に提供する装置及び方法を説明する。
【００３２】
図２に示すように、本発明による自立的に支援されるＧＮＳＳ受信機（２１０）は、ＧＮＳＳ信号を取得する無線周波数（ＲＦ）部（２０２）と、ＣＰＵ（２０７）及び不揮発性メモリ（２０４）を備えたデジタル部（２０５）とを備えている。デジタル部（２０５）は、個人用ナビゲーション装置（例えば、携帯電話またはカーナビゲーション装置）のようなホスト装置と共有することができる。ＣＰＵ（２０７）は、本発明の適用性を制限することなしに、多数の方法で実現することができる。これらは、マイクロプロセッサ、フールドプログラマブル・ゲートアレイ（ＦＰＧＡ：Field Programmable Gate Array）、または特定用途向け集積回路（ＡＳＩＣ：Application Specific Integrated Circuit）を含むが、これらに限定されない。
【００３３】
不揮発性メモリ（２０４）は、次の部分を備えるように構成されている：
１）図３及び図4を参照して後に説明する方法及びルーチンによる軌道及びクロック予測の部分（２０３）。
２）直近に記録した放送暦及びクロック状態を含む部分（２０６）。
３）地球方位値及び軌道の動的パラメータのデータベース（２０８）。
４）図３及び図4を参照して後に説明する方法を適用した結果として得られる予測した軌道状態及びクロック状態を含む部分（２０９）。
【００３４】
図３は、衛星位置の自律的予測の実施例のステップを示すフローチャートである。
【００３５】
ステップ（Ｓ１）に示すように、衛星毎に、図２に示すＧＮＳＳ放送ナビゲーション・メッセージ（１０１）から放送暦を得る。
【００３６】
電力制限または他の動作上の制限により、ＧＮＳＳ受信機が放送暦データを連続的に観測及び記録できないことがある。こうした問題を克服するために、ＧＮＳＳ受信機の動作制限に整合し、放送データの十分な量及び配分を可能にして軌道とクロックの正確な予測を促進するサンプリング方式を考案することができる。
【００３７】
放送暦記録の有用性を最大にするために、これらの記録の有効な２時間間隔（ＧＰＳの場合、他のＧＮＳＳシステムについては、有効な間隔は異なり得る）を利用して、正確性の損失なしに、各放送暦記録を１５分間ずつの、８つの個別記録に拡張することができる。
【００３８】
また、地上波受信機は常時、ＧＮＳＳ衛星群全体のごく一部（任意時点において、その地理的位置から常時見えるＧＰＳ衛星のおよそ３０％）しか見ておらず、ＧＰＳ衛星の地上航跡は毎日繰り返される傾向がある。本発明の実施例によれば、放送暦を１日当たり少なくとも３回サンプリングして、最適なＧＰＳ衛星群のサンプリングを保証することができる。その代わりに、衛星群全体をサンプリングするまで、２４時間の周期性を回避するサンプリング周波数を用いることができる。こうした代案実施例は、２４時間のサンプリング周期の場合に生じるＥＯＰの可観測性の弱まりを克服する。こうした実施例によれば、可能なサンプリング間隔は、７時間間隔または１０時間間隔とすることができる。７時間間隔は、ＧＰＳ衛星群全体をおよそ１日でサンプリングし、７日間は２４時間の周期性がないが、３日後にほぼ２４時間の周期性を有する。１０時間間隔は、ＧＰＳ衛星群全体をおよそ１．５日でサンプリングし、５日間は２４時間の周期性がない。こうした方法はソフトウエアによって実現することができ、不揮発性メモリに記憶することができる。
【００３９】
予想外の受信機操作を機に得られるあらゆる追加的放送記録を用いて、解決法を強化することができることも、当業者にとって明らかである。
【００４０】
図２に示すように、ステップ（Ｓ１）において観測した全ての放送暦及びクロック記録は、受信機（２１０）によって部分（２０６）に保存される。特定日数（例えば１４日間）よりも古い記録は廃棄することができる。
【００４１】
図３のステップ（Ｓ２）に示すように、次に、取得した放送暦を、現在技術において周知でありGPS ICD 200（非特許文献１）に記載された方法により、ケプラー座標からＥＣＥＦ（地球中心−地球固定）デカルト座標に変換し、不揮発性メモリの部分（２０６）に保存する。GPS ICD 200は、その全文を参考文献として本明細書に含める。
【００４２】
ステップ（Ｓ３）は、放送暦の期間の中頃の一時期についての、履歴的記録からの先験的ＥＯＰ（地球方位モデルパラメータ）値の生成を示す。ＥＯＰ値は、ＵＴ１−ＵＴＣ、ＬＯＤ、Ｘ極運動、Ｘ極運動速度、Ｙ極運動、及びＹ極運動速度を含む。
【００４３】
限定的でない一例として、ＥＯＰ値を予測するためのモデルは、国際地球回転観測及び基準座標系事業（ＩＥＲＳ）によって提供されるもののような、測定したＥＯＰ値の履歴的時系列へのパラメータ適合に基づくことができる。パラメータ毎の予測モデルは、履歴的時間記録（ほぼ５０年にわたる）のパワースペクトル分析を実行し、支配的な周波数を長期的な二次多項式と共に選択することによって生成することができる。次に、こうした予測モデルを用いて、将来に向けた履歴的時系列を効率的に推定することができる。このモデルは、図２の受信機の部分（２０３）に記憶されているソフトウエアとして実現することができる。
【００４４】
ステップ（Ｓ３）の予測ＥＯＰ値は後に、これら及び他のパラメータの最小二乗推定法において先験的値として用いることができる。高精度のＥＯＰ予測は、最小二乗推定法が正しい値に収束し、データ処理時間を最小にして電力を節約することを保証する。完全な自律性を保証するために、ＥＯＰ予測モデルは、一般に数年間である受信機の寿命期間中に有効であるべきである。しかし、このモデルは定期的更新により改良することができる。
【００４５】
図３に示す次のステップ（Ｓ４）は、衛星毎に、a) 初期軌道位置及び速度状態x(ｔ0)、x’(ｔ0)、及びb) 初期軌道の動的パラメータｐdを生成するステップである。
【００４６】
この開示で後に示すように、こうしたＥＯＰ（Ｓ３）、軌道状態（Ｓ４）、及び動的パラメータ（Ｓ４）は反復的に推定する。次の段落は、本発明の実施例による、こうした反復推定を実行する方法を説明する。軌道の動的パラメータｐdは、受信機（２１０）の不揮発性メモリ（２０４）の部分（２０８）に記憶されているデータベースから初期化する。
【００４７】
特に、こうした反復推定は、ＥＣＥＦ衛星状態への最小二乗適合を用いて実行することができる。こうした推定法は、受信機（２１０）の部分（２０３）に記憶されるソフトウエアによる方法とすることができる。簡単のため、以下の推定の説明は単一の衛星を参照し、この衛星は、例えば、その一意的な疑似ランダムノイズ（不規則雑音）番号（ＰＲＮ：Pseudorandom Noise Number）によって識別することができる。マトリクス（行列）代数及び記数法を用いる。この推定法は、共にＥＣＩ座標で表した、衛星のモデル化した軌道と、この衛星の放送暦との間の距離の二乗を、衛星の初期状態、動特性パラメータ、及び３つのＥＯＰを調整することによって最小化することを目標とする。
【００４８】
ｔjを、特定の衛星についての、放送軌道暦を取得した時点とする。索引ｊは、暦の時点を数え、ｊ＝1, 2,...Ｎである。ｒejは、時刻ｔjに放送暦から得られた衛星のＥＣＥＦ位置座標を保持する。ｒ’ejは、放送暦から得られた衛星のＥＣＥＦ速度座標を保持する。ｓjは、各時刻ｔjにおける放送暦位置の推定誤差とする。本発明の実施例によれば、全てのｊについてｓj＝１メートルである。ｓ’jは、各時刻ｔjにおける放送暦速度の推定誤差とする。例えば、全てのｊについてｓ’j＝0.001メートル／秒である。Ｔ＝Ｔ(ｔj,ＥＯＰ)は、ＥＣＥＦからＥＣＩへの回転行列とする。Ｔは、時刻の関数、かつ地球方位パラメータ、ＥＯＰの関数である。換言すれば、Ｔは放送暦の各時期における地球方位を表す。
【００４９】
Ｔを計算する好適な方法は、a) 地球方位を、時刻の関数、かつＥＯＰ（ＵＴ１−ＵＴＣ、ＬＯＤ、Ｘ極運動、Ｙ極運動、及びこれらの速度）の関数として計算し、b) 地球方位のＥＯＰに対する感度を計算する手法を用いる方法である。例えば、手法a)及びb)は、地球の歳差運動、章動、極運動、及び地球の回転のモデル化についての２００３年ＩＥＲＳ議定書に従うことができる（IERS Conventions (2003). Dennis D. McCarthy and Gerard Petit. Verlag des Bundesamts fur Katograpgic und Geodasie, 2004（非特許文献２）参照、その全文を参考文献として本明細書に含める）。これらの方法は、受信機（２１０）の部分（２０３）に記憶されているソフトウエェアルーチンとすることができる。
【００５０】
従って、図３のステップ（Ｓ５）に示すように、現在のＥＯＰ値を用い、現在のＥＯＰ値に対して、ＥＣＩ変換した放送軌道暦の偏導関数（感度）を計算することによって、放送軌道暦をＥＣＥＦからＥＣＩに変換することができる。計算した感度はその後、受信機（２１０）内に記憶する。
【００５１】
反復推定法に戻って説明する。ｘ及びｘ’を三次元ベクトルとし、これらのベクトルはそれぞれ、第１時点ｔ1における衛星のＥＣＩ位置及び速度を保持する。最初に、ｘはＴｒe1によって得ることができ、ｘ’はＴｒ’e1＋Ｔ’ｒe1によって得ることができる。本発明による方法の一実施例では、衛星暦からの速度は使用せず、この場合にｘ’は、Ｔｒe1、Ｔｒe2、ｔ1、ｔ2を用いていわゆるガウス問題を解くことによって得ることができる（Bate et al, fundamental of Astrodynamics, Dover, 1971, Chapter 5.2の227〜231ページ参照、 その全文を参考文献として本明細書に含める）。後に、ｘ及びｘ’を推定し、これを反復する。
【００５２】
ｐdを、ｍ個の軌道動特性パラメータの集合とし、これらのパラメータの先験的な値は十分既知ではないが、ＧＮＳＳ衛星の軌道動特性に対して強い影響を有する。本発明の実施例の１つによれば、ｐdは次の７つのパラメータから成る：太陽圧モデルの尺度（一般に「太陽尺度」と称する）、太陽−地球−衛星の平面に直交する方向の定加速度（一般に「Ｙバイアス」称する）、ラジアル（放射）方向及びクロストラック（軌道交差）方向の、軌道周期当たり２回の加速度を表す４つの調和係数、及びトラック（軌道）内の定加速度。
【００５３】
ｐeopを、ｎ個のＥＯＰパラメータの集合とし、これらのパラメータの先験的な値は十分既知ではないが、ＥＣＥＦフレームからＥＣＩフレームへの回転に対して強い影響を有する。本発明の実施例の１つによれば、ｐeopは次の３つのパラメータを備える：ＬＯＤ、Ｘ極運動、及びＹ極運動。
【００５４】
ｒi＝ｒi(ｘ, ｘ’, ｐd)を、時刻ｔjにおける衛星のモデル化したＥＣＩ位置座標とする。ｒiは、初期状態ｘ及びｘ’の関数、かつ動的パラメータｐdの関数であり、時刻ｔ1から時刻ｔjまでの運動方程式の積分によって得られる。ｒiは、ｐeopに弱くしか依存せず、従って、この依存性は明示的に指定しない。同様に、ｒ’i＝ｒ’i(ｘ, ｘ’, ｐd）を、時刻ｔjにおける衛星のモデル化したＥＣＩ速度座標とする。この速度座標が放送暦から抽出されない際は、これらの速度座標はモデル化もされない。
【００５５】
従って、図３のステップ（Ｓ６）に示すように、衛星位置予測モデルは、ｘ、ｘ’、ｐd及びｐeopの現在値から開始して生成される。一旦、このことを行うと、図３のステップ（Ｓ７）に示すように、パラメータｘ、ｘ’、ｐd、及びｐeopは、ＥＣＩ座標で表されるこうした予測モデル及び放送暦に基づいて推定される。
【００５６】
上述した運動／予測モデルの方程式は、例えば、a) 初期状態（所定時刻における位置及び速度）及び特定のモデルパラメータ（例えば地球の重力定数、衛星の質量）の集合からの衛星軌跡の伝播を可能にし、b) 伝播する状態の、こうした初期状態及びモデルパラメータに対する感度を計算するＧＮＳＳ軌道モデルを用いて得ることができる。別の言い方をすれば、ＧＮＳＳ軌道モデルを提供することができ、このモデルは、運動方程式と変分方程式とを十分な忠実性で統合する。このモデルの忠実性は、結果的に予測される軌道の精度に直接影響を与える。本発明の実施例の１つによれば、ＧＮＳＳ軌道モデルは、４〜５次の可変ステップサイズのルンゲクッタ−フェールベルグ（Runge-Kutta-Fehlberg）法の常微分方程式積分器（例えば、Numerical Recipes, Press et al, Cambridge University Press, 1989, Chapter 15.1-15.2の550〜560ページ参照、その全文を参考文献として本明細書に含める）、及び次の軌道動特性のモデルを採用する：地球、太陽、及び月の相対位置についてのＪＰＬ惑星暦を用いた地球、太陽、及び月についての質点引力；次数１２までの地球重力場；Bar-Sever and Kuangの経験的モデル [IPN Progress Report 42-160 JPL, 2004, http://ipnpr.jpl.nasa.gov/progress_report/42-160/title.htm参照、その全文を参考文献として本明細書に含める]のような太陽圧モデル、このモデルはＧＰＳ衛星の軌道予測を改善することが示されている；米国特許第５９６３１６７号明細書（特許文献１）及び第５８２８３３６号明細書（共に、その全文を参考文献として本明細書に含める）に示されているモデルのような地球と月の円錐形の影のモデル；米国特許第６２９５０２１号明細書（特許文献３、その全文を参考文献として本明細書に含める）に示されているＧＰＳ衛星用のモデルのような衛星姿勢モデル；ＲＦスラストモデル；及び相対性モデル。ＧＮＳＳ軌道モデルは、受信機（２１０）の部分（２０３）に記憶されているソフトウェアルーチンとして実現することができる。
【００５７】
反復推定法に戻って説明し。表記目的で、ｄを、ｘ、ｘ’、ｐdから成る複合的な(６＋ｍ)次元ベクトルとする。即ち、ｄt=[ｘt；ｘ’t；ｐd]であり、上付き文字tは移項演算を示す。
は動特性の勾配演算子であり、次式のように定義される：
【数１】
【００５８】
そして、
は、ｘ、ｘ’、及びｐdに対するｒiの感度の３×(６＋ｍ)行列であり、運動（以上で説明したモデル参照）の変分方程式を時刻ｔ1から時刻ｔjまで積分することによって計算される。必要であれば、
を同様に定義し、これは、運動の変分方程式を積分することによって計算される。
をＥＯＰ勾配演算子とし、次式のように定義する。
【数２】
ここに、
は、上述したＴを得る方法によって計算することができる。
【００５９】
ここで、最小二乗問題を次のように定式化する。最小二乗の意味で、ＥＣＩフレーム内で観測される位置Ｔｒejと、ＥＣＩフレーム内でモデル化した位置ｒj（ｊ＝１,...Ｎ）との間の距離を最小にするパラメータｄとｐeopの集合を見出す。この問題を数値的に解くために、この問題を、解となるパラメータｄ及びｐeopの先験的な値付近で線形化し、結果的な線形最小二乗問題を解いて、これらのパラメータについての調整値Δd及びΔｐeopを見出す。次に、この問題を、解となるパラメータの調整値付近で再線形化し、等々によって、解が収束するまで、この最小二乗の解法を反復する。ここで、線形最小二乗問題は次式のようになる：
【数３】
が
【数４】
を近似する。
【００６０】
未知数を左辺に分離することによって、次式のような最小二乗問題の古典的定式化が生じる：
【数５】
が、Ｔ(ｐeop)ｒej−ｒj(ｄ, ｐeop)；ｊ＝１,...Ｎを近似し、関連するデータノイズｓjを伴う。この数値的な最小二乗問題は、平方根情報フィルタ（ＳＲＩＦ：Square Root Information Filter）因数分解法を用いて解くことができる。可能な方法は、Housholderの正規直交三角測量法、上三角行列の逆行列を求める方法、及びHousholderの正規直交変換を用いてＳＲＩＦ配列を結合する方法[例えば、BiermanによるFactorization Methods for Discrete Sequential Estimation, Dover, 2006（非特許文献６）参照、その全文を参考文献として本明細書に含める]を含む。これらのアルゴリズムは不揮発性メモリに記憶することができる。
【００６１】
図３のステップ（Ｓ８）に示すように、ステップＳ５〜Ｓ７を、a) 最終の反復における衛星軌道が、この衛星の放送暦に十分近くなるまで、あるいはb) 推定したパラメータｘ、ｘ’、ｐd及びｐeopの値が、前回の反復から明らかに変化しなくなるまで反復する。換言すれば、解となる「デルタ（Δ）」の値が所定閾値を下回れば、即ち、小さ過ぎて問題にならなければ、あるいは、モデルと観測値との間の距離が所定閾値未満であれば、最小二乗問題が収束したものと見なし、反復プロセスを停止する。これらの条件を満足せず、かつ反復回数が所定閾値を超えなければ、次の新たな仮のパラメータ値：
ｄ＝ｄ＋Δｄ 及び ｐeop＝ｐeop＋Δｐeop
で最小二乗問題を反復する。
【００６２】
本発明の一実施例によれば、この最小二乗パラメータ推定問題は、各軌道面内に１つの衛星について解かれ、ＧＰＳの場合は６つの衛星（６つの軌道面を有する）となる。他の実施例は、衛星の異なる選択についての問題を解くことができ、ＣＰＵ時間と解の精度とのトレードオフになる。最小二乗問題は、衛星毎に大きく異なる動的パラメータｐdの集合を生じさせることが予期される。しかし、ＥＯＰパラメータの推定値ｐeopも異なり得る。ＥＯＰパラメータの推定値ｐeopの精度を向上させるために、図３のステップ（Ｓ９）に示すように、上述したＳＲＩＦ因数分解法を用いて、各衛星からのＥＯＰパラメータの推定値どうしを、それぞれの共分散行列全体と共に組み合わせることができる。このことは、単一の強固なＥＯＰ値の推定を生じさせる。
【００６３】
一旦、不確定のＥＯＰ値ｐeopが強固に推定されると、全ての衛星についての最小二乗問題を形成して解くことができ、ここでは、初期状態パラメータｘ、ｘ’、及び動的パラメータｐdのみを推定し、ＥＯＰは、最終の組み合わせた推定値に固定して保持する。このことも図３のステップ（Ｓ１０）に示す。このプロセスの結果は、全ての衛星についての、初期状態（ｘ, ｘ’）と動的パラメータｐdとの集合となる。図３のステップ（Ｓ１１）〜（Ｓ１３）に示すように、これらの要素は、（以上で説明したモデルと共に）これらの衛星の軌道の将来に向かう伝播を可能にする。これらの値は受信機の不揮発性メモリ（１０４）に記憶する。こうした値は、信号取得及び位置決めを支援するために用いることのできる予測軌道を表す。
【００６４】
ここで図４を参照し、図４は、本発明の実施例による衛星クロック状態を予測する方法のフローチャートを示す。
【００６５】
ステップ（Ｔ１）に示すように、衛星毎の最新の放送クロック記録をメモリからロードする。好適な実施例では、これらの記録は少なくとも1日にわたる。ステップ（Ｔ２）では、衛星毎に、クロックモデルを放送暦値に適合させることによってクロックのモデルパラメータを決定する。
【００６６】
ＧＮＳＳクロックの仮のモデルは、（時間の二次式）＋（１回転に１回の調波）であり、衛星当たり５つのモデルパラメータを生じさせる（うち３つは二次多項式用、２つは１回転に１回の調波のサイン及びコサイン成分用）。これらのモデルパラメータは、放送暦から抽出した数日分の観測クロック値へのモデルの最小二乗適合を用いて推定することができる。図４の反復ステップ（Ｔ３）を参照。図４のステップ（Ｔ４）に示すように、一旦、衛星毎に個別にモデルパラメータを推定すると、衛星毎のモデルを時間的に前向きに、無限に伝播させて、予測クロック値を生成する。衛星毎に推定したモデルパラメータは、受信機（２１０）のメモリ（２０４）に記憶する。これらをいつでも呼び出して、クロック予測を生成することができる。本発明の一実施例によれば、４日分の放送クロック値を、モデルへの最小二乗適合に用いることができる。ＧＮＳＳのオペレータによって、モデルと不整合なクロックジャンプ（クロックの急変）が時として導入され得るので、最小二乗適合誤差（二乗平均平方根：ＲＭＳ：Root-Mean-Square）を所定閾値に対してチェックし、この閾値は仮に１メートルとする。ＲＭＳ適合がこの閾値を超えれば、適合間隔のうち最も適合が見られない１２時間の部分を試しに除去して、モデル不整合のない間隔を生成することができる。このプロセスは、適合のＲＭＳ値が閾値を下回るまで、あるいは適合間隔が１２時間未満である場合に反復することができる。２日間未満の適合間隔については、線形クロックモデルを適用すべきである（衛星当たり２つのパラメータ）。通常のように、こうしたアルゴリズムモデルは受信機（２１０）の部分（２０３）に記憶することができる。
【００６７】
予測プロセスの最後に、伝播させた軌道及び予測したクロックモデルを利用して、将来のＧＮＳＳ信号の取得に当たり受信機を支援することができる。新たな暦データ記録が利用可能になる毎に、受信機の電力管理の考慮及び他の動作上のトレードオフに応じて、この予測プロセスを反復することができる。
【００６８】
上述したように、また図２を参照して述べたように、本発明の実施例の１つによれば、本発明による受信機は不揮発性メモリ部分（２０４）を備え、この不揮発性メモリ部分は次のものを含む：
１）以上で説明した種々のソフトウェアルーチンのオブジェクトコード、例えばi) ＥＯＰ予測モデル；ii) 回転行列を生成する方法；iii) 運動方程式を統合したＧＮＳＳ軌道モデル；iv) ＳＲＩＦアルゴリズム；v) 軌道初期状態、動的パラメータ、及びＥＯＰの反復推定方法；vi) クロック予測方法；及びvii) ＧＮＳＳ放送暦データの最適なサンプリング方法。
２）最新の記録放送暦。
３）初期化目的用の、軌道の動的パラメータｄのデータベース。
４）上記ＥＯＰ予測モデルを用いて生成した先験的ＥＯＰ値。
５）上記１）で述べた方法を適用した結果として得られた、予測した軌道状態及びクロック状態。
【００６９】
要約すれば、本明細書におけるいくつかの実施例によれば、全地球航法衛星システム（ＧＮＳＳ）の軌道状態及びクロック状態の自律的な受信機内予測の方法及び装置を説明している。ＧＮＳＳ放送メッセージのみを使用し、周期的な外部通信情報を必要としない。地球方位情報をＧＮＳＳ放送暦から抽出する。地球方位パラメータの正確な推定により、慣性基準フレーム内で、最良に適合するＧＮＳＳ軌道を、時間的に前向きに伝播させることができる。次に、推定した地球方位パラメータを用いて、予測した軌道を地球中心−地球固定（ＥＣＥＦ）座標に変換して、信号の取得に当たりＧＮＳＳ受信機を支援するために使用する。ＧＮＳＳ衛星のクロック状態も放送暦から抽出して、このデータに、パラメータ化したクロック動作のモデルを適合させる。次に、モデル化したクロックを推定し、時間的に前向きに伝播させて、予測した軌道と共に、より迅速なＧＮＳＳ信号取得を可能にする。
【００７０】
従って、これまで示してきたものは、ＧＮＳＳ暦の自律的な受信機内予測の方法及び装置である。これらの方法及び装置は、特定実施例及びその応用を用いて説明してきたが、本発明の範囲を逸脱することなしに、当業者が多数の変更及び変形を加え得ることは明らかである。従って、特許請求の範囲内で、本明細書で具体的に説明した以外の方法で本発明を実施することができることは明らかである。
【特許請求の範囲】
【請求項１】
衛星の軌道情報を地球中心−地球固定（ＥＣＥＦ）座標で自律的に予測して、衛星受信機を支援する方法において、
i）前記衛星の放送軌道暦を取得するステップと；
ii）収束が得られるまで、前記衛星の軌道位置、速度、及び動的パラメータを、前記衛星の動特性モデルに基づいて反復的に推定し、地球方位パラメータ（ＥＯＰ）を、地球方位モデルに基づいて反復的に推定するステップと；
iii）前記推定した軌道位置、速度、及び動的パラメータを、前記衛星の動特性モデルを用いて時間的に前向きに伝播させ、これにより予測軌道を得るステップと；
iv）前記地球方位モデルを用いて、前記予測軌道を地球中心−地球固定（ＥＣＥＦ）座標に変換するステップと；
v）前記変換した予測軌道を記憶するステップと
を含むことを特徴とする方法。
【請求項２】
ステップi）、ii）、iii）、iv）及びv）を複数の衛星に適用することによって、前記複数の衛星の軌道情報を自律的に予測することを特徴とする請求項１に記載の方法。
【請求項３】
１つの衛星データから推定したＥＯＰ値及びその共分散行例を、他の衛星データからのＥＯＰの推定値及びその共分散行列と最適に組合せることを特徴とする請求項１または２に記載の方法。
【請求項４】
前記ＥＯＰが、Ｘ極運動、Ｙ極運動、及び日長（ＬＯＤ）を含むことを特徴とする請求項１〜３のいずれかに記載の方法。
【請求項５】
前記取得した前記衛星の放送軌道暦を、第１座標系から第２座標系に変換するステップを更に含むことを特徴とする請求項１〜４のいずれかに記載の方法。
【請求項６】
前記反復的な推定中に、前記変換した放送軌道暦のＥＯＰに対する感度を計算するステップを更に含むことを特徴とする請求項５に記載の方法。
【請求項７】
前記第１座標系が地球中心−地球固定（ＥＣＥＦ）座標系であり、前記第２座標系が地球中心慣性（ＥＣＩ）座標系であることを特徴とする請求項５または６に記載の方法。
【請求項８】
前記衛星または前記複数の衛星、及び前記受信機が、全地球航法衛星システム（ＧＮＳＳ）の一部分であることを特徴とする請求項１〜７のいずれかに記載の方法。
【請求項９】
前記衛星の前記放送暦を、サンプリングプロセスを通して取得することを特徴とする請求項１〜８のいずれかに記載の方法。
【請求項１０】
前記サンプリングプロセスを、少なくとも１日３回発生させることを特徴とする請求項９に記載の方法。
【請求項１１】
前記サンプリングプロセスを、１日の周期性を回避するサンプリング周波数で発生させることを特徴とする請求項９または１０に記載の方法。
【請求項１２】
前記サンプリング周波数を、７時間間隔のサンプリング周波数及び１０時間間隔のサンプリング周波数から成るグループから選択することを特徴とする請求項１１に記載の方法。
【請求項１３】
前記地球方位モデルが、測定したＥＯＰ値の履歴的時系列へのパラメータ適合に基づくことを特徴とする請求項３に記載の方法。
【請求項１４】
ステップii)を、衛星状態への最小二乗適合を用いて実行することを特徴とする請求項１〜１３のいずれかに記載の方法。
【請求項１５】
前記衛星状態が、地球中心−地球固定（ＥＣＥＦ）座標の衛星状態であることを特徴とする請求項１４に記載の方法。
【請求項１６】
前記最小二乗適合を、平方根情報フィルタ（ＳＲＩＦ）因数分解法を用いて解くことを特徴とする請求項１５に記載の方法。
【請求項１７】
前記ステップii)における収束は、ａ)現在推定した前記衛星の軌道位置、速度、及び軌道の動的パラメータが、当該衛星の放送暦に十分に近い際に、あるいはｂ)現在推定した前記衛星の軌道位置、速度状態、及び軌道の動的パラメータが、それぞれの前の値に十分近い際に得られることを特徴とする請求項１〜１６のいずれかに記載の方法。
【請求項１８】
前記予測軌道を、前記衛星受信機内に記憶することを特徴とする請求項１〜１７のいずれかに記載の方法。
【請求項１９】
衛星のクロック状態を自律的に予測する方法において、
i）前記衛星についての放送クロック記録を取得するステップと；
ii）収束が得られるまで、前記衛星の将来のクロック状態を、予測モデルに基づいて反復的に推定するステップと；
iii）前記推定したクロック状態を、前記衛星のクロック状態の自律的予測用に記憶するステップと
を含むことを特徴とする方法。
【請求項２０】
ステップi）、ii）及びiii）を複数の衛星に適用することによって、前記複数の衛星のクロック状態を自律的に予測することを特徴とする請求項１９に記載の方法。
【請求項２１】
前記予測モデルが、時間の二次式に、１回転に１回の調波を加えたものであり、衛星当たり５つのモデルパラメータを生じさせることを特徴とする請求項１９または２０に記載の方法。
【請求項２２】
衛星の軌道情報及び衛星のクロック状態を自律的に予測する方法において、
請求項１〜１８のいずれかに記載の方法により、前記衛星の軌道情報を自律的に予測するステップと；
請求項１９〜２１のいずれかに記載の方法により、前記衛星クロック状態を自律的に予測するステップと
を含むことを特徴とする方法。
【請求項２３】
処理装置と、
メモリ装置とを備え、
前記メモリ装置は、
衛星の軌道予測用及び前記衛星のクロック状態予測用のソフトウエアコードを記憶する第１メモリ部分と；
記録した前記衛星の放送軌道及び前記衛星のクロック状態を記憶する第２メモリ部分と、
地球方位パラメータ（ＥＯＰ）値及び軌道の動的パラメータを記憶する第３メモリ部分と、
前記第１メモリ部分に記憶されているソフトウエアコードに従って計算して予測した前記衛星の軌道状態及び前記衛星のクロック状態を記憶する第４メモリ部分と
を備えるように構成されることを特徴とする全地球航法衛星システム（ＧＮＳＳ）受信機。
【請求項２４】
前記衛星の軌道予測用の前記ソフトウエアコードが、請求項１に記載の方法に従って機能し、前記衛星のクロック状態予測用の前記ソフトウエアコードが、請求項２０に記載の方法に従って機能することを特徴とする請求項２３に記載のＧＮＳＳ受信機。
【請求項２５】
請求項２３または２４に記載のＧＮＳＳ受信機を備えた個人用ナビゲーション装置。
【請求項２６】
携帯電話またはカーナビゲーション装置であることを特徴とする請求項２５に記載の個人用ナビゲーション装置。
【請求項１】
衛星の軌道情報を地球中心−地球固定（ＥＣＥＦ）座標で自律的に予測して、衛星受信機を支援する方法において、
i）前記衛星の放送軌道暦を取得するステップと；
ii）収束が得られるまで、前記衛星の軌道位置、速度、及び動的パラメータを、前記衛星の動特性モデルに基づいて反復的に推定し、地球方位パラメータ（ＥＯＰ）を、地球方位モデルに基づいて反復的に推定するステップと；
iii）前記推定した軌道位置、速度、及び動的パラメータを、前記衛星の動特性モデルを用いて時間的に前向きに伝播させ、これにより予測軌道を得るステップと；
iv）前記地球方位モデルを用いて、前記予測軌道を地球中心−地球固定（ＥＣＥＦ）座標に変換するステップと；
v）前記変換した予測軌道を記憶するステップと
を含むことを特徴とする方法。
【請求項２】
ステップi）、ii）、iii）、iv）及びv）を複数の衛星に適用することによって、前記複数の衛星の軌道情報を自律的に予測することを特徴とする請求項１に記載の方法。
【請求項３】
１つの衛星データから推定したＥＯＰ値及びその共分散行例を、他の衛星データからのＥＯＰの推定値及びその共分散行列と最適に組合せることを特徴とする請求項１または２に記載の方法。
【請求項４】
前記ＥＯＰが、Ｘ極運動、Ｙ極運動、及び日長（ＬＯＤ）を含むことを特徴とする請求項１〜３のいずれかに記載の方法。
【請求項５】
前記取得した前記衛星の放送軌道暦を、第１座標系から第２座標系に変換するステップを更に含むことを特徴とする請求項１〜４のいずれかに記載の方法。
【請求項６】
前記反復的な推定中に、前記変換した放送軌道暦のＥＯＰに対する感度を計算するステップを更に含むことを特徴とする請求項５に記載の方法。
【請求項７】
前記第１座標系が地球中心−地球固定（ＥＣＥＦ）座標系であり、前記第２座標系が地球中心慣性（ＥＣＩ）座標系であることを特徴とする請求項５または６に記載の方法。
【請求項８】
前記衛星または前記複数の衛星、及び前記受信機が、全地球航法衛星システム（ＧＮＳＳ）の一部分であることを特徴とする請求項１〜７のいずれかに記載の方法。
【請求項９】
前記衛星の前記放送暦を、サンプリングプロセスを通して取得することを特徴とする請求項１〜８のいずれかに記載の方法。
【請求項１０】
前記サンプリングプロセスを、少なくとも１日３回発生させることを特徴とする請求項９に記載の方法。
【請求項１１】
前記サンプリングプロセスを、１日の周期性を回避するサンプリング周波数で発生させることを特徴とする請求項９または１０に記載の方法。
【請求項１２】
前記サンプリング周波数を、７時間間隔のサンプリング周波数及び１０時間間隔のサンプリング周波数から成るグループから選択することを特徴とする請求項１１に記載の方法。
【請求項１３】
前記地球方位モデルが、測定したＥＯＰ値の履歴的時系列へのパラメータ適合に基づくことを特徴とする請求項３に記載の方法。
【請求項１４】
ステップii)を、衛星状態への最小二乗適合を用いて実行することを特徴とする請求項１〜１３のいずれかに記載の方法。
【請求項１５】
前記衛星状態が、地球中心−地球固定（ＥＣＥＦ）座標の衛星状態であることを特徴とする請求項１４に記載の方法。
【請求項１６】
前記最小二乗適合を、平方根情報フィルタ（ＳＲＩＦ）因数分解法を用いて解くことを特徴とする請求項１５に記載の方法。
【請求項１７】
前記ステップii)における収束は、ａ)現在推定した前記衛星の軌道位置、速度、及び軌道の動的パラメータが、当該衛星の放送暦に十分に近い際に、あるいはｂ)現在推定した前記衛星の軌道位置、速度状態、及び軌道の動的パラメータが、それぞれの前の値に十分近い際に得られることを特徴とする請求項１〜１６のいずれかに記載の方法。
【請求項１８】
前記予測軌道を、前記衛星受信機内に記憶することを特徴とする請求項１〜１７のいずれかに記載の方法。
【請求項１９】
衛星のクロック状態を自律的に予測する方法において、
i）前記衛星についての放送クロック記録を取得するステップと；
ii）収束が得られるまで、前記衛星の将来のクロック状態を、予測モデルに基づいて反復的に推定するステップと；
iii）前記推定したクロック状態を、前記衛星のクロック状態の自律的予測用に記憶するステップと
を含むことを特徴とする方法。
【請求項２０】
ステップi）、ii）及びiii）を複数の衛星に適用することによって、前記複数の衛星のクロック状態を自律的に予測することを特徴とする請求項１９に記載の方法。
【請求項２１】
前記予測モデルが、時間の二次式に、１回転に１回の調波を加えたものであり、衛星当たり５つのモデルパラメータを生じさせることを特徴とする請求項１９または２０に記載の方法。
【請求項２２】
衛星の軌道情報及び衛星のクロック状態を自律的に予測する方法において、
請求項１〜１８のいずれかに記載の方法により、前記衛星の軌道情報を自律的に予測するステップと；
請求項１９〜２１のいずれかに記載の方法により、前記衛星クロック状態を自律的に予測するステップと
を含むことを特徴とする方法。
【請求項２３】
処理装置と、
メモリ装置とを備え、
前記メモリ装置は、
衛星の軌道予測用及び前記衛星のクロック状態予測用のソフトウエアコードを記憶する第１メモリ部分と；
記録した前記衛星の放送軌道及び前記衛星のクロック状態を記憶する第２メモリ部分と、
地球方位パラメータ（ＥＯＰ）値及び軌道の動的パラメータを記憶する第３メモリ部分と、
前記第１メモリ部分に記憶されているソフトウエアコードに従って計算して予測した前記衛星の軌道状態及び前記衛星のクロック状態を記憶する第４メモリ部分と
を備えるように構成されることを特徴とする全地球航法衛星システム（ＧＮＳＳ）受信機。
【請求項２４】
前記衛星の軌道予測用の前記ソフトウエアコードが、請求項１に記載の方法に従って機能し、前記衛星のクロック状態予測用の前記ソフトウエアコードが、請求項２０に記載の方法に従って機能することを特徴とする請求項２３に記載のＧＮＳＳ受信機。
【請求項２５】
請求項２３または２４に記載のＧＮＳＳ受信機を備えた個人用ナビゲーション装置。
【請求項２６】
携帯電話またはカーナビゲーション装置であることを特徴とする請求項２５に記載の個人用ナビゲーション装置。
【図１】
【図２】
【図３】
【図４】
【図２】
【図３】
【図４】
【公表番号】特表２０１２−５０２２９８（Ｐ２０１２−５０２２９８Ａ）
【公表日】平成２４年１月２６日（２０１２．１．２６）
【国際特許分類】
【出願番号】特願２０１１−５２６９００（Ｐ２０１１−５２６９００）
【出願日】平成２１年８月２１日（２００９．８．２１）
【国際出願番号】ＰＣＴ／ＵＳ２００９／０５４６３５
【国際公開番号】ＷＯ２０１０／０４７８７５
【国際公開日】平成２２年４月２９日（２０１０．４．２９）
【出願人】（３０１０５９５７０）カリフォルニア インスティチュート オブ テクノロジー (14)
【Ｆターム（参考）】
【公表日】平成２４年１月２６日（２０１２．１．２６）
【国際特許分類】
【出願日】平成２１年８月２１日（２００９．８．２１）
【国際出願番号】ＰＣＴ／ＵＳ２００９／０５４６３５
【国際公開番号】ＷＯ２０１０／０４７８７５
【国際公開日】平成２２年４月２９日（２０１０．４．２９）
【出願人】（３０１０５９５７０）カリフォルニア インスティチュート オブ テクノロジー (14)
【Ｆターム（参考）】
[ Back to top ]