JP2016127907A - 医用情報処理装置、ｘ線ｃｔ装置、医用情報処理方法、及び医用情報処理プログラム - Google Patents

Abstract

【課題】複数の物質それぞれの投影長を正確に求めること。
【解決手段】実施形態の医用情報処理装置は、取得部と、第１決定部と、推定部と、第２決定部とを備える。取得部は、Ｘ線の光子の入射に応じて信号を出力する光子計数検出器の出力信号を用いて、被検体を透過したＸ線の光子の計数値を、複数のエネルギー帯域ごとに取得する。第１決定部は、複数の物質それぞれの投影長をパラメータとして含むＸ線の減衰モデルに入力する初期投影長を、複数の物質ごとに決定する。推定部は、第１決定部により決定された複数の物質ごとの初期投影長と、減衰モデルと、光子計数検出器のＸ線に対する応答特性を示す検出器モデルとを用いて、複数のエネルギー帯域ごとの計数値を推定する。第２決定部は、推定部が推定した計数値と、取得部が取得した計数値との差に基づくコスト関数を用いて、反復的最適化方法により、被検体内における複数の物質それぞれの投影長を決定する。
【選択図】図１

Description

本発明の実施形態は、医用情報処理装置、Ｘ線ＣＴ装置、医用情報処理方法、及び医用情報処理プログラムに関する。

投影データは、コンピュータ断層撮影、Ｘ線撮影、マンモグラフィ、トモシンセシス等、多くの応用例がある。投影データは、被検体に放射線を透過させ、透過させた放射線について被検体不在時との変化を検出することで、被検体の内部構造を明らかにする。吸収イメージングにおいては、投影データは、放射線による光線に沿った減衰のラドン変換を表す。コンピュータ断層撮影は、一連の投影角度での投影データを使用して、被検体の内部構造の画像を再構成する。

コンピュータ断層撮影（Computed Tomography：ＣＴ）システム、及び、その方法は、広く利用されており、特に医療用の撮影や診断に利用されている。一般的に、ＣＴシステムは、被検体の身体から１つ以上の断面スライスの画像を作成するものである。Ｘ線源等の放射線源が、一方側から身体に放射線を照射する。一般的に、Ｘ線源に隣接するコリメータがＸ線ビームの角度範囲を制限し、それにより、身体に作用する放射線は各身体断面スライスを画定する平面領域（すなわち、Ｘ線投影面）に実質的に制限される。身体の反対側にある少なくとも１つの検出器（一般的には１より更に多い検出器）が、身体の投影面を透過した放射線を受ける。身体を透過した放射線の減衰を、検出器から受信した電気信号を処理することによって測定する。実施形態によっては、マルチスライス検出器の構成が使用され、身体の平面投影画像ではなく３次元的な投影画像が提供される。

通常、Ｘ線源は、身体の長軸を中心に回転するガントリに搭載される。Ｘ線源の反対側で、複数の検出器も同様にガントリに搭載される。一連のガントリ回転角度において投影による減衰測定を行い、ＣＴ再構成アルゴリズムを用いて測定結果のデータを処理することで、単一の身体断面画像が得られる。複数の断面画像又は３次元の画像を得るには、Ｘ線源及び検出器を身体と相対的に平行移動させる必要がある。身体をガントリと相対的に平行移動させ、場合によっては、Ｘ線源を異なる角度及び／又は軸方向の位置にして行った減衰測定を含む複数の像が得られる。ＣＴシステムによっては、身体の平行移動と身体に相対的なガントリ回転の組み合わせによって、Ｘ線源が身体に対して渦巻き状又はらせん状の軌道を通る。複数の像を使用して、スライス又は複数のスライスを示すＣＴ画像を再構成する。

Ｘ線ＣＴは、医療及びその他多様な目的におけるＣＴの最も一般的な形態であり、ＣＴには他の形態（陽電子放出断層撮影、単一光子放出コンピュータ断層撮影等）もあるが、「コンピュータ断層撮影（computed tomography）」という用語のみでＸ線ＣＴを指すことが多い。Ｘ線ＣＴを意味する他の用語には、コンピュータＸ線体軸断層撮影（ＣＡＴ（Computed Axial Tomography）スキャン）及びコンピュータ支援断層撮影（computer-assisted tomography）がある。

Ｘ線ＣＴの一例では、被検体の周囲を回転するＸ線源を使用してＸ線スライスデータを生成するものがあり、Ｘ線源から画像用被検体を挟んで反対側に、複数のＸ線センサが配置される。装置により、固定された被検体の周囲でＸ線源及び検出器を回転させる。１回転が完了すると、被検体は自身の軸に沿って移動され、次の回転が開始される。

比較的新しい装置では、Ｘ線リング内で撮影対象の被検体をゆっくりと滑らかに移動させつつ連続的に回転できるものもある。そのような装置は、ヘリカルＣＴ装置又はスパイラルＣＴ装置と呼ばれている。回転軸に垂直な軸方向に沿って、多数の検出器列を備えるシステムは、Ｘ線ビームの形状からコーンビームＣＴと呼ばれることが多い。

一連の異なる投影角度で身体における投影測定を行うことにより、一方の軸に沿った検出器アレイの空間次元及びもう一方の軸に沿った時間／角度の次元を利用して、投影データからサイノグラムが構成される。身体におけるある特定のボリュームに起因する減衰は、ＣＴシステムの回転軸を中心とする正弦波であって、サイノグラムの空間次元に沿って振動する正弦波を描く。身体のボリュームは、回転の中心から遠いほど、より大きな振幅の正弦波に対応する。サイノグラムにおける各正弦波の位相は、回転軸を中心とする相対的な角度位置に対応する。逆ラドン（Radon）変換（あるいは同等の画像再構成方法）を行うことにより、サイノグラム中の投影データから画像が再構成され、再構成された画像は、身体の１つの断面スライスに対応する。

従来、エネルギー積分検出器は、ＣＴ投影データを測定するために使用されてきた。現在では、最近の技術の進展により光子計数検出器（Photon-Counting Detector：ＰＣＤ）が、エネルギー積分検出器の実用的な代替選択肢となりつつある。ＰＣＤには、スペクトルＣＴの実行能力を含め多くの利点がある。ＰＣＤは、透過されたＸ線データのスペクトル上の性質を得るために、Ｘ線ビームをエネルギー成分すなわちスペクトルビンに分割し、各ビンの光子数をカウントする。

スペクトルＣＴ技術は、物質弁別及びビームハードニング補正の改善をもたらし、多くの臨床上の応用に有益である。更に、半導体に基づくＰＣＤは、従来のスペクトルＣＴ技術（例えば、二重エネルギー源、ｋＶｐスイッチング等）と比較して、提供されるスペクトルの情報が改善されており、スペクトルＣＴにとっての有望な候補である。

スペクトルＣＴ及びスペクトルＸ線イメージング全般の利点の一つとして、異なる原子番号Ｚの原子を持つ物質は、異なるスペクトルプロファイルを示して減衰することが挙げられる。それにより、複数のＸ線エネルギーの減衰を測定すれば、物質が、その構成原子のスペクトル吸収プロファイル（すなわち、物質の有効Ｚ）に基づき識別される。このような物質の識別により、スペクトル領域から物質的領域へのマッピングが可能となる。このマッピングは、従来、物質弁別と呼ばれる。

スペクトルＣＴデータの物質弁別が可能なのは、生体物質内でのＸ線の減衰が、光電効果による散乱及びコンプトン（Compton）散乱という２つの物理的ステップに支配されているからである。したがって、エネルギーの関数である減衰係数は、下記の式（１）の分解によって近似される。

ここで、μ_ＰＥ（Ｅ，ｘ，ｙ）は、光電効果による減衰であり、μ_C（Ｅ，ｘ，ｙ）は、コンプトン減衰である。この減衰係数は、物質１（例えば、骨等の高Ｚ物質）と物質２（例えば、水等の低Ｚ物質）との分解に書き換えることができ、下記の式（２）となる。

ここで、ｃ_１（ｘ，ｙ）及びｃ_２（ｘ，ｙ）は、それぞれ位置（ｘ，ｙ）における物質１及び物質２の密度を表す空間的な関数である。

半導体ベースのＰＣＤは、スペクトルＣＴに利用した場合に固有の利点がある一方で、固有の難点もある。例えば、検出器応答において非線形性及びスペクトルのずれを補正しないままでは、半導体ベースのＰＣＤから再構成された画像の画像品質が良くない場合がある。検出器応答の補正は、パイルアップ、弾道欠損効果（ballistic deficit effect）、極性効果、特性Ｘ線の逃散、及び空間電荷効果に対応して行われる。

検出器応答の補正及び物質弁別を組み合わせると、複雑な問題が生じる。検出器応答を考慮しつつ投影データを投影長と呼ばれる物質成分に分解する反復的方法がある。しかしながら、そのような方法は、問題の定式化における大域的な最小点ではなく、局所的な最小点に対応する誤った解に行き着く可能性がある。したがって、最適な投影長を得るには、投影長について適切に選択された初期推定を用いて反復的方法を開始することが重要である。

米国特許出願第１３／４２６，９０３号明細書 米国特許出願第１３／９０６，１１０号明細書 米国特許第８，１９４，９６１号明細書

