JP6282942B2 - 超音波診断装置、画像処理装置及び画像処理プログラム - Google Patents

Description

本発明の実施形態は、超音波診断装置、画像処理装置及び画像処理プログラムに関する。

従来、超音波診断装置は、生体の血流の観察や診断を行なうために広く用いられている。超音波診断装置は、ドプラ（Doppler）効果に基づくドプラ法により、超音波の反射波から血流情報の生成及び表示を行なう。超音波診断装置により生成表示される血流情報としては、カラードプラ画像や、ドプラ波形（ドプラスペクトラム）等がある。

カラードプラ画像は、カラーフローマッピング（ＣＦＭ：Color Flow Mapping）法により撮像化される超音波画像である。ＣＦＭ法では、超音波の送受信が複数の走査線上で複数回行なわれる。そして、ＣＦＭ法では、同一位置のデータ列に対してＭＴＩ(Moving Target Indicator)フィルタを掛けることで、静止している組織、或いは、動きの遅い組織に由来する信号（クラッタ信号）を抑制して、血流に由来する信号を抽出する。そして、ＣＦＭ法では、この血流信号から血流の速度、血流の分散、血流のパワー等の血流情報を推定し、推定結果の分布を、例えば、２次元でカラー表示した超音波画像（カラードプラ画像）を表示する。

通常、ＭＴＩフィルタとしては、バタワース型のＩＩＲ（Infinite Impulse Response）フィルタや、多項式回帰フィルタ（Polynomial Regression Filter）等、係数が固定されたフィルタが用いられている。一方、入力信号に応じて、係数を変化させる適応型のＭＴＩフィルタも知られている。

一例として、適応型のＭＴＩフィルタは、ＭＴＩフィルタ入力前の信号から組織の速度を求めて、その位相差がキャンセルされた信号を得る。そして、かかるフィルタは、予め用意したＭＴＩフィルタ用の係数から、得られた信号に応じた係数を選択する。また、「Eigenvector Regression Filter」と呼ばれている適応型のＭＴＩフィルタも知られている。この適応型のＭＴＩフィルタは、相関行列から固有ベクトルを計算し、計算した固有ベクトルから、ＭＴＩフィルタに用いる係数を直接計算する方法により、クラッタ成分が抑制された信号を得る。この方法は、主成分分析や、カルーネン・レーベ変換（Karhunen-Loeve transform）、固有空間法で使われている手法を応用したものである。

このような従来の超音波カラードプラ法では、入力データ列は等間隔である場合が多い。ここで、データ列とは、同一走査線において超音波を送信することで受信した反射波に基づいて生成された反射波データの集合を示す。また、データ列が等間隔である場合とは、同一走査線において送信間隔を一定にして超音波を送信して反射波を受信することで、時間的に等間隔の反射波データを得ることを示す。例えば、ある走査線において４回の超音波の送受信を行う場合、１回目と２回目の超音波の送受信の間隔と、２回目と３回目の超音波の送受信の間隔と、３回目と４回目の超音波の送受信の間隔とが同一である。また、超音波カラードプラ法では、データ列が等間隔ではない場合に適用する方法も開示されている。この方法では、例えば、等間隔ではないデータ列に対して多項式の最小２乗法フィッティングを行うことでクラッタを近似し、原信号からこの近似信号を減算することで血流信号を抽出する。なお、「データ列が等間隔ではない」場合とは、同一走査線において送信間隔を一定にせずに超音波を送信して反射波を受信することで、時間的に間隔の異なる反射波データを得ることを示す。なお、「データ列が等間隔ではない」ことを、適宜「データ列が不等間隔である」と記載する。

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

Bjaerum, Torp, Kristoffersen, "Clutter Filters Adapted to Tissue Motion in Ultrasound Color Flow Imaging", IEEE Transactions on Ultrasonics, Ferroelectrics, Frequency Control, vol.49, pp.693-704, June, 2002

本発明が解決しようとする課題は、不等間隔データ列が入力されるカラードプラ法においてモーションアーティファクトを低減することができる超音波診断装置、画像処理装置及び画像処理プログラムを提供することである。

実施形態の超音波診断装置は、相関行列計算部と、フィルタ行列計算部と、抽出部と、推定部と、画像生成部と、制御部とを備える。相関行列計算部は、同一走査線において送信間隔を一定にせずに超音波を送信することで発生した反射波に基づいて生成された反射波データの集合である第１のデータ列を用いて相関行列を計算する。フィルタ行列計算部は、前記相関行列を用いた主成分分析の結果に基づいて、フィルタ係数を計算する。抽出部は、前記第１のデータ列に含まれる、前記同一走査線上に存在した移動体で反射された前記超音波の反射波に由来する反射波データの集合である第２のデータ列を、前記フィルタ係数を用いて抽出する。推定部は、抽出された前記第２のデータ列に基づいて、前記移動体の移動体情報を推定する。画像生成部は、前記移動体情報に基づいて超音波画像データを生成する。制御部は、前記超音波画像データを表示部に表示させる。

図１は、第１の実施形態に係る超音波診断装置の構成例を示すブロック図である。 図２Ａは、第１の実施形態に係る制御部による不等間隔スキャンの一例を示す図である。 図２Ｂは、第１の実施形態に係る制御部による不等間隔スキャンの一例を示す図である。 図２Ｃは、第１の実施形態に係る制御部による不等間隔スキャンの一例を示す図である。 図２Ｄは、第１の実施形態に係る制御部による不等間隔スキャンの一例を示す図である。 図３は、第１の実施形態に係る制御部による同一ラスタ上での不等間隔スキャンの一例を示す図である。 図４は、第１の実施形態に係るランクカット数の決定処理の一例を説明するための図である。 図５は、第１の実施形態に係るランクカット数の決定処理の一例を説明するための図である。 図６は、第１の実施形態に係る超音波画像データの一例を示す図である。 図７は、第１の実施形態に係る超音波画像データの一例を示す図である。 図８は、第１の実施形態に係る超音波診断装置が行なう処理の一例を説明するためのフローチャートである。 図９は、第２の実施形態に係る相関行列計算部が行なう処理の一例を示す図である。 図１０は、第２の実施形態の効果を説明するための図である。 図１１は、第２の実施形態に係る超音波診断装置が行なう処理の一例を説明するためのフローチャートである。

以下、添付図面を参照して、超音波診断装置の実施形態を詳細に説明する。

（第１の実施形態）
まず、第１の実施形態に係る超音波診断装置の構成について説明する。図１は、第１の実施形態に係る超音波診断装置の構成例を示すブロック図である。図１に例示するように、第１の実施形態に係る超音波診断装置は、超音波プローブ１と、モニタ２と、入力装置３と、装置本体１０とを有する。

超音波プローブ１は、超音波の送受信を行なうために、装置本体１０に接続される。超音波プローブ１は、例えば、複数の圧電振動子を有し、これら複数の圧電振動子は、後述する装置本体１０が有する送受信部１１から供給される駆動信号に基づき超音波を発生する。また、超音波プローブ１が有する複数の圧電振動子は、被検体Ｐからの反射波を受信して電気信号に変換する。また、超音波プローブ１は、圧電振動子に設けられる整合層と、圧電振動子から後方への超音波の伝播を防止するバッキング材等を有する。なお、超音波プローブ１は、装置本体１０と着脱自在に接続される。

超音波プローブ１から被検体Ｐに超音波が送信されると、送信された超音波は、被検体Ｐの体内組織における音響インピーダンスの不連続面で次々と反射され、反射波信号として超音波プローブ１が有する複数の圧電振動子にて受信される。受信される反射波信号の振幅は、超音波が反射される不連続面における音響インピーダンスの差に依存する。なお、送信された超音波パルスが、移動している血流や心臓壁等の表面で反射された場合の反射波信号は、ドプラ効果により、移動体の超音波送信方向に対する速度成分に依存して、周波数偏移を受ける。

なお、第１の実施形態は、超音波プローブ１が、被検体Ｐを２次元で走査する１Ｄアレイプローブであっても、被検体Ｐを３次元で走査するメカニカル４Ｄプローブや２Ｄアレイプローブであっても適用可能である。

入力装置３は、マウス、キーボード、ボタン、パネルスイッチ、タッチコマンドスクリーン、フットスイッチ、トラックボール、ジョイスティック等を有する。入力装置３は、超音波診断装置の操作者からの各種設定要求を受け付け、装置本体１０に対して受け付けた各種設定要求を転送する。

モニタ２は、超音波診断装置の操作者が入力装置３を用いて各種設定要求を入力するためのＧＵＩ（Graphical User Interface）を表示したり、装置本体１０において生成された超音波画像データ等を表示したりする。

装置本体１０は、超音波プローブ１が受信した反射波信号に基づいて超音波画像データを生成する装置である。装置本体１０は、図１に例示するように、送受信部１１と、バッファ１２と、Ｂモード処理部１３と、ドプラ処理部１４と、画像生成部１５と、画像メモリ１６と、内部記憶部１７と、制御部１８とを有する。

送受信部１１は、後述する制御部１８の指示に基づいて、超音波プローブ１が行なう超音波送受信を制御する。送受信部１１は、パルス発生器、送信遅延回路、パルサ等を有し、超音波プローブ１に駆動信号を供給する。パルス発生器は、所定の繰り返し周波数（ＰＲＦ：Pulse Repetition Frequency）で送信超音波を形成するためのレートパルスを繰り返し発生する。また、送信遅延回路は、超音波プローブ１から発生される超音波をビーム状に集束し、かつ送信指向性を決定するために必要な圧電振動子ごとの遅延時間を、パルス発生器が発生する各レートパルスに対し与える。また、パルサは、レートパルスに基づくタイミングで、超音波プローブ１に駆動信号（駆動パルス）を印加する。すなわち、送信遅延回路は、各レートパルスに対し与える遅延時間を変化させることで、圧電振動子面から送信される超音波の送信方向を任意に調整する。

