JP2022056084A - 画像処理装置、画像処理方法及びプログラム - Google Patents

Abstract

【課題】エネルギーサブトラクション後の画像に生じ得る不要な構造の消し残りを低減すること。【解決手段】被写体に含まれる第１の物質の画像と第２の物質の画像を生成する画像処理装置は、被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する取得部と、複数の放射線画像からエネルギーサブトラクション処理で用いられる複数の減弱率の画像を生成する生成部と、複数の減弱率の画像を用いたエネルギーサブトラクション処理により実効原子番号の画像と面密度の画像を生成する処理部と、を備える。処理部は、実効原子番号の画像における画素値から第１の物質の減弱情報を推定し、複数の減弱率の画像と減弱情報とを用いたエネルギーサブトラクション処理により、第１の物質の画像と第２の物質の画像とを生成する。【選択図】 図１

Description

本発明は、画像処理装置、画像処理方法及びプログラムに関するものである。より具体的には、医療診断における一般撮影などの静止画撮影や透視撮影などの動画撮影に用いられる放射線撮像装置及び放射線撮像システムで撮影された画像を処理する画像処理装置、画像処理方法及びプログラムに関する。

近年、Ｘ線による医療画像診断や非破壊検査に用いる撮影装置として、平面検出器（Flat Panel Detector、以下、「ＦＰＤ」と略す）を用いた放射線撮像装置が普及している。このような放射線撮像装置は、例えば医療画像診断においては、一般撮影のような静止画撮影や、透視撮影のような動画撮影のデジタル撮像装置として用いられている。

ＦＰＤを用いた撮影方法として、エネルギーサブトラクションがある。エネルギーサブトラクションでは、まず、管電圧の異なるＸ線を照射するなどして、エネルギーの異なる画像を複数枚取得する。それらを演算することで、骨画像と軟部組織画像に分離するなどの処理を行うことができる（特許文献１）。

特開２００８－１６７９４８号公報

エネルギーサブトラクションを用いて骨画像や軟部組織画像を生成したとき、骨画像に軟部組織の構造が残ったり、軟部組織画像に骨の構造が残ったりすることがある。このような不要な構造の消し残りが診断の妨げになるという課題がある。

本発明は、上記の課題に鑑み、エネルギーサブトラクション後の画像に生じ得る不要な構造の消し残りを低減することが可能な画像処理技術を提供する。

本発明の一態様による画像処理装置は以下の構成を備える。すなわち、画像処理装置は、被写体に含まれる第１の物質の画像と第２の物質の画像を生成する画像処理装置であって、
前記被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する取得手段と、
前記複数の放射線画像からエネルギーサブトラクション処理で用いられる複数の減弱率の画像を生成する生成手段と、
前記複数の減弱率の画像を用いた前記エネルギーサブトラクション処理により実効原子番号の画像と面密度の画像を生成する処理手段と、を備え、
前記処理手段は、前記実効原子番号の画像における画素値から前記第１の物質の減弱情報を推定し、
前記複数の減弱率の画像と前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記第１の物質の画像と前記第２の物質の画像とを生成する。

本発明の他の態様による画像処理装置は以下の構成を備える。すなわち、画像処理装置は、被写体に含まれる第１の物質の画像と第２の物質の画像を生成する画像処理装置であって、
前記被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する取得手段と、
前記複数の放射線画像からエネルギーサブトラクション処理で用いられる複数の減弱率の画像を生成する生成手段と、
前記複数の減弱率の画像を用いた前記エネルギーサブトラクション処理により前記第１の物質に含まれる第３の物質の厚みの画像と第４の物質の厚みの画像とに分離した画像を生成する処理手段と、を備え、
前記処理手段は、
前記第３の物質の厚みの画像と前記第４の物質の厚みの画像とから取得した前記第１の物質に含まれる前記第３の物質の比率情報から前記第１の物質の減弱情報を推定し、
前記複数の減弱率の画像と前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記第１の物質の画像と前記第２の物質の画像とを生成する。

本発明によれば、エネルギーサブトラクション後の画像に生じ得る不要な構造の消し残りを低減することが可能になる。

実施形態に係るＸ線撮影システムの構成例を示す図。 実施形態に係るＸ線撮像装置の画素等価回路図。 実施形態に係るＸ線撮像装置のタイミングチャート。 実施形態に係るＸ線撮像装置のタイミングチャート。 実施形態に係る補正処理を説明する図。 実施形態に係る信号処理のブロック図。 実施形態に係る信号処理のブロック図。 実施形態に係るテーブル生成の流れを説明するフローチャート。 実施形態に係るテーブルの模式図。 実施形態に係る画像処理のブロック図。 実施形態に係る画像処理の流れを説明するフローチャート。 実施形態に係る画像処理の模式図。 実際の水の厚みと計算された骨の厚みの関係を示す図。 実施形態に係る画像処理の模式図。 実施形態に係る信号処理のブロック図。 信号処理のブロック図 実際の筋肉の厚みと計算された骨の厚みの関係を示す図。 第２実施形態に係る信号処理のブロック図

以下、添付図面を参照して実施形態を詳しく説明する。なお、以下の実施形態は特許請求の範囲に係る発明を限定するものではない。実施形態には複数の特徴が記載されているが、これらの複数の特徴の全てが発明に必須のものとは限らず、また、複数の特徴は任意に組み合わせられてもよい。さらに、添付図面においては、同一若しくは同様の構成に同一の参照番号を付し、重複した説明は省略する。

なお、本発明における放射線には、放射線崩壊によって放出される粒子（光子を含む）の作るビームであるα線、β線、γ線などの他に、同程度以上のエネルギーを有するビーム、例えばＸ線や粒子線、宇宙線なども、含まれるものとする。以下の実施形態では、放射線の一例としてＸ線を用いた装置を説明する。したがって、以下では、放射線、放射線画像、放射線エネルギー、放射線スペクトル、放射線量、放射線発生装置、放射線撮像装置、放射線撮像システムとして、それぞれＸ線、Ｘ線画像、Ｘ線エネルギー、Ｘ線スペクトル、Ｘ線量、Ｘ線発生装置、Ｘ線撮像装置、Ｘ線撮像システムとして説明する。

（第１実施形態）
図１は、第１実施形態に係る、放射線撮像システムの一例としてのＸ線撮像システムの構成例を示すブロック図である。第１実施形態のＸ線撮像システムは、Ｘ線発生装置１０１、Ｘ線制御装置１０２、撮像制御装置１０３、Ｘ線撮像装置１０４を備える。

Ｘ線発生装置１０１は、Ｘ線を発生し、被写体にＸ線を照射する。Ｘ線制御装置１０２は、Ｘ線発生装置１０１におけるＸ線の発生を制御する。撮像制御装置１０３は、例えば、１つまたは複数のプロセッサー（ＣＰＵ）とメモリを有し、プロセッサーがメモリに格納されたプログラムを実行してＸ線画像の取得及び画像処理を行う。なお、撮像制御装置１０３による画像処理を含む各処理は、専用のハードウエアにより実現されてもよいし、ハードウエアとソフトウエアの協働により実現されてもよい。Ｘ線撮像装置１０４は、Ｘ線を可視光に変換する蛍光体１０５と、可視光を検出する二次元検出器１０６を有する。二次元検出器は、Ｘ線量子を検出する画素２０をＸ列×Ｙ行のアレイ状に配置したセンサであり、画像情報を出力する。

撮像制御装置１０３は、上述したプロセッサーにより放射線画像を処理する画像処理装置として機能する。取得部１３１、補正部１３２、信号処理部１３３、画像処理部１３４は、画像処理装置としての機能構成例を示している。画像処理装置は、放射線画像を処理して被写体に含まれる第１の物質の画像と第２の物質の画像を生成する。ここで、被写体に含まれる第１の物質は、例えば軟部組織であり、第２の物質は、例えば骨である。取得部１３１は、被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する。補正部１３２は、取得部１３１により取得された複数の放射線画像を補正してエネルギーサブトラクション処理で用いられる複数の減弱率の画像（例えば、図５に示すような低エネルギーにおける減弱率の画像Ｌ、高エネルギーにおける減弱率の画像Ｈ）を生成する。

信号処理部１３３は、補正部１３２により生成された複数の減弱率の画像を用いたエネルギーサブトラクション処理により物質特性画像を生成する。物質特性画像とは、例えば、第１の物質（例えば、軟部組織）と第２の物質（例えば、骨）というように物質を分離して表す物質分離画像、実効原子番号とその面密度を表す物質識別画像など、エネルギーサブトラクション処理において取得される画像である。信号処理部１３３の詳細は後述する。画像処理部１３４は、信号処理部１３３により取得された物質特性画像を用いて、ＤＳＡ画像を取得し、エネルギーサブトラクションとＤＳＡを組み合わせた動画像を生成する処理を行う。画像処理部１３４の詳細は後述する。

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

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

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

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

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

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

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

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

Ｘ線撮像装置１０４は、二次元検出器１０６から以上のような画素回路の出力を読み出し、不図示のＡＤ変換器でデジタル値に変換した後、撮像制御装置１０３に画像を転送する。

