Re - ImageJで学ぶ!: 2015-08-09

2015年8月11日火曜日

第8回 Imagejを用いて画像の周波数表現を理解する(その3)で学ぶ!

第6回から「周波数」に着目して記事を書いてきました。その中でフーリエ変換という用語を使っていましたが、今回はこのフーリエ変換と関係する画像のMTFについて解説を行います。
MTFとは、modulation transfer functionの略語で、日本語では「変調伝達関数」と称されています。このMTFが医用画像解析の中でどのように使われるのか?どのように測定するのか?また、それが臨床上、何のために必要なのかを入門レベルとして記述します。

MTF


医用画像分野でもMTFは変調伝達関数と呼ばれる鮮鋭度を測る定量指標です。鮮鋭度とは、画像に写っているものがどれだけはっきりくっきりみえるかという度合いを示すものです。はっきりくっきり見えるものはMTFが高く、ボケて見えるものは低くなります。医療では、眼科学などの光学技術を利用する分野での研究(レンズの集光性など)の歴史が長いです。レントゲン写真、CTやMRIなどの医用画像では、画像の空間分解能の比較(装置間の比較、撮影(像)条件間の比較など)のために1960年代ごろから利用されています。

MTFを測る


MTFの測定は非常に簡単です。よく利用される計測方法として次のようなものがあります。

①矩形もしくは正弦波チャートを用いて、既知の周波数ごとに振幅をプロットする方法
②点像の強度分布(point spread function:PSF)をフーリエ変換する方法
③線像の強度分布(Line Spread Function:LSF)をフーリエ変換する方法
④画像のエッジの強度分布をフーリエ変換する方法

これらの方法は出力を画像処理アプリケーションから定量的に得ますが、この他、人の目で鮮鋭度を採点して評価する方法(Lowry,E.Mら)もあります。

lmageJによるPSFによるMTF測定


ImageJを用いたPSFによるボケを説明するためのシミュレーションツールは、Bob Dougherty氏が、いくつかのプラグインを紹介しているということで、改めてご紹介します。
ここでは、前述の②点像の強度分布をフーリエ変換する方法を簡単に示したいと思います。
まず、Gaussian PSF 3D.classというクラスファイルをダウンロードしてlmageJフォルダ下のPluginフォルダ下に置きます。その後、ImageJを起動し、Plugin/Gaussian_PSF_3Dを選択すれば、下図のようなフレームが現れます。

(512×512、1スライス)

適当なマトリックスを設定(例えば512×512、1スライス)して、各パラメータを変更することによって、オリジナルな三次元PSFが作成できます。スライスを512などに増やすことで、スタック画像も作成できます。このスタック画像はImageJ3Dプラグインを使って3次元処理も可能です。
この画像から、点物体(点光源)の像(点像分布関数)を得て、この関数を1次元フーリエ変換することでMTFを得ます。

簡単に2D画像で試して見ます。画像に線を引いて、中央の線のプロファイルを得てから、1Dフーリエ変換します。

テスト画像

線ROIを作成

1D-FFTによる非線形近似(縦軸がMTF、横軸が空間周波数)


マクロの例
macro "FFT1D Tool - VIS" {
    profile = getProfile();//先頭に余計なものが付いてくる時があるので注意
    fft=Array.fourier(profile,"None");
    s=newArray(lengthOf(fft));
    for (i=0;i<lengthOf(fft);i++) {
        s[i]=i;
    }
    for (i=0;i<lengthOf(fft);i++) {
        print("WANTED: s="+s[i]+" FFT:="+fft[i]);
    }
    Plot.create("Fourier spectrum of the line profile","wavenumber","MTF",s,fft);
    Plot.show();
}

条件の異なるMTFを比較することで、MTFが高い方が、客観的に空間分解能が高いと説明することができます。

また、同ホームページからConvolve 3Dプラグインをダウンロードして、このプラグインを同上の方法で選択すれば、オリジナル画像と三次元点像強度分布を三次元的に重畳することが可能です。

ImageJのMTFの測定プラグインとしては、Jeffrey Kuhn氏によるMeasureMTF_.classがありますが、これはデジタルMTFを測定するツールです。このプラグインについては、次回、MTFの詳細を交えて説明を行う予定です。

Reference
  • 山本修司:ImageJで学ぶ実践医用・バイオ画像処理.INNERVISION(20・8) 2005, p98-100

2015年8月10日月曜日

第7回 Imagejを用いて画像の周波数表現を理解する(その2)で学ぶ!

