4. 2.
MASTER-2009 を⽤いた衝突シミュレーション
本節では，MASTER-2009 を⽤いて衝突データのシミュレーションを⾏う．まず，シミュレーション 条件を表 4.2に⽰す．
表 4.2 シミュレーション条件
観測期間 2007/04/01 - 2009/04/01 [yyyy/mm/dd]
観測衛星の軌道要素
軌道⻑半径 7176.130 [km]
離⼼率 0.0
昇交点⾚経 212.226 [deg]
軌道傾斜⾓ 98.6 [deg]
考慮するデブリサイズ 10 – 0.1 [mm]
観測衛星は⾼度798 kmの太陽同期軌道に投⼊され，2年間の軌道上観測を⾏うことを想定している．
また，観測期間を2007年から2009年としているのは，MASTER-2009の開発時よりも未来のデブリ分 布予測は過去のものと⽐較して正確性が劣るため，それを避けて最も現在に近い2年間を選択したため である．
表 4.2の条件をMASTER-2009に⼊⼒することで，Cell Passage Events (CPE) dump file が出⼒される．
CPEファイルには，観測衛星に対しデブリ衝突の可能性がある全てのセルの位置，セル内での衝突フラ ックス，観測衛星とデブリの相対速度が記録されている．なお，MASTER-2009 の低軌道における計算 ではデブリの昇交点⾚経については考慮せず平均化した計算を⾏っていることに注意が必要である．
CPEファイルに出⼒された衝突フラックスは，1 m2の断⾯積を持つ球体に対する全⽅位からのフラッ クスを計算したものである．そのため，まず式（2.13）を⽤いて⼊射⾓に関する補正を⾏い，センサ⾯
積 0.125 m2を乗じて観測衛星のセンサに対する衝突フラックスを求める．本シミュレーションにおける 観測衛星のセンサ⾯は，図 2.4に⽰すように衛星の進⾏⽅向に対し左右に45°ずつ傾けて2⾯設置され ているとする．
次に，補正されたフラックスから実際に衝突データのシミュレーションを⾏う．シミュレーションで は，衝突の可能性があるそれぞれのセルに対し0から1までの⼀様乱数を計算機により発⽣させ，乱数 の値よりもそのセルにおける衝突フラックスが⼤きい場合に実際に衝突が起きたと判定する．この操作 により，シミュレートされた衝突の時刻および観測衛星の位置が⼀連の観測データとして⽣成される．
なお，このシミュレーションでは乱数を使⽤するため，試⾏する度に結果が変化する．
図 4.2，図 4.3に，以上の⼿法を⽤いて実際に⽣成された観測データ2例を⽰す．図は横軸に観測開 始からの時間，縦軸に衝突時の観測衛星の緯度引数をとっている．2年間の観測シミュレーションにお いて，図 4.2では111個，図 4.3では139個のデブリが観測された．なお，MASTER-2009の予測する 観測衛星に対する総フラックスは年間 61 個程度であり，乱数により衝突回数が増減していることがわ
図 4.2 シミュレートされた観測データ
図 4.3 シミュレートされた観測データ（別例）
0 45 90 135 180 225 270 315 360
0 90 180 270 360 450 540 630 720
Argument of Latitude [deg]
Elapsed Time [day]
0 45 90 135 180 225 270 315 360
0 90 180 270 360 450 540 630 720
Argument of Latitude [deg]
Elapsed Time [day]
また，図 4.2と図 4.3のいずれの例でも，緯度引数90°および270°，つまり極付近における衝突が 多くなっている．これは，低軌道では太陽同期軌道をはじめとする極軌道がミッション軌道として選択 される傾向にあるため⾼い軌道傾斜⾓を持つデブリが多く，また極軌道の物体の衝突は⾼い緯度で発⽣
しやすいためであると考えられる．
4. 3.
モデルの初期分布
衝突データを⽤いて環境推定を⾏うには，まず初期分布となる𝑋Â|Âを設定する必要がある．⼀般的に，
推定の初期分布として真値に近い分布を与えるほど推定結果は向上する．本論⽂で推定対象とする微⼩
デブリの真の環境についての情報は限られているため，何らかの⽅法で妥当な初期分布を与えなければ ならない．
既存の環境モデルであるMASTERおよびORDEMは，これまでの観測データや⼤型のデブリの分布 をもとに微⼩デブリの分布を推定しており，現実に参照できる情報の中では最も真のデブリ分布に近い ものを与えていると考えられる．しかし，4.1 節に⽰したように本論⽂ではこれら既存のモデルを真値 としたシミュレーションにより動的環境推定モデルの性能を検証する．したがって既存のモデルをもと に初期分布を与えることは初めから正解を与えることとなり検証の⽬的にあわないため，本論⽂では既 存のモデルは初期分布としては採⽤しない．なお，実際の観測においては既存のモデルによる微⼩デブ リ分布予測をもとに初期分布を与え，それを観測に合わせ更新していくのが最も真の分布に近い推定結 果を得られると考えられる．
図 4.4 被追跡物体の軌道⾯の分布
i cosΩ
i sin Ω
以上から本論⽂では，軌道が追跡され公開されている⼤型の宇宙機およびデブリの分布をもとに初期 分布の決定を⾏う．微⼩デブリは⼤型の物体から破砕や剥落により放出されるものであるため，⼤型物 体の分布と微⼩デブリの分布に相関性を⾒いだすことは妥当であると考えられる．ただし，燃焼ガスに 含まれる微粒⼦⼤型物体の軌道とは異なる⽅向に⾼速で放出されることとなる存在するため，本節で与 える初期分布は完全に微⼩デブリの分布の傾向を反映したものとはならない．あくまで微⼩デブリの分 布予測使⽤できない状況での便宜的⽅法である．
図 4.4に観測衛星と同⾼度にある被追跡物体の軌道⾯を傾斜⾓ベクトルで⽰す．図 4.4を⾒ると，同
⼼円状の分布が⽬⽴っている．傾斜⾓ベクトルでは周⽅向に昇交点⾚経，半径⽅向に軌道傾斜⾓をとっ ていることから，デブリの分布は昇交点⾚経よりも軌道傾斜⾓に対して強い相関を持つことがわかる．
また，最も外側の同⼼円は軌道傾斜⾓98°程度の太陽同期軌道に多くの衛星とデブリが分布しているこ とを表しており，それよりも軌道傾斜⾓の⼤きい逆⾏軌道にはほとんど物体がないことがわかる．
これらの被追跡物体の分布の特徴から，本論⽂では次の2通りの⽅法により初期分布を与え，Case A, Bとする．
Case A 軌道傾斜⾓100度以下で⼀様に分布
このCaseでは，式（3.24）に定義される状態ベクトルのΩ®, 𝑖® (𝑗 = 1 … 𝑛)について，Ω®は0〜360°，𝑖®は 0〜100°の範囲で⼀様乱数を与える．これは，太陽同期軌道よりも⼤きな軌道傾斜⾓を持つ物体がほと んど存在しないという特徴のみを反映したものである．Case Aにおける初期分布の傾斜⾓ベクトルを図 4.5に⽰す．
i cosΩ
180
90
0
-90
-180
isinΩ
-180 -90 0 90 180
2.0 107
1.5 107
1.0 107
0.5 107
0
Debris Population
Case B ⾮追跡物体の分布をもとに⽣成した分布
このCaseでは，状態ベクトルのΩ®, 𝑖® (𝑗 = 1 … 𝑛)に関して，Ω®は0〜360°の⼀様乱数，𝑖®は図 4.4に⽰
すような被追跡物体の軌道からランダムに選んだ値に正規分布に従うノイズを加えた値とする．この操 作により，軌道傾斜⾓の分布は被追跡物体の分布をノイズによって平滑化したものとなる．この初期分 布は，被追跡物体の分布が昇交点⾚経よりも軌道傾斜⾓に関して⼤きな相関を持つという特徴を反映し たものである．図 4.6にCase B における初期分布の傾斜⾓ベクトルを⽰す．
図 4.6 初期分布 Case B
また，状態ベクトルに含まれる軌道の数であるnは 16とし，それぞれの軌道に含まれるデブリの個 数を表す𝑁°は1.875 × 10Æに正規分布に従うノイズを加えたものとする．これにより，デブリの総数の平 均は3 × 10lÂとなる．デブリの総数に関しては参照すべきデータが他にないため，MASTER-2009の推定 するデブリ総数2.8 × 10lÂと近い値を与えた．
i cosΩ
180
90
0
-90
-180
isinΩ
-180 -90 0 90 180
2.0 107
1.5 107
1.0 107
0.5 107
0
Debris Population
4. 4.
推定結果
本節では，図 4.2に⽰すシミュレートされた観測データを⽤いてCase A, Bそれぞれについて推定を
⾏う．まず，推定モデル内で⽤いるノイズの標準偏差を表 4.3に⽰す．これらのノイズは経験的に推定 の上⼿くいったものを選択しており，実際の観測データを蓄積して運⽤を重ねることでより最適な値を 設定できる可能性がある．また，昇交点⾚経のノイズを⼤きくすることで推定が上⼿くいくのは，𝑑Ω 𝑑𝑡⁄ の値が数度程度と⼤きいこと，モデル化の際無視した離⼼率が𝑑Ω 𝑑𝑡⁄ に影響すること，MASTER内の計 算で昇交点⾚経を平均化して計算していることなどの影響を吸収できるためではないかと考えられる．
表 4.3 ノイズの標準偏差 デブリ数 8 × 10Ç [-/day]
昇交点⾚経 4 [deg/day]
軌道傾斜⾓ 0.01 [deg/day]
SMCフィルタによるモデルの更新操作は，デブリの衝突回数に関わらず5⽇おきに⾏った．そのため，
5⽇のうちの衝突回数が0であったというデータも存在し得るが，3.4で⽰したように，本推定モデルで は衝突が起こらなかったという事実も尤度の計算に組み込まれているため正常に推定を⾏うことが可 能である．
まず，図 4.7，8に1度⽬の更新を⾏ったフィルタ分布，すなわち観測開始から5⽇⽬図 4.7 推定結 果（Case A, Day 5）のデブリ分布推定結果を傾斜⾓ベクトル図で⽰す．いずれの図においても，はっき りと明るい曲線が現れている．観測開始から5⽇までの間に観測衛星は1回のデブリ衝突を経験してお り，この曲線は，その観測データに対応する拘束⽅程式と⼀致する．すなわち，推定の結果として観測 された位置において検出され得る軌道に存在するデブリの個数が増加したことを⽰している．また，特 に図 4.7に顕著であるが，曲線は途中で途切れている．これは，式（2.14）に⽰したセンサ⾯を考慮し た補正においてセンサの裏側から⾶来すると判断される軌道⾯が正しく排除されていることを⽰して いる．更に，曲線中のデブリ数は均等ではなく，軌道により異なった衝突フラックスの計算が⾏われて いることもわかる．以上から，SMCフィルタに組み込まれた拘束⽅程式およびトーラスモデルが正常に 働いたことが確認できた．
続いて，図 4.9-16に，各ステップにおけるデブリ分布の推定結果をCase A，Bそれぞれについて半年 おきに⽰す．いずれのCaseにおいても，図 4.5，6に⽰した初期分布と⽐較して徐々に軌道⾯によって 濃淡の差が⼤きくなっており，観測データに基づき環境の推定が進んでいることがわかる．
続いて図 4.12，16 に⽰す2 年間の観測データを利⽤した推定結果をみると，いずれも同⼼円状の分 布となっており，参考とした⼤型物体の分布と同様に軌道傾斜⾓について⾼い相関を持つ分布となった ことがわかる．また，最も多くのデブリが分布しているのはいずれの結果においても軌道傾斜⾓が80°
〜100°の極軌道である．⼀⽅で，初期分布の違いにより２つの分布は全く同⼀とはなっていない．
MASTER-2009