JP6663094B2

JP6663094B2 - 散乱線強度算出装置及びそれを備えた散乱線補正装置並びに散乱線強度算出方法

Info

Publication number: JP6663094B2
Application number: JP2018137474A
Authority: JP
Inventors: 聖正岡
Original assignee: Asahi Roentgen Industries Co Ltd
Current assignee: Asahi Roentgen Industries Co Ltd
Priority date: 2018-07-23
Filing date: 2018-07-23
Publication date: 2020-03-11
Anticipated expiration: 2038-07-23
Also published as: JP2020014519A; ZA202306810B

Description

本発明は、Ｘ線撮影で得られる画像に含まれる散乱線成分の強度を算出する散乱線強度算出装置及び散乱線強度算出方法に関する。

Ｘ線撮影で得られる画像は、被写体を透過した直接放射線に起因する直接線成分と、直接線成分以外の成分とを含む。直接線成分以外の成分としては、散乱線に起因する散乱線成分と、Ｘ線検出素子の光拡散およびＸ線撮影装置の系から発生する電気的ノイズ等に起因する成分とが存在するが、散乱線成分の方が支配的である。散乱線としては、被写体からの散乱線やフラットパネルディテクタの前面に設置されるカーボンからの散乱線などがある。

特に被写体のＸ線減弱が大きくなる箇所に相当するピクセル位置では、ディテクタに取り込まれるＸ線の大半が散乱放射線となる場合もあり、線形性が大きく崩れ、アーチファクトが現れる要因となる。そして、これがＸ線撮影で得られる画像の画質の劣化を引き起こし、それが誤診の原因となったり、本来治療すべき病変を見逃すことになったりする場合がある。したがって、Ｘ線撮影で得られる画像の画質の劣化を抑えるためには、Ｘ線撮影で得られる画像から散乱線成分を低減した画像を導出する必要がある。

佐久間俊光、外7名、「鉛ディスク法による散乱線測定」、東北大学医療技術短期大学部紀要、東北大学医療技術短期大学部、1995年、第4巻、第2号、p.119―124

非特許文献１では、鉛ディスク法による散乱線測定が提案されている。しかしながら、任意の異なる生体を被写体とした場合でも適用可能な具体的な散乱線強度の算出手法については何ら開示されていない。

本発明は、上記の状況に鑑み、任意の異なる生体を被写体とした場合でも適用可能な散乱線強度算出装置及びそれを備えた散乱線補正装置並びに散乱線強度算出方法を提供することを目的とするものである。

上記目的を達成するために本発明に係る散乱線強度算出装置においては、被写体のＸ線撮影で用いられるディテクタ上の特定のピクセル位置の特定のフレームにおける全放射線強度に対する散乱放射線強度の算出結果の比と、前記特定のフレームでの前記被写体を設けないＸ線撮影で前記特定のピクセル位置で検出される全放射線強度に対する前記特定のフレームでの前記被写体のＸ線撮影で前記特定のピクセル位置で検出される全放射線強度の比との関係を示す曲線式のパラメータを記憶する記憶部と、前記特定のピクセル位置の前記特定のフレームにおいて、前記被写体のＸ線撮影で検出される全放射線強度又は直接放射線強度に対する散乱放射線強度の比を前記曲線式を用いて算出する算出部と、を備える構成（第１の構成）とする。

上記第１の構成の散乱線強度算出装置において、前記特定のピクセル位置が複数設定され、前記特定のフレームが複数設定される構成（第２の構成）としてもよい。また、上記第１又は第２の構成の散乱線強度算出装置において、前記記憶部は、前記ディテクタ上の前記特定以外のピクセル位置における散乱放射線強度を、前記特定のピクセル位置の前記特定のフレームにおける散乱放射線強度の算出結果に基づいて補間する補間式又は前記補間式から算出される補間データを記憶する構成（第３の構成）としてもよい。

上記目的を達成するために本発明に係る散乱線補正装置は、上記第１〜第３いずれかの構成の散乱線強度算出装置と、前記算出部の算出結果に基づいて前記被写体のＸ線撮影で得られる測定画像に対して散乱線補正を行う補正部と、を備える構成（第４の構成）である。

上記目的を達成するために本発明に係る散乱線強度算出方法においては、被写体のＸ線撮影で用いられるディテクタ上の特定のピクセル位置の特定のフレームにおける全放射線強度に対する散乱放射線強度の算出結果の比と、前記特定のフレームでの前記被写体を設けないＸ線撮影で前記特定のピクセル位置で検出される全放射線強度に対する前記特定のフレームでの前記被写体のＸ線撮影で前記特定のピクセル位置で検出される全放射線強度の比との関係を示す曲線式のパラメータを記憶媒体から読み出す読出ステップと、前記特定のピクセル位置の前記特定のフレームにおいて、前記被写体のＸ線撮影で検出される全放射線強度又は直接放射線強度に対する散乱放射線強度の比を前記読出ステップで読み出された前記パラメータによって定まる前記曲線式を用いて算出する算出ステップと、を備える構成（第５の構成）とする。

本発明によると、任意の異なる生体を被写体とした場合でも適用可能な散乱線強度算出装置及びそれを備えた散乱線補正装置並びに散乱線強度算出方法を実現することができる。

鉛円板、Ｘ線源、及び被写体の位置関係を示す図鉛円板の直径を変えた場合の全放射線の輝度値と、鉛円板がない場合の散乱放射線の輝度値とを示すグラフ各近似方法で求めた外挿データの誤差を示す表ディテクタ上の特定のピクセル位置の一例を示す図１０個のプロファイルの位置を示す図１ブロックのピクセルの座標と輝度値に関する模式図既知の輝度値を示す表１ブロック内における輝度値のプロファイルを示す図散乱放射線による輝度値のイメージ画像散乱放射線による輝度値のプロファイルを示す図測定画像上の特定のピクセル位置の一例を示す図鉛円板の直径を変えた場合の全放射線の輝度値と、鉛円板がない場合の散乱放射線の輝度値とを示すグラフアキシャル像を示す図着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフ着目ピクセル領域を示す測定画像および着目ピクセル領域におけるパラメータの計算結果を示すグラフディテクタ上の特定のピクセル位置の一例を示す図全放射線の輝度値と散乱放射線の輝度値とを示すグラフ着目ピクセル領域におけるパラメータの計算結果を示すグラフ４次曲線から計算される散乱放射線による輝度値との差を示すグラフ修正後の全放射線の輝度値と散乱放射線の輝度値とを示すグラフ生体を被写体とした場合のアキシャル像を示す図発明の一実施形態に係る散乱線補正装置の構成を示す図

