JP6452994B2 - 画像処理装置及び磁気共鳴イメージング装置 - Google Patents

Description

本発明の実施形態は、画像処理装置及び磁気共鳴イメージング装置に関する。

磁気共鳴イメージング装置（以下、適宜「ＭＲＩ（Magnetic Resonance Imaging）装置」）は、核磁気共鳴現象を利用して被検体内部の情報を画像化する装置である。ＭＲＩ装置は、コイルを利用して、被検体内部にある特定の原子核（例えば、水素原子の原子核）からの核磁気共鳴信号（ＭＲ信号）をサンプリングすることで、ｋ空間データと呼ばれるデータを収集し、ｋ空間データにフーリエ変換を適用することで、ＭＲ画像を得る。

ＭＲ信号は、１次元データとしてサンプリングされる。ＭＲＩ装置は、２次元のＭＲ画像、又は、３次元のＭＲ画像を得るために、ｋ空間上で１次元のサンプリングを繰り返し行い、ＭＲ画像の生成に必要なデータを収集する。出力対象のＭＲ画像と同じ解像度（フルサンプリング）でｋ空間データがサンプリングされれば、得られたｋ空間データにフーリエ変換を適用することで、ＭＲ画像の生成が可能である。さらに、時系列の画像を撮影することで、被検体の動きを観察することが可能となる。

ＭＲＩ装置で生成される画像は、強度画像と位相画像との２種類が挙げられる。位相画像は流速を画像化する際に用いられ、例えば位相コントラストＭＲ血管撮影（phase contrast MR angiography）では、強度画像の時系列画像と位相画像の時系列画像とが出力される。

ところで、ｋ空間データ収集時やＭＲ画像再構成時にはノイズが生じる場合がある。これらのノイズ成分を除去し、画質を改善する方法の一つとして、フィルタを用いる方法がある。

特許第３１６３５７０号公報

M.Doyle and P.Mansfield, "Real-time image enhancement in NMR", Journal of Physics, vol.19, pp.439-444, 1986.

本発明が解決しようとする課題は、正しい流速データが得られ、かつ、強度画像の画質を向上することができる画像処理装置及び磁気共鳴イメージング装置を提供することである。

実施形態に係る画像処理装置は、変換部と、強度画像生成部と、位相画像生成部とを備える。変換部は、磁気共鳴イメージング装置によって、時系列に沿って収集された時系列ｋ空間データを、空間上の位置を示すｘ空間上の第１の時系列ｘ空間データに変換する変換処理を行う。強度画像生成部は、前記第１の時系列ｘ空間データに第１のフィルタを適用して取得した第２の時系列ｘ空間データから、強度画像を生成する。位相画像生成部は、前記第１の時系列ｘ空間データに、前記第１のフィルタとは異なる第２のフィルタを適用して取得した第３の時系列ｘ空間データから、位相画像を生成する。前記第１のフィルタは、前記第２のフィルタより、データを平滑化する度合いの大きいフィルタである。実施形態に係る磁気共鳴イメージング装置は、収集部と、前記画像処理装置とを備える。収集部は、時系列に沿ったｋ空間データである時系列ｋ空間データを収集する。前記画像処理装置は、時系列ｋ空間データを用いて強度画像及び位相画像を生成する。

第１の実施形態に係るＭＲＩ装置の構成を示す機能ブロック図。 位相コントラスト法の概略を説明するための図。 第１の実施形態における画像生成の概略を説明するための図。 第１の実施形態における画像生成部の処理の流れを説明するための図。 第１の実施形態に係るＭＲＩ装置の処理手順を示すフローチャート。 第１の実施形態に係る画像生成部の処理手順の一例を示すフローチャート。 第１の実施形態におけるローパスフィルタの一例を示す図。 第１の実施形態の第１の変形例におけるフィルタの一例を示す図。 第２の実施形態に係る画像生成部の構成を示す機能ブロック図。 第２の実施形態に係る画像生成部の処理手順の一例を示すフローチャート。

以下、添付図面を参照して、画像処理装置及び磁気共鳴イメージング装置（以下、適宜「ＭＲＩ（Magnetic Resonance Imaging）装置」）の実施形態を詳細に説明する。以下では、画像処理装置の機能を有するＭＲＩ装置を実施形態として説明する。なお、実施形態は、以下の実施形態に限られるものではない。また、各実施形態において説明する内容は、原則として、他の実施形態においても同様に適用することができる。

（第１の実施形態）
図１は、第１の実施形態に係るＭＲＩ装置の構成を示す機能ブロック図である。図１に示すように、ＭＲＩ装置１００は、静磁場磁石１０１と、傾斜磁場コイル１０２と、傾斜磁場電源１０３と、寝台１０４と、寝台制御部１０５と、送信コイル１０６と、送信部１０７と、受信コイル１０８と、受信部１０９と、シーケンス制御部１１０と、計算機システム１２０とを備える。なお、ＭＲＩ装置１００に被検体Ｐ（例えば、人体）は含まれない。また、図１に示す構成は一例に過ぎない。例えば、シーケンス制御部１１０及び計算機システム１２０内の各部は、適宜統合若しくは分離して構成されてもよい。

静磁場磁石１０１は、中空の略円筒形状に形成された磁石であり、内部の空間に一様な静磁場を発生する。静磁場磁石１０１は、例えば、超伝導磁石等であり、図示しない静磁場電源から電流の供給を受けて励磁する。なお、静磁場磁石１０１は、永久磁石でもよく、この場合、ＭＲＩ装置１００は、静磁場電源を備えなくてもよい。

傾斜磁場コイル１０２は、中空の略円筒形状に形成されたコイルであり、静磁場磁石１０１の内側に配置される。傾斜磁場コイル１０２は、互いに直交するＸ，Ｙ，Ｚの各軸に対応する３つのコイルが組み合わされて形成されており、これら３つのコイルは、傾斜磁場電源１０３から個別に電流を受けて、Ｘ，Ｙ，Ｚの各軸に沿って磁場強度が変化する傾斜磁場を発生させる。なお、Ｚ軸方向は、静磁場と同方向とする。

傾斜磁場電源１０３は、傾斜磁場コイル１０２に電流を供給する。ここで、傾斜磁場コイル１０２によって発生するＸ，Ｙ，Ｚ各軸の傾斜磁場は、例えば、スライス選択用傾斜磁場Ｇｓ、位相エンコード用傾斜磁場Ｇｅ、及び、リードアウト用傾斜磁場Ｇｒにそれぞれ対応する。スライス選択用傾斜磁場Ｇｓは、任意に撮像断面を決めるために利用される。位相エンコード用傾斜磁場Ｇｅは、空間的位置に応じて磁気共鳴信号（以下、適宜「ＭＲ信号」）の位相を変化させるために利用される。リードアウト用傾斜磁場Ｇｒは、空間的位置に応じてＭＲ信号の周波数を変化させるために利用される。

寝台１０４は、被検体Ｐが載置される天板１０４ａを備え、寝台制御部１０５による制御のもと、天板１０４ａを、被検体Ｐが載置された状態で傾斜磁場コイル１０２の空洞（撮像口）内へ挿入する。通常、寝台１０４は、長手方向が静磁場磁石１０１の中心軸と平行になるように設置される。寝台制御部１０５は、計算機システム１２０による制御のもと、寝台１０４を駆動して天板１０４ａを長手方向及び上下方向へ移動する。

送信コイル１０６は、傾斜磁場コイル１０２の内側に配置され、送信部１０７からＲＦ（Radio Frequency）パルスの供給を受けて、高周波磁場を発生する。送信部１０７は、対象とする原子核の種類及び磁場の強度で決まるラーモア周波数に対応するＲＦパルスを送信コイル１０６に供給する。

受信コイル１０８は、傾斜磁場コイル１０２の内側に配置され、高周波磁場の影響によって被検体Ｐから発せられるＭＲ信号を受信する。受信コイル１０８は、ＭＲ信号を受信すると、受信したＭＲ信号を受信部１０９へ出力する。なお、第１の実施形態において、受信コイル１０８は、１以上、典型的には複数の受信コイルを有するコイルアレイである。なお、上述した送信コイル１０６及び受信コイル１０８は一例に過ぎない。送信機能のみを備えたコイル、受信機能のみを備えたコイル、若しくは送受信機能を備えたコイルのうち、１つ若しくは複数を組み合わせることによって構成されればよい。

受信部１０９は、受信コイル１０８から出力されるＭＲ信号に基づいてＭＲデータを生成する。具体的には、受信部１０９は、受信コイル１０８から出力されるＭＲ信号をデジタル変換することによってＭＲデータを生成する。また、受信部１０９は、生成したＭＲデータをシーケンス制御部１１０へ送信する。なお、受信部１０９は、静磁場磁石１０１や傾斜磁場コイル１０２等を備える架台装置側に備えられていてもよい。ここで、第１の実施形態において、受信コイル１０８の各コイルエレメントから出力されるＭＲ信号は、適宜分配合成されることで、チャネル等と呼ばれる単位で受信部１０９に出力される。このため、ＭＲデータは、受信部１０９以降の後段の処理においてチャネル毎に取り扱われる。コイルエレメントの総数とチャネルの総数との関係は、同一の場合もあれば、コイルエレメントの総数に対してチャネルの総数が少ない場合、あるいは反対に、コイルエレメントの総数に対してチャネルの総数が多い場合もある。以下において、「チャネル毎」のように表記する場合、その処理が、コイルエレメント毎に行われてもよいし、あるいは、コイルエレメントが分配合成されたチャネル毎に行われてもよいことを示す。なお、分配合成のタイミングは、上述したタイミングに限られるものではない。ＭＲ信号若しくはＭＲデータは、後述する画像生成部１２２による処理の前までに、チャネル単位に分配合成されればよい。

シーケンス制御部１１０は、計算機システム１２０から送信されるシーケンス情報に基づいて、傾斜磁場電源１０３、送信部１０７及び受信部１０９を駆動することによって、被検体Ｐの撮像を行う。ここで、シーケンス情報は、撮像を行うための手順を定義した情報である。シーケンス情報には、傾斜磁場電源１０３が傾斜磁場コイル１０２に供給する電源の強さや電源を供給するタイミング、送信部１０７が送信コイル１０６に送信するＲＦパルスの強さやＲＦパルスを印加するタイミング、受信部１０９がＭＲ信号を検出するタイミング等が定義される。例えば、シーケンス制御部１１０は、ＡＳＩＣ（Application Specific Integrated Circuit）、ＦＰＧＡ（Field Programmable Gate Array）等の集積回路、ＣＰＵ（Central Processing Unit）、ＭＰＵ（Micro Processing Unit）等の電子回路である。

