Rest Term

JavaScriptで2D-FFTによるハイパス/ローパスフィルタ

January 4th, 2012 in tech/study  tags: ,  
Skip to comments (0) ↓
 
このエントリーをはてなブックマークに追加

これまでにJavaScriptとHTML5 Canvas APIでいくつかの画像処理を試してきましたが、今回は二次元離散フーリエ変換(2D-DFT、実装上では2D-FFT)で得られた周波数スペクトルにハイパス/ローパスフィルタ(HPF/LPF)を適用します。

「フーリエ変換」は音声処理でよく耳にする単語かと思います。音声データをフーリエ変換してHTML5 Canvas上でビジュアライズするデモもたくさん公開されています(例: Visualizing an audio spectrum - MDN)。ただし、ここでは"音声"ではなく"画像"に対するフーリエ変換を行います。音声の場合はデータが一次元なので周波数は一つしか持っていませんが、画像の場合は水平/垂直方向の二つの周波数を持つことになります。なので画像データに対する2D-FFTは、

 x軸方向に一次元フーリエ変換 → y軸方向に一次元フーリエ変換

という手順で処理します。

このブログでも2008年にFlashで2D-FFTの検証をしていたことがあるので参考までに。

 * 二次元離散フーリエ変換 – AS3.0
 * ハイパス/ローパスフィルタ
 * 位相画像
 * 位相限定相関法

ブラウザがCanvas APIをサポートするようになって割と早い段階でJavaScriptでも実装していたのですが、当時のブラウザの処理能力だと案の定まともに動作しませんでした。それが今ではどのブラウザも飛躍的に処理能力が向上したので、2D-FFTをWeb Workersを使わずに動かしても大丈夫なようになりました。

FFTおよび空間周波数フィルタのプログラムはGitHubに置きました。全て合わせて8KBちょっとです。外部ライブラリは用いておらず、Canvas APIを直接利用したPure JavaScriptなコードになっています。AS3版ではFFTと逆FFT、ハイパス/ローパスフィルタをそれぞれ一つの関数内で処理していたことで条件分岐が増えて遅くなっていましたが、JS版ではコードが多少冗長になるのを許容し、APIを分けることでコアな部分の処理内で条件分岐を削減しています。

cv/fft at master from wellflat/jslib - GitHub

まずは一次元離散フーリエ変換からテストしてみます。

* 1D-FFT

var x1 = [], y1 = [], x2 = [], y2 = [], x3 = [], y3 = [], N = 32;
for(var i=0; i<N; i++){
  x1[i] = x2[i] = 6*Math.cos(6*Math.PI*i/N) + 4*Math.sin(18*Math.PI*i/N);
  y1[i] = y2[i] = 0.0;
}
FFT.init(N);  // 初期化(三角関数テーブル、ビット逆順テーブルの生成)
FFT.fft1d(x2, y2);  // 一次元FFT
for(var j=0; j<N; j++) {
  x3[j] = x2[j];
  y3[j] = y2[j];
}
FFT.ifft1d(x3, y3);  // 一次元逆FFT
var out = "";
for(var i=0; i<N; i++){
  out += i + ": (" + x1[i] + ", " +y1[i]+ ") (" +
         x2[i] + ", " + y2[i] + ") (" + x3[i] + ", " + y3[i] + ")\n";
}
console.log("N: (元データ) (FFT) (逆FFT)");
console.log(out);

N: (元データ) (FFT) (逆FFT)
0: (6, 0) (1.3100631690576847e-14, 0) (6, 0)
1: (8.911958795428193, 0) (-3.148595663620103e-15, -1.634119351805764e-15) (8.911958795428195, 2.220446049250315e-16)
2: (0.7653668647301803, 0) (-3.582018523007779e-14, 1.0273603512156518e-14) (0.76536686473018, 5.551115123125783e-17)
3: (-4.496420381306951, 0) (96.00000000000003, -4.440892098500626e-14) (-4.496420381306951, -4.440892098500624e-16)
... 省略

浮動小数点数の扱いがASと同じなので、JSでも同じ計算誤差が発生します(上の例でだいたい1.0e-15程度)。

次は本題の二次元離散フーリエ変換による空間周波数フィルタ処理を試してみます。

* 2D-FFTによる空間周波数フィルタ(ハイパス/ローパスフィルタ)

var spectrum = document.querySelector('#Spectrum').getContext('2d'),
    result = document.querySelector('#Result').getContext('2d'),
    image = new Image();