本発明の実施形態について図面を参照して以下に説明する。本発明に係る散乱線強度算出装置及びそれを備えた散乱線補正装置を開発するために検討及び評価した結果について説明し、最後に本発明に係る散乱線強度算出装置及びそれを備えた散乱線補正装置について説明する。なお、本実施形態では、ディテクタの各画素でフレーム毎に計数される値である画素値を輝度値と呼ぶことにする。輝度値には、直接放射線が寄与する分と、散乱放射線が寄与する分とが含まれている。本実施形態では、前者を直接放射線による輝度値又は直接放射線の輝度値と呼び、後者を散乱放射線による輝度値又は散乱放射線の輝度値と呼ぶことにする。直接放射線による輝度値又は直接放射線の輝度値は直接放射線の強度に対応し、散乱放射線による輝度値又は散乱放射線の輝度値は散乱放射線の強度に対応する。

＜１．ディテクタの各ピクセルで取り込まれる散乱放射線の計算＞
＜１−１．鉛円板の位置に対応するピクセルにおける散乱放射線の輝度値＞
ディテクタの有感領域での散乱放射線の輝度値は実験的に求めることができる。その方法としては、直径の異なる数種類の鉛円板１を図１に示すようにＸ線源２と被写体３の間に設置してＣＴ撮影を行い、鉛円板１の中心位置に対応するディテクタ４のピクセル位置の輝度値を比較して外挿する。

鉛円板は十分な厚さがあればＸ線をほぼ遮断できると考えてよい。したがって、鉛円板の中心位置に対応するディテクタのピクセル位置の輝度値は、直接放射線成分を含まず散乱放射線のみが反映されていると考えることができる。この鉛円板の中心位置に対応するディテクタのピクセル位置の輝度値に反映される散乱放射線は、鉛円板の外側を通過したＸ線が被写体に照射した場合に生成する散乱放射線であると考えられる。

しかし、鉛円板に照射するＸ線は吸収されて被写体に照射しないので、その分散乱放射線が減少することになる。そして、その減少分は鉛円板の直径が小さくなればなるほど少なくなる。つまり、鉛円板の直径が小さいほど散乱放射線が多くなる。そして、鉛円板の直径と散乱放射線との関係を外挿して直径が０のときの散乱放射線の輝度値（外挿データ）を求めると、鉛円板を設置せずに通常のＣＴ撮影を行った場合にそのピクセルで反映される散乱放射線の輝度値であるとみなすことができる。

＜１−２．外挿データに関する検討＞
上記１−１．で説明した外挿データを求める具体的な方法を検討するため、頭部ファントムを模した数学ファントムにＸ線管から放出されるＸ線スペクトルのＸ線束を照射したときのＸ線の挙動を、モンテカルロ法を用いて計算することによりディテクタの各ピクセルで取り込まれる直接放射線や散乱放射線の数を計算した。ディテクタは蓄積型であるので厳密には異なるが、ピクセルに取り込まれるＸ線の数と得られる輝度値とはほぼ比例関係にあると考えることができる。

頭部ファントムは、前歯や臼歯、上顎、下顎、頸椎、軟組織等、形状や構成する組成をできるだけ現実のものに近づけるようにしてジオメトリーを作成した。実際のＣＴ撮影では、患者の余分な被曝を抑えるため、ディテクタの有感領域にのみ直接Ｘ線が入射するようＸ線束をＸ線管直下において鉛コリメータでカットしているので、モンテカルロ法を用いた計算でも同様にカットした。Ｘ線が被写体で散乱せずに取り込まれるべきピクセルで検出された場合は直接放射線とみなし、異なるピクセルで取り込まれた場合は散乱放射線とみなして、直接放射線と散乱放射線をそれぞれ計算
した。直接放射線と散乱放射線の和が、実際のＣＴ撮影においてディテクタで検出される放射線に対応する。ここでは、直接放射線と散乱放射線の和を全放射線とする。

Ｘ線源と頭部ファントムとの間に設置した鉛円板の直径を２ｍｍ、４ｍｍ、８ｍｍ、１２ｍｍ、１６ｍｍ、２０ｍｍ、２４ｍｍと変化させてモンテカルロ法で計算した場合の全放射線の輝度値から外挿して算出した０ｍｍの場合、すなわち、鉛円板を設置しない場合の散乱放射線の輝度値と、モンテカルロ法で計算した０ｍｍの場合の散乱放射線の輝度値とを比較する。

図２は、鉛円板の直径を変えた場合の全放射線の輝度値と、鉛円板がない場合の散乱放射線の輝度値とを、正面照射、側面照射、背面照射それぞれにおいて示すグラフである。前者の値から外挿した値を算出するための近似方法として、（I）軸対称な６次曲線で近似した場合、（II）軸対称な６次曲線と、直径２ｍｍ及び４ｍｍの２点を結ぶ直線との間の平均値で近似した場合、（III）１次式の最小二乗法で近似した場合、の３つのパターンを検討した。

図３は、各照射パターンにおける、各近似方法を用いて得られた鉛円板を設置しない場合の散乱放射線の輝度値と、モンテカルロ法で計算した０ｍｍの場合の散乱放射線の輝度値との誤差を示す表である。

全体的に最も安定して誤差が小さいのは近似方法（II）である。したがって、以下の説明では、近似方法（II）で散乱放射線の輝度値を算出することにした。