なお、送受信部１１は、後述する制御部１８の指示に基づいて、所定のスキャンシーケンスを実行するために、送信周波数、送信駆動電圧等を瞬時に変更可能な機能を有している。特に、送信駆動電圧の変更は、瞬間にその値を切り替え可能なリニアアンプ型の発信回路、又は、複数の電源ユニットを電気的に切り替える機構によって実現される。

また、送受信部１１は、アンプ回路、Ａ／Ｄ（Analog／Digital）変換器、受信遅延回路、加算器、直交検波回路等を有し、超音波プローブ１が受信した反射波信号に対して各種処理を行って反射波データを生成する。アンプ回路は、反射波信号をチャンネル毎に増幅してゲイン補正処理を行う。Ａ／Ｄ変換器は、ゲイン補正された反射波信号をＡ／Ｄ変換する。受信遅延回路は、デジタルデータに受信指向性を決定するのに必要な受信遅延時間を与える。加算器は、受信遅延回路により受信遅延時間が与えられた反射波信号の加算処理を行う。加算器の加算処理により、反射波信号の受信指向性に応じた方向からの反射成分が強調される。

そして、直交検波回路は、加算器の出力信号をベースバンド帯域の同相信号（Ｉ信号、Ｉ：In-pahse）と直交信号（Ｑ信号、Ｑ：Quadrature-phase）とに変換する。そして、直交検波回路は、Ｉ信号及びＱ信号（以下、ＩＱ信号と記載する）を反射波データとして、バッファ１２に格納する。なお、直交検波回路は、加算器の出力信号を、ＲＦ（Radio Frequency）信号に変換した上で、バッファ１２に格納しても良い。ＩＱ信号や、ＲＦ信号は、位相情報が含まれる信号（受信信号）となる。以下では、送受信部１１が出力する反射波データを、受信信号と記載する場合がある。

送受信部１１は、被検体Ｐ内の２次元領域を走査（２次元走査）する場合、超音波プローブ１から被検体Ｐ内の２次元領域に超音波ビームを送信させる。そして、送受信部１１は、超音波プローブ１が受信した２次元領域の反射波信号から２次元の反射波データを生成する。また、送受信部１１は、被検体Ｐ内の３次元領域を走査（３次元走査）する場合、超音波プローブ１から被検体Ｐ内の３次元領域に超音波ビームを送信させる。そして、送受信部１１は、超音波プローブ１が受信した３次元領域の反射波信号から３次元の反射波データを生成する。

バッファ１２は、送受信部１１が生成した反射波データ（Ｉ／Ｑ信号）を一時的に記憶するバッファである。具体的には、バッファ１２は、数フレーム分のＩ／Ｑ信号、又は、数ボリューム分のＩ／Ｑ信号を記憶する。例えば、バッファ１２は、ＦＩＦＯ（First-In/First-Out）メモリであり、所定フレーム分のＩ／Ｑ信号を記憶する。そして、例えば、バッファ１２は、新たに１フレーム分のＩ／Ｑ信号が送受信部１１にて生成された場合、生成時間が最も古い１フレーム分のＩ／Ｑ信号を破棄して、新たに生成された１フレーム分のＩ／Ｑ信号を記憶する。

Ｂモード処理部１３及びドプラ処理部１４は、送受信部１１が反射波信号から生成した反射波データに対して、各種の信号処理を行なう信号処理部である。Ｂモード処理部１３は、バッファ１２から読み出した反射波データ（Ｉ／Ｑ信号）に対して、対数増幅、包絡線検波処理、対数圧縮などを行なって、複数のサンプル点それぞれの信号強度が輝度の明るさで表現されるデータ（Ｂモードデータ）を生成する。

なお、Ｂモード処理部１３は、フィルタ処理により、検波周波数を変化させることで、映像化する周波数帯域を変えることができる。このＢモード処理部１３の機能を用いることにより、第１の実施形態に係る超音波診断装置は、コントラストハーモニックイメージング（ＣＨＩ：Contrast Harmonic Imaging）や、ティッシュハーモニックイメージング（ＴＨＩ：Tissue Harmonic Imaging）等のハーモニックイメージングを実行可能である。すなわち、Ｂモード処理部１３は、造影剤が注入された被検体Ｐの反射波データから、造影剤（微小気泡、バブル）を反射源とするハーモニック成分の反射波データ（高調波データ又は分周波データ）と、被検体Ｐ内の組織を反射源とする基本波成分の反射波データ（基本波データ）とを分離する。Ｂモード処理部１３は、ハーモニック成分の反射波データ（受信信号）から、造影画像データを生成するためのＢモードデータを生成することができる。

ドプラ処理部１４は、バッファ１２から読み出した反射波データを周波数解析することで、走査範囲内にある移動体のドプラ効果に基づく運動情報を抽出したデータ（ドプラデータ）を生成する。具体的には、ドプラ処理部１４は、移動体の運動情報として、平均速度、平均分散値、平均パワー値等を、複数のサンプル点それぞれでドプラデータを生成する。ここで、移動体とは、例えば、血流や、心壁等の組織、造影剤である。本実施形態に係るドプラ処理部１４は、血流の運動情報（血流情報）として、血流の平均速度、血流の平均分散値、血流の平均パワー値等を、複数のサンプル点それぞれで推定したドプラデータを生成する。

上記のドプラ処理部１４の機能を用いて、本実施形態に係る超音波診断装置は、カラーフローマッピング法（ＣＦＭ：Color Flow Mapping）とも呼ばれるカラードプラ法を実行可能である。ＣＦＭ法では、超音波の送受信が複数の走査線上で複数回行なわれる。そして、ＣＦＭ法では、同一位置のデータ列に対してＭＴＩ(Moving Target Indicator)フィルタを掛けることで、静止している組織、或いは、動きの遅い組織に由来する信号（クラッタ信号）を抑制して、血流に由来する信号を抽出する。そして、ＣＦＭ法では、この血流信号から血流の速度、血流の分散、血流のパワー等の血流情報を推定する。後述する画像生成部１５は、推定結果の分布を、例えば、２次元でカラー表示した超音波画像データ（カラードプラ画像データ）を生成する。そして、モニタ２は、カラードプラ画像データを表示する。

ＭＴＩフィルタとしては、通常、バタワース型のＩＩＲ（Infinite Impulse Response）フィルタや、多項式回帰フィルタ（Polynomial Regression Filter）等、係数が固定されたフィルタが用いられる。一方、本実施形態に係るドプラ処理部１４は、ＭＴＩフィルタとして、入力信号に応じて係数を変化させる適応型のＭＴＩフィルタを用いる。具体的には、本実施形態に係るドプラ処理部１４は、適応型のＭＴＩフィルタとして、「Eigenvector Regression Filter」と呼ばれているフィルタを用いる。以下、固有ベクトルを用いた適応型ＭＴＩフィルタである「Eigenvector Regression Filter」を、「固有ベクトル型ＭＴＩフィルタ」と記載する。

固有ベクトル型ＭＴＩフィルタは、相関行列から固有ベクトルを計算し、計算した固有ベクトルから、クラッタ成分抑制処理に用いる係数を計算する。この方法は、主成分分析や、カルーネン・レーベ変換（Karhunen-Loeve transform）、固有空間法で使われている手法を応用したものである。

固有ベクトル型ＭＴＩフィルタを用いる第１の実施形態に係るドプラ処理部１４は、図１に例示するように、相関行列計算部１４１と、フィルタ行列計算部１４２と、ＭＴＩフィルタ処理部１４３と、推定部１４４とを有する。相関行列計算部１４１は、同一位置（同一サンプル点）の連続した反射波データのデータ列から、走査範囲の相関行列を計算する。フィルタ行列計算部１４２は、例えば、相関行列の固有値及び当該固有値に対応する固有ベクトルを計算する。そして、フィルタ行列計算部１４２は、例えば、各固有値の大きさに基づいて各固有ベクトルを並べた行列のランクを低減した行列を、クラッタ成分を抑制するフィルタ行列として計算する。

ＭＴＩフィルタ処理部１４３は、フィルタ行列を用いて、同一位置（同一サンプル点）の連続した反射波データのデータ列から、クラッタ成分が抑制され、血流に由来する血流信号が抽出されたデータ列を出力する。推定部１４４は、ＭＴＩフィルタ処理部１４３が出力したデータを用いた自己相関演算等の演算を行なって、血流情報を推定し、推定した血流情報をドプラデータとして出力する。なお、第１の実施形態に係るドプラ処理部１４が行なう具体的な処理については、後に詳述する。

ここで、図１に例示するＢモード処理部１３及びドプラ処理部１４は、２次元の反射波データ及び３次元の反射波データの両方について処理可能である。すなわち、Ｂモード処理部１３は、２次元の反射波データから２次元のＢモードデータを生成し、３次元の反射波データから３次元のＢモードデータを生成する。また、ドプラ処理部１４は、２次元の反射波データから２次元のドプラデータを生成し、３次元の反射波データから３次元のドプラデータを生成する。

