太陽センサの原理 — 太陽方向ベクトルの検出と姿勢への利用

人工衛星が宇宙空間で「自分がどちらを向いているか」を知ることは、ミッション成功の大前提です。通信アンテナを地上局に向ける、太陽電池パネルを太陽に向ける、観測カメラを地表に向ける — いずれも、衛星がまず自分の姿勢を正しく把握していなければ実現できません。では、宇宙空間に浮かぶ衛星はどのようにして方位の基準を得るのでしょうか?

最も直感的な答えの一つが「太陽を見る」ことです。太陽は宇宙空間で圧倒的に明るい天体であり、その方向は地球暦から高精度に計算できます。衛星の機体座標系から見た太陽の方向がわかれば、それは姿勢に関する強力な手がかりとなります。この「太陽がどちらにあるか」を電気信号として出力するのが太陽センサ(Sun Sensor)です。

太陽センサを理解すると、以下のような技術領域が開けます。

  • 小型衛星の姿勢制御: CubeSatでは質量・電力の制約から、太陽センサが主要な姿勢センサとなることが多い
  • 太陽電池パドルの追尾制御: 太陽方向を常に追跡し、発電効率を最大化する
  • 安全モードの姿勢復帰: 異常時に最低限の電力を確保するため、太陽方向を基準にした粗姿勢制御を行う
  • 姿勢決定アルゴリズムへの入力: 後述するTRIADやQUESTなどの姿勢決定法において、基準ベクトルの一つとして利用される

本記事の内容

  • 太陽センサの種類と分類(粗/精、アナログ/デジタル)
  • コサインセンサの原理 — 光電流と入射角の物理的関係
  • 2軸太陽センサによる太陽方向ベクトルの幾何学的導出
  • 視野角(FOV)と死角の問題
  • 地球影(エクリプス)の影響と対策
  • Pythonでの太陽方向ベクトル計算と応答特性の可視化

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

太陽センサの役割と分類

宇宙空間における方位基準

地上では、コンパスや重力の方向など、方位の基準は身の回りにあふれています。しかし宇宙空間では磁場は極めて微弱で、重力の方向を直接感じることもできません。そこで衛星は、外部の天体を「灯台」として利用します。太陽はその中でも最も扱いやすい天体です。なぜなら以下の3つの利点があるからです。

  1. 圧倒的な明るさ: 恒星と比較して約 $10^{10}$ 倍以上明るく、簡素な光学系で検出可能
  2. 方向の既知性: 地球暦(エフェメリス)から太陽の慣性座標系における方向ベクトルを高精度に計算可能
  3. 角直径の小ささ: 太陽の見かけの角直径は約 $0.53°$ であり、点光源に近い扱いが可能

太陽センサが出力するのは「衛星の機体座標系から見た太陽の方向ベクトル」 $\bm{s}_b$ です。一方、暦から計算される慣性座標系での太陽方向ベクトルを $\bm{s}_i$ とすると、姿勢回転行列 $\bm{A}$ を用いて次の関係が成り立ちます。

$$ \bm{s}_b = \bm{A} \bm{s}_i $$

この式は「姿勢行列 $\bm{A}$ が既知のベクトル $\bm{s}_i$ を機体座標系のベクトル $\bm{s}_b$ に変換する」ことを表しています。太陽センサは $\bm{s}_b$ を測定するデバイスであり、 $\bm{s}_i$ は暦から既知なので、両者を組み合わせることで姿勢 $\bm{A}$ に関する情報が得られます。ただし、1本のベクトル観測だけでは太陽方向まわりの回転が決まらないため、姿勢を完全に決定するには地磁気センサやスタートラッカーなど別のセンサとの組み合わせが必要になります。

太陽センサの分類体系

太陽センサは「精度」と「出力形式」の2軸で分類されます。

精度による分類:

分類 精度 用途
粗太陽センサ(CSS: Coarse Sun Sensor) $\pm 2°$ 〜 $\pm 10°$ 安全モード、太陽捕捉、初期姿勢獲得
精太陽センサ(FSS: Fine Sun Sensor) $\pm 0.01°$ 〜 $\pm 0.5°$ 高精度姿勢制御、太陽追尾

粗太陽センサは単純な光電素子で構成され、太陽がおおよそどの方向にあるかを知るために使います。「太陽はだいたいあっちの方角」という情報です。安全モードでは衛星が想定外の姿勢になっている可能性があり、まず太陽を見つけて太陽電池パネルを太陽に向けることが最優先です。このとき、粗い精度で十分なので頑丈で広視野の粗太陽センサが活躍します。

一方、精太陽センサは光学スリットやイメージセンサを用いて太陽方向を高精度に測定します。通常の姿勢制御モードで使用され、精度が1桁〜2桁向上しますが、視野角が狭く、コストと質量も増加します。

出力形式による分類:

  • アナログ型: 光電流のアナログ値から入射角を算出。回路が単純で軽量
  • デジタル型: スリットパターンやAPS(Active Pixel Sensor)で入射角をデジタル符号化。高精度で直線性に優れる

多くの衛星では、全天をカバーする6面の粗太陽センサと、高精度が必要な方向に1〜2台の精太陽センサを搭載する構成をとります。

ここまでで太陽センサの全体像を掴みました。次に、最も基本的な動作原理であるコサインセンサの物理について詳しく見ていきます。