なお、シーケンス制御部１１０は、傾斜磁場電源１０３、送信部１０７及び受信部１０９を駆動して被検体Ｐを撮像した結果、受信部１０９からＭＲデータを受信すると、受信したＭＲデータを計算機システム１２０へ転送する。

計算機システム１２０は、ＭＲＩ装置１００の全体制御や、データ収集、画像生成等を行い、インタフェース部１２１、画像再構成部１２２、記憶部１２３、入力部１２４、表示部１２５及び制御部１２６を有する。また、制御部１２６は、収集部１２６ａを有する。

インタフェース部１２１は、シーケンス情報をシーケンス制御部１１０へ送信し、シーケンス制御部１１０からＭＲデータを受信する。また、インタフェース部１２１は、ＭＲデータを受信すると、受信したＭＲデータを記憶部１２３に格納する。記憶部１２３に格納されたＭＲデータは、制御部１２６が有する収集部１２６ａによってｋ空間に配置される。この結果、記憶部１２３は、ｋ空間データを記憶する。例えば、記憶部１２３は、単一または複数チャネル分の時系列のｋ空間データを記憶する。

画像生成部１２２は、収集部１２６ａで収集されたデータや、記憶部１２３で記憶されたデータを用いた画像処理により、画像の生成を行う画像処理装置として機能する。画像生成部１２２は、図１に示すように、変換部１２２ａ、強度画像生成部１２２ｂ及び位相画像生成部１２２ｃを有する。変換部１２２ａは、収集部１２６ａによって収集されたｋ空間データに対して、例えば、空間方向にフーリエ変換を適用することにより変換後時空間データを得る。強度画像生成部１２２ｂは、上記変換後時空間データから強度画像（時空間強度画像）を生成する。位相画像生成部１２２ｃは、上記変換後時空間データから位相画像（時空間位相画像）を生成する。画像生成部１２２によって得られた画像は、必要に応じて表示部１２５で表示されたり、記憶部１２３に送られて記憶されたりする。また、例えば、位相画像生成部１２２ｃは、時空間位相画像の局所領域から流速値の時間変化を求め、求めた流速値の時間変化を示すグラフを生成する。位相画像生成部１２２ｃが求めた流速値やグラフは、必要に応じて表示部１２５で表示されたり、記憶部１２３に送られて記憶されたりする。なお、画像生成部１２２が行う画像処理については後に詳述する。

記憶部１２３は、インタフェース部１２１によって受信されたＭＲデータや、収集部１２６ａによって収集されたｋ空間データ、画像生成部１２２によって生成された画像データ等を記憶する。例えば、記憶部１２３は、ＲＡＭ（Random Access Memory）、フラッシュメモリ等の半導体メモリ素子、ハードディスク、光ディスク等である。

入力部１２４は、操作者からの各種指示や情報入力を受け付ける。入力部１２４は、例えば、マウスやトラックボール等のポインティングデバイス、あるいはキーボード等の入力デバイスである。表示部１２５は、制御部１２６による制御のもと、画像データ等の各種の情報を表示する。表示部１２５は、例えば、液晶表示器等の表示デバイスである。

制御部１２６は、ＭＲＩ装置１００の全体制御を行う。具体的には、制御部１２６は、上述した収集部１２６ａを備え、収集部１２６ａが、入力部１２４を介してユーザから入力される撮像条件に基づいてシーケンス情報を生成し、生成したシーケンス情報をシーケンス制御部１１０に送信することによって撮像を制御する。また、収集部１２６ａは、ＭＲデータをｋ空間に配置することで、ｋ空間データを収集する。具体的には、収集部１２６ａは、空間方向のサンプリング位置を変化させながら、複数のチャネルについて、時系列に沿ってｋ空間データを収集する。以下では、時系列に沿って収集されたｋ空間データを、時系列ｋ空間データと記載する。また、制御部１２６は、撮像の結果としてシーケンス制御部１１０から送られるＭＲデータを用いて行われる画像生成部１２２による画像の生成を制御したり、表示部１２５による表示を制御したりする。例えば、制御部１２６は、ＡＳＩＣ（Application Specific Integrated Circuit）、ＦＰＧＡ（Field Programmable Gate Array）等の集積回路、ＣＰＵ（Central Processing Unit）、ＭＰＵ（Micro Processing Unit）等の電子回路である。

以上、第１の実施形態に係るＭＲＩ装置１００の全体構成について説明した。かかる構成のもと、画像処理装置としての画像生成部１２２は、収集部１２６ａにより収集された時系列ｋ空間データを取得して、強度画像及び位相画像を生成する。ここで、第１の実施形態に係る画像生成部１２２について説明する前に、以下の実施形態で用いる用語について説明する。

以下の実施形態において、「ｘ空間」とは、空間上の位置を示す「ｘ」の空間であり、便宜上、水素原子以外を撮像の対象とする場合も含めて水素分布画像と呼ぶと、水素分布画像空間のことである。また、「ｘ空間データ」とは、ｘ空間内における信号点の集合のことである。ｘ空間上の異なる信号点は、ｘ空間上の異なる位置の信号点と対応付けられる。例えば、脳の３次元の水素分布画像は、３次元のｘ空間データであり、心臓のある断面を撮像した２次元のｘ空間画像は、２次元のｘ空間データである。

また、「ｋ−ｔ空間データ」とは、時系列ｋ空間データである。ｋ−ｔ空間データに対して、ｋ空間データを、ｘ空間データに変換する変換、例えばフーリエ変換を適用すると、時系列ｘ空間データ（ｘ−ｔ空間データ）が得られる。

また、「ｔ−ｆ変換」とは、時間軸方向の軸ｔにおけるデータを、軸ｔと異なる軸ｆにおけるデータに変換する変換処理である。このようなｔ−ｆ変換の具体例としては、例えば、フーリエ変換が挙げられる。ｔ−ｆ変換の他の例としては、例えば、フーリエ変換の場合とは異なる正規直交関数で展開する場合や、例えば球面調和関数や球ベッセル関数で展開する場合が考えられる。なお、「ｆ」は、周波数（frequency）に対応し、フーリエ変換の場合、対応する正弦波の周波数と関連づけられる量である。

上記の時系列ｘ空間データ（ｘ−ｔ空間データ）に対して時間軸方向にｔ−ｆ変換を適用すると、ｘ−ｆ空間データが得られる。ｔ−ｆ変換として、フーリエ変換を採用した場合、ｘ−ｆ空間データは、ｘ−ｔ空間データがフーリエ変換後の係数値で表現されたものとなる。ｘ−ｆ空間データ内における各信号点は、ｘ軸位置とｆ軸位置とが定まれば、その位置が定まる。以下の実施形態では、ｔ−ｆ変換の例として、フーリエ変換を採用した場合を用いて説明する。

次に、強度画像及び位相画像について簡単に説明する。複素数であるＭＲ信号は、複素数の絶対値と、複素数平面での偏角とが存在する。前者を強度（振幅）、後者を位相と呼ぶ。ＭＲ信号の強度をｘ空間の位置ごとに求めた画像を強度画像と呼び、ＭＲ信号の位相をｘ空間の位置ごとに求めた画像を位相画像と呼ぶ。

次に、位相画像から流速を求める手順について簡単に説明する。位相画像から流速を求める手法の一つに、位相コントラスト（Ｐｈａｓｅ ｃｏｎｔｒａｓｔ）法がある。

図２は、位相コントラスト法を用いて流速を求める手順について簡単に説明した図である。位相コントラスト法において、重要な概念として、双極傾斜磁場の印加が挙げられる。双極傾斜磁場とは、空間的に場所によってその強さが変化する磁場（傾斜磁場）であって、その印加される方向が、時間によって、正方向・負方向に逆転するような磁場を言う。典型的には、正方向に印加される時間と負方向に印加される時間が等しく、また、印加される磁場の強度の絶対値が互いに等しいような双極傾斜磁場が使用される。

図２（ａ）は、位相コントラスト法において印加される「第１の双極傾斜磁場」の具体例を示した図である。この例では、グラフ６０１のように、時間Ａから時間Ｂまでの間は、正方向に傾斜磁場が印加される。また、グラフ６０２のように、時間Ｂから時間Ｃまでの間は、負方向に傾斜磁場が印加される。

ところで、ここで印加される磁場が、傾斜磁場であるから、空間的に磁場の強度は異なる。図２（ｃ）は、そのことを模式的に示した図である。横軸は、空間的位置を、縦軸は、その位置にかかるトータルの磁場の大きさを示す。例えば、グラフ６１０のように、正の傾斜磁場がかかっている場合には、位置が「Ｐ」「Ｑ」「Ｒ」「Ｓ」「Ｔ」と変化するにしたがって、磁場の大きさは大きくなる。逆に、グラフ６１１のように、負の傾斜磁場がかかっている場合には、位置が、「Ｐ」「Ｑ」「Ｒ」「Ｓ」「Ｔ」と変化するにしたがって、磁場の大きさは小さくなる。例えば、図２（ａ）の、時間Ａから時間Ｂまでの間では、グラフ６１０のような磁場が印加され、図２（ａ）の、時間Ｂから時間Ｃまでの間では、グラフ６１１のような磁場が印加される。

次に、静止した流体に、双極傾斜磁場を印加した場合のＭＲ信号について説明する。空間的な位置ｘにおけるＭＲ信号の位相θ（ｘ）は、点ｘにおける局所磁場をＢ（ｘ）、定数をγ、時間をｔとして、θ（ｘ）＝∫γＢ（ｘ）ｄｔで表される。従って、ＭＲ信号の位相θ（ｘ）は、局所磁場Ｂ（ｘ）に比例する。

次に、図２（ｃ）の、位置Ｐにおける静止した流体の位相変化について考える。位置Ｐで静止した流体には、時間Ａから時間Ｂの間にはグラフ６１０のように正の磁場がかかり、時間Ｂから時間Ｃの間には、グラフ６１１のように（一様磁場を基準にすると）負の磁場がかかる。しかし、時刻Ａから時刻Ｃの間で磁場の強度を積分すると「０」であるから、位置Ｐにおける静止した流体の位相変化はない。

次に、動いている流体の位相変化について考える。まず、図２（ｃ）において、「流速低」である状況を考える。この流体は、時間Ａにおいて位置Ｐにあり、時間Ｂにおいては位置Ｑに移動し、時間Ｃにおいては位置Ｒに移動している。ここで、時間Ａから時間Ｂにおいては、グラフ６１０のように、正の磁場がかかっており、時間Ｂから時間Ｃにおいては、グラフ６１１のように、負の磁場がかかっているが、傾斜磁場の大きさは、位置Ｐから離れたところほど大きい。したがって、この流体にとっては、時間Ｂから時間Ｃにおいて印加された磁場の効果のほうが、時間Ａから時間Ｂにおいて印加された磁場の効果より大きいため、位相変化は０にならず、有限の位相変化を持つことになる。

