JP7496243B2 - 画像処理装置及び画像処理方法。 - Google Patents

Description

本発明は、画像処理装置及び画像処理方法に関する。

超音波撮像装置やＭＲＩ（Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ｉｍａｇｉｎｇ）装置、Ｘ線ＣＴ（Ｃｏｍｐｕｔｅｄ Ｔｏｍｏｇｒａｐｈｙ）装置に代表される医療用の画像表示装置は、目視できない生体内の情報を数値又は画像の形態で提示する装置として広く利用されている。特に、超音波撮像装置は、他の装置と比較して高い時間分解能を備えており、拍動下の心臓を滲みなく画像化できる性能を有する。

非特許文献１に記載の技術は、微細な血管網を画像情報として表示するための技術であり、特異値分解（Ｓｉｎｇｕｌａｒ Ｖａｌｕｅ Ｄｅｃｏｍｐｏｓｉｔｉｏｎ）に基づく血流映像化技術を開示している。非特許文献１に記載の技術は、複数枚分の画像のエコー信号から２次元行列（Ｃａｓｏｒａｔｉ行列）を生成し、２次元行列に対して特異値分解を施して、空間に関する特異ベクトルと時間に関する特異ベクトルと、を取得する。また、非特許文献１に記載は、特異ベクトルを周波数方向に対するフィルタリングをする処理を開示している。

また、特許文献１は、画素ごとの相関行列を平均化してから固有値分解を適用し、一定範囲の固有ベクトルを排除する方法を開示している。さらに、特許文献１に記載の技術は、処理範囲を小領域に分割し、各小領域にて相関行列の平均を取る処理を開示している。

特開２０１４－１５８６９８号公報

Ｃ. Demene, T. Deffieux, M. Pernot, B. Osmanski, V. Biran, J. Gennisson, L. Sieu, A. Bergel, S. Franqui, J. Correas, I. Cohen, O. Baud, M. Tanter, IEEE Trans. Med. Img., vol. 34, no.11, p.2271, 2015年.

非特許文献１には、雑音成分を抑制するための閾値を設定することについては記載されておらず、当該閾値の必要性についても記載されてない。特許文献１には、雑音成分について記載されていない。つまり、非特許文献１に記載の技術及び特許文献１に記載の技術は、雑音成分が除去された画像を生成することができないおそれがある。そこで、本発明の一態様は、雑音成分が高精度に除去された画像を生成することを目的とする。

上記課題を解決するため、本発明の一態様は以下の構成を採用する。画像処理装置は、プロセッサとメモリとを含み、前記メモリは、対象の第１の複数の画像の時系列を示す行列を保持し、前記プロセッサは、前記行列に対して、画像生成処理を実行して、前記対象の第１の複数の画像を生成し、前記画像生成処理において、前記行列を特異値分解して左行列、特異値行列、及び右行列を取得し、前記右行列を各列が示す周波数に変換し、前記変換した周波数に基づいて、前記右行列に含まれる列群の列それぞれが示す周波数のパワー分布を算出し、前記パワー分布に基づいて、前記列群の列それぞれが示す周波数パワーを示すカーブを算出し、所定の条件に基づいて、前記算出したカーブが示す周波数パワーが一定となる列番号を上側閾値に決定し、前記特異値行列において前記上側閾値以上の行番号の行に位置する特異値を含む特異値群の値を０に変換し、前記左行列、前記変換した特異値行列、及び前記右行列と、に基づいて、前記対象の第１の複数の画像を生成する。

本発明によれば、雑音成分が高精度に除去された画像を生成することができる。

上記した以外の課題、構成及び効果は、以下の実施形態の説明により明らかにされる。

実施例１における画像処理装置のブロック図の一例を示すブロック図である。 実施例１における超音波画像処理の一例を示すフローチャートである。 実施例１における超音波画像処理に含まれる閾値設定処理の一例を示す説明図である。 実施例２における治療効果判定処理の一例を示すフローチャートである。 実施例２における治療効果判定処理の一例を示す説明図である。 実施例２における位置合わせ処理の一例を示す説明図である。 実施例３における移動補正を含む血管画像処理の一例を示すフローチャートである。 実施例３における移動補正を含む血管画像処理の一例を示す説明図である。 実施例４における周波数フィルタ最適化処理の一例を示すフローチャートである。 実施例４における周波数フィルタ最適化処理の一例を示す説明図である。

以下、本発明の実施形態を図面に基づいて詳細に説明する。本実施形態において、同一の構成には原則として同一の符号を付け、繰り返しの説明は省略する。なお、本実施形態は本発明を実現するための一例に過ぎず、本発明の技術的範囲を限定するものではないことに注意すべきである。

本実施形態は、画像処理装置について説明する。微細血管網の映像情報は、腫瘍鑑別及び炎症検査への有効性が高いと期待される情報である。乳腺組織、腹部臓器、及び筋肉等、生体内の様々な組織が検査対象として想定される。検査対象が置かれている環境により体動の影響は異なり、血管径や流速も多様である。特に深部に位置する低速血流が検査対象である場合、エコー信号に含まれる血流成分は微弱であり、当該エコー信号に含まれる雑音成分との区別が難しくなる。

従って、画像処理装置が特異値分解型フィルタにより血流情報を取得するためには、エコー信号から組織成分や体動成分などのクラッタ成分を排除する他、画像処理装置内に備える増幅器等に起因する雑音成分を排除する必要がある。そして、実用性の観点から、画像処理装置は、クラッタ成分および雑音成分を抑制するために抽出する特異ベクトルのランク範囲を指定するための閾値（後述するランク閾値）は、自動的に設定する。特に、検査者が探触子を手動で操作しながら検査対象を探索する、超音波検査に特有の扱いを想定すれば、撮像面の変化に応じてランク閾値が自動的に切り替わる機能の必要性は高い。

本実施形態の画像処理装置は、エコー信号から、自動的に組織成分や体動成分などのクラッタ成分を排除するのみならず、雑音成分も排除することにより、血管情報を高精度に取得することができる。なお、本実施形態では、画像処理装置が、超音波画像に含まれる血管情報を取得する例を主に説明しているが、任意の画像（例えば、ＣＴ（Ｃｏｍｐｕｔｅｄ Ｔｏｍｏｇｒａｐｈｙ）画像やＭＲＩ（Ｍａｇｎｅｔｉｃ Ｒｅｓｏｎａｎｃｅ Ｉｍａｇｉｎｇ）画像）の任意の情報（例えば組織情報やリンパ情報）を取得する場合にも本実施形態の画像処理装置を適用することができる。

