JP2018029656A - 磁気共鳴イメージング装置および画像再構成方法 - Google Patents

Abstract

【課題】ＣＳを適用したＭＲＩ装置において、検査終了後画像を表示するまでの時間を短縮する。
【解決手段】ＭＲＩ装置の画像処理部２２０は、並列して処理を行う第一演算部２２１と第二演算部２２２とを備える。画像処理部２２０は、疎な計測データを用いて画像再構成する際に、疎な計測データを復元する処理と、それに続く画像再構成や濃淡補正等の処理とを実施する。その際、第一演算部にて、復元する処理を行い、それと並行して、第二演算部で、復元処理以外の画像再構成として行われる処理の一部を実施し、第一演算部による復元処理後に残りの画像再構成に係る処理を実施する。
【選択図】図５

Description

本発明は、圧縮センシング技術を適用した磁気共鳴イメージング装置（以下、ＭＲＩ装置という）に関し、特に、その画像再構成手法に関する。

ＭＲＩ装置は、主に水素原子核の核磁気共鳴現象を利用した医用画像診断装置であるが、その特性上、検査時間が長いという問題がある。一般的に一検体当たり数十分を要し、患者への負担も大きいことから、撮像の高速化が望まれている。そこで近年、圧縮センシング（Ｃｏｍｐｒｅｓｓｅｄ Ｓｅｎｓｉｎｇ：ＣＳ）と呼ばれる手法のＭＲＩ装置への応用が進められている（例えば、特許文献１、２）。ＣＳは信号の持つ疎性を利用し、疎な観測結果から本来の信号を高精度に復元することが可能である。

撮像の高速化を図る技術として、パラレルイメージングがある。パラレルイメージングは、感度分布が異なる複数の受信コイルを用いて、通常画像再構成に必要なデータよりも少ないデータを計測し（アンダーサンプリング）、それにより画像に含まれることになる折り返しを受信コイルの感度分布を用いた演算によって除去する技術である。このパラレルイメージングとＣＳを組み合わせることも提案されている（例えば、特許文献３）。

米国特許第７６４６９２４号明細書 特開２０１５−２０５０３７号公報 米国特許公開第２０１０／０３０８８２４号明細書

通常、操作者は検査終了後、なるべく早く画像を確認する必要がある。また次の検査に早く移ることで、より多くの検査を実施することが可能となる。このような理由により、撮像から画像が表示されるまでの時間がなるべく短いことが望ましい。しかしＣＳは、信号復元の処理量が膨大であることが多く、その結果、検査完了後に画像が表示されるまでの時間が長いという問題がある。さらに、ＭＲＩ装置では、ＣＳによって本来の信号を高精度に復元した後で、画像化するための再構成処理も実施する必要がある。例えば２ＤＦＦＴ、パラレルイメージング演算や、コイルの感度差を補正するための処理等が挙げられる。これらの再構成処理は、信号復元後のデータに対して行われるため、順番に実施され、画像表示まで時間はＣＳに係る処理と再構成処理との単純な加算となり延長する。

本発明は、ＣＳを適用したＭＲＩ装置において、検査終了後画像を表示するまでの時間を短縮することを課題とする。

本発明のＭＲＩ装置は、ＣＳ演算以外の再構成に関わる演算のうち、ＣＳ演算の完了を待たずに開始可能な演算について、ＣＳ演算と並列して処理する手段を設けることにより上記課題を解決する。

具体的には、本発明のＭＲＩ装置は、被検体から核磁気共鳴信号を計測し、疎の計測データを取得する計測部と、前記計測部が計測した疎の計測データを用いて画像再構成に必要な処理を行う画像処理部と、前記計測部及び前記画像処理部の動作を制御する制御部と、を備え、前記画像処理部は、並行して動作する第一演算部と第二演算部とを有し、前記制御部は、前記疎の計測データを用いて、計測されていないデータの復元処理を前記第一演算部が行い、前記復元処理以外の画像再構成に必要な演算の少なくとも一部を前記第二演算部が行うように前記画像処理部の動作を制御することを特徴とする。

本発明によれば、ＣＳを適用した撮像において検査完了から画像が表示されるまでの時間を短縮することができる。

（Ａ）はＭＲＩ装置１００の全体構成を示すブロック図、（Ｂ）は受信部を示す図。 ＣＳにおける処理を示す図。 ＣＳにおける観測パターンの例を示す図。 ＣＳの復元処理を概念的に示す図。 実施形態の画像処理部２２０の構成例を示すブロック図。 実施形態の画像処理部２２０の処理の流れを示す図。 処理時間の短縮効果を示す図で、（Ａ）は従来手法の処理時間、（Ｂ）は実施形態の処理時間を示す。 従来手法による画像処理部の処理の流れを示す図。 第１の実施形態の画像処理部の処理の流れを示す図。 第１の実施形態の処理時間の短縮効果を示す図で、（Ａ）は従来手法の処理時間、（Ｂ）は第１の実施形態の処理時間を示す。 第２の実施形態の画像処理部の処理の流れを示す図。 第３の実施形態の画像処理部の処理の流れを示す図。 第４の実施形態の画像処理部の処理の流れを示す図。 第５の実施形態の画像処理部の処理の流れを示す図。 第６の実施形態の画像処理部の処理の流れを示す図。 第６の実施形態の処理時間の短縮効果を示す図で、（Ａ）は従来手法の処理時間、（Ｂ）は第６の実施形態の処理時間を示す。 第７の実施形態の画像処理部の処理の流れを示す図。 第７の処理時間の短縮効果を示す図で、（Ａ）は従来手法の処理時間、（Ｂ）は第７の実施形態の処理時間を示す。

以下、本発明を適用した実施形態を、図面を参照して説明する。

まず、本実施形態のＭＲＩ装置の全体構成について説明する。図１（Ａ）は、本実施形態のＭＲＩ装置１００の全体構成を示すブロック図である。本実施形態のＭＲＩ装置１００は、大きく分けて、計測部１１０と、信号処理部１２０と、操作部１３０とを備える。