次に、図２（ｃ）において、「流速高」である状況を考える。この流体は、時間Ａにおいて、位置Ｐにあり、時間Ｂのとき、位置Ｒにあり、時間Ｃの時、位置Ｔにある。より遠くまで流体が進むほど、その間の傾斜磁場の差は大きくなるので、「流速低」の状況に比べて、時間Ａから時間Ｂまでの間に印加された磁場の効果と、時間Ｂから時間Ｃまでの印加された磁場の効果の差が、大きくなっている。したがって、「流速低」の状況に比べて、「流速高」の場合は、位相変化が大きくなる。

すなわち、双極傾斜磁場中で流体の位置が移動すると、その流体の移動速度に応じたＭＲ信号の位相変化を観測することができる。

実際には、例えば、図２（ａ）のように、第１の双極傾斜磁場を印加して第１のデータを作成し、図２（ｂ）のように、第１の双極傾斜磁場の極性を反転させた第２の双極傾斜磁場を印加して第２のデータを作成し、第１のデータと第２のデータとの差分データを作成してもよい。この差分データは、流速０のときデータ値０になり、また、流速が大きいほどデータ値が大きくなることから、この位相画像の差分データから、流速を求めることができる。そして、位相コントラスト法では、「第１の双極傾斜磁場」や第１の双極傾斜磁場の印加中に収集された時系列ｋ空間データから強度画像を生成することができる。

位相コントラスト法を適用することのできる部位の例としては、例えば大動脈（心臓）があげられるが、その他の部位についても適用することができる。例えば、位相コントラスト法は、肺動脈（肺）、頸動脈（頸部）などに適用することもできる。

しかしながら、例えば、位相コントラスト法で、強度画像の時系列画像と位相画像の時系列画像との両方を生成する場合、強度画像の時間変化ノイズを除去する程度の強いローパスフィルタを適用すると、流速を正しく求めることができない。また、ローパスフィルタが弱いと、強度画像の時間変化ノイズを除去できない。

そこで、第１の実施形態では、正しい流速データが得られ、かつ、強度画像の画質を向上するために、以下の処理が画像生成部１２２により行われる。変換部１２２ａは、時系列ｋ空間データを、ｘ空間上の時系列データである第１の時系列ｘ空間データに変換する変換処理を行う。そして、強度画像生成部１２２ｂは、第１の時系列ｘ空間データに第１のフィルタを適用して取得した第２の時系列ｘ空間データから、強度画像を生成する。そして、位相画像生成部１２２ｃは、第１の時系列ｘ空間データに、第１のフィルタとは異なる第２のフィルタを適用して取得した第３の時系列ｘ空間データから、位相画像を生成する。

第１の実施形態では、第１のフィルタ及び第２のフィルタそれぞれは、時間軸方向の軸ｔと異なる軸ｆの方向に関するフィルタである。強度画像生成部１２２ｂは、変換部１２２ａが出力した第１の時系列ｘ空間データに対して、ｔ−ｆ変換を適用してｘ−ｆ空間データを取得し、ｘ−ｆ空間データに対して、軸ｆの方向に第１のフィルタを適用した後に、ｔ−ｆ変換の逆変換を行うことで、第２の時系列ｘ空間データを取得する。また、位相画像生成部１２２ｃは、変換部１２２ａが出力した第１の時系列ｘ空間データに対して、ｔ−ｆ変換を適用してｘ−ｆ空間データを取得し、ｘ−ｆ空間データに対して、軸ｆの方向に第２のフィルタを適用した後に、ｔ−ｆ変換の逆変換を行うことで、第２の時系列ｘ空間データを取得する。

ここで、第１のフィルタ及び第２のフィルタそれぞれは、出力画像の種類に応じて設定される。例えば、強度画像生成部１２２ｂは、変換後時空間データである第１の時系列ｘ空間データに、時間変化の周波数成分の直流成分を含む第１の範囲に限定する第１のフィルタをかけ、時空間強度画像を得る。また、例えば、位相画像生成部１２２ｃは、変換後時空間データである第１の時系列ｘ空間データに、時間変化の周波数成分の直流成分を含む第１の範囲より高周波成分まで含む第２の範囲に限定する第２のフィルタをかけ、時空間位相画像を得る。

以下、第１の実施形態に係る画像生成部１２２が行う処理について、図３及び図４を用いて簡単に概略を説明する。以下では、「ｔ−ｆ変換」としてフーリエ変換を用いた後に、第１のフィルタ及び第２のフィルタとしてローパスフィルタを使う例を説明するが、後述するように実施形態はこれに限られない。

図３は、第１の実施形態における画像生成の概略を説明するための図であり、図４は、第１の実施形態における画像生成部１２２が行う処理の流れを説明するための図である。図４は、データの入出力形態に合わせて、図１に示す変換部１２２ａ、強度画像生成部１２２ｂ及び位相画像生成部１２２ｃの配置を変更した図である。

収集部１２６ａは、シーケンス制御部１１０からＭＲデータを取得して、時系列データｋ空間データ４００を収集する。変換部１２２ａは、収集部１２６ａから時系列データｋ空間データ４００を取得する。変換部１２２ａは、図３に示すように、時系列データｋ空間データ４００に対して、ｋ方向に逆フーリエ変換(ＩＦＦＴ：Inverse Fast Fourier Transform)を行うことにより、第１の時系列ｘ空間データ（第１のｘ−ｔ空間データ）を得る。この第１のｘ−ｔ空間データは複素数の値を持つ。したがって、この時系列データをもとに、強度画像４０１及び位相画像４０６を生成することができる。

ただし、図３に示すように、第１のｘ−ｔ空間データから生成した場合の強度画像４０１及び位相画像４０６にはノイズが含まれているので、これらのノイズを取り除く必要がある。強度画像生成部１２２ｂは、変換部１２２ａが変換した第１の時系列ｘ空間データ（第１のｘ−ｔ空間データ）をもとに、第１のフィルタを用いたフィルタ処理を行い、ノイズを除去した強度画像を生成する。位相画像生成部１２２ｃは、変換部１２２ａが変換した第１の時系列ｘ空間データ（第１のｘ−ｔ空間データ）をもとに、第１のフィルタとは異なるフィルタである第２のフィルタを用いてフィルタ処理を行い、ノイズを除去した位相画像を生成する。

はじめに、強度画像生成部１２２ｂ及び位相画像生成部１２２ｃは、第１の時系列ｘ空間データ（第１のｘ−ｔ空間データ）に対して、ｔに関してｔ−ｆ変換（逆フーリエ変換：ＩＦＦＴ）を行い、第１のｘ−ｆ空間データを得る。画像生成部１２２は、この第１のｘ−ｆ空間データをもとに、強度成分及び位相成分を算出することにより、強度画像４０２及び位相画像４０７を生成することができる。ただし、図３に示すように、第１のｘ−ｆ空間データから生成した強度画像４０２及び位相画像４０７にも、ノイズが含まれている。

次に、強度画像生成部１２２ｂは、第１のｘ−ｆ空間データに対して、軸ｆの方向に、第１のフィルタを適用する。位相画像生成部１２２ｃは、第１のｘ−ｆ空間データに対して、軸ｆの方向に、第１のフィルタとは異なるフィルタである第２のフィルタを適用する。第１のフィルタ及び第２のフィルタそれぞれは、時間方向の軸ｔと異なる軸ｆの方向に関するフィルタである。

流速は時々刻々と変化しているため、位相画像には、様々な周波数のモードが、複雑に寄与している。したがって、強すぎるフィルタをかけてしまうと、流速を正確に求めることができない。これに対して、強度画像においては、静的な情報（直流成分）及び、心臓の周期や呼吸の周期などに連動する、特定の周波数付近の情報が得られれば十分なことも多いので、強いフィルタをかけることにより、信号雑音強度比の良いデータを得ることができる。

このように、強度画像と位相画像とで、かけるべきフィルタの強さが違うため、強度画像と位相画像で、同一のフィルタを用いたのでは、精度のよい結果が得られない。したがって、本実施形態に係る画像処理装置としての画像生成部１２２は、強度画像と位相画像とで、異なるフィルタを用いる。特に、流速の時間変化を求めるためには、位相画像の周波数依存性を正確に算出することが必要になるため、例えば、位相画像生成に用いる第２のフィルタを、強度画像生成に用いる第１のフィルタに比べて弱いフィルタとすることができる。

例えば、第１のフィルタは、第２のフィルタより、データを平滑化する度合いの大きいフィルタとする。すなわち、第１のフィルタの例としては、例えば、図３に示すローパスフィルタ４０３のように、ｆの小さい領域のみ残し、ｆのほとんどの領域をカットするローパスフィルタが好適である。また、第２のフィルタとしては、例えばローパスフィルタ４０８のように、ｆのほとんどの領域を残し、ｆの高周波の領域のみをカットするローパスフィルタが好適である。

強度画像生成部１２２ｂは、第１のｘ−ｆ空間データに、例えばローパスフィルタ４０３を適用することで、強度画像用ｘ−ｆ空間データを生成する。強度画像用ｘ−ｆ空間データは複素数であるので、強度成分を算出することにより、強度画像４０４を生成することができる。また、位相画像生成部１２２ｃは、第１のｘ−ｆ空間のデータに、ローパスフィルタ４０８を適用することで、位相画像用ｘ−ｆ空間のデータを生成することができる。位相画像用ｘ−ｆ空間データは複素数であるので、位相成分を算出することにより、位相画像４０９を生成することができる。

そして、強度画像生成部１２２ｂは、強度画像用ｘ−ｆ空間データに対して、ｆ方向に関して前述のｔ−ｆ変換の逆変換（フーリエ変換、例えばＦＦＴ）を行うことで、ノイズが除去された第２のｘ−ｔ空間データを得ることができる。強度画像生成部１２２ｂは、この第２のｘ−ｔ空間データの強度成分を算出することにより、図３に示すように、ノイズ除去後の強度画像４０５を生成する。

また、位相画像生成部１２２ｃは、位相画像用ｘ−ｆ空間データに対して、ｆ方向に関して前述のｔ−ｆ変換の逆変換（フーリエ変換）を行うことで、ノイズが除去された第３のｘ−ｔ空間データを得ることができる。位相画像生成部１２２ｃは、第３のｘ−ｔ空間データの位相成分を算出することにより、ノイズ除去後の位相画像４１０を生成する。さらに、位相コントラスト法では、位相画像４１０と流速フローの間に関係があるので、位相画像生成部１２２ｃは、ノイズ除去後の位相画像４１０を用いて、図３に示すように、時間ｔの関数としての流速グラフ４１１を算出する。