image.src = '/path/to/image';
image.addEventListener('load', function(e) {
  var w = image.width,
      h = image.height, // w == h
      re = [],
      im = [];

  // 初期化、次数Nは基数2
  FFT.init(w);
  FrequencyFilter.init(w);
  // 初期化、表示用の CanvasRenderingContext2D を渡す
  SpectrumViewer.init(spectrum);
  spectrum.drawImage(image, 0, 0);
  var src = spectrum.getImageData(0, 0, w, h),
      data = src.data,
      radius = 30,
      i, val, p;
  // 実数部、虚数部を格納する配列を生成
  // 画像データは実数部配列に詰める
  for(var y=0; y<h; y++) {
    i = y*w;
    for(var x=0; x<w; x++) {
      re[i + x] = data[(i << 2) + (x << 2)];
      im[i + x] = 0.0;
    }
  }
  FFT.fft2d(re, im);  // 二次元FFT
  FrequencyFilter.swap(re, im)  // 象限入れ替え
  FrequencyFilter.HPF(re, im, radius); // ハイパスフィルタ
  //FrequencyFilter.LPF(re, im, radius);  // ローパスフィルタ
  //FrequencyFilter.BPF(re, im, radius, radius/2);  // バンドパスフィルタ
  SpectrumViewer.render(re, im);  // スペクトル画像を描画
  FrequencyFilter.swap(re, im);  // 象限入れ替え
  FFT.ifft2d(re, im);  // 二次元逆FFT
  // 実数部配列のデータを表示用の CanvasPixelArray に詰める
  for(var y=0; y<h; y++) {
    i = y*w;
    for(var x=0; x<w; x++) {
      val = re[i + x];
      p = (i << 2) + (x << 2);
      data[p] = data[p + 1] = data[p + 2] = val;
    }
  }
  result.putImageData(src, 0, 0);
}, false);

Demo: Spatial Frequency Filtering by 2D-FFT
(IE9.0, Firefox8.0, Chrome16.0, Safari5.0, Opera11.60で動作確認済み)

このデモでは周波数スペクトルを見やすいように緑色に着色しています。人間の目は「低周波成分には敏感だが、高周波成分には鈍感である」という性質があるので、画像の高周波成分を間引いても元画像との違いを判別することは難しいです。JPEG圧縮アルゴリズムはこの性質を利用してデータ量を減らしています。デモでローパスフィルタを使ってぜひ実感してみてください。

今回書いたモジュールは基数2以外の次数だと例外を投げるようにしています。基数4と余裕があれば8のやつも実装しようかと思ったのですが、式も図も書かずにバタフライ演算の部分を実装してたら混乱してきたので止めておきました。。C言語ならともかくJavaScriptのレベルだとCPUのキャッシュライン等を考えて実装したりはしないので、落とし所としてはこのあたりが妥協点かなと思います。とりあえずはこれで位相限定相関法など多くのアルゴリズムが記述できるので、適切にWeb Workersを併用することでより高度な画像解析も行うことができそうです。

クロスブラウザ対応について

FFTの処理自体はECMAScriptベースの言語ならどれも同じ結果が得られるはずですが、Canvas APIの実装は各ブラウザで異なっているのでこちらはケアしなければなりません。Firefox, Chrome, Safari, IE9では特に問題なく動作してコードの改修も必要なかったのですが、Operaの場合はサチュレーション(飽和演算)を手抜きしてると動作しないようです(2012/01現在)。Opera以外のブラウザだと内部でクリッピングしてくれてるのだと思います。

min(max(value,0),255) の処理を忘れずに。OpenCVでいうところの saturate_cast<typename> ですね。

JavaScriptエンジンの進化

ちょうどFirefoxが怒濤のラピッドリリースに突入した1年くらい前から(Firefox3.6の頃)、ブラウザのJavaScript処理能力が急速に上がっていったようで、もちろんFirefoxに限らずChromeやOperaも大きく進化しています。悪評高いInternet Explorerもversion9になってからはいろいろな面でまともに動作するようになっているようです。

これも昔Flashで作ったサンプルのCanvasへの移植なのですが、クラインの壺をCanvas上で回してみました。お絵かき系の処理能力はまだまだといった印象ですが、こちらも急速に進化していくのではないかと思います。

Demo: Klein Bottle

やっぱり激化するブラウザ競争の中で生き残るためにどのベンダーも必死にがんばってるんでしょうね。こういうオープンな「競争」はとても刺激があって良いことだと思います。僕は普段仕事でJavaScriptを使うことがないので、もうちょっとフロントエンド寄りのプログラマやデザイナーと交流していろいろ学んでいきたいです。

合わせて読む:

 

Tags: ,

Comments

No comments so far.

  • Leave a Reply
     
    Your gravatar
    Your Name