4.WSの測定 

【ノイズ特性】
ノイズ特性の復習・確認.増感紙−フィルム系の場合,一様にX線を照射して得られたX線フィルム上の黒化銀の分布は一様ではない.画像を注意深く眺めると,黒化銀の集まり(mottle:モトル)が不規則に変化する模様として観察できる.その一例として,図1にビーズ玉のX線写真を示す.


 図1 異なるノイズレベルで撮影されたビーズ玉のX線写真.
 左:ノイズレベルが低い(粒状性が良い.図2の@の状態).右:ノイズレベルが高い(粒状性が悪い.図2のAの状態).
 画像診断では,小さくて低コントラストな信号も容易に識別できるような画質をもった画像を提供することが大切である.

 図2 微小で低コントラストな信号と画像のノイズレベル

このような画像のザラツキを粒状といい,粒状の示す性質を粒状性あるいはノイズ特性という.ノイズ特性は,RMS粒状度(root mean square)ウィーナースペクトル(WS:Wiener spectrum)[ノイズパワースペクトル(NPS:noise power spectrum)とも呼ばれる]で評価することができる.粒状の主たる原因は,X線が増感紙−フィルム系で吸収される現象がランダムなためである.X線光子のランダムな現象の結果として,X線光子が吸収される位置や密度は統計的なゆらぎをもつ.X線光子のゆらぎの統計的な分布はPoisson分布に従い,

で表わされる.Nは単位面積当たりに吸収された平均のX線光子の個数,σはX線光子の変動を示す標準偏差である.例えば,増感紙−フィルム系で写真濃度1.0を得るとき,単位面積()当たりに吸収された平均のX線光子数が10000個であったと仮定すると,0.1の非常に小さな面積には平均100個のX線光子が吸収されたことになる.しかし,この小さな面積に吸収された平均のX線光子数に対する標準偏差を計算すると,

となり,統計的には10%もゆらいでいることがわかる.つまり,0.1の面積には,常に100個のX線光子がぴったり吸収されるのではなく,約90〜110個のX線光子が吸収されることを意味する.こうしたX線光子の数や分布が統計的にゆらぐことに起因するノイズはX線量子モトル(quantum mottle)とよばれ,増感紙−フィルム系のノイズ=X線写真モトル(radiographic mottle)のなかで最も大きな割合(60%以上)を占めている.図3に増感紙−フィルム系によるX線写真のノイズの構成を示す.


 図3 増感紙−フィルム系によるX線写真のノイズ構成

増感紙の構造モトル(structure mottle)は,増感紙の蛍光体の構造の不均一性に関係するものである.フィルムの粒状(film graininess)は,X線フィルム自体がもともともっている粒状である.また,X線量子モトルと増感紙の構造モトルは,増感紙が発光し,X線フィルムに記録される過程でX線写真モトルに影響するものであり,これらを併せて増感紙モトル(screen mottle)という.ただし,増感紙モトルに占める増感紙の構造モトルの割合は,X線量子モトルと比べて非常に少ない.これは先にも述べたように,X線量子モトルがX線写真モトルの中で最も大きな割合を占めているからである.また,低い空間周波数領域では,X線量子モトルが最も支配的なノイズである一方,X線フィルムの粒状性は空間周波数に依存せずほぼ一定の値を示し,高い空間周波数領域で最も支配的なノイズ因子となる.
図4に,増感紙とX線フィルムを密着させて撮影したX線写真,および,増感紙とX線フィルムの間隔を5段階に変化させて撮影したX線写真を示す.すべての写真は同じX線光子数で撮影されている.

 図4 (a)増感紙とX線フィルムを密着させて得たX線写真,(b)−(f)増感紙とX線フィルムの間にスペースを設けて撮影したX線写真(b:スペース=0.02mm,c:0.05mm,d:0.12mm,e:0.23mm,f:13.0mm)

これらのX線写真を観察すると,増感紙とX線フィルムが密着したとき(aの写真)には粗くて鮮明な粒状の様子が見える.しかし,増感紙とX線フィルムの間隔が大きくなるに従って,粒状の模様がぼかされ,やがて,高い空間周波数成分をもった鮮鋭で細かなX線フィルムの粒状だけが観察できる(fの写真).この結果より,X線写真モトルは,増感紙でX線が吸収され,発光した光がX線フィルムに露光される過程で,光の拡散による影響,すなわち,増感紙によるボケ(MTF)の影響を受けることを示している.
X線写真モトルのなかで最も大きな割合を占めるX線量子モトルのウィーナースペクトル(WSq)は,MTFのほかにも,フィルム特性曲線のグラジェント(G)と,単位面積当たりに吸収した平均のX線光子数(n)の影響を受けることを考慮すると,次式のように近次できることが知られている.

