画像速度推定の改善
類似度を用いて連続する画像のブロック間の類似度を計算するブロックマッチング技法を用いる画像処理における画像速度推定方法。類似度を用いて候補速度の確率密度関数を計算する。この計算は、値がフレーム内の位置に依存しないパラメータを類似度に掛ける類似度の指数関数に基づく。候補速度を閾値処理して確率の低いものを除外する。全てのフレームを最初のフレームに同時位置合わせすること、位置合わせ誤差を計算すること、および、パラメータと閾値の値を変化させて位置合わせ誤差を最小化することによってパラメータおよび閾値の値を最適化する。類似度は、例えば比較されるブロック中の画像サンプルの数で割ることによってブロックのサイズに対して正規化される。用いられる類似度は、類似度を計算する前に、比較される２ブロックの平均および標準偏差が同じになるように調節されるＣＤ２−ｂｉｓ類似度であってもよい。これにより、類似度は特に超音波画像に適したものとなる。さらに、第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、ブロックの強度を比較し、第２フレームのブロックに最も良く一致する第３フレーム中のブロック、およびそのブロックの第１フレームにおける対応位置を見つけることによって、シーケンスの３フレームにわたってブロックマッチングを行ってもよい。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、画像処理に関し、特に一連の画像フレームにおける画像速度（image velocity）推定の改善に関する。
【０００２】
画像中の被写体が動いており、フレーム間で被写体の動きを追跡または測定することが望ましいイメージングの状況は多い。この動きはオプティカルフローまたは画像速度として知られている。このような画像速度の推定または測定は例えば、画像の符号化効率を改善するため、あるいは画像の何らかの特定の追跡される部分の動きの表示または測定の向上を可能にし、その画像を解釈しようとする見る人を助けるために行うことができる。画像速度推定には多くの技法が提案され使用されており、基本的な技法の１つがブロックマッチングとして知られている。ブロックマッチングでは、第１フレームにおいて画素ブロックが定義され、そうした上での目的は、次の第２フレームにおいてそれらのブロックの位置を特定することである。１つの手法は、差の二乗和のような類似度（similarity measure）を用いて、第１フレーム中のブロックの画素強度(intensities of the pixels)を、第２フレーム中のこれに続く変位した候補ブロックと比較することである。第２フレームにおいて、差の二乗和が最小となる（あるいは選択されたいずれかの類似度によって最適一致をもたらす）ブロックを、被写体の動きによって変位した同一ブロックとして解釈する。第１画像フレーム中の連続したブロックについてこのプロセスを繰り返すことにより、画像中の各位置における被写体の動き（画像速度フィールド）の推定値が得られる。
【０００３】
図１はこの概念を概略的に示す。２つのフレーム、すなわちフレーム１およびフレーム２が示される。これらは、シーケンス中の連続するフレームであってもよいが、そうである必要はない。フレーム１は、辺の長さがそれぞれ（２ｎ＋１）画素、すなわち中心画素（ｘ，ｙ）を中心として−ｎから＋ｎまでである正方形の画素ブロックに分割される。図１には１ブロックＷｃを示す。第２フレームには、第２フレーム中の対応する中心画素（ｘ，ｙ）の位置の周囲に探索窓Ｗｓが規定される。図１に示されるように、これは、辺の長さが（２Ｎ＋１）画素である正方形の探索領域である。次に、フレーム１の画素ブロックＷｃの強度を探索窓Ｗｓ中のそのブロックの全ての可能な位置において比較する。よって例えば、探索窓Ｗｓの左上の角の対応する（２ｎ＋１）×（２ｎ＋１）のブロックを用いて最初の比較を行い、次にかかるブロックを右へ１画素変位したものを用いて比較を行い、次にブロックを右へ２画素変位したものを用いて比較を行う、ということを探索窓の端に達するまで続ける。次にこの手順を、探索窓において１行目から下へ１画素変位した候補ブロックの行について繰り返す、ということを探索窓の下端に達するまで続ける。類似度は例えば次のような差の二乗和であってもよい。−Ｎ≦ｕ，ｖ≦Ｎである（ｕ，ｖ）の各値について
【０００４】
【数１】
【０００５】
であり、ここで、ｉおよびｊは、画素（ｘ，ｙ）を中心とするブロックＷｃをそれぞれｘ方向およびｙ方向においてインデックス付けし、ｕおよびｖは、探索窓Ｗｓをインデックス付けする異なる変位値である。ｕおよびｖの値は、フレーム間の時間差が与えられると、速度とみなすことができる。これにより、各推定変位についてＥｃの値が得られる。Ｅｃが最小となる推定変位は多くの場合に、ブロックの実際の変位とみなされる。これを、速度フィールドを得るためにフレーム１の全ての位置について、次にシーケンス中の全てのフレームについて繰り返す。当然、異なる類似度を用いてもよい。また、シーケンスの全てのフレームに対して、または各フレーム中の全ての画素またはブロックについてブロックマッチングを行うことが常に必要であるわけではない。ブロックＷｃは、フレーム中の画素をサブサンプリングすることができ、候補変位ｕおよびｖを２つ以上の画素によってインデックス付けすることができる。したがって、異なる分解能およびスケールで探索を行うことができる。時として、マルチスケールおよび／または多分解能手法を用いて、ブロックマッチングを先ず粗い分解能または大きなスケールで、その後、次の、より細かい分解能で、前に計算した速度値を用いて行い、より細かい分解能で必要とされる探索量を低減することもできる。
【０００６】
医用画像は、通常高レベルのノイズが見られるため、画像処理に多くの問題を呈する。例えば、心臓の超音波画像において心臓壁(cardiac walls)を追跡することは、超音波画像に高レベルのノイズが見られるため、また、心臓の動きの特徴のために困難である。特に、心臓の動きは心周期(cardiac cycle)中に変化する。超音波心臓検査図において心臓壁を識別および追跡する様々な方法がＷＯ０１／１６８８６およびＷＯ０２／４３００４において提案されているが、これは改良の余地がある難しい課題である。
【０００７】
上述のようなブロックマッチング技法の発展がA. Singh著「Image-flow computation: An estimation-theoretic framework and a unified perspective」（CVGIP: Image understanding, vol. 65, no.2, pp. 152-177, 1992）（参照により本明細書中に援用される）に提案されている。この手法では、例えば上述のようなブロックマッチング技法からの保存情報(conservation information)と、近傍情報(neighbourhood information)（すなわち、周辺画素の速度を見る）との両方が、それらに伴う誤差の推定値に基づく重みをもって組み合わされる。したがって、保存情報に基づく第１のステップでは、確率質量関数において類似度の値Ｅｃを用いて応答Ｒｃを計算する。探索窓中の各位置におけるＲｃの値は、対応する変位の確率を表す。確率質量関数は次式によって与えられる。
【０００８】
【数２】
【０００９】
ここで、Ｚは、全ての確率の総和が１となるように定義される。すなわち、
【００１０】
【数３】
【００１１】
応答に関する関数では、パラメータｋは各位置において、計算上の理由から探索窓中の最大応答が１(unity)に近くなる（正規化前は０．９５）ように選択される。その場合、速度の期待値は、各候補値にその確率を掛け、その結果を総和することによって次のように求められる。
【００１２】
【数４】
【００１３】
近傍情報を用いて別の速度推定値を得ることもできる。言い換えれば、各画素における速度がその近傍の速度と完全に無関係である可能性は低い。したがって、小さな近傍Ｗｐ中の各画素の速度が推定されたと仮定すると、その近傍画素の速度を使用することによって各画素の速度推定値を精緻化することができる。明らかに、近い近傍の速度ほど、遠く離れた画素よりも関連性が高い可能性が高い。したがって、近傍画素について計算した速度に重みを割り当て、中心画素からの距離が大きくなるほど重みを小さくする（サイズ（２ｗ＋１）（２ｗ＋１）の窓Ｗｐの２Ｄガウスマスク(Gaussian mask)を用いる）。これらの重みは、確率質量関数Ｒｎ＝（ｕｉ，ｖｉ）として解釈することができ、ここで、ｕｖ空間において、
【００１４】
【数５】
【００１５】
である（ｉはＷｐ中の画素のインデックスである）。次に、近傍情報を用いて、中心画素の速度推定値
【００１６】
【数６】
【００１７】
を、次のように計算することができる。
【００１８】
【数７】
【００１９】
Singh手法の重要な側面は、保存情報に基づく推定値と近傍情報に基づく推定値との両方について、各速度推定値の共分散行列が計算されることである。これらの共分散行列を用いて誤差を計算することができ、この誤差を、２つの推定値を合成して融合した最適推定値を得る際に重みとして用いる。
【００２０】
推定値Ｕｃｃに対応する共分散行列は次式によって与えられる。
【００２１】
【数８】
【００２２】
近傍推定値
【００２３】
【数９】
【００２４】
に対応する共分散行列は次の通りである。
【００２５】
【数１０】
【００２６】
したがって、これらのステップにより、速度の２つの推定値、すなわちＵｃｃおよび
【００２７】
【数１１】
【００２８】
がそれぞれ保存情報および近傍情報から得られ、それぞれがその誤差を表す共分散行列を有する。すると、保存情報と近傍両方の両方を考慮に入れる速度の推定値Ｕを計算することができる。対応する共分散行列によって適切に重み付けされた
【００２９】
【数１２】
【００３０】
からのこの新たな推定値の距離は、近傍情報を満たす際の誤差を表す。これを近傍誤差と呼ぶことができる。同様に、その共分散行列によって重み付けされたＵｃｃからのこの推定値の距離は、保存情報を満たす際の誤差を表す。これを保存誤差と呼ぶことができる。近傍誤差と保存誤差の和は、次のような融合された速度推定値の二乗誤差を表す。
【００３１】
【数１３】
【００３２】
速度の最適値は、この誤差を最小化する値であり、Ｕに対する誤差の勾配を０に等しくなるように設定することによって得ることができ、次式が得られる。
【００３３】
【数１４】
【００３４】
【数１５】
【００３５】
およびＳｎの値は、近傍の各画素の速度が事前に分かっているという仮定に基づいて導出されるため、式（１１）は実際には、各画素における速度の初期値を保存情報のみから取得して（ガウス−ザイデル緩和法により）反復法で解かれる。したがって、次式が成り立つ。
【００３６】
【数１６】
【００３７】
ここで、上付き文字ｍは反復数を示す。反復は、Ｕｏｐの２つの連続する値の差が所定の値よりも小さくなるまで続く。
【００３８】
この技法は、保存情報と近傍情報を有用に、かつ確率論的に合成するが、実際には、特に一般に医用画像および超音波イメージングに見られるようなタイプのノイズの多い画像を用いる場合、常に良好に機能するとは限らない。
【００３９】
マッチング技法を用いた画像速度の推定における別の共通の問題は、（例えば既知の開口問題、または特に探索窓のサイズが大きい場合の不一致による）マルチモーダル応答(multi-modal response)として知られる。マルチモーダル応答の影響を低減するための一般的な方法は、上記で説明したような２フレームだけでなく３フレームの強度を比較することである。よって図面の図２に示すように、２フレームｘｉおよびｙｉのブロックＷｃ間、ならびにｙｉおよびｚｉの２ブロックＷｃ間の類似度を求める。２フレームの比較では、時刻ｔにおける一方のフレームｘｉのブロックＷｃの強度を、探索窓Ｗｓ中の（ｕ，ｖ）の全ての値について、時刻ｔ＋１における次のフレームｙｉの候補速度（ｕ，ｖ）だけ変位した対応するブロックの強度と比較する。３フレーム手法では、ブロックＷｃの強度を、同じく探索窓Ｗｓ中の（ｕ，ｖ）の値について、時刻ｔ＋２における１つおいた次のフレームｚｉの（２ｕ，２ｖ）だけ変位したブロックの強度とも比較する。差の二乗和を類似度として用いる場合、これは次のように書き表すことができる。
【００４０】
【数１７】
【００４１】
ここで、第１項は、変位（ｕ，ｖ）だけ離れたｔおよびｔ＋１におけるフレームのブロックを比較しており、第２項は、その２倍、すなわち（２ｕ，２ｖ）だけ離れたｔおよびｔ＋２におけるフレームのブロックを比較している。これは、シーケンスの３フレームにわたって速度が一定であると仮定することに等しい。言い換えれば、時刻ｔ、ｔ＋１およびｔ＋２における３フレームについて、ｔおよびｔ＋１間の変位がｔ＋１およびｔ＋２間の変位と同じであると仮定する。この仮定は、高フレームレートのシーケンスには適当であるが、いくつかの超音波イメージング様式を含むいくつかの医用画像技法において遭遇するような低フレームレートのシーケンスには不十分である。
【００４２】
本発明は、医用画像、特に本質的にノイズの多い超音波画像に特に有効なブロックマッチングの改良に関する。
【００４３】
本発明の第１の態様は、画像フレームシーケンスを処理して、シーケンスの画像速度を推定する方法を提供する。本方法は、シーケンスの２フレーム中の画像ブロックの強度を比較し、ブロック間の類似度をブロックの強度に基づいて計算することによって、類似度を用いてブロックマッチングすることと、類似度から、比較される２ブロックが同じである確率測度を計算することと、確率測度に基づいて画像速度を推定することとを含み、確率測度は、画像フレーム中の位置に依存しない類似度のパラメトリック関数を用いて計算される。
【００４４】
好ましくは、パラメトリック関数のパラメータは、画像フレーム中の位置に依存しない。関数は類似度の単調関数、例えば指数関数であってもよく、類似度には、位置に対して不変のパラメータが掛けられる。パラメータは、計算された画像速度に基づいてシーケンス中のフレームを同時位置合わせする（coregister）こと、位置合わせ誤差(registration error)を計算すること、および、位置合わせ誤差を最小化するようにパラメータの少なくとも１つを変化させることによって最適化されることができる。位置合わせ誤差は、同時位置合わせされたフレームの強度差、例えば差の二乗和から計算することができる。
【００４５】
したがって、Singhによって提案される手法の具体例では、各位置について（探索窓中の最大応答が１に近くなるように）パラメータｋの値が設定される。これは、ｋがフレームにわたって位置毎に異なることを意味する。しかし、本発明のこの態様では、ｋの値はフレームにわたって一定であり、フレーム内の位置によって変化しない。高度に非線形な（指数）関数では応答（確率）を計算する際にｋが用いられるため、速度および誤差推定値は一様でないことに留意すべきである。これは、ｋの値の変化が大きな影響を持つためである。一方、本発明のこの態様では、画像中の全ての画素についてｋが一定であるため、処理は画像にわたって、またフレーム間で一様である。
【００４６】
上記のように、ｋの値は、例えばシーケンス中の全てのフレームを１番目のフレームに位置合わせする、すなわち計算された画像速度を用いて画像の位置を調整し、動きを相殺すること（動き補正が完全であれば、各フレームの画像が完全に位置合わせされる）、および、例えば強度差の二乗和を計算することによって位置合わせ誤差を計算することによって最適化することができる。最小の位置合わせ誤差を生じるｋの値が選択される。
【００４７】
計算された類似度は、ブロック中の画素数、またはブロック中で用いられる画像サンプル数（画像をサブサンプリングしている場合）で割ることによって正規化されることができる。
【００４８】
したがって、ここでも上記の特定の例において、Ｒｃに対する上記式（２）中のｋの値をｋ／（２ｎ＋１）２で置き換えることができる。これは、ブロックサイズが変化してもｋの値を変化させる必要がないことを意味する。特に、ｋの値は再最適化する必要がないため、１つのスケールおよび分解能の１フレームシーケンスを用いて、ある用途（例えば胸部超音波）について一旦最適化されれば、同じ用途について、他のスケールおよび分解能の他のシーケンスにも同じｋの値を用いることができる。
【００４９】
確率測度は、確率が特定の閾値よりも低い画像速度における動きが無視されるように閾値処理されることができる。閾値は、上記パラメータｋの最適化に用いたのと同じプロセスによって、すなわち、計算された画像速度に基づいてシーケンス中のフレームを同時位置合わせすること、位置合わせ誤差を計算すること、および、位置合わせ誤差を最小化するように閾値を変化させることによって最適化されることができる。閾値は位置に依存しなくてもよい。
【００５０】
本発明の第２の態様は、画像速度推定において用いられる類似度に関し、類似度が計算される前に、比較されるフレーム中のブロックＷｃの強度が、同じ平均および標準偏差を持つように正規化されることを提供する。類似度は、（Singhの差の二乗和ではなく）ＣＤ２類似度であってもよく、これは特に超音波画像に適している（B. CohenおよびI. Dinstein著「New maximum likelihood motion estimation schemes for noisy ultrasound images」（Pattern Recognition 35 (2002）, pp. 455-463を参照）。
【００５１】
本発明の第３の態様は、３フレーム間に等速を仮定するのではなく、観測される移動組織が時間を通して（少なくとも連続する３〜４フレームの間）その統計学的挙動を保存すると仮定することによって、マルチモーダル応答を回避するようにSinghの手法を変更する。
【００５２】
本発明のこの態様は、３フレームのうち、第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、ブロックの強度を比較すること、および、この比較した強度に基づいて類似度を計算することによって、シーケンスの３フレームにわたるブロックマッチングを行う。
【００５３】
第１フレームおよび第２フレームのブロックは、前の画像速度推定値（すなわち、先行するフレームを処理することから得られる画像速度推定値）に基づいて互いに対応するものとして計算されるブロックであることが好ましい。
【００５４】
したがって、本方法は、第２フレームの各ブロックについて、第３フレームのいくつかのブロックを包含する探索窓を規定すること、および、第２フレーム中の上記ブロックおよび上記ブロックの第１フレームにおける対応位置（前の画像速度推定値から推定される）に対する探索窓中の各ブロックの類似度を計算することを含んでもよい。したがって、これにより、３フレームにわたって速度が変わらないと仮定することを回避する。したがって本方法は、等速仮定が当てはまらない傾向があるフレームレートの比較的低い画像フレームシーケンスに適している。
【００５５】
本発明の異なる態様は、例えばSinghのそれに似た全体的な枠組みに有利に組み合わせることができる。したがって、Singhの手法と同様に、上記の技法を用いて推定される画像速度は、探索窓にわたって各候補変位にその対応する確率測度を掛けた値を総和することによって得ることができる。さらに、推定値は、周辺位置の推定画像速度（いわゆる近傍情報）を用いて変更することによって精緻化することができる。
【００５６】
本発明の技法は、医用画像、特に超音波画像のようなノイズの多い画像シーケンスに特に適している。
【００５７】
本発明はまた、上記の方法に従って画像を処理する装置を提供する。本発明は、例えば記憶媒体上に符号化され、適切にプログラムされたコンピュータ上で実行されると上記方法を実行するコンピュータプログラムとして具現することができる。
【００５８】
本発明がさらに、添付図面を参照して例として説明される。
【００５９】
画像速度を計算することが望まれる画像フレームシーケンスが与えられると、本発明の第１の態様は、用いられる類似度、すなわちＥｃ（ｕ，ｖ）の計算に関する。Singhが提案している画像処理アルゴリズムは、類似度として差の二乗和を用いるが、ＣＤ２および正規化相互相関（ＮＣＣ）のような他の類似度も知られている。本実施形態では、ＣＤ２類似度の修正版を用いる。ＣＤ２類似度を用いた場合、速度の最尤値（most likely value）は次のように定義される。
【００６０】
【数１８】
【００６１】
ここで、ｉはブロックを指し、ｊはブロック中の画素をインデックス付けし、ブロック中には２ｎ＋１個の画素があり、ｘｉｊおよびｙｉｊは比較する２ブロックの強度である。
【００６２】
この類似度は超音波画像の場合に、超音波画像中のノイズが乗算レイリーノイズ(multiplicative Rayleigh noise)であること、および表示される超音波画像が対数圧縮されていることを考慮に入れるため、差の二乗和や正規化相互相関のような他の類似度よりも優れている。しかしこの類似度は、両方のブロックＷｃのノイズ分布が同じであると仮定しており、この仮定は超音波画像については正しくない。超音波の減衰により均質組織の画像に不均質性が生じる。組織に依存せず、一般に取得中の所与の位置について一定である、時間利得および横方向利得の補償（深い組織ほど暗く見えるという効果およびビーム強度の変化をそれぞれ補償する）は、この減衰を完全には補償しない。したがって、本発明のこの実施形態では、ＣＤ２類似度を計算する前に強度の正規化を行う。これは、２つのデータブロックＷｃの少なくとも平均および分散が確実に同じになるようにすることによって達成される。より詳細には、上記の元の強度値ｘｉｊおよびｙｉｊに対して、ブロックの強度値の平均値を減算し、標準偏差（分散の平方根）で割ることによって、
【００６３】
【数１９】
【００６４】
および
【００６５】
【数２０】
【００６６】
という新たな値に変換する。これによって、ＣＤ２−ｂｉｓとして示すことができる新たな類似度が次のように得られる。
【００６７】
【数２１】
【００６８】
これは、本実施形態において、使用されるＥｃｃの値を計算するために用いられる類似度である。
【００６９】
マルチモーダル応答を回避するために、類似度を連続する３フレームにわたって計算することができる。しかし、上記の図２に関連して説明したような、時刻ｔにおける第１フレームと時刻ｔ＋１における次のフレームとの、および、時刻ｔ＋２における第３フレームとの比較に基づく類似度を結果として生じる通常の等速仮定を行う代わりに、時刻ｔ−１における前のフレームと時刻ｔにおける現フレームとの間の速度を計算した結果を用いる。時刻ｔ＋１におけるフレームｙｉ中の探索窓Ｗｓ内のブロックと比較される時刻ｔにおけるフレームｘｉ中のブロックを与えられると、時刻ｔ−１における前のフレーム（ｏｉとして示す）中のそのブロックの位置は既に計算されているため、前のフレーム中のその位置を（ｘ−ｕ０，ｙ−ｖ０）として示すことができ、ここで、（ｕ０，ｖ０）は前の速度（画像速度）の計算の結果であった。したがって、本発明のこの実施形態における３フレームの手法では、探索窓Ｗｓ内の各候補ブロックの強度は、時刻ｔにおけるフレームｘｉ中の（ｘ，ｙ）にあるブロックの強度、および時刻ｔ−１におけるフレームｏｉ中のそのブロックの計算された位置（ｘ−ｕ０，ｙ−ｖ０）と比較される。（ｘｉとｙｉとの、および、ｏｉとｙｉとの）比較毎にＥの値が計算され、これらの値が総和される。
【００７０】
これを図３に概略的に示す。この手法は、類似度を用いて強度を比較するものなら何にでも適用することができる。差の二乗和の場合、新たな類似度は次のようになる。
【００７１】
【数２２】
【００７２】
ここで、第１の項はフレームｏｉとｙｉの、すなわち時刻ｔ−１およびｔ＋１における強度を比較し、第２の項はフレームｘｉおよびｙｉの、すなわち時刻ｔおよびｔ＋１における強度を比較する。
【００７３】
上で定義したＣＤ２−ｂｉｓ類似度の場合、３フレームにわたるＥの計算は次のようになる。
【００７４】
【数２３】
【００７５】
あるいはより詳細には次のようになる。
【００７６】
【数２４】
【００７７】
ここで、
【００７８】
【数２５】
【００７９】
は、上記で詳述したように（ただし当然ながら画像全体ではなく考慮中のブロック内のみで）変換された強度データＩを表す。
【００８０】
これにより、速度が３フレームにわたって同じであるという仮定が回避される。その代わりに、これは、ｘｉ中のブロックおよび計算された（ｏｉ中の）ブロックの前の位置に対するフレームｙｉにおける最適一致を探す。これによりマッチングプロセスは、特に例えば２０〜３０Ｈｚの低いフレームレートにおいて改善される。よってこのマッチングプロセスは、通常フレームレートが低いコントラスト超音波心臓検査図、腹部画像、組織ドップラーおよびリアルタイム３Ｄ画像の場合に特に適したものとなる。
【００８１】
新たな類似度を確立したら、次の段階は、類似度から確率質量関数を計算することである。Singh手法では、これは上記式（２）によるものであった。上述のように、それは、フレーム中の各位置についてｋの値を、探索窓中の最大応答が１に近くなるように設定することを伴った。しかし、本発明のこの実施形態では、ｋの値は代わりに、フレーム中の全位置およびシーケンス中の全フレームについて同じになるように設定される。ｋの値は、後述の最適化手法において求められる。ｋの値が与えられると、本実施形態の確率質量関数は次式によって与えられる。
【００８２】
【数２６】
【００８３】
ここで、ｍは、数値の不安定性を回避するためにＥｃ（ｕ，ｖ）から演繹される探索窓Ｗｓ中の（すなわち−Ｎ≦ｕ，ｖ≦Ｎの）類似度の最大値である。
【００８４】
したがって、類似度は、ｋの値をブロックＷｃのサイズで割ることによって変更されることが分かる。これは、１画像シーケンスについて計算されるｋの最適化値をそのシーケンスに対し全てのスケールおよび分解能で（すなわち選択されるブロックＷｃのサイズに関係なく）用いることができるようにするために必要である。
【００８５】
次に、この式を用いて計算される応答Ｒｃの値を用いて、上記の式（４）、（５）および（８）により速度の期待値（ｕｃｃ，ｖｃｃ）および対応する共分散行列を計算する。しかし、本実施形態では、式（４）および（５）の速度推定値において特定の閾値αよりも高い確率を有する候補速度のみを使用することによって速度（ｕｃｃ，ｖｃｃ）の計算をさらに変更するが、共分散の計算では全ての候補速度を用いる。したがって、本実施形態では、速度推定値を次のように計算する。
【００８６】
【数２７】
【００８７】
ここで、
【００８８】
【数２８】
【００８９】
である。
【００９０】
閾値αは次のように定義される。
【００９１】
【数２９】
【００９２】
ここで、
【００９３】
【数３０】
【００９４】
および
【００９５】
【数３１】
【００９６】
はそれぞれ確率質量関数Ｒｃ（ｕ，ｖ）の最大値および最小値である。
【００９７】
したがって、Ｔが１に設定された場合、閾値αはＲｃの最小値となる、すなわち候補速度の全ての値が計算に用いられ、この計算はSingh手法の計算に等しくなることが分かる。一方、Ｔ＝０である場合、閾値αは応答の最大値となり、最大の確率を有する候補速度のみが推定速度とみなされる。したがって、推定値は完全に優勢なモード(predominant mode)にバイアスされる。実際に、以下で説明するｋに用いられるのと同じ最適化プロセスにおいて最適化されるＴの値は、事実上０と１との間となる。
【００９８】
上述の反復法における速度および共分散行列の推定を近傍情報とともに用いて、上記の式（１２）に従って最適化された速度値を計算する。
【００９９】
図４は、２Ｄ空間においてｋおよびＴの両方の値を最適化する方法を概略的に示す。ステップ４０において画像シーケンスを取得し、ステップ４１においてｋおよびＴの値を初期化する。次にステップ４２において、ｋおよびＴの初期値を用いて画像速度を推定する。これらの初期値は、イメージング機器のタイプおよびイメージングシーケンスの被写体に基づく経験から選択することができる。このプロセスは、ｋおよびＴの選択について比較的頑強であるため、例えば、Ｔ＝０．５およびｋ＝０．５という初期値が超音波イメージングシーケンスに適している可能性がある。ステップ４２において画像速度を計算したら、ステップ４３において以後のフレームをすべて最初のフレームと位置合わせすることが可能である。フレームの「位置合わせ」は、画像を互いに重ね合わせ、それらの相対位置を調整して最適一致を得ることに等しい。実際には、このプロセスは、計算された画像速度を用いて以後のフレームの動きを修正することを伴う。フレームを位置合わせしたら、ステップ４４において誤差関数を用いて位置合わせ誤差ξを計算する。一例として、誤差関数は、フレームの強度における差の二乗和であってもよい。画像速度の推定が完全であった場合、（動き補正が完全であるため）強度に差はなく、誤差関数はゼロになる。当然、実際に、誤差関数は非ゼロであり、ステップ４５においてｋおよびＴの値を変化させて誤差関数ξを最小化する。これは、パウェル（Powell）アルゴリズムのような多次元最小化アルゴリズムを用いて達成することができる（William H. Press等著「Numerical recipes in C: The art of scientific computing」（Cambridge University Press）を参照）。最適化プロセスは、誤差関数ξの値の変化が所定の閾値を下回るまで続けることができる。胸部圧縮シーケンスの歪みを補償する１つの実験において、最適値はＴ＝０．６６０およびｋ＝０．２３７であることが分かった。図６は、超音波胸部データに対して行った実験の結果を示す。図示の誤差は位置合わせ誤差ξである。以下の２つの観測を行うことができる。
【０１００】
１．Ｔ＝０の場合、速度推定は、確率の最大値の項（argument）をとることに等しい。したがって、理論的には、パラメータｋは結果に何ら影響を及ぼさない。これは、この実験において容易に観測することができ、最大誤差に対応する。この場合、画像速度は画像の画素分解能によって定量化されるため、画像速度の誤差は画素分解能のオーダーである。さらに、この手法はノイズに対して頑強でない。これは、速度推定に対するこの高い誤差を説明する。
【０１０１】
２．Ｔ＝１の場合（Singhと同様）、速度推定値は確率の平均をとることに等しい。この結果はＴ＝０の場合よりも良好であるが、最適値に対応しない。この結果は、確率の平均を速度の推定値としてとることは、あまり正確でなく、ｐｄｆがモノモーダルでないか、ピークをあまり示さない（non-well-peaked）ｐｄｆである場合にバイアスされた推定値をもたらす可能性があるということによって説明することができる。２つのパラメータ（Ｔおよびｋ）間の予期される関数依存性も認める。この最後の点は、Ｔおよびｋの最適値の探索が２Ｄ空間で行われるべきであることを示す。この実験において、最適値はＴ＝０．６６０およびｋ＝０．２３７である。したがって、Singhの結果とは明確な違いを示す。
【０１０２】
上記の改良は疎密（coarse-to-fine）方法、すなわち速度を初めに低分解能で推定し、次のより細かい分解能において、これらの推定値を推定プロセスにおける初期推定値（first guess）として用いて精緻化する多分解能手法において用いることができることに留意すべきである。これは、第２フレームの（ｘ，ｙ）周辺で窓内を探索する代わりに、（ｘ＋ｕｅｓｔ，ｙ＋ｖｅｓｔ）周辺を探索することができることを意味し、ここで、ｕｅｓｔおよびｖｅｓｔは、より粗いレベルから伝播された速度推定値である。この手法は計算上効率が良い。さらに、画像速度の推定は、画像全体にわたって行うのではなく、動いている画像領域に集中することができる。
【０１０３】
図５は、ｋおよびＴの値（例えば、これらが所与の用途の最初の推定値である場合は初期値、または最適値）を与えられた場合の画像速度推定のプロセスフローを示す。先ず、ステップ５０において画像フレームシーケンスを取得する。次にステップ５１において、ＣＤ２−ｂｉｓ類似度を用いて、すなわち所望のスケールおよび分解能における式（１８）を用いてシーケンスの３フレームセットにわたる類似度を計算する。「分解能」は、全ての画素をサンプリングしているのか、あるいはブロックＷｃ中の特定の画素のみをサンプリングしているのかを意味し、「スケール」は、探索窓Ｗｓにおいてブロックがどれだけ変位したか（例えば１画素、または数画素）を指す。類似度の値を計算したら、ステップ５２において式（１９）を用いて応答Ｒｃの値を計算することができる。次にステップ５３において式（２０）を用いてＵｃｃの値を、式（８）を用いて対応する共分散行列Ｓｃｃを計算する。ステップ５４において、式（６）、（７）および（９）を用いて、近傍推定値の共分散および
【０１０４】
【数３２】
【０１０５】
の値を計算する。次にステップ５５において式（１２）の反復法を用いて保存情報および近傍情報を融合し、最適化された速度推定値Ｕｏｐを得る。
【０１０６】
ステップ５６によって示されるように、このプロセスは、より細かいスケールおよび分解能で繰り返すことができ、その際、既に取得した画像速度推定値を活用することによって計算的な負荷が緩和される。
【０１０７】
ブロックマッチング技法における上記の改良は特に、超音波心臓検査図シーケンス中の心臓境界画素の追跡を可能にするという点で成功する。ブロックマッチングステップは、心臓境界を画定する輪郭の周囲のリボン(ribbon)（帯(band)）に集中して、計算負荷を低減することができる。しかし、本技法は、超音波イメージングの、心臓以外の用途にも適用することができる。
【図面の簡単な説明】
【０１０８】
【図１】ブロックマッチングプロセスの概略図である。
【図２】３フレームについて等速仮定を用いた類似度計算の概略図である。
【図３】３フレームについて移動組織の統計的な保存の仮定を用いた類似度計算の図である。
【図４】本発明の１実施形態において用いられる最適化プロセスのフロー図である。
【図５】本発明の１実施形態の全体的なプロセスを示す図である。
【図６】胸部超音波画像シーケンスのｋおよびＴの最適化を示す図である。
【技術分野】
【０００１】
本発明は、画像処理に関し、特に一連の画像フレームにおける画像速度（image velocity）推定の改善に関する。
【０００２】
画像中の被写体が動いており、フレーム間で被写体の動きを追跡または測定することが望ましいイメージングの状況は多い。この動きはオプティカルフローまたは画像速度として知られている。このような画像速度の推定または測定は例えば、画像の符号化効率を改善するため、あるいは画像の何らかの特定の追跡される部分の動きの表示または測定の向上を可能にし、その画像を解釈しようとする見る人を助けるために行うことができる。画像速度推定には多くの技法が提案され使用されており、基本的な技法の１つがブロックマッチングとして知られている。ブロックマッチングでは、第１フレームにおいて画素ブロックが定義され、そうした上での目的は、次の第２フレームにおいてそれらのブロックの位置を特定することである。１つの手法は、差の二乗和のような類似度（similarity measure）を用いて、第１フレーム中のブロックの画素強度(intensities of the pixels)を、第２フレーム中のこれに続く変位した候補ブロックと比較することである。第２フレームにおいて、差の二乗和が最小となる（あるいは選択されたいずれかの類似度によって最適一致をもたらす）ブロックを、被写体の動きによって変位した同一ブロックとして解釈する。第１画像フレーム中の連続したブロックについてこのプロセスを繰り返すことにより、画像中の各位置における被写体の動き（画像速度フィールド）の推定値が得られる。
【０００３】
図１はこの概念を概略的に示す。２つのフレーム、すなわちフレーム１およびフレーム２が示される。これらは、シーケンス中の連続するフレームであってもよいが、そうである必要はない。フレーム１は、辺の長さがそれぞれ（２ｎ＋１）画素、すなわち中心画素（ｘ，ｙ）を中心として−ｎから＋ｎまでである正方形の画素ブロックに分割される。図１には１ブロックＷｃを示す。第２フレームには、第２フレーム中の対応する中心画素（ｘ，ｙ）の位置の周囲に探索窓Ｗｓが規定される。図１に示されるように、これは、辺の長さが（２Ｎ＋１）画素である正方形の探索領域である。次に、フレーム１の画素ブロックＷｃの強度を探索窓Ｗｓ中のそのブロックの全ての可能な位置において比較する。よって例えば、探索窓Ｗｓの左上の角の対応する（２ｎ＋１）×（２ｎ＋１）のブロックを用いて最初の比較を行い、次にかかるブロックを右へ１画素変位したものを用いて比較を行い、次にブロックを右へ２画素変位したものを用いて比較を行う、ということを探索窓の端に達するまで続ける。次にこの手順を、探索窓において１行目から下へ１画素変位した候補ブロックの行について繰り返す、ということを探索窓の下端に達するまで続ける。類似度は例えば次のような差の二乗和であってもよい。−Ｎ≦ｕ，ｖ≦Ｎである（ｕ，ｖ）の各値について
【０００４】
【数１】
【０００５】
であり、ここで、ｉおよびｊは、画素（ｘ，ｙ）を中心とするブロックＷｃをそれぞれｘ方向およびｙ方向においてインデックス付けし、ｕおよびｖは、探索窓Ｗｓをインデックス付けする異なる変位値である。ｕおよびｖの値は、フレーム間の時間差が与えられると、速度とみなすことができる。これにより、各推定変位についてＥｃの値が得られる。Ｅｃが最小となる推定変位は多くの場合に、ブロックの実際の変位とみなされる。これを、速度フィールドを得るためにフレーム１の全ての位置について、次にシーケンス中の全てのフレームについて繰り返す。当然、異なる類似度を用いてもよい。また、シーケンスの全てのフレームに対して、または各フレーム中の全ての画素またはブロックについてブロックマッチングを行うことが常に必要であるわけではない。ブロックＷｃは、フレーム中の画素をサブサンプリングすることができ、候補変位ｕおよびｖを２つ以上の画素によってインデックス付けすることができる。したがって、異なる分解能およびスケールで探索を行うことができる。時として、マルチスケールおよび／または多分解能手法を用いて、ブロックマッチングを先ず粗い分解能または大きなスケールで、その後、次の、より細かい分解能で、前に計算した速度値を用いて行い、より細かい分解能で必要とされる探索量を低減することもできる。
【０００６】
医用画像は、通常高レベルのノイズが見られるため、画像処理に多くの問題を呈する。例えば、心臓の超音波画像において心臓壁(cardiac walls)を追跡することは、超音波画像に高レベルのノイズが見られるため、また、心臓の動きの特徴のために困難である。特に、心臓の動きは心周期(cardiac cycle)中に変化する。超音波心臓検査図において心臓壁を識別および追跡する様々な方法がＷＯ０１／１６８８６およびＷＯ０２／４３００４において提案されているが、これは改良の余地がある難しい課題である。
【０００７】
上述のようなブロックマッチング技法の発展がA. Singh著「Image-flow computation: An estimation-theoretic framework and a unified perspective」（CVGIP: Image understanding, vol. 65, no.2, pp. 152-177, 1992）（参照により本明細書中に援用される）に提案されている。この手法では、例えば上述のようなブロックマッチング技法からの保存情報(conservation information)と、近傍情報(neighbourhood information)（すなわち、周辺画素の速度を見る）との両方が、それらに伴う誤差の推定値に基づく重みをもって組み合わされる。したがって、保存情報に基づく第１のステップでは、確率質量関数において類似度の値Ｅｃを用いて応答Ｒｃを計算する。探索窓中の各位置におけるＲｃの値は、対応する変位の確率を表す。確率質量関数は次式によって与えられる。
【０００８】
【数２】
【０００９】
ここで、Ｚは、全ての確率の総和が１となるように定義される。すなわち、
【００１０】
【数３】
【００１１】
応答に関する関数では、パラメータｋは各位置において、計算上の理由から探索窓中の最大応答が１(unity)に近くなる（正規化前は０．９５）ように選択される。その場合、速度の期待値は、各候補値にその確率を掛け、その結果を総和することによって次のように求められる。
【００１２】
【数４】
【００１３】
近傍情報を用いて別の速度推定値を得ることもできる。言い換えれば、各画素における速度がその近傍の速度と完全に無関係である可能性は低い。したがって、小さな近傍Ｗｐ中の各画素の速度が推定されたと仮定すると、その近傍画素の速度を使用することによって各画素の速度推定値を精緻化することができる。明らかに、近い近傍の速度ほど、遠く離れた画素よりも関連性が高い可能性が高い。したがって、近傍画素について計算した速度に重みを割り当て、中心画素からの距離が大きくなるほど重みを小さくする（サイズ（２ｗ＋１）（２ｗ＋１）の窓Ｗｐの２Ｄガウスマスク(Gaussian mask)を用いる）。これらの重みは、確率質量関数Ｒｎ＝（ｕｉ，ｖｉ）として解釈することができ、ここで、ｕｖ空間において、
【００１４】
【数５】
【００１５】
である（ｉはＷｐ中の画素のインデックスである）。次に、近傍情報を用いて、中心画素の速度推定値
【００１６】
【数６】
【００１７】
を、次のように計算することができる。
【００１８】
【数７】
【００１９】
Singh手法の重要な側面は、保存情報に基づく推定値と近傍情報に基づく推定値との両方について、各速度推定値の共分散行列が計算されることである。これらの共分散行列を用いて誤差を計算することができ、この誤差を、２つの推定値を合成して融合した最適推定値を得る際に重みとして用いる。
【００２０】
推定値Ｕｃｃに対応する共分散行列は次式によって与えられる。
【００２１】
【数８】
【００２２】
近傍推定値
【００２３】
【数９】
【００２４】
に対応する共分散行列は次の通りである。
【００２５】
【数１０】
【００２６】
したがって、これらのステップにより、速度の２つの推定値、すなわちＵｃｃおよび
【００２７】
【数１１】
【００２８】
がそれぞれ保存情報および近傍情報から得られ、それぞれがその誤差を表す共分散行列を有する。すると、保存情報と近傍両方の両方を考慮に入れる速度の推定値Ｕを計算することができる。対応する共分散行列によって適切に重み付けされた
【００２９】
【数１２】
【００３０】
からのこの新たな推定値の距離は、近傍情報を満たす際の誤差を表す。これを近傍誤差と呼ぶことができる。同様に、その共分散行列によって重み付けされたＵｃｃからのこの推定値の距離は、保存情報を満たす際の誤差を表す。これを保存誤差と呼ぶことができる。近傍誤差と保存誤差の和は、次のような融合された速度推定値の二乗誤差を表す。
【００３１】
【数１３】
【００３２】
速度の最適値は、この誤差を最小化する値であり、Ｕに対する誤差の勾配を０に等しくなるように設定することによって得ることができ、次式が得られる。
【００３３】
【数１４】
【００３４】
【数１５】
【００３５】
およびＳｎの値は、近傍の各画素の速度が事前に分かっているという仮定に基づいて導出されるため、式（１１）は実際には、各画素における速度の初期値を保存情報のみから取得して（ガウス−ザイデル緩和法により）反復法で解かれる。したがって、次式が成り立つ。
【００３６】
【数１６】
【００３７】
ここで、上付き文字ｍは反復数を示す。反復は、Ｕｏｐの２つの連続する値の差が所定の値よりも小さくなるまで続く。
【００３８】
この技法は、保存情報と近傍情報を有用に、かつ確率論的に合成するが、実際には、特に一般に医用画像および超音波イメージングに見られるようなタイプのノイズの多い画像を用いる場合、常に良好に機能するとは限らない。
【００３９】
マッチング技法を用いた画像速度の推定における別の共通の問題は、（例えば既知の開口問題、または特に探索窓のサイズが大きい場合の不一致による）マルチモーダル応答(multi-modal response)として知られる。マルチモーダル応答の影響を低減するための一般的な方法は、上記で説明したような２フレームだけでなく３フレームの強度を比較することである。よって図面の図２に示すように、２フレームｘｉおよびｙｉのブロックＷｃ間、ならびにｙｉおよびｚｉの２ブロックＷｃ間の類似度を求める。２フレームの比較では、時刻ｔにおける一方のフレームｘｉのブロックＷｃの強度を、探索窓Ｗｓ中の（ｕ，ｖ）の全ての値について、時刻ｔ＋１における次のフレームｙｉの候補速度（ｕ，ｖ）だけ変位した対応するブロックの強度と比較する。３フレーム手法では、ブロックＷｃの強度を、同じく探索窓Ｗｓ中の（ｕ，ｖ）の値について、時刻ｔ＋２における１つおいた次のフレームｚｉの（２ｕ，２ｖ）だけ変位したブロックの強度とも比較する。差の二乗和を類似度として用いる場合、これは次のように書き表すことができる。
【００４０】
【数１７】
【００４１】
ここで、第１項は、変位（ｕ，ｖ）だけ離れたｔおよびｔ＋１におけるフレームのブロックを比較しており、第２項は、その２倍、すなわち（２ｕ，２ｖ）だけ離れたｔおよびｔ＋２におけるフレームのブロックを比較している。これは、シーケンスの３フレームにわたって速度が一定であると仮定することに等しい。言い換えれば、時刻ｔ、ｔ＋１およびｔ＋２における３フレームについて、ｔおよびｔ＋１間の変位がｔ＋１およびｔ＋２間の変位と同じであると仮定する。この仮定は、高フレームレートのシーケンスには適当であるが、いくつかの超音波イメージング様式を含むいくつかの医用画像技法において遭遇するような低フレームレートのシーケンスには不十分である。
【００４２】
本発明は、医用画像、特に本質的にノイズの多い超音波画像に特に有効なブロックマッチングの改良に関する。
【００４３】
本発明の第１の態様は、画像フレームシーケンスを処理して、シーケンスの画像速度を推定する方法を提供する。本方法は、シーケンスの２フレーム中の画像ブロックの強度を比較し、ブロック間の類似度をブロックの強度に基づいて計算することによって、類似度を用いてブロックマッチングすることと、類似度から、比較される２ブロックが同じである確率測度を計算することと、確率測度に基づいて画像速度を推定することとを含み、確率測度は、画像フレーム中の位置に依存しない類似度のパラメトリック関数を用いて計算される。
【００４４】
好ましくは、パラメトリック関数のパラメータは、画像フレーム中の位置に依存しない。関数は類似度の単調関数、例えば指数関数であってもよく、類似度には、位置に対して不変のパラメータが掛けられる。パラメータは、計算された画像速度に基づいてシーケンス中のフレームを同時位置合わせする（coregister）こと、位置合わせ誤差(registration error)を計算すること、および、位置合わせ誤差を最小化するようにパラメータの少なくとも１つを変化させることによって最適化されることができる。位置合わせ誤差は、同時位置合わせされたフレームの強度差、例えば差の二乗和から計算することができる。
【００４５】
したがって、Singhによって提案される手法の具体例では、各位置について（探索窓中の最大応答が１に近くなるように）パラメータｋの値が設定される。これは、ｋがフレームにわたって位置毎に異なることを意味する。しかし、本発明のこの態様では、ｋの値はフレームにわたって一定であり、フレーム内の位置によって変化しない。高度に非線形な（指数）関数では応答（確率）を計算する際にｋが用いられるため、速度および誤差推定値は一様でないことに留意すべきである。これは、ｋの値の変化が大きな影響を持つためである。一方、本発明のこの態様では、画像中の全ての画素についてｋが一定であるため、処理は画像にわたって、またフレーム間で一様である。
【００４６】
上記のように、ｋの値は、例えばシーケンス中の全てのフレームを１番目のフレームに位置合わせする、すなわち計算された画像速度を用いて画像の位置を調整し、動きを相殺すること（動き補正が完全であれば、各フレームの画像が完全に位置合わせされる）、および、例えば強度差の二乗和を計算することによって位置合わせ誤差を計算することによって最適化することができる。最小の位置合わせ誤差を生じるｋの値が選択される。
【００４７】
計算された類似度は、ブロック中の画素数、またはブロック中で用いられる画像サンプル数（画像をサブサンプリングしている場合）で割ることによって正規化されることができる。
【００４８】
したがって、ここでも上記の特定の例において、Ｒｃに対する上記式（２）中のｋの値をｋ／（２ｎ＋１）２で置き換えることができる。これは、ブロックサイズが変化してもｋの値を変化させる必要がないことを意味する。特に、ｋの値は再最適化する必要がないため、１つのスケールおよび分解能の１フレームシーケンスを用いて、ある用途（例えば胸部超音波）について一旦最適化されれば、同じ用途について、他のスケールおよび分解能の他のシーケンスにも同じｋの値を用いることができる。
【００４９】
確率測度は、確率が特定の閾値よりも低い画像速度における動きが無視されるように閾値処理されることができる。閾値は、上記パラメータｋの最適化に用いたのと同じプロセスによって、すなわち、計算された画像速度に基づいてシーケンス中のフレームを同時位置合わせすること、位置合わせ誤差を計算すること、および、位置合わせ誤差を最小化するように閾値を変化させることによって最適化されることができる。閾値は位置に依存しなくてもよい。
【００５０】
本発明の第２の態様は、画像速度推定において用いられる類似度に関し、類似度が計算される前に、比較されるフレーム中のブロックＷｃの強度が、同じ平均および標準偏差を持つように正規化されることを提供する。類似度は、（Singhの差の二乗和ではなく）ＣＤ２類似度であってもよく、これは特に超音波画像に適している（B. CohenおよびI. Dinstein著「New maximum likelihood motion estimation schemes for noisy ultrasound images」（Pattern Recognition 35 (2002）, pp. 455-463を参照）。
【００５１】
本発明の第３の態様は、３フレーム間に等速を仮定するのではなく、観測される移動組織が時間を通して（少なくとも連続する３〜４フレームの間）その統計学的挙動を保存すると仮定することによって、マルチモーダル応答を回避するようにSinghの手法を変更する。
【００５２】
本発明のこの態様は、３フレームのうち、第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、ブロックの強度を比較すること、および、この比較した強度に基づいて類似度を計算することによって、シーケンスの３フレームにわたるブロックマッチングを行う。
【００５３】
第１フレームおよび第２フレームのブロックは、前の画像速度推定値（すなわち、先行するフレームを処理することから得られる画像速度推定値）に基づいて互いに対応するものとして計算されるブロックであることが好ましい。
【００５４】
したがって、本方法は、第２フレームの各ブロックについて、第３フレームのいくつかのブロックを包含する探索窓を規定すること、および、第２フレーム中の上記ブロックおよび上記ブロックの第１フレームにおける対応位置（前の画像速度推定値から推定される）に対する探索窓中の各ブロックの類似度を計算することを含んでもよい。したがって、これにより、３フレームにわたって速度が変わらないと仮定することを回避する。したがって本方法は、等速仮定が当てはまらない傾向があるフレームレートの比較的低い画像フレームシーケンスに適している。
【００５５】
本発明の異なる態様は、例えばSinghのそれに似た全体的な枠組みに有利に組み合わせることができる。したがって、Singhの手法と同様に、上記の技法を用いて推定される画像速度は、探索窓にわたって各候補変位にその対応する確率測度を掛けた値を総和することによって得ることができる。さらに、推定値は、周辺位置の推定画像速度（いわゆる近傍情報）を用いて変更することによって精緻化することができる。
【００５６】
本発明の技法は、医用画像、特に超音波画像のようなノイズの多い画像シーケンスに特に適している。
【００５７】
本発明はまた、上記の方法に従って画像を処理する装置を提供する。本発明は、例えば記憶媒体上に符号化され、適切にプログラムされたコンピュータ上で実行されると上記方法を実行するコンピュータプログラムとして具現することができる。
【００５８】
本発明がさらに、添付図面を参照して例として説明される。
【００５９】
画像速度を計算することが望まれる画像フレームシーケンスが与えられると、本発明の第１の態様は、用いられる類似度、すなわちＥｃ（ｕ，ｖ）の計算に関する。Singhが提案している画像処理アルゴリズムは、類似度として差の二乗和を用いるが、ＣＤ２および正規化相互相関（ＮＣＣ）のような他の類似度も知られている。本実施形態では、ＣＤ２類似度の修正版を用いる。ＣＤ２類似度を用いた場合、速度の最尤値（most likely value）は次のように定義される。
【００６０】
【数１８】
【００６１】
ここで、ｉはブロックを指し、ｊはブロック中の画素をインデックス付けし、ブロック中には２ｎ＋１個の画素があり、ｘｉｊおよびｙｉｊは比較する２ブロックの強度である。
【００６２】
この類似度は超音波画像の場合に、超音波画像中のノイズが乗算レイリーノイズ(multiplicative Rayleigh noise)であること、および表示される超音波画像が対数圧縮されていることを考慮に入れるため、差の二乗和や正規化相互相関のような他の類似度よりも優れている。しかしこの類似度は、両方のブロックＷｃのノイズ分布が同じであると仮定しており、この仮定は超音波画像については正しくない。超音波の減衰により均質組織の画像に不均質性が生じる。組織に依存せず、一般に取得中の所与の位置について一定である、時間利得および横方向利得の補償（深い組織ほど暗く見えるという効果およびビーム強度の変化をそれぞれ補償する）は、この減衰を完全には補償しない。したがって、本発明のこの実施形態では、ＣＤ２類似度を計算する前に強度の正規化を行う。これは、２つのデータブロックＷｃの少なくとも平均および分散が確実に同じになるようにすることによって達成される。より詳細には、上記の元の強度値ｘｉｊおよびｙｉｊに対して、ブロックの強度値の平均値を減算し、標準偏差（分散の平方根）で割ることによって、
【００６３】
【数１９】
【００６４】
および
【００６５】
【数２０】
【００６６】
という新たな値に変換する。これによって、ＣＤ２−ｂｉｓとして示すことができる新たな類似度が次のように得られる。
【００６７】
【数２１】
【００６８】
これは、本実施形態において、使用されるＥｃｃの値を計算するために用いられる類似度である。
【００６９】
マルチモーダル応答を回避するために、類似度を連続する３フレームにわたって計算することができる。しかし、上記の図２に関連して説明したような、時刻ｔにおける第１フレームと時刻ｔ＋１における次のフレームとの、および、時刻ｔ＋２における第３フレームとの比較に基づく類似度を結果として生じる通常の等速仮定を行う代わりに、時刻ｔ−１における前のフレームと時刻ｔにおける現フレームとの間の速度を計算した結果を用いる。時刻ｔ＋１におけるフレームｙｉ中の探索窓Ｗｓ内のブロックと比較される時刻ｔにおけるフレームｘｉ中のブロックを与えられると、時刻ｔ−１における前のフレーム（ｏｉとして示す）中のそのブロックの位置は既に計算されているため、前のフレーム中のその位置を（ｘ−ｕ０，ｙ−ｖ０）として示すことができ、ここで、（ｕ０，ｖ０）は前の速度（画像速度）の計算の結果であった。したがって、本発明のこの実施形態における３フレームの手法では、探索窓Ｗｓ内の各候補ブロックの強度は、時刻ｔにおけるフレームｘｉ中の（ｘ，ｙ）にあるブロックの強度、および時刻ｔ−１におけるフレームｏｉ中のそのブロックの計算された位置（ｘ−ｕ０，ｙ−ｖ０）と比較される。（ｘｉとｙｉとの、および、ｏｉとｙｉとの）比較毎にＥの値が計算され、これらの値が総和される。
【００７０】
これを図３に概略的に示す。この手法は、類似度を用いて強度を比較するものなら何にでも適用することができる。差の二乗和の場合、新たな類似度は次のようになる。
【００７１】
【数２２】
【００７２】
ここで、第１の項はフレームｏｉとｙｉの、すなわち時刻ｔ−１およびｔ＋１における強度を比較し、第２の項はフレームｘｉおよびｙｉの、すなわち時刻ｔおよびｔ＋１における強度を比較する。
【００７３】
上で定義したＣＤ２−ｂｉｓ類似度の場合、３フレームにわたるＥの計算は次のようになる。
【００７４】
【数２３】
【００７５】
あるいはより詳細には次のようになる。
【００７６】
【数２４】
【００７７】
ここで、
【００７８】
【数２５】
【００７９】
は、上記で詳述したように（ただし当然ながら画像全体ではなく考慮中のブロック内のみで）変換された強度データＩを表す。
【００８０】
これにより、速度が３フレームにわたって同じであるという仮定が回避される。その代わりに、これは、ｘｉ中のブロックおよび計算された（ｏｉ中の）ブロックの前の位置に対するフレームｙｉにおける最適一致を探す。これによりマッチングプロセスは、特に例えば２０〜３０Ｈｚの低いフレームレートにおいて改善される。よってこのマッチングプロセスは、通常フレームレートが低いコントラスト超音波心臓検査図、腹部画像、組織ドップラーおよびリアルタイム３Ｄ画像の場合に特に適したものとなる。
【００８１】
新たな類似度を確立したら、次の段階は、類似度から確率質量関数を計算することである。Singh手法では、これは上記式（２）によるものであった。上述のように、それは、フレーム中の各位置についてｋの値を、探索窓中の最大応答が１に近くなるように設定することを伴った。しかし、本発明のこの実施形態では、ｋの値は代わりに、フレーム中の全位置およびシーケンス中の全フレームについて同じになるように設定される。ｋの値は、後述の最適化手法において求められる。ｋの値が与えられると、本実施形態の確率質量関数は次式によって与えられる。
【００８２】
【数２６】
【００８３】
ここで、ｍは、数値の不安定性を回避するためにＥｃ（ｕ，ｖ）から演繹される探索窓Ｗｓ中の（すなわち−Ｎ≦ｕ，ｖ≦Ｎの）類似度の最大値である。
【００８４】
したがって、類似度は、ｋの値をブロックＷｃのサイズで割ることによって変更されることが分かる。これは、１画像シーケンスについて計算されるｋの最適化値をそのシーケンスに対し全てのスケールおよび分解能で（すなわち選択されるブロックＷｃのサイズに関係なく）用いることができるようにするために必要である。
【００８５】
次に、この式を用いて計算される応答Ｒｃの値を用いて、上記の式（４）、（５）および（８）により速度の期待値（ｕｃｃ，ｖｃｃ）および対応する共分散行列を計算する。しかし、本実施形態では、式（４）および（５）の速度推定値において特定の閾値αよりも高い確率を有する候補速度のみを使用することによって速度（ｕｃｃ，ｖｃｃ）の計算をさらに変更するが、共分散の計算では全ての候補速度を用いる。したがって、本実施形態では、速度推定値を次のように計算する。
【００８６】
【数２７】
【００８７】
ここで、
【００８８】
【数２８】
【００８９】
である。
【００９０】
閾値αは次のように定義される。
【００９１】
【数２９】
【００９２】
ここで、
【００９３】
【数３０】
【００９４】
および
【００９５】
【数３１】
【００９６】
はそれぞれ確率質量関数Ｒｃ（ｕ，ｖ）の最大値および最小値である。
【００９７】
したがって、Ｔが１に設定された場合、閾値αはＲｃの最小値となる、すなわち候補速度の全ての値が計算に用いられ、この計算はSingh手法の計算に等しくなることが分かる。一方、Ｔ＝０である場合、閾値αは応答の最大値となり、最大の確率を有する候補速度のみが推定速度とみなされる。したがって、推定値は完全に優勢なモード(predominant mode)にバイアスされる。実際に、以下で説明するｋに用いられるのと同じ最適化プロセスにおいて最適化されるＴの値は、事実上０と１との間となる。
【００９８】
上述の反復法における速度および共分散行列の推定を近傍情報とともに用いて、上記の式（１２）に従って最適化された速度値を計算する。
【００９９】
図４は、２Ｄ空間においてｋおよびＴの両方の値を最適化する方法を概略的に示す。ステップ４０において画像シーケンスを取得し、ステップ４１においてｋおよびＴの値を初期化する。次にステップ４２において、ｋおよびＴの初期値を用いて画像速度を推定する。これらの初期値は、イメージング機器のタイプおよびイメージングシーケンスの被写体に基づく経験から選択することができる。このプロセスは、ｋおよびＴの選択について比較的頑強であるため、例えば、Ｔ＝０．５およびｋ＝０．５という初期値が超音波イメージングシーケンスに適している可能性がある。ステップ４２において画像速度を計算したら、ステップ４３において以後のフレームをすべて最初のフレームと位置合わせすることが可能である。フレームの「位置合わせ」は、画像を互いに重ね合わせ、それらの相対位置を調整して最適一致を得ることに等しい。実際には、このプロセスは、計算された画像速度を用いて以後のフレームの動きを修正することを伴う。フレームを位置合わせしたら、ステップ４４において誤差関数を用いて位置合わせ誤差ξを計算する。一例として、誤差関数は、フレームの強度における差の二乗和であってもよい。画像速度の推定が完全であった場合、（動き補正が完全であるため）強度に差はなく、誤差関数はゼロになる。当然、実際に、誤差関数は非ゼロであり、ステップ４５においてｋおよびＴの値を変化させて誤差関数ξを最小化する。これは、パウェル（Powell）アルゴリズムのような多次元最小化アルゴリズムを用いて達成することができる（William H. Press等著「Numerical recipes in C: The art of scientific computing」（Cambridge University Press）を参照）。最適化プロセスは、誤差関数ξの値の変化が所定の閾値を下回るまで続けることができる。胸部圧縮シーケンスの歪みを補償する１つの実験において、最適値はＴ＝０．６６０およびｋ＝０．２３７であることが分かった。図６は、超音波胸部データに対して行った実験の結果を示す。図示の誤差は位置合わせ誤差ξである。以下の２つの観測を行うことができる。
【０１００】
１．Ｔ＝０の場合、速度推定は、確率の最大値の項（argument）をとることに等しい。したがって、理論的には、パラメータｋは結果に何ら影響を及ぼさない。これは、この実験において容易に観測することができ、最大誤差に対応する。この場合、画像速度は画像の画素分解能によって定量化されるため、画像速度の誤差は画素分解能のオーダーである。さらに、この手法はノイズに対して頑強でない。これは、速度推定に対するこの高い誤差を説明する。
【０１０１】
２．Ｔ＝１の場合（Singhと同様）、速度推定値は確率の平均をとることに等しい。この結果はＴ＝０の場合よりも良好であるが、最適値に対応しない。この結果は、確率の平均を速度の推定値としてとることは、あまり正確でなく、ｐｄｆがモノモーダルでないか、ピークをあまり示さない（non-well-peaked）ｐｄｆである場合にバイアスされた推定値をもたらす可能性があるということによって説明することができる。２つのパラメータ（Ｔおよびｋ）間の予期される関数依存性も認める。この最後の点は、Ｔおよびｋの最適値の探索が２Ｄ空間で行われるべきであることを示す。この実験において、最適値はＴ＝０．６６０およびｋ＝０．２３７である。したがって、Singhの結果とは明確な違いを示す。
【０１０２】
上記の改良は疎密（coarse-to-fine）方法、すなわち速度を初めに低分解能で推定し、次のより細かい分解能において、これらの推定値を推定プロセスにおける初期推定値（first guess）として用いて精緻化する多分解能手法において用いることができることに留意すべきである。これは、第２フレームの（ｘ，ｙ）周辺で窓内を探索する代わりに、（ｘ＋ｕｅｓｔ，ｙ＋ｖｅｓｔ）周辺を探索することができることを意味し、ここで、ｕｅｓｔおよびｖｅｓｔは、より粗いレベルから伝播された速度推定値である。この手法は計算上効率が良い。さらに、画像速度の推定は、画像全体にわたって行うのではなく、動いている画像領域に集中することができる。
【０１０３】
図５は、ｋおよびＴの値（例えば、これらが所与の用途の最初の推定値である場合は初期値、または最適値）を与えられた場合の画像速度推定のプロセスフローを示す。先ず、ステップ５０において画像フレームシーケンスを取得する。次にステップ５１において、ＣＤ２−ｂｉｓ類似度を用いて、すなわち所望のスケールおよび分解能における式（１８）を用いてシーケンスの３フレームセットにわたる類似度を計算する。「分解能」は、全ての画素をサンプリングしているのか、あるいはブロックＷｃ中の特定の画素のみをサンプリングしているのかを意味し、「スケール」は、探索窓Ｗｓにおいてブロックがどれだけ変位したか（例えば１画素、または数画素）を指す。類似度の値を計算したら、ステップ５２において式（１９）を用いて応答Ｒｃの値を計算することができる。次にステップ５３において式（２０）を用いてＵｃｃの値を、式（８）を用いて対応する共分散行列Ｓｃｃを計算する。ステップ５４において、式（６）、（７）および（９）を用いて、近傍推定値の共分散および
【０１０４】
【数３２】
【０１０５】
の値を計算する。次にステップ５５において式（１２）の反復法を用いて保存情報および近傍情報を融合し、最適化された速度推定値Ｕｏｐを得る。
【０１０６】
ステップ５６によって示されるように、このプロセスは、より細かいスケールおよび分解能で繰り返すことができ、その際、既に取得した画像速度推定値を活用することによって計算的な負荷が緩和される。
【０１０７】
ブロックマッチング技法における上記の改良は特に、超音波心臓検査図シーケンス中の心臓境界画素の追跡を可能にするという点で成功する。ブロックマッチングステップは、心臓境界を画定する輪郭の周囲のリボン(ribbon)（帯(band)）に集中して、計算負荷を低減することができる。しかし、本技法は、超音波イメージングの、心臓以外の用途にも適用することができる。
【図面の簡単な説明】
【０１０８】
【図１】ブロックマッチングプロセスの概略図である。
【図２】３フレームについて等速仮定を用いた類似度計算の概略図である。
【図３】３フレームについて移動組織の統計的な保存の仮定を用いた類似度計算の図である。
【図４】本発明の１実施形態において用いられる最適化プロセスのフロー図である。
【図５】本発明の１実施形態の全体的なプロセスを示す図である。
【図６】胸部超音波画像シーケンスのｋおよびＴの最適化を示す図である。
【特許請求の範囲】
【請求項１】
画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法であって、
前記シーケンスの２フレーム中の画像ブロックの強度を比較し、前記ブロック間の類似度を前記ブロックの強度に基づいて計算することによって、前記類似度を用いてブロックマッチングすることと、
前記類似度から、前記比較される２ブロックが同じである確率測度を計算することと、
前記確率測度に基づいて前記画像速度を推定することと
を含み、
前記確率測度は、前記画像フレーム中の位置に依存しない前記類似度のパラメトリック関数を用いて計算される、画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法。
【請求項２】
前記パラメトリック関数のパラメータは、前記画像フレーム中の位置に依存しない請求項１に記載の方法。
【請求項３】
前記パラメータの少なくとも１つは、
前記計算された画像速度に基づいて前記シーケンス中の前記フレームを同時位置合わせすることと、
位置合わせ誤差を計算することと、
前記位置合わせ誤差を最小化するように前記パラメータのうちの少なくとも１つを変化させることと
によって最適化される請求項２に記載の方法。
【請求項４】
前記位置合わせ誤差は、前記位置合わせしたフレームの強度差から計算される請求項３に記載の方法。
【請求項５】
前記位置合わせ誤差は、前記同時位置合わせしたフレームの強度差の二乗和から計算される請求項４に記載の方法。
【請求項６】
前記計算された類似度を前記ブロックのサイズに対して正規化し、前記正規化した類似度に基づいて前記確率測度を計算するステップをさらに含む請求項１〜５のいずれか一項に記載の方法。
【請求項７】
前記計算された類似度は、前記ブロック中の画像サンプル数で割ることによって正規化される請求項６に記載の方法。
【請求項８】
前記計算された類似度は、前記ブロック中の画素数で割ることによって正規化される請求項６に記載の方法。
【請求項９】
前記確率測度は前記類似度の単調関数である請求項１〜８のいずれか一項に記載の方法。
【請求項１０】
前記確率測度は、確率が閾値と所定の関係を有する前記画像速度における動きが無視されるように閾値処理される請求項１〜９のいずれか一項に記載の方法。
【請求項１１】
前記閾値は、
前記計算された画像速度に基づいて前記シーケンス中の前記フレームを同時位置合わせすることと、
位置合わせ誤差を計算することと、
前記位置合わせ誤差を最小化するように前記閾値を変化させることと
によって最適化される請求項１０に記載の方法。
【請求項１２】
前記閾値は位置に依存しない請求項１０または１１に記載の方法。
【請求項１３】
前記閾値およびパラメータは一緒に最適化される請求項１０、１１または１２に記載の方法。
【請求項１４】
前記類似度を計算する前に、前記２ブロック中の強度を、同じ平均および標準偏差を持つように正規化することをさらに含む請求項１〜１３のいずれか一項に記載の方法。
【請求項１５】
前記類似度はＣＤ２−ｂｉｓ類似度である請求項１〜１４のいずれか一項に記載の方法。
【請求項１６】
前記ブロックマッチングは、
前記シーケンスの３フレームのうち第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、ブロックの強度を比較することと、
前記比較した強度からの類似度を計算することと
によって、前記３フレームにわたって行われる請求項１〜１５のいずれか一項に記載の方法。
【請求項１７】
前記第１フレームおよび第２フレームのブロックは、前の画像速度推定値に基づいて互いに対応するものとして計算されるブロックである請求項１６に記載の方法。
【請求項１８】
画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法であって、
前記シーケンスの３フレームのうち第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、画像ブロックの強度を比較することによって前記３フレーム中の前記ブロックの強度を比較することと、
前記強度に基づいて前記ブロック間の類似度を計算することと
によって類似度を用いてブロックマッチングすることを含む、画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法。
【請求項１９】
前記第１フレームおよび前記第２フレームのブロックは、前の画像速度推定値に基づいて互いに対応するものとして計算されるブロックである請求項１８に記載の方法。
【請求項２０】
前記第２フレーム中の各ブロックについて、前記第３フレームのいくつかのブロックを包含する探索窓を規定することと、
前記前の画像速度推定値に基づいて、前記第２フレームの前記ブロックに対する、および、前記ブロックの前記第１フレームにおける対応位置に対する、前記探索窓中の各ブロックの類似度を計算することと
を含む請求項１９に記載の方法。
【請求項２１】
画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法であって、
前記シーケンスの２フレームの画像ブロックの強度を比較すること、および、前記強度に基づいて前記ブロック間の類似度を計算することによって類似度を用いてブロックマッチングすること
を含み、さらに、
前記類似度を計算する前に、前記２ブロックの強度を、同じ平均および標準偏差を持つように正規化すること
を含む、画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法。
【請求項２２】
前記類似度はＣＤ２−ｂｉｓ類似度である請求項２１に記載の方法。
【請求項２３】
前記ブロックマッチングは、
前記シーケンスの３フレームのうち、第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、ブロックの強度を比較することと、
前記比較した強度からの類似度を計算することと
によって前記３フレームにわたって行われる請求項２１または２２に記載の方法。
【請求項２４】
前記第１フレームおよび第２フレームのブロックは、前の画像速度推定値に基づいて互いに対応するものとして計算されるブロックである請求項２３に記載の方法。
【請求項２５】
前記画像速度推定値は、周辺位置の前記推定される画像速度を用いて前記画像中の各位置における前記画像速度推定値を変更することによって精緻化される請求項１〜２４のいずれか一項に記載の方法。
【請求項２６】
前記画像は医用画像である請求項１〜２５のいずれか一項に記載の方法。
【請求項２７】
前記画像は超音波画像である請求項１〜２６のいずれか一項に記載の方法。
【請求項２８】
請求項１〜２７のいずれか一項に記載の方法に従って画像速度を推定するようになっている画像速度推定器を備える画像処理装置。
【請求項２９】
プログラムされたコンピュータ上で請求項１〜２７のいずれか一項に記載の方法を実行するためのプログラムコード手段を含むコンピュータプログラム。
【請求項３０】
請求項２９に記載のコンピュータプログラムを記憶するコンピュータ読み取り可能な記憶媒体。
【請求項１】
画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法であって、
前記シーケンスの２フレーム中の画像ブロックの強度を比較し、前記ブロック間の類似度を前記ブロックの強度に基づいて計算することによって、前記類似度を用いてブロックマッチングすることと、
前記類似度から、前記比較される２ブロックが同じである確率測度を計算することと、
前記確率測度に基づいて前記画像速度を推定することと
を含み、
前記確率測度は、前記画像フレーム中の位置に依存しない前記類似度のパラメトリック関数を用いて計算される、画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法。
【請求項２】
前記パラメトリック関数のパラメータは、前記画像フレーム中の位置に依存しない請求項１に記載の方法。
【請求項３】
前記パラメータの少なくとも１つは、
前記計算された画像速度に基づいて前記シーケンス中の前記フレームを同時位置合わせすることと、
位置合わせ誤差を計算することと、
前記位置合わせ誤差を最小化するように前記パラメータのうちの少なくとも１つを変化させることと
によって最適化される請求項２に記載の方法。
【請求項４】
前記位置合わせ誤差は、前記位置合わせしたフレームの強度差から計算される請求項３に記載の方法。
【請求項５】
前記位置合わせ誤差は、前記同時位置合わせしたフレームの強度差の二乗和から計算される請求項４に記載の方法。
【請求項６】
前記計算された類似度を前記ブロックのサイズに対して正規化し、前記正規化した類似度に基づいて前記確率測度を計算するステップをさらに含む請求項１〜５のいずれか一項に記載の方法。
【請求項７】
前記計算された類似度は、前記ブロック中の画像サンプル数で割ることによって正規化される請求項６に記載の方法。
【請求項８】
前記計算された類似度は、前記ブロック中の画素数で割ることによって正規化される請求項６に記載の方法。
【請求項９】
前記確率測度は前記類似度の単調関数である請求項１〜８のいずれか一項に記載の方法。
【請求項１０】
前記確率測度は、確率が閾値と所定の関係を有する前記画像速度における動きが無視されるように閾値処理される請求項１〜９のいずれか一項に記載の方法。
【請求項１１】
前記閾値は、
前記計算された画像速度に基づいて前記シーケンス中の前記フレームを同時位置合わせすることと、
位置合わせ誤差を計算することと、
前記位置合わせ誤差を最小化するように前記閾値を変化させることと
によって最適化される請求項１０に記載の方法。
【請求項１２】
前記閾値は位置に依存しない請求項１０または１１に記載の方法。
【請求項１３】
前記閾値およびパラメータは一緒に最適化される請求項１０、１１または１２に記載の方法。
【請求項１４】
前記類似度を計算する前に、前記２ブロック中の強度を、同じ平均および標準偏差を持つように正規化することをさらに含む請求項１〜１３のいずれか一項に記載の方法。
【請求項１５】
前記類似度はＣＤ２−ｂｉｓ類似度である請求項１〜１４のいずれか一項に記載の方法。
【請求項１６】
前記ブロックマッチングは、
前記シーケンスの３フレームのうち第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、ブロックの強度を比較することと、
前記比較した強度からの類似度を計算することと
によって、前記３フレームにわたって行われる請求項１〜１５のいずれか一項に記載の方法。
【請求項１７】
前記第１フレームおよび第２フレームのブロックは、前の画像速度推定値に基づいて互いに対応するものとして計算されるブロックである請求項１６に記載の方法。
【請求項１８】
画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法であって、
前記シーケンスの３フレームのうち第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、画像ブロックの強度を比較することによって前記３フレーム中の前記ブロックの強度を比較することと、
前記強度に基づいて前記ブロック間の類似度を計算することと
によって類似度を用いてブロックマッチングすることを含む、画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法。
【請求項１９】
前記第１フレームおよび前記第２フレームのブロックは、前の画像速度推定値に基づいて互いに対応するものとして計算されるブロックである請求項１８に記載の方法。
【請求項２０】
前記第２フレーム中の各ブロックについて、前記第３フレームのいくつかのブロックを包含する探索窓を規定することと、
前記前の画像速度推定値に基づいて、前記第２フレームの前記ブロックに対する、および、前記ブロックの前記第１フレームにおける対応位置に対する、前記探索窓中の各ブロックの類似度を計算することと
を含む請求項１９に記載の方法。
【請求項２１】
画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法であって、
前記シーケンスの２フレームの画像ブロックの強度を比較すること、および、前記強度に基づいて前記ブロック間の類似度を計算することによって類似度を用いてブロックマッチングすること
を含み、さらに、
前記類似度を計算する前に、前記２ブロックの強度を、同じ平均および標準偏差を持つように正規化すること
を含む、画像フレームシーケンスを処理して前記シーケンスの画像速度を推定する方法。
【請求項２２】
前記類似度はＣＤ２−ｂｉｓ類似度である請求項２１に記載の方法。
【請求項２３】
前記ブロックマッチングは、
前記シーケンスの３フレームのうち、第１フレームと第３フレームとにおける、および、第２フレームと第３フレームとにおける、ブロックの強度を比較することと、
前記比較した強度からの類似度を計算することと
によって前記３フレームにわたって行われる請求項２１または２２に記載の方法。
【請求項２４】
前記第１フレームおよび第２フレームのブロックは、前の画像速度推定値に基づいて互いに対応するものとして計算されるブロックである請求項２３に記載の方法。
【請求項２５】
前記画像速度推定値は、周辺位置の前記推定される画像速度を用いて前記画像中の各位置における前記画像速度推定値を変更することによって精緻化される請求項１〜２４のいずれか一項に記載の方法。
【請求項２６】
前記画像は医用画像である請求項１〜２５のいずれか一項に記載の方法。
【請求項２７】
前記画像は超音波画像である請求項１〜２６のいずれか一項に記載の方法。
【請求項２８】
請求項１〜２７のいずれか一項に記載の方法に従って画像速度を推定するようになっている画像速度推定器を備える画像処理装置。
【請求項２９】
プログラムされたコンピュータ上で請求項１〜２７のいずれか一項に記載の方法を実行するためのプログラムコード手段を含むコンピュータプログラム。
【請求項３０】
請求項２９に記載のコンピュータプログラムを記憶するコンピュータ読み取り可能な記憶媒体。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図２】
【図３】
【図４】
【図５】
【図６】
【公表番号】特表２００６−５０８７２３（Ｐ２００６−５０８７２３Ａ）
【公表日】平成１８年３月１６日（２００６．３．１６）
【国際特許分類】
【出願番号】特願２００４−５５６４７３（Ｐ２００４−５５６４７３）
【出願日】平成１５年１１月１９日（２００３．１１．１９）
【国際出願番号】ＰＣＴ／ＧＢ２００３／００５０４７
【国際公開番号】ＷＯ２００４／０５２０１６
【国際公開日】平成１６年６月１７日（２００４．６．１７）
【出願人】（５００２０５２２０）アイシス イノヴェイション リミテッド (10)
【Ｆターム（参考）】
【公表日】平成１８年３月１６日（２００６．３．１６）
【国際特許分類】
【出願日】平成１５年１１月１９日（２００３．１１．１９）
【国際出願番号】ＰＣＴ／ＧＢ２００３／００５０４７
【国際公開番号】ＷＯ２００４／０５２０１６
【国際公開日】平成１６年６月１７日（２００４．６．１７）
【出願人】（５００２０５２２０）アイシス イノヴェイション リミテッド (10)
【Ｆターム（参考）】
[ Back to top ]