JP2008513775A - ナビゲーション用途のための、改良されたｇｐｓ累積デルタ距離処理方法 - Google Patents

Abstract

移動物体の位置と速度の測定に使用されるＧＰＳナビゲーション技術において、擬似距離（ＰＲ）測定値および累積デルタ距離（ＡＤＲ）測定値が、物体において受信ＧＰＳ信号から生成される。ナビゲーションパラメータ（例えば位置、速度、クロック）はＰＲ測定値、およびＡＤＲ測定値の間の差から推定される。本明細書で述べるＡＤＲ測定式は、よりいっそう精密な方法で公式化され、ＡＤＲ測定値の間の差を計算する時間間隔は、現在のＧＰＳ用途のＡＤＲ差分手法で使用されるものよりもはるかに大きい。故にＡＤＲの差はより正確であり、非常に精密なナビゲーションソリューションを生じる。ＡＤＲ差分手法はカルマンフィルター処理の収束時間短縮に寄与し、宇宙飛行体のナビゲーション精度を向上する。これらの高精度なＡＤＲ処理アルゴリズムを統合型ＧＰＳ／ＩＭＵナビゲーション用途に拡げるための技術も提供される。

Description

（関連出願の引用）
本出願は、米国特許仮出願第６０／６１０，６０９号、名称「ＧＰＳ ＡＤＲ Ｐｒｏｃｅｓｓｉｎｇ ｆｏｒ Ｎａｖｉｇａｔｉｏｎｏｆ ａ Ｓｐａｃｅｃｒａｆｔ」（２００４年９月１７日出願）に対する優先権を主張し、該仮出願の全体が本明細書において参考として援用される。

（技術分野）
本発明は、グローバル測位システム（ＧＰＳ）ナビゲーション技術に関し、特に、より精密なナビゲーション結果を得るための、累積デルタ距離（ＡＤＲ）測定値を使用するナビゲーションアルゴリズムに関する。

ＧＰＳは、宇宙飛行体、航空機だけでなく、地上の自動車または個人のナビゲーションの手段として十分に確立されている。基本的に、使用者は、ＧＰＳ衛星から信号を受信し、ＧＰＳ信号に含まれるデータからの擬似距離（ＰＲ）測定値を使用して、位置、速度および他のパラメータを決定する。用語「利用者」とは、その位置がＧＰＳまたはそれが受信する他の測距信号に基づいて計算される、物体を意味する。用途が高精度のナビゲーションを必要とする場合、ＰＲ測定値は、ポイントソリューションアルゴリズムによるよりも、むしろカルマンフィルターによって処理される。カルマンフィルターは、ＰＲ測定値および利用者の位置、速度および他の状態ベクトルパラメータを推定するための複雑な伝播モデルを使用する。伝播モデルは、時間Ｎ−１での状態が既知である場合に、時間Ｎでの利用者の状態を計算するように設計されている。これらは、カルマンフィルターが現在の推定を計算するために、利用者の状態の事前情報を使用することを可能とし、よって精度が向上する。

残念なことに、伝播モデルはかなり精密であるが、伝播された状態の精度を保証しない。伝播モデルは、インプットとして以前の基準時点での状態ベクトル推定値を使用し、後者の精度は限定要因となり得る。伝播された状態の精度は、特に、物体速度推定値の誤差を受けやすい。速度は、位置との相関のために、カルマンフィルターによって間接的に推定される。

高精度なナビゲーションの達成は、正確な推定速度と密接に連結する。これは、位置の正確さと速度の正確さの間に円環状の依存関係を作り出す。実施において、推定処理は段階的であり、速度および位置の精度の向上は、ある時間にわたって相互を促進する。全体の処理は、カルマンフィルターの収束とみなされる。典型的な例として、宇宙利用におけるフィルター収束は、数時間から１日またはそれ以上かかる。これは、フィルターを再収束しなければならないときには、特に宇宙飛行体の操縦にとって問題である。

生成されるナビゲーションソリューションおよびカルマンフィルター計算の収束時間の双方において、ＧＰＳナビゲーション処理の精度向上の余地が存在する。

簡単に説明すると、移動物体の位置および速度を測定するために使用されるＧＰＳナビゲーションの技術を提供する。擬似距離（ＰＲ）および累積デルタ距離（ＡＤＲ）測定値は、物体において、受信ＧＰＳ信号から測定される。連続するＡＤＲ測定値の間の時間間隔よりも大きな時間間隔によって分離される、ＡＤＲ測定値の間の差が計算される。ナビゲーションパラメータ（例えば、位置、速度およびクロック）は、ＰＲ測定値およびＡＤＲ差から推定される。

本明細書で説明される技術の一つの類似は、以下のようである。車で移動している個人が、車が移動し始めるときにストップウォッチを始動する。ある程度の時間間隔の後、個人は、ストップウォッチを停止し、車が移動した距離を測定する。ストップウォッチが始動および停止されるごとに、誤差が生じる。短時間で距離を測定する場合、車の速度の推定は、非常に正確であるとは限らない。逆に、車が、一定の速度で走行する場合、および個人が長時間測定をする場合、速度の推定は、非常に正確となる。二番目の選択が好ましいが、一定速度の仮定は妥当でないことがあり得るので、それが常に実施可能とは限らない。さらに比較的長時間間隔の速度測定は、この時間間隔に対する平均速度を算出する。ほとんどの用途が現在の状態の情報を必要とするために、平均速度はあまり実用性がない。

本明細書で説明されるこの技術は、比較的長時間間隔で物体速度を測定することの利点を、そのような測定を処理時間における瞬間速度に正確にマップする能力と組み合わせる。これは、異なったＡＤＲ測定値を処理することと、そのような差を利用者の現在の状態ベクトルに正確にマップする数式を適用することによって達成される。

本明細書で述べるＡＤＲ測定式は、ＡＤＲ差を計算するために使用されるＡＤＲ測定値の間の時間間隔が、ＧＰＳナビゲーションの用途において、現在のＡＤＲ差分技術で使用される時間間隔よりも非常に大きいために、非常に精密な方法で公式化されている。従って、ＡＤＲ差はより正確であり、これは、より正確なナビゲーションソリューションを生じる。さらに、ＡＤＲ差分技術は、カルマンフィルター処理の収束時間の短縮に貢献するために、ナビゲーションの精度を向上する。

ＡＤＲ測定値の使用が困難なのは、カルマンフィルターが瞬間測定値を必要とするのに対し、これらがある程度の時間間隔に対して累積／積分された量に関係することである。この困難は、ＡＤＲ測定値を瞬間測定値に正確にマップする数式により本明細書に説明される技術を用いて克服される。さらに、このような瞬間測定値は、第一義的に速度に対して提供され、よって、物体の速度を推定するための直接的な手段を提供する。

本明細書で説明される数式を用いたＡＤＲ測定値の使用には多くの利点がある。伝播モデルと併用した速度推定値は、直ちにより正確な伝播状態をもたらし、それはカルマンフィルターの収束時間の大幅な削減を可能にする。ＰＲに加えてＡＤＲ測定値の使用は、処理に使用される測定値の合計数を約２倍にする。さらに、ＡＤＲ測定値は、ＰＲ測定値より非常に低いノイズを有する。より多くの測定値とより低いノイズは、全体的により良い精度をもたらす。

信頼のおける伝播モデルが特定の移動物体に対して利用不可能であるナビゲーション用途に対して、これらの高精度なＡＤＲ処理アルゴリズムを拡げる技術も提供する。例えば、精密な伝播モデルは、ある種の航空機、地上のビークル、個人等に対して利用不可能である。しかし、これらの物体が、加速度および回転速度などの慣性計測値を測定する能力を有する場合には、慣性計測値は、ＡＤＲ差分処理を補助するために使用され得る。

添付図面と組み合わせて以下の説明を参照するときには、他の目的および利点がより容易に明らかとなる。

図１は、本明細書で説明するナビゲーション技術が使用可能な環境を示す。物体１０が示されており、その位置は、複数のグローバル測位システム（ＧＰＳ）衛星２０（１）〜２０（Ｎ）から受信した信号から計算される。物体１０は、衛星または宇宙船などの宇宙飛行体、地球の地上または空中のビークル（例えば、車、トラック、タンク、航空機、ヘリコプター、船など）であり得る。物体１０は、ＧＰＳ衛星２０（１）〜２０（Ｎ）から信号を受信できる受信機１２、受信信号から測定値を測定して物体の位置、速度、方位などを計算する処理装置１４と、所望するなら、個人にナビゲーション情報を伝達するための物体インターフェース１６（キーボードおよびディスプレイなど）を備える。さらに、物体は、物体の加速度および回転速度の測定値データをアウトプットする、慣性計測ユニット（ＩＭＵ）１８を含み得る。処理装置１４は、受信したＧＰＳ信号から得た測定値で計算を実行する。

受信機１４は、受信したＧＰＳ信号から擬似距離（ＰＲ）測定値および累積デルタ距離（ＡＤＲ）測定値を測定する。処理装置１４は、それから、測定式を使用してＰＲおよびＡＤＲ測定値のカルマンフィルター計算を実行し、伝播モデル計算を使用して、またはＩＭＵデータを使用して状態ベクトルを更新する。

図１は、ＧＰＳ衛星２０（１）〜２０（Ｎ）が一般的に位置を決定するために使用される一方、測距信号を送信する他の測距源が存在することも示す。参考数字２２は、地上システムまたは屋内位置追跡システムなどの他の測距源を表すことを意図する。ＧＰＳ信号は、本明細書で説明する技術が使用する測距信号のタイプの一例にすぎないことが、意図される。

カルマンフィルターは、様々なタイプ（例えば、ＰＲおよびＡＤＲ測定値）の複数の測定値から状態ベクトル（例えば、位置、速度およびクロック）を推定する最適な線形推定装置である。ＰＲおよびＡＤＲ測定値は反復基準で受信ＧＰＳ信号から測定され、各測定イベントが「基準時点（ｅｐｏｃｈ）」と呼ばれる。カルマンフィルターは、最小二乗の意味での、最小分散推定値を生成する。ナビゲーション状態に加えて、カルマンフィルターは、一般的に全体的なナビゲーション状態の誤差を推定する。フィルターは、誤差状態ベクトル推定値の精度の測定をアウトプットする。この精度の測定は、状態推定値の誤差の二次中心モーメントのマトリックスである、共分散マトリックスによって示される。

一般的に、カルマンフィルターは、アルゴリズムが利用可能な推定値の処理（測定式により示される測定値モデルを用いて）と状態ベクトルのソリューションの次の基準時点への伝播（伝播式により示される測定値モデルを用いて）との間で交替する、周期で適用される。伝播モデルは、どのように誤差状態ベクトルが時間で変化するかを説明するプロセスである。測定値モデルは、誤差状態ベクトルとフィルターにより処理されたあらゆる測定値間と間の関係を定義する。ソリューションの精度は、測定値の質および量、測定式の精度、および伝播ルーチンの精度に依存する。衛星状態の推定などの非線形問題に対して、測定式は一般的に線形化される。

線形化は、状態ベクトルの成分に関して測定された量の偏導関数の計算を含む。ＰＲ測定値におけるこのような偏導関数の計算は単純で、従来のナビゲーションアルゴリズムの基礎である。しかし、ＡＤＲ測定値における偏導関数の計算は、より複雑である。第一の問題は、時間の不一致である。状態ベクトルは、ある特定の時間での衛星状態を示すが、ＡＤＲ測定値は、ある時間間隔に対する位相測定値の累積を示す。一つの技術は、ＡＤＲ測定値を時間間隔から実時間へマップすることである。このようなマッピングは必然的に近似であり、測定値ノイズからの誤差に加え、ナビゲーションシステムに誤差を導入する。

概念的に、カルマンフィルターは、測定値および状態ベクトルの動的挙動の相対寄与に対して重み付けをする。測定値および状態ベクトルは、それらの相対共分散によって重み付けされる。測定値が状態ベクトル推定値と比較して不正確（大きい分散値）な場合、フィルターは、測定値を重視しない。測定値が、状態ベクトル推定値と比較して非常に正確（小さい分散値）な場合、フィルターは、非常に測定値を重視するため、以前に計算された状態推定値は、最新の状態推定値にほとんど寄与しない。

