JP4134187B2 - 拡散強調画像処理装置及び画像処理プログラム - Google Patents

Description

本発明は、被検体内の原子核スピンの磁気共鳴現象を利用した磁気共鳴イメージング（ＭＲＩ）に係り、特に、原子核スピンの位相拡散を強調した画像を得るディフュージョンイメージングに関する。

近年、静磁場中に置かれた被検体の原子核スピンをラーモア周波数の高周波信号で磁気的に励起し、この励起に伴って発生するＭＲ信号に基づいて画像を再構成したり、スペクトルデータを得るＭＲＩ装置にかかる技術が盛んに研究開発されている。

このＭＲＩ装置で実施される撮像法の中に、拡散強調撮像法（Diffusion Weighted Imaging）がある。この撮像法は、被検体の組織内部にある水分子の微細でランダムないわゆるブラウン運動の特性を特に強調して画像化するイメージング法であり、従来のＴ1強調画像・Ｔ２強調画像と比較して、水分子の状態が把握しやすいため、脳梗塞巣、脳神経線維や末梢神経繊維の走行や腫瘍の描画等に用いられている。

この拡散強調撮像法によって撮像された拡散強調画像（Diffusion Weighted Image）は、撮像される際における勾配磁場パルスの方向によって、出力される画像が歪むという特性を持っている。

このような問題に対し、例えば、非特許文献１には基準となる画像を用いて、拡散強調画像の歪みを補正する技術が開示されている。
"Eddy-Current Distortion Correction and Robust Tensor Estimation for MR Diffusion Imaging" J.F. Mangin, C. Poupon, C. Clark, D. Le Bihan, I. Bloch 10 Oct. 2002

また、拡散強調画像では撮像される際における勾配磁場パルスの方向によって、強調される水分子の拡散方向に偏りを生じることから、異方性を持って特徴が強調されるという特性がある。このため、従来は複数の勾配磁場パルスの方向を変えた拡散強調画像を撮影し、それを診断の時に並べて比較するということが行われた。

ところが、複数の勾配磁場パルスの方向を変えた拡散強調画像を並べて比較する場合は、この複数の画像の違いの把握は観察する人間が行う必要があった。ここで、異なる勾配磁場パルスの下で撮像された複数の拡散強調画像を重ね合わせる（フュージョン）ことにより、より精密な診断を行なうことが期待できる。

ところで、上述した拡散強調撮像法は、特に撮像時間が長いため患者の負担が大きい。そのため、できる限り患者の負担にならないよう撮像時間の短縮が要求されるところ、上述した非特許文献１に記載の手法では、歪み補正時に歪みの含まれていない参照用の画像を必要とする。その為に、画像種の異なる基準画像を複数撮像しなければならず、撮像時間の長時間化を招くという問題がある。

更に、非特許文献１に記載の手法であっては、画像に発生する剪断歪みを過度に単純化しているので不正確である。

そこで、本発明は上記課題に鑑みてなされたものであり、拡散強調画像の歪みを補正する拡散強調画像処理法において、被検体への負担を強いることなく、拡散強調画像の歪みを補正することができる拡散強調画像処理装置等を提供することができる。

上記課題を解決すべく、請求項１に記載の発明は、撮像された際の勾配磁場パルス方向がそれぞれ異なる複数の拡散強調画像を、撮像された際の前記勾配磁場パルス方向に沿って複数の帯状領域に分割する第１分割手段と、前記各拡散強調画像を互いに参照し、かつ前記各帯状領域に基づいて前記勾配磁場パルスによって発生する前記勾配磁場パルス方向の歪みを相対的に補正する補正手段を有することを特徴とする。

これによれば、複数の拡散強調画像を重ね合わせる際に、当該複数の拡散強調画像が勾配磁場パルスによって歪んでいる場合であっても、該勾配磁場パルス方向に沿って分割された帯状領域に基づいて歪みを補正することができるため、被検体への負担を強いることなく、かつ、精度よく歪みを補正することができる。

上記課題を解決すべく、請求項７に記載の発明は、コンピュータを、撮像された際の勾配磁場パルス方向がそれぞれ異なる複数の拡散強調画像を、撮像された際の前記勾配磁場パルス方向に沿って複数の帯状領域に分割する第１分割手段、前記各拡散強調画像を互いに参照し、かつ前記各帯状領域に基づいて前記勾配磁場パルスによって発生する前記勾配磁場パルス方向の歪みを相対的に補正する補正手段として機能させることを特徴とする。

これによれば、複数の拡散強調画像を重ね合わせる際に、当該複数の拡散強調画像が勾配磁場パルスによって歪んでいる場合であっても、該勾配磁場パルス方向に沿って分割された帯状領域に基づいて歪みを補正することができるため、被検体への負担を強いることなく、かつ、精度よく歪みを補正することができる。

本発明によれば、被検体への負担を強いることなく、拡散強調画像の歪みを補正することができる拡散強調画像処理装置等を提供することを目的とする。

以下、本発明の好適な実施の形態を添付図面に基づいて説明する。

本実施形態にかかる拡散強調画像処理装置は、寝台上などに寝かされた被検体から発生するＭＲ信号を検出する拡散強調撮像法によって取得された勾配磁場パルスがそれぞれ異なる複数の拡散強調画像（例えば、２枚〜ｎ枚）を装置内部に取り込んで、これら複数の拡散強調画像を互いに参照して各拡散強調画像の歪みを相対的に補正した後に、重ね合わせ（フュージョン）て合成画像を作成するよう構成される。

図１は、本発明に係る拡散強調画像処理装置の構成を示すブロック図である。同図に示す如く拡散強調画像処理装置１は、演算機能を有するコンピュータとしてのＣＰＵ（Central Processing Unit）、作業用ＲＡＭ（Random Access Memory）、各種データ及びプログラムを記憶するＲＯＭ（Read Only Memory）等から構成された制御部１１と、図示しない撮像装置に対して各種制御信号を送信すると共に、当該撮像装置にて撮像された拡散強調画像（以下の説明において、「拡散強調画像」は、単に「画像」と言う。）を画像データとして装置内部に取り込むための外部機器接続部１２と、当該外部機器接続部１２を介して装置内部に取り込んだ画像を記憶するための画像記憶部１３Ａを含んで構成された記憶部１３、拡散強調画像処理装置１によって作成された合成画像等を表示する表示装置としての表示部１４、各画像を補正する際に各種指示を与えたり、或いは当該拡散強調画像処理装置１を介して図示しない撮像装置に対して各種入力指示を行なう入力部１５を備えて構成され、制御部１１、外部機器接続部１２、記憶部１３、表示部１４及び入力部１５は、バス１６を介して相互に接続されている。

制御部１１は、図示しないＣＰＵ、作業用ＲＡＭ、本発明の画像処理プログラム等を含む各種制御プログラムやデータ等を記憶するＲＯＭ及び発振回路等を備えて構成されており、図示しない操作部からの操作信号に基づいて、当該操作信号に含まれている操作情報に対応する動作を実現すべく上記各構成部材を制御するための制御情報を生成し、バス１６を介して当該制御情報を該当する構成部材に出力して当該各構成部材の動作を統轄制御する。

また、制御部１１は、ＲＯＭ等に記憶された画像処理プログラムを実行することにより、他の構成部材と協動して本発明の補正手段、合成手段、第１分割手段、検索手段、第２分割手段、特徴点決定手段、対応点取得手段、特徴点歪み量取得手段として機能するようになっている。

また、外部機器接続部１２は、シリアル方式、ＵＳＢ方式、ＩＥＥＥ１３９４、ＬＡＮ方式、或いはその他の適宜な方式で図示しない撮像装置へ指示信号を送出する。

なお、本実施形態では表示部１４を拡散強調画像処理装置１内に内蔵したが、外部接続したモニタ等を表示装置として用いてもよい。この場合には、拡散強調画像処理装置１に内蔵されたビデオカード、及び、ＶＧＡケーブル、ＤＶＩケーブル、ＢＮＣケーブルなどを介してモニタ等を表示装置へ表示制御指示信号を送出するよう構成する。