計測部１１０は、被検体ここでは人１５０が配置される空間に静磁場を発生する静磁場発生部１０１、同空間に傾斜磁場を発生する傾斜磁場コイル１０２、検査対象１５０に高周波磁場パルスを照射する送信コイル１０３、及び、被検体１５０が発生する核磁気共鳴信号を受信する受信コイル１０４を備える。

静磁場発生部１０１は、永久磁石方式、常電導方式あるいは超電導方式の静磁場発生源により構成され、垂直磁場方式であれば、被検体１の周りの空間にその体軸と直交する方向に、水平磁場方式であれば、被検体１の体軸方向に均一な静磁場を発生させる。

傾斜磁場コイル１０２は、ＭＲＩ装置１００の座標系（静止座標系）であるＸ、Ｙ、Ｚの３軸方向に巻かれた傾斜磁場コイルから成り、それぞれ傾斜磁場電源１１２に接続されている。傾斜磁場コイル１０２は傾斜磁場電源１１２により駆動されて、Ｘ、Ｙ、Ｚの３軸方向に傾斜磁場Ｇｘ，Ｇｙ，Ｇｚを発生し、その組み合わせで任意の方向の傾斜磁場を発生することができる。例えば、撮影時には、スライス面（撮影断面）に直交する方向にスライス方向傾斜磁場パルス（Ｇｓ）を印加して被検体１に対するスライス面を設定し、そのスライス面に直交して且つ互いに直交する残りの２つの方向に位相エンコード方向傾斜磁場パルス（Ｇｐ）と周波数エンコード方向傾斜磁場パルス（Ｇｆ）とを印加して、エコー信号にそれぞれの方向の位置情報をエンコードする。

送信コイル１０３は、被検体１５０に近接して配置された高周波コイルからなり、図示しない高周波発振器、変調器、及び高周波増幅器を備えた送信部１１３に接続されている。高周波発振器から出力されたＲＦパルスを変調器により振幅変調した後、高周波増幅器で増幅し、送信コイルから被検体に照射することで、被検体１５０の生体組織を構成する原子の原子核スピンに核磁気共鳴を起こさせる。

受信コイル１０４は、被検体１５０から核磁気共鳴により放出されるエコー信号（ＮＭＲ信号）を検出する高周波コイルで、図示しない信号増幅器、直交位相検波器、及びＡ／Ｄ変換器などを備えた受信部１１４に接続されている。送信コイル１０３から照射されたＲＦパルスによって誘起された応答のＮＭＲ信号は、被検体１５０に近接して配置された受信コイル１０４で検出され、信号増幅器で増幅され、直交位相検波器により直交する二系統の信号に分割され、それぞれがＡ／Ｄ変換器でディジタル量に変換されて、計測データとして信号処理部１２０に送られる。

受信コイル１０４は、パラレルイメージングを行う場合、図１（Ｂ）に示すように、複数の小型ＲＦコイル１０４ａ〜１０４ｄを組み合わせたマルチプルアレイコイル（ＭＡＣ）を用いる。その場合、上述した受信部１１４はチャンネルごと（マルチプルアレイコイルを構成する小型ＲＦコイル毎）に設けられている。なお図では４チャンネルのマルチプルアレイコイルの例を示しているが、小型ＲＦコイル（チャンネル）の数は図示するものに限らず任意である。

信号処理部１２０は、装置全体の制御、各種データ処理、処理結果の表示及び保存等を行うもので、制御部２１０と、画像処理部２２０と、記憶部２３０とを備える。制御部２１０は、主として、計測部１１０を制御する機能を有し、例えば、所定のパルスシーケンスに従って、ＲＦパルスと傾斜磁場パルスとを繰り返し印加するよう制御するとともに、計測データ収集に必要な種々の命令を送信部１１３、傾斜磁場電源１１２、および受信部１１４に送る。本実施形態では、ＣＳを適用した計測データを収集するように撮像部１１０を制御する。ＣＳを適用した計測データについては後述する。

画像処理部２２０は、計測したデータに対し画像再構成に必要な種々の演算を行う。本実施形態では、ＣＳ演算、画像の補正、またパラレルイメージングの場合にはそのための演算などを行う。このような機能を実現するため、信号処理部１２０は、具体的には、ＣＰＵ、ＲＯＭ、ＲＡＭ等の記憶装置や、光ディスク、磁気ディスク等の外部記憶装置を備え、制御部２１０による撮像制御や画像処理部２２０による処理は、予め記憶装置に格納されたプログラムや演算のアルゴリズムを読み込むことにより実現される。また演算の一部は、ＡＳＩＣやＦＰＧＡのようなハードウェアで実現することも可能である。

操作部１３０は、ＭＲＩ装置１００自体の各種制御情報および信号処理部１２０で行う処理の各種制御情報の入力を受け付けるもので、トラックボール又はマウス、キーボード等の入力装置からなる入力部２４０と、画像やＧＵＩ等を表示す表示部２５０とを備える。操作者は、表示部２５０を見ながら入力部２４０を介してインタラクティブにＭＲＩ装置１００の各種処理に必要な情報を入力する。

信号処理部１２０は、受信部１１４からの計測データが入力されると、画像処理部２２０が信号処理、画像再構成等の処理を実行し、その結果である被検体１５０の断層画像を表示部２５０に表示すると共に、記憶装置に記録する。

次に本実施形態のＭＲＩ装置に適用されるＣＳについて説明する。
図２に、ＣＳにおける基本的な処理を示す。図示するように、ＣＳでは、疎なＫ空間データの取得Ｓ１０１、信号或いはデータの復元処理Ｓ１０２、復元処理以外の画像再構成処理Ｓ１０３が実行される。以下、各処理を詳述する。