次に、上述した構成を備えた第１実施形態のＸ線撮像システムの動作について説明する。図３は、第１実施形態に係るＸ線撮像システムにおいてエネルギーサブトラクションに提供するための、互いにエネルギーの異なる複数のＸ線画像を得る場合のＸ線撮像装置１０４の駆動タイミングを示す。図３中の波形は横軸を時間として、Ｘ線の曝射、同期信号、光電変換素子２０１のリセット、サンプルホールド回路２０７、信号線２１からの画像の読み出しのタイミングを示している。

リセット信号により光電変換素子２０１のリセットが行われてからＸ線が曝射される。Ｘ線の管電圧は理想的には矩形波となるが、管電圧の立ち上がりと立下りには有限の時間がかかる。特に、パルスＸ線で曝射時間が短い場合は、管電圧はもはや矩形波とはみなせず、Ｘ線３０１～３０３に示すような波形となる。立ち上がり期のＸ線３０１、安定期のＸ線３０２、立下り期のＸ線３０３ではそれぞれＸ線のエネルギーが異なる。したがって、サンプルホールドによって区切られる期間の放射線に対応したＸ線画像を得ることにより、互いにエネルギーが異なる複数種類のＸ線画像が得られる。

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

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

サンプルホールド回路２０７及び光電変換素子２０１のリセットを行うタイミングは、Ｘ線発生装置１０１からＸ線の曝射が開始されたことを示す同期信号３０７を用いて決定される。Ｘ線の曝射開始を検出する方法としては、Ｘ線発生装置１０１の管電流を測定し、電流値が予め設定された閾値を上回るか否かを判定する構成を用いることができるがこれに限られるものではない。例えば、光電変換素子２０１のリセットが完了した後、画素２０を繰り返して読み出し、画素値が予め設定された閾値を上回るか否かを判定することによりＸ線の曝射開始を検出する構成が用いられてもよい。

あるいは、例えば、Ｘ線撮像装置１０４に二次元検出器１０６とは異なるＸ線検出器を内蔵し、その測定値が予め設定された閾値を上回るか否かを判定することによりＸ線の曝射開始を検出する構成が用いられてもよい。いずれの方式の場合も、Ｘ線の曝射開始を示す同期信号３０７の入力から予め指定した時間が経過した後に、信号サンプルホールド回路２０７Ｓのサンプリング、ノイズサンプルホールド回路２０７Ｎのサンプリング、光電変換素子２０１のリセットが行われる。

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

図４は、第１実施形態に係るＸ線撮像システムにおいてエネルギーサブトラクションに提供するための、互いにエネルギーの異なる複数のＸ線画像を得る、図３とは異なるＸ線撮像装置１０４の駆動タイミングを示す。図３とは、Ｘ線発生装置１０１の管電圧を能動的に切り替えている点で異なる。

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

次に、Ｘ線撮像装置１０４は、低エネルギーのＸ線４０３の曝射と、画像４０４の読み出しとが完了してから、再び信号サンプルホールド回路２０７Ｓでサンプリングを行う。その後、Ｘ線撮像装置１０４は、光電変換素子２０１のリセットを行い、再びノイズサンプルホールド回路２０７Ｎでサンプリングを行い、信号線２１Ｎと信号線２１Ｓの差分を画像として読み出す。このとき、ノイズサンプルホールド回路２０７ＮにはＸ線が曝射されていない状態の信号が保持され、信号サンプルホールド回路２０７Ｓには低エネルギーのＸ線４０１の信号と高エネルギーのＸ線４０２と低エネルギーのＸ線４０３の信号（Ｒ_２）の和（Ｒ_１＋Ｂ＋Ｒ_２）が保持されている。従って、低エネルギーのＸ線４０１の信号と高エネルギーのＸ線４０２の信号と低エネルギーのＸ線４０３の信号に対応した画像４０６が読み出される。

その後、画像４０６と画像４０４の差分を計算することで、低エネルギーのＸ線４０１と低エネルギーのＸ線４０３の和に対応した画像４０５が得られる。この計算は、Ｘ線撮像装置１０４で行われてもよいし、撮像制御装置１０３で行われてもよい。同期信号４０７については、図３と同様である。このように、管電圧を能動的に切り替えながら画像を取得することで、図３の方法に比べて低エネルギーと高エネルギーの放射線画像の間のエネルギー差をより大きくすることが出来る。

次に、撮像制御装置１０３によるエネルギーサブトラクション処理について説明する。第１実施形態におけるエネルギーサブトラクション処理は、補正部１３２による補正処理、信号処理部１３３による信号処理、画像処理部１３４による画像処理の３段階に分かれている。以下、それぞれの処理について説明する。

（補正処理の説明）
補正処理は、Ｘ線撮像装置１０４から取得された複数の放射線画像を処理してエネルギーサブトラクション処理における後述の信号処理で用いられる複数の画像を生成する処理である。図５に、実施形態に係るエネルギーサブトラクション処理のための補正処理を示す。まず、取得部１３１は、Ｘ線撮像装置１０４にＸ線を曝射しない状態での撮像を行わせ、図３または図４に示した駆動で画像を取得する。この駆動により２枚の画像が読み出される。以下、１枚目の画像（画像３０４または画像４０４）をＦ＿ＯＤＤ、２枚目の画像（画像３０６または画像４０６）をＦ＿ＥＶＥＮとする。Ｆ＿ＯＤＤとＦ＿ＥＶＥＮは、Ｘ線撮像装置１０４の固定パターンノイズ（ＦＰＮ：Fixed Pattern Noise）に対応する画像である。

次に、取得部１３１は、被写体がない状態でＸ線撮像装置１０４にＸ線を曝射して撮像を行わせ、図３又は図４に示した駆動によりＸ線撮像装置１０４から出力されるゲイン補正用の画像を取得する。この駆動により、上記と同様に２枚の画像が読み出される。以下、１枚目のゲイン補正用の画像（画像３０４または画像４０４）をＷ＿ＯＤＤ、２枚目のゲイン補正用の画像（画像３０６または画像４０６）をＷ＿ＥＶＥＮとする。Ｗ＿ＯＤＤとＷ＿ＥＶＥＮは、Ｘ線撮像装置１０４のＦＰＮとＸ線による信号の和に対応する画像である。補正部１３２は、Ｗ＿ＯＤＤからＦ＿ＯＤＤを、Ｗ＿ＥＶＥＮからＦ＿ＥＶＥＮを減算することで、Ｘ線撮像装置１０４のＦＰＮが除去された画像ＷＦ＿ＯＤＤとＷＦ＿ＥＶＥＮを得る。これをオフセット補正と呼ぶ。

ＷＦ＿ＯＤＤは安定期のＸ線３０２に対応する画像であり、ＷＦ＿ＥＶＥＮは立ち上がり期のＸ線３０１、安定期のＸ線３０２、立下り期のＸ線３０３の和に対応する画像である。従って、補正部１３２は、ＷＦ＿ＥＶＥＮからＷＦ＿ＯＤＤを減算することで、立ち上がり期のＸ線３０１と立下り期のＸ線３０３の和に対応する画像を得る。このように、複数枚の画像の減算により、サンプルホールドによって区切られる特定の期間のＸ線に対応した画像を得る処理を色補正と呼ぶ。立ち上がり期のＸ線３０１と立下り期のＸ線３０３のエネルギーは、安定期のＸ線３０２のエネルギーに比べて低い。従って、色補正により、ＷＦ＿ＥＶＥＮからＷＦ＿ＯＤＤを減算することで、被写体がない場合の低エネルギー画像Ｗ＿Ｌｏｗが得られる。また、ＷＦ＿ＯＤＤから、被写体がない場合の高エネルギー画像Ｗ＿Ｈｉｇｈが得られる。

次に、取得部１３１は、被写体がある状態でＸ線撮像装置１０４にＸ線を曝射して撮像を行わせ、図３または図４に示した駆動によりＸ線撮像装置１０４から出力される画像を取得する。このとき２枚の画像が読み出される。以下、１枚目の画像（画像３０４または画像４０４）をＸ＿ＯＤＤ、２枚目の画像（画像３０６または画像４０６）をＸ＿ＥＶＥＮとする。補正部１３２は、被写体がない場合と同様のオフセット補正および色補正を行うことで、被写体がある場合の低エネルギー画像Ｘ＿Ｌｏｗと、被写体がある場合の高エネルギー画像Ｘ＿Ｈｉｇｈを得る。

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

［数１］式を変形すると、以下の［数２］式が得られる。［数２］式の右辺は被写体の減弱率を示す。被写体の減弱率は０～１の間の実数である。

従って、補正部１３２は、被写体がある場合の低エネルギー画像Ｘ＿Ｌｏｗを、被写体がない場合の低エネルギー画像Ｗ＿Ｌｏｗで除算することで、低エネルギーにおける減弱率の画像Ｌを得る。同様に、補正部１３２は、被写体がある場合の高エネルギー画像Ｘ＿Ｈｉｇｈを、被写体がない場合の高エネルギー画像Ｗ＿Ｈｉｇｈで除算することで、高エネルギーにおける減弱率の画像Ｈを得る。このように、被写体ありの状態で得られた放射線画像に基づいて得られた画像を被写体なしの状態で得られた放射線画像に基づいて得られた画像で除算することにより減弱率の画像を取得する処理をゲイン補正と呼ぶ。以上が、第１実施形態の補正部１３２による補正処理の説明である。

