(57)【特許請求の範囲】
【請求項１】
コンピュータ断層像撮影手段で撮影した被検体の断層像を入力する断層像入力手段と、前記断層像における各画素ごとに画素値の時間変化情報を取り出す手段と、各画素ごとに取り出した時間変化情報からノイズを除去するノイズ除去手段を有し、各画素ごとにノイズ除去した時間変化情報に基づいて血流動態を解析する演算手段と、からなる血流動態解析装置において、
前記ノイズ除去手段は、各画素ごとに取り出された画素値の時間変化情報と、該画素値の時間変化情報をカーブフィッティングした関数との誤差に基づく値であるパラメータを算出するパラメータ算出手段と、前記算出したパラメータに対応するマトリクスサイズを有する画像フィルタを選択するフィルタ選択手段と、を有し、画素ごとに画像フィルタを選択してスムージング処理を行うことを特徴とする血流動態解析装置。
【請求項２】
前記フィルタ選択手段は、前記パラメータと、前記画像フィルタとを対応付けて記憶するルックアップテーブルを有し、前記算出したパラメータに基づいて前記ルックアップテーブルから対応する画像フィルタを選択することを特徴とする請求項１の血流動態解析装置。
【請求項３】
前記パラメータは、前記画素値の時間変化情報と、該画素値の時間変化情報をカーブフィッティングした関数との誤差の絶対値の標準偏差値であることを特徴とする請求項１又は２に記載の血流動態解析装置。
【請求項４】
コンピュータ断層像撮影手段で撮影した被検体の断層像を入力する断層像入力手段と、前記断層像における各画素ごとに画素値の時間変化情報を取り出す手段と、各画素ごとに取り出した時間変化情報の周波数空間でのノイズを除去するノイズ除去手段を有し、各画素ごとにノイズ除去した時間変化情報に基づいて血流動態を解析する演算手段と、からなる血流動態解析装置において、
前記ノイズ除去手段は、各画素ごとに取り出された時間変化情報に基づいて高周波数領域制御フィルタの形状を決定するフィルタ形状決定手段を有し、画素ごとに高周波数領域制御フィルタを決定してフィルタリング処理を行うことを特徴とする血流動態解析装置。
【請求項５】
コンピュータ断層像撮影手段で撮影した被検体の断層像の入力を受付ける工程と、
前記断層像における各画素ごとに画素値の時間変化情報を取り出す工程と、
各画素ごとに取り出された画素値の時間変化情報と、該画素値の時間変化情報をカーブフィッティングした関数との誤差に基づく値であるパラメータを算出する工程と、
前記算出したパラメータに対応するマトリクスサイズを有する画像フィルタを選択するフィルタ選択工程と、
前記フィルタ選択工程により選択された画素ごとの画像フィルタに基づいてスムージング処理を行なうことにより、各画素ごとに取り出した時間変化情報からノイズを除去する工程と、
前記ノイズを除去する工程により各画素ごとにノイズ除去した時間変化情報に基づいて血流動態を解析する演算工程と、
を含むことを特徴とする血流動態解析方法。
【請求項６】
コンピュータ断層像撮影手段で撮影した被検体の断層像の入力を受付ける断層像入力工程と、
前記断層像における各画素ごとに画素値の時間変化情報を取り出す工程と、
各画素ごとに取り出された時間変化情報に基づいて高周波数領域制御フィルタの形状を決定するフィルタ形状決定工程と、
画素ごとに高周波数領域制御フィルタを決定してフィルタリング処理を行うことにより各画素ごとに取り出した時間変化情報の周波数空間でのノイズを除去する工程と、
前記ノイズを除去する工程により各画素ごとにノイズ除去した時間変化情報に基づいて血流動態を解析する演算工程と、
を含むことを特徴とする血流動態解析方法。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は、血流動態解析装置及び方法に係り、特にＸ線や電磁波を用いて被検体の断層像から血流動態の情報を得る血流動態解析装置及び方法に関する。
【０００２】
【従来の技術】
血流動態を計測するコンピュータ断層診断装置には、陽電子放射断層撮影装置（Positron Emission Tomography：PET ）や単光子放射ＣＴ（Single Photon Emission CT ：SPECT ）などの核医学装置がある。核医学装置による血流動態計測では、放射性核種を被検体内に注入し、核種から放射された放射線をシンチレーションカメラで計測することで、横断面についての放射性核種分布像が得られる。この放射性核種分布像を解析することにより、血流量や平均通過時間に代表される血流動態の情報を得ていた。
【０００３】
血流動態を計測する他のコンピュータ断層診断装置には、Ｘ線ＣＴ装置がある。単純Ｘ線ＣＴ像により病変の超早期相を診断することが困難な場合に、Ｘ線ＣＴ装置による血流動態計測では、コントラストを付けるために主にヨード系の造影剤を用いていた。造影剤の濃度とＣＴ値とは比例関係にあるため、造影剤を被検体に注入しダイナミックスキャンを行うことで、撮影断面の時間変化の情報が得られる。この撮影断面の時間変化から、各組織のＣＴ値の時間変化を示す時間−濃度曲線が得られ、この時間−濃度曲線を解析することにより、血流量や平均通過時間に代表される血流動態の情報を得ていた。
【０００４】
血流動態を計測する更に他のコンピュータ断層診断装置には、ＭＲＩ（Magnetic Resonance Imaging）装置がある。高磁場ＭＲＩ装置により非造影で診断する方法もあるが、高磁場ＭＲＩ装置によりコントラストの良好な画像を得る場合、低磁場ＭＲＩ装置により診断する場合においては、主にガドリニウム系の造影剤を用いていた。ＭＲＩ装置についても、撮影断面の時間変化の情報から、各組織のＭＲ値の時間変化を示す時間−濃度曲線が得られ、この時間−濃度曲線を解析することにより、血流量や平均通過時間に代表される血流動態の情報を得ていた。
【０００５】
Ｘ線ＣＴ装置やＭＲＩ装置等の画像診断装置により、脳などの器官の血流動態を調べる場合に、コントラストの良好な投影像を得るために造影剤を用いるが、この造影剤は被検体内に急速に注入する必要がある。造影剤の注入量には限度があるため、組織のＣＴ値やＭＲ信号値の変化は、ノイズを打ち消す程度に大きくならない。このため、ＳＮ比（signal-to-noise ratio : 信号と雑音との振幅の比）を上げるためにＣＴ画像やＭＲＩ画像に対して一定のマトリクスサイズでスムージング処理を行っていた。
【０００６】
高磁場ＭＲＩ装置により脳などの器官の血流動態を調べる場合、造影剤を用いずに検査することも可能である。この場合、撮影時に印加できる傾斜磁場の大きさには限度があるため、ＭＲ信号値の変化はノイズを打ち消す程度に大きくならない。このため、ＳＮ比を上げるためにＭＲＩ画像に対して一定のマトリクスサイズでスムージング処理を行っていた。
【０００７】
従来、血流動態を解析するアルゴリズムには、ガンマ関数解析法がある。この方法は、ダイナミックスキャンによって得られ、ＣＴ値やＭＲ信号値の時間変化の情報を示す時間−濃度曲線に対して、ガンマ関数でカーブフィッティングを行い、フィッティング曲線の曲線下面積やピーク値から血流量を算出していた。
【０００８】
血流動態を解析する他のアルゴリズムには、最大勾配（maximum slope ）法がある。この方法は、時間−濃度曲線の最大勾配から血流量を算出していた。最大勾配法は、ガンマ関数解析法に比べて定量性が改善される。
【０００９】
血流動態を解析する更に他のアルゴリズムには、逆たたみこみ（deconvolution ）法がある。この方法は、動脈の時間−濃度曲線と各画素の時間−濃度曲線とをdeconvolution 演算することでインパルス応答を算出し、このインパルス応答の曲線下面積やピーク値から血流量を算出していた。deconvolution 法は、ガンマ関数解析法や最大勾配法に比べて、定量性や造影レートについての問題が改善される。
【００１０】
【発明が解決しようとする課題】
ところで、機能解析を行う器官について、断層像でのコントラストの付き具合は画素毎に異なるため、時間−濃度曲線のＳＮ比も画素毎に異なる。従来、断層像のすべての画素に対して同じマトリクスサイズでスムージングを行っていた。このため、時間−濃度曲線にてＳＮ比の低い画素でのノイズを打ち消すため、スムージングのマトリクスサイズを大きくすると、ＳＮ比の高い画素での分解能が必要以上に低下するという問題があった。
【００１１】
分解能が必要以上に低下するのを防止するためには、スムージングのマトリクスサイズを小さくすればよい。しかし、マトリクスサイズを小さくした場合には、時間−濃度曲線のＳＮ比が低い組織についてノイズを打ち消すことができないという問題があった。
【００１２】
機能解析する器官の断層像についてスムージングする際、最適なマトリクスサイズは撮影条件や撮影部位等により異なる。このため、最適なマトリクスサイズを選択するには、オペレータの経験と習熟とが必要とされるという問題があった。
【００１３】
最適なマトリクスサイズが選択された場合でも、最適にスムージングを行うための重み付けは、撮影対象の性質により異なる。そのため、最適にスムージングを行うための重み付けを決定するには、オペレータの経験と習熟とが必要とされるという問題があった。
【００１４】
deconvolution 演算は、大別して行列式を用いる方法と周波数空間に変換する方法とがある。周波数空間に変換する方法を用いた場合、周波数空間でのdeconvolution 演算は、高周波数成分を必要以上に強調するため、ノイズをできる限り排除する必要がある。このため、周波数空間でdeconvolution 演算する場合には、一般に高周波数領域を抑制するフィルタ（高周波数領域制御フィルタ）によりフィルタリングを行う。フィルタリングを行う際に、最適な高周波数領域制御フィルタの形状は、撮影条件や撮影部位等に依存して決定されるため、オペレータの経験と習熟とが必要とされるという問題があった。
【００１５】
本発明はこのような事情に鑑みてなされたもので、断層像から血流動態の情報を得る場合に、機能画像の分解能を必要以上に低下させることなく、かつ良好にノイズを除去することができる血流動態解析装置及び方法を提供することを目的とする。
【００１６】
【課題を解決するための手段】
請求項１に係る血流動態解析装置は、コンピュータ断層像撮影手段で撮影した被検体の断層像を入力する断層像入力手段と、前記断層像における各画素ごとに画素値の時間変化情報を取り出す手段と、各画素ごとに取り出した時間変化情報からノイズを除去するノイズ除去手段を有し、各画素ごとにノイズ除去した時間変化情報に基づいて血流動態を解析する演算手段と、からなる血流動態解析装置において、前記ノイズ除去手段は、各画素ごとに取り出された画素値の時間変化情報と、該画素値の時間変化情報をカーブフィッティングした関数との誤差を示す値であるパラメータを算出するパラメータ算出手段と、前記算出したパラメータに対応するマトリクスサイズを有する画像フィルタを選択するフィルタ選択手段と、を有し、画素ごとに画像フィルタを選択してスムージング処理を行うことを特徴としている。
請求項２に係る血流動態解析装置は、前記フィルタ選択手段は、前記パラメータと、前記画像フィルタとを対応付けて記憶するルックアップテーブルを有し、前記算出したパラメータに基づいて前記ルックアップテーブルから対応する画像フィルタを選択することを特徴とする。
請求項３に係る血流動態解析装置は、前記パラメータは、前記画素値の時間変化情報と、該画素値の時間変化情報をカーブフィッティングした関数との誤差の絶対値の標準偏差値であることを特徴とする。
【００１７】
請求項４に係る血流動態解析装置は、コンピュータ断層像撮影手段で撮影した被検体の断層像を入力する断層像入力手段と、前記断層像における各画素ごとに画素値の時間変化情報を取り出す手段と、各画素ごとに取り出した時間変化情報の周波数空間でのノイズを除去するノイズ除去手段を有し、各画素ごとにノイズ除去した時間変化情報に基づいて血流動態を解析する演算手段と、からなる血流動態解析装置において、前記ノイズ除去手段は、各画素ごとに取り出された時間変化情報に基づいて高周波数領域制御フィルタの形状を決定するフィルタ形状決定手段を有し、画素ごとに高周波数領域制御フィルタを決定してフィルタリング処理を行うことを特徴とする。
【００１８】
請求項５に記載の血流動態解析方法は、コンピュータ断層像撮影手段で撮影した被検体の断層像の入力を受付ける工程と、前記断層像における各画素ごとに画素値の時間変化情報を取り出す工程と、各画素ごとに取り出された画素値の時間変化情報と、該画素値の時間変化情報をカーブフィッティングした関数との誤差に基づく値であるパラメータを算出する工程と、前記算出したパラメータに対応するマトリクスサイズを有する画像フィルタを選択するフィルタ選択工程と、前記フィルタ選択工程により選択された画素ごとの画像フィルタに基づいてスムージング処理を行なうことにより、各画素ごとに取り出した時間変化情報からノイズを除去する工程と、前記ノイズを除去する工程により各画素ごとにノイズ除去した時間変化情報に基づいて血流動態を解析する演算工程と、を含むことを特徴とする。
また、請求項６に係る血流動態解析方法は、コンピュータ断層像撮影手段で撮影した被検体の断層像の入力を受付ける断層像入力工程と、前記断層像における各画素ごとに画素値の時間変化情報を取り出す工程と、各画素ごとに取り出された時間変化情報に基づいて高周波数領域制御フィルタの形状を決定するフィルタ形状決定工程と、画素ごとに高周波数領域制御フィルタを決定してフィルタリング処理を行うことにより各画素ごとに取り出した時間変化情報の周波数空間でのノイズを除去する工程と、前記ノイズを除去する工程により各画素ごとにノイズ除去した時間変化情報に基づいて血流動態を解析する演算工程と、を含むことを特徴とする。
【００１９】
本発明によれば、断層像から血流動態の解析を行う場合、入力画像から血流動態の情報を解析することにより、オペレータの経験と習熟がなくても、撮影条件や撮影部位等によって異なる最適な画像フィルタ（最適なマトリクスサイズと重み付け等を有するフィルタ）を、画素毎に自動で決定することができる。
【００２０】
また、周波数空間でノイズ除去（高周波数領域制御）を行う場合、画素ごとに入力画像から血流動態の情報を解析することにより、オペレータの経験と習熟がなくても、撮影条件や撮影部位等によって異なる最適な高周波数領域制御フィルタの形状を自動で決定することができる。
【００２１】
【発明の実施の形態】
以下添付図面に従って本発明に係る血流動態解析装置及び方法の好ましい実施の形態について詳説する。
【００２２】
図１は、本発明に係る血流動態解析装置のハードウェア構成例を示すブロック図である。同図に示すように、本発明に係る血流動態解析装置は、主として断層像撮影手段で撮影した断層像を入力する断層像入力装置１と、血流動態の解析などの諸演算を行う演算装置２とから構成され、機能画像の分解能の低下を防止する。
【００２３】
断層像入力装置１と演算装置２は必ずしも別個の装置である必要はなく、一体化していてもよい。また、断層像入力装置１と演算装置２は、必ずしも断層像撮影手段と別個の装置である必要はなく、一体化していてもよい。
【００２４】
演算装置２は、断層像入力装置１から断層像の取り込みを行うインターフェイス（Ｉ／Ｆ）３と、断層像や演算結果を一時的に格納するメモリ４と、血流動態の解析などの諸演算を行う中央処理装置（ＣＰＵ）５と、演算結果や機能画像を記録するハードディスクなどの記録媒体６と、マウスやキーボードなどの外部入力装置７と、処理結果を表示する表示装置８と、各構成要素を接続する共通バス１０とから構成される。
【００２５】
図２は、本発明に係る血流動態解析装置における断層像の入力から機能画像の表示までの処理を示すフローチャートである。図２に示すように、まず、断層像入力装置１から血流解析を行う被検体の断層像を選択し、選択された断層像を演算装置２に入力する（ステップＳ２０１）。
【００２６】
入力された断層像はメモリ４に一時格納されるか、あるいは記録媒体６に保存される。また、断層像入力装置１と演算装置２とが一体化して構成される場合には、血流解析を行う被検体の断層像を外部入力装置７により選択し、選択された断層像を記録媒体６から読み出し、メモリ４に一時格納する。続いて、中央処理装置５により格納された断層像から、画素値が時間変化する情報を表す時間−濃度曲線を画素毎に取り出す（ステップＳ２０２）。次に、取り出された時間−濃度曲線から、Ｐｅａｋ値やＳＤ値（standard deviation：標準偏差）などのマトリクスサイズ決定パラメータを算出する（ステップＳ２０３）。尚、マトリクスサイズ決定パラメータの算出方法については後述する。
【００２７】
次に、中央処理装置５により、算出したマトリクスサイズ決定パラメータに対応して、スムージング処理の際のマトリクスサイズを決定する（ステップＳ２０４）。任意のマトリクスサイズ決定パラメータに対応するマトリクスサイズを決定するため、予めルックアップテーブル等の対応表を作成しておく。尚、対応表の作成方法については後述する。
【００２８】
次に、各断層像における局所領域のＳＤ値等の重み付け決定パラメータを算出する (ステップＳ２０５) 。尚、重み付け決定パラメータの算出方法については後述する。次に、中央処理装置５により、算出した重み付け決定パラメータから、最適な重み付けを決定する (ステップＳ２０６) 。尚、重み付けの決定方法については後述する。
【００２９】
次に、決定したマトリクスサイズおよび重み付けの画像フィルタにより、画素毎にスムージング処理を行う（ステップＳ２０７）。ステップＳ２０２からステップＳ２０７までの処理を断層像の全画素に対して行うことにより、スムージング処理を施した断層像が得られる。次に、中央処理装置５により、スムージング処理が施された断層像から、画素毎に時間−濃度曲線を取り出す（ステップＳ２０８）。次に、中央処理装置５により、この時間−濃度曲線に対して画素毎に機能情報を算出する（ステップＳ２０９）。この機能情報を解析する際のアルゴリズムは、ガンマ関数解析法、最大勾配法、deconvolution 法等の公知のアルゴリズムでよい。ステップＳ２０８とステップＳ２０９の処理は、スムージング処理後の断層像の全画素に対して行う。次に、断層像の全画素について解析された機能情報をマッピングすることで機能画像が得られ（ステップＳ２１０）、この機能画像は表示装置８に表示される（ステップＳ２１１）。
【００３０】
次に、時間−濃度曲線のＰｅａｋ値やＳＤ値などのマトリクスサイズ決定パラメータの算出方法について説明する。図３（ａ）、（ｂ）は、時間−濃度曲線のＳＤ値を算出する方法の説明図である。図３（ａ）に示すように、断層像の任意の画素から時間−濃度曲線を取り出し、この時間−濃度曲線の最大値をＰｅａｋ値とする。続いて、取り出した時間−濃度曲線に対して、ガンマ関数等の公知の関数でカーブフィッティングを行う（図３（ｂ））。尚、カーブフィッティングを行う際のアルゴリズムは、最小二乗法等の公知のアルゴリズムで良い。また、フィッティング関数により得られた最大値をＰｅａｋ値としても良い。
【００３１】
ここで、濃度時間曲線をTDC(t)、フィッティング関数をfit(t)とすると、両者の誤差の絶対値g(t)は次式で表される。
【００３２】
【数１】
0004193168.tif 000002
g(t)のＳＤ値は次式で表される。
【００３３】
【数２】
0004193168.tif 000003
ここで、ｇａｖ（ｔ）は、ｇ（ｔ）の平均値を表す。[数２]式により算出されるＳＤ値を、時間−濃度曲線のＳＤ値と定義する。
【００３４】
ここでは、マトリクスサイズ決定パラメータの代表例として、Ｐｅａｋ値とＳＤ値とを取り上げたが、マトリクスサイズ決定パラメータはこれらに限定されない。例えば［数１］式に示す誤差の絶対値g(t)又は時間−濃度曲線の分散、あるいはＰｅａｋ値とＳＤ値との比など、時間−濃度曲線が解析精度を確保するのに十分なＳＮ比を持つように、スムージング処理を行う際のマトリクスサイズを決定するパラメータならば、いかなるものでもよい。
【００３５】
次に、マトリクスサイズ決定パラメータとマトリクスサイズとの対応表を作成する方法について図４（ａ）〜（ｅ）を用いて説明する。
【００３６】
まず、図４（ａ）に示すように、血流動態を解析する器官と同等の大きさの水ファントム４０１中に任意の径の管４０２を配置したファントム４００を用意する。このファントム４００外壁の断層像の形状は、血流動態を解析する器官の断層像の形状と類似していることが望ましい。例えば、人間の脳の血流動態を解析する場合には、直径１６０mm〜１９０mm程度の円筒形のファントム等が適切である。
【００３７】
ファントム４００内部の模擬血管を示す管４０２に、任意の造影条件（濃度、量、造影レートなど）で造影剤を流しながら任意の撮影条件（管電流、管電圧、スライス厚、スキャン時間など）で撮影し、断層像を作成する。この断層像から時間−濃度曲線を取得し、この時間−濃度曲線に対して算出したＳＤ値をａとする（図４（ｂ））。ここでＳＤ値ａが閾値ＴＨよりも大きい場合、Ｎ×Ｎのマトリクスでスムージング処理を行う（図４（ｃ））。Ｎは、スムージング処理後の断層像のＳＤ値ｂが閾値ＴＨ以下になるような最小の整数である。閾値ＴＨは、十分な解析精度が得られる範囲内で任意に設定すればよい。これにより、マトリクスサイズ決定パラメータとマトリクスサイズとの関係の一例が得られる。同様にして、設定可能な複数の造影条件、撮影条件でファントムを撮影することにより、様々なマトリクスサイズ決定パラメータにおけるマトリクスサイズの関係が得られる（図４（ｄ））。これらの結果から、マトリクスサイズはマトリクスサイズ決定パラメータの何らかの関数で表すことができる。よって、任意のマトリクスサイズ決定パラメータに対応するマトリクスサイズを求めることができ、図４（ｅ）に示す対応表を作成することができる。この対応表を用いれば、画素毎に最適なマトリクスサイズを自動決定することが可能である。尚、本実施の形態では、便宜上、ファントム４００内部に一本の管４０２のみを配置したが、管の径や本数は任意でよい。
【００３８】
次に、重み付け決定パラメータの算出方法および重み付けの決定方法について図５（ａ）〜（ｅ）を用いて説明する。
【００３９】
例えば、図５（ａ）に示すように、断層像５００に対してＮ×Ｎのマトリクスでスムージング処理を行う場合において、スムージングを行う画素の座標を（i,j ）とおく。ここで、i-(N-1)/2 ≦i ≦i ，j-(N-1)/2 ≦j ≦j の領域（図５（ｂ））の標準偏差をＳＤ１、スムージングの重み付けをｗ１とし、i ≦i ≦i+(N-1)/2 ，j-(N-1)/2 ≦j ≦j の領域（図５（ｃ））の標準偏差をＳＤ２、スムージングの重み付けをｗ２とし、i-(N-1)/2 ≦i ≦i ，j ≦j ≦j+(N-1)/2 の領域（図５（ｄ））の標準偏差をＳＤ３、スムージング処理の重み付けをｗ３とし、i ≦i ≦i+(N-1)/2 ，j ≦j ≦j+(N-1)/2 の領域（図５（ｅ））の標準偏差をＳＤ４、スムージング処理の重み付けをｗ４とする。ここで必要以上に分解能を低下させることなくスムージング処理を行うためには、標準偏差の小さい領域の重み付けを大きく、標準偏差の大きい領域の重み付けを小さく設定すればよく、例えば次式のように、スムージング処理の重み付けを標準偏差の逆数の比に等しくなるように設定すればよい。
【００４０】
【数３】
0004193168.tif 000004
［数３］式のように重み付けを算出することで、必要以上に分解能を低下させることなく、オペレータの経験や習熟に依らずにスムージング処理の重み付けを自動決定することが可能である。
【００４１】
また、重み付け方法は［数３］式に限らず、例えば最も標準偏差の小さい領域の重み付けを１、その他の領域の重み付けを０に決定する等の、必要以上に分解能が低下しないならば任意の重み付け方法を用いてもよい。
【００４２】
更に、本実施の形態では、便宜上、重み付け決定パラメータに標準偏差を用いたが、重み付け決定パラメータはこれに限らず、分散やスムージング処理を行う画素とマトリクス内のその他の画素との差の絶対値和等の、マトリクス内の画素値のばらつきを判定できるパラメータであれば任意のパラメータでよい。
【００４３】
更に、本実施の形態では、便宜上、標準偏差を算出する領域を図５（ｂ）から図５（ｅ）までの４つの領域に分割したが、分割の方向や形状、個数等の分割方法はこれに限らず、必要以上に分解能が低下しないならば任意の分割方法でよい。
【００４４】
次に、ステップＳ２０９にて機能情報を解析するためのアルゴリズムにdeconvolution 法を用いた場合に、高周波数領域制御フィルタの形状を自動決定する方法について図６（ａ）〜（ｃ）を用いて説明する。
【００４５】
前述のように、deconvolution 演算を周波数空間で行う場合、一般に、図７に示す形状の高周波数領域制御フィルタでフィルタリングすることにより、無用なノイズ強調を抑制する。例えば、対象組織への動脈入力関数、即ち流入動脈における時間−濃度曲線が、図６（ａ）に示すように実測されたとする。この時間−濃度曲線の平均通過時間（ＭＴＴ）は、次式で求められる。
【００４６】
【数４】
0004193168.tif 000005
ここで、ｔは時間、ｎは計測点数、Ca(t) は時間−濃度曲線を表す。流入動脈における時間−濃度曲線のピークタイムＰＴ（図６（ａ））と、［数４］式から求められる平均通過時間から、図６（ｂ）に示すようなピークタイムＰＴを中心に平均通過時間の幅を持つ理想的な矩形の入力関数ideal ＿Ca(t) が求められる。このideal ＿Ca(t) のフーリエ変換をＦ｛ideal ＿Ca(t) ｝とすると、図６（ｃ）に示すようになる。
【００４７】
このＦ｛ideal ＿Ca(t) ｝に対して、例えば、周波数０から２π／ｋまでの成分のみをそのまま残し、２π／ｋから４π／ｋまでの成分を減弱し、４π／ｋよりも高い周波数成分はカットするようにフィルタの形状を決定する。この高周波数領域制御フィルタは、次式で与えられる。
【００４８】
【数５】
0004193168.tif 000006
ここで、ωs は周波数制御開始周波数、ωe は周波数制御終了周波数であり、ideal ＿Ca(t) の幅、すなわちCa(t) の平均通過時間をｋとすれば、ωs ＝２π／ｋ，ωe ＝４π／ｋである。よって、流入動脈における時間−濃度曲線から、高周波領域制御フィルタの形状を決定することができる。
【００４９】
本実施の形態では、便宜上、ωs ＝２π／ｋ，ωe ＝４π／ｋとして説明したが、ωs とωe はこれに限らず、ωs ＝ｍπ／ｋ，ωe ＝ｎπ／ｋ（ｍ≦ｎ）として任意に設定してもよい。更に、本実施の形態では、流入動脈における時間−濃度曲線からωs とωe を算出したが、ωs とωe の算出方法はこれに限らず、図７のように最大周波数から任意に設定してもよい。
【００５０】
【発明の効果】
以上説明したように本発明に係る血流動態解析装置及び方法によれば、血流動態の情報を得る場合に、機能画像の分解能を必要以上に低下させることなく、画素毎に最適な画像フィルタを選択してスムージング処理を行なうことができる。
【００５１】
更に、周波数空間で高周波数領域の制御を行う場合に、最適な高周波数領域制御フィルタの形状を画素ごとに自動で決定してフィルタリング処理を行うことができる。
【図面の簡単な説明】
【図１】本発明に係る血流動態解析装置のハードウェア構成例を示すブロック図
【図２】本発明に係る血流動態解析装置における処理手順を示すフローチャート
【図３】時間−濃度曲線のＳＤ値を算出する方法を説明する図
【図４】マトリクスサイズ決定パラメータとマトリクスサイズとの対応表を作成する方法を説明する図
【図５】重み付け決定パラメータと重み付けの決定方法を説明する図
【図６】高周波数領域制御フィルタの作成方法を説明する図
【図７】高周波数領域制御フィルタの例を示す図
【符号の説明】
１…断層像入力装置、２…演算装置、３…インターフェイス（Ｉ／Ｆ）、４…メモリ、５…中央処理装置（ＣＰＵ）、６…記録媒体（ハードディスク）、７…外部入力装置、８…表示装置
【図１】
0004193168.tif 000007
【図２】
0004193168.tif 000008
【図３】
0004193168.tif 000009
【図４】
0004193168.tif 000010
【図５】
0004193168.tif 000011
【図６】
0004193168.tif 000012
【図７】
0004193168.tif 000013