まず、疎なＫ空間データの取得Ｓ１０１では、計測部１１０が、傾斜磁場の印加パターンを調整し、Ｋ空間データのデータ点より少ないデータ点数で且つ所定のパターン（観測パターン）でＫ空間データを収集する。一般的にＣＳでは、信号成分を偏りなく観測するため観測対象の空間をランダムに観測するのがよいとされる。しかしながら、Ｋ空間においては、原点付近の低周波成分を密に観測しないと、再構成される画像のコントラストや明るさに大きな変動が発生してしまうため、原点付近の低周波成分に関してはすべての要素を観測する。このような観測方法による観測パターンとしては、図３に示すような平行線状３１、放射状３２（例えば、特許文献２参照）、ランダム状３３、更には図示を省略したらせん状等が考えられる。本実施形態ではそのいずれを採用してもよい。

次に復元処理Ｓ１０２は、信号処理部１２０がＣＳアルゴリズムによる演算を行う。ＣＳアルゴリズムは、計測したＫ空間データをフーリエ変換して画像空間にしたのち、スパース化変換して、Ｌ１ノルム最小化等の最適化問題を反復アルゴリズムで解くもので、スパース化変換にはウェーブレット変換、それとＣｕｒｖｅｌｅｔ変換の組み合わせ、ＴＶ（Ｔｏｔａｌ Ｖａｒｉａｔｉｏｎ）、Ｒｉｄｇｅｌｅｔ変換などがある。このような圧縮センシングのアルゴリズムは、公知であり、例えば、Ｌ１−ＳＰＩＲｉＴ、Ｌ１−ＥＳＰＩＲｉＴ、ＳＡＫＥ−Ｌ１ＥＳＰＩＲｉＴ、ＴＶＧ（Ｔｏｔａｌ Ｇｅｎｅｒａｌｉｚｅｄ Ｖａｒｉａｔｉｏｎ：ＴＶの改良版）など、ソフトウェアのパッケージとして入手可能である。ＣＳのアルゴリズムはＭＲＩ装置の記憶装置或いは他の記憶装置にプログラムとして格納されており、画像処理部２２０が、そのプログラムをアップロードして実行する。

図４に、ＣＳのＬ１ノルム最小化を用いた復元処理の一例を概念的に示す。観測データ４１、ここでは計測部１１０により計測された疎なＫ空間データ、をフーリエ変換し実空間データ４２とする（Ｓ２０１）。次いで実空間データ４２をスパース変換演算子で圧縮し、スパース画像４３を得る（Ｓ２０２）。このスパース処理は、ＭＲ画像では理想的には同じ組織の信号強度は等しく組織の境界でのみ信号強度が変化するということを利用し、画像を境界とそれ以外の情報に分離してスパース性の高い（すなわちゼロ成分が多い）画像に変換する処理であり、スパース変換演算子としては、ウェーブレット変換、離散余弦変換等がある。このスパース画像４３に対し、閾値によるマスク処理を実施して、閾値マスク後のスパース画像４４を得る。（Ｌ１ノルム最小化では、反復するにつれて、この閾値を徐々に小さくしていくことで、最小化を計る。）次いでスパース画像４４の逆変換を行い、ＣＳ後の実空間データ４５とし（Ｓ２０４）、これを更に逆フーリエ変換し（Ｓ２０５）、復元後のＫ空間データ４６を得る。

初期値である疎な観測データ４１と復元後のＫ空間データ４６との差分（例えば平均二乗誤差）を評価する評価関数を設定し、評価関数の値が所定の閾値に達していない場合には、直前の演算で得られた復元後のＫ空間データ４６と疎な観測データ４１を合わせ（疎な観測データ４１の中で、観測していない部分を復元後の観測データ４６で埋めて）、それを初期値４１として、上述した処理Ｓ２０１〜Ｓ２０５を繰り返す（Ｓ２０６）。

なお図４では、復元処理にＬ１ノルム最小化を用いたが、復元アルゴリズムとしてはＬ１ノルム最小化以外にも種々のものが知られており、Ｌ１ノルム最小化に限定されず採用することができる。

復元処理Ｓ１０２の繰り返し演算は、受信コイル１０４が図１（Ｂ）に示すようなマルチプルアレイコイルの場合、Ｋ空間データはマルチプルアレイコイルを構成する小型ＲＦコイルごとに収集されるので、各Ｋ空間データに対し行う。

復元処理以外の画像再構成処理Ｓ１０３には、例えば、パラレルイメージングの場合の折り返し除去演算（パラレルイメージング演算）、マルチプルアレイコイルを用いた場合のシェーディング（画像濃淡）補正などがある。パラレルイメージング演算には、大きく分けて、計測空間（Ｋ空間）上で収集していないデータを作成する方法（ＳＭＡＳＨやＧＲＡＰＰＡ等）と、実空間上で折り返しを展開して除去する方法（ＳＥＮＳＥなど）がある。前者は、上述した復元処理Ｓ１０２において、ＣＳと複合化したアルゴリズム（例えば特許文献３）で処理することも可能である。ここでは、実空間上で折り返しを展開する方法を説明する。

まずパラレルイメージングでは、計測部１１０が、マルチプルアレイコイルを用いて、位相エンコードの繰り返し回数を一定の割合で間引いて計測する。位相エンコード数を等間隔に１／Ｎ間引いて計測を行った場合には、２次元フーリエ変換後の画像にはＮ回の折り返しが発生する。ここで、Ｎ以上のＭ個の受信コイルを用いて１／Ｎ間引いた計測を行い、Ｎ回折り返しの発生したＭ個の画像（受信コイル毎の画像）を取得する。そして、信号処理部１２０が、Ｍ個の受信コイルの感度分布を用いてＮ回の折り返し画像を展開し、折り返しを除去した１枚の画像を取得する。

