ImageMagick 使用例 -- フーリエ変換

ImageMagick 使用例 前書きと索引

はじめに

フーリエ変換

• 画像は波である

• 画像内の2次元の波

  ImageMagick での FFT/IFT

• FFT の使い方

• 振幅のみ・位相のみの画像
• FFT スペクトル画像
• HDRI FFT 画像

• 実部・虚部成分としての FFT

  フーリエ変換の性質

• 定数画像の FFT

• DC カラーの効果
• サイン波画像のスペクトル
• FFT 画像を直接生成する
• 矩形パターン画像のスペクトル
• 平坦な円形パターン画像のスペクトル
• ガウス円形パターン画像のスペクトル

• グリッドパターン画像のスペクトル

  実用例

• 画像のコントラスト変更 - 係数のべき乗

• 画像をぼかす - ローパスフィルタ
• 画像のエッジ検出 - ハイパスフィルタ
• 画像のシャープ化 - ハイブーストフィルタ

• ノイズ除去 - ノッチフィルタ

  発展的な応用

  FFT 乗算と除算 (低レベルの例 - サブページ)

はじめに

画像処理で最も理解しにくい概念の一つがフーリエ変換です。理由は二つあります。第一に数学的に高度であること、第二に得られる画像が元の画像に似ておらず、解釈が難しいことです。それでもフーリエ変換を活用すれば、明るさやコントラストの強調、ぼかし、シャープ化、ノイズ除去といった馴染みのある処理を新しい方法で行えます。さらに、通常の画像領域ではできない新たな能力も得られます。たとえば、モーションブラーやレンズのピントぼけといった典型的なカメラの歪みのデコンボリューション（デブラーとも呼ばれます）や、正規化相互相関を用いた画像マッチングなどです。このページの目的は、フーリエ変換の背景と簡略化した数学を説明し、フーリエ変換を使って行える処理の例を示すことです。難しすぎると感じたら、読み飛ばして性質と例に集中してもかまいません。ImageMagick での FFT/IFT から始めてください。興味のある方には、光学とのアナロジーを含む別の分かりやすい解説が An Intuitive Explanation of Fourier Theory にあります。ヴァンダービルト大学工学部の講義ノートも、より数学的な思考を好む方にとって非常に有益です: 1 & 2 Dimensional Fourier Transforms と Frequency Filtering。その他の数学的な参考資料には、Wikipedia の Fourier Transform、Discrete Fourier Transform、Fast Fourier Transform のページや、Complex Numbers があります。元のデモのコーディングを行った Sean Burke 氏と、それを ImageMagick に統合した ImageMagick の作者に感謝します。どちらも英雄的な働きでした。 多くの例では、変換後の画像の精度を保つために必要な HDRI 版 ImageMagick を使用します。これらの技法を最大限に活用したい場合は、自分用に HDRI 版をコンパイルすることをお勧めします。

フーリエ変換

通常、画像は「ピクセル」の配列で構成され、各ピクセルは赤・緑・青、ときには透明度といった値の集合で定義されます。しかしここでの目的のためには透明度は無視します。したがって赤・緑・青の各「チャンネル」は、それぞれ「強度」または「グレースケール」の値の集合を含みます。これは「空間領域」のラスター画像として知られています。これは要するに、画像が各「位置」あるいは「空間内の位置」で持つ「強度値」によって定義されている、という洒落た言い方にすぎません。しかし画像は別の方法でも表現でき、それは画像の「周波数領域」として知られています。この領域では、各画像チャンネルが正弦波（サイン波）によって表現されます。このような「周波数領域」では、各チャンネルは「振幅」値を持ち、それらは X,Y の「空間」座標ではなく X,Y の「周波数」に基づく位置に格納されます。これはデジタル表現なので、周波数は「最小」あるいは単位周波数の倍数となり、ピクセル座標はこの単位周波数の指数すなわち整数倍を表します。これは「適切に振る舞う任意の関数は、正弦波の重ね合わせ（組み合わせまたは和）によって表現できる」という原理から導かれます。言い換えれば、「周波数領域」表現は「空間領域」の画像を格納し再現するもう一つの方法にすぎません。しかし画像をどうやって「波」として表現できるのでしょうか?

画像は波である

任意の 画像から1行または1列のピクセルを取り出してグラフ化すると（"gnuplot" を用い、スクリプト "[im_profile](../static/img/scripts/im_profile)" で生成）、それはかなり波のように見えることが分かります。

  magick holocaust_tn.gif -colorspace gray miff:- |\
    im_profile -s - image_profile.gif

もし変動の間隔や振幅がもっと規則的であれば、次のような波のパターンに近いものが得られます。

  magick -size 20x150 gradient: -rotate 90 \
          -function sinusoid 3.5,0,.4   wave.gif
  im_profile -s wave.gif wave_profile.gif

しかし、この規則的な波のパターンは上に示した画像のプロファイルにある程度似ているものの、規則的すぎます。しかし、より多くの波を足し合わせれば、画像から得られたパターンにさらに近いパターンを作ることができます。

  magick -size 1x150 gradient: -rotate 90 \
          -function sinusoid 3.5,0,.25,.25     wave_1.png
  magick -size 1x150 gradient: -rotate 90 \
          -function sinusoid 1.5,-90,.13,.15   wave_2.png
  magick -size 1x150 gradient: -rotate 90 \
          -function sinusoid 0.6,-90,.07,.1    wave_3.png

  magick wave_1.png wave_2.png wave_3.png \
          -evaluate-sequence add added_waves.png

上記の代替例については Adding Biased Gradients も参照してください。この「波の重ね合わせ」（波の加算）はずっと近づきましたが、それでも画像のパターンと正確には一致しません。しかし、このように波を追加し調整し続ければ、得られる合成波は元の画像の実際のプロファイルにどんどん近づいていきます。最終的に十分な数の波を足し合わせれば、画像の元のプロファイルを正確に再現できます。これが数学者 ジョゼフ・フーリエ によってなされた発見でした。その現代的な解釈は「適切に振る舞う任意の関数は、正弦波の重ね合わせによって表現できる」というものです。言い換えれば、ちょうど適切な周波数と振幅を持つ十分な数のサイン波を足し合わせれば、どんな変動するパターンも再現できるのです。したがって、「周波数領域」表現は「空間領域」の画像を格納し再現するもう一つの方法にすぎません。「フーリエ変換」とは、上の例で行ったように、画像がどの「波」で構成されているかを割り出す処理のことです。

画像内の2次元の波

