JP7161297B2 - 放射線撮影システム、画像処理装置、撮影制御装置及びその作動方法、プログラム - Google Patents

Description

本発明は、放射線撮影システム、画像処理装置、撮影制御装置及びその作動方法、プログラムに関する。

現在、放射線による医療画像診断や非破壊検査に用いる撮影装置として、半導体材料によって形成された平面検出器（Ｆｌａｔ Ｐａｎｅｌ Ｄｅｔｅｃｔｏｒ、以下ＦＰＤと略す）を用いた放射線撮影装置が普及している。このような放射線撮影装置は、例えば医療画像診断においては、一般撮影のような静止画撮影や、透視撮影のような動画撮影のデジタル撮影装置として用いられている。

ＦＰＤを用いた動画撮影として、造影剤を注入しながら撮影を行う血管造影や、カテーテルやステント等の医療デバイスを用いた血管内治療が盛んに行われている。この造影の方法の一つに、デジタル差分造影法（以下、ＤＳＡ（digital subtraction angiography））がある。ＤＳＡでは、造影剤の注入前に撮影を行ってマスク画像を取得した後、造影剤の注入後に撮影を行ってライブ画像を取得する。マスク画像とライブ画像を用いた演算により、造影剤のみを画像化して血管の診断を行うことができる。

一方、ＦＰＤを用いた撮影方法のひとつに、エネルギーサブトラクションがある。エネルギーサブトラクションでは、まず、管電圧の異なる放射線を複数回照射するなどして、エネルギーの異なる画像を複数枚取得する。それらを演算することで、骨画像と軟部組織画像に分離するなどの処理を行うことができる。これにより造影剤のみを画像化することで、動く対象の血管を診断する方法が提案されている（特許文献１）。

特開昭５８－２２１５８０号公報

特許文献１で生成された画像を用いてＤＳＡを適用すると、マスク画像に動く物体（ステント、石灰化など）がある場合、ライブ画像とマスク画像の差分を取った際にアーチファクトが生じてしまう。

本発明は、差分による画像生成においてアーチファクトの発生をより低減することを目的とする。

本発明の一態様による放射線撮影装置は以下の構成を備える。すなわち、
被写体を透過した互いに異なる放射線エネルギーによる複数の画像をそれぞれ動画撮影により取得する取得手段と、
前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内である第１の時間において得られた複数の画像を用いて第１のエネルギーサブトラクション画像を生成し、前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内であり且つ前記第１の時間の後である第２の時間において得られた複数の画像を用いて第２のエネルギーサブトラクション画像を生成する第１の生成手段であって、前記第１のエネルギーサブトラクション画像と前記第２のエネルギーサブトラクション画像とを含む複数のエネルギーサブトラクション画像であって、前記第１の時間から前記第２の時間までの期間において得られた前記複数のエネルギーサブトラクション画像を前記動画撮影の期間内に生成する前記第１の生成手段と、
前記複数のエネルギーサブトラクション画像に基づくマスク画像と、前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内であり且つ前記第１の時間及び前記第２の時間の後である第３の時間において得られた複数の画像を用いて生成された第３のエネルギーサブトラクション画像に基づくライブ画像との差分画像を前記動画撮影の期間内に生成する第２の生成手段と、を備える。

本発明によれば、差分による画像生成においてアーチファクトの発生をより低減することができる。

実施形態による放射線撮影システムの構成例を示すブロック図。 放射線撮影装置へ供給される信号および制御用コンピュータの機能構成を示す図。 実施形態の放射線撮影装置における画素の等価回路図。 放射線エネルギーの異なる画像の時分割による取得処理を示すタイミングチャート。 放射線エネルギーの異なる画像の時分割による取得処理を示すタイミングチャート。 放射線画像の補正処理を説明する図。 （ａ）第１実施形態による物質分離に係る信号処理を示すブロック図、（ｂ）第２実施形態による物質識別に係る信号処理を示すブロック図、（ｃ）第２実施形態に係るカラー画像生成の画像処理を示すブロック図。 第１実施形態によるデジタル差分血管造影法による処理を示すフローチャート。 第１実施形態によるデジタル差分血管造影法による処理を説明する模式図。 第２実施形態によるデジタル差分血管造影法による処理を説明する模式図。 第３実施形態によるデジタル差分血管造影法による処理を示すフローチャート。 第２実施形態によるデジタル差分血管造影法による処理を説明する模式図。 第４実施形態に係るテーブル生成の処理を示すフローチャート。 第４実施形態に係るテーブル参照を説明する模式図。

以下、本発明の実施形態について図面を参照しながら詳細に説明する。なお、以下の実施形態において、放射線という用語は、Ｘ線の他、例えば、α線、β線、γ線粒子線、宇宙線などを含み得る。

＜第１実施形態＞
図１Ａは、第１実施形態に係る放射線撮影システム１００の構成例を示すブロック図である。第１実施形態の放射線撮影システム１００は、放射線発生装置１０１、放射線制御装置１０２、制御用コンピュータ１０３、放射線撮影装置１０４を備える。

放射線発生装置１０１は放射線制御装置１０２の制御下で放射線を発生する。制御用コンピュータ１０３は、放射線制御装置１０２を制御して放射線を曝射させるとともに、放射線撮影装置１０４を制御して画像情報を取得する撮影制御装置として機能する。制御用コンピュータ１０３は、例えば、不図示の１つまたは複数のプロセッサ（ＣＰＵ）、メモリを有し、１つまたは複数のプロセッサがメモリに格納されたプログラムを実行することで各種処理を実行する。放射線撮影装置１０４は、放射線を可視光に変換する蛍光体１０５と、可視光を検出する二次元検出器１０６を備えている。二次元検出器１０６は、放射線量子を検出する画素２０をＸ列×Ｙ行のアレイ状に配置したセンサであり、画像情報を出力する。

図２に、第１実施形態に係る画素２０の等価回路図を示す。画素２０は、光電変換素子２０１と、出力回路部２０２とを含む。光電変換素子２０１は、典型的にはフォトダイオードである。出力回路部２０２は、増幅回路部２０４、クランプ回路部２０６、サンプルホールド回路部２０７、選択回路部２０８を含む。

光電変換素子２０１は、電荷蓄積部を含む。この電荷蓄積部は、増幅回路部２０４のＭＯＳトランジスタ２０４ａのゲートに接続されている。ＭＯＳトランジスタ２０４ａのソースは、ＭＯＳトランジスタ２０４ｂを介して電流源２０４ｃに接続されている。ＭＯＳトランジスタ２０４ａと電流源２０４ｃとによってソースフォロア回路が構成されている。ＭＯＳトランジスタ２０４ｂは、そのゲートに供給されるイネーブル信号ＥＮがアクティブレベルになるとオンしてソースフォロア回路を動作状態にするイネーブルスイッチである。

図２に示す例では、光電変換素子２０１の電荷蓄積部およびＭＯＳトランジスタ２０４ａのゲートが共通のノードを構成していて、このノードは、電荷蓄積部に蓄積された電荷を電圧に変換する電荷電圧変換部として機能する。即ち、電荷電圧変換部には、電荷蓄積部に蓄積された電荷Ｑと電荷電圧変換部が有する容量値Ｃとによって定まる電圧Ｖ（＝Ｑ／Ｃ）が現れる。電荷電圧変換部は、リセットスイッチ２０３を介してリセット電位Ｖｒｅｓに接続されている。リセット信号ＰＲＥＳがアクティブレベルになると、リセットスイッチ２０３がオンして、電荷電圧変換部の電位がリセット電位Ｖｒｅｓにリセットされる。

クランプ回路部２０６は、リセットした電荷電圧変換部の電位に応じて増幅回路部２０４によって出力されるノイズをクランプ容量２０６ａによってクランプする。つまり、クランプ回路部２０６は、光電変換素子２０１で光電変換により発生した電荷に応じてソースフォロア回路から出力された信号から、ノイズをキャンセルするための回路である。このノイズはリセット時のｋＴＣノイズを含む。クランプは、クランプ信号ＰＣＬをアクティブレベルにしてＭＯＳトランジスタ２０６ｂをオン状態にした後に、クランプ信号ＰＣＬを非アクティブレベルにしてＭＯＳトランジスタ２０６ｂをオフ状態にすることによってなされる。クランプ容量２０６ａの出力側は、ＭＯＳトランジスタ２０６ｃのゲートに接続されている。ＭＯＳトランジスタ２０６ｃのソースは、ＭＯＳトランジスタ２０６ｄを介して電流源２０６ｅに接続されている。ＭＯＳトランジスタ２０６ｃと電流源２０６ｅとによってソースフォロア回路が構成されている。ＭＯＳトランジスタ２０６ｄは、そのゲートに供給されるイネーブル信号ＥＮ_０がアクティブレベルになるとオンしてソースフォロア回路を動作状態にするイネーブルスイッチである。

光電変換素子２０１で光電変換により発生した電荷に応じてクランプ回路部２０６から出力される信号は、光信号として、光信号サンプリング信号ＴＳがアクティブレベルになることによってスイッチ２０７Ｓａを介して容量２０７Ｓｂに書き込まれる。電荷電圧変換部の電位をリセットした直後にＭＯＳトランジスタ２０６ｂをオン状態とした際にクランプ回路部２０６から出力される信号は、クランプ電圧である。ノイズ信号は、ノイズサンプリング信号ＴＮがアクティブレベルになることによってスイッチ２０７Ｎａを介して容量２０７Ｎｂに書き込まれる。このノイズ信号には、クランプ回路部２０６のオフセット成分が含まれる。スイッチ２０７Ｓａと容量２０７Ｓｂによって信号サンプルホールド回路２０７Ｓが構成され、スイッチ２０７Ｎａと容量２０７Ｎｂによってノイズサンプルホールド回路２０７Ｎが構成される。以上のように、サンプルホールド回路部２０７は、信号サンプルホールド回路２０７Ｓとノイズサンプルホールド回路２０７Ｎとを含む。