ＡＤＲデータの処理に対する単純なアプローチでさえ、根底にあるＡＤＲ測定値よりも、より正確な測定式を与える（つまり、測定式の誤差は、測定値ノイズより小さい）。しかし、ナビゲーション精度における測定式の一見非有意な誤差の作用は、かなり大きい。この理由は、この誤差が非常に時間に相関があり、よって、フィルターにより除去されないからである。

測定式の長期相関誤差による危険は、伝播における同様の作用を調査することにより認識することができる。測定値処理および伝播は、カルマンフィルターアルゴリズムの二つの工程であり、各周期で実行される。伝播モデルにおける小さい計算外の加速度誤差は、時間とともに累積され、ナビゲーション精度に影響を与えることは周知である。この作用は、複雑で正確な伝播モデルの改善を強いる。

同様の作用が測定式においても当てはまることは、直感的に明らかである。実際、測定式における誤差は、処理の更新工程（データ処理）での各基準時点で物体（例えば、宇宙飛行体）を「引く」傾向があり、伝播工程でモデル化されない加速度の作用と同様の作用をもつ。ナビゲーションソリューションに対するこのような誤差からの作用は、大きい可能性がある。その場合、長期間誤差を伴う測定値の使用は、ＡＤＲ測定値をまったく使用しないよりも悪い可能性がある。よって、ナビゲーションでのＡＤＲ測定値の使用は、長期相関誤差がほとんどまたは全くなく状態ベクトルに測定値をマップする能力に非常に依存する。この精度条件は、以下のように定量化することができる。測定式の相関誤差は、伝播式の誤差以下またはほぼ同等であるべきである。この場合、ＡＤＲに対して測定式を使用することは、精度をより良くする可能性がある。

これは、相対ナビゲーション、例えば、衛星編隊飛行に対しては明確な可能性を有する。差分および二倍差分はバイアスを除去し、相対ナビゲーションソリューションにおける特別な精度を可能にする。しかし、絶対的ナビゲーション用途に対してＡＤＲ測定値を適用するために、どのように精度なアルゴリズムを向上するかは、あまり明らかではない。これが、多くの以前の用途で宇宙飛行体の絶対的ナビゲーションに対してＡＤＲ測定値が使用されなかった理由であり得る。

本明細書で説明される技術の重要な側面は、宇宙飛行体の絶対的ナビゲーションに対して正確なＡＤＲ測定値処理を可能にする、数式である。このＡＤＲ処理は、従来のＰＲ処理と相補的であり、併用して使用される。説明される以下のアルゴリズムの主な区別点は、その精度である。測定式は、計算外の誤差を含むが、これらの誤差の規模は、モデル化されていない伝播誤差と同様の程度に低減される。

図２を参照すると、高レベルのナビゲーション処理の画像的表示が示されている。移動物体は、図３と併用して、以下に説明される技術を使用してＰＲ測定値およびＡＤＲ差（時間により分けられた二つのＡＤＲ測定値からの）を導き出すために受信ＧＰＳ信号を使用する。ＰＲ測定値およびＡＤＲ差は、ブロック７０により示されるカルマンフィルターに供給される。また、物体１０および／または物体１０がＩＭＵを含む場合のＩＭＵデータに対する伝播モデルは、既知の場合には、カルマンフィルターに供給される。ＩＭＵデータは、物体１０に対して信頼できる伝播モデルが知られていないときなど、統合型ＧＰＳ−ＩＭＵナビゲーション計算が実行される環境において役立つ。

図３に戻ると、ＧＰＳのみまたはＧＰＳ−ＩＭＵナビゲーションアルゴリズムの基本ステップが、説明されている。ステップ８０において、物体１０は、ＧＰＳ衛星からのＧＰＳナビゲーション信号または図１に示す他の距離源からの測距信号など、測距信号を受信する。次に、ステップ８２において、物体は、ステップ８０で受信した距測信号からＰＲおよびＡＤＲ測定値を測定する。ステップ８４において、差が、連続ＡＤＲ測定値間の時間間隔よりも大きな時間間隔によって分離された二つの基準時点（例えば、実時間）で、二つのＡＤＲ測定値の間で計算される。このＡＤＲ差分手法は、図４と併用してより詳細に説明される。ステップ８６において、ＡＤＲ差は、カルマンフィルターの修正偏導関数方程式を使用して処理され、状態ベクトルは、伝播モデル（およびＰＲ測定値）を使用して、または物体により得られたＩＭＵデータを使用して状態ベクトルを伝播するかのいずれかにより、伝播される。修正偏導関数方程式は、図５と併用して、後により詳細に説明する。ＩＭＵデータを使用して状態ベクトルを伝播する技術は、図６と併用して説明する。また、図６の技術は、測定される物体の位置が予測可能な加速度をもたず、よって、信頼できる伝播モデルが知られていない場合には、ナビゲーションの用途において役立つ。

（Ｉ．ＧＰＳのみおよびＧＰＳ−ＩＭＵ統合ナビゲーションに対するＡＤＲ計測処理方程式）
（Ａ．バイアスの相殺）
ＡＤＲデータは、未知の一定のバイアス（周期のずれがない場合）を有する。このバイアスのため、ＡＤＲ測定値は、ＰＲ測定値と同様の方法で適用することができない。バイアスを処理する周知の方法の一つは、例えば、ＰＲ測定値を平滑化するためにＡＤＲを使用することによって、バイアスを推定することである。平滑化よりも優れた技術が、ＡＤＲ測定値の処理およびバイアス除去のために、本明細書で説明される。

バイアスを除去するための一般的な方法は、以下の形で測定値を使用することである。

式中、Ｐ（ｔ）は、ＡＤＲ測定値であり、重み関数Ｗ（ｔ）は、

となるように選択される。方程式（１．１）が、ＡＤＲ測定値Ｐ（ｔ）に含まれる一定のバイアスを相殺することは、直接代入することによって検証され得る。重み関数Ｗ（ｔ）は、測定値の結果

における誤差を最小限にするため、およびカルマンフィルターにより処理された結果を改善するために選択され得る。

方程式（１．１）を使用する利点を例示するために、二つの基準時点ｔ_ｎ−１、ｔ_ｎでのＡＤＲ値の差

を単純にもたらす、重み関数Ｗ（ｔ）の特別な例が考慮される。これは恐らくバイアスを相殺する最も簡単な（最適ではないが）式であり、以下の測定値を生成する。

式中、

は、物体の位置ベクトル、

は、ＧＰＳ衛星の位置ベクトル、

は、物体クロックドリフト、および文字ｎ、ｎ−１は、二つの基準時点（連続とは限らない）を表す。方程式（１．２）の性質はある時間間隔にわたって積分された、物体とＧＰＳ衛星と間の相対速度の性質である。公知のＡＤＲ処理技術において、ＡＤＲ値が差分される二つの実時間ｔ_ｎ−１、ｔ_ｎは、二つの連続する基準時点に対応する。一方、本明細書で説明されるＡＤＲ処理技術においては、ＡＤＲ測定値の質は、これらの二つの実時間がより長い時間間隔によって分離される場合には、大いに改善され得る。より一般的な公式においては、これは、重み関数Ｗ（ｔ）がゼロでない、ｔのより広い範囲と対応する。測定式（１．２）を使用する簡単な例は、実質的な利点をもたらす。さらなる改善は、式（１．１）の形で測定値を使用することにより可能である。

（Ｂ．ＡＤＲ差計算のためのＡＤＲ測定値の間の時間間隔の選択）
方程式（１．２）は、ＡＤＲ差計算が実行される、二つの基準時点の間の時間間隔の期間を特定しない。この時間間隔の期間の選択は、以下に提案する方法の可能性に完全に到達するために、とても重要である。

上に示される質的考慮に基づく、二つのＡＤＲ測定値から導き出される速度推定値の基本的で概略的な形は、

である。

以下の説明は、カルマンフィルターで使用される測定式を特定する。この説明の目的のために、伝播モデルは、カルマンフィルターを公式化するために使用され、伝播モデルの誤差は、ＡＤＲ測定値の測定値ノイズとともに速度推定値の誤差に寄与する。よって、ＡＤＲ差を計算するためのＡＤＲ測定値の間の時間間隔を拡げることは、伝播モデルの精密さに依存する。より精密な伝播モデルは、ＡＤＲ差を計算するためのＡＤＲ測定値の間のよりも長い時間間隔を可能にする。同様に、伝播目的のためのより正確なＩＭＵデータおよびＩＭＵデータの使用は（図６と併用して後に説明されるように）、差計算のＡＤＲ測定値の間のより長い時間間隔を可能にする。

伝播モデルの誤差

は、モデル化されていない加速度のためとみなされ、以下のように推定される。

式中、δａは、モデル化されていない加速度である。二つのＡＤＲ測定値の測定値ノイズが考慮され、相関関係を有しないとみなされる。これらの分散は、等しいとみなされ、s^２ _Ｐにより表される。

速度推定値誤差の合計分散は、以下のように概算される。

時間間隔の期間に関してこの分散の最小を見つけるのは、簡単である。これは、本明細書で説明される技術を適用するのに使用される最適な期間を定義する。

伝播モデルの質が、モデル化されていない加速度が｜δａ｜≒１０^−６ｍ／ｓ^２の程度であり、ＡＤＲノイズが１０^−２ｍの程度である場合には、ＡＤＲ差の計算が実行される二つの基準時点の間の間隔の最適な期間は１００ｓの程度であり、各基準時点は約１秒である。

図４を参照にして、どのようにＡＤＲ差が計算されるかを例示するための例を示す。現在のＧＰＳナビゲーション技術が連続するＡＤＲ測定値の間の差を計算する一方で、ＡＤＲ測定値の間の時間間隔よりも大きい、特に実際に大きい、時間間隔によって分離されるＡＤＲ測定値の間の差の計算は、より正確なナビゲーションソリューションを提供し、カルマンフィルターのより速い収束を補助できることが解った。例えば、図４は、ＡＤＲ差は、基準時点Ｅ_１でのＡＤＲ測定値と基準時点Ｅ_１００（各基準時点が約１秒間隔である）でのＡＤＲ測定値との間、基準時点Ｅ_２でのＡＤＲ測定値と基準時点Ｅ_１０１でのＡＤＲ測定値との間、基準時点Ｅ_３でのＡＤＲ測定値と基準時点Ｅ_１０２でのＡＤＲ測定値との間、等で計算されることを示す。従って、この例において、ＡＤＲ差は、各新しいＡＤＲ測定値が受信されると、その新しいＡＤＲ測定値に対して計算される。

ＡＤＲ差は、多数の基準時点がそれに先行する以前のＡＤＲ測定値に関して、各新しいＡＤＲ測定値に対して計算される必要はない。ＡＤＲ差は、測定された各および全てのＡＤＲ測定値を必ずしも使用しない適切な時間間隔（例えば、基準時点の数）によって分離されるＡＤＲ測定値の間で計算され得る。例えば、ＡＤＲ差は、基準時点Ｅ_１でのＡＤＲ測定値と基準時点Ｅ_１００でのＡＤＲ測定値との間、基準時点Ｅ_１０１でのＡＤＲ測定値と基準時点Ｅ_２００でのＡＤＲ測定値との間、等で計算され得る。また、他の可能性としては、基準時点Ｅ_１と基準時点Ｅ_１００でのＡＤＲ測定値の間、基準時点Ｅ_１０と基準時点Ｅ_１１１でのＡＤＲ測定値の間、基準時点Ｅ_２と基準時点Ｅ_１２１でのＡＤＲ測定値の間などでＡＤＲ差を計算することである。これらのＡＤＲ差は、以下に説明する修正された偏導関数方程式を使用して処理される。

短時間（例えば１秒）または非常に長時間（例えば１０^４秒）の時間間隔の双方は、速度推定値において大規模な誤差を生じる恐れがあることを留意する。最適な期間に対して、速度推定値の誤差は、最適値を方程式（２）の期間代入することによって推定され得、

を得る。

誤差に対するこの推定値は、他の推定速度の方法と比較して小さい。これは、信号測定値の処理結果として、特に印象的である。もちろん、このような高精度の測定値を得る可能性は、同程度の精度をもつ測定式を使用することによってのみ実現され得る。以下のセクションは、上述のＡＤＲ差の処理において有用な測定式の導出を示す。