コサインセンサの原理

光電流と入射角の関係

太陽センサの動作原理を最も単純な形で理解するには、1枚のフォトダイオード(光電変換素子)を考えるのが最適です。日常的な例で言えば、懐中電灯を壁にまっすぐ当てると明るい丸い光が見えますが、壁を斜めに傾けると同じ光が広い面積に広がって暗くなります。フォトダイオードでも全く同じことが起きます。太陽光がまっすぐ当たれば強い電流が流れ、斜めに当たれば弱くなります。

この関係を定量化しましょう。フォトダイオードの受光面の法線ベクトルを $\hat{\bm{n}}$、太陽方向の単位ベクトルを $\hat{\bm{s}}$ とします。太陽光の入射角 $\theta$ は、法線と太陽方向のなす角です。

$$ \cos\theta = \hat{\bm{n}} \cdot \hat{\bm{s}} $$

受光面に入射する光束(単位面積あたりのエネルギー)は、ランベルトのコサイン則に従います。受光面の面積を $A$、太陽放射照度を $E_{\text{sun}}$ とすると、受光面に入射する光パワーは次のようになります。

$$ P = E_{\text{sun}} \cdot A \cdot \cos\theta $$

フォトダイオードの光電流 $I$ は入射光パワーに比例するので、量子効率と分光感度をまとめた比例定数を $\eta$ とすると次の式が成り立ちます。

$$ \begin{equation} I(\theta) = I_0 \cos\theta \end{equation} $$

ここで $I_0 = \eta E_{\text{sun}} A$ は太陽光が法線方向から入射したとき( $\theta = 0$ )の最大光電流です。この関係から、太陽センサは「コサインセンサ」(cosine sensor)とも呼ばれます。

式(1)の物理的な意味を図で考えると、太陽光線に垂直な断面積が $A\cos\theta$ になることに対応しています。入射角が大きくなるほど有効受光面積が減少し、光電流も減少するのです。

コサインセンサの限界

式(1)は美しい関係ですが、実用上いくつかの課題があります。

1. 角度の一意性がない: $\cos\theta$ だけでは入射角 $\theta$ はわかりますが、太陽がどの方位から来ているかはわかりません。$\theta = 30°$ と言われても、360°のうちどの方向の30°なのか区別できないのです。これは1軸の情報しかないためです。

2. 非線形性: $\cos\theta$ は $\theta = 0$ の付近では変化が緩やかですが、$\theta = 90°$ に近づくと急激に変化します。数学的に言えば次のとおりです。

$$ \frac{dI}{d\theta} = -I_0 \sin\theta $$

$\theta = 0$ 付近では $\sin\theta \approx 0$ なので感度がほぼゼロとなり、微小な角度変化を検出しにくくなります。逆に $\theta = 90°$ 付近では光電流そのものがゼロに近づくため、S/N比が劣化します。

3. 絶対値の校正: $I_0$ は太陽放射照度やフォトダイオードの経年劣化に依存するため、事前校正が必要です。宇宙環境では放射線による素子劣化も問題になります。

これらの限界を克服するために、実際のセンサではフォトダイオードを複数組み合わせて「差動方式」を用いたり、光学スリットを導入して入射方向を直接検出する方式が採用されます。

次に、2つのフォトダイオードを組み合わせた差動検出の原理を見ていきましょう。

差動方式による1軸角度検出

2素子差動検出の原理

コサインセンサの「絶対値の校正が必要」という問題を解消する巧妙な方法があります。それは2つのフォトダイオードを互いに傾けて配置し、両者の出力の比を取る差動検出方式です。

日常のアナロジーで言えば、これは「2つの耳で音源の方向を判定する」ことに似ています。左右の耳に届く音量の差から、音源が左にあるか右にあるかを判断できるのと同じ原理です。

2つのフォトダイオードを法線方向が半角 $\alpha$ だけ対称に傾くように配置したとします。太陽方向からの入射角を基準軸(2素子の対称軸)からの角度 $\phi$ で表すと、各フォトダイオードの光電流は次のようになります。

$$ \begin{align} I_1 &= I_0 \cos(\phi – \alpha) \\ I_2 &= I_0 \cos(\phi + \alpha) \end{align} $$

ここで差と和を計算します。まず三角関数の和差の公式を使います。

$$ \cos(\phi – \alpha) – \cos(\phi + \alpha) = 2\sin\phi\sin\alpha $$

同様に和の公式は次のとおりです。

$$ \cos(\phi – \alpha) + \cos(\phi + \alpha) = 2\cos\phi\cos\alpha $$

したがって、差と和の比を取ると $I_0$ がキャンセルされます。

$$ \begin{equation} \frac{I_1 – I_2}{I_1 + I_2} = \frac{2\sin\phi\sin\alpha}{2\cos\phi\cos\alpha} = \tan\phi \cdot \tan\alpha \end{equation} $$

$\alpha$ はセンサの設計パラメータとして既知なので、この比から入射角 $\phi$ が求まります。

$$ \begin{equation} \phi = \arctan\left(\frac{1}{\tan\alpha} \cdot \frac{I_1 – I_2}{I_1 + I_2}\right) \end{equation} $$