上では、画像の1行のプロファイルを複数のサイン波でどのように近似できるかの一例を示しました。しかし画像は2次元であり、そのため画像を「周波数領域」で表現するために使う波も2次元である必要があります。次にそのような2次元の波の一例を示します。この波にはいくつかの構成要素があります。 画像の例

ImageMagick での FFT/IFT の使用

実装上の注意

ImageMagick は FFTW（離散フーリエ変換ライブラリ） を利用しています。これは画像を浮動小数点値（複素数）に変換し、また元に戻すことを必要とし、IM バージョン 6.5.4-3 で初めて実装されました。一般に画像に対して人々が期待する動作をさせるため、正方形でない画像や奇数の寸法を持つ画像は、画像の幅と高さの最大値の正方形になるよう（仮想ピクセル を用いて）パディングされます。「FFT の原点」を画像の中心に正しく配置できるようにするため、寸法は偶数（2の倍数）にも強制されます。この結果として、逆フーリエ変換を適用した後、画像はパディングを取り除くために元の寸法に切り戻す必要があります。フーリエ変換は「複素数」で構成されているため、変換の結果を直接可視化することはできません。そのため、複素変換は次の2つの形式のいずれかで2つの成分画像に分離されます。
複素数
実部/虚部

実部と虚部

「複素数」の通常の数学的・数値的表現は、「実部」(a) と「虚部」(b) からなる浮動小数点値のペアです。残念ながらこれらの2つの数は負の値を含むことがあるため、表示可能な画像を形成しません。そのため、この表現は通常版の ImageMagick では使用できず、そのような画像はクリッピングされてしまいます（その結果の効果の例は後述）。しかし HDRI 版 ImageMagick を使えば、フーリエ変換された画像のこの表現を生成し、使用し、さらには保存することができます。それらは画像としては有用でないどころか表示すらできないかもしれませんが、多くの数学的演算を適用することは可能です。この表現を生成するには、演算子の「プラス」形式である "[+fft](https://usage.imagemagick.org/fourier/option_link.cgi#fft)" と "[+ift](https://usage.imagemagick.org/fourier/option_link.cgi#ift)" を使用します。これは後で 実部・虚部成分としての FFT で詳しく見ていきます。
複素極形式
振幅/位相

振幅と位相

「複素数」の直接的な数値表現は画像処理にはあまり役立ちません。しかしその値を2次元平面上にプロットすれば、値を「振幅」(r) と「位相」(θ) 成分からなる 極形式表現 に変換できます。この形式は画像処理に非常に役立ち、特に振幅成分は本質的に画像を構成するすべての周波数を指定します。「振幅」成分は正の値のみを含み、そのまま画像の値に直接マッピングされます。値に固定された範囲はありませんが、DC すなわちゼロ周波数の色を除いて、値は一般にかなり小さくなります。この結果として、振幅画像は一般に非常に暗く（ほとんど真っ黒に）見えます。視覚的なディテールを引き出すには、通常、振幅をスケーリングしてその強度値に対数変換を適用する必要があります。得られた「対数変換された」振幅画像は、画像の「スペクトル」として知られています。ただし、逆変換に使用すべきなのは「スペクトル」画像ではなく「振幅」画像であることを忘れないでください。画像の中央の「原点」に現れる DC（「直流 (Direct Current)」の略）すなわち「ゼロ周波数」の色は、画像全体の平均的な色の値になります。また入力画像は「虚部」成分を含まないため、DC の位相値も常にゼロ位相となり、純粋なグレー色を生み出します。一方、「位相」成分は -π から +π の範囲です。これはまず 0 から 2π の範囲にバイアスされ、次に 0 から QuantumRange（コンパイル時のメモリ品質 によって決まる）までの実際の画像値にスケーリングされます。この結果、ゼロ位相は（各チャンネルにとって適切な）純粋なグレー値となり、負の位相は純粋な黒（'0'）の値になります。なお、純粋な白（'_QuantumRange_'）はほとんど同じですが完全に同じではないことに注意してください。画像の振幅と位相による FFT 表現は、通常の FFT 演算子 "[+fft](https://usage.imagemagick.org/fourier/option_link.cgi#fft)" と "[+ift](https://usage.imagemagick.org/fourier/option_link.cgi#ift)" を使って生成されます。これは最初に FFT 画像とその逆変換の生成 で見ていきます。

FFT 画像とその逆変換の生成

(振幅と位相)

では、単純に Lena 画像でフーリエ変換の往復を試してみましょう。つまり、順変換を行ってから直ちに逆変換を適用し、元の画像を取り戻します。そして結果を比較し、得られる品質のレベルを見てみます。

  time magick lena.png -fft -ift lena_roundtrip.png


  echo -n "RMSE = "

  magick compare -metric RMSE lena.png lena_roundtrip.png null:

  echo -n "PAE = "

  magick compare -metric PAE  lena.png lena_roundtrip.png null:

上記の "[compare](basics.html#compare)" プログラムは、2つの画像がどれだけ異なるかの尺度を返します。この場合、全体的な差は非常に小さく、約 0.22% であることが分かります。少なくとも1つのピクセルにおけるピーク値の差（"PAE"、Peak Absolute Error）はわずか約 1% です。これは HDRI 版 ImageMagick を使うことで改善できます。（後述の HDRI を使った FFT を参照）。上記の往復で生成された FFT 画像をもう少し詳しく見てみましょう。

  magick lena.png -fft    +depth +adjoin lena_fft_%d.png

元画像 | |
振幅 |
位相
---|---|---|---

John M. Brayer はフーリエ変換について次のように述べています... 我々は通常 PHASE（位相）画像を表示しない。なぜなら、それを見た人のほとんどが、まもなく幻覚剤に手を出すか、チベットの僧院に行き着いてしまうからだ。 "[-fft](https://usage.imagemagick.org/fourier/option_link.cgi#fft)" 演算子は2つの画像を生成し、最初の画像は「振幅」成分（そう、これはほとんど真っ黒で、中央に1つの色付きの点があるだけです）であり、2番目のほとんどランダムに見える画像は「位相」成分を含むことに注意してください。PNG 画像は1ファイルにつき1画像しか保存できないため、"[+adjoin](https://usage.imagemagick.org/fourier/option_link.cgi#adjoin)" も出力ファイル名の '%d' も実際には必要ありませんでした。IM がこれを処理してくれるからです。しかし、上記では完全性のためにこれらのオプションを含めており、1つではなく2つの別々の画像ファイルを生成していることを明確にしています。詳細は 複数画像のシーケンスを書き出す を参照してください。2つの画像が生成されるので、振幅画像（最初すなわち0番目の画像）は "lena_fft_0.png" に、位相画像（2番目の画像）は "lena_fft_1.png" に保存されます。 | _FFT 画像を保存することで生じるおそれのある歪みを防ぐため、それらをディスクに保存せず、画像を処理する間メモリ内に保持しておくのが最善です。