本発明が解決しようとする課題は、複数の物質それぞれの投影長を正確に求めることができる医用情報処理装置、Ｘ線ＣＴ装置、医用情報処理方法、及び医用情報処理プログラムを提供することである。

実施形態の医用情報処理装置は、取得部と、第１決定部と、推定部と、第２決定部とを備える。取得部は、Ｘ線の光子の入射に応じて信号を出力する光子計数検出器の出力信号を用いて、被検体を透過したＸ線の光子の計数値を、複数のエネルギー帯域ごとに取得する。第１決定部は、複数の物質それぞれの投影長をパラメータとして含むＸ線の減衰モデルに入力する初期投影長を、前記複数の物質ごとに決定する。推定部は、前記第１決定部により決定された前記複数の物質ごとの初期投影長と、前記減衰モデルと、前記光子計数検出器のＸ線に対する応答特性を示す検出器モデルとを用いて、前記複数のエネルギー帯域ごとの計数値を推定する。第２決定部は、前記推定部が推定した計数値と、前記取得部が取得した計数値との差に基づくコスト関数を用いて、反復的最適化方法により、前記被検体内における前記複数の物質それぞれの投影長を決定する。

図１は、実施形態に係るＣＴスキャナの構成例を示すブロック図である。 図２は、実施形態に係るＣＴスキャナによる処理手順の一例を示すフローチャートである。 図３は、コスト関数と投影長との関係の一例を示す図である。 図４は、実施形態に係るＣＴスキャナによる処理手順の一例を示すフローチャートである。 図５は、実施形態に係るＣＴスキャナによる処理手順の一例を示すフローチャートである。 図６は、コスト関数と投影長との関係の一例を示す図である。 図７は、コスト関数と投影長との関係の一例を示す図である。 図８は、実施形態に係るＣＴスキャナによる処理手順の一例を示すフローチャートである。 図９Ａは、エネルギー積分検出器からの非スペクトルＣＴデータを用いた、物質１に関する推定上の物質弁別の再構成画像を示す図である。 図９Ｂは、スペクトルＣＴデータを用いた、物質１に関する物質弁別の再構成画像を示す図である。 図１０Ａは、図９Ａに示す推定上の物質１再構成画像上の順投影から得られるサイノグラムの図である。 図１０Ｂは、図９Ｂに示す物質１画像を再構成する際に使用されるスペクトルＣＴデータからの物質１の最適投影長のサイノグラムの図である。 図１１は、図１０Ａのサイノグラムのラインアウトを図１０Ｂのサイノグラムの対応するラインアウトと比較するグラフである。 図１２Ａは、エネルギー積分検出器からの非スペクトルＣＴデータを用いた、物質２に関する推定上の物質弁別の再構成画像を示す図である。 図１２Ｂは、スペクトルＣＴデータを用いた、物質２に関する物質弁別の再構成画像を示す図である。 図１３Ａは、図１２Ａに示す推定上の物質２再構成画像上の順投影から得られるサイノグラムの図である。 図１３Ｂは、図１２Ｂに示す物質２画像を再構成する際に使用されるスペクトルＣＴデータからの物質２の最適投影長のサイノグラムの図である。 図１４は、図１３Ａのサイノグラムのラインアウトを図１３Ｂのサイノグラムの対応するラインアウトと比較するグラフである。 図１５は、物質弁別投影長の推定を取得及び／又は改善するランダムサンプリング方法の実施形態のフロー図である。

本開示は、添付の図面とともに以下の詳細な説明を参照することによって、より完全に理解されるであろう。

光子計数検出器（ＰＣＤ）を用いたスペクトルＣＴでは、画像再構成は、検出器応答及び物質弁別の補正ステップを含む複数の前処理ステップの後に行われる。実施形態のなかには、そのような前処理ステップのいくつかを組み合わせて行うことが有効となることもある。

例えば、検出器応答及び物質弁別の補正は、下記の式（３）に示す検出器応答の等式を最も良く満たす物質１及び物質２の投影長を求めることで同時に解決できる。

ここで、Ｎ_ｍは、ＰＣＤのｍ個目のエネルギービンにおける計数（計数値）であり、Ｔは、積分時間であり、Ｗ_ｍ（Ｅ）は、ｍ個目のエネルギービンの検出スペクトルを示すウィンドウ関数（すなわち、Ｗ_ｍ（Ｅ）は、エネルギーＥの検出事象がｍ個目のエネルギービンに記録される確率密度を与える）である。Ｓ_Ｌｉｎ（Ｅ）は、ＰＣＤの線形反応であり、下記の式（４）で与えられる。

また、Ｓ_{Ｎｏｎｌｉｎ}（Ｅ）は、ＰＣＤの非線形な検出器応答であり、下記の式（５）で与えられる。

ここで、省略記号は、高次の項を表し、それらは高い線束密度において意味が大きくなる。また、Ｒ_０は、線形な検出器応答の関数であり、Ｒ_１は、１次のパイルアップ検出器応答の関数であり、τは、検出器の不感時間である。被積分関数中のスペクトル項Ｓ（Ｅ）は、Ｘ線の光吸収の後に起こる量子効率、空間電荷効果、弾道欠損効果、パイルアップ、特性Ｘ線の逃散等の効果を含まない検出器への入射スペクトルである。このスペクトル項は、下記の式（６）で与えられる。

ここで、μ_１及びμ_２は、物質弁別についての基本物質の減衰係数であり、Ｌ_１及びＬ_２は、各投影長である。各検出器についてのＸ線束は、下記の式（７）で与えられる。

ｎ_ａｉｒ＝Ａ×Ｉ_ｒｅｆは、Ｘ線源の参照強度測定Ｉ_ｒｅｆに基づく計算上の線束であり、Ａは較正係数である。係数Ｓ_ａｉｒは、撮像対象の被検体ＯＢＪが不在の場合のＸ線スペクトルを表す（すなわち、Ｓ_ａｉｒは、撮像対象の被検体ＯＢＪを空気で置き換えた場合のＳと等しい）。同様に、係数ｎ_ａｉｒは、撮像対象の被検体ＯＢＪが不在の場合の検出器におけるＸ線束を表す。各投影測定において、参照強度Ｉ_ｒｅｆは、被検体ＯＢＪより手前のＸ線源付近の参照検出器によって測定される。参照強度は、非スペクトル型の参照検出器（例えば、エネルギー積分検出器）で測定される。

上記の検出器応答の等式を解き投影長Ｌ_１及びＬ_２を求める方法は、少なくとも２つある。一実施形態としては、検出器モデルに基づく反復的方法において、投影長の推定を使用して検出器の信号を補正し、その補正された信号を用いて反復的に投影長の推定を更新する。検出器モデルに基づく反復的な方法は、検出器応答のモデルを仮定して、投影長Ｌ_１及びＬ_２の第１の推定を計算する。次に、この第１の推定を用いて、例えば２つのエネルギー物質への分解による物質弁別を行うことにより、投影長Ｌ_１及びＬ_２の第２の推定を計算する。第２の推定は、次に、仮定された検出器応答のモデルを用いて第３の推定を計算するのに用いられ、更に、投影長Ｌ_１及びＬ_２の推定が収束するまで同様の計算が行われる。

上記の検出器応答の等式を解き投影長Ｌ_１及びＬ_２を求める別の方法では、コスト関数を利用する。コスト関数による方法では、検出器応答の等式を使用して、測定上の計数Ｎ′_ｍと計算上の計数Ｎ_ｍとの差を示すコスト関数の大域的な最適点を求めることにより、投影長の値を求める。いくつかの異なるコスト関数φ（Ｌ_１、Ｌ_２）が考えられる。一実施形態では、コスト関数は、測定上の計数Ｎ′_ｍと計算上の計数Ｎ_ｍとの差の最小二乗、すなわち、下記の式（８）である。

一実施形態では、コスト関数は、測定上の計数Ｎ′_ｍと計算上の計数Ｎ_ｍとの差の重み付き最小二乗、すなわち、下記の式（９）である。なお、σ_ｍは、Ｎ′_ｍの標準偏差である。

一実施形態では、コスト関数はポアソン尤度関数、すなわち、下記の式（１０）である。

この投影長検出器応答の問題への解は、どちらの方法の場合も、投影長の初期推定から始まる反復的方法に基づく。それらの反復的方法は、投影長の初期推定が適切に選択されて開始された場合に最も効果が高い。ここで、投影長の初期推定を選択する方法を説明する。

初期推定方法の一実施形態では、減衰画像は、非スペクトル型のエネルギー積分検出器からの投影データを用いて再構成される。この再構成画像は、次に、第１の物質成分と第２の物質成分にそれぞれ対応する領域に更に分割される。第１の物質成分の投影長は、第１の物質に対応する再構成画像の領域に対する順投影により得られ、第２の物質成分の投影長は、第２の物質に対応する再構成画像の領域に対する順投影により得られる。この実施形態は、エネルギー積分検出器及びＰＣＤをともに備えるハイブリッド型ＣＴシステムに特に適用可能であるが、例えば、スペクトルＣＴデータのスペクトル次元にわたって総和を求めるというような他のやり方で非スペクトルＣＴデータが得られる場合にも適用可能である。すなわち、取得部としての後述の処理部１７０は、被検体を透過したＸ線の強度を検出するエネルギー積分検出器の出力信号を取得する。そして、第１決定部としての処理部１７０は、エネルギー積分検出器の出力信号から再構成される再構成画像に基づいて、初期投影長を決定する。例えば、第１決定部としての処理部１７０は、再構成画像を第１物質成分及び第２物質成分に物質弁別を行い、第１物質成分に対して順投影を行うことで第１物質成分に対応する第１初期投影長を決定し、第２物質成分に対して順投影を行うことで第２物質成分に対応する第２初期投影長を決定する。