図１は、画像処理装置のブロック図の一例を示すブロック図である。画像処理装置１００は、例えば圧電素子によって実現される探触子２００に接続され、探触子２００によって得られた超音波画像に対する画像処理を行う。画像処理装置１００は、例えば、ＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）１１０、メモリ１２０、補助記憶装置１３０、通信装置１４０、入力装置１５０、及び出力装置１６０を有する計算機によって構成される。

ＣＰＵ１１０は、プロセッサを含み、メモリ１２０に格納されたプログラムを実行する。メモリ１２０は、不揮発性の記憶素子であるＲＯＭ（Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）及び揮発性の記憶素子であるＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）を含む。ＲＯＭは、不変のプログラム（例えば、ＢＩＯＳ（Ｂａｓｉｃ Ｉｎｐｕｔ／Ｏｕｔｐｕｔ Ｓｙｓｔｅｍ））などを格納する。ＲＡＭは、ＤＲＡＭ（Ｄｙｎａｍｉｃ Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）のような高速かつ揮発性の記憶素子であり、ＣＰＵ１１０が実行するプログラム及びプログラムの実行時に使用されるデータを一時的に格納する。

補助記憶装置１３０は、例えば、磁気記憶装置（ＨＤＤ（Ｈａｒｄ Ｄｉｓｋ Ｄｒｉｖｅ））、フラッシュメモリ（ＳＳＤ（Ｓｏｌｉｄ Ｓｔａｔｅ Ｄｒｉｖｅ））等の大容量かつ不揮発性の記憶装置であり、ＣＰＵ１１０が実行するプログラム及びプログラムの実行時に使用されるデータを格納する。すなわち、プログラムは、補助記憶装置１３０から読み出されて、メモリ１２０にロードされて、ＣＰＵ１１０によって実行される。

なお、本実施形態において、画像処理装置１００が使用する情報は、データ構造に依存せずどのようなデータ構造で表現されていてもよい。例えば、リスト、テーブル、データベース、又はキューから適切に選択したデータ構造体が、情報を格納することができる。

入力装置１５０は、マウスやキーボード等の、オペレータからの入力を受けつけるため装置である。出力装置１６０は、ディスプレイ装置やプリンタ等の、プログラムの実行結果をオペレータが視認可能な形式で出力する装置である。なお、入力装置１５０や出力装置１６０が画像処理装置１００に含まれずに、画像処理装置１００が有する入力インターフェースに入力装置１５０が接続され、画像処理装置１００が有する出力インターフェースに出力装置１６０が接続されていてもよい。

通信装置１４０は、所定のプロトコルに従って、他の装置との通信を制御するネットワークインターフェース装置である。また、通信装置１４０は、例えば、ＵＳＢ（Ｕｎｉｖｗｓａｌ Ｓｅｒｉａｌ Ｂｕｓ）等のシリアルインターフェースを含む。なお、画像処理装置１００は、通信装置１４０を含まなくてもよい。

ＣＰＵ１１０が実行するプログラムは、リムーバブルメディア（ＣＤ－ＲＯＭ、フラッシュメモリなど）又はネットワークを介して画像処理装置１００に提供され、非一時的記憶媒体である不揮発性の補助記憶装置１３０に格納される。このため、画像処理装置１００は、リムーバブルメディアからデータを読み込むインターフェースを有するとよい。

画像処理装置１００は、物理的に一つの計算機上で、又は、論理的又は物理的に構成された複数の計算機上で構成される計算機システムであり、同一の計算機上で別個のスレッドで動作してもよく、複数の物理的計算機資源上に構築された仮想計算機上で動作してもよい。

ＣＰＵ１１０は、例えば、送受信制御部１１１及び信号処理部１１２を含む。送受信制御部１１１は、送信ビームフォーマ１１３及び受信ビームフォーマ１１４を含み、これらによる処理を制御する。送信ビームフォーマ１１３は、探触子２００が検査対象に向けて送信する超音波信号を所定の時間遅延を素子毎に与えて生成する。受信ビームフォーマ１１４は、当該時間遅延の補正を含めて検査対象からの受信信号を生成し、当該受信信号からＩ／Ｑデータを生成する。

信号処理部１１２は、行列変換部１１５、特異値解析部１１６、特異ベクトル解析部１１７、フィルタ解析部１１８、画像再構成部１１９、及び画質調整部１１１０を含む。行列変換部１１５は、送受信制御部１１１が取得した複数枚の超音波画像のＩ／Ｑデータを行列（例えば行方向が空間を、列方向が時間を示す行列）に変換する。特異値解析部１１６は、行列変換部１１５が生成した行列に対して特異値分解を施し、特異値分解によって得られた特異値を解析し、第１のランク閾値（下側閾値）を決定する。ランク閾値の詳細については後述する。

特異ベクトル解析部１１７は、特異値解析部１１６による特異値分解によって得られた時間的な特異ベクトルを周波数解析し、第２のランク閾値（上側閾値）を決定する。フィルタ解析部１１８は、第１のランク閾値と第２のランク閾値とによって決定される範囲の特異ベクトルを抽出する。画像再構成部１１９は、抽出された特異ベクトルを、超音波画像の取得枚数分のＩ／Ｑデータに再変換する。画質調整部１１１０は、エッジ強調、背景除去、及びスムージング効果等の画像処理を行い、画像処理が行われた画像を出力装置１６０に表示する。

例えば、ＣＰＵ１１０は、メモリ１２０にロードされた送受信制御プログラムに従って動作することで、送受信制御部１１１として機能し、メモリ１２０にロードされた信号処理プログラムに従って動作することで、信号処理部１１２として機能する。ＣＰＵ１１０に含まれる他の機能部（即ち、送受信制御部１１１及び信号処理部１１２それぞれに含まれる機能部）についても、プログラムと機能部の関係は同様である。

なお、ＣＰＵ１１０に含まれる機能部による機能の一部又は全部が、例えば、ＡＳＩＣ（Ａｐｐｌｉｃａｔｉｏｎ Ｓｐｅｃｉｆｉｃ Ｉｎｔｅｇｒａｔｅｄ Ｃｉｒｃｕｉｔ）やＦＰＧＡ（Ｆｉｅｌｄ－Ｐｒｏｇｒａｍｍａｂｌｅ Ｇａｔｅ Ａｒｒａｙ）等のハードウェアによって実現されてもよい。例えば、送受信制御部１１１は、ＦＰＧＡによって実現されることが望ましい。

