ImageMagick 例 --
フーリエ乗算/除算

索引

ImageMagick 例前書きと索引

フーリエ変換

DIY FFT 数学入門

FFT 乗算

振幅/位相画像、IM Q16 非HDRI を使用
実数/虚数画像、IM Q16 HDRI を使用

FFT 除算

振幅/位相画像、IM Q16 非HDRI を使用
実数/虚数画像、IM Q16 HDRI を使用

DIY FFT 数学入門

以下は、2種類の高速フーリエ変換画像による乗算と除算を行うためのいくつかの技法の例です。

数学的手法には基本的に2つの方法があります。FX演算と合成演算です。

FX 演算: FX、DIY 演算子を使用すると、式を直接適用できます。ただし、解釈されるため非常に遅くなります。しかし、式はすべて単一の式で適用されるため、中間画像が不要になります。IMの非HDRIバージョンでは、量子化丸め誤差が減少します。
合成演算: これは画像合成を使用して算術演算を実行するもので、「FX 演算」よりもはるかに高速です。ただし、一度に1つの算術演算しか実行できないため、各演算後に画像を別のメモリ内画像に保存する必要があります。非HDRIバージョンでは、演算間で追加の「丸め誤差」が発生します。

上記で触れられているもう一つの側面は、HDRIの使用です。

ImageMagick の HDRI バージョンは、量子範囲でスケールされた整数（品質を参照）ではなく、浮動小数点数としてカラー値を保存します。そのため、FFT画像の生成時や中間結果を新しい画像に保存する際に、「丸め誤差」が発生しません。特に小さな値（FFT画像では一般的）が関与する場合に有効です。

また、実数/虚数FFT画像のペアには負の数を使用する必要があるため、これらのタイプのFFT画像を処理する際には、ImageMagickのHDRIバージョンを使用する必要があります。

ただし、振幅/位相FFT画像では負の数は問題にならないため、HDRIを使用せずにそのような画像を処理できますが、それでも各画像処理ステップごとに強い「丸め誤差」が発生します。

このため、使用しているFFT画像の種類に関係なく、高速フーリエ変換はHDRIバージョンのImageMagickで使用するのが最適です。HDRIと非HDRIのImageMagickバージョンで、まったく同じコマンドを使用します。

中間HDRI画像をディスクファイルに保存する場合は、非常に少ない数の浮動小数点画像ファイル形式のいずれかを使用する必要があります。これらの形式には、NetPBM PFMファイル形式、TIFF、または特別な設定「-define quantum:format=floating-point」を伴うMIFFファイル形式が含まれます。

使用する畳み込みカーネル画像を以下に示します。


  # motion blur kernel
  #magick -size 128x128 xc:black -fill white -draw 'line 60,64 68,64' \
  #        -alpha off convolve_kernel.png

  # disk blur kernel (Fred Wienhaus's version)
  magick -size 128x128 xc:black -fill white -draw 'circle 64,64 54,64' \
          -alpha off convolve_kernel.png

畳み込みカーネルは、ピクセル0,0が畳み込みの中心を含むようにロールする必要があります。この場合、カーネル画像の読み込み直後に「-roll -64-64」を適用する必要があります。

カーネル画像は、FFT乗算または除算の対象となる画像の強度を維持するために、平均値で除算する必要もあります。これは、FFT乗算と除算における大きな問題です。

非HDRI画像はこのように正規化できないため、これは一般的にFFT変換が適用された後、乗算の前に実行されます。

IMのHDRIバージョンでは、いつでもこれを実行できます。

FFT 乗算( )

振幅/位相画像のFFT乗算、IM Q16を使用

    R = A ⊗ B       ( FFT Multiply )

    Rm = Am × Bm
    Rp = mod( Ap + Bp +0.5, 1.0)

    mean = the DC value (center pixel) of the magnitude image