画像を、間引きをしないで計測した場合には、画像Ｓは、式（１）に示すように、被検体の密度分布ｐと、コイルの感度分布Ｃの積で表される。ここで、ｘ，ｙ方向の画像マトリスクをそれぞれＸ、Ｙとした場合、ｉ（ただし、１≦ｉ≦Ｍ）番目の小型受信コイルの画像内の座標（ｘ，ｙ）（ｘ：１≦ｘ≦Ｘ，ｙ：１≦ｙ≦Ｙ）の画素値をＳｉ（ｘ，ｙ）、感度分布をＣｉ（ｘ，ｙ）、同じく被検体の密度分布をｐ（ｘ，ｙ）とする。

位相エンコード方向を１／Ｎに間引きして計測した場合には、画像の位相エンコード方向のマトリクスは、Δｙ＝Ｙ／Ｎとなり、折り返し画像は、

となる。

簡単のために、式（３）によって式（１）を置き換えると、式（４）となる。

この式（４）は、Ｍ個のＲＦ受信コイルについて同様に成り立つので、行列で表わすと、式（５）のようになる。

したがって、被検体の密度分布Ｐ、即ち、折り返しの無い画像は、以下のように感度分布の逆行列を求めることで、得ることができる。

次に、復元処理以外の画像再構成処理Ｓ１０３の別の例として、シェーディング補正について説明する。マルチプルアレイコイルを用いた計測では、各局所コイルの感度分布が局在するため、画像にシェーディングが生じる。シェーディング補正では、まずそのシェーディングを補正するための補正係数を求め（処理１）、次いでこの補正係数を計測データ（画像データ）に乗じて補正後の画像データを得る（処理２）。これによりシェーティングが補正された計測データが得られる。

係数を求める処理１では、最初にスライスごとに、独立してフィッティング曲面Ｆ（ｘ、ｙ）を算出し、ＭＡＣ（合成したコイル全体）の感度分布を近似する。この時、算出するＦ（ｘ、ｙ）は、ａ〜ｆの６つのフィッティング係数を用いて式（８）で表せる。

なお式（８）において、ｘ、ｙはそれぞれ周波数エンコード方向、位相エンコード方向の位置を示す。

次に、スライスごとに式（８）で算出したフィッティング曲面に対し、スライス方向にも式（９）で表す２次関数でフィッティングする。マルチプルアレイコイルの感度は３次元空間に滑らかに分布しているため、スライス方向にもフィッティングを行うことで、スライス方向にも連続性を持たせることができる。

なお式（９）において、ｚはスライス方向の位置を示す。

上述した式（８）、（９）を合わせて、フィッティング曲面を式（１０）で求める。

さらに、このフィッティング曲面の逆数を計算し、補正係数（Ｇ（ｘ，ｙ，ｚ））に変換する。
［数１１］
Ｇ（ｘ，ｙ，ｚ）＝１／Ｆ’ （ｘ，ｙ，ｚ） （１１）

シェーディング補正（処理２）は、上述のように求めた補正係数Ｇを画像データに乗じることで、受信コイル感度だけでなく照射不均一に起因する濃淡も含めたシェーディングが補正された画像データが得られる。
［数１２］
Ｓｉ（ｘ，ｙ，ｚ）corrected＝Ｓｉ（ｘ，ｙ，ｚ）×Ｇ（ｘ，ｙ，ｚ）
（１２）

以上、ＣＳを適用したＭＲＩ装置の画像処理部２２０が主として行う処理について説明した。本実施形態のＭＲＩ装置は、これらの処理のうち、特に演算部１２０が行う処理を短縮化するために、処理の内容に応じた並列化を行うことが特徴である。すなわち処理の内容によっては、シーケンシャルであることが必須である（前後関係に拘束がある）が、一つの処理を適切に分割することで分割した処理を並列化することができ、処理時間の短縮化を図ることができる。

このような処理時間の短縮化を図る画像処理部２２０の構成例を、図５に示す。図示するように画像処理部２２０は、複数の演算を並列処理するために第一演算部２２１と第二演算部２２２とを備え、これら第一演算部２２１及び第二演算部２２２の動作を制御する主制御部２２３を更に備えている。画像処理部２２０は、前述したようにＣＰＵを備えた信号処理部１２０の一機能部であり、これら第一演算部２２１及び第二演算部２２２は、そのＣＰＵ上に構築することができるが、並列処理を同一のＣＰＵ上で実施した場合は、それぞれの処理時間が延長してしまい、効果が十分得られない可能性もある。そのため、画像処理部２２０は、複数のＣＰＵを割り当てる、もしくはマルチコアＣＰＵやＧＰＵを用いるなど、同時に複数の処理を実施できる構成であることが望ましい。なお主制御部２２３は、計測部１１０を制御する制御部２１０がその機能を兼ねてもよい。

このような構成における画像処理部２２０の処理の流れを、図６を参照して説明する。
まず計測部１１０が計測した疎なＫ空間データを取り込む（Ｓ１０１）。次いで第一演算部２２１において、疎なＫ空間データの復元処理を開始する（Ｓ１０２）。復元処理は、例えば図４に示すような、スパース変換と閾値処理を含む繰り返し演算であり、複数の受信コイルがある場合には、受信コイル毎に実行される。受信コイル毎の処理は第一演算部において並列処理することができる。第一演算部２２１における処理と並行して、第二演算部２２２において、取得した疎なＫ空間データを用いた処理が開始する（Ｓ２０１）。この処理は、復元後に行われる画像再構成処理（Ｓ２０２）の一部であり、復元前のデータを用いて行う処理、例えば比較的密である低周波領域のデータを利用した処理で、マップの作成、補正係数の算出などである。

第一演算部２２１及び第二演算部２２２における処理が終了後、画像再構成処理Ｓ２０２を実行する。画像再構成処理Ｓ２０２は、２ＤＦＦＴ、パラレルイメージング演算、各受信コイルの画像データの合成（ＭＡＣ合成）、受信コイルの感度マップを用いたシェーディング補正などであり、第一演算部２２１及び第二演算部２２２のいずれが行ってもよいし、或いは第一演算部２２１及び第二演算部２２２とは別に設けた第三の演算部が実行してもよい。この復元処理後に行われる画像再構成処理Ｓ２０２は、処理に必要なマップや係数などが、既に、第二演算部２２２による処理Ｓ２０１で得られているので、それを用いて短時間で行うことができる。