信号線２１は、信号線２１Ｓと信号線２１Ｎを含む。駆動回路部が行選択信号をアクティブレベルに駆動すると、容量２０７Ｓｂに保持された信号（光信号）がＭＯＳトランジスタ２０８Ｓａおよび行選択スイッチ２０８Ｓｂを介して信号線２１Ｓに出力される。また、同時に、容量２０７Ｎｂに保持された信号（ノイズ）がＭＯＳトランジスタ２０８Ｎａおよび行選択スイッチ２０８Ｎｂを介して信号線２１Ｎに出力される。ＭＯＳトランジスタ２０８Ｓａは、信号線２１Ｓに設けられた不図示の定電流源とソースフォロア回路を構成する。同様に、ＭＯＳトランジスタ２０８Ｎａは、信号線２１Ｎに設けられた不図示の定電流源とソースフォロア回路を構成する。ＭＯＳトランジスタ２０８Ｓａと行選択スイッチ２０８Ｓｂによって信号用選択回路部２０８Ｓが構成され、ＭＯＳトランジスタ２０８Ｎａと行選択スイッチ２０８Ｎｂによってノイズ用選択回路部２０８Ｎが構成される。選択回路部２０８は、信号用選択回路部２０８Ｓとノイズ用選択回路部２０８Ｎとを含む。

画素２０は、隣接する複数の画素２０の光信号を加算する加算スイッチ２０９Ｓを有してもよい。加算モード時には、加算モード信号ＡＤＤがアクティブレベルになり、加算スイッチ２０９Ｓがオン状態になる。これにより、隣接する画素２０の容量２０７Ｓｂが加算スイッチ２０９Ｓによって相互に接続されて、光信号が平均化される。同様に、画素２０は、隣接する複数の画素２０のノイズを加算する加算スイッチ２０９Ｎを有してもよい。加算スイッチ２０９Ｎがオン状態になると、隣接する画素２０の容量２０７Ｎｂが加算スイッチ２０９Ｎによって相互に接続されて、ノイズが平均化される。加算部２０９は、加算スイッチ２０９Ｓと加算スイッチ２０９Ｎを含む。

画素２０は、感度を変更するための感度変更部２０５を有してもよい。画素２０は、例えば、第１感度変更スイッチ２０５ａおよび第２感度変更スイッチ２０５'ａ、並びにそれらに付随する回路素子を含みうる。第１変更信号ＷＩＤＥがアクティブレベルになると、第１感度変更スイッチ２０５ａがオンして、電荷電圧変換部の容量値に第１付加容量２０５ｂの容量値が追加される。これによって画素２０の感度が低下する。第２変更信号ＷＩＤＥ２がアクティブレベルになると、第２感度変更スイッチ２０５'ａがオンして、電荷電圧変換部の容量値に第２付加容量２０５'ｂの容量値が追加される。これによって画素２０の感度が更に低下する。このように画素２０の感度を低下させる機能を追加することによって、より大きな光量を受光することが可能となり、ダイナミックレンジを広げることができる。第１変更信号ＷＩＤＥがアクティブレベルになる場合には、イネーブル信号ＥＮ_ｎをアクティブレベルにして、ＭＯＳトランジスタ２０４ａに変えてＭＯＳトランジスタ２０４'ａをソースフォロア動作させてもよい。ＭＯＳトランジスタ２０４'ｂは、そのゲートに供給されるイネーブル信号ＥＮ_ｎがアクティブレベルになるとオンしてソースフォロア回路を動作状態にするイネーブルスイッチである。

放射線撮影装置１０４は、以上のような画素回路の出力を、不図示のＡＤ変換器でデジタル値に変換した後、制御用コンピュータ１０３に画像として転送する。撮影制御装置としての制御用コンピュータ１０３は、図１Ｂに示されるように、画像取得部１３１、画像処理部１３２、ＤＳＡ処理部１３３を有する。これらは、プロセッサがプログラムを実行することで実現されてもよいし、１部または全部が専用のハードウエアにより実現されてもよい。画像取得部１３１、画像処理部１３２、ＤＳＡ処理部１３３の詳細については、のちの説明により明らかになる。タイミングコントローラ１２０は、制御用コンピュータ１０３からの同期信号に従って、ＲＥＳＥＴ、ＳＨ＿Ｓ、ＳＨ＿Ｎを放射線撮影装置１０４に提供する。なお、ＲＥＳＥＴ信号は図２のリセット信号ＰＲＥＳに対応し、ＳＨ＿Ｓ信号とＳＨ＿Ｎ信号はそれぞれ図２の光サンプリング信号ＴＳとノイズサンプリング信号ＴＮに対応する。

次に第１実施形態による放射線撮影システム１００の駆動について説明する。制御用コンピュータ１０３の画像取得部１３１は、１ショット間にエネルギーが変化する放射線が被写体に照射され、当該１ショット間に前記被写体を透過した前記放射線を複数回検出して得られた、互いに異なる放射線エネルギーによる複数の画像を取得する。本実施形態では、連続する２つのリセットに挟まれた期間において二次元検出器１０６で少なくとも２回のサンプリングを行い、１フレームの撮影期間において互いに異なる放射線エネルギーによる第１の画像と第２の画像を取得する。図３は、第１実施形態に係る放射線撮影システムにおいてエネルギーサブトラクションを行う場合の二次元検出器１０６の駆動タイミング（サンプリングと読出しのタイミング）を示すタイミングチャートである。図中の波形は横軸を時間として、放射線の曝射、同期信号、光電変換素子２０１のリセット、サンプルホールド回路２０７、信号線２１からの画像の読み出しのタイミングを示している。図３に示される各信号（Ｒｅｓｅｔ、ＳＨ＿Ｓ、ＳＨ＿Ｎ）は例えば、図１Ｂに示されるようなタイミングコントローラ１２０から供給され得る。なお、タイミングコントローラ１２０は、制御用コンピュータ１０３に組み込まれていてもいし、放射線撮影装置１０４に組み込まれていてもよいし、制御用コンピュータ１０３と放射線撮影装置１０４の間に介在する独立した装置であってもよい。二次元検出器１０６からは、信号サンプルホールド回路２０７Ｓとノイズサンプルホールド回路２０７Ｎでホールドされた信号の差分（Ｓ－Ｎ）が読み出される。

まず、光電変換素子２０１のリセットを行ってから、放射線を曝射する。放射線の管電圧は理想的には矩形波であるが、管電圧の立ち上がりと立下りには有限の時間がかかる。特に、パルス放射線で曝射時間が短い場合は、管電圧はもはや矩形波とはみなせず、図３に示すような波形となる。すなわち、放射線の立ち上がり期、安定期、立下り期で放射線のエネルギーが異なる。

そこで、立ち上がり期の放射線３０１が曝射された後に、ノイズサンプルホールド回路２０７Ｎでサンプリングを行い、さらに安定期の放射線３０２が曝射された後に信号サンプルホールド回路２０７Ｓでサンプリングを行う。その後、信号線２１Ｎと信号線２１Ｓの差分を画像として読み出す。このとき、ノイズサンプルホールド回路２０７Ｎには立ち上がり期の放射線３０１の信号（Ｒ_１）が保持されている。また、信号サンプルホールド回路２０７Ｓには立ち上がり期の放射線３０１の信号（Ｒ_１）と安定期の放射線３０２の信号（Ｂ）の和となる信号（Ｒ_１＋Ｂ）が保持されている。従って、信号（Ｒ_１＋Ｂ）から信号（Ｒ_１）が差し引かれた画像、すなわち安定期の放射線３０２の信号に対応した画像３０４（Ｂ）が読み出される。

次に、立下り期の放射線３０３の曝射と、画像３０４の読み出しとが完了してから、再び信号サンプルホールド回路２０７Ｓでサンプリングを行う。このサンプリングにより画像（Ｒ_１＋Ｂ＋Ｒ_２）が取得される。その後、光電変換素子２０１のリセットを行い、再びノイズサンプルホールド回路２０７Ｎでサンプリングを行い、信号線２１Ｎと信号線２１Ｓの差分を画像として読み出す。このとき、ノイズサンプルホールド回路２０７Ｎには放射線が曝射されていない状態の信号（０）が保持される。また、信号サンプルホールド回路２０７Ｓには立ち上がり期の放射線３０１の信号と安定期の放射線３０２と立下り期の放射線３０３の信号の和（Ｒ_１＋Ｂ＋Ｒ_２）が保持されている。従って、立ち上がり期の放射線３０１の信号と安定期の放射線３０２の信号と立下り期の放射線３０３の信号に対応した画像３０６（Ｒ_１＋Ｂ＋Ｒ_２）が読み出される。その後、画像３０６と画像３０４の差分を計算することで、立ち上がり期の放射線３０１と立下り期の放射線３０３の和に対応した画像３０５（Ｒ_１＋Ｒ_２）が得られる。

タイミングコントローラ１２０は、サンプルホールド回路２０７におけるサンプリング及び光電変換素子２０１のリセットを行うタイミングを、放射線発生装置１０１から放射線の曝射が開始されたことを示す同期信号３０７を用いて決定する。放射線の曝射開始を検出する方法としては、例えば、放射線発生装置１０１の管電流を測定し、測定された管電流（電流値）が予め設定された閾値を上回るか否かを判定することにより曝射開始を検出する方法を用いることができる。或いは、光電変換素子２０１のリセットが完了した後、画素２０を繰り返して読み出し、画素値が予め設定された閾値を上回るか否かを判定することにより曝射開始を検出する方法を用いることもできる。さらには、放射線撮影装置１０４に二次元検出器１０６とは異なる放射線検出器を内蔵し、その測定値が予め設定された閾値を上回るか否かを判定することにより曝射開始を検出する方法も用いられ得る。いずれの方式の場合も、同期信号３０７の入力から予め指定した時間が経過した後に、信号サンプルホールド回路２０７Ｓのサンプリング、ノイズサンプルホールド回路２０７Ｎのサンプリング、光電変換素子２０１のリセットが行われる。

以上のようにして、パルス放射線の安定期に対応した画像３０４と、立ち上がり期と立下り期の和に対応した画像３０５とが制御用コンピュータ１０３により得られる。これら２枚の画像を形成する際に曝射された放射線のエネルギーが異なるため、これらの画像を用いてエネルギーサブトラクション処理を行うことができる。

図４に、第１実施形態に係る放射線撮影システム１００においてエネルギーサブトラクションを行う場合の駆動タイミングの他の例を示す。図３に示される駆動タイミングとは、放射線の管電圧を能動的に切り替えている点で異なる。

まず、光電変換素子２０１のリセットを行ってから、低エネルギーの放射線４０１を曝射する。その後、ノイズサンプルホールド回路２０７Ｎでサンプリングを行ってから、管電圧を切り替えて高エネルギーの放射線４０２が曝射された後に、信号サンプルホールド回路２０７Ｓでサンプリングを行う。その後、管電圧を切り替えて低エネルギーの放射線４０３の曝射を行う。さらに、信号線２１Ｎと信号線２１Ｓの差分を画像として読み出す。このとき、ノイズサンプルホールド回路２０７Ｎには低エネルギーの放射線４０１の信号（Ｒ_１）が保持され、信号サンプルホールド回路２０７Ｓには低エネルギーの放射線４０１の信号と高エネルギーの放射線４０２の信号の和（Ｒ_１＋Ｂ）が保持されている。従って、高エネルギーの放射線４０２の信号に対応した画像４０４（Ｂ）が読み出される。