［画像処理手順］
続いて、図を参照して複数の画像の重ね合わせ（フュージョン）にかかる画像処理の具体的手法について説明する。

本発明に関わる手法は、画像の歪みが勾配磁場パルスの方向によって、一方向に平行移動、剪断及び拡大縮小するような歪みであるという特徴を利用する。詳しくは、その一方向に沿った直線を規定した場合に、歪みによって直線は変化しないことを利用する。また、直線上の点は歪みによって直線上を移動するが、直線は全体として平行移動及び拡大縮小のみするので、移動後の直線上の点が２点が定まれば他の点に関しても一意に定まることを利用する。なお、剪断は複数の一方向に沿った直線がそれぞれ異なる歪み方をすることにより表現できる。

１.画像処理の実施形態１
図２（Ａ）は、画像処理の実施形態１にかかる画像処理を示すフローチャート、（Ｂ）は、ステップＳ３における補正処理を示すフローチャートである。なお、当該処理は拡散強調画像処理装置１の制御部１１が画像処理プログラムを実行することによって行われる。

先ず、制御部１１は、画像記憶部１３Ａから重ね合わせ（フュージョン）を所望する画像ｉ、画像ｉ２を取得する（ステップＳ１）。例えば、画像記憶部１３Ａに記憶されている複数の画像を表示部１４に表示させて、医師等の操作者が入力部１５を操作することによって、表示された画像の中から、重ね合わせ（フュージョン）を所望する画像ｉ、画像ｉ２を指定するよう構成する。

続いて、取得した画像ｉ、画像ｉ２の夫々から対象物体の輪郭Ｃ［ｉ］及びＣ［ｉ２］を抽出する。輪郭の抽出には領域拡張法(Region Grow)やLevel Set Method等の既存の方法を用いる事ができる（ステップＳ２）。

そして、制御部１１は補正手段として機能し、画像ｉ、画像ｉ２を相互に参照して各画像に対して補正処理を行なう（ステップＳ３）。

続いて、図２（Ｂ）の補正処理を示すフローチャートに移行し、ステップＳ３の補正処理について説明する。なお、ここでは画像ｉ２を用いて画像ｉを補正し、補正後の画像ｉ’を取得する場合の例に従って説明する。

先ず、画像ｉ上の対象物体の輪郭Ｃ１に基づいて、特徴点Ｚ１［ｉ］、Ｚ２［ｉ］、Ｚ３［ｉ］、、、を決定する（ステップＳ３０１）。

ここで、図３（Ａ）に示す特徴点の決定の説明図を用いて具体的に説明する。

先ず、制御部１１が第１分割手段として機能し、画像ｉを撮像の際における勾配磁場パルスの方向に沿って複数の帯状領域に分割する。そして、制御部１１が特徴点決定手段として機能し、各帯状領域に含まれる輪郭Ｃ［ｉ］線上の点からそれそれ特徴点Ｚを２点ずつ決定していく。従って、帯状領域の数の２倍の数の特徴点Ｚが決定されることとなる。決定の際には、制御部１１が、対象物体の輪郭が急激に湾曲していたり、尖鋭形状である、或いは曲率が高い、曲率の変動が大きいなど特徴のある形状を成す箇所を選択する。

続いて、制御部１１は対応点取得手段として機能し、決定された各特徴点Ｚに対応する画像ｉ２上の対応点を検索し、取得する（ステップＳ３０２）。

この際、特徴点Ｚ１［ｉ］、Ｚ２［ｉ］、Ｚ３［ｉ］、、、の周囲の領域（以下、「局所領域ＲＯＩ」と言う。）と最も類似する場所（対応領域）を探すことにより、特徴点Ｚ１［ｉ］、Ｚ２［ｉ］、Ｚ３［ｉ］、、、に対応する画像ｉ２上の対応点Ｕ１［ｉ２］、Ｕ２［ｉ２］、Ｕ３［ｉ２］、、、を検索する。

検索する際には、特徴点Ｚ１［ｉ］上の輪郭線の向きや曲率を、対応すると考えられるＵ１［ｉ２］上の輪郭線の向きや曲率と比較すればよい。また、特徴点Ｚ１［ｉ］の周辺の輪郭線の形状を、対応すると考えられるＵ１［ｉ２］周辺の輪郭線の形状とパターンマッチングしても良い。尚、パターンマッチングとして画像の類似性を評価しても良い。また、フーリエ変換やウェーブレット変換を行った上でパターンマッチングを行っても良い。

このように、全体の歪みを無視できる程度局所的に画像を捉えることで、２つの画像上の夫々対応する点同士を結びつけることができるのである。なお、検索する際には特徴点Ｚの有する位置座標と同一の位置座標を有する画像ｉ２上の場所から検索を始めると効率がよい。

続いて、検索された画像ｉ２上の点Ｕ１［ｉ２］、Ｕ２［ｉ２］、Ｕ３［ｉ２］、、、を基準として、対応する画像ｉ上の特徴点Ｚ１［ｉ］、Ｚ２［ｉ］、Ｚ３［ｉ］、、、を補正する（ステップＳ３０３）。

図３（Ｂ）に示すように、画像ｉ上の特徴点Ｚ［ｉ］の補正後の位置Ｚｆ［ｉ］は、各画像ｉ、ｉ２の勾配磁場パルスの方向に従って求めることができる。

すなわち、画像ｉ上の特徴点Ｚ［ｉ］の補正後の位置Ｚｆ［ｉ］は、画像ｉ上の特徴点Ｚ［ｉ］を通り、画像ｉの勾配磁場パルスの方向に平行な直線Ｌ［ｉ］と、当該特徴点Ｚ［ｉ］に対応する画像ｉ２上の点Ｕ［ｉ２］を通り、画像ｉ２の勾配磁場パルスの方向に平行な直線Ｌ［ｉ２］の交点として求めることができる。

このように、画像ｉ２上の点Ｕ１［ｉ２］、Ｕ２［ｉ２］、Ｕ３［ｉ２］、、、を基準として、対応する画像ｉ上の特徴点Ｚ１［ｉ］、Ｚ２［ｉ］、Ｚ３［ｉ］、、、を補正する。

これによって、画像ｉの勾配磁場パルスに沿って分割した帯状領域毎に２つずつ決定した特徴点Ｚ１［ｉ］、Ｚ２［ｉ］、Ｚ３［ｉ］、、、を、画像ｉ２を用いて補正して、補正後の特徴点Ｚｆ１［ｉ］、Ｚｆ２［ｉ］、Ｚｆ３［ｉ］、、、に基づいて、画像全体の歪みを補正する（ステップＳ３０４）。

画像ｉの歪みは勾配磁場パルスの方向に沿っていることから、勾配磁場パルスに沿って分割した帯状領域上の特徴点であって、補正後の２つの特徴点に基づいて、各帯状領域における歪みベクトルが、一様に分布するよう帯状領域の特徴点以外の箇所を補間演算することにより、画像ｉ全体の歪みを補正し、補正後の画像を画像ｉ’として取得することができる。

そして、同様にしてステップＳ３０５からステップＳ３０８の処理を行なって、画像ｉを用いて画像ｉ２を補正し、補正後の画像ｉ２’を取得する。なお、ステップＳ３０５からステップＳ３０８の処理は上述したステップＳ３０１乃至ステップＳ３０４の処理と同様であるので、説明を省略する。

そして、制御部１１は合成手段として機能し、図２（Ａ）のフローチャートのステップＳ４に移行して、補正後の画像ｉ’と画像ｉ２’とを夫々重ね合わせ（フュージョン）て、合成画像を作成して処理を終了する（ステップＳ４）。

なお、上述した実施形態では、画像の補正処理を行なう際に、各帯状領域中に２つの特徴点を決定し、他の箇所は補間演算することにより全体の歪みを補正したが、本発明はこれに限らず、２つ以上であれば特徴点は何点あってもよい。