図２は、超音波画像処理の一例を示すフローチャートであり、図３は、超音波画像処理に含まれる閾値設定処理（後述するＳ２０１～Ｓ２０７）の一例を示す説明図である。送受信制御部１１１は、探触子２００から検査対象のＮフレーム分の画像（必ずしも連続するフレームの画像である必要はないが複数の異なる時刻に撮影された画像）それぞれのＩ／Ｑデータを取得し、メモリ１２０及び／又は補助記憶装置１３０に格納する（Ｓ２０１）。各画像の空間サイズは（Ｚ，Ｘ）である、つまり送受信制御部１１１は、検査対象のＺ×Ｘ画素のＮフレーム分の画像を取得する。

行列変換部１１５は、空間サイズ（Ｚ，Ｘ）のＮフレーム分の画像を、Ｚ行Ｘ列に変換する（Ｓ２０２）。具体的には、例えば、行列変換部１１５は、Ｎフレーム分のＺ×Ｘ行列それぞれについて、全ての要素を縦に並べた列ベクトルを生成し、当該生成した列ベクトルをフレームの順序で列方向に並べた（例えば時間が古いフレームほど左側に位置するように並べた）ＺＸ行Ｎ列のＣａｓｏｒａｔｉ行列を生成する。つまり、このＺＸ×Ｎ行列は、検査対象のＮフレームの画像の時系列を示し、このＺＸ×Ｎ行列を構成する各列ベクトルは、画像情報の１フレーム分のＩ／Ｑデータを示す。

特異値解析部１１６は、このＺＸ×Ｎ行列を特異値分解することで、ＶΛＷを得る（Ｓ２０３）。ここで、Ｖは空間に関する情報を示すＺＸ行ＺＸ列の左特異ベクトル行列であり、Ｗは時間に関する情報を示すＮ行Ｎ列の右特異ベクトル行列である。

ΛはＺＸ行Ｎ列の特異値行列である。なお、特異値行列は、各特異値を対角成分とする対角行列が行方向（下側）及び／又は列方向（右側）に拡張されたＺＸ行Ｎ列の行列であり、当該拡張された部分の成分は全て０である。特異値行列の対角成分をλ_１，・・・，λ_ｑ（ｑ＝ｍｉｎ｛ＺＸ，Ｎ｝）としたとき、λ_１≧λ_２≧・・・≧λ_ｑ－１≧λ_ｑ≧０であるとする。また、特異値行列のランクは（０でない）特異値の個数と等しい。以下、説明の便宜のために、λ_ｋ（１≦ｋ≦ｑ）をランクｋの特異値とも呼び、行列の積であるＶΛＷを算出する際にλ_ｋ（１≦ｋ≦ｑ）を含む成分に掛けられる行をランクｋの行とも呼び、ＶΛＷを算出する際にλ_ｋ（１≦ｋ≦ｑ）を含む成分に掛ける列をランクｋの列とも呼ぶ。

特異ベクトル解析部１１７は、特異値行列Λの各対角成分を、特異値行列Λの対角成分の総和で除算した商（σ_ｉ＝λ_ｉ／Σλ_ｉ（Σの範囲はｉ＝１からｉ＝ｑ））を算出する（Ｓ２０４）ことで、各ランクの成分が全てのＩ／Ｑデータに占める情報の割合を算出することができる。

一般的に、生体からのエコー信号のうち、クラッタ成分と血流成分とには、約４０ｄＢ（１００対１）程度の強度差がある。従って、したがって、特異ベクトル解析部１１７は、情報量（前述した商）が０．０１（予め定められた値）に最も近い値を示した特異値に対応するランクを、クラッタ成分と血流成分とを判別する第１のランク閾値ｋ_ｃｒ（単に第１閾値とも呼ぶ）として設定する（Ｓ２０５）。なお、送受信制御部１１１が取得する画像の種類や分析の目的等によって、前述した０．０１という値を、画像処理装置１００のユーザが変更可能であることが望ましい。

なお、特異ベクトル解析部１１７は、図３に示す、情報量とランクとの関係を示すグラフ、並びに当該グラフ上に示された血流閾値（上述した予め定められた値）及び第１のランク閾値を、出力装置１６０に出力してもよい。

なお、特異ベクトル解析部１１７は、情報量（前述した商）が０．０１（予め定められた値）に最も近い値を示した特異値に対応する複数のランクがあると判定した場合には、例えば、当該複数のランクのうち特異値行列において最も左上に位置する特異値に対応するランクを第１のランク閾値として設定してもよいし、最も右下に位置する特異値に対応するランクを第１のランク閾値として設定してもよいし、これらの中央に位置する特異値に対応するランクを第１のランク閾値として設定してもよい。

続いて、特異ベクトル解析部１１７は、右特異ベクトル行列Ｗを行方向にフーリエ変換することで右特異ベクトル行列Ｗを周波数変換する（Ｓ２０６）。特異ベクトル解析部１１７は、フーリエ変換後の右特異ベクトル行列Ｗに対して周波数解析を行って、第２のランク閾値ｋ_ｎｒ（単に第２閾値とも呼ぶ）を算出する（Ｓ２０７）。ステップＳ２０７の具体例について説明する。

高周波域かつ特異値が低い高次側（即ち特異値行列においてより右下に位置する対角成分に対応するランク）のランク範囲（つまり当該ランク範囲の特異値に掛ける右特異ベクトル行列Ｗの列）が除去すべき雑音成分の支配域である。従って、例えば、特異ベクトル解析部１１７は、高周波域（例えば１ｋＨ（予め定められた値）以上の範囲、全範囲の周波数について算出してもよい）のランクごとの（即ち右特異ベクトル行列Ｗの各列に対応する）周波数のパワー分布を算出し、算出したパワー分布から、ランクごとの周波数パワーの平均値を計算することで、エネルギーカーブを算出する。

エネルギーカーブが示す周波数パワーは、クラッタや血流の成分を含むランク範囲より高いランクにおいてほぼ一定、つまりクラッタや血流の成分を含むランク範囲より高いランクにおいては雑音成分が支配的である。つまりエネルギーカーブが示す周波数パワーが一定となる起点付近が、血流成分と雑音成分とを判別する第２のランク閾値として適当である。