初期推定方法の別の実施形態では、投影長の初期推定は、サンプル空間内で無作為に選択した１組の投影長についてコスト関数を計算することで得られる。初期推定は、次に、コスト関数の最小値に対応する投影長を選択することにより、それらの無作為に選択した投影長のなかから選択される。コスト関数は、検出器応答の等式を用いて、測定上の計数Ｎ′_ｍと計算上の計数Ｎ_ｍとの差を用いて表される。この実施形態は、スペクトルＣＴデータを収集できるＣＴシステムであれば適用可能である。すなわち、第１決定部としての処理部１７０は、探索領域内の投影長に対応するサンプル点を無作為に複数選択し、選択した複数のサンプル点のうちコスト関数が最小値となるサンプル点の投影長を初期投影長として決定する。例えば、取得部としての処理部１７０は、更に、被検体を透過したＸ線の強度を検出するエネルギー積分検出器の出力信号を取得する。そして、第１決定部としての処理部１７０は、エネルギー積分検出器の出力信号から再構成される再構成画像を第１物質成分及び第２物質成分に分解し、第１物質成分に対応する第１投影長と、第２物質成分に対応する第２投影長との組み合わせのパターンの中から複数の組み合わせを無作為に選択し、選択した組み合わせのうちコスト関数が最小値となる組み合わせの投影長を初期投影長として決定する。

初期推定方法の更に別の実施形態では、初期推定の選択は、コスト関数におけるより急な勾配に対応する初期推定及び／又は大域的な最小点である解に至る可能性がより高い初期推定を優先的に選択するために、入射線束密度ｎ_ａｉｒに応じて調整される。より大きな勾配に対応する初期推定を選択すれば、勾配降下に基づく反復的方法は、高速に収束する可能性が高くなる。すなわち、第１決定部としての処理部１７０は、光子計数検出器に入射されるＸ線の線束密度に応じて、初期投影長を調整する。

コスト関数の一実施形態では、低線束密度の場合のコスト関数の勾配は、一般的に、最適な投影長より短い投影長ほど大きい。それに対し、高線束密度の場合のコスト関数の勾配は、一般的に、最適な投影長より長い投影長ほど大きい。高線束密度の場合、最適な投影長より長い投影長の初期推定ほど、極小点が避けられる。

上記の実施形態は、種々の変更を加えて組み合わせることができる。すなわち、取得部としての処理部１７０は、複数の物質それぞれの投影長をパラメータとして含むＸ線の減衰モデルに入力する初期投影長を、前記複数の物質ごとに決定する。

以降、図面を参照するが、複数の図面間において同一又は対応する要素には同様の参照符号を付す。図１は、複数のエネルギー積分検出器が第３世代の配列で配置され、複数の光子計数検出器（ＰＣＤ）が第４世代の配列で配置されたハイブリッド型ＣＴスキャナシステムの概略図を示す。その開示内容全体が参照として本明細書に組み入れられる特許文献１に記載されるように、ＣＴスキャナは、ＰＣＤリングの内側にＸ線源１１２を備える結合リング型トポロジーで配置され、Ｘ線検出器ユニット１０３はＰＣＤリングの外側に位置する。

図１は、実施形態に係るＣＴスキャナの構成例を示すブロック図である。図１は、ＣＴスキャナにおいて、所定の第３世代の配列で配置された検出器ユニット１０３と組み合わせて、所定の第４世代の配列でＰＣＤを配置する実施形態を示す。この図は、寝台１１６に横になっているスキャン対象である被検体ＯＢＪ、Ｘ線源１１２、コリメータ又はフィルタ１１４、Ｘ線検出器１０３、及び光子計数検出器ＰＣＤ１〜ＰＣＤＮの相対位置を示す。これらのＰＣＤは、前面が被検体ＯＢＪに対向しており、背面は被検体ＯＢＪに対向していない。被検体ＯＢＪを透過するＸ線は、ＰＣＤ（の前面）において検出されるか、疎に配置されたＰＣＤ間の空間を通過して、Ｘ線検出器１０３に密に実装されたエネルギー積分検出器に検出されるかのいずれかである。

図１は、Ｘ線投影データを収集し、記憶し、処理し、提供するための回路及びハードウェアも示す。これらの回路及びハードウェアは、処理部１７０、ネットワーク制御部１７４、メモリ１７８、及びデータ収集システム１７６を備える。

一実施形態においては、Ｘ線源１１２及びコリメータ又はフィルタ１１４は、ガントリ１４０に対して回転可能に取り付けられた回転部材１１０に固定して取り付けられる。Ｘ線検出器１０３も同様に、ガントリ１４０に対して回転可能に取り付けられた回転部材１３０に固定して取り付けられる。ＰＣＤ１〜ＰＣＤＮは、ガントリ１４０に固定して取り付けられた円形部材１２０に固定して取り付けられる。ガントリ１４０には、ＣＴスキャナの多くの構成要素が収容される。

ＣＴスキャナのガントリは、開放されたボア（開口部）１１５を備えており、寝台１１６上に載せられた被検体ＯＢＪを、Ｘ線源から出射され各ＰＣＤ及び検出器ユニット１０３まで進むＸ線の投影面に位置付けることができる。「投影面」とは、Ｘ線がＸ線源１１２から各ＰＣＤ及び検出器ユニット１０３を含む検出器まで透過するボリュームである。「被検体空間」とは、投影面とガントリの開放されたボア１１５との交差領域である。「画像空間」は、Ｘ線源１１２がガントリのボア１１５の周囲を回転する際の、Ｘ線源１１２のすべての投影角度に対応する投影面の和集合である。

被検体ＯＢＪが被検体空間に置かれ、Ｘ線源が一連の投影角度を通過しながら回転して、ＣＴスキャナが各投影角度で被検体ＯＢＪを通るＸ線の透過や減衰の投影データを収集することで、スキャンが行われる。

一般には、各光子計数検出器ＰＣＤ１〜ＰＣＤＮは、所定数のエネルギービンのそれぞれについて光子の数を出力する。図１に示す実施形態では、第４世代の配列で配置された光子計数検出器ＰＣＤ１〜ＰＣＤＮに加え、エネルギー積分検出器が従来の第３世代の配列で配置された検出器ユニット１０３を備える。検出器ユニット１０３の検出器素子は、検出器ユニットの表面に、光子計数検出器よりも高い密度で配置できる。

一実施形態において、光子計数検出器は、例えば円形のような所定の配列で被検体ＯＢＪの周囲に疎に配置される。例えば、光子計数検出器ＰＣＤ１〜ＰＣＤＮは、ガントリ内の所定の第２の円形部材１２０に固定して配置される。一実施形態において、光子計数検出器ＰＣＤ１〜ＰＣＤＮは、円形部材１２０上の等間隔の所定の位置に固定して配置される。別の実施形態においては、光子計数検出器ＰＣＤ１〜ＰＣＤＮは、円形部材１２０上の間隔が異なる所定の位置に固定して配置される。円形部材１２０は、データ収集中には被検体ＯＢＪに対して固定された状態のまま回転しない。

光子計数検出器ＰＣＤ１〜ＰＣＤＮは、被検体ＯＢＪに対して固定されるが、Ｘ線源１１２、コリメータ１１４（例えば、線量補償フィルタ）、及び検出器ユニット１０３は、被検体ＯＢＪの周囲を回転する。一実施形態においては、Ｘ線源１１２は、疎に配置された光子計数検出器ＰＣＤ１〜ＰＣＤＮの外側で被検体ＯＢＪの周囲を回転しながら、被検体ＯＢＪに向けて所定の線源ファンビーム角度θＡでＸ線照射を行う。更に、検出器ユニット１０３は、被検体ＯＢＪを挟んでＸ線源１１２と１８０度反対の位置に設置され、光子計数検出器ＰＣＤ１〜ＰＣＤＮが所定の疎な配置に固定された円形部材１２０の外側を回転する。

一実施形態においては、Ｘ線源１１２が、被検体ＯＢＪに対して相対的にらせん状の軌道で動くことを選択可能にする。その態様では、回転部材１１０によりＸ線源１１２及び検出器ユニット１０３を回転部材１１０の回転面で回転させながら、寝台１１６は被検体ＯＢＪを回転面と直交する所定の方向に直線的に移動させる。

被検体ＯＢＪの周囲での回転部材１１０の動きは、動作制御システムによって制御される。動作制御システムは、データ収集システムと統合されてもよく、又は、分離され、回転部材１１０の角度位置及び寝台１１６の直線上の位置についての一方的な情報を提供してもよい。動作制御システムは、位置エンコーダを備え、回転部材１１０及び寝台１１６の位置を制御するフィードバックを行うようにしてよい。動作制御システムは、開ループ系であっても、閉ループ系であっても、開ループ系と閉ループ系の組み合わせであってもよい。動作制御システムは、回転部材１１０の位置及び寝台１１６の位置についてのフィードバックに直線エンコーダ及び回転エンコーダを使用できる。動作制御システムは、回転部材１１０の動き及び寝台１１６の動きの駆動にアクチュエータを使用できる。これらのポジショナ及びアクチュエータは、例えば、ステッパモータ、ＤＣモータ、ウォーム駆動、ベルト駆動、及び当業者既知のその他のアクチュエータである。