次に、低エネルギーの放射線４０３の曝射と、画像４０４の読み出しとが完了してから、再び信号サンプルホールド回路２０７Ｓでサンプリングを行う。その後、光電変換素子２０１のリセットを行い、再びノイズサンプルホールド回路２０７Ｎでサンプリングを行い、信号線２１Ｎと信号線２１Ｓの差分を画像として読み出す。このとき、ノイズサンプルホールド回路２０７Ｎには放射線が曝射されていない状態の信号（０）が保持されている。また、信号サンプルホールド回路２０７Ｓには低エネルギーの放射線４０１の信号と高エネルギーの放射線４０２と低エネルギーの放射線４０３の信号の和（Ｒ_１＋Ｂ＋Ｒ_２）が保持されている。従って、低エネルギーの放射線４０１の信号と高エネルギーの放射線４０２の信号と低エネルギーの放射線４０３の信号に対応した画像４０６（Ｒ_１＋Ｂ＋Ｒ_２）が読み出される。その後、画像４０６と画像４０４の差分を計算することで、低エネルギーの放射線４０１と低エネルギーの放射線４０３の和に対応した画像４０５（Ｒ_１＋Ｒ_２）が得られる。同期信号４０７については、図３と同様である。このように、管電圧を能動的に切り替えながら画像を取得することで、図３の方法に比べて低エネルギーと高エネルギーの画像の間のエネルギー差をより大きくすることが出来る。

次に、エネルギーサブトラクション処理について説明する。以下に説明するように、画像取得部１３１と画像処理部１３２は、複数の画像を用いてエネルギーサブトラクション処理を行い、エネルギーサブトラクション画像を生成する。本実施形態では、画像取得部１３１が、連続する２つのリセットに挟まれた期間に、１つの二次元検出器１０６において少なくとも２回のサンプリングが行われる。これにより、画像取得部１３１は、互いに異なる放射線エネルギーの２つの画像を１フレームの撮影期間において取得する。画像取得部１３１は、さらに、取得した２つの画像について図５に示すような補正処理を施し、減弱率の画像である第１の画像と第２の画像を得る。制御用コンピュータ１０３の画像処理部１３２は、画像取得部１３１により得られた第１の画像と第２の画像についてエネルギーサブトラクション処理（図６（ａ）に示される信号処理）を行ってエネルギーサブトラクション画像（以下、処理画像）を得る。

図５は、第１実施形態によるエネルギーサブトラクション処理において、画像取得部１３１が行う補正処理を説明する図である。まず、放射線撮影装置１０４に放射線を曝射せずに撮影を行い、図３または図４に示した駆動で画像を取得する。このとき２枚の画像が読み出されるが、１枚目の画像をF_ODD、２枚目の画像をF_EVENとする。F_ODDとF_EVENは、放射線撮影装置１０４の固定パターンノイズ（ＦＰＮ）に対応する画像である。次に、被写体がない状態で放射線撮影装置１０４に放射線を曝射して撮影を行い、図３または図４に示した駆動で画像を取得する。このとき２枚の画像が読み出されるが、１枚目の画像をW_ODD、２枚目の画像をW_EVENとする。W_ODDとW_EVENは、放射線撮影装置１０４のＦＰＮと放射線による信号の和に対応する画像である。

従って、W_ODDからF_ODDを、W_EVENからF_EVENを減算することで、放射線撮影装置１０４のＦＰＮが除去された画像WF_ODDとWF_EVENが得られる。WF_ODDは、例えば図３における安定期の放射線３０２に対応する画像であり、WF_EVENは立ち上がり期の放射線３０１と安定期の放射線３０２と立下り期の放射線３０３との和に対応する画像である。従って、WF_EVENからWF_ODDを減算することで、立ち上がり期と立下り期の放射線の和に対応する画像が得られる。放射線の立ち上がり期と立下り期のエネルギーは、放射線の安定期のエネルギーに比べて低い。従って、WF_EVENからWF_ODDを減算することで、被写体がない場合の低エネルギー画像W_Lowが得られる。また、WF_ODDから、被写体がない場合の高エネルギー画像W_Highが得られる。

次に、被写体がある状態で放射線撮影装置１０４に放射線を曝射して撮影を行い、図３、または図４に示した駆動で画像を取得する。このとき１回の曝射に対して２枚の画像が読み出されるが、１枚目の画像をX_ODD、２枚目の画像をX_EVENとする。被写体がない場合と同様の減算を行うことで、被写体がある場合の低エネルギー画像X_Lowと、被写体がある場合の高エネルギー画像X_Highが得られる。

ここで、被写体の厚みをｄ、被写体の線減弱係数をμ、被写体がない場合の画素２０の出力をＩ_０、被写体がある場合の画素２０の出力をＩとすると、以下の式が成り立つ。

［数１］を変形すると、以下の式が得られる。

［数２］の右辺は被写体の減弱率を示す。すなわち、Ｉ／Ｉ_０は、減弱率の画像である。なお、被写体の減弱率は０～１の間の実数である。

従って、被写体がある場合の低エネルギー画像X_Lowを、被写体がない場合の低エネルギー画像W_Lowで除算することで、低エネルギーにおける減弱率の画像Ｌが得られる。同様に、被写体がある場合の高エネルギー画像X_Highを、被写体がない場合の高エネルギー画像W_Highで除算することで、高エネルギーにおける減弱率の画像Ｈが得られる。

図６（ａ）に、画像処理部１３２におけるエネルギーサブトラクション処理の信号処理のブロック図を示す。本実施形態のエネルギーサブトラクション処理では、２種類の放射線エネルギーに対する減弱率を示す第１の画像と第２の画像から骨と軟部組織の画像を処理画像として取得する。すなわち、本実施形態における信号処理では、図５に示した補正によって得られた低エネルギーにおける減弱率の画像Ｌと高エネルギーにおける減弱率の画像Ｈから、骨の厚さの画像Ｂと軟部組織の厚さの画像Ｓを求める。なお、以下では、画像Ｌを低エネルギーの画像、画像Ｈを高エネルギーの画像ともいう。

放射線フォトンのエネルギーをＥ、エネルギーＥにおけるフォトン数をＮ（Ｅ）、骨の厚さをＢ、軟部組織の厚さをＳ、エネルギーＥおける骨の線減弱係数をμ_Ｂ（Ｅ）、エネルギーＥおける軟部組織の線減弱係数をμ_Ｓ（Ｅ）、減弱率をＩ／Ｉ_０とする。すると、以下の式が成り立つ。

エネルギーＥにおけるフォトン数Ｎ（Ｅ）は、放射線のスペクトルである。放射線のスペクトルは、シミュレーション又は実測により得られる。また、エネルギーＥおける骨の線減弱係数μ_Ｂ（Ｅ）とエネルギーＥおける軟部組織の線減弱係数μ_Ｓ（Ｅ）は、ＮＩＳＴなどのデータベースから得られる。すなわち、任意の骨の厚さＢ、軟部組織の厚さＳ、放射線のスペクトルＮ（Ｅ）における減弱率Ｉ／Ｉ_０を計算することが可能である。

ここで、低エネルギーの放射線におけるスペクトルをＮ_Ｌ（Ｅ）、高エネルギーの放射線におけるスペクトルをＮ_Ｈ（Ｅ）とすると、以下の式が成り立つ。

また、高エネルギーの放射線における平均エネルギーＥ_Ｈと低エネルギーの放射線における平均エネルギーＥ_Ｌを、以下の式で表すものとする。

ここで、放射線のスペクトルＮ（Ｅ）を、平均エネルギーにおいてピークを持つ単色放射線として近似すると、［数４］は以下のように近似できる。

［数６］の両辺の対数をとると、［数７］となる。

［数７］の係数を行列としたときの行列式をｄｅｔとすると、骨の厚さの画像Ｂ、軟部組織の厚さの画像Ｓは以下の式で表される。

エネルギーＥにおける線減弱係数μ（Ｅ）は既知の値である。従って、全ての画素について［数８］を計算することで、低エネルギーにおける減弱率の画像Ｌと高エネルギーにおける減弱率の画像Ｈから、骨の厚さの画像Ｂ、軟部組織の厚さの画像Ｓを得ることができる。具体的には、低エネルギーにおける減弱率の画像Ｌの対数の画像ｌｎ（Ｌ）と、高エネルギーにおける減弱率の画像Ｈの対数の画像ｌｎ（Ｈ）に対して、線減弱係数μ（Ｅ）によって決まる係数を掛けてから差分する。これにより、骨の厚さの画像Ｂと軟部組織の厚さの画像Ｓが求められる。

なお、本実施形態では［数８］における係数を平均エネルギーＥ_ＨとＥ_Ｌにおける線減弱係数から求めたが、本実施形態はこれに限定されない。例えば、平均エネルギー以外のエネルギーにおける線減弱係数から［数８］の係数を求めてもよい。また、［数８］において使用する線弱係数をユーザーが指定してもよいし、画像の統計量から推定するなどしてもよい。また、第１実施形態では説明を簡略化するために、骨の厚さＢと軟部組織の厚さＳを算出していたが、このような形態に限定されるものではない。例えば、水の厚さＷと造影剤の厚さＩを算出してもよい。要するに第１実施形態における信号処理とは、任意の二種類の物質の厚さに分解する処理であると言える。

以上のように、第１の画像（Ｈ）と第２の画像（Ｌ）についてエネルギーサブトラクションを行って、第１の物質（骨）の厚さに関わる第１の処理画像（Ｂ）と第１の物質とは異なる第２の物質（軟部組織）の厚さに関わる第２の処理画像（Ｓ）が取得される。より具体的には、第１の放射線エネルギーにおける第１の物質と第２の物質のそれぞれの線減弱係数と、第２の放射線エネルギーにおける第１の物質と第２の物質のそれぞれの線減弱係数とに基づいて第１の処理画像と第２の処理画像が取得される。