この方式の優れた点は3つあります。

  1. $I_0$ に依存しない: 太陽放射照度の変動や素子劣化の影響を受けにくい
  2. $\phi = 0$ 付近で高い線形性: $\phi$ が小さいとき $\tan\phi \approx \phi$ なので、比は $\phi$ にほぼ比例する
  3. 符号から方向がわかる: $\phi > 0$ なら $I_1 > I_2$、$\phi < 0$ なら $I_1 < I_2$ なので、入射方向の正負を判別できる

線形領域の近似

入射角 $\phi$ が十分小さい($|\phi| \ll 1$)場合、$\tan\phi \approx \phi$(ラジアン)と近似でき、式(2)は次のように単純化されます。

$$ \frac{I_1 – I_2}{I_1 + I_2} \approx \phi \cdot \tan\alpha $$

これは入射角 $\phi$ に対する線形応答であり、信号処理が格段に容易になります。$\alpha$ を大きくすると感度(角度あたりの出力変化)が上がりますが、同時に線形範囲が狭くなるというトレードオフがあります。典型的には $\alpha = 30°$ 〜 $60°$ 程度が採用されます。

しかし、この1軸の差動検出では太陽方向の1成分しか測れません。太陽の方向ベクトルを完全に求めるには、直交する2軸の情報が必要です。次に、2軸太陽センサの幾何学を見ていきましょう。

2軸太陽センサによる太陽方向ベクトルの算出

2軸検出の幾何学

2軸太陽センサは、直交する2つの差動検出軸をもちます。これにより太陽方向の2つの角度成分を同時に測定し、3次元の太陽方向ベクトルを求めることができます。

センサの局所座標系を設定しましょう。センサ面の法線方向(ボアサイト)を $z$ 軸、2つの検出軸をそれぞれ $x$ 軸と $y$ 軸とします。太陽方向の単位ベクトル $\hat{\bm{s}}$ を次のように表します。

$$ \hat{\bm{s}} = \begin{pmatrix} s_x \\ s_y \\ s_z \end{pmatrix} $$

$x$ 軸方向の差動検出系は、太陽光が $xz$ 平面となす角度 $\phi_x$ を測定します。$y$ 軸方向の差動検出系は、太陽光が $yz$ 平面となす角度 $\phi_y$ を測定します。

これらの角度と太陽方向ベクトルの成分の関係は、幾何学的に次のように導出できます。太陽方向ベクトルを $xz$ 平面に投影すると、 $z$ 軸からの角度 $\phi_x$ は次のとおりです。

$$ \tan\phi_x = \frac{s_x}{s_z} $$

同様に $yz$ 平面への投影から次の関係を得ます。

$$ \tan\phi_y = \frac{s_y}{s_z} $$

$\hat{\bm{s}}$ は単位ベクトル($|\hat{\bm{s}}| = 1$)なので、$s_x^2 + s_y^2 + s_z^2 = 1$ の拘束条件があります。$\tan\phi_x$ と $\tan\phi_y$ の2つの測定値からベクトルの3成分を求めましょう。

$t_x = \tan\phi_x$、$t_y = \tan\phi_y$ とおくと、$s_x = s_z t_x$、$s_y = s_z t_y$ です。単位ベクトルの条件に代入します。

$$ s_z^2 t_x^2 + s_z^2 t_y^2 + s_z^2 = 1 $$

$s_z^2$ でくくると次のようになります。

$$ s_z^2 (t_x^2 + t_y^2 + 1) = 1 $$

太陽がセンサの表側($s_z > 0$)にある条件のもとで $s_z$ について解きます。

$$ \begin{equation} s_z = \frac{1}{\sqrt{1 + \tan^2\phi_x + \tan^2\phi_y}} \end{equation} $$

したがって、太陽方向の単位ベクトルは次のように得られます。

$$ \begin{equation} \hat{\bm{s}} = \frac{1}{\sqrt{1 + \tan^2\phi_x + \tan^2\phi_y}} \begin{pmatrix} \tan\phi_x \\ \tan\phi_y \\ 1 \end{pmatrix} \end{equation} $$

この式は、2つの角度測定値 $(\phi_x, \phi_y)$ から太陽方向の3次元単位ベクトルを一意に決定できることを示しています。

精太陽センサ(デジタル型)の原理

精太陽センサでは、光学スリットとラインセンサ(1次元イメージセンサ)を組み合わせて入射角を直接測定する方式が広く用いられます。

原理は単純です。薄いスリット(細い開口部)を通った太陽光が、その奥に置かれたラインセンサ上にスポットを作ります。スリットとラインセンサの間隔を $d$、スポットの中心位置を $p$ とすると、入射角は次のとおりです。

$$ \tan\phi = \frac{p}{d} $$

この方式では、角度がスポットの位置という幾何学的な量に変換されるため、光の強度変動に依存せず、高い精度と線形性が得られます。典型的な精太陽センサの精度は $0.01°$ 〜 $0.1°$ 程度であり、粗太陽センサの100倍以上の精度を実現します。

直交する2組のスリット+ラインセンサを配置すれば、前述の $(\phi_x, \phi_y)$ が直接得られ、式(5)から太陽方向ベクトルを計算できます。

太陽センサの検出原理がわかったところで、次にセンサの視野角(FOV)の制約と、それに伴う「死角」の問題を考えます。

視野角(FOV)と死角の問題

視野角の定義と制約