しかしながら、近似方法（I）、近似方法（III）、或いは近似方法（I）〜（III）以外の近似方法を採用してよい。なお、近似方法（II）や近似方法（I）のように６次式を用いることで誤差が小さくなるので、６次式を用いた近似を行うことが望ましい。

ディテクタ上の特定のピクセル位置における散乱放射線による輝度値は上記のようにして求めることができる。以下の説明では、図１において設置する鉛円板の位置を、それに対応するディテクタ上のピクセル位置が図４に示す５×５＝２５箇所の黒丸部分のいずれかになるように調整して各撮影を行い、得られる実験データから散乱放射線を算出する。

＜１−３．ディテクタ上の全ピクセルにおける散乱放射線の輝度値の計算＞
ディテクタ上において、各ピクセル位置での散乱放射線はピクセル座標が変化してもそれほど急激に変化せず、緩やかに変化する。図５の２５箇所の散乱放射線の輝度値が分かっているので、これが滑らかに変化するように残りのピクセルにおける散乱放射線の輝度値を計算することが可能である。撮影実験に先立って、計算によりこれを実証したので以下に示す。

上記１−２．と同様にモンテカルロ法を用いて、鉛円板を取り除いたジオメトリーでＸ線を照射した場合にディテクタの各ピクセルで取り込まれる直接放射線と散乱放射線の数を計算する。このとき、横軸または縦軸方向に平行な１列のピクセル上で取り込まれる散乱放射線の輝度値のピクセルナンバーに対するプロファイルを考えると、これらは緩やかな曲線を描く。このプロファイルは、その内５点の値が既知であれば、最小二乗法により４次曲線で近似できる。図４に示す２５箇所の黒丸部分で考えると、縦方向と横方向で計１０個のプロファイルをそれぞれ近似することができる。すると、図５に示す実線上のピクセルにおける散乱放射線の輝度値は４次曲線での近似により算出することができる。これらの輝度値から、その他のピクセルにおける散乱放射線の輝度値を算出することになる。その手順を以下に示す。

ディテクタ上のピクセルは図５に示す実線によって区切られる１６個のブロックのいずれかに必ず属するので、その属するブロック内の座標(x,y)のピクセルにおける散乱放射線の輝度値Bv(x,y)を算出することを考える。この１ブロックのピクセルの座標と輝度値に関する模式図を図６に示す。

ブロックの角の座標をそれぞれ(x₀,y₀)、(x₁,y₀)、(x₀,y₁)、(x₁,y₁)とする。すると、ブロックの角４点の輝度値Bv(x₀,y₀)、Bv(x₁,y₀)、Bv(x₀,y₁)、Bv(x₁,y₁)および、着目ピクセル(x,y)を通り横軸方向および縦軸方向に平行な直線とブロックの境界位置のピクセルにおける輝度値Bv(x,y₀)、Bv(x₁,y)、Bv(x,y₁)、Bv(x₀,y)は、全て上述した４次曲線での近似により輝度値を算出することができるピクセルにおける輝度値であるので、すでに既知である。

ピクセル位置を(x,y₀)から(x,y)、 (x,y₁)へと移動するとき、Bv(x,y₀)、Bv(x,y)、 Bv(x,y₁)への輝度値の変化は、Bv(x₀,y₀)、Bv(x₀,y)、Bv(x₀,y₁)への輝度値の変化、および、Bv(x₁,y₀)、Bv(x₁,y)、Bv(x₁,y₁)への輝度値の変化の影響を受けることになる。このとき、輝度値は緩やかに変化するので、輝度値の比Bv(x,y)/Bv(x₀,y)は、Bv(x,y₀)/Bv(x₀,y₀)とBv(x,y₁)/Bv(x₀,y₁)の間の値を緩やかに変化すると考えるのが妥当である。

ここで、y 座標の変化に伴って輝度値が単調に変化すると仮定すると、輝度値の比Bv(x,y)/Bv(x₀,y)を求めることができ、Bv(x₀,y)が既知なので、Bv(x,y)を計算することができる。これをBvs1(x,y) とする。一方、Bv(x,y)/Bv(x₁,y)はBv(x,y₀)/Bv(x₁,y₀)とBv(x,y₁)/Bv(x₁,y₁)の間を緩やかに変化すると考えることもできるので、このことからもBv(x,y)を計算することができ、これをBvl1(x,y) とする。一般にBvs1(x,y) とBvl1(x,y)は一致しないので、これをBv(x,y)の仮の値としてBv1(x,y)とした。

しかし、これをブロック内のすべてのピクセルに適用しようとすると、ブロックの境界部分（x= x₀およびx= x₁ のとき）で輝度値の不整合が生じる。そこで、x= x₀およびx= x₁のときのBv(x,y)の値をそれぞれBv00(y), Bv01(y) とする。本来、Bv1(x₀,y) はBv00(y)に、Bv1(x₁,y) はBv01(y)に一致しなければならないので、x= x₀ およびx= x₁ でそれぞれ一致するよう以下のようにBv1(x,y) を補正して、その補正した値をBv(x,y)とする。

ここで、５１ピクセル×５１ピクセルのディテクタで２次元のピクセル座標(x,y)にそれぞれ０〜５０のピクセル番号を与える。検出される輝度値の内、x=0,x=50およびy=0, y=50のピクセルライン上のピクセルについてはすべて既知であり、図７に当てはめて考えると、x₀=0, x₁=50, y₀=0, y₁=50である。例えば、x=25, y=20 として、既知の輝度値を図７に示す値とする。

図７に示す既知の輝度値からBv1 を計算したところ、1139 であった。そして、式（２）を適用してBv(25,20)の値を算出したところ、1135となった。また、ピクセルラインで囲まれた内部のピクセルの輝度値を同様に計算し、ピクセルの位置を水平な平面上に、輝度値を高さ方向に取って３次元で表わすと、プロファイルは図８のようになった。図８に子示すプロファイルの位置による変化は違和感のないものとなっている。