FX演算を使用
(注: 値 v.p{64x64} は、畳み込みカーネルのDC、つまり平均です)


  # non-HDRI
  magick convolve_kernel.png -roll -64-64 -fft \
          \( cameraman_sm.png -fft \) \
          \
          \( -clone 0,2 -fx 'u*v / p{64,64} ' \) \
          \( -clone 1,3 -fx 'mod(u + v + 0.5, 1.0)' \) \
          \
          -delete 0-3 -ift  cameraman_convolve_1.png

合成演算を使用


  # non-HDRI
  magick convolve_kernel.png -roll -64-64 -fft \
          \( -clone 0 -crop 1x1+64+64 +repage -scale 128x128 \
             -clone 0 -compose divide -composite \) -swap 0 +delete \
          \( cameraman_sm.png -fft \) \
          \
          \( -clone 0,2 -compose multiply -composite \) \
          \( -clone 1,3 -compose add -background gray50 -flatten \) \
           -delete 0-3 -ift  cameraman_convolve_2.png

「add」はモジュロ加算合成であり、50％バイアス値を加算して50％バイアス値を生成します。

代替案
元の画像の平均で除算


  # non-HDRI
  magick convolve_kernel.png \( -clone 0 -roll -64-64 -fft \) \
          \( -clone 0 -scale 1x1 -scale 128x128 \
             -clone 1 -compose divide -composite \) \
          -delete 0,1 +swap \
          \( cameraman_sm.png -fft \) \
          \
          \( -clone 0,2 -compose multiply -composite \) \
          \( -clone 1,3 -compose add -background gray50 -flatten \) \
           -delete 0-3 -ift  cameraman_convolve_2b.png

代替案
DC値/DC値 == 1.0に注意してください。

DC値は振幅スペクトルの中で最大の値であるため、正規化しないのはなぜですか！例：-auto-level"

これは振幅画像に対してのみ実行でき、実数/虚数ペアでは同じ量だけ引き伸ばす必要があるため機能しません。


  # non-HDRI
  magick convolve_kernel.png -roll -64-64 -fft \
          \( -clone 0 -auto-level \) -swap 0 +delete \
          \( cameraman_sm.png -fft \) \
          \
          \( -clone 0,2 -compose multiply -composite \) \
          \( -clone 1,3 -compose add -background gray50 -flatten \) \
          \
          -delete 0-3 -ift  cameraman_convolve_2c.png

実数/虚数画像のFFT乗算、IM HDRIを使用

    R = A ⊗ B       ( FFT Multiply )

    Rr = Ar×Br - Ai×Bi
    Ri = Ar×Bi + Ai×Br

    mean = the DC value (center pixel) of the real image

FXを使用
5つの画像が関与し、5番目の画像がスケールされた平均です - FXが除算を実行します


  # HDRI
  magick convolve_kernel.png -roll -64-64 \( +clone +fft \) \
          \( cameraman_sm.png +fft \) +matte \
          \( -clone 0 -scale 1x1 \) -delete 0 \
          \
          \( -clone 0-4 -fx '( u[0]*u[2] - u[1]*u[3] ) / u[4].p{0,0}' \) \
          \( -clone 0-4 -fx '( u[0]*u[3] + u[1]*u[2] ) / u[4].p{0,0}' \) \
          -delete 0-4 +ift  cameraman_convolve_3.png

FXを使用
平均値にDC値を使用（5番目の画像はなし）


  # HDRI
  magick convolve_kernel.png -roll -64-64 +fft \
          \( cameraman_sm.png +fft \) \
          \
          \( -clone 0-3 -fx '( u[0]*u[2] - u[1]*u[3] ) / u[0].p{64,64}' \) \
          \( -clone 0-3 -fx '( u[0]*u[3] + u[1]*u[2] ) / u[0].p{64,64}' \) \
          -delete 0-3 +ift  cameraman_convolve_3b.png