保存しなければならない場合は、画像を最高品質（ビット深度）で保持するために Magick ファイル形式 "[MIFF](https://usage.imagemagick.org/files/#miff")" を使うのが最善です。この形式は複数の画像を1ファイルに保存することもできます。スクリプト作業では、冗長な "[TXT](files.html#txt)" 列挙ピクセル形式も使えます。

"[JPEG](formats.html#jpg)"、"[GIF](formats.html#gif)" 画像形式を使って保存しては いけません。

これらの画像を、Web ブラウザでの実際の表示などのためにファイルに保存しなければならない場合は、"[+depth](https://usage.imagemagick.org/fourier/option_link.cgi#depth)" を内部デフォルトにリセットした "[PNG](formats.html#png)" 画像形式を使ってください（これらの例で行っているように）。ただし1ファイルにつき1画像しか保存できません。

"[TIFF](formats.html#tiff)" ファイル形式も使えますが、Web ブラウザにはあまり適していません。とはいえ1ファイルに複数画像を保存できます。

_
---|---
複数の中間画像を1つのファイルに保存する最善の方法は、"[MIFF](https://usage.imagemagick.org/files/#miff")" ファイル形式を使うことです...

  magick lena.png -fft  +depth lena_fft.miff

あるいは、"[-write](https://usage.imagemagick.org/fourier/option_link.cgi#write)" を使って完全に別々のファイル名に保存することもできます（画像を書き出す を参照）...

  magick lena.png -fft  +depth \
          \( -clone 0 -write lena_magnitude.png +delete \) \
          \( -clone 1 -write lena_phase.png +delete \) \
          null:

上記では、特別な "[NULL:](files.html#null)" 画像形式を使って2つの画像を捨てていますが、これらはさらなる処理のためメモリ内に保持されていることに注意してください。そして最後に、再び2つの画像を読み込んで、通常の「空間」画像に変換し直します...

  magick lena_magnitude.png lena_phase.png -ift lena_restored.png

FFT 処理によって生成される両方の画像は変更に非常に敏感で、わずかな変更でも結果が大きく歪んでしまうことがあります。そのため、それらの値を歪めるおそれのある画像形式で決して保存しないことが重要です。周波数領域から画像を復元する際には両方の画像が必要であることを覚えておくのも重要です。したがって、画像の再構成に使うつもりなら、一方の画像を保存して他方を捨てるのはよくありません。

振幅のみ・位相のみの画像

最後に、振幅成分のみ、または位相成分のみから画像を再構成してみましょう。

  magick lena_fft_0.png  -size 128x128 xc:'gray(50%)' \
                                                  -ift lena_magitude_only.png

  magick -size 128x128 xc:gray1  lena_fft_1.png  -ift lena_phase_only.png

振幅のみ |
位相のみ
---|---

これから分かるのは、画像の位置情報のほとんどを実際に含んでいるのは位相画像であり、一方で振幅は色情報の多くを保持しているということです。情報に多少の重なりがあるため厳密ではありませんが、概してそうなります。「振幅のみ」の画像は、一定の 50% 位相画像を使ったため、常に白い角を持ちます。これらの白い斑点は、ランダム化した位相画像を使うことで取り除けます。ただし、中央ピクセルの位相が完璧な 50% グレーであることを確かめてください。さもないと画像全体が暗くなってしまいます。「位相のみ」の画像は、変換のために一定の 1% グレー（ほぼ純粋な黒）の振幅画像を使いました。この一定の振幅でも、特にエッジに沿って非常に強いピクセルの斑点が生じます。元の画像を再構成するには両方の画像が必要であることを覚えておく必要があります。

周波数スペクトル画像

振幅画像（最初すなわち0番目の画像）がほとんど完全に黒く見えることにお気づきでしょう。実際にはそうではないのですが、私たちの目にはすべての値が非常に非常に小さく見えるのです。このような画像は研究のために眺めてもあまり面白くないので、対数変換で結果を強調して「周波数スペクトル」画像を生成しましょう。これは 正規化 した「振幅」画像に強い Evaluate 対数変換 を適用することで行います。

  magick lena_fft_0.png -auto-level -evaluate log 10000 \
          lena_spectrum.png

これで振幅画像のスペクトル版でディテールが見えるようになりました。スペクトル画像にいくつかの特定の色が見えるかもしれませんが、一般にそうした色はスペクトル画像では重要ではありません。各周波数の全体的な強度と、それらが生み出すパターンの方がはるかに重要です。そのため、強調後にスペクトル画像をグレースケール化したい場合もあるでしょう。どれだけの対数強調を使う必要があるかは画像によるので、画像の周波数スペクトルのパターンをはっきり見るのに必要なディテールが得られるまで調整するとよいでしょう。 あるいは、次の小さなシェルスクリプトを使って、特定の振幅画像に使う 対数スケーリング係数 を計算することもできます。 |

  scale=`magick lena_fft_0.png -auto-level \
          -format "%[fx:exp(log(mean)/log(0.5))]" info:`
  magick lena_fft_0.png -auto-level \
          -evaluate log $scale    lena_spectrum_auto.png

ただし、スペクトル画像は逆変換 "[-ift](https://usage.imagemagick.org/fourier/option_link.cgi#ift)" には使えないことを覚えておいてください。過度に明るい画像が生成されてしまうからです。 |

  magick lena_spectrum.png lena_fft_1.png -ift lena_roundtrip_fail.png

基本的に「振幅」画像を強調したことで、得られる画像も同じように強調され、示されているようなひどく「クリッピングされた」結果になってしまいます。

HDRI FFT 画像

フーリエ変換の結果を画像表現にマッピングする際、私たちは値を浮動小数点の「複素数」から整数の画像値へとスケーリングし変換しました。これは当然 丸め誤差やその他の「量子化」効果 を生み、特に小さな低周波の振幅で顕著です。画像処理において精度が重要な場合は、ビット品質（Q32 や Q64 ビット版 ImageMagick など）を使うか、さらに良いのは値が浮動小数点数として格納される HDRI 版 ImageMagick を使う必要があります。HDRI 版 IM をフーリエ変換の 振幅と位相 表現とともに使う場合、振幅成分は依然としてすべて正の値となり、上で示したように使えますが、はるかに正確になります。一方、位相成分は依然として前述のとおりバイアスとスケーリングが行われます。言い換えれば、HDRI での振幅と位相の表現はまったく同じですが、はるかに正確です。
たとえば、ここでは HDRI 版 ImageMagick を使って、別の画像の「往復」変換を生成します。 |

  # HDRI version of IM used
  time magick lena.png -fft -ift lena_roundtrip_hdri.png


  echo -n "RMSE = "

  magick compare -metric RMSE lena.png lena_roundtrip_hdri.png null:

  echo -n "PAE = "

  magick compare -metric PAE  lena.png lena_roundtrip_hdri.png null:

上記の結果を、先ほどの非 HDRI の比較と見比べると...

実部・虚部成分としての FFT

ここまでは、フーリエ変換された画像の「振幅」と「位相」表現のみを見てきました。しかし HDRI 版 IM をコンパイルしていれば、浮動小数点の「実部」と「虚部」成分を使って画像を処理することもできます。これはオプション "[+fft](https://usage.imagemagick.org/fourier/option_link.cgi#fft)" と "[+ift](https://usage.imagemagick.org/fourier/option_link.cgi#ift)" の「プラス」版を使うことで行います。たとえば、ここでは HDRI 版 IM を使って画像の「往復」FFT も行いますが、今回は実部/虚部画像を生成します。 |

  # HDRI version of IM used
  time magick lena.png   +fft +ift   lena_roundtrip_ri.png


  echo -n "RMSE = "

  magick compare -metric RMSE lena.png lena_roundtrip_ri.png null:

  echo -n "PAE = "

  magick compare -metric PAE  lena.png lena_roundtrip_ri.png null:

プラス形式を使って実部/虚部 FFT 画像を生成するときは、必ず HDRI 版を使わなければなりません。そうしないと、約半分の値がゼロになり、「汚れて」見える画像になってしまいます。たとえば... |

  # non-HDRI Q16 version of IM used  -- THIS IS BAD
  magick lena.png   +fft +ift   lena_roundtrip_ri_bad.png

もう一つ覚えておくべきことは、どの形式の FFT 画像を生成するかによって、その FFT 画像に適用したいすべての画像処理演算も影響を受けるということです。それらは非常に異なる画像であり、そのため非常に異なる方法で、異なる数学的演算によって処理しなければなりません。また以前と同様に、最終的な画像を復元するには実部と虚部の両方の成分画像が必要です。たとえば、成分の一方を「黒」画像で置き換えるとどうなるかを示します。

  # HDRI version of IM used
  magick lena.png +fft -delete 1 \
          -size 128x128 xc:black +ift lena_real_only.png
  magick lena.png +fft -delete 0 \
          -size 128x128 xc:black +ift lena_imaginary_only.png

実部のみ |
虚部のみ
---|---

これから、実部/虚部の両方の FFT 画像が元の画像についての重要な情報をほぼ等しく含んでいることが分かります。両者の最大の違いは、特別な DC すなわち「平均色」には虚部成分がないため、振幅画像にしか存在しないという点です。両方の画像に見られる対角線方向の鏡像（実際には 180 度回転）効果は、もう一方の成分に含まれる「符号」情報の損失によって生じます。もう一方の成分がないと、波は 180 度位相がずれていると考えることもでき、この奇妙な見た目を生み出します。この情報損失は2種類の画像で等しくなっています。

フーリエ変換の性質

定数画像の FFT

これらの性質のいくつかを実演してみましょう。まず単純に定数色の画像を取り、その振幅を求めます。

  magick -size 128x128 xc:gold constant.png
  magick constant.png -fft +delete constant_magnitude.png

この場合の振幅画像は、画像のちょうど中央、ピクセル位置 width/2, height/2 にある1つの色付きピクセルを除いて、本当に純粋な黒であることに注意してください。このピクセルは画像のゼロ周波数すなわち DC（「直流 (Direct Current)」）の値であり、サイン波を表さない唯一のピクセルです。言い換えれば、この値はまさに FFT の定数成分なのです! この単一ピクセルをよりはっきり見るために、その画像領域を拡大してみましょう... |

  magick constant_magnitude.png -gravity center -extent 5x5 \
           -scale 2000% constant_dc_zoom.gif

DC 点の色が元画像と同じであることに注意してください。実際、見えているものは3つの値であることを覚えておくとよいでしょう。つまり、生成された画像は実際には3つの別々の高速フーリエ変換です。赤・緑・青の各画像チャンネルごとに1つの FFT があります。FFT 自体は色について実際の知識を持たず、色の値すなわち「グレーレベル」だけを扱います。実際、FFT 変換はほとんどどんな色空間にも適用できます。なぜなら本当に... それは気にしないからです! フーリエ変換にとって画像は単なる値の配列であり、それだけです。 | DC 値の「位相」は重要ではありませんが、常に「ゼロ」角（50% グレーの位相色値）であるべきです。50% グレーに設定されていないと、DC 値は「非実数」成分を持ち、その値が与えられた角度によって変調されてしまいます。
---|---

DC カラーの効果

より典型的な非定数の画像では、DC 値は画像の平均色です。これは、画像を完全にぼかしたり、平均化したり、単一のピクセルや色までリサイズした場合に一般に得られる色です。たとえば「Lena」画像の FFT から DC ピクセルを抽出してみましょう。

  magick lena.png -fft +delete lena_magnitude.png
  magick lena_magnitude.png -gravity center -extent 1x1 \
          -scale 60x60   lena_dc_zoom.gif

ご覧のとおり、画像の平均色は一種の「暗いピンク」色です。この特別なピクセルについてのもう一つの考え方は、それが中心の「バイアス」レベルを表し、その周りで他のすべてのサイン波が画像の色を変調する、というものです。 たとえば、その「暗いピンク」の DC ピクセルを、よりオレンジ色の「tomato」のような別の色で置き換えてみましょう... |

  magick lena.png -fft \
          \( -clone 0  -draw "fill tomato  color 64,64 point" \) \
          -swap 0 +delete -ift lena_dc_replace.png

実際に起きているのは、FFT 画像の DC 値を変えることで、画像全体を同じように変えているということです。実際、DC 値のあらゆる変化（その差分）は、結果として得られる画像のすべてのピクセルに加算（または減算）されます。これはまさに、元の画像の各ピクセルに何らかの定数を加算しているのと同じです。そのため、再構成された画像の最終的なピクセル色は、最大値（白）または最小値（黒）の限界によってクリッピングされることもあります。したがって、これは画像に色を付ける方法としては推奨されません。これは画像全体のすべてのピクセルを変更するよりも適用が簡単ですが、FFT の往復によって、全体として色付けの技法としてははるかに遅くなります。