（信号処理の説明）
図６に、本実施形態に係るエネルギーサブトラクション処理の信号処理のブロック図を示す。信号処理部１３３は、補正部１３２から得られる複数の画像を用いて物質特性画像を生成する。以下では、骨の厚さの画像Ｂと軟部組織の厚さの画像Ｓからなる物質分離画像の生成を説明する。信号処理部１３３は、以下の処理により図５に示した補正によって得られた低エネルギーにおける減弱率の画像Ｌと高エネルギーにおける減弱率の画像Ｈから、骨の厚さの画像Ｂと軟部組織の厚さの画像Ｓを求める。

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

エネルギーＥにおけるフォトン数Ｎ（Ｅ）は、Ｘ線のスペクトルである。Ｘ線のスペクトルは、シミュレーション又は実測により得られる。また、エネルギーＥおける骨の線減弱係数μ_Ｂ（Ｅ）とエネルギーＥおける軟部組織の線減弱係数μ_Ｓ（Ｅ）は、それぞれＮＩＳＴ（National Institute of Standards and Technology）などのデータベースから得られる。したがって、［数３］によれば、任意の骨の厚さＢ、軟部組織の厚さＳ、Ｘ線のスペクトルＮ（Ｅ）における減弱率Ｉ／Ｉ_０を計算することが可能である。

ここで、低エネルギーのＸ線におけるスペクトルをＮ_Ｌ（Ｅ）、高エネルギーのＸ線におけるスペクトルをＮ_Ｈ（Ｅ）とすると、以下の［数４］式の各式が成り立つ。なお、Ｌは低エネルギーの減弱率の画像における画素値、Ｈは高エネルギーの減弱率の画像における画素値である。

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

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

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

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

従って、［数７］式に［数８］式を代入すると、以下の[数９]式が求まる。

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

なお、第１実施形態では、骨の厚みの画像Ｂと軟部組織の厚みの画像Ｓを算出していたが、本発明はこのような形態に限定されない。例えば、水の厚みＷと造影剤の厚みＩを算出してもよい。すなわち、任意の二種類の物質の厚さに分解してもよい。

また、本実施形態では、ニュートンラフソン法を用いて非線形連立方程式を解いていた。しかしながら本発明はこのような形態に限定されない。例えば、最小二乗法や二分法などの反復解法を用いてもよい。また、第１実施形態では非線形連立方程式を反復解法で解いていたが、本発明はこのような形態に限定されない。様々な組み合わせの高エネルギーの減弱率Ｈと低エネルギーの減弱率Ｌに対する骨の厚みＢや軟部組織の厚みＳを事前に求めてテーブルを生成し、このテーブルを参照することで骨の厚みＢや軟部組織の厚みＳを高速に求める構成を用いても良い。

図７に、本実施形態に係るエネルギーサブトラクション処理の信号処理のブロック図を示す。信号処理部１３３は、補正部１３２から得られる複数の画像を用いて物質特性画像を生成する。以下では、実効原子番号Ｚと面密度Ｄからなる画像の生成を説明する。信号処理部１３３は、以下の処理により図５に示した補正処理によって得られた低エネルギーにおける減弱率の画像Ｌと高エネルギーにおける減弱率の画像Ｈから、実効原子番号Ｚの画像と面密度Ｄの画像を求める。実効原子番号Ｚとは混合物の等価的な原子番号のことであり、面密度Ｄ［ｇ／ｃｍ^２］とは被写体の密度［ｇ／ｃｍ^３］と被写体の厚み［ｃｍ］の積である。

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

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

ここで、低エネルギーのＸ線のスペクトルをＮ_Ｌ（Ｅ）、高エネルギーのＸ線のスペクトルをＮ_Ｈ（Ｅ）とすると、以下の[数１１]式が成り立つ。なお、Ｌは低エネルギーの減弱率の画像における画素値、Ｈは高エネルギーの減弱率の画像における画素値である。

［数１１］式の非線形連立方程式を解くことで、実効原子番号Ｚと面密度Ｄが求まる。非線形連立方程式を解く代表的な方法として、ここではニュートンラフソン法を用いた場合について説明する。まず、ニュートンラフソン法の反復回数をｍ、ｍ回目の反復後の実効原子番号をＺ^ｍ、ｍ回目の反復後の面密度をＤ^ｍとしたとき、ｍ回目の反復後の高エネルギーの減弱率Ｈ^ｍ、ｍ回目の反復後の低エネルギーの減弱率Ｌ^ｍは以下の［数１２］式で表される。

また、厚みが微小に変化したときの減弱率の変化率を、以下の［数６］式で表す。また、実効原子番号と面密度が微小に変化したときの減弱率の変化率を、以下の[数１３]式で表す。

このとき、ｍ＋１回目の反復後の実効原子番号Ｚ^ｍ＋１と面密度Ｄ^ｍ＋１を、高エネルギーの減弱率Ｈと低エネルギーの減弱率Ｌを用いて、以下の［数１４］式で表す。

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

従って、［数１４］式に［数１５］式を代入すると、以下の[数１６]式が求まる。

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

本実施形態では、ニュートンラフソン法を用いて非線形連立方程式を解くことを説明したが、このような実施形態に限定されない。例えば、最小二乗法や二分法などの反復解法を用いてもよい。

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

なお、テーブルは、骨の厚みＢや軟部組織の厚みＳに関するものの他、実効原子番号Ｚや面密度Ｄに関するものでも同様である。例えば、様々な組み合わせの高エネルギーの減弱率Ｈと低エネルギーの減弱率Ｌに対する実効原子番号Ｚや面密度Ｄを事前に求めてテーブルを生成し、このテーブルを参照することで実効原子番号Ｚや面密度Ｄを高速に求めるようにしてもよい。

図８に、本実施形態に係るテーブル生成の流れを説明するフローチャートを示す。テーブルの分割数をＭ、テーブルの低エネルギーの座標をｌ、テーブルの高エネルギーの座標をｈとする。なお、テーブルの座標ｌ及びｈは整数である。まず、撮像制御装置１０３は、低エネルギーの座標をｌ＝０として初期化してから（Ｓ８０１）、高エネルギーの座標をｈ＝０として初期化する（Ｓ８０２）。

次に、撮像制御装置１０３は、座標ｌにおける低エネルギーの減弱率Ｌ［ｌ］と座標ｈにおける高エネルギーの減弱率Ｈ［ｈ］を以下の[数１７]式により求める。

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

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

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

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

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

［数２０］式によって座標から減弱率を求めた場合は、以下の[数２１]式を用いて座標を求めることができる。

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

また、本実施形態において、テーブルの生成と参照を行うとき、解が存在しない高エネルギーの減弱率Ｈと低エネルギーの減弱率Ｌの組み合わせが生じうる。例えば、高エネルギーの減弱率Ｈは、通常、低エネルギーの減弱率Ｌよりも大きくなる。従って、［数１７］式や［数２０］式を用いて生成されたテーブルでは、Ｈ＜Ｌとなる領域の解が求められない。

本実施形態においては、テーブル上で解が求められない領域を減らすように、座標を選択してもよい。例えば、縦軸の座標をｌｎ（Ｌ）／ｌｎ（Ｈ）、横軸の座標を－ｌｎ（Ｈ）などにしてもよい。あるいは、縦軸の座標として、単色放射線で近似して解いた値を用いるなどしてもよい。また、本実施形態においてテーブルの参照を行う際に、テーブルの範囲外の座標が指定される可能性や、解が求められない領域が参照される可能性がある。このような場合は、指定された座標の周囲にある、解が存在する領域の値を用いる構成が用いられる。

図１０に、本実施形態に係るエネルギーサブトラクション処理の画像処理のブロック図を示す。本実施形態における画像処理では、例えば、図６に示した信号処理によって得られた骨画像Ｂや軟部組織画像Ｓに階調処理を施して、表示する構成が用いられる。骨画像Ｂからは軟部組織が除去されているため、エネルギーサブトラクション処理を行っていない画像と比べて、医療用デバイスや石灰化などの視認性が向上する可能性がある。

また、骨画像Ｂや軟部組織画像Ｓに対して、リカーシブフィルタ等の時間方向のフィルタや、ガウシアンフィルタ等の空間方向のフィルタを適用し、ノイズ低減を行ってから階調処理を施して表示する構成も用いられる。しかしながら、本実施形態はそのような形態に限定されない。例えば、図１０に示すように、骨画像Ｂと軟部組織画像Ｓから、合計の厚み画像Ｔを得る構成が好適に用いられる。

人体において、骨は軟部組織に包まれている。すなわち、骨がある部分の軟部組織の厚みは、骨の厚みの分だけ減少している。従って、軟部組織画像Ｓに階調処理を施して表示すると、骨のコントラストが現れてしまう。一方、軟部組織画像Ｓと骨画像Ｂを加算すると、軟部組織画像における骨の厚みの減少が埋め戻される。すなわち、厚み画像Ｔは骨のコントラストが除去された画像となる。このような画像は、下肢におけるＩＶＲのように、骨のコントラストが造影剤の視認性を低下させる状況で用いられる。当然のことながら、厚み画像Ｔに対してリカーシブフィルタ等の時間方向のフィルタや、ガウシアンフィルタ等の空間方向のフィルタを適用し、ノイズ低減を行ってから階調処理を施して表示する構成としてもよい。また、造影剤のＫ吸収端を利用できるよう、例えば、信号処理部１３３の信号処理において、低エネルギーのＸ線の管電圧を４０ｋＶ～５０ｋＶに設定する構成が好適に用いられる。