具体的には、特異ベクトル解析部１１７は、例えば、エネルギーカーブの微分係数が０となる始点のランクを第２のランク閾値に決定する。また、特異ベクトル解析部１１７は、例えば、エネルギーカーブの始点と終点を結ぶ直線が示す周波数パワーと、エネルギーカーブが示す周波数パワーと、の差分が最大となるランクを第２のランク閾値に決定してもよい。また、特異ベクトル解析部１１７は、例えば、エネルギーカーブの変曲点を通過する直線グラフ（例えばエネルギーカーブに対するフィッテイング処理を用いて算出される）と、最終ランク（即ち０でない特異値のうち最も小さい特異値に対応するランク）を通過する水平の直線グラフと、の交点のような、エネルギーカーブに対する解析結果から得られる起点を推定し、当該起点が示すランクによって第２のランク閾値を設定する。

なお、特異ベクトル解析部１１７は、第２のランク閾値として複数の値の候補があると判定した場合には、例えば、当該複数のランクのうち特異値行列において最も左上に位置する特異値に対応するランクを第２のランク閾値として設定してもよいし、最も右下に位置する特異値に対応するランクを第２のランク閾値として設定してもよいし、これらの中央に位置する特異値に対応するランクを第２のランク閾値として設定してもよい。

また、特異ベクトル解析部１１７は、第１のランク閾値＋１≧第２のランク閾値、であると判定した場合、例えば、第１のランク閾値及び／又は第２のランク閾値を再設定するとよい。具体的には、例えば、特異ベクトル解析部１１７は、第１のランク閾値に２を加えた値以上のランクから、ステップＳ２０７と同様の方法で第２のランク閾値を再設定する。また、特異ベクトル解析部１１７は、単に第１のランク閾値又は第２のランク閾値の一方を固定して、第２のランク閾値－第１のランク閾値＝２以上の所定値（例えば２）となるように、他方を再設定するようにしてもよい。

なお、特異ベクトル解析部１１７は、図３に示すパワー分布、高周波域の平均パワー推移（第２閾値を含む）、及び差分グラフ（第２閾値を含む）を出力装置１６０に出力してもよい。

続いて、フィルタ解析部１１８は、特異値行列Λの、第１のランク閾値より小さい、又は第２のランク閾値よりも大きいランクに対応する特異値を全て０に変換した行列である血流抽出フィルタΛ^（ｋ）を生成する（Ｓ２０８）。フィルタ解析部１１８は、ＺＸ×Ｎ行列であるＳ＝Ｖ（Λ^（ｋ））Ｗ^Ｈ（但し、Ｗ^ＨはＷの共役転置行列）を算出する（Ｓ２０９）。画像再構成部１１９は、ＺＸ×Ｎ行列のＳを、Ｎ個のＺ×Ｘ行列へと分解する（ステップＳ２０２で実行された変換の逆変換をする）ことで、Ｚ×Ｘ画素のＮフレーム分の２次元画像を再構成する（Ｓ２１０）。

これにより、フィルタ解析部１１８は、第１のランク閾値と第２のランク閾値との間で規定されるランク範囲によって血流成分を抽出し、元のＩ／Ｑデータと同じサイズのデータを生成することができる。

画質調整部１１１０は、ステップＳ２１１～ステップＳ２１４において、血管画像としての視認性を向上する画像フィルタを適用する。なお、ステップＳ２１１～ステップＳ２１４に示す処理はあくまで画像フィルタ適用の例であり、ステップＳ２１１～ステップＳ２１４の少なくとも一部の処理が行われなくてもよいし、さらに追加の画質調整処理が行われてもよい。

画質調整部１１１０は、再構成されたＺ×Ｘ画素のＮフレーム分の画像を時間軸方向に加算平均する（Ｓ２１１）ことで、画像の背景雑音を抑制することができる。さらに、画質調整部１１１０は、再構成されたＺ×Ｘ画素のＮフレーム分の画像の背景除去を行う（Ｓ２１２）ことで、血管領域の信号をより強調することができる。画像全体の平均値を各画素の輝度から一律に差分する方法は、背景除去の一例である。

画質調整部１１１０は、再構成されたＺ×Ｘ画素のＮフレーム分の画像に対し、移動平均とメディアンフィルタを適用する（Ｓ２１３）。画質調整部１１１０は、さらに構造物の輪郭を強調する画像処理を適用する。例えば、画質調整部１１１０は、特定方向の輪郭を強調する機能を備えるモルフォロジフィルタを、入力画像に対して複数の角度について適用し、モルフォロジフィルタ適用後の画像に対して加算平均を行う等の多重解像度処理を行う（Ｓ２１４）。これにより、画質調整部１１１０は、画像内の血管部の繋がりを改善し、輪郭が明瞭化された画像を得ることができる。

続いて、画質調整部１１１０は、画質調整後の画像を出力装置１６０に出力する（Ｓ２１５）。以上の超音波画像処理により、術者は探触子２００を検査対象に向けて設置するだけで、明瞭な血管網の情報を出力装置１６０に表示される映像情報として取得することができる。

なお、本実施例の画像処理装置１００は、超音波画像から血流成分を得るために第１のランク閾値と第２のランク閾値とを算出したが、組織成分を主とするクラッタ成分のみを得たいのであれば第１のランク閾値のみを算出して第１のランク閾値以上のランクに対応する特異値を０に変換してもよいし、雑音のみを除去したいのであれば第２のランク閾値のみを算出して第２のランク閾値以上のランクに対応する特異値を０に変換してもよい。

本実施例では、治療効果判定処理について説明する。実施例１と同様の構成や処理については説明を一部省略する。図４は、治療効果判定処理の一例を示すフローチャートであり、図５は、治療効果判定処理の一例を示す説明図である。図６は、後述するステップＳ４０７の位置合わせ処理の一例を示す説明図である。

例えば癌の化学療法に用いられる薬剤は、抗血管新生薬などの血管を縮退させる作用を持つものが少なくない。癌に対する放射線治療において効果がある場合も、血管網の縮退が認められる。図４乃至図６に示す治療効果判定処理は、治療前後の血管画像に対する画像処理を行うことで、治療効果判定の定量化する。

送受信制御部１１１は、ステップＳ２０１と同様の処理により、治療前の検査対象（血管を含む対象）の複数の画像と、治療後の当該検査対象の複数の画像と、を収集する（Ｓ４０１、Ｓ４０２）。行列変換部１１５は、ステップＳ２０２と同様の処理によって、治療前の検査対象の複数の画像に対して行列変換を施し、治療後の検査対象の複数の画像に対しても行列変換を施し、特異値解析部１１６は、行列変換部１１５が生成した行列それぞれを特異値分解する（Ｓ４０３）。