サイン波画像のスペクトル

次に、画像の幅にわたって 4 サイクルを持つ単一のサイン（またはコサイン）波画像のスペクトルを見てみましょう。

  magick -size 128x129 gradient: -chop 0x1 -rotate 90 -evaluate sine 4 \
          sine4.png
  magick sine4.png -fft +delete \
          -auto-level -evaluate log 100  sine4_spectrum.png

| _上記でグラデーション画像をこのように変わった方法で作成しているのは、得られるサイン波画像が画像全体で完璧にタイル状に並ぶようにするために必要だからです。

通常の "[gradient:](canvas.html#gradient)" 画像は完璧にはタイル状に並ばないため、そこから生成されるサイン波もそうなりません。そのような不完全なタイルの FFT 変換は、フーリエ変換スペクトル内の単一の「点」ではなく、望ましくない高調波の配列を生み出します。

この問題の詳細については 完璧なグラデーションの生成 を参照してください。

_
---|---
上記のスペクトル画像（強調された振幅画像）には、3つの点があるのが分かります。中央の点は以前と同じく平均 DC 値です。他の2つの点は、フーリエ演算子が画像内で見つけた完璧なサイン波を表します。画像の幅にわたる周波数がちょうど 4 サイクルなので、結果として2つの周波数ピクセルは中央の DC 値からちょうど 4 ピクセル離れています。しかしなぜ2つのピクセルなのでしょう? それは、単一のサイン波は2つの完全に異なる方法で記述できる（一方は負の方向と位相を持つ）からです。どちらの記述も数学的に正しく、フーリエ変換はそれらを区別しません。これを 16 サイクルのサイン波で繰り返すと、やはり3つの点を持つことが分かりますが、点同士はより離れています。この場合、両側の点は中央の点の左右に 16 ピクセル離れて配置されます。

  magick -size 128x129 gradient: -chop 0x1 -rotate 90 -evaluate sine 16 \
          -write sine16.png -fft -delete 1 \
          -auto-level -evaluate log 100 sine16_spectrum.png

これから、完璧なサイン波は適切な位置の2つの点で単純に表現されることが分かります。この位置が中央の DC 値からどれだけ離れているかが、サイン波の周波数を決めます。波長が短いほど周波数が高くなるので、点は DC 値からより遠くに離れます。実際、画像のサイズを周波数（中央からの点の距離）で割ると、波の波長（山と山の間の距離）が得られます。上記の場合: 128 ピクセルを 16 サイクルで割ると、各「帯」の間の波長は 8 ピクセルになります。これはFFT 変換の最も重要な際立った特徴の一つです。元の画像の小さな特徴のパターンは小さな波長を必要とし、したがって大きな周波数を必要とします。それが周波数領域では大規模な効果をもたらします。同様に、大きな特徴はより小さな周波数を使うので、特に中央近くで小規模なパターンを生み出します。フーリエ変換では...

小さいものは大きくなり、大きいものは小さくなります。

これはフーリエ変換を扱う際に覚えておくべき最も重要な側面の一つです。これは画像の全体的なより大きな側面を保ちながら、画像からノイズ（小さな特徴）を除去する鍵だからです。 これら3つの「周波数」をより詳しく見るために、元の振幅（対数スペクトルではない）をプロットしてみましょう。 |

  magick sine16.png -fft -delete 1  miff:- |\
     im_profile - sine16_magnitude_pf.png

DC 値（画像の平均またはバイアス）が 1/2 の値を持つことに注目してください。これは予想どおりです（画像の平均値は完璧な 50% グレーです）。しかしフーリエ変換が見つけた2つの 16 サイクルのサイン波の実際の振幅は、最大値のわずか 1/4 です。元のサイン波の振幅は実際には 1/2 ですが、フーリエ変換はその振幅を2つに分割し、プロットされた両方の周波数波に結果を分け合わせるので、2つの成分のそれぞれは 1/4 の振幅しか持ちません。それがフーリエ変換の通常の一部分です。 FFT 画像におけるこの正の周波数と負の周波数の二重性が、すべての FFT 画像スペクトル（左に再掲した Lena のスペクトルなど）が常に中心に関して対称である理由を説明します。画像の片側にあるどの点についても、画像の中心を挟んで回転対称の位置に同様の「点」が必ず得られます。同じことが FFT 画像ペアの「位相」成分でも起きますが、値が 180 度シフト（負の位相）します。つまり各画像の半分は実際にはもう半分の複製ですが、元の画像を再現するには両方の画像が必要です。言い換えれば、2つの画像は依然としてまったく同じ量の情報を含んでおり、半分は一方の画像に、半分はもう一方にあります。両者が合わさって全体を作ります。 | _生成時、FFT アルゴリズムは画像の左半分だけを生成します。残りの半分は、生成されたデータの回転と複製によって生成されます。

周波数領域の画像を空間領域の画像に戻すときも、アルゴリズムは再び画像の左半分だけを見ます。右半分は複製にすぎないため、完全に無視されます。

そのため、（後の例で）FFT 振幅画像を「ノッチフィルタ」する際には、実際には振幅画像の左半分だけをフィルタすればよいのです。右半分も無視することで作業を省けます。ただし、分かりやすさのため、私は両方の半分を「ノッチ」します。

_
---|---

FFT 画像を直接生成する

上記の情報を使って、実際にサイン波の画像を生成できます。必要なのは、黒の画像と 50% グレーの画像のペアを作り、適切な振幅と位相を持つ「点」を追加することだけです。たとえば...

  magick -size 128x128  xc:black \
          -draw 'fill gray(50%)  color 64,64 point' \
          -draw 'fill gray(50%)  color 50,68 point' \
          -draw 'fill gray(25%)  color 78,60 point' \
          generated_magnitude.png
  magick generated_magnitude.png \
          -auto-level -evaluate log 3  generated_spectrum.png
  magick -size 128x128  xc:gray50  generated_phase.png
  magick generated_magnitude.png generated_phase.png \
          -ift  generated_wave.png

_
---|---
FUTURE: Perlin Noise Generator using FFT

垂直線のスペクトル

細い線と太い線の FFT スペクトルを示す 画像の FFT において、小さな特徴が「大きく」なり、大きな特徴が「小さく」なることを実演する。それを、単一の高調波を持つ「線」とみなせるサイン波に関連付ける。 線を回転させる

矩形パターン画像のスペクトル

次に、黒い背景の中に幅 8、高さ 16 の白い矩形があるスペクトルを見てみましょう。

  magick -size 8x16 xc:white -gravity center \
          -gravity center -background black -extent 128x128 rectangle.png
  magick rectangle.png -fft +delete \
          -auto-level -evaluate log 100 rect_spectrum.png

