天体写真のキャリブレーションとノイズの伝搬(ver.0.90)

著者:だいこもん

概要:天体写真でのキャリブレーションでは、異なる画像の加減算や除算を行う。それによって画像のノイズがどのように変化するかを計算する方法を考える。公式を頭ごなしに与えるのではなく、すべて基本原理から導出することを心がけた。読者として高校レベルの微積分を理解していることを前提としている。

準備

画像の輝度値は,各フレームごとに揺らぐ確率変数である。ここでは,確率変数の平均,
揺らぎから決まるノイズの大きさ,SN比,それぞれの定義を確認する。

輝度の平均

以下では特に断らない限り,「画像」といったら,その画像とある特定のピクセルのことを言っていることにする。
ある画像の輝度を \( x \) とする。いま \( n \) 枚の画像を撮影したとして,それぞれの画像に番号 \(k=1,2,3\cdots n\) を振っておく。それぞれの画像の輝度の値を
$$
x=x_1, x_2, x_3,\cdots,x_n
$$
と書くと,これらの値は同じではなく,ある平均値の周りで揺らいでいる。このような \(x\) は確率変数と呼ばれる。\(n\) 枚の画像についての輝度値の平均は
$$
\frac{ x_1+x_2+x_3+\cdots+x_n }{n}
$$
である。平均は、輝度の「確率分布」を使って次のように書き換えることもできる。つまり,ピクセルの輝度が \(x\) という値を持つ確率を \(P(x)\) とすれば,平均値は
$$
\int_{-\infty}^{\infty} xP(x)dx~~~~~~~~~~~~(1)
$$
とも書ける。これは期待値と呼ばれるが,期待値と平均値は同じものと思って良いのか場合によっ
ては違うのか,著者はよく理解できていない。ここでは同じものとする。以下では,確率変数 \(x\) に対してその平均を記号 \(\langle x\rangle\) でらわす。確率変数を \(“\langle”\) と \(“\rangle”\) で挟んだら,それは \(x\) に \(P(x)\) を掛けて積分せよ,という意味で使う。

ノイズの大きさ

画像の輝度の値はその平均値周りで揺らいでいるから,そのノイズの大きさを議論する場合は
画像の輝度の平均値からのズレ,つまり
$$
x-\langle x\rangle
$$
が問題になる。しかしこの値は正負両方の値を持つので,扱いづらい。そこで同じ量を
2乗して平均をとった量,つまり
$$
\langle \left(x-\langle x\rangle \right)^2\rangle
$$
を考える。これは分散と呼ばれ、今後は記号\(\sigma^2\) で書き表す。平均値の定義(1)式をつかって分散を計算すると
$$
\sigma^2=\langle x^2\rangle -\langle x\rangle ^2~~~~~~~~~(2)
$$
と変形でき,こちらの表現の方が計算する時に便利。分散の平方根
$$
\sigma=\sqrt{\langle x^2\rangle -\langle x\rangle ^2}
$$
は標準偏差と呼ばれる。この標準偏差をノイズの大きさとして定義する。

SN比

一枚の画像の輝度を\(x\),その輝度のノイズの大きさを \(\sigma\) として,比
$$
SNR\equiv\frac{\langle x\rangle}{\sigma}
$$
はSN比と呼ばれる。このSN比に対して,常用対数をとって10倍した値
$$
10\log_{10}\frac{\langle x\rangle}{\sigma}
$$
がよく聞く「デシベル」の定義。たとえば \(\frac{x}{\sigma}=1\) なら0デシベル,ノイズが小さくて \(\frac{x}{\sigma}=10\) なら10デシベルと言った具合。

誤差の伝搬

加算した画像のノイズの大きさ