画像生成部１５は、Ｂモード処理部１３及びドプラ処理部１４が生成したデータから超音波画像データを生成する。画像生成部１５は、Ｂモード処理部１３が生成した２次元のＢモードデータから反射波の強度を輝度で表した２次元Ｂモード画像データを生成する。また、画像生成部１５は、ドプラ処理部１４が生成した２次元のドプラデータから血流情報が映像化された２次元ドプラ画像データを生成する。２次元ドプラ画像データは、速度画像データ、分散画像データ、パワー画像データ、又は、これらを組み合わせた画像データである。画像生成部１５は、ドプラ画像データとして、血流情報がカラーで表示されるカラードプラ画像データを生成したり、１つの血流情報がグレースケールで表示されるドプラ画像データを生成したりする。

ここで、画像生成部１５は、一般的には、超音波走査の走査線信号列を、テレビ等に代表されるビデオフォーマットの走査線信号列に変換（スキャンコンバート）し、表示用の超音波画像データを生成する。具体的には、画像生成部１５は、超音波プローブ１による超音波の走査形態に応じて座標変換を行なうことで、表示用の超音波画像データを生成する。また、画像生成部１５は、スキャンコンバート以外に、種々の画像処理として、例えば、スキャンコンバート後の複数の画像フレームを用いて、輝度の平均値画像を再生成する画像処理（平滑化処理）や、画像内で微分フィルタを用いる画像処理（エッジ強調処理）等を行なう。また、画像生成部１５は、超音波画像データに、種々のパラメータの文字情報、目盛り、ボディーマーク等を合成する。

すなわち、Ｂモードデータ及びドプラデータは、スキャンコンバート処理前の超音波画像データであり、画像生成部１５が生成するデータは、スキャンコンバート処理後の表示用の超音波画像データである。なお、Ｂモードデータ及びドプラデータは、生データ（Raw Data）とも呼ばれる。画像生成部１５は、スキャンコンバート処理前の２次元超音波画像データから、表示用の２次元超音波画像データを生成する。

更に、画像生成部１５は、Ｂモード処理部１３が生成した３次元のＢモードデータに対して座標変換を行なうことで、３次元Ｂモード画像データを生成する。また、画像生成部１５は、ドプラ処理部１４が生成した３次元のドプラデータに対して座標変換を行なうことで、３次元ドプラ画像データを生成する。画像生成部１５は、「３次元のＢモード画像データや３次元ドプラ画像データ」を「３次元超音波画像データ（ボリュームデータ）」として生成する。

更に、画像生成部１５は、ボリュームデータをモニタ２にて表示するための各種の２次元画像データを生成するために、ボリュームデータに対してレンダリング処理を行なう。画像生成部１５が行なうレンダリング処理としては、例えば、断面再構成法（ＭＰＲ：Multi Planer Reconstruction）を行なってボリュームデータからＭＰＲ画像データを生成する処理がある。また、画像生成部１５が行なうレンダリング処理としては、例えば、３次元の情報を反映した２次元画像データを生成するボリュームレンダリング（ＶＲ：Volume Rendering）処理がある。

画像メモリ１６は、画像生成部１５が生成した表示用の画像データを記憶するメモリである。また、画像メモリ１６は、Ｂモード処理部１３やドプラ処理部１４が生成したデータを記憶することも可能である。画像メモリ１６が記憶するＢモードデータやドプラデータは、例えば、診断の後に操作者が呼び出すことが可能となっており、画像生成部１５を経由して表示用の超音波画像データとなる。また、画像メモリ１６は、送受信部１１が出力した反射波データを記憶することも可能である。

内部記憶部１７は、超音波送受信、画像処理及び表示処理を行なうための制御プログラムや、診断情報（例えば、患者ＩＤ、医師の所見等）や、診断プロトコルや各種ボディーマーク等の各種データを記憶する。また、内部記憶部１７は、必要に応じて、画像メモリ１６が記憶する画像データの保管等にも使用される。また、内部記憶部１７が記憶するデータは、図示しないインターフェースを経由して、外部装置へ転送することができる。また、内部記憶部１７は、外部装置から図示しないインターフェースを経由して転送されたデータを記憶することも可能である。

制御部１８は、超音波診断装置の処理全体を制御する。具体的には、制御部１８は、入力装置３を介して操作者から入力された各種設定要求や、内部記憶部１７から読込んだ各種制御プログラム及び各種データに基づき、送受信部１１、Ｂモード処理部１３、ドプラ処理部１４及び画像生成部１５の処理を制御する。例えば、制御部１８は、送受信部１１を介して超音波プローブ１を制御することで、超音波走査の制御を行なう。通常、ＣＦＭ法では、血流像データであるカラードプラ画像データとともに、組織像データであるＢモード画像データを表示する。かかる表示を行なうため、制御部１８は、第１走査範囲内の血流情報を取得する第１超音波走査を超音波プローブ１に実行させる。第１超音波走査は、例えば、ドプラモードでカラードプラ画像データを収集するための超音波走査である。また、制御部１８は、第１超音波走査とともに、第２走査範囲内の組織形状の情報を取得する第２超音波走査を超音波プローブ１に実行させる。第２超音波走査は、例えば、ＢモードでＢモード画像データを収集するための超音波走査である。

制御部１８は、送受信部１１を介して超音波プローブ１を制御することで、第１超音波走査及び第２超音波走査を実行させる。なお、第１走査範囲と第２走査範囲は、同じ範囲であっても、第１走査範囲が第２走査範囲より小さい範囲であっても、第２走査範囲が第１走査範囲より小さい範囲であっても良い。

また、制御部１８は、画像メモリ１６や内部記憶部１７が記憶する表示用の超音波画像データをモニタ２にて表示するように制御する。なお、装置本体１０に内蔵される送受信部１１等は、集積回路などのハードウェアで構成されることもあるが、ソフトウェア的にモジュール化されたプログラムである場合もある。

また、制御部１８は、カラードプラのデータ列が不等間隔になるようなスキャンを行うように送受信部１１を介して超音波プローブ１を制御する。図２Ａ〜図２Ｄは、第１の実施形態に係る制御部１８による不等間隔スキャンの一例を示す図である。例えば、図２Ａ〜図２Ｄに示す例において、横方向に示す各升目はラスタ方向を示し、縦方向の各升目は送受信順序を示す。言い換えると、横方向の各升目は走査線の位置を示し、縦方向の各升目は同一の位置における時間経過を示す。なお、各升目の時間単位は同一であり、例えば「Ｔ」とする。

ここで、図２Ａ〜図２Ｄにおいて、ドット表示の升目は超音波を送受信していることを示し、ハッチング表示の升目は超音波を送受信しているが演算には使用しないダミーデータであることを示す。このダミーデータとは、残留多重エコーの条件をそろえるためのものである。すなわち、不等間隔スキャンにおいて、ドット表示の升目の信号のみが演算に使用される。

説明の便宜上、送受信している場合を「○」、ダミーデータである場合を「△」、送受信していない場合を「−」で表すと、図２Ａに示す例では、左端のラスタにおいて、「△○−−△○−−△○○○○−−△○−−△○」であるパターンの不等間隔スキャンが行われることを示す。また、図２Ａに示す例では、左端のラスタにおいて１３番目の送受信後に、隣の（左から２番目の）ラスタにおいて、左端のラスタと同じパターンの不等間隔スキャンが開始することを示す。

また、図２Ｂに示す例では、左端のラスタにおいて、「−○−−−−−○−−−−−○−○−○−○−−−−−○−−−−−○」であるパターンの不等間隔スキャンが行われることを示す。また、図２Ｂに示す例では、左から２番目のラスタにおいて、「△−○−−−△−○−−−△−○−○−○−○−−−△−○−−−△−○」であるパターンの不等間隔スキャンが行われることを示す。

また、図２Ｃに示す例では、左端のラスタにおいて、「△○○○−−−−△○○○−−−−△○○○」であるパターンの不等間隔スキャンが行われることを示す。また、図２Ｃに示す例では、左端のラスタにおいて４番目の送受信後に、隣の（左から２番目の）ラスタにおいて、左端のラスタと同じパターンの不等間隔スキャンが開始することを示す。

また、図２Ｄに示す例では、左端のラスタにおいて、「−○−○−○−−−−−−−−−○−○−○−−−−−−−−−○−○−○」であるパターンの不等間隔スキャンが行われることを示す。また、図２Ｂに示す例では、左から２番目のラスタにおいて、「△−○−○−○−−−−−−△−−○−○−○−−−−−−−△−○−○−○」であるパターンの不等間隔スキャンが行われることを示す。なお、制御部１８による不等間隔スキャンのパターンは、図２Ａ〜図２Ｄに図示したものに限定されるものではなく、変更可能である。

図３は、第１の実施形態に係る制御部１８による同一ラスタ上での不等間隔スキャンの一例を示す図である。図３では、図２Ａ〜図２Ｄに示すスキャンパターンそれぞれについて、同一のラスタ上での受信信号の順序を数字で示している。図３の（Ａ）は、図２Ａのスキャンパターンに対応し、図３の（Ｂ）は、図２Ｂのスキャンパターンに対応し、図３の（Ｃ）は、図２Ｃのスキャンパターンに対応し、図３の（Ｄ）は、図２Ｄのスキャンパターンに対応する。なお、図３に示す例では、ダミーデータを図示せず、演算に使用される受信信号のみを図示する。言い換えると、図３は、１つのラスタ上の同一地点でのデータ列を示す。また、図２Ａ〜図２Ｄと同様に、各升目の時間単位は同一であり、例えば「Ｔ」とする。