また、特徴に乏しい帯状領域であっては近傍の帯状領域より補正量を補間して求めることもできる。それは、画像の歪みは勾配磁場パルスの方向に沿っていることが自明であることから、少なくとも２点あればステップＳ３０４で全体の歪みを補正することができるからである。

以上、説明した実施形態によれば、被検体への負担を強いることなく、勾配磁場パルスによって発生する画像（拡散強調画像）の歪みを補正して所望の数の画像（拡散強調画像）からなる合成画像を取得することができる。この際、各画像（拡散強調画像）の撮像時における勾配磁場パルスの方向が異なる場合であっても、勾配磁場パルスの方向の画像の歪みを確実に補正して、各画像の勾配磁場パルスの方向に起因する異方性を加味した合成画像を作成することができる。

また、当該実施形態の変形例として、制御部１１が第２分割手段として機能し、画像ｉ、ｉ２を夫々帯状領域に分割した後に、更に各帯状領域を複数の分割領域に分割して、夫々の分割領域毎に同一性を判断すべくパターンマッチングを行なうことにより、各分割領域単位で歪みを補正するよう構成してもよい。

制御部１１が検索手段として機能し、各画像の各帯状領域に含まれる少なくとも２つの分割領域について、他の画像上で対応する対応領域を夫々検索して、上記実施形態１と同様に、画像ｉ、ｉ２の勾配磁場パルスの方向に基づいて各分割領域について補正を行なってから、補正後の各分割領域を用いて、当該分割領域が属する帯状領域の他の領域についても補正を行なうよう構成すればよい。

分割領域の大きさは特に限定されないが、例えば、画像ｉ、画像ｉ２が夫々２５６×２５６ピクセル程度の大きさである場合には、各画像を３２×３２ピクセル程度の分割領域に分割すればよい。

また、各画像に対してメジアンフィルタリング処理及びガウシャンフィルタリング処理を施して、特徴のある分割領域のみを選抜してパターンマッチングを実行するよう構成してもよい。これにより、パターンマッチングの処理演算時間を短縮することができる。なお、複数の分割領域のうち、例えば対象物体の平坦な領域でなく、輪郭部分等の特徴のある箇所を有する分割領域のみを選択して、選択された分割領域のみパターンマッチングを行なうよう構成してもよい。

なお、上記実施形態では、２つの画像を相互に参照して各画像を補正し、補正後の画像の合成画像を作成する場合について説明したが、３以上の画像に対しても適用可能である。

２.画像処理の実施形態２
次に、画像処理の実施形態２について説明する。本実施形態は、輪郭抽出を用いず処理する方法を示す。すなわち、実施形態1では輪郭抽出の結果得られた輪郭の形状を用いて対応点を検索したが、本実施例であっては画像を小領域に分解したそれぞれの領域に対応する領域をパターンマッチングで検索する。本手法は輪郭抽出が難しい画像や、輪郭形状に特徴の少ない画像に対して効果的である。また、画像処理の実施形態１では２枚の画像を合成したが、本実施形態では、３以上の画像を合成する例を示す。

図４は、画像処理の実施形態２にかかる画像処理を示すフローチャートである。なお、当該処理は拡散強調画像処理装置１の制御部１１が画像処理プログラムを実行することによって行われる。

先ず、制御部１１は、画像記憶部１３Ａから重ね合わせ（フュージョン）を所望するｎ枚の画像を取得する（ステップＳ１１）。例えば、画像記憶部１３Ａに記憶されている複数の画像を表示部１４に表示させて、医師等の操作者が入力部１５を操作することによって、表示された画像の中から、重ね合わせ（フュージョン）を所望するｎ枚の画像を指定したり、或いは、同時に撮像されたｎ枚の画像の全てを指定して取得する。

続いて、画像について歪み量を表現する配列Ｗａｒｐ［p,q］［ｉ］を制御部１１の図示しないＲＡＭ等に確保する（ステップＳ１２）。なお、ｉは画像を一意に識別する固有の画像IDであって、以下、画像を識別するために変数ｉ（０≦ｉ＜ｎ）として用いる。例えば、ステップＳ１１において画像［ｉ］（ｉ＝0、１）の４枚の画像を取得した場合（ｎ＝４）には、各画像毎に歪み量を表現する配列Ｗａｒｐ［p,q］［0］、配列Ｗａｒｐ［p,q］［１］、配列Ｗａｒｐ［p,q］［２］、配列Ｗａｒｐ［p,q］［３］を確保する。また、配列Ｗａｒｐ［p,q］［ｉ］のそれぞれの要素Ｗａｒｐ［p,q］もまたp,q を添え字とした２次元配列であり、p,qは画像上の当該帯状領域を格子形状に区切った領域を識別するための値である。そしてＷａｒｐ［p,q］の値はその示す格子の歪み量を表す。ここで、本発明の場合は歪み量は勾配磁場パルス方向に拘束されるので歪み量は方向を含んだベクトルではなく勾配磁場パルス方向のスカラー値として表現できる。

そして、各画像の撮像条件であり、補正すべき歪みにも起因する勾配磁場パルスの方向Ａｎｇｌｅ［ｉ］を取得する（ステップＳ１３）。

次に、変数ｉを０として初期化し（ステップＳ１４）、画像［ｉ］のその他の各画像に対する「歪み量Ｗａｒｐ［ｉ］取得処理」を行なう（ステップＳ１５）。すなわち、ある画像［ｉ］を（ｎ−１）枚の他の画像と重ね合わせる前に、他の画像を参照して当該他の画像に対する相対的な歪み量を取得して、当該歪み量に基づいて画像［ｉ］を補正する。

画像全体を補正する前に、ステップＳ１５にて取得した歪み量Ｗａｒｐ［ｉ］の異常検出（ステップＳ１６）、歪み量Ｗａｒｐ［ｉ］の補間（ステップＳ１７）を行ない、取得した歪み量Ｗａｒｐ［ｉ］の不自然な部分を異常個所として除去し、また、取得した歪み量Ｗａｒｐ［ｉ］を用いて勾配磁場パルスの方向に補間演算処理を施す。

続いて、変数ｉを1増加させ（ステップＳ１８）、変数ｉがステップＳ１１にて取得した画像の数ｎよりも多いか否かを判定し（ステップＳ１９）、取得したｎ枚の画像全てに対してステップＳ１５乃至ステップＳ１７の処理を行ない各画像の他の画像に対する歪み量Ｗａｒｐ［ｉ］（０≦ｉ＜ｎ）を取得した後に、制御部１１が合成手段として機能し、「ｎ枚の画像の合成処理」を行ない処理を終了する（ステップＳ２０）。

なお、ステップＳ１５における「歪み量Ｗａｒｐ［ｉ］取得処理」、及びステップＳ１６、ステップＳ１７の異常検出及び補間、ステップＳ２０における「ｎ枚の画像の合成処理」については後にフローチャート等の図面を用いて詳細に説明する。

２−１．歪み量Ｗａｒｐ［ｉ］取得処理
続いて、上述したステップＳ１５における「歪み量Ｗａｒｐ［ｉ］取得処理」について、図５に示すフローチャートを用いて説明する。ここでは、制御部１１が第１分割手段、特徴点決定手段、対応点取得手段、特徴点歪み量取得手段として機能することとなる。

先ず、画像［ｉ］の大きさをｘ_ｍａｘ、ｙ_ｍａｘ、画像［ｉ］を勾配磁場パルスの方向に沿って複数の帯状領域に分割して、更に当該帯状領域を格子形状に分割する際の格子幅をＧＷ、画像［ｉ］を局所領域毎に対応領域を検索すべく、パターンマッチングを行なう際に当該マッチングに対する寄与度の最小許容値を最小特徴量ＣＶ_ＭＩＮとして定義する（ステップＳ３０）。