複数の画像を加算した結果のノイズの大きさと,もとの画像のノイズの大きさの関係を考えたい。まず簡単のために2枚の画像 \(x_a,~x_b\) を考える。加算した画像のノイズの大きさは, \(x_a+x_b\) の標準偏差で与えられるから,(2)式を用いて、
\begin{align}
(\text{ノイズの大きさ})^2=&\langle(x_a+x_b)^2\rangle -\langle x_a+x_b\rangle ^2\\ =&\langle{x_a}^2+2x_ax_b+{x_b}^2\rangle -\langle x_a+x_b\rangle ^2
\end{align}
ここで「和をとってから平均すること」と「平均してから和を取ること」は等しいので,うえの2行目の式は
$$
\langle{x_a}^2\rangle +2\langle x_ax_b\rangle +\langle {x_b}^2\rangle -\left({\langle x_a\rangle}^2 + 2\langle x_a\rangle\langle x_b\rangle + {\langle x_b\rangle}^2\right)
$$
と変形できる。さらに,2枚の画像に含まれるノイズがお互いに独立(*1)であるなら,
$$
\langle x_ax_b\rangle=\langle x_a\rangle\langle x_b\rangle
$$
であることが(1)式から証明できるので,結局
\begin{align}
\text{(ノイズの大きさ)}^2 &=\left(\langle{x_a}^2\rangle-\langle{x_a}\rangle^2\right)+\left(\langle{x_b}^2\rangle-\langle{x_b}\rangle^2\right)\\
&\equiv{\sigma_a}^2+{\sigma_b}^2
\end{align}
を得る。

*1 ノイズが独立とは?
ある画像のノイズ値と別の画像のノイズ値が全く関係なく決まる状態のことを呼ぶ。反対に、それぞれの画像のノイズ値がお互いに影響を及ぼし合うような場合、ノイズが相関をもつという。もしx,yの二つの画像のノイズに相関があるなら、(1)式に表れた画像がxという輝度を持つ確率は P(x,y)の形でxとyの両方に依存する。ノイズが独立であれば、一方の画像が輝度xをもち、もう一方の画像がyを持つ確率はP(x)P(y)となる。

以上から,2枚の画像を加算したとき,そのノイズの大きさは
$$
\text{(ノイズの大きさ)}=\sqrt{ {\sigma_a}^2+{\sigma_b}^2 }
$$
と,それぞれの画像の分散の和の平方根で与えられることがわかった。
一般的に \(n\) 枚の画像があり,それぞれの分散が \({\sigma_1}^2, {\sigma_2}^2,\cdots{\sigma_n}^2\) である場合に,
これらを全て加算して得られる画像のノイズの大きさは
$$
\sqrt{{\sigma_1}^2+{\sigma_2}^2+\cdots+{\sigma_n}^2}
$$
となる。特に単純な場合として,全ての画像の分散が等しく、値 \({\sigma}^2\) を持つ場合は,加算した画像のノイズの大きさは \(\sqrt{n}\sigma\) である。つまり加算によってノイズは増える

練習問題:減算した画像のノイズの大きさを,\(x_a-x_b\) の標準偏差を計算することで求め,結果が加算と同じであることを確かめよ。

誤差の伝搬法則

上では加算における誤差を考えたが、もっと一般的な演算によって得られたデータの誤差が、元のデータでどのようにあらわされるかを考える。つまり、二枚の画像 \(x,~y\) があって、その平均輝度 \(\langle x \rangle,\langle y\rangle\)と、分散 \(\sigma_x,\sigma_y\) がわかっているとき
$$
z=f(x,y)
$$
と表される \(z\) の分散がどうなるかを求めたい。