ＤＳＡ処理部１３３は、以上のようにして画像処理部１３２で生成された複数のエネルギーサブトラクション画像を用いて差分画像を生成する。本実施形態では、差分画像の生成の一例としてデジタル差分血管造影法（以下、ＤＳＡ（digital subtraction angiography））が用いられる。ＤＳＡ処理部１３３は、エネルギーサブトラクション処理により得られた処理画像を用いてＤＳＡによる差分画像を生成する。ＤＳＡ画像の生成では、造影剤の注入前後の第１の処理画像（Ｂ）が用いられる。すなわち、第１実施形態では、造影剤の注入前にエネルギーサブトラクション処理により得られる処理画像を用いてマスク画像を生成する。そして、造影剤の注入後にエネルギーサブトラクション処理により得られる処理画像を用いてライブ画像を生成し、マスク画像とライブ画像の差分画像（ＤＳＡ画像）を生成する。

図７に、第１実施形態に係るフローチャートを示す。まず、画像取得部１３１は、造影剤を注入する前に図３又は図４に示したタイミングで撮影を行い、図５に示した補正を行うことで、高エネルギーのマスク画像Ｈ_Ｍと、低エネルギーのマスク画像Ｌ_Ｍを取得する（Ｓ７０１）。造影剤が注入された後（Ｓ７０２）、画像取得部１３１は、マスク画像の取得時と同様に撮影と補正を行い、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌを取得する（Ｓ７０３）。その後、画像処理部１３２は、高エネルギーのマスク画像Ｈ_Ｍと低エネルギーのマスク画像Ｌ_Ｍに対して図６（ａ）に示した信号処理を行い、骨のマスク画像Ｂ_Ｍと軟部組織のマスク画像Ｓ_Ｍに分離する。さらに、画像処理部１３２は、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌに対しても同様の処理を行い、骨のライブ画像Ｂ_Ｌと軟部組織のライブ画像Ｓ_Ｌに分離する（Ｓ７０４）。そして、ＤＳＡ処理部１３３は、骨のライブ画像Ｂ_Ｌから骨のマスク画像Ｂ_Ｍを引くことで骨のＤＳＡ画像Ｂ_ＤＳＡを、軟部組織のライブ画像Ｓ_Ｌから軟部組織のマスク画像Ｓ_Ｍを引くことで軟部組織のＤＳＡ画像Ｓ_ＤＳＡを得る（Ｓ７０５）。処理の終了が指示されるまでの間（Ｓ７０６でＮＯ）、上述したライブ画像の撮影と画像処理を繰り返し行うことで、エネルギーサブトラクションとＤＳＡを組み合わせた動画像を得ることができる。終了が指示されると（Ｓ７０６でＹＥＳ）、本処理は終了する。

図８に、第１実施形態に係る画像処理の模式図を示す。図８では、骨のマスク画像Ｂ_Ｍと骨のライブ画像Ｂ_ＬとからＤＳＡ画像Ｂ_ＤＳＡが生成される様子が示されている。高エネルギーのマスク画像Ｈ_Ｍと、低エネルギーのマスク画像Ｌ_Ｍは、造影剤を注入する前に撮影した画像であるため、骨と軟部組織の情報のみを含んでいる。これらのマスク画像に対してエネルギーサブトラクションを行って得た骨のマスク画像Ｂ_Ｍには骨の情報が、軟部組織のマスク画像Ｓ_Ｍ（不図示）には軟部組織の情報が含まれる。一方、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌは、造影剤を注入した後に撮影した画像であるため、骨と軟部組織と造影剤の情報を含んでいる。これらのライブ画像に対してエネルギーサブトラクションを行うと、三種類目の物質である造影剤の情報が骨のライブ画像Ｂ_Ｌに現れる。骨の線弱係数と造影剤の線弱係数が近いためである。すなわち、骨のライブ画像Ｂ_Ｌには骨と造影剤の情報が、軟部組織のライブ画像Ｓ_Ｌ（不図示）には軟部組織の情報が含まれる。

マスク画像とライブ画像の間で、心臓の拍動や呼吸などによって臓器が動くと、軟部組織の画像には大きな変化が現れる。しかし骨の位置は大きく動かない。すなわち、骨のライブ画像Ｂ_Ｌと骨のマスク画像Ｂ_Ｍに含まれる骨の情報はほぼ同じである。従って、骨のＤＳＡ画像Ｂ_ＤＳＡでは骨の情報が相殺され、ほぼ造影剤の情報のみが残る。このようにして、心臓の拍動や呼吸などによる動きがあっても、造影剤のみを画像化して、血管の診断を行うことができる。

また、第１実施形態では図３や図４で示したように、サンプルホールド回路を用いることで高エネルギーと低エネルギーの放射線画像を取得していたが、このような形態に限定されるものではない。例えば、高エネルギーの放射線を曝射して画像の読み出しを行ってから、低エネルギーの放射線を曝射して画像の読み出しを行い、高エネルギーと低エネルギーの放射線画像を取得してもよい。しかしながら、本実施形態では心臓の拍動や呼吸などにより被写体が動く状況下で用いることを想定している。被写体が動いた場合、高エネルギーと低エネルギーの放射線画像の間で時間差があると、エネルギーサブトラクション後の画像にモーションアーチファクトを生じるという課題がある。このような理由から、本実施形態では図３や図４で示したように、１つのフレームを撮影する間に（次のリセットがかかる前に）高エネルギーと低エネルギーの放射線画像を撮影する。このような撮影は、二次元検出器１０６の画素２０が有するサンプルホールド回路部２０７を用いることで実現され得る。このような撮影を時分割駆動と称する。時分割駆動は、高エネルギーと低エネルギーの放射線画像の撮影の時間差を短縮することができるため、ＤＳＡ画像の生成にとってより望ましい形態である。すなわち、異なるエネルギーの画像を得るための撮影の間隔を短くすることができ、高精度でアーチファクトのないＤＳＡ画像を得ることができる。また、サンプルホールド回路を用いることで高エネルギーと中エネルギーと低エネルギーの放射線画像（３以上の放射線画像）を取得し、エネルギーサブトラクション処理を行ってもよい。

本実施形態においては、互いにエネルギーの異なる放射線画像を取得するために、ＦＰＤを２枚重ねることで、放射線照射側から見て前面のＦＰＤで低エネルギーの放射線画像を、背面のＦＰＤで高エネルギーの放射線画像を取得する方法を用いてもよい。この方法では、高エネルギーと低エネルギーの放射線画像の間で時間差が生じないため、モーションアーチファクトが生じない。しかしながらＦＰＤを積層した場合、前面のＦＰＤを放射線が透過することになる。このＦＰＤがアモルファスシリコンであれば問題ないが、結晶シリコンの場合、放射線によって画素回路の劣化が進んで、画質が悪化するという問題がある。血管造影や血管内治療では、長時間にわたって放射線を照射し続けるため、画素回路の劣化が進行しやすい。また、結晶シリコンのＦＰＤはアモルファスシリコンと比較して高価である。すなわち、結晶シリコンのＦＰＤを積層することは、画質やコストの点から望ましくない。一方、アモルファスシリコンのＦＰＤは、結晶シリコンと比較してオフセットや残像が生じやすい。特に動画の場合は、より大きなオフセットや残像が発生し、エネルギーサブトラクション後の画像にアーチファクトを生じるという課題がある。このような理由からも、本実施形態の時分割駆動による高エネルギーと低エネルギーの放射線画像の撮影は、より望ましい形態である。

また、互いに異なる放射線エネルギーの第１の画像と第２の画像を取得するための駆動タイミングは図３、図４に示した駆動タイミングに限られない。互いに異なる放射線エネルギーの第１の画像と第２の画像を１フレームの撮影期間におい取得できる駆動タイミングであればよい。

＜第２実施形態＞
第１実施形態では、高エネルギーにおける減弱率の放射線画像（Ｈ）と低エネルギーにおける減弱率の放射線画像（Ｌ）とから、二種類の物質の厚さの画像を取得した（例えば、図６（ａ）の軟部組織（Ｓ）と骨（Ｂ）の画像）。第２実施形態では、高エネルギーの放射線画像（Ｈ）と低エネルギーの放射線画像（Ｌ）とから、実効原子番号（Ｚ）の画像と面密度（Ｄ）の画像を取得する（図６（ｂ））。なお、第２実施形態の放射線撮影システム１００の構成及び高／低エネルギーの放射線画像を取得する動作は、第１実施形態（図１Ａ、１Ｂ、２～５）と同様である。

図６（ｂ）は、第２実施形態の画像処理部１３２によるエネルギーサブトラクション処理の信号処理のブロック図である。上述したように、第２実施形態における信号処理では、図５に示した補正処理によって得られた低エネルギーにおける減弱率の画像Ｌと高エネルギーにおける減弱率の画像Ｈから、実効原子番号Ｚの画像と面密度Ｄの画像を求める。なお、実効原子番号Ｚとは混合物の等価的な原子番号のことである。また、面密度Ｄとは、被写体の密度［ｇ／ｃｍ^３］と被写体の厚み［ｃｍ］の積であり、次元は［ｇ／ｃｍ^２］となる。

まず、放射線フォトンのエネルギーをＥ、エネルギーＥにおけるフォトン数をＮ（Ｅ）、実効原子番号をＺ、面密度をＤ、実効原子番号ＺかつエネルギーＥおける質量減弱係数をμ（Ｚ，Ｅ）、減弱率をＩ／Ｉ_０とすると、以下の式が成り立つ。

エネルギーＥにおけるフォトン数Ｎ（Ｅ）は、放射線のスペクトルである。放射線のスペクトルは、シミュレーション又は実測により得られる。また、実効原子番号ＺかつエネルギーＥおける質量減弱係数μ（Ｚ，Ｅ）は、ＮＩＳＴなどのデータベースから得られる。すなわち、［数９］により、任意の実効原子番号Ｚ、面密度Ｄ、放射線のスペクトルＮ（Ｅ）における減弱率Ｉ／Ｉ_０を計算することが可能である。

ここで、低エネルギーの放射線のスペクトルをＮ_Ｌ（Ｅ）、高エネルギーの放射線のスペクトルをＮ_Ｈ（Ｅ）とすると、以下の式が成り立つ。

ここで第１実施形態と同様に、放射線のスペクトルＮ（Ｅ）を、平均エネルギーにおいてピークを持つ単色放射線として近似すると、［数１０］は以下のように近似できる。

［数１１］の両辺の対数をとると、［数１２］となる。

さらに、［数１２］の比をとることで面密度の項を相殺すると、以下の式が得られる。

すなわち、低エネルギーにおける減弱率Ｌと高エネルギーにおける減弱率Ｈの対数の比は、実効原子番号から一意に求められる。従って逆関数を求めれば、低エネルギーにおける減弱率のＬと高エネルギーにおける減弱率のＨの対数の比から実効原子番号Ｚを求めることができる。また、実行原子番号Ｚが求まると、［数１２］から面密度Ｄを求めることができる。