低エネルギーに基づいた放射線画像を取得するための放射線のスペクトルの平均エネルギーが、造影剤のＫ吸収端よりも低いエネルギーとなるように線質を選択することで、厚み画像Ｔにおける造影剤のコントラストを強調することができる。

図１１に、本実施形態に係る画像処理の流れを説明するフローチャートを示す。まず、取得部１３１は、造影剤を注入する前に図３又は図４に示したタイミングで撮影を行い、補正部１３２は、図５に示した補正を行うことで、高エネルギーのマスク画像Ｈ_Ｍと、低エネルギーのマスク画像Ｌ_Ｍを取得する（Ｓ１１０１）。次に、造影剤が注入された後（Ｓ１１０２）、取得部１３１及び補正部１３２は、同様に撮影と補正を行い、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌを取得する（Ｓ１１０３）。

その後、信号処理部１３３は、高エネルギーのマスク画像Ｈ_Ｍと低エネルギーのマスク画像Ｌ_Ｍに対して図６に示した信号処理（エネルギーサブトラクション（ＥＳ）処理）を行い、骨のマスク画像Ｂ_Ｍと軟部組織のマスク画像Ｓ_Ｍに分離する。さらに、信号処理部１３３は、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌに対しても同様の処理を行い、骨のライブ画像Ｂ_Ｌと軟部組織のライブ画像Ｓ_Ｌに分離する（Ｓ１１０４）。

最後に、画像処理部１３４は、骨のライブ画像Ｂ_Ｌから骨のマスク画像Ｂ_Ｍを引くことで骨のＤＳＡ画像Ｂ_ＤＳＡを取得し、軟部組織のライブ画像Ｓ_Ｌから軟部組織のマスク画像Ｓ_Ｍを引くことで軟部組織のＤＳＡ画像Ｓ_ＤＳＡを取得する。このようなライブ画像の撮影と画像処理を繰り返し行うことで、エネルギーサブトラクションとＤＳＡを組み合わせた動画像を得る。処理の終了が指示されるまでの間（Ｓ１１０６でＮＯ）、上述した処理を繰り返し行うことで、エネルギーサブトラクションとＤＳＡを組み合わせた動画像を得ることができる。終了が指示されると（Ｓ１１０６でＹＥＳ）、本処理は終了する。

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

一方、高エネルギーのライブ画像Ｈ_Ｌと低エネルギーのライブ画像Ｌ_Ｌは、造影剤を注入した後に撮影した画像であるため、骨と軟部組織と造影剤の情報を含んでいる。これらのマスク画像に対してエネルギーサブトラクションを行うと、三種類目の物質である造影剤の情報が骨のライブ画像Ｂ_Ｌに現れる。すなわち、骨のライブ画像Ｂ_Ｌには骨と造影剤の情報が、軟部組織のライブ画像Ｓ_Ｌ（不図示）には軟部組織の情報が含まれる。

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

以上に示したように、本実施形態におけるエネルギーサブトラクション処理を行うことで、軟部組織が除去された骨画像や、骨が除去された厚み画像、軟部組織と骨が除去された造影剤画像を得ることが可能である。いずれの場合も、不要な構造が除去されることで、視認性の向上が期待される。

本願発明者が更に検討を行ったところ、図６の信号処理を行った骨画像から軟部組織が除去しきれない場合があることが判明した。また、図１０の画像処理を行った厚み画像から骨が除去しきれない場合や、図１２の画像処理を行った造影剤画像から軟部組織が除去しきれない場合があることも判明した。不要な構造が除去しきれない、すなわち消し残りがあると、視認性が低減するという課題が生じ得る。

本願発明者が検討を行った結果、散乱線や画素２０の非直線性によって、[数３]式の左辺である減弱率Ｉ／Ｉ_０に誤差が生じた場合、消し残りが生じうることが判明した。また、[数４]式の右辺である、低エネルギーのＸ線におけるスペクトルＮ_Ｌ（Ｅ）、高エネルギーのＸ線におけるスペクトルＮ_Ｈ（Ｅ）が実際のスペクトルと異なる場合にも、消し残りが生じうることが判明した。

同様に、骨の線減弱係数μ_Ｂ（Ｅ）、軟部組織の線減弱係数μ_Ｓ（Ｅ）が、実際の線減弱係数と異なる場合にも、消し残りが生じうる。説明を簡略化するため、散乱線、画素２０の非直線性、スペクトルの誤差は除去または補正済であるものとする。本実施形態においては、線減弱係数の誤差を低減する方法を提案する。

図１３に、実際の水の厚みと計算された骨の厚みの関係を示す。ここでは、被写体の構成として水を設定する。図３又は図４に示したタイミングで取得部１３１は撮影を行い、補正部１３２が図５に示した補正を行い、信号処理部１３３が図６に示した信号処理を行うと、骨画像Ｂ及び軟部組織画像Ｓが得られる。被写体の構成は水のみなので、骨画像Ｂの出力、すなわち骨の厚みはゼロであることが望ましい。[数４]式に入力する軟部組織の線減弱係数μ_Ｓ（Ｅ）として水の線減弱係数を用いた場合、図１３の破線１３０１に示すように、実際の水の厚みが変化しても、計算された骨の厚みは変化しない。一方、軟部組織の線減弱係数μ_Ｓ（Ｅ）として脂肪の線減弱係数を用いた場合、図１３の実線１３０２に示すように、実際の水の厚みが変化すると、計算された骨の厚みが変化してしまう。図１４は実施形態に係る画像処理の模式図であり、骨画像１４０１は軟部組織として脂肪の線減弱係数を用いた場合を示し、骨画像１４０２は軟部組織として水の線減弱係数を用いた場合を示している。この場合、図１４の実施形態に係る画像処理の模式図に示すように、骨画像１４０１では、軟部組織の厚みの変化が骨画像に現れて、消し残りが生じてしまう。

以上より本願発明者は、実際の軟部組織の線減弱係数と、図６に示した信号処理に入力する線減弱係数が異なると、消し残りが生じることを見出した。また、実際の軟部組織の線減弱係数を推定し、図６に示した信号処理に入力することで、消し残りを低減できることを見出した。線減弱係数を推定する方法としては、実効原子番号が好適に用いられる。図７の信号処理で説明したとおり、実効原子番号Ｚとは混合物の等価的な原子番号のことであり、面密度Ｄ［ｇ／ｃｍ^２］とは被写体の密度［ｇ／ｃｍ^３］と被写体の厚み［ｃｍ］の積である。すなわち実効原子番号は厚みによらないパラメータである。水の実効原子番号は７．４２、脂肪は６．３３、骨は１３．１８、アクリルは６．４７である。また、本願発明者が生物の組織片を測定したところ、臓器によって実効原子番号が異なることが判明した。すなわち、実際の軟部組織の線減弱係数は、撮影する部位によって異なると言える。従って、撮影中に軟部組織の実効原子番号を測定し、線減弱係数を推定する構成が好適に用いられる。

図１５に、本実施形態に係るエネルギーサブトラクション処理の信号処理のブロック図を示す。信号処理部１３３は、複数の減弱率の画像（Ｈ、Ｌ）を用いたエネルギーサブトラクション処理により実効原子番号Ｚの画像と面密度Ｄの画像を生成する。まず、信号処理ブロック１５１において、信号処理部１３３は、図５で示した補正処理で得た高エネルギーにおける減弱率の画像Ｈと低エネルギーにおける減弱率の画像Ｌに対して、図７で示した信号処理を行い、実効原子番号Ｚの画像と面密度Ｄの画像を求める。

信号処理部１３３は、実効原子番号の画像における画素値から第１の物質（例えば、軟部組織）の減弱情報を推定する。信号処理部１３３は信号処理ブロック１５２において、信号処理ブロック１５１で得た実効原子番号Ｚの画像における複数画素の画素値を平均化して、軟部組織の実効原子番号Ｚ_Ｓを推定する。次に、信号処理部１３３は信号処理ブロック１５３において、信号処理ブロック１５２で推定した軟部組織の実効原子番号Ｚ_Ｓに基づいて、軟部組織の線減弱係数μ_Ｓ（Ｅ）を推定する。

信号処理部１３３は、複数の減弱率の画像（Ｈ，Ｌ）と、減弱情報とを用いたエネルギーサブトラクション処理により、第１の物質の厚みの画像と第２の物質の厚みの画像とに分離した画像を生成する。

ここで、信号処理部１３３は、減弱情報として、放射線のエネルギーにおける第１の物質の線減弱係数を推定すると、信号処理部１３３は、複数の減弱率の画像（Ｈ，Ｌ）と、減弱情報（第１の物質の線減弱係数）とを用いたエネルギーサブトラクション処理により、第１の物質の厚みの画像（軟部組織画像Ｓ）と第２の物質の厚みの画像（骨画像Ｂ）とに分離した画像を生成する。

すなわち、信号処理部１３３は、信号処理ブロック１５４において、信号処理ブロック１５３で推定した軟部組織の線減弱係数μ_Ｓ（Ｅ）を用いて、図６で示した信号処理を行い、骨画像Ｂと軟部組織画像Ｓを取得する。

このような信号処理を行うことで、実際の軟部組織の線減弱係数と信号処理において入力する軟部組織の線減弱係数μ_Ｓ（Ｅ）との誤差が小さくなり、骨画像における軟部組織の消し残りを低減することができる。