（Ｃ．物体クロック誤差の軽減）
以下の導出は、優れた伝播モデルが、一般的に衛星物体の場合である物体位置および速度（つまり、物体１０の位置および速度）に対して利用可能であると仮定する。残念なことに、クロックに対する伝播モデルは、（あまり高価でないクロックが使用される場合）あまり精度がよくない可能性があり、従って、方程式（１．２）の最後の項は、正確に推定するのが難しい。方程式（１．２）の右辺に対して偏微分が計算される場合には、クロックドリフトによるものが優勢となり、それは状態ベクトルの他の成分（主に速度）の正確な推定に対して問題を生じる可能性がある。

ＡＤＲ測定値のカルマンフィルター処理において、クロックドリフトの推定は、速度の推定から分離される。完全な分離は、実行不可能であり、また必要でもない。目的は、不安定なクロックの存在下で速度が推定できることである。方程式（１．２）の最後の項は視野にある全ての衛星に対して同じ値を有し、よって、測定値を差分することによって相殺され得ることが観察される。従って、カルマンフィルターで使用される衛星ペア（ｐ，ｑ）の測定値残余は、以下のようになる。

この残余は、利用可能なＡＤＲ測定値およびＧＰＳ衛星および物体の既知または推定される位置から容易に計算され得る。

式（１．７）の形での測定値残余は、速度を推定するために使用され得る。クロックドリフトに対しては、方程式（１．２）により与えられる形の単一測定値の任意の一つが使用される。

カルマンフィルターの測定値を使用するために、偏微分マトリックスを計算しなければならない。偏微分マトリックス計算の問題は、瞬間測定値モデルが必用とされるが、しかし方程式（１．２）がある時間の期間にわたって積分された速度に対するものであることである。この問題は、図５のフローチャートを併用して以下に説明する処理方法論により克服される。図５に示される、および以下に説明する処理ステップの目的は、状態ベクトルの速度成分を更新するために使用される瞬間速度推定値を提供するために、状態ベクトルの現在の時間での一つ以上の成分（例えば、状態ベクトル速度および位置成分）に関する、ＡＤＲ差の導関数を計算することである。ＰＲ測定値も、位置成分を更新するためにカルマンフィルターで処理され、既知の伝播モデルまたはＩＭＵデータ（図６と併用して後に説明する）は、状態ベクトルを伝播するために使用される。これらの技術は、後に明らかになるが、ＡＤＲ測定値の間の時間分離に関して、従来のＧＰＳナビゲーション技術よりも高精度を有する。例えば、誤差は〜χ^４で、このχは、従来の当該分野ナビゲーション技術の誤差〜ｘと比較して、測定値の間の時間差に比例する小パラメータである。

本方法は、現在の基準時点での状態ベクトルに関する、過去の基準時点での物体位置の偏導関数の計算を含む。本計算は、伝播モデルを根底にあるエンジンとして使用する、基準時点の間の時間差に関する高次近似（例えば、４次の）である。

（Ｄ．慣性座標系の位置および速度に関する偏微分）
方程式（１．２）は、スカラのため、任意の座標系を使用して公式化される。本サブセクションにおいて、慣性座標系が使用される。地球固定座標系（ＥＣＥＦ）における偏微分を計算するために、偏微分は、基準時点ｎにおけるＥＣＥＦ座標系と同じ軸方向を有する慣性座標系において最初に計算される。結果は、ＥＣＥＦ座標系で偏微分を計算するために後で使用される。

方程式をより扱いやすくするために、ペアの衛星の一つと対応する方程式（１．２）を使用して偏微分が得られる、つまり、

を計算する。

式中、

は、物体の位置、速度、クロックおよびクロックドリフトを含む物体状態ベクトルである。

よって、ステップ１００において、状態ベクトルは、ＥＣＥＦ座標系から慣性座標系に変換され得る。計算１００は、任意である。ＥＣＥＦ座標系におよびからの変換は任意である。

式（１．８）の偏導関数は、基準時点ｎにおける状態ベクトルの位置および速度に関してのみ計算される。状態ベクトルのこの区分は、

として表される。クロックおよびクロックドリフトに関する偏微分は、別に計算される。

物体位置

は、式（１．８）で明示的に使用されるために、それは物体位置に関する偏微分のみを計算する必要があるという印象を生む。しかし、物体位置

は、基準時点ｎ−１からの状態ベクトル伝播の結果である。従って、

は、

の関数である、よって、

は、

の（逆）関数であり、

（

の一部）は、

の関数である。これは、導関数のより複雑な計算を生じる。

式中、

およびχは、１に等しい主対角線の最初の３成分を除いて、全てゼロを含む６×６マトリックスである。方程式（１．９）の右辺は、測定値処理においてカルマンフィルターで使用される偏微分マトリックスの区分である。

方程式（１．９）において、マトリックスLが使用される。上で説明したように、本マトリックスは、状態ベクトルの伝播のためである、二つの基準時点における状態ベクトルの間の依存関係のためである。

（ＩＩ．さらなるＡＤＲ処理）
（Ａ．ＧＰＳのみのナビゲーション）
状態ベクトルを伝播するために、状態ベクトルの全成分が必要である。速度成分が不足している場合、伝播は行うことができない。マトリックスLは、基準時点ｎ−１で位置のみを、基準時点ｎで位置および速度の双方を使用することに留意する。これは、逆方向に、つまり基準時点ｎから基準時点ｎ−１に向かって、言い換えれば、差が計算されるＡＤＲ測定値を分離する時間間隔にわたって、伝播することを必要とする。

従って、ステップ１１０において、状態ベクトルは、ルンゲクッタ法などの伝播ルーチンを使用して、一ステップ戻って伝播される。伝播は、通常の微分方程式によって支配される。多くの伝播ルーチンは、高精度を確実にするために、いわゆるルンゲクッタ法を使用する。マトリックスLを計算する時、高精度が、根底にある伝播アルゴリズムの第４次ルンゲクッタ法を使用することにより、維持される。本伝播アルゴリズムは、状態ベクトルの伝播に使用されるアルゴリズムと正確にマッチする必要はない（いくつかの同様に優れたルンゲクッタ法がある）。有用なルンゲクッタ法の一例は、以下の一連の方程式により示される。

式中、ｈは伝播のステップである。一連の方程式（１．１１）は、標準形式（つまり、方程式が基準時点ｎから基準時点ｎ＋１に伝播されるとき）で書かれる。反対方向への伝播が必要であるために、一連の方程式（１．１１）でｈは（−ｈ）と置換される。

マトリックスLを計算するために、偏導関数

が計算される。

式中、δ_ｉｊは、クロネッカーのデルタ関数であり、

簡略化のために、以下の表記法が導入される。

一連の方程式（１．１４）を一連の方程式（１．１３）に代入し、偏導関数の計算により以下を得る。

一連の方程式（１．１５）の最初の方程式からの

を二番目の方程式に代入し、このような連鎖代入を最後の方程式まで続け、χ_ｎに関するｋの導関数が得られる。これは、マトリックスΨ_１−４を通して表示される。次に、これらの方程式の全てを方程式（１．１２）に代入し、以下を得る。

これは、６×６マトリックスの複数積を含む、面倒な公式である。幸いにも、Ψマトリックスの特別な形式により、かなりの簡略化が行われ得る。

方程式（１．１６）のχ_ｎ＋１、χ_ｎは、方程式（１．９）のｒ_ｎ-１、χ_ｎ，ｊと対応することに留意する。ベクトル

は、慣性座標系で示される。さらに、軌道における物体の運動は、第一義的に重力ポテンシャル場Ｕ（ｔ，χ）のためであると仮定する。これは、一連の（１．１１）の最初の方程式の関数

が、以下の形式を有することを意味する。

つまり、

の最初の３成分（位置の導関数に要因がある）は、単純に速度の成分であり、最後の３成分（速度の導関数に要因がある）は、位置のみの関数である。この事実は、方程式（１６）の右辺を簡潔化するために使用される。

の導関数（一連の方程式（１．１４）に必要とされる）は、以下である。

式中、

はゼロの３×３マトリックス、Ｉは、恒等３×３マトリックス、および

は、以下のように定義される３×３マトリックスである。

ここにΨマトリックスの積は、簡略化され得る。

式中、添字ａ、ｂ、ｃ、ｄは、値１から４の値を採る。方程式（１．２０）を方程式（１．１６）に代入して以下の式を得る。

式中、Φマトリックスは、３×３の次元を有し、以下のように定義される。

方程式（１．２１）のマトリックス

は、最高４つまでの６×６マトリックスの積を含む方程式（１．１６）の元の表式と比較すると、二つまでの３×３マトリックスの積を含む。マトリックスLの次元は、６×３であり、よって、４つの３×３マトリックスのうちの２つのみが方程式（１．２１）の右辺に必要であることに留意する。

前のセクションの一つで示された粗い推定値は、ＡＤＲ差計算に関連する２つの基準時点の間の時間間隔がかなり長い可能性があることを示す。偏導関数（１．２１）がｈに関する４次の精度で計算されても、結果の精度は、長い時間間隔に対しては充分ではないことがあり得る。精度を向上するための自然なソリューションは、いくつかのサブステップを採用することと、これらのサブステップに対応する、マトリックスの積として偏導関数（１．２１）を計算することである。

マトリックスは、方程式（１．２１）を計算するための基本である。これらは、Ｊ２重力ポテンシャルの場合に計算され得る。

式中、ｎ＝｛ｎ_ｉ｝は、ｚ方向における単位ベクトルであり、

であり、Ｒ_Ｅは、地球半径である。

伝播モデルは、しばしば、地球の重力場以外の作用を含む。これらの作用のいくつかは、宇宙飛行体の位置のみの関数であるが、（物体）宇宙飛行体の加速、例えば、太陽と月による、および太陽の輻射圧による加速を生じさせる。この場合、ＡＤＲ処理のための偏微分の導出は、上の導出と類似する。その他の作用（例えば、大気抵抗）において、加速度は、位置および速度の双方の結果である。この場合、方程式の簡略化を行わず完全な形式（１．１６）を使用するか、またはそのような作用を低オーダーの精度を伴う揺動として扱う。

従って、ステップ１２０において、

マトリックスは、方程式（１．２３）または以下に別に提供される、Ｊ２を超える高次の重力項に対する公式を使用して、各ルンゲクッタサブステップに対して計算される。

マトリックスが計算された後、ステップ１３０において、これらは、一連の方程式（１．２２）に代入され、方程式（１．２２）の結果は、マトリックスLを計算するために、一連の方程式（１．２１）で使用される。

マトリックスLは、偏微分方程式（１．９）に代入される。

（１．地上座標系における位置および速度に関する偏微分）
前のサブセクションは、慣性座標系における偏微分を表す。衛星の位置および高次の重力項は、ＥＣＥＦにおいて利用可能であるため、本座標系において偏微分を計算する必要がある。方程式（１．９）は、慣性座標系の物体の位置および速度に関して偏微分を計算するために直接適用できるため、ステップ１４０において、方程式（１．９）の計算結果は、位置および速度（対の測定値）にする偏微分を計算するための以下の変換方程式（１．２４）を使用して、ＥＣＥＦに変換される。

結果をＥＣＥＦに変換するために、以下の公式が使用される。

式中、

は、基準時点ｎにおけるＥＣＥＦ座標系の位置および速度ベクトルのｉ番目の成分であり、

は、慣性座標系の同様のベクトルのｊ番目の成分である。ヤコビアンマトリックスＪにおける要素

は明示的に計算され、その結果は、

であり、式中、

であり、

は地球角速度のベクトルであり、

は

から形成される３×３非対称マトリックスである。

偏微分は、最初に慣性座標系で計算されることに留意する。これは、方程式（１．１８）のΨマトリックスの特別な形式によるものであり、大幅に結果を簡略化し、最初に慣性座標系で偏微分を計算してから方程式（１．２４）を使用してＥＣＥＦに変換することをより経済的にする。しかし、これは必要条件ではなく、便利なためである。

（２．クロックおよびクロックドリフトに関する偏微分）
一部の宇宙飛行体の実装は、比較的低い安定性のクロックを使用する。従って、クロックの誤差、好ましくは、クロックドリフトの推定が必要である。ステップ１５０において、上で計算された位置および速度に関する偏微分は、ここで説明するクロックおよびクロックドリフトに関する偏微分により増加される。

これは、クロックドリフトの推定によって開始する。方程式（１．２）により与えられる形式のＡＤＲ測定値の一つが使用され、クロック誤差およびクロックドリフトに関するその導関数が得られる。ドリフトは、基準時点ｎ−１とｎとの間の時間間隔中一定であると仮定される。これは、以下の偏微分を生成する。