なお、画像［ｉ］の大きさｘ_ｍａｘ、ｙ_ｍａｘは、特に限定されないが、例えば２５６×２５６ピクセル程度とし、また、格子幅ＧＷも特に限定されないが、例えば、６４×６４ピクセル程度とする。また、以下の説明において、画像［ｉ］の勾配磁場パルスの方向は画像のｙ軸方向と一致するものとする。

そして、画像［ｉ］上の格子形状領域を特定する変数ｐ及びｑを夫々０として初期化し（ステップＳ３１、Ｓ３２）、（p*ｇｗ、ｑ*ｇｗ）の周辺領域を局所領域ＲＯＩ［p,q］として設定する（ステップＳ３３）。変数ｐ及びｑは、画像［ｉ］の大きさが２５６×２５６ピクセル、格子幅ＧＷが６４×６４ピクセルとすると、格子形状領域は４×４＝１６領域となり、変数ｐ及びｑは夫々０から３の値となる。

次に、当該「局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］［ｉ］算出処理」を行ない（ステップＳ３４）、後に説明するように画像［ｉ］上の領域を特定する変数ｐ及びｑを徐々に増加させながら画像［ｉ］上の局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］［ｉ］を順次求めていく。

図６に示す例を用いて説明すると、画像［ｉ］の左下隅の点（ｐ、ｑ）を有する点を中心とする局所領域ＲＯＩ［p,q］（一点鎖線枠にて図示。）の特徴量ＣＶ［p,q］［ｉ］が算出される。しかし、この場所には対象物体が映し出されていないことから、パターンマッチングに対する寄与度が無く、特徴量ＣＶ［p,q］［ｉ］として非常に小さな値が算出されることとなる。

特徴量ＣＶは、後のパターンマッチングにおいて、画像の特徴を表す箇所を構成する有効な画素が含まれているかを表すパラメータであって、このような対象物体でない箇所（空気）や、対象物体であっても平坦な箇所など、画像の特徴を表していない個所は、パターンマッチングに対する寄与度が低く、比較的小さな値の特徴量ＣＶが算出されることとなる。

一方、図６中点線枠にて示した局所領域ＲＯＩ［p,q］のように、対象物体の輪郭に相当する箇所など、画像の特徴を表している箇所は、比較的大きな値の特徴量ＣＶが算出されることとなる。なお、当該「局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］算出処理」については後にフローチャートを用いて具体的に説明する。

続いて、領域を特定する変数ｑを１増加させ（ステップＳ３５）、変数ｑが、ステップＳ３０にて定義した画像の大きさｙ_ｍａｘを格子幅ＧＷで除算した値よりも小さいか否かが判定され（ステップＳ３６）、小さい場合（ステップＳ３６：Ｙｅｓ）には、ステップＳ３３に移行して再び画像［ｉ］上の新たな局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］［ｉ］を算出しつづける。

一方、ステップＳ３６の判定において領域を特定する変数ｐが、ステップＳ３０にて定義した画像の大きさｘ_ｍａｘを格子幅ＧＷで除算した値以上である場合（ステップＳ３６：Ｎｏ）には、ステップＳ３７に移行する。

続いて、領域を特定する変数ｐを１増加させ（ステップＳ３７）、変数ｐが、ステップＳ３０にて定義した画像の大きさｘ_ｍａｘを格子幅ＧＷで除算した値よりも小さいか否かが判定され（ステップＳ３８）、小さい場合（ステップＳ３８：Ｙｅｓ）には、ステップＳ３２に移行して再び画像［ｉ］上の新たな局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］［ｉ］を算出しつづける。

一方、ステップＳ３８の判定において、変数ｐが、ステップＳ３０にて定義した画像の大きさｘ_ｍａｘを格子幅ＧＷで除算した値以上である場合（ステップＳ３８：Ｎｏ）には、ステップＳ３９に移行する。

そして、領域を特定する変数ｐを０として初期化し（ステップＳ３９）、最大特徴量ＣＶ_ＭＡＸ及び第２位特徴量ＣＶ_ＳｅｃｏｎｄをステップＳ３０にて定義した最小特徴量ＣＶ_ＭＩＮとし、最大特徴量ＣＶ_ＭＡＸ及び第２位特徴量ＣＶ_Ｓｅｃｏｎｄを有する点の変数ｑ（以下、夫々「ＣＶ_ＭＡＸＱ」「ＣＶ_ＳｅｃｏｎｄＱ」と言う。）を夫々「−１」として定義する（ステップＳ４０）。

最大特徴量ＣＶ_ＭＡＸ及び第２位特徴量ＣＶ_Ｓｅｃｏｎｄとは、ステップＳ３４の「局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］算出処理」において算出された特徴量ＣＶ［p,q］の最大値と、これに続く第２位の値である。そして、これら最大特徴量ＣＶ_ＭＡＸ及び第２位特徴量ＣＶ_Ｓｅｃｏｎｄを有する点が特徴点Ｚ１及びＺ２となる。

続いて、領域を特定する変数ｑを０として初期化し（ステップＳ４１）、画像［ｉ］上の局所領域ＲＯＩ［p,q］における特徴量ＣＶ［p,q］［ｉ］が最大特徴量ＣＶ_ＭＡＸ以上であるか否かを判定し（ステップＳ４２）、最大特徴量ＣＶ_ＭＡＸ以上である場合（ステップＳ４２：Ｙｅｓ）には、ステップＳ４３へ、最大特徴量ＣＶ_ＭＡＸより小さい場合（ステップＳ４２：Ｎｏ）には、ステップＳ４４へ夫々移行する。

そして、ステップＳ４３では、第２位特徴量ＣＶ_Ｓｅｃｏｎｄを最大特徴量ＣＶ_ＭＡＸとし、ＣＶ_ＳｅｃｏｎｄＱをＣＶ_ＭＡＸＱとし、最大特徴量ＣＶ_ＭＡＸを画像［ｉ］上の局所領域ＲＯＩ［p,q］における特徴量ＣＶ［p,q］［ｉ］とし、ＣＶ_ＭＡＸＱを現在注目している領域を特定する変数ｑとして定義する（ステップＳ４３）。

ステップＳ４４では、領域を特定する変数ｑを１増加させた（ステップＳ４４）後に、変数ｑが、ステップＳ３０にて定義した画像の大きさｙ_ｍａｘを格子幅ＧＷで除算した値よりも小さいか否かが判定され（ステップＳ４５）、小さい場合（ステップＳ４５：Ｙｅｓ）には、ステップＳ４２に移行して画像［ｉ］上の領域（ｐ、ｑ）における特徴量ＣＶ［p,q］［ｉ］が最大特徴量ＣＶ_ＭＡＸ以下であるか否かが判定される。

一方、ステップＳ４５の判定において、変数ｑが、ステップＳ３０にて定義した画像の大きさｙ_ｍａｘを格子幅ＧＷで除算した値以上である場合（ステップＳ４５：Ｎｏ）には、ステップＳ４６に移行して、最大特徴量ＣＶ_ＭＡＸが最小特徴量ＣＶ_ＭＩＮより小さいか否かが判定される（ステップＳ４６）。

判定の結果、最大特徴量ＣＶ_ＭＡＸが最小特徴量ＣＶ_ＭＩＮと等しい場合（ステップＳ４６：Ｙｅｓ）には、当該特徴点Ｚ１は許容範囲よりも小さいので、歪み量Ｗａｒｐを取得することなくステップＳ５０へ移行する。なお、これに対応する特徴点Ｚ２は、その第２位特徴量ＣＶ_Ｓｅｃｏｎｄが最大特徴量ＣＶ_ＭＡＸ以下なことは明らかなので、ステップＳ４６の判定を行うは必要はない。一方、最大特徴量ＣＶ_ＭＡＸが最小特徴量ＣＶ_ＭＩＮと等しくない場合（ステップＳ４６：Ｎｏ）、すなわち、最大特徴量ＣＶ_ＭＡＸが最小特徴量ＣＶ_ＭＩＮより大きい場合には、「点Ｚ１（p、ＣＶ_ＭＡＸＱ）における各画像に対する歪み量Ｗａｒｐ［p、ＣＶ_ＭＡＸＱ］［ｉ］の取得処理」を行なって（ステップＳ４７）、特徴点Ｚ１の歪み量を取得する。