太陽センサには検出可能な角度範囲、すなわち視野角(FOV: Field of View)が存在します。これはカメラのレンズに画角があるのと同じです。FOVの外に太陽がある場合、センサは有効な出力を返すことができません。

粗太陽センサの場合、コサインセンサの原理上、入射角が $90°$ に近づくと光電流がゼロになるため、実用的な FOV は半角 $\pm 60°$ 〜 $\pm 80°$ 程度です。つまり1つのセンサで全天の約 $1/3$ 程度しかカバーできません。

精太陽センサの場合、スリット方式の幾何学的制約から FOV はさらに狭く、半角 $\pm 16°$ 〜 $\pm 64°$ 程度が一般的です。

全天カバーの設計

衛星がどのような姿勢をとっても太陽を検出できるようにするには、複数のセンサを衛星の異なる面に配置して全天をカバーする必要があります。

最も典型的な配置は、衛星の6面体(直方体)の各面に1つずつ粗太陽センサを搭載する方式です。各センサの法線は $\pm x$、$\pm y$、$\pm z$ の6方向を向き、半角 $60°$ 以上の FOV があれば全天をカバーできます。

$N$ 個のセンサで全天をカバーするために必要な最小 FOV の半角 $\theta_{\min}$ は、センサの配置幾何に依存します。立方体の6面配置の場合は次のとおりです。

$$ \theta_{\min} = \arctan(\sqrt{2}) \approx 54.7° $$

これは立方体の面の中心と辺の中点を結ぶ角度に対応します。この値は、隣接する2つのセンサの FOV が辺の方向で接する条件から求まります。

複数センサの統合

複数のセンサが太陽を同時に検出する場合、どのセンサの値を使うかの選択が必要です。一般的な方法は次の2つです。

方法1: 最大出力センサの選択 — 光電流(またはそれに相当する信号振幅)が最大のセンサを選ぶ。最大出力は太陽が最も正面に近いセンサを意味するため、精度が最も高い測定値を利用できます。

方法2: 重み付き平均 — 各センサの出力を信号強度に応じた重みで統合する。FOV の重なり領域で平滑化が可能ですが、計算コストが増加します。

衛星設計では、安全モード用の6面粗太陽センサに加えて、通常運用モード用の精太陽センサを太陽パネル面の方向に配置するのが典型的な構成です。

ここまでで太陽センサのハードウェア的な特性を理解しました。次に、もう一つの重要な制約 — 衛星が地球の影に入る「エクリプス」の問題を考えます。

地球影(エクリプス)の影響

エクリプスの幾何

地球の周りを周回する衛星は、軌道の一部で地球の影に入ります。この期間をエクリプス(eclipse)と呼びます。エクリプス中は太陽光が遮断されるため、太陽センサは機能しません。

エクリプスの幾何学を考えましょう。地球を半径 $R_E \approx 6371$ km の球とし、太陽光を平行光線と近似します(太陽は衛星から約1億5千万km離れているため、十分良い近似です)。

衛星の軌道高度を $h$、軌道半径を $r = R_E + h$ とすると、円軌道の場合のエクリプス最大半角 $\beta_{\max}$ は次のように求まります。

地球の影の円錐を考えると、衛星軌道面上でエクリプスが生じる半角は次のとおりです。

$$ \begin{equation} \beta_{\max} = \arcsin\left(\frac{R_E}{r}\right) \end{equation} $$

したがって、エクリプスの最大持続時間は軌道周期 $T$ に対して次のようになります。

$$ \Delta t_{\text{eclipse}} = \frac{2\beta_{\max}}{2\pi} T = \frac{\beta_{\max}}{\pi} T $$

低軌道(LEO)の例として、軌道高度 $h = 500$ km の場合を計算しましょう。

$r = 6371 + 500 = 6871$ km なので、

$$ \beta_{\max} = \arcsin\left(\frac{6371}{6871}\right) \approx \arcsin(0.9273) \approx 68.1° $$

軌道周期は $T \approx 94.5$ 分(ケプラーの第3法則から)なので、

$$ \Delta t_{\text{eclipse}} \approx \frac{68.1°}{180°} \times 94.5 \approx 35.8 \text{ 分} $$

つまり、低軌道衛星は1周回の約 $38\%$ の時間でエクリプスに入り、太陽センサが使えなくなります。これは無視できない大きな割合です。

エクリプス中の姿勢維持

エクリプス中に太陽センサが使えない問題に対する対策は、主に以下の3つです。

1. ジャイロスコープによる伝搬: エクリプス突入前の最後の姿勢推定値を初期値として、ジャイロ(角速度センサ)の出力を積分して姿勢を伝搬します。ただし、ジャイロにはドリフト(バイアス)があるため、時間とともに推定誤差が蓄積します。30分のエクリプスで角度にして数度の誤差が生じることもあります。

2. 地磁気センサの併用: 地球の磁場方向は地磁気モデル(IGRFなど)から予測可能であり、地磁気センサで測定した磁場方向と比較することで姿勢情報を補完できます。太陽センサほど精度は高くありませんが、エクリプス中の粗姿勢維持には有効です。

3. スタートラッカーの使用: 恒星は太陽の有無に関係なく(むしろ太陽光が遮られた方が)観測できるため、スタートラッカーはエクリプス中でも動作可能です。高精度な姿勢決定が可能ですが、コストと質量が増加します。