図３の（Ａ）に示す例では、８つの受信信号を得る場合の一例を示す。言い換えると、１フレームのドプラデータ（血流情報）の推定計算に使用するデータ長が８である場合を示す。図３の（Ａ）に示すように、１番目の受信信号と２番目の受信信号とのデータ間隔、２番目の受信信号と３番目の受信信号とのデータ間隔、６番目の受信信号と７番目の受信信号とのデータ間隔、及び７番目の受信信号と８番目の受信信号とのデータ間隔が同一であり、４Ｔとなる。一方で、３番目の受信信号と４番目の受信信号とのデータ間隔、４番目の受信信号と５番目の受信信号とのデータ間隔、及び５番目の受信信号と６番目の受信信号とのデータ間隔が同一であり、Ｔとなる。

また、図３の（Ｂ）に示す例では、８つの受信信号を得る場合の他の一例を示す。図３の（Ｂ）に示すように、１番目の受信信号と２番目の受信信号とのデータ間隔、２番目の受信信号と３番目の受信信号とのデータ間隔、６番目の受信信号と７番目の受信信号とのデータ間隔、及び７番目の受信信号と８番目の受信信号とのデータ間隔が同一であり、６Ｔとなる。一方で、３番目の受信信号と４番目の受信信号とのデータ間隔、４番目の受信信号と５番目の受信信号とのデータ間隔、及び５番目の受信信号と６番目の受信信号とのデータ間隔が同一であり、２Ｔとなる。

図３の（Ｃ）に示す例では、９つの受信信号を得る場合の一例を示す。図３の（Ｃ）に示すように、１番目の受信信号と２番目の受信信号とのデータ間隔、２番目の受信信号と３番目の受信信号とのデータ間隔、４番目の受信信号と５番目の受信信号とのデータ間隔、５番目の受信信号と６番目の受信信号とのデータ間隔、７番目の受信信号と８番目の受信信号とのデータ間隔、及び８番目の受信信号と９番目の受信信号とのデータ間隔が同一であり、Ｔとなる。一方で、３番目の受信信号と４番目の受信信号とのデータ間隔、及び６番目の受信信号と７番目の受信信号とのデータ間隔が同一であり、６Ｔとなる。

図３の（Ｄ）に示す例では、９つの受信信号を得る場合の他の一例を示す。図３の（Ｄ）に示すように、１番目の受信信号と２番目の受信信号とのデータ間隔、２番目の受信信号と３番目の受信信号とのデータ間隔、４番目の受信信号と５番目の受信信号とのデータ間隔、５番目の受信信号と６番目の受信信号とのデータ間隔、７番目の受信信号と８番目の受信信号とのデータ間隔、及び８番目の受信信号と９番目の受信信号とのデータ間隔が同一であり、２Ｔとなる。一方で、３番目の受信信号と４番目の受信信号とのデータ間隔、及び６番目の受信信号と７番目の受信信号とのデータ間隔が同一であり、１０Ｔとなる。以下では、ドプラ処理部１４において、この不等間隔なデータ列に対して固有ベクトル型ＭＴＩフィルタを掛ける方法を説明する。

以上、第１の実施形態に係る超音波診断装置の全体構成について説明した。かかる構成のもと、第１の実施形態に係る超音波診断装置は、不等間隔のデータ列に対して、ＭＴＩフィルタを用いて推定された血流情報（ドプラデータ）によりＣＦＭ法を行なう。

ところで、従来技術に係る超音波診断装置は、不等間隔のデータ列に対して多項式の最小２乗法フィッティングを行うことでクラッタを近似し、原信号からこの近似信号を減算することで血流信号を抽出する。多項式の最小２乗法フィッティングでは、データ間隔の短いパルスペア間でラグ１の自己相関演算を行って速度を計算することにより高い折り返し速度にできると同時に、少ないデータ数で長い観測時間の信号が得られるのでよりＭＴＩフィルタのカットオフ周波数を下げることが可能になり低流速の血流まで組織から分離して抽出することが可能である。また、不等間隔データ列を等間隔に分割することなく一度に処理するために、ブラインド周波数によるノイズも発生しない。つまり、低流速から高流速までの広い範囲の血流速度を折り返しなく表示することが可能になる。

しかしながら、多項式の最小２乗法フィッティングでは、フィルタ係数は事前に決定されており、組織の動きが大きい場合にモーションアーティファクトが発生することを避けることができなかった。なお、モーションアーティファクトを効果的に低減する方法としては、前述の入力信号に応じて適応的にＭＴＩフィルタ特性を変化させる方法（Adaptive MTI Filter）が有効であるが、不等間隔なデータ列にAdaptive MTI Filterを適応した例は知られていない。

このようなことから、第１の実施形態に係る超音波診断装置は、不等間隔データ列が入力されるカラードプラ法において適応的にＭＴＩフィルタ係数を変化させることによりモーションアーティファクトを低減する。例えば、第１の実施形態に係る超音波診断装置は、同一走査線において送信間隔を一定にせずに超音波を送信することで発生した反射波に基づいて生成された反射波データの集合である第１のデータ列を用いて相関行列を計算する。また、第１の実施形態に係る超音波診断装置は、相関行列を用いた主成分分析の結果に基づいて、フィルタ係数を計算する。続いて、第１の実施形態に係る超音波診断装置は、第１のデータ列に含まれる、同一走査線上に存在した移動体で反射された超音波の反射波に由来する反射波データの集合である第２のデータ列を、フィルタ係数を用いて抽出する。また、第１の実施形態に係る超音波診断装置は、抽出された第２のデータ列に基づいて、移動体の移動体情報を推定する。続いて、第１の実施形態に係る超音波診断装置は、移動体情報に基づいて超音波画像データを生成する。そして、第１の実施形態に係る超音波診断装置は、超音波画像データをモニタ２に表示させる。以下では、第１の実施形態に係る超音波診断装置の詳細について説明する。なお、以下では、移動体情報の一例として、血流情報を説明する。

ドプラ処理部１４は、不等間隔データ列の相関行列を計算してその固有値、固有ベクトルを求める。固有値の大きい固有ベクトルでクラッタを近似して原信号から減算することで血流信号を得る。固有値、固有ベクトルは入力信号によって変化するので、ＭＴＩフィルタの特性は入力信号によって変化する。これにより、モーションアーティファクトのような大きな組織の動きでも抑圧できる。なお、この方法は、主成分分析或いはカルーネン・レーベ変換と呼ばれる手法でクラッタを近似すると考えることができる。この手法では、エネルギーが上位固有値に集中する。カラードプラのデータ列の中でエネルギーが大きいものはクラッタである。従って固有値の大きい固有ベクトルで近似した信号をクラッタとみなすことができる。そしてこの手法ではクラッタの発生源となる組織の動きが変化しても分離して抽出できるので、モーションアーティファクトを除去することが可能になる。

以下、ドプラ処理部１４が有する、相関行列計算部１４１、フィルタ行列計算部１４２、ＭＴＩフィルタ処理部１４３、及び推定部１４４における処理の詳細について説明する。まず、第１の実施形態に係る相関行列計算部１４１は、複数の走査線で形成される走査範囲（第１走査範囲）での超音波送受信により複数フレームにわたって収集された同一位置の反射波データのデータ列から、走査範囲（第１走査範囲）の相関行列を計算する。具体的には、第１の実施形態に係る相関行列計算部１４１は、走査範囲（第１走査範囲）にわたり、１回の超音波送受信を行なうごとに、走査線を切り替える走査形態を繰り返すことで収集された同一位置の連続した反射波データのデータ列を、走査範囲の相関行列を計算するデータ列として用いる。より具体的には、第１の実施形態に係る相関行列計算部１４１は、走査範囲（第１走査範囲）での超音波送受信を各走査線で１回とする走査形態を繰り返すことで収集された同一位置の連続した反射波データのデータ列を、走査範囲の相関行列を計算するデータ列として用いる。ここで、相関行列計算部１４１は、同一走査線において送信間隔を一定にせずに超音波を送信することで発生した反射波に基づいて生成された反射波データの集合である第１のデータ列を用いて相関行列を計算する。

具体的には、相関行列計算部１４１は、式（１）に示すように、１つのラスタ上の同一地点でのデータ列を列ベクトルｘで表す。言い換えると、列ベクトル「ｘ」は、ある走査位置におけるデータ列を示す。この列ベクトルｘには、要素ｘ_１，ｘ_２，・・・ｘ_Ｌが含まれる。ここで列ベクトルｘの各要素に付された添字の数字は、図３に示す升目の中の数字に対応する。例えば、図３の（Ａ）の場合、Ｌ＝８である。「Ｌ」は、列ベクトル「ｘ_ｍ」の長さを示し、１フレームのドプラデータ（血流情報）の推定計算に使用するデータ長である。

次に、相関行列計算部１４１は、以下に示す式（２）により相関行列「Ｒ_ｘｘ」を計算する。すなわち、相関行列計算部１４１は、式（２）により、空間的に別の位置にあるデータ列を加算平均することでアンサンブル平均を取る。言い換えると、相関行列計算部１４１は、式（２）により、複数のサンプル点それぞれで、データ列の自己相関行列を計算し、複数のサンプル点それぞれの自己相関行列の平均を計算する。ここで、ｘ_ｍの添字ｍは空間の位置を表し、位置ｍの総数をＭとする。なお、位置「ｍ」は、２次元スキャンの場合は、２次元座標系で示され、３次元スキャンの場合は、３次元座標系で示される。Ｈは行列の各要素の複素共役を取った行列の転置（エルミート転置）を表す。相関行列「Ｒ_ｘｘ」は、式（２）により、Ｌ行Ｌ列の行列となる。上述したように、相関行列が計算されるデータ列のデータ長「Ｌ」は、任意に変更可能である。