そして、第２位特徴量ＣＶ_Ｓｅｃｏｎｄが最小特徴量ＣＶ_ＭＩＮと等しいか否かが判定され（ステップＳ４８）、第２位特徴量ＣＶ_Ｓｅｃｏｎｄが最小特徴量ＣＶ_ＭＩＮと等しい場合（ステップＳ４８：Ｙｅｓ）には、当該特徴点Ｚ２は許容範囲よりも小さいので、歪み量Ｗａｒｐを取得することなくステップＳ５０へ移行する。

一方、第２位特徴量ＣＶ_Ｓｅｃｏｎｄが最小特徴量ＣＶ_ＭＩＮと等しくない場合（ステップＳ４８：Ｎｏ）、すなわち、第２位特徴量ＣＶ_Ｓｅｃｏｎｄが最小特徴量ＣＶ_ＭＩＮより大きい場合には、「点Ｚ２（ｐ、ＣＶ_ＳｅｃｏｎｄＱ）における各画像に対する歪み量Ｗａｒｐ［ｐ、ＣＶ_ＳｅｃｏｎｄＱ］［ｉ］の取得処理」を行なって（ステップＳ４９）、特徴点Ｚ２の歪み量を取得する。なお、ステップＳ４７及びＳ４９における「点Ｚ１及びＺ２の歪み量Ｗａｒｐ取得処理」は、後にフローチャートを用いて詳細に説明する。

そして、領域を特定する変数ｐを１増加させた（ステップＳ５０）後に、変数ｐが、ステップＳ３０にて定義した画像の大きさｘ_ｍａｘを格子幅ＧＷで除算した値よりも小さいか否かが判定され（ステップＳ５１）、小さい場合（ステップＳ５１：Ｙｅｓ）には、ステップＳ４０に移行して再び特徴点Ｚ１及びＺ２の歪み量Ｗａｒｐを取得する処理が行なわれる。

一方、ステップＳ５１の判定において、変数ｐが、ステップＳ３０にて定義した画像の大きさｘ_ｍａｘを格子幅ＧＷで除算した値以上である場合（ステップＳ５１：Ｎｏ）には、各帯状領域における特徴点Ｚ１及びＺ２であって、最小特徴量ＣＶ_ＭＩＮより大きい特徴量を有する特徴点Ｚ１及びＺ２の歪み量Ｗａｒｐを求めることができたので、処理を終了する。

２−１−１「局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］算出処理」
続いて、ステップＳ３４における局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］算出処理について図７を参照して説明する。図７（Ａ）は局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］［ｉ］算出処理を示すフローチャートである。

先ず、局所領域ＲＯＩ［p,q］に対してメジアンフィルター処理及びガウシャンフィルター処理を行なう（ステップＳ６１、Ｓ６２）。図７（Ｂ１）は例えば臓器及びその周辺の組織を含むために対象物体の輪郭部分など比較的特徴を有する局所領域ＲＯＩ［p,q］の各フィルター処理過程を示す説明図であり、（Ｂ２）及び（Ｂ３）は、特徴部分がない局所領域ＲＯＩ［p,q］の各フィルター処理過程を示す説明図であって、（Ｂ２）は例えば全体的に一様な筋肉組織であるために対象物体の平坦な領域等であり、（Ｂ３）は例えば対象物体でない箇所（空気）等である。

同図に示すように、局所領域ＲＯＩ［p,q］に対してメジアンフィルター処理を施すことにより、雑音を除去し、その後にガウシャンフィルター処理を施すことで、繊維などの細かい構造物を除去し、局所領域ＲＯＩ［p,q］を分散値取得に適した画像データとすることができる。その為に単純にフィルターを用いずに分散値を取得した場合と比較して、筋繊維を含むために細かい構造物が存在する領域で不当に高い特徴量を示してしまうことを防ぐことができる。

そして、局所領域ＲＯＩ［p,q］の分散値を算出し、これを特徴量ＣＶ［p,q］［ｉ］として処理を終了する（ステップＳ６３）。

２−１−２「特徴点Ｚの各画像に対する歪み量Ｗａｒｐ［ｚ］取得処理」
続いて、ステップＳ４７及びステップＳ４９における、特徴点Ｚ１及びＺ２の各画像に対する歪み量Ｗａｒｐ［ｚ］［ｉ］取得処理について説明する。

はじめに、特徴点について図８及び図９を用いて説明する。図８は拡散強調画像の歪みの特徴を示す説明図、図９は特徴点の説明図である。

図８（Ａ）は、ステップＳ１２において確保された歪み量を表現する配列Ｗａｒｐ［p,q］［ｉ］を模式的に示し、格子幅ＧＷで区切られた状態を示す。

歪みは勾配磁場パルスの方向に沿って起こることから、例えば、図８（Ｂ１）に図示するように画像［ｉ］のある部分が勾配磁場パルスの方向に沿って並行移動したり、図８（Ｂ２）に図示するように画像［ｉ］のある部分が勾配磁場パルスの方向に沿って拡大したり、縮小したりという歪みが発生する。

つまり、画像の歪みは勾配磁場パルスの方向に限定されるものであり、具体的には、勾配磁場パルスの方向に沿って分割した帯状領域（図８（Ｃ）参照。）内の全体が略均一に歪み、例えば、図８（Ｄ）に示すように帯状領域内の特定の箇所で拡大縮小率が異なるように歪むということは起こらない。従って、帯状領域の中で少なくとも2箇所の歪み量がわかれば、この歪み量を用いて、その帯状領域全体の歪み量が一意に定まる。そこで、帯状領域全体の歪み量を用いて補正することにより、帯状領域の全体の歪みを解消することができる。そこで、各帯状領域に２つの特徴点Ｚを求め（図９（Ａ）参照。）、当該特徴点Ｚの歪み量Ｗａｒｐを夫々求める（図９（Ｂ）参照。）。

以下、歪み量Ｗａｒｐの取得について図１０に示すフローチャートを参照して説明する。図１０は特徴点Ｚ１及びＺ２の各画像に対する歪み量Ｗａｒｐ［ｚ］［ｉ］取得処理を示すフローチャートである。特徴点Ｚ１に対する当該処理（ステップＳ４７からの処理）も、特徴点Ｚ２に対する当該処理（ステップＳ４９からの処理）も同様であるので、以下の説明において特徴点Ｚ１及びＺ２を単に「特徴点Ｚ」として説明する。

先ず、特徴点Ｚの周囲に局所領域ＲＯＩ［ｚ］を設定する（ステップＳ７１）（図１１（Ａ）参照）。

続いて、相対的な歪み量取得の対象となる画像を示す変数ｉ２の値をｉ+１とする（ステップＳ７２）。この時、ｉ＝ｎ―１である場合にはｉ２の値は０とする。そして、画像［ｉ］と画像［ｉ２］の対象物体の方向を一致させるべく、画像［ｉ２］を、画像［ｉ２］の勾配磁場パルスの方向Ａｎｇｌｅ［ｉ２］から画像［ｉ１］の勾配磁場パルスの方向Ａｎｇｌｅ［ｉ］を減算した値分だけ回転させる（ステップＳ７３）（図１１（Ａ）及び（Ｂ）参照）。図１１（Ａ）及び（Ｂ）に示すように、画像に写っている対象物体の向きを揃えると、勾配磁場パルスの方向は異なる。

続いて、各画像上で局所領域ＲＯＩ［ｚ］にマッチングする領域を探し（図１１（Ｃ）参照）、歪み量Ｗａｒｐ［ｚ］［ｉ］を取得する（ステップＳ７４）。