本実施形態による時間短縮効果を図７に示す。図７（Ａ）は、従来の処理時間を、（Ｂ）は本実施形態の処理時間を模式的に示したものである。図示するように、従来では復元処理に続いて画像再構成処理を行っていたのに対し、本実施形態では、時間を要する復元処理の間に、画像再構成処理の一部を行うことで、復元処理後の画像再構成処理の時間を大幅に短縮することができる。

以下、画像再構成処理の具体的な内容に即して、本実施形態をさらに詳述する。

＜第１の実施形態＞
最初に説明する実施形態では、主制御部２２３（或いは制御部２１０）は、第一演算部２２１において、復元処理とパラレルイメージングに基づく計算とを連続して行い、第一演算部２２１の処理に並行して、第二演算部２２２において、パラレルイメージングに基づく計算と補正処理に用いる補正値の算出とを行うよう画像処理部２２０を制御する。第一演算部２２１と第二演算部２２２による処理終了後に、第二演算部が算出した補正値を用いた補正処理を行う。補正処理は、例えば、シェーディング補正であり、第二演算部は、復元前の疎なＫ空間データを使ってシェーディング補正に用いる補正係数の演算を行う。

一例として、復元処理の後に、パラレルイメージングと濃淡の補正を行う場合の処理フローを説明する。比較のために従来の処理フローを図８に示し、本実施形態の処理フローを図９に示す。両図において同じ処理は同一の符号でしめす。また図８及び図９では、処理を示すブロックの横に、その処理を実施した時点でのデータ（頭部撮像した際のイメージ図）を併せて示している（以下の処理フローでも同様）。なお以下の説明で「画像」はＫ空間データに対しフーリエ変換を施し実空間の画像としたものを云うものとする。

まず従来手法の画像化処理フローについて図８を用いて説明する。ステップＳ８０１で、ＣＳにより計測した疎なＫ空間データを取得する。この時点では、このＫ空間データを再構成した画像８１は、輪郭がぼやけた画像になっている。

ステップＳ８０２では、取得した疎なＫ空間データに対して復元処理を実施し、Ｋ空間データを復元する。Ｋ空間データを復元したことで、輪郭がはっきりした画像８２になる。

ステップＳ８０３では、画像（８２）データに対してパラレルイメージング演算を実施して画像に発生している折り返しを展開する。その結果、折り返しがない画像８３になる。

ステップＳ８０４では、パラレルイメージング演算後の画像に対して、コイル感度差や照射不均一などによって発生する濃淡を補正するための値（補正係数）を計算する。この計算は、上述の式（８）〜（１１）で説明した補正係数を算出する処理である。補正係数８４は、信号値が高くなっている部分は信号値を低く補正し、低くなっている部分は信号値を高く補正するような値になる。

ステップＳ８０５では、ステップＳ８０４で作成した補正係数を使って補正処理を実施する（式（１２））。その結果、濃淡がなくなった画像８５が完成する。

このように従来の画像化処理フローでは、各ステップを順番に実施する必要があるので、疎なＫ空間データ取得後、画像が完成するまでの時間は、図１０（Ａ）に示すように、各ステップの処理時間の合計となる。図１０に模式的に示したように各ステップの処理時間は、処理内容によって異なり、復元処理Ｓ８０２及びフィッティング曲面を求める補正係数計算Ｓ８０４は、パラレルイメージング演算Ｓ８０３や補正処理Ｓ８０５に比べ長時間を要する。

続いて、本実施形態の処理フローについて図９を用いて説明する。
ステップＳ８０１で、ＣＳにより計測した疎なＫ空間データを取得した後、長時間を要する２つの処理、復元処理Ｓ８０２と画像濃淡の補正係数計算Ｓ８０４とを並列処理する。例えば第一演算部２２１が復元処理Ｓ８０２と、それに続いてパラレルイメージング演算Ｓ８０３を行い、それと並行して第二演算部２２２がパラレルイメージング演算Ｓ９０１、それに続いて補正係数計算Ｓ９０２とを行う。Ｓ８０２とＳ８０３は、従来の処理と内容が同じであるので説明を省略する。

ステップＳ９０１では、まず、疎なＫ空間データに対してパラレルイメージング処理を実施し、折り返しがない輪郭がぼやけた画像９１を作成する。そして、ステップＳ９０２では、Ｓ９０１で作成された画像データを使って濃淡を補正するための補正係数を作成する。この処理で用いる画像データ（画像９１）は、従来手法のステップＳ８０４で用いたデータ８４と違い、元にするデータの輪郭がぼやけているが、濃淡は特徴的に表れているので、従来と変わらない補正係数を求めることができる。
並列して行われた処理が終了した時点で、ステップＳ９０２で作成した補正係数を使って、ステップＳ８０５の補正処理を実施する。その結果、従来手法と変わらない濃淡がなくなった画像８５が完成する。

本実施形態によれば、比較的長時間の演算を含む「復元処理」と「画像濃淡の補正値計算」とを並列で実施した結果、図１０（Ｂ）に示すように、画像化の時間は従来手法と比べて短くなる。

＜第２の実施形態＞
本実施形態は、画像処理部２２０が、復元処理と並行して、疎なＫ空間データの中心部分を使って、画像の濃淡を補正する値を計算する。以下、図１１を参照して、本実施形態の処理を説明する。図１１において、図９に示す処理と同じ内容の処理は同一符号で示し重複する説明は省略する。

第１の実施形態では、第二演算部２２２がまずパラレルイメージング演算を行った後、折り返し除去後のデータを用いて補正係数計算を行ったが、本実施形態では、疎なＫ空間データのうち中心データを用いて補正係数計算を行う点が異なる。