なお、相関行列計算部１４１は、空間的に別の位置にあるデータとしては、カラードプラのスキャン範囲全体のデータを使用して１つの共分散行列を計算しても良いし、スキャン範囲をブロックに分割してそれぞれに共分散行列を計算しても良い。

次に、フィルタ行列計算部１４２は、相関行列を用いた主成分分析の結果に基づいて、フィルタ係数を計算する。言い換えると、フィルタ行列計算部１４２は、相関行列を用いて主成分分析を行ない、クラッタ成分を主成分として近似し低減する行列演算を行なうことで組織からのクラッタを抑圧するフィルタ係数を計算する。例えば、フィルタ行列計算部１４２は、Ｌ行Ｌ列の行列Ｒ_ｘｘの固有値と固有ベクトルを計算する。固有値の大きい順に固有ベクトルを列ベクトルとして左から並べた行列をＶとする。信号を上位Ｋ個の主成分で近似すると以下のような式になる。ここで、式（３）では、中央の行列は対角行列で、対角要素は左上から１がＫ個あり残りは０である場合を示す。ここで、式（３）に示す「Ｖ^Ｈ」は、「Ｖ」の複素共役転置行列である。また、式（３）の右辺において、「Ｖ」と「Ｖ^Ｈ」との間の行列は、Ｌ行Ｌ列の対角行列である。

原信号から式（３）を減算すると、式（４）となる。但し、Ｉは単位行列を表す。つまり、式（５）に示す行列Ｗを入力データ列ｘに掛ければＭＴＩフィルタが構成できることになる。ＭＴＩフィルタ行列「Ｗ」は、式（５）により、Ｌ行Ｌ列の行列となる。

式（５）の右辺における対角行列において、０の個数をＫとする。これは行列のランクを低減する処理となり、Ｋをランクカット数と呼ぶことにする。ランクカット数Ｋは固定でも、組織が動くとその固有値は大きくなるので、モーションアーティファクトをそのランクカット数を用いて最大限に除去する。しかし、組織の動きが大きい場合には組織の動きに対応する固有ベクトルの数が大きくなってくるので、ランクカット数を大きくすることが望ましい。装置上のスイッチ等によって操作者がランクカット数を変更しても良いが、固有値の大きさから適応的にランクカット数を変化させる方法が好適である。

例えば、フィルタ行列計算部１４２は、予め設定された値、或いは、操作者が指定した値により、低減される主成分の数、すなわち、ランクカット数の値を決定する。しかし、心臓や血管等、拍動により移動速度が時間により変化する組織が走査範囲内に含まれる場合、ランクカット数の値は、固有値の大きさから適応的に決定されることが好適である、すなわち、フィルタ行列計算部１４２は、相関行列の固有値の大きさに応じて、低減する主成分の数を変更する。本実施形態では、フィルタ行列計算部１４２は、固有値の大きさに応じて、低減するランク数を変更する。

固有値の大きさからランクカット数を適応的に決定する論理は、超音波走査が行なわれる部位に応じて最適化する必要があるが、例えば、フィルタ行列計算部１４２は、図４に示す閾値と、図５に示すアルゴリズムとにより、ランクカット数を決定する。図５に示すアルゴリズムは、２番目に大きい固有値を１番大きい固有値で割った値に基づいてランクカット数を決定するアルゴリズムである。図４及び図５は、第１の実施形態に係るランクカット数の決定処理の一例を説明するための図である。

まず、フィルタ行列計算部１４２は、相関行列「Ｒ_ｘｘ」の固有値を大きい順に並べた配列の中で、ｋ番目の固有値を「ｅｉｇ（ｋ）」と定義する。「ｋ」は、「１≦ｋ≦Ｌ」の整数となる。そして、フィルタ行列計算部１４２は、ｋ番目に大きい固有値を１番大きい固有値「ｅｉｇ（１）」で割った値をｄＢ単位で表した値「ｅｉｇｄＢ（ｋ）」を、以下の式（６）により計算する。

なお、式（６）において、「ａｂｓ」は絶対値を演算する関数となる。式（６）において、「ｋ＝２」とした「ｅｉｇｄＢ（２）」は、２番目に大きい固有値を１番大きい固有値「ｅｉｇ（１）」で割った値をｄＢ単位で表した値となる。

また、フィルタ行列計算部１４２は、Ｌ個の固有値が得られることから、ランクカット数決定用のＬ個の閾値（ＴＨ_ｉ、１≦ｉ≦Ｌ）を用いる。ＴＨ_ｉは、ｉの値によって変化する値が設定される。例えば、「Ｌ＝８」の場合、８個の閾値「ＴＨ_１〜ＴＨ_８」は、図４に示すように設定される。図４では、ＴＨ_１及びＴＨ_２は「１００００００ｄＢ」に設定される。また、図４では、ＴＨ_３は「２０ｄＢ」に設定され、ＴＨ_４は「１５ｄＢ」に設定される。また、図４では、ＴＨ_５は「１０ｄＢ」に設定され、ＴＨ_６は「５ｄＢ」に設定される。また、図４では、ＴＨ_７及びＴＨ_８は「−１ｄＢ」に設定される。図４に例示する閾値を用いた場合、ランクカット数は、以下に説明する図５のアルゴリズムにより、２以上６以下の値となる。図５では、ランクカット数を「RankCut」と示している。

まず、フィルタ行列計算部１４２は、「RankCut＝０」に設定し（ステップＳ１）、「ｉ＝１」に設定する（ステップＳ２）。そして、フィルタ行列計算部１４２は、「ｉ」が「Ｌ」より大きいか、又は、「−ｅｉｇｄＢ（２）」が「ＴＨ_ｉ」より大きいかを判定する（ステップＳ３）。ここで、「ｉ」が「Ｌ」以下であり、かつ、「−ｅｉｇｄＢ（２）」が「ＴＨ_ｉ」以下である場合（ステップＳ３否定）、フィルタ行列計算部１４２は、ランクカット数をインクリメントして「RankCut＝RankCut＋１」とする（ステップＳ４）。

そして、フィルタ行列計算部１４２は、「ｉ＝ｉ＋１」に設定して（ステップＳ５）、ステップＳ３の判定処理を行なう。なお、例えば、最初のステップＳ５の処理の後にステップＳ３の判定処理で用いられる「−ｅｉｇｄＢ（２）」は、Ｌ個の固有値を大きい順に並べた配列から１番大きい固有値を除外した配列において、２番目に大きい固有値を１番大きい固有値「ｅｉｇ（１）」で割った値をｄＢ単位で表した値に「−１」を乗算した値となる。

一方、「ｉ」が「Ｌ」より大きい場合、又は、「−ｅｉｇｄＢ（２）」が「ＴＨ_ｉ」より大きい場合（ステップＳ３肯定）、フィルタ行列計算部１４２は、最新の「RankCut」を、低減するランク数として決定する（ステップＳ６）。例えば、フィルタ行列計算部１４２は、図４の条件下において、ｅｉｇｄＢ（２）＝−１２ｄＢの場合、ランクカット数を４に決定する。

なお、ランクカット数を固有値の大きさから適応的に決定するアルゴリズムは、上記のアルゴリズム以外にも、様々なアルゴリズムにより行なうことができる。これらのアルゴリズムは、例えば、撮影部位に応じて、選択可能である。

フィルタ行列計算部１４２は、図５に例示したアルゴリズムにより、各表示フレームでランクカット数を決定して、ＭＴＩフィルタ行列「Ｗ」を計算する。

そして、ＭＴＩフィルタ処理部１４３は、第１のデータ列に含まれる、同一走査線上に存在した移動体で反射された超音波の反射波に由来する反射波データの集合である第２のデータ列を、フィルタ係数を用いて抽出する。言い換えると、ＭＴＩフィルタ処理部１４３は、フィルタ係数を用いて、同一位置（同一サンプル点）の連続した反射波データのデータ列から、クラッタ成分が抑制され、血流に由来する血流信号が抽出されたデータ列を出力する。本実施形態では、ＭＴＩフィルタ処理部１４３は、フィルタ行列を用いて、同一位置（同一サンプル点）の連続した反射波データのデータ列から、クラッタ成分が抑制され、血流に由来する血流信号が抽出されたデータ列を出力する。具体的には、ＭＴＩフィルタ処理部１４３は、位置「ｍ」の列ベクトル「ｘ_ｍ」を入力データとし、入力データとＭＴＩフィルタ行列「Ｗ」とから、位置「ｍ」の出力データである列ベクトル「ｙ_ｍ」を以下の式（７）により計算する。列ベクトル「ｙ_ｍ」の長さは、「Ｌ」となる。

ＭＴＩフィルタ処理部１４３は、「Ｍ」個のサンプル点それぞれにおいて、式（７）の計算を行なう。これにより、ＭＴＩフィルタ処理部１４３は、「Ｍ」個のサンプル点それぞれの出力データを、推定部１４４に出力する。

推定部１４４は、抽出された移動体からの信号のデータ列に基づいて、移動体の移動体情報を推定する。例えば、推定部１４４は、位置「ｍ」の出力データである列ベクトル「ｙ_ｍ」から、自己相関演算処理及び速度・分散・パワー推定処理を行なって、位置「ｍ」の血流情報を推定する。まず、推定部１４４は、列ベクトル「ｙ_ｍ」から、ラグ０及びラグ１の自己相関値の演算を行なう。ラグ０を「Ｃ_０」とし、ラグ１を「Ｃ_１」とすると、推定部１４４は、以下の式（８）によりラグ０「Ｃ_０」を計算する。