\(x, y\) がもつノイズを \(\Delta x, \Delta y\) とするとこれらの間には
$$
x=\langle x \rangle + \Delta x\\
y=\langle y \rangle + \Delta y
$$
という関係がある。これを \(z=f(x,y)\) に代入し,ノイズが小さいとしてテイラー展開の1次までを実行すれば
\begin{align}
z&=f\big(\langle x \rangle, \langle y \rangle\big)+\frac{\partial f\big(\langle x \rangle, \langle y \rangle\big)}{\partial \langle x \rangle}\Delta x+\frac{\partial f\big(\langle x \rangle,\langle y \rangle\big)}{\partial \langle y \rangle}\Delta y
\end{align}
ここで \(f\big(\langle x \rangle, \langle y \rangle\big)=\langle z\rangle\)であることに留意すれば,上の式は次のように整理できる。
$$
z-\langle z\rangle=\frac{\partial f\big(\langle x \rangle, \langle y \rangle\big)}{\partial \langle x \rangle}\big(x-\langle x\rangle\big)+\frac{\partial f\big(\langle x \rangle,\langle y \rangle\big)}{\partial \langle y \rangle}\big(y-\langle y\rangle\big)
$$
最後に,この両辺を2乗して全体の平均を取る(\(\langle\)と\(\rangle\)で挟む)。クロスタームの \( \langle(x-\langle x \rangle)(y-\langle y \rangle)\rangle=0\) になるので,
$$
{\sigma_z}^2 = \left(\frac{\partial f\big(\langle x \rangle, \langle y \rangle\big)}{\partial \langle x \rangle}\right)^2{\sigma_x}^2+\left(\frac{\partial f\big(\langle x \rangle,\langle y \rangle\big)}{\partial \langle y \rangle}\right)^2{\sigma_y}^2
$$
これは誤差の伝搬法則と呼ばれる。誤差が平均輝度に対して小さい場合のみ使える式であることに注意。たとえばダーク画像に対してこの式を使うのはキケンである。

この関係は,あとでフラット除算のノイズを評価する時に使う。

キャリブレーション結果のSN比

加算した画像のSN比

\(n\) 枚の画像を加算した時,その輝度は
$$
x_1+x_2+\cdots+x_n
$$
であり,ノイズの大きさは
$$
\sqrt{{\sigma_1}^2+{\sigma_2}^2+\cdots+{\sigma_n}^2}
$$
であるから,加算した画像のSN比は
$$
(SNR)_n\equiv\frac{x_1+x_2+\cdots+x_n}{\sqrt{{\sigma_1}^2+{\sigma_2}^2+\cdots+{\sigma_n}^2}}
$$
ここで単純な場合として,それぞれの画像の撮影条件が同じで分散は全て同じ値 \(\sigma^2\) を持つとすれば,
$$
(SNR)_n=\sqrt{n}\frac{\langle x\rangle}{\sigma}
$$
つまり加算して得られた画像のSN比は,1枚の画像の \(\sqrt{n}\) 倍になる。加算ではなくて,加算平均を考える場合は,分母と分子が両方とも \(n\) で割られることになるので,結果は同じ。

ライト画像のノイズは輝度のショットノイズが支配的
ライト画像には由来の異なる複数のノイズが含まれていて、ダーク減算で除去できる固定ノイズを除けば、それらは
・輝度信号のショットノイズ
・電子の熱揺らぎに起因するランダムノイズ
・センサーの読み出しエラーに起因するランダムノイズ
がある。このなかで特に、輝度のショットノイズは飛来する光子の頻度と関係し、これはポアソン分布に従う。ポアソン分布は、分散が平均輝度に等しいという単純な性質を持つので、輝度が大きいほどノイズも大きくなり、極端な短時間露光をおこなっていない限りは、ライト画像のなかで上にあげた残りの二つのノイズに比べて大きいノイズになる。加算処理の最も大きな目的は、このショットノイズを小さくすることにある。

ダーク減算した画像のSN比

ライト画像を \( x \),ダーク画像を \(y\) として画像 \(x-y\) のノイズの大きさを考える。減算後の画像の輝度は\( \langle x-y\rangle=\langle x\rangle-\langle y\rangle\) で、ノイズの大きさは
$$
\sqrt{{\sigma_x}^2+{\sigma_y}^2}
$$
だから、減算後のSN比は
$$
(SNR)_{dark}=\frac{\langle x\rangle-\langle y\rangle}{ \sqrt{{\sigma_x}^2+{\sigma_y}^2} }
$$
ダーク画像はライト画像よりも十分に暗いので\(\langle y\rangle\)は無視して良い。よってダーク減算後のSN比は次のようになる
$$
(SNR)_{dark}=\frac{\langle x\rangle}{ \sqrt{{\sigma_x}^2+{\sigma_y}^2} }
$$
元々の画像 \(x\) のSN比は \(\langle x\rangle/\sigma_x\)だから、ダーク減算によってSN比は必ず小さく(悪く)なる。