また、信号処理部１３３は、減弱情報として、放射線のエネルギーにおける第１の物質の質量減弱係数を推定すると、信号処理部１３３は、複数の減弱率の画像（Ｈ，Ｌ）と、減弱情報（第１の物質の質量減弱係数）とを用いたエネルギーサブトラクション処理により、第１の物質の面密度の画像と第２の物質の面密度の画像を生成する。

図１５の信号処理ブロック１５１では、図７で示した信号処理を行い、実効原子番号Ｚの画像と面密度Ｄの画像を取得する例を説明したが、本実施形態は、このような形態に限定されない。例えば、図８で示した方法で実効原子番号と面密度のテーブルを生成し、図９で示した方法でテーブルを参照することで、図７と同等の信号処理を高速に行う構成を用いることも可能である。

また、図１５の信号処理ブロック１５２では、信号処理部１３３は、信号処理ブロック１５１で得た実効原子番号Ｚの画像における複数画素の画素値を平均化して、軟部組織の実効原子番号Ｚ_Ｓを推定する。信号処理部１３３は信号処理ブロック１５３において、信号処理ブロック１５２で推定した軟部組織の実効原子番号Ｚ_Ｓに基づいて、第１の物質の減弱情報（軟部組織の線減弱係数μ_Ｓ（Ｅ））を推定する。そして、信号処理部１３３は、複数の減弱率の画像（Ｈ、Ｌ）と、推定した第１の物質の減弱情報とを用いたエネルギーサブトラクション処理により、第１の物質の画像と第２の物質の画像とに分離した画像を生成する。

ここで、複数画素の画素値の具体例として、実効原子番号の画像Ｚの全ての画素の画素値を平均化する構成が考えられる。しかしながら本実施形態はこのような形態に限定されない。例えば、骨がある領域における画素の画素値を平均化の計算に含めると、第１の物質（例えば、軟部組織）の実効原子番号Ｚ_Ｓの推定精度が低下するという課題が生じ得る。そのため、信号処理部１３３は、実効原子番号が予め定めた閾値（実効原子番号の閾値）を上回る画素を、実効原子番号の画像Ｚにおいて第２の物質（例えば、骨）が存在する領域の画素と判定し、第２の物質が存在する領域の画素と判定された画素値を平均化の計算に含めないように処理することも可能である。

また、被写体が存在しない領域における画素の画素値を平均化の計算に含めた場合に、第１の物質（例えば、軟部組織）の実効原子番号Ｚ_Ｓの推定精度が低下するという課題も生じ得る。そのため、信号処理部１３３は、面密度が予め定めた閾値（面密度の閾値（面密度閾値）を下回る画素を、被写体が存在しない領域における画素と判定し、被写体が存在しない領域の画素と判定された画素値を平均化の計算に含めないように処理することも可能である。

さらには、撮影する範囲が広い場合、異なる臓器が存在するため、画面内の領域ごとに軟部組織の実効原子番号Ｚ_Ｓが異なるという課題も生じ得る。そのため、信号処理部１３３は、実効原子番号の画像（画面全体）を複数の領域に分割し、分割した領域ごとに実効原子番号の画像の画素値を平均化して、各領域における第１の物質（例えば、軟部組織）の実効原子番号Ｚ_Ｓを推定することも可能である。そして、信号処理部１３３は推定した第１の物質の実効原子番号Ｚ_Ｓに基づいて、第１の物質の減弱情報（例えば、線減弱係数μ_Ｓ（Ｅ））を推定することも可能である。
この場合において、分割した領域ごとに実効原子番号画像Ｚの画素値を平均化する際に、骨がある領域と被写体が存在しない領域の値を平均化の計算に含めないように処理することも可能である。

図１５の信号処理ブロック１５３では、信号処理部１３３は、信号処理ブロック１５２で推定した軟部組織の実効原子番号Ｚ_Ｓに基づいて、軟部組織の線減弱係数μ_Ｓ（Ｅ）を求めている。例えば、軟部組織を水と脂肪の混合物と仮定し、水の実効原子番号をＺ_Ｗ、水の面密度をＤ_Ｗ、脂肪の実効原子番号をＺ_Ａ、脂肪の面密度をＤ_Ａとすると、軟部組織の面密度Ｄ_Ｓ及び実効原子番号Ｚ_Ｓは[数２２]式で表される。

また、軟部組織に含まれる水の比率αは、水の密度をρ_Ｗ、脂肪の密度をρ_Ａとすると、[数２３]式で表される。

[数２２]式と[数２３]式より、以下の[数２４]式が導かれる。

水の密度ρ_Ｗ、と脂肪の密度ρ_Ａは、ＮＩＳＴなどのデータベースから得られる。より精度を高めるために、実測した水の密度ρ_Ｗ、と脂肪の密度ρ_Ａを使用する構成が好適に用いられる。

なお、水の密度ρ_Ｗ、と脂肪の密度ρ_Ａが等しいと近似した場合、水の比率αは[数２５]で表される。

以上より、軟部組織の線減弱係数μ_Ｓ（Ｅ）は、以下の[数２６]式で表される。ここで、μ_Ｗ（Ｅ）は水の線減弱係数であり、μ_Ａ（Ｅ）は脂肪の線減弱係数である。例えば、軟部組織を水と脂肪の混合物と仮定すると、[数２６]の上式において、左辺の軟部組織は、右辺において、水と脂肪の関係式の和として示される。[数２６]の下式は軟部組織の線減弱係数μ_Ｓ（Ｅ）を、[数２３]式の水の比率αと、水の線減弱係数μ_Ｗ（Ｅ）と、脂肪の線減弱係数μ_Ａ（Ｅ）と、を用いて表した式である。

[数２６]式の水の線減弱係数μ_Ｗ（Ｅ）と脂肪の線減弱係数μ_Ａ（Ｅ）は、ＮＩＳＴなどのデータベースから得た値が用いられる。しかしながら、より精度を高めるために、実測した線減弱係数を用いてもよい。信号処理部１３３は、減弱情報として、放射線のエネルギーにおける第１の物質（例えば、軟部組織）の線減弱係数を推定する。ここで、信号処理部１３３は、第３の物質（例えば、水）の線減弱係数と、第３の物質の密度と、第４の物質（例えば、脂肪）の線減弱係数と、第４の物質の密度と、を用いて、第１の物質の線減弱係数を推定する。信号処理部１３３は、第３の物質の線減弱係数と第４の物質の線減弱係数、又は、第３の物質の密度と第４の物質の密度として、実測値を用いて第１の物質（例えば、軟部組織）の線減弱係数を推定することが可能である。

図１５の信号処理ブロック１５４では、信号処理部１３３は、信号処理ブロック１５３で推定した軟部組織の線減弱係数μ_Ｓ（Ｅ）を用いて、図６で示した信号処理を行い、骨画像Ｂと軟部組織画像Ｓを取得する例を説明したが、本実施形態はこのような形態に限定されない。例えば、図８で示した方法で軟部組織の厚みと骨の厚みのテーブルを生成し、図９で示した方法でテーブルを参照することで、図６と同等の信号処理を高速に行う構成を用いることも可能である。

ただし、軟部組織に含まれる水の比率αが変化すると、軟部組織の線減弱係数μ_Ｓ（Ｅ）も変化するため、軟部組織の厚みと骨の厚みのテーブルの再生成が必要になるという課題が生じ得る。

図１６に、本実施形態に係るエネルギーサブトラクション処理の信号処理のブロック図を示す。まず図８で示した方法で、水の厚み（Ｗ）と骨の厚み（Ｂ_Ｗ）のテーブルと、脂肪の厚み（Ａ）と骨の厚み（Ｂ_Ａ）のテーブルを生成する。信号処理部１３３は、複数の減弱率の組合せに応じて第２の物質の厚み（例えば、骨の厚み（Ｂ_Ｗ））と第３の物質の厚み（例えば、水の厚み（Ｗ））が登録された第１のテーブルと、複数の減弱率の組合せに応じて第２の物質の厚み（例えば、骨の厚み（Ｂ_Ａ））と第４の物質の厚み（例えば、脂肪の厚み（Ａ））が登録された第２のテーブルと、を生成する。

そして、信号処理ブロック１６１において、信号処理部１３３は、水の厚みと骨の厚みのテーブルを用いて、水の厚みＷと骨の厚みＢ_Ｗを求める。

さらに信号処理ブロック１６２において、信号処理部１３３は、脂肪の厚みと骨の厚みのテーブルを用いて、脂肪の厚みＡと骨の厚みＢ_Ａを求める。

次に信号処理ブロック１６３において、信号処理部１３３は、信号処理ブロック１６１の処理で求めた骨の厚みＢ_Ｗ、信号処理ブロック１６２の処理で求めた骨の厚みＢ_Ａ、軟部組織に含まれる水の比率αを用いて、骨の厚みＢを求める。すなわち、信号処理部１３３は、信号処理ブロック１６１で使用した第１のテーブルを用いて取得した第２の物質の厚み（骨の厚みＢ_Ｗ）と、信号処理ブロック１６２で使用した第２のテーブルを用いて取得した第２の物質の厚み（骨の厚みＢ_Ａ）と、比率情報（水の比率α）とに基づいて、第２の物質の厚み画像（骨の厚みＢ）を取得する。