なお、図３では、説明上、第１のｘ−ｔ空間データから生成した場合の強度画像４０１及び位相画像４０６と、第１のｘ−ｆ空間データから生成した場合の強度画像４０２及び位相画像４０７と、第２のｘ−ｆ空間データから生成した場合の強度画像４０４及び第３のｘ−ｆ空間データから生成した場合の位相画像４０９を示しているが、これら画像は、実際には生成されても、されなくてもよい。ただし、本実施形態は、例えば、強度画像４０５及び位相画像４１０とともに、強度画像４０１及び位相画像４０６を表示することで、第１のフィルタ及び第２のフィルタによるノイズ除去効果を操作者に確認させてもよい。同様に、本実施形態は、例えば、流速グラフ４１１とともに、位相画像４０６を用いて算出した流速グラフを表示することで、第１のフィルタ及び第２のフィルタによるノイズ除去効果を操作者に確認させてもよい。

続いて、図５を用いて、第１の実施形態に係るＭＲＩ装置１００が実行する処理について説明する。図５は、第１の実施形態に係るＭＲＩ装置の処理手順を示すフローチャートである。

図５に示すように、まず、制御部１２６が、操作者から撮像条件の入力を受け付ける（ステップＳ１０１）。続いて、収集部１２６ａが、操作者から入力された撮像条件に基づいてシーケンス情報を生成し、生成したシーケンス情報をシーケンス制御部１１０に送信することで、準備スキャンの実行を制御する（ステップＳ１０２）。ここで、準備スキャンには、例えば、位置決め用の画像を収集するスキャンや、静磁場の不均一性を補正するシミングスキャン、感度分布情報を収集するスキャン等が含まれる。

準備スキャンの終了後、続いて、収集部１２６ａは、操作者から入力された撮像条件に基づいてシーケンス情報を生成し、生成したシーケンス情報をシーケンス制御部１１０に送信することで、出力画像（例えば、診断用に出力される画像）を収集するイメージングスキャンの実行を制御する（ステップＳ１０３）。第１の実施形態において、収集部１２６ａは、イメージングスキャンにおいて、単一または複数チャネルについて、時系列ｋ空間データを収集する。例えば、第１の実施形態においては、収集部１２６ａは、データを間引かずにフルサンプリングで時系列ｋ空間データを収集する。

そして、画像生成部１２２は、ステップＳ１０３において収集され、記憶部１２３に格納された単一または複数チャネル分の時系列のｋ空間データを用いて、画像の生成を行う（ステップＳ１０４）。生成された画像は、必要に応じて、記憶部１２３に格納されたり、表示部１２５に表示されたりする。

次に、図５のステップＳ１０４において、第１の実施形態に係る画像生成部１２２が実行する処理の一例について説明する。図６は、第１の実施形態に係る画像生成部の処理手順の一例を示すフローチャートである。

図６に示すように、まず、変換部１２２ａが、単一または複数チャネル分の時系列のｋ空間データを記憶部１２３から読み出し、各チャネルの時系列のｋ空間データに対して、例えば、空間方向にフーリエ変換を適用することで、各チャネルのｘ−ｔ空間データを得る（ステップＳ２０１）。変換部１２２ａは、得られた各チャネルのｘ−ｔ空間データを、強度画像生成部１２２ｂと位相画像生成部１２２ｃに受け渡す。

具体的には、変換部１２２ａが、記憶部１２３から読みだす単一または複数チャネル分の時系列ｋ空間データは、例えば、ｋ空間のＲＯ（Read Out）方向、ｋ空間のＰＥ（Phase Encoding）方向及び時間それぞれが軸となる時系列データになっており、各々がチャネルの個数分だけ用意されている。この時系列データに対して、変換部１２２ａは、ｋ空間のＲＯ方向、ｋ空間のＰＥ方向の２軸についてフーリエ変換を実行して、水素分布画像のｘ方向、水素分布画像のｙ方向及び時間を各軸とする多次元データを出力する。

そして、強度画像生成部１２２ｂは、ステップＳ２０１で得られた各チャネルのｘ−ｔ空間データに対して、時間方向に例えばローパスフィルタをかけ、ｘ−ｔ空間強度画像を生成する（ステップＳ２０２）。具体的には、強度画像生成部１２２ｂは、ステップＳ２０１で得られた各チャネルのｘ−ｔ空間データ（第１のｘ−ｔ空間データ）の時間方向にフーリエ変換を適用することで、各チャネルのｘ−ｆ空間データを得る。続いて、強度画像生成部１２２ｂは、ｆ軸に対して範囲１のローパスフィルタをかけ、範囲１以外の信号点の信号値を、「ゼロ（０）」で置き換える。強度画像生成部１２２ｂは、ローパスフィルタ適用後のｘ−ｆ空間データに、逆フーリエ変換を適用し、ｘ−ｔ空間データ（第２のｘ−ｔ空間データ）を得る。各ｘ−ｔ空間データの強度を、各ｘ−ｔ空間強度画像の信号値として設定し、ｘ−ｔ空間強度画像を得る。

そして、ステップＳ２０２の後、又は、ステップＳ２０２と並行して、位相画像生成部１２２ｃは、ステップＳ２０１で得られた各チャネルのｘ−ｔ空間データに対して時間方向に例えばローパスフィルタをかけ、ｘ−ｔ空間位相画像を生成する（ステップＳ２０３）。具体的には、位相画像生成部１２２ｃは、ステップＳ２０１で得られた各チャネルのｘ−ｔ空間データ（第１のｘ−ｔ空間データ）の時間方向にフーリエ変換を適用することで各チャネルのｘ−ｆ空間データを得る。

続いて、位相画像生成部１２２ｃは、ｆ軸に対して範囲１より広い範囲２のローパスフィルタをかけ、範囲２以外の信号点の信号値を、「ゼロ（０）」で置き換える。ローパスフィルタ適用後のｘ−ｆ空間データに、逆フーリエ変換を適用し、ｘ−ｔ空間データ（第３のｘ−ｔ空間データ）を得る。各ｘ−ｔ空間データの位相を、各ｘ−ｔ空間位相画像の信号値として設定し、ｘ−ｔ空間位相画像を得る。

上記の範囲１及び範囲２との関係は、図３に示すローパスフィルタ４０３及びローパスフィルタ４０８にも示されているが、改めて、図７を用いて説明する。図７は、第１の実施形態におけるローパスフィルタの一例を示す図である。ｘ−ｔ空間強度画像は、被検体Ｐの動きを観察する目的で取得される画像であり、例えば、心臓の周期的な拍動の様子などを観察する際に取得される。そのため、範囲１は、図７に示すように、直流成分と被検体Ｐの動きがもつ周波数を含む範囲に設定する必要がある。一方、ｘ−ｔ空間位相画像は、血流などの流速を観察する目的で取得される画像であり、例えば大動脈起始部など画像の局所領域内の信号の時間変化から流速値の時間変化を求めるために利用される場合もある。流速値には高周波成分も含まれるため、範囲２は流速値に影響を与えない程度の範囲に設定する必要がある。このため、範囲２は、図７に示すように、範囲１＜範囲２となるように設定される。

なお、本実施形態は、範囲１＜範囲２となる範囲内で、範囲１及び範囲２は、被検体ごとに変化させてもよい。これについては、後に詳述する。また、範囲２については全周波数帯域を設定してもよく、その場合にはステップＳ２０３でローパスフィルタは適用しない場合と等価である。

なお、本実施形態は、前述した範囲１及び範囲２のようにカットオフ周波数が設定されたローパスフィルタを用いる場合について説明したが、本実施形態は、部分的な範囲を指定したバンドパスフィルタが第１のフィルタ及び第２のフィルタとして設定される場合であってもよい。その場合にも、範囲２は範囲１よりも高周波成分まで含む必要がある。

上述してきたように、第１の実施形態によれば、出力画像の種類に応じて、位相画像生成用の第２のフィルタに設定する範囲２を強度画像生成用の第１のフィルタに設定する範囲１よりも広く設定しておき、ｘ−ｔ空間データに対して、第１のフィルタ及び第２のフィルタを適用することで、強度画像及び位相画像が生成される。これにより、第１の実施形態では、正しい流速データが得られ、かつ、強度画像の画質を向上することができる。

（第１の実施形態の第１の変形例）
上記の第１の実施形態では、一例として、フーリエ変換及びローパスフィルタを用いたフィルタ処理が画像生成部１２２により行われる場合について説明した。すなわち、上記の第１の実施形態では、画像生成部１２２は、ｋ−ｔ空間データを変換した第１のｘ−ｔ空間データ（第１の時系列ｘ空間データ）をフーリエ変換して第１のｘ−ｆ空間データを取得する。そして、画像生成部１２２は、第１のｘ−ｆ空間データに対して強度画像用のローパスフィルタを掛けて強度画像用ｘ−ｆ空間データを取得し、第１のｘ−ｆ空間データに対して位相画像用のローパスフィルタを掛けて位相画像用ｘ−ｆ空間データを取得する。そして、画像生成部１２２は、強度画像用ｘ−ｆ空間データをフーリエ（逆）変換して、第２のｘ−ｔ空間データ（第２の時系列ｘ空間データ）を取得し、位相画像用ｘ−ｆ空間データをフーリエ（逆）変換して、第３のｘ−ｔ空間データ（第３の時系列ｘ空間データ）を取得する。そして、画像生成部１２２は、第２のｘ−ｔ空間データからノイズが除去された時系列強度画像を生成し、第３のｘ−ｔ空間データからノイズが除去された時系列位相画像を生成する。しかし、第１の実施形態は、これに限られない。例えば、画像生成部１２２が用いる第１のフィルタ及び第２のフィルタそれぞれは、時間軸方向に平滑化を行うフィルタ（空間平滑化フィルタ）を適用することもできる。

具体的には、第１の実施形態の第１の変形例では、強度画像生成部１２２ｂは、第１の時系列ｘ空間データに対して、時間方向に第１のフィルタを適用して、第２の時系列ｘ空間データを取得する。また、第１の実施形態の第１の変形例では、位相画像生成部１２２ｃは、第１の時系列ｘ空間データに対して、時間方向に第２のフィルタを適用して、第３の時系列ｘ空間データを取得する。

第１の実施形態の第１の変形例においては、ステップＳ２０２及びステップＳ２０３において、ｘ−ｆ空間データへの変換を行わず、第１のｘ−ｔ空間データに対して、時間方向に２種類の空間平滑化フィルタを適用することで、直接、ノイズ除去後のｘ−ｔ空間データ（第２のｘ−ｔ空間データ及び第３のｘ−ｔ空間データ）を生成する。