式中、ｈは基準時点の間の時間差であり、一番大きい期間である。従って、本測定値の処理は、クロックドリフトの値を基本的に更新する。本更新は、優れた伝播モデルがクロックにおいて利用可能でない場合（一つのモデルにおいて、ドリフトは、基準時点間において一定であると仮定される）には、あまり正確ではない。クロックドリフトの推定における不正確は、速度の推定において、この測定値をほぼ無益なものとする。速度を推定するために、差分測定値が必要である。

差分測定値から物体速度を推定するための偏微分計算が上に示されたが、しかしながら、全ての偏微分は計算されていない。つまり、物体クロックおよびクロックドリフトに関する差分測定値の偏微分が導出されていない。

方程式（１．８）の分子は、明示的にクロック誤差またはドリフトを含んでいないが、間接的な依存性がある。基準時点ｎにおけるクロック誤差は、τと仮定され、クロックドリフトは、

である。つまり、ｎ番目の基準時点測定値は、整数秒と−τだけ異なる時間において測定され、（ｎ−１）番目の基準時点測定値は、整数秒と

だけ異なる時間において測定される。τおよび

の任意の変化は、スナップショットの時間に影響するために、これらは、測定されたＡＤＲにも影響する。対応する偏導関数の計算から以下の公式を得る。

方程式（１．２７）の右辺において、公式はスカラおよび双方に対して有効なため、位置および速度の値が慣性またはＥＣＥＦ座標系のいずれかで使用される。

方程式（１．２７）の結果を使用して、偏微分マトリックスが作成され、ＡＤＲ差から瞬間速度推定値を提供するためのカルマンフィルター測定式を形成するためにステップ１６０で使用される。上述の偏微分方程式は、カルマンフィルターがＡＤＲ差から瞬間速度を提供するように公式化される。よって、ステップ１７０において、この瞬間速度は、ステップ１７０で示すようにカルマンフィルターの測定位相中、状態ベクトルの速度成分を更新するために使用される。

よって、カルマンフィルターは、対応する微分マトリックスを使用して、方程式（１．７）の形式で対の測距信号源（例えば、ＧＰＳ衛星）のＡＤＲ測定値残余を処理する。偏微分マトリックスは、２つの区分からなる。第一の区分は、（物体、例えば、宇宙飛行体）位置および速度のためで、各項が測距信号源（例えば、ＧＰＳ衛星）に対応し、方程式（１．９）の右辺により定義される、２つの項の差によって定義される。第二の区分は、（物体、例えば、宇宙飛行体）クロックおよびクロックドラフトのためで、方程式（１．２７）の右辺により定義される。宇宙飛行体の状態ベクトルが追加成分（例えば、クロック加速度、太陽輻射反射係数等）を含む場合には、対応する偏微分が同様の方法で導出され使用される。

前段は、ＡＤＲデータ処理の測定式の導出を示す。ＰＲ処理の測定式は、当該分野において周知であり、本明細書では説明されていない。カルマンフィルターはまたステップ１７０において、状態ベクトルの位置成分を更新するために、ＰＲ測定値を処理する。

ＡＤＲ差を計算するために使用されるＡＤＲ測定値の間の時間間隔は、ＧＰＳナビゲーション用途の現在のＡＤＲ差分手法において使用される時間間隔よりも非常に大きいために、上述のＡＤＲ測定式は、よりいっそう精密な方法で公式化される。本手法の利点を明確にするためのアナロジーを以下に記す。

オフィスビルなどの大きな建物の寸法（長さ、幅または高さ）測定の作業を想定する。これらの測定値を測定する方法の一つは、比較的短い物差しを使用する。物差しが短いため、建物の一つの寸法の測定を終えるのに、多数の測定値が物差しを用いて必要とされる。測定値誤差は、短い物差しで測定された複数の測定値に累積される。一方、これらの測定値を測定するもう一つの方法は、非常に長い測定テープを使用することである。事実、測定テープが十分長ければ、たった１回の測定が建物の一寸法に必要とされ、よりいっそう正確な測定値が得られる。従来のＡＤＲ ＧＰＳナビゲーション技術は、非常に長い距離を測定するために短い物差しを使用するのに似ている。各ＡＤＲ測定の間の時間間隔は短く（例えば１秒）、したがって、「短い物差し」問題に類似するため、測定値誤差は各ＡＤＲ差計算に累積する。このようなＡＤＲ測定値誤差は、その累積作用のため、ナビゲーション計算の全体的な精度を低下させる。

カルマンフィルターで使用される伝播モデルが改善され、いっそう精密になる一方で、ＡＤＲデータに使用される方程式に対しては、これは該当しない。上で示されたＡＤＲ測定式は、数多くの既知の伝播モデルにより達成された精度と同程度の、測定値の精度を達成する。よって、新しい測定値データを用いて状態ベクトルを更新（例えば、伝播）する時、ＡＤＲ測定式は、伝播式に高精度な測定値を提供し、それは全体としてナビゲーション計算の精度をいっそう向上する。

別のアナロジーは、ストップウォッチを開始して、車で移動し、ある時間間隔後、ストップウォッチを停止し、車が移動した距離を測定することである。ストップウォッチが開始および停止するごとに、誤差が生じる。短い時間を使用する場合、推定速度は、非常に正確な可能性がある。逆に、車が一定の速度で移動する場合、および長期間測定を実行する場合、推定速度は、よりいっそう正確になる。二番目のオプションが好まれるが、一定速度の仮定が有効でない可能性があるため、常に実行可能ではない。この場合、長期間の速度の測定は、この移動に対しての平均速度のみを提供する。多くの用途は現在の状態の情報を必要とするため、平均速度はほとんど実用性がない。

これらの技術は、比較的長い時間間隔に対して物体の速度を測定することの利点を、そのような測定値を処理時間における瞬間速度に正確にマップする能力と組み合わせる。これは、差分ＡＤＲ測定値を処理することと、このような差を現在の利用者の状態ベクトルに正確にマップする数式を適用することによって達成される。

本発明の方法を説明する方法は、以下である。車の速度アナロジーにおいて、移動した距離（当面の問題におけるＡＤＲ差に類似する）は、方向的に車の速度に比例する。よって、車の速度に関する距離の導関数は、測定値の時間範囲に等しい。カルマンフィルターは、測定値そのもの（車のアナロジーにおいて、移動した距離）および未知数に関するその偏導関数（車のアナロジーにおいて、この導関数は、測定値の時間範囲である）を使用して、速度を形式的に推定する。本発明において、高精度の方法が、カルマンフィルターにより必要とされる偏導関数を計算するために考案され、各測定値の長い時間範囲と組み合わされて、物体に対する高精度な速度推定値を提供する。

（３．伝播）
カルマンフィルターにより使用される伝播式を以下に述べる。

状態ベクトルは強制モデルおよび積分ルーチン（例えば、ルンゲクッタ法）を使用して伝播され得る。また、これは当該分野において周知であり、本明細書では説明されない。共分散マトリックスの伝播は、より困難な課題である。出発点は、基準時点ｎにおける状態ベクトルに関する、基準時点（ｎ＋１）における状態ベクトルの偏導関数である。これらの偏導関数が計算される場合には、共分散マトリックスは、標準式を使用して伝播され得る。必要とされる導関数は、

である。

状態ベクトルは、位置、速度およびクロック／ドリフトの区分を有する。よって、マトリックス

は、速度に関する位置の、クロックに関する位置の、位置に関する速度の等の導関数の区分を有する。これらの区分は、以下で別に考慮される。

（ａ．クロック成分を含まない導関数）
方程式（１．１０）との比較は、本マトリックスがマトリックスLと類似しており、Lの計算は、

の計算を補助すると考えられることを示す。しかし、方程式（１．２８）の状態ベクトルの成分はＥＣＥＦで表示されており、これは、単にLを適用すること以上のことを必要とする。

方程式（１．２９）において、方程式（１．２５）に定義されるように

、および方程式（１．２１）により与えられる

の表示を伴う

、を認識する。上付き「＋」は、方程式（１．１０）のL_ｋｍの定義における逆方向伝播に対して、この場合、位置および速度が順方向に伝播されることを示す。（これらの２つの場合は、方程式（１．２１）のｈに対する正および負の値の使用に対応する。）
残された、計算されるべき量は

である。この量において、分子および分母の双方は、同様のタイムスタンプを有し、伝播は含まれない。全ての導関数は、座標系の間で変換される。

方程式（１．２９）において、慣性系は、基準時点ｎにおいてＥＣＥＦと一致する軸を有すると仮定される。それは基準時点（ｎ＋１）において同じ方向を保つ。従って、

から

への変換は、基準時点ｎにおける慣性座標系からＥＣＥＦ座標系への変換と、次に基準時点ｎにおけるＥＣＥＦ座標系から基準時点（ｎ＋１）におけるＥＣＥＦへの変換の２つの操作を含むとみなされる。この変換は以下のように表示される。

式中、

であり、Ｍは基準時点ｎにおけるＥＣＥＦから基準時点（ｎ＋１）におけるＥＣＥＦへの３×３回転マトリックスである。Ｊ^−１は、

を

に変化することによりＪから形成され得る。全マトリックスを（１．２９）に代入することにより、クロック成分を含まない、

区分に対する以下の表示を得る。

（ｂ．クロック／ドリフトに関する位置／速度の導関数）
本サブセクションは、導関数

および

の計算を説明する。これらの導関数を計算するとき、基準時点ｎ＋１における位置／速度は、基準時点ｎにおける位置／速度に依存し、後者は、同様に、クロック誤差およびクロックドリフトに依存することに留意する。従って、基準時点ｎ＋１における位置／速度は、基準時点ｎにおける位置／速度を通して、基準時点ｎにおけるクロック／ドリフトに間接的に依存する。しかし、偏導関数の計算が必要とされ、この間接的な依存は、除外されなければならない。結果は以下である。

保存的方程式（ＥＣＥＦ座標系における地球の重力場の衛星力学の場合のような）に対しては、クロック誤差への依存性はない（すなわち、導関数はゼロである）が、クロックドリフトに対するわずかな依存性がある。それは以下の公式により提供される。

式中、

は、ＥＣＥＦ座標系の基準時点ｎ＋１における物体加速度である（加速度は強制モデルを使用して計算される）。

（ｃ．他の成分に関するクロック／ドリフトの導関数）
これらは、あまり重要ではなく、以下によって与えられる。

（４．高次の重力項）
本セクションは、高次の重力項におけるマトリックス

の計算を示す。高次の項の重力ポテンシャルは、一般的に、極座標（Ｍ．Ｋａｐｌａｎ，Ｍｏｄｅｒｎ Ｓｐａｃｅｃｒａｆｔ Ｄｙｎａｍｉｃｓ＆Ｃｏｎｔｒｏｌ，Ｊｏｈｎ Ｗｉｌｅｙ＆Ｓｏｎｓ，１９７６を参照）において与えられる。

公式（１．３４）は、以下のように、Ｋａｐｌａｎの公式（７．１６）と異なる。それは緯度にθ（極角ではないことに留意する）、経度にφ、および−Ｊ_ｋの位置に

を使用する。

極座標に関するＵの一次および二次の偏導関数は、

である。

デカルト座標に関する偏導関数は、以下のように計算される。

式中、ベクトル

は、

のように定義される。（Ａ３）は、デカルト座標に関する極座標の一次および二次の導関数を含むことに留意する。これらの導関数は、

により与えられる。

要約すれば、ＡＤＲ差データを処理するためにカルマンフィルターにおいて使用される偏微分は、以下の方法論を使用して計算される。
１．状態ベクトルをＥＣＥＦから、瞬間的にＥＣＥＦと整合される（回転ではない）慣性座標系に変換する。
２．第４次のルンゲクッタ法を使用して、状態ベクトルを一ステップ戻して伝播する。
３．各サブステップにおいて、方程式（１．２３）または高次重力項の公式を使用して、

要素マトリックスを計算する。
４．これらのマトリックスを一連の（１．２２）に代入し、結果を一連の（１．２１）に代入する。これはLを計算する。
５．Lを方程式（１．９）に代入し、次に、位置および速度に関する偏微分を計算するための変換式（１．２４）に適用する。これらを対の測定値のために使用する。
６．方程式（１．２６）を単一の測定値に対して使用し、方程式（１．２７）を対に対して使用して、クロックおよびクロックドリフト成分に関する偏微分を増加させる。