以上のようにして、全ての画素について［数１４］を計算することで、低エネルギーにおける減弱率の画像Ｌと高エネルギーにおける減弱率の画像Ｈから、実効原子番号Ｚの画像と面密度Ｄの画像を得ることができる。

図６（ｃ）に第２実施形態に係るエネルギーサブトラクション処理の画像処理のブロック図を示す。第２実施形態における画像処理では、図６（ｂ）に示した信号処理によって得られた実効原子番号Ｚの画像と面密度Ｄの画像から、カラーの透視画像を得る。このときの変換には、ＨＬＳ色空間が好適に用いられる。ＨＬＳ色空間とは、色相Ｈ（Ｈｕｅ）、彩度Ｓ（Ｓａｔｕｒａｔｉｏｎ）、輝度Ｌ（Ｌｉｇｈｔｎｅｓｓ）の３つ成分からなる色空間である。ＨＬＳ色空間を利用する方法としては、実効原子番号Ｚにより色相Ｈを決定し、面密度Ｄによって輝度Ｌを決定する構成が用いられ得る。このような画像処理によって生成されたカラーの放射線画像は、空港における手荷物検査などでも利用されている。なお、彩度Ｓについては、ユーザーが指定できるようにしてもよいし、固定値が用いられてもよい。

ＨＬＳ色空間の色相Ｈ、彩度Ｓ、輝度ＬからＲＧＢ値に変換してカラー画像を生成する方法として、円柱モデルや円錐モデルを用いる方法がある。例えば円錐モデルにおいては、以下の式を用いる。

彩度Ｓは０～１の間の実数である。例えば彩度Ｓ＝０の場合はＭａｘ＝Ｌ，Ｍｉｎ＝Ｌとなることから、色相Ｈがどのような値であっても、ＲＧＢ値は（Ｌ，Ｌ，Ｌ）となる。すなわち、輝度Ｌのみが反映されたモノクロの画像となる。これに対して彩度Ｓが１に近づくほど色相Ｈの情報が強く反映されるため、被写体を構成する材料を色によって識別しやすくなる。このため、彩度Ｓを手動で調整する形態が好ましい。

第２実施形態に係るフローチャートは、図７と同様である。ただし、エネルギーサブトラクション（Ｓ７０４）とＤＳＡの処理方法（Ｓ７０５）が異なる。まず、画像取得部１３１は、造影剤を注入する前に図３又は図４に示したタイミングで撮影を行い、図５に示した補正を行うことで、高エネルギーのマスク画像Ｈ_Ｍと、低エネルギーのマスク画像Ｌ_Ｍを取得する。次に、画像取得部１３１は、造影剤が注入された後、上記と同様に撮影と補正を行い、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌを取得する。その後、画像処理部１３２は、高エネルギーのマスク画像Ｈ_Ｍと低エネルギーのマスク画像Ｌ_Ｍに対して図６（ｂ）に示した信号処理を行い、実効原子番号のマスク画像Ｚ_Ｍと面密度のマスク画像Ｄ_Ｍに分離する。さらに、画像処理部１３２は、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌに対しても同様の処理を行い、実効原子番号のライブ画像Ｚ_Ｌと面密度のマスク画像Ｄ_Ｌに分離する。このとき、各画像の関係は以下の［数１７］で表される。

ここで、ｎは２．５以上、３以下の実数であり、例えば、ｎ＝２．９４である。これを変形すると、以下の［数１８］が得られる。

ＤＳＡ処理部１３３は、［数１８］を用いて、面密度のＤＳＡ画像Ｄ_ＤＳＡと実効原子番号のＤＳＡ画像Ｚ_ＤＳＡを求める。すなわち、実行原子番号のＤＳＡ画像は、面密度のマスク画像およびライブ画像（Ｄ_ＭおよびＤ_Ｌ）と実行原子番号のマスク画像およびライブ画像（Ｚ_ＭおよびＺ_Ｌ）に基づいて取得される。なお、実効原子番号のライブ画像Ｚ_Ｌから実効原子番号のマスク画像Ｚ_Ｍを単純に差分するなどして実効原子番号のＤＳＡ画像Ｚ_ＤＳＡを求めてもよい。最後に、面密度のＤＳＡ画像Ｄ_ＤＳＡと実効原子番号のＤＳＡ画像Ｚ_ＤＳＡに対して図１０の画像処理を行うことで、カラーのＤＳＡ画像Ｃ_ＤＳＡを得る。このようなライブ画像の撮影と画像処理を繰り返し行うことで、エネルギーサブトラクションとＤＳＡを組み合わせたカラーの動画像を得ることができる。

図９に、第２実施形態に係る画像処理の模式図を示す。図８と同様、マスク画像とライブ画像の間で、心臓の拍動や呼吸などによって臓器は動くが、骨の位置は大きく動かない。従って、臓器が動いた部分に関して、面密度のＤＳＡ画像Ｄ_ＤＳＡは変化するが、実効原子番号のＤＳＡ画像Ｚ_ＤＳＡはあまり変化しない。一方、造影剤が注入された部分は、面密度のＤＳＡ画像Ｄ_ＤＳＡと実効原子番号のＤＳＡ画像Ｚ_ＤＳＡが変化する。従って、実効原子番号のＤＳＡ画像Ｚ_ＤＳＡには造影剤の情報のみが残る。これをそのまま表示してもよい。

或いは、実効原子番号のＤＳＡ画像Ｚ_ＤＳＡと面密度のＤＳＡ画像Ｄ_ＤＳＡとを用いて、新たなＤＳＡ画像を生成してもよい。この場合、例えば、実効原子番号のＤＳＡ画像Ｚ_ＤＳＡｓと面密度のＤＳＡ画像Ｄ_ＤＳＡの画素値に対してそれぞれ閾値を設定し、閾値を下回る場合に画素値を０にクリップした後で両者を乗算するなどの処理が好適に用いられる。さらに、図６（ｃ）の画像処理を行って得たカラーのＤＳＡ画像Ｃ_ＤＳＡを生成してもよい。いずれの場合も、面密度の変化と実効原子番号の変化を両方用いるため、より正確に造影剤のみを画像化することができる。このようにして、心臓の拍動や呼吸などによる動きがあっても、造影剤のみを画像化して、血管の診断を行うことができる。

また、第２実施形態における実効原子番号のＤＳＡ画像Ｚ_ＤＳＡでは、マスク画像を引くことで軟部組織や骨によりもたらされる実効原子番号の変化が除去されている。従って、マスク画像の取得後に挿入された物質の実効原子番号（種類）を測定することができる。すなわち、実効原子番号を用いて、造影剤（ヨウ素）とカテーテルやステント（鉄）を識別することができるという利点がある。また、呼吸などによって骨の位置が多少ずれてしまった部分は、マスク画像の取得後に骨が挿入されたものとみなせる。従って、実効原子番号を用いて造影剤（ヨウ素）と骨（カルシウム）を識別し、造影剤のみを画像化するなどの処理も可能となる。

なお、実効原子番号と面密度の画像を用いたＤＳＡ画像の生成において、２つの異なる放射線エネルギーによる減弱率の画像を得るために時分割駆動を用いることは好ましい実施形態であるが、これに限られるものではない。例えば、第１実施形態で上述したように、時分割駆動以外の方法で２つの異なる放射線エネルギーによる減弱率の画像が得られてもよい。

＜第３実施形態＞
第１、第２実施形態では、造影剤の注入前に撮影を行ってマスク画像を生成した。第３実施形態では、エネルギーサブトラクション処理により得られる処理画像のうち、現在のフレームよりも過去の複数のフレームの処理画像に基づいてマスク画像を生成する。そして、現在のフレームの処理画像をライブ画像として用いてマスク画像との差分画像であるＤＳＡ画像を生成する。より具体的には、画像処理部１３２が、異なる放射線エネルギーによる複数の画像を用いてエネルギーサブトラクション処理を行い、エネルギーサブトラクション画像を生成する。ＤＳＡ処理部１３３は、過去に生成された複数のエネルギーサブトラクション画像に基づくマスク画像と現在のエネルギーサブトラクション画像に基づくライブ画像との差分画像を生成する。なお、第３実施形態の放射線撮影システム１００の構成及び動作は、第１実施形態～第２実施形態と同様である。以下、第３実施形態について詳細に説明する。

図１０に、第３実施形態に係るフローチャートを示す。第１実施形態とは異なり、マスク画像の撮影は行わず、造影剤の注入の有無は問わない。まず、画像取得部１３１は、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌを取得する（Ｓ１００１）。その後、画像処理部１３２は、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌに対して図６に示した信号処理を行い、骨のライブ画像Ｂ_Ｌと軟部組織のライブ画像Ｓ_Ｌに分離する（Ｓ１００２）。このようなライブ画像の撮影とエネルギーサブトラクション処理を繰り返して行うことで動画像が得られる。

ここで、現在のフレームの骨のライブ画像をＢ_Ｌ［０］、ｎフレーム前の骨のライブ画像をＢ_Ｌ［－ｎ］とする。なお、ｎは２以上の整数である。ＤＳＡ処理部１３３は、ｎフレーム前の骨のライブ画像Ｂ_Ｌ［－ｎ］から１フレーム前の骨のライブ画像Ｂ_Ｌ［－１］までのｎ枚の画像を用いて、骨のマスク画像Ｂ_Ｍを生成する（Ｓ１００３）。最後に、ＤＳＡ処理部１３３は、現在のフレームの骨のライブ画像Ｂ_Ｌ［０］を骨のライブ画像Ｂ_Ｌとして、骨のライブ画像Ｂ_Ｌから骨のマスク画像Ｂ_Ｍを引くことで骨のＤＳＡ画像Ｂ_ＤＳＡを得る（Ｓ１００４）。軟部組織についても同様にして、軟部組織のライブ画像Ｓ_Ｌから軟部組織のマスク画像Ｓ_Ｍを引くことで軟部組織のＤＳＡ画像Ｓ_ＤＳＡを得ることができる。このようなライブ画像の撮影と画像処理を繰り返し行うことで、エネルギーサブトラクションとＤＳＡを組み合わせた動画像を得る。

図１１に、第３実施形態に係る画像処理の模式図を示す。説明を単純にするため、心臓の冠動脈に造影剤を注入しながらライブ画像を撮影した場合を想定する。高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌは、骨と軟部組織と造影剤の情報を含んでいる。これらのライブ画像に対してエネルギーサブトラクションを行うと、三種類目の物質である造影剤の情報が骨のライブ画像Ｂ_Ｌに現れる。すなわち、骨のライブ画像Ｂ_Ｌには骨と造影剤の情報が含まれ得る。