また、推定部１４４は、ラグ１「Ｃ_１」の計算には不等間隔データ列で最も短いパルスペアを使用する。例えば、図３（Ａ）に示す場合、３番目の受信信号と４番目の受信信号とのデータ間隔、４番目の受信信号と５番目の受信信号とのデータ間隔、及び５番目の受信信号と６番目の受信信号とのデータ間隔がＴであり最も短いパルスペアとなる。かかる場合、推定部１４４は、以下の式（９）でラグ１を計算する。

また、図３（Ｂ）に示す場合は、図３（Ｂ）の最も短いパルスペアは図３（Ａ）の最も短いパルスペアと同じなので、上記と同じ式（９）になる。なお、ＭＴＩフィルタの特性により速度の精度は短いパルスペアの端部で落ちる。このため、推定部１４４は、式（９）においてｉ＝４である１つのパルスペアだけを用いることでＭＴＩフィルタの影響を小さくしてもよい。なお、かかる場合、加算による安定性は落ちる。

また、例えば、推定部１４４は、図３（Ｃ）に示す場合、１番目の受信信号と２番目の受信信号とのデータ間隔、２番目の受信信号と３番目の受信信号とのデータ間隔、４番目の受信信号と５番目の受信信号とのデータ間隔、５番目の受信信号と６番目の受信信号とのデータ間隔、７番目の受信信号と８番目の受信信号とのデータ間隔、及び８番目の受信信号と９番目の受信信号とのデータ間隔がＴであり最も短いパルスペアとなる。かかる場合、推定部１４４は、以下の式（１０）でラグ１を計算する。

同様に、推定部１４４は、図３（Ｄ）の最も短いパルスペアは図３（Ｃ）の最も短いパルスペアと同じなので、上記と同じ式（１０）になる。なお、推定部１４４は、式（１０）のｉ＝４，５である２つのパルスペアを用いることでＭＴＩフィルタの影響を小さくするようにしてもよい。

なお、式（８）〜式（１０）では、列ベクトル「ｙ_ｍ」の位置を表す添字「ｍ」を省略し、代わりに列ベクトル「ｙ_ｍ」の要素「ｉ」を添字で表している。また、式（９）及び式（１０）において、上付きのアスタリスク「＊」は、複素共役を示す。推定部１４４は、「Ｍ」個のサンプル点それぞれの「Ｃ_０」及び「Ｃ_１」を計算する。

そして、推定部１４４は、「Ｃ_０」及び「Ｃ_１」から速度「Ｖ」を以下の式（１１）により計算し、「Ｃ_０」及び「Ｃ_１」から分散「Ｔ」を以下の式（１２）により計算し、「Ｃ_０」からパワー「Ｐ」を以下の式（１３）により計算する。なお、推定部１４４は、移動体情報として、移動体からの信号のデータ列の少なくとも一部を用いて移動体のパワーを計算してもよい。

なお、式（１１）の「ａｔａｎ２」は、「−π〜＋π」までの角度を出力する「アークタンジェント関数」であり、「ｉｍａｇ」は、複素数から虚数部のみを出力する関数であり、「ｒｅａｌ」は、複素数から実数部のみを出力する関数である。パワーは、自己相関演算によりラグ０として求められ、速度及び分散は、自己相関演算の結果を周波数解析することで求められる。式（１１）で得られる速度は−π〜＋πまでの角度で表現されているが、これをｍ／ｓの単位の速度Ｖ’に直すのには、不等間隔データ列の最も短いパルスペアの間隔をＴとすると、式（１４）となる。また、Ｃは超音波の生体での音速、ｆ_０は超音波の中心周波数である。折り返し速度Ｖ’_ｍａｘは、式（１５）となる。

ここで、例えば、図３（Ａ）ですべてのデータ間隔が１番目と２番目のデータ間隔であったとすると、周期は４Ｔとなる。かかる場合、式（１５）より折り返し速度は１／４になってしまう。ここで、図３（Ａ）に示す３番目と４番目のデータ間隔、４番目と５番目のデータ間隔、５番目と６番目のデータ間隔ように、周期Ｔのデータを挿入することで、折り返し速度を周期４Ｔの４倍に向上することができる。

また、データ間隔がＴのままで図３（Ａ）のデータ長を維持するには、図３（Ａ）に示す升目の数に相当する２３個のデータが必要となる。図３（Ａ）の実際のデータ数は８個なので、同じ観測時間において１／２．８７５のデータ数に減らすことができる。このように、データ数を少なくすることは、フレームレート向上に寄与する。

また、観測時間が長いほど、横軸をｍ／ｓの速度とした場合に急峻な肩特性のＭＴＩフィルタが構成できる。すなわち、よりクラッタを抑圧可能になる。つまり、例えば、図２Ａのような不等間隔データ列のスキャンを行うことで、高フレームレート、高流速検出能、低流速検出能という３つの相反する要求を満足することが可能になる。

推定部１４４は、「Ｍ」個のサンプル点それぞれの「Ｖ、Ｔ、Ｐ」を計算する。そして、推定部１４４は、「Ｍ」個のサンプル点それぞれの「Ｖ、Ｔ、Ｐ」を１フレーム分のドプラデータとして、画像生成部１５に出力する。

そして、画像生成部１５は、フィルタ係数を用いて推定された血流情報（ドプラデータ）から、超音波画像データ（カラードプラ画像データ）を生成する。本実施形態では、画像生成部１５は、ＭＴＩフィルタ行列を用いて推定された血流情報（ドプラデータ）から、超音波画像データ（カラードプラ画像データ）を生成する。そして、制御部１８は、超音波画像データ（カラードプラ画像データ）をモニタ２に表示させる。

図６及び図７は、第１の実施形態に係る超音波画像データの一例を示す図である。図６及び図７では、腎臓血流の一例を示す。また、図６では、腎臓の動きが少ない場合の画像を示し、図７では、腎臓の動きが大きい場合の画像を示す。

図６及び図７において、上段にパワー画像を示し、下段に速度画像を示す。また、列の左側から「オリジナルデータ」、「間引きデータ」、「多項式フィッテイングデータ」、「固有ベクトルデータ」の順に示す。ここで、「オリジナルデータ」は、３２データの等間隔データ列である。また、「間引きデータ」は、等間隔の３２データのデータ列を４データ毎に等間隔に間引いて８データを使用して演算した画像である。

「多項式フィッテイングデータ」は、不等間隔データ列を３次多項式で近似して原信号から減算した画像である。「固有ベクトルデータ」は、第１の実施形態に係る固有ベクトル型ＭＴＩフィルタで、データ長が８である不等間隔データ列を処理した画像である。ここで、「多項式フィッテイングデータ」及び「固有ベクトルデータ」では、中央だけ密にして間引いたデータ長が８である図３（Ａ）に示すスキャンパターンの不等間隔データ列を使用するものとする。なお、固有ベクトル型ＭＴＩフィルタのランクカット数は４に固定される。

図６に示すように、「オリジナルデータ」の速度画像では、血流の方向によって血管の分岐が表示されている。なお、速度画像のおける血流方向と速度とをスケールで示す。一方、「間引きデータ」では、折り返し速度が「オリジナルデータ」の１／４になるので速度画像で折り返しが発生している。このため、血流の方向が正確に再現できず、血管の分岐が表示されていない。

また、「多項式フィッテイングデータ」及び「固有ベクトルデータ」の速度画像は、ほぼ同等であり、「間引きデータ」の速度画像に比べると折り返しが大幅に少ない画像になっている。また、「多項式フィッテイングデータ」及び「固有ベクトルデータ」の速度画像は、「オリジナルデータ」の速度画像にかなり近い画像が表示されている。

また、図７に示すように、「オリジナルデータ」、「間引きデータ」、及び「多項式フィッテイングデータ」のパワー画像には、あるブロックにモーションアーティファクトが混入している。一方、「固有ベクトルデータ」のパワー画像には、モーションアーティファクトが混入しておらず、高精細な血流画像が表示されている。同様に、「オリジナルデータ」、「間引きデータ」、及び「多項式フィッテイングデータ」の速度画像にも、あるブロックにモーションアーティファクトが混入しているが、「固有ベクトルデータ」の速度画像には、モーションアーティファクトが混入していない。なお、動きが少ないときと同様に、「固有ベクトルデータ」の速度画像は、「間引きデータ」の速度画像に比べると折り返しが大幅に少ない画像になっている。

次に、図８を用いて、第１の実施形態に係る超音波診断装置が行なう処理の一例について説明する。図８は、第１の実施形態に係る超音波診断装置が行なう処理の一例を説明するためのフローチャートである。なお、図８に例示するフローチャートは、第１の実施形態に係る超音波診断装置が行なうドプラ画像データの生成表示処理を説明したフローチャートである。

図８に示すように、第１の実施形態に係る超音波診断装置の制御部１８は、まず１フレーム分の不等間隔データ列が入力されたか否かを判定する（ステップＳ１０１）。ここで、制御部１８は、１フレーム分の不等間隔データ列が入力されたと判定した場合（ステップＳ１０１、Ｙｅｓ）、ステップＳ１０２に移行する。一方、制御部１８は、１フレーム分の不等間隔データ列が入力されたと判定しなかった場合（ステップＳ１０１、Ｎｏ）、ステップＳ１０１を繰り返し実行する。

次に、相関行列計算部１４１は、走査範囲の相関行列を計算し（ステップＳ１０２）、フィルタ行列計算部１４２は、相関行列からＬ組の固有値及び固有ベクトルを計算する（ステップＳ１０３）。