ＣＴスキャナは、光子計数検出器及び検出器ユニット１０３からの投影測定結果を、データ収集システム１７６、処理部１７０、メモリ１７８、ネットワーク制御部１７４へ送信するデータチャネルを備える。データ収集システム１７６は、検出器からの投影データの収集、デジタル化、及び経路指定を制御する。データ収集システム１７６は、更に、環状の回転フレーム１１０及び１３０の回転を制御するＸ線撮影制御回路を備える。一実施形態において、データ収集システム１７６は、更に、寝台１１６の動き、Ｘ線源１１２の作動、及びＸ線検出器１０３の作動を制御する。データ収集システム１７６は、集中型のシステムであってよいが、分散型のシステムであってもよい。一実施形態においては、データ収集システム１７６は処理部１７０と一体化されていてもよい。処理部１７０は、投影データからの画像再構成、投影データの再構成前処理、及び画像データの再構成後処理等の機能を実行する。

投影データの再構成前処理には、較正、検出器非線形性、極性効果、ノイズバランシング、及び物質弁別についての補正を含めることができる。

再構成後処理には、必要に応じて、画像のフィルタリング及び平滑化、ボリュームレンダリング処理、及び画像の差分処理を含めることができる。画像再構成ステップは、フィルタ補正逆投影、反復的画像再構成法、又は確率的画像再構成法を用いて行うことができる。処理部１７０及びデータ収集システム１７６のいずれもが、例えば、投影データ、再構成した画像、較正のデータ及びパラメータ、及びコンピュータプログラムを記憶させるためにメモリ１７８を使用できる。

処理部１７０は、離散論理ゲートとして、特定用途向け集積回路（Application Specific Integrated Circuit：ＡＳＩＣ）、フィールドプログラマブルゲートアレイ（Field Programmable Gate Array：ＦＰＧＡ）、又はその他の複合プログラマブル論理デバイス（Complex Programmable Logic Device：ＣＰＬＤ）としての実装が可能な中央演算処理装置（Central Processing Unit：ＣＰＵ）を備えることができる。ＦＰＧＡ又はＣＰＬＤとしての実施形態は、超高速集積回路設計用ハードウェア記述言語（VHSIC（Very High Speed Integrated Circuit) Hardware Description Language：ＶＨＤＬ）、Ｖｅｒｉｌｏｇ、又は他のどのようなハードウェア記述言語でプログラムされてもよく、そのプログラムコードはＦＰＧＡ又はＣＰＬＤ内部の電子メモリに直接記憶されてもよいし、別の電子メモリに記憶されてもよい。更に、メモリは、ＲＯＭ（Read Only Memory）、ＥＰＲＯＭ（Erasable Programmable Read Only Memory）、ＥＥＰＲＯＭ（Electrically Erasable Programmable Read Only Memory）、フラッシュメモリ等のように不揮発性であってよい。メモリは、スタティックＲＡＭ（Random Access Memory）、ダイナミックＲＡＭ等のように揮発性とすることもでき、その場合、電子メモリだけでなくＦＰＧＡ又はＣＰＬＤとメモリとの間の連携を管理するマイクロコントローラ、マイクロプロセッサ等の処理部を設けてもよい。

また、再構成処理部におけるＣＰＵは、本明細書に記載の機能を実施するコンピュータで読み取り可能な命令の集合を含むコンピュータプログラムを実行してもよく、そのプログラムは、上述の非一時的な電子メモリ及び／又はハードディスクドライブ、ＣＤ（Compact Disc）、ＤＶＤ（Digital Versatile Disc）、フラッシュドライブ、又は他の任意の既知の記憶媒体に記憶されている。更に、それらのコンピュータで読み取り可能な命令は、ユーティリティアプリケーション、バックグラウンドデーモン、オペレーティングシステムの一部、又はそれらの組み合わせとして提供され、米国インテル社のＸｅｎｏｎプロセッサ又は米国ＡＭＤ社のＯｐｔｅｒｏｎプロセッサ等のプロセッサ、及びＭｉｃｒｏｓｏｆｔ ＶＩＳＴＡ、ＵＮＩＸ（登録商標）、Ｓｏｌａｒｉｓ、ＬＩＮＵＸ（登録商標）、Ａｐｐｌｅ ＭＡＣ−ＯＳ、及び当業者既知のその他のオペレーティングシステム等のオペレーティングシステムと連携して実行される。更に、ＣＰＵは、命令を実行するために並行して協調して動作する複数のプロセッサとして実装することもできる。

一実施形態においては、再構成した画像は、ディスプレイに表示することができる。ディスプレイは、ＬＣＤ（Liquid Crystal Display）、ＣＲＴ（Cathode-Ray Tube）ディスプレイ、プラズマディスプレイ、ＯＬＥＤ（Organic Light Emitting Diode）、ＬＥＤ（Light Emitting Diode）、又は当業者が知る他の任意のディスプレイでよい。

メモリ１７８は、ハードディスクドライブ、ＣＤ−ＲＯＭ（Compact Disc Read Only Memory）ドライブ、ＤＶＤドライブ、フラッシュドライブ、ＲＡＭ、ＲＯＭ、又は当業者が知る他の任意の電子記憶装置でよい。

米国インテル社のネットワークインターフェースカードであるＩｎｔｅｌ Ｅｔｈｅｒｎｅｔ（登録商標） ＰＲＯ等のネットワーク制御部１７４により、ＣＴスキャナの各部間のインターフェースをとることができる。また、ネットワーク制御部１７４により、外部ネットワークとのインターフェースをとることができる。当然のことながら、外部ネットワークは、インターネット等の公共ネットワーク、ＬＡＮ又はＷＡＮ等の専用ネットワーク、又はその組み合わせであってよく、ＰＳＴＮ又はＩＳＤＮのサブネットワークを含んでもよい。また、外部ネットワークは、イーサネット（Ethernet、登録商標）ネットワークのように有線であってもよいし、ＥＤＧＥ（Enhanced Data rates for Global Evolution）、３Ｇ、及び４Ｇ等の無線携帯電話通信システムを含む携帯電話通信ネットワークのように無線であってもよい。更に、無線ネットワークは、ＷｉＦｉ、Ｂｌｕｅｔｏｏｔｈ（登録商標）、又は知られている他の任意の無線通信形態であってよい。

一実施形態においては、Ｘ線源１１２は、単一のエネルギー源とする選択ができる。別の一実施形態においては、Ｘ線源１１２は、所定の高レベルエネルギー及び所定の低レベルエネルギーでＸ線を曝射する管電圧切り替え機能を実行するように構成される。更に、別の実施形態においては、Ｘ線源１１２は、広範囲なＸ線エネルギースペクトルを曝射する単一の線源である。更に、他の実施形態においては、Ｘ線源１１２は、複数のＸ線エミッタを備え、各エミッタは空間的にもスペクトル的にも相異なる。

検出器ユニット１０３は、光電子増倍管又はアバランシェフォトダイオードを用いたシンチレータ素子等のエネルギー積分検出器を使用して、Ｘ線照射とシンチレータ素子とが相互作用して起こるシンチレーション現象の結果生じるシンチレーション光子を検出することができる。シンチレータ素子は、結晶構造のもの（例えば、ＮａＩ（Ｔｌ）、ＣｓＩ（Ｔｌ）、ＣｓＩ（Ｎａ）、ＣｓＩ（純粋）、ＣｓＦ、ＫＩ（Ｔｌ）、ＬｉＩ（Ｅｕ）、ＢａＦ_２、ＣａＦ_２（Ｅｕ）、ＺｎＳ（Ａｇ）、ＣａＷＯ_４、ＣｄＷＯ_４、ＹＡＧ（Ｃｅ）、Ｙ_３Ａｌ_５Ｏ_１２（Ｃｅ）、ＧＳＯ、ＬＳＯ、ＬａＣｌ_３（Ｃｅ）、ＬａＢｒ_３（Ｃｅ）、ＬＹＳＯ、ＢＧＯ、Ｃ_１４Ｈ_１０、Ｃ_１４Ｈ_１２、及びＣ_１０Ｈ_８）、有機液体（例えば、ｐ−テルフェニル（Ｃ_１８Ｈ_１４）、ＰＢＤ（Ｃ_２０Ｈ_１４Ｎ_２Ｏ）、ブチルＰＢＤ（Ｃ_２４Ｈ_２２Ｎ_２Ｏ）、又はＰＰＯ（Ｃ_１５Ｈ_１１ＮＯ）等の蛍光体を含む有機溶液）、可塑性物質（例えば、固体ポリマーマトリックス中に懸濁した蛍光体を含む）、又は既存の他の任意のシンチレータであってよい。