すなわちステップＳ９１１では疎なＫ空間データの中心部分を抜き出す。中心部分とは、Ｋ空間原点を含む低周波数領域を意味し、図３に示したように、データ点を密に計測した領域のデータであり、ＣＳ計測パターン（３１〜３３）によって決まる。ステップＳ８０１で取得するデータは、画像８１で表したように、疎なＫ空間データであるが、Ｋ空間中心部分は密にデータが存在しているため、折り返しが発生していない。従ってパラレルイメージングを実施せずに、中心部分のＫ空間データを再構成することで折り返しのない画像９３を得ることができる。

そして、ステップＳ９１２では、Ｓ９１１で作成された画像データを使って、フィッティングにより濃淡を補正するための補正係数を算出する（式（８）〜（１１））。
その後、ステップＳ８０５で、ステップＳ９１２で作成した補正値（補正係数）９４を使って補正処理を実施することは第１の実施形態と同様であり、従来と変わらず、濃淡がなくなった画像８５が完成する。

本実施形態でも、図１０に示したように、処理時間の短縮化を図ることができる。

＜第３の実施形態＞
本実施形態は、画像処理部２２０が、復元処理に先立ってパラレルイメージング演算を実施し、その後、復元処理と、シェーディング補正のための補正係数計算とを並列して行うことが特徴である。以下、図１２を参照して、本実施形態の処理を説明する。図１１において、図９に示す処理と同じ内容の処理は同一符号で示し重複する説明は省略する。

疎なＫ空間データを取得後（Ｓ８０１）、最初に、ステップＳ９２１にて、パラレルイメージング演算を行う。この演算は、第一演算部２２１及び第二演算部２２２のいずれが行ってもよいし、或いは第一演算部２２１及び第二演算部２２２とは別に設けた第三の演算部が実行してもよい。パラレルイメージング演算は、Ｋ空間データに対する演算（例えばＳＭＡＳＨ）でもよいし、画像データに対する演算（例えばＳＥＮＳＥ）のいずれでもよい。この処理で、図９の画像９１と同様に、折り返しがない輪郭がぼやけた画像９１が得られる。

次いで、ステップＳ９２２では、第二演算部２２２が、第１の実施形態と同様に、ステップＳ９２１で得られた画像データを用いて、濃淡を補正するための補正係数を作成する。この第二演算部２２２の処理と並行して、ステップＳ９２３で、第一演算部２２１が復元処理を実施する。

復元後の画像データと、ステップＳ９２２で算出した補正係数を用いて補正処理Ｓ８０５を行うことは、第１の実施形態と同様である。また本実施形態においても第１の実施形態と同様の処理時間の短縮効果が得られる。

＜第４の実施形態＞
本実施形態は、第一演算部２２１が、復元処理とパラレルイメージング演算とを複合化アルゴリズムを用いて実施することが特徴である。図１３に、本実施形態の処理フローを示す。図１３において、図９に示す処理と同じ内容の処理は同一符号で示している。

図９に示す処理との比較からわかるように、第一演算部２２１が行う処理（ステップＳ８０２、Ｓ８０３）がステップＳ９３１に置き換わった以外は、第１の実施形態と同様である。

このステップＳ９３１では、例えば、Ｗａｖｅｌｅｔ変換とＧＲＡＰＰＡなどのパラレルイメージング演算とを組み合わせたＣＳアルゴリズムによる復元処理を実施する。このようなＣＳアルゴリズムは、例えば、特許文献３にも開示され、また前掲のソフトウェアパッケージ（Ｌ１−ＳＰＩＲｉＴ、Ｌ１−ＥＳＰＩＲｉＴ）としても入手可能である。

その他のステップは、第１の実施形態と同様であり、従来手法より短い処理時間で、折り返しがなく、濃淡がなくなった画像８５が完成する。

＜第５の実施形態＞
本実施形態は、第一演算部２２１で復元処理を行いながら、それと並行して、第二演算部２２２で、パラレルイメージング用のマップ作成と、シェーディング用の補正係数計算とを実施することが特徴である。即ち、本実施形態では、第２の実施形態における第二演算部２２２の処理として、パラレルイメージング用のマップ作成が追加され、この処理の終了後に、パラレルイメージング演算のステップＳ８０３が実施される。

以下、図１４を参照して、本実施形態の処理を説明する。図１４において、図１１に示す処理と同じ内容の処理は同一符号で示し重複する説明は省略する。

疎なＫ空間データを取得した後（Ｓ８０１）、ステップＳ８０２で第一演算部２２１は、所定のＣＳアルゴリズムを用いて復元処理を行う。一方、第二演算部２２２は、ステップＳ９１１にて、小型コイル毎に、疎なＫ空間データの中心部分を抜き出す。中心部分は第２の実施形態で説明したとおり、密にデータが存在しているため、折り返しが発生しないデータである。

ステップＳ９４１では、この中心部分のデータを再構成した画像９３からパラレルイメージング用のマップを作成する。パラレルイメージング演算は、前述したように、各小型コイルの感度マップ（式（１）或いは（４）のＣｉ（ｘ，ｙ）、ｉ＝１〜Ｍの整数、Ｍはコイルの数）を用いる。マップＣｉ（ｘ，ｙ）は、各小型コイルの画像データＳｉ（ｘ，ｙ）を次のように算出することができる。
［数１３］
Ｃｉ（ｘ，ｙ）＝Ｓｉ（ｘ，ｙ）／√（ΣＣｉ（ｘ，ｙ）^２） （１３）

次のステップＳ９１２では、第二演算部２２２が、第２の実施形態と同様に、ステップＳ９４１で得られた画像データを用いて、濃淡を補正するための補正係数を作成する。

第一演算部２２１における処理（復元処理）と、上述した第二演算部２２２における処理が完了した後、ステップＳ８０３において、ステップＳ９４１で作成したマップを用いてパラレルイメージング演算を行い、次いで、ステップＳ８０５において、ステップＳ９１２で算出した補正係数９４を用いてシェーディング補正を行う。その結果、他の実施形態と同様に短時間で、濃淡がなくなった画像８５が完成する。

