(57)【特許請求の範囲】
【請求項１】
被検体の所定の撮像部位のサブトラクション像を得るためのディジタルサブトラクション装置であって、（ａ）前記撮像部位にＸ線を照射し、その部位のＸ線透過像を撮像するＸ線透視手段と、（ｂ）前記Ｘ線透過像をディジタルデータに変換するデータ変換手段と、（ｃ）前記ディジタルデータに変換されたＸ線透過像（以下、基本画像という）から所定のしきい値周波数以上の周波数成分を除去した画像（以下、高周波数成分除去画像）を得る周波数特性変換手段と、（ｄ）前記基本画像と高周波数成分除去画像とのサブトラクションを行い、その撮像部位のサブトラクション像を求める演算手段と、（ｅ）前記演算手段への基本画像の供給を遅らせて、前記基本画像と高周波数成分除去画像とを前記演算手段に同期して供給させる遅延手段と、（ｆ）サブトラクション像を表示する表示手段と、（ｇ）与えられた変更指示に従って前記のしきい値周波数を変更するしきい値周波数変更手段とを備え、
前記周波数特性変換手段は、実空間上で直線形状で表されるテンプレートフィルタでもって前記基本画像の高周波数成分を実空間上で除去する実空間フィルタリング部を備え、
前記実空間フィルタリング部のテンプレートフィルタ形状は、互いに直交するＸ，Ｙ，Ｚ軸で表される、実空間上の３次元直交座標系において、Ｘ，Ｙ軸方向を前記テンプレートフィルタの矩形形状の底面大きさとし、Ｚ軸方向を前記テンプレートのゲインとし、Ｚ軸方向に先細りとなる立体形状であって、かつ、前記テンプレートフィルタのＺ方向の少なくとも中央から先細りの先端部にかけて角錐形状で表されるものであり、
前記実空間フィルタリング部は、前記基本画像の行方向および列方向に三角形フィルタをかけることにより、前記テンプレートフィルタを実現することを特徴とするディジタルサブトラクション装置。
【請求項２】
請求項１に記載のディジタルサブトラクション装置において、前記基本画像の高周波数成分を強調した画像（以下、高周波数成分強調画像という）を得る第２の周波数特性変換手段を備え、前記演算手段は、前記基本画像と高周波数成分除去画像とのサブトラクションを行うのではなく、前記高周波数成分強調画像と前記高周波数成分除去画像とのサブトラクションを行うことで、その撮像部位のサブトラクション像を求めることを特徴とするディジタルサブトラクション装置。
【請求項３】
請求項１に記載のディジタルサブトラクション装置において、前記ディジタルデータに変換されたＸ線透過像（以下、基本画像という）を順次取り込み、複数回数分の撮像で得られた基本画像の平均画像を求める積分処理手段を備え、前記周波数特性変換手段は、前記平均画像から所定のしきい値周波数以上の周波数成分を除去した画像（以下、高周波数成分除去画像）を得るようにし、前記演算手段は、前記平均画像と前記高周波数成分除去画像とのサブトラクションを行い、撮像部位のサブトラクション像を求めることを特徴とするディジタルサブトラクション装置。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
この発明は、医療分野、工業分野などに用いられる、被検体の撮像部位をＸ線透視撮像してその撮像部位のサブトラクション像を得るディジタルサブトラクション装置に係り、特に、被検体の撮像部位のサブトラクション像を好適に得る技術に関する。
【０００２】
【従来技術】
従来のディジタルサブトラクション装置としては、例えば、医療分野で用いられる、被検体の所定の撮像部位のサブトラクション像を得るディジタルアンギオグラフィ装置がある。このディジタルアンギオグラフィ装置としては、例えば、被検体の撮像部位を１回Ｘ線透視撮像するだけでその撮像部位のサブトラクション像を得れるものがある。以下に、このディジタルアンギオグラフィ装置により、被検体の所定の撮像部位のサブトラクション像を得る動作について、説明する。
【０００３】
まず、造影剤が投与された被検体の所定の撮像部位を、Ｘ線透視撮像装置（Ｘ線管と、イメージインテンシファイアとテレビカメラ、ＦＰＤ（Ｆｌａｔ  Ｐａｎｅｌ  Ｄｅｔｅｃｔｏｒ）などで構成される撮像系とを備えたもの）でもってＸ線透過像として撮像し、この撮像したＸ線透視像をディジタルデータに変換して基本画像を取得する。この基本画像は、造影剤が投与された血管像などの高周波数成分が残っている画像であり、これをライブ像として用いる。一方、周波数特性変換回路によって、前記基本画像を構成する周波数成分を空間／周波数変換処理で取り出し、所定のしきい値周波数以上の周波数成分を除去し、これに周波数／空間変換処理を施して、前記基本画像から前記しきい値周波数以上の周波数成分を除去した高周波数成分除去画像を得ている。この高周波数成分除去画像は、造影剤が投与された血管像などの高周波数成分が除去された画像であり、これをマスク像として用いる。前記の空間／周波数変換処理としては、ＦＦＴ（高速フーリエ変換）、カルーネン・レーベ変換、ＤＣＴ（離散コサイン変換）、アダマール変換などの各変換方式がある。また、前記の周波数／空間変換処理としては、前記空間／周波数変換処理の逆変換（逆ＦＦＴ、逆カルーネン・レーベ変換、逆ＤＣＴ、逆アダマール変換など）がある。
【０００４】
次に、遅延回路によって、後段の演算器へのライブ像（基本画像）の供給を前記周波数特性変換回路の処理時間分だけ遅延させることで、ライブ像（基本画像）とマスク像（高周波数成分除去画像）とを後段の演算器に同期して供給している。そして、演算器は、ライブ像（基本画像）からマスク像（高周波数成分除去画像）をサブトラクションしてサブトラクション像（サブトラクション画像）を求めている。なお、上述の所定のしきい値周波数は、サブトラクション像に残したい関心物（血管像など）を好適に除去し得る周波数値のことであり、このしきい値周波数を予め理論的あるいは実験的に求めておいて設定することで、サブトラクション像を得るようにしている。
【０００５】
このように、造影剤が投与された被検体の所定の撮像部位を１回だけＸ線透視撮像するだけで、その撮像部位のサブトラクション像を得ることができるので、造影剤の投与前後の２回にわたってＸ線透視撮像する場合に比べて、被検体へのＸ線曝射線量が軽減できるとともに、ライブ像（基本画像）からマスク像（高周波数成分除去画像）を生成していることから、造影剤投与前後の２回の撮像における被検体の体動などに起因するマスク像とライブ像の画像ずれを完全に無くすことができる。
【０００６】
【発明が解決しようとする課題】
しかしながら、このような構成を有する従来例の場合には、次のような問題がある。
すなわち、上述した従来例では、関心物（例えば、血管像など）を好適に除去するための、予め実験的に求めたしきい値周波数を用いて、ライブ像（基本画像）からそのしきい値周波数以上の周波数成分を除去したマスク像（高周波数成分除去画像）を得て、ライブ像（基本画像）からマスク像（高周波数成分除去画像）をサブトラクションして、関心物（例えば、血管像など）を残したサブトラクション像を得るようにしているので、例えば、Ｘ線透視撮像毎に、関心物の空間的形状（大きさや形など）の傾向が異なるような場合において、全てに渡って適正なしきい値を用意しておくことは事実上不可能であり、標準値から外れた関心物、すなわち、設定されたしきい値では十分に除去されない関心物については、好適なサブトラクション像を得ることができないという問題がある。
【０００７】
また、一回のＸ線透視撮像において、空間的形状の異なる複数種類の関心物が存在する場合、つまり、撮像部位内に空間的形状の異なる複数種類の関心物が存在する場合にも、全ての関心物に対して最適となるようなしきい値を求めることは不可能であり、一部の関心物については好適なサブトラクション像を得ることができず、好適に観察できなくなる。例えば、関心物を例えば血管とした場合であっても、動脈や静脈などで臓器に近い部分の太い血管と、動脈や静脈の末梢に位置する毛細血管などの細い血管とでは、それらの大きさの違いからわかるように、上述のしきい値周波数は異なる。
【０００８】
この発明は、このような事情に鑑みてなされたものであって、被検体の撮像部位のサブトラクション像を好適に得ることができるディジタルサブトラクション装置を提供することを目的とする。
【０００９】
【課題を解決するための手段】
この発明は、このような目的を達成するために、次のような構成をとる。
すなわち、請求項１に記載のディジタルサブトラクション装置は、被検体の所定の撮像部位のサブトラクション像を得るためのディジタルサブトラクション装置であって、（ａ）前記撮像部位にＸ線を照射し、その部位のＸ線透過像を撮像するＸ線透視手段と、（ｂ）前記Ｘ線透過像をディジタルデータに変換するデータ変換手段と、（ｃ）前記ディジタルデータに変換されたＸ線透過像（以下、基本画像という）から所定のしきい値周波数以上の周波数成分を除去した画像（以下、高周波数成分除去画像）を得る周波数特性変換手段と、（ｄ）前記基本画像と高周波数成分除去画像とのサブトラクションを行い、その撮像部位のサブトラクション像を求める演算手段と、（ｅ）前記演算手段への基本画像の供給を遅らせて、前記基本画像と高周波数成分除去画像とを前記演算手段に同期して供給させる遅延手段と、（ｆ）サブトラクション像を表示する表示手段と、（ｇ）与えられた変更指示に従って前記のしきい値周波数を変更するしきい値周波数変更手段とを備え、前記周波数特性変換手段は、実空間上で直線形状で表されるテンプレートフィルタでもって前記基本画像の高周波数成分を実空間上で除去する実空間フィルタリング部を備え、前記実空間フィルタリング部のテンプレートフィルタ形状は、互いに直交するＸ，Ｙ，Ｚ軸で表される、実空間上の３次元直交座標系において、Ｘ，Ｙ軸方向を前記テンプレートフィルタの矩形形状の底面大きさとし、Ｚ軸方向を前記テンプレートのゲインとし、Ｚ軸方向に先細りとなる立体形状であって、かつ、前記テンプレートフィルタのＺ方向の少なくとも中央から先細りの先端部にかけて角錐形状で表されるものであり、前記実空間フィルタリング部は、前記基本画像の行方向および列方向に三角形フィルタをかけることにより、前記テンプレートフィルタを実現することを特徴とするものである。
【００１０】
また、請求項２に記載のディジタルサブトラクション装置は、請求項１に記載のディジタルサブトラクション装置において、前記基本画像の高周波数成分を強調した画像（以下、高周波数成分強調画像という）を得る第２の周波数特性変換手段を備え、前記演算手段は、前記基本画像と高周波数成分除去画像とのサブトラクションを行うのではなく、前記高周波数成分強調画像と前記高周波数成分除去画像とのサブトラクションを行うことで、その撮像部位のサブトラクション像を求めることを特徴とするものである。
【００１１】
また、請求項３に記載のディジタルサブトラクション装置は、請求項１に記載のディジタルサブトラクション装置において、前記ディジタルデータに変換されたＸ線透過像（以下、基本画像という）を順次取り込み、複数回数分の撮像で得られた基本画像の平均画像を求める積分処理手段を備え、前記周波数特性変換手段は、前記平均画像から所定のしきい値周波数以上の周波数成分を除去した画像（以下、高周波数成分除去画像）を得るようにし、前記演算手段は、前記平均画像と前記高周波数成分除去画像とのサブトラクションを行い、撮像部位のサブトラクション像を求めることを特徴とするものである。
【００１５】
【作用】
この発明の作用は次の通りである。
すなわち、請求項１に記載の発明によれば、被検体の所定の撮像部位に対して、Ｘ線透視手段でＸ線透過像を撮像し、データ変換手段でディジタルデータに変換して基本画像を得る。この基本画像は、関心物の高周波数成分が残っている画像であり、これをライブ像として用いる。一方、周波数特性変換手段は、前記基本画像から所定のしきい値周波数以上の周波数成分を除去した高周波数成分除去画像を得る。この高周波数成分除去画像は、関心物の高周波数成分が除去された画像であり、これをマスク像として用いる。遅延手段は、後段の演算手段へのライブ像（基本画像）の供給を前記周波数特性変換手段での所定の処理時間分だけ遅延させることで、ライブ像（基本画像）とマスク像（高周波数成分除去画像）とを後段の演算手段に同期して供給している。そして、演算手段は、ライブ像（基本画像）からマスク像（高周波数成分除去画像）をサブトラクションしてサブトラクション像を求める。表示手段は、サブトラクション像を表示する。しきい値周波数変更手段は、与えられた変更指示に従って前記のしきい値周波数を変更する。したがって、しきい値周波数の変更指示により、しきい値周波数変更に係る新たなサブトラクション像が表示手段に表示されるので、適正なしきい値を対話的に設定することができ、表示手段に表示されるサブトラクション像をモニタリングしながら好適なサブトラクション像が得れる。
【００１６】
また、請求項２に記載の発明によれば、請求項１に記載の発明と同様の作用で基本画像と高周波数成分除去画像とが得られる。そして、第２の周波数特性変換手段は、得られた基本画像の高周波数成分を強調した高周波数成分強調画像を得る。演算手段は、この高周波数成分強調画像と高周波数成分除去画像とのサブトラクションを行うことで、その撮像部位の関心物をより強調したサブトラクション像を求める。したがって、関心物をより強調したサブトラクション像を得る場合であっても、しきい値周波数の変更により、好適なサブトラクション像が得られる。
【００１７】
また、請求項３に記載の発明によれば、積分処理手段は、複数回数分の撮像で得られた基本画像の平均画像を求める。この平均画像をライブ像として用い、このライブ像（平均画像）から請求項１に記載の発明と同様の処理で高周波数成分除去画像を得ているので、ライブ像及びマスク像のＳ／Ｎ比が向上する。演算手段は、ライブ像（平均画像）からマスク像（平均画像から得られた高周波数成分除去画像）をサブトラクションして、撮像部位のサブトラクション像を求める。したがって、高画質なサブトラクション像を得る場合であっても、しきい値周波数の変更により、好適なサブトラクション像が得られる。
【００１８】
また、本発明によれば、上述の周波数特性変換手段での高周波数成分除去画像の生成は、実空間フィルタリング部によって行われる。実空間フィルタリング部は、実空間上で直線形状で表されるテンプレートフィルタでもって前記基本画像の高周波数成分を実空間上で除去する。従来のガウス関数形状のテンプレートフィルタでは、フィルタサイズが例えば５１（ＴＡＰ数）×５１（ＴＡＰ数）のように大きくなることでフィルタリングの演算量が膨大になり、リアルタイムに近い短時間での処理は不可能であったが、この請求項１に記載の発明では、実空間上で直線形状で表されるテンプレートフィルタでもってフィルタリングするので、フィルタ移動によって更新される部分についての加減算のみを行なうだけで良く、フィルタサイズに依らずにそのフィルタリングの演算量が一定量に低減され、リアルタイムに近い短時間での処理が可能となり、しきい値周波数の変更により、リアルタイムに好適なサブトラクション像が得られる。
【００２０】
また、本発明によれば、実空間フィルタリング部のテンプレートフィルタ形状は、互いに直交するＸ，Ｙ，Ｚ軸で表される、実空間上の３次元直交座標系において、Ｘ，Ｙ軸方向をテンプレートフィルタの矩形形状の底面大きさとし、Ｚ軸方向を前記テンプレートのゲインとし、Ｚ軸方向に先細りとなる立体形状であって、かつ、テンプレートフィルタのＺ方向の少なくとも中央から先細りの先端部にかけて角錐形状で表されるものとしている。したがって、フィルタ移動によって更新される部分についての加減算のみを行なうだけで良く、フィルタサイズに依らずにそのフィルタリングの演算量が一定量に低減され、リアルタイムに近い短時間での処理が可能となる。また、しきい値周波数の変更により、リアルタイムに好適なサブトラクション像が得られる。さらに、テンプレートフィルタ形状が四角柱形状の場合よりも、高周波数成分を除去する周波数特性に優れる。
【００２１】
【発明の実施の形態】
以下、図面を参照してこの発明のディジタルサブトラクション装置に係る一実施例としてのＸ線ディジタルアンギオグラフィ装置について、図面を参照しながら説明する。
【００２２】
＜第１実施例＞
図１は、この発明の第１実施例に係るＸ線ディジタルアンギオグラフィ装置の全体構成を示す正面図であり、図２は、Ｘ線透視装置を側面から見た図であり、図３は、画像処理部の構成を示すブロック図である。
【００２３】
この第１実施例のＸ線ディジタルアンギオグラフィ装置は、ベッド１、Ｘ線透視手段としてのＸ線透視装置２、画像処理部３、表示手段としてのモニタ４、制御部５、操作盤６などを備えて構成されている。
【００２４】
ベッド１は、床面に設置されたベッド基台１１と天板１２を備えている。被検体Ｍは天板１２上に載置される。この天板１２はモータ１３の駆動で水平移動可能であり、天板１２上の被検体ＭとＸ線透視装置２との相対的な位置関係を被検体Ｍの体軸方向に変位することができる。モータ１３の駆動制御は、制御部５により行われる。
【００２５】
Ｘ線透視装置２は、Ｘ線管２１、撮像系２２を支持するＣ型アーム２３がベッド１の近傍に定置された装置基台２４の上部に支持されて構成されている。Ｃ型アーム２３は、モータ２５の駆動で図２の矢印方向に変位可能に装置基台２４に支持されており、Ｘ線管２１、撮像系２２を被検体Ｍの体軸回りに変位可能に構成し、Ｘ線透過像の撮像方向の調整が可能となっている。モータ２５の駆動制御は制御部５により行われる。
【００２６】
Ｘ線管２１と撮像系２２とはＣ型アーム２３の両端部に取り付けられており、天板１２上の被検体Ｍを挟み込んだ状態で対向配置されている。Ｘ線管２１から被検体Ｍの任意の撮像部位に向けて照射され、被検体Ｍを透過したＸ線は、撮像系２２で受像され、その部位のＸ線透過像が撮像される。Ｘ線管２１からのＸ線の照射は、Ｘ線高電圧発生装置２６から所定の電力（Ｘ線管電圧およびＸ線管電流）がＸ線管２１に供給されて行われる。Ｘ線高電圧発生装置２６からＸ線管２１への所定の電力の供給は制御部５に制御されて行われる。撮像系２２は、イメージインテンシファイアやテレビカメラなどで構成されている。撮像されたＸ線透過像は画像処理部３に与えられる。
【００２７】
画像処理部３は、図３に示すように、データ変換手段としてのＡ／Ｄ（アナログtoディジタルデータ）変換器３１、高周波除去手段としての周波数特性変換回路３２、遅延回路３３、演算手段としての演算器３４、階調変換回路３５、Ｄ／Ａ（ディジタルtoアナログ）変換器３６で構成されている。
【００２８】
造影剤が投与された被検体Ｍの所定の撮像部位のＸ線透過像が撮像されると、撮像系２２からのその像の画像信号（アナログ信号）は、Ａ／Ｄ変換器３１でディジタルデータに変換され基本画像が得られる。この基本画像は、骨格などの低周波数成分や、造影剤が投与された血管像などの高周波数成分を含んだ画像であり、これをライブ像として用いる。この基本画像（ライブ像）は、周波数特性変換回路３２と遅延回路３３とに与えられる。周波数特性変換回路３２は、後述する処理によって基本画像から血管像などの高周波数成分を除去しマスク像を得る。演算器３４では、遅延回路３３を経て供給されるライブ像と、周波数特性変換回路３２を経て供給されるマスク像とのサブトラクションを行いサブトラクション像を求めて階調変換回路３５に与える。階調変換回路３５では、サブトラクション像をモニタ５に表示したとき見やすい画像にするために、サブトラクション像を構成する各画素の濃度を調整（全画素を対象に、各画素の濃度に所定濃度を加算したり減算する）する。階調変換されたサブトラクション像はＤ／Ａ変換器３６に与えられ、そこでＤ／Ａ変換されてモニタ４に表示される。遅延回路３３は、周波数特性変換回路３２での処理時間による時間的な遅れを補償するために設けており、これにより、ライブ像とマスク像とが同期されて演算器３４に供給される。なお、画像処理部３を構成する各部の動作制御は、制御部５により行われる。
【００２９】
周波数特性変換回路３２は、基本画像（ライブ像）から血管像などの高周波数成分を除去して高周波数成分除去画像（マスク像）を得ることを、実空間データのままで高速処理することを目的とするものである。また、この周波数特性変換回路３２は、図３に示すように、しきい値変更手段７に接続されており、しきい値変更手段７から変更指示（しきい値周波数を変更させるための指示）が与えられるようになっており、この変更指示に従ってしきい値周波数が変更される。このしきい値変更手段７は、制御部５の機能うちで画像処理部３の周波数特性変換回路３２を制御する機能と、操作盤６の機能うちでしきい値周波数の変更指示を制御部５に与える機能とを備えたものであり、制御部５および操作盤６の構成うちの一部の構成に属するものである。
【００３０】
周波数特性変換回路３２は、実空間上で直線形状で表されるテンプレートフィルタでもって前記基本画像の高周波数成分を実空間上で除去する実空間フィルタリング部３２ａを備えている。この実空間フィルタリング部３２ａのテンプレートフィルタ形状は、例えば、図７に示すように、互いに直交するＸ，Ｙ，Ｚ軸で表される、実空間上の３次元直交座標系において、Ｘ，Ｙ軸方向をテンプレートフィルタの矩形形状である底面大きさとし、Ｚ軸方向をこのテンプレートのゲインとした四角柱形状で表されるものとしている。この実空間上で四角柱形状のテンプレートフィルタ（以下、適宜に「ボックスフィルタ」と呼ぶ）による処理方法は、図８に示すように、入力画像（入力される基本画像）における画素（ｉ，ｊ）の近傍（Ｎ×Ｎ点の正方形分）の平均濃度を出力画像の画素（ｉ，ｊ）の値とし、この平均処理を全画素について行うことで、基本画像から所定のしきい値周波数以上の周波数成分（高周波数成分）を除去した高周波数成分除去画像（マスク像）を得るものであり、移動平均フィルタ法と適宜に呼ぶこととする。
【００３１】
なお、図８は、図７に示したボックスフィルタを、フィルタ処理しようとする基本画像の所定の複数個の画素上に位置させた状態をＺ軸方向から見下ろした図である。図８では、ボックスフィルタのフィルタサイズを、説明の便宜上、３×３の正方形として図示し、画素（ｉ，ｊ）とその周囲の８画素とからなる９画素（９点）の平均値を画素（ｉ，ｊ）の値としているが、このフィルタサイズは、次に説明するように所定のしきい値周波数以上の周波数成分（高周波数成分）を除去する大きさに設定される。すなわち、このボックスフィルタは、そのフィルタサイズに応じてしきい値周波数が変更されることになる。つまり、フィルタサイズを大きくすればしきい値周波数を下げることになりぼかしの程度が大きくなるし、フィルタサイズを小さくすればしきい値周波数を上げることになりぼかしの程度が小さくなるのである。ボックスフィルタのフィルタサイズは、図７に示すように、例えば２１×２１点（ＴＡＰ数×ＴＡＰ数）の正方形としている。このＴＡＰ数はカーネルの点数のことでもある。また、標準的な空間的形状の血管像を好適に除去し得るしきい値周波数は、実験的あるいは理論的に求めておくことができるので、標準的な空間的形状の血管像を好適に除去し得るボックスフィルタのフィルタサイズも同様に実験的あるいは理論的に求めておいても良いし、そのフィルタサイズは初期値として設定しておいても良い。この第１実施例では、しきい値変更手段７からの変更指示は、フィルタサイズを変更するという形式で行われ、操作盤６の例えばダイヤルや操作ボタンなどを操作することで、変更量を連続的に変化させるように入力したり、操作盤６の例えばテンキーなどで特定の数値（変更後の値）を指定入力したりするなど、各種の入力方法でもって使用者により操作される。使用者は、モニタ４に表示されるサブトラクション像をモニタリングしながら、必要に応じて、対話的に上記の変更を行うことになる。
【００３２】
なお、上述の移動平均フィルタ法は、次に説明するようにすることで、高速に演算処理ができるように改良されたものである。例えば、フィルタサイズがＮ点×Ｎ点であるボックスフィルタを用いて移動平均を求めるには、ボックスフィルタが位置するＮ×Ｎ点の画素を正直に加算した合計値をボックスフィルタの点数分（Ｎ×Ｎ点）で平均することで、ボックスフィルタの中心の画素の値を算出し、この演算を全画素について個別に行なっていたのでは、演算量が膨大になり高速に演算処理することができないし、フィルタサイズが大きくなるにつれてその演算量は増加することになる。
【００３３】
そこで、次に説明するようにして、上述の移動平均フィルタ法を高速に演算処理可能に改良している。すなわち、図９に示すように、Ｎ×Ｎ点の移動平均を求めるには、各画素についてＮ×Ｎ点の加算を正直に計算することなく、前回の結果に対して、更新される部分についての加減算を行うことで、簡単に計算でき、演算時間が短縮できる。例えば、１次元状の「Ａ」〜「Ｆ」の画素に対して、１次元状の３点のボックスフィルタで移動平均する場合を例に挙げて説明する。「Ａ」〜「Ｃ」の画素にボックスフィルタがある場合は、「Ａ＋Ｂ＋Ｃ」の３点平均値が「Ｂ」画素の値となる。そして、ボックスフィルタを次の位置（「Ｂ」〜「Ｄ」の画素の位置）に移動させたときには、前回の値（「Ａ＋Ｂ＋Ｃ」）から更新される部分（「Ｄ」が加わり、「Ａ」が外される）についての加減算（−Ａ＋Ｄ）を行うだけでよく、前回と今回とで重複する加算を行なうという無駄が排除できるし、フィルタサイズが大きくなってもその演算量（演算の総量）は低減したまま一定である。この例では、説明の便宜上、ボックスフィルタを３点としているが、ボックスフィルタの点数が大きくなればなる程、全点について正直に順次加算していく演算量と、更新される部分の加減算のみを行なう演算量との差は大きくなり、効果的に改善されることがわかる。なお、図７，図８に示すようにＮ×Ｎ点のボックスフィルタの２次元移動平均の場合には、縦方向のＮ行累算値を、図９と同様の方法で求め、さらに横方向に図９と同様の方法を繰り返していけばよい。すなわち、Ｎ×Ｎ点のボックスフィルタの２次元移動平均の場合であっても、更新される部分の加減算は、１画素の加算および１画素の減算のみである。上述したように演算の総量を低減しているので、汎用のチップを用いて実空間フィルタリング部３２ａを構成することができ、基本画像を撮像してから極めて短時間で（リアルタイムに近く）サブトラクション像をモニタ４に表示させることができる。
【００３４】
図１，図２に戻って、制御部５は、操作盤６からの各種の指示などによって、各装置、各部の駆動制御や動作制御を行う。例えば、制御部５は、操作盤６から、ボックスフィルタのフィルタサイズを変更する変更指示を受けると、実空間フィルタリング部３２ａのボックスフィルタのフィルタサイズを変更指示に応じたサイズに変更する。この制御部５は、例えば、後述する動作を実現するプログラムを遂行するＣＰＵ（中央処理装置）で構成されている。
【００３５】
操作盤６は、撮像部位や条件の設定、処理開始指示、フィルタサイズを変更する変更指示などを、操作者が行うためのものである。操作盤６における、フィルタサイズを変更する構成としては、例えば、マウス、トラックボール、キーボード、ジョイスティックなどの入力装置が挙げられる。
【００３６】
上記構成を有する実施例装置の動作を以下に説明する。
まず、被検体Ｍのある１箇所の撮像部位（例えば胸部）のサブトラクション像を得る場合の動作を説明する。
【００３７】
この場合、まず、操作者により操作盤６から設定された撮像部位や条件（撮像方向など）に従って、制御部５は、モータ１３を駆動制御して被検体Ｍを載置した天板１２を水平移動させ、設定された撮像部位（胸部とする）を、Ｘ線管２１、撮像系２２の間の撮像位置に位置させ、モータ２５を駆動制御してＸ線管２１、撮像系２２を被検体Ｍ（の撮像部位）の体軸回りに変位させ、撮像方向を調節する。この状態を図４に示す。図４では、撮像部位ＳＢの下方からＸ線を照射してＸ線透過像を撮像するように撮像方向が調節されている。
【００３８】
次に、被検体Ｍに造影剤を投与する。なお、造影剤を投与してから上記位置合わせ動作などを行ってもよい。いずれにしても、以下の撮像動作の前に、被検体Ｍに造影剤を投与しておき、撮像部位ＳＢに造影剤が拡散した状態で操作者が操作盤６から処理開始を指示し、以下の撮像動作が実行される。
【００３９】
処理開始が指示されると、制御部５はＸ線高電圧発生装置２６を制御して、Ｘ線管２１に所定の電力を供給させてＸ線を照射させ、造影剤が拡散された撮像部位ＳＢのＸ線透過像を撮像させる。そして、制御部５は、画像処理部３の各部を制御して、基本画像（ライブ像）を得るとともに、その基本画像からマスク像を求め、ライブ像とマスク像とのサブトラクションを行わせサブトラクション像をモニタ４に表示させる。
【００４０】
したがって、ある撮像部位のサブトラクション像を得るための被検体へのＸ線照射は１回でよく、１フレームの基本画像からマスク像とライブ像を得ているので、被検体の体動によるマスク像とライブ像の画像のずれも完全に無い。
【００４１】
ここで仮に、モニタ４に表示されたサブトラクション像が好適なものでない場合には、使用者は操作盤６によってフィルタサイズの変更指示を入力する。制御部５は、実空間フィルタリング部３２ａのボックスフィルタのフィルタサイズを、操作盤６からの変更指示に応じたサイズに変更する。そして、Ｘ線透視撮像を行って基本画像（ライブ像）を得るか、または、Ａ／Ｄ変換器３１の出力側に設けられたメモリ（図示省略）に記憶しておいたディジタルデータとしての基本画像（ライブ像）を読み出すようにする。実空間フィルタリング部３２ａは、前記の変更されたフィルタサイズでもって基本画像（ライブ像）をフィルタリングして、マスク像（高周波数成分除去画像）を生成する。減算器３４により、基本画像（ライブ像）からマスク像（高周波数成分除去画像）をサブトラクションして得られたサブトラクション像は、階調変換回路３６やＤ／Ａ変換器３６を介して、フィルタサイズ変更に係る新たなサブトラクション像として、モニタ４にリアルタイムに表示される。このように、操作盤６への変更指示入力から、フィルタサイズ変更に係る新たなサブトラクション像のモニタ４への表示までは、即時に行われる。したがって、フィルタサイズ変更に係る新たなサブトラクション像がモニタ４にリアルタイムに表示されるので、適正なしきい値（適正なフィルタサイズ）を対話的にリアルタイムに設定することができ、モニタ４に表示されるサブトラクション像をモニタリングしながら好適なサブトラクション像を表示させることができる。
【００４２】
次に、被検体ＭとＸ線透視装置２との相対的な位置関係を被検体Ｍの体軸方向に変位させて、例えば、図５に示すように、胸部から腹部にかけての領域ＳＲ内の複数の撮像部位のサブトラクション像を得る場合の動作を説明する。なお、この実施例では、Ｘ線透視装置２を固定し、これに対して被検体Ｍを載置した天板１２を水平移動するように構成しているが、撮像状況をわかり易くするために、図５では、天板１２上の被検体Ｍを固定し、これに対してＸ線透視装置２（Ｘ線管２１、撮像系２２）が変位しているように描いている。
【００４３】
この場合、制御部５は、最初の撮像部位（図では、撮像領域ＳＲの左端部側の撮像部位）を撮像位置に位置させ、撮像方向を調節する。そして、以下の撮像動作の前に被検体Ｍに造影剤を投与しておく。
【００４４】
被検体Ｍの各撮像部位（領域ＳＲ）に造影剤が拡散し、処理開始が指示されると、上記１箇所の撮像領域ＳＢのサブトラクション像を求めた手順と同様の手順で、最初の撮像部位のサブトラクション像を求め、天板１２を図５の左方向に定速で移動させながら、以降の各撮像部位が撮像位置に位置するごとに、その撮像部位のサブトラクション像を順次求めていく。
【００４５】
したがって、被検体とＸ線透視装置との相対的な位置関係を被検体の体軸方向に変位させて、複数の連続する撮像部位に対するサブトラクション像を得る場合であっても、被検体へのＸ線曝射線量を必要最小限とし、被検体の体動や撮像する部位のずれなどに起因する、各撮像部位ごとの対となるマスク像とライブ像との画像のずれを完全に無くしている。
【００４６】
ここで仮に、上述の領域ＳＲを撮像している途中で好適なサブトラクション像が得られなくなった場合について説明する。好適なサブトラクション像が得られなくなると、使用者は操作盤６にフィルタサイズの変更指示を入力する。制御部５は、実空間フィルタリング部３２ａのボックスフィルタのフィルタサイズを、操作盤６からの変更指示に応じたサイズに変更する。操作盤６からのフィルタサイズの変更指示に従って、フィルタサイズ変更に係る新たなサブトラクション像がモニタ４にリアルタイムに表示されるので、適正なしきい値（適正なフィルタサイズ）を対話的にリアルタイムに設定することができ、モニタ４に表示されるサブトラクション像をモニタリングしながら好適なサブトラクション像を表示させることができる。
【００４７】
なお、天板１２を固定し、Ｘ線透視装置２を天板１２上の被検体Ｍの体軸方向に移動させることで、被検体とＸ線透視装置との相対的な位置関係を被検体の体軸方向に変位させるように構成してもよい。
【００４８】
また、被検体Ｍのある部位（例えば、胸部）を撮像位置に位置させた状態で、図６に示すように、Ｘ線管２１、撮像系２２をその部位の回り（体軸回り）に回転変位させながら、各撮像方向からのサブトラクション像を求めることもあるが、このような場合であっても、上記各動作と同様に、被検体ＭへのＸ線曝射線量を必要最小限とし、各撮像方向からのサブトラクション像を求めるための一対のマスク像とライブ像とに画像のずれが生じることがない。
【００４９】
また、実空間フィルタリング部３２ａは、実空間上で直線形状で表されるテンプレートフィルタでもって、基本画像の高周波数成分を実空間上で除去しており、次に説明するような効果がある。
【００５０】
例えば、従来のガウス関数形状のテンプレートフィルタでは、図１６に示すそのガウス曲線形状からもわかるように、基本画像を単純な加減算のみで処理することできないし、フィルタサイズ（Ｎ×Ｎ点）が例えば５１（ＴＡＰ数）×５１（ＴＡＰ数）点のように大きくなることでフィルタリングの演算量が膨大になり、従来例の場合における演算量は次に示す式（１）で表されて、リアルタイムに近い短時間での処理は不可能であった。なお、式（１）中のＰは、１回の乗算と１回の加算との２回であり、式（１）中のＧ2 は、基本画像の全画素数である。
例えば、参考までに式（１）に具体的な数値を代入してみる。
演算量＝Ｎ（ＴＡＰ数）×Ｎ（ＴＡＰ数）×Ｐ×Ｇ2     …  （１）
＝５１×５１×２×Ｇ2
【００５１】
これに対して、上述した第１実施例では、実空間フィルタリング部３２ａは、実空間上で直線形状で表される、四角柱形状のテンプレートフィルタでもって、基本画像をフィルタリングするので、基本画像を単純な加減算のみで処理することできるし、さらに、フィルタリング処理における重複する演算（加算）を繰り返し実行することがないように演算量を軽減することができる。この第１実施例の場合におけるフィルタリングの演算量は次に示す式（２）で表される。なお、式（２）中のＱは、１回の加算と１回の減算との２回であり、式（２）中のＧ2 は、前述の式（１）と同様に、基本画像の全画素数である。例えば、参考までに式（２）に具体的な数値を代入してみる。
演算量＝Ｑ×Ｇ2     …  （２）
＝２×Ｇ2
【００５２】
式（１），（２）を比してわかるように、第１実施例の場合におけるフィルタリングの演算量は、フィルタ移動によって更新される部分についての単純な加減算を実行するだけで良く、式（２）中にＴＡＰ数が含まれていないことから、フィルタサイズに依らずに一定量でしかも式（１）に比べて大幅に低減され、リアルタイムに近い短時間での処理が可能となることがわかる。したがって、上述のボックスフィルタを用いた場合であっても、しきい値変更手段７からしきい値周波数の変更指示（ボックスフィルタのフィルタサイズ変更指示）が入力されると、フィルタサイズ変更に係る新たなサブトラクション像をモニタ４にリアルタイムに表示させることができる。
【００５３】
＜第２実施例＞
続いて、この発明のＸ線ディジタルサブトラクション装置の第２実施例に係るＸ線ディジタルアンギオグラフィ装置について説明する。図１０は、この発明の第２実施例に係る実空間フィルタリング部４０のテンプレートフィルタ形状を示す模式図である。図１１は、この第２実施例に係る実空間フィルタリング部４０のブロック図である。この第２実施例のＸ線ディジタルアンギオグラフィ装置は、第１実施例の実空間フィルタリング部３２ａに替えて、図１１に示す実空間フィルタリング部４０を採用し、第２実施例の実空間フィルタリング部４０のテンプレートフィルタ形状を、四角柱形状から後述するような八角錐形状に替えた点以外については、前述の第１実施例と同様であるので、特にテンプレートフィルタ形状が実空間上で後述するような八角錐形状である実空間フィルタリング部４０の構成および機能について詳細に説明するものとする。
【００５４】
まず、実空間フィルタリング部４０のテンプレートフィルタ形状について説明する。実空間フィルタリング部４０のテンプレートフィルタ形状は、図１０に示すように、互いに直交するＸ，Ｙ，Ｚ軸で表される、実空間上の３次元直交座標系において、Ｘ，Ｙ軸方向をテンプレートフィルタの矩形形状である底面大きさとし、Ｚ軸方向をこのテンプレートのゲインとし、Ｚ軸方向に先細りとなる立体形状であって、かつ、テンプレートフィルタのＺ方向の少なくとも中央から先細りの先端部にかけて八角錐形状で表されるものとしている。なお、このような実空間フィルタリング部４０のテンプレートフィルタ形状を、説明の便宜上、略八角錐形状とも呼ぶこととする。
【００５５】
実空間フィルタリング部４０は、図１１に示すように、入力される基本画像の行方向に後述する所定のフィルタをかける行方向三角形フィルタ４１と、この行方向三角形フィルタ４１でフィルタ処理されたデータを記憶する記憶手段４２と、この記憶手段４２に記憶された行方向処理後の基本画像の列方向に後述する所定のフィルタをかける列方向三角形フィルタ４３と、行方向三角形フィルタ４１と記憶手段４２と列方向三角形フィルタ４３とを制御する制御手段４４とを備えている。
【００５６】
行方向三角形フィルタ４１は、図１２に示すように、更新される部分のうちで新たに加えられる、後述する三角形状の右側の下がり勾配部分の画素の値（図１４（ｂ）参照）を生成する第１の１次元フィルタブロック５１と、更新される部分のうちで削除される、後述する三角形状の左側の上がり勾配部分の画素の値（図１４（ｂ）参照）を生成する第２の１次元フィルタブロック５２と、これらの第１の１次元フィルタブロック５１と第２の１次元フィルタブロック５２とからの画像データを減算する減算器５３と、この減算器５３からの出力を累算して出力する累算回路５４と、入力される基本画像を第１の１次元フィルタブロック５１よりも所定量遅延させて第２の１次元フィルタブロック５２に入力するための遅延回路５５とを備えている。
【００５７】
第１の１次元フィルタブロック５１は、一方の入力ポートに入力される画像データ（基本画像の行方向についての一連の画像データ）から、他方の入力ポートから入力される画像データ（基本画像の行方向についての一連の画像データ）を減算して出力する減算器６１と、入力される基本画像を減算器６１の一方の入力ポートよりも所定量遅延させて減算器６１の他方の入力ポートに入力するための遅延回路６２と、この減算器６１からの出力を累算して出力する累算回路６３とを備えている。なお、第２の１次元フィルタブロック５２も、第１の１次元フィルタブロック５１と同様に構成されており、減算器６１と累算回路６２と遅延回路６３とを備えている。
【００５８】
なお、列方向三角形フィルタ４３は、前述の行方向三角形フィルタ４１と同様に構成されており、行方向三角形フィルタ４１で処理されて記憶手段４２に記憶された基本画像の行方向のデータを、列方向にフィルタリング処理するものである。
【００５９】
ここで、テンプレートフィルタ形状が実空間上で略八角錐形状である実空間フィルタリング部４０について、そのフィルタリング機能について説明する。この実空間フィルタリング部４０は、実空間の１次元形状が三角形をしている図１３（ａ）に示す１次元フィルタを、図１０に示すように２次元に拡張したものと等価である。図１３（ｃ）に示すように、行列（縦横）分離型のフィルタで構成することで、高速化に有利となる。なお、フィルタ係数の一例を図１３（ｂ），（ｃ）に示す。例えば、図１３（ａ）に示すように実空間形状が三角形である１次元フィルタを５点のものとすると、図１３（ｂ）に示すように左から順に「１，２，３，２，１」の係数となる。図１３（ｂ）に示す行（横）方向の１次元フィルタを列（縦）方向に並べると、図１３（ｃ）の右上側に示すような２次元フィルタとなる。また、図１３（ｂ）に示す行（横）方向の１次元フィルタを縦方向にしたもの、つまり上から順に「１，２，３，２，１」の係数としたものを、行（横）方向に並べると、図１３（ｃ）の左上側に示すような２次元フィルタとなる。そして、これらの２次元フィルタをその同一点同士の係数をかけることで、図１３（ｃ）の下側に示すような、行列（縦横）分離型のフィルタ、つまり、実空間上で略八角錐形状（略ピラミッド形状）であるテンプレートフィルタと等価な効果が得られる。
【００６０】
なお、図１３に示してきた係数（「１，２，３，４，６，９」）は、図１４に説明するように、重み度を示すものである。すなわち、図１４（ａ）に示すように、図１３（ｂ）に示す１次元フィルタが基本画像の「Ａ」〜「Ｅ」の画素に位置する場合に、それらの画素の値を何倍するかを示している。具体的に、画素「Ａ」は図１３（ｂ）に示す１次元フィルタの係数「１」にあるので「Ａ」を１倍したもの、つまり「Ａ」のままとなり、画素「Ｂ」は図１３（ｂ）に示す１次元フィルタの係数「２」にあるので「Ｂ」を２倍したもの、つまり「２Ｂ」となり、画素「Ｃ」は図１３（ｂ）に示す１次元フィルタの係数「３」にあるので「Ｃ」を３倍したもの、つまり「３Ｃ」となり、画素「Ｄ」は図１３（ｂ）に示す１次元フィルタの係数「２」にあるので「Ｄ」を２倍したもの、つまり「２Ｄ」となり、画素「Ｅ」は図１３（ｂ）に示す１次元フィルタの係数「１」にあるので「Ｅ」を１倍したもの、つまり「Ｅ」のままとなる。
【００６１】
図１３（ｃ）の下側に示すような、行列（縦横）分離型のフィルタ、つまり、実空間上で略八角錐形状（略ピラミッド形状）であるテンプレートフィルタは、列（縦）方向と行（横）方向とに分離すると、図１４に示すように、前述の第１実施例での移動平均のアルゴリズムの延長で計算でき、フィルタリング処理における重複する演算（加算）を繰り返し実行することがないように演算量を軽減することができる。
【００６２】
すなわち、第１の１次元フィルタブロック５１は、図１４（ｂ）に示すように、更新される部分のうちで新たに加えられる、三角形状の下がり勾配部分の画素の値を生成する。図１４（ｂ）では、更新される部分のうちで新たに加えられる、三角形状の下がり勾配部分の画素の値は、「Ｄ＋Ｅ＋Ｆ」である。
【００６３】
また、第２の１次元フィルタブロック５２は、図１４（ｂ）に示すように、更新される部分のうちで削除される、三角形状の上がり勾配部分の画素の値を生成する。図１４（ｂ）では、更新される部分のうちで削除される、三角形状の上がり勾配部分の画素の値は、「Ａ＋Ｂ＋Ｃ」である。
【００６４】
図１３（ｂ）に示す１次元フィルタが基本画像の「Ａ」〜「Ｅ」の画素に位置する場合には、図１４（ａ）に示すように、重み付けされたもの「Ａ＋２Ｂ＋３Ｃ＋２Ｄ＋Ｅ」が得られる。次に、図１３（ｂ）に示す１次元フィルタが画素１つ分移動して基本画像の「Ｂ」〜「Ｆ」の画素に位置する場合には、図１４（ａ）に示したものから「Ａ＋Ｂ＋Ｃ」を引くとともに「Ｄ＋Ｅ＋Ｆ」を加えることで、図１４（ｂ）に示すように、重み付けされたもの「Ｂ＋２Ｃ＋３Ｄ＋２Ｅ＋Ｆ」が得られる。したがって、図１３（ｂ）に示す１次元フィルタが画素１つ分移動して基本画像の「Ｂ」〜「Ｆ」の画素に位置する場合に、最初から正直に演算して、つまり、前回（図１４（ａ）に示すもの）と今回（図１４（ｂ）に示すもの）とで重複する演算を行なうことで、図１４（ｂ）に示す重み付けされたもの「Ｂ＋２Ｃ＋３Ｄ＋２Ｅ＋Ｆ」を得る必要はない。次に、図１３（ｂ）に示す１次元フィルタがさらに画素１つ分移動して基本画像の「Ｃ」〜「Ｇ」の画素に位置する場合には、同様に、図１４（ｂ）に示したもの「Ｂ＋２Ｃ＋３Ｄ＋２Ｅ＋Ｆ」から「Ｂ＋Ｃ＋Ｄ」を引くとともに「Ｅ＋Ｆ＋Ｇ」を加えて重み付けされたもの「Ｃ＋２Ｄ＋３Ｅ＋２Ｆ＋Ｇ」が得られる。
【００６５】
このように、フィルタ移動によって更新される部分の加減算は、３画素分の加算および３画素の減算のみであるので、この更新される部分の加減算を行なうだけで良い。
【００６６】
次に、行方向三角形フィルタ４１によって、図１４（ｂ）に示す「Ｂ＋２Ｃ＋３Ｄ＋２Ｅ＋Ｆ」を得る動作を、具体的に説明する。なおこの図１４では、説明の便宜上、フィルタサイズを５点（ＴＡＰ数＝５点）としているので、第１の１次元フィルタブロック５１の遅延回路６２と、第２の１次元フィルタブロック５２の遅延回路６２と遅延回路５５との遅延量は、次に示す式（３）により、それぞれ「３」に設定されている。
遅延量＝（ＴＡＰ数＋１）／２    …  （３）
＝（５＋１）／２＝３
【００６７】
これらの遅延回路のパラメータ（遅延量）などがデジタルフィルタで言うところのＴＡＰ数に相当する。ＴＡＰ数を変えることで、しきい値周波数の値を変更できる。なお、図１４では、説明の便宜上、フィルタサイズを５点としていたが、この第２実施例では、図１０に示すようにＴＡＰ数を５１×５１点とした場合には、第１の１次元フィルタブロック５１の遅延回路６２と、第２の１次元フィルタブロック５２の遅延回路６２と遅延回路５５との遅延量は、「２６」に設定される。制御手段４４は、第１の１次元フィルタブロック５１の累算回路６３と、第２の１次元フィルタブロック５２の累算回路６３と、累算回路５４との累算値の初期化や、第１の１次元フィルタブロック５１の遅延回路６２と、第２の１次元フィルタブロック５２の遅延回路６２と遅延回路５５との遅延量の設定を行なう。
【００６８】
例えば、基本画像の「Ａ」〜「Ｆ」が行方向三角形フィルタ４１に入力された時点での、この行方向三角形フィルタ４１の各構成での処理状況について見てみる。第１の１次元フィルタブロック５１の減算器６１には「Ｆ」と「Ｃ」とが入力されて「Ｆ−Ｃ」が第１の１次元フィルタブロック５１の累算回路６３に入力される。第１の１次元フィルタブロック５１の累算回路６３は、（「Ｃ＋Ｄ＋Ｅ」＋「Ｆ−Ｃ」）により「Ｄ＋Ｅ＋Ｆ」が出力される。また、遅延回路５５の遅延量が「３」であるので、第２の１次元フィルタブロック５２の減算器６１には未だ「Ｃ」のみが入力されるだけで「Ｃ」が第２の１次元フィルタブロック５２の累算回路６３に入力される。第２の１次元フィルタブロック５２の累算回路６３は、（「Ａ＋Ｂ」＋「Ｃ」）により「Ａ＋Ｂ＋Ｃ」が出力される。減算器５３では、第１の１次元フィルタブロック５１の累算回路６３から出力された「Ｄ＋Ｅ＋Ｆ」と、第２の１次元フィルタブロック５２の累算回路６３から出力された「Ａ＋Ｂ＋Ｃ」との減算が行なわれ、（「Ｄ＋Ｅ＋Ｆ」−「Ａ＋Ｂ＋Ｃ」）が出力される。累算回路５４は、前回の値「Ａ＋２Ｂ＋３Ｃ＋２Ｄ＋Ｅ」と、減算器５３からの（「Ｄ＋Ｅ＋Ｆ」−「Ａ＋Ｂ＋Ｃ」）とを累算し、「Ｂ＋２Ｃ＋３Ｄ＋２Ｅ＋Ｆ」を記憶手段４２に出力する。
【００６９】
以上、上述した第２実施例では、実空間フィルタリング部４０は、実空間上で直線形状で表され、略八角錐形状のテンプレートフィルタでもって、基本画像の高周波数成分を実空間上で除去しており、更新される部分の加減算のみを行なえば良いので、前述の第１実施例の場合と同等に、フィルタリングの演算量は、式（２）で表され、この式（２）中にＴＡＰ数が含まれていないことから、フィルタサイズに依らずに一定量でしかも、従来例の式（１）に比べて大幅に低減され、リアルタイムに近い短時間での処理が可能となる効果を有している。さらに、この第２実施例では、次に説明する点で、前述の第１実施例よりもさらに優れている。
【００７０】
すなわち、第２実施例のテンプレートフィルタ形状（図１０参照）は、前述の第１実施例のテンプレートフィルタ形状（図７参照）に比べて、理想的なテンプレートフィルタ形状（図１６参照）に近似しているので、図１５に示すように、前述の第１実施例の場合よりも、周波数特性に優れている。具体的には、前述の第１実施例のテンプレートフィルタでは０．２５ｌｐ／ｍｍ（ラインペア／ミリ）付近でリンギング（ＭＡＸが０．１５程度の山）が生じているが、このようなリンギングは第２実施例のテンプレートフィルタではかなり小さくなっていて僅かに存在しているだけであることから、前述の第１実施例に比べて、所定のしきい値周波数以上の周波数成分を良好に除去できることがわかる。
【００７１】
また、しきい値周波数の値の変更が、遅延回路の遅延量の調整で実現でき、演算量は増減しないため、しきい値周波数の値に依らず、低減した一定時間で処理できる。したがって、所定の処理時間を確保したまま、しきい値周波数の値の変更の自由度が高いディジタルサブトラクション装置が提供できる。
【００７２】
したがって、上述の略八角錐形状のテンプレートフィルタを用いた場合であっても、しきい値変更手段７からしきい値周波数の変更指示（略八角錐形状のテンプレートフィルタのフィルタサイズ変更指示）が入力されると、フィルタサイズ変更に係る新たなサブトラクション像をモニタ４にリアルタイムに表示させることができる。
【００７３】
なお、この第２実施例では、実空間フィルタリング部４０のテンプレートフィルタ形状を、図１０に示すように、テンプレートフィルタのＺ方向の少なくとも中央から先細りの先端部にかけて八角錐形状で表されるものとした略八角錐形状としているが、底面から先細りの先端部までの全部を八角錐形状としても良い。
【００７４】
＜第３実施例＞
図１７は、第３実施例装置の画像処理部の構成を示すブロック図である。この第３実施例は、基本画像（ライブ像）の高周波数成分を強調するための第２の周波数特性変換回路７１を付設したことを特徴とする。また、この第３実施例装置の周波数特性変換回路３２は、前述の第１，２実施例のような実空間フィルタリング部３２ａ，４０を有するものではなく、前記基本画像を構成する周波数成分を空間／周波数変換処理で取り出し、所定のしきい値周波数以上の周波数成分を除去し、これに周波数／空間変換処理を施して、前記基本画像から前記しきい値周波数以上の周波数成分を除去した高周波数成分除去画像を得るものである。この高周波数成分除去画像は、造影剤が投与された血管像などの高周波数成分が除去された画像であり、これをマスク像として用いる。前記の空間／周波数変換処理としては、ＦＦＴ（高速フーリエ変換）、カルーネン・レーベ変換、ＤＣＴ（離散コサイン変換）、アダマール変換などの各変換方式がある。また、前記の周波数／空間変換処理としては、前記空間／周波数変換処理の逆変換（逆ＦＦＴ、逆カルーネン・レーベ変換、逆ＤＣＴ、逆アダマール変換など）がある。
【００７５】
基本画像の高周波数成分を強調する処理は、例えば、空間／周波数変換して得られた各基本画像の周波数成分のうちの高周波数成分部分を増分させ、それを周波数／空間変換することで実現できる。第２の周波数特性変換回路７１ではこのような処理を実現するように構成される。
【００７６】
これにより、関心物（例えば血管像）がより強調されたサブトラクション像を得ることができる。また、関心物（例えば血管像）がより強調されたサブトラクション像を得る場合であっても、しきい値変更手段７からしきい値周波数の変更指示を入力することにより、しきい値周波数変更に係る新たなサブトラクション像をモニタ４にリアルタイムに表示させることができる。なお、この第３実施例に係る特徴部分は、以下の第４実施例にも同様に適用することができる。
【００７７】
＜第４実施例＞
図１８は、第４実施例装置の画像処理部の構成を示すブロック図である。この第４実施例は、複数回分の撮像で得られた基本画像の平均画像を求める積分処理部８１を付設したことを特徴とする。
【００７８】
積分処理部８１は、図１９（ａ）に示すように、第２の階調変換回路９１と加算器９２とメモリ９３とで構成してもよいし、図１９（ｂ）に示すように加算器９２とメモリ９３と除算器９４とで構成してもよい。
【００７９】
例えば、Ｎ回分の撮像で得られた基本画像の平均画像を求める場合、図１９（ａ）の構成では、順次与えられる基本画像の濃度（基本画像を構成する各画素の濃度）を第２の階調変換回路９１で１／Ｎにするように階調変換させ、加算器９２で、メモリ９３に記憶されている加算画像に加算していく。ただし、最初の基本画像が与えられるときには、メモリ９３には何も記憶されていない。
【００８０】
また、図１９（ｂ）の構成では、加算器９２で、順次与えられる基本画像を、メモリ９３に記憶されている加算画像（最初の基本画像が与えられるときには、メモリ９３には何も記憶されていない）に順次加算していき、最後に、メモリ９３に記憶されたＮ回分の基本画像の加算結果を除算器９４で１／Ｎにする。
【００８１】
このように、複数回分の撮像で得られた基本画像の平均画像を求めることにより、基本画像のＳ／Ｎ比を向上させることができ、この発明では、基本画像からマスク像とライブ像を求めるので、マスク像とライブ像のＳ／Ｎ比をまとめて向上させることができ、高画質なサブトラクション像を得ることができる。また、高画質なサブトラクション像を得る場合であっても、しきい値変更手段７からしきい値周波数の変更指示を入力することにより、しきい値周波数変更に係る新たなサブトラクション像をモニタ４にリアルタイムに表示させることができる。なお、周波数特性変換回路３２に前述の第１，２実施例の実空間フィルタリング部３２ａを備えた場合には、しきい値変更手段７からしきい値周波数の変更指示としてのフィルタサイズの変更指示を入力することにより、フィルタサイズ変更に係る新たなサブトラクション像がモニタ４にリアルタイムに表示される。
【００８２】
なお、例えば、ある撮像部位に対する１撮像方向からのサブトラクション像を得る場合には、その撮像部位に１方向からのみＸ線照射して基本画像を撮像する動作をＮ回繰り返すことになる。従来装置においても、マスク像やライブ像のＳ／Ｎ比を向上させるために、マスク像とライブ像をそれぞれ複数（Ｎ）回分撮像し、各画像の平均画像を求めることがあるが、この場合、被検体へのＸ線照射は、マスク像とライブ像を別々に撮像するので、２×Ｎ回行うことになる。これに対してこの実施例では被検体へのＸ線照射は従来装置の半分でよい。
【００８３】
また、被検体とＸ線透視装置との相対的な位置関係を被検体の体軸方向に変位させながら、複数の撮像部位のサブトラクション像を得る場合には、対象となる撮像部位の相前後する複数の撮像部位で撮像された複数回分の基本画像からその対象となる撮像部位の平均画像を求めるようにしてもよい。ある部位の周回方向にＸ線管、撮像系を回転変位させながら複数の撮像方向からのサブトラクション像を得る場合も同様に、対象となる撮像方向の相前後する複数の撮像方向から撮像された複数回分の基本画像からその対象となる撮像方向からの平均画像を求めるようにしてもよい。これらの場合でも、被検体へのＸ線曝射線量は従来の場合（従来装置で同様の手順で平均画像を求める場合）の半分になる。
【００８４】
この発明は、上記の各実施例に限られるものではなく、下記のように変形実施することができる。
【００８５】
（１）上述の各実施例のディジタルサブトラクション装置における実空間フィルタリング３２ａ，４０を汎用チップでもってソフトウエア的に実現することもできる。
【００８６】
（２）上述の各実施例のディジタルサブトラクション装置は、上述の各実施例のように、被検体Ｍを人体などとして医療用に用いることもできるし、被検体ＭをＢＧＡ（Ball Grid Array）基板やプリント配線基板など各種の電子部品などとして非破壊検査用に用いることもできる。
【００８７】
【発明の効果】
以上の説明から明らかなように、請求項１に記載のディジタルサブトラクション装置によれば、サブトラクション像を表示する表示手段と、与えられた変更指示に従ってしきい値周波数を変更するしきい値周波数変更手段とを備えているので、しきい値周波数の変更指示により、しきい値周波数変更に係る新たなサブトラクション像が表示手段に表示されるので、適正なしきい値を対話的に設定することができ、表示手段に表示されるサブトラクション像をモニタリングしながら好適なサブトラクション像を得ることができる。
【００８８】
また、請求項２に記載のディジタルサブトラクション装置によれば、基本画像の高周波数成分を強調した高周波数成分強調画像を得て、この高周波数成分強調画像と高周波数成分除去画像とのサブトラクションを行ってサブトラクション像を求めるので、関心物（例えば血管像）がより強調されたサブトラクション像を得ることができる。また、関心物がより強調されたサブトラクション像を得る場合であっても、しきい値周波数の変更により、好適なサブトラクション像を得ることができる。
【００８９】
また、請求項３に記載のディジタルサブトラクション装置によれば、複数回数分の撮像で得られた基本画像の平均画像を求め、この平均画像をライブ像として用いるとともに、このライブ像（平均画像）から高周波数成分除去画像を得て、ライブ像（平均画像）からマスク像（平均画像から得られた高周波数成分除去画像）をサブトラクションしてサブトラクション像を求めるので、ライブ像及びマスク像のＳ／Ｎ比を向上させることができ、高画質なサブトラクション像を得ることができる。また、高画質なサブトラクション像を得る場合であっても、しきい値周波数の変更により、好適なサブトラクション像を得ることができる。
【００９０】
また、本発明のディジタルサブトラクション装置によれば、実空間上で直線形状で表されるテンプレートフィルタでもってフィルタリングするので、フィルタ移動によって更新される部分についての加減算のみを行なうだけで良く、フィルタサイズに依らずにそのフィルタリングの演算量を一定量に低減でき、リアルタイムに近い短時間での処理ができ、しきい値周波数の変更により、リアルタイムに好適なサブトラクション像が得られる。
【００９２】
また、本発明のディジタルサブトラクション装置によれば、実空間フィルタリング部のテンプレートフィルタ形状は、互いに直交するＸ，Ｙ，Ｚ軸で表される、実空間上の３次元直交座標系において、Ｘ，Ｙ軸方向をテンプレートフィルタの矩形形状の底面大きさとし、Ｚ軸方向をテンプレートのゲインとし、Ｚ軸方向に先細りとなる立体形状であって、かつ、テンプレートフィルタのＺ方向の少なくとも中央から先細りの先端部にかけて角錐形状で表されるものとしているので、フィルタ移動によって更新される部分についての加減算のみを行なうだけで良く、フィルタサイズに依らずにそのフィルタリングの演算量を一定量に低減でき、リアルタイムに近い短時間での処理が可能となる。また、しきい値周波数の変更により、リアルタイムに好適なサブトラクション像が得られる。さらに、テンプレートフィルタ形状が四角柱形状の場合よりも、高周波数成分を除去する周波数特性に優れる。
【図面の簡単な説明】
【図１】  この発明の第１実施例に係るディジタルアンギオグラフィ装置の全体構成を示す正面図である。
【図２】  Ｘ線透視装置を側面から見た図である。
【図３】  第１実施例装置に備えられた画像処理部の構成を示すブロック図である。
【図４】  ある１箇所の撮像部位のサブトラクション像を得る場合の動作を説明するための図である。
【図５】  被検体の体軸方向の複数の撮像部位のサブトラクション像を得る場合の動作を説明するための図である。
【図６】  ある部位に対する複数の撮像方向からのサブトラクション像を得る場合の動作を説明するための図である。
【図７】  第１実施例のテンプレートフィルタ形状を示す模式図である。
【図８】  第１実施例の移動平均フィルタ法を説明するための模式図である。
【図９】  第１実施例の移動平均フィルタの計算アルゴリズムを説明するための模式図である。
【図１０】  第２実施例のテンプレートフィルタ形状を示す模式図である。
【図１１】  この発明の第２実施例に係る実空間フィルタリング部のブロック図である。
【図１２】  この発明の第２実施例に係る行方向三角形フィルタの構成を示すブロック図である。
【図１３】  （ａ）〜（ｃ）は第２実施例のピラミッド形フィルタを説明するための模式図である。
【図１４】  （ａ）、（ｂ）は第２実施例のピラミッド形フィルタの計算アルゴリズムを説明するための模式図である。
【図１５】  各種フィルタの周波数特性を示す特性図である。
【図１６】  従来のガウス関数形状で表される理想的なテンプレートフィルタ形状を示す模式図である。
【図１７】  第３実施例装置の画像処理部の構成を示すブロック図である。
【図１８】  第４実施例装置の画像処理部の構成を示すブロック図である。
【図１９】  （ａ），（ｂ）は第４実施例の積分処理部の構成を示すブロック図である。
【符号の説明】
２  …  Ｘ線透視装置
３  …  画像処理部
４  …  モニタ
７  …  しきい値変更手段
３１  …  Ａ／Ｄ変換器
３２  …  周波数特性変換回路
３３  …  遅延回路
３４  …  演算器
３２ａ…  実空間フィルタリング部
４０  …  実空間フィルタリング部
７１  …  第２の周波数特性変換回路
８１  …  積分処理部
Ｍ  …  被検体
ＳＢ  …  撮像部位
【図１】
0004626086.tif 000002
【図２】
0004626086.tif 000003
【図３】
0004626086.tif 000004
【図４】
0004626086.tif 000005
【図５】
0004626086.tif 000006
【図６】
0004626086.tif 000007
【図７】
0004626086.tif 000008
【図８】
0004626086.tif 000009
【図９】
0004626086.tif 000010
【図１０】
0004626086.tif 000011
【図１１】
0004626086.tif 000012
【図１２】
0004626086.tif 000013
【図１３】
0004626086.tif 000014
【図１４】
0004626086.tif 000015
【図１５】
0004626086.tif 000016
【図１６】
0004626086.tif 000017
【図１７】
0004626086.tif 000018
【図１８】
0004626086.tif 000019
【図１９】
0004626086.tif 000020