ODOGカーネルのフーリエ変換

ここではOriented Difference of Gaussians (ODOG)という視覚モデル^[1]のフーリエ変換を手計算して確認してみる． ODOGは，明度知覚に関する様々な効果をうまく説明できるモデルとして広く知られている． ODOGモデルは1枚の画像が与えっれたとき，人間が知覚する明度を画素値として持つ画像に変換しフィルタとして機能する数式モデル． このモデルは，初等関数の簡単な組み合わせでできているので，簡単にフーリエ変換できる． フーリエ変換するメリットは，このフィルタを入力画像に畳み込む計算を，フーリエ変換すると単なる掛け算になって，高速化できること． とりあえずそうした背景は一切忘れて，ODOGのカーネルをフーリエ変換してみる．

ODOGのカーネル関数って？

その肝心のODOGのモデルだけど，全ての数式を説明するっだけでそれなりのスペースが必要になるので，ここではODOGのカーネル関数だけをフーリエ変換してみる． ODOGのカーネル関数には，角度とスケールのパラメータがあり，様々な角度とスケールのカーネル関数でフィルタした結果を足し合わせて最終的な画像を出力する仕組みになっている． そのカーネル関数$f(y_1, y_2)$を数式で表すと次のようになる．

$$$$ $$\begin{array}{rcl}f(y_1,y_2)& =& \frac{1}{2\pi {\sigma }_1^2}{e}^{-\frac{y_1^2+y_2^2}{2{\sigma }_1^2}}-\frac{1}{2\pi {\sigma }_1{\sigma }_2}{e}^{-\frac{1}{2}(\frac{y_1^2}{{\sigma }_2^2}+\frac{y_2^2}{{\sigma }_1^2})}\end{array}$$

ただし，画素の位置を表す変数は$(x,y)$や$(u,v)$ではなくて，論文[2]に合わせて$(y_1,y_2)$にした．

この式は，異なる2つのガウス関数の差に見える． 実際に$y_1$-軸方向の断面を見るために$y_2=0$を代入すると次のように2つの異なるガウス関数の差になっている．

$$\begin{array}{rcl}f(y_1,0)& =& \frac{1}{2\pi {\sigma }_1^2}{e}^{-\frac{y_1^2}{2{\sigma }_1^2}}-\frac{{\sigma }_2}{{\sigma }_1}\frac{1}{2\pi {\sigma }_2^2}{e}^{-\frac{y_1^2}{2{\sigma }_2^2}}\end{array}$$

一方，$y_2$-軸方向の断面を見るために$y_1=0$を代入すると次のように1つの異なるガウス関数になっている．

$$\begin{array}{rcl}f(0,y_2)& =& \frac{1}{2\pi {\sigma }_1^2}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}-\frac{1}{2\pi {\sigma }_1{\sigma }_2}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}=(1-\frac{{\sigma }_1}{{\sigma }_2})\frac{1}{2\pi {\sigma }_1^2}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}\end{array}$$

このように，普通のDoG^[3]と異なり，方向をもつ構造になっている．

画像に対する二次元フーリエ変換

教科書でよく出てくるのは一次元の信号$f(t)$なんかのフーリエ変換． $f$のフーリエ変換を$\mathcal{F}[f]$と表すと，周波数$\omega$を用いて次のように定義されている．

$$\mathcal{F}[f](\omega )={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}f(t)e^{-i\omega t}d\omega$$

画像に対する二次元フーリエ変換を計算するには，$x$方向と$y$方向（正確にはここでは$y_1$方向と$y_2$方向）の一次元のフーリエ変換を二重にかければよい． フーリエ変換は積分なので二重積分になる． 周波数空間($s_1,s_2$)でのカーネル関数を$F(s_1,s_2)$とすると次のようになる．

$$\begin{array}{r}F(s_1,s_2)={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}({\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}f(y_1,y_2)e^{-i{y}_1{s}_1}d{y}_1)e^{-i{y}_2{s}_2}d{y}_2\end{array}$$

ここでも論文[2]に合わせて周波数の変数は$(s_1,s_2)$とした．

足し算のフーリエ変換は分けることができる