特異ベクトル解析部１１７は、特異値分解された行列それぞれの特異値行列及び右特異ベクトル行列に対して、ステップＳ２０４～ステップＳ２０７と同様の処理によって、第１のランク閾値及び第２のランク閾値を決定する（Ｓ４０４）。フィルタ解析部１１８、画像再構成部１１９、及び画質調整部１１１０は、ステップＳ２０８～ステップＳ２１４の処理によって、治療前の検査対象及び治療後の検査対象の血管が明瞭化された画像を得る（Ｓ４０５）。これにより、図５に示すように、治療前と治療後とそれぞれについて、血管５０３（第１の対象の一例）を含む腫瘍領域５０２（第２の対象の一例）の画像データ５０１が得られる。

画質調整部１１１０は、治療前と治療後それぞれの画像を小領域に分割し（例えば5ｍｍ四方）、各小領域について輝値を加算し血管値として算出し、各小領域の面積に対する血管値を当該小領域の血管密度として算出することで、血管密度画像を生成する（Ｓ４０６）。画質調整部１１１０は、このようにして血管密度を算出することにより、僅かな治療効果が反映される血管密度画像を算出することができる。

なお、画質調整部１１１０は、血管密度を正確に試算するために、一定値以上の輝度の画素を血管領域と定義し、血管領域に該当する画素の総数を血管値に決定することが望ましい。画質調整部１１１０は、このように血管密度を算出することにより、流量や流速の影響を含みにくく、血管の有無だけを判定しやすい血管密度画像を生成することができる。

続いて、画質調整部１１１０は、治療前の血管密度画像と治療後の血管密度画像と、を血管密度が高い領域同士で位置合わせする（Ｓ４０７）。例えば、画質調整部１１１０は、治療前後の血管密度分布画像を重畳させる。具体的には、例えば、画質調整部１１１０は、最も簡便な方法として、血管密度が最大となる画素を位置整合の基準として治療前後の血管密度分布画像を重畳させる。

また、画質調整部１１１０は、図６に示すように、治療前後の血管密度分布画像それぞれに対して重心を定義し、当該重心を基準に位置整合を行なってもよい。一般的に重心は質点系における基準点からの距離と質量から算出されるが、ここでは、画質調整部１１１０は、質量を血管密度に置き換えて重心として定義する。

具体的には、画質調整部１１１０は、治療前の血管密度分布画像の小領域Ｐｉそれぞれについて、ｒ_ｉ＝ＯＰ_ｉ（ＯＰ_ｉはＰｉの位置ベクトル）としたとき、治療前の血管密度分布画像の重心ベクトルｒ_Ｇ＝Σｍ_ｉｒ_ｉ／Σｍ_ｉ（ｍ_ｉはＰ_ｉの血管密度、Σの範囲はいずれもｉ＝１からｉ＝小領域の個数）を算出して、重心Ｐ_Ｐｒｅを決定する。画質調整部１１１０は、同様の方法で、治療後の血管密度分布画像の重心Ｐ_Ｐｏｓｔを決定する。画質調整部１１１０は、治療前の血管密度画像と治療後の血管密度画像との重心の位置が一致するように、これら２つの血管密度画像の位置合わせをする。画質調整部１１１０は、このように位置合わせを行うことで、血管密度画像において局所的に血管密度が高い領域があっても、その影響を軽減した正確な位置合わせを行うことができる。

続いて、画質調整部１１１０は、位置合わせをした治療前と治療後の血管密度画像を用いて、効果判定指標を算出し（Ｓ４０８）、位置合わせをした血管密度画像と効果判定指標とを出力装置１６０に出力する（Ｓ４０９）。

具体的には、例えば、画質調整部１１１０は、血管密度が一定値（予め定められた値又は画像処理装置１００のユーザが設定可能な値）以上の範囲を定義して、ユーザの関心領域を可視化することができる。例えば、画質調整部１１１０は、治療前後の関心領域の面積を算出し、関心領域の面積変化率を効果判定指標として算出する。

また、画質調整部１１１０は、関心領域の血管密度変化率を効果判定指標として算出してもよい。また、画質調整部１１１０は、腫瘍領域を球体と仮定して、関心領域の体積変化率を効果判定指標として算出してもよい。なお、画質調整部１１１０は、関心領域ではなく腫瘍領域５０２全体の面積変化率、血管密度変化率、及び体積変化率を効果判定指標として算出してもよい。

本実施例の画像処理装置１００により、各種治療の効果判定の定量化が実現し、ひいては、画像処理装置１００のユーザは、客観的情報に基づく治療最適化（治療方法や薬剤選定の最適化）を迅速に行うことができる。

本実施例では、移動補正を含む血管画像処理について説明する。実施例１と同様の構成や処理については説明を一部省略する。図７は、移動補正を含む血管画像処理の一例を示すフローチャートであり、図８は、移動補正を含む血管画像処理の一例を示す説明図である。

送受信制御部１１１は、ステップＳ２０１と同様の処理により、検査対象（血管を含む対象）の複数の画像を収集する（Ｓ７０１）。行列変換部１１５は、ステップＳ２０２と同様の処理によって、検査対象の複数の画像に対して行列変換を施し、特異値解析部１１６は、行列変換部１１５が生成した行列を特異値分解する（Ｓ７０２）。特異ベクトル解析部１１７は、特異値分解された行列の特異値行列及び右特異ベクトル行列に対して、ステップＳ２０４～ステップＳ２０７と同様の処理によって、第１のランク閾値及び第２のランク閾値を決定する（Ｓ７０３）。

フィルタ解析部１１８は、特異値行列Λの、第１のランク閾値以上のランクに対応する特異値を全て０に変換した行列である組織抽出フィルタΛ^（ｊ）を生成し、ＺＸ×Ｎ行列であるＳ’＝Ｖ（Λ^（ｊ））Ｗ^Ｈ（但し、Ｗ^ＨはＷの共役転置行列）を算出する（Ｓ７０４）。組織抽出フィルタにより、血流成分を含まないクラッタ成分を抽出することができ、組織の移動を推定することができる。画像再構成部１１９は、ＺＸ×Ｎ行列のＳ’を、Ｎ個のＺ×Ｘ行列へと分解する（ステップＳ２０２で実行された変換の逆変換をする）ことで、Ｚ×Ｘ画素のＮフレーム分の２次元組織画像を再構成する（Ｓ７０５）。