このとき、心臓の拍動に伴って血管が動くため、造影剤も血管と一緒に移動する。一方で、骨は心臓の拍動によっては移動しない。すなわち、骨のライブ画像Ｂ_Ｌのある画素の時系列に注目すると、骨の情報のみを含むフレームと、骨と造影剤の情報を含むフレームが存在する。具体的には、造影剤が注入された血管がある画素を横切ったとき、造影剤の分だけ骨のライブ画像Ｂ_Ｌの画素値が上昇する。従って、ｎフレーム前の骨のライブ画像Ｂ_Ｌ［－ｎ］から１フレーム前の骨のライブ画像Ｂ_Ｌ［－１］までのｎ枚の画像において、画素ごとに最小値をとると、動く造影剤の情報が除去されて、動かない骨の情報のみが残る。これを骨のマスク画像Ｂ_Ｍとして、現在のフレームの骨のライブ画像Ｂ_Ｌから引くと、骨のＤＳＡ画像Ｂ_ＤＳＡが得られる。骨のＤＳＡ画像Ｂ_ＤＳＡでは骨の情報が相殺されて、造影剤の情報のみが残る。このようにして、心臓の拍動による動きがあっても、造影剤のみを画像化して、血管の診断を行うことができる。

なお、第３実施形態の図１１の説明においては、心臓の冠動脈に造影剤を注入しながらライブ画像を撮影した場合を想定したが、本実施形態はそのような場合に限定されない。例えば、心臓の冠動脈にガイドワイヤやステント等を挿入しながら撮影した場合が考えられる。このとき、造影剤の場合と同様、心臓の拍動に伴って血管が動くため、血管内を通って挿入されたガイドワイヤやステントも移動する。従って、撮影により得られたライブ画像に対してエネルギーサブトラクションを行うと、三種類目の物質であるガイドワイヤやステント（鉄）の情報が骨のライブ画像Ｂ_Ｌに現れる。この場合も図１１の説明と同様、ガイドワイヤやステントが挿入された血管がある画素を横切ったとき、ガイドワイヤやステントの分だけ骨のライブ画像Ｂ_Ｌの画素値が上昇する。従って、ｎフレーム前の骨のライブ画像Ｂ_Ｌ［－ｎ］から１フレーム前の骨のライブ画像Ｂ_Ｌ［－１］までのｎ枚の画像において、画素ごとに最小値をとると、動くガイドワイヤやステントの情報が除去されて、動かない骨の情報のみが残る。従って骨のＤＳＡ画像Ｂ_ＤＳＡでは、骨の情報が相殺されて、ガイドワイヤやステントの情報のみが残る。すなわち、ガイドワイヤやステントを選択的に画像化することができる。同様の処理によって、血管内手術が完了して冠動脈内にステントが留置された状態でステントを選択的に画像化する形態や、心臓の冠動脈内に石灰化（カルシウム）が沈着しているときに、石灰化を選択的に画像化する形態が好適に用いられる。また、第３実施形態では、心臓の冠動脈を撮影する例を示したが、このような形態に限定されない。例えば、腹部大動脈や腎臓等の造影を行う際に、呼吸に伴う横隔膜の上下や、腸の蠕動などによって、ＤＳＡ画像にアーチファクトが生じるという課題がある。本実施形態によれば、そのような状況に適用することで横隔膜や腸が動いても、造影剤やステントや石灰化を選択的に画像化することができる。

第３実施形態の図１１では、ｎフレーム前の骨のライブ画像Ｂ_Ｌ［－ｎ］から１フレーム前の骨のライブ画像Ｂ_Ｌ［－１］までのｎ枚の画像（ｎ≧２）、すなわち過去の複数フレームの画像において、画素ごとに最小値をとることで、骨のマスク画像Ｂ_Ｍを生成した。しかしながら本実施形態はこのような形態に限定されない。例えば、過去の複数フレームの骨のライブ画像Ｂ_Ｌがあるとき、画素ごとに最大値や平均値、中央値などの統計量を計算することで、骨のマスク画像Ｂ_Ｍを生成してもよい。また、画素ごとに異常値を除去することで、骨のマスク画像Ｂ_Ｍを生成してもよい。例えば、骨のライブ画像Ｂ_Ｌでは、骨の画像にはほとんど移動がないものの造影剤が通過する血管が移動するため、造影剤による画素値の変動を異常値として判定することで、骨の画像を残すことができる。さらには、骨のライブ画像Ｂ_Ｌに対してリカーシブフィルタを適用することで、骨のマスク画像Ｂ_Ｍを生成してもよい。以上のように、過去の複数フレームの骨のライブ画像Ｂ_Ｌから、任意の方法で骨のマスク画像Ｂ_Ｍを生成してよい。

また、上記の説明では骨のライブ画像Ｂ_Ｌについて説明したが、本実施形態はこのような形態に限定されない。軟部組織のマスク画像Ｓ_Ｍを、過去の複数フレームの軟部組織のライブ画像Ｓ_Ｌから生成してもよい。

また、第２実施形態のように、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌとから実効原子番号の画像Ｚと面密度の画像Ｄに分離する場合にも、第３実施形態の上記処理を適用することができる。この場合は、過去の複数フレームの実効原子番号のライブ画像Ｚ_Ｌから実効原子番号のマスク画像Ｚ_Ｍを、過去の複数フレームの面密度のライブ画像Ｄ_Ｌから面密度のマスク画像Ｄ_Ｍを生成するようにすればよい。

以上のように、第３実施形態は、互いにエネルギーが異なる複数の放射線画像からエネルギーサブトラクションによって得た画像をエネルギーサブトラクション後のライブ画像とする。そして、過去の複数フレームのエネルギーサブトラクション後のライブ画像から、エネルギーサブトラクション後のマスク画像を生成する。そして、現在のフレームのエネルギーサブトラクション後のライブ画像から、上記のようにして生成されたエネルギーサブトラクション後のマスク画像を引くことで、エネルギーサブトラクション後のＤＳＡ画像を生成する。したがって、例えば、ステントや石灰化などの、軟部組織以外の動く物体がマスク画像に存在しても、ＤＳＡ画像におけるアーチファクトの発生を低減することができる。

＜第４実施形態＞
次に、第４実施形態を説明する。第１実施形態においては、図６（ａ）で示した信号処理における［数４］を、単色放射線で近似することで解いていた。しかしながら、一般に、放射線が物質を通過するとき、低エネルギーの方がより多く吸収されるため、結果的にエネルギーピークが高い側にシフトする、ビームハードニングという現象が発生する。したがって、被写体の厚みが増すと、このビームハードニングにより、放射線の平均エネルギーがより高エネルギー側にシフトする。これによって、骨の厚みＢと軟部組織の厚みＳに誤差が生じる可能性がある。そこで第４実施形態では、第１の物質である骨と第２の物質である軟部組織の厚みの変化により生じる減弱率の変化を考慮して、高エネルギー画像と低エネルギー画像から第１の物質と第２の物質の厚みを得るようにする。すなわち、以下で説明するようにして［数４］の非線形連立方程式を解くことで、骨の厚みＢと軟部組織の厚みＳを求める。第４実施形態の放射線撮影システム１００の構成及び動作は、第１実施形態～第３実施形態と同様である。

非線形連立方程式を解く代表的な方法として、ニュートンラフソン法を用いた場合について説明する。まず、ニュートンラフソン法の反復回数をｍ、ｍ回目の反復後の骨の厚みをＢ^ｍ、ｍ回目の反復後の軟部組織の厚みをＳ^ｍとしたとき、ｍ回目の反復後の高エネルギーの減弱率Ｈ^ｍ、および、ｍ回目の反復後の低エネルギーの減弱率Ｌ^ｍを、以下の式で表す。

また、厚みが微小に変化したときの減弱率の変化率を、以下の式で表す。

このとき、ｍ＋１回目の反復後の骨の厚みＢ^ｍ＋１と軟部組織の厚みＳ^ｍ＋１を、高エネルギーの減弱率Ｈと低エネルギーの減弱率Ｌを用いて、以下の式で表す。

２ｘ２の行列の逆行列は、行列式をｄｅｔとすると、クラメルの公式より以下の式で表される。

従って、［数２１］に［数２２］を代入すると、以下の式が求まる。

このような計算を繰り返すことで、ｍ回目の反復後の高エネルギーの減弱率Ｈ^ｍと実測した高エネルギーの減弱率Ｈの差分が限りなく０に近づいていく。低エネルギーの減弱率Ｌについても同様である。これによって、ｍ回目の反復後の骨の厚みＢ^ｍが骨の厚みＢに収束し、ｍ回目の軟部組織の厚みＳ^ｍが軟部組織の厚みＳに収束する。このようにして、［数４］に示した非線形連立方程式を解くことができる。

また、第２実施形態で示したような実効原子番号Ｚと面密度Ｄを求める計算についても、同様の解法が成り立つ。すなわち、実効原子番号と面密度の変化により生じる減弱率の変化を考慮して、高エネルギー画像と低エネルギー画像から実効原子番号の画像と面密度の画像を得ることができる。例えば、［数１０］の非線形連立方程式をニュートンラフソン法により解いた場合について説明する。まず、ニュートンラフソン法の反復回数をｍ、ｍ回目の反復後の実効原子番号をＺ^ｍ、ｍ回目の反復後の面密度をＤ^ｍとしたとき、ｍ回目の反復後の高エネルギーの減弱率をＨ^ｍ、ｍ回目の反復後の低エネルギーの減弱率Ｌ^ｍを、以下の式で表す。

また、実効原子番号と面密度が微小に変化したときの減弱率の変化率を、以下の式で表す。

このとき、ｍ＋１回目の反復後の骨の厚みＢ^ｍ＋１と軟部組織の厚みＳ^ｍ＋１を、高エネルギーの減弱率Ｈと低エネルギーの減弱率Ｌを用いて、以下の式で表す。

２ｘ２の行列の逆行列は、行列式をｄｅｔとすると、クラメルの公式より以下の式で表される。

従って、［数２６］に［数２７］を代入すると、以下の式が求まる。

このような計算を繰り返すことで、ｍ回目の反復後の高エネルギーの減弱率Ｈ^ｍと実測した高エネルギーの減弱率Ｈの差分が限りなく０に近づいていく。低エネルギーの減弱率Ｌについても同様である。これによって、ｍ回目の反復後の実効原子番号Ｚ^ｍが実効原子番号Ｚに収束し、ｍ回目の反復後の面密度Ｄ^ｍが面密度Ｄに収束する。以上のようにして、［数１０］に示した非線形連立方程式を解くことができる。