合成演算を使用
スケール画像平均を使用！


  # HDRI
  magick convolve_kernel.png -roll -64-64 \
          \( -clone 0 -scale 1x1 -scale 128x128 \) \
          \( -clone 0 +fft \) \
          \( -clone 1,2 -compose divide -composite \) \
          \( -clone 1,3 -compose divide -composite \) \
          -delete 0--3 \
          \
          \( cameraman_sm.png  +fft \) \
          \
          \( -clone 0,2 -compose multiply -composite \) \
          \( -clone 1,3 -compose multiply -composite \) \
          \( -clone 0,3 -compose multiply -composite \) \
          \( -clone 1,2 -compose multiply -composite \) \
          \
          \( -clone 4,5 +swap +matte -compose minus -composite \) \
          \( -clone 6,7 -compose plus -composite \) \
          \
          -delete 0--3 +ift  cameraman_convolve_4.png

最適化された合成演算
画像スケーリングからの平均


  # HDRI
  magick convolve_kernel.png -roll -64-64 \
          \( -clone 0 -scale 1x1 -scale 128x128 \) \
          \( -clone 0 +fft    null: +insert    -clone 1 +insert \
             -compose divide -layers composite \) -delete 0,1  \
          \
          \( cameraman_sm.png +fft \) \
          \
          \( -clone 0,1,0,1 null: -clone 2,3,3,2 \
                -compose multiply -layers composite \) \
          \( -clone 4,5 +swap +matte -compose minus -composite \) \
          \( -clone 6,7 -compose plus -composite \) \
          \
          -delete 0--3 +ift  cameraman_convolve_4b.png

代替案
平均値にDC値を使用


  # HDRI
  magick convolve_kernel.png -roll -64-64 +fft \
          \( -clone 0 -crop 1x1+64+64 +repage -scale 128x128 \) \
          null: +insert +insert   -compose divide -layers composite \
          \
          \( cameraman_sm.png +fft \) \
          \
          \( -clone 0,1,0,1 null: -clone 2,3,3,2 \
                -compose multiply -layers composite \) \
          \( -clone 4,5 +swap +matte -compose minus -composite \) \
          \( -clone 6,7 -compose plus -composite \) \
          \
          -delete 0--3 +ift  cameraman_convolve_4c.png

代替案
FFT乗算の後、DC値平均を適用


  # HDRI
  magick convolve_kernel.png -roll -64-64 +fft \
          \( cameraman_sm.png  +fft \) \
          \
          \( -clone 0,1,0,1 null: -clone 2,3,3,2 \
                -compose multiply -layers composite \) \
          \( -clone 4,5 +swap -compose minus -composite \) \
          \( -clone 6,7 -compose plus -composite \) \
          \
          \( -clone 0 -crop 1x1+64+64 +repage -scale 128x128 \
             null: -clone -2,-1  -compose divide -layers composite \) \
          \
          -delete 0--3 +ift  cameraman_convolve_4d.png

FFT 除算( )

ここでは、除算を使用して、上記の画像に追加されたぼかしを除去またはデコンボリューションします（より良い言葉が見つからないため）。基本的にまったく同じですが、「正規化された」畳み込みカーネルがメイン画像から除算されます。

各例では、上記と同じ手法を使用して生成された画像を使用します。

振幅/位相のFFT除算、IM Q16を使用

    R = B ø A       ( FFT Divide )

    Rm = Bm / Am
    Rp = mod( -Ap + Bp +1.5, 1.0)

    mean = the DC value (center pixel) of the magnitude image

注：ゼロによる完全な除算を回避するために、分母に量子化スケールの最小値が追加されます。これは後で画像からノイズを除去するためにも使用されます。