式（$\text{???}$）の（）の中の一次元フーリエ変換をまず計算してみる． そのために式（$\text{???}$）を$f(y_1, y_2)$に代入する． 式（$\text{???}$）は，$a-b$の形になっている．フーリエ変換の大事な性質は線形性． つまり$a$と$b$それぞれのフーリエ変換を$A$と$B$だとすると，$a-b$のフーリエ変換は単純に$A-B$になるので， 式（$\text{???}$）のフーリエ変換は次のように各項単独のフーリエ変換の差になる．

$$\begin{array}{rcl}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}f(y_1,y_2)e^{-i{y}_1{s}_1}d{y}_1& =& {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{1}{2\pi {\sigma }_1^2}{e}^{-\frac{y_1^2+y_2^2}{2{\sigma }_1^2}}{e}^{-i{y}_1{s}_1}d{y}_1\\ & -& {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{1}{2\pi {\sigma }_1{\sigma }_2}{e}^{-\frac{1}{2}(\frac{y_1^2}{{\sigma }_2^2}+\frac{y_2^2}{{\sigma }_1^2})}{e}^{-i{y}_1{s}_1}d{y}_1\\ & =& \frac{1}{2\pi {\sigma }_1^2}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{e}^{-\frac{y_1^2}{2{\sigma }_1^2}}{e}^{-i{y}_1{s}_1}d{y}_1\\ & -& \frac{1}{2\pi {\sigma }_1{\sigma }_2}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{e}^{-\frac{y_1^2}{2{\sigma }_2^2}}{e}^{-i{y}_1{s}_1}d{y}_1\end{array}$$

積分変数$y_1$を含まない係数は，積分の外に出して単純にした． そうすると積分しないといけない，ヤバそうな部分は結構シンプルになる． 式（$\text{???}$）と（$\text{???}$）の2つあるけど，よく見ると$\sigma_1$と$\sigma_2$が違うだけで全く同じ形．

積分の部分を計算する

その肝心の積分変数は，「ガウス関数のフーリエ変換はガウス関数になる」という公式が使えて実際には積分計算をする必要がない．

$$\begin{array}{r}\mathcal{F}\left[\frac{1}{\sqrt{2\pi }\sigma }{e}^{-t^2/2{\sigma }^2}\right](\omega )=\frac{1}{\sqrt{2\pi }\sigma }{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{e}^{-t^2/2{\sigma }^2}{e}^{-i\omega t}dt=e^{-{\sigma }^2{\omega }^2/2}\end{array}$$

この公式の変数を置換するだけで，式（$\text{???}$）と（$\text{???}$）の積分の部分を計算できる．

$$\begin{array}{r}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{e}^{-\frac{y_1^2}{2{\sigma }_1^2}}{e}^{-i{y}_1{s}_1}d{y}_1=\sqrt{2\pi }{\sigma }_1{e}^{-\frac{{\sigma }_1^2{s}_1^2}{2}}\\ {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{e}^{-\frac{y_1^2}{2{\sigma }_2^2}}{e}^{-i{y}_1{s}_1}d{y}_1=\sqrt{2\pi }{\sigma }_2{e}^{-\frac{{\sigma }_2^2{s}_1^2}{2}}\end{array}$$

一重目のフーリエ変換を計算

これで， 式（$\text{???}$）に上の結果を代入すると， 式（$\text{???}$）の内側のフーリエ変換はちゃんと計算できる．

$$\begin{array}{rcl}{\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}f(y_1,y_2)e^{-i{y}_1{s}_1}d{y}_1& =& \frac{1}{2\pi {\sigma }_1^2}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}\sqrt{2\pi }{\sigma }_1{e}^{-\frac{{\sigma }_1^2{s}_1^2}{2}}\\ & -& \frac{1}{2\pi {\sigma }_1{\sigma }_2}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}\sqrt{2\pi }{\sigma }_2{e}^{-\frac{{\sigma }_2^2{s}_1^2}{2}}\\ & =& \frac{1}{\sqrt{2\pi }{\sigma }_1}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}{e}^{-\frac{{\sigma }_1^2{s}_1^2}{2}}\\ & -& \frac{1}{\sqrt{2\pi }{\sigma }_1}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}{e}^{-\frac{{\sigma }_2^2{s}_1^2}{2}}\\ & =& \frac{1}{\sqrt{2\pi }{\sigma }_1}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}(e^{-\frac{{\sigma }_1^2{s}_1^2}{2}}-e^{-\frac{{\sigma }_2^2{s}_1^2}{2}})\end{array}$$