ＰＣＤは、テルル化カドミウム（ＣｄＴｅ）、テルル化カドミウム亜鉛（ＣＺＴ）、シリコン（Ｓｉ）、ヨウ化水銀（ＨｇＩ_２）、及びヒ化ガリウム（ＧａＡｓ）等の半導体をベースにした直接的なＸ線照射検出器であってよい。一般的に、半導体ベースの直接的なＸ線検出器は、シンチレータ検出器等の間接的な検出器より時間反応がはるかに速い。直接的な検出器では反応速度が速いため、Ｘ線検出事象を個別に解決することができる。しかしながら、Ｘ線の臨床上の適用例で通常使用される高エネルギーＸ線束では、Ｘ線検出事象のパイルアップが起こる。検出されたＸ線のエネルギーは、直接的な検出器で生成された信号と正比例するので、検出事象はエネルギービンに整理することができ、スペクトルＣＴ用にスペクトル的に分解されたＸ線データが生成される。すなわち、取得部としての処理部１７０は、Ｘ線の光子の入射に応じて信号を出力する光子計数検出器の出力信号を用いて、被検体を透過したＸ線の光子の計数値を、複数のエネルギー帯域ごとに取得する。

スペクトルＣＴ投影データの取得後、処理部１７０を利用したスペクトルＣＴイメージングシステムは、スペクトルＣＴ投影データに対して物質弁別を行い、第１の物質及び第２の物質（例えば、骨等の高Ｚ物質及び水等の低Ｚ物質）にそれぞれ対応する投影長Ｌ_１及びＬ_２を求める。この物質弁別は、ＰＣＤ反応の非線形性、特性Ｘ線の逃散、Ｘ線エネルギーの関数である吸収係数の変動によるビームハードニングのために複雑である。

そのような複雑さの結果、それぞれが物質弁別の問題の解に思えるような、複数の投影長の組が存在する可能性がある（例えば、図３に示す投影長のコスト関数における多数の極小点）。しかしながら、投影長の１組のみが、物質弁別の問題の本当の解である（例えば、図３に示す投影長コスト関数の大域的な最小点）。物質弁別の問題と検出器応答の問題との結合問題に対する解は存在しないが、このような結合問題の解は、勾配降下法あるいは関連する探索法を用いる反復的方法で、適切に選択された投影長Ｌ_１及びＬ_２の初期推定から開始すれば求められる。

図２及び図４に、物質弁別の問題を解くための２つの反復的方法を示す。どちらの方法も、投影長Ｌ_１及びＬ_２の初期推定が適切に選択されるかどうかに左右される。図２に示す方法は、コスト関数手法を利用したものである。図４に示す方法は、物質弁別と検出器応答の補正が別々のステップで行われる反復的手法を利用したものである。

検出器応答は、一般的に、少なくとも極性効果が原因でＰＣＤ毎及び投影角度毎に異なる。したがって、較正係数は、一般的には、イメージングスキャンの投影角度毎に行われる較正測定によって得られる。

図２は、実施形態に係るＣＴスキャナによる処理手順の一例を示すフローチャートである。図２は、コスト関数φ（Ｌ_１，Ｌ_２）の最小化により反復的に最適投影長Ｌ^{Ｆｉｎａｌ}＝（Ｌ_１，Ｌ_２）の値を求める処理２００を示す。このコスト関数は、測定上の投影データＮ′_ｍを、前述の検出器モデルから得られた対応する計算値Ｎ_ｍと組み合わせたものである。

いくつかの異なるコスト関数φ（Ｌ_１，Ｌ_２）が考えられる。すなわち、推定部としての処理部１７０は、第１決定部により決定された複数の物質ごとの初期投影長と、減衰モデルと、光子計数検出器のＸ線に対する応答特性を示す検出器モデルとを用いて、複数のエネルギー帯域ごとの計数値を推定する。一実施形態では、コスト関数は測定上の計数Ｎ′_ｍと計算上の計数Ｎ_ｍとの差の最小二乗、すなわち、上記の式（８）である。

一実施形態では、コスト関数は測定上の計数Ｎ′_ｍと計算上の計数Ｎ_ｍとの差の重み付き最小二乗、すなわち、上記の式（９）であり、σ_ｍは、検出器のｍ番目のエネルギービンにおける測定の不確実性の指標である。

一実施形態では、コスト関数はポアソン尤度関数、すなわち、上記の式（１０）である。

図２では、処理２００は、Ｌ^（０）＝（Ｌ_１ ^（０），Ｌ_２ ^（０））の初期推定として乱数値が選ばれるステップから始まる。例えば、ステップＳ２０６として、初期投影長推定Ｌ^（０）が選択される。続いて、ステップＳ２１０において、ループ変数ｎがインクリメントされる。

ステップＳ２１０に続いて、処理２００は、ステップＳ２２０に進み、新しいサンプル点Ｌ′が、現在の投影長の組Ｌ^{（ｎ−１）}＝（Ｌ_１ ^{（ｎ−１）}，Ｌ_２ ^{（ｎ−１）}）の周囲にあるサンプル空間から無作為に選ばれる。

ステップＳ２３０に進むと、処理２００では、コスト関数φ（Ｌ^{（ｎ−１）}）の値とφ（Ｌ′）の値のどちらが小さいかが調べられる。ステップＳ２４０及びＳ２５０において、コスト関数の小さい方の値に対応する引数が、次のループ反復のための次の投影長の組Ｌ^（ｎ）＝（Ｌ_１ ^（ｎ），Ｌ_２ ^（ｎ））に指定される。

処理２００のステップＳ２６０では、ループ停止基準が満たされたかどうかを判定する。異なる停止基準を使用することもできるが、図２には、ループ反復の最大数ｎ_ｍａｘに到達した場合、あるいはコスト関数が所定の閾値εを下回った場合に、ループが停止する実施形態を示す。処理２００では、停止基準が満たされた場合、ステップＳ２６０においてループが終了し、ステップＳ２７０として現在の投影長Ｌ^（ｎ）＝（Ｌ_１ ^（ｎ），Ｌ_２ ^（ｎ））が最終的な投影長Ｌ^{（Ｆｉｎａｌ）}として通知される。すなわち、第２決定部としての処理部１７０は、推定した計数値と、取得した計数値との差に基づくコスト関数を用いて、反復的最適化方法により、被検体内における複数の物質それぞれの投影長を決定する。停止基準が満たされない場合、ステップＳ２６０からステップＳ２１０に戻り、ループ処理が継続される。

図３は、コスト関数と投影長との関係の一例を示す図である。図３は、極小点と大域的な最小点を有する最小二乗によるコスト関数の例を示す。コスト関数の大域的な最小点を見つける方法では、反復的探索から得られる投影長が、コスト関数の極小点ではなく、大域的な最小点に対応することを確実にする必要がある。一般的に、パイルアップに伴う非線形性及び多色Ｘ線ビームに起因するビームハードニングのため、すべての場合に収束を保証するようなアルゴリズムはない。しかしながら、投影長に適切に選択された初期値を使用することで、求めた投影長が大域的な最小点に対応する可能性が高くなる。

図４は、実施形態に係るＣＴスキャナによる処理手順の一例を示すフローチャートである。図４は、最適投影長Ｌ^{Ｆｉｎａｌ}＝（Ｌ_１，Ｌ_２）の値を反復的に求める処理４００を示す。処理４００と処理２００との相違点は、処理２００では物質弁別及び検出器応答補正がコスト関数に内在しており単一のステップに対応するのに対し、処理４００では、物質弁別が検出器応答補正ステップ（すなわち、ステップＳ４３０）と異なるステップ（すなわち、ステップＳ４５０）で行われることである。

図４に示すように、処理４００は、反復ループの最初のステップＳ４２０に投影長Ｌ_１及びＬ_２の初期推定が入力されて開始される。ステップＳ４２０からステップＳ４６０を含む大ループの１回目の反復の後、投影長Ｌ_１及びＬ_２が、ステップＳ４５０における物質弁別から大反復ループの最初のステップＳ４２０にフィードバックされる。ステップＳ４５０の物質弁別で計算される投影長を使用する２回目以降のループ反復とは異なり、１回目の反復では、ステップＳ４１０による投影長の初期推定が使用される。

投影長推定の受け取り後、大ループは、ステップＳ４２０で、上記の式（７）によるＸ線束率及び前述の非線形スペクトル項Ｓ_{Ｎｏｎｌｉｎ．}（Ｅ）を計算することから開始される。図４に示す実施形態では、非線形スペクトルには下記の式（１１）で与えられる１次のパイルアップしか含まれない。

ここで、Ｓ_ｉｎ（Ｅ）は、下記の式（１２）である。

続いて、処理４００のステップＳ４３０で、パイルアップを補正した補正後の検出器スペクトルが下記の式（１３）により計算される。

ここで、Ｓ_Ｒａｗ（Ｅ）は、検出器応答の補正前の未処理の測定スペクトルである。補正後のエネルギーの計数は、下記の式（１４）で与えられる。Ｔは、積分時間である。

一実施形態では、ＰＣＤのｍ番目のエネルギービンについての補正後計数は、下記の式（１５）で与えられる。

ここで、Ｎ_ｍ ^{Ｃｏｒｒ．}は、ＰＣＤのｍ番目のエネルギービンについての補正後の計数値であり、Ｎ_ｍ ^Ｒａｗは、検出器で記録された未処理（生の）の計数値であり、Ｎ_ｍ ^{Ｎｏｎｌｉｎ．}は、非線形な検出器応答からの計算上の計数である。この非線形な計数値Ｎ_ｍ ^{Ｎｏｎｌｉｎ．}は、下記の式（１６）により計算される。