FX演算を使用


  # non-HDRI
  noise=QuantumScale
  magick convolve_kernel.png -roll -64-64 -fft \
          \( cameraman_convolve_1.png -fft \) \
          \
          \( -clone 0,2 -fx "v/(u/p{64,64}+$noise)" \) \
          \( -clone 1,3 -fx 'mod( -u + v + 1.5, 1.0)' \) \
          \
          -delete 0--3 -ift cameraman_deconvolve_1.png


  echo -n "Peak = "
  magick compare -metric PAE cameraman_sm.png cameraman_deconvolve_1.png null:
  echo -n "Avg = "
  magick compare -metric RMSE cameraman_sm.png cameraman_deconvolve_1.png null:

IMのHDRIバージョンを使用して、ラウンドトリップMP畳み込みFFT（FX演算）を実行

合成演算を使用


  # non-HDRI
  noise=1
  magick convolve_kernel.png -roll -64-64 -fft \
           \( -clone 0 -crop 1x1+64+64 +repage -scale 128x128 \
              -clone 0 -compose divide -composite \) -swap 0 +delete \
          \( cameraman_convolve_2.png -fft \) \
          \
          \( -clone  0  -evaluate add $noise \
             -clone  2  -compose divide -composite \) \
          \( -clone 1,3 -compose subtract -background gray50 -flatten \) \
          \
          -delete 0--3 -ift cameraman_deconvolve_2.png


  echo -n "Peak = "
  magick compare -metric PAE cameraman_sm.png cameraman_deconvolve_2.png null:
  echo -n "Avg = "
  magick compare -metric RMSE cameraman_sm.png cameraman_deconvolve_2.png null:

IMのHDRIバージョンを使用して、ラウンドトリップMP畳み込みFFT（合成演算）を実行

実数/虚数のFFT除算、IMのHDRIバージョンを使用

    R = B ø A       ( FFT Divide )

    Denom = Ar×Ar + Ai×Ai + noise
    Rr = ( Ar×Br + Ai×Bi ) / Denom
    Ri = ( Ar×Bi - Ai×Br ) / Denom

    mean = the DC value (center pixel) of the real image

FXを使用


  # HDRI
  #noise=QuantumScale
  noise=Epsilon
  magick convolve_kernel.png -roll -64-64 +fft \
          \( -clone 0 -fx "u/p{64,64}" \) \( -clone 0,1 -fx "v/u.p{64,64}" \) \
          -delete 0,1 \
          \( cameraman_convolve_3.png +fft \) \
          \
          \( -clone 0-3 -fx "u[0]*u[0] + u[1]*u[1] + $noise" \) \
          \( -clone 0-3 -fx 'u[0]*u[2] + u[1]*u[3]' \) \
          \( -clone 0-3 -fx 'u[0]*u[3] - u[1]*u[2]' \) \
          \( -clone 4,5 -fx 'u==0?0:v/u' \) \
          \( -clone 4,6 -fx 'u==0?0:v/u' \) \
          \
          -delete 0--3 +ift cameraman_deconvolve_3.png


  echo -n "Peak = "
  magick compare -metric PAE cameraman_sm.png cameraman_deconvolve_3.png null:
  echo -n "Avg = "
  magick compare -metric RMSE cameraman_sm.png cameraman_deconvolve_3.png null:

上記で使用した-fx除算は、合成除算ほど優れていません。おそらく、関与する追加のFXスケーリングによるものです。


  # HDRI
  #noise=QuantumScale
  noise=Epsilon
  magick convolve_kernel.png -roll -64-64 +fft \
          \( -clone 0 -fx "u/p{64,64}" \) \( -clone 0,1 -fx "v/u.p{64,64}" \) \
          -delete 0,1 \
          \( cameraman_convolve_3.png +fft \) \
          \
          \( -clone 0-3 -fx "u[0]*u[0] + u[1]*u[1] + $noise" \) \
          \( -clone 0-3 -fx 'u[0]*u[2] + u[1]*u[3]' \) \
          \( -clone 0-3 -fx 'u[0]*u[3] - u[1]*u[2]' \) \
          \( -clone 4,5 -compose divide -composite \) \
          \( -clone 4,6 -compose divide -composite \) \
          \
          -delete 0--3 +ift cameraman_deconvolve_3b.png


  echo -n "Peak = "
  magick compare -metric PAE cameraman_sm.png cameraman_deconvolve_3b.png null:
  echo -n "Avg = "
  magick compare -metric RMSE cameraman_sm.png cameraman_deconvolve_3b.png null:

合成演算を使用（詳細）


  # HDRI
  noise=0.00000001  # Epsilon
  magick \( convolve_kernel.png -roll -64-64 +fft \
             \( -clone 0 -crop 1x1+64+64 +repage -scale 128x128 \) \
             \( -clone 2,0 -compose divide -composite \) \
             \( -clone 2,1 -compose divide -composite \) \
             -delete 0-2 \) \
          \( cameraman_convolve_4.png +fft \) \
          \
          \( -clone 0,0 -compose multiply -composite \) \
          \( -clone 1,1 -compose multiply -composite \) \
          \( -clone 0,2 -compose multiply -composite \) \
          \( -clone 1,3 -compose multiply -composite \) \
          \( -clone 0,3 -compose multiply -composite \) \
          \( -clone 1,2 -compose multiply -composite \) \
          \
          \( -clone 4,5 -compose plus  -composite -evaluate add $noise \) \
          \( -clone 6,7 -compose plus -composite \) \
          \( -clone 8,9 +swap -compose minus -composite \) \
          \
          \( -clone 10,11 -compose divide -composite \) \
          \( -clone 10,12 -compose divide -composite \) \
          \
          -delete 0--3 +ift cameraman_deconvolve_4.png


  echo -n "Peak = "
  magick compare -metric PAE cameraman_sm.png cameraman_deconvolve_4.png null:
  echo -n "Avg = "
  magick compare -metric RMSE cameraman_sm.png cameraman_deconvolve_4.png null:

最適化された合成演算


  # HDRI
  noise=0.00000001  # Epsilon
  magick \( convolve_kernel.png -roll -64-64 +fft \
             \( -clone 0 -gravity center -extent 1x1 -scale 128x128 \) \
             null: +insert +insert -compose divide -layers composite \) \
          \( cameraman_convolve_4.png +fft \) \
          \
          \( -clone 0,1,0,1,0,1 null: -clone 0,1,2,3,3,2 \
                -compose multiply -layers composite \) \
          \( -clone 4,5 -compose plus  -composite -evaluate add $noise \) \
          \( -clone 6,7 -compose plus  -composite \) \
          \( -clone 9,8 -compose minus -composite \) \
          \( -clone 10 null: -clone 11,12 \
                -compose divide -layers composite \) \
          \
          -delete 0--3 +ift cameraman_deconvolve_4b.png


  echo -n "Peak = "
  magick compare -metric PAE cameraman_sm.png cameraman_deconvolve_4b.png null:
  echo -n "Avg = "
  magick compare -metric RMSE cameraman_sm.png cameraman_deconvolve_4b.png null:

更新：上記の色「gray50」の使用は、「gray(50%)」に変更して、より正確な50％グレイ値を生成する必要があります。「名前付き」の色は実際には8ビットカラーのみであり、そのためあまり正確ではありません。
また、新しいカラー空間処理では、gray50はsRGBカラー空間にあるのに対し、gray(50%)は線形カラー空間にあるため、上記の計算で使用すべきです。
例は、線形グレースケールカラー空間の使用についても確認する必要があります。

ImageMagick 例 -- フーリエ乗算/除算

DIY FFT 数学入門

FFT 乗算( )

振幅/位相画像のFFT乗算、IM Q16を使用

実数/虚数画像のFFT乗算、IM HDRIを使用

FFT 除算( )

振幅/位相のFFT除算、IM Q16を使用

実数/虚数のFFT除算、IMのHDRIバージョンを使用

ImageMagick 例 --
フーリエ乗算/除算