第１の実施形態の第１の変形例で用いる空間平滑化フィルタの一例としては、例えば移動平均法による平滑化フィルタが挙げられる。図８は、第１の実施形態の第１の変形例における空間ローパスフィルタの一例を示す図である。図８では、空間平滑化フィルタの一例として、移動平均法による平滑化フィルタの一例を示している。図８の左図は、２次元空間上において、前後左右の９マスの値を基に移動平均を取る平滑化フィルタの例を表す。各正方形は、データ点を表し、各正方形の中の数字は、それぞれのデータ点にかけられる重みを表す。例えば、図８の左のフィルタでは、ｘ座標及びｙ座標が（ｘ、ｙ）のデータ点の平滑化後の値は、ｘ座標及びｙ座標が（ｘ−１、ｙ−１）、（ｘ−１、ｙ）、（ｘ−１、ｙ＋１）、（ｘ、ｙ−１）、（ｘ、ｙ）、（ｘ、ｙ＋１）、（ｘ＋１、ｙ−１）、（ｘ、ｙ）、（ｘ、ｙ＋１）の９個のデータ点を基に平滑化される。具体的には、それら９個の点それぞれのデータ点の値に、それらの重みである「１／９」を乗じて９個の点について和をとったものが、平滑化後の（ｘ、ｙ）の値になる。

図８の右図は、前後左右２５マスについて、同様に平滑化するフィルタを示したものである。例えば、図８の左のフィルタでは、（ｘ、ｙ）のデータ点の平滑化後の値は、（ｘ、ｙ）のデータ点と、（ｘ、ｙ）のデータ点の周囲の２４個のデータ点とを基に平滑化される。具体的には、それら２５個の点それぞれのデータ点の値に、それらの重みである「１／２５」を乗じて２５個の点について和をとったものが、平滑化後の（ｘ、ｙ）の値になる。ここで、図８の左のフィルタと右のフィルタの違いを考察する。今、平滑化の対象となる点が（ｘ、ｙ）であり、その他のデータ点の値がすべて「０」であったとすると、点（ｘ、ｙ）の、平滑化後の値は、左のフィルタを適用したあとにおいては元のデータの１／９であり、右のフィルタを適用したあとにおいては、元のデータの１／２５となる。このことから、右のフィルタのほうが、データ点が強く潰れてしまうので、左のフィルタに比べて、強い平滑化フィルタ（データを平滑化する度合いの大きいフィルタ）である。

また、このフィルタは、実際には時間方向（１次元）に対して適用するので、図８の左図のフィルタの重みは、例えば、「１／３，１／３，１／３」、図８の右の図のフィルタの重みは、例えば、「１／５，１／５，１／５，１／５，１／５」となる。

例えば、第１の実施形態の第１の変形例では、強度画像生成用には強い平滑化フィルタである図８の右のフィルタを適用し、位相画像生成用には弱い平滑化フィルタである図８の左のフィルタを適用することで、正しい流速データが得られ、かつ、強度画像の画質を向上することができる。なお、第１の実施形態の第１の変形例では、近接する画素からの寄与の重みを大きく取る平滑化フィルタ、例えばガウシアンフィルタを適用してもよい。

また、スパイクノイズ（周辺の値から大きな飛びがあるようなノイズ）に対しては、移動平均フィルタや、ガウシアンフィルタ、ローパスフィルタに比べて、非線形フィルタを適用することが有利になる。このような場合においては、例えば、メディアンフィルタのような非線形フィルタを適用してもよい。

上述してきたように、第１の実施形態の第１の変形例では、時系列ｘ空間データの時間方向に直接フィルタを適用することにより、正しい流速データが得られ、かつ、強度画像の画質を向上することができる。また、第１の実施形態の第１の変形例では、例えば、非線形フィルタを用いることにより、ローパスフィルタが得意でないある種の雑音、例えば周波数によらず一定の強度をもつ白色雑音や、ｘ空間上で特異的な値をもつスパイクノイズの除去を効果的に行うことができる。なお、第１の実施形態の第１変形例では、第１のｘ−ｔ空間データに対して、位相画像生成用には平滑化フィルタを適用せずに、第１のｘ−ｔ空間データから時系列位相画像が生成する場合であってもよい。

（第１の実施形態の第２の変形例）
これまでの実施形態では、データを間引かずにフルサンプリングで収集された時系列ｋ空間データを基に処理を行う場合を説明したが、第１の実施形態はこれに限られない。すなわち、画像生成部１２２は、シーケンス制御部１１０は、データを間引きながらサンプリングされた時系列ｋ空間データを基に処理を行う場合でも、以下に説明する処理により、正しい流速データが得られ、かつ、強度画像の画質を向上することができる。

ここで、データを間引きながらサンプリングを行う場合としては、パラレルイメージングのように複数のコイルを用いてデータを間引きながら収集を行う場合であってもよいし、パラレルイメージングを用いず、一つのコイルを用いてデータを間引きながら収集を行う場合であってもよい。

間引かれたデータを取り扱う場合、大きくわけて二つの方法が存在する。第１の方法としては、間引かれたデータをもとにデータを、強度画像及び位相画像での再構成手法で時系列ｘ空間データを生成し、それに対して、第１の実施形態で説明した画像処理方法を適用し、強度画像及び位相画像で異なるフィルタを用いる方法である。第２の方法としては、間引かれたデータを再構成する処理の中で、強度画像及び位相画像で異なる再構成処理を用いる方法である。この場合、具体的には、強度画像及び位相画像で、異なる再構成領域で再構成を行う。本変形例では、前者の方法について説明し、後者の方法については、第２の実施形態で詳細に説明する。

具体的には、第１の実施形態の第２の変形例では、収集部１２６ａが収集する時系列ｋ空間データ（ｋ−ｔ空間データ）は、ｋ空間上で第１の間隔で時系列に沿って収集されたｋ空間データとなる。そして、第１の実施形態の第２の変形例では、変換部１２２ａは、ｋ空間上で第１の間隔で時系列に沿って収集された時系列ｋ空間データに対してフーリエ変換を行って、第４の時系列ｘ空間データを生成する。そして、変換部１２２ａは、第４の時系列ｘ空間データを所定の方法で再構成することで、第１の間隔より狭い第２の間隔でｋ空間上で収集された場合における時空間ｋ空間データに対応する第５の時系列ｘ空間データを生成する。すなわち、第５の時系列ｘ空間データは、アンフォールド処理により第４の時系列ｘ空間データから再構成される。そして、変換部１２２ａは、第５の時系列ｘ空間データを第１の時系列ｘ空間データとして取得し、強度画像生成部１２２ｂ及び位相画像生成部１２２ｃに出力する。強度画像生成部１２２ｂは、第１の時系列ｘ空間データに対して、第１の実施形態、又は、第１の実施形態の第１の変形例で説明したフィルタ処理を行なって、強度画像を生成する。また、位相画像生成部１２２ｃは、第１の時系列ｘ空間データに対して、第１の実施形態、又は、第１の実施形態の第１の変形例で説明したフィルタ処理を行なって、位相画像を生成する。

つまり、第１の実施形態の第２の変形例で用いられる時系列ｋ空間データ（ｋ−ｔ空間データ）は、「Reduction factor R」に従って、データを間引きながら収集されたｋ空間データである。すなわち、第１の実施形態の第２の変形例で用いられる時系列ｋ空間データの個数は、フルサンプリングに比べて１／Ｒとなる。

変換部１２２ａは、各チャネルのｋ−ｔ空間データに対して、空間方向にフーリエ変換を適用することで、各チャネルのｘ−ｔ空間データ（第４の時系列ｘ空間データ）を得る。ここで得られるｘ−ｔ空間データ（第４の時系列ｘ空間データ）は、折り返しを含むデータであるため、所定の再構成方法として、例えばパラレルイメージングの一種であるＳＥＮＳＥ（Ｓｅｎｓｉｔｉｖｉｔｙ ｅｎｃｏｄｉｎｇ）によるアンフォールド処理により、アンフォールドされたｘ−ｔ空間データ（第５の時系列ｘ空間データ）を再構成し、第５の時系列ｘ空間データを第１の時系列ｘ空間データとして取得する。例えば、ＳＥＮＳＥにおいては、準備スキャン等によって収集された各チャネルごとの感度分布情報と、折り返しを含むデータを基に、「折り返しを含むデータは、折り返しのないデータに対して、チャネルの感度分布情報と、折り返しを含む各チャネルにおけるデータとの積和演算を施して得られた画像である」と考えて、折り返しのないデータを推定する。また、所定の再構成方法として用いるアンフォールド処理の方法としては、ＳＥＮＳＥ以外のアンフォールド手法でもよい。また、パラレルイメージングのアンフォールド手法でなくともよく、例えば１チャネルで再構成する方法でもよい。また、当該方法をパラレルイメージングと組み合わせた方法でもよい。

また、第１の実施形態の第２の変形例は、ｋ空間上で第１の間隔で時系列に沿って収集された時系列ｋ空間データに対して、以下に説明する処理により、第１の時系列ｘ空間データを取得してもよい。例えば、第１の実施形態の第２の変形例では、データを間引きながら収集されたｋ―ｔ空間データから折り返しを含むｘ−ｆ空間データを生成し、折り返しを含むｘ−ｆ空間データに対して、例えば、ＳＥＮＳＥを拡張した手法の１種であるｋ−ｔ ＳＥＮＳＥ（ｋ−ｔ ｓｅｎｓｉｔｉｖｉｔｙ ｅｎｃｏｄｉｎｇ）などのアンフォールド処理を用いた再構成を行ってもよい。ｋ−ｔ ＳＥＮＳＥでは、時系列ｋ空間データに対して時間軸方向にフーリエ変換を適用した空間上で、ＳＥＮＳＥによるアンフォールド処理を実行する。

この場合、第１の実施形態の第２の変形例では、変換部１２２ａは、ｋ空間上で第１の間隔で時系列に沿って収集された時系列ｋ空間データに対してフーリエ変換を行うとともに、上述した「ｔ−ｆ変換」を適用することで、第１のｘ−ｆ空間データを生成する。そして、変換部１２２ａは、第１のｘ−ｆ空間データを所定の方法で再構成することで、第１の間隔より狭い第２の間隔でｋ空間上で収集された場合における時系列ｋ空間データに対応する第２の時系列ｘ−ｆ空間データを生成する。そして、変換部１２２ａは、第２のｘ−ｆ空間データにｔ−ｆ変換の逆変換を適用することで、第１の時系列ｘ空間データを取得し、強度画像生成部１２２ｂ及び位相画像生成部１２２ｃに出力する。