実際の衛星では、これらの手法を組み合わせた統合的な姿勢推定が行われます。これは後の記事で扱うカルマンフィルタの重要な応用の一つです。

ここまでで太陽センサの理論的な基盤をすべてカバーしました。次に、これらの知識をPythonで実装し、実際に太陽方向ベクトルの計算やセンサ応答の特性を確認していきましょう。

Pythonでの実装: 太陽方向ベクトルの計算

コサインセンサの応答特性

まず、単一のコサインセンサの応答特性と差動検出方式の応答を比較するコードを実装します。入射角に対する各方式の出力を可視化し、差動方式の線形性の優位性を確認しましょう。

import numpy as np
import matplotlib.pyplot as plt

# --- パラメータ設定 ---
theta = np.linspace(-90, 90, 500)  # 入射角 [deg]
theta_rad = np.radians(theta)

# コサインセンサの応答
I0 = 1.0  # 最大光電流(正規化)
I_cosine = I0 * np.cos(theta_rad)
I_cosine = np.maximum(I_cosine, 0)  # 負の値はクリップ

# 差動検出の応答(半角 alpha = 30°, 45°, 60°)
alphas = [30, 45, 60]
colors = ['#00bcd4', '#ff9800', '#e91e63']

fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# 左: コサインセンサの応答
axes[0].plot(theta, I_cosine, 'w', linewidth=2, label='$I(\\theta) = I_0 \\cos\\theta$')
axes[0].set_xlabel('Incident angle $\\theta$ [deg]', fontsize=12)
axes[0].set_ylabel('Photocurrent $I / I_0$', fontsize=12)
axes[0].set_title('Cosine Sensor Response', fontsize=13)
axes[0].legend(fontsize=11)
axes[0].grid(True, alpha=0.3)
axes[0].set_xlim(-90, 90)
axes[0].set_ylim(-0.1, 1.1)
axes[0].set_facecolor('#1a1a2e')

# 右: 差動検出の応答
for alpha_deg, color in zip(alphas, colors):
    alpha_rad = np.radians(alpha_deg)
    I1 = I0 * np.cos(theta_rad - alpha_rad)
    I2 = I0 * np.cos(theta_rad + alpha_rad)
    I1 = np.maximum(I1, 0)
    I2 = np.maximum(I2, 0)

    # 差動比
    diff_ratio = np.where(
        (I1 + I2) > 0.01,
        (I1 - I2) / (I1 + I2),
        np.nan
    )
    axes[1].plot(theta, diff_ratio, color=color, linewidth=2,
                 label=f'$\\alpha = {alpha_deg}°$')

# 理想的な線形応答(参考線)
axes[1].plot(theta, np.radians(theta) * np.tan(np.radians(45)),
             '--', color='gray', alpha=0.5, linewidth=1,
             label='Ideal linear ($\\alpha=45°$)')

axes[1].set_xlabel('Incident angle $\\phi$ [deg]', fontsize=12)
axes[1].set_ylabel('Differential ratio $(I_1 - I_2)/(I_1 + I_2)$', fontsize=12)
axes[1].set_title('Differential Detection Response', fontsize=13)
axes[1].legend(fontsize=10)
axes[1].grid(True, alpha=0.3)
axes[1].set_xlim(-90, 90)
axes[1].set_ylim(-2.5, 2.5)
axes[1].set_facecolor('#1a1a2e')

fig.patch.set_facecolor('#0f0f23')
for ax in axes:
    ax.tick_params(colors='white')
    ax.xaxis.label.set_color('white')
    ax.yaxis.label.set_color('white')
    ax.title.set_color('white')
    for spine in ax.spines.values():
        spine.set_color('white')

plt.tight_layout()
plt.savefig('sun_sensor_response.png', dpi=150, bbox_inches='tight',
            facecolor='#0f0f23')
plt.show()

左のグラフから、コサインセンサの応答は $\theta = 0°$ 付近で変化が緩やか(感度が低い)であることが確認できます。これは $d(\cos\theta)/d\theta = -\sin\theta$ が $\theta = 0$ でゼロになるためです。一方、右のグラフでは差動方式が $\phi = 0$ 付近で良好な線形性を示していることがわかります。半角 $\alpha$ が大きいほど傾き(感度)が急峻になりますが、線形領域が狭くなるトレードオフが明確に見て取れます。$\alpha = 45°$ の場合、$\pm 30°$ 程度の範囲で差動比はほぼ線形に変化しており、信号処理が容易であることがわかります。

2軸太陽センサによる太陽方向ベクトルの計算

次に、2軸太陽センサの出力から太陽方向ベクトルを計算する実装を行います。式(5)に基づく計算を行い、さまざまな入射方向に対するベクトルの復元精度を確認しましょう。

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

def sun_vector_from_angles(phi_x_deg, phi_y_deg):
    """2軸太陽センサの角度測定値から太陽方向単位ベクトルを計算"""
    tx = np.tan(np.radians(phi_x_deg))
    ty = np.tan(np.radians(phi_y_deg))
    sz = 1.0 / np.sqrt(1.0 + tx**2 + ty**2)
    sx = sz * tx
    sy = sz * ty
    return np.array([sx, sy, sz])