最後に信号処理ブロック１６４において、信号処理部１３３は、信号処理ブロック１６１の処理で求めた水の厚みＷ、信号処理ブロック１６２の処理で求めた脂肪の厚みＡ、軟部組織に含まれる水の比率αを用いて、軟部組織の厚みＳを求める。すなわち、信号処理部１３３は、信号処理ブロック１６１で使用した第１のテーブルを用いて取得した第３の物質の厚み（水の厚み（Ｗ））と、信号処理ブロック１６２で使用した第２のテーブルを用いて取得した第４の物質の厚み（脂肪の厚み（Ａ））と、比率情報（水の比率α）とに基づいて、第１の物質の厚み画像（軟部組織の厚みＳ）を取得する。信号処理ブロック１６３で骨の厚みＢを求める方法と、信号処理ブロック１６４で軟部組織の厚みＳを求める方法としては、[数２７]式が用いられる。

以上のようにして、軟部組織に含まれる水の比率αが変化しても、テーブルを再生成することなく、骨画像Ｂと軟部組織画像Ｓを求めることができる。

図１６の説明では、水の厚み（Ｗ）と骨の厚み（Ｂ_Ｗ）のテーブルと、脂肪の厚み（Ａ）と骨の厚み（Ｂ_Ａ）のテーブルを用いていた。しかしながら本実施形態はこのような形態に限定されない。例えば、水と脂肪以外の物質を用いてテーブルを生成してもよい。また、水と脂肪の比率を変化させて２個以上のテーブルを作成し、それら２個以上のテーブルと軟部組織に含まれる水の比率αから、骨の厚みＢや軟部組織の厚みＳを求めてもよい。また、被写体の軟部組織が、水と脂肪の混合物以外のもので構成されていてもよい。

図１７に、実際の筋肉の厚みと計算された骨の厚みの関係を示す。ここでは、被写体の構成として筋肉を設定する。図３又は図４に示したタイミングで取得部１３１は撮影を行い、補正部１３２が図５に示した補正を行い、信号処理部１３３が図６に示した信号処理を行うと、骨画像Ｂ及び軟部組織画像Ｓが得られる。被写体の構成は筋肉のみなので、骨画像Ｂの出力、すなわち骨の厚みはゼロであることが望ましい。[数４]式に入力する軟部組織の線減弱係数μ_Ｓ（Ｅ）として脂肪の線減弱係数を用いた場合、図１７の実線１７０１に示すように、実際の筋肉の厚みが変化すると、計算された骨の厚みが変化してしまう。また、軟部組織の線減弱係数μ_Ｓ（Ｅ）として水の線減弱係数を用いた場合も、図１７の一点鎖線１７０２に示すように、実際の筋肉の厚みが変化すると、計算された骨の厚みが変化してしまう。筋肉の実効原子番号は７．４５であり、水の実効原子番号である７．４２よりも高いためである。このような状況で、図１６で示した方法で軟部組織に含まれる水の比率αを計算するとα＞１となり、軟部組織に含まれる脂肪の比率は（１－α）＜０となってしまう。これは物理的には起こりえない比率である。

しかしながら本願発明者が検討を行ったところ、[数４]式に入力する軟部組織の線減弱係数μ_Ｓ（Ｅ）として、軟部組織に含まれる水の比率α（α＞１）から計算した線減弱係数を用いることが可能であることが判明した。また、そのような線減弱係数を用いると、図１７の破線１７０３に示すように、実際の筋肉の厚みが変化しても、計算された骨の厚みが変化しないことが判明した。すなわち、軟部組織の実効原子番号が水以上の場合や脂肪以下の場合であっても、軟部組織に含まれる水の比率αを推定することで、消し残りを低減することが可能であることが判明した。要するに、軟部組織に含まれる水の比率αや、軟部組織に含まれる脂肪の比率（１－α）は、１以上や０以下となってもよい。当然のことながら、被写体の軟部組織に含まれる物質として、水と脂肪以外の物質である物質１と物質２の線減弱係数を用いてテーブルの生成や比率αの計算を行ったときも、物質１の比率αや物質２の比率（１－α）は、１以上や０以下となってよい。

また、図１５で説明した信号処理ブロック１５４では、推定された軟部組織の線減弱係数μ_Ｓ（Ｅ）を用いて、軟部組織が除去された骨画像や、骨が除去された厚み画像、軟部組織と骨が除去された造影剤画像における不要な構造の消し残りを低減していた。

しかしながら本実施形態はこのような形態に限定されない。例えば、信号処理ブロック１５４で線減弱係数（単位は１／ｃｍ）を用いた場合、出力される画像の単位はｃｍ、すなわち長さである。これに対して質量減弱係数（単位はｃｍ^２／ｇ）を用いる構成としてもよい。質量減弱係数を用いる場合には、線減弱係数を密度で割る構成が好適に用いられる。この場合、出力される画像の単位はｇ／ｃｍ^２、すなわち面密度となる。

また、本実施形態の信号処理ブロック１５４では、軟部組織画像Ｓと骨画像Ｂを出力していた。しかしながら本実施形態はこのような形態に限定されない。任意の二物質の画像を出力してもよい。以上より本実施形態では、互い異なるエネルギーで取得した複数の画像から、第１の物質の減弱係数を推定し、推定された第１の物質の減弱係数と、予め定めた第２の物質の減弱係数を用いて、物質１の画像又は物質２の画像を生成することも可能である。

（第２実施形態）
第１実施形態では、軟部組織の線減弱係数を推定することで、不要な構造の消し残りを低減していた。この際、図７で示した信号処理を用いて計算した軟部組織の実効原子番号から、軟部組織の線減弱係数を推定していた。これに対して第２実施形態では、図６で示した信号処理を用いて計算した水の厚みと脂肪の厚みから、軟部組織の線減弱係数を推定する構成について説明する。

図１８に、第２実施形態に係る信号処理のブロック図を示す。信号処理部１３３は、複数の減弱率の画像（Ｈ、Ｌ）を用いたエネルギーサブトラクション処理により第１の物質（例えば、軟部組織）に含まれる第３の物質（例えば、軟部組織に含まれる水）の厚みの画像と第４の物質（例えば、軟部組織に含まれる脂肪）の厚みの画像とに分離した画像を生成する。まず、信号処理ブロック１８１において、信号処理部１３３は、図５で示した補正処理で得た高エネルギーにおける減弱率の画像Ｈと低エネルギーにおける減弱率の画像Ｌに対して、図６で示した信号処理を行う。ただし、骨の線減弱係数μ_Ｂ（Ｅ）と軟部組織の線減弱係数μ_Ｓ（Ｅ）ではなく、水の線減弱係数μ_Ｗ（Ｅ）と脂肪の線減弱係数μ_Ａ（Ｅ）を用いることで、第３の物質の厚みの画像（水画像Ｗ）と第４の物質の厚みの画像（脂肪画像Ａ）を求める。

信号処理部１３３は、第３の物質の厚みの画像と第４の物質の厚みの画像とから取得した第１の物質に含まれる第３の物質の比率情報から第１の物質の減弱情報を推定する。まず、信号処理ブロック１８２において、信号処理部１３３は、信号処理ブロック１８１で得た第３の物質の厚みの画像（水画像Ｗ）と第４の物質の厚みの画像（脂肪画像Ａ）から、第１の物質（例えば、軟部組織）に含まれる第３の物質の比率情報（水の比率α）を推定する。

次に信号処理ブロック１８３において、信号処理部１３３は、信号処理ブロック１８２で推定した第１の物質（例えば、軟部組織）に含まれる第３の物質の比率情報（水の比率α）から、第１の物質の減弱情報（例えば、線減弱係数μ_Ｓ（Ｅ））を推定する。

信号処理部１３３は、複数の減弱率の画像（Ｈ、Ｌ）と、第１の物質の減弱情報とを用いたエネルギーサブトラクション処理により、第１の物質の画像と第２の物質の画像とに分離した画像を生成する。

ここで、信号処理部１３３は、減弱情報として、放射線のエネルギーにおける第１の物質の線減弱係数を推定すると、信号処理部１３３は、複数の減弱率の画像（Ｈ，Ｌ）と、減弱情報（第１の物質の線減弱係数）とを用いたエネルギーサブトラクション処理により、第１の物質の厚みの画像（軟部組織画像Ｓ）と第２の物質の厚みの画像（骨画像Ｂ）とに分離した画像を生成する。

すなわち、信号処理ブロック１８４において、信号処理部１３３は、信号処理ブロック１８３で推定した軟部組織の線減弱係数μ_Ｓ（Ｅ）を用いて、図６で示した信号処理を行い、骨画像Ｂと軟部組織画像Ｓを取得する。

このような信号処理を行うことで、実際の軟部組織の線減弱係数と信号処理に入力する軟部組織の線減弱係数μ_Ｓ（Ｅ）の誤差が小さくなり、骨画像における軟部組織の消し残りを低減することができる。

また、信号処理部１３３は、減弱情報として、放射線のエネルギーにおける第１の物質の質量減弱係数を推定すると、信号処理部１３３は、複数の減弱率の画像（Ｈ，Ｌ）と、減弱情報（第１の物質の質量減弱係数）とを用いたエネルギーサブトラクション処理により、第１の物質（例えば、軟部組織）の面密度の画像と第２の物質（例えば、骨）の面密度の画像を生成する。