つまり、変換部１２２ａは、各チャネルの時系列ｋ−ｔ空間データに対して、空間方向にフーリエ変換を適用するとともに時間軸方向にフーリエ変換を適用することで、各チャネルのｘ−ｆ空間データ（第１のｘ−ｆ空間データ）を取得する。そして、変換部１２２ａは、第１のｘ−ｆ空間データから、所定の再構成方法として、例えばｋ−ｔＳＥＮＳＥによるアンフォールド処理により、アンフォールドされたｘ−ｆ空間データ（第２のｘ−ｆ空間データ）を再構成し、時間軸方向に逆フーリエ変換を適用してｘ−ｔ空間データ（第１の時系列ｘ空間データ）を取得する。

なお、第１の実施形態の第２の変形例で用いられる再構成方法は、ＳＥＮＳＥやｋ−ｔＳＥＮＳＥに限定されるものではなく、例えば、事前知識を利用した正則化（Regularization）や圧縮センシング（Compressed Sensing）といったその他の手法を用いてｘ−ｔ空間データを再構成してもよい。

以上のように、第１の実施形態の第２の変形例によれば、撮像時間を短縮するためにデータを間引きながらサンプリングされた時系列ｋ空間データが収集される場合であっても、正しい流速データが得られ、かつ、強度画像の画質を向上することができる。

（第１の実施形態の第３の変形例）
第１の実施形態の第３の変形例では、撮像時間を短縮するためにデータを間引きながらサンプリングされた時系列ｋ空間データが収集される場合に、第２の変形例とは異なる方法で、間引かれていない時系列ｘ空間データ（第１の時系列ｘ空間データ）を取得する方法について説明する。

第１の実施形態の第３の変形例では、変換部１２２ａは、時系列ｋ空間データがｋ空間上で第１の間隔で時系列に沿って収集されたｋ空間データである場合、第１の間隔より狭い第２の間隔でｋ空間上で収集された場合における時系列ｋ空間データを推定し、当該推定した時系列ｋ空間データから、第１の時系列ｘ空間データを取得する。

すなわち、第１の実施形態の第３の変形例では、例えば、間引かれたデータを補間することで、間引かれていない時系列ｘ空間データを単純に推定し、推定した時系列ｘ空間データを、第１の時系列ｘ空間データとして、強度画像生成部１２２ｂ及び位相画像生成部１２２ｃに出力する。

例えば、第１の実施形態の第３の変形例では、変換部１２２ａは、各チャネルのサンプリングしたｋ−ｔ空間データから、各チャネルのサンプリングしていないｋ−ｔ空間データを推定し、フルサンプリング相当のｋ−ｔ空間データを得る。そして、変換部１２２ａは、例えば、フルサンプリング相当のｋ−ｔ空間データに対して空間方向にフーリエ変換を適用することで、各チャネルのｘ−ｔ空間データを得る。

以上のように、第１の実施形態の第３の変形例でも、撮像時間を短縮するためにデータを間引きながらサンプリングされた時系列ｋ空間データが収集される場合であっても、正しい流速データが得られ、かつ、強度画像の画質を向上することができる。

（第１の実施形態の第４の変形例）
第１の実施形態及び第１の実施形態の第１の変形例〜第３の変形例では、例えば、ユーザや、ＭＲＩ装置１００のユーザにより、撮像前に範囲１及び範囲２が設定されることで、特性が決定された第１のフィルタ及び第２のフィルタが用いられることを前提として説明した。かかる手動設定は、予め初期設定されたパラメータ値が用いられる場合でも、被検体Ｐごとに行われる場合であってもよい。しかし、第１のフィルタの特性と、第２のフィルタの特性とは、以下に説明する様々な方法により、自動的に決定される場合であってもよい。

第１の実施形態の第４の変形例では、制御部１２６は、所定の情報を取得し、取得した所定の情報に基づいて、第１のフィルタの特性を定める第１のパラメータ及び第２のフィルタの特性を定める第２のパラメータを決定する。そして、制御部１２６は、第１のパラメータを強度画像生成部１２２ｂに通知し、第２のパラメータを位相画像生成部１２２ｃに通知する。

例えば、制御部１２６は、上記の所定の情報として、撮像条件に関する情報及び撮像される被検体Ｐに関する情報の少なくとも一方の情報を取得して、第１のパラメータ及び第２のパラメータを決定する。ここで、撮像条件に関する情報とは、例えば印加するパルスの種類、パルスシーケンス、撮像部位及び撮像部位内における関心領域の設定などの情報である。また、撮像される被検体Ｐに関する情報とは、例えば、被検体Ｐの撮影直前に測定された血圧や、心拍数、性別、身長、体重などの情報である。

例えば、ＭＲＩ装置１００のユーザは、撮像部位（撮像対象の臓器）に応じて、強度画像及び位相画像それぞれの画質を向上させることができる第１のパラメータの値及び第２のパラメータの値を経験的に求めておき、求めた値を、撮像部位に対応付けて、予め、記憶部１２３に記憶させておく。そして、制御部１２６は、例えば、検査情報として入力された撮像部位に対応する第１のパラメータの値及び第２のパラメータの値を、記憶部１２３から取得し、強度画像生成部１２２ｂ及び位相画像生成部１２２ｃに通知する。

また、記憶部１２３に格納される第１のパラメータ及び第２のパラメータの値は、上記の一例の他に、例えば、「印加パルスの種類」、又は、「パルスシーケンス」、又は、「撮像部位及び印加パルスの種類」、又は、「撮像部位、パルスシーケンス」、又は、「撮像部位、印加パルスの種類及びパルスシーケンス」等、撮像条件に関する情報に含まれる複数の情報から選択された任意数の情報に応じて、予め設定される場合であってもよい。

また、記憶部１２３に格納される第１のパラメータ及び第２のパラメータの値は、上記の一例の他に、被検体Ｐの撮影直前に測定された各種測定値から選択された任意数の情報に応じて、予め設定される場合であってもよい。また、記憶部１２３に格納される第１のパラメータ及び第２のパラメータの値は、撮像条件に関する情報に含まれる複数の情報と、被検体Ｐの撮影直前に測定された各種測定値とを組み合わせた任意数の情報に応じて、予め設定される場合であってもよい。

或いは、第１の実施形態の第４の変形例では、撮像部位に応じてフィルタのパラメータの値を変化させる場合、例えば、撮像部位（撮像対象の臓器）を流動する流体（例えば、血液）のピークフロー（最大流速）に応じて、第１のパラメータ及び第２のパラメータを決定してもよい。例えば、撮像部位が、大動脈（心臓）であるか、肺動脈（肺）であるか、頸動脈（頸部）であるかによって、血管の太さや心臓からの距離が異なるので、時間の関数としての流速の最大値であるピークフロー（最大流速）の大きさが異なる。そこで、制御部１２６は、撮像部位のピークフローに応じて、第１のフィルタ及び第２のフィルタの強さを変化させる。例えば、流速や最大流速が大きい場合には、範囲２を狭めることで、強めの第２のフィルタにし、逆に流速や最大流速が小さい場合には、範囲２を広げることで、弱めの第２のフィルタとすることが好適である。なお、第１のフィルタについては、撮像部位のピークフローに応じて、フィルタ特性を変化させてもよいし、変化させなくてもよい。

例えば、撮像部位と、撮像部位の平均的なピークフローの値と、第１のパラメータの値及び第２のパラメータの値とが対応付けられた設定情報を記憶部１２３に記憶させておく。制御部１２６は、被検体Ｐの撮像部位に対応する第１のパラメータの値及び第２のパラメータの値を取得することで、第１のフィルタ及び第２のフィルタの特性を決定する。

また、上記では、予め設定されたフィルタ特性の設定情報を用いて、制御部１２６がフィルタ特性の自動決定を行う場合について説明したが、本変形例は、制御部１２６が以下に例示する方法に基づいて、第１のパラメータの値及び第２のパラメータの値を算出する場合であってもよい。

例えば、制御部１２６は、ドップラー計測で計測された被検体Ｐの撮像部位の平均流速の値を取得する。或いは、例えば、制御部１２６は、血圧及び心拍数の情報をもとに、被検体Ｐの撮像部位の平均流速の値を推定する。

例えば、「血圧＝心拍数×１回拍出量×血管抵抗」という公式が成りたつので、制御部１２６は、被検体Ｐの血圧と心拍数とから、１回拍出量を算出する。「血量＝平均流速×血管の断面積」であるので、制御部１２６は、１回拍出量から求めた血量と、撮像部位の平均的な血管の断面積とから、平均流速を推定する。ここで、位相画像生成において、残す必要があるデータの周波数の値と、平均流速との関係式は、経験的に得ることができる。制御部１２６は、かかる関係式と、被検体Ｐの撮像部位の平均流速とから、第２のフィルタの特性を定める第２のパラメータの値を算出することができる。

或いは、例えば、制御部１２６は、第１のフィルタ処理を行なわずに生成された強度画像（以下、第１の強度画像）及び第２のフィルタ処理を行なわずに生成された位相画像（以下、第１の位相画像）を元に、第１のパラメータの値と第２のパラメータの値とを算出してもよい。かかる場合、第１段階として、強度画像生成部１２２ｂは、第１の時系列ｘ空間データから強度情報を生成する。また、第１段階として、位相画像生成部１２２ｃは、第１の時系列ｘ空間データから位相情報を生成する。そして、制御部１２６は、上記の所定の情報として、強度情報及び位相情報を取得して、第１のパラメータ及び第２のパラメータを決定する。

例えば、制御部１２６は、第１の強度画像及び第１の位相画像それぞれのノイズの大きさを推定する。そして、例えば、制御部１２６は、第１の強度画像のノイズの大きさから第１のパラメータの値を算出し、第１の位相画像のノイズの大きさから第２のパラメータの値を算出する。なお、本変形例は、制御部１２６の指示により、位相画像生成部１２２ｃが、第１の位相画像から流速の時系列データを生成し、制御部１２６が流速の時系列データのノイズの大きさから、第２のパラメータの値を算出してもよい。

ここで、制御部１２６は、強度画像に関しては、例えばＭＲ信号が小さい場所におけるノイズの大きさから、ノイズの大きさを推定することができる。また、制御部１２６は、位相画像に関しては、例えば流速の時系列データにおいて、時間の関数として滑らかに変化する成分を除去したあとのノイズの大きさから、ノイズの大きさを推定することができる。

上述してきたように、第１の実施形態の第４の変形例では、正しい流速データが得られ、かつ、強度画像の画質を向上することができるフィルタのパラメータの値を自動的に算出することにより、ユーザの利便性を向上することができる。