一部の実施形態では、非線形スペクトル補正に１次のパイルアップしか含まないが、非線形スペクトル補正に２次以上のパイルアップ項を含む実施形態もある。

処理４００は、次にステップＳ４４０に進む。ステップＳ４４０では、物質弁別の準備として、検出器からの計数を高エネルギー成分及び低エネルギー成分に分けることによりノイズバランシングが行われる。計数を高エネルギー成分及び低エネルギー成分に分配するノイズバランシング処理は、その開示内容全体が参照として本明細書に組み入れられる特許文献２に記載されている。高エネルギー成分計数及び低エネルギー成分計数のいずれか一方の信号対雑音比が他方のものより大幅に低い（すなわち、高エネルギー成分計数及び低エネルギー成分計数の信号対雑音比が均衡しない）ことで再構成画像の画像品質が低下する状況を回避するために、ノイズバランシングにより、高スペクトルと低スペクトル間の線引きが検出器毎に異なるようにできる。ノイズバランシングは、ＰＣＤ毎に異なってよい。

ステップＳ４４０におけるノイズバランシングは、各エネルギービンからの計数を、下記の式（１７）及び式（１８）によって高エネルギー成分及び低エネルギー成分に分ける。

ここで、下記の式（１９）及び式（２０）が成り立つ。また、ａ_ｍ ^（Ｈ）及びａ_ｍ ^（Ｌ）は、ノイズバランシング処理によって決まる。

続いて、処理４００は、ステップＳ４５０に進む。ステップＳ４５０では、物質弁別が行われ、投影長Ｌ_１及びＬ_２の新しい値が計算される。

最後に、ステップＳ４６０で、停止基準が満たされたかどうかの判定を行う。停止基準は、投影長Ｌ_１及びＬ_２の収束、及びループ反復の最大数に到達したかを基にしてよい。

また、物質弁別のステップＳ４５０は、それ自体に停止基準を設けた反復的なステップ（図５に示すようなステップ）であってよい。このような反復的なステップは、ステップＳ４２０からステップＳ４６０を含む大ループと対比して、小ループと呼ぶことが多い。ステップＳ４５０に対応する小ループの停止基準は、大ループの現在の反復回数を基にしてよい。例えば、一実施形態では、大ループの反復変数が小さい場合、小ループの停止基準は、収束閾値を比較的緩くして反復の最大数を下げ、推定投影長が大ループの反復毎にまだ大幅に変化している間は小ループが過度に最適化されないようにする。

一実施形態では、ステップＳ４５０は、その開示内容全体が参照として本明細書に組み入れられる特許文献３に記載されているように、図５に示す方法で行われる。図５は、実施形態に係るＣＴスキャナによる処理手順の一例を示すフローチャートである。高エネルギー成分及び低エネルギー成分についての検出器による計数と同様に、高スペクトル及び低スペクトルは、下記の式（２１）及び式（２２）で与えられる。

ここで、Ｓ_Ｈ（Ｅ）及びＳ_Ｌ（Ｅ）は被検体ＯＢＪ不在（すなわち、被検体ＯＢＪが空気である）時の高エネルギースペクトル及び低エネルギースペクトルであり、Ｓ_Ｈ（Ｅ）及びＳ_Ｌ（Ｅ）は、下記の式（２３）のように正規化されている。

更に、検出器の計数の自然対数をとることで、下記の式（２４）及び式（２５）のような対数投影データが得られる。

一実施形態では、図５に示すように、Ｌ_１及びＬ_２は、減衰係数μ_１（Ｅ）及びμ_２（Ｅ）の平均に対する変動を摂動として扱うことにより摂動論を利用して求めることができる。まず、高エネルギースペクトル及び低エネルギースペクトルの平均減衰は、下記の式（２６）で与えられる。

ここで、平均に対する変動は、下記の式（２７）で与えられる。

したがって、対数投影データは、下記の式（２８）及び式（２９）で表すことができる。

これらの数式表現を簡略化すると、対数投影データは、下記の式（３０）及び式（３１）のようになる。

ここで、ビームハードニング摂動は、下記の式（３２）である。

ステップＳ４５０の最初のステップＳ５１２では、反復変数がｎ＝０に初期化され、投影長Ｌ_１及びＬ_２の値も初期化される。一実施形態では、投影長の初期値は、ステップＳ４２０の検出器応答の補正計算で使用されるものと同じ値である。別の実施形態では、投影長の初期値は、行列方程式である下記の式（３３）を解いて、下記の式（３４）を求めることで計算されるゼロ次の摂動値である。

ここで、Ｄは、行列式である下記の式（３５）である。

ステップＳ４５０の２番目のステップＳ５１４である、反復ループの最初のステップでは、上記の式（３２）の中のｎ次の摂動を使用してビームハードニング摂動値が更新される。

ステップＳ４５０の３番目のステップＳ５１６は、行列方程式である下記の式（３６）を解いて、下記の式（３７）を求めることで（ｎ＋１）番目の摂動値を求めることにより、Ｌ_１及びＬ_２の値を更新するステップである。

ステップＳ５１６の後、ステップＳ４５０のステップＳ５２０において、停止基準が満たされたかどうかを判定する。一実施形態では、値Ｌ_１及び値Ｌ_２が所定の収束基準を満たした場合に停止基準が満たされたとする。所定の収束基準は、Ｌ_１及びＬ_２のそれぞれの現在の値と直前の値との差が所定の閾値より小さいかどうか等である。停止基準は、反復の最大数に到達したかどうかで条件付けしてもよい。停止基準が満たされていなかった場合、ループ変数ｎは、ステップＳ５１８でインクリメントされ、ループはステップＳ５１４から再び開始される。

処理２００及び処理４００のいずれにおいても、最適投影長の反復的探索の効率、更に、多くの場合はその探索の結果も、投影長の初期推定次第である。適切に選択された初期推定は、適切な初期推定方法を使用して求められる。初期推定方法は、対象となっているＣＴスキャンのパラメータに合わせて調整することもできる。例えば、図６及び図７に、Ｘ線源においてＸ線束密度が高い場合と低い場合に計算されたコスト関数間の違いを示す。Ｘ線束密度によっては、投影長の初期推定を長めもしくは短めに偏らせると都合がよい場合がある。図６及び図７に示す実施形態では、物質弁別における２つの物質成分として水及び骨を示すが、当業者には、他の物質成分も使用可能であることが理解されよう。

図６及び図７は、コスト関数と投影長との関係の一例を示す図である。図６及び図７に、それぞれ、第１の物質及び第２の物質の関数としてグラフ化された最小二乗のコスト関数のグラフを示す。図６は、Ｘ線源において線束密度が低い場合のコスト関数を示し、図７は、線束密度が高い場合のコスト関数を示す。図６及び図７のいずれにおいても、最適な投影長はＬ_１＝１０ｍｍ、Ｌ_２＝４０ｍｍである。

図６及び図７に、初期投影長の推定を偏らせることに関連する係数を示す。各ＰＣＤへの入射線束が低い場合は、図６のように、コスト関数のより大きな勾配が探索アルゴリズムに利用されるように、初期投影長をより短くする方向に偏らせる必要がある。これに対して参照強度Ｉ_ｒｅｆが高い場合は、図７のように、コスト関数のより大きな勾配が探索アルゴリズムに利用されるように、更に、短い投影長においては極小点でしかない点が避けられるように、初期投影長をより長くする方向に偏らせる必要がある。

初期推定を偏らせる根拠は、最適な投影長の探索に勾配降下法が利用される場合に特に説得力がある。処理４００において最適な投影長の探索に勾配降下法が利用される場合、収束速度はコスト関数の勾配に左右される。したがって、入射線束密度が低い場合は、図６に示すように、初期推定をより低い方に偏らせることで収束速度を上げることができる。

これに対して、入射線束密度が高い場合は、図７に示すように、投影長の初期推定が高い方に偏るほど、勾配が大きくなる。更に、入射線束密度が高い場合には、短い投影長においてコスト関数の極小点でしかない点が発生している。投影長の初期推定を最適な投影長より大きくなるように偏らせることで、この望ましくない極小点を避けることができる。

要約すれば、コスト関数のより大きな勾配に対応する初期投影長推定を選ぶことで、勾配降下法は収束が速くなる。低線束密度では、投影長の推定が全体的に最適投影長より短い（すなわち、最適投影長より短い方に偏っている）と、一般的には、コスト関数の勾配はより大きくなる。一方で、高線束密度では、投影長の推定が全体的に最適投影長より長い（すなわち、最適投影長より長い方に偏っている）と、一般的には、コスト関数の勾配はより大きくなる。更に、高線束密度では、大域的最小点とは異なる極小点は、全体的に最適投影長より長い初期投影長推定であれば避けられる可能性が高くなる。すなわち、第１決定部としての処理部１７０は、線束密度が第１閾値より高い場合、初期投影長が長くなるように調整し、線束密度が第１閾値より低い第２閾値より低い場合、初期投影長が短くなるように調整する。例えば、取得部としての処理部１７０は、更に、被検体を透過したＸ線の強度を検出するエネルギー積分検出器の出力信号を取得する。そして、第１決定部としての処理部１７０は、エネルギー積分検出器の出力信号から再構成される再構成画像を第１物質成分及び第２物質成分に分解し、第１物質成分に対して順投影を行うことで第１物質成分に対応する第１初期投影長を決定し、第２物質成分に対して順投影を行うことで第２物質成分に対応する第２初期投影長を決定し、決定した第１初期投影長及び第２初期投影長を、線束密度に応じて調整する。