そして、フィルタ行列計算部１４２は、Ｌ個の固有値の大小関係から、低減するランク数を決定し（ステップＳ１０４）、ＭＴＩフィルタ行列を計算する（ステップＳ１０５）。そして、ＭＴＩフィルタ処理部１４３は、ＭＴＩフィルタ処理を行ない（ステップＳ１０６）、推定部１４４は、ＭＴＩフィルタ処理で出力された出力データを用いて、自己相関演算処理を行なう（ステップＳ１０７）。そして、推定部１４４は、自己相関演算処理の結果から血流情報を推定する（ステップＳ１０８）。

そして、画像生成部１５は、血流情報からカラードプラ画像データを生成し（ステップＳ１０９）、制御部１８の制御により、モニタ２は、カラードプラ画像データを表示し（ステップＳ１１０）、処理を終了する。

上述したように、第１の実施形態では、固有ベクトル型ＭＴＩフィルタを用いた血流情報の推定処理を、不等間隔データ列スキャンによる高フレームレート超音波走査と組み合わせて行なう。これにより、第１の実施形態では、走査範囲全体で１つの相関行列を計算して、画像全体に同一の固有ベクトル型ＭＴＩフィルタを掛けることができる。これにより第１の実施形態では、不等間隔データ列に対して固有ベクトルを用いた適応型ＭＴＩフィルタによりモーションアーティファクトを大幅に低減した画像を提供することができる。そして、不等間隔データ列のデータ間隔の最も短いパルスペアのラグ１の自己相関関数から速度を計算することにより、折り返し速度が高い速度を表示することができる。

また、不等間隔データ列を使用することで、少ないデータ数で低いカットオフ周波数で肩特性の急峻なＭＴＩフィルタを構成できる、すなわち、高フレームレート（＝少ないデータ数）で低流速まで検出可能（＝低いカットオフ周波数で肩特性の急峻なＭＴＩフィルタ）なシステムを実現することができる。

また、第１の実施形態では、固有値の大きさから適応的に決定したランクカット数により計算された固有ベクトル型ＭＴＩフィルタを用いることで、ＭＴＩフィルタの特性を、組織の動きに応じて変化させることができる。その結果、第１の実施形態では、モーションアーティファクトを更に大幅に抑制することができる。

（第２の実施形態）
第１の実施形態では１フレーム内で１つの相関行列を計算したが、領域に分割してその領域毎に相関行列を計算してＭＴＩフィルタ行列を計算して処理を行っても良い。このようなことから、第２の実施形態では、領域に分割してその領域毎に相関行列を計算してＭＴＩフィルタ行列を計算して処理を行うことで、固有ベクトル型ＭＴＩフィルタを用いたドプラ画像データの画質を更に向上される方法について説明する。

第２の実施形態に係る超音波診断装置は、図１を用いて説明した第１の実施形態に係る超音波診断装置と同様の構成となる。ただし、第２の実施形態に係るドプラ処理部１４は、以下に説明する処理により、ＭＴＩフィルタ行列を計算する。

第１の実施形態では、血流を表示する走査範囲全体で１つの相関行列を計算し、画像全体に同一のＭＴＩフィルタを掛けている。一方、血流を表示する走査範囲において、場所により組織の動きが大きく異なる場合、当該走査範囲を複数の処理ブロックに分割し、各処理ブロックで固有ベクトル型ＭＴＩフィルタを計算することで、各処理ブロックで最適なＭＴＩフィルタを掛けることが望ましい。しかし、かかる場合、フィルタの特性差により、ドプラ画像データには、処理ブロック間で不連続な境界が発生してしまう。

そこで、第２の実施形態では、各処理ブロックで計算した相関行列を空間的に補間することで、固有ベクトル型ＭＴＩフィルタの係数を計算する。

まず、第２の実施形態に係る相関行列計算部１４１は、超音波送受信が行われた走査範囲を分割した複数の範囲（複数の処理ブロック）それぞれで、相関行列を計算する。そして、相関行列計算部１４１は、複数の処理ブロックそれぞれより細分化された複数の細分化範囲（複数の細分化処理ブロック）それぞれの相関行列を、複数の処理ブロックそれぞれの相関行列を用いた補間処理により計算する。

そして、第２の実施形態に係るフィルタ行列計算部１４２は、複数の細分化処理ブロックそれぞれの相関行列から、各細分化処理ブロックのフィルタ係数を計算する。具体的には、フィルタ行列計算部１４２は、複数の細分化処理ブロックそれぞれの相関行列から、各細分化処理ブロックのフィルタ行列を計算する。そして、第２の実施形態に係るＭＴＩフィルタ処理部１４３は、各細分化処理ブロックのフィルタ係数を用いて、各サンプル点の血流信号を抽出する。具体的には、ＭＴＩフィルタ処理部１４３は、各細分化処理ブロックのフィルタ行列を用いて、各サンプル点の血流信号を抽出する。そして、第２の実施形態に係る推定部１４４は、複数の細分化範囲それぞれの移動体からの信号のデータ列に基づいて、移動体情報を推定する。例えば、第２の実施形態に係る推定部１４４は、各サンプル点の血流情報を推定する。これにより、画像生成部１５は、複数の細分化処理ブロックそれぞれのフィルタ係数を用いて推定された血流情報から、超音波画像データ（ドプラ画像データ）を生成する。具体的には、画像生成部１５は、複数の細分化処理ブロックそれぞれのフィルタ行列を用いて推定された血流情報から、超音波画像データ（ドプラ画像データ）を生成する。

図９は、第２の実施形態に係る相関行列計算部が行なう処理の一例を示す図である。例えば、走査範囲は、図９に例示するように、４つの処理ブロック（Ｒ１、Ｒ２、Ｒ３、Ｒ４）に分割される。なお、走査範囲の分割数や各処理ブロックは、任意に設定可能である。相関行列計算部１４１は、第１の実施形態で説明した式（２）を用いた計算処理を、処理ブロックＲ１〜Ｒ４それぞれで行なう。これにより、相関行列計算部１４１は、処理ブロックＲ１の相関行列「Ｒ_１」と、処理ブロックＲ２の相関行列「Ｒ_２」と、処理ブロックＲ３の相関行列「Ｒ_３」と、処理ブロックＲ４の相関行列「Ｒ_４」とを計算する。ここで、図９に示すように、処理ブロックＲ１の中心座標を（ｘ_１，ｙ_１）とし、処理ブロックＲ２の中心座標を（ｘ_２，ｙ_１）とし、処理ブロックＲ３の中心座標を（ｘ_１，ｙ_２）とし、処理ブロックＲ４の中心座標を（ｘ_２，ｙ_２）とする。

例えば、相関行列計算部１４１は、複数の細分化処理ブロックそれぞれが、走査範囲内の個々のサンプル点として設定された場合、サンプル点Ｒ（ｘ，ｙ）の相関行列を、以下の式（１６）を用いたバイリニア補間により計算する。

そして、フィルタ行列計算部１４２は、式（５）により、Ｍ個のサンプル点それぞれの相関行列を用いて、Ｍ個のサンプル点それぞれのＭＴＩフィルタ行列を計算する。なお、フィルタ行列計算部１４２は、固有値の大きさに応じて、Ｍ個のサンプル点それぞれで、ランクカット数を決定する。そして、フィルタ処理部１４３は、各サンプル点の血流信号を、該当するＭＴＩフィルタ行列を用いて抽出し、推定部１４４は、各サンプル点の血流情報を推定する。

図１０は、第２の実施形態の効果を説明するための図である。図１０に例示する左図は、血流表示が行なわれる走査範囲を１６個の処理ブロックに分割し、各処理ブロックでＭＴＩフィルタ行列を計算した場合に表示されるドプラ画像データである。一方、図１０に例示する右図は、血流表示が行なわれる走査範囲を１６個の処理ブロックに分割し、各サンプル点の相関行列を補間処理で計算し、更に、各サンプル点でＭＴＩフィルタ行列を計算した場合に表示されるドプラ画像データである。

図１０の左図に例示するように、上記の補間処理を行なわない「補間無し」のドプラ画像データには、１６個の処理ブロック間に、不連続な境界が発生している。一方、図１０の右図に例示するように、上記の補間処理を行なった「補間有り」のドプラ画像データには、１６個の処理ブロック間に、不連続な境界が発生していない。

なお、上記では、細分化処理ブロックが１個のサンプル点で形成される場合について説明した。しかし、第２の実施形態は、処理負荷を軽減するために、例えば、１０個のサンプル点で走査範囲を分割した複数の細分化処理ブロックそれぞれでＭＴＩフィルタ行列を計算しても良い。

次に、図１１を用いて、第２の実施形態に係る超音波診断装置が行なう処理の一例について説明する。図１１は、第２の実施形態に係る超音波診断装置が行なう処理の一例を説明するためのフローチャートである。なお、図１１に例示するフローチャートは、第２の実施形態に係る超音波診断装置が行なうドプラ画像データの生成表示処理を説明したフローチャートである。また、図１１では、細分化処理ブロックが１つのサンプル点として設定された場合を例示している。

図１１に示すように、第２の実施形態に係る超音波診断装置の制御部１８は、まず１フレーム分の不等間隔データ列が入力されたか否かを判定する（ステップＳ２０１）。ここで、制御部１８は、１フレーム分の不等間隔データ列が入力されたと判定した場合（ステップＳ２０１、Ｙｅｓ）、ステップＳ２０２に移行する。一方、制御部１８は、１フレーム分の不等間隔データ列が入力されたと判定しなかった場合（ステップＳ２０１、Ｎｏ）、ステップＳ２０１を繰り返し実行する。