この近次式は,ノイズ特性(この場合はX線量子モトル)が,解像特性や入出力特性とどのような関係をもつのかを理解するのに役立つ.たとえば,X線に対する感度は同じであるが,異なった解像特性(MTF)をもった2種類の増感紙と1種類のX線フィルムを組み合わせた2つのシステムを考える.このときnとGは同じとすると,X線量子モトルのウィーナースペクトルは,解像特性の優れた(MTFが高い)システムのほうが値が高くなることが理解できる.同様に,MTFとnが変わらないと仮定した場合,フィルムコントラストの高い(Gの値が大きい)システムのほうがノイズレベルが高くなることがわかる.

ノイズ特性の評価には,RMS粒状度(root mean square)やウィーナースペクトル(Wiener spectrum)[ノイズパワースペクトル(noise power spectrum)とも呼ばれる]が用いられている.求め方を簡単にまとめると次のようになる.

RMS粒状度は具体的には次の式で求めることができる.Nはデータ数,Diは個々の写真濃度データ,は写真濃度データの平均値である.

この式は,標準偏差の式そのものである.標準偏差は,平均値に対して値がどの程度変動するのかを示すもので,平均値±1σ(σ=RMS)の範囲内に,約68%のデータが含まれることを意味する.RMSの値が大きいほど,変動が大きいことを示すので,粒状性は悪い.RMS粒状度は,計算が簡単で有用ではあるけれど,ノイズ特性の一部を表わしているに過ぎないため,詳細なノイズ解析を行うことはできない.一方,ウィーナースペクトルは,空間周波数の視点からノイズレベルを測定することができ,より詳細なノイズ解析が可能である.ウィーナースペクトルを求める方法は,写真濃度データから自己相関関数(ACF:auto-correlation function)を求めそれをフーリエ変換する方法と,写真濃度データを直接フーリエ変換する方法がある.自己相関関数は次式で表わされる.

自己相関関数は,ξ=0で最大の値をとり,相関距離が大きくなるにつれて小さな値をとる.自己相関関数の式は,同じ関数の重畳積分で表わされていることから,自己相関関数をフーリエ変換すると,2つの同じ関数f(x)のフーリエ変換の掛け算,すなわち,が求まる.これをノイズパワースペクトル(NPS)といい,画像では,ウィーナースペクトルと呼んでいる.したがって,逆に,ウィーナースペクトルを逆フーリエ変換すれば自己相関関数が求まる.もう一方の写真濃度データを直接フーリエ変換する方法についてはあとで詳細を述べるとして,ここではまずはウィーナースペクトルの見方やその解釈について述べる.図5にウィーナースペクトルのグラフの例を示す.


 図5 2つの異なるシステムのウィーナースペクトル

横軸は空間周波数(cycles/mm)であり,縦軸がウィーナースペクトルの値()を示している.ウィーナースペクトルのグラフは,両対数でプロットすることが多い.図5も縦軸,横軸ともに対数で表示されている.ここで注意してもらいたいのは,ウィーナースペクトルの単位はmmの二乗で表わされる面積であるということ!それでは,ある空間周波数においてウィーナースペクトルの値の大小は何を意味するのか.ウィーナースペクトルの値が大きいということは,ノイズのレベルが高い(粒状性が悪い=ノイズが多い)ことを意味している.図5の例でみると,1.0cycles/mmあたりまでの低空間周波数域では,点線のシステムのほうが実線のシステムに比べてノイズ成分が少ないといえる.さらに理解を深めるために,図6に4つの異なる写真濃度の分布と,それらのウィーナースペクトルの模式図を示す.図6のこの関係図は必ず覚えておこう!


 図6 ウィーナースペクトルによるノイズ成分の模式図

(a)と(c)は比較的低い空間周波数成分をもち,写真濃度の変動が大きな(a)のほうがウィーナースペクトルの値が高い.また(b)と(d)は,(a)と(c)より高い空間周波数の成分をもつことを示し,写真濃度の変動が大きな(b)のほうがウィーナースペクトルの値が高い.ここで,RMS粒状度とウィーナースペクトルWS(u)の間には,

の関係があり,ウィーナースペクトルの面積が,RMS粒状度の二乗に対応することが知られている.この関係は,たとえRMS粒状度が同じでも,ウィーナースペクトルは異なった様子を示すことがあることを意味する.例えば,図6において,(a)と(d)でウィーナースペクトルの面積が同じであるとすると,異なるウィーナースペクトルを示すにもかかわらず,RMS粒状度はまったく同じ値と計算されてしまう.これは,RMS粒状度では,ノイズ特性のすべて表わすことはできないことを明示的に示している.したがって,画像のノイズ特性を詳細に調べるためには,RMS粒状度だけでは十分でなく,ウィーナースペクトルの測定が必要かつ重要となる.