図８は、実施形態に係るＣＴスキャナによる処理手順の一例を示すフローチャートである。図８に、投影長の初期推定を求める別の方法を示す。この方法では、非スペクトル投影データによる再構成画像を物質成分に分解してから、物質成分の画像を順投影して投影長の初期推定を求める。

図８に、非スペクトル投影データから再構成したＣＴ画像のハンスフィールド単位（Hounsfield Unit：ＨＵ）に基づく投影長推定を導出する処理８００を示す。非限定的かつ例示的な実施形態では、この方法は、再構成画像の減衰が非常に大きい領域が骨に対応し、それより減衰が小さな領域は水に対応すると想定する。また、大きい減衰と小さい減衰の中間である減衰領域については、水と骨とが重なり合った混合領域とする。

処理８００の最初のステップＳ８１０は、非スペクトル投影データを用いて画像を再構成するものである。画像再構成の方法は、フィルタ補正逆投影法、反復的画像再構成法、及び確率的画像再構成法のいずれでもよい。画像スケールは、関数ＨＵ（ｘ，ｙ）によりＨＵで与えられる。

図８に示す非限定的な実施形態では、投影データが、第３世代の幾何形状で配置されたエネルギー積分検出器から導出されることを示す。一実施形態では、再構成画像は、任意のスペクトルＣＴ投影データから導出してよい。また、一実施形態では、各検出器のスペクトル次元及び投影角度にわたってデータの総和を求めて、スペクトルＣＴ投影データから非スペクトル投影データを導出してよい。スペクトルデータから非スペクトルデータが取得されると、この合算された非スペクトル投影データから画像を再構成できる。ただし、このスペクトル合算による実施形態では、各ＰＣＤの非線形な検出器応答について補正するために、スペクトル投影データに最初に前処理を行う必要がある場合がある。

ステップＳ８１０で減衰画像が再構成されると、処理８００は、ステップＳ８２０に進み、再構成画像の領域が第１の物質及び第２の物質に対応する物質成分に分解される。このステップにおいて、物質１及び物質２それぞれの密度を説明する空間的な関数は、下記の式（３８）及び式（３９）で与えられる。

定数ＨＵ_１及びＨＵ_２は、それぞれ、完全な第１の物質から混合物質への移行、及び混合物質から完全な第２の物質への移行を表す閾値である。

ステップＳ８２０で再構成画像が物質成分に分解されると、処理８００は、ステップＳ８３０に進む。ステップＳ８３０では、空間的な関数ｃ_１（ｘ，ｙ）及びｃ_２（ｘ，ｙ）に順投影（例えば、検出されたＸ線の軌跡に沿って線積分するラドン変換）を行うことで、投影長が求められる。すなわち、第１決定部としての処理部１７０は、再構成画像のある１点における減衰が第１閾値未満である場合、その１点における第１物質成分の減衰はその１点における減衰と等しく、その１点における第２物質成分の減衰は０であり、再構成画像の１点における減衰が第１閾値より大きい第２閾値を超える場合、その１点における第１物質成分の減衰は０であり、その１点における前記第２物質成分の減衰は１点における減衰と等しく、再構成画像の１点における減衰が第１閾値以上かつ第２閾値未満である場合、その１点における第１物質成分の減衰と第２の物質成分の減衰との重畳で得られる減衰は、その１点における減衰と等しくなるように、再構成画像を第１物質成分及び第２物質成分に物質弁別を行うことで、第１初期投影長及び第２初期投影長を決定する。

順投影問題でもある画像再構成問題は、行列方程式である下記の式（４０）として定式化される。

ここで、ｇは、被検体ＯＢＪが存在する被検体空間を透過したＸ線の投影測定値であり、Ａは、被検体空間を通るＸ線の離散線積分（すなわち、ラドン変換）を説明するシステム行列であり、ｆは、被検体ＯＢＪの画像（すなわち、システム行列等式を解くことで求められる数量）である。すなわち、第１決定部としての処理部１７０は、第１物質成分及び第２物質成分それぞれに対する順投影としてラドン変換を用いて、第１投影長及び第２投影長を決定する。画像ｆは、位置の関数である減衰のマップである。これにより、投影推定は下記の式（４１）及び式（４２）で与えられる。

ここで、ｃ_１及びｃ_２は、それぞれ、物質密度ｃ_１（ｘ，ｙ）及びｃ_２（ｘ，ｙ）を表す列ベクトル形式である。

ステップＳ８３０からの投影長の出力が、処理２００のステップＳ２０６において、又は処理４００のステップＳ４１０において初期推定として使用できる。

図９Ａ及び図９Ｂに、ファントムの再構成画像の一例を示す。図９Ａに、処理８００においてエネルギー積分検出器による投影データ及びＨＵ物質弁別の等式を用いた、第１の物質ｃ_１（ｘ，ｙ）についての推定画像を示す。図９Ｂに、スペクトル投影データを用いた、第１の物質の実際のファントムについての再構成画像を示す。これらの図からわかるように、図９Ａ及び図９Ｂの画像は大部分が一致するが、相違点は主に遷移や境界に相当する領域に見られる。

図１０Ａ及び図１０Ｂに、それぞれ図９Ａ及び図９Ｂに対応するサイノグラムを示す。図１０Ａは、第１の物質の投影長Ｌ_１の推定を示す。図１０Ｂは、図９Ｂに示す画像に対応する、第１の物質の実際のファントムにおける投影長のサイノグラムを示す。図９Ａと図９Ｂとの間の類似性と一致して、図１０Ａには、図１０Ｂとの強い類似性が見られる。このことから、処理８００により、適切に選択された投影長の初期推定が得られることが示されている。これら２つのサイノグラムの底辺に沿って表示されているのは、これらサイノグラムから得られ、比較のため図１１に示すラインアウトに対応する線である。

図１１は、図１０Ａ及び１０Ｂのサイノグラム上に示す線に対応する推定投影長及び実投影長のラインアウト間の比較である。図１１は、縦軸に投影長の値、横軸に画素番号を示す。このグラフでも、推定投影長と実投影長との強い相関が見られ、処理８００により、図２及び４の反復的方法に対して適切な初期推定が得られることが示されている。

図１２Ａ及び図１２Ｂに、ファントムの再構成画像を示す。図１２Ａに、処理８００においてエネルギー積分検出器による投影データ及びＨＵ物質弁別等式を用いた、第２の物質ｃ_２（ｘ，ｙ）についての推定画像を示す。図１２Ｂに、スペクトル投影データを用いた、第２の物質の実際のファントムについての再構成画像を示す。これらの図からわかるように、図１２Ａ及び図１２Ｂの画像は一致度が極めて高い。

図１３Ａ及び図１３Ｂに、それぞれ図１２Ａ及び図１２Ｂに対応するサイノグラムを示す。図１３Ａは、第２の物質の投影長Ｌ_２の推定を示す。図１３Ｂは、図１２Ｂに示す画像に対応する、第２の物質の実際のファントムにおける投影長のサイノグラムを示す。図１２Ａと図１２Ｂとの間の類似性と同様に、図１３Ａには、図１３Ｂとの強い類似性が見られる。このことから、処理８００により、適切に選択された投影長の初期推定が得られることが示されている。これら２つのサイノグラムの底辺に沿って表示されているのは、これらサイノグラムから得られ、比較のため図１４に示すラインアウトに対応する線である。

図１４は、図１３Ａ及び１３Ｂのサイノグラム上に示す線に対応する推定投影長及び実投影長のラインアウト間の比較である。図１４は、縦軸に投影長の値、横軸に画素番号を示す。このグラフでも、推定投影長と実投影長との強い相関が見られ、処理８００により、図２及び４の反復的方法に対して適切な初期推定が得られることが示されている。

図１５は、投影長の初期推定を求める処理１５００を示す。この処理１５００は、所定のサンプル空間を受け付けることで開始される。例えば、所定のサンプル空間は、処理４００を用いて求められ初期推定の所定の距離内の、投影長の一群として与えられる。また、所定のサンプル空間は他の方法でも求められる。例えば、サンプル空間は、結果的に検出された吸収になる投影長の組み合わせすべてを含むようにしてもよい。

処理１５００のステップＳ１５１０では、任意の投影長対の一群が、所定のサンプル空間から無作為に選択される。

続いて、ステップＳ１５２０で、その投影長対の一群における各投影長対Ｌ^（ｉ）について、コスト関数φ（Ｌ^（ｉ））が計算される。

続いて、ステップＳ１５３０で、最も低い値のコスト関数に対応する投影長が、投影長推定として選択される。この投影長推定が、処理１５００の出力である。

処理１５００は、投影長の初期推定を求めるのに単独で使用できる。処理１５００は、偏りを持たせる方法と処理８００のいずれかと組み合わせて使用することや、偏りを持たせる方法及び処理８００の両方と組み合わせて使用することが可能である。例えば、一実施形態においては、処理８００を使用して第１の初期推定を求めることもでき、その第１の初期推定を使用して処理１５００における所定のサンプル空間を定義できる。次に、処理１５００では、処理８００による所定のサンプル空間を使用して第２の初期推定を求めることができる。更に、この第２の初期推定を、偏りを持たせる方法によって高くあるいは低く偏らせて、第３の初期推定を求めることができる。これら３つの初期推定のいずれか、又はこれら３つの初期推定の重み付き平均を、反復的方法で初期推定として使用して、最適投影長を求めることができる。あるいは、処理８００を使用して求めた第１の初期推定を、偏りを持たせる方法によって高くあるいは低く偏らせて、第４の推定を求めることができる。当業者には理解されるであろうが、初期推定を導出するこれらの初期推定方法の組み合わせ方及びその一部の組み合わせ方には多くの態様がある。