このようにして求めた輝度値に関しても、頭部ファントムを模した数学ファントムを被写体としてモンテカルロ法で計算した散乱放射線による輝度値との誤差を見積もった。２５箇所のピクセル領域で、モンテカルロ法で計算した散乱放射線による輝度値を利用して、縦軸方向・横軸方向それぞれ５本ずつ、計１０本のピクセルラインの輝度値のプロファイルを４次曲線で近似し、それ以外のピクセルの輝度値を上記で述べた方法で式（２）より算出する。

そして、上述した１０本のピクセルライン以外のピクセルの輝度値を式（２）で算出したプロファイル（以下、「第１プロファイル」と称す）を、すでにモンテカルロ法で計算している全ピクセルの輝度値のプロファイル（以下、「第２プロファイル」と称す）と比較する。図９（ａ）は第２プロファイルから得られる正面照射の場合の散乱放射線による輝度値のイメージ画像であり、図９（ｂ）は第１プロファイルから得られる正面照射の場合の散乱放射線による輝度値のイメージ画像である。両者を比較すると、図９（ｂ）は図９（ａ）をよく再現していることが分かる。

第１プロファイル及び第２プロファイルに対して、正面照射、側面照射、背面照射の場合のY=350における散乱放射線による輝度値のプロファイルを調べたところそれぞれ図１０（ａ）、図１０（ｂ）、図１０（ｃ）のようになった。図１０には、第２プロファイルをスムージングした第３プロファイルも合わせて示している。第２プロファイルと第１プロファイルの平均誤差をそれぞれ計算したところ、4.09%, 4.02%, 6.44%となった。一方、第３プロファイルと第１プロファイルの平均誤差はそれぞれ0.86%, 2.76%, 3.03%となり、モンテカルロ法で計算した輝度値のピクセル間のバラツキが誤差に影響を与えていたと思われる。Y=300, 400についても同様に誤差を計算したところ、表１のようになり、同様の傾向を示している。ここで、表１中のMCは第２プロファイル、Smは第３プロファイル、F(1)は第１プロファイルを示している。Y=400の背面照射は誤差が大きいのは歯列部分にあまり影響を与えない両端部分のずれが大きいためである。中央付近の350ピクセルのみで誤差を計算したところ1.42%であったので、ＣＴ再構成画像への影響は少ない。

＜２．実験方法＞
Ｘ線撮影装置は朝日レントゲン工業株式会社製のAUJE SOLIO を使用した。歯科用Ｘ線ＣＴ装置では、被写体の周りをＸ線管とディテクタがそれぞれ両端に設置されたアームが３６０度回転する。フレーム数は５１０枚、フレームレートは３０ｆｐｓである。ディテクタの有感領域は１５６．９６ｍｍ×１５９．３６ｍｍ（６５４ピクセル×６６４ピクセル）であり、ディテクタの素子はＣｓＩ、厚さは７００μｍである。焦点からディテクタまでの距離が６００ｍｍ、焦点から回転中心までの距離が３７５ｍｍとなる。Ｘ線管の管電圧は８５ｋＶ、管電流は４ｍＡとした。

撮影は、頭部ファントムを最適な位置に設置した場合、横方向に−５ｍｍ、−３ｍｍ、＋３ｍｍ、＋５ｍｍずらした場合、前後方向に−５ｍｍ、＋５ｍｍずらした場合、上下方向に−５ｍｍ、＋５ｍｍずらした場合の９パターン行った。実際はファントムを固定してアームの中心位置および高さを調整することで対応した。

それぞれのパターンについて鉛円板を設置しない通常の撮影、図４に示す２５箇所に対応する位置に中心位置がくるように異なる直径の鉛円板をそれぞれ1 つずつ設置した場合の撮影を行った。鉛円板の直径は、２．３ｍｍ、４．７ｍｍ、８．２ｍｍ、１１．８ｍｍ、１５．５ｍｍの５種類である。したがって、撮影は１パターンで５×２５＋１＝１２６回、全部で９×１２６＝１１３４回行った。それぞれの撮影に対して、５１０枚（フレーム数分）の測定画像のデータがある。

＜３．実験結果＞
上述した＜１．ディテクタの各ピクセルで取り込まれる散乱放射線の計算＞で説明した方法に基づいて、各パターン、各鉛円板位置、各フレームの撮影データについて、鉛円板の直径の変化に対する、鉛円板位置に対応するピクセル領域の平均輝度値をグラフ化して、直径が０の場合の輝度値を外挿することによって、有感領域内の２５箇所の散乱放射線による輝度値を算出した。これらの内、頭部ファントムを最適な位置に設置したパターンにおいて、図１１(i)〜図１１(iv)に示すフレームNo.1（正面照射）, 128（側面照射）, 256（背面照射）, 384（側面照射）の測定画像の４箇所のピクセル領域ａ〜ｄについて、鉛円板の直径を横軸、輝度値を縦軸に取った場合のグラフを図１２（ａ）〜図１２（ｄ）に示す。ただし、直径０ｍｍのときの輝度値は、直径２．３ｍｍから１５．５ｍｍのときの５点の輝度値を用いて外挿により算出した値である。また、図１１（i）はｃの位置、図１１（ii）はｄの位置、図１１（iii）はａの位置、図１１（iv）はｂの位置に鉛円板が該当する場合の測定画像を例として示している。なお、暗電流による輝度値は差し引いている。

図１２（ａ）〜図１２（ｄ）に示すグラフでは、被写体の厚みが薄く測定画像の輝度値が高い箇所では散乱放射線の輝度値が大きく、被写体が厚い部分では逆に小さくなっており、散乱放射線の特徴が現れている。また、いずれも外挿した０ｍｍのときの輝度値と撮影で得られた他の輝度値とは滑らかな曲線を描いており、この方法で妥当な散乱放射線の輝度値が算出されていると考えてよい。そして、頭部ファントムの位置がずれた他のパターンについても同様にディテクタの有感領域内における散乱放射線の輝度値を算出することができる。

＜４．散乱線補正の効果＞
上述した方法で算出した散乱放射線の輝度値を用いて散乱線補正を行うことによって、ＣＴ再構成画像の画質がどの程度向上するのかを確認する必要がある。