ご覧のとおり、得られた画像には非常に特徴的なパターンがあり、多くの高調波周波数を持っています。また、矩形が 90 度回転したように見えることも分かります。それは間違いで、見えているのは前に述べたのと同じ規則です... 大きな特徴は小さくなり、小さな特徴は大きくなります。そのため、矩形の小さい方の寸法が大きくなり、大きい方の寸法が小さくなったのです。では、矩形を 45 度回転させてみましょう。スペクトルも同じ方向に 45 度回転することが分かります。

  magick rectangle.png -rotate 45 -gravity center -extent 128x128 \
          -write rect_rot45.png -fft -delete 1 \
          -auto-level -evaluate log 100 rect_rot45_spectrum.png

ご覧のとおり、周波数領域でも同じ回転が起きています。つまり、何らかの回転した物体の効果は、そのフーリエ変換でも回転します。しかし、今度は矩形を移動させると...

  magick rectangle.png -rotate 45  -geometry +30+20 -extent 128x128 \
          -write rect_rot45off.png -fft -delete 1 \
          -auto-level -evaluate log 100 rect_rot45off_spectrum.png

周波数パターンは移動しませんでした。それは、すべての位置情報が位相画像に含まれているからです。周波数パターン（振幅またはそのスペクトル）は、物体が移動しても変化しません。この位置の分離は、フーリエ変換を非常に重要なものにしている主要な特徴の一つです。これにより、そのフーリエスペクトルパターンを生み出した物体の位置に関係なく、より大きな画像の中から特定の画像パターンを探すことができるようになります。

平坦な円形パターン画像のスペクトル

次に、白く平坦な円形パターンを持つ画像のスペクトルを見てみましょう。一方は直径 12（半径 6）、もう一方は直径 24（半径 12）です。

  magick -size 128x128 xc:black -fill white  \
          -draw "circle 64,64 64,70" -write circle6.png -fft -delete 1 \
          -auto-level -evaluate log 100 circle6_spectrum.png

  magick -size 128x128 xc:black -fill white  \
          -draw "circle 64,64 64,76" -write circle12.png -fft -delete 1 \
          -auto-level -evaluate log 100 circle12_spectrum.png

最初の画像は、上のほうで jinc の例として生成したものに非常に近いことに注意してください。ただし少し崩れています。これらのアーティファクトは円が小さいサイズであるために生じます。デジタルで表現されているため、その周は完璧には円形になっていません。ここでも、小さなディテールが変換後の周波数空間で大きくなることが分かります。一方、より大きな円の変換のほうが良好です。その周がより真円に近い近似になっているからです。したがって、平坦な円形の変換は確かに jinc 関数であり、より小さい直径の円を含む画像は、より広がった幅広い変換の特徴を生み出すと結論づけられます。フーリエ変換の数学的性質によれば、中心からスペクトル内の最初の暗い環の中央までの距離は 1.22N/d になります。円の直径が d=12 のとき、距離は 1.22128/12=13 になります。同様に、円の直径が d=24 のとき、距離は 1.22*128/24=6.5 になります。

ガウスパターン画像のスペクトル

次に、白いガウス円形パターンを持つ2つの画像のスペクトルを見てみましょう。シグマはそれぞれ 8 と 16 です。

  magick -size 128x128 xc:black -fill white \
          -draw "point 64,64" -gaussian-blur 0x8 -auto-level \
          -write gaus8.png -fft -delete 1 \
          -auto-level -evaluate log 1000 gaus8_spectrum.png

  im_profile -s gaus8.png gaus8_pf.gif
  im_profile -s gaus8_spectrum.png gaus8_spectrum_pf.gif

  magick -size 128x128 xc:black -fill white \
          -draw "point 64,64" -gaussian-blur 0x16 -auto-level \
          -write gaus16.png -fft -delete 1 \
          -auto-level -evaluate log 1000 gaus16_spectrum.png

  im_profile -s gaus16.png gaus16_pf.gif
  im_profile -s gaus16_spectrum.png gaus16_spectrum_pf.gif

パターンの矩形配列によって生じるノイズを除けば、ガウスパターンはほぼ同一のガウス周波数パターンを生み出すという結果になりました。さらに重要なのは、そのパターンの見た目がかなりきれいだったことです。もちろんサイズの差はあり、これもやはり同じ規則、すなわち大きいものは小さくなり、小さいものは大きくなる、に従います。数学的性質から、スペクトル内のシグマはちょうど N/(2*sigma) になります。ここで sigma は元画像のものです。したがって、サイズ N=128、sigma=8 の画像では、スペクトル内のシグマは 128/16=8 になります。同様に、画像のシグマが 16 なら、スペクトル内のシグマは 128/32=4 になります。これは「大きいものは小さくなり、その逆も成り立つ」という規則の数学的な関係であり、知っておくと役立ちます。

グリッドパターン画像のスペクトル

次に、16x8 ピクセル間隔で配置された一連のグリッド線だけを含む画像を変換してみましょう。

  magick -size 16x8 xc:white -fill black \
          -draw "line 0,0 15,0" -draw "line 0,0 0,7" \
          -write mpr:tile +delete \
          -size 128x128 tile:mpr:tile \
          -write grid16x8.png -fft -delete 1 \
          -auto-level -evaluate log 100000 grid16x8_spectrum.png

得られたスペクトルは単なる点の配列で、より密に配置されたグリッド線はより離れた点を生み出し、その逆も成り立ちます。上記の性質によれば、グリッド線は 16x8 ピクセル間隔で配置されているので、点は N/a=128/16=8 と M/b=128/8=16 の間隔で配置されるはずで、これはこの画像で測定されたとおりです。このパターンは特に重要です。なぜなら、フーリエ変換と画像内の規則的なタイリングパターンとの関係を教えてくれるからです。そのようなタイリングパターンは、そのフーリエ変換内に非常に強い非中心のグリッドパターンを生み出します。ここでの要点は、形状情報は中心にあるが、タイリング情報はそのフーリエ変換の中心から離れたグリッド状の配列にある、ということです。

スペクトルに関するさらなる情報

スペクトル画像とその性質についてもっと知りたい方のために、いくつかのリンクを挙げます。

• Wikipedia: Fourier Transform
• Fred Weinhaus, Properities of a Fourier Transform
• Wolfram MathWorld: Fourier Transform

実用例