特定の実施形態を説明したが、これらの実施形態は、単なる事例として提示したものであり、本開示の教示を限定するものではない。実際、本明細書で説明した新規の方法、装置、及びシステムは、様々な別の態様で具体化が可能である。更には、本開示の要旨を逸脱することなく、本明細書で説明した方法、装置、及びシステムの態様において様々な省略、置換、及び変更が可能である。

以上説明した少なくともひとつの実施形態によれば、複数の物質それぞれの投影長を正確に求めることができる。

本発明のいくつかの実施形態を説明したが、これらの実施形態は、例として提示したものであり、発明の範囲を限定することは意図していない。これら実施形態は、その他の様々な形態で実施されることが可能であり、発明の要旨を逸脱しない範囲で、種々の省略、置き換え、変更を行うことができる。これら実施形態やその変形は、発明の範囲や要旨に含まれると同様に、特許請求の範囲に記載された発明とその均等の範囲に含まれるものである。

１７０ 処理部

Claims (14)

1. Ｘ線の光子の入射に応じて信号を出力する光子計数検出器の出力信号を用いて、被検体を透過したＸ線の光子の計数値を、複数のエネルギー帯域ごとに取得する取得部と、
   複数の物質それぞれの投影長をパラメータとして含むＸ線の減衰モデルに入力する初期投影長を、前記複数の物質ごとに決定する第１決定部と、
   前記第１決定部により決定された前記複数の物質ごとの初期投影長と、前記減衰モデルと、前記光子計数検出器のＸ線に対する応答特性を示す検出器モデルとを用いて、前記複数のエネルギー帯域ごとの計数値を推定する推定部と、
   前記推定部が推定した計数値と、前記取得部が取得した計数値との差に基づくコスト関数を用いて、反復的最適化方法により、前記被検体内における前記複数の物質それぞれの投影長を決定する第２決定部と、
   を備える、医用情報処理装置。
2. 前記取得部は、更に、前記被検体を透過したＸ線の強度を検出するエネルギー積分検出器の出力信号を取得し、
   前記第１決定部は、前記エネルギー積分検出器の出力信号から再構成される再構成画像に基づいて、前記初期投影長を決定する、
   請求項１に記載の医用情報処理装置。
3. 前記第１決定部は、前記再構成画像を第１物質成分及び第２物質成分に物質弁別を行い、前記第１物質成分に対して順投影を行うことで第１物質成分に対応する第１初期投影長を決定し、前記第２物質成分に対して順投影を行うことで第２物質成分に対応する第２初期投影長を決定する、
   請求項２に記載の医用情報処理装置。
4. 前記第１決定部は、
   前記再構成画像の１点における減衰が第１閾値未満である場合、前記１点における前記第１物質成分の減衰は前記１点における減衰と等しく、前記１点における前記第２物質成分の減衰は０であり、
   前記再構成画像の前記１点における前記減衰が前記第１閾値より大きい第２閾値以上である場合、前記１点における前記第１物質成分の減衰は０であり、前記１点における前記第２物質成分の減衰は前記１点における前記減衰と等しく、
   前記再構成画像の前記１点における前記減衰が前記第１閾値以上かつ前記第２閾値未満である場合、前記１点における前記第１物質成分の減衰と前記第２物質成分の減衰との重畳で得られる減衰は、前記１点における前記減衰と等しくなるように、
   前記再構成画像を前記第１物質成分及び前記第２物質成分に物質弁別を行うことで、前記第１初期投影長及び前記第２初期投影長を決定する、
   請求項３に記載の医用情報処理装置。
5. 前記第１決定部は、探索領域内の前記投影長に対応するサンプル点を無作為に複数選択し、選択した複数のサンプル点のうち前記コスト関数が最小値となるサンプル点の投影長を前記初期投影長として決定する、
   請求項１に記載の医用情報処理装置。
6. 前記取得部は、更に、前記被検体を透過したＸ線の強度を検出するエネルギー積分検出器の出力信号を取得し、
   前記第１決定部は、前記被検体を透過したＸ線の強度を検出するエネルギー積分検出器の出力信号から再構成される再構成画像を第１物質成分及び第２物質成分に分解し、前記第１物質成分に対応する第１投影長と、前記第２物質成分に対応する第２投影長との組み合わせのパターンの中から複数の組み合わせを無作為に選択し、選択した組み合わせのうち前記コスト関数が最小値となる組み合わせの投影長を前記初期投影長として決定する、
   請求項５に記載の医用情報処理装置。
7. 前記第１決定部は、前記第１物質成分及び前記第２物質成分それぞれに対する順投影としてラドン変換を用いて、前記第１投影長及び前記第２投影長を決定する、
   請求項６に記載の医用情報処理装置。
8. 前記第１決定部は、前記光子計数検出器に入射されるＸ線の線束密度に応じて、前記初期投影長を調整する、
   請求項１に記載の医用情報処理装置。
9. 前記第１決定部は、
   前記線束密度が第１閾値より高い場合、前記初期投影長が長くなるように調整し、
   前記線束密度が前記第１閾値より低い第２閾値より低い場合、前記初期投影長が短くなるように調整する、
   請求項８に記載の医用情報処理装置。
10. 前記取得部は、更に、前記被検体を透過したＸ線の強度を検出するエネルギー積分検出器の出力信号を取得し、
    前記第１決定部は、前記被検体を透過したＸ線の強度を検出するエネルギー積分検出器の出力信号から再構成される再構成画像を第１物質成分及び第２物質成分に分解し、前記第１物質成分に対して順投影を行うことで第１物質成分に対応する第１初期投影長を決定し、前記第２物質成分に対して順投影を行うことで第２物質成分に対応する第２初期投影長を決定し、決定した前記第１初期投影長及び前記第２初期投影長を、前記線束密度に応じて調整する、
    請求項９に記載の医用情報処理装置。
11. 前記取得部は、前記出力信号の総和を算出し、
    前記第１決定部は、前記取得部により算出された総和から再構成される再構成画像を第１物質成分及び第２物質成分に分解し、前記第１物質成分に対して順投影を行うことで第１物質成分に対応する第１初期投影長を決定し、前記第２物質成分に対して順投影を行うことで第２物質成分に対応する第２初期投影長を決定する、
    請求項１に記載の医用情報処理装置。
12. Ｘ線を発生させるＸ線源と、
    Ｘ線の光子の入射に応じて信号を出力する光子計数検出器と、
    前記光子計数検出器の出力信号を用いて、被検体を透過したＸ線の光子の計数値を、複数のエネルギー帯域ごとに取得する取得部と、
    複数の物質それぞれの投影長をパラメータとして含むＸ線の減衰モデルに入力する初期投影長を、前記複数の物質ごとに決定する第１決定部と、
    前記第１決定部により決定された前記複数の物質ごとの初期投影長と、前記減衰モデルと、前記光子計数検出器のＸ線に対する応答特性を示す検出器モデルとを用いて、前記複数のエネルギー帯域ごとの計数値を推定する推定部と、
    前記推定部が推定した計数値と、前記取得部が取得した計数値との差に基づくコスト関数を用いて、反復的最適化方法により、前記被検体内における前記複数の物質それぞれの投影長を決定する第２決定部と、
    を備える、Ｘ線ＣＴ装置。
13. Ｘ線の光子の入射に応じて信号を出力する光子計数検出器の出力信号を用いて、被検体を透過したＸ線の光子の計数値を、複数のエネルギー帯域ごとに取得し、
    複数の物質それぞれの投影長をパラメータとして含むＸ線の減衰モデルに入力する初期投影長を、前記複数の物質ごとに決定し、
    決定された前記複数の物質ごとの初期投影長と、前記減衰モデルと、前記光子計数検出器のＸ線に対する応答特性を示す検出器モデルとを用いて、前記複数のエネルギー帯域ごとの計数値を推定し、
    推定した計数値と、取得した計数値との差に基づくコスト関数を用いて、反復的最適化方法により、前記被検体内における前記複数の物質それぞれの投影長を決定する
    ことを含む、医用情報処理方法。
14. Ｘ線の光子の入射に応じて信号を出力する光子計数検出器の出力信号を用いて、被検体を透過したＸ線の光子の計数値を、複数のエネルギー帯域ごとに取得し、
    複数の物質それぞれの投影長をパラメータとして含むＸ線の減衰モデルに入力する初期投影長を、前記複数の物質ごとに決定し、
    決定された前記複数の物質ごとの初期投影長と、前記減衰モデルと、前記光子計数検出器のＸ線に対する応答特性を示す検出器モデルとを用いて、前記複数のエネルギー帯域ごとの計数値を推定し、
    推定した計数値と、取得した計数値との差に基づくコスト関数を用いて、反復的最適化方法により、前記被検体内における前記複数の物質それぞれの投影長を決定する
    処理をコンピュータに実行させる、医用情報処理プログラム。