なお、第４実施形態においては、ニュートンラフソン法を用いて非線形連立方程式を解くことを説明したが、このような形態に限定されない。例えば、最小二乗法や二分法などの反復解法を用いてもよい。また、反復解法を用いるとき、０回目の反復後の骨の厚みＢ^０と軟部組織の厚みＳ^０、すなわち骨の厚みと軟部組織の厚みの初期値が必要となる。初期値としては任意の定数を与えてもよい。また、第１実施形態で説明したように、放射線のスペクトルを単色放射線で近似した［数８］を用いて骨の厚みと軟部組織の厚みを求め、これらを初期値とする構成が好適に用いられる。第２実施形態で示した実効原子番号Ｚと面密度Ｄの場合についても同様に、放射線のスペクトルを単色放射線で近似した［数１４］を用いて実効原子番号と面密度を求め、これらを初期値とする構成が好適に用いられる。以上のような計算により、骨の厚みＢや軟部組織の厚みＳ、実効原子番号Ｚや面密度Ｄを高精度に求めることができる。

また、上記では、［数４］や［数１０］で示した非線形連立方程式を反復解法で解いていた。しかしながら、この過程で［数１９］、［数２０］、［数２４］、［数２５］を計算するとき、数値積分が必要となる。しかも、ｍ回の反復計算を行うたびに再計算が必要である。さらには、全ての画素においてこのような演算を行うことになる。従って、図６（ａ）や図６（ｂ）で示したエネルギーサブトラクションの信号処理に、時間がかかるという課題がある。特に本実施形態では動画撮影を想定しているため、信号処理に許容される時間は１フレーム以下となる。例えばフレームレートが２０ｆｐｓの場合、５０ｍｓ以下で信号処理と画像処理を行わなければ、処理が追い付かなくなる。そこで様々な組み合わせの高エネルギーの減弱率Ｈと低エネルギーの減弱率Ｌに対する骨の厚みＢや軟部組織の厚みＳを事前に求めてテーブルを生成しておき、このテーブルを参照することで骨の厚みＢや軟部組織の厚みＳを高速に求めるようにしてもよい。

図１２に、第４実施形態におけるテーブル生成のフローチャートを示す。テーブルの分割数をＭ、テーブルの低エネルギーの座標をｌ、テーブルの高エネルギーの座標をｈとする。なお、テーブルの座標ｌ及びｈは整数である。まず、制御用コンピュータ１０３は、低エネルギーの座標をｌ＝０として初期化してから（Ｓ１２０１）、高エネルギーの座標をｈ＝０として初期化する（Ｓ１２０２）。

次に、制御用コンピュータ１０３は、座標ｌにおける低エネルギーの減弱率Ｌ［ｌ］と座標ｈにおける高エネルギーの減弱率Ｈ［ｈ］を以下の式により求める。

このようにして求めた低エネルギーの減弱率Ｌ［ｌ］と高エネルギーの減弱率をＨ［ｈ］に対して、［数４］で示した非線形連立方程式を解き、骨の厚みＢと軟部組織の厚みＳを求める。この結果を、骨の厚みＢに対するテーブルＢ［ｌ，ｈ］及び軟部組織の厚みＳに対するテーブルＳ［ｌ，ｈ］に格納する（Ｓ１２０３）。その後、制御用コンピュータ１０３は、ｈ＝ｈ＋１とする（Ｓ１２０４）。高エネルギーの座標ｈがテーブルの分割数Ｍを越えていない場合（Ｓ１２０５でＮＯ）は、制御用コンピュータ１０３は、Ｓ１２０３からの処理を繰り返す。高エネルギーの座標ｈがテーブルの分割数Ｍを越えた場合（Ｓ１２０５ＹＥＳ）は、制御用コンピュータ１０３は、ｌ＝ｌ＋１とする（Ｓ１２０６）。そして、低エネルギーの座標ｌがテーブルの分割数Ｍを越えていない場合（Ｓ１２０７でＮＯ）は、制御用コンピュータ１０３は、ｈ＝０として（Ｓ１２０２）、Ｓ１２０３～Ｓ１２０５の処理を繰り返す。低エネルギーの座標ｌがテーブルの分割数Ｍを越えた場合（Ｓ１２０７でＹＥＳ）は、テーブルの生成が終了する。以上のようにして、あらゆるｌとｈの組み合わせについて骨の厚みＢと軟部組織の厚みＳを求めて、テーブルに格納することができる。

図１３に、第４実施形態におけるテーブル参照の模式図を示す。このテーブルは、高エネルギー画像が示す減弱率と低エネルギー画像が示す減弱率の組合せに応じて骨と軟部組織の厚みが登録された２次元のテーブルである。ある画素において低エネルギーの減弱率がＬ、高エネルギーの減弱率がＨであるとする。このとき、以下の式によって減弱率をテーブルの座標ｌ'及びｈ'に変換する。

この座標を用いて、骨の厚みＢに対するテーブルＢ［ｌ，ｈ］を参照することで、骨の厚みＢを求めることができる。軟部組織の厚みについても同様である。ただし、テーブルの座標ｌ'及びｈ'は小数となるが、テーブルは配列に格納されているため、整数でなければ参照できない。従って、座標ｌ'及びｈ'を整数に変換してから、補間によって骨の厚みＢ及び軟部組織の厚みＳを求める構成が好適に用いられる。例えば、座標ｌ'の小数点を切り捨てて整数に変換した値をｌ、座標ｈ'の小数点を切り捨てて整数に変換した値をｈとし、バイリニア補間で骨の厚みＢと軟部組織の厚みＳを求める場合は、以下の式を用いる。

従って、事前にテーブルを生成しておけば、非線形連立方程式を解くのに比べて格段に少ない計算量で、骨の厚みＢと軟部組織の厚みＳを求めることができる。このようなテーブルは、放射線のスペクトルＮ（Ｅ）が変化しない限り有効である。一般的に、動画撮影中に放射線のスペクトルＮ（Ｅ）が変化することはないため、撮影前に一回テーブルを生成すれば十分である。当然のことながら、図１２と図１３で示したテーブルの生成と参照は、任意の二物質の厚みを求める計算についても同様に用いることができる。また、実効原子番号Ｚと面密度Ｄを求める計算に対しても、同様に用いることができる。この場合、高エネルギー画像が示す減弱率と低エネルギー画像が示す減弱率の組合せに応じて実効原子番号と面密度が登録された２次元のテーブルが用いられることになる。このようにして、図６（ａ）や図６（ｂ）で示したエネルギーサブトラクションの信号処理を高精度かつ高速に行うことができる。

なお、第４実施形態では、式（２９）を用いて、座標ｌにおける低エネルギーの減弱率Ｌ［ｌ］と座標ｈにおける高エネルギーの減弱率Ｈ［ｈ］を求めていた。このようにして生成されたテーブルは、縦軸の座標が高エネルギーの減弱率Ｈ、横軸の座標が低エネルギーの減弱率Ｌとなるが、テーブルは０～１の間の減弱率を等間隔に分割したものになる。しかしながら、人体の組成と厚みに対する減弱率は０付近の値となることが多い。このため、テーブルの分割数Ｍが小さいと、テーブルの参照と補間によって得た値と、非線形連立方程式を解いて得た値との誤差が大きくなるという欠点がある。そこで、座標の範囲を決定する定数をｋ（０＜ｋ）としたとき、以下の式を用いて減弱率を求める構成が好適に用いられる。

［数３２］によって座標から減弱率を求めた場合は、以下の式を用いて座標を求める。

テーブルの参照と補間については、［数３１］を用いる。このようにして生成されたテーブルは、縦軸の座標が－ｌｎ（Ｈ）、横軸の座標が－ｌｎ（Ｌ）となっている。従って、減弱率の値が０付近になっても、テーブルの参照と補間によって得た値と、非線形連立方程式を解いて得た値との誤差を小さくできる。

また、第４実施形態において、テーブルの生成と参照を行うとき、解が存在しない高エネルギーの減弱率Ｈと低エネルギーの減弱率Ｌの組み合わせが生じうる。例えば、高エネルギーの減弱率Ｈは、通常、低エネルギーの減弱率Ｌよりも大きくなる。従って、［数２９］や［数３２］を用いて生成されたテーブルでは、Ｈ＜Ｌとなる領域の解が求められない。本実施形態においては、テーブル上で解が求められない領域を減らすように、座標を選択してもよい。例えば、縦軸の座標をｌｎ（Ｌ）／ｌｎ（Ｈ）、横軸の座標を－ｌｎ（Ｈ）などにしてもよい。あるいは、縦軸の座標として、［数８］や［数１４］のように、単色放射線で近似して解いた値を用いるなどしてもよい。また、第４実施形態においてテーブルの参照を行う際に、テーブルの範囲外の座標が指定される可能性や、解が求められない領域が参照される可能性がある。このような場合は、指定された座標の周囲にある、解が存在する領域の値を用いる構成が好適に用いられる。

なお、第１実施形態では、図７のようにエネルギーサブトラクションを行ってからＤＳＡを行う順序で画像処理を行っていた。しかしながら、［数８］のように、放射線のスペクトルＮ（Ｅ）を単色放射線で近似した場合は、ＤＳＡを行ってからエネルギーサブトラクションを行う構成としてもよい。例えば、［数８］をマスク画像とライブ画像の両方に適用すると、厚みＢと軟部組織の厚みＳは、以下の式で表される。

［数３４］より、ＤＳＡ画像は以下の式で表される。

すなわち、エネルギーサブトラクションを行った後にライブ画像からマスク画像を減算してＤＳＡを行った結果と、ライブ画像をマスク画像で除算することでＤＳＡを行った後にエネルギーサブトラクションを行った結果が一致する。しかしながら、第４実施形態のように、非線形連立方程式を画素ごと解く構成、又は事前に解いてテーブルを生成してから参照する構成では、このような演算の順不同性は保証されない。すなわち、ＤＳＡを行ってからエネルギーサブトラクションを行うと、誤差が大きくなるという課題がある。一方で、非線形方程式を画素ごとに解いてエネルギーサブトラクション行った後にＤＳＡを行う構成にすると、誤差は小さくなるが計算量が多くなるという課題がある。このような理由から、第４実施形態で示したように、非線形連立方程式を解くことで、互いに異なる複数のエネルギーにおける画素値とエネルギーサブトラクション後の画素値との対応関係を示すテーブルを生成して置くことが望ましい。画像処理部１３２は、こうして生成されたテーブルを参照することでエネルギーサブトラクション後の画像を生成し、ＤＳＡ処理部１３３は、エネルギーサブトラクション後の画像を用いてＤＳＡを行う。