ステップＳ７４における歪み量Ｗａｒｐ［ｚ］［ｉ］の取得手順について、図１１（Ｃ）及び図１２を参照して具体的に説明する。図１１及び図１２は、画像［ｉ２］に対する相対的な歪み量Ｗａｒｐ［ｚ］［ｉ］を求める場合の例である。

図１１（Ｃ）に示すように、画像［ｉ］上の特徴点Ｚの局所領域ＲＯＩ［ｚ］と最も類似する箇所を、パターンマッチングにより画像［ｉ２］上で探すことにより、特徴点Ｚに対応する画像［ｉ２］上の対応点Ｕを検索し、取得する。この場合、特徴点Ｚの有する領域と同一の領域を特定する変数ｐ及びｑを有する画像［ｉ２］上の場所から検索を始めると効率がよい。そうすると図１２に示すように、２つの画像上の夫々対応する点同士を結びつけることができる。

続いて、検索された画像［ｉ２］上の点Ｕに対する画像［ｉ］の特徴点Ｚの歪み量Ｗａｒｐ［ｚ］［ｉ］を求める。

図１２を参照して説明すると、先ず、特徴点Ｚを通り、画像［ｉ］の勾配磁場パルスの方向に平行な直線Ｌ［ｉ］と、当該特徴点Ｚに対応する画像［ｉ２］上の点Ｕを通り、画像［ｉ２］の勾配磁場パルスの方向に平行な直線Ｌ［ｉ２］と、の交点を求める。そして、この交点と特徴点Ｚの距離を、特徴点Ｚの歪み量Ｗａｒｐ［ｚ］［ｉ］として取得する。

以上説明したように、ステップＳ４７では、画像［ｉ］上の全ての特徴点Ｚ１について、当該特徴点Ｚ１の各画像に対する歪み量Ｗａｒｐ［ｚ］［ｉ］を取得して当該ステップＳ４７の処理を終了し、一方、ステップＳ４９では、画像［ｉ］上の全ての特徴点Ｚ２について、当該特徴点Ｚ２の各画像に対する歪み量Ｗａｒｐ［Ｚ２］［ｉ］を取得して当該ステップＳ４９の処理を終了する。

２−２．歪み量Ｗａｒｐ［ｉ］の異常検出及び補間
続いて、上述したステップＳ１６及びステップＳ１７における歪み量Ｗａｒｐ［ｉ］の異常検出と補間について、図１３及び図１４を用いて詳細に説明する。

図１３（Ａ）はステップＳ１２において確保された歪み量を表現する配列Ｗａｒｐ［p,q］［ｉ］を模式的に示した図であり、（Ｂ）は、ステップＳ１５にて取得した各特徴点Ｚの歪み量Ｗａｒｐ［ｚ］［ｉ］を当該配列表現によって模式的に示した図である。

図１３（Ｂ）に示すように、画像［ｉ］のある部分が勾配磁場パルスの方向に沿って縮小した歪みが発生している場合に、特徴点Ｚの歪み量Ｗａｒｐ［ｚ］［ｉ］の配列が交差しているなど、不自然な箇所がある場合には、この特徴点Ｚは異常であるとして削除し、削除された特徴点Ｚの歪み量Ｗａｒｐ［ｚ］［ｉ］と、特徴点Ｚ以外の点の歪み量Ｗａｒｐ［ｚ］［ｉ］に対して補間処理を行なう（図１３（Ｃ）参照。）。

当該ステップＳ１７における補間処理について、図１４を参照して説明する。図１４は、異常検出によって削除された特徴点Ｚの歪み量Ｗａｒｐ［ｚ］［ｉ］と、特徴点Ｚ以外の領域（ｐ、ｑ）の歪み量Ｗａｒｐ［p,q］［ｉ］の補間を示す説明図である。

図１４（Ａ）に示すように、歪み量Ｗａｒｐ［p,q］は各特徴点Ｚに対応する要素のみ計算されているので疎である。各特徴点Ｚを“■”で表し、当該特徴点Ｚの添え字はマッチングされた歪み量Ｗａｒｐ［p,q］の添え字ｑ（０〜３）である。そして、歪み量Ｗａｒｐ［p,q］の値は勾配磁場パルス方向への歪み量を表す。各帯状領域毎に２つの特徴点Ｚがあることになるが、異常検出によって削除された箇所はステップＳ１６の処理で削除されている。

ここで、帯状領域の中で少なくとも２箇所の歪み量がわかれば、この歪み量を用いて、その帯状領域全体の歪み量が一意に定まる。よって、図１４（Ｂ）に示すように、勾配磁場パルス方向には情報の得られている2箇所を用いて勾配磁場パルス方向の歪み量が定まる。更に、各配列を同一q座標値同士で繋ぐよう構成することにより、歪み量Ｗａｒｐ［p,q］の値を持たない要素を補間することができる。これによって、画像［ｉ］の全ての領域について、他の画像、厳密には、ステップＳ１４又はＳ１８にて変数ｉによって定義される他の画像（例えば、画像［ｉ＋１］）に対する歪み量Ｗａｒｐ［p,q］［ｉ］が取得できる。

そして、上述したステップＳ１８乃至ステップＳ１９を経て、ステップＳ１１にて取得した全ての画像［ｉ］（０≦ｉ＜ｎ）について、他の全ての画像に対する歪み量Ｗａｒｐ［ｉ］を取得した後に、ステップＳ２０にて合成処理を行なうこととなる。

２−３．歪み量Ｗａｒｐ［ｉ］の異常検出及び補間
続いて、ステップＳ２０における「ｎ枚の画像の合成処理」について、図１５に示すフローチャートを用いて説明する。図１５は、「ｎ枚の画像の合成処理」を示すフローチャートである。

先ず、変数ｉを０として初期化（ステップＳ８１）し、ステップＳ１３で取得した画像［ｉ］の勾配磁場パルスの方向Ａｎｇｌｅ［ｉ］と、ステップＳ１５乃至Ｓ１９にて取得した画像［ｉ］の他の画像に対する歪み量Ｗａｒｐ［ｉ］に基づいて、画像［ｉ］を補正する（ステップＳ８２）。そして、変数ｉを１増加（ステップＳ８３）させ、当該変数ｉが最初にステップＳ１１において取得した画像の枚数ｎよりも小さいか否かを判定（ステップＳ８４）し、小さい場合（ステップＳ８４：Ｙｅｓ）には、まだ補正処理を行なうべき画像があるので、ステップＳ８２へ戻り、補正処理を次の画像に対して行ない、変数ｉがｎ以上である場合（ステップＳ８４：Ｎｏ）にはステップＳ８５へ移行する。

そして、補正処理後のｎ枚の画像を重ね合わせ（フュージョン）て、合成画像を作成して処理を終了する（ステップＳ８５）。なお、作成した合成画像は表示部１４に表示して医師等の操作者に提示したり、或いは図示しないプリンタ等の出力部にて画像出力することが可能である。

３.画像処理の実施形態３
次に、画像処理の実施形態３について説明する。

上述した画像処理の実施形態２では、１の画像の各局所領域に対応する他の画像上の対応領域を、パターンマッチングにより順番に検索したが、例えば、図１６に示すように、複数の画像（同図では４枚の画像（ｎ＝４））について互いに補正を行なったとしても、１の画像は他の画像の内の1つ（画像ｉは画像ｉ＋１）を参照している場合にすぎなかった。

本実施形態では、１の画像は他の画像全てを参照することによって、より歪み補正の精度を上げようとするものである。すなわち、本実施形態では歪み量Ｗａｒｐは、画像ｉは画像ｉ以外の全ての画像を検索対象とし、歪み量Ｗａｒｐ２［p, q］［ｉ］［ｉ２］（ｉ≠ｉ２）と表現する。これは画像ｉが画像ｉ２を参照したときの歪み量である。歪み量Ｗａｒｐ［p, q］［ｉ］は画像ｉの歪みの無い画像と比較しての歪みの量を表現しているが、この歪み量の算出は別の歪みのある画像ｉ２を用いて作成されたものである。その為に比較対象とする画像ｉ２が複数存在する場合には機器の誤差等により、計算される画像ｉの歪み量が異なることになる。本実施形態では比較対象によって計算される歪み量の差異を吸収することを目的としている。画像処理の実施形態２において、歪み量Ｗａｒｐ２［p, q］［ｉ］［ｉ＋１］のみを計算に用いたことに対応する。なお、歪み量Ｗａｒｐ２の値の取得手順等、拡散強調画像処理装置１の構成及び機能は、歪み量Ｗａｒｐの値を取得する構成に対して存在する２枚の画像毎に複数回計算する点を除いて、画像処理の実施形態２と同様であるので詳細な説明を省略する。