画質調整部１１１０は、再構成したＮフレーム分の２次元組織画像に対して、フレーム間に生じるＺ方向及びＸ方向の移動量を計測する等の移動解析を行う（Ｓ７０６）。具体的には、例えば、画質調整部１１１０は、一般的に知られる有効な手法であるトラッキングアルゴリズム（特に例えば、正規化相互相関演算）によって移動量を計測する。

フィルタ解析部１１８は、ステップＳ２０８と同様の処理によって、～ステップＳ２１４の処理と同様の処理によって血管画像を明瞭化して抽出するための血流抽出フィルタを生成し、ステップＳ２０９と同様の処理によってＺＸ×Ｎ行列であるＳを算出する（Ｓ７０７）。画像再構成部１１９は、ステップ２１０と同様の処理によって、Ｚ×Ｘ画素のＮフレーム分の２次元血管画像を再構成する（Ｓ７０８）。

画質調整部１１１０は、各フレームの２次元血管画像に対して、対応するフレームのステップＳ７０６において計測した２次元組織画像の移動量に応じた移動補正（例えば移動量をそのまま差し引く）を行う（Ｓ７０９）。画質調整部１１１０は、移動補正後のＮフレームの２次元血管画像に対して加算平均を行う等（加算平均に加えてステップＳ２１２～ステップＳ２１４の処理を行ってもよい）の画像解析処理行う（Ｓ７１０）。

なお、画質調整部１１１０は、ステップＳ７０９における加算平均処理に代えて又は加えて、移動補正後の２次元血管画像の画素ごとの輝度値の時間変化を算出し、最大の輝度値を選択して、最大輝度値保持画像を取得してもよいし、移動補正後の２次元血管画像から血管内の血流の流量の変化値を示す画像を生成してもよいし、移動補正後の２次元血管画像から血流の流速を示す画像を取得してもよい。

画質調整部１１１０は、ステップＳ７１０で生成した画像を出力装置１６０に出力する（Ｓ７１１）。なお、画質調整部１１１０は、ステップＳ７１０で得られた解析結果をさらに出力装置１６０に出力してもよい。本実施例の画像処理装置１００は、組織画像から組織の移動量を推定し、推定した移動量を血管画像に適用することで、フレーム間で大きな体動があった場合に平均加算等を行っても血管画像に滲みが発生しづらくなる上に、他の血流や血管の変化値画像についても正確な画像を生成することができる。なお、画質調整部１１１０は、血管画像から血管の移動量を推定し、推定した移動量を組織画像に適用してもよい。

本実施例の画像処理装置１００は、特異値分解によって得られるクラッタ成分の特性を分析し、分析結果を周波数フィルタの設計に適応的に反映させる周波数フィルタ最適化処理を行う。実施例１と同様の構成や処理については説明を一部省略する。図９は、周波数フィルタ最適化処理の一例を示すフローチャートであり、図１０は、周波数フィルタ最適化処理の一例を示す説明図である。

送受信制御部１１１は、ステップＳ２０１と同様の処理により、検査対象（血管を含む対象）の複数の画像を収集する（Ｓ９０１）。行列変換部１１５は、ステップＳ２０２と同様の処理によって、検査対象の複数の画像に対して行列変換を施し、特異値解析部１１６は、行列変換部１１５が生成した行列を特異値分解する（Ｓ９０２）。

特異ベクトル解析部１１７は、ステップＳ２０４及びステップＳ２０５と同様の処理によって、第１のランク閾値を決定する（Ｓ９０３）。特異ベクトル解析部１１７は、ステップＳ２０６及びステップＳ２０７で説明した処理と同様の処理によって、全範囲の周波数のパワー分布を算出し、算出したパワー分布からランクごとの周波数パワーの平均値を計算してエネルギーカーブを算出する（Ｓ９０４）。

特異ベクトル解析部１１７は、図１０に示す、情報量とランクとの関係を示すグラフ、当該グラフ上に示された血流閾値（上述した予め定められた値）及び第１のランク閾値、平均勾配、並びにパワー分布を出力装置１６０に出力してもよい。

特異ベクトル解析部１１７は、第１のランク閾値以下のランクに対応する周波数の範囲である組織周波数（組織成分を主とするクラッタ成分の周波数）解析を行う（Ｓ９０５）。具体的には、例えば、特異ベクトル解析部１１７は、組織周波数における周波数ごとの平均パワーを算出する。例えば、第１のランク閾値以下のランクがランク１からランク５、組織周波数範囲が±３０００Ｈｚである場合、特異ベクトル解析部１１７は、５つのデータ平均によって得られる０～３０００Ｈｚの平均周波数パワーを算出結果として得る。特異ベクトル解析部１１７は、この平均周波数パワーを微分することにより、平衡周波数ｂ（平均周波数パワーが一定、つまり微分曲線が０に近接する（又は０になる）周波数）や、ピーク周波数ａ（平均周波数の最大勾配、つまり微分曲線の極値となる周波数）、を算出することができる。

また、例えば、特異ベクトル解析部１１７は、最小のランクであるランク１の周波数である０Ｈｚを起点に、各ランクの最大パワーを選択して得られる周波数ピークグラフを算出し、当該ピークグラフの勾配（角度ｃ）を算出する。また、例えば、特異ベクトル解析部１１７は、組織周波数の範囲において、周波数毎の平均値を算出する。

なお、特異ベクトル解析部１１７は、図１０に示すクラッタ範囲の平均周波数パワー（ピーク周波数ａ及び平衡周波数ｂを含む）及び組織領域の周波数ピークグラフ（角度ｃを含む）を出力装置１６０に出力してもよい。

続いて、特異ベクトル解析部１１７は、ステップＳ９０６で算出した情報に基づく解析を行う（Ｓ９０６）。具体的には、例えば、特異ベクトル解析部１１７は、クラッタのパワー、内在する周波数レンジ、及び速度を推定する。ピーク周波数ａや平衡周波数ｂは、クラッタ範囲に含まれる周波数範囲と強さを示す指標であり、角度ｃは、この勾配が大きいほど、クラッタの速度は高くなる傾向があるため、ｃ＝Δｆ／ΔＲ(ランク範囲)が速度関連因子であるため、特異ベクトル解析部１１７は、これらの指標や関連因子を用いてクラッタ解析を行うことができる。なお、特異ベクトル解析部１１７は、クラッタ解析結果を出力装置１６０に出力してもよい。