そこで、頭部ファントムが最適な位置に設置されているパターンで撮影して得られるデータから、各フレームについて２５箇所のピクセル領域の散乱放射線の輝度値を算出した。そして、フレームごとに、式（２）を用いてディテクタのピクセル全体の散乱放射線の輝度値を算出して、撮影で得られた測定画像の各ピクセルにおける放射線（直接放射線＋散乱放射線）の輝度値からその散乱放射線の輝度値を差し引くことによって、散乱線補正を実行した。

そして、これら処理後の測定画像を再構成した三次元ＣＴ画像と、特開２０１５−６５９９３号公報で提案されている散乱線補正を適用した画像とを比較した。その結果を以下に示す。

図１３（ａ）は測定画像を、散乱線補正を行わずに画像再構成した場合の、アキシャル画像を１箇所ピックアップして示したものである。

図１３（ｂ）は測定画像に対して、特開２０１５−６５９９３号公報で提案されている散乱線補正を適用した後に画像再構成した場合の、アキシャル画像を１箇所ピックアップして示したものである。

図１３（ｃ）は測定画像に対して、上述した方法で算出した散乱放射線の輝度値を用いて散乱線補正を行った後に画像再構成した場合の、アキシャル画像を１箇所ピックアップして示したものである。

図１３（ｂ）では歯列の内側や歯間に沿った方向や臼歯間にアーチファクトが現れており、図１３（ａ）と比較してそれほど変化がなくあまり画質が改善されていないが、図１３（ｃ）ではアーチファクトが解消されている。すなわち、上述した方法で算出した散乱放射線の輝度値を用いた散乱線補正が有効であることが分かる。

＜５．ＣＴ撮影における頭部ファントムの異なった位置づけで得られる散乱放射線の輝度値の割合の曲線近似＞
＜５−１．近似曲線の算出方法＞
上述した散乱線補正を、頭部ファントムの任意の位置づけ、さらには、任意の生体データにも適用できるようにするためには、鉛円板の位置に対応する25箇所のピクセル領域のフレームごとに得られる散乱放射線の輝度値について、計算式を用いて算出できるようにする必要がある。本節では、頭部ファントムの９種類の位置づけパターンにおいて得られるデータを用いて曲線の計算式を導き、２５箇所のピクセル領域のフレームごとのパラメータを算出したので、その方法を以下に示す。なお、以下で述べるピクセルの輝度値について、暗電流による影響はすでに取り除いているものとする。

曲線の横軸は、被写体を設置しないで撮影を行った場合（このとき得られる測定画像を白画像と呼ぶことにする）に各ピクセルで得られる輝度値に対する、頭部ファントムを設置して撮影した場合（直接放射線も散乱放射線も含んでいるので、このとき得られる放射線を全放射線と呼ぶことにする。また、鉛円板は設置しない）の輝度値の比とし、それをＸとする。また、曲線の縦軸は、全放射線の輝度値に対する、上記＜３．実験結果＞
において算出した散乱放射線の輝度値の比とし、それをＹとする。

すると、ＸとＹはともに１より小さい値をもち、被写体が薄い経路を通るときはＸが１に近づき、逆に厚い場合は０に近づく。散乱放射線の割合が大きいときはＹが１に近づき、小さいときは０に近づく。頭部ファントムのような大きい被写体においては、被写体が厚いと散乱放射線の割合が大きく、薄いと小さくなる傾向がある。以上の傾向に加え、グラフ化した場合のこれらの分布を考慮し、ＸとＹの関係を表す曲線の式を、パラメータａおよびパラメータｂを用いた下記の式（３）とした。

前述したように、これを各ピクセル領域の各フレームで得られる9 つのデータ分布に対して適用する。式（３）は、両辺に対数をとってＹ_log＝logＹ、Ｘ_log＝logＸとすることにより、下記の式（４）に変形でき、Ｙ_logをＸ_logの１次式で表わすことができる。

このため、最小二乗法を適用することによって、各ピクセル領域のフレームごとに式（３）で近似して、そのパラメータａの値およびパラメータｂの値を算出することができる。したがって、Ｘの値とＹの値が式（３）で一意的に決まるので、任意のＸの値に対してＹの値を求めることができる。すなわち、特定のピクセル領域の特定のフレームにおいて、撮影によって得られる任意の輝度値に対して散乱放射線による輝度値を一意的に算出することができる。

＜５−２．最小二乗法を用いたパラメータの計算＞
撮影した９パターンの頭部ファントムの位置づけに対して得られたデータ群すべてを式（３）で近似して、それぞれについてパラメータａの値およびパラメータｂの値を決め、データ群と近似曲線とを比較した。これらの内、特徴的な曲線のグラフについて取り上げ、撮影で得られたデータと曲線の間の誤差を計算したので以下に示す。

図１４〜図２５の（ａ）は、着目したフレームの測定画像の内頭部ファントムの位置づけが最適な位置となっているパターンの測定画像を示し、着目したピクセル領域を四角のマークで示している。

図１４〜図２５の（ｂ）は、その着目したフレームの着目したピクセル領域におけるＸとＹの分布と曲線の近似式をグラフ化したものである。

図１４から図１８までは、歯列付近のピクセル領域の曲線近似を表している。図１５、図１８では、頭部ファントムの位置がずれることで透過Ｘ線の経路上の骨の厚さが変化するためＸ線減弱が大きく変化しており、Ｘの変化量が大きい。しかし、これらのデータ群の分布は式（３）で表わされる曲線で近似することができ、よく合っていることがグラフからも分かる。その誤差は５％程度以内となっている。誤差の精度は、データ数が増えればさらによくなっていくことが期待される。

図１６、図１７の着目ピクセル領域は、頭部ファントムの厚みが大きい部分であり、この部分はファントムの位置づけがずれてもＸ線の透過経路におけるファントムの厚みはほとんど変わらない。このため、ＸもＹも変化量が小さい。したがって、散乱放射線による輝度値もほとんど変化しない。