以下、本実施形態における画像ｉの、画像ｉ以外の全ての画像に対する歪み量Ｗａｒｐ［p, q］［ｉ］の取得処理について説明する。拡散強調画像は、一方向に平行移動、剪断及び拡大縮小する歪みであるため、勾配磁場パルスの方向の平行移動量と、勾配磁場パルスの方向の拡大量を把握することにより、画像ｉの、画像ｉ以外の全ての画像に対する歪み量Ｗａｒｐ［p, q］［ｉ］を求めることができる。なお、ここでの勾配磁場パルスによって発生する歪みの方向をｑ方向とする。

先ず、画像ｉの画像ｉ２（０≦ｉ２＜０）に対する歪み量Ｗａｒｐ２［p, q］［ｉ］［ｉ２］から、勾配磁場パルスの方向の平行移動量Translation１［p］［ｉ］［ｉ２］と、勾配磁場パルスの方向の拡大量Expansion１［p］［ｉ］［ｉ２］を取得する（ステップＳ９１）。歪み量Ｗａｒｐ２［p, q］［ｉ］［ｉ２］と平行移動量Translation１［p］［ｉ］［ｉ２］と、拡大量Expansion１［p］［ｉ］［ｉ２］の関係を、以下の式1に示す。

次に、画像ｉの、他の画像ｉ２（０≦ｉ２＜０）に対する平行移動量Translation１［p］［ｉ］［ｉ２（０≦ｉ２＜０）］の平均値と、拡大量Expansion１［p］［ｉ］［ｉ２（０≦ｉ２＜０）］の平均値を、夫々、平行移動量Translation２［p］［ｉ］及び拡大量Expansion２［p］［ｉ］として取得する（ステップＳ９２、Ｓ９３）。

そして、取得した平行移動量Translation２［p］［ｉ］、拡大量Expansion２［p］［ｉ］から、以下の式２に従って、画像ｉの、画像ｉ以外の全ての画像に対する歪み量Ｗａｒｐ［p, q］［ｉ］を取得して処理を終了する（ステップＳ９４）。

上記ステップＳ９１〜Ｓ９４の処理を、画像ｉ以外の他の画像の全てに対して繰り返し行なうことにより、画像ｉの、画像ｉ以外の全ての画像に対する歪み量Ｗａｒｐ［p, q］［ｉ］を求めることができる。

４．画像処理の実施形態４
次に、画像処理の実施形態４について説明する。上述した画像処理の各実施形態１〜３を用いても、取得された画像の品質が乏しい場合や、歪みの程度が著しい場合には安定した補正が行なえない場合がある。本実施形態では、歪みの自由度を制限することによってより安定した歪み補正を行う事を目的とする。

上述した画像処理の各実施形態１〜３で補正対象としている歪みは、一方向に平行移動、剪断及び拡大縮小するような歪みであるから、歪みの自由度は勾配磁場パルスの方向への平行移動及び拡大縮小の２自由度しかない。その為にＷａｒｐ［p, q］［ｉ］の形態で表現される歪みは勾配磁場パルスの方向の平行移動量Translation3［p］［ｉ］及び、勾配磁場パルスの方向の拡大量Expansion3［p］［ｉ］の組み合わせとして表現することができる。

平行移動量Translation3［p］［ｉ］及び拡大量Expansion3［p］［ｉ］は、いずれもｐを引数とした区分連続関数としての平行移動量Translation４［p］［ｉ］及び拡大量Expansion４［p］［ｉ］とも解釈できる。そこで、Expansion4（ｐ）［ｉ］及び、Expansion4（ｐ）［ｉ］をそれぞれ３次関数Ｆ（ａ×ｐ^３＋ｂ×ｐ^２＋ｃ×ｐ＋ｄ）［ｉ］に近似する。なお、近似は例えば最小二乗法を用いればよい。また、近似する関数は３次関数に限定されず、容易に表現できる関数であればよく、例えば、スプライン曲線や２次関数等でもよい。

更に、区分連続関数から３次関数Ｆへの変換時に、区分連続関数と３次関数Ｆの相違度を最小化するのではなく、ｉ〜（ｎ−１）の歪み補正が行なわれた夫々の画像の相違度を最小化しても良い。かかる相違度は例えば差分画像の夫々のピクセルの二乗和で求めることができる。

以上、説明した実施形態によれば、被検体への負担を強いることなく、勾配磁場パルスによって発生する画像（拡散強調画像）の歪みをを補正して所望の数の画像（拡散強調画像）からなる合成画像を取得することができる。
なお、上述した実施形態において、各画像の撮像条件に起因する特質は、勾配磁場パルス方向に依存する水分子の拡散に基づく画像強調としたが、本発明はこれに限定されるものでは無く、撮像条件により異方性を含んだ歪みが発生する種々の画像の特質を、各画像を重ね合わせて合成画像により観察（解析）を所望する場合のいずれにおいても適用可能である。

更に、上述した図２、図４、図５、図７（Ａ）、図１０、図１５及び図１７に示すフローチャートに対応するプログラムを、フレキシブルディスク又はハードディスク等の情報記録媒体に記録しておき、又はインターネット等を介して取得して記録しておき、これらを汎用のコンピュータで読み出して実行することにより、当該コンピュータを実施形態に係る制御部１１として機能させることも可能である。

以上夫々説明したように、本発明は拡散強調撮像法によって取得された複数の拡散強調画像を重ね合わせ（フュージョン）た合成画像を用いて行なう画像診断、画像解析処理の分野に利用することが可能であり、特に勾配磁場パルスがそれぞれ異なる複数の拡散強調画像の合成画像解析処理の分野に適用すれば特に顕著な効果が得られる。

本発明に係る拡散強調画像処理装置の構成を示すブロック図である。 画像処理の実施形態１にかかる画像処理を示すフローチャートである。 特徴点決定手法を示す説明図である。 画像処理の実施形態２にかかる画像処理を示すフローチャートである。 歪み量Ｗａｒｐ［ｉ］取得処理を示すフローチャートである。 局所領域ＲＯＩ［p,q］における特徴量ＣＶ［p,q］［ｉ］算出の説明図である。 （Ａ）は、局所領域ＲＯＩ［p,q］の特徴量ＣＶ［p,q］［ｉ］算出処理を示すフローチャートである。（Ｂ１）、（Ｂ２）及び（Ｂ３）は種々の局所領域ＲＯＩにかかる特徴量算出の説明図である。 拡散強調画像の歪みの特徴を示す説明図である。 特徴点の説明図である。 特徴点Ｚ１及びＺ２の各画像に対する歪み量Ｗａｒｐ［ｚ］［ｉ］取得処理を示すフローチャートである。 特徴点Ｚ１及びＺ２の各画像に対する歪み量Ｗａｒｐ［ｚ］［ｉ］取得処理の説明図である。 特徴点Ｚ１及びＺ２の各画像に対する歪み量Ｗａｒｐ［ｚ］［ｉ］取得処理の説明図である。 各特徴点Ｚの歪み量Ｗａｒｐ［ｚ］［ｉ］を当該配列表現によって模式的に示した図である。 歪み量Ｗａｒｐ［p,q］［ｉ］の補間を示す説明図である。 ｎ枚の画像の合成処理を示すフローチャートである。 ４枚の画像の画像の相互補正を示す説明図である。 画像処理の実施形態３にかかる歪み量Ｗａｒｐ取得処理を示すフローチャートである。