なお、第１～第４実施形態では、放射線撮影装置１０４は蛍光体を用いた間接型の放射線センサとしたが、このような形態に限定されない。例えばＣｄＴｅ等の直接変換材料を用いた直接型の放射線センサを用いてもよい。また、第１～第４実施形態では放射線発生装置１０１の受動的な管電圧変化を利用するか、能動的に管電圧を切り替えるなどしていたが、このような形態に限定されない。放射線発生装置１０１のフィルタを時間的に切り替えるなどして、放射線撮影装置１０４に曝射される放射線のエネルギーを変化させてもよい。

さらに、第１～第４実施形態では放射線撮影装置１０４に曝射される放射線のエネルギーを変化させることで、エネルギーサブトラクションを行っていたがこのような形態に限定されない。例えば、第１実施形態でも述べたが、センサを２枚積層することで、前面のセンサと背面のセンサで検出する放射線のスペクトルを変化させる方式が用いられてもよい。また、放射線量子の個数をエネルギー別にカウントする、フォトンカウンティング方式のセンサを用いることで、互いにエネルギーが異なる複数の画像を取得するなどしてもよい。

また、第１と第２実施形態では、放射線撮影システムの制御用コンピュータ１０３を用いてエネルギーサブトラクション処理を行っていたが、このような形態に限定されない。制御用コンピュータ１０３で取得した画像を別のコンピュータに転送して、エネルギーサブトラクション処理を行ってもよい。例えば、取得した画像を医療用のＰＡＣＳを介して別のパソコンに転送し、エネルギーサブトラクション処理を行ってから表示する構成が好適に用いられる。すなわち、上記各実施形態では、互いにエネルギーが異なる放射線画像をエネルギーサブトラクション処理に提供できればよく、互いにエネルギーが異なる放射線画像を取得するための方法は、上記実施形態に限定されるものではない。

本発明は、上述の実施形態の１以上の機能を実現するプログラムを、ネットワーク又は記憶媒体を介してシステム又は装置に供給し、そのシステム又は装置のコンピュータにおける１つ以上のプロセッサーがプログラムを読出し実行する処理でも実現可能である。また、１以上の機能を実現する回路（例えば、ＡＳＩＣ）によっても実現可能である。

１０１：放射線発生装置、１０２：放射線制御装置、１０３：制御用コンピュータ、１０４放射線撮影装置、１０５：蛍光体、１０６：二次元検出器

Claims (19)

1. 被写体を透過した互いに異なる放射線エネルギーによる複数の画像をそれぞれ動画撮影により取得する取得手段と、
   前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内である第１の時間において得られた複数の画像を用いて第１のエネルギーサブトラクション画像を生成し、前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内であり且つ前記第１の時間の後である第２の時間において得られた複数の画像を用いて第２のエネルギーサブトラクション画像を生成する第１の生成手段であって、前記第１のエネルギーサブトラクション画像と前記第２のエネルギーサブトラクション画像とを含む複数のエネルギーサブトラクション画像であって、前記第１の時間から前記第２の時間までの期間において得られた前記複数のエネルギーサブトラクション画像を前記動画撮影の期間内に生成する前記第１の生成手段と、
   前記複数のエネルギーサブトラクション画像に基づくマスク画像と、前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内であり且つ前記第１の時間及び前記第２の時間の後である第３の時間において得られた複数の画像を用いて生成された第３のエネルギーサブトラクション画像に基づくライブ画像との差分画像を前記動画撮影の期間内に生成する第２の生成手段と、を備えることを特徴とする撮影制御装置。
2. 前記取得手段は、１ショット間にエネルギーが変化する放射線が被写体に照射され、当該１ショット間に前記被写体を透過した前記放射線を複数回検出して得られた、互いに異なるエネルギーを有する前記複数の画像を取得することを特徴とする請求項１に記載の撮影制御装置。
3. 前記取得手段は、連続する２つのリセットに挟まれた期間に、１つの二次元検出器において少なくとも２回のサンプリングを行うことにより、前記複数の画像を１フレームの撮影期間において取得することを特徴とする請求項１に記載の撮影制御装置。
4. 前記第１の生成手段は、前記複数の画像を用いてエネルギーサブトラクション処理を行って、第１の物質の厚さに関わるエネルギーサブトラクション画像と、前記第１の物質とは異なる第２の物質の厚さに関わる第２のエネルギーサブトラクション画像を生成し、
   前記第２の生成手段は、前記第１の物質の厚さに関わるエネルギーサブトラクション画像を用いて前記差分画像を生成することを特徴とする請求項１乃至３のいずれか１項に記載の撮影制御装置。
5. 前記第１の物質は骨であり、前記第２の物質は軟部組織であることを特徴とする請求項４に記載の撮影制御装置。
6. 前記第１の生成手段は、前記複数の画像についてエネルギーサブトラクション処理を行って、実効原子番号の画像と面密度の画像を前記エネルギーサブトラクション画像として取得し、
   前記第２の生成手段は、前記実効原子番号の画像を用いて前記差分画像を生成することを特徴とする請求項１乃至３のいずれか１項に記載の撮影制御装置。
7. 前記第１の生成手段は、前記複数の画像を用いてエネルギーサブトラクション処理を行って、実効原子番号の画像と面密度の画像を前記エネルギーサブトラクション画像として取得し、
   前記第２の生成手段は、前記実効原子番号の画像と面密度の画像を用いて前記差分画像を生成することを特徴とする請求項１または２に記載の撮影制御装置。
8. 前記第２の生成手段は、前記実効原子番号の画像と前記面密度の画像に基づいて生成されたマスク画像をそれぞれＺＭ、ＤＭとし、前記実効原子番号の画像と前記面密度の画像に基づいて生成されたライブ画像をそれぞれＺＬ、ＤＬとした場合に、面密度の差分画像であるＤＤＳＡと実行原子番号の差分画像であるＺＤＳＡを、
   

   

   （ただし、ｎは２．５以上、３以下の実数）
   により生成することを特徴とする請求項７に記載の撮影制御装置。
9. 前記第２の生成手段は、過去に生成された前記複数のエネルギーサブトラクション画像の画素ごとの画素値の統計量を計算し、前記マスク画像の画素値とすることを特徴とする請求項１乃至８のいずれか１項に記載の撮影制御装置。
10. 前記統計量は、最小値、最大値、平均値、中央値のいずれかであることを特徴とする請求項９に記載の撮影制御装置。
11. 前記第２の生成手段は、過去に生成された前記複数のエネルギーサブトラクション画像に基づいて、画素ごとに異常値を除去することで前記マスク画像を生成することを特徴とする請求項１乃至８のいずれか１項に記載の撮影制御装置。
12. 前記第２の生成手段は、過去に生成された前記複数のエネルギーサブトラクション画像についてリカーシブフィルタを適用することで前記マスク画像を生成することを特徴とする請求項１乃至８のいずれか１項に記載の撮影制御装置。
13. 前記複数の画像は、第１の放射線エネルギーによる第１の画像と前記第１の放射線エネルギーとは異なる第２の放射線エネルギーによる第２の画像を含み、
    前記第１の生成手段は、前記エネルギーサブトラクション画像を構成する画素の値の変化により生じる減弱率の変化を考慮して、前記第１の画像と前記第２の画像から前記エネルギーサブトラクション画像を生成することを特徴とする請求項１乃至１２のいずれか１項に記載の撮影制御装置。
14. 前記第１の生成手段は、前記第１の画像と前記第２の画像が示す減弱率の組合せに応じて、前記エネルギーサブトラクション画像における画素の値が登録された２次元のテーブルを参照することを特徴とする請求項１３に記載の撮影制御装置。
15. 前記複数のエネルギーサブトラクション画像は骨の厚さに関わるエネルギーサブトラクション画像であって、
    前記第１の生成手段は、前記複数のエネルギーサブトラクション画像において、画素ごとに最小値を用いて前記マスク画像を生成することを特徴とする請求項１乃至１４のいずれか1項に記載の撮影制御装置。
16. 被写体を透過した互いに異なる放射線エネルギーにより動画撮影の期間である第１の時間において得られた複数の画像を用いて生成された第１のエネルギーサブトラクション画像と、前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内であり且つ前記第１の時間の後である第２の時間において得られた複数の画像を用いて生成された第２のエネルギーサブトラクション画像とを含む複数のエネルギーサブトラクション画像であって、前記第１の時間から前記第２の時間までの期間において得られた前記複数のエネルギーサブトラクション画像を前記動画撮影の期間内に生成する生成手段であって、前記複数のエネルギーサブトラクション画像に基づくマスク画像と、前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内であり且つ前記第１の時間及び前記第２の時間の後である第３の時間において得られた複数の画像を用いて生成された第３のエネルギーサブトラクション画像に基づくライブ画像との差分画像を前記動画撮影の期間内に生成する生成手段を備えることを特徴とする画像処理装置。
17. 放射線撮影装置のための撮影制御装置の作動方法であって、
    被写体を透過した互いに異なる放射線エネルギーによる複数の画像をそれぞれ動画撮影により取得する取得工程と、
    前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内である第１の時間において得られた複数の画像を用いて第１のエネルギーサブトラクション画像を生成し、前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内であり且つ前記第１の時間の後である第２の時間において得られた複数の画像を用いて第２のエネルギーサブトラクション画像を生成する第１の生成工程であって、前記第１のエネルギーサブトラクション画像と前記第２のエネルギーサブトラクション画像とを含む複数のエネルギーサブトラクション画像であって、前記第１の時間から前記第２の時間までの期間において得られた前記複数のエネルギーサブトラクション画像を前記動画撮影の期間内に生成する前記第１の生成工程と、
    前記複数のエネルギーサブトラクション画像に基づくマスク画像と、前記被写体を透過した互いに異なる放射線エネルギーにより前記動画撮影の期間内であり且つ前記第１の時間及び前記第２の時間の後である第３の時間において得られた複数の画像を用いて生成された第３のエネルギーサブトラクション画像に基づくライブ画像との差分画像を前記動画撮影の期間内に生成する第２の生成工程と、を備えることを特徴とする撮影制御装置の作動方法。
18. 請求項１乃至１５のいずれか１項に記載された撮影制御装置と、
    二次元検出器を含む放射線撮影装置と、
    放射線を発生する放射線発生装置と、を備えることを特徴とする放射線撮影システム。
19. コンピュータを、請求項１乃至１５のいずれか１項に記載された撮影制御装置の各手段として機能させるためのプログラム。

Images