図１４の（ｂ）は、ちょうど歯列が重なっている部分とその周辺に位置する場合の曲線近似であり、重なり位置となった場合は散乱放射線の割合が非常に高くなる。しかし、この場合でもデータ群は、曲線に沿った形で分布している。

図１９〜図２２の（ｂ）は、歯列付近から離れたピクセル領域の曲線近似を表している。これらについても、いずれのデータ群も曲線に沿った分布をしており、誤差は歯列付近よりもよい値が得られている。ピクセルラインに沿った輝度値の変化が歯列付近よりも緩やかであるため、散乱放射線の割合の変化も少ないと考えられる。

図２３〜図２５の（ｂ）は、ディテクタの端の、頭部ファントムの厚みが非常に薄い領域の曲線近似を表している。これらについても、いずれのデータ群も曲線に沿った分布をしている。

＜６．臼歯の重なり部分等、Ｘ線減弱が大きい部分の取扱い＞
歯列部分、特に歯の重なり部分などＸ線減弱の大きい部分に対応するディテクタ上のピクセルにおいては、Ｘ線減弱が小さい場合よりも散乱放射線による輝度値が小さい。すると、任意の生体データに適用する場合、実際の投影画像においては、歯の重なり部分が２５箇所のピクセル領域の間に位置する場合などは、上記の＜１−３．ディテクタ上の全ピクセルにおける散乱放射線の輝度値の計算＞で求めた４次曲線によるプロファイルでは散乱放射線による輝度値が実際よりも大きい値で計算される場合がある。そうすると、散乱線補正後の歯列部分の輝度値が少なくなり、ＣＴ再構成画像に少なからず影響を与えるので、散乱放射線の輝度値を修正する必要がある。この修正方法の一例を以下に示す。

図２６に１フレームの投影画像と２５箇所の領域が示されている。歯列部分に相当するのは点７と点１２の間、点８と点１３の間、点９と点１４の間である。図２７は図２６中の点２、点７、点１２、点１７、点２２を通るピクセルライン上における全放射線の輝度値と４次曲線から計算した散乱放射線の輝度値のプロファイルを図示したものである。

図２７より、ピクセルNo.200 から300 の間（この部分は点７と点１２の間に該当する）で双方の輝度値がほぼ等しくなっているピクセルがあることが分かる。これは直接放射線が検出されておらず、Ｘ線が透過していないことになり、金属などＸ線減弱係数の大きい物質がなければこのようなことは起こらないので矛盾する。

したがって、この部分の散乱線補正を厳密に行うため、上記の＜１−３．ディテクタ上の全ピクセルにおける散乱放射線の輝度値の計算＞に基づいて縦軸方向および横軸方向５本ずつ計１０本に関して４次曲線の計算を行った後に、上記２点間の輝度値のプロファイルの計算を追加でそれぞれ行い、散乱放射線の輝度値を修正する必要がある。

この各２点間のピクセルにおいても、鉛円板を用いた実験を行うことは可能である。しかし、矢印部分の領域においては、上記の式（３）のパラメータはあまり変化しないので、このことを利用する。例えば、フレームNo. 40 における点７と点１２の両方について、白画像に対する全放射線の輝度値の比を横軸に、全放射線に対する散乱放射線の輝度値の比を縦軸にとって図示したグラフ（図１４〜図２５のグラフと同等）を図２８に示す。これらはいずれも１本の曲線に沿って分布しており、式（３）のパラメータは共通の値を適用できることが分かる。この曲線で算出される値を利用して、図２７の散乱放射線による輝度値のプロファイルを修正する。

点７と点１２の間における散乱放射線と全放射線の輝度値の関係がこの曲線で表わされると仮定して、この２点間を４等分する３点を取り、３点それぞれについて、散乱放射線による輝度値を算出する。そして、この算出した輝度値が正しいとする。

一方、上記の＜１−３．ディテクタ上の全ピクセルにおける散乱放射線の輝度値の計算＞に基づいて計算される４次曲線からすでに散乱放射線による輝度値を計算しているので、両者の値を比較してその差を求め、３点から得られる輝度値差のプロファイルを４次曲線から計算した散乱放射線による輝度値のプロファイルから差し引くことによって修正する。

まず、３点における散乱放射線の輝度値差の大小関係は、
（１）ピクセルナンバー順に小さくなる、
（２）ピクセルナンバー順に大きくなる、
（３）真ん中のピクセルが最大となる、
（４）真ん中のピクセルが最小となる、
の４パターンに大別できる。ピクセルナンバーに対するこの輝度値の差のプロファイルが滑らかになるようにするため、Sinカーブを採用する。それを図示したものを図２９（ａ）〜（ｄ）に示す。横軸がピクセルナンバーで縦軸が上記２種類の方法で算出した散乱放射線の輝度値の差である。図のピクセルナンバー0が点７に該当し、ピクセルナンバー165 が点１２に該当する。これら２つのピクセルにおける輝度値差は０であり、プロファイルの両端は０からゆっくり大きくなるので、点７および点１２におけるSin カーブの位相は−π/2＋2nπ（n; 整数)である。点７と点１２を４等分した各領域を順に領域I、領域II、領域III、領域IVとする（図２９（ａ）参照）。

パターン（１）については、領域Iで輝度値差が０から最大値に達するので、その２分の１を振幅としてSin カーブの位相が−π/2 から＋π/2 まで変化するようにする。領域II〜IVは、各領域で位相が＋π/2＋2nπから(3/2)π＋2nπに変化するよう位相の定数のパラメータを調整する。振幅は各領域の輝度値差の２分の１である。

パターン（２）については、ピクセルに対するプロファイルがパターン（１）と左右逆になるが、上記と同様に修正できる。

パターン（３）については、領域I及びIIで−π/2＋2nπからπ/2＋2nπまで、領域III及びIVではπ/2＋2nπから(3/2)π＋2nπまで位相が変化するようパラメータを決める。

パターン（４）については、領域Iで−π/2 からπ/2、領域IIでπ/2 から(3/2)π、領域IIIで(3/2)πから(5/2)π、領域IVで(5/2)πから(7/2)πに位相を変化させ、同様に修正する。