nFAQ:じゃあなんでダーク減算するの?
たしかにダーク減算によってノイズは必ず増加する。しかしそれは画像に含まれるランダムノイズだけに着目した場合にそう見えるだけであって、本来のダーク減算の目的は固定ノイズの除去である。固定ノイズ除去のメリットが、ランダムノイズ増加のマイナスを上回るなら、ダーク減算は常に行う価値がある。またほとんどの場合、信号の輝度揺らぎに由来するライト画像のショットノイズは、ダーク画像のもつ熱ノイズよりもずっと大きい。ダーク減算によるランダムノイズの増加は、実は問題にならないくらい小さい。

フラット除算した画像のSN比

ライトフレーム画像を \( x \),フラットフレーム画像を \(y\) として \(z=x/y\) とした \(z\) のノイズの大きさを考える。誤差の伝搬法則を用いて計算すると
$$
{\sigma_z}^2=\frac{1}{\langle y\rangle^2}{\sigma_x}^2+\frac{\langle x \rangle^2}{\langle y\rangle^4}{\sigma_y}^2
$$
となる。ただし,実際にはフラット除算の前後で平均輝度が一定になるようにノーマライズが行われる。それは大雑把には上の結果に\(\langle y\rangle\)をかけてあげればいいので,実際のフラット除算後のノイズの大きさは
$$
\langle y\rangle\sigma_z=\langle x \rangle\sqrt{\frac{{\sigma_x}^2}{\langle x\rangle^2}+\frac{{\sigma_y}^2}{\langle y\rangle^2}}
$$
となる。よってフラット除算後の画像のSN比 \((SNR)_{flat}\) は、
$$
(SNR)_{flat}=\frac{\langle x\rangle}{\langle y\rangle\sigma_z}=\dfrac{1}{\sqrt{\frac{{\sigma_x}^2}{\langle x\rangle^2}+\frac{{\sigma_y}^2}{\langle y\rangle^2}}}
$$
見通しを良くするためにさらにちょっと変形すると,フラット除算前のSN比 \(\frac{\langle x\rangle}{\sigma_x}\) に対して,フラット除算後のSN比は
\begin{align}
(SNR)_{flat}&=\frac{\langle x\rangle}{{\sigma_x}}\frac{1}{\sqrt{1+\frac{{\sigma_y}^2}{{\sigma_x}^2}\frac{\langle x\rangle^2}{\langle y\rangle^2}}}\\
&=(SNR)_x\frac{1}{\sqrt{1+\frac{{(SNR)_x}^2}{{(SNR)_y}^2}}}
\end{align}
と書くことができる。つまり、フラット除算によってノイズは必ず増える。\((SNR)_y\) が \((SNR)_x\) に比較して大きければ(フラット画像が良いSN比を持てば)上の平方根の中身が1に近づくので、フラット除算によるノイズ増加は抑えられる(これはアタリマエ)。反対に、\((SNR)_y\) が \((SNR)_x\) に比較して小さければ(ライト画像のSN比が悪い)、マスターフラットの質を高めてもあまり意味がない(ほんとか?)。

具体例で計算してみる(未完)

具体的な画像データをもとに、ダーク減算やフラット除算が、キャリブレーション後のライト画像にどのような影響を与えるか評価してみよう(近日いつか追記予定)。

謝辞

文中の「誤差の伝搬法則」については、Twitterにてお付き合いのあるRamb氏にご教示いただいた。またBlog上の\(\rm\TeX \) 環境はぐらすのすち氏が用意してくれた。若いお二人に感謝いたします。

コメント

タイトルとURLをコピーしました