（第２の実施形態）
第２の実施形態では、第１の実施形態の第２の変形例において記載した第２の方法（間引かれたデータを再構成する処理の中で、強度画像及び位相画像で異なる再構成処理を用いる方法）について、図９等を用いて説明する。図９は、第２の実施形態に係る画像生成部の構成を示す機能ブロック図である。

第２の実施形態に係るＭＲＩ装置１００は、図１に示す第１の実施形態に係るＭＲＩ装置１００と比較して、画像生成部１２２が、図９に示すように、変換部１２２ａに換えて変換部１２２ｄを有し、強度画像生成部１２２ｂに換えて強度画像再構成部１２２ｅを有し、位相画像生成部１２２ｃに換えて位相画像再構成部１２２ｆを有する点が異なる。

変換部１２２ｄは、時系列ｋ空間データ（間引かれたｋ空間で収集され、時系列に沿って収集されたｋ空間データ）に対して、空間上の位置を示すｘの空間であるｘ空間上の時系列データである時系列ｘ空間データに変換する変換処理、及び、ｔ−ｆ変換を行う変換処理を適用して、第１のｘ−ｆ空間データを生成する。例えば、上記の時系列ｋ空間データは、「Reduction factor R」に従って、データを間引きながら時系列に沿って収集されたｋ空間データである。

そして、強度画像再構成部１２２ｅは、軸ｆの方向において第１の再構成範囲にある第１のｘ−ｆ空間データの各データ点を再構成の対象範囲として、所定の画像再構成手法により、第１のｘ−ｆ空間データがアンフォールディングされた第２のｘ−ｆ空間データを再構成する。そして、強度画像再構成部１２２ｅは、第２のｘ−ｆ空間データにｔ−ｆ変換の逆変換を適用して得たｘ−ｔ空間データから、強度画像を生成する。

また、位相画像再構成部１２２ｆは、軸ｆの方向において第１の再構成範囲と異なる再構成範囲である第２の再構成範囲にある第１のｘ−ｆ空間データの各データ点を再構成の対象範囲として、所定の画像再構成手法により、第１のｘ−ｆ空間データがアンフォールディングされた第３のｘ−ｆ空間データを再構成する。そして、位相画像再構成部１２２ｆは、第３のｘ−ｆ空間データにｔ−ｆ変換の逆変換を適用して得たｘ−ｔ空間データから、位相画像を生成する。

まず、所定の再構成手法として、ｋ−ｔ ＳＥＮＳＥを用いた場合を例にとり、強度画像及び位相画像で異なる再構成処理を用いる点について簡単に説明する。もっとも、再構成手法はｋ−ｔ ＳＥＮＳＥに限られるものではない。また、複数のコイルを用いてデータを収集するパラレルイメージングに限られるものでもなく、一つのコイルを用いてデータを間引いて収集したものに対して、所定の再構成手法を用いてもよい。

高速撮像を行うため、複数のコイルでデータを収集する代わりに、ｋ空間上のデータを間引いて収集して撮像を行うことがある。この場合、ｋ空間上のデータ点がそのままでは不足しているが、複数のコイルでデータ収集を行っているため、ＳＥＮＳＥなどの所定の再構成手法を用いることで、ｋ空間上のデータを間引かずに収集した場合におけるデータに対応するデータを再構成することができる。時系列ｋ空間データに対しても、例えばｘ−ｆ空間で再構成を行うｋ−ｔＳＥＮＳＥと呼ばれる再構成手法がある。

ところで、画像の再構成を行うにあたっては、一般に再構成誤差が生じる。例えば、信号強度が低く、ノイズが大きいデータ点を再構成に使用すると再構成誤差が生じる。また、系のｇ−ｆａｃｔｏｒ（geometry factor）が大きい場合にも再構成誤差が生じる。このように、再構成においては一般に再構成誤差を伴うので、全てのデータ点を用いて画像再構成を行うよりも、一部のデータ点のみを用いて画像再構成を行った方が、結果として出力結果の画像が、ノイズの少なくなる場合が存在する。

実際、ｘ−ｆ空間でｋ−ｔ ＳＥＮＳＥにより画像再構成を行う場合、周波数ｆの小さい範囲のみを、画像再構成の対象範囲とすることで、画像再構成誤差を小さく抑えることができる場合がある。

ここで、どのようなｆの領域を、画像再構成の対象範囲にするのが良いかは、画像再構成対象の性質によって異なる。画像再構成の対象範囲が狭すぎると、画像再構成の対象から外されることによって新たに生じる誤差が大きくなる。一方、画像再構成の対象範囲が広すぎると、画像再構成誤差が大きくなる。この両者を考慮して好適な画像再構成範囲が決定されるが、当該好適な画像再構成範囲は、画像再構成となる画像の性質及び利用目的に応じて異なる。

従って、第２の実施形態では、強度画像生成のための画像再構成処理において、第１の画像再構成範囲（以下、範囲３とも記載する）のデータに対して画像再構成を行い、位相画像生成のための画像再構成処理において、第１の画像再構成範囲とは異なる画像再構成範囲である第２の画像再構成範囲（以下、範囲４とも記載する）のデータに対して画像再構成を行う。第２の実施形態では、第２の画像再構成範囲（範囲４）を、第１の画像再構成範囲（範囲３）より広くとるのが好適である。また、例えば、第２の実施形態では、第１の画像再構成範囲及び第２の画像再構成範囲が、ｆに関して最初の折り返しが生じる手前の範囲にしてもよい。また、例えば、第２の画像再構成範囲を、ｆに関して全範囲としてもよい。

次に、図１０のフローチャートを用いて、第２の実施形態に係る画像生成部の処理内容の一例について、説明する。図１０は、第２の実施形態に係る画像生成部の処理手順の一例を示すフローチャートである。なお、図１０に示すフローチャートは、図５に示すステップＳ１０４において、第２の実施形態に係る画像生成部１２２が行う処理に対応する。

図１０に示すように、まず、変換部１２２ｄが、単一または複数チャネル分の時系列ｋ空間データを記憶部１２３から読み出し、各チャネルの時系列ｋ空間データに対して、時空間方向にフーリエ変換を適用することで、各チャネルのｘ−ｆ空間データ（第１のｘ−ｆ空間データ）を得る（ステップＳ３０１）。すなわち、変換部１２２ｄは、空間方向にフーリエ変換を適用するとともに、時間軸方向にフーリエ変換を適用することで、各チャネルのｘ−ｆ空間データ（第１のｘ−ｆ空間データ）を得る。変換部１２２ｄは、得られた各チャネルのｘ−ｆ空間データ（第１のｘ−ｆ空間データ）を、強度画像再構成部１２２ｅと位相画像再構成部１２２ｆとに受け渡す。

強度画像再構成部１２２ｅは、ステップＳ３０１で得られた各チャネルのｘ−ｆ空間データ（第１のｘ−ｆ空間データ）に対して、ｆ軸の周波数帯域の範囲３に限定した再構成処理を行い、ｘ−ｔ空間強度画像を生成する（ステップＳ３０２）。

具体的には、強度画像再構成部１２２ｅは、ステップＳ３０１で得られた各チャネルのｘ−ｆ空間データ（第１のｘ−ｆ空間データ）に対して、ｋ−ｔ ＳＥＮＳＥなどの時空間アンフォールド処理を適用する。このとき、強度画像再構成部１２２ｅは、時間変化の周波数成分であるｆ軸で再構成する範囲を、直流成分を含む範囲３（第１の再構成範囲）に限定する。そして、強度画像再構成部１２２ｅは、アンフォールドされたｘ−ｆ空間データ（第２のｘ−ｆ空間データ）に対し、時間軸方向に逆フーリエ変換を適用してｘ−ｔ空間強度画像を得る。

そして、ステップＳ３０２の後、又は、ステップＳ３０２と並行して、位相画像再構成部１２２ｆは、ステップＳ３０１で得られた各チャネルのｘ−ｆ空間データ（第１のｘ−ｆ空間データ）に対して、時間軸方向（軸ｔの方向）の周波数帯域をｆ軸の周波数帯域の範囲４に限定した再構成処理を行い、ｘ−ｔ空間位相画像を生成する（ステップＳ３０３）。

具体的には、位相画像再構成部１２２ｆは、ステップＳ３０１で得られた各チャネルのｘ−ｆ空間データ（第１のｘ−ｆ空間データ）に対して、ｋ−ｔＳＥＮＳＥなどの時空間アンフォールド処理を適用する。このとき、位相画像再構成部１２２ｆは、時間変換の周波数成分であるｆ軸で再構成する範囲を、直流成分を含み範囲３（第１の再構成範囲）より広い範囲４（第２の再構成範囲）に限定する。そして、位相画像再構成部１２２ｆは、アンフォールドされたｘ−ｆ空間データ（第３のｘ−ｆ空間データ）に対し、時間軸方向に逆フーリエ変換を適用してｘ−ｔ空間位相画像を得る。

範囲３（第１の再構成範囲）と範囲４（第２の再構成範囲）との関係は、第１の実施形態における範囲１と範囲２との関係と同様であり、両者の関係は範囲３＜範囲４とするのが好適である。また、範囲４については全周波数帯域を設定してもよく、その場合にはステップＳ３０３ですべての周波数帯域について再構成処理を行うことになる。なお、範囲３及び範囲４それぞれは、被検体Ｐに応じて、設定値を変える必要があるが、初期設定値が用いられる場合であってもよい。また、範囲３及び範囲４それぞれは、例えば、第１の実施形態の第４の変形例で説明した範囲１と範囲２との自動設定のように、制御部１２６により自動的に、設定値が決定される場合でもよい。

また、前述した範囲３および範囲４は、直流成分を含む１つの範囲でもよいし、直流成分を含み複数の範囲に分割された範囲であってもよい。その場合にも、範囲４は範囲３よりも高周波成分まで含むのが好適である。

上述してきたように、第２の実施形態によれば、位相画像再構成時と強度画像再構成時で、時間変化を表す周波数軸の再構成範囲を別々に設定することで、流速の時間変化を正しく求め、かつ、再構成ノイズを抑制し、強度画像（時系列強度画像）の画質を改善することができる。

（第２の実施形態の第１の変形例）
第２の実施形態では、範囲３（第１の再構成範囲）及び範囲４（第２の再構成範囲）は、固定値として取り扱ったが、再構成の間に動的に変更してもよい。第２の実施形態の第１の変形例では、画像再構成手法として、反復的（iterative）な演算を複数のステップで繰り返し行うことにより最終的な出力画像を算出する方法であり、複数のステップそれぞれでの演算の結果に応じて、前述の第１の再構成範囲（範囲３）及び第２の再構成範囲（範囲４）を、複数のステップそれぞれで動的に変更する手法を用いる。