３点の輝度値の差の大小関係から上記パターンの内のいずれに該当するかを判断し、輝度値差のプロファイルを計算して、４次曲線からすでに計算している散乱放射線による輝度値からこの輝度値差を差し引くことによって、より理想に近い散乱放射線の輝度値を計算することができる。

このように、散乱放射線の輝度値の計算に処理を追加した結果、図２７のプロファイルは図３０のようになった。これにより歯列部分で散乱線補正後の輝度値が極端に小さくなるということは回避され、散乱放射線のプロファイルも滑らかに変化していることが分かる。

＜７．生体を被写体とした場合における散乱線補正の効果＞
図３１（ａ）は生体を被写体とした測定画像を、散乱線補正を行わずに画像再構成した場合の、アキシャル画像を１箇所ピックアップして示したものである。

図３１（ｂ）は生体を被写体とした測定画像に対して、特開２０１５−６５９９３号公報で提案されている散乱線補正を適用した後に画像再構成した場合の、アキシャル画像を１箇所ピックアップして示したものである。

図３１（ｃ）は生体を被写体とした測定画像に対して、上述した＜１．ディテクタの各ピクセルで取り込まれる散乱放射線の計算＞で説明した方法で算出した散乱放射線の輝度値を用いて散乱線補正を行った後に画像再構成した場合の、アキシャル画像を１箇所ピックアップして示したものである。

図３１（ｂ）では歯列の内側や歯間に沿った方向や臼歯間にアーチファクトが現れており、図３１（ａ）と比較してそれほど変化がなくあまり画質が改善されていないが、図３１（ｃ）ではアーチファクトが解消されている。すなわち、上述した＜１．ディテクタの各ピクセルで取り込まれる散乱放射線の計算＞で説明した方法で算出した散乱放射線の輝度値を用いた散乱線補正は、生体を被写体とした場合でも有効であることが分かる。

＜８．散乱線補正装置＞
図３２は、本発明の一実施形態に係る散乱線補正装置１００の構成を示す図である。

散乱線補正装置１００は、ＲＯＭ１２やＨＤＤ１７に格納されているプログラムに従って散乱線補正装置１００全体を制御するＣＰＵ１１と、固定的なプログラムやデータを記録するＲＯＭ１２と、作業メモリを提供するＲＡＭ１３と、外部との通信を行うための通信インターフェース部１４と、画像データを一時的に記憶するＶＲＡＭ１５と、ＶＲＡＭ１５に記憶された画像データに基づいて画像を表示する表示部１６と、詳細は後述するＨＤＤ１７と、キーボード、ポインティングデバイス等の操作部１８とを備えている。

散乱線補正装置１００の通信インターフェース部１４と外部との通信方法は、有線通信でもよく、無線通信でもよく、有線と無線を組み合わせた通信であってもよい。散乱線補正装置１００としては、例えば、パーソナルコンピュータを挙げることができる。

ＨＤＤ１７は、画像再構成処理プログラム、散乱線補正プログラム等の各種プログラム、及び、歯科用Ｘ線撮影装置によるＣＴ撮影の撮影データ（測定画像）、上記画像再構成処理プログラムによって生成された再構成ボリュームデータ、各種プログラムを実行する際に用いられる各種パラメータの設定値等の各種データを記憶する。

また、散乱線補正装置１００は、ディテクタからの情報を取得して測定画像を生成し、その生成した測定画像をＨＤＤ１７に記憶してもよく、通信インターフェース部１４を用いて外部から測定画像を取得し、その取得した測定画像をＨＤＤ１７に記憶してもよい。外部の装置がＣＴ再構成をおこなってもよい。この場合、本実施形態とは異なりＨＤＤ１７が画像再構成処理プログラムを記憶しない構成にすることが可能であり、散乱線補正装置１００が散乱線補正後の測定画像を、ＣＴ再構成を行う外部の装置に送信すればよい。

画像再構成処理プログラムは、歯科用Ｘ線撮影装置によるＣＴ撮影の撮影データを再構成して再構成ボリュームデータを生成するためのプログラムである。散乱線補正プログラムは、散乱線補正を行うためのプログラムである。なお、散乱線補正プログラムは、散乱線強度を算出する散乱線強度算出プログラムを含む。

ＨＤＤ１７に記憶されている各プログラムは、散乱線補正装置１００にプリインストールされていてもよく、光ディスク等の記憶媒体に格納された形態で流通されて散乱線補正装置１００にインストールされてもよく、ネットワークを介して流通されて散乱線補正装置１００にインストールされてもよい。ＨＤＤ１７に記憶されている各プログラム及び各データの一部を、ＨＤＤ１７ではなくＲＯＭ１２に記憶するようにしてもよい。

ＨＤＤ１７は、被写体のＸ線撮影で用いられるディテクタ上の特定のピクセル位置（例えば図４に示す２５箇所の黒丸部分）の特定のフレーム（例えば上述した５１０個のフレーム）における全放射線強度に対する散乱放射線強度の算出結果（例えば上述した外挿データ）の比と、特定のフレームでの被写体を設けないＸ線撮影で特定のピクセル位置で検出される全放射線強度に対する特定のフレームでの被写体のＸ線撮影で前記特定のピクセル位置で検出される全放射線強度の比との関係を示す曲線式のパラメータ（例えば上述したパラメータａ及びパラメータｂ）を記憶する。ＨＤＤ１７に記憶されている散乱線補正プログラムは、ＣＰＵ１１で実行されることで散乱線補正装置１００のハードウェアを、特定のピクセル位置の特定のフレームにおいて、被写体のＸ線撮影で検出される全放射線強度又は直接放射線強度に対する散乱放射線強度の比を算出する算出部、及び、当該記算出部の算出結果に基づいて測定画像に対して散乱線補正を行う補正部、として機能させる。さらに、ＨＤＤ１７は、例えば上述した式（２）又は当該式（２）から算出される補間データを記憶してもよい。なお、ＨＤＤ１７に記憶されている散乱線強度算出プログラムは、ＣＰＵ１１で実行されることで散乱線補正装置１００のハードウェアを、上記の算出部として機能させる。すなわち、散乱線補正プログラムのうち散乱線強度算出プログラムのみを実行した場合、散乱線補正装置１００は、散乱線強度を算出する散乱線強度算出装置となる。