フィルタ解析部１１８は、ステップＳ９０６のクラッタ解析によって得られたクラッタ成分の推定結果から、ＩＩＲ（Ｉｎｆｉｎｉｔｅ Ｉｍｐｌｕｓｅ Ｒｅｓｐｏｎｓｅ）型フィルタの種類やカットオフ周波数等のパラメータを決定することで、適応型ＩＩＲフィルタを設計する（Ｓ９０７）。

フィルタ解析部１１８は、元のＩ／Ｑデータに自己相関演算を適用して流速解析の結果を取得し（Ｓ９０８）、画像再構成部１１９は、流速解析の結果に対して適応型ＩＩＲフィルタを適用することにより、クラッタ成分を抑制した血流画像（及び流速画像）を生成する（Ｓ９０９）。なお、画質調整部１１１０が、例えば、ステップＳ２１１～ステップＳ２１４等の画質調整処理を、ステップＳ９０８で得られた血流画像（及び流速画像）に対して行ってもよい。画質調整部１１１０は、クラッタ成分を抑制した血流画像（及び流速画像）を出力装置１６０に出力する（Ｓ９１１）。

本実施例の画像処理装置１００は、フレームの取り込みと同時にＩＩＲフィルタも再設計できるため、検査者が探触子２００を操作していても、常に適切な血流画像を表示することが可能であり、ひいては高精度な超音波診断が実現する。

なお、本発明は上記した実施例に限定されるものではなく、様々な変形例が含まれる。例えば、上記した実施例は本発明を分かりやすく説明するために詳細に説明したものであり、必ずしも説明した全ての構成を備えるものに限定されるものではない。また、ある実施例の構成の一部を他の実施例の構成に置き換えることも可能であり、また、ある実施例の構成に他の実施例の構成を加えることも可能である。また、各実施例の構成の一部について、他の構成の追加・削除・置換をすることが可能である。

また、上記の各構成、機能、処理部、処理手段等は、それらの一部又は全部を、例えば集積回路で設計する等によりハードウェアで実現してもよい。また、上記の各構成、機能等は、プロセッサがそれぞれの機能を実現するプログラムを解釈し、実行することによりソフトウェアで実現してもよい。各機能を実現するプログラム、テーブル、ファイル等の情報は、メモリや、ハードディスク、ＳＳＤ（Ｓｏｌｉｄ Ｓｔａｔｅ Ｄｒｉｖｅ）等の記録装置、または、ＩＣカード、ＳＤカード、ＤＶＤ等の記録媒体に置くことができる。

また、制御線や情報線は説明上必要と考えられるものを示しており、製品上必ずしも全ての制御線や情報線を示しているとは限らない。実際には殆ど全ての構成が相互に接続されていると考えてもよい。

１００ 画像処理装置、１１０ ＣＰＵ、１１１ 送受信制御部、１１２ 信号処理部、１１３ 送信ビームフォーマ、１１４ 受信ビームフォーマ、１１５ 行列変換部、１１６ 特異値解析部、１１７ 特異ベクトル解析部、１１８ フィルタ解析部、１１９ 画像再構成部、１２０ メモリ、１３０ 補助記憶装置、１４０ 通信装置、１５０ 入力装置、１６０ 出力装置、２００ 探触子、１１１０ 画質調整部

Claims (14)

1. 画像処理装置であって、
   プロセッサとメモリとを含み、
   前記メモリは、対象の第１の複数の画像の時系列を示す行列を保持し、
   前記プロセッサは、
   前記行列に対して、画像生成処理を実行して、前記対象の第１の複数の画像を生成し、
   前記画像生成処理において、
   前記行列を特異値分解して左行列、特異値行列、及び右行列を取得し、
   前記右行列を行方向にフーリエ変換することで周波数に変換し、
   前記変換した周波数に基づいて、前記右行列の各列に対応する周波数のパワー分布を算出し、
   前記パワー分布に基づいて、前記右行列の列ごとの周波数パワーを示すカーブを算出し、
   所定の条件に基づいて、前記算出したカーブが示す周波数パワーが一定となる列番号を上側閾値に決定し、
   前記特異値行列において前記上側閾値以上の行番号の行に位置する特異値を含む特異値群の値を０に変換し、
   前記左行列、前記変換した特異値行列、及び前記右行列と、に基づいて、前記対象の第１の複数の画像を生成する、画像処理装置。
2. 請求項１に記載の画像処理装置であって、
   前記プロセッサは、
   前記画像生成処理において、
   前記特異値行列の各特異値が示す情報量を算出し、
   前記算出した情報量が所定値と最も近い特異値の前記特異値行列における行番号を下側閾値に決定し、
   前記上側閾値が前記下側閾値に１を加えた値以下であると判定した場合には、前記上側閾値と前記下側閾値との少なくとも一方を再設定し、
   前記特異値行列において前記上側閾値以上の行番号の行に含まれる特異値と、前記下側閾値以下の行番号の行に含まれる特異値と、を０に変換する、画像処理装置。
3. 請求項１又は２に記載の画像処理装置であって、
   前記第１の複数の画像それぞれは、血管を含む組織が撮像された超音波画像であり、
   前記左行列は、前記第１の複数の画像の空間に関する情報を示し、
   前記右行列は、前記時系列における時間に関する情報を示す、画像処理装置。
4. 請求項２に記載の画像処理装置であって、
   前記対象は、第１の対象と第２の対象とを含み、
   前記メモリは、前記対象の第２の複数の画像の時系列を示す行列を保持し、
   前記プロセッサは、
   前記第２の複数の画像の時系列を示す行列に対して、前記画像生成処理を実行して、前記対象の第２の複数の画像を生成し、
   前記生成した第１の複数の画像から前記対象の第１の画像を生成し、
   前記生成した第２の複数の画像から前記対象の第２の画像を生成し、
   所定の条件に基づいて、前記生成した第１の画像に含まれる前記第１の対象の高密度領域と、前記生成した第２の画像に含まれる前記第２の対象の高密度領域と、を特定し、
   前記生成した第１の画像に含まれる前記第１の対象の高密度領域と、前記第２の画像に含まれる前記生成した第２の対象の高密度領域と、で前記生成した第１の画像と、前記生成した第２の画像と、を位置合わせして、前記対象の第３の画像を生成する、画像処理装置。