さて、基礎を説明したところで、フーリエ変換を使った実用的な応用にはどんなものがあるでしょうか? 行えることには次のようなものがあります: 1) 画像のコントラストを上げる・下げる、2) ぼかし、3) シャープ化、4) エッジ検出、5) ノイズ除去。

画像のコントラスト変更 - 係数のべき乗

画像のコントラストは、順フーリエ変換を行い、振幅画像をべき乗し、それを位相とともに逆フーリエ変換に使うことで調整できます。コントラストを上げるには 1 よりわずかに小さい指数を使い、コントラストを下げるには 1 よりわずかに大きい指数を使います。では、まず指数 0.9 を使って Lena 画像のコントラストを上げ、次に指数 1.1 を使ってコントラストを下げてみましょう。

  magick lena.png -fft \
          \( -clone 0 -evaluate pow 0.9 \) -delete 0 \
          +swap -ift lena_plus_contrast.png

  magick lena.png -fft \
          \( -clone 0 -evaluate pow 1.1 \) -delete 0 \
          +swap -ift lena_minus_contrast.png

しかし、これを元の画像に対して行うのは、元の画像に対して同じことをするのと同じ効果になります。つまり、振幅をグローバルに変更することは、元の画像をグローバルに変更したのと同じ効果を持ちます。

画像をぼかす - ローパスフィルタ

フーリエ変換の最も重要な性質の一つは、空間領域での畳み込みが周波数領域での単純な乗算と等価であるということです。空間領域では、小さな正方形の単純な畳み込みフィルタ（カーネル）を使い、-convole オプションで画像をぼかします。これはローパスフィルタと呼ばれます。最も単純なフィルタは、均等な重みを持つ正方形の配列です。つまり、すべての値が 1 であり、畳み込みを適用する前にその合計で割って正規化されます。これは局所的すなわち近傍の平均と等価です。もう一つのローパスフィルタは、-gaussian-blur または -blur によって提供される、ガウス重み付けの円形フィルタです。周波数領域では、ローパスのぼかしフィルタの一種は、黒に囲まれた一定の強度の白い円にすぎません。これは空間領域での円形の平均化畳み込みフィルタに似ています。しかし、空間領域での畳み込みは周波数領域での乗算と等価なので、必要なのは順フーリエ変換を行い、次にフィルタを振幅画像と乗算し、最後に逆フーリエ変換を行うことだけです。小さいサイズの畳み込みフィルタは周波数領域での大きな円に対応することに注意してください。乗算は -compose の multiply 設定とともに -composite を介して実行されます。では、2 つのサイズの円形フィルタ、一方は直径 40（半径 20）、もう一方は直径 28（半径 14）でこれを試してみましょう。

  magick -size 128x128 xc:black -fill white \
          -draw "circle 64,64 44,64" circle_r20.png
  magick lena.png -fft \
       \( -clone 0 circle_r20.png -compose multiply -composite \) \
       \( +clone -evaluate log 10000 -write lena_blur_r20_spec.png +delete \) \
       -swap 0 +delete -ift lena_blur_r20.png

  magick -size 128x128 xc:black -fill white \
          -draw "circle 64,64 50,64" circle_r14.png
  magick lena.png -fft \
       \( -clone 0 circle_r14.png -compose multiply -composite \) \
       \( +clone -evaluate log 10000 -write lena_blur_r14_spec.png +delete \) \
       -swap 0 +delete -ift lena_blur_r14.png

小さい直径のフィルタを使った画像のほうが、よりぼかしが強くなることが分かります。また、得られた画像のエッジ付近に「リンギング」または「波紋」効果があることにも気づきます。これは、先ほど見たように、円のフーリエ変換が jinc 関数だからです。jinc 関数は中心から外側に進むにつれて減衰する振動を持ちます。ただしここでは、jinc 関数とその振動は、先ほど上で示したように周波数領域ではなく空間領域にあります。では、これについて何ができるでしょうか? 最も単純なのは、さまざまな 窓関数 を使って円のエッジを先細りにすることです。あるいは、定義上すでに先細りになっているガウス形状のようなフィルタを使うこともできます。では後者を行い、2 つのガウスぼかし円を使って厳しい「リンギング」効果のほとんどを取り除いてみましょう。

  magick circle_r20.png -blur 0x4 -auto-level gaussian_r20.png
  magick lena.png -fft \
       \( -clone 0 gaussian_r20.png -compose multiply -composite \) \
       \( +clone -evaluate log 10000 -write lena_gblur_r20_spec.png +delete \) \
       -swap 0 +delete -ift lena_gblur_r20.png

  magick circle_r14.png -blur 0x4 -auto-level gaussian_r14.png
  magick lena.png -fft \
       \( -clone 0 gaussian_r14.png -compose multiply -composite \) \
       \( +clone -evaluate log 10000 -write lena_gblur_r14_spec.png +delete \) \
       -swap 0 +delete -ift lena_gblur_r14.png

これはもちろんずっと良くなっています。理想的なローパスフィルタは円をまったくぼかすのではなく、実際には 半径 ではなく シグマ の適切なガウス曲線を使うことです。もちろんこの例では、ぼかしを行うために結局ぼかしを行うことになりました! しかし FFT 振幅画像に対して乗算するぼかしパターンは固定されており、実際には事前に生成したキャッシュから取得できます。また、乗算する画像は元の画像の全サイズである必要はなく、より小さい画像を使えます。そのため、上記の方法は大きな画像や多数の画像を扱う場合にはずっと高速になりえます。より重要な点は、強い大きなぼかしの場合、周波数領域の画像は小さく、元の画像の各ピクセルごとに多数のピクセルを平均化するのではなく、単一の乗算を行うだけだということです。小さいサイズのぼかしには、より直接的な畳み込みぼかしのほうが良いかもしれません。

画像のエッジ検出 - ハイパスフィルタ

空間領域では、画像からエッジを抽出するハイパスフィルタは、しばしば合計がゼロになるような正と負の重みを持つ畳み込みとして実装されます。周波数領域では物事ははるかに単純です。ここでハイパスフィルタは、ローパスフィルタを単に反転させたものにすぎません。つまり、ローパスフィルタが明るいところでハイパスフィルタは暗く、その逆も成り立ちます。そのため ImageMagick では、必要なのはローパスフィルタ画像を -negate するだけです。では、円画像を使って Lena 画像にハイパスフィルタを適用してみましょう。そしてもう一度、純粋なガウス曲線を使ってみます。

  magick circle_r14.png -negate circle_r14i.png
  magick lena.png -fft \
       \( -clone 0 circle_r14i.png -compose multiply -composite \) \
       \( +clone -evaluate log 10000 -write lena_edge_r14_spec.png +delete \) \
       -delete 0 +swap -ift -normalize lena_edge_r14.png

  magick -size 128x128 xc: -draw "point 64,64" -blur 0x14 \
          -auto-level   gaussian_s14i.png
  magick lena.png -fft \
       \( -clone 0 gaussian_s14i.png -compose multiply -composite \) \
       \( +clone -evaluate log 10000 -write lena_edge_s14_spec.png +delete \) \
       -delete 0 +swap -ift -normalize lena_edge_s14.png