（Ｂ．ＧＰＳ−ＩＭＵ統合型ナビゲーションのためのＩＭＵデータを使用して状態ベクトルを伝播する）
以下は、統合型ＧＰＳ−ＩＭＵシステムのＡＤＲデータ処理のためのナビゲーションアルゴリズムを説明する。ＩＭＵは、物体の位置または速度のいかなる絶対的な測定値を提供しない、よって、主としてナビゲーションの補助として使用され得る。また、位置および速度の実際の推定は、ＧＰＳ測定値の処理結果である。しかし、ＩＭＵデータは、特に、精密な伝播モデルが対象物体に対してまだ利用不可能な場合には、高ナビゲーション機能を得るために重要な、精密な伝播手段として重要である。

ＩＭＵ測定値の重要性を認識するために、これらの役割を宇宙用の伝播モデルの役割と比較することは有用である。宇宙飛行体ナビゲーション用途における優れた精度を得るために、精巧な伝播モデルが不可欠である。その役割は、長期間にわたり収集された測定値も間の一貫性を確実にすることである。優れた伝播モデルは、処理ノイズを相応的に低レベルに設定することが可能であり、それは、最近の記録からの各過去の測定値を、任意の現在の測定値とほぼ等しい現在の推定値に対して暗黙的に寄与させる。結果は、任意の与えられた時間において、ソリューションを得るためのさらに多数の測定値があるかのようである。（もちろん、過去の測定値は現在の時間に直接には適用されないが、以前の推定値は実時間に伝播され、ソリューションにおいて大きく重み付けされる。）例えば、伝播モデルが、１０００秒間「有効である」場合、一基準時点で測定されるよりも約１０００倍の測定値の使用を助長し、これは、

倍の割合でナビゲーション誤差を減少させる。

別の重要な考慮すべきことがある。伝播モデルの誤差は、測定値のランダムな誤差とは逆に、一般的に、時間に高い相関（つまり系統的）を有する。系統的な誤差は、フィルター平均化によって除去されず、小さな系統的誤差は、より大きなランダム誤差よりもより有害となる可能性がある。これは、優れた伝播モデルの必要性をさらに強調する。伝播モデルの誤差源は探索され、軽減される必要がある。これを行うに際して、数学的モデルのコンピュータ実装が追加の誤差、しばしば系統的な誤差を導入する可能性があることに留意する。例えば、宇宙飛行体伝播における運動の微分方程式の積分は、ルンゲクッタ法などの高次の方法を使用すべきである。低次の方法（例えば、オイラー）は、離散型の微分方程式表示の中の誤差によるモデル化されていない加速度を有することと類似している。

良好に導出された伝播モデルにおいて興味のある用途は、ＡＤＲ測定値を介する速度の推定方法を統合型ＧＰＳ／ＩＭＵナビゲーション処理に拡大することである。より優れた伝播モデルが統合型ＧＰＳ／ＩＭＵナビゲーションシステムに対して開発された場合には、ＡＤＲ処理をこのような統合型ＧＰＳ／ＩＭＵナビゲーションシステムに拡大することが可能となる。

統合型ＧＰＳ／ＩＭＵシステムにおいて、システムのＩＭＵ部分は、ある程度、ナビゲーションソリューションの伝播手段として役立つ。以下は、ＩＭＵデータを使用してソリューションを伝播するための高精度のモデルを示す。

（１．観測可能な量および状態ベクトル）
ＩＭＵは、加速度および回転速度を測定する。これらの量は、速度および姿勢パラメータの時間導関数であり、状態ベクトルに含まれ得る。しかし、これに関係する著しい処理コストがあり得る。ＩＭＵは、基本的に１００Ｈｚから１０００Ｈｚでデータをアウトプットし、このような高速度でフィルター処理を実行するのは過剰である。代わりに、ＩＭＵによりアウトプットされた加速度および回転速度データが、伝播ツールとして厳密に使用される。

ＩＭＵデータが伝播ツールとして使用されるときでさえ、姿勢、回転速度、加速度、速度および位置の間でのカップリングが存在する。このカップリングは、共分散マトリックスの非対角要素のためであり、後に詳細に検討される。これは、位置測定値から速度を推定することに類似する。位置／範囲が測定量であっても、速度は、共分散マトリックスの非対角要素を介して位置とカップリングされるために、推定され得る。

本アプローチは、状態ベクトルの加速度および回転速度の使用を避けることを可能にする。状態ベクトルに含まれる候補パラメータは、以下である。
１．位置の３成分
２．速度の３成分
３．物体のクロックおよびクロック速度
４．姿勢
５．バイアス、不整合、尺度係数および他の系統的ＩＭＵ誤差
状態ベクトルにおける姿勢情報の特定の表示は後に説明される。ここにおいては、姿勢は３つのパラメータによって完全に説明されることが注目される。

（２．伝播に対するＩＭＵ測定誤差の影響）
本サブセクションは、ランダムＩＭＵ誤差源を説明する。

伝播は、ステップサイズｈ（例えば、ｈ＝０．０１秒）において利用可能な加速度測定値ａ_ｉを使用して、ある期間Ｔ（例えば、１秒）にわたり考慮される。基準時点ごとに、合計で

個の加速度測定値がある。

加速測定値が誤差

を有する場合、結果は、

として推定される速度の誤差である。この速度の誤差は、基準時点の終わりまで維持される時間のあいだ、つまり、ｔ_ｉ＝ｉｈである時間間隔（Ｔ−ｔ_ｉ）のあいだ、物体位置の伝播に寄与する。この特定の測定値からの位置において得られる誤差は、次に、以下のように推定される。

この誤差の分散は、方程式（２．１）を二乗することにより推定され得る。全サブ基準時点からの分散は、基準時点の終わりに合計誤差分散を形成するために加算される。

合計の計算により

を得る。

同様に、速度に対する誤差分散の計算により以下を得る。

加速度計誤差は、物体姿勢における誤差に寄与しない。

ジャイロ測定値誤差は、物体位置の誤差、速度の誤差および姿勢の誤差に寄与する。回転の速度

が、誤差

を伴って各基準時点で測定される場合には、姿勢、位置および速度の誤差は、

により与えられる。

式中、ａは、この基準時点のあいだの物体加速度（時間により変動するが、誤差分散計算の目的のためには、基準時点の任意の、例えば最初の、ポイントにおける加速度が使用され得る）であり、位置および速度の誤差は、加速度に垂直方向においてゼロでない。

上述のように、これはランダム誤差であり、処理ノイズを伝播位置、速度および姿勢に適用することによって説明される。バイアス、不整合、尺度係数および他の系統的誤差は別に扱われる。

（３．測定値の更新）
測定値更新方程式は、ＰＲ ＧＰＳ測定値および任意でＡＤＲ測定値を処理する。ＰＲ測定値更新方程式は、当該分野において周知であり、よって、本明細書において説明されない。ＡＤＲ測定値更新方程式は、ＧＰＳのみのナビゲーションに対して上で示したものと同様である。測定値マトリックスのゼロでない成分のみが、位置およびクロックバイアスに対応するゼロでない成分である。

ＡＤＲ測定値において、差分また適応方法が速度の推定のために使用される。本明細書で説明される方法において、ＡＤＲ測定値は、バイアスを相殺するために時間において差分される。可能性を最大限に達成するために、ＡＤＲ差分が、ＧＰＳのみのナビゲーションで上述したように、多数の基準時間を含む時間間隔に対して実行される。

カルマンフィルターによるように時間差分ＡＤＲ測定値を処理するために、状態ベクトル成分に関する測定値の偏導関数を含む、観測マトリックス（Ｈ）が計算される。

式中、

は、位置、速度、姿勢およびＩＭＵ誤差状態を含む物体状態ベクトルである。

方程式（８）および（９）に関連する上述のように、物体位置

は、方程式（２．８）において明示的に使用されるため、それに関する偏微分を計算しなければならないという印象を受ける。しかしながら、物体状態

は、基準時点ｎ−Ｍからの状態ベクトル伝播の結果である。従って、

は、

の関数である。よって、

は、

の（逆）関数であり、

（

の一部）は、

の関数である。これは、導関数のより複雑な計算を生じる。

式中、

であり、χは、１と等しい主対角線の最初の３成分を除き、全てゼロを含む。

（４．伝播）
（ａ．順方向伝播）
順方向伝播アルゴリズムは、状態ベクトルおよび共分散マトリックスの双方に必要である。図６を参照に、ＰＲ測定値およびＡＤＲ差から計算された状態ベクトルを伝播するためのＩＭＵデータを使用する方法が説明される。ステップ２００において、観測マトリックスＨは、方程式２．１０で示されるように、状態ベクトルから計算される。

（ｉ．姿勢の伝播）
次に、ステップ２１０において、姿勢成分は、ベクトルをＩＭＵ体軸からＥＣＩ座標系に変換する回転マトリックスＣ（ｔ）の伝播によって伝播される。目的は、ジャイロセンサーデータを使用する本マトリックスの伝播アルゴリズムを導出することである。

各サブ基準時点において、マトリックスＣ（ｔ）は、わずかな変化を受ける。この変化は、このサブ基準時点中のＩＭＵの回転によるものであり、変化の速度は、体座標系のジャイロにより測定される。従って、マトリックスＣ（ｔ＋ｄｔ）が、時間（ｔ＋ｄｔ）においてベクトルを体座標系からＥＣＩに変換することを必要とする場合、それは２つの連続的ステップを適用することにより行われる。
１．時間（ｔ＋ｄｔ）における体座標系から時間ｔにおける体座標系に回転する。
２．マトリックスＣ（ｔ）を適用することによって、時間ｔにおける体座標系からＥＣＩに回転する。

最初の回転は小さく、よって、全ての３つの角度の回転は可換である。それは時間ｔにおける体座標系から時間（ｔ＋ｄｔ）における体座標系への回転の逆であり、これは、

の形式であり、
式中、マトリックス

は、ジャイロにより測定された回転速度の３成分ベクトルから形成された歪対称（ｓｋｅｗ−ｓｙｍｍｅｔｒｉｃ）マトリックスである。本マトリックスの逆マトリックスは、単に同様のマトリックスだが、

後にマイナス表示をもつ。よって、マトリックスＣ（ｔ）は、

により記述される。

本公式は、伝播が極短時間間隔ｄｔに対して必要とされる場合には、非常に優れている。しかしながら、１基準時点（例えば、１秒）わたって公式を伝播するために修正される。単純な一般化により、

を得る。

式中、マトリックス乗算の順序は、後のサブ基準時点が積の一番右の乗数に対応する順序である。（この表式は、ボルテラ乗法積分と密接に関係する。）本伝播スキームは、各ステップにおいて一次の精度（すなわち、ｏ（ｄｔ^２）オーダーの誤差）と、全基準時点に対する伝播後のｏ（ｄｔ）オーダーの誤差を有する。

高精度が必要とされない用途において、この単純な一般化は、十分であり得るが、さらなる改善の余地がある。第一に、二次の伝播が各ステップで求められる。

（回転マトリックスに対する二次の伝播スキーム）
第一ステップは、伝播スキームを中心化することである。（実際、乗法積分よりも加算積分が伝播される場合には、これが必要とされることである。）簡潔に、下付き添字０は、時間ｔを意味し、下付き添字１は時間（ｔ＋ｄｔ）を意味する。そして、伝播スキームは、以下の形式で求められる。

式中、Ｃ_１は、時間（ｔ＋ｄｔ）における回転マトリックスの正確な値である。正確な伝播とこの公式との間の全差は、ｏ（ｄｔ^２）項によって説明される。

はじめに、２つの連続時間間隔に対する伝播を考える。これを行うのに、２ｄｔ長の１ステップと、ｄｔそれぞれの２ステップの２つの方法がある。厳密な回転マトリックスが計算されるので、結果はマッチしなければならない。

最初の場合、誤差ｏ（ｄｔ^２）は、単一ｄｔ式と比較して、４倍に増加することに留意する。これは、ｏ（ｄｔ^２）がｄｔに対して二次従属性を有し、ｄｔはここで２倍大きいためである。二番目の場合、誤差は、単一ｄｔ式と比較して、２倍に増加する。これは、同じ誤差が伝播を２回適用することにより２回得られるためである。

に対してテイラー展開を使用した後、（２．１５）の二番目の方程式が、以下のように近似される。

