今回はTIGREというCTなどに利用可能な再構成ソフトを紹介します。
TIGREはCERNが公開している再構成ソフトウェアで、Tomographic Iterative GPU-based Reconstruction Toolbox（TIGRE）といいます。ソフトウェアはgithubにあります。名前にもある通りGPUを使って高速に再構成することが可能です。
GPUとCUDAの準備など手間だなと思ってしまいますが、Google colab.で実行可能でしたので、その記事を書いていきます。
Google colabratoryの準備
Google colab.を使うことによってほとんど準備を必要とせずにpythonを使える環境が手に入ります。
ネットに膨大な情報がありますので、ここでgoogle colab.の使い方を紹介する必要はないでしょう。このサイトを見て、準備してください。
Googleアカウントさえあればすぐに使えます。
ランタイムの変更
Google colabは初期設定ではCPUを使うようになっていますので、GPUを使えるようにランタイムを変更します。
T4 GPUというのを使用します。
たったこれだけです。
ただし、GPUの利用には制限がありますので限度を超えると一定期間は使えなくなります。
その際には違うアカウントでgoogle colab.を使えばまた使えます。
作業手順
# 必要なパッケージをインストール
!pip install numpy scipy matplotlib scikit-image
# TIGREのリポジトリをクローン
!git clone https://github.com/CERN/TIGRE.git最初のパッケージは念のためで、必要ないかもしれません。2つ目のコマンドでgithubのコードをコピーします。
実行した結果は以下の通りです。
次にクローンしてできたフォルダの中に入り、インストールします。
# 作業フォルダに移動
%cd TIGRE/Python
# インストール
!pip install .
とりあえず、errorは出ていないようなのでOKです。
インストールが終わったらexample.pyを実行してもよいです。動作確認になります。
ただし、実行しただけでは画像が表示されないので、画像を保存するためにexample.pyの一番最後の行にsavefigを追記します。
上記のようにexample.pyを選択して（ダブルクリック）、表示されたファイル内容の最後に以下を加えます。
fig.savefig('reconstruction_results.png')そうしたら、以下のコマンドで実行します。
# example.pyの実行
!python example.py
reconstruction_results.pngがPythonフォルダに保存されているので、自分のPC（ローカル）にダウンロードしましょう。
Total variation（TV）正則化再構成
では、最後に一つ例題をやろうと思います。
いくつかサンプルコードがありますが、今回はTotal variation（TV）を利用したコードを使ってみます。TVはMRIなどでよく使われますが、スパースな情報から再構成するのに適した手法です。また、ノイズの低減にも有用です。TV正則化は、画像内の変化の合計（Total Variation）を最小化するように画像を更新することによって、平滑な部分を滑らかにしつつ、エッジのような重要な部分を保ちます。
先ほどと同様に編集します。
まず、以下の関数を定義します。これは画像を保存するための関数です。
パッケージをimportする文のすぐ後にでも書いておいてください。
# 各再構成画像の全スライスを TIFF 形式で保存
def save_all_slices_as_tiff(image_array, filename):
# 各スライスをそのままPillowに渡し、無圧縮TIFFで保存
tiff_images = []
for i in range(image_array.shape[0]):
# 画像データを8bitに変換せず、そのままPillowのImageに変換
tiff_image = Image.fromarray(image_array[i, :, :])
tiff_images.append(tiff_image)
# TIFFファイルとして無圧縮で保存
tiff_images[0].save(filename, save_all=True, append_images=tiff_images
, compression=None)
次に、途中でprojectionデータが出てくるので、一応サイノグラムを保存しておこうと思います。projectionデータとsinogramは軸が違うだけで、同じものです。
以下のコードをnoise_projections = …の後ろに追加してください。
from PIL import Image
# サイノグラム保存用のlist
tiff_images = []
for i in range(noise_projections.shape):
# 各スライスのデータ（全角度からのサイノグラム）を抽出
sinogram_slice = noise_projections[:, i, :]
tiff_image = Image.fromarray(sinogram_slice)
tiff_images.append(tiff_image)
# マルチページのTIFF画像として保存
tiff_images[0].save("sinogram_all_slices.tiff", save_all=True, append_images=tiff_images, compression=None)
ちなみに保存したサイノグラムは以下のようなものです。
軸を変えて表示すると・・・projection画像も見れますね。
最後に、再構成画像を保存するために一番最後に以下のコードを追加してください。
# 各アルゴリズムの再構成結果を TIFF に保存
save_all_slices_as_tiff(imgAWASDPOCS, "AWASDPOCS.tiff")
save_all_slices_as_tiff(imgOSASDPOCS, "OSASDPOCS.tiff")
save_all_slices_as_tiff(imgASDPOCS, "ASDPOCS.tiff")
save_all_slices_as_tiff(imgOSSART, "OSSART.tiff")
#save_all_slices_as_tiff(imgCGLSTV, "CGLSTV.tiff")
再構成された画像です。
左から順にOSSART、ASDPOCS、OSASDPOCS、AWASDPOCSです。
OSSARTは代数的再構成法のSARTのorder-subset版です。MLEMとOSEMの関係のような感じですね。
それ以外の3つはTV再構成になります。
ASD-POCSが基本、OS-ASD-POCSが高速化、AW-ASD-POCSがエッジ保存平滑化という関係になります。どれもOSSARTよりも画質が良いのが分かります。
ちなみにd09_Algorithms04.pyのgeoを以下のように変更すると（False→True）高解像度の結果が得られますが、かなり時間がかかります。
#geo = tigre.geometry_default(high_resolution=False)
geo = tigre.geometry_default(high_resolution=True)high_resolution=Trueはどのような条件なのかを確認したい場合には、
TIGRE/Python/tigre/utilities/geometry_default.pyの中身を確認しましょう。
※ちょっと書き換えていますので、オリジナルとは少し違います
自分が使用したいCT装置で実験したい場合には、Geometryを正しく設定すればファントムから投影データを作成、再構成ができます。ファントムは
head = sample_loader.load_head_phantom(geo.nVoxel)ここで指定しているようにサンプルファントムデータを読み込んでいますので、これを自分の使いたいファントムに変更可能です。
noise_projections = CTnoise.add(projections, Poisson=1e5, Gaussian=np.array([0, 10]))もしプロジェクションデータやサイノグラムがあるのであれば、noise_projectionsに入れれば良いです。
ファントムもプロジェクションデータもndarray形式です。
他にもいろいろなサンプルコードがありますが、機会があればいつか説明します。