次に、相関行列計算部１４１は、各処理ブロックで相関行列を計算し（ステップＳ２０２）、各サンプル点の相関行列を補間処理により計算する（ステップＳ２０３）。そして、フィルタ行列計算部１４２は、サンプル点ごとに、相関行列からＬ組の固有値及び固有ベクトルを計算する（ステップＳ２０４）。

そして、フィルタ行列計算部１４２は、サンプル点ごとに、低減するランク数を決定し（ステップＳ２０５）、サンプル点ごとに、ＭＴＩフィルタ行列を計算する（ステップＳ２０６）。そして、ＭＴＩフィルタ処理部１４３は、ＭＴＩフィルタ処理を行ない（ステップＳ２０７）、推定部１４４は、ＭＴＩフィルタ処理で出力された出力データを用いて、自己相関演算処理を行なう（ステップＳ２０８）。そして、推定部１４４は、自己相関演算処理の結果から血流情報を推定する（ステップＳ２０９）。

そして、画像生成部１５は、血流情報からカラードプラ画像データを生成し（ステップＳ２１０）、制御部１８の制御により、モニタ２は、カラードプラ画像データを表示し（ステップＳ２１１）、処理を終了する。

上述したように、第２の実施形態では、各処理ブロックの相関行列を用いた補間処理により各細分化処理ブロックの相関関数を計算することで、各細分化処理ブロックで最適化されたＭＴＩフィルタ行列を計算することができる。その結果、第２の実施形態では、場所により組織の動きが大きく異なる場合であっても、モーションアーティファクトが大幅に抑制されたドプラ画像データを生成表示することができる。

なお、第２の実施形態で説明した補間処理が行なわれる場合、ドプラ用の第１超音波走査は、高フレームレート用超音波走査以外の走査形態により実行されても良い。すなわち、第２の実施形態は、血流情報の推定が可能な反射波データのデータ列が収集可能であれば、任意の走査形態を適用可能である。例えば、交互スキャンを実行する場合でも、第２の実施形態で説明した補間処理を行なうことで、走査ブロック及び処理ブロックに起因するアーティファクトが発生することを大幅に抑制することができる。

かかる場合、相関行列計算部１４１は、任意の走査形態で超音波送受信が行なわれた走査範囲を複数の処理ブロックに分割し、各範囲で収集された同一位置の反射波データのデータ列から、複数の処理ブロックそれぞれの相関行列を計算する。そして、相関行列計算部１４１は、更に、複数の細分化処理ブロックそれぞれの相関行列を、複数の処理ブロックそれぞれの相関行列を用いた補間処理により計算する。そして、フィルタ行列計算部１４２は、各細分化範囲において該当する相関行列を用いて主成分分析を行ない、クラッタ成分を主成分として近似し低減する行列演算を行なうことで組織からのクラッタを抑圧するフィルタ係数を、前記複数の細分化範囲それぞれで計算する。具体的には、フィルタ行列計算部１４２は、各細分化処理ブロックにおいて該当する相関行列の固有値及び当該固有値に対応する固有ベクトルを計算し、各固有値の大きさに基づいて各固有ベクトルを並べた行列のランクを低減した行列を、クラッタ成分を抑制するＭＴＩフィルタ行列として計算する。これにより、フィルタ行列計算部１４２は、複数の細分化処理ブロックそれぞれのフィルタ行列を計算する。

なお、上記では、第１〜第２の実施形態で説明した画像処理方法が、超音波診断装置で実行される場合について説明した。しかし、第１〜第２の実施形態で説明した画像処理方法は、送受信部１１が出力した反射波データ（ＩＱ信号）を取得可能な画像処理装置において実行される場合であっても良い。

また、上記の実施形態において、図示した各装置の各構成要素は機能概念的なものであり、必ずしも物理的に図示の如く構成されていることを要しない。すなわち、各装置の分散・統合の具体的形態は図示のものに限られず、その全部または一部を、各種の負荷や使用状況などに応じて、任意の単位で機能的または物理的に分散・統合して構成することができる。更に、各装置にて行なわれる各処理機能は、その全部または任意の一部が、ＣＰＵおよび当該ＣＰＵにて解析実行されるプログラムにて実現され、或いは、ワイヤードロジックによるハードウェアとして実現され得る。

また、第１の実施形態〜第２の実施形態で説明した画像処理方法は、予め用意された画像処理プログラムをパーソナルコンピュータやワークステーションなどのコンピュータで実行することによって実現することができる。この画像処理プログラムは、インターネットなどのネットワークを介して配布することができる。また、この画像処理プログラムは、ハードディスク、フレキシブルディスク（ＦＤ）、ＣＤ−ＲＯＭ、ＭＯ、ＤＶＤ、ＵＳＢメモリ及びＳＤカードメモリ等のＦｌａｓｈメモリ等、コンピュータで読み取り可能な非一時的な記録媒体に記録され、コンピュータによって非一時的な記録媒体から読み出されることによって実行することもできる。

以上説明した少なくともひとつの実施形態によれば、不等間隔データ列が入力されるカラードプラ法においてモーションアーティファクトを低減することができる。

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

１０ 装置本体
１４ ドプラ処理部
１４１ 相関行列計算部
１４２ フィルタ行列計算部
１４３ ＭＩＴフィルタ処理部
１５ 画像生成部
１８ 制御部

Claims (7)

1. 同一走査線において送信間隔を一定にせずに超音波を送信することで発生した反射波に基づいて生成された反射波データの集合である第１のデータ列を用いて相関行列を計算する相関行列計算部と、
   前記相関行列を用いた主成分分析の結果に基づいて、フィルタ係数を計算するフィルタ行列計算部と、
   前記第１のデータ列に含まれる、前記同一走査線上に存在した移動体で反射された前記超音波の反射波に由来する反射波データの集合である第２のデータ列を、前記フィルタ係数を用いて抽出する抽出部と、
   抽出された前記第２のデータ列に基づいて、前記移動体の移動体情報を推定する推定部と、
   前記移動体情報に基づいて超音波画像データを生成する画像生成部と、
   前記超音波画像データを表示部に表示させる制御部と、
   を備えたことを特徴とする超音波診断装置。
2. 前記推定部は、前記移動体情報として、前記第２のデータ列の少なくとも一部を用いて前記移動体のパワーを計算することを特徴とする請求項１に記載の超音波診断装置。
3. 前記推定部は、前記移動体情報として、前記第２のデータ列の間隔の短いパルスペアから前記移動体の速度を計算することを特徴とする請求項１又は２に記載の超音波診断装置。
4. 前記フィルタ行列計算部は、前記主成分分析において抑圧する主成分の数を、前記相関行列の固有値の大きさに基づいて決定することを特徴とする請求項１〜３のいずれか１つに記載の超音波診断装置。
5. 前記相関行列計算部は、超音波送受信が行われた走査範囲を分割した複数の範囲それぞれで、相関行列を計算し、更に、前記複数の範囲それぞれより細分化された複数の細分化範囲それぞれの相関行列を、前記複数の範囲それぞれの相関行列を用いた補間処理により計算し、
   前記フィルタ行列計算部は、前記複数の細分化範囲それぞれの相関行列から、各細分化範囲のフィルタ係数を計算し、
   前記抽出部は、前記複数の細分化範囲それぞれのフィルタ係数を用いて、前記複数の細分化範囲それぞれの前記第２のデータ列を抽出し、
   前記推定部は、前記複数の細分化範囲それぞれの前記第２のデータ列に基づいて、前記移動体情報を推定し、
   前記画像生成部は、前記移動体情報から、超音波画像データを生成する
   ことを特徴とする請求項１〜４のいずれか１つに記載の超音波診断装置。
6. 同一走査線において送信間隔を一定にせずに超音波を送信することで発生した反射波に基づいて生成された反射波データの集合である第１のデータ列を用いて相関行列を計算する相関行列計算部と、
   前記相関行列を用いた主成分分析の結果に基づいて、フィルタ係数を計算するフィルタ行列計算部と、
   前記第１のデータ列に含まれる、前記同一走査線上に存在した移動体で反射された前記超音波の反射波に由来する反射波データの集合である第２のデータ列を、前記フィルタ係数を用いて抽出する抽出部と、
   抽出された前記第２のデータ列に基づいて、前記移動体の移動体情報を推定する推定部と、
   前記移動体情報に基づいて超音波画像データを生成する画像生成部と、
   前記超音波画像データを表示部に表示させる制御部と、
   を備えたことを特徴とする画像処理装置。
7. コンピュータに
   同一走査線において送信間隔を一定にせずに超音波を送信することで発生した反射波に基づいて生成された反射波データの集合である第１のデータ列を用いて相関行列を計算する相関行列計算ステップと、
   前記相関行列を用いた主成分分析の結果に基づいて、フィルタ係数を計算するフィルタ行列計算ステップと、
   前記第１のデータ列に含まれる、前記同一走査線上に存在した移動体で反射された前記超音波の反射波に由来する反射波データの集合である第２のデータ列を、前記フィルタ係数を用いて抽出する抽出ステップと、
   抽出された前記第２のデータ列に基づいて、前記移動体の移動体情報を推定する推定ステップと、
   前記移動体情報に基づいて超音波画像データを生成する画像生成ステップと、
   前記超音波画像データを表示部に表示させる制御ステップと、
   を実行させることを特徴とする画像処理プログラム。

Images