(19)【発行国】日本国特許庁(JP)
(12)【公報種別】特許公報(B2)
(11)【特許番号】
(24)【登録日】2024-10-11
(45)【発行日】2024-10-22
(54)【発明の名称】撮影装置
(51)【国際特許分類】
A61B 5/055 20060101AFI20241015BHJP
G01N 24/00 20060101ALI20241015BHJP
G01R 33/56 20060101ALI20241015BHJP
【ＦＩ】
A61B5/055 382
A61B5/055 390
G01N24/00 510Y
G01R33/56
【請求項の数】
2
(21)【出願番号】P 2021113419
(22)【出願日】2021-07-08
(65)【公開番号】P2023009819
(43)【公開日】2023-01-20
【審査請求日】2024-04-05
(73)【特許権者】
【識別番号】305027401
【氏名又は名称】東京都公立大学法人
(74)【代理人】
【識別番号】100137752
【弁理士】
【氏名又は名称】亀井 岳行
(72)【発明者】
【氏名】沼野 智一
(72)【発明者】
【氏名】伊東 大輝
【審査官】宮川 数正
(56)【参考文献】
【文献】特開２０２１－９７８２２（ＪＰ，Ａ）
【文献】特開２０１９－１２２５５４（ＪＰ，Ａ）
【文献】特開２００９－１８０７９（ＪＰ，Ａ）
【文献】中国特許出願公開第１１２２４４８１３（ＣＮ，Ａ）
(58)【調査した分野】(Int.Cl.，ＤＢ名)
Ａ６１Ｂ ５／０５５
Ｇ０１Ｎ ２４／００－２４／１４
Ｇ０１Ｒ ３３／００－３３／６４
(57)【特許請求の範囲】
【請求項1】
静磁場と、位置に応じて磁場が変化する傾斜磁場と、プロトンの磁気共鳴条件に基づいて予め設定された交番磁場と、を被検者の被検査部に対して発生させる磁場発生装置と、被検査部に振動を付与する振動付与部材と、
前記傾斜磁場を印加する時期に応じて取得された電磁波に基づいて、電磁波信号の位相に応じたＭＲ位相画像を取得するＭＲ位相画像の取得手段であって、予め定められた時間間隔をあけて複数の前記ＭＲ位相画像を取得する前記ＭＲ位相画像の取得手段と、
ＭＲ位相画像の全画像領域の中から区分けされた判別対象領域において、伝播する振動の波長が長くなると、前記判別対象領域のサイズを大きくする領域サイズ調整手段と、
前記領域サイズ調整手段で調整された判別対象領域に基づいて、伝播する振動の主方向が揃っているか否かを判別する判別手段と、
前記判別手段で振動の主方向が揃っていない判別対象領域の波が除去されたＭＲ位相画像に基づいて、ＭＲエラストグラフィ画像を作成するＭＲＥ画像作成手段と、
を備えたことを特徴とする撮影装置。
【請求項2】
前記振動付与部材で付与される振動の伝播波長が大きくなるにつれて、前記判別対象領域の大きさを大きくする前記領域サイズ調整手段、を備えたことを特徴とする請求項１に記載の撮影装置。
【発明の詳細な説明】
【技術分野】
【0001】
本発明は、撮影装置に関し、特に、組織の硬さを測定するＭＲエラストグラフィ（Magnetic Resonance Elastography : ＭＲＥ）や超音波エラストグラフィの撮影が可能な撮影装置に関する。
【背景技術】
【0002】
医療現場において、患者に対して放射線被曝等の影響がないＭＲＩ装置（Magnetic Resonance Imaging装置：磁気共鳴画像装置）が使用されている。ＭＲＩ装置では、人体の各細胞に含まれる水素原子核（プロトン）に対して、プロトンのスピンに応じた高周波磁場を印加して励起し、励起したプロトンが元の状態に戻る（緩和）際に発する電磁波に基づいて、プロトンの密度が異なる部分（例えば、水と脂肪）を濃淡で表した画像を得ることが可能である。
ＭＲＥ（Magnetic Resonance Elastography：磁気共鳴エラストグラフィ）は、対象物に振動（体幹部の場合には、５０Ｈｚ～１００Ｈｚ程度）を加えながら、ＭＲＩ装置で撮像することで、対象部内部の「硬さ」の違いによる振動波の伝播の違いを利用し、硬さを画像化する撮像法である（特許文献１，２参照）。
ＭＲＥ（Magnetic Resonance Elastography：磁気共鳴エラストグラフィ）は、対象物に振動（体幹部の場合には、５０Ｈｚ～１００Ｈｚ程度）を加えながら、ＭＲＩ装置で撮像することで、対象部内部の「硬さ」の違いによる振動波の伝播の違いを利用し、硬さを画像化する撮像法である（特許文献１，２参照）。
【0003】
ＭＲエラストグラフィにおいて、体内を振動波が伝播する場合、振動源から付与される振動の進行方向に対して、硬さの異なる臓器や骨、筋肉、体表等で伝播波が反射したり散乱したりして進行方向が変わることがある。したがって、計測される伝播波は、さまざまな方向の成分を含むこととなる。主方向の伝播波と、反射波、散乱波が干渉すると、伝播波画像（Wave Image：ＷＩ）の局所波長が変化する、すなわち、誤差が発生する可能性がある。この変化（誤差）が弾性率マップ（elastogram）の精度を低下させる恐れがある。
【0004】
このような反射波、散乱波に対して、これまでにも、Directional Filter（ＤＦ）と呼ばれる画像フィルタを使用して、反射等で進行方向が変わった伝播波を除去することが行われている（非特許文献１参照）。非特許文献１に記載の技術では、主方向の伝播波（非特許文献１では上部から下部に向かう伝播波）に対してＤＦを適用することで、反射波（非特許文献１では下部から上部に向かう波、もしくは波紋状に上方向に広がる波）を除去している。
【先行技術文献】
【特許文献】
【0005】
【文献】特表２００５－５０７６９１号公報（「０００８」、「００１４」～「００２３」）
【文献】特開２０１１－９８１５８号公報（「０００３」～「００１８」、「００２６」～「００２９」）
【非特許文献】
【0006】
【文献】A.Manduca et al.,”Spatio-temporal directional filtering forimproved inversion of MR elastography images”, Medical Image Analysis 7 (2003)465-473
【発明の概要】
【発明が解決しようとする課題】
【0007】
（従来技術の問題点）
非特許文献１に記載のＤＦ（Directional Filter）では、特定の一方向の波以外を除去する技術であるため、波が特定の一方向に伝播することが前提となっている。ファントム（性能評価用の試料）のような理想的な状態では、非特許文献１に記載の技術は適用可能であるが、生体内は、臓器や骨等のため、伝播する波の主方向が変化することがほとんどである。逆に言えば、生体内で伝播波が一方向だけに進むことは、ほぼあり得ない。
したがって、非特許文献１に記載の技術では、生体内を伝播する波に対して、途中で主方向が変化した波まで除去してしまう恐れが高い問題がある。主方向の波が除去されると、ＭＲＥ画像の精度が低下したり、ＭＲＥ画像が形成できない問題がある。
非特許文献１に記載のＤＦ（Directional Filter）では、特定の一方向の波以外を除去する技術であるため、波が特定の一方向に伝播することが前提となっている。ファントム（性能評価用の試料）のような理想的な状態では、非特許文献１に記載の技術は適用可能であるが、生体内は、臓器や骨等のため、伝播する波の主方向が変化することがほとんどである。逆に言えば、生体内で伝播波が一方向だけに進むことは、ほぼあり得ない。
したがって、非特許文献１に記載の技術では、生体内を伝播する波に対して、途中で主方向が変化した波まで除去してしまう恐れが高い問題がある。主方向の波が除去されると、ＭＲＥ画像の精度が低下したり、ＭＲＥ画像が形成できない問題がある。
【0008】
なお、ＭＲＥではなくＭＲＩ画像やＣＴ画像において、撮像パラメータや画像処理パラメータなどに国際的な基準がある（ＱＩＢＡ：Quantitative Imaging Biomarkers Alliance）。ＭＲＥにおいても同様な枠組みが形成されつつあり、「伝播波に反射や散乱が含まれていない一方向に伝播する領域（コヒーレントな伝播波領域）を選択して硬さを測定する」というものがテーマとして挙げられている。現状では、コヒーレントな伝播波領域の判断は、観察者の主観によって決定されており、客観的な領域設定がなされていない。
【0009】
本発明は、利用者の主観で領域設定を行う場合に比べて、領域設定の精度および再現性の低下を解決することを技術的課題とする。
【課題を解決するための手段】
【0010】
前記技術的課題を解決するために、請求項１に記載の発明の撮影装置は、
静磁場と、位置に応じて磁場が変化する傾斜磁場と、プロトンの磁気共鳴条件に基づいて予め設定された交番磁場と、を被検者の被検査部に対して発生させる磁場発生装置と、
被検査部に振動を付与する振動付与部材と、
前記傾斜磁場を印加する時期に応じて取得された電磁波に基づいて、電磁波信号の位相に応じたＭＲ位相画像を取得するＭＲ位相画像の取得手段であって、予め定められた時間間隔をあけて複数の前記ＭＲ位相画像を取得する前記ＭＲ位相画像の取得手段と、
ＭＲ位相画像の全画像領域の中から区分けされた判別対象領域において、伝播する振動の波長が長くなると、前記判別対象領域のサイズを大きくする領域サイズ調整手段と、
前記領域サイズ調整手段で調整された判別対象領域に基づいて、伝播する振動の主方向が揃っているか否かを判別する判別手段と、
前記判別手段で振動の主方向が揃っていない判別対象領域の波が除去されたＭＲ位相画像に基づいて、ＭＲエラストグラフィ画像を作成するＭＲＥ画像作成手段と、
を備えたことを特徴とする。
静磁場と、位置に応じて磁場が変化する傾斜磁場と、プロトンの磁気共鳴条件に基づいて予め設定された交番磁場と、を被検者の被検査部に対して発生させる磁場発生装置と、
被検査部に振動を付与する振動付与部材と、
前記傾斜磁場を印加する時期に応じて取得された電磁波に基づいて、電磁波信号の位相に応じたＭＲ位相画像を取得するＭＲ位相画像の取得手段であって、予め定められた時間間隔をあけて複数の前記ＭＲ位相画像を取得する前記ＭＲ位相画像の取得手段と、
ＭＲ位相画像の全画像領域の中から区分けされた判別対象領域において、伝播する振動の波長が長くなると、前記判別対象領域のサイズを大きくする領域サイズ調整手段と、
前記領域サイズ調整手段で調整された判別対象領域に基づいて、伝播する振動の主方向が揃っているか否かを判別する判別手段と、
前記判別手段で振動の主方向が揃っていない判別対象領域の波が除去されたＭＲ位相画像に基づいて、ＭＲエラストグラフィ画像を作成するＭＲＥ画像作成手段と、
を備えたことを特徴とする。
【0011】
請求項２に記載の発明は、請求項１に記載の撮影装置において、
前記振動付与部材で付与される振動の伝播波長が大きくなるにつれて、前記判別対象領域の大きさを大きくする前記領域サイズ調整手段、
を備えたことを特徴とする。
前記振動付与部材で付与される振動の伝播波長が大きくなるにつれて、前記判別対象領域の大きさを大きくする前記領域サイズ調整手段、
を備えたことを特徴とする。
【発明の効果】
【0012】
請求項１に記載の発明によれば、利用者の主観でコヒーレント領域を判別する場合に比べて精度よく判別可能であり、弾性率計測の精度および再現性の向上を可能にすることができる。
請求項２に記載の発明によれば、振動付与部材で付与される振動の伝播波長に応じて、判別対象領域の大きさを自動的に調整できる。
請求項２に記載の発明によれば、振動付与部材で付与される振動の伝播波長に応じて、判別対象領域の大きさを自動的に調整できる。
【図面の簡単な説明】
【0013】
【図6】図６は実施例１のコヒーレント領域を判別するために領域を区分けする説明図であり、図６Ａは小領域（カーネル）に区分けする説明図、図６Ｂは画素（ボクセル）とカーネルの説明図、図６Ｃは図６Ａの全域において各画素の局所波長を画素値で示した一例の図である。
【図7】図７はWave-direction MapとCoherence Mapであり、図７ＡはWave-direction Map、図７Ｂは図６Ｃに対応する局所波長Map、図７ＣはCoherence Mapである。
【図8】図８は比較例の説明図であり、固定したカーネルサイズを用いて抽出したコヒーレント領域をＭＲ位相画像上に重ね合わせて表示した図である。図８Ａは解像度が低い場合のＭＲ位相画像に重ね合わせた図、図８Ｂは図８Ａの場合において解像度が高い場合にＭＲ位相画像に重ね合わせた図、図８Ｃは図８Ｂと同じ解像度の場合において付与される振動の周波数が低い（波長が長い）場合のＭＲ位相画像に重ね合わせた図である。
【図9】図９は比較例の説明図であり、固定したカーネルサイズを用いて抽出したコヒーレント領域を弾性率マップ上に重ね合わせて表示した図である。図９Ａは図８Ａと同条件の場合で弾性率マップ上に重ね合わせた図、図９Ｂは図８Ｂと同条件の場合で弾性率マップ上に重ね合わせた図、図９Ｃは図８Ｃと同条件の場合で弾性率マップ上に重ね合わせた図である。
【発明を実施するための形態】
【0014】
次に図面を参照しながら、本発明の実施の形態の具体例（以下、実施例と記載する）を説明するが、本発明は以下の実施例に限定されるものではない。
なお、以下の図面を使用した説明において、理解の容易のために説明に必要な部材以外の図示は適宜省略されている。
なお、以下の図面を使用した説明において、理解の容易のために説明に必要な部材以外の図示は適宜省略されている。
【実施例1】
【0015】
図１は本発明の実施例１の磁気共鳴撮影装置の説明図である。
図１において、本発明の撮影装置の一例としての実施例１の磁気共鳴撮影装置１は、磁場発生装置の一例としての磁石部２を有する。磁石部２には、内部を水平方向に貫通する貫通孔３が形成されている。貫通孔３には、寝た状態の被検者４が支持される寝台６が貫通可能である。
磁石部２は、静磁場印加部材の一例としての静磁場発生磁石１１を有する。なお、静磁場発生磁石として、超電導電磁石や永久磁石を使用することが可能である。静磁場発生磁石１１の内側には、傾斜磁場印加部材の一例としての傾斜磁場発生コイル１２が配置されている。傾斜磁場発生コイル１２の内側には、励起磁場印加部材の一例としての高周波磁場発生コイル１３が配置されている。高周波磁場発生コイル１３の内側には、受信部の一例として、電磁波を受信する受信コイル１４が配置されている。
図１において、本発明の撮影装置の一例としての実施例１の磁気共鳴撮影装置１は、磁場発生装置の一例としての磁石部２を有する。磁石部２には、内部を水平方向に貫通する貫通孔３が形成されている。貫通孔３には、寝た状態の被検者４が支持される寝台６が貫通可能である。
磁石部２は、静磁場印加部材の一例としての静磁場発生磁石１１を有する。なお、静磁場発生磁石として、超電導電磁石や永久磁石を使用することが可能である。静磁場発生磁石１１の内側には、傾斜磁場印加部材の一例としての傾斜磁場発生コイル１２が配置されている。傾斜磁場発生コイル１２の内側には、励起磁場印加部材の一例としての高周波磁場発生コイル１３が配置されている。高周波磁場発生コイル１３の内側には、受信部の一例として、電磁波を受信する受信コイル１４が配置されている。
【0016】
実施例１では、被検者４には、被検査部の一例としての大腰筋のＭＲＥ測定をするために、背骨側の体表面の部分に、振動付与部材の一例としての振動パッド１６が支持されている。
【0017】
前記磁石部２には、情報処理装置の一例としてのコンピュータ装置２１がケーブルＣｂを介して電気的に接続されている。したがって、コンピュータ装置２１は、磁石部２との間で、静磁場発生磁石１１等の制御信号や受信コイル１４での検知信号等が送受信可能に構成されている。コンピュータ装置２１は、コンピュータ本体２２と、表示部の一例としてのディスプレイ２３と、入力部の一例としてのキーボード２４およびマウス２５と、を有する。なお、実施例１では、コンピュータ装置２１と磁石部２とをケーブルＣｂで接続する構成を例示したが、これに限定されず、携帯電話回線やBluetooth（登録商標）、無線ＬＡＮ等、任意の無線通信方式で情報の送受信を行うことも可能である。
【0018】
（実施例１のコンピュータ本体２２の制御部の説明）
図２は実施例１の磁気共鳴撮影装置におけるコンピュータ本体の機能ブロック図である。
図２において、実施例１のコンピュータ本体２２の制御部４１は、外部との信号の入出力および入出力信号レベルの調節等を行うＩ／Ｏ（入出力インターフェース）、必要な起動処理を行うためのプログラムおよびデータ等が記憶されたＲＯＭ（リードオンリーメモリ）、必要なデータ及びプログラムを一時的に記憶するためのＲＡＭ（ランダムアクセスメモリ）、ＲＯＭ等に記憶された起動プログラムに応じた処理を行うＣＰＵ（中央演算処理装置）ならびにクロック発振器等を有するコンピュータ装置により構成されており、前記ＲＯＭ及びＲＡＭ等に記憶されたプログラムを実行することにより種々の機能を実現することができる。
制御部４１には、基本動作を制御する基本ソフト、いわゆる、オペレーティングシステムＯＳ、アプリケーションプログラムの一例としての撮影装置制御プログラムＡＰ１、その他の図示しないソフトウェアが記憶されている。
図２は実施例１の磁気共鳴撮影装置におけるコンピュータ本体の機能ブロック図である。
図２において、実施例１のコンピュータ本体２２の制御部４１は、外部との信号の入出力および入出力信号レベルの調節等を行うＩ／Ｏ（入出力インターフェース）、必要な起動処理を行うためのプログラムおよびデータ等が記憶されたＲＯＭ（リードオンリーメモリ）、必要なデータ及びプログラムを一時的に記憶するためのＲＡＭ（ランダムアクセスメモリ）、ＲＯＭ等に記憶された起動プログラムに応じた処理を行うＣＰＵ（中央演算処理装置）ならびにクロック発振器等を有するコンピュータ装置により構成されており、前記ＲＯＭ及びＲＡＭ等に記憶されたプログラムを実行することにより種々の機能を実現することができる。
制御部４１には、基本動作を制御する基本ソフト、いわゆる、オペレーティングシステムＯＳ、アプリケーションプログラムの一例としての撮影装置制御プログラムＡＰ１、その他の図示しないソフトウェアが記憶されている。
【0019】
（実施例１の制御部４１に接続された要素）
制御部４１には、キーボード２４やマウス２５、受信コイル１４等の信号出力要素からの出力信号が入力されている。
また、実施例１の制御部４１は、ディスプレイ２３、静磁場発生磁石１１、傾斜磁場発生コイル１２、高周波磁場発生コイル１３等の被制御要素へ制御信号を出力している。
制御部４１には、キーボード２４やマウス２５、受信コイル１４等の信号出力要素からの出力信号が入力されている。
また、実施例１の制御部４１は、ディスプレイ２３、静磁場発生磁石１１、傾斜磁場発生コイル１２、高周波磁場発生コイル１３等の被制御要素へ制御信号を出力している。
【0020】
（制御部４１の機能）
実施例１の制御部４１の撮影装置制御プログラムＡＰ１は、下記の機能手段（プログラムモジュール）５１～５５を有する。
実施例１の制御部４１の撮影装置制御プログラムＡＰ１は、下記の機能手段（プログラムモジュール）５１～５５を有する。
【0021】
磁場制御手段５１は、磁石部２を制御して、被検者４の被検査部をＭＲ撮影するための磁場を制御する。実施例１の磁場制御手段５１は、繰り返し時間記憶手段５１ａと、エコー時間記憶手段５１ｂと、静磁場印加手段５１ｃと、傾斜磁場印加手段５１ｄと、交番磁場の印加手段の一例としての高周波磁場印加手段５１ｅと、を有する。
【0022】
図３は実施例１の磁場の印加および振動の付与の説明図であり、横軸に時間を取ったグラフである。
繰り返し時間記憶手段５１ａは、被検者４の被検査部に含まれるプロトンを励起するために印加される交番磁場の一例としての高周波磁場を印加する間隔である繰り返し時間ＴＲを記憶する。
繰り返し時間記憶手段５１ａは、被検者４の被検査部に含まれるプロトンを励起するために印加される交番磁場の一例としての高周波磁場を印加する間隔である繰り返し時間ＴＲを記憶する。
【0023】
エコー時間記憶手段５１ｂは、高周波磁場が印加されてから、励起されたプロトンが元の状態に戻る（緩和する）際に発する電磁波を取得するまでの間隔であるエコー時間ＴＥを記憶する。なお、実施例１では、繰り返し時間ＴＲおよびエコー時間ＴＥは、予め設定されているが、磁気共鳴撮影装置１の利用者が手動で入力して、設定、変更が可能に構成することも可能である。
【0024】
静磁場印加手段５１ｃは、静磁場発生磁石１１を制御して、静磁場を発生させる。実施例１の静磁場印加手段５１ｃは、一例として、３［Ｔ］の静磁場を発生させる。
傾斜磁場印加手段５１ｄは、傾斜磁場発生コイル１２を制御して、位置に応じて磁場が変化するＭＲＥ用の傾斜磁場（勾配磁場）を発生させる。従って、傾斜磁場が振動検出傾斜磁場ＭＥＧ（motion encoding gradient）と呼ばれる磁場である。実施例１の傾斜磁場印加手段５１ｄは、図３に示すように、互いに直交するスライス(slice)方向、リードアウト(read out)方向およびフェーズ(phase)方向の３軸方向において、スライス方向（スライス軸）に傾斜磁場を発生させる。
高周波磁場印加手段５１ｅは、高周波磁場発生コイル１３を制御して、プロトンを励起する周波数に対応する交番磁場である高周波磁場を発生させる。
傾斜磁場印加手段５１ｄは、傾斜磁場発生コイル１２を制御して、位置に応じて磁場が変化するＭＲＥ用の傾斜磁場（勾配磁場）を発生させる。従って、傾斜磁場が振動検出傾斜磁場ＭＥＧ（motion encoding gradient）と呼ばれる磁場である。実施例１の傾斜磁場印加手段５１ｄは、図３に示すように、互いに直交するスライス(slice)方向、リードアウト(read out)方向およびフェーズ(phase)方向の３軸方向において、スライス方向（スライス軸）に傾斜磁場を発生させる。
高周波磁場印加手段５１ｅは、高周波磁場発生コイル１３を制御して、プロトンを励起する周波数に対応する交番磁場である高周波磁場を発生させる。
【0025】
振動付与制御手段５２は、振動パッド１６により被検査部に振動を付与する。
受信手段の一例としての信号取得手段５３は、エコー時間ＴＥの時期に、受信コイル１４を介して被検者４のプロトンが緩和する際に発生する電磁波信号を取得する。したがって、実施例１では、図３に示す振動が振動付与制御手段５２で付与された状態で、エコー時間ＴＥにおいて受信コイル１４で信号を測定することで、振動（Vibration）を付与しながらのＭＲ画像を撮像する。
受信手段の一例としての信号取得手段５３は、エコー時間ＴＥの時期に、受信コイル１４を介して被検者４のプロトンが緩和する際に発生する電磁波信号を取得する。したがって、実施例１では、図３に示す振動が振動付与制御手段５２で付与された状態で、エコー時間ＴＥにおいて受信コイル１４で信号を測定することで、振動（Vibration）を付与しながらのＭＲ画像を撮像する。
【0026】
ＭＲ画像取得手段（ＭＲＥ画像作成手段）５４は、信号取得手段５３が取得した電磁波信号に基づいて、ＭＲ画像を作成する。
信号処理手段５４ａは、受信した電磁波信号において信号処理をする。実施例１の信号処理手段５４ａは、信号取得手段５３が取得した信号を実数部ｒ（real part）とし、信号取得手段５３が取得した信号をπ／２位相を遅らせた信号を虚数部ｉ（imaginary part）とする。すなわち、受信した電磁波信号に基づいた複素数、いわゆる、ＭＲＩの技術分野におけるｋ空間（周波数空間）の信号を生成する。そして、実数部ｒと虚数部ｉに対して、フーリエ逆変換（実施例では高速フーリエ逆変換）を行って、実空間の信号Ｒ，Ｉに変換する。そして、実空間における実数部Ｒと虚数部Ｉとに基づいて、複素平面における強度Ｍ＝（Ｒ2＋Ｉ2）1/2と、位相φ＝ｔａｎ-1（Ｉ／Ｒ）とを演算する。なお、この演算は、従来のＭＲＩ装置において導入されており、公知であるため、これ以上の詳細な説明は省略する。
信号処理手段５４ａは、受信した電磁波信号において信号処理をする。実施例１の信号処理手段５４ａは、信号取得手段５３が取得した信号を実数部ｒ（real part）とし、信号取得手段５３が取得した信号をπ／２位相を遅らせた信号を虚数部ｉ（imaginary part）とする。すなわち、受信した電磁波信号に基づいた複素数、いわゆる、ＭＲＩの技術分野におけるｋ空間（周波数空間）の信号を生成する。そして、実数部ｒと虚数部ｉに対して、フーリエ逆変換（実施例では高速フーリエ逆変換）を行って、実空間の信号Ｒ，Ｉに変換する。そして、実空間における実数部Ｒと虚数部Ｉとに基づいて、複素平面における強度Ｍ＝（Ｒ2＋Ｉ2）1/2と、位相φ＝ｔａｎ-1（Ｉ／Ｒ）とを演算する。なお、この演算は、従来のＭＲＩ装置において導入されており、公知であるため、これ以上の詳細な説明は省略する。
【0027】
ＭＲ強度画像の作成手段（ＭＲ強度画像の取得手段）５４ｂは、信号処理手段５４ａで算出された強度Ｍに基づいて、ＭＲ強度画像を作成する。なお、ＭＲ強度画像は、一般的にＭＲＩ画像として、診断に使用される画像である。
【0028】
ＭＲ位相画像の作成手段（ＭＲ位相画像の取得手段）５４ｃは、電磁波信号から算出された位相φに基づいて、ＭＲ位相画像（ＭＲＥではWave Imageとして利用される）を作成する。ここで、実施例１のＭＲ位相画像の作成手段５４ｃでは、エコー時間ＴＥで測定された電磁波信号から算出された位相φに基づいて、ＭＲ位相画像を作成する。
【0029】
判別方向設定手段５４ｄは、各ボクセル（最小画素）における主方向を判別するために、主方向の候補としての判別方向を設定する。実施例１では、判別方向は、一例として、８方向（図４Ｂにおける上から下、下から上、右から左、左から右、右上から左下、左下から右上、左上から右下、右下から左上）が順に設定される。したがって、まず、上から下の方向が設定されて、主方向の判別が行われると、次の方向が設定されて判別が行われるという処理が８つの方向に対して順に行われる。なお、実施例１では、一例として判別方向が８つの場合を例示したが、これに限定されず、１度ずつずらした３６０度分の方向とする等、要求される精度や計算処理速度等に応じて任意に変更可能である。
【0030】
図４は判別方向と変位の読取方向の説明図であり、図４Ａは判別方向と読取方向が同一の場合の説明図、図４Ｂは判別方向と読取方向が逆方向の場合の説明図である。
読取方向設定手段５４ｅは、各ボクセルにおいて設定された判別方向１１１に沿って並んだ複数の位置での振動の変位を取得する際に、変位を読み取る読取方向１１２（図４Ｂ参照）を設定する。実施例１では、図４Ｂに示すように、判別方向１１１とは逆方向に読取方向１１２が設定されている。一例として、判別方向１１１が、図４Ｂにおける「上から下」の場合、判別方向に沿って並んだ３つの点１０３ａ～１０３ｃの各位置での振動の変位を取得する際に、読取方向１１２は「下から上」となり、１０３ｃ，１０３ｂ，１０３ａの順に振動の変位が読み取られるように、読取方向１１２が設定される。
読取方向設定手段５４ｅは、各ボクセルにおいて設定された判別方向１１１に沿って並んだ複数の位置での振動の変位を取得する際に、変位を読み取る読取方向１１２（図４Ｂ参照）を設定する。実施例１では、図４Ｂに示すように、判別方向１１１とは逆方向に読取方向１１２が設定されている。一例として、判別方向１１１が、図４Ｂにおける「上から下」の場合、判別方向に沿って並んだ３つの点１０３ａ～１０３ｃの各位置での振動の変位を取得する際に、読取方向１１２は「下から上」となり、１０３ｃ，１０３ｂ，１０３ａの順に振動の変位が読み取られるように、読取方向１１２が設定される。
【0031】
変位取得手段５４ｆは、ＭＲ位相画像の各ボクセルにおいて、設定された判別方向に沿って３つの点１０３ａ～１０３ｃを抽出して、各点（各位置）における振動の変位を取得する。この時、変位の取得は、読取方向に沿って順に取得する。また、実施例１の変位取得手段５４ｆでは、時間ごとに変化していく各時期のＭＲ位相画像に対して、振動パッド１６の振動の周期Ｔの１／４の時間間隔をあけて、各点１０３ａ～１０３ｃの変位を取得するように設定されている。具体的には、実施例１では、時刻（位相）０，Ｔ／４，Ｔ／２，３Ｔ／４における各点１０３ａ～１０３ｃの変位について、点１０３ｃの時刻０の変位、点１０３ｂの時刻０の変位、点１０３ａの時刻０の変位、点１０３ｃの時刻Ｔ／４の変位、点１０３ｂの時刻Ｔ／４の変位、点１０３ａの時刻Ｔ／４の変位、点１０３ｃの時刻Ｔ／２の変位、点１０３ｂの時刻Ｔ／２の変位、点１０３ａの時刻Ｔ／２の変位、点１０３ｃの時刻３Ｔ／４の変位、点１０３ｂの時刻３Ｔ／４の変位、点１０３ａの時刻３Ｔ／４の変位の順に取得する。なお、実施例１では、各点１０３ａ～１０３ｃを等間隔の点としたが、これに限定されない。すなわち、後述する正弦波への近さの判定では、各点１０３ａ～１０３ｃの間隔が等間隔である方が誤差は少なくなりやすいが、判定自体については等間隔でなくても可能であり、等間隔に制限されない。また、３点に限定されず、４点以上とすることも可能である。
【0032】
波形生成手段５４ｇは、変位取得手段５４ｆで取得された変位に基づいて、判別波形１１３（図５参照）を生成する。
【0033】
図５は残差の演算方法の説明図である。
残差算出手段５４ｈは、判別波形１１３が正弦波に対する近さを導出するために、振動パッド１６で付与される振動である正弦波との差である残差を算出する。なお、実施例１では、振動パッド１６で正弦波状の振動が付与されるが、付与される振動が正弦波状ではない場合は、付与される振動（ｆ（ｔ））との残差を導出することが可能である。
実施例１では、正弦波との残差Ｒを算出するために、波の変位量をｕとし、振幅（最大の変位量）をｕ0とし、波数ベクトルをｋとし、位置ベクトルをｒとし、振動の各振動数をωとし、時間をｔとし、初期位相をθとし、磁場の不均一性に起因する成分をＣとした場合に、以下の式（１）を考える。
ｕ（ｒ，ｔ）＝ｕ0・ｓｉｎ（ｋｒ＋ωｔ＋θ）＋Ｃ …式（１）
残差算出手段５４ｈは、判別波形１１３が正弦波に対する近さを導出するために、振動パッド１６で付与される振動である正弦波との差である残差を算出する。なお、実施例１では、振動パッド１６で正弦波状の振動が付与されるが、付与される振動が正弦波状ではない場合は、付与される振動（ｆ（ｔ））との残差を導出することが可能である。
実施例１では、正弦波との残差Ｒを算出するために、波の変位量をｕとし、振幅（最大の変位量）をｕ0とし、波数ベクトルをｋとし、位置ベクトルをｒとし、振動の各振動数をωとし、時間をｔとし、初期位相をθとし、磁場の不均一性に起因する成分をＣとした場合に、以下の式（１）を考える。
ｕ（ｒ，ｔ）＝ｕ0・ｓｉｎ（ｋｒ＋ωｔ＋θ）＋Ｃ …式（１）
【0034】
この式において、未知の定数は、ｕ0、ｋ、θ、Ｃである。ｋの値を決めると、他の変数は正弦波の線形最小二乗近似で推定可能であり、ｕ（ｒ，ｔ）の近似解ｕ′（ｒ，ｔ）が求まる。ｕ（ｒ，ｔ）とｕ′（ｒ，ｔ）との残差が最小となるｋを推定するために以下の手順（１）～（３）の処理を行う。
（１）適当なｋの初期値を与える。
（２）最小二乗近似により近似解ｕ′（ｒ，ｔ）を算出し、実測値ｕ（ｒ，ｔ）、すなわち、判別波形１１３との残差Ｒを算出する。
（３）最急降下法により残差Ｒが最小になるまでｋを更新する。
（１）～（３）の手順で得られた最小の残差Ｒを、８つの判別方向１１１毎に算出する。
（１）適当なｋの初期値を与える。
（２）最小二乗近似により近似解ｕ′（ｒ，ｔ）を算出し、実測値ｕ（ｒ，ｔ）、すなわち、判別波形１１３との残差Ｒを算出する。
（３）最急降下法により残差Ｒが最小になるまでｋを更新する。
（１）～（３）の手順で得られた最小の残差Ｒを、８つの判別方向１１１毎に算出する。
【0035】
主方向判別手段５４ｉは、８つの判別方向１１１における判別波形１１３の中で最も正弦波に近い判別波形１１３が得られた判別方向１１１を主方向と判別する。実施例１では、残差算出手段５４ｈで得られた残差Ｒを、プロファイル方向（判別方向）を１８０度変えた残差Ｒと比較し、残差Ｒが小さい方向を記録する。これを、全判別方向（３６０度分、実施例１では８方向分）行う。よって、８つの方向で行ったとすると、１８０度変える前と後の２方向の比較が４セット繰り返されることとなる。それぞれで残差Ｒが小さい方向が記録されるので、４方向のベクトルとなる。それぞれのベクトルの大きさを同じ（例えば「１」）に正規化して、ベクトルの和を算出する。この合成ベクトルが各ボクセルの波の進行方向（主方向）となる。
【0036】
図６は実施例１のコヒーレント領域を判別するために領域を区分けする説明図であり、図６Ａは小領域（カーネル）に区分けする説明図、図６Ｂは画素（ボクセル）とカーネルの説明図、図６Ｃは図６Ａの全域において各画素の局所波長を画素値で示した一例の図である。
局所波長推定手段５４ｊは、ＭＲ位相画像の作成手段５４ｃで作成されたＭＲ位相画像に基づいて、ＭＲ位相画像内の各所における伝播波長の分布を導出する。前記主方向判別手段５４ｉで判別された各ボクセルの主方向にＤＦを適用し、反射や散乱の影響を取り除いた状態で、伝播波の局所位相変化から局所波長を推定する。この、ＤＦの適用と局所波長の推定を、全ボクセル分繰り返す。
ただし、局所波長の推定方法は上記に限定されず、推定精度が担保されている方法であれば代用することが可能である。
局所波長推定手段５４ｊは、ＭＲ位相画像の作成手段５４ｃで作成されたＭＲ位相画像に基づいて、ＭＲ位相画像内の各所における伝播波長の分布を導出する。前記主方向判別手段５４ｉで判別された各ボクセルの主方向にＤＦを適用し、反射や散乱の影響を取り除いた状態で、伝播波の局所位相変化から局所波長を推定する。この、ＤＦの適用と局所波長の推定を、全ボクセル分繰り返す。
ただし、局所波長の推定方法は上記に限定されず、推定精度が担保されている方法であれば代用することが可能である。
【0037】
領域区分け手段５４ｋは、ＭＲ位相画像の全域１０１を、予め定められた広さの領域の一例として小領域（カーネル）１０３に区分けする。図６Ａにおいて、小領域（カーネル）１０３のサイズは、ＭＲ位相画像の全域１０１の大きさや分解能等の、仕様に応じて任意に変更可能であるが、少なくとも３画素×３画素以上の広さが必要である。
【0038】
領域サイズ調整手段５４ｍは、各ボクセル１０２において、コヒーレント領域の判別を行うための判別対象領域（カーネル）１０３の大きさを調整する。図６Ｂ、図６Ｃにおいて、実施例１の領域サイズ調整手段５４ｍは、カーネル１０３のサイズを、ボクセル１０２を伝播する振動の波長に応じて調整する。実施例１では、波長が大きくなるほど、カーネル１０３のサイズが大きくなるように調整する。具体的には、振動付与制御手段５２で付与される振動の伝播波長（局所波長）が、予め定められた閾値よりも高い場合（図６Ｃで黒色の領域１０１ｂ、弾性率の高い領域）は、カーネル１０３のサイズを大きくし、前記閾値よりも低い場合（図６Ｃで白色の領域１０１ａ、弾性率の低い領域）はカーネル１０３のサイズを小さくする。
一例として、長波長の場合はカーネル１０３を８画素×８画素とし、低波長の場合はカーネル１０３を４画素×４画素とすることが可能であるが、具体的な数値は例示した数値に限定されず、要求される精度や仕様、利用者の設定等に応じて任意に変更可能である。
また、伝播波長と閾値を使用して判別を行う構成を例示したがこれに限定されない。例えば、伝播波長に限らず、弾性率のように「伝播波長に関連するパラメータ」を使用することも可能である。
一例として、長波長の場合はカーネル１０３を８画素×８画素とし、低波長の場合はカーネル１０３を４画素×４画素とすることが可能であるが、具体的な数値は例示した数値に限定されず、要求される精度や仕様、利用者の設定等に応じて任意に変更可能である。
また、伝播波長と閾値を使用して判別を行う構成を例示したがこれに限定されない。例えば、伝播波長に限らず、弾性率のように「伝播波長に関連するパラメータ」を使用することも可能である。
【0039】
図７はWave-direction MapとCoherence Mapであり、図７ＡはWave-direction Map、図７Ｂは図６Ｃに対応する局所波長Map、図７ＣはCoherence Mapである。
コヒーレント領域判別手段５４ｎは、領域サイズ調整手段５４ｍで調整されたカーネル１０３に基づいて、各カーネル１０３がコヒーレントな領域であるか否かを判別する。実施例１のコヒーレント領域判別手段５４ｎは、まず、ＭＲ位相画像（Wave Image画像）から主方向を算出した画像（Wave-direction map：図７Ａ参照）に対し、ボクセル毎に局所波長を算出して（図７Ｂ参照）、カーネル１０３のサイズを調整し、カーネル１０３内のベクトル平均（主方向のベクトルの大きさを１と仮定）を算出する（Coherence map：図７Ｃ参照）。そして、コヒーレント領域判別手段５４ｎは、予め定められた閾値に基づいて、各カーネル１０３のベクトル平均の値が１に近ければ（閾値以上であれば）コヒーレント、０に近ければ（閾値未満であれば）インコヒーレントと判定する。すなわち、１に近ければ、その領域を伝播する振動はコヒーレント、すなわち、振動の位相（振動の進行方向）が揃っていて、０に近ければ、その領域を伝播する振動はインコヒーレント、すなわち、振動の進行方向が揃っていない、と判別される。
コヒーレント領域判別手段５４ｎは、領域サイズ調整手段５４ｍで調整されたカーネル１０３に基づいて、各カーネル１０３がコヒーレントな領域であるか否かを判別する。実施例１のコヒーレント領域判別手段５４ｎは、まず、ＭＲ位相画像（Wave Image画像）から主方向を算出した画像（Wave-direction map：図７Ａ参照）に対し、ボクセル毎に局所波長を算出して（図７Ｂ参照）、カーネル１０３のサイズを調整し、カーネル１０３内のベクトル平均（主方向のベクトルの大きさを１と仮定）を算出する（Coherence map：図７Ｃ参照）。そして、コヒーレント領域判別手段５４ｎは、予め定められた閾値に基づいて、各カーネル１０３のベクトル平均の値が１に近ければ（閾値以上であれば）コヒーレント、０に近ければ（閾値未満であれば）インコヒーレントと判定する。すなわち、１に近ければ、その領域を伝播する振動はコヒーレント、すなわち、振動の位相（振動の進行方向）が揃っていて、０に近ければ、その領域を伝播する振動はインコヒーレント、すなわち、振動の進行方向が揃っていない、と判別される。
【0040】
除去手段５４ｐは、主方向以外の振動を除去する。実施例１では、コヒーレント領域判別手段５４ｎで判別されたコヒーレントな領域内に関して、ボクセル毎に、主方向判別手段５４ｉで判別された主方向以外の方向についてＤＦ：Directional Filterで振動の成分を除去する。これにより、弾性率算出に誤差を引き起こす散乱波や反射波の除去が可能となる。なお、ＤＦ自体は、公知であるため、詳細な説明は省略する。
したがって、実施例１のＭＲ画像取得手段（ＭＲＥ画像作成手段）５４は、ＭＲ位相画像の作成手段５４ｃで作成されたＭＲ位相画像から、除去手段５４ｐでボクセル毎に主方向以外の振動が除去されて弾性率が導出された画像を、最終的なＭＲＥ画像として出力（作成）する。
したがって、実施例１のＭＲ画像取得手段（ＭＲＥ画像作成手段）５４は、ＭＲ位相画像の作成手段５４ｃで作成されたＭＲ位相画像から、除去手段５４ｐでボクセル毎に主方向以外の振動が除去されて弾性率が導出された画像を、最終的なＭＲＥ画像として出力（作成）する。
【0041】
画像表示手段５５は、ＭＲ画像取得手段５４で作成された画像を、ディスプレイ２３に表示する。すなわち、被検査部の断面画像であるＭＲ強度画像（通常の診断で利用する画像）と、ＭＲＥ画像とがディスプレイ２３に表示される。
【0042】
（実施例１の作用）
前記構成を備えた実施例１の磁気共鳴撮影装置１では、ＭＲ位相画像において、領域１０１をボクセル１０２に区分けして、各ボクセル１０２において伝播する振動の波長に応じてサイズが調整されたカーネル１０３毎に、そのカーネル内がコヒーレントな領域であるかが判別される。生体内では音響インピーダンスが異なる組織が含まれるため反射や散乱が生じやすく、たとえコヒーレントな領域であっても硬さの算出に誤差が生じる場合がある。この問題に対し、判別された主方向にあわせたＤＦを適用することで、反射や散乱を抑制すると共に、主方向以外の振動成分を除去できる。したがって、反射波や散乱波等の主方向以外の波が除去されない場合に比べて、ＭＲＥ画像において各部位の硬さ等の算出を精度よく行うことができ、診断の精度向上にも寄与する。
前記構成を備えた実施例１の磁気共鳴撮影装置１では、ＭＲ位相画像において、領域１０１をボクセル１０２に区分けして、各ボクセル１０２において伝播する振動の波長に応じてサイズが調整されたカーネル１０３毎に、そのカーネル内がコヒーレントな領域であるかが判別される。生体内では音響インピーダンスが異なる組織が含まれるため反射や散乱が生じやすく、たとえコヒーレントな領域であっても硬さの算出に誤差が生じる場合がある。この問題に対し、判別された主方向にあわせたＤＦを適用することで、反射や散乱を抑制すると共に、主方向以外の振動成分を除去できる。したがって、反射波や散乱波等の主方向以外の波が除去されない場合に比べて、ＭＲＥ画像において各部位の硬さ等の算出を精度よく行うことができ、診断の精度向上にも寄与する。
【0043】
また、カーネル１０３のサイズが大きくなると、対象となる波の分布領域が大きく（さまざまな主方向をもつ波がカーネル内に含まれるように）なり、潜在的にコヒーレンス度が低くなる傾向にある。一方、カーネルサイズが大きくなると、対象となる波の分布領域が小さくなるため、コヒーレンス度が高くなる傾向にある。ＭＲＥではこれと同様の問題が弾性率マップの空間分解能に影響する。具体的には、波長が短いと、ある領域内に含まれる波の数（ピーク、波の極大と極小の数）が多くなるため、弾性率マップの空間分解能が向上する。一方、波長が長いと、弾性率マップの空間分解能は低下する。反射や散乱などの波の干渉は弾性率算出にエラーを引き起こすが、空間分解能が低い場合、干渉領域だけでなくその周囲においても弾性率算出にエラーが生じてしまう。
【0044】
したがって、カーネル１０３のサイズが波長に関わらず固定の場合は、付与される振動の波長や、振動が伝播する部位の弾性率によっては、十分な精度でコヒーレント領域を判定できなかったり、エラー領域が含まれたりして、弾性率計測の精度が低下する恐れがある。
これに対して、実施例１では、振動の波長に応じてカーネル１０３のサイズが調整されている。したがって、伝播する振動の波長が長い場合に、１つのカーネル１０３に含まれる波の数を増やすことが可能であり、カーネル１０３のサイズが固定の場合に比べて、コヒーレント領域の抽出精度を向上し、弾性率計測の精度を向上させることができる。
これに対して、実施例１では、振動の波長に応じてカーネル１０３のサイズが調整されている。したがって、伝播する振動の波長が長い場合に、１つのカーネル１０３に含まれる波の数を増やすことが可能であり、カーネル１０３のサイズが固定の場合に比べて、コヒーレント領域の抽出精度を向上し、弾性率計測の精度を向上させることができる。
【0045】
（実験例）
図８は比較例の説明図であり、固定したカーネルサイズを用いて抽出したコヒーレント領域をＭＲ位相画像上に重ね合わせて表示した図である。図８Ａは解像度が低い場合のＭＲ位相画像に重ね合わせた図、図８Ｂは図８Ａの場合において解像度が高い場合にＭＲ位相画像に重ね合わせた図、図８Ｃは図８Ｂと同じ解像度の場合において付与される振動の周波数が低い（波長が長い）場合のＭＲ位相画像に重ね合わせた図である。
図９は比較例の説明図であり、固定したカーネルサイズを用いて抽出したコヒーレント領域を弾性率マップ上に重ね合わせて表示した図である。図９Ａは図８Ａと同条件の場合で弾性率マップ上に重ね合わせた図、図９Ｂは図８Ｂと同条件の場合で弾性率マップ上に重ね合わせた図、図９Ｃは図８Ｃと同条件の場合で弾性率マップ上に重ね合わせた図である。
図１０は実験例の説明図であり、カーネルサイズを可変にした実施例１で抽出したコヒーレント領域をＭＲ位相画像上に重ね合わせて表示している。図１０Ａは解像度が図７Ａと同条件の場合でＭＲ位相画像上に重ね合わせた図、図１０Ｂは図７Ｂと同条件の場合でＭＲ位相画像上に重ね合わせた図、図１０Ｃは図７Ｃと同条件の場合でＭＲ位相画像上に重ね合わせた図である。
図１１は実験例の説明図であり、カーネルサイズを可変にした実施例１で抽出したコヒーレント領域を弾性率マップ上に重ね合わせて表示した図である。図１１Ａは図１０Ａと同条件の場合で弾性率マップ上に重ね合わせた図、図１１Ｂは図１０Ｂと同条件の場合で弾性率マップ上に重ね合わせた図、図１１Ｃは図１０Ｃと同条件の場合で弾性率マップ上に重ね合わせた図である。
図８は比較例の説明図であり、固定したカーネルサイズを用いて抽出したコヒーレント領域をＭＲ位相画像上に重ね合わせて表示した図である。図８Ａは解像度が低い場合のＭＲ位相画像に重ね合わせた図、図８Ｂは図８Ａの場合において解像度が高い場合にＭＲ位相画像に重ね合わせた図、図８Ｃは図８Ｂと同じ解像度の場合において付与される振動の周波数が低い（波長が長い）場合のＭＲ位相画像に重ね合わせた図である。
図９は比較例の説明図であり、固定したカーネルサイズを用いて抽出したコヒーレント領域を弾性率マップ上に重ね合わせて表示した図である。図９Ａは図８Ａと同条件の場合で弾性率マップ上に重ね合わせた図、図９Ｂは図８Ｂと同条件の場合で弾性率マップ上に重ね合わせた図、図９Ｃは図８Ｃと同条件の場合で弾性率マップ上に重ね合わせた図である。
図１０は実験例の説明図であり、カーネルサイズを可変にした実施例１で抽出したコヒーレント領域をＭＲ位相画像上に重ね合わせて表示している。図１０Ａは解像度が図７Ａと同条件の場合でＭＲ位相画像上に重ね合わせた図、図１０Ｂは図７Ｂと同条件の場合でＭＲ位相画像上に重ね合わせた図、図１０Ｃは図７Ｃと同条件の場合でＭＲ位相画像上に重ね合わせた図である。
図１１は実験例の説明図であり、カーネルサイズを可変にした実施例１で抽出したコヒーレント領域を弾性率マップ上に重ね合わせて表示した図である。図１１Ａは図１０Ａと同条件の場合で弾性率マップ上に重ね合わせた図、図１１Ｂは図１０Ｂと同条件の場合で弾性率マップ上に重ね合わせた図、図１１Ｃは図１０Ｃと同条件の場合で弾性率マップ上に重ね合わせた図である。
【0046】
次に、本発明の効果を確認するために実験を行った。
なお、実施例の除去手段５４ｐでＤＦを適用すると、弾性率の誤差が小さくなり、比較例と実験例との差が分かりにくくなる（弾性率の誤差がごまかされてしまう）。よって、図９、図１１では、比較例と実験例との差をわかりやすくするために、除去手段５４ｐでのＤＦを適用する前の弾性率マップを示している。
図８Ａでは、付与された振動の周波数が１００Ｈｚ、マトリクスサイズが１２８×１２８画素のＭＲ位相画像に対し、１３×１３画素に固定したカーネルサイズを用いてコヒーレント領域を抽出した。図９Ａでは、図８ＡのＭＲ位相画像から作成した弾性率マップ上に抽出したコヒーレント領域を重ね合わせて表示している。図１０Ａでは、図８ＡのＭＲ位相画像に対し、カーネル１０３のサイズが可変の実施例１の構成でコヒーレント領域の抽出を行った。図１１Ａでは、弾性率マップ上に可変のカーネルサイズで抽出したコヒーレント領域を重ね合わせて表示している。図９Ａ（図８Ａ）、図１１Ａ（図１０Ａ）からわかるように、図９Ａでは、波長に対してカーネル１０３のサイズが大きすぎるために一部（特に画像の左上部分）のコヒーレント領域が抽出できていなかった部分についても、図１１Ａに示す実施例１では十分な精度で抽出できることがわかる。
なお、実施例の除去手段５４ｐでＤＦを適用すると、弾性率の誤差が小さくなり、比較例と実験例との差が分かりにくくなる（弾性率の誤差がごまかされてしまう）。よって、図９、図１１では、比較例と実験例との差をわかりやすくするために、除去手段５４ｐでのＤＦを適用する前の弾性率マップを示している。
図８Ａでは、付与された振動の周波数が１００Ｈｚ、マトリクスサイズが１２８×１２８画素のＭＲ位相画像に対し、１３×１３画素に固定したカーネルサイズを用いてコヒーレント領域を抽出した。図９Ａでは、図８ＡのＭＲ位相画像から作成した弾性率マップ上に抽出したコヒーレント領域を重ね合わせて表示している。図１０Ａでは、図８ＡのＭＲ位相画像に対し、カーネル１０３のサイズが可変の実施例１の構成でコヒーレント領域の抽出を行った。図１１Ａでは、弾性率マップ上に可変のカーネルサイズで抽出したコヒーレント領域を重ね合わせて表示している。図９Ａ（図８Ａ）、図１１Ａ（図１０Ａ）からわかるように、図９Ａでは、波長に対してカーネル１０３のサイズが大きすぎるために一部（特に画像の左上部分）のコヒーレント領域が抽出できていなかった部分についても、図１１Ａに示す実施例１では十分な精度で抽出できることがわかる。
【0047】
図８Ｂでは、付与された振動の周波数が１００Ｈｚ、マトリクスサイズが２５６×２５６画素で実験を行った。なお、図８Ｂおよび図９Ｂで使用されたカーネル１０３のサイズは、１３×１３画素で固定されている。図９Ｂ（図８Ｂ）、図１１Ｂ（図１０Ｂ）において、図８Ｂの条件では、波長に対してカーネル１０３のサイズが適切であったため、実験結果に差異はほとんどなかった。
図８Ｃでは、付与された振動の周波数が５０Ｈｚ、マトリクスサイズが２５６×２５６画素で実験を行った。なお、図８Ｃおよび図９Ｃで使用されたカーネル１０３のサイズは、１３×１３画素で固定されている。図９Ｃ（図８Ｃ）、図１１Ｃ（図１０Ｃ）において、図８Ｃの条件では、波長に対してカーネル１０３のサイズが小さすぎたため、エラー領域（特に、図９Ｃの右上の部分）がコヒーレント領域として認識されてしまったが、図１１Ｃではエラー領域がコヒーレント領域内に含まれなくなっている。
図８Ｃでは、付与された振動の周波数が５０Ｈｚ、マトリクスサイズが２５６×２５６画素で実験を行った。なお、図８Ｃおよび図９Ｃで使用されたカーネル１０３のサイズは、１３×１３画素で固定されている。図９Ｃ（図８Ｃ）、図１１Ｃ（図１０Ｃ）において、図８Ｃの条件では、波長に対してカーネル１０３のサイズが小さすぎたため、エラー領域（特に、図９Ｃの右上の部分）がコヒーレント領域として認識されてしまったが、図１１Ｃではエラー領域がコヒーレント領域内に含まれなくなっている。
【0048】
したがって、図８～図１１の実験結果から、図８Ｂ（図９Ｂ）のようにたまたまカーネル１０３のサイズが適切であれば、精度が十分なコヒーレント領域の抽出を行うことも可能であるが、カーネル１０３の適切なサイズは、利用者の経験や技量に依存するところがあり、安定して高精度の弾性率計測を行うことが困難である。これに対して、実施例１では、伝播する振動の波長に応じて自動的にカーネル１０３のサイズが適切なサイズに調整されており、利用者の経験に頼ることなく、安定して高精度な弾性率計測を行うことが可能である。
【0049】
また、実施例１では、振動波が伝播中に伝播の主方向が変化していっても、ボクセル毎に主方向が判別されて、主方向以外の振動波が除去される。よって、特定の一方向以外の振動波が除去される従来の構成に比べて、振動の主方向の波を誤って除去したり減衰させたりすることが低減され、ＭＲＥ画像の精度を向上させることができる。
また、振動波の主方向を作業者（測定者）が手動で選択する従来構成では、作業者の技量や経験によって主方向以外の波を除去する精度が大きく異なってくる。これに対して、実施例１では、自動的に各ボクセルの主方向が判定されており、作業者の技量等に依存せず、安定して振動の主方向以外の波を除去することができる。
また、振動波の主方向を作業者（測定者）が手動で選択する従来構成では、作業者の技量や経験によって主方向以外の波を除去する精度が大きく異なってくる。これに対して、実施例１では、自動的に各ボクセルの主方向が判定されており、作業者の技量等に依存せず、安定して振動の主方向以外の波を除去することができる。
【0050】
（変更例）
以上、本発明の実施例を詳述したが、本発明は、前記実施例に限定されるものではなく、特許請求の範囲に記載された本発明の要旨の範囲内で、種々の変更を行うことが可能である。本発明の変更例（Ｈ01）～（Ｈ07）を下記に例示する。
（Ｈ01）前記実施例において、磁石部２がリング状、いわゆる、トンネル型の磁気共鳴撮影装置を例示したが、これに限定されない。例えば、磁石部２がコの字型、いわゆる、オープン型の磁気共鳴撮影装置にも適用可能である。
（Ｈ02）前記実施例において、例示した具体的な数値は、設計や仕様等に応じて、任意に変更可能である。
以上、本発明の実施例を詳述したが、本発明は、前記実施例に限定されるものではなく、特許請求の範囲に記載された本発明の要旨の範囲内で、種々の変更を行うことが可能である。本発明の変更例（Ｈ01）～（Ｈ07）を下記に例示する。
（Ｈ01）前記実施例において、磁石部２がリング状、いわゆる、トンネル型の磁気共鳴撮影装置を例示したが、これに限定されない。例えば、磁石部２がコの字型、いわゆる、オープン型の磁気共鳴撮影装置にも適用可能である。
（Ｈ02）前記実施例において、例示した具体的な数値は、設計や仕様等に応じて、任意に変更可能である。
【0051】
（Ｈ03）前記実施例において、振動パッド１６の具体的な形状等は、設計や仕様等に応じて変更可能である。
（Ｈ04）前記実施例において、波長の判別は、閾値を超える／超えない、の２段階で判別を行う場合を例示したがこれに限定されない。２つ以上の閾値を使用して、３段階以上の判別を行うことも可能である。これに伴って、カーネル１０３のサイズも３段階以上とすることも可能である。
（Ｈ05）前記実施例において、カーネル１０３の形状が四角形状の場合を例示したがこれに限定されない。円形状や三角形状、六角形状等、任意の形状とすることが可能である。
（Ｈ06）前記実施例において、磁場を利用するＭＲエラストグラフィ技術に適用する場合を例示したが、これに限定されない。例えば、超音波を利用して弾性率を導出する超音波エラストグラフィにも適用可能である。
（Ｈ04）前記実施例において、波長の判別は、閾値を超える／超えない、の２段階で判別を行う場合を例示したがこれに限定されない。２つ以上の閾値を使用して、３段階以上の判別を行うことも可能である。これに伴って、カーネル１０３のサイズも３段階以上とすることも可能である。
（Ｈ05）前記実施例において、カーネル１０３の形状が四角形状の場合を例示したがこれに限定されない。円形状や三角形状、六角形状等、任意の形状とすることが可能である。
（Ｈ06）前記実施例において、磁場を利用するＭＲエラストグラフィ技術に適用する場合を例示したが、これに限定されない。例えば、超音波を利用して弾性率を導出する超音波エラストグラフィにも適用可能である。
【0052】
（Ｈ07）前記実施例において、各手段５４ｄ～５４ｉで説明したように３つの点１０３ａ～１０３ｃの各位置での振動の変位を使用して主方向の算出を行う場合を例示したがこれに限定されず、主方向が導出可能な任意の方法を採用可能である。例えば、主方向の導出方法の他の例として、以下のProfile Gradient法を使用することが可能である。
まず、各ボクセルにおいて離散フーリエ変換を行って、位相画像を導出し、位相画像で判別方向（ＰＡ：Profile Angle、プロファイル角度φ）に沿った位相の勾配（ＰＧ: Profile Gradient、振動初期位相の空間勾配、ＰＧ（φ））を取得する。この処理を、各ＰＡに対して行う。
そして、ＰＧ（φ）に対して、ＰＧの最大値（振幅）をＰＧ0とし、プロファイル角度に対するＰＧ（φ）の初期位相をθとした場合、ＰＧ（φ）＝ＰＧ0・ｓｉｎ（φ＋θ）が成立する。
この式において、未知の定数はＰＧ0とθであり、実測値ＰＧ（φ）と推定値ＰＧ′（φ）との残差が最小となるＰＧ0、θを最小二乗法で算出できる。
求まったＰＧ0、θにおける推定値ＰＧ′（φ）が最小となるプロファイル角度φが、振動の伝播の主方向となる。
まず、各ボクセルにおいて離散フーリエ変換を行って、位相画像を導出し、位相画像で判別方向（ＰＡ：Profile Angle、プロファイル角度φ）に沿った位相の勾配（ＰＧ: Profile Gradient、振動初期位相の空間勾配、ＰＧ（φ））を取得する。この処理を、各ＰＡに対して行う。
そして、ＰＧ（φ）に対して、ＰＧの最大値（振幅）をＰＧ0とし、プロファイル角度に対するＰＧ（φ）の初期位相をθとした場合、ＰＧ（φ）＝ＰＧ0・ｓｉｎ（φ＋θ）が成立する。
この式において、未知の定数はＰＧ0とθであり、実測値ＰＧ（φ）と推定値ＰＧ′（φ）との残差が最小となるＰＧ0、θを最小二乗法で算出できる。
求まったＰＧ0、θにおける推定値ＰＧ′（φ）が最小となるプロファイル角度φが、振動の伝播の主方向となる。
【符号の説明】
【0053】
１…撮影装置、
２…磁場発生装置、
４…被検者、
１６…振動付与部材、
５４…ＭＲＥ画像作成手段、
５４ｃ…ＭＲ位相画像の取得手段、
５４ｆ…変位取得手段、
５４ｉ…主方向判別手段、
５４ｍ…領域サイズ調整手段、
５４ｐ…除去手段、
１０１…ＭＲ位相画像の全画像領域、
１０３…判別対象領域、
１１１…判別方向、
１１３…判別波形、
ｕ（ｒ、ｔ）…変位。
２…磁場発生装置、
４…被検者、
１６…振動付与部材、
５４…ＭＲＥ画像作成手段、
５４ｃ…ＭＲ位相画像の取得手段、
５４ｆ…変位取得手段、
５４ｉ…主方向判別手段、
５４ｍ…領域サイズ調整手段、
５４ｐ…除去手段、
１０１…ＭＲ位相画像の全画像領域、
１０３…判別対象領域、
１１１…判別方向、
１１３…判別波形、
ｕ（ｒ、ｔ）…変位。
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
【図10】
【図11】