これら 2 つの結果を注意深く調べると、単純な円はガウスほど良くないことが分かります。「リンギング」アーティファクトがあり、それほど鮮明でもないからです。

画像のシャープ化 - ハイブーストフィルタ

画像をシャープ化する最も単純な方法は、（正規化のストレッチなしで）ハイパスフィルタをかけ、それを元の画像とブレンドすることです。

  magick lena.png -fft \
       \( -size 128x128 xc: -draw "point 64,64" -blur 0x14 -auto-level \
          -clone 0 -compose multiply -composite \) \
       -delete 0 +swap -ift \
       lena.png -compose blend -set option:compose:args 100x100 -composite \
       lena_sharp14.png

ここでは、ハイパスフィルタが周波数領域で行われ、その結果が空間領域に変換し戻され、そこで元の画像とブレンドされて画像のエッジが強調されます。

ノイズ除去 - ノッチフィルタ

多くのノイズの多い画像には、何らかのパターン化されたノイズが含まれています。この種のノイズは、パターンがいくつかの点または線のパターンとして現れるため、周波数領域で簡単に除去できます。単純なサイン波は繰り返しパターンであり、スペクトル内ではわずか 3 つの点として現れることを思い出してください。このノイズを除去するには、残念ながら、振幅画像内の点または線を手作業でマスク（またはノッチ）するしかありません。これは、周波数領域に変換し、スペクトルのグレースケール版を作り、点または線をマスクし、しきい値処理し、二値マスク画像を振幅画像と乗算し、それから空間領域に変換し戻すことで行います。これを、斜めの縞模様のディザのようなパターンを含む クラウン画像 で試してみましょう。まずクラウン画像を変換して、その振幅画像と位相画像を作ります。

  magick clown_orig.jpg -fft \
          \( +clone  -write clown_phase.png +delete \) +delete \
          -write clown_magnitude.png  -colorspace gray \
          -auto-level -evaluate log 100000  clown_spectrum.png

元画像 | |
スペクトル |
位相
---|---|---|---

スペクトルには 4 つの明るい星のような点が、各象限に 1 つずつ含まれていることが分かります。これらの異常な点は、取り除きたい画像内のパターンを表します。画像の中央にある明るい点と線は気にする必要はありません。それらは DC（平均画像色）と画像のエッジの効果を表すもので、変更すべきではないからです。スペクトル画像を生成する際、得られる画像を純粋なグレースケール画像にするよう強制したことに注意してください。これにより、画像をエディタに読み込み、（赤などの）任意の非グレー色を使って、それら 4 つの星のようなパターンの領域をマスクできます。編集を終えたら、未編集版との差分画像を抽出することで、色を付けた領域を抽出できます。このように...

  magick clown_spectrum_edited.png clown_spectrum.png \
          -compose difference -composite \
          -threshold 0 -negate clown_spectrum_mask.png

あとは単純にマスクを振幅と乗算し、その結果を元の位相画像とともに使って空間領域に変換し戻します。比較のため、その横に元の画像を表示します。

  magick clown_magnitude.png clown_spectrum_mask.png \
          -compose multiply -composite \
          clown_phase.png -ift clown_filtered.png

非常に良い結果です。しかし、さらに良くすることができます。前の例で見たように、単純な「円」は FFT 画像にとって特に好ましくないので、マスクを少しぼかしてみましょう... |

  magick clown_spectrum_mask.png \
          -blur 0x5 -level 50x100%  clown_mask_blurred.png

そしてクラウンをフィルタします。今回は FFT 画像をメモリ内で再生成します。 |

  magick clown_orig.jpg -fft \
          \( -clone 0 clown_mask_blurred.png -compose multiply -composite \) \
          -swap 0 +delete -ift clown_filtered_2.png

実に見事な結果です! そして、そのマスクを「星」の形によりよく合うように調整すれば、さらに改善できる可能性があります。 元画像と結果の差分を取って、ノイズが除去された領域の画像を作ることもできます。 |

  magick clown_orig.jpg clown_filtered_2.png -compose difference \
          -composite -normalize clown_noise.png

これを別の例で試してみましょう。今回は RoboRealm のウェブサイトで見つけた「Twigs」画像で、水平・垂直の縞模様の不規則なパターンを含んでいます。やはり、前と同じようにグレースケールのスペクトル画像を抽出します。

  magick twigs.jpg -fft +delete -colorspace gray \
          -auto-level -evaluate log 100000 twigs_spectrum.png

この場合、画像内のノイズが水平・垂直方向に向いているため、これは中心線に沿った太い水平・垂直の帯として現れますが、画像の実際の中心には現れません。やはり画像エディタを使って各部分をマスクします。今回は「青」色を使います（実際にはどの色を使っても問題ありません）...

  magick twigs_spectrum_edited.png twigs_spectrum.png \
          -compose difference -composite \
          -threshold 0 -negate twigs_spectrum_mask.png

ここで再びマスクを FFT 振幅画像と乗算し、画像を再構成します。

  magick twigs.jpg -fft \
          \( -clone 0 twigs_spectrum_mask.png -compose multiply -composite \) \
          -swap 0 +delete  -ift twigs_filtered.png

そして元画像と結果の差分を取って、ノイズが除去された領域の画像を作ることができます。

  magick twigs.jpg twigs_filtered.png -compose difference -composite \
          -normalize twigs_noise.png

マスクに少しぼかしを加えれば、結果をさらに改善できるでしょう。練習として、画像から紐を取り除いてみてください。ヒントとして、実際の画像内の線の効果が FFT では 90 度回転することを思い出してください。これを間違えると、おそらく代わりに小枝を取り除いてしまうでしょう。

発展的な応用

フーリエ変換を使ったその他のより発展的な応用には、次のようなものがあります: 1) モーションブラーやピントぼけの画像のデコンボリューション（デブラー）、2) より大きな画像の中で小さな画像が最もよく一致する場所を見つけるための正規化相互相関。FFT 乗算と除算（デコンボリューション）の例は、より正式に定義された画像処理演算子を待っている状態のため、サブディレクトリ に移されました。