例えば、ステップＳ３０２およびステップＳ３０３で、事前知識を利用した正則化（Regularization）や圧縮センシング（Compressed Sensing）といった手法を用いて繰り返し処理によりｘ−ｔ空間データ（第２のｘ−ｆ空間データ及び第３のｘ−ｆ空間データ）の再構成を行う。その際、範囲３および範囲４は繰り返し処理ごとに変更する。

上述した第２の実施形態の第１の変形例によれば、再構成範囲を動的に変更することで、流速の時間変化を正しく求め、かつ、再構成ノイズを抑制し、時系列強度画像の画質を改善するとともに、強度画像及び位相画像を適切に再構成することができる。

（その他の実施形態）
なお、実施形態は、上述した実施形態に限られるものではない。

（画像処理装置）
上述した実施形態においては、医用画像診断装置であるＭＲＩ装置１００が各種処理を実行する場合を説明したが、実施形態はこれに限られるものではない。例えば、ＭＲＩ装置１００に替わり、画像処理装置や、ＭＲＩ装置１００と画像処理装置とを含む画像処理システムが、上述した各種処理を実行してもよい。ここで、画像処理装置とは、例えば、ワークステーション、ＰＡＣＳ（Picture Archiving and Communication System）の画像保管装置（画像サーバ）やビューワ、電子カルテシステムの各種装置等である。この場合、例えば、画像処理装置は、ＭＲＩ装置１００によって収集されたｋ空間データを、ＭＲＩ装置１００から、若しくは、画像サーバからネットワーク経由で受信することで、あるいは、記録媒体を介して操作者から入力されること等で、受け付けて、記憶部に記憶する。そして、画像処理装置は、記憶部に記憶したこのｋ空間データを対象として、上述した各種処理（例えば、変換部１２２ａ、強度画像生成部１２２ｂ、位相画像生成部１２２ｃによる処理）を実行すれば良い。

（プログラム）
上述した実施形態の中で示した処理手順に示された指示は、ソフトウェアであるプログラムに基づいて実行されることが可能である。汎用の計算機システムが、このプログラムを予め記憶しておき、このプログラムを読み込むことにより、上述した実施形態のＭＲＩ装置や画像処理装置による効果と同様な効果を得ることも可能である。上述した実施形態で記述された指示は、コンピュータに実行させることのできるプログラムとして、磁気ディスク（フレキシブルディスク、ハードディスク等）、光ディスク（ＣＤ−ＲＯＭ、ＣＤ−Ｒ、ＣＤ−ＲＷ、ＤＶＤ−ＲＯＭ、ＤＶＤ±Ｒ、ＤＶＤ±ＲＷ等）、半導体メモリ、又はこれに類する記録媒体に記録される。コンピュータ又は組み込みシステムが読み取り可能な記憶媒体であれば、その記憶形式は何れの形態であってもよい。コンピュータは、この記録媒体からプログラムを読み込み、このプログラムに基づいてプログラムに記述されている指示をＣＰＵで実行させれば、上述した実施形態のＭＲＩ装置や画像処理装置と同様な動作を実現することができる。もちろん、コンピュータがプログラムを取得する場合又は読み込む場合はネットワークを通じて取得又は読み込んでもよい。

また、記憶媒体からコンピュータや組み込みシステムにインストールされたプログラムの指示に基づきコンピュータ上で稼働しているＯＳ（オペレーティングシステム）や、データベース管理ソフト、ネットワーク等のＭＷ（ミドルウェア）等が、上述した実施形態を実現するための各処理の一部を実行してもよい。

更に、記憶媒体は、コンピュータあるいは組み込みシステムと独立した媒体に限らず、ＬＡＮ（Local Area Network）やインターネット等により伝達されたプログラムをダウンロードして記憶又は一時記憶した記憶媒体も含まれる。

また、記憶媒体は１つに限られず、複数の媒体から、上述した実施形態における処理が実行される場合も、実施形態における記憶媒体に含まれ、媒体の構成は何れの構成であってもよい。

なお、実施形態におけるコンピュータ又は組み込みシステムは、記憶媒体に記憶されたプログラムに基づき、上述した実施形態における各処理を実行するためのものであって、パソコン、マイコン等の１つからなる装置、複数の装置がネットワーク接続されたシステム等の何れの構成であってもよい。

また、実施形態におけるコンピュータとは、パソコンに限らず、情報処理機器に含まれる演算処理装置、マイコン等も含み、プログラムによって実施形態における機能を実現することが可能な機器、装置を総称している。

また、当業者には自明であるが、これまですべての実施形態において、フーリエ変換及び逆フーリエ変換は、互いに逆フーリエ変換及びフーリエ変換で置き換えても良い。また、フーリエ変換の代わりに、フーリエ級数展開をしてもよい。

以上述べた少なくとも一つの実施形態によれば、正しい流速データが得られ、かつ、強度画像の画質を向上することができる。

本発明のいくつかの実施形態を説明したが、これらの実施形態は、例として提示したものであり、発明の範囲を限定することは意図していない。これら実施形態は、その他の様々な形態で実施されることが可能であり、発明の要旨を逸脱しない範囲で、種々の省略、置き換え、変更を行うことができる。これら実施形態やその変形は、発明の範囲や要旨に含まれると同様に、特許請求の範囲に記載された発明とその均等の範囲に含まれるものである。

１００ ＭＲＩ装置
１２２ 画像生成部
１２２ａ 変換部
１２２ｂ 強度画像生成部
１２２ｃ 位相画像生成部

Claims (8)

1. 磁気共鳴イメージング装置によって、時系列に沿って収集された時系列ｋ空間データを、空間上の位置を示すｘ空間上の第１の時系列ｘ空間データに変換する変換部と、
   前記第１の時系列ｘ空間データに第１のフィルタを適用して取得した第２の時系列ｘ空間データから、強度画像を生成する強度画像生成部と、
   前記第１の時系列ｘ空間データに、前記第１のフィルタとは異なる第２のフィルタを適用して取得した第３の時系列ｘ空間データから、位相画像を生成する位相画像生成部と
   を備え、
   前記第１のフィルタは、前記第２のフィルタより、データを平滑化する度合いの大きいフィルタである、画像処理装置。
2. 前記第１のフィルタ及び前記第２のフィルタは、時間軸方向の軸ｔと異なる軸ｆの方向に関するフィルタであり、
   前記強度画像生成部は、前記第１の時系列ｘ空間データに対して、前記軸ｔを前記軸ｆに変換するｔ−ｆ変換を適用してｘ−ｆ空間データを取得し、前記ｘ−ｆ空間データに対して、前記第１のフィルタを適用した後に、前記ｔ−ｆ変換の逆変換を行うことで、前記第２の時系列ｘ空間データを取得し、
   前記位相画像生成部は、前記第１の時系列ｘ空間データに対して、前記ｔ−ｆ変換を適用して前記ｘ−ｆ空間データを取得し、前記ｘ−ｆ空間データに対して、前記第２のフィルタを適用した後に、前記ｔ−ｆ変換の逆変換を行うことで、前記第３の時系列ｘ空間データを取得する請求項１に記載の画像処理装置。
3. 前記第１のフィルタ及び前記第２のフィルタは、時間軸方向に平滑化を行うフィルタであり、
   前記強度画像生成部は、前記第１の時系列ｘ空間データに対して、前記時間軸方向に前記第１のフィルタを適用して、前記第２の時系列ｘ空間データを取得し、
   前記位相画像生成部は、前記第１の時系列ｘ空間データに対して、前記時間軸方向に前記第２のフィルタを適用して、前記第３の時系列ｘ空間データを取得する請求項１に記載の画像処理装置。
4. 撮像条件に関する情報及び撮像される被検体に関する情報のうち少なくとも一方の情報を記憶部から取得し、前記少なくとも一方の情報に基づいて、前記第１のフィルタの特性を定める第１のパラメータ及び前記第２のフィルタの特性を定める第２のパラメータを決定し、前記第１のパラメータを前記強度画像生成部に通知し、前記第２のパラメータを前記位相画像生成部に通知する制御部を更に備えた請求項１〜３のいずれか１項に記載の画像処理装置。
5. 時系列に沿って収集されたｋ空間データである時系列ｋ空間データに対して、空間上の位置を示すｘの空間であるｘ空間上の時系列データである時系列ｘ空間データに変換する変換処理、及び、時間軸方向の軸ｔにおけるデータを、前記軸ｔと異なる軸ｆにおけるデータに変換するｔ−ｆ変換を行う変換処理を適用して、第１のｘ−ｆ空間データを生成する変換部と、
   前記軸ｆの方向において第１の再構成範囲にある前記第１のｘ−ｆ空間データの各データ点を再構成の対象範囲として、所定の画像再構成手法により、前記第１のｘ−ｆ空間データがアンフォールディングされた第２のｘ−ｆ空間データを再構成し、前記第２のｘ−ｆ空間データに前記ｔ−ｆ変換の逆変換を適用して得たｘ−ｔ空間データから、強度画像を生成する強度画像再構成部と、
   前記軸ｆの方向において前記第１の再構成範囲より広い再構成範囲である第２の再構成範囲にある前記第１のｘ−ｆ空間データの各データ点を再構成の対象範囲として、前記所定の画像再構成手法により、前記第１のｘ−ｆ空間データがアンフォールディングされた第３のｘ−ｆ空間データを再構成し、前記第３のｘ−ｆ空間データに前記ｔ−ｆ変換の逆変換を適用して得たｘ−ｔ空間データから、位相画像を生成する位相画像再構成部と、
   を備えた画像処理装置。
6. 前記所定の画像再構成手法は、反復的な演算を複数のステップで繰り返し行うことにより最終的な出力画像を算出する方法であり、前記複数のステップそれぞれでの前記演算の結果に応じて、前記第１の再構成範囲及び前記第２の再構成範囲を、前記複数のステップそれぞれで変更する手法である、請求項５に記載の画像処理装置。
7. 時系列に沿ったｋ空間データである時系列ｋ空間データを収集する収集部と、
   前記時系列ｋ空間データを用いて強度画像及び位相画像を生成する請求項１〜６のいずれか１項に記載の画像処理装置と、
   を備える磁気共鳴イメージング装置。
8. 前記収集部は、位相コントラスト法に基づいて所定の双極傾斜磁場を印加するパルスシーケンスを用いて、前記時系列ｋ空間データを収集し、
   前記画像処理装置は、前記位相画像を基に、流速を算出する、
   請求項７に記載の磁気共鳴イメージング装置。