前回は画像における周波数表現の基本的な概念を説明しました。
今回はImageJのプラグインやマクロ言語の作り方を簡単に説明し、画像の周波数表現の解説を交えていきたいと思います。

プラグインとは、自分自身で作ったオリジナルのプログラムをImageJに付け加えることです。
例えば、JavaやImageJの機能を簡単に組み合わせるためのマクロ言語を使ってプログラミングすることによって実行可能になります。
ImageJを通じてコンピュータを使う技術者としての後押しになれば嬉しいです。

どのようにしてImageJのプラグイン機能を使うか?


ImageJに自分自身で作ったプログラムを実装(プラグイン)する場合、画像を必要としないプラグインインターフェースと画像入力を必要とするPlugInFilterインターフェースを使い分けることができます。
インターフェースが持つメソッド(プログラムで作る関数のようなもの)は、前者が、文字列が引数となるrunメソッドvoid run(java.lang.String arg)を宣言するのに対し、後者は、画像処理ベースのrunメソッドvoid run (ImageProccessor ip)を用います。画像処理の場合はint setup(java.langString arg,ImagePlus imp)のメソッドで初期化を行いますが、その際に、取扱う画像の種類に合わせて変数を宣言します。その変数の種類はこちらのField Detailを参照して下さい。

http://rsb.info.nih.gov/ij/developer/api/ij/plugin/filter/PlugInFilter.html

早速、読み込んだ画像の濃淡値を反転させるプログラムを用いて、ImageJプラグインの例を解説します。プラグインコードは次の通りです。

==ここから==
/*
 *ImageJ>Plugins>New>Plugin Filterで、自動生成されるプラグインプログラム。
 */

import ij.*;
import ij.process.*;
import ij.gui.*;
import java.awt.*;
import ij.plugin.filter.*;

public class Filter_Plugin_TEST implements PlugInFilter {
    //画像オブジェクトを定義
    ImagePlus imp;

    public int setup(String arg, ImagePlus imp) {
        //ImageJに読み込んだ画像を指定
        this.imp = imp;
        //ImageJ読み込み可能な全画像形式を対象にする
        return DOES_ALL;
    }

    public void run(ImageProcessor ip) {
        //画像を反転
        ip.invert();
    }

}
==ここまで==

基本的には、Java言語によるクラスファイルの作成方法と同じであるため、Java言語のプログラミング入門書の前半までは読んでいるという読者を対象に話を進めます。
前述したsetupメソッドの中の変数として、return DOES_ALLと記述していますが、これは、ImageJに読み込んだ画像に対して、このプログラムが適応できるようにする設定です。ROIのピクセル値取得コードを含む、プラグインの基本については、次の参考資料を参考にしてみてください。

http://www.gm.fh-koeln.de/~konen/WPF-BV/tutorial-ImageJ_V1.71.pdf

余談:ROIのピクセル値取得コード(8bitグレースケール画像の場合

Rectangle r == ip.getRoiO;
int offset, i;
for (int y = r.y; y〈(r.y+r.height); y++)i
offset = y*width;
for (int x = r.x; x〈(r.x+r.width); x++)l
 i= offset + xl
 pixels[i] = (byte)(255-pixels[i]);
} }

ImageJで取り扱う座標は画像処理では一般的な左上隅が(0,0)の座標値です。Javaは二次元画像でのピクセル配列であっても、一次元配列にデータが納められるので、上記のようにoffset値はスキャンラインごとに取り出すことができるように、offset = y*width;となっています。"i"は常に、スキャンラインの先頭位置から始まり、x方向の幅分だけpixel[i]に格納されていくことになります。このときに扱う画像が8bitグレースケールなので、ピクセル値はbyteでキャストされています。


ちなみに,ここではプラグインに慣れてもらうために,このようなコードの説明を細かく述べましたが,ImageJのメニューのPlugins/New/からPlugin Filterを選択すると,自動的に濃淡反転するコードが現れ,これを.javaの拡張子でPluginsフォルダに保存して、そのエディタのメニューのFile/Compile and Runを選択すると,このクラスのメソッドが実行されて、表示した画像の濃淡が反転します。

ここで,当然,さまざまなクラスやそのメソッドの処理機能を追加するために,ImageJのクラスのリファレンスが必要となります。このAPIドキュメントが次のURLから参照できます。

http://rsb.info.nih.gov/ij/developer/api/overview-summary.html

このImageJの多くのクラスは、システマティックにリンクしているため、クラスの継承などを理解するのに、公開されているクラスダイヤグラムが役立ちます。