符号の説明

１ 拡散強調画像処理装置
１１ 制御部
１２ 外部機器接続部
１３ 記憶部
１４ 表示部
１５ 入力部
１６ バス
Ｚ 特徴点
Ｕ 特徴点Ｚに対応する点
Ｗａｒｐ 歪み量
ＲＯＩ 局所領域
ＣＶ 特徴量
ＣＶ_ＭＡＸ 最大特徴量
ＣＶ_Ｓｅｃｏｎｄ 第２位特徴量
ＣＶ_ＭＩＮ 最小特徴量

Claims (12)

1. 撮像された際の勾配磁場パルス方向がそれぞれ異なる複数の拡散強調画像を、撮像された際の前記勾配磁場パルス方向に沿って複数の帯状領域に分割する第１分割手段と、
   前記各拡散強調画像を互いに参照し、かつ前記各帯状領域に基づいて前記勾配磁場パルスによって発生する前記勾配磁場パルス方向の歪みを相対的に補正する補正手段を有することを特徴とする拡散強調画像処理装置。
2. 請求項１に記載の拡散強調画像処理装置において、
   前記帯状領域を更に複数の分割領域に分割する第２分割手段と、
   前記帯状領域毎に少なくとも２以上の分割領域と、当該各分割領域と対応する他の前記拡散強調画像の対応領域を、パターンマッチングにより、それぞれ検索する検索手段と、を有し、
   前記補正手段は、各前記対応領域を参照し、かつ、前記一の拡散強調画像と前記他の拡散強調画像のそれぞれの前記勾配磁場パルスの方向に基づいて、前記一の拡散強調画像の前記各分割領域について補正を行なうことを特徴とする拡散強調画像処理装置。
3. 請求項１に記載の拡散強調画像処理装置において、
   前記一の拡散強調画像のうち、画像の特徴を表す箇所を、前記帯状領域毎に少なくとも２点決定する特徴点決定手段と、
   各前記特徴点の周囲を局所領域として、当該局所領域に対応する他の前記拡散強調画像の対応領域をパターンマッチングによりそれぞれ検索し、検索された各前記対応領域上の点を、前記各特徴点に対応する対応点としてそれぞれ取得する対応点取得手段と、を有し、
   前記補正手段は、各前記対応点の位置座標と、前記一の拡散強調画像と前記他の拡散強調画像のそれぞれの前記勾配磁場パルスの方向と、に基づいて、前記各特徴点の補正を行なうと共に、補正後の前記各特徴点を用いて、前記一の拡散強調画像の前記各特徴点が属する前記帯状領域の他の領域について補正を行なうことを特徴とする拡散強調画像処理装置。
4. 請求項１に記載の拡散強調画像処理装置において、
   前記一の拡散強調画像のうち、画像の特徴を表す箇所を、前記帯状領域毎に少なくとも２点決定する特徴点決定手段と、
   各前記特徴点の周囲を局所領域として、当該局所領域に対応する他の前記拡散強調画像の対応領域をパターンマッチングによりそれぞれ検索し、検索された各前記対応領域上の点を、前記各特徴点に対応する対応点としてそれぞれ取得する対応点取得手段と、
   各前記対応点の位置座標と、前記一の拡散強調画像と前記他の拡散強調画像のそれぞれの前記勾配磁場パルスの方向と、に基づいて、前記特徴点毎に歪み量をそれぞれ取得する特徴点歪み量取得手段と、を有し、
   前記補正手段は、前記帯状領域毎に当該各帯状領域に属する前記各特徴点の前記歪み量を用いて各前記拡散強調画像の補正を行なうことを特徴とする拡散強調画像処理装置。
5. 請求項１乃至４のいずれか一項に記載の拡散強調画像処理装置において、
   前記補正手段によって補正された前記各拡散強調画像を重ね合わせて、合成画像を作成する合成手段を有することを特徴とする拡散強調画像処理装置。
6. 請求項１乃至５のいずれか一項に記載の拡散強調画像処理装置において、
   前記勾配磁場パルスに起因する特質は、前記拡散強調画像を撮影する際の前記勾配磁場パルス方向に依存する水分子の拡散の異方性に基づく画像強調であることを特徴とする拡散強調画像処理装置。
7. コンピュータを、撮像された際の勾配磁場パルス方向がそれぞれ異なる複数の拡散強調画像を、撮像された際の前記勾配磁場パルス方向に沿って複数の帯状領域に分割する第１分割手段、
   前記各拡散強調画像を互いに参照し、かつ前記各帯状領域に基づいて、前記勾配磁場パルスによって発生する前記勾配磁場パルス方向の歪みを相対的に補正する補正手段として機能させることを特徴とする画像処理プログラム。
8. 請求項７に記載の画像処理プログラムにおいて、
   前記コンピュータを、
   前記帯状領域を更に複数の分割領域に分割する第２分割手段、
   前記帯状領域毎に少なくとも２以上の分割領域と、当該各分割領域と対応する他の前記拡散強調画像の対応領域を、パターンマッチングによりそれぞれ検索する検索手段として更に機能させ、
   前記補正手段を、各前記対応領域を参照し、かつ、前記一の拡散強調画像と前記他の拡散強調画像のそれぞれの前記勾配磁場パルスの方向に基づいて、前記一の拡散強調画像の前記各分割領域について補正を行なうと共に、補正後の前記各分割領域を用いて、前記一の拡散強調画像の前記各分割領域が属する前記帯状領域の他の領域について補正を行なうよう機能させることを特徴とする画像処理プログラム。
9. 請求項７に記載の画像処理プログラムにおいて、
   前記コンピュータを、
   前記一の拡散強調画像のうち、画像の特徴を表す箇所を、前記帯状領域毎に少なくとも２点決定する特徴点決定手段、
   各前記特徴点の周囲を局所領域として、当該局所領域に対応する他の前記拡散強調画像の対応領域をパターンマッチングによりそれぞれ検索し、検索された各前記対応領域の中心点を、前記各特徴点に対応する対応点としてそれぞれ取得する対応点取得手段としてさらに機能させ、
   前記補正手段を、各前記対応点の位置座標と、前記一の拡散強調画像と前記他の拡散強調画像のそれぞれの前記勾配磁場パルスの方向と、に基づいて、前記各特徴点の補正を行なうよう機能させることを特徴とする画像処理プログラム。
10. 請求項７に記載の画像処理プログラムにおいて、
    前記コンピュータを、
    前記一の拡散強調画像のうち、画像の特徴を表す箇所を、前記帯状領域毎に少なくとも２点決定する特徴点決定手段と、
    各前記特徴点の周囲を局所領域として、当該局所領域に対応する他の前記拡散強調画像の対応領域をパターンマッチングによりそれぞれ検索し、検索された各前記対応領域上の点を、前記各特徴点に対応する対応点としてそれぞれ取得する対応点取得手段、
    各前記対応点の位置座標と、前記一の拡散強調画像と前記他の拡散強調画像のそれぞれの前記勾配磁場パルスの方向と、に基づいて、前記特徴点毎に歪み量をそれぞれ取得する特徴点歪み量取得手段として更に機能させ、 前記補正手段を、前記帯状領域毎に当該各帯状領域に属する前記各特徴点の前記歪み量を用いて補正を行なうよう機能させることを特徴とする画像処理プログラム。
11. 請求項７乃至１０のいずれか一項に記載の画像処理プログラムにおいて、
    前記コンピュータを、補正された前記各拡散強調画像を重ね合わせて、合成画像を作成する合成手段として更に機能させることを特徴とする画像処理プログラム。
12. 請求項７乃至１１のいずれか一項に記載の画像処理プログラムにおいて、
    前記勾配磁場パルスに起因する特質は、前記拡散強調画像を撮影する際の前記勾配磁場パルス方向に依存する水分子の拡散の異方性に基づく画像強調であることを特徴とする画像処理プログラム。

Images