5. 請求項２に記載の画像処理装置であって、
   前記対象は、第１の対象を含み、
   前記プロセッサは、
   前記左行列と、前記下側閾値以上の行番号の行に位置する特異値を０に変換した前記特異値行列と、及び前記右行列と、に基づいて、前記対象の第２の複数の画像を生成し、
   前記生成した複数の画像に含まれる前記対象の移動解析をし、
   前記移動解析の結果に基づいて、前記生成した第１の複数の画像それぞれに含まれる前記第１の対象を移動補正し、
   前記移動補正した第１の対象を含む画像を生成する、画像処理装置。
6. 請求項２に記載の画像処理装置であって、
   前記メモリは、前記第１の複数の画像を保持し、
   前記プロセッサは、
   前記カーブに基づいて、前記右行列の前記下側閾値以下の列番号の列それぞれが示す周波数の範囲における、平衡周波数及びピーク周波数の少なくとも一方を含む指標を算出し、
   前記指標に基づいて、ＩＩＲフィルタを生成し、
   前記第１の複数の画像に対して、流速解析を実行し、
   前記流速解析の結果に前記ＩＩＲフィルタを適用して、前記第１の複数の画像を生成する、画像処理装置。
7. 請求項１に記載の画像処理装置であって、
   出力装置を含み、
   前記プロセッサは、前記パワー分布又は前記カーブの少なくとも一方を、前記出力装置に出力する、画像処理装置。
8. 画像処理装置による画像処理方法であって、
   前記画像処理装置は、プロセッサとメモリとを含み、
   前記メモリは、対象の第１の複数の画像の時系列を示す行列を保持し、
   前記画像処理方法は、
   前記プロセッサが、前記行列に対して、画像生成処理を実行して、前記対象の第１の複数の画像を生成し、
   前記画像生成処理において、
   前記プロセッサが、前記行列を特異値分解して左行列、特異値行列、及び右行列を取得し、
   前記プロセッサが、前記右行列を行方向にフーリエ変換することで周波数に変換し、
   前記プロセッサが、前記変換した周波数に基づいて、前記右行列の各列に対応する周波数のパワー分布を算出し、
   前記プロセッサが、前記パワー分布に基づいて、前記右行列の列ごとの周波数パワーを示すカーブを算出し、
   前記プロセッサが、所定の条件に基づいて、前記算出したカーブが示す周波数パワーが一定となる列番号を上側閾値に決定し、
   前記プロセッサが、前記特異値行列において前記上側閾値以上の行番号の行に位置する特異値を含む特異値群の値を０に変換し、
   前記プロセッサが、前記左行列、前記変換した特異値行列、及び前記右行列と、に基づいて、前記対象の第１の複数の画像を生成する、画像処理方法。
9. 請求項８に記載の画像処理方法であって、
   前記画像生成処理において、
   前記プロセッサが、前記特異値行列の各特異値が示す情報量を算出し、
   前記プロセッサが、前記算出した情報量が所定値と最も近い特異値の前記特異値行列における行番号を下側閾値に決定し、
   前記上側閾値が前記下側閾値に１を加えた値以下であると判定した場合には、前記上側閾値と前記下側閾値との少なくとも一方を再設定し、
   前記プロセッサが、前記特異値行列において前記上側閾値以上の行番号の行に含まれる特異値と、前記下側閾値以下の行番号の行に含まれる特異値と、を０に変換する、画像処理方法。
10. 請求項８又は９に記載の画像処理方法であって、
    前記第１の複数の画像それぞれは、血管を含む組織が撮像された超音波画像であり、
    前記左行列は、前記第１の複数の画像の空間に関する情報を示し、
    前記右行列は、前記時系列における時間に関する情報を示す、画像処理方法。
11. 請求項９に記載の画像処理方法であって、
    前記対象は、第１の対象と第２の対象とを含み、
    前記メモリは、前記対象の第２の複数の画像の時系列を示す行列を保持し、
    前記画像処理方法は、
    前記プロセッサが、前記第２の複数の画像の時系列を示す行列に対して、前記画像生成処理を実行して、前記対象の第２の複数の画像を生成し、
    前記プロセッサが、前記生成した第１の複数の画像から前記対象の第１の画像を生成し、
    前記プロセッサが、前記生成した第２の複数の画像から前記対象の第２の画像を生成し、
    前記プロセッサが、所定の条件に基づいて、前記生成した第１の画像に含まれる前記第１の対象の高密度領域と、前記生成した第２の画像に含まれる前記第２の対象の高密度領域と、を特定し、
    前記プロセッサが、前記生成した第１の画像に含まれる前記第１の対象の高密度領域と、前記第２の画像に含まれる前記生成した第２の対象の高密度領域と、で前記生成した第１の画像と、前記生成した第２の画像と、を位置合わせして、前記対象の第３の画像を生成する、画像処理方法。
12. 請求項９に記載の画像処理方法であって、
    前記対象は、第１の対象を含み、
    前記画像処理方法は、
    前記プロセッサが、前記左行列と、前記下側閾値以上の行番号の行に位置する特異値を０に変換した前記特異値行列と、及び前記右行列と、に基づいて、前記対象の第２の複数の画像を生成し、
    前記プロセッサが、前記生成した複数の画像に含まれる前記対象の移動解析をし、
    前記プロセッサが、前記移動解析の結果に基づいて、前記生成した第１の複数の画像それぞれに含まれる前記第１の対象を移動補正し、
    前記プロセッサが、前記移動補正した第１の対象を含む画像を生成する、画像処理方法。
13. 請求項９に記載の画像処理方法であって、
    前記メモリは、前記第１の複数の画像を保持し、
    前記画像処理方法は、
    前記プロセッサが、前記カーブに基づいて、前記右行列の前記下側閾値以下の列番号の列それぞれが示す周波数の範囲における、平衡周波数及びピーク周波数の少なくとも一方を含む指標を算出し、
    前記プロセッサが、前記指標に基づいて、ＩＩＲフィルタを生成し、
    前記プロセッサが、前記第１の複数の画像に対して、流速解析を実行し、
    前記プロセッサが、前記流速解析の結果に前記ＩＩＲフィルタを適用して、前記第１の複数の画像を生成する、画像処理方法。
14. 請求項８に記載の画像処理方法であって、
    前記画像処理装置は、出力装置を含み、
    前記画像処理方法は、前記プロセッサが、前記パワー分布又は前記カーブの少なくとも一方を、前記出力装置に出力する、画像処理方法。

Images