http://rsb.info.nih.gov/ij/developer/diagram.html

ここまでは、2次元画像についてお話してきましたが、ImageJではスタック画像を用いた3次元解析ももちろん可能で、そのプラグインも自分で開発することができます。その一例として、ImageJ3DViewerなど(http://3dviewer.neurofly.de/)があります。

ImageJのマクロ機能について


ImageJのマクロ言語を使用すると、前述のプラグインコードを書くよりも非常に簡単にオリジナルの画像処理や画像解析が可能になります。サンプルコードも豊富に紹介されています。
次のURLにマクロ言語の説明が記されています。

http://rsbweb.nih.gov/ij/docs/macro_reference_guide.pdf

ImageJのマクロ言語の作法は非常に簡単です。お決まりの“Hello World”という言葉をImageJで表示する場合,Plugins>New>Macroを選択すると,ブランクのエディタが起動するので、そこに print("Hello world");を記述し、マクロウィンドウのMacrosからRunを選択します。すると別ウィンドウで実行結果が表示されます。マクロはテキスト形式や".ijm"の拡張子で保存可能です。

次に、サンプルとして周波数フィルタのマクロのサンプルを起動させることにします。
まず,ImageJのメニューのPlugins/Macros/Runを選択すると,ファイルダイアログが現れるの
で,ImageJフオルダの中のmacrosフォルダから,FFTCustomFilterDemo.ijmを選択すると,自動的にBridge(174K)というサンプルファイルが呼び込まれ,自動でローパスフィルタとしてのカスタムフィルターを作成し、ガウシアンフィルタ処理(平滑化処理)を行います。そして、このフィルタを用いたローパスフィルタ処理画像を表示します。その後、このローパスフィルタを反転したフィルタをハイパスフィルタとしてハイパスフィルタ処理画像を表示します。
これを少し変更して、マクロの中の画像へのパスを変更してPET画像を対象として起動してみます。(run ("Bridge(174K)")→任意の画像ファイルへのパスに変更、Windowsの場合、パスの区切りはバックスラッシュ2つが必要です。)

画像に含まれる周波数成分のうち低周波成分を残し,高周波成分は除去するようなフィルタをローパスフィルタと言い。その逆をハイパスフィルタと言います。ガウス分布に従って,周波数の遮断を滑らかに移行することによって,画像も自然に近い平滑化がなされます。

ここではPET画像を示していますが、これに限らず、微少石灰化などの高周波成分を解析するときには、ハイパスフィルタなどが強調表示に効力を発揮します。

カスタムフィルター(ガウシアン処理後)




カスタムフィルター プロットプロファイル

原画像


ローパスフィルタ―処理後



 ハイパスフィルタ―処理後


プラグインやマクロはすでに多くのサンプルが公開されているので、これらを複製して、カスタマイズするのが、ImageJを自分なりに使いこなす近道になります。

マクロ言語を用いた周波数空間フィルタについて

前回の記事で画像の「位置と濃淡情報」の世界を,単位距離あたりの濃淡周波数のスペクトルとして表現するための変換方法の1つにフーリエ変換という方法があると説明しました。ここで,たまたまImageJにおけるマクロ言語の説明の流れで,フィルタリングの話になりましたが、フーリエ変換の利用分野は,直交変換の一般的特徴を持つと同時に、次のような応用が試みられています。

①医用画像の特徴抽出
②畳み込み積分の高速計算
③画像の圧縮処理
④CT, MRIやPETにおける再構成
⑤画像問のマッチング処理(レジストレーション)など

画像処理の話の中で,よく「直交変換」という言葉が使われますが、画像というのは、数値が二次元に並んだ行列と考えれば、そのフーリエ変換も実は行列式で説明できるということです。行列式を使うに当たり、ある値(ピクセル値や信号)を、ある位置情報をデカルト座標(x,y, z)で表すと、その基底ベクトル(長さが1で互いに直交している)がフーリエ変換後も直交しているということで、「直交変換」と称されています。このような直交変換は、ほかにもカルーネン・レーべ変換,離散コサイン・サイン変換,ウォルシュ・アダマール変換,スラント変換やバール変換などがあるそうです。

デジタル画像のフーリエ変換の数式については、参考記事やWiki linkをご参照ください。

次回は画像の特徴量解析に話を進めます!

Reference
  • 山本修司:ImageJで学ぶ実践医用・バイオ画像処理.INNERVISION(20・7) 2005, p86-91