def true_sun_angles(s_vec):
    """太陽方向ベクトルから2軸角度を逆算(検証用)"""
    phi_x = np.degrees(np.arctan2(s_vec[0], s_vec[2]))
    phi_y = np.degrees(np.arctan2(s_vec[1], s_vec[2]))
    return phi_x, phi_y

# --- テスト: 様々な太陽方向に対するベクトル復元 ---
np.random.seed(42)
n_samples = 200

# ランダムな太陽方向を生成(FOV内: 半角60°以内)
phi_x_true = np.random.uniform(-55, 55, n_samples)
phi_y_true = np.random.uniform(-55, 55, n_samples)

# ノイズを加える(精太陽センサ: σ = 0.05°)
noise_sigma = 0.05  # deg
phi_x_meas = phi_x_true + np.random.normal(0, noise_sigma, n_samples)
phi_y_meas = phi_y_true + np.random.normal(0, noise_sigma, n_samples)

# 真のベクトルと測定ベクトルを計算
errors = []
for i in range(n_samples):
    s_true = sun_vector_from_angles(phi_x_true[i], phi_y_true[i])
    s_meas = sun_vector_from_angles(phi_x_meas[i], phi_y_meas[i])

    # 角度誤差 [deg]
    cos_err = np.clip(np.dot(s_true, s_meas), -1, 1)
    err_deg = np.degrees(np.arccos(cos_err))
    errors.append(err_deg)

errors = np.array(errors)

# --- 可視化 ---
fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# 左: 入射角に対する復元誤差
incidence = np.sqrt(phi_x_true**2 + phi_y_true**2)
scatter = axes[0].scatter(incidence, errors * 3600,  # arcsecに変換
                          c=incidence, cmap='cool', s=15, alpha=0.7)
axes[0].set_xlabel('Total incidence angle [deg]', fontsize=12)
axes[0].set_ylabel('Sun vector error [arcsec]', fontsize=12)
axes[0].set_title('Sun Vector Recovery Error vs Incidence Angle', fontsize=13)
axes[0].grid(True, alpha=0.3)
axes[0].set_facecolor('#1a1a2e')
plt.colorbar(scatter, ax=axes[0], label='Incidence [deg]')

# 右: 誤差のヒストグラム
axes[1].hist(errors * 3600, bins=30, color='#00bcd4', edgecolor='white',
             alpha=0.8)
axes[1].axvline(np.mean(errors) * 3600, color='#ff9800', linewidth=2,
                linestyle='--', label=f'Mean: {np.mean(errors)*3600:.1f} arcsec')
axes[1].axvline(np.std(errors) * 3600 * 3, color='#e91e63', linewidth=2,
                linestyle='--', label=f'3$\\sigma$: {np.std(errors)*3600*3:.1f} arcsec')
axes[1].set_xlabel('Sun vector error [arcsec]', fontsize=12)
axes[1].set_ylabel('Count', fontsize=12)
axes[1].set_title('Error Distribution ($\\sigma_{noise}$ = 0.05°)', fontsize=13)
axes[1].legend(fontsize=10)
axes[1].grid(True, alpha=0.3)
axes[1].set_facecolor('#1a1a2e')

fig.patch.set_facecolor('#0f0f23')
for ax in axes:
    ax.tick_params(colors='white')
    ax.xaxis.label.set_color('white')
    ax.yaxis.label.set_color('white')
    ax.title.set_color('white')
    for spine in ax.spines.values():
        spine.set_color('white')

plt.tight_layout()
plt.savefig('sun_vector_error.png', dpi=150, bbox_inches='tight',
            facecolor='#0f0f23')
plt.show()

print(f"角度測定ノイズ: σ = {noise_sigma}°")
print(f"太陽方向ベクトル誤差(平均): {np.mean(errors)*3600:.1f} arcsec")
print(f"太陽方向ベクトル誤差(3σ): {np.std(errors)*3600*3:.1f} arcsec")

左のグラフから、太陽方向ベクトルの復元誤差は入射角が大きくなるにつれて増大する傾向が見て取れます。これは $\tan\phi$ の非線形性が大きい入射角で顕著になるためです。特に入射角が $50°$ を超えると誤差の増大が目立ちます。右のヒストグラムでは、$\sigma = 0.05°$ の角度測定ノイズに対して、太陽方向ベクトルの誤差が平均で約 $180$ arcsec(約 $0.05°$)となっていることが確認できます。これは精太陽センサの典型的な性能仕様と整合しています。

エクリプス判定と太陽センサの可用性

最後に、低軌道衛星の1周回にわたるエクリプス判定と太陽センサの可用性をシミュレーションします。

import numpy as np
import matplotlib.pyplot as plt

# --- 軌道パラメータ ---
R_E = 6371.0       # 地球半径 [km]
h = 500.0           # 軌道高度 [km]
r = R_E + h         # 軌道半径 [km]
mu = 398600.4418    # 地球の重力定数 [km^3/s^2]

# 軌道周期
T = 2 * np.pi * np.sqrt(r**3 / mu)  # [s]
T_min = T / 60.0
print(f"軌道高度: {h} km")
print(f"軌道半径: {r} km")
print(f"軌道周期: {T_min:.1f} 分")

# エクリプス半角
beta_max_rad = np.arcsin(R_E / r)
beta_max_deg = np.degrees(beta_max_rad)
print(f"エクリプス最大半角: {beta_max_deg:.1f}°")