図１８の信号処理ブロック１８１では、水の線減弱係数μ_Ｗ（Ｅ）と脂肪の線減弱係数μ_Ａ（Ｅ）を用いて図６で示した信号処理を行い、水画像Ｗと脂肪画像Ａを取得する例を説明したが、本実施形態は、このような形態に限定されない。例えば、図８で示した方法で水と脂肪のテーブルを生成し、図９で示した方法でテーブルを参照することで、図６と同等の信号処理を高速に行う構成を用いることも可能である。

また、信号処理部１３３は、第３の物質の厚みの画像（例えば、水画像Ｗ）と第４の物質の厚みの画像（例えば、脂肪画像Ａ）の複数の画素の画素値を平均化した比率情報（例えば、第１の物質に含まれる第３の物質の比率α＝Ｗ／（Ｗ＋Ａ））から第１の物質の減弱情報（例えば、軟部組織の線減弱係数μ_Ｓ（Ｅ））を推定する。

図１８の信号処理ブロック１８２では、信号処理部１３３は、信号処理ブロック１８１で得た水画像Ｗと脂肪画像Ａから、画素毎に、水の比率α（Ｗ／（Ｗ＋Ａ））を計算する。さらに、信号処理部１３３は、複数の画素の画素値における水の比率Ｗ／（Ｗ＋Ａ）を平均化して、軟部組織に含まれる水の比率αを推定している。

複数の画素の具体例として、全ての画素の画素値を平均化する構成が考えられる。しかしながら本実施形態はこのような形態に限定されない。例えば、骨がある領域における画素の画素値を平均化の計算に含めると、軟部組織に含まれる水の比率αの推定精度が低下するという課題が生じ得る。

そのため、信号処理部１３３は、比率情報（第１の物質に含まれる第３の物質の比率α）が比率閾値を上回る画素を、第３の物質の厚みの画像（例えば、水画像Ｗ）と第４の物質の厚みの画像（例えば、脂肪画像Ａ）において第２の物質（例えば、骨）が存在する領域の画素と判定し、当該領域の画素と判定された画素値を平均化の計算に含めないように処理する。すなわち、信号処理部１３３は、水の比率Ｗ／（Ｗ＋Ａ）が予め定めた閾値（水の比率閾値）を上回る画素を、骨が存在する領域の画素と判定し、骨が存在する領域の画素と判定された画素の画素値を平均化の計算に含めないように処理することも可能である。

また、被写体が存在しない領域における画素の画素値を平均化の計算に含めた場合に、軟部組織に含まれる水の比率αが低下するという課題も生じ得る。そのため、信号処理部１３３は、第３の物質（例えば、水）と第４の物質（例えば、脂肪）との合計の厚みが厚み閾値を下回る画素を、被写体が存在しない領域における画素と判定し、当該領域の画素と判定された画素値を平均化の計算に含めないように処理する。すなわち、信号処理部１３３は、水画像Ｗにおける水の厚みと脂肪画像Ａにおける脂肪の厚みとの合計の厚み（Ｗ＋Ａ）が予め定めた閾値（厚みの閾値）を下回る画素を、被写体が存在しない領域における画素と判定し、被写体が存在しない領域における画素と判定された画素値を平均化の計算に含めないように処理することも可能である。

さらには、撮影する範囲が広い場合、異なる臓器が存在するため、画面内の領域ごとに軟部組織に含まれる水の比率αが異なるという課題も生じ得る。そのため、信号処理部１３３は、第３の物質の厚みの画像（例えば、水画像Ｗ）と第４の物質の厚みの画像（例えば、脂肪画像Ａ）を複数の領域に分割し、領域ごとに比率情報（第１の物質に含まれる第３の物質の比率α）を平均化して、各領域における第１の物質の減弱情報（例えば、軟部組織の線減弱係数μ_Ｓ（Ｅ））を推定する。

すなわち、信号処理部１３３は、画像（画面全体）を複数の領域に分割し、領域ごとに、水の比率Ｗ／（Ｗ＋Ａ）を平均化して、各領域における軟部組織に含まれる水の比率αを推定することも可能である。この場合において、分割した領域ごとに水の比率Ｗ／（Ｗ＋Ａ）を平均化する際に、骨がある領域と被写体が存在しない領域の値を平均化の計算に含めないように処理することも可能である。

図１８の信号処理ブロック１８３では、信号処理部１３３は、信号処理ブロック１８２で推定した軟部組織に含まれる水の比率αに基づいて、軟部組織の線減弱係数μ_Ｓ（Ｅ）を求めている。このとき、図１５の信号処理ブロック１５３と同様に、信号処理部１３３は、[数２６]式を用いて、軟部組織の線減弱係数μ_Ｓ（Ｅ）を求めることができる。[数２６]式の水の線減弱係数μ_Ｗ（Ｅ）と脂肪の線減弱係数μ_Ａ（Ｅ）は、ＮＩＳＴなどのデータベースから得た値が用いられる。しかしながら、より精度を高めるために、実測した線減弱係数として、第３の物質（例えば、水）の線減弱係数の実測値と、第４の物質（例えば、脂肪）の線減弱係数の実測値とを用いてもよい。実測した第３の物質の密度（例えば、水の密度ρ_Ｗ）と第４の物質（例えば、脂肪の密度ρ_Ａ）とを用いて、データベースから得た線減弱係数を補正することも可能である。

図１８の信号処理ブロック１８４では、信号処理部１３３は、図１５の信号処理ブロック１５４と同様に、信号処理ブロック１８３で推定した軟部組織の線減弱係数μ_Ｓ（Ｅ）を用いて、図６で示した信号処理を行い、骨画像Ｂと軟部組織画像Ｓを取得する。また、図８で示した方法で軟部組織の厚みと骨の厚みのテーブルを生成し、図９で示した方法でテーブルを参照することで、図６と同等の信号処理を高速に行う構成を用いることも可能である。

また、図１６と同様に、水の厚み（Ｗ）と骨の厚み（Ｂ_Ｗ）のテーブルと、脂肪の厚み（Ａ）と骨の厚み（Ｂ_Ａ）のテーブルを生成し、軟部組織に含まれる水の比率αを用いて、骨の厚みＢ及び軟部組織の厚みＳを求めてもよい。また、信号処理部１３３は、複数の減弱率の画像（Ｈ，Ｌ）を用いたエネルギーサブトラクション処理により物質特性画像（物質分離画像、物質識別画像を含む）を生成する。信号処理部１３３は、物質特性画像における画素値（例えば、実効原子番号の画像における画素値）、または、物質特性画像から取得した被写体に含まれる物質の比率（例えば、物質分離画像とから取得した軟部組織に含まれる物質の比率）から、被写体に含まれる分離対象の物質の減弱情報（線減弱係数、または質量減弱係数）を推定する。

そして、信号処理部１３３は、複数の減弱率の画像（Ｈ、Ｌ）と推定した減弱情報とを用いたエネルギーサブトラクション処理により、分離対象の物質（例えば、被写体に含まれる軟部組織や骨）を分離した画像とを生成することも可能である。この場合、信号処理部１３３は、分離対象の物質の減弱情報として線減弱係数を推定した場合、分離対象の物質に関する物質分離画像を生成する。また、信号処理部１３３は、減弱情報として質量弱係数を推定した場合、分離対象の物質に関する物質識別画像を生成する。

また、上記の第１実施形態と第２実施形態では、Ｘ線撮像装置１０４として蛍光体を用いた間接型のＸ線センサを用いた。しかしながら本発明はこのような形態に限定されない。例えばＣｄＴｅ等の直接変換材料を用いた直接型のＸ線センサを用いてもよい。すなわち、Ｘ線センサは、間接型・直接型のどちらでもよい。

また、第１、第２実施形態では、例えば、図４の動作においてＸ線発生装置１０１の管電圧を変化させていた。しかしながら本発明はこのような形態に限定されない。Ｘ線発生装置１０１のフィルタを時間的に切り替えるなどして、Ｘ線撮像装置１０４に曝射されるＸ線のエネルギーを変化させるようにしてもよい。すなわち、Ｘ線撮像装置１０４へ曝射されるＸ線のエネルギーを変更させる方法は、どのような方法は、何等限定されない。

また、第１実施形態及び第２実施形態では、Ｘ線のエネルギーを変化させることで異なるエネルギーの画像を得ていたが、本発明はこのような形態に限定されない。例えば、複数の蛍光体１０５および二次元検出器１０６を重ねることで、Ｘ線の入射方向に対して前面の二次元検出器と背面の二次元検出器から異なるエネルギーの画像を得る、積層型の構成としてもよい。

また、第１実施形態及び第２実施形態では、Ｘ線撮影システムの撮像制御装置１０３を用いてエネルギーサブトラクション処理を行っていた。しながらこの本発明はこのような形態に限定されない。例えば、撮像制御装置１０３で取得した画像を別のコンピュータに転送して、エネルギーサブトラクション処理を行ってもよい。例えば、取得した画像を医療用のＰＡＣＳを介して別のパソコンに転送し、エネルギーサブトラクション処理を行ってから表示する構成としてもよい。すなわち、上記実施形態で説明した補正処理を行う装置は、撮影装置とセットでなくてよい（画像ビューアでもよい）。

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

発明は上記実施形態に制限されるものではなく、発明の精神及び範囲から離脱することなく、様々な変更及び変形が可能である。従って、発明の範囲を公にするために請求項を添付する。