二重目のフーリエ変換を計算

ようやく計算し終えた一重目のフーリ変換（式（eq:15））を元の二次元フーリ変換の式（$\text{???}$）に代入する．

$$\begin{array}{rcl}F(s_1,s_2)& =& {\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}\frac{1}{\sqrt{2\pi }{\sigma }_1}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}(e^{-\frac{{\sigma }_1^2{s}_1^2}{2}}-e^{-\frac{{\sigma }_2^2{s}_1^2}{2}})e^{-i{y}_2{s}_2}d{y}_2\\ & =& \frac{1}{\sqrt{2\pi }{\sigma }_1}(e^{-\frac{{\sigma }_1^2{s}_1^2}{2}}-e^{-\frac{{\sigma }_2^2{s}_1^2}{2}}){\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{e}^{-\frac{y_2^2}{2{\sigma }_1^2}}{e}^{-i{y}_2{s}_2}d{y}_2\end{array}$$

また，積分の部分がガウス関数のフーリ変換の形になっているので，一重目と同様にして公式を使うと計算できる．

$$\begin{array}{rcl}F(s_1,s_2)& =& \frac{1}{\sqrt{2\pi }{\sigma }_1}(e^{-\frac{{\sigma }_1^2{s}_1^2}{2}}-e^{-\frac{{\sigma }_2^2{s}_1^2}{2}})\sqrt{2\pi }{\sigma }_1{e}^{-\frac{{\sigma }_1^2{s}_2^2}{2}}\\ & =& e^{-\frac{{\sigma }_1^2{s}_2^2}{2}}(e^{-\frac{{\sigma }_1^2{s}_1^2}{2}}-e^{-\frac{{\sigma }_2^2{s}_1^2}{2}})\end{array}$$

論文中のフーリエ変換の式と比較

論文[2]の式は次のようになっていて，上の式と微妙に違うので，正しいの確認してみる．

$$\begin{array}{rcl}F(s_1,s_2)& =& e^{-2{\pi }^2{\sigma }_1^2{t}_2^2}(e^{-2{\pi }^2{\sigma }_1^2{t}_1^2}-e^{-2{\pi }^2{\sigma }_2^2{t}_1^2})\end{array}$$

まず$(t_1, t_2)$というのが出てくるけど，実際の論文では回転を考慮しているので，$(s_1, s_2)$を回転させた$(t_1,t_2)$というのが出てくる． 予め回転させた画像を与えると，ここは関係ないので，これは問題ない． あと，$2\pi^2$というなぞの係数が登場するけど，これはフーリエ変換の定義で周波数を$\omega=2\pi k$のように， 単位画素あたり位相$\omega$ではなてく，波数$k$で定義すると現れる係数．要するに単位の違い． 実際に式（$\text{???}$）に$s_1=2\pi k_1$と$s_2 =2\pi k_2$を代入すると，式（$\text{???}$）と同じ形になる．

参考文献

1. Barbara Blakeslee and Mark E. McCourt: "A multiscale spatial filtering account of the White effect, simultaneous brightness contrast and grating induction," Vision Research, vol. 39, issue 26, pp. 4361-4377, 1999.
   doi: 10.1016/S0042-6989(99)00119-4
2. Barbara Blakeslee, Davis Cope, and Mark E. McCourt1: "The Oriented Difference of Gaussians (ODOG) Model of Brightness Perception: Overview and Executable," Behavior Research Methods, pp. 306–312, vol. 48, no. 1, 2016.
   doi: 10.3758/s13428-015-0573-4
3. "Difference of Gaussians" in Wikipedia
   https://en.wikipedia.org/wiki/Difference_of_Gaussians