ここで、方程式（２．１６）の積が計算され得、結果が（２．１５）の最初の方程式と等値される。その結果（ｄｔに関する２次の精度に対して）は、以下の通りである

この推定は、最終結果のために方程式（２．１４）に代入され得、ｄｔに関する２次の精度を有する。

（回転マトリックスに対する３次の伝播スキーム）
同様の技術が、ｄｔに関する３次に有効な式を得るために繰り返され得る。２ｄｔ長の１ステップ後と、ｄｔそれぞれの２ステップ後の回転マトリックスの方程式は、以下である。

同様の代数計算の後、

を得る。

従って、回転マトリックスを伝播するための３次のスキームは、以下である。

式中、

である。

（ｉｉ．位置および速度の伝播）
ステップ２２０において、上述の伝播アルゴリズムは、位置および速度に拡大される。伝播は、高次積分スキームを使用して、姿勢マトリックスを含む式の時間に対する積分を含む。

積分は、有限の和として示され、各ステップで２次の精度に有効であるスキームが導出される。これは、姿勢マトリックスを各ステップに対して３次まで計算することを必要とする。実際に、現在のステップ以前の数多くのステップは、ステップサイズに反比例する。よって、姿勢マトリックスにおける誤差は、ステップサイズに関して２次になるまで累積し得、それはスキームの残りと一致する。しかし、この考えは、２次のスキームが適用され得る、位置の公式を導く目的のためのステップサイズ内の姿勢マトリックスの近似に対しては適用されない。

以下は、位置伝播の導出である。速度に対する式は、多少少ない手間により、同様に導出され得る。

基準時点（Ｎ個のサブ基準時点を含む）の終わりにおける位置は、

により提供される。

式中、

および

は、基準時点の初めにおける位置および速度であり、

は、加速時計により測定される体座標における加速度である。部分積分が、上の式を単一積分に簡略化する。

目的は、この積分を計算することである。台形法は、有限の和により積分を近似するために使用され得る。各ステップは、以下の寄与を行う。

これは、厳密な公式であることに留意する。台形法の全ての誤差は、ｏ（ｄｔ^２）に組み入れられる。同様の式が速度伝播に対して導かれ得るが、乗数（Ｔ−ｔ_ｎ−１）および（Ｔ−ｔ_ｎ）がないことのみが（２．２５）と異なる。

は、３次のスキーム（ｎ−１のサブ基準時点に対して第２レベル誤差を生成する可能性がある）を使用して、時間ｔ_ｎ−１にすでに伝播されると仮定される。この導出の目的のみのために、時間ｔ_ｎにおける姿勢マトリックスへの更新が、２次のスキームを使用して時間ｔ_ｎ−１における姿勢マトリックスから推定される。

残りの導出は、以下である。方程式（２．２６）は、方程式（２．２５）に代入され、公式の結果は、１ステップに対する単一伝播または２連続の半ステップ伝播の２つの場合に適用される。結果は、（式（２．２５）は厳密なため）同じはずである。代数計算の後、これはｏ（ｄｔ^２）に対する方程式を生成し、以下が得られる。

式中、高次の項における全パラメータの値は、（ｔ_ｎ−１，ｔ_ｎ）の間隔内（スキームは、まだ次数を保持する）の任意の時間点において計算され、

は、先に導出される。位置における全体の増加は、式（２．２７）の右辺のｎに関する合計である。

速度伝播に対する同様の方程式は、以下である。

新しい表記法が、簡略化のために導入される。

式中、

である。

（ｉｉｉ．不整合）
今までのところ、ＩＭＵ誤差の取り扱いは、ジャイロセンサーにおける不整合を考慮することに限られていた。括弧内に示されるように、これは、モデルにとって扱いにくい誤差である。ステップ２３０において、以下に述べるように、ベクトル

は、ＩＭＵセンサーにおける不整合およびバイアスをカルマンフィルター状態ベクトルに特徴づけるように定義され、対応する導関数が計算される。

ＩＭＵジャイロセンサーが不整合の場合、測定された回転速度は、体座標系の実際の回転速度と異なる。この作用を説明するモデルは、以下である。

式中、

は体座標系の回転速度の実際の回転ベクトル、

は測定された（アウトプット）ベクトル、

は、

を

に回転するマトリックスである。以前のように、マトリックス

は、

によって定義される。マトリックス

は、回転性であり、恒等マトリックスに近く、よって、

の形式で示すことができる。
式中、

であり、ベクトル

は、３つの不整合パラメータを特徴づける。これらのパラメータは、一般的に、ナビゲーション用途の初めには不明であり、よって推定されるべきである。

このモデルをナビゲーションアルゴリズムに組み入れるために、以下が実行される。
１．カルマンフィルター状態ベクトルにおいてベクトル

を含み、対応する偏導関数を計算する。
２．伝播式の

を使用する前に、公式（２．３１）および（２．３２）を適用して修正されたジャイロ速度により状態ベクトルを伝播する。
３．不整合パラメータを推定する。これらは、一定値であると仮定される。

項目（３）に達するために、不整合パラメータに関する位置および速度の偏導関数を含む、伝播マトリックスが計算される。また、同様の計算が、これらのパラメータに関する測定値の偏導関数を計算するために適用される。

偏導関数の計算は、以下である。定義は、

であり、

に対する一般式が導かれ、任意の

、

等が付け加えられる。次に、これらの構成単位は、所望の偏導関数を探すために使用される。

の導出は後に示され、それは以下の結果をもたらす。

公式（２．３４）の右辺は、不整合に対して修正されていない（よって、下付き添字ｔをもたない）回転速度

を含むことに留意する。これは、公式（２．３４）の導出が明示的に不整合の修正を取扱うためである。

姿勢マトリックスの導関数に対して、以下の回帰公式が提供される（以下の全ての

およびそれらの導関数が、不整合に対してすでに修正された値を表す）。

または、同等に、

式中、

位置および速度の導関数に対して、

式中、

である。

（ｉｖ．順方向伝播マトリックス）
次のステップ２４０は、伝播マトリックスを導くためである。状態ベクトルの偏導関数は、基準時点の初めにおける状態ベクトルに関して基準時点の終わりに計算される。

状態ベクトルの姿勢の表示における選択は、もはや延期できない。姿勢マトリックスは、９つの要素を有するために適切ではないが、たった３つの独立のパラメータにより定義される。オイラー角は、特定の値において縮退を示す不利を有する。最後に、４元数表示は、二つの不利を有する。それは、４つの要素（３つのみが独立である）を含み、各ステップで姿勢マトリックスの計算を必要とする。２つの４元数回転の後の用途は重要である。これは、姿勢表示のための以下の選択へと導く。各基準時点において、小角度の回転は、現在の姿勢マトリックスに関して３つの小オイラー角を用いて定義される。よって、姿勢マトリックスは、

の形式で示される。式中、

は、姿勢の現在の推定であり、

は、

の形式における小回転である。

角度

は、小さく、対応する回転は可換である。これらの角度は、状態ベクトルの姿勢区分を形成する。

このアプローチの綿密な調査は、２つの問題の可能性を明らかにする。
１．ある一定時間後、誤差の累積は、

がもはや小さくない結果を生じる。この問題におけるソリューションは、姿勢マトリックスおよび

値を周期的に（例えば、基準時点ごとに）再調整することである。各基準時点で、

に対する回転

が適用され、よって、現在の姿勢マトリックスを更新し、

をゼロにリセットする。対応する変換が、

の共分散にも適用されるべきである。これは、新しい、更新された推定位置に関する各ステップにおける疑似距離測定の線形化に類似する。唯一の違いは、マトリックスの「線形化」が、加算的よりもむしろ乗算的であることである。
２．小回転の連続した適用は、もはや純粋に回転性でないマトリックス（例えば、正規直交条件

が、不正確になる）を生成し得る。この問題のソリューションは、周期的に姿勢マトリックスを正規直交形式に再投入することである。

偏微分が慣性座標形における１基準時点に対する伝播のために導出される。

（補助計算：

）

不整合が存在する場合には、最終姿勢マトリックスは公称値と異なる。

一方、姿勢におけるこの変化は、パラメータ

によりモデル化される。

（２．４４）および（２．４５）の比較により、

を得る。（（２．４６）の右辺は歪対称マトリックスである。これは

から得られる。）
位置、速度および不整合を含む状態ベクトルに対して、偏導関数は以下に提供される。

全ての他の成分はゼロである。クロックおよびドリフトに関する導関数は示されていないが、それらは、ＧＰＳのみをベースのナビゲーションに対する公式と変わらない。

これらの公式は、扱いにくく見えるが、計算をより効率的にする多くの機会がある。例えば、３×３マトリックス

の複数の積は、簡略化され、ある計算は再利される、などである。

本伝播マトリックスは、共分散を伝播するために使用し得る。状態ベクトルは、位置に対する方程式（２．２９）および（２．３０）を使用して伝播されなければならないが、姿勢は方程式（２．２１〜２．２２）を使用して伝播される。

（ｂ．逆方向伝播）
次に、ステップ２５０において、逆方向伝播アルゴリズムが、複数基準時点による状態ベクトルを逆方向伝播マトリックスを用いて逆方向に伝播するためと、現在の状態の成分に関する過去の状態の偏導関数を計算するために必要である。逆方向伝播は、ＡＤＲ処理において使用される。

状態および共分散マトリックスが逆方向に伝播されても、これは時間の逆方向にステップすることによって計算されるとは限らない。これは、逆方向伝播の間、ＩＭＵデータを積み重ねて保管することを必要とする。伝播時間が数十秒で、ＩＭＵ速度が数百ヘルツの場合、保管および処理の要求が過剰になる可能性がある。よって、順方向に行われるような、逆伝播ルーチンが考案される。言い換えれば、逆方向伝播は、順方向伝播に基づく（および、順方向伝播から導かれる）べきである。この導出は本サブセクションの主題である。逆方向伝播は、順方向伝播式（２．４２）に基づき、ただし

および

、

の値が、複数基準時点伝播間隔の終わりにおける状態の関数として計算される。

各基準時点において、（２．４２）の

および

、

の値は、複数基準時点逆方向伝播間隔の終わりにおいて有効である、不整合パラメータ

の値を使用して計算される。残念なことに、前の基準時点中、この未来値は不明である。

のこのような未来値を取扱うために、固定された

（例えば、

）における

および

、

の値と

に関する同様の導関数を推定するべきである。

の値が伝播間隔の終わりに利用可能になるとき、

および

、

の値は、

における値と

に関する導関数を使用して、二項テイラー展開により近似される。

（ｉ．逆方向伝播マトリックス）
ＡＤＲ測定値を最適に使用するために、逆方向伝播マトリックスが複数の基準時点に対して計算される。よって、このような量

が計算され、ここでＭは逆方向に進む基準時点の数である。

Ｍ−１からＭに行くために、修正された方程式（２．４２）が開始点として使用され、

基準時点Ｍにおける状態ベクトルの成分に関するこれらの公式の偏導関数を計算する。状態ベクトルは、４つの区分からなる。

であり、ここで恒等式

である。姿勢マトリックス

は、基準時点Ｍ−１においてタイムスタンプされる、よって、導関数計算の目的のために、基準時点ＭからＭ−１への姿勢マトリックスの逆方向伝播の結果を考慮するべきである。逆方向伝播マトリックスの公式を得る前に、２つの補助公式が以下で導かれる。

（ｉｉ．補助計算：姿勢マトリックスの導関数）
姿勢マトリックスに対して、以下が使用される。

式中、

は、（２．２１）により提供される（（２．２１）は、１サブ基準時点に対する伝播を記述し、１基準時点に対する伝播を記述する

を得るために、複数回繰り返されなければならないことに留意する）。これにより、

および

を得る。

（ｉｉｉ．補助計算：

）
これは、以前の

の計算に類似するが、ここで、導関数は、

を使用して逆方向伝播に対して計算される。

不整合が存在する場合には、最終姿勢マトリックスは、公称値とは異なる。

一方、姿勢におけるこの変化は、パラメータ

によりモデル化される。

（２．５３）および（２．５４）の比較により、

を得る。

ここで、逆伝播マトリックスに対する導関数は、以下のように書かれる。

全ての他の成分はゼロである。

これらの公式を用いて、

に対する回帰公式が書かれる。

（ｃ．状態ベクトルの逆方向伝播）
回帰アルゴリズムは

のように定義される。式中、

であり、

である。積

は、個別のマトリックス

の

積以外のなにものでもない。