# エクリプス時間
t_eclipse = (beta_max_rad / np.pi) * T
t_eclipse_min = t_eclipse / 60.0
eclipse_fraction = t_eclipse / T * 100
print(f"エクリプス時間: {t_eclipse_min:.1f} 分 ({eclipse_fraction:.1f}%)")

# --- 1周回のシミュレーション ---
# 簡単化: 衛星の軌道角度 theta を 0 → 360° で回す
# 太陽方向は固定(beta_sun = 0°: 太陽は軌道面内)
t = np.linspace(0, T, 1000)
orbit_angle = 2 * np.pi * t / T  # 軌道角 [rad]

# 太陽方向角(軌道面内で太陽が 0° 方向にいると仮定)
sun_angle = 0.0  # rad

# 衛星-太陽の角度
sat_sun_angle = orbit_angle - sun_angle  # 衛星軌道角から太陽方向への角度

# エクリプス判定: 太陽の反対側(π ± beta_max)にいるとき
# 衛星が地球の影に入る条件
in_eclipse = np.abs(np.pi - np.mod(sat_sun_angle, 2*np.pi)) < beta_max_rad

# 太陽センサの出力(コサインセンサ、ボアサイトが太陽パネル法線方向)
# 太陽入射角はセンサ面法線と太陽方向の角度
# 簡単のため、センサは +z 方向を向いていると仮定
cos_incidence = np.cos(sat_sun_angle)
cos_incidence[cos_incidence < 0] = 0  # 裏側は検出不可
cos_incidence[in_eclipse] = 0  # エクリプス中は出力ゼロ

# --- 可視化 ---
fig, axes = plt.subplots(3, 1, figsize=(12, 10), sharex=True)

# 1: 軌道角度
axes[0].plot(t / 60, np.degrees(orbit_angle), color='#00bcd4', linewidth=2)
axes[0].fill_between(t / 60, 0, 360, where=in_eclipse,
                     color='#444444', alpha=0.4, label='Eclipse')
axes[0].set_ylabel('Orbit angle [deg]', fontsize=12)
axes[0].set_title('Low Earth Orbit: Eclipse and Sun Sensor Availability', fontsize=13)
axes[0].legend(fontsize=10)
axes[0].grid(True, alpha=0.3)
axes[0].set_facecolor('#1a1a2e')

# 2: エクリプスフラグ
axes[1].fill_between(t / 60, 0, 1, where=~in_eclipse,
                     color='#ff9800', alpha=0.6, label='Sunlit')
axes[1].fill_between(t / 60, 0, 1, where=in_eclipse,
                     color='#333333', alpha=0.6, label='Eclipse')
axes[1].set_ylabel('Illumination', fontsize=12)
axes[1].set_yticks([0, 1])
axes[1].set_yticklabels(['Shadow', 'Sunlit'])
axes[1].legend(fontsize=10, loc='upper right')
axes[1].grid(True, alpha=0.3)
axes[1].set_facecolor('#1a1a2e')

# 3: 太陽センサ出力
axes[2].plot(t / 60, cos_incidence, color='#ff9800', linewidth=2,
             label='Sun sensor output')
axes[2].fill_between(t / 60, 0, cos_incidence, color='#ff9800', alpha=0.2)
axes[2].fill_between(t / 60, 0, 1, where=in_eclipse,
                     color='#333333', alpha=0.4)
axes[2].set_xlabel('Time [min]', fontsize=12)
axes[2].set_ylabel('Sensor output (normalized)', fontsize=12)
axes[2].legend(fontsize=10)
axes[2].grid(True, alpha=0.3)
axes[2].set_facecolor('#1a1a2e')

fig.patch.set_facecolor('#0f0f23')
for ax in axes:
    ax.tick_params(colors='white')
    ax.xaxis.label.set_color('white')
    ax.yaxis.label.set_color('white')
    ax.title.set_color('white')
    for spine in ax.spines.values():
        spine.set_color('white')

plt.tight_layout()
plt.savefig('eclipse_simulation.png', dpi=150, bbox_inches='tight',
            facecolor='#0f0f23')
plt.show()

3つのグラフから、低軌道衛星の1周回における太陽センサの可用性が視覚的に理解できます。第1グラフでは軌道角度の進行とエクリプス領域(灰色帯)の対応が示されています。第2グラフでは、約35分間のエクリプス(軌道周期の約38%)が存在し、その間は完全に太陽センサが利用不可であることが明確です。第3グラフでは、太陽センサの出力がコサイン特性に従い、太陽方向に正対するとき最大になることが確認できます。さらに、エクリプス直後の復帰時には太陽を再捕捉する必要があり、このタイミングで粗太陽センサの広いFOVが重要になることが理解できます。

6面粗太陽センサによる全天カバーの可視化

最後に、衛星の6面に配置した粗太陽センサがどのように全天をカバーするかを可視化します。

import numpy as np
import matplotlib.pyplot as plt

def css_output(sun_vec, normal_vec):
    """粗太陽センサの出力(コサインセンサ)"""
    cos_angle = np.dot(sun_vec, normal_vec)
    return max(0.0, cos_angle)

# 6面のセンサ法線ベクトル
normals = {
    '+X': np.array([1, 0, 0]),
    '-X': np.array([-1, 0, 0]),
    '+Y': np.array([0, 1, 0]),
    '-Y': np.array([0, -1, 0]),
    '+Z': np.array([0, 0, 1]),
    '-Z': np.array([0, 0, -1]),
}