１０１：Ｘ線発生装置、１０２：Ｘ線制御装置、１０３：撮像制御装置、１０４：Ｘ線撮像装置、１０５：蛍光体、１０６：二次元検出器、１３１：取得部、１３２：補正部、１３３：信号処理部、１３４：画像処理部

Claims (22)

1. 被写体に含まれる第１の物質の画像と第２の物質の画像を生成する画像処理装置であって、
   前記被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する取得手段と、
   前記複数の放射線画像からエネルギーサブトラクション処理で用いられる複数の減弱率の画像を生成する生成手段と、
   前記複数の減弱率の画像を用いた前記エネルギーサブトラクション処理により実効原子番号の画像と面密度の画像を生成する処理手段と、を備え、
   前記処理手段は、前記実効原子番号の画像における画素値から前記第１の物質の減弱情報を推定し、
   前記複数の減弱率の画像と前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記第１の物質の画像と前記第２の物質の画像とを生成することを特徴とする画像処理装置。
2. 被写体に含まれる第１の物質の画像と第２の物質の画像を生成する画像処理装置であって、
   前記被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する取得手段と、
   前記複数の放射線画像からエネルギーサブトラクション処理で用いられる複数の減弱率の画像を生成する生成手段と、
   前記複数の減弱率の画像を用いた前記エネルギーサブトラクション処理により前記第１の物質に含まれる第３の物質の厚みの画像と第４の物質の厚みの画像とに分離した画像を生成する処理手段と、を備え、
   前記処理手段は、
   前記第３の物質の厚みの画像と前記第４の物質の厚みの画像とから取得した前記第１の物質に含まれる前記第３の物質の比率情報から前記第１の物質の減弱情報を推定し、
   前記複数の減弱率の画像と前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記第１の物質の画像と前記第２の物質の画像とを生成することを特徴とする画像処理装置。
3. 前記第３の物質は、前記被写体の軟部組織に含まれる水であり、前記第４の物質は、前記軟部組織に含まれる脂肪であることを特徴とする請求項２に記載の画像処理装置。
4. 前記第１の物質は、前記被写体に含まれる軟部組織であり、前記第２の物質は、前記被写体に含まれる骨であることを特徴とする請求項１乃至３のいずれか１項に記載の画像処理装置。
5. 前記処理手段は、
   前記減弱情報として、前記放射線のエネルギーにおける前記第１の物質の線減弱係数を推定し、
   前記複数の減弱率の画像と、前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記第１の物質の厚みの画像と前記第２の物質の厚みの画像とに分離した画像を生成することを特徴とする請求項１乃至４のいずれか１項に記載の画像処理装置。
6. 前記処理手段は、
   前記減弱情報として、前記放射線のエネルギーにおける前記第１の物質の線減弱係数を推定し、
   前記第３の物質の線減弱係数と、前記第３の物質の密度と、前記第４の物質の線減弱係数と、前記第４の物質の密度と、を用いて、前記第１の物質の線減弱係数を推定することを特徴とする請求項２または３に記載の画像処理装置。
7. 前記処理手段は、前記第３の物質の線減弱係数と前記第４の物質の線減弱係数、又は、前記第３の物質の密度と前記第４の物質の密度として、実測値を用いることを特徴とする請求項６に記載の画像処理装置。
8. 前記処理手段は、
   前記減弱情報として、前記放射線のエネルギーにおける前記第１の物質の質量減弱係数を推定し、
   前記複数の減弱率の画像と、前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記第１の物質の面密度の画像と前記第２の物質の面密度の画像を生成することを特徴とする請求項１乃至４のいずれか１項に記載の画像処理装置。
9. 前記処理手段は、
   前記複数の減弱率の組合せに応じて前記第２の物質の厚みと前記第３の物質の厚みが登録された第１のテーブルと、
   前記複数の減弱率の組合せに応じて前記第２の物質の厚みと前記第４の物質の厚みが登録された第２のテーブルと、を有することを特徴とする請求項２または３に記載の画像処理装置。
10. 前記処理手段は、
    前記第１のテーブルを用いて取得した前記第２の物質の厚みと、前記第２のテーブルを用いて取得した前記第２の物質の厚みと、前記比率情報とに基づいて、前記第２の物質の厚み画像を取得し、
    前記第１のテーブルを用いて取得した前記第３の物質の厚みと、前記第２のテーブルを用いて取得した前記第４の物質の厚みと、前記比率情報とに基づいて、前記第１の物質の厚み画像を取得する
    ことを特徴とする請求項９に記載の画像処理装置。
11. 前記処理手段は、前記実効原子番号の画像における画素値を平均化して前記第１の物質の減弱情報を推定することを特徴とする請求項１に記載の画像処理装置。
12. 前記処理手段は、前記実効原子番号の画像を複数の領域に分割し、分割した領域ごとに前記実効原子番号の画像の画素値を平均化して、各領域における前記第１の物質の減弱情報を推定することを特徴とする請求項１１に記載の画像処理装置。
13. 前記処理手段は、
    前記実効原子番号が閾値を上回る画素を、前記実効原子番号の画像において前記第２の物質が存在する領域の画素と判定し、当該領域の画素と判定された画素値を前記平均化の計算に含めないように処理する
    ことを特徴とする請求項１１または１２に記載の画像処理装置。
14. 前記処理手段は、
    前記面密度が面密度閾値を下回る画素を、前記被写体が存在しない領域における画素と判定し、当該領域の画素と判定された前記画素値を平均化の計算に含めないように処理する
    ことを特徴とする請求項１１乃至１３のいずれか１項に記載の画像処理装置。
15. 前記処理手段は、前記第３の物質の厚みの画像と前記第４の物質の厚みの画像の複数の画素の画素値を平均化した前記比率情報から前記第１の物質の減弱情報を推定することを特徴とする請求項２または３に記載の画像処理装置。
16. 前記処理手段は、前記第３の物質の厚みの画像と前記第４の物質の厚みの画像を複数の領域に分割し、領域ごとに前記比率情報を平均化して、各領域における前記第１の物質の減弱情報を推定することを特徴とする請求項１５に記載の画像処理装置。
17. 前記処理手段は、
    前記比率情報が比率閾値を上回る画素を、前記第３の物質の厚みの画像と前記第４の物質の厚みの画像において前記第２の物質が存在する領域の画素と判定し、当該領域の画素と判定された画素値を前記平均化の計算に含めないように処理する
    ことを特徴とする請求項１５または１６に記載の画像処理装置。
18. 前記処理手段は、
    前記第３の物質と前記第４の物質との合計の厚みが厚み閾値を下回る画素を、前記被写体が存在しない領域における画素と判定し、当該領域の画素と判定された前記画素値を平均化の計算に含めないように処理することを特徴とする請求項１５乃至１７のいずれか１項に記載の画像処理装置。
19. 被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する取得手段と、
    前記複数の放射線画像からエネルギーサブトラクション処理で用いられる複数の減弱率の画像を生成する生成手段と、
    前記複数の減弱率の画像を用いた前記エネルギーサブトラクション処理により物質特性画像を生成する処理手段と、を備え、
    前記処理手段は、前記物質特性画像における画素値、または、前記物質特性画像から取得した前記被写体に含まれる物質の比率から、前記被写体に含まれる分離対象の物質の減弱情報を推定し、
    前記複数の減弱率の画像と前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記分離対象の物質を分離した画像とを生成することを特徴とする画像処理装置。
20. 被写体に含まれる第１の物質の画像と第２の物質の画像を生成する画像処理装置における画像処理方法であって、
    前記被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する取得工程と、
    前記複数の放射線画像からエネルギーサブトラクション処理で用いられる複数の減弱率の画像を生成する生成工程と、
    前記複数の減弱率の画像を用いた前記エネルギーサブトラクション処理により実効原子番号の画像と面密度の画像を生成する処理工程と、を有し、
    前記処理工程では、前記実効原子番号の画像における画素値から前記第１の物質の減弱情報を推定し、
    前記複数の減弱率の画像と前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記第１の物質の画像と前記第２の物質の画像とを生成することを特徴とする画像処理方法。
21. 被写体に含まれる第１の物質の画像と第２の物質の画像を生成する画像処理装置における画像処理方法であって、
    前記被写体に放射線を照射して撮影を行うことで得られた互いに異なる複数の放射線エネルギーに対応した複数の放射線画像を取得する取得工程と、
    前記複数の放射線画像からエネルギーサブトラクション処理で用いられる複数の減弱率の画像を生成する生成工程と、
    前記複数の減弱率の画像を用いた前記エネルギーサブトラクション処理により前記第１の物質に含まれる第３の物質の厚みの画像と第４の物質の厚みの画像とに分離した画像を生成する処理工程と、を有し、
    前記処理工程では、
    前記第３の物質の厚みの画像と前記第４の物質の厚みの画像とから取得した前記第１の物質に含まれる前記第３の物質の比率情報から前記第１の物質の減弱情報を推定し、
    前記複数の減弱率の画像と前記減弱情報とを用いた前記エネルギーサブトラクション処理により、前記第１の物質の画像と前記第２の物質の画像とを生成することを特徴とする画像処理方法。
22. コンピュータを、請求項１乃至１９のいずれか１項に記載された画像処理装置の各手段として機能させるためのプログラム。

Images