に対する（２．５８）の最初の方程式を解き、結果を（２．４８）の最初の方程式に代入することにより、以下を得る。

Ｍに関する（２．５９）の総和は、

に対する所望の逆方向伝播式を提供する。

同様に、

に対する（２．５８）の２番目の方程式を解き、結果を（２．４８）の２番目の方程式に代入する。

次に方程式（２．６０）を使用し、下付き添字をＭからＭ−１に変更し、得られた方程式を

について解き、（２．６１）に代入すると、

を得る。

Ｍに関する（２．６２）の総和は、

に対する所望の逆方向伝播式を提供する。

、

が（２．５８）を適用することにより計算される、公式（２．６０）および（２．６３）は、逆方向伝播アルゴリズムの基礎として利用される。

よって、ステップ２６０に示すように、逆方向伝播マトリックスは、カルマンフィルター測定式を形成するために使用され、ステップ２７０において、カルマンフィルターは、ＰＲ測定値およびＡＤＲ差を処理し、状態ベクトルを更新する。

（Ｃ．ナビゲーション性能の説明）
図７および図８は、時間の関数としてのナビゲーション誤差（メートル、ｍによる位置の誤差、メートル／秒、ｍ／ｓによる速度の誤差）のシュミュレーション結果を、図５を併用して上述されたＰＲ＋ＡＤＲ処理技術と比較して、従来のＧＰＳアルゴリズム（ＰＲ測定値のみを使用）に対して示す。本シュミュレーションの環境は以下である。

静止軌道の衛星。

Ｓ／Ｎ_０＞２７ｄＢ・Ｈｚで受信したＧＰＳ信号。これは、ＧＰＳ送信アンテナの主ローブおよび第一のサイドローブのピークにおける追跡信号に対応する。

コード（ＰＲ）追跡誤差は、正規分布、非相関、およびＧＰＳアンテナの照準の信号に対して１ｍの標準偏差をもつように較正される。弱い信号に対しては、追跡誤差は、１／（Ｓ／Ｎ_０の平方根）によって倍率される。

天文歴誤差は、正規分布、時間に高い相関、擬似距離測定値に対して１．４ｍの標準偏差を有する。

キャリア（ＡＤＲ）追跡誤差は、正規分布、非相関、およびコードに対する標準偏差の１／１００の標準偏差を有するように較正される。

利用者クロックは、ドリフトおよび１５ｃｍ／ｓの標準偏差のランダムな運動誤差を有する。

追跡１は、従来のＰＲ処理に対する位置誤差である。追跡２は、ＰＲ＋ＡＤＲ処理を使用した位置誤差である。追跡３は、従来のＰＲ処理に対する速度誤差である。追跡４は、ＰＲ＋ＡＤＲ処理を使用した速度誤差である。従来のＰＲ測定値のみを使用するＧＰＳナビゲーションに対する、本明細書で説明したＰＲ＋ＡＤＲ処理技術のナビゲーション誤差における改善は、ほぼ１桁の規模である。

図８は、図７で示すプロットの最初の２分間の拡大部を示す。図８において、ＡＤＲ処理は、ｔ＝５で開始する。ＡＤＲ処理（追跡４）の開始後すぐの速度ナビゲーション誤差の急な低下は、非常に明白である。

本明細書で説明したＧＰＳ−ＩＭＵナビゲーションアルゴリズムに関して、シミュレーションが以下の条件で実行された。

ＩＭＵは、理想的なデータ（ノイズおよびバイアスはゼロである）を生成する。

擬似距離測定値は、１ｍの標準偏差の正規分布誤差を有する。

ＡＤＲ測定値は、１ｃｍの標準偏差の正規分布誤差を有する。

図９〜図１０はこれらのシミュレーションの結果を示し、二乗和平方根（ＲＳＳ）ナビゲーション誤差が、位置（メートル）、速度（ｍ／ｓ）および姿勢（姿勢マトリックスの成分、ほぼラジアンに対応する）に対して時間関数として示されている。図９は、ＰＲ測定値のみが処理され、伝播アルゴリズムが１次に限られるときの、これらのパラメータに対する誤差を示す。図１０は、ＰＲ測定値のみを使用し２次伝播アルゴリズムを使用したこれらのパラメータに対する誤差を示す。図１１は、ＰＲおよびＡＤＲ差分測定値に対する３次伝播アルゴリズムを使用したこれらのパラメータに対する誤差を示す。これらの図は、精密な伝播アルゴリズムを併用してＡＤＲ処理を使用する、姿勢および速度の推定における２桁規模の改善、および位置推定における１桁規模の改善を示す。

正確な伝播アルゴリズムは宇宙用途におけるＡＤＲ処理の新しい方法を可能にする重要な要因である。本明細書で説明する方法は、宇宙飛行体の速度の測定において、特別な精度を達成する。ＧＥＯ宇宙飛行体に対する推定は、各個別の測定値に対して０．１ｍｍ／ｓの程度の測定精度の達成を約束する。

統合型ＧＰＳ／ＩＭＵシステムに対する伝播式は、同様の方法で使用され得る。より優れた質の伝播は、よりよい精度のために多くの基準時点間隔を測定したＡＤＲ測定値の差の使用を可能にする。空間伝播とＩＭＵ誘導伝播との間の差は、後者において、ＡＤＲ測定値が速度および姿勢の推定に適用されることである。

（ＩＩＩ．その他の導出）
（Ａ．Ｇに対する式の導出）
この導出をとおして、以下の恒等式が使用される。

（１．式１）
最初に、

が計算され、式中、

は任意のベクトルである。故に、

ベクトル

は、座標系のｋ番目のユニットベクトルと等しくセットされる。次いで、導関数が、以下のように計算される。

（２．式２）
次のステップは、以下の式を計算することである。

同様に、

は、座標系のｋ番目のユニットベクトルと等しくセットされる。

故に、

（３．最終式）
これらの構成単位を用いて、以下が計算される。

（Ｂ．ｘ、Ｖ伝播式の導出）

から開始する。

この積分は、位置および速度伝播式の双方を導くために使用され得る。位置においては、最終結果においてα＝１である。速度においては、α＝０、Ｔ＝１と置く。台形法を使用する。

次に、同じ量が推定されるが、２つの半ステップを使用する。

式中、

は、ｔ_ｎ−１からｔ_{ｎ−１／２}への姿勢マトリックスの更新である。第二項Ｊが表記され、それを二つの部分に分ける。

Ｊに対して、

の２つの実現値に対する２つの別々の式（これ以降、厳密な下付き添字が目的の精度限界において問題がないとき、項の下付き添字ｎ、ｎ−１／２、ｎ−１を省く）が代入される。

に対する２つの別の式を代入する。最初の式は、

の一般式からである。

に対する２番目の式は、恒等式から得られる。

次に、

に対して、以下が使用される。

全てのこれらの代入により

を得る。

ｄｔ^２の項を保持する。

主要な項、およびｄｔに関する２次までの推定の全ての他の項を抽出する。本計算に対して、

、等が使用される。

さらにもう一度、ｄｔに関する２次（またはそれ以下）の項のみが保持される。

この式を、Ｉ_ｎ−１に対する第２の方程式（２つの半ステップ公式より得た方程式）に代入することにより、

を得る。

１ステップ式との比較により、

を得る。これは、最終式、

を生成する。

要約すると、測距信号源から受信した信号から物体における擬似距離（ＰＲ）測定値を測定するステップと、測距信号源から受信した信号から物体における累積デルタ距離（ＡＤＲ）測定値を測定するステップと、連続ＡＤＲ測定値間の時間間隔よりも大きな時間間隔によって分離されたＡＤＲ測定値の間のＡＤＲ差を計算するステップと、ＰＲ測定値およびＡＤＲ差から移動物体の少なくとも１つのナビゲーションパラメータを推定するステップとを包含する、移動物体のナビゲーションパラメータを決定する方法が提供される。

同様に、測距信号源から受信した信号から物体における擬似距離（ＰＲ）測定値を測定するステップと、測距信号源から受信した信号から物体における累積デルタ距離（ＡＤＲ）測定値を測定するステップと、連続ＡＤＲ測定値の間の時間間隔よりも大きな時間間隔によって分離されたＡＤＲ測定値の間のＡＤＲ差を計算するステップと、ＰＲ測定値およびＡＤＲ差から現在の状態における物体の少なくとも位置および速度と以前の状態における物体の状態ベクトルを含む、状態ベクトルの成分を推定するステップとを包含する、移動物体の位置および速度を決定する方法が提供される。

さらに、本明細書で説明される方法論は、処理装置またはコンピュータにより実行されるときには、処理装置またはコンピュータが、測距信号源から受信した信号から物体における擬似距離（ＰＲ）測定値を得、測距信号源から受信した信号から物体における累積デルタ距離（ＡＤＲ）測定値を得、連続ＡＤＲ測定値の間の時間間隔よりも大きな時間間隔によって分離されたＡＤＲ測定値の間のＡＤＲ差を計算し、ＰＲ測定値およびＡＤＲ差から移動物体の少なくとも１つのナビゲーションパラメータを推定する、命令（例えば、ソフトウェア）を記憶する、処理装置またはコンピュータで読み出し可能な媒体において具体化され得る。

本明細書で説明されたシステムおよび方法は、その精神または本質的特徴から逸脱することなく、他の特定の形式で具体化され得る。前述の実施形態は、よって、全ての観点において例示的であると考えられ、制限されるべきものではない。

移動物体に搭載され、ＧＰＳ信号のみに基づき、またはＧＰＳ信号とＩＭＵデータとに基づきナビゲーションパラメータを計算するためのＰＲおよびＡＤＲ測定値を使用する、ＧＰＳ衛星およびナビゲーションユニットを含むナビゲーションシステムのブロック図である。 ＰＲおよびＡＤＲ測定値に基づくナビゲーション計算のために、カルマンフィルターにインプットされる測定値およびパラメータを示す図である。 ＰＲおよびＡＤＲ測定値に基づく、全体的なナビゲーション計算処理を示す、流れ図である。 複数の基準時点によって時間で分離されれた二つのＡＤＲ測定値から、ＡＤＲ差が如何にして計算されるかを絵図的に示す図である。 ＰＲおよびＡＤＲ測定値に基づくＧＰＳのみのナビゲーション計算の処理ステップを示す流れ図である。 統合型ＧＰＳ−ＭＵナビゲーション計算の状態ベクトルを伝播するために、ＩＭＵデータを使用するための処理ステップを示す流れ図である。 従来のＧＰＳ処理と比較した、図の流れ図によって示されるＧＰＳのみのナビゲーション計算の、ナビゲーション誤差対時間のプロットを含む図である。 図７に示すプロットから拡大した、ＡＤＲ処理の最初の２分間のナビゲーション誤差対時間のプロットを含む図である。 ＰＲ測定値のみにおけるＩＭＵデータの１次伝播を使用した、ナビゲーション誤差の二乗和平方根の、時間の関数としてのプロットを含む図である。 ＰＲ測定値のみにおけるＩＭＵデータの３次伝播を使用した、ナビゲーション誤差の二乗和平方根の、時間の関数としてのプロットを含む図である。 ＰＲおよびＡＤＲ測定値におけるＩＭＵデータの３次伝播を使用した、ナビゲーション誤差の二乗和平方根の、時間の関数としてのプロットを含む図である。

Claims (41)

1. 移動物体のナビゲーションパラメータを決定する方法であって、
   （ａ）測距信号源から受信する信号から該記物体での擬似距離（ＰＲ）測定値を測定するステップと、
   （ｂ）該記測距信号源から受信する信号から該物体での累積デルタ距離（ＡＤＲ）測定値を測定するステップと、
   （ｃ）連続するＡＤＲ測定値の間の時間間隔よりも大きな時間間隔によって分離される、ＡＤＲ測定値の間のＡＤＲ差を計算するステップと、
   （ｄ）該ＰＲ測定値および該ＡＤＲ差から該移動物体の少なくとも一つのナビゲーションパラメータを推定するステップと
   を包含する、方法。