さて,ここからはディジタルX線画像システムのノイズ特性評価について述べる.
ディジタルX線画像システムにおけるノイズの評価方法も,増感紙−フィルム系と同様に,RMS粒状度による簡便な方法と,ウィーナースペクトルを用いた方法の両方をそのまま使うことができる.ただ,システムで発生するノイズの構成などディジタル系ならではものがあるので,まずそれらについて簡単に述べた後に,具体的にディジタルX線画像システムでのウィーナースペクトルの測定手順に進む.図7は,ディジタルX線画像システムとしてCRシステムのおもなノイズの要因をまとめたものである.

 図7 CRシステムのおもなノイズ構成

ディジタルX線画像システムの構成は増感紙−フィルム系と比べて複雑であり,ノイズの主たる原因であるX線量子モトルのほかにも,X線検出器の構造ノイズ(この例ではイメージングプレートの構造ノイズ),輝尽発光の光量子ノイズ電気系ノイズ量子化ノイズなど多くの因子がある.これらの中で,構造ノイズ,電気系ノイズ,量子化ノイズは,撮影線量に依存せず一定の値を示すため固定ノイズと呼ばれている.さらにシステム全体のノイズ特性には画像表示部のノイズも付加される.ディジタルX線システムでは,被曝線量軽減の観点から撮影線量をコントロールして使用するため,ディジタルX線システムでもX線量子モトルが支配的であることは増感紙−フィルム系と変わらない.一方,高線量域ではX線検出器の構造ノイズが支配的であり,この領域では撮影線量を増加させてもノイズ特性は改善されないので注意が必要である.このほかにも,画像処理によってノイズ特性が改善したり,逆にノイズが増加することもある.次式は,ディジタルX線画像システムの総合的な二次元のウィーナースペクトル[オーバーオールウィーナースペクトル,WSoverall(u,v)]を示したものである.

WSA(u,v)はアナログ成分のウィーナースペクトルで,通常,X線量子モトルが大きな割合を示す.OTFS(u,v),OTFF(u,v),OTFD(u,v)は,それぞれサンプリングアパーチャ,画像処理フィルタ,ディスプレイアパーチャのレスポンス関数である.WSE(u,v)は,画像出力部の付加ノイズである.{ }内がディジタルウィーナースペクトルであり,[ ]内がプリサンプリングウィーナースペクトルである.エリアシングエラーの効果が無視できる場合には,プリサンプリングウィーナースペクトルとディジタルウィーナースペクトルは等しくなる.実際には,WSA(u,v)にOTFS(u,v)の二乗が掛かっているため,これらのスペクトルはほぼ等しくなると考えられ,ディジタル系のノイズ解析にはディジタルウィーナースペクトルが有用となる.さらにオーバーオールウィーナースペクトルには,OTFF(u,v)の二乗とOTFD(u,v)の二乗が掛かってくるので,高い空間周波数域におけるエリアシングの影響はほとんど無視できるくらいに小さくなる.このような場合には,MTFとは異なり,オーバーオールウィーナースペクトルもノイズの解析に有効となる.

ディジタルウィーナースペクトルの測定には均一にX線照射を行って得られたディジタルX線画像を利用する.被写体を置かずに撮影し,そのときの線量を記録しておく.こうした画像データからディジタルウィーナースペクトルを測定する方法として,仮想スリット法二次元フーリエ変換法がある.図8に仮想スリット法の測定手順,図9に二次元フーリエ変換法の測定手順を示す.

  図8 仮想スリット法によるウィーナースペクトル測定の概略

 


  図9 二次元フーリエ変換法によるウィーナースペクトル測定の概略

 

実験1: 仮想スリット法でディジタルウィーナースペクトルを測定してみよう

【目的】 仮想スリット法によるWS測定法を習得する.

【手順】

1.スリットトレースデータの取得

⇒まず一様曝射で撮影したノイズサンプルの画像データ(Routine1436.raw:1800×1800ピクセル,2byteデータ,サンプリングピッチ100μm,X線装置:島津製作所 UD150L−30,画像読取装置:FCR XG−1N,管電圧:42 kV,管電流時間積:2.8mAs,焦点サイズ:小焦点,FFD:60cm,付加フィルタ:0.5mm Al + 0.1mm Cu)をダウンロードする.以降の作業はImageJで行う.

I.画像(Routine1436.raw)を開く.ImageJのメニューFile→Import→RawでRoutine1436.rawを選択.

←このように入力してOKをクリック

II.スリット長30ピクセル,トレース長512ピクセルのスリットトレースを実行


