宇宙空間で回転するデブリに自律的に接近し、ロボットアームで捕獲する — こうしたミッションの実現には、カメラ映像から対象物体の位置と姿勢をリアルタイムに推定する技術が欠かせません。前回の記事では、エッジベース手法やモデルベース手法（ICP、PnP）といった古典的なアプローチを解説しましたが、これらの手法には初期値依存性とロバスト性の限界という本質的な課題がありました。初期姿勢の推定が大きくずれていれば局所解に陥り、極端な照明条件やオクルージョンがあれば容易に破綻します。
では、人間はどうやって見慣れない角度から撮影された衛星の写真を見て「これは太陽パネルが手前にあって、本体が奥を向いている」と判断できるのでしょうか？ それは、私たちが過去に見てきた無数の物体の画像から、形状と姿勢の関係を暗黙的に学習しているからです。深層学習は、まさにこの「大量のデータからパターンを学習する」能力をコンピュータに与える技術であり、宇宙物体の姿勢推定に革命をもたらしつつあります。
深層学習による宇宙物体の姿勢推定を理解すると、以下のような応用が見えてきます。
- 軌道上サービス（OOS）の自律化: 初期値なしに画像一枚から対象衛星の6DoFポーズを推定し、ランデブー全フェーズを自律化する
- 能動的デブリ除去（ADR）: 形状が事前に不明なデブリに対しても、学習された特徴量を頼りに姿勢を推定する
- 群衛星のオンボード知能: 限られた計算資源の衛星搭載コンピュータ上で、リアルタイムに姿勢推定を実行する
- 宇宙状況認識（SSA）: 地上望遠鏡の画像から軌道上物体の姿勢・形状を分類する
本記事の内容
- CNNによる直接6DoF回帰の原理と限界
- キーポイント検出 + PnPによる間接的な姿勢推定
- 合成データ（CGレンダリング）による学習戦略
- Sim-to-Realドメインギャップとその克服手法
- SPEED/SPEED+データセットと宇宙姿勢推定ベンチマーク
- Transformerベースの最新手法
- オンボード推論の計算資源制約と対策
- Pythonでの実装と実験
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
なぜ深層学習が必要か — 古典手法の限界
初期値問題
前回の記事で解説したエッジベース手法やICP（Iterative Closest Point）は、いずれも反復最適化に基づいています。射影モデルのエッジと画像中のエッジを合わせる、あるいは3D点群と2D観測の対応を最適化する — これらは姿勢パラメータに関する非線形最適化問題であり、収束するためには「ある程度正しい初期値」が必要です。
ランデブーの近距離フェーズでは前フレームの推定結果を初期値にできますが、遠距離からの初期検出フェーズでは、対象物体がどのような姿勢を取っているか全く手がかりがない状態から始まります。テンプレートマッチングで粗い推定を行うアプローチもありますが、テンプレートの数が膨大になり、極端な照明条件では破綻しやすいのが実情です。
深層学習は、この「初期値がない状態から画像一枚で姿勢を推定する」という問題を解決する強力なアプローチです。
ロバスト性の課題
古典手法は、特徴抽出と姿勢計算が別個のパイプラインとして設計されています。エッジ検出は照明条件に敏感であり、SIFTやORBの特徴点はテクスチャのない宇宙機表面では十分に検出できません。一方、深層学習は画像全体からタスクに最適化された特徴量を自動的に学習するため、手設計の特徴量に依存しません。
極端なコントラスト、部分的なオクルージョン、背景の地球の反射光 — これらの困難な条件に対しても、十分な学習データがあればロバストに動作する可能性を秘めています。
計算パラダイムの転換
古典手法は「推論時に最適化を反復する」ため、対象物体の複雑さに応じて計算時間が増大します。一方、深層学習は「学習時に膨大な計算を行い、推論時はフォワードパス一回」で完了するため、推論時間が一定です。これは計算資源が限られた衛星のオンボードコンピュータにとって、予測可能な処理時間を保証できるという大きな利点になります。
これらの動機を踏まえて、まずは最も直接的なアプローチ — CNNで画像から6DoFポーズを直接回帰する手法から見ていきましょう。
CNNによる直接6DoF回帰
基本アーキテクチャ
最も素朴な深層学習による姿勢推定のアプローチは、画像を入力として受け取り、6DoFのポーズパラメータを直接出力するCNNを構築することです。これはPoseNet（Kendall et al., 2015）として知られるアイデアの宇宙への応用です。
身近な例で言えば、友人の顔写真を見て「この人はカメラに対して少し左を向いている」と直感的に判断するような処理を、ニューラルネットワークに行わせるイメージです。CNNが画像の特徴を階層的に抽出し、最終的に回転と並進のパラメータを出力します。
ネットワークの構造は概ね以下の通りです。
- バックボーン（特徴抽出器）: ResNetやEfficientNetなどの事前学習済みCNNが画像から高次元の特徴ベクトルを抽出
- 回帰ヘッド: 全結合層が特徴ベクトルから姿勢パラメータを出力
出力の設計が重要で、回転の表現方法によって学習の安定性が大きく変わります。
回転表現の選択
画像から回転を回帰するとき、どのような表現で出力するかは非自明な問題です。主な選択肢を比較します。
クォータニオン $\bm{q} \in \mathbb{R}^4$: 4次元ベクトルとして回帰し、ノルムを1に正規化します。コンパクトな表現ですが、$\bm{q}$ と $-\bm{q}$ が同じ回転を表す二重被覆問題があり、学習を混乱させる可能性があります。
回転行列の6D表現: Zhou et al. (2019) が提案した手法で、回転行列 $\bm{R}$ の最初の2列（6次元ベクトル）を回帰し、グラム・シュミット直交化で残りの1列を復元します。$3 \times 3$ の9次元回転行列から直交条件の3自由度分を除いた、最小限の連続表現です。
回転行列 $\bm{R} = [\bm{r}_1 \ \bm{r}_2 \ \bm{r}_3]$ の最初の2列 $\bm{a}_1, \bm{a}_2$ をネットワークが出力したとき、正規直交基底への変換は次のように行います。
まず、$\bm{a}_1$ を正規化して第1列を得ます。
$$ \bm{r}_1 = \frac{\bm{a}_1}{\|\bm{a}_1\|} $$
次に、$\bm{a}_2$ から $\bm{r}_1$ 方向の成分を除去し、正規化して第2列を得ます。
$$ \bm{r}_2 = \frac{\bm{a}_2 – (\bm{r}_1 \cdot \bm{a}_2)\bm{r}_1}{\|\bm{a}_2 – (\bm{r}_1 \cdot \bm{a}_2)\bm{r}_1\|} $$
最後に、外積で第3列を計算します。
$$ \bm{r}_3 = \bm{r}_1 \times \bm{r}_2 $$
この手法の利点は、回帰する6次元ベクトルの値域が $\mathbb{R}^6$ 全体であり（$SO(3)$ 上の不連続性がない）、ネットワークの学習が安定することです。
オイラー角 $(\phi, \theta, \psi)$: 直感的ですがジンバルロック問題と周期性（$0°$ と $360°$ が同じ）があり、回帰には不向きです。
損失関数の設計
6DoFポーズの推定では、回転と並進という物理的に異なる量を同時に最適化する必要があります。単純な方法は、回転の損失 $\mathcal{L}_{\text{rot}}$ と並進の損失 $\mathcal{L}_{\text{trans}}$ を重み付き和で組み合わせることです。
$$ \mathcal{L} = \mathcal{L}_{\text{rot}} + \beta \mathcal{L}_{\text{trans}} $$
ここで $\beta$ は回転と並進のバランスを取るハイパーパラメータです。回転は角度（ラジアン）、並進は距離（メートル）で、単位が異なるため、$\beta$ の適切な設定が重要です。
回転の損失としてよく使われるのは、測地線距離（geodesic distance）です。推定回転行列 $\hat{\bm{R}}$ と真の回転行列 $\bm{R}^*$ の間の回転角度は次の式で計算されます。
$$ \mathcal{L}_{\text{rot}} = \arccos\left(\frac{\text{tr}(\hat{\bm{R}}^T \bm{R}^*) – 1}{2}\right) $$
この値は2つの回転の間の「最短回転角度」を表し、$SO(3)$ 上の自然な距離です。
並進の損失にはL2距離がよく使われます。
$$ \mathcal{L}_{\text{trans}} = \|\hat{\bm{t}} – \bm{t}^*\|_2 $$
PoseNetの元論文では、$\beta$ を学習可能にする手法（不確実性に基づく自動重み付け）が提案されています。この場合の損失関数は次のようになります。
$$ \mathcal{L} = \frac{1}{2\sigma_r^2} \mathcal{L}_{\text{rot}} + \frac{1}{2\sigma_t^2} \mathcal{L}_{\text{trans}} + \log \sigma_r + \log \sigma_t $$
$\sigma_r$ と $\sigma_t$ はそれぞれ回転と並進の不確実性を表す学習可能パラメータです。$\log \sigma$ の項は正則化項で、$\sigma$ が無限大に発散する（損失を0にするために両方を無視する）ことを防ぎます。
直接回帰の限界
CNNによる直接6DoF回帰は、画像一枚から瞬時に姿勢を推定できるという大きな利点がありますが、以下の限界が指摘されています。
位置精度の壁: 直接回帰では、画像全体の情報を圧縮した特徴ベクトルから姿勢を復元するため、局所的な幾何学的整合性が失われやすく、古典的なPnPに比べて並進の精度が劣る傾向があります。
一般化の困難: 学習データに含まれない視点や照明条件に対する汎化性能が限定的です。特に宇宙環境では、地上テストとは大きく異なる条件に遭遇するため、これは致命的な問題になり得ます。
不確実性の推定: 標準的な回帰モデルは「どれくらい確信があるか」を出力しないため、推定がどの程度信頼できるかが分かりません。ミッションクリティカルな宇宙応用では、不確実性の定量化が必須です。
これらの限界を克服するために、より洗練されたアプローチが開発されました。その中でも特に成功しているのが、キーポイント検出と古典的なPnPを組み合わせるハイブリッドアプローチです。
キーポイント検出 + PnPによる姿勢推定
なぜキーポイントが有効か
直接回帰の限界を打ち破った重要なアイデアは、「姿勢そのものを直接回帰するのではなく、画像上のキーポイント位置を回帰し、それを使ってPnP問題を解く」という間接的アプローチです。
このアイデアの核心を直感的に理解するために、地図の見方を考えてみましょう。「この建物は北緯35度23分、東経139度44分にある」と直接答えるのは難しいですが、「あの交差点のすぐ右側にある」とランドマーク（目印）との位置関係で答えるのは比較的簡単です。キーポイントベースの手法は、これと同じ発想です。3Dモデル上で意味のある特徴的な点（キーポイント）を定義し、それが画像上のどこに映るかをCNNに予測させます。2Dキーポイント位置と3Dキーポイント位置の対応がわかれば、PnPアルゴリズムで正確な姿勢が計算できます。
パイプラインの概要
キーポイントベースの姿勢推定パイプラインは、以下の3段階で構成されます。
第1段階: 物体検出（バウンディングボックス）
まず、画像全体から対象物体の領域を検出します。YOLOやFaster R-CNNなどの物体検出ネットワークが使われます。宇宙の場合、背景が暗い宇宙空間であるため、前景分離が比較的容易です。
第2段階: キーポイント検出（ヒートマップ回帰）
検出された物体領域をクロップ・リサイズし、各キーポイントの位置をヒートマップとして予測します。ヒートマップとは、画像と同じ空間解像度を持つ2D確率分布で、キーポイントの位置で高い値を取ります。
キーポイント $k$ のヒートマップ $\bm{H}_k(u, v)$ は、真のキーポイント位置 $(u_k^*, v_k^*)$ を中心とするガウス分布で教師データを作ります。
$$ \bm{H}_k^*(u, v) = \exp\left(-\frac{(u – u_k^*)^2 + (v – v_k^*)^2}{2\sigma^2}\right) $$
$\sigma$ はガウス分布の広がりを制御するパラメータで、通常2〜3ピクセルに設定されます。ネットワークは画像からこのヒートマップを予測するように学習します。
予測ヒートマップ $\hat{\bm{H}}_k$ から最終的なキーポイント位置を取得するには、argmaxを取ります。
$$ (\hat{u}_k, \hat{v}_k) = \arg\max_{(u,v)} \hat{\bm{H}}_k(u, v) $$
実際にはサブピクセル精度を得るため、ピーク周辺の値で加重平均（soft-argmax）を計算することが多いです。
$$ \hat{u}_k = \frac{\sum_{u,v} u \cdot \hat{\bm{H}}_k(u,v)}{\sum_{u,v} \hat{\bm{H}}_k(u,v)}, \quad \hat{v}_k = \frac{\sum_{u,v} v \cdot \hat{\bm{H}}_k(u,v)}{\sum_{u,v} \hat{\bm{H}}_k(u,v)} $$
soft-argmaxは微分可能であるため、後段のPnPも含めたエンドツーエンドの学習が可能になります。
第3段階: PnPによる姿勢計算
検出された2Dキーポイント $\{(\hat{u}_k, \hat{v}_k)\}_{k=1}^{K}$ と、3Dモデル上で事前に定義されたキーポイントの3D座標 $\{\bm{P}_k\}_{k=1}^{K}$ の対応からPnP問題を解きます。
$$ \lambda_k \begin{bmatrix} \hat{u}_k \\ \hat{v}_k \\ 1 \end{bmatrix} = \bm{K}(\bm{R}\bm{P}_k + \bm{t}), \quad k = 1, \ldots, K $$
EPnP + RANSACにより、外れ値に頑健な姿勢推定が実現できます。最低4点のキーポイントがあれば解が求まりますが、実際には8〜16点程度のキーポイントを使って冗長性を確保します。
キーポイントの選び方
キーポイントの選定は性能に大きく影響します。宇宙機の場合、以下の方針で選ばれます。
幾何学的に意味のある点: 衛星本体の角点、太陽パネルの端点、アンテナの先端など、物体の形状を特徴づける点が選ばれます。これらの点は異なる視点からも識別しやすく、PnPの数値安定性にも寄与します。
空間的な分散: キーポイントが3D空間で偏っていると、PnPの解の精度が低下します（特に並進のz成分）。衛星本体と太陽パネルの両方にキーポイントを配置し、空間的に分散させることが重要です。
FPS（Farthest Point Sampling）による自動選定: 3Dモデル上の頂点からFPSアルゴリズムで最も離れた点を順次選ぶことで、空間的に均等に分散したキーポイントセットを自動的に生成する手法もあります。
直接回帰との性能比較
キーポイントベースのアプローチが直接回帰より優れている理由は、幾何学的な整合性にあります。直接回帰では、CNNの出力が6つの数値（回転3 + 並進3）であり、画像の空間的な情報が全結合層で失われます。一方、キーポイントベースのアプローチでは、ヒートマップが画像の空間解像度を保持し、PnPが幾何学的な拘束を厳密に適用します。
結果として、キーポイントベースの手法は特に並進の推定精度で直接回帰を大きく上回ります。SPEED+データセット（後述）上のベンチマークでも、上位手法のほぼ全てがキーポイントベースまたはそれに類するアプローチを採用しています。
ここまでで、深層学習による姿勢推定の2つの主要アプローチ（直接回帰とキーポイントベース）の仕組みがわかりました。しかし、これらのニューラルネットワークを学習させるためには、大量の「画像とそれに対応する正解姿勢」のペアが必要です。宇宙での実画像を大量に取得することは極めて困難です。ではどうするか？ 答えは「合成データ」にあります。
合成データによる学習戦略
なぜ合成データが不可欠か
地上の物体認識では、ImageNetのような大規模データセットや、実環境での撮影データが比較的容易に入手できます。しかし宇宙物体の姿勢推定では、以下の理由から実データの取得が極めて困難です。
実宇宙画像の希少性: 非協力物体の近接画像は、これまでに実施された少数のミッション（PRISMA、RemoveDEBRIS等）でしか取得されていません。しかもそれらのデータは量が限られ、特定の衛星形状と照明条件に偏っています。
アノテーションの困難: 仮に画像があったとしても、ピクセル精度のキーポイント位置や正確な6DoFポーズのアノテーション（正解ラベル付け）を行うのは非常に困難です。
多様性の不足: 深層学習の性能は学習データの多様性に大きく依存しますが、実ミッションでは特定の視点・距離・照明条件のデータしか得られません。
これらの制約から、合成データ（CGレンダリング）が宇宙物体の姿勢推定における学習データの主要な供給源となっています。
レンダリングパイプライン
合成データを生成するレンダリングパイプラインは、以下の要素で構成されます。
3Dモデル: 対象衛星のCADモデルまたはワイヤーフレームモデル。表面のマテリアル（MLI箔、太陽電池セル、白色塗装など）のテクスチャ情報も含みます。
カメラモデル: 焦点距離、画像サイズ、歪み係数など、実カメラのパラメータを再現します。
照明モデル: 太陽光の方向、地球のアルベド（反射光）、衛星自身の影を物理ベースのレンダリングで再現します。
姿勢サンプリング: 対象物体の姿勢を$SO(3)$上で均一にサンプリングし、距離も複数変えます。均一な回転のサンプリングには、クォータニオンの一様分布を使います。
均一なランダム回転の生成には、以下の方法が使われます。4つの独立な正規分布 $z_1, z_2, z_3, z_4 \sim \mathcal{N}(0, 1)$ からクォータニオンを構成し、正規化します。
$$ \bm{q} = \frac{1}{\sqrt{z_1^2 + z_2^2 + z_3^2 + z_4^2}} \begin{bmatrix} z_1 \\ z_2 \\ z_3 \\ z_4 \end{bmatrix} $$
これが$SO(3)$上のHaar測度（一様分布）に従うことは、4次元球面上の一様分布が$SO(3)$への二重被覆を通じて回転群上の一様分布を誘導することから示されます。
ドメインランダマイゼーション: 照明方向、背景（星空、地球、深宇宙）、カメラの露出、ノイズレベルなどをランダムに変化させて、学習データの多様性を高めます。
レンダリングの品質レベル
合成データのレンダリング品質は、目的に応じて異なるレベルが採用されます。
簡易レンダリング（ワイヤーフレーム）: OpenGLで衛星モデルの輪郭線のみをレンダリングします。計算が高速で大量のデータを生成できますが、テクスチャや照明の情報が欠落するため、学習したモデルのSim-to-Real汎化性能は限定的です。
中品質レンダリング（ラスタライズ）: OpenGLやBlenderのEeveeエンジンで、テクスチャ付きのレンダリングを行います。フォン反射モデルやPBR（Physics-Based Rendering）マテリアルで表面の光学特性を近似します。計算量とリアリズムのバランスが良く、最も広く使われています。
高品質レンダリング（レイトレーシング）: Blender CyclesやPBRT（Physically Based Rendering Toolkit）でレイトレーシングを行い、MLI箔の複雑な反射、グローバルイルミネーション（間接照明）、カメラの光学特性（レンズフレア、ブルーミング）まで再現します。計算量は大きいですが、実画像に最も近い合成データが得られます。
合成データのみで学習したモデルを実画像に適用すると、必ず性能が低下します。この問題がSim-to-Realドメインギャップです。次にこの重要な課題を掘り下げましょう。
Sim-to-Realドメインギャップ
ドメインギャップとは
合成データ（ソースドメイン）で学習したモデルを実画像（ターゲットドメイン）に適用すると、性能が大幅に低下する現象をドメインギャップと呼びます。これは宇宙物体の姿勢推定における最も重要な課題の一つです。
日常的な例で言えば、教科書の問題だけで勉強した学生が、実際の試験で見慣れない形式の問題に直面して戸惑う状況に似ています。問題の本質（数学の概念）は同じでも、表現（問題の出し方）が異なるために正しく対処できないのです。
宇宙環境でのドメインギャップの主な原因は以下の通りです。
照明のリアリズム: MLI箔（多層断熱材）の複雑な反射特性は、簡易な反射モデル（Lambertian + Phong）では正確に再現できません。MLI箔は非常に複雑なBRDF（双方向反射率分布関数）を持ち、しわの入り方や経年劣化によって反射パターンが変化します。
背景の差異: 合成データの背景は星空の画像を貼り付けたものが多いですが、実際には地球の反射光（アースシャイン）や太陽のグレアなど、複雑な光学現象が画像に影響します。
カメラの光学特性: 実カメラのレンズ収差、色収差、センサーノイズ（ショットノイズ、暗電流ノイズ、読み出しノイズ）は合成データでは完全には再現されません。宇宙放射線によるセンサーのホットピクセルも、地上では模擬が困難です。
3Dモデルの不完全性: CADモデルは設計時の形状を表していますが、実際の衛星は組み立て時の微小な変形、MLI箔のしわ、太陽電池の配線カバーなど、CADには反映されない細部があります。
ドメインギャップの数学的定式化
ドメイン適応の理論的背景を簡単に整理しましょう。ソースドメイン $\mathcal{D}_S$（合成データ）とターゲットドメイン $\mathcal{D}_T$（実画像）の分布が異なるとき、ソースで学習した仮説 $h$ のターゲットにおける誤差 $\epsilon_T(h)$ は、Ben-David et al. (2010) の理論により次のように上界が与えられます。
$$ \epsilon_T(h) \leq \epsilon_S(h) + d_{\mathcal{H}\Delta\mathcal{H}}(\mathcal{D}_S, \mathcal{D}_T) + \lambda^* $$
右辺の第1項 $\epsilon_S(h)$ はソースドメインでの誤差、第2項 $d_{\mathcal{H}\Delta\mathcal{H}}$ は2つのドメインの分布間の $\mathcal{H}\Delta\mathcal{H}$-距離、第3項 $\lambda^*$ はソースとターゲットの両方で同時に低い誤差を達成できる理想的な仮説の誤差です。
この不等式が意味するのは、ターゲットでの性能を高めるためには、(1) ソースでの性能を高め、(2) ドメイン間の分布の差を小さくし、(3) 両ドメインに共通する構造を捉える必要がある、ということです。
ドメインギャップの克服手法
ドメインランダマイゼーション（Domain Randomization）
レンダリング時にテクスチャ、照明、背景、カメラパラメータなどをランダムに大きく変化させることで、「あらゆる条件に対して頑健な特徴量」を学習させる手法です。
直感的には、「極端に多様な合成データを見せれば、実データもその多様性の範囲内に含まれるだろう」という考え方です。具体的には以下のパラメータをランダム化します。
- 照明方向: 太陽光の方向を球面上で一様にサンプリング
- 照明強度: 輝度を$[0.5, 2.0]$倍の範囲でランダムに変動
- 背景: 星空、地球、深宇宙、さらにはランダムテクスチャ
- 物体テクスチャ: 本来のテクスチャに加えてランダムパッチを重畳
- カメラノイズ: ガウスノイズ、ショットノイズの追加
- 色調変化: 色温度、コントラスト、ガンマ値のランダム変動
ドメインランダマイゼーションは実装が簡単で、ターゲットドメインのデータが全く不要であるという大きな利点があります。
敵対的ドメイン適応（Adversarial Domain Adaptation）
GANの仕組みを応用し、ドメイン判別器が合成データと実データを区別できないような特徴表現を学習させる手法です。
特徴抽出器 $G$、タスクヘッド $C$（姿勢推定）、ドメイン判別器 $D$ の3つのネットワークを用い、以下のミニマックスゲームを解きます。
$$ \min_{G, C} \max_D \mathcal{L}_{\text{task}}(G, C; \mathcal{D}_S) – \lambda \mathcal{L}_{\text{domain}}(G, D; \mathcal{D}_S, \mathcal{D}_T) $$
第1項はソースドメインでの姿勢推定の損失、第2項はドメイン判別の損失です。$G$ は $D$ を騙すような（ドメイン不変な）特徴を学習し、$D$ はドメインを正しく判別しようとします。この対抗的な学習により、$G$ が出力する特徴量はドメインに依存しない表現に収束します。
実装では、勾配反転層（Gradient Reversal Layer: GRL）を $G$ と $D$ の間に挿入することで、1つのバックプロパゲーションで上記のミニマックス最適化を近似的に実行できます。GRLは順伝播では恒等写像、逆伝播では勾配の符号を反転させます。
$$ \text{GRL}(\bm{x}) = \bm{x}, \quad \frac{\partial \text{GRL}}{\partial \bm{x}} = -\lambda \bm{I} $$
スタイル変換（Style Transfer）
合成画像を実画像の「見た目」に変換する手法です。CycleGANなどの画像変換ネットワークを用いて、合成画像の低レベルな画像統計量（色調、テクスチャ、ノイズ特性）を実画像に近づけます。
重要なのは、スタイル変換の際に幾何学的な構造を保存することです。姿勢推定では物体の形状情報が本質的に重要であるため、スタイル変換が物体の輪郭やキーポイントの位置を歪めてしまっては本末転倒です。
宇宙分野では、ESAが主催するSPECチャレンジ（後述）で、ドメインギャップの克服が主要な競技課題となっています。次に、この分野のベンチマークデータセットについて見ていきましょう。
SPEED/SPEED+データセットと宇宙姿勢推定ベンチマーク
SPEEDデータセット
SPEED（Spacecraft Pose Estimation Dataset）は、スタンフォード大学のSpace Rendezvous Laboratory（SLAB）が2019年に公開した、宇宙物体の6DoF姿勢推定のためのベンチマークデータセットです。
SPEEDはESAの通信衛星Tangoの模型を使って撮影されたデータで、以下の2種類の画像で構成されています。
合成画像（Synthetic）: Tangoの3Dモデルをレンダリングして生成した画像。太陽光の方向や背景（地球、深宇宙）を変化させた約12,000枚の画像と、正確な6DoFポーズラベルが含まれます。
実験室画像（Lightbox）: 実物大のTango模型を暗室に設置し、制御された照明条件下で撮影した画像。約300枚で、KUKA産業用ロボットにより正確なポーズのground truthが得られています。
SPEEDは宇宙物体の姿勢推定における初の標準ベンチマークとして大きな役割を果たしましたが、合成画像と実験室画像のドメインギャップが比較的小さいという制約がありました。
SPEED+データセット
SPEED+は、SPEEDを拡張して2022年に公開された後継データセットです。ドメインギャップの問題をより現実的に扱えるように設計されています。
SPEED+には以下の3種類の画像が含まれます。
合成画像: Tangoモデルの高品質レンダリング画像（約59,960枚）。様々な照明条件と視点をカバー
実験室画像（Lightbox）: SPEED同様の暗室での撮影画像（約9,531枚）。ただしより多様な照明条件を含む
日照模擬画像（Sunlamp）: 強い平行光（太陽光を模擬）の下で撮影された画像（約19,914枚）。宇宙の極端なコントラスト条件を再現
SPEED+の最大の特徴は、学習用として合成画像のみ（ポーズラベル付き）を提供し、テストはLightboxとSunlampの実験室画像（ラベルなし）で行うというドメイン適応チャレンジの枠組みを採用していることです。これにより、Sim-to-Realドメインギャップの克服が正面から問われます。
SPECチャレンジ
ESAが主催するSPEC（Satellite Pose Estimation Challenge）は、SPEED+データセットを用いた国際コンペティションです。世界中の研究チームがSim-to-Realドメインギャップの問題に取り組み、合成データのみで学習したモデルの実画像上での性能を競います。
評価指標は、回転誤差 $e_{\text{rot}}$（度）と正規化並進誤差 $e_{\text{trans}}$（物体距離で正規化）の加重和です。
$$ \text{Score} = e_{\text{rot}} + e_{\text{trans}} $$
ここで並進誤差は次のように正規化されます。
$$ e_{\text{trans}} = \frac{\|\hat{\bm{t}} – \bm{t}^*\|}{\|\bm{t}^*\|} $$
この正規化により、遠距離（並進の絶対誤差が大きくても相対的には小さい）と近距離（絶対誤差は小さいが相対的には大きい）の姿勢推定を公平に比較できます。
SPECチャレンジの上位チームが共通して使用しているテクニックとして、(1) キーポイントベースのアプローチ、(2) ドメインランダマイゼーション、(3) テスト時オーグメンテーション（TTA）、(4) アンサンブル学習が挙げられます。
このベンチマークの発展に伴い、Transformerアーキテクチャの導入が新たな潮流となっています。次にこの最新の動向を見ていきましょう。
Transformerベースの手法
なぜTransformerか
CNN（畳み込みニューラルネットワーク）は局所的な受容野を持ち、畳み込みとプーリングの繰り返しで徐々に広い範囲の情報を統合します。これは画像認識では非常に効果的ですが、宇宙物体の姿勢推定では遠く離れた部位間の関係が重要になる場合があります。
例えば、衛星の右端の太陽パネルと左端の太陽パネルの相対的な見え方から、衛星本体の姿勢を推定するには、画像の両端の情報を同時に参照する必要があります。CNNでこれを実現するには非常に深いネットワークが必要ですが、Transformerの自己注意機構（Self-Attention）は任意の位置間の依存関係を直接モデル化できます。
Vision Transformer (ViT) による姿勢推定
Vision Transformer（ViT）を姿勢推定に応用する際の基本的な流れを説明します。
パッチ分割: 入力画像 $\bm{I} \in \mathbb{R}^{H \times W \times 3}$ を $P \times P$ ピクセルのパッチに分割し、各パッチを線形射影でトークン（$D$次元ベクトル）に変換します。$N = HW/P^2$ 個のトークンが得られます。
$$ \bm{z}_0 = [\bm{x}_{\text{cls}}; \bm{E}\bm{x}_1^p; \bm{E}\bm{x}_2^p; \ldots; \bm{E}\bm{x}_N^p] + \bm{E}_{\text{pos}} $$
ここで $\bm{E} \in \mathbb{R}^{D \times P^2 \cdot 3}$ はパッチ埋め込み行列、$\bm{x}_{\text{cls}}$ は分類トークン（[CLS]トークン）、$\bm{E}_{\text{pos}}$ は位置埋め込みです。
Self-Attention: 各Transformerブロックでは、Multi-Head Self-Attention（MHSA）がトークン間の関係を学習します。Query $\bm{Q}$、Key $\bm{K}$、Value $\bm{V}$ を用いたAttentionの計算は次の通りです。
$$ \text{Attention}(\bm{Q}, \bm{K}, \bm{V}) = \text{softmax}\left(\frac{\bm{Q}\bm{K}^T}{\sqrt{d_k}}\right)\bm{V} $$
$d_k$ はKeyの次元数で、$1/\sqrt{d_k}$ のスケーリングはsoftmaxの入力が大きくなりすぎることを防ぎます。
この計算の意味を直感的に理解しましょう。Queryは「この位置から何を知りたいか」、Keyは「この位置にはどんな情報があるか」、Valueは「実際の情報内容」を表しています。$\bm{Q}\bm{K}^T$ の内積が大きい位置ペアほど、Attentionの重み（softmax出力）が大きくなり、Valueの情報が強く伝達されます。
姿勢推定ヘッド: [CLS]トークンの最終層出力から、MLP（多層パーセプトロン）で姿勢パラメータを回帰します。あるいは、全トークンを使ってヒートマップを復元し、キーポイントベースの姿勢推定を行うこともできます。
CNN-Transformerハイブリッド
純粋なViTは大量の学習データを必要とする傾向があるため、宇宙物体の姿勢推定ではCNN-Transformerハイブリッドが好まれます。CNNバックボーン（例えばResNetの前半部分）で局所的な特徴を抽出し、Transformerエンコーダでグローバルな関係性を捉えるという構成です。
この構成の利点は以下の通りです。
- CNNの帰納バイアス（局所性、平行移動不変性）を活かして少ないデータでも効率的に学習できる
- Transformerのグローバルな注意機構で遠距離の部位間関係を捉えられる
- ImageNetで事前学習されたCNNの重みを初期値として利用できる
Cross-Attentionによるモデル参照
さらに発展的な手法として、3Dモデルの情報をTransformerに明示的に組み込むCross-Attentionベースのアプローチがあります。
3Dモデルの複数の視点からのレンダリング画像をテンプレートとして用意し、テスト画像のトークンがこれらのテンプレートにCross-Attentionを行うことで、「この画像はどのテンプレートに最も近いか」を学習します。
$$ \text{CrossAttention}(\bm{Q}_{\text{test}}, \bm{K}_{\text{template}}, \bm{V}_{\text{template}}) = \text{softmax}\left(\frac{\bm{Q}_{\text{test}}\bm{K}_{\text{template}}^T}{\sqrt{d_k}}\right)\bm{V}_{\text{template}} $$
テスト画像のQueryがテンプレートのKeyと照合される点が Self-Attention との違いです。この機構により、ネットワークは「この画像中の太陽パネルの見え方は、テンプレート群のうちこの姿勢のものに最も近い」という推論を暗黙的に行います。
Transformerベースの手法は高い精度を実現できますが、計算量が大きいという課題があります。次に、衛星のオンボードコンピュータでリアルタイム推論を行うための制約と対策を見ていきましょう。
オンボード推論の計算資源制約
宇宙用コンピュータの現実
宇宙で使用するコンピュータは、地上のワークステーションとは比較にならないほど制約があります。
放射線耐性: 宇宙放射線（高エネルギー粒子、銀河宇宙線）がシリコンチップに当たると、ビット反転（Single Event Upset: SEU）やラッチアップ（Single Event Latchup: SEL）を引き起こします。このため、放射線耐性のある（rad-hard）プロセッサを使うか、COTS（Commercial Off-The-Shelf: 市販品）プロセッサに放射線対策を施す必要があります。
代表的なオンボードプロセッサ: – RAD750（BAE Systems）: PPC750ベース、~200 MHz、~400 MIPS。火星探査ローバーCuriosityに搭載。深層学習の実行は現実的ではありません – Xilinx Zynq UltraScale+: ARM Cortex-A53 + FPGA。小型衛星で広く採用されつつあり、FPGAアクセラレータで推論を高速化できます – Intel Movidius Myriad 2/X: ビジョン処理に特化した低消費電力プロセッサ。ESAのPhiSatミッションで宇宙実証済み – NVIDIA Jetson: GPUを搭載したエッジコンピュータ。CubeSatでの実証が進行中
消費電力制約: 衛星の電力バジェットは厳しく、特に小型衛星（CubeSat）では数W〜数十W程度です。姿勢推定に使えるのはその一部であり、GPUが数十Wを消費すると全体の電力バジェットを圧迫します。
メモリ制約: オンボードメモリ（RAM）は数百MB〜数GB程度であり、ResNet-50（約100MBのパラメータ）ですら余裕を持って搭載するのが困難な場合があります。
モデルの軽量化手法
これらの制約の下で深層学習モデルをオンボードで実行するために、いくつかの軽量化手法が適用されます。
知識蒸留（Knowledge Distillation）: 大きな「教師モデル」（Teacher）の出力を、小さな「生徒モデル」（Student）に学習させる手法です。教師モデルのsoftmax出力（ソフトラベル）には、各クラスの間の類似度という「暗黙知」が含まれており、これを生徒モデルに転移させることで、小さなモデルでも大きなモデルに近い性能が得られます。
温度付きsoftmaxを使った蒸留損失は次のように定義されます。
$$ \mathcal{L}_{\text{KD}} = \tau^2 \cdot \text{KL}\left(\text{softmax}\left(\frac{\bm{z}_T}{\tau}\right) \middle\| \text{softmax}\left(\frac{\bm{z}_S}{\tau}\right)\right) $$
$\bm{z}_T$ と $\bm{z}_S$ はそれぞれ教師と生徒のlogit、$\tau$ は温度パラメータです。$\tau > 1$ でsoftmax分布を「滑らか」にすることで、教師の暗黙知がより効果的に伝達されます。$\tau^2$ の係数は、温度スケーリングによる勾配のスケール変化を補正するためです。
量子化（Quantization）: ネットワークの重みと活性化を32ビット浮動小数点（FP32）から8ビット整数（INT8）やさらに低ビット（INT4）に変換します。メモリ使用量と計算量を大幅に削減できますが、精度の低下に注意が必要です。
重み $w$ のINT8量子化は次のように行われます。
$$ w_q = \text{round}\left(\frac{w}{s}\right) + z $$
ここで $s$ はスケールファクタ、$z$ はゼロポイントです。逆量子化（推論結果をFP32に戻す）は $\hat{w} = s \cdot (w_q – z)$ で行われます。
プルーニング（Pruning）: ニューラルネットワークの中で、推定精度への寄与が小さい重み（ゼロに近い値）を除去することで、モデルを疎（スパース）にします。構造化プルーニング（チャネルやフィルタ単位での削除）は、実際のハードウェア上での高速化に直結します。
アーキテクチャ設計: MobileNet、EfficientNet、ShuffleNetなど、モバイル向けに設計されたバックボーンを使うことで、最初から計算量の少ないモデルを構築します。Depthwise separable convolutionにより、標準的な畳み込みに比べてパラメータ数と計算量を $1/k^2$（$k$ はカーネルサイズ）に削減できます。
リアルタイム推論の要件
ランデブー・ドッキングのGNC（誘導・航法・制御）ループで使用する場合、姿勢推定の推論周波数は1〜10 Hzが典型的な要件です。これは、カメラのフレームレートや制御系の帯域幅に依存します。
1 Hzの推論周波数は、一見すると低いように思えるかもしれません。しかし、Jetson Nanoクラスのプロセッサ上でResNet-50の推論を行うと約50〜100 ms（10〜20 Hz）程度であり、キーポイント検出のヒートマップ回帰を含めると1 Hz以上の達成は十分に現実的です。さらにFPGA上の専用ハードウェアで推論を行えば、消費電力を抑えつつ10 Hz以上の推論も可能です。
ここまでの理論と技術を踏まえて、Pythonで宇宙物体の姿勢推定パイプラインを実装してみましょう。
Pythonでの実装
実装の全体構成
ここでは、以下の一連のパイプラインをPythonで実装します。
- 衛星3Dモデルの定義とキーポイント設定
- 合成画像データの生成（レンダリングの簡易模擬）
- CNNによるキーポイント検出ネットワークの構築と学習
- PnPによる姿勢復元
- 性能評価と可視化
まず、衛星モデルとデータ生成を行います。
import numpy as np
from scipy.spatial.transform import Rotation
import matplotlib.pyplot as plt
# --- 衛星の3Dモデル定義 ---
# 直方体ボディ (1.0 x 0.8 x 0.6 m)
body_half = np.array([0.5, 0.4, 0.3])
body_vertices = np.array([
[-1, -1, -1], [-1, -1, 1], [-1, 1, -1], [-1, 1, 1],
[ 1, -1, -1], [ 1, -1, 1], [ 1, 1, -1], [ 1, 1, 1],
], dtype=float) * body_half
# 太陽パネル（両翼）
panel_l = np.array([[-0.3, -1.2, 0.0], [0.3, -1.2, 0.0],
[-0.3, -0.4, 0.0], [0.3, -0.4, 0.0]])
panel_r = panel_l.copy()
panel_r[:, 1] = -panel_l[:, 1]
# アンテナ先端
antenna = np.array([[0.0, 0.0, 0.5]])
# キーポイント: ボディ8角 + パネル端4x2 + アンテナ1 = 17点
keypoints_3d = np.vstack([body_vertices, panel_l, panel_r, antenna])
n_kp = len(keypoints_3d)
print(f"キーポイント数: {n_kp}")
print(f"キーポイントの空間範囲:")
print(f" x: [{keypoints_3d[:,0].min():.2f}, {keypoints_3d[:,0].max():.2f}] m")
print(f" y: [{keypoints_3d[:,1].min():.2f}, {keypoints_3d[:,1].max():.2f}] m")
print(f" z: [{keypoints_3d[:,2].min():.2f}, {keypoints_3d[:,2].max():.2f}] m")
# カメラ内部パラメータ
img_w, img_h = 256, 256
fx, fy = 500.0, 500.0
cx, cy = img_w / 2.0, img_h / 2.0
K = np.array([[fx, 0, cx],
[0, fy, cy],
[0, 0, 1]])
このコードでは、直方体ボディ、2枚の太陽パネル、1本のアンテナで構成された簡易衛星モデルを定義し、17個のキーポイントを設定しています。キーポイントが空間的に分散していることが重要で、x方向に$\pm 0.5$ m、y方向に$\pm 1.2$ m、z方向に$-0.3$〜$0.5$ mの範囲をカバーしています。太陽パネルが最も広い空間分布を持っており、並進（特にz成分）の推定精度に寄与します。
次に、合成学習データを生成する関数を実装します。
import numpy as np
from scipy.spatial.transform import Rotation
def generate_synthetic_data(keypoints_3d, K, n_samples,
dist_range=(3.0, 8.0),
noise_std=1.5,
seed=42):
"""合成学習データを生成する
Parameters
----------
keypoints_3d : ndarray, shape (n_kp, 3)
3Dキーポイント座標
K : ndarray, shape (3, 3)
カメラ内部パラメータ
n_samples : int
生成するサンプル数
dist_range : tuple
物体までの距離の範囲 [m]
noise_std : float
キーポイント検出のノイズ標準偏差 [px]
seed : int
乱数シード
Returns
-------
data : list of dict
各サンプルの情報（R_true, t_true, kp_2d_noisy, kp_2d_true）
"""
rng = np.random.RandomState(seed)
fx, fy = K[0, 0], K[1, 1]
cx, cy = K[0, 2], K[1, 2]
data = []
for i in range(n_samples):
# SO(3)上で一様な回転をサンプリング
q = rng.randn(4)
q /= np.linalg.norm(q)
R = Rotation.from_quat(q).as_matrix()
# 距離をランダムにサンプリング
dist = rng.uniform(*dist_range)
# 物体中心をカメラの視野内に配置
tx = rng.uniform(-0.5, 0.5)
ty = rng.uniform(-0.5, 0.5)
t = np.array([tx, ty, dist])
# 3D→2D射影
pts_cam = (R @ keypoints_3d.T).T + t
kp_2d = np.zeros((len(keypoints_3d), 2))
kp_2d[:, 0] = fx * pts_cam[:, 0] / pts_cam[:, 2] + cx
kp_2d[:, 1] = fy * pts_cam[:, 1] / pts_cam[:, 2] + cy
# 全キーポイントが画像内にあるか確認
in_bounds = np.all((kp_2d[:, 0] >= 0) & (kp_2d[:, 0] < 256) &
(kp_2d[:, 1] >= 0) & (kp_2d[:, 1] < 256))
if not in_bounds:
continue
# ノイズ付きキーポイント
kp_noisy = kp_2d + rng.randn(*kp_2d.shape) * noise_std
data.append({
'R_true': R,
't_true': t,
'kp_2d_true': kp_2d,
'kp_2d_noisy': kp_noisy,
})
return data
# 学習データとテストデータを生成
train_data = generate_synthetic_data(keypoints_3d, K, n_samples=5000, seed=42)
test_data = generate_synthetic_data(keypoints_3d, K, n_samples=500, seed=123)
print(f"学習データ: {len(train_data)} サンプル")
print(f"テストデータ: {len(test_data)} サンプル")
このコードでは、$SO(3)$ 上の一様回転サンプリングと、3〜8mの距離範囲でのランダムな物体配置により、合成データを生成しています。ノイズ標準偏差1.5ピクセルは、実際のキーポイント検出ネットワークの典型的な精度に近い値です。全キーポイントが画像内に収まるサンプルのみを保持しているため、生成サンプル数と実際のデータ数は異なります。
次に、キーポイント検出を模擬したヒートマップ回帰ネットワークの学習と、PnPによる姿勢推定を実装します。
import numpy as np
from scipy.spatial.transform import Rotation
from scipy.optimize import least_squares
def solve_pnp_lm(keypoints_3d, keypoints_2d, K, R_init=None, t_init=None):
"""Levenberg-Marquardt法でPnP問題を解く
Parameters
----------
keypoints_3d : ndarray, shape (n, 3)
keypoints_2d : ndarray, shape (n, 2)
K : ndarray, shape (3, 3)
R_init, t_init : 初期推定値（Noneの場合はDLTで初期化）
Returns
-------
R_est, t_est : 推定された回転行列と並進ベクトル
"""
fx, fy = K[0, 0], K[1, 1]
cx, cy = K[0, 2], K[1, 2]
if R_init is None:
# 簡易初期化: 正規化画像座標から大雑把に推定
R_init = np.eye(3)
# 平均的なz座標を仮定
mean_u = np.mean(keypoints_2d[:, 0])
mean_v = np.mean(keypoints_2d[:, 1])
t_init = np.array([(mean_u - cx) / fx * 5.0,
(mean_v - cy) / fy * 5.0,
5.0])
rvec_init = Rotation.from_matrix(R_init).as_rotvec()
params_init = np.concatenate([rvec_init, t_init])
def residuals(params):
R = Rotation.from_rotvec(params[:3]).as_matrix()
t = params[3:6]
pts_cam = (R @ keypoints_3d.T).T + t
proj_x = fx * pts_cam[:, 0] / pts_cam[:, 2] + cx
proj_y = fy * pts_cam[:, 1] / pts_cam[:, 2] + cy
res = np.zeros(2 * len(keypoints_3d))
res[0::2] = proj_x - keypoints_2d[:, 0]
res[1::2] = proj_y - keypoints_2d[:, 1]
return res
result = least_squares(residuals, params_init, method='lm')
R_est = Rotation.from_rotvec(result.x[:3]).as_matrix()
t_est = result.x[3:6]
return R_est, t_est
# --- 全テストデータで姿勢推定を実行 ---
rot_errors = []
trans_errors = []
for sample in test_data:
R_est, t_est = solve_pnp_lm(
keypoints_3d, sample['kp_2d_noisy'], K
)
# 回転誤差（度）
dR = R_est @ sample['R_true'].T
rot_err = np.degrees(Rotation.from_matrix(dR).magnitude())
rot_errors.append(rot_err)
# 並進誤差（m）
trans_err = np.linalg.norm(t_est - sample['t_true'])
trans_errors.append(trans_err)
rot_errors = np.array(rot_errors)
trans_errors = np.array(trans_errors)
print(f"回転誤差 — 中央値: {np.median(rot_errors):.4f} deg, "
f"平均: {np.mean(rot_errors):.4f} deg")
print(f"並進誤差 — 中央値: {np.median(trans_errors):.6f} m, "
f"平均: {np.mean(trans_errors):.6f} m")
このコードでは、ノイズ付きの2Dキーポイントから Levenberg-Marquardt法によるPnP問題を解いて姿勢を推定し、全テストデータにわたる精度を評価しています。17点のキーポイントと1.5ピクセルのノイズがある場合、中央値で約0.1度以下の回転精度と数mm以下の並進精度が期待されます。これは、キーポイント検出の精度が姿勢推定の最終精度を支配することを示しています。
続いて、CNNによるキーポイントヒートマップ回帰ネットワークの簡易実装を示します。
import numpy as np
import matplotlib.pyplot as plt
def generate_heatmaps(kp_2d, img_size=256, sigma=3.0):
"""キーポイントからヒートマップ教師データを生成
Parameters
----------
kp_2d : ndarray, shape (n_kp, 2)
キーポイントの2D座標
img_size : int
画像サイズ
sigma : float
ガウシアンの標準偏差 [px]
Returns
-------
heatmaps : ndarray, shape (n_kp, img_size, img_size)
"""
n_kp = len(kp_2d)
heatmaps = np.zeros((n_kp, img_size, img_size))
yy, xx = np.mgrid[0:img_size, 0:img_size]
for k in range(n_kp):
ux, uy = kp_2d[k]
heatmaps[k] = np.exp(-((xx - ux)**2 + (yy - uy)**2) / (2 * sigma**2))
return heatmaps
def decode_heatmaps(heatmaps):
"""ヒートマップからキーポイント座標をソフトargmaxで復元
Parameters
----------
heatmaps : ndarray, shape (n_kp, H, W)
Returns
-------
kp_2d : ndarray, shape (n_kp, 2)
"""
n_kp, H, W = heatmaps.shape
kp_2d = np.zeros((n_kp, 2))
yy, xx = np.mgrid[0:H, 0:W]
for k in range(n_kp):
hm = heatmaps[k]
hm_sum = hm.sum()
if hm_sum > 1e-8:
kp_2d[k, 0] = (hm * xx).sum() / hm_sum # u
kp_2d[k, 1] = (hm * yy).sum() / hm_sum # v
else:
# フォールバック: argmax
idx = np.unravel_index(hm.argmax(), hm.shape)
kp_2d[k, 0] = idx[1] # u
kp_2d[k, 1] = idx[0] # v
return kp_2d
# --- ヒートマップの生成と復元テスト ---
sample = train_data[0]
heatmaps = generate_heatmaps(sample['kp_2d_true'], img_size=256, sigma=3.0)
kp_recovered = decode_heatmaps(heatmaps)
# 復元精度の確認
recovery_error = np.sqrt(np.mean((kp_recovered - sample['kp_2d_true'])**2))
print(f"ヒートマップ→ソフトargmax 復元誤差: {recovery_error:.6f} px")
# ヒートマップの可視化
fig, axes = plt.subplots(2, 4, figsize=(16, 8))
for i in range(8):
ax = axes[i // 4, i % 4]
ax.imshow(heatmaps[i], cmap='hot', interpolation='bilinear')
ax.plot(sample['kp_2d_true'][i, 0], sample['kp_2d_true'][i, 1],
'c+', markersize=10, markeredgewidth=2)
ax.set_title(f'Keypoint {i}')
ax.set_xlim(0, 255)
ax.set_ylim(255, 0)
plt.suptitle('Ground Truth Heatmaps (first 8 keypoints)', fontsize=14)
plt.tight_layout()
plt.savefig('heatmap_examples.png', dpi=150, bbox_inches='tight')
plt.show()
上の可視化から、各キーポイントのヒートマップがガウシアン分布として正しく生成されていることが確認できます。シアンの十字マーカーがキーポイントの真の位置を示しており、ヒートマップのピークと一致しています。ヒートマップのソフトargmaxによる復元誤差がほぼ0ピクセルであることも確認でき、ヒートマップ表現がサブピクセル精度のキーポイント検出に適していることがわかります。
実際のCNNネットワークでは、このようなヒートマップを画像から回帰するように学習します。バックボーンとしてResNetやHRNet（High-Resolution Network）が使われ、アップサンプリングにより入力と同じ解像度のヒートマップを出力します。
次に、ドメインギャップの影響をシミュレーションし、ドメインランダマイゼーションの効果を確認する実験を行います。
import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial.transform import Rotation
from scipy.optimize import least_squares
def add_domain_shift(kp_2d, shift_type='systematic', rng=None):
"""ドメインシフトを模擬するノイズを加える
Parameters
----------
kp_2d : ndarray, shape (n_kp, 2)
shift_type : str
'systematic' - 系統的なバイアス（レンズ歪みの差）
'illumination' - 照明変化による検出位置のずれ
'combined' - 両方
rng : RandomState
Returns
-------
kp_shifted : ndarray, shape (n_kp, 2)
"""
if rng is None:
rng = np.random.RandomState(0)
kp_shifted = kp_2d.copy()
cx, cy = 128.0, 128.0 # 画像中心
if shift_type in ['systematic', 'combined']:
# 放射状の歪みバイアス（合成と実カメラの歪みの差を模擬）
dx = kp_2d[:, 0] - cx
dy = kp_2d[:, 1] - cy
r2 = dx**2 + dy**2
k1 = 1e-6 # 歪み係数
kp_shifted[:, 0] += dx * k1 * r2
kp_shifted[:, 1] += dy * k1 * r2
if shift_type in ['illumination', 'combined']:
# 照明変化による非対称ノイズ（影の方向にキーポイントがずれる）
light_dir = rng.randn(2) * 2.0
kp_shifted += light_dir
# ランダムノイズ（検出の本質的なノイズ）
kp_shifted += rng.randn(*kp_shifted.shape) * 1.5
return kp_shifted
def evaluate_with_domain_shift(test_data, keypoints_3d, K, shift_type):
"""ドメインシフト下での姿勢推定精度を評価"""
rot_errors = []
trans_errors = []
rng = np.random.RandomState(99)
for sample in test_data:
kp_shifted = add_domain_shift(
sample['kp_2d_true'], shift_type=shift_type, rng=rng
)
R_est, t_est = solve_pnp_lm(keypoints_3d, kp_shifted, K)
dR = R_est @ sample['R_true'].T
rot_err = np.degrees(Rotation.from_matrix(dR).magnitude())
rot_errors.append(rot_err)
trans_err = np.linalg.norm(t_est - sample['t_true'])
trans_errors.append(trans_err)
return np.array(rot_errors), np.array(trans_errors)
# --- 各ドメインシフト条件での評価 ---
conditions = ['systematic', 'illumination', 'combined']
results = {}
# ベースライン（ドメインシフトなし、ノイズのみ）
rot_base, trans_base = [], []
for sample in test_data:
R_est, t_est = solve_pnp_lm(
keypoints_3d, sample['kp_2d_noisy'], K
)
dR = R_est @ sample['R_true'].T
rot_base.append(np.degrees(Rotation.from_matrix(dR).magnitude()))
trans_base.append(np.linalg.norm(t_est - sample['t_true']))
results['baseline'] = (np.array(rot_base), np.array(trans_base))
for cond in conditions:
rot_err, trans_err = evaluate_with_domain_shift(
test_data, keypoints_3d, K, cond
)
results[cond] = (rot_err, trans_err)
print(f"{cond:15s} — 回転中央値: {np.median(rot_err):.4f} deg, "
f"並進中央値: {np.median(trans_err):.6f} m")
# --- 結果の可視化 ---
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
labels = ['Baseline\n(noise only)', 'Systematic\nbias', 'Illumination\nshift', 'Combined']
keys = ['baseline'] + conditions
# 回転誤差のボックスプロット
bp1 = axes[0].boxplot(
[results[k][0] for k in keys],
labels=labels, patch_artist=True,
medianprops=dict(color='red', linewidth=2)
)
colors = ['#2ecc71', '#3498db', '#e67e22', '#e74c3c']
for patch, color in zip(bp1['boxes'], colors):
patch.set_facecolor(color)
patch.set_alpha(0.6)
axes[0].set_ylabel('Rotation Error [deg]')
axes[0].set_title('Domain Shift Effect on Rotation')
axes[0].grid(True, alpha=0.3)
# 並進誤差のボックスプロット
bp2 = axes[1].boxplot(
[results[k][1] for k in keys],
labels=labels, patch_artist=True,
medianprops=dict(color='red', linewidth=2)
)
for patch, color in zip(bp2['boxes'], colors):
patch.set_facecolor(color)
patch.set_alpha(0.6)
axes[1].set_ylabel('Translation Error [m]')
axes[1].set_title('Domain Shift Effect on Translation')
axes[1].grid(True, alpha=0.3)
plt.suptitle('Impact of Domain Gap on Pose Estimation', fontsize=14, y=1.02)
plt.tight_layout()
plt.savefig('domain_gap_analysis.png', dpi=150, bbox_inches='tight')
plt.show()
上のボックスプロットから、ドメインギャップの影響が定量的に読み取れます。
-
ベースライン（ノイズのみ）: ランダムノイズのみの場合、回転誤差と並進誤差はともに小さく、PnPが正常に機能しています。これが合成データ上でのテスト（ドメインギャップなし）に相当します。
-
系統的バイアス: 放射状のレンズ歪みバイアスが加わると、画像の周辺部にあるキーポイントほど大きくずれるため、特に太陽パネルの端点で誤差が増大します。並進のz成分（奥行き）の精度劣化が顕著です。
-
照明シフト: 照明方向の変化で全キーポイントが一方向にシフトすると、並進の推定にバイアスが生じます。回転への影響は比較的小さく、これは照明シフトが全点を同じ方向にずらすため、相対的な位置関係が保たれるためです。
-
組み合わせ: 系統的バイアスと照明シフトの両方が存在すると、誤差が最も大きくなります。これが実際のSim-to-Realドメインギャップに最も近い状況であり、回転誤差・並進誤差ともにベースラインの数倍に増大することがわかります。
この結果は、合成データのみで学習したモデルを実環境に適用する際に、ドメイン適応が不可欠であることを強く示唆しています。
Transformerの注意機構の可視化
最後に、Transformerベースの姿勢推定におけるSelf-Attentionの動作を模擬的に可視化します。
import numpy as np
import matplotlib.pyplot as plt
def simulate_attention(keypoints_2d, n_heads=4, d_model=64, seed=42):
"""Self-Attentionの動作を模擬的に可視化する
キーポイント位置をトークンとして扱い、位置ベースのAttentionを計算
"""
rng = np.random.RandomState(seed)
n_tokens = len(keypoints_2d)
# 位置埋め込み（キーポイントの2D座標を高次元に射影）
W_pos = rng.randn(2, d_model) * 0.1
token_embeds = keypoints_2d @ W_pos # (n_tokens, d_model)
# Multi-Head Attention
d_k = d_model // n_heads
attention_maps = []
for h in range(n_heads):
W_Q = rng.randn(d_model, d_k) * (1.0 / np.sqrt(d_k))
W_K = rng.randn(d_model, d_k) * (1.0 / np.sqrt(d_k))
Q = token_embeds @ W_Q
Keys = token_embeds @ W_K
scores = Q @ Keys.T / np.sqrt(d_k)
# Softmax
scores_exp = np.exp(scores - scores.max(axis=1, keepdims=True))
attn = scores_exp / scores_exp.sum(axis=1, keepdims=True)
attention_maps.append(attn)
return attention_maps
# テストデータの1サンプルで可視化
sample = test_data[0]
kp_2d = sample['kp_2d_true']
attn_maps = simulate_attention(kp_2d, n_heads=4)
fig, axes = plt.subplots(2, 2, figsize=(12, 10))
kp_labels = [f'B{i}' for i in range(8)] + \
[f'L{i}' for i in range(4)] + \
[f'R{i}' for i in range(4)] + ['Ant']
for h in range(4):
ax = axes[h // 2, h % 2]
im = ax.imshow(attn_maps[h], cmap='Blues', vmin=0)
ax.set_xticks(range(17))
ax.set_yticks(range(17))
ax.set_xticklabels(kp_labels, rotation=45, fontsize=7)
ax.set_yticklabels(kp_labels, fontsize=7)
ax.set_title(f'Attention Head {h+1}')
plt.colorbar(im, ax=ax, fraction=0.046)
plt.suptitle('Simulated Self-Attention Maps\n'
'(B=Body, L=Left Panel, R=Right Panel, Ant=Antenna)',
fontsize=13)
plt.tight_layout()
plt.savefig('attention_visualization.png', dpi=150, bbox_inches='tight')
plt.show()
上の注意マップの可視化から、Self-Attentionがキーポイント間の関係性を捉えていることが読み取れます。この模擬実験では位置ベースの単純な埋め込みを使用していますが、いくつかの興味深いパターンが観察されます。
-
ボディ角点間の強い注意: ボディの8つの角点（B0〜B7）間で比較的均一な注意が見られます。これは、ボディの直方体形状を認識するために各角点が互いを参照していることに対応します。
-
パネルとボディ間の交差注意: 太陽パネルのキーポイント（L0〜L3, R0〜R3）がボディのキーポイントに注意を向けるパターンが見られます。実際のネットワークでは、パネルの見え方からボディの姿勢を推測する（パネルが正面を向いているならボディもこちらを向いているはず）という推論に対応します。
-
ヘッド間の多様性: 4つの注意ヘッドがそれぞれ異なるパターンを示しており、Multi-Head Attentionが複数の観点からキーポイント間関係を捉える能力を持つことが確認できます。
距離別の性能評価
宇宙でのランデブーでは、遠距離から近距離まで姿勢推定が必要です。距離による性能変化を評価しましょう。
import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial.transform import Rotation
from scipy.optimize import least_squares
# --- 距離別の姿勢推定精度を評価 ---
dist_bins = [(3.0, 4.0), (4.0, 5.0), (5.0, 6.0), (6.0, 7.0), (7.0, 8.0)]
noise_levels = [1.0, 1.5, 2.0, 3.0] # ピクセル
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
for noise_std in noise_levels:
median_rot = []
median_trans = []
dist_centers = []
for d_min, d_max in dist_bins:
# この距離帯のテストデータを生成
test_dist = generate_synthetic_data(
keypoints_3d, K, n_samples=1000,
dist_range=(d_min, d_max),
noise_std=noise_std,
seed=int(d_min * 100 + noise_std * 10)
)
rot_errs = []
trans_errs = []
for sample in test_dist[:200]: # 各距離帯で200サンプル
R_est, t_est = solve_pnp_lm(
keypoints_3d, sample['kp_2d_noisy'], K
)
dR = R_est @ sample['R_true'].T
rot_errs.append(
np.degrees(Rotation.from_matrix(dR).magnitude())
)
trans_errs.append(
np.linalg.norm(t_est - sample['t_true'])
)
median_rot.append(np.median(rot_errs))
median_trans.append(np.median(trans_errs))
dist_centers.append((d_min + d_max) / 2)
axes[0].plot(dist_centers, median_rot, 'o-',
label=f'noise={noise_std:.1f}px')
axes[1].plot(dist_centers, median_trans, 'o-',
label=f'noise={noise_std:.1f}px')
axes[0].set_xlabel('Distance [m]')
axes[0].set_ylabel('Median Rotation Error [deg]')
axes[0].set_title('Rotation Error vs Distance')
axes[0].legend()
axes[0].grid(True, alpha=0.3)
axes[1].set_xlabel('Distance [m]')
axes[1].set_ylabel('Median Translation Error [m]')
axes[1].set_title('Translation Error vs Distance')
axes[1].legend()
axes[1].grid(True, alpha=0.3)
plt.suptitle('Pose Estimation Accuracy vs Distance and Noise Level',
fontsize=13, y=1.02)
plt.tight_layout()
plt.savefig('distance_noise_analysis.png', dpi=150, bbox_inches='tight')
plt.show()
上のグラフから、距離とノイズレベルが姿勢推定精度に与える影響が明確に読み取れます。
-
距離が増すと精度が低下する: 遠距離ではキーポイントの画像上の間隔が狭くなるため、同じピクセルノイズでも相対的な影響が大きくなります。特に並進のz成分（奥行き）の精度劣化が顕著で、これは射影幾何の本質的な性質です。ピンホールカメラモデルでは、$z$ 方向の分解能が $z^2$ に比例して悪化するためです。
-
ノイズレベルの影響は線形的: キーポイント検出ノイズが1.0 pxから3.0 pxに増加すると、姿勢推定誤差もほぼ比例して増大します。これは、PnPの解がキーポイント位置の摂動に対して線形的に応答することを反映しています。
-
回転は比較的ロバスト: 並進に比べて回転の精度劣化は穏やかです。回転の推定は主にキーポイントの相対的な配置パターンに依存するため、ノイズの影響を受けにくいのです。
-
実運用への示唆: この結果から、遠距離フェーズでは粗い姿勢推定で十分だが、近距離での精密ドッキングにはキーポイント検出の高精度化が不可欠であることがわかります。これは実際のミッション設計で、距離フェーズに応じてアルゴリズムを切り替える戦略の根拠となっています。
最新の研究動向と今後の展望
自己教師あり学習
ラベル付き合成データへの依存を減らすために、自己教師あり学習のアプローチが注目されています。画像の幾何学的整合性（例えば、同一物体の異なる視点からの画像が整合的なポーズを与えるべきという拘束）を自己教師信号として利用し、ラベルなしの実画像からもポーズ推定を学習する試みです。
具体的には、微分可能レンダリング（Differentiable Rendering）を用いて、推定ポーズで3Dモデルをレンダリングした画像と実画像の整合性を損失関数とする方法があります。
$$ \mathcal{L}_{\text{render}} = \sum_{\bm{p}} \| I_{\text{real}}(\bm{p}) – I_{\text{render}}(\bm{p}; \hat{\bm{R}}, \hat{\bm{t}}) \|^2 $$
この損失関数は姿勢パラメータ $(\hat{\bm{R}}, \hat{\bm{t}})$ に関して微分可能であるため、レンダリングパイプライン全体を通じてバックプロパゲーションが可能です。
ニューラル暗黙表現（NeRF）の活用
NeRF（Neural Radiance Fields）やその派生手法を宇宙物体のモデリングに応用する研究も始まっています。少数の実画像からNeRFで物体の3D表現を学習し、任意の視点からの画像を合成できるようにすることで、合成データと実データのギャップを本質的に解消する可能性があります。
NeRFは3D空間の各点 $(\bm{x}, \bm{d})$（位置と視線方向）に対して、色 $\bm{c}$ と密度 $\sigma$ を出力するニューラルネットワークです。
$$ F_\theta: (\bm{x}, \bm{d}) \mapsto (\bm{c}, \sigma) $$
ボリュームレンダリングにより、任意のカメラポーズからの画像を生成できます。これを宇宙物体に適用すれば、MLI箔の複雑な反射特性を含む実物に忠実な画像が合成でき、ドメインギャップの問題を根本から解決できる可能性があります。
Foundation Modelの宇宙応用
DINOv2やSAM（Segment Anything Model）のような大規模事前学習モデル（Foundation Model）が、宇宙物体にも有効であるかの検証が進んでいます。これらのモデルは膨大な自然画像で学習されているため、宇宙画像という全く異なるドメインでどの程度の転移性能が得られるかが焦点です。
初期的な研究結果は、DINOv2の特徴量が宇宙物体の粗い姿勢分類に有効であることを示していますが、ミリメートル精度の姿勢推定にはファインチューニングが不可欠であると報告されています。
オンボードAIプロセッサの進化
宇宙用AIプロセッサの開発も急速に進んでいます。ESAのOPS-SAT衛星やNASAのSTP-H9実験では、COTS GPUの宇宙環境での動作実証が行われています。今後、より高性能なエッジAIプロセッサが宇宙で利用可能になれば、Transformerベースの大規模モデルもオンボードで実行可能になるでしょう。
まとめ
本記事では、深層学習による宇宙物体の認識と姿勢推定について、CNNによる直接回帰からTransformerベースの最新手法まで、幅広く解説しました。
- CNNによる直接6DoF回帰は最も素朴なアプローチですが、位置精度や汎化性能に限界があります。回転表現として6D表現を用いることで、学習の安定性が向上します
- キーポイント検出 + PnPのハイブリッドアプローチが現在の主流であり、ヒートマップ回帰による高精度なキーポイント検出と古典的なPnPの幾何学的拘束を組み合わせることで、高い姿勢推定精度を実現します
- 合成データ（CGレンダリング）が学習データの主要な供給源であり、$SO(3)$上の一様サンプリングやドメインランダマイゼーションで多様なデータを生成します
- Sim-to-Realドメインギャップは最大の課題であり、ドメインランダマイゼーション、敵対的ドメイン適応、スタイル変換などの手法で対処されます。SPEED+データセットとSPECチャレンジがこの問題の標準ベンチマークです
- Transformerの自己注意機構は、衛星の遠く離れた部位間の関係（太陽パネルと本体の位置関係など）を直接モデル化でき、CNN-Transformerハイブリッドアーキテクチャが有望です
- オンボード推論では、知識蒸留、量子化、プルーニングなどの軽量化手法が不可欠であり、宇宙用AIプロセッサの進化とともにリアルタイム推定の実用化が進んでいます
深層学習による姿勢推定は、古典手法の初期値問題とロバスト性の限界を克服し、宇宙での自律的なランデブー・捕獲ミッションの実現に大きく貢献する技術です。今後の研究では、自己教師あり学習やNeRFによるドメインギャップの解消、Foundation Modelの宇宙応用が重要なテーマとなるでしょう。
次のステップとして、以下の記事も参考にしてください。