＜第６の実施形態＞
第１〜第５の実施形態では、第一演算部２２１の処理と第二演算部２２２の処理は、独立して並列に行われる例であるが、本実施形態は、第一演算部２２１の復元途中のデータを第二演算部２２２に渡し、第二演算部２２２は復元途中のデータを用いて補正係数計算を行うことが特徴である。

以下、図１５を参照して、本実施形態の処理を説明する。図１５において、図９に示す処理と同じ内容の処理は同一符号で示し重複する説明は省略する。

ステップＳ８０１で取得した疎なＫ空間データを、第一演算部２２１が取込み、第１の実施形態と同様に、復元処理Ｓ８０２、続いてパラレルイメージング演算Ｓ８０３を実施する。また第一演算部２２１は、復元処理の途中で、ある程度復元されたデータを第二演算部２２２に渡す。第二演算部２２２に渡すデータの復元度は、第二演算部２２における処理を進める上で問題にならない程度であることが好ましく、特に限定されるものではないが、例えば、閾値処理に用いられた閾値（ソフト閾値）の、当初の閾値に対する割合や繰り返し回数で定めることができる。これら割合或いは繰り返し回数が所定の値に達した時点で、第二演算部２２２は途中データを取り込み（ステップＳ９５１）、ステップＳ９５２でこの途中データに対しパラレルイメージング演算を行い、折り返しが除去された画像９５を得る。その後、ステップＳ９５３において、補正係数計算を行い、補正係数９６を得ることは第１の実施形態と同様である。

一般に、ＣＳ計測において、Ｋ空間を疎にする度合が大きい場合は、上述した実施形態（第１〜第５の実施形態）のように全く復元しないデータをそのまま使うと、復元後の再構成処理を途中まで進めることができなかったり、進めることはできても、画質が劣化してしまったりする可能性がある。本実施形態では、ある程度、復元が進んだデータを使うことで、より正確に、シェーディング補正用の補正係数を計算することができる。

本実施形態における処理時間の短縮効果を、図１６（Ｂ）に示す。図１６（Ａ）は従来法の処理時間である。図示するように、本実施形態では、「復元処理」を途中まで実施した段階から「画像濃淡の補正係数計算」を並列で実施する。従って第１の実施形態と比べると「画像濃淡の補正値計算」を開始するタイミングが遅くなるが、「復元処理」に必要な時間は長いため、「画像濃淡の補正値計算」の方が早く完了し、第１の実施形態と同様の時間短縮効果が得られる。また、仮に「復元処理」が先に完了した場合でも、「画像濃淡の補正値計算」が完了するのを待ってから補正処理を実施することで画像化の時間を短くすることができる。

＜第７の実施形態＞
上述した第６の実施形態は、第一演算部２２１の処理途中でデータを第二演算部２２２に渡し、第二演算部２２２が途中データを用いて処理を行う例であるが、本実施形態では、第一演算部２２１と第二演算部２２２との間でデータのやりとりを行い、第一演算部２２１が第二演算部２２２の処理結果を利用して復元の繰り返し演算を制御することが特徴であり、これにより復元処理にかかる時間を短縮する。

以下、図１７を参照して、本実施形態の処理を説明する。本実施形態では、一例として、第二演算部２２２が疎なＫ空間データを使って作成したパラレルイメージング後のデータを復元処理に使う場合を説明する。なお図１７において、図９に示す処理と同じ内容の処理は同一符号で示し重複する説明は省略する。

疎なＫ空間データ取得（Ｓ８０１）後に第一演算部２２１が行う処理は、第１の実施形態と同様に、ステップＳ８０２、ステップＳ８０３と続く処理で、ステップＳ８０２でＫ空間データが復元され、ステップＳ８０３でパラレルイメージング演算によって折り返しが展開される。従来手法と同じ処理フローである。

第二演算部２２２が行う処理は、まずステップＳ９６１で疎なＫ空間データに対してパラレルイメージング処理を実施し、折り返しがない輪郭がぼやけた画像９１を作成する。次にステップＳ９６２では、ステップＳ８０２の復元処理途中の第一演算部２２１との間でデータをやり取りして、復元処理に利用する。例えば、パラレルイメージング演算（Ｓ９６１）で作成した画像９１と、復元処理Ｓ８０２で途中まで復元した画像を比較し、その差の大きさで復元度合を判断し、ソフト閾値の値を決める係数の値を変更したりする。これにより判断した復元度合いに応じて計算量を減らすことができる。

その後、ステップＳ９６３で、復元前のデータ（ステップＳ９６１でパラレルイメージング演算を行った後のデータ）を使って濃淡を補正するための値を作成し、補正係数を求める。また復元前のデータを使う代わりに、第６の実施形態と同様に、ステップＳ９６２でやり取りした復元途中のデータを用いてもよい。

最後にステップＳ８０５でシェーディング補正を行い、折り返しが除去され且つ濃淡が補正された画像８５を得る。
本実施形態によれば、復元処理をより最適な精度で実施することができ、無駄な処理を省くことで処理時間を短縮することができる。

本実施形態における処理時間の短縮効果を、図１８（Ｂ）に示す。図１８（Ａ）は従来法の処理時間である。図示するように、本実施形態では、「復元処理」と「画像濃淡の補正値計算」を並列で実施するため、画像化の時間は従来手法と比べて短くなる。さらに、疎なＫ空間データを使って作成したデータを復元処理に利用することで、復元処理の処理時間も短くなる。

以上、本実施形態のＭＲＩ装置の、特に画像処理部の処理の各実施形態を説明したが、本発明はこれら実施形態に限定されることなく、種々の変更が可能である。例えば、上記実施形態では、復元処理と並行して実施される処理の主なものとして、シェーディング補正の補正係数の計算を例示したが、画像の再構成に関わる演算であって、疎なＫ空間データ或いは復元途中のデータを用いた演算であれば、適用することが可能である。