↑このように方形ROIを設定して,メニューのAnalyze→Plot Plofileを実行する.


↑これが1×30ピクセルの仮想スリットでトレースしたデータ(トレース長512)となる.Copyボタンを押すと数値データがメモリーにコピーされるのでそれをExcelに貼り付けて以降の計算で使用する.

2.以降,トレンド除去処理,トレースデータのフーリエ変換,ピクセル値のウィーナースペクトルWSΔp算出,X線強度のウィーナースペクトルWSΔE/Eへの変換を順にExcel上で行う.

A列: トレースの位置 →トレース長は512ピクセルなので0〜511.
B列: 仮想スリットのトレースデータ →手順1のデータ
C列: トレースデータ(B列)を多項式で近次した値(=多項式近次によって平滑化したのと同義).
    →よくフィッティングするように次数を選ぶこと.
D列: トレンド除去の実施 →オリジナルのトレースデータから多項式近次したデータを減算する(B列−C列).
    (ここまでの処理結果をグラフ化して,トレンド除去の成否を確認しておくこと!)
E列: 測定に必要な諸条件 →ピクセル寸法,スリット長,トレース長は,各自で計算して入力すること.
    E2→ピクセル寸法(mm).サンプリング間隔が100μmなので... 
    E4→スリット長(mm).30ピクセルだったのでこれをミリメートルに変換すると...
    E5→トーレス長(mm).512ピクセルだったのでこれをミリメートルに変換すると... 
    E6→特性曲線傾きG.別実験で実測済の4038.8をそのまま使用する. 
    E8→グレア含有率k.別実験で実測済の0.065をそのまま使用する. 
    E10→log10e.計算式の通り.
F列: トレンド除去後のトレースデータ(D列)を高速フーリエ変換(FFT)する.
    →Excelのメニュー:ツール→分析ツール→フーリエ解析を選択.
    
    入力範囲はトレンド除去後のトレースデータ(D列),出力先はF列に指定する.
G列:複素数で出力されるフーリエ変換後データの絶対値を求める.→複素数の絶対値を求めるIMABS( )関数を使用する.
H列:ピクセル値のウィーナースペクトルWSΔp算出 →算出式は図8中の説明を参照すること.
I 列:ピクセル値のウィーナースペクトルWSΔpをX線強度のウィーナースペクトルWSΔE/Eへ変換 →図8の定義式そのまま.
J列:ウィーナースペクトルのグラフの横軸に対応する空間周波数(cycles/mm)の計算
   →1/トレース長(mm)×n (nはA列のトレース位置.なぜこの式になるのか考え理解した上で使用すること!)
K列:X線強度のウィーナースペクトルWSΔE/E(I列)に対して,変動成分を取り除くために移動平均法で平滑化を行う.
   →20区間で移動平均を行う.例:K2の値は,I2からI21までの平均値.K3の値は,I3からI22までの平均値....

*)複数のスリットトレースデータから求めた複数本のWSを平均することで測定精度を高めることができる.

グラフの作成:
@空間周波数(J列)を横軸に最終的なウィーナースペクトル値(K列)を縦軸にプロットする.
A縦軸,横軸ともに対数表示とする.目盛りの表示形式は,縦軸は指数で,横軸は標準にする.
B空間周波数0およびナイキスト周波数以上は表示しない.

【考察】
1.今回の実験画像のナイキスト周波数はいくつか?
2.グラフ作成途中@の段階でナイキスト周波数を境にグラフはどのようになっているか確認し,その理由を考える.


実験2: 解像特性の異なるノイズサンプルでウィーナースペクトルを測定してみよう.

【目的】 解像特性とノイズ特性の関係を理解する.

【手順】 実験1とは解像特性の異なるノイズサンプルの画像データ(Routine1435.raw:1800×1800ピクセル,2byteデータ,サンプリングピッチ100μm,X線装置:島津製作所 UD150L−30,画像読取装置:FCR XG−1N,管電圧:42 kV,管電流時間積:2.8mAs,焦点サイズ:小焦点,FFD:60cm,付加フィルタ:0.5mm Al + 0.1mm Cu)を使用する.測定手順は,実験1と同じ.

【考察】
実験1と実験2のノイズサンプルの画像データは,同じ撮影条件と同じ装置で撮られた画像である.異なるのは,検出器(イメージングプレート)の種類だけである.一方は,標準タイプ(ST-VI)のイメージングプレートで撮影された画像,もう一方は,高鮮鋭度タイプ(HR-V)のイメージングプレートで撮影された画像である.高鮮鋭度タイプは,標準タイプに比べて,解像特性に優れている(MTFが高い).実験1と実験2の測定結果から,どちらが標準タイプでどちらが高鮮鋭度タイプであるか考えてみよう.