2. 前記ＡＤＲ差を計算するために使用される、ＡＤＲ測定値の間の最適な時間間隔を決定することをさらに包含する、請求項１に記載の方法。
3. ＡＤＲ測定値の間の前記最適な時間間隔を決定することは、二つのＡＤＲ測定値に対する速度推定値誤差の分散の最小値を計算することを包含する、請求項２に記載の方法。
4. 前記最適な時間間隔を決定することは伝播モデルの精密さに依存し、より精密な伝播モデルが使用されるときには、該最適な時間間隔はより長くなり得る、請求項３に記載の方法。
5. 推定することは、前記物体の瞬間速度推定値を導き出すために、該物体の位置および速度を示す状態ベクトルの少なくとも速度成分に関する、前記ＡＤＲ差の偏導関数を計算することを包含する、請求項１に記載の方法。
6. 推定することは、前記物体の瞬間速度推定値を導き出すために、該物体の位置および速度を示す状態ベクトルの少なくとも速度成分に関する、前記ＡＤＲ差の偏導関数を計算することと、時間Ｎ−１での該状態ベクトルに基づいて時間Ｎにおける該物体に対する該状態ベクトルを計算する伝播方程式を用いて、該状態ベクトルを伝播することを包含する、請求項１に記載の方法。
7. 偏導関数を計算することは、前記状態ベクトルの位置および速度に関する前記ＡＤＲ差の偏導関数を計算することと、前記物体のクロックおよびクロックドリフトに関する該ＡＤＲ差の偏導関数を計算することと、ＡＤＲ差を処理するためのカルマンフィルターにおいてこれらの偏導関数を使用することによって、前記物体の瞬間速度推定値を計算することを包含する、請求項６に記載の方法。
8. 前記移動物体の位置を推定することは、前記ＰＲ測定値および前記推定瞬間速度測定値を使用して、前記状態ベクトルの位置および速度成分を計算することを包含する、請求項７に記載の方法。
9. 前記ＡＤＲ差の偏導関数を計算することは、現在の基準時点での状態ベクトルに関する、過去の基準時点での前記物体の位置の導関数を計算することを包含し、前記状態ベクトルは、ある基準時点での該記物体の少なくとも位置および速度の成分を備える、請求項８に記載の方法。
10. 現在の基準時点での前記状態ベクトルに関する、過去の基準時点における前記物体の位置の導関数を計算することは、連続する基準時点の間の時間差に関する高次近似を計算することを包含する、請求項９に記載の方法。
11. 高次近似を計算することは、複数のサブステップを有する伝播ルーチンを使用して、ＡＤＲ測定値を分離する時間間隔上で逆方向に前記状態ベクトルを伝播することと、該伝播ルーチンの該サブステップのそれぞれに対して、Λマトリックスの要素に対応する
    

    

    マトリックスを計算することと、それぞれの該
    

    

    マトリックスから該Λマトリックを計算することと、該Λマトリックスを使用して、該状態ベクトルの位置および速度に関する、過去の基準時点での前記物体位置の偏微分を計算することを包含する、請求項１０に記載の方法。
12. 強制モデルおよび統合型ルーチンを使用して、前記状態ベクトルを伝播することと、現在の基準時点の該状態ベクトルに関する次の基準時点の該状態ベクトルの偏導関数を計算することによって、共分散マトリックスを伝播することをさらに包含する、請求項１１に記載の方法。
13. 前記状態ベクトルは、前記物体の位置、速度、クロックおよびクロックドリフトに対する区分を有し、該状態ベクトルを伝播することは、速度に関する位置の導関数と、クロックに関する位置の導関数と、クロックドリフトに関する位置の導関数と、位置に関する速度の導関数と、クロックに関する速度の導関数と、クロックドリフトに関する速度の導関数とに対する区分を有するマトリックスを計算することを包含する、請求項１２に記載の方法。
14. 計算するステップ（ｃ）は、ＡＤＲ測定値が測定される連続する基準時点の間の期間よりも大きな前記期間によって分離される二つの基準時点での、ＡＤＲ測定値の間のＡＤＲ差を計算することを包含する、請求項１に記載の方法。
15. 計算するステップ（ｃ）は、ＡＤＲ測定値が測定される連続する基準時点の間の期間よりも数段階も規模の大きな前記期間によって分離される二つの基準時点での、ＡＤＲ測定値の間のＡＤＲ差を計算することを包含する、請求項１４に記載の方法。
16. 推定することは、前記物体に関連する慣性測定ユニット（ＩＭＵ）によって生成される慣性測定データに基づいて、前記状態ベクトルを伝播することをさらに包含する、請求項１に記載の方法。
17. 伝播することは、前記物体の状態ベクトルに関する前記ＡＤＲ差の偏導関数を備える観測マトリックスを計算することを包含し、該状態ベクトルは、該物体に関連する慣性測定ユニット（ＩＭＵ）によって生成される慣性測定データの位置、速度、姿勢および誤差状態の推定を備える、請求項１６に記載の方法。
18. 前記観測マトリックスを計算することは、前記ＩＭＵの体軸におけるジャイロセンサーデータを示すベクトルを、地球中心慣性（ＥＣＩ）座標系に変換する回転マトリックスを計算することと、前記状態ベクトルの前記姿勢成分を伝播するために基準時点に対して該回転マトリックスを伝播することを包含する、請求項１７に記載の方法。
19. 前記観測マトリックスを計算することは、前記状態ベクトルの前記位置および速度成分を伝播することを包含する、請求項１８に記載の方法。
20. カルマンフィルター状態ベクトルの不整合ベクトルを含む、前記ＩＭＵのジャイロセンサーの不整合パラメータを特徴づける不整合ベクトルを計算することをさらに包含する、請求項１９に記載の方法。
21. 順方向伝播マトリックスを提供するために、前記基準時点の開始の前記状態ベクトルに関する、各基準時点の終わりの該状態ベクトルの偏導関数を計算することをさらに包含する方法であって、各基準時点において、該状態ベクトルの前記姿勢区分は、該姿勢の現在の推定に関するオイラー角によって示される、請求項２０に記載の方法。
22. 複数の基準時点にわたり逆方向伝播マトリックスを計算することをさらに包含する、請求項２１に記載の方法。
23. 回帰的アルゴリズムを使用して、前記逆方向伝播マトリックスを用いて前記状態ベクトルを逆方向に伝播することをさらに包含する、請求項２２に記載の方法。
24. 測定するステップ（ａ）は、グローバル測位衛星（ＧＰＳ）から受信する信号から前記物体でのＰＲ測定値を測定することを包含する、請求項１に記載の方法。
25. 移動物体の位置および速度を決定する方法であって、
    （ａ）測距信号源から受信する信号から、該物体での擬似距離（ＰＲ）測定値を測定するステップと、
    （ｂ）該測距信号源から受信する信号から、該物体での累積デルタ距離（ＡＤＲ）測定値を測定するステップと、
    （ｃ）連続するＡＤＲ測定値の間の時間間隔よりも大きな時間間隔によって分離される、ＡＤＲ測定値の間のＡＤＲ差を計算するステップと、
    （ｄ）該ＰＲ測定値および該ＡＤＲ差および以前の状態での状態ベクトルから、現在の状態での前記物体の少なくとも位置および速度を備える状態ベクトルの成分を推定するステップと
    を包含する、方法。
26. 計算するステップ（ｃ）は、前記時間間隔によって分離されるＡＤＲ測定値に対するＡＤＲ差を計算することを包含し、該時間間隔は、二つのＡＤＲ測定値に対する速度推定値誤差の分散の最小値を計算することによって決定される、請求項２５に記載の方法。
27. 推定することは、前記物体の瞬間速度推定値を導き出すために、前記状態ベクトルの少なくとも速度成分に関する、前記ＡＤＲ差の偏導関数を計算することを包含する、請求項２５に記載の方法。
28. 前記ＡＤＲ差の偏導関数を計算することは、現在の基準時点での前記状態ベクトルに関する、過去の基準時点での前記物体の位置の導関数を計算することを包含する、請求項２７に記載の方法。
29. 現在の基準時点での前記状態ベクトルに関する、過去の基準時点での前記物体の位置の導関数を計算することは、連続する基準時点の間の時間差に関する高次近似を計算することを包含する、請求項２８に記載の方法。
30. 高次近似を計算することは、複数のサブステップを有する伝播ルーチンを使用して、ＡＤＲ測定値を分離する時間間隔上で逆方向に前記状態ベクトルを伝播することと、該伝播ルーチンの該サブステップのそれぞれに対して、Λマトリックスの要素に対応するjマトリックスを計算することと、それぞれの該jマトリックスから該Λマトリックスを計算することと、該Λマトリックスを使用して、該状態ベクトルの位置および速度に関する、過去の基準時点での該物体位置の偏微分を計算することを包含する、請求項２９に記載の方法。
31. 強制モデルおよび統合型ルーチンを使用して、前記状態ベクトルを伝播することと、現在の基準時点の該状態ベクトルに関する、次の基準時点の該状態ベクトルの偏導関数を計算することによって、共分散マトリックスを伝播することをさらに包含する、請求項３０に記載の方法。
32. 推定することは、前記物体の瞬間速度推定値を導き出すために、前記状態ベクトルの少なくとも速度成分に関する、前記ＡＤＲ差の偏導関数を計算することと、時間Ｎ−１における該状態ベクトルに基づいて時間Ｎにおける該物体に対する該状態ベクトルを計算する、伝播方程式を用いて該状態ベクトルを伝播することを包含する、請求項２５に記載の方法。
33. 推定することは、各項が対応する測距信号源に対する、二つの項の差によって定義される前記物体の位置および速度に関する第一区分と、該物体のクロックおよびクロックドリフトに関する第二区分とを備える、偏微分マトリックスを含む一対の測距信号源に対するＡＤＲ測定値残差を、カルマンフィルターを用いて処理することを包含する、請求項２５に記載の方法。
34. 計算するステップ（ｃ）は、ＡＤＲ測定値が測定される連続する基準時点の間の期間よりも数段階規模の大きな前記期間によって分離される二つの基準時点での、ＡＤＲ測定値の間のＡＤＲ差を計算することを包含する、請求項２５に記載の方法。
35. 推定することは、前記物体に関連する慣性計測ユニット（ＩＭＵ）によって提供される慣性計測データに基づいて、前記状態ベクトルを伝播することをさらに包含する、請求項２５に記載の方法。
36. 伝播することは、前記物体の状態ベクトルに関する、前記ＡＤＲ差の偏導関数を備える観測マトリックスを計算することを包含し、該状態ベクトルは、前記物体に関連する慣性計測ユニット（ＩＭＵ）によって生成される慣性計測値データの位置、速度、姿勢および誤差状態の推定を備える、請求項３５に記載の方法。
37. 前記観測マトリックスを計算することは、前記ＩＭＵの体軸におけるジャイロセンサーデータを示すベクトルを、地球中心慣性（ＥＣＩ）座標系に変換する回転マトリックスを計算することと、該状態ベクトルの前記姿勢成分を伝播するために、基準時点に対して該回転マトリックスを伝播することを包含する、請求項３６に記載の方法。
38. 前記観測マトリックスを計算することは、前記状態ベクトルの位置および速度成分を伝播することを包含し、カルマンフィルター状態ベクトルの中の不整合ベクトルを含む、前記ＩＭＵの前記ジャイロセンサーの不整合パラメータを特徴づける不整合ベクトルを計算することをさら包含する、請求項３７に記載の方法。
39. 順方向伝播マトリックスを作成するために、前記基準時点の開始の前記状態ベクトルに関する、各基準時点の終わりの該状態ベクトルの偏導関数を計算することをさらに包含する方法であって、各基準時点において、該状態ベクトルの前記姿勢区分は、該姿勢の現在の推定に関するオイラー角によって示される、請求項３８に記載の方法。
40. 複数の基準時点に対する逆方向伝播マトリックスを計算することと、回帰法を使用して、該逆方向伝播マトリックスを用いて前記状態ベクトルを逆方向に伝播することをさらに包含する、請求項３９に記載の方法。
41. 命令を記憶する、処理装置で読み取り可能な媒体であって、該命令は、処理装置によって実行されるときには、該処理装置に、
    （ａ）測距信号源から受信する信号から、物体での擬似距離（ＰＲ）測定値を取得させ、
    （ｂ）該測距信号源から受信する信号から、該物体での累積デルタ距離（ＡＤＲ）測定値を取得させ、
    （ｃ）連続するＡＤＲ測定値の間の時間間隔よりも大きな時間間隔によって分離される、ＡＤＲ測定値の間のＡＤＲ差を計算させ、
    （ｄ）該ＰＲ測定値および該ＡＤＲ差から移動物体の少なくとも一つのナビゲーションパラメータを推定させる、
    媒体。

Images