以上、本発明の一実施形態について説明したが、本発明の範囲はこれに限定されるものではなく、発明の主旨を逸脱しない範囲で種々の変更を加えて実施することができる。

例えば、鉛円板の代わりに、円形状以外の形状（例えば矩形形状）である鉛板材を用いてもよい。また、鉛以外のＸ線吸収材を用いてもよい。また、上述した実施形態のように散乱放射線の外挿データを求めるのではなく、微小な貫通孔を有するＸ線吸収材を用いて特定のピクセル位置の直接放射線を求めるようにしてもよい。

例えば、上述した実施形態では、頭部ファントムにインプラント等の金属を含めていないが、頭部ファントムにインプラント等の金属を含めてもよい。

１００散乱線補正装置
１１ＣＰＵ
１２ＲＯＭ
１３ＲＡＭ
１４通信インターフェース部
１５ＶＲＡＭ
１６表示部
１７ＨＤＤ
１８操作部

Claims

被写体のＸ線撮影で用いられるディテクタ上の特定のピクセル位置の特定のフレームにおける全放射線強度に対する散乱放射線強度の算出結果の比である第１の比と、前記特定のフレームでの前記被写体を設けないＸ線撮影で前記特定のピクセル位置で検出される全放射線強度に対する前記特定のフレームでの前記被写体のＸ線撮影で前記特定のピクセル位置で検出される全放射線強度の比である第２の比との関係を示す曲線式のパラメータを記憶する記憶部と、
前記特定のピクセル位置の前記特定のフレームにおいて、前記被写体のＸ線撮影で検出される全放射線強度又は直接放射線強度に対する散乱放射線強度の比である第３の比を前記曲線式を用いて算出する算出部と、
を備え、
前記第１の比は、前記特定のピクセル位置の前記特定のフレームにおける散乱放射線の輝度値の外挿データを用いて算出される、散乱線強度算出装置。
被写体のＸ線撮影で用いられるディテクタ上の特定のピクセル位置の特定のフレームにおける全放射線強度に対する散乱放射線強度の算出結果の比である第１の比と、前記特定のフレームでの前記被写体を設けないＸ線撮影で前記特定のピクセル位置で検出される全放射線強度に対する前記特定のフレームでの前記被写体のＸ線撮影で前記特定のピクセル位置で検出される全放射線強度の比である第２の比との関係を示す曲線式のパラメータを記憶する記憶部と、
前記特定のピクセル位置の前記特定のフレームにおいて、前記被写体のＸ線撮影で検出される全放射線強度又は直接放射線強度に対する散乱放射線強度の比である第３の比を前記曲線式を用いて算出する算出部と、
を備え、
前記曲線式はＹ＝ａＸ ^ｂ（ただし、前記Ｙは前記第１の比であり、前記Ｘは前記第２の比であり、前記ａ及び前記ｂはそれぞれ前記パラメータである）で表される、散乱線強度算出装置。
前記特定のピクセル位置が複数設定され、前記特定のフレームが複数設定される、請求項１又は請求項２に記載の散乱線強度算出装置。
前記記憶部は、前記ディテクタ上の前記特定以外のピクセル位置における散乱放射線強度を、前記特定のピクセル位置の前記特定のフレームにおける散乱放射線強度の算出結果に基づいて補間する補間式又は前記補間式から算出される補間データを記憶する、請求項１〜３のいずれか一項に記載の散乱線強度算出装置。
請求項１〜４のいずれか一項に記載の散乱線強度算出装置と、
前記算出部の算出結果に基づいて前記被写体のＸ線撮影で得られる測定画像に対して散乱線補正を行う補正部と、
を備える、散乱線補正装置。
被写体のＸ線撮影で用いられるディテクタ上の特定のピクセル位置の特定のフレームにおける全放射線強度に対する散乱放射線強度の算出結果の比である第１の比と、前記特定のフレームでの前記被写体を設けないＸ線撮影で前記特定のピクセル位置で検出される全放射線強度に対する前記特定のフレームでの前記被写体のＸ線撮影で前記特定のピクセル位置で検出される全放射線強度の比である第２の比との関係を示す曲線式のパラメータを記憶媒体から読み出す読出ステップと、
前記特定のピクセル位置の前記特定のフレームにおいて、前記被写体のＸ線撮影で検出される全放射線強度又は直接放射線強度に対する散乱放射線強度の比である第３の比を前記読出ステップで読み出された前記パラメータによって定まる前記曲線式を用いて算出する算出ステップと、
を備え、
前記第１の比は、前記特定のピクセル位置の前記特定のフレームにおける散乱放射線の輝度値の外挿データを用いて算出される、散乱線強度算出方法。
被写体のＸ線撮影で用いられるディテクタ上の特定のピクセル位置の特定のフレームにおける全放射線強度に対する散乱放射線強度の算出結果の比である第１の比と、前記特定のフレームでの前記被写体を設けないＸ線撮影で前記特定のピクセル位置で検出される全放射線強度に対する前記特定のフレームでの前記被写体のＸ線撮影で前記特定のピクセル位置で検出される全放射線強度の比である第２の比との関係を示す曲線式のパラメータを記憶媒体から読み出す読出ステップと、
前記特定のピクセル位置の前記特定のフレームにおいて、前記被写体のＸ線撮影で検出される全放射線強度又は直接放射線強度に対する散乱放射線強度の比である第３の比を前記読出ステップで読み出された前記パラメータによって定まる前記曲線式を用いて算出する算出ステップと、
を備え、
前記曲線式はＹ＝ａＸ ^ｂ（ただし、前記Ｙは前記第１の比であり、前記Ｘは前記第２の比であり、前記ａ及び前記ｂはそれぞれ前記パラメータである）で表される、散乱線強度算出方法。