(57)【要約】【課題】衛星の姿勢を推定する推定誤差分散が小さく、推定の計算による計算機負荷が小さな人工衛星の姿勢決定装置を得る。【解決手段】人工衛星の姿勢決定装置は、２台以上のスタートラッカと、各スタートラッカで取得した画像から星像を検出するスター検出部と、星像とスターカタログに格納された星像との比較照合および同定を行うスター同定部と、同定結果から各スタートラッカの姿勢を推定する姿勢推定部と、姿勢推定部で推定した各スタートラッカの姿勢から衛星姿勢を推定する最尤推定部とを備えた人工衛星の姿勢決定装置において、最尤推定部は、先に姿勢推定部で推定した各スタートラッカの姿勢から衛星姿勢の第１近似解を求め、後に各スタートラッカの姿勢推定誤差統計値と衛星姿勢の第１近似解とを用いて衛星姿勢の最尤推定値を求め、最尤推定値を衛星姿勢決定値とする。【選択図】        図１
2004224171.tif 000002
【特許請求の範囲】
【請求項１】
複数台のスタートラッカと、
上記各スタートラッカで取得した画像から星像を検出するスター検出部と、
上記スター検出部で検出した星像とスターカタログに格納された星像との比較照合および同定を行うスター同定部と、
上記スター同定部における同定結果から上記各スタートラッカの姿勢を推定する姿勢推定部と、
上記姿勢推定部で推定した上記各スタートラッカの姿勢から衛星姿勢を推定する最尤推定部と
を備えた人工衛星の姿勢決定装置において、
上記最尤推定部は、上記姿勢推定部で推定した上記各スタートラッカの姿勢から上記衛星姿勢の第１近似解を求め、各スタートラッカの姿勢推定誤差統計値と上記衛星姿勢の上記第１近似解とを用いて上記衛星姿勢の最尤推定値を求め、
上記最尤推定値を衛星姿勢決定値とすることを特徴とする人工衛星の姿勢決定装置。
【請求項２】
上記衛星姿勢がクォータニオン形式で表されたことを特徴とする請求項１に記載の人工衛星の姿勢決定装置。
【請求項３】
上記スター同定部で同定した星像の位置から上記各スタートラッカの姿勢推定誤差統計値を算出する姿勢推定誤差統計値算出部を備え、
上記姿勢推定部で推定した各スタートラッカの姿勢と上記姿勢推定誤差統計値算出部で求めた各スタートラッカの姿勢推定誤差統計値とから衛星姿勢を推定することを特徴とする請求項１乃至２のいずれか一項に記載の人工衛星の姿勢決定装置。
【請求項４】
人工衛星の姿勢角速度を検出する慣性基準装置と、
上記姿勢角速度から求まる衛星姿勢を伝播する姿勢伝播部および上記最尤推定部で求められた推定値を姿勢基準値として、上記衛星姿勢伝播値を更新して姿勢更新値を求める姿勢更新部とを備えた時間前向きフィルタ部とを備え、
上記姿勢更新部で更新された姿勢更新値を衛星姿勢の姿勢決定値とすることを特徴とする請求項１乃至３のいずれか一項に記載の人工衛星の姿勢決定装置。
【請求項５】
上記時間前向きフィルタ部における姿勢決定値を初期値として、時間反転させて姿勢決定値を得る時間後向きフィルタ部と、
上記時間前向きフィルタ部による姿勢決定値と上記時間後向きフィルタ部による姿勢決定値とを平均化する平均化処理部とを備え、
上記平均化処理部は、上記時間前向きフィルタ部による姿勢決定値と上記時間後向きフィルタ部による姿勢決定値の重み付き平均を行い、得られた姿勢を衛星姿勢の姿勢決定値とすることを特徴とする請求項４に記載の人工衛星の姿勢決定装置。
【請求項６】
上記時間前向きフィルタ部における姿勢決定値を用いて、時間後向きに姿勢決定値を得る時間後向きスムーザ部を備え、
上記時間後向きスムーザ部は、上記時間前向きフィルタ部より得られた姿勢伝播値と上記時間後向きスムーザ部より求めた姿勢決定値との偏差を求め、上記偏差を用いて上記時間前向きフィルタ部より得られた姿勢決定値に対する更新処理を行う姿勢更新部を備え、
上記時間後向きスムーザ部で得られた姿勢更新値を、衛星姿勢の姿勢決定値とすることを特徴とする請求項４に記載の人工衛星の姿勢決定装置。
【請求項７】
上記時間後向きフィルタ部における姿勢決定値を用いて、時間前向きに姿勢決定値を得る時間前向きスムーザ部を備え、
上記時間前向きスムーザ部は、上記時間後向きフィルタ部より得られた姿勢逆伝播値と上記時間前向きスムーザ部より求めた姿勢決定値から偏差を求め、これを用いて上記時間後向きフィルタ部より得られた姿勢決定値に対する更新処理を行う姿勢更新部を備え、
上記時間前向きスムーザ部で得られた姿勢更新値を、衛星姿勢の姿勢決定値とすることを特徴とする請求項５に記載の人工衛星の姿勢決定装置。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
この発明は、スタートラッカを用いて人工衛星の姿勢を決定する装置に関する。
【０００２】
【従来の技術】
人工衛星の姿勢決定装置は複数台のスタートラッカと１台の地球センサと、中央処理部から構成される。スタートラッカでは、視野内の星像を検出、追尾するとともに、検出した星像とスターカタログとの照合を行い、スター同定を行う。スタートラッカからは、検出した星像のスタートラッカ固定座標系における位置ベクトルＳｉｊｍと、スターカタログで照合した対応する恒星の慣性座標系での位置ベクトルｒｉｊｃが出力される。また地球センサからは、衛星から見た地球指向角φＥＬ、θＡＺが出力される。中央処理部では、各スタートラッカ毎に、検出した星像の位置ベクトルと、スターカタログにおける対応する恒星の位置ベクトルに関して平均化処理が行われ、各スタートラッカにつき１組の位置ベクトルが得られる。また地球センサ出力からは、衛星から見た地球の位置ベクトルが求められる。これらの複数の位置ベクトルと１つの位置ベクトルを用いて、最小二乗法により、衛星の姿勢をオイラー角を用いて決定している（例えば特許文献１参照。）。
【０００３】
【特許文献１】
米国特許第５４１２５７４号明細書
【０００４】
【発明が解決しようとする課題】
上記先行例では、衛星の姿勢変動が、基準とする軌道座標系に対して微少である場合を想定している。そのため、軌道座標系に対する衛星姿勢角をオイラー角で表現し、オイラー角に関して、軌道座標系まわりに線形近似して姿勢決定を行っている。このような手法は、宇宙サービス機のように大きな姿勢変動を伴う衛星には適用できない。
【０００５】
一般に、１台のスタートラッカにより求められる姿勢（慣性座標系に対するスタートラッカ固定座標系の姿勢）に関しては、そのｂｏｒｅｓｉｇｈｔ軸まわりの姿勢推定誤差は大きく、それに直交する軸まわりの姿勢推定誤差は小さいという誤差特性を持っている。しかし上述の従来技術では、複数台のスタートラッカから得られた位置ベクトルに対して、このような誤差特性を考慮しない、すなわち重み付けを行わない最小二乗法により姿勢角を決定している。そのため決定される衛星姿勢角は、その分散が最小となるような最尤推定値ではない。
【０００６】
この発明の目的は、衛星の姿勢を推定する推定誤差分散が小さく、推定の計算による計算機負荷が小さな人工衛星の姿勢決定装置を得ることである。
【０００７】
【課題を解決するための手段】
この発明に係わる人工衛星の姿勢決定装置は、複数台のスタートラッカと、各スタートラッカで取得した画像から星像を検出するスター検出部と、スター検出部で検出した星像とスターカタログに格納された星像との比較照合および同定を行うスター同定部と、スター同定部における同定結果から各スタートラッカの姿勢を推定する姿勢推定部と、姿勢推定部で推定した各スタートラッカの姿勢から衛星姿勢を推定する最尤推定部とを備えた人工衛星の姿勢決定装置において、最尤推定部は、先に姿勢推定部で推定した各スタートラッカの姿勢から衛星姿勢の第１近似解を求め、次に各スタートラッカの姿勢推定誤差統計値と衛星姿勢の第１近似解とを用いて衛星姿勢の最尤推定値を求め、最尤推定値を衛星姿勢決定値とする。
【０００８】
【発明の実施の形態】
実施の形態１．
図１は、この発明の実施の形態１に係る人工衛星の姿勢決定装置の構成を示すブロック図である。図２は、図１の最尤推定部の詳細な構成を示すブロック図である。
【０００９】
人工衛星の姿勢決定装置は、人工衛星に搭載された複数台（ここではＮ台とする。）のスタートラッカ１と、スタートラッカ１で取得した画像から星像を検出するスター検出部２と、恒星情報を格納したスターカタログ３と、スター検出部２で検出した星像とスターカタログ３に格納されている恒星情報との比較照合および同定を行うスター同定部４と、スター同定部４で求められた同定結果から各スタートラッカ１の姿勢を求める姿勢推定部５と、姿勢推定部５で推定された各スタートラッカ１の姿勢から衛星姿勢の最尤推定値を求める最尤推定部６とを備えている。
【００１０】
最尤推定部６は、姿勢推定部５で推定した各スタートラッカ１の姿勢から衛星姿勢の初期推定値を求める第１近似解算出部７と、各スタートラッカ１に対する姿勢推定誤差共分散行列Ｗｓｊを記憶している姿勢推定誤差統計値記憶部８と、各スタートラッカ１に対する姿勢推定誤差共分散行列Ｗｓｊと衛星姿勢の初期推定値とを用いて衛星姿勢の最尤推定値を求め、求めた最尤推定値を衛星姿勢決定値とする衛星姿勢推定部９である。
【００１１】
次にこの人工衛星の姿勢決定装置の動作について図１を参照しながら説明する。スタートラッカ１において天球の画像が取得されると、スター検出部２において取得画像からスター検出部２の視野内にある星像を抽出し、各スタートラッカ１に固定された座標系で表した星像の中心位置や輝度等の星像データが検出される。スターカタログ３には、慣性座標系における恒星の位置ベクトルおよび輝度等の恒星情報が記憶されている。スター同定部４では、スター検出部２で検出された星像のスタートラッカに固定された座標系の中心位置から求めた星像間の距離および輝度等をスターカタログ３に格納されている恒星情報と比較照合して、検出した星像とスターカタログ３に格納されている恒星とを同定する。
【００１２】
スター検出部２で検出された星像の位置ベクトルをｓｉｊｍで表す。ここでｓｉｊｍは第ｊ番目（ｊ＝１、２、…、Ｎ）のスタートラッカ１で検出された第ｉ番目（ｉ＝１、２、・・・、ξｊ）の星位置ベクトルを表す。上添え字のｍは観測値を示す。この星像位置ベクトルｓｉｊｍは、各スタートラッカ１に固定されたスタートラッカ固定座標系で記述されている。
【００１３】
一方、スターカタログ３においてｓｉｊｍと同定された恒星の星位置ベクトルをｃｉｊで表す。この星位置ベクトルｃｉｊは慣性座標系で記述されている。
【００１４】
姿勢推定部５では、スタートラッカ固定座標系における星像位置ベクトルｓｉｊｍと、慣性座標系における星位置ベクトルｃｉｊの対応関係から、慣性座標系における第ｊ番目のスタートラッカの姿勢を推定して、クォータニオンｑｓｊｍで記述する。これらのξｊ組の同定された星位置ベクトルから姿勢を推定する際には、例えば特開平１０−７２０００号公報、人工衛星の姿勢決定装置に示された手法を用いて推定できる。こうして推定された第ｊ番目のスタートラッカの姿勢をｑｓｊｍで記述する。
【００１５】
最尤推定部６では、姿勢推定部５で推定された各スタートラッカの姿勢ｑｓｊｍ（ｊ＝１、２、…、Ｎ）から衛星姿勢の最尤推定値ｑｂ（ハット）を求める。衛星姿勢はクォータニオン形式で表されている。ここで下述の説明のために、以下の記号を定義する。ｑｓｊは慣性座標系においてクォータニオン形式で表された第ｊ番目のスタートラッカの姿勢真値、ｑｓｊｍは慣性座標系においてクォータニオン形式で表された第ｊ番目のスタートラッカの姿勢推定値、ｑｓｊｂは衛星固定座標系におけるクォータニオン形式で表された第ｊ番目のスタートラッカの姿勢、ｑｂは慣性座標系におけるクォータニオン形式で表された衛星の姿勢、Ｗｓｊは第ｊ番目のスタートラッカにおける姿勢推定誤差共分散行列、Ｉｍはｍ次元単位行列である。
【００１６】
この姿勢推定誤差共分散行列Ｗｓｊは、仕様値として予めスタートラッカ毎に与えられ、姿勢推定誤差統計値記憶部８に記憶されている。一般的に推定されるスタートラッカの姿勢の姿勢推定誤差は、スタートラッカのｂｏｒｅｓｉｇｈｔ軸まわりには大きく、それに直交する軸まわりには小さいという誤差特性を持っている。
【００１７】
任意のクォータニオン形式で表された値ｑに関して、以下の記号を定義する。Ｖ（ｑ）はクォータニオン形式で表された値ｑのベクトル部、Ｓ（ｑ）はクォータニオン形式で表された値ｑのスカラー部、数１はクォータニオン形式で表された値ｑのベクトル部Ｖ（ｑ）の符号を反転したものである。
【００１８】
【数１】
2004224171.tif 000003
【００１９】
任意の２つのクォータニオン形式で表された値ｑ１とｑ２との積を数２で定義する。
【００２０】
【数２】
2004224171.tif 000004
【００２１】
ここでベクトル部Ｖ（ｑ１）に対する上付き添字の波（〜）は、Ｖ（ｑ１）の外積演算Ｖ（ｑ１）×を示している。上付き添えのＴは転置マトリックスを示している。
【００２２】
次にこの最尤推定部６の動作について説明する。第ｊ番目のスタートラッカにおける姿勢推定誤差をδｅｊで記述すると、定義よりδｅｊは次のように表される。
【００２３】
【数３】
2004224171.tif 000005
【００２４】
従って、衛星姿勢ｑｂの最尤推定値を求める問題は、クォータニオン形式で表された衛星姿勢に関する正規化条件が課された次の重み付き最小二乗問題となる。拘束条件すなわち数４で示される正規化条件の下で、数５で示される評価関数Ｊが最小となるｑｂを求めることである。
【００２５】
【数４】
2004224171.tif 000006
【００２６】
【数５】
2004224171.tif 000007
【００２７】
ここで第ｊ番目のスタートラッカにおける姿勢推定誤差共分散行列Ｗｓｊは、姿勢推定誤差統計値記憶部８から得ることができる。
【００２８】
次に拘束条件付きの最小二乗問題を解くことであるが、これを直接解くには計算負荷が過大になるので近似解法により、衛星姿勢ｑｂの最尤推定値を求める。
【００２９】
ｑｂ０は求める衛星姿勢ｑｂの第１近似解であり、初期推定値になり、ｑｂを第１近似解ｑｂ０からの偏差ｑδを用いて数６のように表わすことができる。
【００３０】
【数６】
2004224171.tif 000008
【００３１】
ここで偏差ｑδはクォータニオン形式で表された微少オイラー角から構成されている。また第１近似解ｑｂ０は、第１近似解算出部７において、次のようにして求めることができる。基準とするスタートラッカ（ここでは仮に第ｋ番目のスタートラッカを基準とするものと仮定する。）に対して推定された姿勢ｑｓｋｍを用いて数７を用いてｑｂ０を求める。
【００３２】
【数７】
2004224171.tif 000009
【００３３】
また次のようなクォータニオン形式で表された微少オイラー角ｑεｊを定義する。
【００３４】
【数８】
2004224171.tif 000010
【００３５】
数７と数８を用いて、衛星姿勢推定誤差δｅｊは数９のように書き直される。
【００３６】
【数９】
2004224171.tif 000011
【００３７】
ここでＣｓｊｍは、慣性座標系に対する第ｊ番目のスタートラッカ固定座標系への座標変換行列で、クォータニオン形式で表された第ｊ番目のスタートラッカの姿勢推定値ｑｓｊｍを用いて構成される。
【００３８】
数９を数５に示される評価関数Ｊに代入すると、上述の最小二乗問題は、数１０に示すような評価関数Ｊを最小とするＶ（ｑδ）を求めることになる。
【００３９】
【数１０】
2004224171.tif 000012
【００４０】
この最小二乗問題を解くと、数１１に示すようにＶ（ｑδ）を求めることができる。
【００４１】
【数１１】
2004224171.tif 000013
【００４２】
数１１から求められた解Ｖ（ｑδ）およびＳ（ｑδ）＝１を用いて偏差ｑδを構成し、それを用いて衛星姿勢ｑｂの最尤推定値ｑｂ（ハット）を数６より求める。このようにして求められた衛星姿勢ｑｂの最尤推定値ｑｂ（ハット）を衛星姿勢決定値ｑｂ（バー）として出力する。
【００４３】
このような人工衛星の姿勢決定装置は、各スタートラッカに対してあらかじめ決められた各スタートラッカの姿勢推定誤差統計値を利用し衛星姿勢を推定するため、得られる姿勢推定値は、その分散が小さく高精度な姿勢決定値を得ることができる。
【００４４】
また、最尤推定部において衛星姿勢の最尤推定値を求める際に、クォータニオンに関して拘束条件（正規化条件）のもと重み付き最小二乗問題を近似解法を採用することにより、拘束条件無しの重み付き最小二乗問題に変換され、計算負荷を小さくすることができ、衛星に搭載されたオンボード計算機でも十分実行することができる。
【００４５】
また、最尤推定部６では、衛星姿勢をクォータニオン形式で表現しているため、大きな姿勢角変動を伴う衛星等にも適用することができる。
【００４６】
実施の形態２．
図３は、この発明の実施の形態２に係る人工衛星の姿勢決定装置の構成を示すブロック図である。なお、図中、同一符号は図１と同様であり、同じ部分の説明は省略する。
【００４７】
図３において、スター同定部４において、スターカタログ３との同定された星像の位置から、姿勢推定部５で推定される姿勢推定値に関する姿勢推定誤差共分散の下界を求める姿勢推定誤差統計値算出部１０と、姿勢推定部５で推定された各スタートラッカの姿勢と、姿勢推定誤差統計値算出部１０で求めた姿勢推定誤差共分散の下界から衛星姿勢の最尤推定値を求め、求めた最尤推定値を衛星姿勢決定値とする最尤推定部１１とである。
【００４８】
次に図３の人工衛星の姿勢決定装置の動作について説明する。一般的にスタートラッカにより推定される姿勢の姿勢推定誤差に関しては、そのｂｏｒｅｓｉｇｈｔ軸まわりには大きく、それに直交する軸まわりには小さいという誤差特性を持っている。ただし、この誤差特性は、後述するように、スター同定部４で同定された星像数や、同定された星像の位置に依存している。そこで、これらの情報を利用することにより、姿勢推定部５において推定される各スタートラッカの姿勢の姿勢推定誤差をより詳細に評価することができる。ここではその一例としてクラメル・ラオの下界を用いて姿勢推定部５において推定される姿勢の姿勢推定誤差共分散を評価し、姿勢推定誤差統計値算出部１０において、その下界を姿勢推定誤差共分散行列Ｗｓｊとして出力する。最尤推定部１１では、この姿勢推定誤差共分散行列Ｗｓｊを用いて衛星姿勢の最尤推定値ｑｂ（ハット）を求める。
【００４９】
実施の形態１で説明したのと同様に、スター検出部２で検出された星像の位置ベクトルをｓｉｊｍで表す。また第ｊ番目のスタートラッカに関して、以下の記号を定義する。
Θｊは視野角、Ｍｊは画素数、ｆｊは焦点距離である。
ｓｉｊｍ、ｓｉｊ、ｃｉｊをそれぞれ次のように成分表現する。
【００５０】
【数１２】
2004224171.tif 000014
【００５１】
【数１３】
2004224171.tif 000015
【００５２】
【数１４】
2004224171.tif 000016
【００５３】
ここでα、β、γはスタートラッカに固定された座標系で表された星像の座標値（観測値）、Ｘ、Ｙ、Ｚはスタートラッカに固定された座標系で表された星像の座標値（真値）、ｘ、ｙ、ｚは慣性座標系で表されたスターカタログに記憶されている恒星の座標値をそれぞれ表す。
【００５４】
画素の配置が等方的であると仮定し、スター検出部２における星像位置の検出誤差がホワイトノイズであると仮定すれば、第ｊ番目のスタートラッカにおいて、スター検出部２での星像位置検出に関する条件付き確率密度関数は、数１５のように与えられる。衛星姿勢角を表すオイラー角をφ、θ、ψで記述する。
【００５５】
【数１５】
2004224171.tif 000017
【００５６】
ここでσｊ２は第ｊ番目のスタートラッカに対するスター検出部２における星像位置の検出誤差の分散を表している。例えば検出誤差として量子化誤差を想定すると、σｊ２は画素の大きさにより規定される。検出される星像が全て独立であると仮定すると、第ｊ番目のスタートラッカに対するスター検出部２における検出値全体の確率密度関数は数１６である。
【００５７】
【数１６】
2004224171.tif 000018
【００５８】
これより第ｊ番目のスタートラッカ１における星像位置検出に関する対数尤度ηｊおよびＦｉｓｈｅｒの情報行列Ｊｊが、それぞれ数１７と数１８で表される。
【００５９】
【数１７】
2004224171.tif 000019
【００６０】
【数１８】
2004224171.tif 000020
【００６１】
ここでＥ［  ］は期待値を表す。
【００６２】
クラメル・ラオの下界は、パラメータ推定における分散の最小値を与え、推定値の良さを評価する一つの尺度となる。このクラメル・ラオの下界を用いると、姿勢推定誤差共分散は次のように評価することができる。
（ロール角φ推定誤差共分散）≧  Ｊｒｊ−１
（ピッチ角θ推定誤差共分散）≧  Ｊｐｊ−１
（ヨー角ψ推定誤差共分散）  ≧  Ｊｙｊ−１
【００６３】
ここでＪｒｊ−１、Ｊｐｊ−１、Ｊｙｊ−１がクラメル・ラオの下界を示し、これは前述のＦｉｓｈｅｒの情報行列Ｊｊの逆行列の対角成分で与えられる。具体的にφ＝θ＝ψ＝０の場合に、Ｆｉｓｈｅｒの情報行列Ｊｊの逆行列を求めると、数１９のようになる。
【００６４】
【数１９】
2004224171.tif 000021
【００６５】
数１９の対角成分に着目すると、一般に第３成分は、第１、２成分に比べて小さいため、スタートラッカにより推定される姿勢の姿勢推定誤差は、そのｂｏｒｅｓｉｇｈｔ軸まわりには大きく、それに直交する軸まわりには小さい誤差特性を持つ。この解析結果を用いて、姿勢推定誤差統計値算出部１０では、姿勢推定部５で推定された第ｊ番目のスタートラッカの姿勢に対する姿勢推定誤差共分散行列を数２０のように構成する。
【００６６】
【数２０】
2004224171.tif 000022
【００６７】
ここでＪｒｊ−１、Ｊｐｊ−１、Ｊｙｊ−１は、スター同定部４で同定した星像の位置ベクトルを用いて構成される、次に示すＦｉｓｈｅｒの情報行列Ｊｊの逆行列の対角成分である。
【００６８】
【数２１】
2004224171.tif 000023
【００６９】
ここでσｊは、スター検出部２における検出誤差として、１画素分の量子化誤差を考えると、星像位置の検出誤差の３σ値が１画素の大きさに相当するものと考えて、３σｊは次のように与えることができる。
【００７０】
【数２２】
2004224171.tif 000024
【００７１】
こうして姿勢推定部５で推定される各スタートラッカの姿勢に関する推定誤差共分散行列Ｗｓｊを、数２０、２１、２２より求めることができる。姿勢推定誤差統計値算出部１０では、こうして求めた推定誤差共分散の下界を姿勢推定誤差共分散行列Ｗｓｊとして出力する。
【００７２】
最尤推定部１１では、姿勢推定部５で推定した各スタートラッカの姿勢と、姿勢推定誤差統計値算出部１０で求めた姿勢推定誤差共分散行列Ｗｓｊとから、衛星姿勢の最尤推定値ｑｂ（ハット）を求め、求めた最尤推定値を衛星姿勢決定値ｑｂ（バー）とする。
【００７３】
このような人工衛星の姿勢決定装置を用いると、各スタートラッカに対して推定される姿勢に関して、姿勢推定誤差統計値算出部においてその姿勢推定誤差を評価し、姿勢推定誤差共分散行列を構成し、それを利用して衛星姿勢を推定するため、スター同定部で同定された星像数や、同定された星像の位置に依存している姿勢推定誤差を反映しているので、得られる姿勢推定値は、実施の形態１よりさらに、その分散が小さくなり、より高精度な姿勢決定値を得ることができる。
【００７４】
実施の形態３．
図４は、この発明の実施の形態３に係る人工衛星の姿勢決定装置の構成を示すブロック図である。なお、図中、同一符号は図１と同様であり、同じ部分の説明は省略する。
【００７５】
図４に示すように、図１の人工衛星の姿勢決定装置に、衛星の姿勢角速度を検出する慣性基準装置（以下ＩＲＵと称す。）１２と、ＩＲＵ１２により検出された姿勢角速度と衛星姿勢決定値とを用いて、衛星姿勢を伝播させる姿勢伝播部１３と、姿勢伝播部１３で伝播させた衛星姿勢と、最尤推定部６で推定した衛星姿勢の最尤推定値との間の偏差を用いて衛星姿勢を更新し、更新した姿勢を衛星姿勢決定値とする姿勢更新部１４とがさらに備えられている。ここで姿勢伝播部１３と姿勢更新部１４とから時間前向きフィルタ部１５を構成している。時間前向きとは時刻的に先の測定値をもとに、後の各時刻の状態変数を推定することを意味する。
【００７６】
次に図４の人工衛星の姿勢決定装置の動作について説明する。なお、最尤推定部６で推定された衛星姿勢の最尤推定値ｑｂ（ハット）を、ここでは姿勢基準値ｑｍ（ハット）と名前を呼び変えて説明するが、同一のものである。
【００７７】
衛星には、姿勢角速度を検出するためのＩＲＵ１２には姿勢角速度センサとしてジャイロが搭載されており、これらのセンサ信号を用いてカルマンフィルタにより衛星姿勢が決定される。
【００７８】
カルマンフィルタの状態変数を、姿勢角推定誤差とジャイロドリフト推定誤差の６変数から構成する。任意の変数ｘに対して、その推定値を上付き添字の＾（ハット）で表す。また時刻ｔまでの観測値ｘｔが与えられている場合に、この観測値ｘｔに基づき時刻ｔ＋ｍにおける条件付き期待値をｘｔ＋ｍ／ｔで表す。
【００７９】
また、ｄはジャイロドリフト、Ｐはカルマンフィルタにおける状態推定誤差共分散行列、Ｑは状態方程式におけるシステムノイズ共分散行列、Ｈは観測方程式における観測行列、Ｒは観測方程式における観測ノイズ共分散行列を示す。
【００８０】
姿勢伝播部１３では、時刻ｔ−１においてＩＲＵ１２で観測された姿勢角速度ωｔ−１ｍを用いて、時刻ｔ−１における衛星姿勢決定値ｑｂ，ｔ−１／ｔ−１（ハット）とジャイロドリフト決定値ｄｔ−１／ｔ−１（ハット）とを伝播させ、それぞれの時刻ｔにおける伝播値を求める。なお、姿勢角速度決定値ωｔ−１／ｔ−１（ハット）は数２３で求められる。
【００８１】
【数２３】
2004224171.tif 000025
【００８２】
【数２４】
2004224171.tif 000026
【００８３】
【数２５】
2004224171.tif 000027
【００８４】
また数２４の右辺第２項におけるΔｔは離散時間を表す。
【００８５】
同様に、状態推定誤差共分散行列に関しても、時刻ｔ−１における値を伝播させ、時刻ｔにおける値を数２６に基づいて求める。
【００８６】
【数２６】
2004224171.tif 000028
【００８７】
ここで行列Ｆｔ−１は、時刻ｔ−１における状態遷移行列である。
【００８８】
次にここで求めた時刻ｔにおける衛星姿勢伝播値ｑｂ，ｔ／ｔ−１（ハット）と、最尤推定部６で求めた時刻ｔにおける衛星姿勢の姿勢基準値ｑｍ（ハット）との間の偏差νｔを数２７で求める。
【００８９】
【数２７】
2004224171.tif 000029
【００９０】
ここで、ｑｍ（ハット）は最尤推定部６で求めた衛星姿勢の最尤推定値ｑｂ（ハット）を姿勢基準値と定義したものである。
【００９１】
姿勢更新部１４では、ここで求めた偏差νｔを用いて時刻ｔ−１における衛星姿勢ｑｂ，ｔ／ｔ−１（ハット）とジャイロドリフトｄｔ／ｔ−１（ハット）の更新処理を行う。時刻ｔにおけるカルマンゲインＫｔが数２８によって求められる。
【００９２】
【数２８】
2004224171.tif 000030
【００９３】
数２９に示すように、姿勢偏差νｔにカルマンゲインＫｔを乗じて、更新量Δλｔを得る。
【００９４】
【数２９】
2004224171.tif 000031
【００９５】
ここで求めた更新量Δλｔを用いて時刻ｔにおける衛星姿勢ｑｂ，ｔ／ｔ（ハット）とジャイロドリフトｄｔ／ｔ（ハット）をそれぞれ数３０と数３１で更新する。
【００９６】
【数３０】
2004224171.tif 000032
【００９７】
【数３１】
2004224171.tif 000033
【００９８】
ここでΔλｑ、Δλｄはそれぞれ、前述の時刻ｔにおける更新量Δλｔの衛星姿勢成分とジャイロドリフト成分とである。
【００９９】
また状態推定誤差共分散行列も数３２で更新する。
【０１００】
【数３２】
2004224171.tif 000034
【０１０１】
姿勢更新部１４では、こうして得られた姿勢更新値ｑｂ，ｔ／ｔ（ハット）を姿勢決定値ｑｂ（バー）とする。
【０１０２】
このような人工衛星の姿勢決定装置は、最尤推定部で推定された衛星姿勢推定値と衛星角速度とを使用して、カルマンフィルタを用いたフィルタリング処理を行って衛星姿勢を推定するので、衛星の時系列の角速度とその統計情報を反映した衛星姿勢決定値を得られ、実施の形態２よりさらに高精度な姿勢決定値を得ることができる。
【０１０３】
実施の形態４．
図５は、この発明の実施の形態４に係る人工衛星の姿勢決定装置の構成を示すブロック図である。なお、図中、同一符号は図４と同様であり、同じ部分の説明は省略する。
【０１０４】
図５に示すように、図４の姿勢決定装置に、図４の時間前向きフィルタ部１５で求まる姿勢更新値と時間後向きフィルタ部１６で求められた時間後向き姿勢推定値とを平均化して衛星姿勢推定値を求める平均化処理部１７とをさらに備えている。時間後向きフィルタ部１６は、ＩＲＵ１２により検出された姿勢角速度と姿勢決定値とを用いて、衛星姿勢を時間後向き方向に伝播させる姿勢逆伝播部１８と、姿勢逆伝播部１８で伝播させた衛星姿勢と、最尤推定部６で推定した衛星姿勢の最尤推定値ｑｍ（ハット）との間の偏差を用いて衛星姿勢を更新し、更新した姿勢を衛星姿勢決定値とする姿勢更新部１９とを備えている。
【０１０５】
時間後向きフィルタ部１６は地上でのオフライン処理により人工衛星の姿勢決定を行う場合を想定している。
【０１０６】
次に図５の人工衛星の姿勢決定装置の動作について説明する。時間前向きフィルタ部１５では、図４と同様な手法に従って人工衛星の姿勢決定値を得る。所望の時刻範囲に含まれる衛星姿勢決定値とジャイロドリフト決定値と状態推定誤差共分散行列とを時間前向きフィルタ部１５からオフラインで収集し、地上に送信し、それら値を時間後向きフィルタ部１６における初期値とする。
【０１０７】
姿勢逆伝播部１８では、時刻ｔ＋１においてＩＲＵ１２で観測された姿勢角速度ωｔ＋１ｍを用いて、初期値とされた時刻ｔ＋１における衛星姿勢決定値とジャイロドリフト決定値とを時間後向きに伝播させ、それぞれの時刻ｔにおける推定値を数３３と数３４で求める。ここで姿勢角速度決定値ωｔ＋１／ｔ＋１（ハット）は数３５で求まる。
【０１０８】
【数３３】
2004224171.tif 000035
【０１０９】
【数３４】
2004224171.tif 000036
【０１１０】
【数３５】
2004224171.tif 000037
【０１１１】
ここで数３３の右辺第２項におけるΔｔは離散時間を表す。
【０１１２】
同様に、状態推定誤差共分散行列に関しても、時刻ｔ＋１における値を用いて次のように伝播させ、時刻ｔにおける値を数３６に基づいて求める。
【０１１３】
【数３６】
2004224171.tif 000038
【０１１４】
次にここで求めた時刻ｔにおける衛星姿勢伝播値ｑｂ，ｔ／ｔ＋１（ハット）と、最尤推定部６で求めた時刻ｔにおける衛星姿勢の姿勢基準値ｑｍ（ハット）との間の偏差νｔを数３７で求める。
【０１１５】
【数３７】
2004224171.tif 000039
【０１１６】
姿勢更新部１９では、ここで求めた偏差νｔを用いて時刻ｔ＋１における衛星姿勢ｑｂ，ｔ／ｔ＋１（ハット）とジャイロドリフトｄｔ／ｔ＋１（ハット）の更新処理を行う。時刻ｔにおけるカルマンゲインＫｔが数３８によって求められる。
【０１１７】
【数３８】
2004224171.tif 000040
【０１１８】
数３９に示すように、姿勢偏差νｔにカルマンゲインＫｔを乗じて、更新量Δλｔを得る。
【０１１９】
【数３９】
2004224171.tif 000041
【０１２０】
ここで求めた更新量Δλｔを用いて時刻ｔにおける衛星姿勢ｑｂ，ｔ／ｔ（ハット）とジャイロドリフトｄｔ／ｔ（ハット）をそれぞれ数４０と数４１で更新する。
【０１２１】
【数４０】
2004224171.tif 000042
【０１２２】
【数４１】
2004224171.tif 000043
【０１２３】
ここでΔλｑ、Δλｄはそれぞれ、前述の時刻ｔにおける更新量Δλｔの衛星姿勢成分とジャイロドリフト成分とである。
【０１２４】
また状態推定誤差共分散行列も数４２で更新する。
【０１２５】
【数４２】
2004224171.tif 000044
【０１２６】
姿勢更新部１９では、こうして得られた姿勢更新値ｑｂ，ｔ／ｔ（ハット）を姿勢決定値ｑｂ（バー）とする。
【０１２７】
このように、時間前向きフィルタ部１５から時間後向きフィルタ部１６の順にカルマンフィルタを実行し、時間前向きフィルタ部１５より得られた時刻ｔにおける衛星の姿勢決定値を上付き添字ｆで表す。また時間後向きフィルタ部１６より得られた時刻ｔにおける衛星の姿勢決定値を上付き添字ｂで表す。平均化処理部１７は、数４３に示すようにこれらの姿勢決定値を用いて重み付け平均を行い、得られた姿勢を衛星の姿勢決定値ｑｂ（バー）とする。
【０１２８】
【数４３】
2004224171.tif 000045
【０１２９】
ここでスカラー係数ｃは重み係数で、０≦ｃ≦１の値をとるものとする。また重み付け平均を行った後には正規化を行って、クォータニオンの正規化条件を満足させる。
【０１３０】
重み係数ｃに関しては、例えば、時間前向きフィルタ部１５による姿勢決定値と、時間後向きフィルタ部１６による姿勢決定値が、同程度に確かな場合には、ｃを０．５とすればよい。
【０１３１】
このような人工衛星の姿勢決定装置は、衛星姿勢基準値として時間前向きフィルタ部を構成し、また各スタートラッカに対してあらかじめ定められた誤差特性を先見情報として用いて衛星姿勢の最尤推定値を得、これを衛星姿勢基準値として時間後向きフィルタ部を構成し、これらのフィルタにより決定された姿勢決定値の重み付け平均値を衛星の姿勢決定値とするため、高精度な姿勢決定値を得ることができる。
【０１３２】
実施の形態５．
図６は、この発明の実施の形態５に係る人工衛星の姿勢決定装置の構成を示すブロック図である。なお、図中、同一符号は図４と同様であり、同じ部分の説明は省略する。
【０１３３】
図６に示すように、図４の人工衛星の姿勢決定装置にさらに時間後向きスムーザ部２０を備えている。時間後向きスムーザ部２０は、姿勢伝播部１３から出力される姿勢伝播値と姿勢決定値との偏差を求め、その偏差と時間前向きフィルタ部１５において決定した姿勢決定値とを用いて衛星姿勢を更新し、更新した姿勢を衛星姿勢決定値とする姿勢更新部２１を備えている。
【０１３４】
次に図６の人工衛星の姿勢決定装置の動作について説明する。図６は、図５と同様に地上でのオフライン処理により人工衛星の姿勢決定を行う場合を想定している。
【０１３５】
このオフライン処理のように、ある時刻区間［０，Ｎ］の観測データが与えられた時、その区間での全ての時点における状態の最適推定値を求める問題は、固定区間スムージングと言われている。
【０１３６】
時間前向きフィルタ部１５では、図４に記載した手法に従って人工衛星の姿勢決定値を得る。時間前向きフィルタ部１５で求めた各時刻における姿勢伝播値、姿勢決定値、ジャイロドリフト伝播値、ジャイロドリフト決定値、状態推定誤差共分散行列伝播値、状態推定誤差共分散行列更新値、状態遷移行列を記憶しておく。これら値は人工衛星から送信され、地上の記憶装置に記憶されている。
【０１３７】
時間後向きスムーザ部２０では、時間前向きフィルタ部１５におけるこれらの計算結果を用いてスムージング処理が行われる。
【０１３８】
時間前向きフィルタ部１５における時刻ｔ＋１における衛星姿勢伝播値、ジャイロドリフト伝播値と、時間後向きスムーザ部２０における時刻ｔ＋１での姿勢決定値、ジャイロドリフト決定値の間の偏差を数４４で求める。
【０１３９】
【数４４】
2004224171.tif 000046
【０１４０】
姿勢更新部２１では、ここで求めた偏差を用いて、時間前向きフィルタ部１５で決定した時刻ｔにおける衛星姿勢とジャイロドリフトの更新処理を行う。時刻ｔにおける平滑ゲインが数４５から求められる。
【０１４１】
【数４５】
2004224171.tif 000047
【０１４２】
前述の偏差に平滑ゲインを乗じて、更新量を得る。
【０１４３】
【数４６】
2004224171.tif 000048
【０１４４】
ここで求めた更新量を用いて、時刻ｔにおける衛星姿勢とジャイロドリフトを更新する。
【０１４５】
【数４７】
2004224171.tif 000049
【０１４６】
【数４８】
2004224171.tif 000050
【０１４７】
ここでΔλｑ、Δλｄはそれぞれ、前述の時刻ｔにおける更新量Δλｔの衛星姿勢成分、ジャイロドリフト成分である。姿勢更新部２１では、こうして得られた姿勢更新値を衛星姿勢決定値とする。
【０１４８】
このように各スタートラッカに対してあらかじめ定められた誤差特性を先見情報として用いて衛星姿勢の最尤推定値を得、これを衛星姿勢基準値として時間前向きフィルタ部を構成し、時間前向きフィルタ部における姿勢決定値、状態遷移行列、状態推定誤差共分散行列を用いて時間後向きスムーザ部を構成し、時間後向きスムーザ部で更新した姿勢を衛星の姿勢決定値とするため、高精度な姿勢決定値を得ることができる。
【０１４９】
実施の形態６．
地上でのオフライン処理により人工衛星の姿勢決定を行う場合には、実施の形態４で説明した時間前向きフィルタ部および時間後向きフィルタ部と、実施の形態５で説明した時間前向きフィルタ部および固定区間スムーザを組み合わせて衛星の姿勢決定値を得る手法以外も考えられる。
【０１５０】
図７は、この発明の実施の形態６に係る人工衛星の姿勢決定装置の構成を示すブロック図である。図７は図５の姿勢決定装置に、時間前向きスムーザ部２２とを備えている。
【０１５１】
時間前向きスムーザ部２２は、時間後向きフィルタ部１６の姿勢逆伝播部１８から出力される姿勢逆伝播値と姿勢決定値との偏差を求め、その偏差と時間後向きフィルタ部１６の姿勢決定値とを用いて姿勢決定値を更新する姿勢更新部２３を備えている。
【０１５２】
次に図７の人工衛星の姿勢決定装置の動作について説明する。図５の姿勢決定装置で説明したように、時間前向きフィルタ部１５から姿勢角速度と姿勢伝播値と姿勢決定値とジャイロドリフト伝播値、ジャイロドリフト決定値、状態推定誤差共分散行列伝播値、状態推定誤差共分散行列更新値を衛星から送信して、地上の記憶装置に記憶する。さらに図５の姿勢決定装置で説明したように、これら値と入力として時間後向きフィルタ部１６でそれぞれの逆伝播したのち、決定値を更新して、時間後向きフィルタ部１６で求めた各時刻における姿勢逆伝播値、姿勢決定値、ジャイロドリフト逆伝播値、ジャイロドリフト決定値、状態推定誤差共分散行列逆伝播値、状態推定誤差共分散行列更新値、状態遷移行列を記憶しておく。
【０１５３】
次に、時間前向きスムーザ部２２では、時間後向きフィルタ部１６におけるこれらの計算結果を用いてスムージング処理が行われる。ここでの処理は、時間の向きが前向きになったことを除いて、本質的に実施の形態５で記載したのと同様である。姿勢更新部２３において、こうして得られた姿勢更新値を衛星姿勢決定値とする。
【０１５４】
このように各スタートラッカに対してあらかじめ定められた誤差特性を先見情報として用いて衛星姿勢の最尤推定値を得、これを衛星姿勢基準値として時間前向きフィルタ部、時間後向きフィルタ部を構成し、時間後向きフィルタ部における姿勢決定値、状態遷移行列、状態推定誤差共分散行列を用いて時間前向きスムーザ部を構成し、時間前向きスムーザ部で更新した姿勢を衛星の姿勢決定値とするため、高精度な姿勢決定値を得ることができる。
【０１５５】
【発明の効果】
この発明に係わる人工衛星の姿勢決定装置による効果は、複数台のスタートラッカと、各スタートラッカで取得した画像から星像を検出するスター検出部と、スター検出部で検出した星像とスターカタログに格納された星像との比較照合および同定を行うスター同定部と、スター同定部における同定結果から各スタートラッカの姿勢を推定する姿勢推定部と、姿勢推定部で推定した各スタートラッカの姿勢から衛星姿勢を推定する最尤推定部とを備えた人工衛星の姿勢決定装置において、最尤推定部は、先に姿勢推定部で推定した各スタートラッカの姿勢から衛星姿勢の第１近似解を求め、後に各スタートラッカの姿勢推定誤差統計値と衛星姿勢の第１近似解とを用いて衛星姿勢の最尤推定値を求め、最尤推定値を衛星姿勢決定値とするので、各スタートラッカに対してあらかじめ決められた各スタートラッカの姿勢推定誤差統計値を利用し衛星姿勢を推定するため、得られる姿勢推定値は、その分散が小さく高精度な姿勢決定値を得ることができる。
【０１５６】
また、最尤推定部において衛星姿勢の最尤推定値を求める際に、クォータニオンに関して拘束条件（正規化条件）のもと重み付き最小二乗問題を近似解法を採用することにより、拘束条件無しの重み付き最小二乗問題に変換され、計算負荷を小さくすることができ、衛星に搭載されたオンボード計算機でも十分実行することができる。
【図面の簡単な説明】
【図１】この発明の実施の形態１に係る人工衛星の姿勢決定装置の構成を示すブロック図である。
【図２】図１の最尤推定部の詳細な構成を示すブロック図である。
【図３】この発明の実施の形態２に係る人工衛星の姿勢決定装置の構成を示すブロック図である。
【図４】この発明の実施の形態３に係る人工衛星の姿勢決定装置の構成を示すブロック図である。
【図５】この発明の実施の形態４に係る人工衛星の姿勢決定装置の構成を示すブロック図である。
【図６】この発明の実施の形態５に係る人工衛星の姿勢決定装置の構成を示すブロック図である。
【図７】この発明の実施の形態６に係る人工衛星の姿勢決定装置の構成を示すブロック図である。
【符号の説明】
１  スタートラッカ（ＳＴＴ）、２  スター検出部、３  スターカタログ、４スター同定部、５  姿勢推定部、６、１１  最尤推定部、７  第１近似解算出部、８  姿勢推定誤差統計値記憶部、９  衛星姿勢推定部、１０  姿勢推定誤差統計値算出部、１２  慣性基準装置（ＩＲＵ）、１３  姿勢伝播部、１４、１９、２１、２３  姿勢更新部、１５  時間前向きフィルタ部、１６  時間後向きフィルタ部、１７  平均化処理部、１８  姿勢逆伝播部、２０  時間後向きスムーザ部、２２  時間前向きスムーザ部。
【図１】
2004224171.tif 000051
【図２】
2004224171.tif 000052
【図３】
2004224171.tif 000053
【図４】
2004224171.tif 000054
【図５】
2004224171.tif 000055
【図６】
2004224171.tif 000056
【図７】
2004224171.tif 000057