JP2005043104A - Method for calibrating radiation position detector - Google Patents
Method for calibrating radiation position detector Download PDFInfo
- Publication number
- JP2005043104A JP2005043104A JP2003200616A JP2003200616A JP2005043104A JP 2005043104 A JP2005043104 A JP 2005043104A JP 2003200616 A JP2003200616 A JP 2003200616A JP 2003200616 A JP2003200616 A JP 2003200616A JP 2005043104 A JP2005043104 A JP 2005043104A
- Authority
- JP
- Japan
- Prior art keywords
- radiation
- scintillator
- position detector
- dimensionally
- clustering
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Images
Landscapes
- Measurement Of Radiation (AREA)
- Nuclear Medicine (AREA)
Abstract
Description
【0001】
【発明の属する技術分野】
本発明は、放射線位置検出器の校正方法に係り、特に、3次元的に配置された多数のシンチレータ素子と、それらに光結合された光センサから構成された放射線位置検出器を校正する際に用いるのに好適な、放射線の位置検出機能や放射線吸収エネルギ選別機能を持たせた放射線位置検出器の校正方法に関する。
【0002】
【従来の技術】
被検体に放射性同位元素(RI)が投与されると、消滅γ線と呼ばれる2個の光子が正反対方向に同時に放出されるイベントが発生する。陽電子放出断層装置(PET)は、この2個の光子を同時計数で検出することにより、被検体内のRIの分布像を得る装置である。
【0003】
このようなPET装置に用いられる放射線3次元位置検出器としては、複数のシンチレータ素子を、屈折率の異なる透明板を挟みながら光センサ上に積層してシンチレータブロックを構成し、光センサに到達する光の透過率を各シンチレータ素子毎に異ならせることにより、光センサにおける受光量の差に基づいて、放射線が入射して蛍光を発したシンチレータ素子を同定するようにした放射線検出器が、特許文献1に記載されている。
【0004】
又、複数のシンチレータ素子を、その中心位置が光位置検出器の受光面に平行な方向に偏倚するように光位置検出器上に積層し、光位置検出器からの出力光の空間分布の重心位置を各シンチレータセル毎に異ならせることにより、重心位置演算に基づいて、放射線が入射して蛍光を発したシンチレータ素子を同定するようにした放射線検出器が、特許文献2に記載されている。
【0005】
又、図1に示す如く、シンチレータブロック12と、光センサとして4つの光電子増倍管(PMT)14A、14B、14C、14Dから構成されるブロック検出器10の場合、シンチレータブロック12中のγ線との相互作用が起こった蛍光発光シンチレータ素子(結晶素子とも称する)の同定を、図2(簡単のため1層の例で説明する)に示すフロントエンド回路の如く、PMT14A、14B、14C、14Dの出力信号(A,B,C,D)を用いて次式(1)、(2)に示すような重心演算を行い、更に後続回路(図示せず)により図3に示す9×10領域に区分けした2次元位置マップ(位置弁別テーブルとも称する)等に得られた重心位置を投影し、続いてこの作業を繰り返すことで一シンチレータ素子に対応する統計的なシンチレータ素子対応域を作成し、更に他のシンチレータ素子についても同様に作成し、最終的に全シンチレータ素子に対する対応域を合成した2次元マップを作成することも行われている。
【0006】
X=(A+B)/(A+B+C+D) ・・・(1)
Y=(A+C)/(A+B+C+D) ・・・(2)
【0007】
但し、2次元位置マップ上での各シンチレータ素子に対応する対応域は、検出器光子数の統計的ばらつきなどによってある広がりを持つ上、PMT14A〜Dの各ゲインの違いやシンチレータブロック12の作り方によって検出器毎に異なる。従って、予め作成した位置ルックアップテーブル(位置LUTと称する)28に従って直線性補正を行う。又、A+B+C+Dはシンチレータ内の全ての発光光子数に相当するので、γ線のエネルギに比例する。エネルギの波高値も、シンチレータ素子毎の効率の違いによりシンチレータブロック12内のどのシンチレータ素子で蛍光を発したか(発光したシンチレータ素子位置)で異なるため、同じようにエネルギルックアップテーブル(エネルギLUTと称する)30を作成して蛍光発光シンチレータ素子位置により補正する。
【0008】
図2において、20は、アンプ22A、22B、22C、22D、加算回路24A、24B、24C、A/D変換回路26A、26B、前記位置LUT28、前記エネルギLUT30、エネルギウィンドウ32を含んで構成される位置演算回路である。
【0009】
LUT作成のためには従来から、均一な並行ビームのγ線を照射し、PMT14A〜Dの出力を用いて重心演算することにより、各イベントに対する2次元位置マップ上の対応位置(X,Y)を算出し、図3に示したような2次元位置マップに描画することが行われている。これらは各蛍光発光シンチレータ素子(結晶素子)位置に応じたピークを持つ分布として描かれ、シンチレータ対応域を示すため、従来法では、それらピーク間の谷間を境界線として、クラスタに分割し、ルックアップテーブル(LUT)を作成していた。実際の測定では、入射した各イベント毎に重心演算を行い、LUTを参照しながら、各クラスタに対応したシンチレータ素子位置を出力する。
【0010】
【特許文献1】
特開昭63−47686号公報
【特許文献2】
特開平5−75940号公報
【0011】
【発明が解決しようとする課題】
しかしながら、特許文献1,2は、いずれも放射線位置検出器自体に関するもので、位置弁別のための校正方法については述べられていない。
【0012】
従って、上述したような従来の技術においては、放射線が入射して蛍光を発したシンチレータ素子の同定を正確に行うために、2次元位置マップを直線や円形で区分けする人の主観に左右される方法が用いられてきた。
【0013】
又、図1乃至図3に示した例においては、図2及び図3は簡単のため1層で説明したが、シンチレータ素子を多層に配置した3次元放射線位置検出器の場合、層毎に位置マップへ反映されるパターンが異なり、これらの重ね合わせとして投影される。そのため、規則的に位置マップのピークが配列されず、単純にピーク間の谷から境界を決めるのが難しい。更に各シンチレータ素子の感度やシンチレータ素子対応域の広がりが異なる場合などでは、ピークの谷をクラスタ間の境界と設定しただけでは、位置弁別の精度が低下することがある。
【0014】
又、ピーク判別のために多くのカウント数が必要で、LUT作成に時間がかかる上、ピークとピークが接近した場合、境界が判別しにくくなり、全自動で境界を判別することが困難である。
【0015】
更に比較的低い原子番号・密度のシンチレータの場合、各シンチレータ素子を構成する結晶に入射したγ線の一部は、光電吸収されずにコンプトン散乱によって結晶外に逃げてしまうが、その際エネルギの一部を結晶に付与するので、このコンプトン散乱がPMTの出力信号の統計誤差を大きくする。エネルギ弁別によって散乱成分の大部分は取り除くことができるが、多重散乱成分は残るため、それらが位置マップのピークのバックグランド成分として現れる。この多重散乱成分のため従来のLUT作成方法では、これらのイベントを取り除くことが困難である。例えば分解能を優先させる場合、多重散乱成分は位置弁別については不確定領域なので取り除く必要があるが、取り除きすぎて利用できるイベント数が減少してしまうことは、放射線位置検出器の感度を低下させることになるので好ましくない。
【0016】
本発明は、前記従来の問題点を解決するべくなされたもので、2次元位置マップを作成するのに統計的クラスタリング処理を用い、位置弁別の精度を高めることを課題とする。
【0017】
【課題を解決するための手段】
本発明は、1次元、2次元、あるいは3次元的に配置された多数のシンチレータ素子と、それらに光結合された光センサから構成された放射線位置検出器を校正するに際して、前記光センサから得られる複数の出力信号から統計的クラスタリング処理によって、放射線の検出位置に属する領域弁別を行なうようにして、前記課題を解決したものである。
【0018】
又、前記クラスタリング処理に際して、まず、1段階目のクラスタリングで各シンチレータ素子の中心位置を求め、次いで、全データの中心を加えた2段階目のクラスタリングで各シンチレータ素子の広がりを求めるようにしたものである。
【0019】
又、前記光センサに吸収された放射線エネルギの選別情報も得るようにしたものである。
【0020】
本発明は、又、1次元、2次元、あるいは3次元的に配置された多数のシンチレータ素子と、それらに光結合された光センサから構成された放射線位置検出器を校正するに際して、前記光センサから得られる複数の出力信号から統計的クラスタリング処理によって、該光センサに吸収された放射線エネルギの選別情報を得るようにしたものである。
【0021】
本発明は、又、前記の校正方法を用いて作成されたルックアップテーブル(位置ルックアップテーブルやエネルギルックアップテーブル)を備えたことを特徴とする放射線位置検出器を提供するものである。
【0022】
本発明は、上記の目的を達成するため、複数の光センサの出力がシンチレータ素子数に応じたクラスタに分類できることを利用し、混合ガウスモデル(Mixture Gaussian Model: MGM)に代表される統計的クラスタリングによってシンチレータ素子の弁別を行う。MGMにおいては、実測データから各クラスタを構成するガウス分布の平均、分散および各クラスタの存在確率(事前確率とも呼ばれる)を推定する。
【0023】
予め各クラスタに付番しておき、p(x|j)をjクラスタから値xが出力される条件付き確率、P(j)をjクラスタが存在する確率、Mをクラスタ数とすると、xが出力される確率p(x)は、混合モデルを用いて、
【数1】
【数2】
となる。ここで、μとΣは、それぞれ各クラスタ固有の平均と共分散であり、Dは実測データの次元である。例えば、一様照射のデータから(1)(2)式で得られる2次元位置マップを利用する場合、D=2となる。
【0024】
実測データからP(j)、p(x|j)を計算すると、次式で示したBayes定理を適用することで、各ピクセルに対応する出力xがjクラスタに属する確率(Bayes定理でいう事後確率)を各々求めることができる(P(x)=1と仮定)。
【0025】
【数3】
【0026】
これによって2次元位置マップ上における全てのクラスタでの事後確率を取得する。その後、出力xをデジタル化した各ピクセルについて事後確率が最大のクラスタ番号で位置マップ全体をラベリングし、テーブルを作成する。テーブル参照方式を用いることで、テーブル作成以降は、従来のPET装置の回路系をそのまま利用することができる。
【0027】
このようにして、事後確率p(j|x)を使うことで、jクラスタと対応するγ線が入射するシンチレータ素子位置を確率論的に求めることができるため、手動あるいはピーク間の谷を区切るだけの単純な自動化手法に比べて、正確に位置弁別を行うことができる。また、各シンチレータ素子対応域に対応するクラスタと複数のシンチレータ素子にまたがる多重散乱成分のクラスタを分けることで、多重散乱成分を除いて、各シンチレータ素子対応域のみを選別することができる。更にクラスタの広がりを閾値によって調整することで、分解能(位置弁別精度)と感度を最適化することができる。
【0028】
すなわち、シンチレータ素子対応域を広げていくと多重散乱を取り込んでしまうため位置弁別精度は低下するが、利用できるイベント数は増えて感度を向上できる。一方でシンチレータ素子対応域を狭めると感度は低下するが、位置弁別精度は向上する。
【0029】
各シンチレータ素子対応域の範囲をσのn倍で変動させた時、多重散乱がシンチレータ素子対応域に含まれるイベントと損失するイベントの全イベントに対する割合を比較した結果を図4に示す。図の横軸は、各クラスタの分散σの倍数である。図から明らかなように、感度を重視するか位置分解能を重視するかを利用目的に応じて選択することが可能である。
【0030】
またガウス分布の確率密度関数は、他の確率分布の近似関数としての性質も良好であり、P(j)、p(x|j)の計算に期待値最大化(Expectation Maximization:EM)アルゴリズムのように収束性が強固な推定アルゴリズムが使用可能である。
そのため、比較的少ないカウント数でもLUTの作成が可能である上に、近接したピークでも比較的安定してピークを抽出することが可能となる。
【0031】
また、説明に用いる次元マップは、回路内あるいはコンピュータ内などでの処理概念を視覚的に説明するため示したものであり、確認などのため実際に図5のように表示することもあるが、本発明の放射線位置検出器がPET装置などに組み込まれた場合は、処理速度が優先されるので表示させる必要はない。
【0032】
【発明の実施の形態】
以下、図面により本発明を詳しく説明する。
【0033】
図6(A)(平面図)及び(B)(正面図)に示す如く、3次元結晶配列による検出器の深さ方向の位置情報を取得するDepth−of−interaction(DOI)検出器を例にとって説明する。
【0034】
検出器モジュール40は、図7に示すように2×2の各2.9mm×2.9mm×7.5mmのGSOシンチレータ素子のアレイ44A、44B、44C、44Dを4段積み重ねたものを1ブロックとし、これを8×8個並べたシンチレータブロック42と、その最下段の16×16個のシンチレータ素子と各アノードが光結合された16×16のマルチアノードを持つ52mm角のフラットパネルPS−PMT(位置弁別型光電子増倍管)46からなる。図において、48は素子間に配設された反射板であり、46は1アノードを示す。
【0035】
前記PS−PMT46からの出力は抵抗チェーンを用いて4ヶ所にまとめられ、重心演算(アンガー方式)によって2次元位置マップ上にパターン表示される。
【0036】
前記検出器モジュール40によって得られる、図8に示すような位置マップを用いて、γ線が相互作用したシンチレータ素子位置を求める。ここで奇数層(1層44A、3層44C)と偶数層(2層44B、4層44D)とでセリウム濃度の異なるシンチレータ素子を用いて、減衰時間の違いから奇数層マップと偶数層マップに分けている。つまり1回の計測で1つの検出器当たり2つのマップが作成されることになる。従来は、これらのマップを直線や円形で区分けする人の主観に左右される方法が用いられてきた。
【0037】
本実施形態では、137Cs線源を検出器上方から一様照射して従来の2次元位置マップを得た後、MGMによるクラスタリングを行い、マップ上に投影される全ての範囲での事後確率を取得する。その後、各ピクセルについて事後確率が最大のクラスタ番号でマップ全体をラベリングし、テーブルを作成する(本実施例では4096×4096ピクセルに分割している)。又、テーブルに閾値を設けることでシンチレータ素子対応域と多重散乱領域を分け、シンチレータ素子対応域のみを選別する。
【0038】
PS−PMT46の4つのアノード出力A、B、C、Dから作成された、図9(1)の2次元位置マップに対して2段階のMGMによるクラスタリングを行う場合、まず図10に示す如く、ステップ100の1段階目のクラスタリングにおいて、クラスタ数をシンチレータ素子数(ここでは8個、即ち1層目44Aと3層目44C、あるいは2層目44Bと4層目44Dの各2層の合計のシンチレータ数)より少し多め(例えば10)に取ってクラスタリングを行い、シンチレータ素子位置を得る。
【0039】
ステップ110で、8つのシンチレータ素子対応域が全て取得できたと判定された場合、1段階目のMGMで得られたガウス分布の平均をシンチレータ素子対応域の確定位置として固定すると共に、ステップ120で、多重散乱成分として全データの重心を中心とする9番目のクラスタを追加する。そして、ステップ130で、共分散及び事後確率のみの推定を2段階目のMGMによるクラスタリングで行い、クラスタの範囲を得る。
【0040】
このようにして、ステップ140で、クラスタリング結果から、それぞれのクラスタの平均と分散が得られる。
【0041】
クラスタリングによって得られた結果を図9(2)に示す。2段階のMGMを行うことによって、8つのシンチレータ素子対応域および多重散乱成分を弁別できている。図において、楕円はMGMにより推定された共分散楕円の3倍として示した。
【0042】
この結果から、ステップ150で、事後確率が最大になるように、それぞれのクラスタでのシンチレータ素子位置を得て、LUTを作成すると、図5(1)のようになる。この場合、1から8番までのシンチレータ素子に対応する各クラスタにおいて事後確率が最大のクラスタの中で3σ以内の領域のクラスタを利用し、これ以外の領域に来たイベントは捨てることなる。
【0043】
一方、従来のテーブル作成手法を拡張したものでは、全データの中心に対して十字に領域を切った後、奇数層マップ(図5(1)(b))と偶数層マップ(図5(1)(a))についてそれぞれ谷の部分を算出し領域分けすることにより、図5(2)のようなLUTが得られる。こうしたLUTを図2に示すフロントエンド回路上に作成して、ステップ160でLUTを参照し、ステップ170で、実際の測定時の各入射γ線イベントに対して位置弁別を行い、相互作用したシンチレータ素子位置を得ていた。
【0044】
光センサに吸収された放射線エネルギの選別情報は、LUTを基にエネルギがエネルギウインドウ32に入っているか否かを選別する。エネルギウインドウ32は、図2における和信号24Cの波高幅を予め設定したものである。
【0045】
本実施形態においては、奇数の1、3層及び偶数の2、4層にセリウム濃度の異なるシンチレータ素子を用いて、減衰時間の違いから奇数層マップと偶数層マップに分けているので、2つのマップが容易に作成できる。なお、マップを作り分ける方法はこれに限定されない。
【0046】
又、前記実施形態においては、GSOシンチレータ素子のアレイを4段積み重ねたものを1ブロックとし、これを8×8個並べたシンチレータブロック42と、16×16のマルチアノードを持つ52mm角のフラットパネルPS−PMT46からなる検出器モジュール40が用いられていたが、シンチレータの種類や形状、構成、光センサの種類は、これに限定されず、例えばシンチレータとしてBGOやLSOを用いたり、光センサとしてAPDを用いることも可能である。
【0047】
又、前記実施形態においては、光センサの出力A、B、C、Dより重心演算されていたが、光センサの出力A、B、C、Dを直接クラスタリングすることも可能である。又、主成分分析等により、それらを低次元のデータに変換したものをクラスタリングしてもよい。
【0048】
又、クラスタリングの手法は混合モデルが好適であるが、他のクラスタリング手法を用いることも可能である。
【0049】
本発明によって特定される放射線の入射位置を基に校正される物理量としては、位置の他、γ線エネルギ、入射時間(タイミング)等が考えられる。
【0050】
本発明は、特許文献1や2に記載された放射線位置検出器にも適用可能であるが、適用対象は、これらに限定されない。
【0051】
【発明の効果】
本発明では、2次元位置マップを統計的にクラスタリングすることで、位置弁別の精度を増すことができる。実際、前記実施形態に対して、シミュレーションによって得られた誤答率を比較すると、図11のように、従来法に比較して誤答率が大幅に改善できている。
【0052】
また、本発明による手法は、図4に示す如くシンチレータ素子対応域を狭めれば、多重散乱成分を可能な限り取り込まないようにできる。このとき感度は犠牲になるものの、位置弁別の精度は向上させることが可能である。反対に、シンチレータ素子対応域を広げることで、多重散乱成分が増すため位置弁別精度は低下するものの、利用できるイベント数をできるだけ取り込むようにすることもできる。このように、感度を重視するか位置分解能を重視するかを利用目的に応じて選択することが可能である。
【0053】
更に、イベント数が少ない場合でも比較的安定に処理が可能であるため、LUTを作成するための校正時間を短縮することができる。更に任意の形状の2次元マップに対してもLUTの自動作成が可能になる。
【図面の簡単な説明】
【図1】従来のブロック検出器の構成例を示す斜視図
【図2】同じく光センサと位置弁別回路の構成例を示すブロック図
【図3】同じく位置弁別テーブルの例を示す平面図
【図4】本発明の実施形態による結晶領域を広げた場合の位置弁別精度と感度を示す線図
【図5】図9のクラスタリング結果の最大事後確率から得られる位置弁別テーブルと従来法で得られたテーブルを比較して示す線図
【図6】本発明の実施形態で用いられる深さ弁別型(DOI)検出器の構成を示す平面図及び正面図
【図7】前記DOI検出器のモジュール構造を示す斜視図
【図8】前記実施形態で得られる2次元位置マップの例を示す平面図
【図9】一様照射による2次元位置マップの例とそのクラスタリング例を示す平面図
【図10】本発明の実施形態の処理手順を示す流れ図
【図11】本発明と従来法による位置弁別の精度を比較して示す線図
【符号の説明】
28…位置ルックアップテーブル(LUT)
30…エネルギルックアップテーブル(LUT)
40…検出器モジュール
42…シンチレータブロック
46…位置弁別型光電子増倍管(PS−PMT)[0001]
BACKGROUND OF THE INVENTION
The present invention relates to a calibration method for a radiation position detector, and in particular, when calibrating a radiation position detector composed of a large number of scintillator elements arranged three-dimensionally and an optical sensor optically coupled thereto. The present invention relates to a method for calibrating a radiation position detector having a radiation position detection function and a radiation absorption energy selection function suitable for use.
[0002]
[Prior art]
When a radioisotope (RI) is administered to a subject, an event occurs in which two photons called annihilation gamma rays are emitted simultaneously in opposite directions. The positron emission tomography apparatus (PET) is an apparatus that obtains an RI distribution image in a subject by detecting these two photons by coincidence counting.
[0003]
As a radiation three-dimensional position detector used in such a PET apparatus, a plurality of scintillator elements are stacked on an optical sensor while sandwiching transparent plates having different refractive indexes to form a scintillator block, and reach the optical sensor. A radiation detector that identifies a scintillator element that emits fluorescence by entering radiation based on the difference in the amount of light received by the optical sensor by making the light transmittance different for each scintillator element is disclosed in
[0004]
In addition, a plurality of scintillator elements are stacked on the optical position detector so that the center position thereof is deviated in a direction parallel to the light receiving surface of the optical position detector, and the center of gravity of the spatial distribution of the output light from the optical position
[0005]
Further, as shown in FIG. 1, in the case of a
[0006]
X = (A + B) / (A + B + C + D) (1)
Y = (A + C) / (A + B + C + D) (2)
[0007]
However, the corresponding area corresponding to each scintillator element on the two-dimensional position map has a certain spread due to statistical variations in the number of detector photons and the like, depending on differences in gains of the
[0008]
In FIG. 2, 20 includes
[0009]
Conventionally, in order to create an LUT, a uniform parallel beam of γ rays is irradiated, and the center of gravity is calculated using the outputs of the
[0010]
[Patent Document 1]
JP 63-47686 A [Patent Document 2]
Japanese Patent Laid-Open No. 5-75940
[Problems to be solved by the invention]
However,
[0012]
Therefore, in the conventional technology as described above, in order to accurately identify the scintillator element that emits fluorescence when radiation is incident, it depends on the subjectivity of the person who divides the two-dimensional position map into straight lines and circles. A method has been used.
[0013]
In the example shown in FIGS. 1 to 3, FIGS. 2 and 3 have been described with one layer for the sake of simplicity. However, in the case of a three-dimensional radiation position detector in which scintillator elements are arranged in multiple layers, the position is determined for each layer. The pattern reflected on the map is different and projected as a superposition of these. Therefore, the peaks of the position map are not regularly arranged, and it is difficult to simply determine the boundary from the valley between the peaks. Further, when the sensitivity of each scintillator element and the spread of the scintillator element corresponding region are different, the accuracy of position discrimination may be reduced only by setting the peak valley as the boundary between clusters.
[0014]
In addition, a large number of counts are required for peak discrimination, and it takes time to create an LUT. In addition, when a peak and a peak approach each other, it becomes difficult to discriminate the boundary, and it is difficult to discriminate the boundary automatically. .
[0015]
Furthermore, in the case of a scintillator with a relatively low atomic number and density, a part of the γ rays incident on the crystal constituting each scintillator element escapes from the crystal by Compton scattering without being photoelectrically absorbed. Since a part is imparted to the crystal, this Compton scattering increases the statistical error of the output signal of the PMT. Although most of the scattered components can be removed by energy discrimination, multiple scattered components remain, so they appear as background components at the peak of the position map. Because of this multiple scattering component, it is difficult to remove these events with the conventional LUT creation method. For example, when priority is given to resolution, multiple scattered components need to be removed because they are indeterminate for position discrimination, but reducing the number of events that can be used excessively reduces the sensitivity of the radiation position detector. This is not preferable.
[0016]
The present invention has been made to solve the above-described conventional problems, and an object of the present invention is to improve the accuracy of position discrimination by using a statistical clustering process to create a two-dimensional position map.
[0017]
[Means for Solving the Problems]
The present invention obtains a radiation position detector composed of a large number of scintillator elements arranged one-dimensionally, two-dimensionally or three-dimensionally and a photosensor optically coupled to the scintillator elements from the photosensor. The above-described problems are solved by discriminating regions belonging to the radiation detection position by statistical clustering processing from a plurality of output signals.
[0018]
In the clustering process, first, the center position of each scintillator element is obtained by the first stage clustering, and then the spread of each scintillator element is obtained by the second stage clustering including the centers of all data. It is.
[0019]
In addition, selection information of radiation energy absorbed by the optical sensor is also obtained.
[0020]
The present invention also provides a method for calibrating a radiation position detector comprising a large number of scintillator elements arranged one-dimensionally, two-dimensionally or three-dimensionally, and an optical sensor optically coupled thereto. The selection information of the radiation energy absorbed by the optical sensor is obtained from a plurality of output signals obtained from the above by a statistical clustering process.
[0021]
The present invention also provides a radiation position detector comprising a look-up table (position look-up table or energy look-up table) created using the calibration method.
[0022]
In order to achieve the above object, the present invention utilizes the fact that the outputs of a plurality of photosensors can be classified into clusters corresponding to the number of scintillator elements, and statistical clustering represented by a mixed Gaussian model (MGM). To discriminate the scintillator elements. In the MGM, the average and variance of the Gaussian distribution constituting each cluster and the existence probability (also referred to as prior probability) of each cluster are estimated from actually measured data.
[0023]
If each cluster is numbered in advance, p (x | j) is a conditional probability that the value x is output from the j cluster, P (j) is the probability that the j cluster exists, and M is the number of clusters. The probability p (x) that is output using the mixed model is
[Expression 1]
It becomes. Since p (x | j) assumes a Gaussian distribution,
[Expression 2]
It becomes. Here, μ and Σ are the mean and covariance specific to each cluster, respectively, and D is the dimension of the measured data. For example, when using a two-dimensional position map obtained from the uniform irradiation data (1) and (2), D = 2.
[0024]
When P (j) and p (x | j) are calculated from the measured data, the probability that the output x corresponding to each pixel belongs to the j cluster by applying the Bayes theorem shown by the following equation (the posterior in the Bayes theorem) Probability) can be obtained (assuming P (x) = 1).
[0025]
[Equation 3]
[0026]
As a result, the posterior probabilities in all clusters on the two-dimensional position map are acquired. Thereafter, the entire position map is labeled with the cluster number having the maximum posterior probability for each pixel obtained by digitizing the output x, and a table is created. By using the table reference method, the circuit system of the conventional PET apparatus can be used as it is after the table creation.
[0027]
In this way, by using the posterior probability p (j | x), the scintillator element position where the γ ray corresponding to the j cluster is incident can be obtained probabilistically. Compared to a simple automated method, position discrimination can be performed accurately. In addition, by dividing the cluster corresponding to each scintillator element corresponding region and the cluster of multiple scattering components extending over a plurality of scintillator elements, only the scintillator element corresponding region can be selected except for the multiple scattering component. Furthermore, the resolution (position discrimination accuracy) and sensitivity can be optimized by adjusting the spread of the cluster according to the threshold value.
[0028]
That is, if the scintillator element corresponding area is expanded, multiple scattering is taken in, so that the position discrimination accuracy is lowered, but the number of available events increases and the sensitivity can be improved. On the other hand, when the scintillator element corresponding area is narrowed, the sensitivity is lowered, but the position discrimination accuracy is improved.
[0029]
When the range of each scintillator element corresponding region is fluctuated by n times σ, FIG. 4 shows the result of comparing the ratios of the events included in the scintillator element corresponding region and the lost events to all the events. The horizontal axis in the figure is a multiple of the variance σ of each cluster. As is apparent from the figure, it is possible to select whether to place importance on sensitivity or position resolution depending on the purpose of use.
[0030]
Further, the probability density function of the Gaussian distribution also has good properties as an approximation function of other probability distributions, and the expectation maximization (Expectation Maximization: EM) algorithm is used for calculating P (j) and p (x | j). Thus, it is possible to use an estimation algorithm with strong convergence.
For this reason, it is possible to create an LUT with a relatively small number of counts, and it is possible to extract a peak relatively stably even with adjacent peaks.
[0031]
In addition, the dimension map used for explanation is shown for visually explaining the concept of processing in the circuit or in the computer, and may be actually displayed as shown in FIG. When the radiation position detector of the present invention is incorporated in a PET apparatus or the like, there is no need to display it because processing speed is given priority.
[0032]
DETAILED DESCRIPTION OF THE INVENTION
Hereinafter, the present invention will be described in detail with reference to the drawings.
[0033]
As shown in FIGS. 6A (plan view) and (B) (front view), a depth-of-interaction (DOI) detector that acquires position information in the depth direction of the detector by a three-dimensional crystal array is taken as an example. I will explain to you.
[0034]
As shown in FIG. 7, the detector module 40 is one block in which 2 × 2 arrays of 2.9 mm × 2.9 mm × 7.5 mm
[0035]
Outputs from the PS-
[0036]
A position map as shown in FIG. 8 obtained by the detector module 40 is used to determine the scintillator element position where the γ rays interact. Here, scintillator elements having different cerium concentrations are used for the odd layer (1
[0037]
In this embodiment, a 137Cs radiation source is uniformly irradiated from above the detector to obtain a conventional two-dimensional position map, and then clustering by MGM is performed to obtain posterior probabilities in the entire range projected on the map. To do. Thereafter, the entire map is labeled with the cluster number having the maximum posterior probability for each pixel, and a table is created (in this embodiment, it is divided into 4096 × 4096 pixels). Further, by providing a threshold value in the table, the scintillator element corresponding area and the multiple scattering area are separated, and only the scintillator element corresponding area is selected.
[0038]
When two-stage MGM clustering is performed on the two-dimensional position map of FIG. 9A created from the four anode outputs A, B, C, and D of the PS-
[0039]
When it is determined in
[0040]
Thus, at
[0041]
The result obtained by clustering is shown in FIG. By performing the two-stage MGM, it is possible to discriminate eight scintillator element corresponding areas and multiple scattering components. In the figure, the ellipse is shown as three times the covariance ellipse estimated by MGM.
[0042]
From this result, when the scintillator element positions in the respective clusters are obtained in
[0043]
On the other hand, in the extension of the conventional table creation method, the odd layer map (FIG. 5 (1) (b)) and the even layer map (FIG. ) (A)), the valley portion is calculated and divided into regions to obtain the LUT as shown in FIG. Such an LUT is created on the front-end circuit shown in FIG. 2, and the LUT is referred to in
[0044]
The selection information of the radiation energy absorbed by the optical sensor is used to select whether the energy is in the
[0045]
In this embodiment, scintillator elements having different cerium concentrations are used for odd 1, 3 layers and even 2, 4 layers, and are divided into an odd layer map and an even layer map from the difference in attenuation time. Maps can be created easily. Note that the method of creating a map is not limited to this.
[0046]
Further, in the above-described embodiment, a 4-mm stacked array of GSO scintillator elements is used as one block, and a 8 × 8 array of scintillator blocks 42 and a 52 mm square flat panel having 16 × 16 multi-anodes. Although the detector module 40 comprising the PS-
[0047]
In the above embodiment, the center of gravity is calculated from the outputs A, B, C, and D of the optical sensor. However, the outputs A, B, C, and D of the optical sensor can be directly clustered. Further, data obtained by converting them into low-dimensional data by principal component analysis or the like may be clustered.
[0048]
A clustering method is preferably a mixed model, but other clustering methods can also be used.
[0049]
As physical quantities to be calibrated based on the incident position of radiation specified by the present invention, γ-ray energy, incident time (timing), etc. can be considered in addition to the position.
[0050]
The present invention can also be applied to the radiation position detector described in
[0051]
【The invention's effect】
In the present invention, the accuracy of position discrimination can be increased by statistically clustering the two-dimensional position map. In fact, when the error rate obtained by the simulation is compared with the embodiment, the error rate can be greatly improved as compared with the conventional method as shown in FIG.
[0052]
In addition, the technique according to the present invention can prevent multiple scattering components from being taken in as much as possible by narrowing the scintillator element corresponding area as shown in FIG. At this time, although the sensitivity is sacrificed, the accuracy of position discrimination can be improved. On the contrary, by expanding the area corresponding to the scintillator element, the number of available events can be taken in as much as possible although the position discrimination accuracy is lowered because the multiple scattering components are increased. In this way, it is possible to select whether to place importance on sensitivity or position resolution depending on the purpose of use.
[0053]
Furthermore, since the process can be performed relatively stably even when the number of events is small, the calibration time for creating the LUT can be shortened. Furthermore, an LUT can be automatically created for a two-dimensional map having an arbitrary shape.
[Brief description of the drawings]
FIG. 1 is a perspective view showing a configuration example of a conventional block detector. FIG. 2 is a block diagram showing a configuration example of an optical sensor and a position discrimination circuit. FIG. 3 is a plan view showing an example of a position discrimination table. 4 is a diagram showing position discrimination accuracy and sensitivity when the crystal region is expanded according to an embodiment of the present invention. FIG. 5 is a position discrimination table obtained from the maximum posterior probability of the clustering result of FIG. FIG. 6 is a diagram showing a comparison of tables. FIG. 6 is a plan view and a front view showing a configuration of a depth discrimination type (DOI) detector used in the embodiment of the present invention. FIG. 7 is a diagram showing a module structure of the DOI detector. FIG. 8 is a plan view showing an example of a two-dimensional position map obtained in the embodiment. FIG. 9 is a plan view showing an example of a two-dimensional position map by uniform irradiation and an example of clustering thereof. Processing hand of embodiment of invention Diagram showing a comparison of the accuracy of the position discrimination by the flow diagram 11 the present invention and the conventional method shown an [Description of symbols]
28 ... Position Lookup Table (LUT)
30 ... Energy Lookup Table (LUT)
40 ...
Claims (6)
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2003200616A JP4367903B2 (en) | 2003-07-23 | 2003-07-23 | Radiation position detector calibration method and radiation position detector |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2003200616A JP4367903B2 (en) | 2003-07-23 | 2003-07-23 | Radiation position detector calibration method and radiation position detector |
Publications (2)
Publication Number | Publication Date |
---|---|
JP2005043104A true JP2005043104A (en) | 2005-02-17 |
JP4367903B2 JP4367903B2 (en) | 2009-11-18 |
Family
ID=34260965
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2003200616A Expired - Lifetime JP4367903B2 (en) | 2003-07-23 | 2003-07-23 | Radiation position detector calibration method and radiation position detector |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP4367903B2 (en) |
Cited By (15)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2008023451A1 (en) * | 2006-08-25 | 2008-02-28 | National Institute Of Radiological Sciences | Positron emission tomography scanner and radiation detector |
JP2009085683A (en) * | 2007-09-28 | 2009-04-23 | Shimadzu Corp | Nuclear medicine equipment |
WO2009116174A1 (en) | 2008-03-21 | 2009-09-24 | 株式会社島津製作所 | Method of correcting two-dimensional position map and radiation detector |
WO2009141861A1 (en) | 2008-05-19 | 2009-11-26 | 株式会社島津製作所 | Method for calibrating two-dimensional position map, and radiation detector |
WO2010074400A2 (en) * | 2008-12-26 | 2010-07-01 | 연세대학교 산학협력단 | Multi-planar detector for measuring gamma-ray image and method for detecting 3d position |
JP2010160045A (en) * | 2009-01-08 | 2010-07-22 | Shimadzu Corp | Method for correcting light amount output of pet detector |
KR101025513B1 (en) | 2008-06-05 | 2011-04-04 | 연세대학교 산학협력단 | A qusai-monolothic detector module for depth of interaction measurement for PET |
JP2011185840A (en) * | 2010-03-10 | 2011-09-22 | Shimadzu Corp | Radiation signal-processing unit and radiation detector provided with the same |
JP2012037523A (en) * | 2005-09-14 | 2012-02-23 | Dainippon Printing Co Ltd | Image processing apparatus, and image processing program |
CN102933983A (en) * | 2010-06-10 | 2013-02-13 | 株式会社岛津制作所 | Two-dimensional position map correction method |
JP2013113790A (en) * | 2011-11-30 | 2013-06-10 | Hamamatsu Photonics Kk | Radiation detection device |
WO2014001926A1 (en) * | 2012-06-27 | 2014-01-03 | Koninklijke Philips N.V. | Digital positron emission tomography (dpet) energy calibration method |
JP2015530597A (en) * | 2012-10-08 | 2015-10-15 | コミッサリア ア レネルジー アトミーク エ オ ゼネルジ ザルタナテイヴ | Method and apparatus for detecting ionizing radiation using a pixel-type photodetector |
CN107643538A (en) * | 2017-10-18 | 2018-01-30 | 中国检验检疫科学研究院 | A kind of device for being used to verify scintillator detector performance |
JP2020190446A (en) * | 2019-05-21 | 2020-11-26 | 株式会社島津製作所 | Calibration method of two-dimensional position map of radiation detector and radiation detector |
-
2003
- 2003-07-23 JP JP2003200616A patent/JP4367903B2/en not_active Expired - Lifetime
Cited By (26)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JP2012037523A (en) * | 2005-09-14 | 2012-02-23 | Dainippon Printing Co Ltd | Image processing apparatus, and image processing program |
JP2008051701A (en) * | 2006-08-25 | 2008-03-06 | Natl Inst Of Radiological Sciences | Positron emission tomographic equipment and radiation detector |
US7671339B2 (en) | 2006-08-25 | 2010-03-02 | National Institute Of Radiological Sciences | Positron emission tomography scanner and radiation detector |
WO2008023451A1 (en) * | 2006-08-25 | 2008-02-28 | National Institute Of Radiological Sciences | Positron emission tomography scanner and radiation detector |
JP2009085683A (en) * | 2007-09-28 | 2009-04-23 | Shimadzu Corp | Nuclear medicine equipment |
WO2009116174A1 (en) | 2008-03-21 | 2009-09-24 | 株式会社島津製作所 | Method of correcting two-dimensional position map and radiation detector |
US8278625B2 (en) | 2008-03-21 | 2012-10-02 | Shimadzu Corporation | Two-dimensional position map correcting method |
WO2009141861A1 (en) | 2008-05-19 | 2009-11-26 | 株式会社島津製作所 | Method for calibrating two-dimensional position map, and radiation detector |
US8471211B2 (en) | 2008-05-19 | 2013-06-25 | Shimadzu Corporation | Two-dimensional position map correcting method, and radiation detecting apparatus |
JP4983981B2 (en) * | 2008-05-19 | 2012-07-25 | 株式会社島津製作所 | Two-dimensional position map calibration method and radiation detection apparatus |
KR101025513B1 (en) | 2008-06-05 | 2011-04-04 | 연세대학교 산학협력단 | A qusai-monolothic detector module for depth of interaction measurement for PET |
WO2010074400A2 (en) * | 2008-12-26 | 2010-07-01 | 연세대학교 산학협력단 | Multi-planar detector for measuring gamma-ray image and method for detecting 3d position |
WO2010074400A3 (en) * | 2008-12-26 | 2010-08-12 | 연세대학교 산학협력단 | Multi-planar detector for measuring gamma-ray image and method for detecting 3d position |
KR101111011B1 (en) * | 2008-12-26 | 2012-02-15 | 연세대학교 산학협력단 | Multi-layer scintillation detector with 3D positioning capability for pinhole gamma camera |
JP2010160045A (en) * | 2009-01-08 | 2010-07-22 | Shimadzu Corp | Method for correcting light amount output of pet detector |
JP2011185840A (en) * | 2010-03-10 | 2011-09-22 | Shimadzu Corp | Radiation signal-processing unit and radiation detector provided with the same |
US9045097B2 (en) | 2010-06-10 | 2015-06-02 | Shimadzu Corporation | Two-dimensional position map correcting method |
CN102933983A (en) * | 2010-06-10 | 2013-02-13 | 株式会社岛津制作所 | Two-dimensional position map correction method |
JP2013113790A (en) * | 2011-11-30 | 2013-06-10 | Hamamatsu Photonics Kk | Radiation detection device |
WO2014001926A1 (en) * | 2012-06-27 | 2014-01-03 | Koninklijke Philips N.V. | Digital positron emission tomography (dpet) energy calibration method |
RU2620863C2 (en) * | 2012-06-27 | 2017-05-30 | Конинклейке Филипс Н.В. | Method for energy calibrating digital positron emission tomography (dpet) |
US10371836B2 (en) | 2012-06-27 | 2019-08-06 | Koninklijke Philips N.V. | Digital positron emission tomography (DPET) energy calibration method |
JP2015530597A (en) * | 2012-10-08 | 2015-10-15 | コミッサリア ア レネルジー アトミーク エ オ ゼネルジ ザルタナテイヴ | Method and apparatus for detecting ionizing radiation using a pixel-type photodetector |
CN107643538A (en) * | 2017-10-18 | 2018-01-30 | 中国检验检疫科学研究院 | A kind of device for being used to verify scintillator detector performance |
JP2020190446A (en) * | 2019-05-21 | 2020-11-26 | 株式会社島津製作所 | Calibration method of two-dimensional position map of radiation detector and radiation detector |
JP7247745B2 (en) | 2019-05-21 | 2023-03-29 | 株式会社島津製作所 | Method for calibrating two-dimensional position map of radiation detection device and radiation detection device |
Also Published As
Publication number | Publication date |
---|---|
JP4367903B2 (en) | 2009-11-18 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Gonzalez-Montoro et al. | Evolution of PET detectors and event positioning algorithms using monolithic scintillation crystals | |
Mohammadi et al. | Minimization of parallax error in positron emission tomography using depth of interaction capable detectors: methods and apparatus | |
US7601963B2 (en) | High-resolution depth-of-interaction PET detector | |
Borghi et al. | A 32 mm× 32 mm× 22 mm monolithic LYSO: Ce detector with dual-sided digital photon counter readout for ultrahigh-performance TOF-PET and TOF-PET/MRI | |
Ling et al. | Depth of interaction decoding of a continuous crystal detector module | |
US6943355B2 (en) | SPECT gamma camera | |
US9335425B2 (en) | Modelling of ToF-DOI detector arrays | |
US20100044571A1 (en) | Method for determining the three-dimensional position of a scintillation event | |
JP5641930B2 (en) | Diagnostic imaging system, time stamp calculation method, processor for executing the method, and computer-readable medium programmed with the method | |
US7397035B2 (en) | Scatter correction for time-of-flight positron emission tomography data | |
JP4367903B2 (en) | Radiation position detector calibration method and radiation position detector | |
US7227149B2 (en) | Method and system for positron emission tomography image reconstruction | |
US20090224164A1 (en) | Scintillation detector for positron emission tomography | |
US20110150181A1 (en) | Apparatus and methods for detector scatter recovery for nuclear medicine imaging systems | |
CN107735694A (en) | Apparatus and method for the interaction depth position emissron tomography detector using two points of sensings | |
US9091771B2 (en) | System and method for improving detection of gamma interactions in a positron emission tomography system | |
CN107850677A (en) | Compton camera system and method for detecting gamma radiation | |
Binder et al. | Performance evaluation of a staggered three-layer DOI PET detector using a 1 mm LYSO pitch with PETsys TOFPET2 ASIC: comparison of HAMAMATSU and KETEK SiPMs | |
EP0892284B1 (en) | Method and apparatus for diagnostic imaging | |
JP7453359B2 (en) | Method for position and energy determination in scintillation detectors | |
Yoshida et al. | Calibration procedure for a DOI detector of high resolution PET through a Gaussian mixture model | |
Jin et al. | Experimental evaluation of a 3-D CZT imaging spectrometer for potential use in Compton-Enhanced PET imaging | |
Astakhov et al. | Effect of depth of interaction decoding on resolution in PET: A simulation study | |
LeBlanc et al. | A novel PET detector block with three dimensional hit position encoding | |
Kitamura et al. | Preliminary design studies of a high sensitivity small animal DOI-PET scanner: jPET-RD |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20050908 |
|
A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20071228 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20080805 |
|
A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20081003 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20090512 |
|
A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20090708 |
|
TRDD | Decision of grant or rejection written | ||
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20090728 |
|
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 |
|
A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20090824 |
|
R150 | Certificate of patent or registration of utility model |
Ref document number: 4367903 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20120904 Year of fee payment: 3 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20140904 Year of fee payment: 5 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
S533 | Written request for registration of change of name |
Free format text: JAPANESE INTERMEDIATE CODE: R313533 |
|
R350 | Written notification of registration of transfer |
Free format text: JAPANESE INTERMEDIATE CODE: R350 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
S533 | Written request for registration of change of name |
Free format text: JAPANESE INTERMEDIATE CODE: R313533 |
|
R350 | Written notification of registration of transfer |
Free format text: JAPANESE INTERMEDIATE CODE: R350 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
S111 | Request for change of ownership or part of ownership |
Free format text: JAPANESE INTERMEDIATE CODE: R313117 |
|
R350 | Written notification of registration of transfer |
Free format text: JAPANESE INTERMEDIATE CODE: R350 |
|
EXPY | Cancellation because of completion of term |