また復元処理や補正係数の計算等、実施形態中で説明した手法や式は単なる例示であって、それ以外の手法を用いることも本発明に包含される。さらに必須ではない要素を追加或いは削除した構成も本発明に包含される。

１００・・・ＭＲＩ装置、１１０・・・計測部、１２０・・・信号処理部、１３０・・・操作部、２１０・・・制御部、２２０・・・画像処理部、２２１・・・第一演算部、２２２・・・第二演算部、２２３・・・主制御部、２３０・・・記憶部、２４０・・・表示部、２５０・・・入力部。

Claims (14)

1. 被検体から核磁気共鳴信号を計測し、疎の計測データを取得する計測部と、
   前記計測部が計測した疎の計測データを用いて、復元処理を含む、画像再構成に必要な処理を行う画像処理部と、
   前記計測部及び前記画像処理部の動作を制御する制御部と、を備え、
   前記画像処理部は、並行して動作する第一演算部と第二演算部とを備え、
   前記制御部は、前記復元処理を前記第一演算部が行い、前記復元処理以外の画像再構成に必要な演算の少なくとも一部を前記第二演算部が行うように前記画像処理部の動作を制御することを特徴とする磁気共鳴イメージング装置。
2. 請求項１に記載の磁気共鳴イメージング装置であって、
   前記画像処理部が行う処理は、前記復元処理後の計測データに対する補正処理を含み、
   前記第二演算部は、前記第一演算部の処理と並行して、前記補正処理に用いる補正値の算出を行うことを特徴とする磁気共鳴イメージング装置。
3. 請求項１に記載の磁気共鳴イメージング装置であって、
   前記計測部は、核磁気共鳴信号を計測する受信部として２以上の受信コイルを有し、
   前記画像処理部が行う処理は、前記２以上の受信コイルの感度分布を用いた画像再構成を含み、
   前記第二演算部は、前記第一演算部の処理と並行して、前記感度分布を用いた画像再構成に必要なデータの算出を行うことを特徴とする磁気共鳴イメージング装置。
4. 請求項３に記載の磁気共鳴イメージング装置であって、
   前記２以上の受信コイルの感度分布を用いた画像再構成処理は、パラレルイメージングに基づく計算を含むことを特徴とする磁気共鳴イメージング装置。
5. 請求項１に記載の磁気共鳴イメージング装置であって、
   前記計測部は、核磁気共鳴信号を計測する受信部として２以上の受信コイルを有し、
   前記画像処理部が行う処理は、前記２以上の受信コイルの感度分布を用いたパラレルイメージングに基づく計算と、前記復元処理後の計測データに対する補正処理を含み、第二演算部は、前記第一演算部の処理と並行して、前記補正処理に用いる補正値の算出を行うことを特徴とする磁気共鳴イメージング装置。
6. 請求項５に記載の磁気共鳴イメージング装置であって、
   前記第一演算部は、前記復元処理と前記パラレルイメージングに基づく計算とを連続して行い、前記第二演算部は、前記第一演算部の処理に並行して、パラレルイメージングに基づく計算と前記補正処理に用いる補正値の算出とを行うことを特徴とする磁気共鳴イメージング装置。
7. 請求項６に記載の磁気共鳴イメージング装置であって、
   前記第一演算部は、前記パラレルイメージングに基づく計算と前記復元処理とを複合化したアルゴリズムに基づく処理を実行することを特徴とする磁気共鳴イメージング装置。
8. 請求項５に記載の磁気共鳴イメージング装置であって、
   前記第一演算部は、前記復元処理と前記パラレルイメージングに基づく計算とを連続して行い、前記第二演算部は、前記第一演算部の処理に並行して、前記疎の計測データから低周波帯域データを抽出する処理と、当該抽出した低周波帯域データを用いた、前記補正処理に用いる補正値の算出と、を行うことを特徴とする磁気共鳴イメージング装置。
9. 請求項５に記載の磁気共鳴イメージング装置であって、
   前記画像処理部は、前記第一演算部による復元処理後に、前記パラレルイメージングに基づく計算を行い、
   前記第二演算部は、前記第一演算部による復元処理と並行して、前記パラレルイメージングに基づく計算に用いるマップの作成を行うことを特徴とする磁気共鳴イメージング装置。
10. 請求項５に記載の磁気共鳴イメージング装置であって、
    前記第一演算部は、復元処理の途中データを前記第二演算部に渡し、
    前記第二演算部は、前記復元処理の途中データを用いて、前記補正処理に用いる補正値の算出を行うことを特徴とする磁気共鳴イメージング装置。
11. 請求項５に記載の磁気共鳴イメージング装置であって、
    前記第一演算部は、復元処理の途中データを前記第二演算部に渡し、
    前記第二演算部は、前記パラレルイメージングに基づく演算を行うとともに、その演算結果と前記第一演算部から受け取った前記復元処理の途中データとを用いて、前記第一演算部における前記復元処理に用いる閾値を変更して前記第一演算部に渡すことを特徴とする磁気共鳴イメージング装置。
12. 請求項１に記載の磁気共鳴イメージング装置であって、
    前記第一演算部は、復元処理の途中データを前記第二演算部に渡し、
    前記第二演算部は、前記演算として、前記第一演算部における前記復元処理に用いる閾値の算出を行い、前記第一演算部に渡すことを特徴とする磁気共鳴イメージング装置。
13. 請求項１に記載の磁気共鳴イメージング装置であって、
    前記画像処理部は、ＧＰＵ（Ｇｒａｐｈｉｃｓ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）またはＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）を複数備えることを特徴とする磁気共鳴イメージング装置。
14. 磁気共鳴イメージング装置で計測した、疎な計測データを用いて画像再構成する方法であって、
    疎な計測データを復元する処理と、それに続く画像濃淡補正処理とを含み、
    前記復元する処理と並行して、前記画像濃淡補正処理に用いる補正値の算出を行うことを特徴とする画像再構成方法。
    
    

Images