# 全天の太陽方向を球面座標で走査
n_az = 360
n_el = 180
azimuth = np.linspace(0, 2*np.pi, n_az)
elevation = np.linspace(-np.pi/2, np.pi/2, n_el)
AZ, EL = np.meshgrid(azimuth, elevation)

# 各方向での最大センサ出力と最適センサ番号
max_output = np.zeros_like(AZ)
best_sensor = np.zeros_like(AZ, dtype=int)

sensor_names = list(normals.keys())
sensor_normals = list(normals.values())

for i in range(n_el):
    for j in range(n_az):
        # 太陽方向の単位ベクトル
        sun = np.array([
            np.cos(EL[i, j]) * np.cos(AZ[i, j]),
            np.cos(EL[i, j]) * np.sin(AZ[i, j]),
            np.sin(EL[i, j])
        ])

        outputs = [css_output(sun, n) for n in sensor_normals]
        max_output[i, j] = max(outputs)
        best_sensor[i, j] = np.argmax(outputs)

# --- 可視化 ---
fig, axes = plt.subplots(1, 2, figsize=(16, 5),
                         subplot_kw={'projection': 'mollweide'})

# 左: 最大センサ出力
im1 = axes[0].pcolormesh(AZ - np.pi, EL, max_output,
                         cmap='inferno', shading='auto')
axes[0].set_title('Maximum CSS Output Over Full Sky', fontsize=13,
                  color='white', pad=15)
axes[0].grid(True, alpha=0.3)
plt.colorbar(im1, ax=axes[0], label='Normalized output',
             orientation='horizontal', pad=0.08, shrink=0.8)

# 右: 最適センサ番号
colors_map = plt.cm.Set2(np.linspace(0, 1, 6))
from matplotlib.colors import ListedColormap, BoundaryNorm
cmap_discrete = ListedColormap(colors_map)
bounds = np.arange(-0.5, 6.5, 1)
norm = BoundaryNorm(bounds, cmap_discrete.N)

im2 = axes[1].pcolormesh(AZ - np.pi, EL, best_sensor,
                         cmap=cmap_discrete, norm=norm, shading='auto')
axes[1].set_title('Best Sensor Selection (6-face CSS)', fontsize=13,
                  color='white', pad=15)
axes[1].grid(True, alpha=0.3)
cbar = plt.colorbar(im2, ax=axes[1], orientation='horizontal',
                    pad=0.08, shrink=0.8, ticks=range(6))
cbar.ax.set_xticklabels(sensor_names, fontsize=9)
cbar.set_label('Active sensor')

fig.patch.set_facecolor('#0f0f23')
for ax in axes:
    ax.set_facecolor('#1a1a2e')
    ax.tick_params(colors='white')
    ax.title.set_color('white')

plt.tight_layout()
plt.savefig('css_coverage.png', dpi=150, bbox_inches='tight',
            facecolor='#0f0f23')
plt.show()

# 出力がゼロの方向があるか確認
zero_fraction = np.sum(max_output < 0.01) / max_output.size * 100
print(f"出力がほぼゼロの方向の割合: {zero_fraction:.1f}%")
print("→ 6面CSSで全天カバーが達成されています" if zero_fraction < 0.1
      else f"→ {zero_fraction:.1f}%の死角が存在します")

左のモルワイデ投影図から、6面の粗太陽センサにより全天がカバーされていることが確認できます。ただし、センサの切り替わり境界(立方体の辺や頂点に対応する方向)では出力が低下し、最小値は $\cos(54.7°) \approx 0.577$ 程度です。右の図では、太陽の入射方向に応じてどのセンサが最大出力を返すかが色分けされており、6つの領域がほぼ等しい立体角を持つことがわかります。実際の衛星では、これらの境界付近での精度劣化を考慮して、FOVのオーバーラップを持たせた設計が行われます。

まとめ

本記事では、太陽センサの動作原理について、基礎物理から2軸ベクトル算出、運用上の制約まで体系的に解説しました。

  • コサインセンサの原理: 光電流が入射角のコサインに比例する物理法則(ランベルトのコサイン則)に基づいている
  • 差動検出方式: 2素子の出力比を取ることで、絶対校正に依存しない頑健な角度検出が可能になる
  • 2軸太陽センサ: 直交する2軸の角度測定 $(\phi_x, \phi_y)$ から太陽方向の3次元単位ベクトルを一意に復元できる
  • 視野角とカバレッジ: 6面に粗太陽センサを配置することで全天をカバーし、安全モードでの太陽捕捉を実現する
  • エクリプスの影響: 低軌道では軌道周期の約38%でエクリプスに入り、太陽センサが利用不可になる。他のセンサとの併用が必須

太陽センサは太陽方向の1本のベクトルしか提供しないため、これだけでは3軸姿勢を完全には決定できません。完全な姿勢決定には、もう1本の独立なベクトル観測が必要です。次の記事では、恒星のパターンマッチングにより高精度な姿勢決定を実現するスタートラッカーを解説します。スタートラッカーは太陽センサとは異なるアプローチで姿勢情報を提供し、両者を組み合わせることで強力な姿勢決定が可能になります。

次のステップとして、以下の記事も参考にしてください。