特開 2022-44393 画像分析方法および画像分析装置

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】公開特許公報(A)

(11)【公開番号】P2022044393

(43)【公開日】2022-03-17

(54)【発明の名称】画像分析方法および画像分析装置

(51)【国際特許分類】

   G01N 21/17 20060101AFI20220310BHJP

   G06T 7/246 20170101ALI20220310BHJP

   G01N 21/45 20060101ALI20220310BHJP

【ＦＩ】

G01N21/17 A

G06T7/246

G01N21/45 Z

【審査請求】未請求

【請求項の数】10

【出願形態】ＯＬ

(21)【出願番号】P 2020149993

(22)【出願日】2020-09-07

【新規性喪失の例外の表示】特許法第３０条第２項適用申請有り 開催年月日：令和２年８月３日 集会名、開催場所： 第２３回 画像の認識・理解シンポジウム ＭＩＲＵ２０２０ オンライン開催 ｈｔｔｐｓ：／／ｓｉｔｅｓ．ｇｏｏｇｌｅ．ｃｏｍ／ｖｉｅｗ／ｍｉｒｕ２０２０／％Ｅ３％８３％９７％Ｅ３％８３％ＡＤ％Ｅ３％８２％Ｂ０％Ｅ３％８３％Ａ９％Ｅ３％８３％Ａ０＃ｈ．ｐ＿ｖｄＮｊｍＪ１ｆＹ＿ＮＸ

(71)【出願人】

【識別番号】301021533

【氏名又は名称】国立研究開発法人産業技術総合研究所

(71)【出願人】

【識別番号】504145342

【氏名又は名称】国立大学法人九州大学

(71)【出願人】

【識別番号】510108951

【氏名又は名称】公立大学法人広島市立大学

(74)【代理人】

【識別番号】110001243

【氏名又は名称】特許業務法人 谷・阿部特許事務所

(72)【発明者】

【氏名】佐川 立昌

(72)【発明者】

【氏名】川崎 洋

(72)【発明者】

【氏名】樋口 雄介

(72)【発明者】

【氏名】古川 亮

【テーマコード（参考）】

2G059

5L096

【Ｆターム（参考）】

2G059AA05

2G059EE02

2G059EE09

2G059FF01

2G059FF04

2G059GG01

2G059KK03

2G059KK04

2G059MM01

5L096AA06

5L096BA18

5L096CA04

5L096CA05

5L096DA02

5L096FA66

(57)【要約】

【課題】撮影対象の動きを判定する。
【解決手段】制御装置は、第１の撮影信号に基づいて、第１のスペックルパターンを含む第１の画像を生成し、撮像装置が第２の時間に撮影対象を撮影することによって生成された第２の撮影信号に基づいて、第２のスペックルパターンを含む第２の画像を生成し、第１の画像内の第１の画素集合の間で第１の特徴ベクトルの変位を近似させ、第２の画像内の第２の画素集合の間で第２の特徴ベクトルの変位を近似させ、第１の特徴ベクトルと第２の特徴ベクトルとの間で、画素ごとに画素間の距離を計算することによって、第１のスペックルパターンと第２のスペックルパターンとの間の変化を判定する。
【選択図】図５

【特許請求の範囲】

【請求項1】

撮影対象の動きを判定するコンピュータデバイスであって、
撮像装置が第１の時間に前記撮影対象を撮影することによって生成された第１の撮影信号に基づいて、第１のスペックルパターンを含む第１の画像を生成し、
前記撮像装置が第２の時間に前記撮影対象を撮影することによって生成された第２の撮影信号に基づいて、第２のスペックルパターンを含む第２の画像を生成し、
前記第１の画像内の第１の画素集合における輝度から第１の特徴ベクトルを計算し、前記第１の画素集合の間で前記第１の特徴ベクトルの変位を近似させ、
前記第２の画像内の第２の画素集合における輝度から第２の特徴ベクトルを計算し、前記第２の画素集合の間で前記第２の特徴ベクトルの変位を近似させ、
前記近似させた第１の特徴ベクトルと前記近似させた第２の特徴ベクトルとの間で、画素ごとに画素間の距離を計算することによって、前記第１のスペックルパターンと前記第２のスペックルパターンとの間の変化を判定する、
ように構成された制御装置を備えたことを特徴とするコンピュータデバイス。

【請求項2】

前記制御装置は、
第１の次元数で前記第１の特徴ベクトルおよび前記第２の特徴ベクトルを計算し、
前記第１の特徴ベクトルに基づいて、前記第１の次元数よりも少ない第２の次元数で前記第１の特徴ベクトルを計算し、
前記第２の特徴ベクトルに基づいて、前記第２の次元数で前記第２の特徴ベクトルを計算する、
ように更に構成されている、ことを特徴とする請求項１に記載のコンピュータデバイス。

【請求項3】

前記制御装置は、拡散マップによって、前記第１の特徴ベクトルおよび前記第２の特徴ベクトルを計算するように更に構成されている、ことを特徴とする請求項２に記載のコンピュータデバイス。

【請求項4】

前記第２の次元数は、３～２０のいずれかである、ことを特徴とする請求項２または３に記載のコンピュータデバイス。

【請求項5】

前記制御装置は、
前記第１の画像から、特定の領域を判定することによって第１の局所平面を抽出し、
前記第２の画像から、前記第１の局所平面に対応する領域を判定することによって第２の局所平面を抽出する、ように更に構成され、
前記第１の画素集合は、前記第１の局所平面内の画素集合であり、
前記第２の画素集合は、前記第２の局所平面内の画素集合である、
ことを特徴とする請求項１乃至３のいずれか一項に記載のコンピュータデバイス。

【請求項6】

前記制御装置は、
前記第１の特徴ベクトルに線形変換を適用することによって、前記第１の画素集合の間で前記第１の特徴ベクトルの変位を近似させ、
前記第２の特徴ベクトルに線形変換を適用することによって、前記第２の画素集合の間で前記第２の特徴ベクトルの変位を近似させる、
ように更に構成されている、ことを特徴とする請求項１に記載のコンピュータデバイス。

【請求項7】

前記制御装置は、前記第１の画素集合と前記第２の画素集合との間の画素ごとの画素間の距離が近似するという制約の下、前記第１の特徴ベクトルおよび前記第２の特徴ベクトルに前記線形変換を適用するように更に構成されている、ことを特徴とする請求項６に記載のコンピュータデバイス。

【請求項8】

前記制御装置は、前記第１のスペックルパターンと前記第２のスペックルパターンとの間の前記変化に基づいて、前記撮影対象の動きを表示装置に表示する、ことを特徴とする請求項１乃至７のいずれか一項に記載のコンピュータデバイス。

【請求項9】

撮影対象の動きを判定するコンピュータデバイスによって実行される方法であって、
撮像装置が第１の時間に前記撮影対象を撮影することによって生成された第１の撮影信号に基づいて、第１のスペックルパターンを含む第１の画像を生成するステップと、
前記撮像装置が第２の時間に前記撮影対象を撮影することによって生成された第２の撮影信号に基づいて、第２のスペックルパターンを含む第２の画像を生成するステップと、
前記第１の画像内の第１の画素集合における輝度から第１の特徴ベクトルを計算し、前記第１の画素集合の間で前記第１の特徴ベクトルの変位を近似させるステップと、
前記第２の画像内の第２の画素集合における輝度から第２の特徴ベクトルを計算し、前記第２の画素集合の間で前記第２の特徴ベクトルの変位を近似させるステップと、
前記近似させた第１の特徴ベクトルと前記近似させた第２の特徴ベクトルとの間で、画素ごとに画素間の距離を計算することによって、前記第１のスペックルパターンと前記第２のスペックルパターンとの間の変化を判定するステップと、
を備えたことを特徴とする方法。

【請求項10】

コンピュータ実行可能命令を含むコンピュータプログラムであって、前記コンピュータ実行可能命令は、コンピュータデバイスによって実行されるとき、前記コンピュータデバイスに、請求項９に記載の方法を実行させることを特徴とするコンピュータプログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本開示は、画像分析方法および画像分析装置に関し、特に、撮影対象から反射した光の干渉から生じるスペックルに基づいて、撮影対象のわずかな動きを判定する画像分析方法および画像分析装置に関する。

【背景技術】

【0002】

物体の表面にレーザ光などの光を照射すると、物体の照射面の凹凸形状などによって反射した光が散乱および干渉し、照射面上に斑点状のパターンが出現する。このパターンをスペックルと称する。コヒーレント光を照射した物体から生成された画像を観察すると、画像内にスペックルが観測される。そして、物体が動くと、その動きに応じて、スペックルパターンにも変化が生じる。

【先行技術文献】

【特許文献】

【0003】

【特許文献1】特表２０１９－５１９２６９号公報

【発明の概要】

【発明が解決しようとする課題】

【0004】

画像内のスペックルパターンを分析することによって、光が照射された物体の動きを判定する技術が従来から提案されている。特許文献１は、二次的スペックルパターンに関連した画像データ片を収集することによって、生物学的検体に関連した機械的および／または音響学的振動のような、機械的および／または音響学的信号をモニタリングする技術を開示している。

【0005】

特許文献1に開示された技術は、スペックルパターンの間の空間的相関を示す相関関数を決定し、周波数で生ずる相関の変化を選択するために、決定された相関関数にフィルタリングを適用している。

【0006】

例えば、カメラと物体との間の距離が短い場合、物体がわずかに動いた場合でさえ、スペックルパターンが急激に変化することがある。このような場合、スペックルパターンを連続して観察した場合、変化の前後でパターンが無相関になり、物体の連続したわずかな動きを判定することができない。特許文献１は、上述した物体のわずかな動きを判定するための解決手段を何ら提示していない。

【課題を解決するための手段】

【0007】

一実施形態に係るコンピュータデバイスは、撮影対象の動きを判定するコンピュータデバイスであって、撮像装置が第１の時間に前記撮影対象を撮影することによって生成された第１の撮影信号に基づいて、第１のスペックルパターンを含む第１の画像を生成し、前記撮像装置が第２の時間に前記撮影対象を撮影することによって生成された第２の撮影信号に基づいて、第２のスペックルパターンを含む第２の画像を生成し、前記第１の画像内の第１の画素集合における輝度から第１の特徴ベクトルを計算し、前記第１の画素集合の間で前記第１の特徴ベクトルの変位を近似させ、前記第２の画像内の第２の画素集合における輝度から第２の特徴ベクトルを計算し、前記第２の画素集合の間で前記第２の特徴ベクトルの変位を近似させ、前記近似させた第１の特徴ベクトルと前記近似させた第２の特徴ベクトルとの間で、画素ごとに画素間の距離を計算することによって、前記第１のスペックルパターンと前記第２のスペックルパターンとの間の変化を判定する、ように構成された制御装置を含む。

【0008】

また、別の実施形態に係る方法は、撮影対象の動きを判定するコンピュータデバイスによって実行される方法であって、撮像装置が第１の時間に前記撮影対象を撮影することによって生成された第１の撮影信号に基づいて、第１のスペックルパターンを含む第１の画像を生成するステップと、前記撮像装置が第２の時間に前記撮影対象を撮影することによって生成された第２の撮影信号に基づいて、第２のスペックルパターンを含む第２の画像を生成するステップと、前記第１の画像内の第１の画素集合における輝度から第１の特徴ベクトルを計算し、前記第１の画素集合の間で前記第１の特徴ベクトルの変位を近似させるステップと、前記第２の画像内の第２の画素集合における輝度から第２の特徴ベクトルを計算し、前記第２の画素集合の間で前記第２の特徴ベクトルの変位を近似させるステップと、前記近似させた第１の特徴ベクトルと前記近似させた第２の特徴ベクトルとの間で、画素ごとに画素間の距離を計算することによって、前記第１のスペックルパターンと前記第２のスペックルパターンとの間の変化を判定するステップと、を含む。

【発明の効果】

【0009】

実施形態に係る画像分析システムによれば、スペックルパターンの変化の前後で相関性を維持することができ、物体のわずかな動きを精度よく判定することができる。

【図面の簡単な説明】

【0010】

【図1】画像分析システムの構成の例を示すブロック図である。

【図2】コンピュータデバイスの構成の例を示すブロック図である。

【図3】撮影対象を撮影する例を示す図である。

【図4】局所平面の例を示す図である。

【図5】画像分析システムが実行する処理の例を示すフローチャートである。

【図6】局所平面の例を示す図である。

【図7】布を指で押した状態の例を示す図である。

【図8】振動の変化の例を示すグラフである。

【発明を実施するための形態】

【0011】

以下、添付した図面を参照して、一実施形態に係る画像分析方法および画像分析装置を詳細に説明する。本実施形態では、画像分析方法および画像分析装置は、画像分析システムにおいて実装される。また、本実施形態では、光が照射され、動きなどを判定される対象となるいずれかの物体を「撮影対象」と称する。

【0012】

まず、図１を参照して、画像分析システム１００の構成の例を説明する。本実施形態では、画像分析システム１００は、コンピュータデバイス１、撮像装置２ａ、２ｂ、および２ｃ（総称して、撮像装置２）、投光装置３ａ、３ｂ、および３ｃ（総称して、投光装置３）、ならびに表示装置４を含む。コンピュータデバイス１は、少なくとも撮像装置２および表示装置４と、バスまたはネットワークなどを介して相互に結合され、任意選択で、投光装置３と結合されている。

【0013】

コンピュータデバイス１は、少なくとも演算機能を含むいずれかの情報処理装置である。コンピュータデバイス１は、撮像装置２から受信した撮影信号に基づいて画像を生成し、画像内のスペックルパターンに基づいて撮影対象Ｏの動きを判定する。

【0014】

撮像装置２は、１つまたは複数のＣＣＤイメージセンサまたはＣＭＯＳイメージセンサなどを含み、撮影対象Ｏを撮影し、撮影信号をコンピュータデバイス１に送信する。本実施形態では、撮像装置２は、カメラによって実装される。撮像装置２は、高画素の画像を生成するために、より多くのイメージセンサが配列されることが望ましい。

【0015】

なお、本実施形態では、撮影対象Ｏの連続した動きを判定するために、複数の異なる位置にある撮影対象Ｏを撮影する複数の撮像装置２（撮像装置２ａ、２ｂ、および２ｃ）が設けられるが、撮像装置２の数は、１以上のいずれかの数であってもよい。例えば、１台の撮像装置２が設けられる場合、撮像装置２は、あらゆる方向において撮影対象Ｏを撮影することができるような回転式および／または可動式光源であってもよい。

【0016】

投光装置３は、撮影対象Ｏに向けてコヒーレント光を照射する光源である。コヒーレント光は、光束内の任意の二点における光波の位相関係が時間的に不変で一定であり、任意の方法で光束を分割した後、大きな光路差を与えて再び重ねあわせても完全な干渉性を示す光をいう。本実施形態では、投光装置３は、コヒーレント光としてレーザ光を放射する。

【0017】

なお、撮像装置２と同様に、撮影対象Ｏの連続した動きを判定するために、複数の異なる位置にある撮影対象Ｏにレーザ光を照射する複数の投光装置３（投光装置３ａ、３ｂ、および３ｃ）が設けられるが、投光装置３の数は、１以上のいずれかの数であってもよい。例えば、１台の投光装置３が設けられる場合、投光装置３は、あらゆる方向にレーザ光を放射することができるような回転式および／または可動式光源であってもよい。

【0018】

表示装置４は、コンピュータデバイス１が物体の動きを判定した結果、その動きを表示するディスプレイである。表示装置４は、本実施形態での必須の構成要素ではないが、本実施形態が、例えば、物体に伝わる音の振動を観察するために適用される場合、その音の振動を視覚的に表現することができる。

【0019】

なお、本実施形態では、コンピュータデバイス１、撮像装置２、投光装置３、および表示装置４がそれぞれ独立した装置によって実装されるが、そのような構成に限定されない。それらが全て統合された単一の装置が設けられてもよく、またはそれらのうちのいずれかが組み合わされた装置が存在してもよい。

【0020】

次に、図２を参照して、コンピュータデバイス１の詳細な構成要素を説明する。コンピュータデバイス１は、制御装置１１、メモリ１２、記憶装置１３、入力ドライバ１４、および出力ドライバ１５を含む。メモリ１２、記憶装置１３、入力ドライバ１４、および出力ドライバ１５はそれぞれ、内部バスを通じて制御装置１１に結合され、制御装置１１によって制御される。

【0021】

制御装置１１は、プロセッサとも称され、中央処理装置（ＣＰＵ）およびグラフィックプロセシングユニット（ＧＰＵ）などを含む。制御装置１１は、撮像装置２から受信した撮影信号に基づいて画像を生成し、画像内のスペックルパターンに基づいて、撮影対象Ｏの動きを判定するなどの演算を実行する。なお、制御装置１１は、コンピュータデバイス１に結合された撮像装置２および投光装置３を制御してもよい。この場合、例えば、制御装置１１は、投光装置３が撮影対象Ｏにレーザ光を照射するタイミング、および撮像装置２が撮影対象Ｏを撮影するタイミングなどを制御してもよい（撮像装置２および投光装置３に制御信号を送信することによって）。

【0022】

メモリ１２は、制御装置１１が処理する、コンピュータ実行可能命令、および当該命令による演算処理後のデータなどを記憶した揮発性データ記憶装置である。メモリ１２は、ＲＡＭ（ランダムアクセスメモリ）（例えば、ＳＲＡＭ（スタティックＲＡＭ）およびＤＲＡＭ（ダイナミックＲＡＭ））などで実装されてもよい。

【0023】

記憶装置１３は、上述したコンピュータ実行可能命令を含むプログラムなどを記憶した不揮発性データ記憶装置である。記憶装置１３は、ＲＯＭ（リードオンリメモリ）などの不揮発性半導体メモリ、磁気記憶装置（ハードディスクドライブなど）、および光ディスクなどで実装されてもよい。なお、プログラムなどのデータは、記憶装置１３に加えまたはその代わりに、ＮＡＳ（Network Attached Storage）および／またはＳＡＮ（Storage Area Network）などに記憶されてもよい。

【0024】

入力ドライバ１４は、撮像装置２から撮影信号を受信するなど、コンピュータデバイス１と結合された外部機器からの入力を受信し、その入力を制御装置１１に送信する。出力ドライバ１５は、表示装置４に画像信号を送信するなど、コンピュータデバイス１と結合された外部機器への出力を制御装置から受信し、その出力を外部機器に送信する。

【0025】

次に、図３を参照して、撮像装置２が撮影対象Ｏを撮影することによって生成された画像内のスペックルパターンを説明する。上述したように、撮影対象Ｏが動くと、その動きに応じて投光装置３が照射したレーザ光が散乱および干渉し、スペックルパターンにも変化が生じる。

【0026】

図３は、破線矢印Ａに示す方向に動いている撮影対象Ｏを３つの撮影ポイントで撮影し、対応して生成された３つの画像を示す。撮像装置２ａは、撮影対象Ｏを撮影ポイントＳ１において撮影し、撮像装置２ｂは、撮影対象Ｏを撮影ポイントＳ２において撮影し、撮像装置２ｃは、撮影対象Ｏを撮影ポイントＳ３において撮影する。本実施形態では、撮影対象Ｏは、３つの異なる位置である撮影ポイントＳ１～Ｓ３のそれぞれにおいて撮影されるが、撮影対象Ｏは、同一の位置において撮影されてもよい。つまり、撮影対象Ｏは、異なる時間点において、異なる位置または同一の位置において複数回撮影される。

【0027】

撮像装置２ａが撮影対象Ｏを撮影することによって（制御装置１１により）生成された画像は、画像Ｉ１に対応する。同様に、撮像装置２ｂが撮影対象Ｏを撮影することによって生成された画像は、画像Ｉ２に対応し、撮像装置２ｃが撮影対象Ｏを撮影することによって生成された画像は、画像Ｉ３に対応する。

【0028】

上述したように、撮影対象Ｏの連続した動きを、スペックルパターンを分析することによって判定する場合、複数の撮影ポイントにおいて撮影対象Ｏを撮影することになる。それぞれの撮影ポイントにおいて、投光装置３から照射したレーザ光が散乱および干渉することになり、スペックルパターンに変化を生じさせる。図３に示すように、３つの画像Ｉ１～Ｉ３内のスペックルパターンはそれぞれ異なる。

【0029】

本実施形態では、スペックルパターンを含む画像内の特定の領域を平面と見なし、平面内の特定の画素（ｘ，ｙ）の周りの３Ｄ平面パラメータを推定する。本実施形態では、この特定の平面領域を「局所平面」と称する。図４を参照して、動いている撮影対象Ｏを（３回）撮影することによって生成された３つの画像内の局所平面を説明する。

【0030】

図４は、３つの局所平面ＬＰ１～ＬＰ３を示している。局所平面ＬＰ１は、図３に示した画像Ｉ１内の特定の局所平面に相当する。局所平面ＬＰ２は、図３に示した画像Ｉ２内の、局所平面ＬＰ１に対応する局所平面に相当する。局所平面ＬＰ３は、図３に示した画像Ｉ３内の、局所平面ＬＰ１および局所平面ＬＰ２に対応する局所平面に相当する。

【0031】

図３で説明したように、画像Ｉ１～Ｉ３は、動いている撮影対象Ｏを撮影することによって生成された画像であるので、局所平面ＬＰ２は、局所平面ＬＰ１に対して画像内の座標位置において変位している。局所平面ＬＰ３も同様に、局所平面ＬＰ１および局所平面ＬＰ２に対して画像内の座標位置において変位している。

【0032】

図４は、局所平面ＬＰ１を基準に、局所平面ＬＰ２および局所平面ＬＰ３の変位を示している。変位は、２つの局所平面内の画素間の距離によって表される。図４に示す例では、局所平面ＬＰ１内の特定の画素Ｐ１ａおよび画素Ｐ１ａの周りの隣接した画素Ｐ１ｂを基準とする。

【0033】

まず、局所平面をフレームｉとみなし、フレームｉ内の画素Ｐ１ａ（ｘ，ｙ）の周りの３Ｄ平面パラメータを、式（１）に従って推定する。

【0034】

【数1】

【0035】

【数2】

は、３次元ベクトル、すなわち、（ｘ，ｙ）の各点で局所平面を定義する３次元パラメータを表す関数であり、

【0036】

【数3】

は、（ｘ，ｙ）の近隣の画素Ｐ１ｂを表す。

【0037】

図４に示す通り、局所平面ＬＰ２内の画素Ｐ２ａの周りの３Ｄ平面パラメータは、ψ_i（ｘ，ｙ）として表し、局所平面ＬＰ３内の画素Ｐ３ｂの周りの３Ｄ平面パラメータは、

【0038】

【数4】

【0039】

として表す。

【0040】

上述したように、図３に示した画像Ｉ１～Ｉ３は、撮影対象Ｏの動きを表すので、局所平面ＬＰ１内の画素Ｐ１ａと局所平面ＬＰ２内の画素Ｐ２ａと間のオフセットは、変位ｄ１（ｄ（ｘ，ｙ））として表される。変位ｄ１は、図３に示した、撮影ポイントＳ１と撮影ポイントＳ２との間の撮影対象Ｏの移動に対応する。

【0041】

変位ｄ１は、式（２）に従って計算される。

【0042】

【数5】

【0043】

同様に、局所平面ＬＰ１内の画素Ｐ１ｂと局所平面ＬＰ３内の画素Ｐ２ｂと間のオフセットは、変位ｄ２（

【0044】

【数6】

【0045】

）として表される。変位ｄ２は、図３に示した、撮影ポイントＳ１と撮影ポイントＳ３との間の撮影対象Ｏの移動に対応する。

【0046】

上述した変位ｄ１およびｄ２を計算することによって、画像内のスペックルパターンの変化を分析することができる。スペックルパターンは、撮像装置２と撮影対象Ｏとの間の相対位置に影響を受けやすい。よって、撮影対象Ｏを複数回撮影することにより生成されたそれぞれの画像内のスペックルパターンの変化を分析することによって、撮影対象Ｏの動きを判定することができる。

【0047】

スペックルパターンが急激に変化すると、変化の前後でスペックルパターンが無相関になる。このことは、撮影対象Ｏの連続した動きを判定することができず、わずかな変化を正確に判定することができないことにつながる。本実施形態に係る画像分析システム１００は、スペックルパターンの変化に関わらず、スペックルパターンの一貫性を維持しながら、スペックルパターンに基づいて撮影対象Ｏの動きを判定する。

【0048】

図５のフローチャートを参照して、画像分析システム１００が実行する処理を説明する。図５において説明する処理は、動いている撮影対象Ｏを撮影ポイントＳ１およびＳ２において撮影することによって生成された画像Ｉ１およびＩ２内のスペックルパターンの変化を分析する。

【0049】

まず、投光装置３ａおよび投光装置３ｂがそれぞれ、撮影対象Ｏにレーザ光を照射する（ステップＳ５０１）。次に、撮像装置２ａが撮影ポイントＳ１において撮影対象Ｏを撮影し、撮影信号をコンピュータデバイス１に送信し、コンピュータデバイス１の制御装置１１が撮影信号に基づいて画像Ｉ１を生成する（ステップＳ５０２）。撮像装置２ｂが撮影ポイントＳ２において撮影対象Ｏを撮影し、撮影信号をコンピュータデバイス１に送信し、制御装置１１が撮影信号に基づいて画像Ｉ２を生成する（ステップＳ５０３）。

【0050】

次に、制御装置１１は、画像Ｉ１から局所平面ＬＰ１を抽出する（ステップＳ５０４）。局所平面ＬＰ１は、例えば、ハフ変換などの技術を使用して、特定の領域を特定することによって抽出される。なお、以下では、局所平面ＬＰ１をフレームｉとみなす。本実施形態では、局所平面ＬＰ１は、矩形領域に相当するが、その領域の形状は任意の形状を有してもよい。

【0051】

また、制御装置１１は、局所平面ＬＰ１に基づいて、局所平面ＬＰ１（フレームｉ）に対応する局所平面ＬＰ２を抽出する（ステップＳ５０５）。局所平面ＬＰ２は、例えば、射影行列を使用した射影変換を使用して、局所平面ＬＰ１に対応する領域を特定することによって抽出される。なお、以下では、局所平面ＬＰ１をフレームｊと見なす。ステップＳ５０５の処理によって、フレームｉ内の画素Ｐ１ａおよびＰ１ｂに対応する、フレームｊ内の画素Ｐ２ａおよびＰ２ｂも特定される。

【0052】

次に、制御装置１１は、フレームｉとフレームｊとの間の変位を求めるよう、それぞれのフレーム内の画素ごとに、画素の間の距離を計算する（ステップＳ５０６）。ここで、フレームｉ内の特定の画素Ｐ１ａ（ｘ，ｙ）についての輝度をＩ_i（ｘ，ｙ）とし、画素集合（ｘ_k，ｙ_k）についての輝度集合をＩ_i（ｘ_k，ｙ_k）とする。ｋは、フレーム内の画素の次元数を表す。本実施形態では、ｋ＝１１²とする。

【0053】

また、フレームｊ内の画素Ｐ１ａに対応する画素Ｐ２ａ（ｘ，ｙ）についての輝度をＩ_j（ｘ，ｙ）とし、画素集合（ｘ_k，ｙ_k）についての輝度集合をＩ_j（ｘ_k，ｙ_k）とする。ステップＳ５０６における画素間の距離は、上述した特徴ベクトルの間のユークリッド距離によって計算される。

【0054】

【数7】

【0055】

輝度集合Ｉ_i（ｘ_k，ｙ_k）は、フレームｉの状態を示す特徴ベクトルＦ_i＝Ｉ_i（ｘ₁，ｙ₁），…，Ｉ_i（ｘ_k，ｙ_k）を形成する。また、輝度集合Ｉ_j（ｘ_k，ｙ_k）は、フレームｊの状態を示す特徴ベクトルＦ_j＝Ｉ_j（ｘ₁，ｙ₁），…，Ｉ_j（ｘ_k，ｙ_k）を形成する。

【0056】

次に、制御装置１１は、フレームｉの３Ｄ平面パラメータψ_i（ｘ，ｙ）に基づいて、次元数ｋを次元数ｌ（ｌ＜ｋ）に削減する次元削減技術によって、フレームｉについての低次元ベクトルΨ_iを計算する（ステップＳ５０７）。

【0057】

上述した次元数削減技術は、例えば、拡散マップ（Diffusion map）およびIsomapなどを含む。本実施形態では、拡散マップによって、式（４）に従って、削減した次元数ｌでのフレームｉについての特徴ベクトルである低次元ベクトルΨ_iが計算される。

【0058】

【数8】

【0059】

次元削減の入力は、特徴ベクトルの集合Ｆ_i（ｘ，ｙ）（ｉ＝１，…，Ｎ）である。Ｎは、フレームｉ内の画素数である。

【0060】

また、制御装置１１は、フレームｊの３Ｄ平面パラメータψ_j（ｘ，ｙ）に基づいて、次元数ｋを次元数ｌに削減する次元削減技術（拡散マップ）によって、フレームｊについての低次元ベクトルΨ_jを計算する（ステップＳ５０８）。

【0061】

上述した次元数削減技術も同様に、拡散マップおよびIsomapなどを含む。本実施形態では、拡散マップによって、式（５）に従って、削減した次元数ｌでのフレームｉについての特徴ベクトルである低次元ベクトルΨ_jが計算される。なお、ステップＳ５０７およびステップＳ５０８において削減した次元数であるｌは、例えば、３～２０のうちのいずれかの値であることが望ましい。

【0062】

【数9】

【0063】

次元削減の入力は、特徴ベクトルの集合Ｆ_j（ｘ，ｙ）（ｊ＝１，…，Ｎ’）である。Ｎ’は、フレームｊ内の画素数である。

【0064】

削減した次元数ｌによって、式（３）に従って距離を計算することは、式（３）の代わりに、式（６）に従って計算される。

【0065】

【数10】

【0066】

ステップＳ５０７およびステップＳ５０８における処理は、次元数ｋで計算したフレームｉについての特徴ベクトルとフレームｊについての特徴ベクトルとの間の変位（ユークリッド距離）に基づいて、削減した次元数ｌでフレームｉについての特徴ベクトルとフレームｊについての特徴ベクトルとの間の変位を計算することになる。

【0067】

ステップＳ５０７およびステップＳ５０８における次元数の削減によって、後述するように削減した次元数ｌにおいて距離が計算される。例えば、物体の動きが撮像装置２に対して奥行方向の動きがなく、撮像装置２に対して平行な動きのみを含む場合、高い次元数で変位を計算すると、計算においてノイズが発生する。削減した次元数ｌにおいて距離（式（６）に従って）を計算することによって、このノイズを除去することができる。

【0068】

また、削減した次元数ｌにおいて変位を計算することによって、フレームｉ内およびフレームｊ内のそれぞれの画素の間の特徴ベクトルの変位が小さくなるので、スペックルパターンの線形性を維持することができる。

【0069】

更に、削減した次元数ｌにおいて距離を計算することによって、上述し距離を計算する際の計算量を削減することができる。特に、式（３）では、画素ごとに距離を計算するので、次元数を削減しないで距離を計算することは、莫大な計算量を要することになる。

【0070】

次に、制御装置１１は、式（７）に従って、フレームｉについての低次元ベクトルΨ_iに線形変換を適用する（ステップＳ５０９）。また、制御装置１１は、式（８）に従って、フレームｉについての低次元ベクトルΨ_iに線形変換を適用する（ステップＳ５１０）。

【0071】

ステップＳ５０９の処理は、フレームｉ内の輝度集合をＩ_i（ｘ_l，ｙ_l）（ｌは、ステップＳ５０７において次元数削減が適用された次元数を表す）のそれぞれの画素の間で、特徴ベクトルの変位を近似させる。同様に、ステップＳ５１０の処理は、フレームｊ内の輝度集合をＩ_j（ｘ_l，ｙ_l）のそれぞれの画素の間で、特徴ベクトルの変位を近似させる。

【0072】

【数11】

【0073】

【数12】

【0074】

式（７）および式（８）におけるＭ（ｘ，ｙ）は、３×ｌの変換行列であり、Ψ_i（ｘ，ｙ）は、画素Ｐ１ａ（ｘ，ｙ）についてのｌ×１列ベクトルである。

【0075】

式（７）および式（８）において、フレーム内で、特定の画素とその画素の特定の範囲内にある画素との間で特徴ベクトルの変位を近似させるように、式（９）および（１０）に従った制約Ｅ_Sが課される。この制約は、例えば、図６に示すように、局所平面ＬＰ２（フレームｊ）内の画素Ｐ２ａ（ｘ₁，ｙ₁）と局所平面ＬＰ１（フレームｉ）内の画素Ｐ１ａとの間の距離（画素間の距離（変位）ｄ１）および画素Ｐ２ｂ（ｘ₂，ｙ₂）と画素Ｐ１ｂとの間の距離（画素間の距離（変位）ｄ２）を近似させる制約である。この例では、画素間の距離ｄ１および画素間の距離ｄ２は，撮像装置２から見て同一画素に対応する実在の点までの距離の変化（変位）に対応するが、画素間の距離を表す別の定義であってもよい。

【0076】

【数13】

【0077】

【数14】

【0078】

フレームｊ内の画素Ｐ２ａ（ｘ₁，ｙ₁）および画素Ｐ２ｂ（ｘ₂，ｙ₂）のペアについての制約は、式（１１）に従った制約が課される。

【0079】

【数15】

【0080】

【数16】

【0081】

は、要素単位の積である。

【0082】

式（１０）では、全ての画素のベクトルについて線形行列Ｍ（ｘ，ｙ）を適用することによって得られる列ベクトルをｍとする。制約がＡｍ＝０によって表現され、Ａは、式（１１）から計算された係数行列であるので、ｍは、Ａ^TＡの最小固有値と関連付けられた固有ベクトルとして与えられる。線形行列Ｍ（ｘ，ｙ）が計算されると、式（７）および式（８）に従って計算された、線形変換が適用された、低次元ベクトルψ_i（ｘ，ｙ）および低次元ベクトルψ_j（ｘ，ｙ）は、空間的に一貫性が保持される（つまり、低次元ベクトルによって表される、フレーム内の画素の特徴ベクトルの変位が一定の範囲内にある）。

【0083】

固有ベクトルを計算する際の問題の１つは、例えば、フレームｉ内の画素についての低次元ベクトルΨ_i（ｘ，ｙ）の時間的分布が異なることである。画素（ｘ，ｙ）のΨ_i（ｘ，ｙ）がほぼゼロである場合、Ｍ（ｘ，ｙ）の大きさが固有ベクトルにおいて支配的であってもΨ_i（ｘ，ｙ）が誤差に影響を与えないため、Ｍ（ｘ，ｙ）の大きさは、１にほぼ等しく、他の画素についてのパラメータは、ほぼゼロになる。したがって、式（１２）に従って、時間軸に沿って係数行列Ａを計算する前にΨ_i（ｘ，ｙ）の標準偏差を正規化する。

【0084】

【数17】

【0085】

Ψ_ikは、Ψ_iのｋ番目の成分であり、

【0086】

【数18】

【0087】

は、Ψ_ikの平均値であり、εは、ゼロ割り算を回避するための小さい数である。

【0088】

例えば、高解像度を有する撮像装置２によって撮影して生成された画像から固有ベクトルを計算する負荷が高いので、画素をサブサンプリングした後に変換行列Ｍ（ｘ，ｙ）が計算される。

【0089】

図５の説明に戻ると、制御装置１１は、ステップＳ５０９において線形変換が適用された低次元ベクトルψ_i（ｘ，ｙ）およびステップＳ５１０において線形変換が適用された低次元ベクトルψ_j（ｘ，ｙ）に基づいて、フレームｉとフレームｊとの間の変位を画素ごとに計算する（ステップＳ５１１）。

【0090】

上記説明した変換行列Ｍ（ｘ，ｙ）の計算は、式（３）によって計算された変位を考慮しない。フレームｉについて、上述した低次元ベクトルψ_i（ｘ，ｙ）が特徴ベクトルＦ_i（ｘ，ｙ）の間の距離を保存するように変換行列Ｍ（ｘ，ｙ）が最適化される。

【0091】

フレームｉとフレームｊとの距離（変位）が非常に短いと想定することによって、低次元ベクトルψ_i（ｘ，ｙ）によって示される距離が特徴ベクトルＦ_i（ｘ，ｙ）によって示される距離を保存するような制約が、式（１３）および式（１４）に従って、最小化されることになる誤差関数に追加される。

【0092】

【数19】

【0093】

【数20】

【0094】

Ｅ_sは、式（１０）によって定義された空間的制約であり、λはその重みである。変換行列Ｍ（ｘ，ｙ）が最適化されると、平面パラメータψ_i（ｘ，ｙ）が最適化され、画素ごとの変位は、ψ_i（ｘ，ｙ）の第３の要素として与えられる。

【0095】

式（１３）の非線形最小化のための初期推測値ψ_i ^'は、式（１５）、式（１６）、および式（１７）に従って、画素ごとの大きさを正規化した後に、上記計算されたψ_i（ｘ，ｙ）によって与えられる。

【0096】

【数21】

【0097】

【数22】

【0098】

【数23】

【0099】

計算された変換行列Ｍ（ｘ，ｙ）がサブサンプリングされるので、変換行列Ｍ（ｘ，ｙ）は、元の解像度の画像についての変位を取得するよう、パラメータを最適化する前に補間される。補間は、初期推測値ψ_i ^'を使用することによって計算された、サブサンプリングされた特徴ベクトルの変位を補間することによって行われる。元の解像度の画素Ｐ１ａ（ｘ，ｙ）および変換行列Ｍ（ｘ，ｙ）は、画素Ｐ１ｂ

【0100】

【数24】

に対して式（１８）を満たす必要がある。

【0101】

【数25】

【0102】

この変換行列Ｍ（ｘ，ｙ）の線形式が局所パッチ内の全ての画素に対して取得されるので、変換行列Ｍ（ｘ，ｙ）は最小２乗法による解法として取得される。変換行列Ｍ（ｘ，ｙ）が全ての画素に対して取得されると、元の解像度の画像に最適化が適用される。

【0103】

図５に示した処理によって、フレームｉとフレームｊとの間の変位を計算することができる。図５で説明した例では、２つのフレームの間の変位のみを計算しているが、物体Ｏの連続した動きを判定する際は、ステップＳ５０１乃至ステップＳ５１１の処理が繰り返される。つまり、動いている物体を任意の回数で撮影し、対応する数の画像の局所平面の間の変位を計算する。

【0104】

上述したように、物体Ｏの動きがわずかであっても、物体Ｏと撮像装置２との間の距離によっては、スペックルパターンが大きく変化することがある。本実施形態によれば、それぞれの局所平面内の画素に対し、特徴ベクトルの変位を近似させるので、局所平面内での特徴ベクトルの変化が小さくなる。そして、後続の局所平面内の画素に対しても、特徴ベクトルの変位を近似させるので、局所平面の間でも特徴ベクトルの変化が小さくなる。このことは、スペックルパターンの変化の前後で相関性を維持することにつながる。

【0105】

また、変位を計算する際に、次元数を削減するので（ｋからｌに、ｌは３～２０の任意の数が望ましい）、高い次元数で変位を計算することによって生じるノイズを削減することができ、計算量を削減することができる。

【0106】

本実施形態に係る処理によって、変化するスペックルパターンの変化に従って、例えば、物体の動きを視覚的に表現することができる。図７および図８を参照して、その例を説明する。図７は、キャンバス性の布の表面を指で押している状態を示す。布の表面を指で押すと、指の力によって、押した箇所から振動が伝達する。

【0107】

このような状態で、振動を伝達している布の任意の箇所（同一または異なる箇所）において複数回に布の表面を撮影し、図５で説明した処理に従ってスペックルパターンの変化を分析することによって、例えば、図８に示すように、振動の変化を表現することができる。図８に示すグラフの横軸は、ｘ軸（つまり、振動による波が横方向に進んだ距離）を示し、縦軸は、変位（つまり、２つの画像についての局所平面の間の変位）を示す。

【0108】

図８に示すグラフは、制御装置１１による制御の下、計算した変位に基づいて、表示装置４に表示される。変位は、振動の振幅に対応し、ｘ軸は、変位を加算することによって導出される。このようにして、動いている物体から生成された画像内のスペックルパターンを分析することによって、物体のわずから動きの変化を視覚的に表現することができる。

【0109】

以上のように、本実施形態に係る画像分析システム１００を説明した。本実施形態によれば、物体のわずかな動きを精度よく判定することができる。また、スペックルパターンの変化に応じて、動いている撮影対象Ｏ（物体）を表現することができる。例えば、判定したスペックルパターンの変化に応じて、対応する撮影対象Ｏの動きを表示装置４において動的に表示してもよい。

【0110】

上記実施形態で説明したハードウェアの構成要素は例示的なものにすぎず、その他の構成も可能であることに留意されたい。また、上記実施形態で説明した処理の順序は、必ずしも説明した順序で実行される必要がなく、任意の順序で実行されてもよい。更に、本発明の基本的な概念から逸脱することなく、追加のステップが新たに加えられてもよい。

【0111】

また、本発明の一実施形態に係る画像判定、コンピュータデバイス１によって実行されるコンピュータプログラムによって実装されるが、当該コンピュータプログラムは、非一時的記憶媒体に記憶されてもよい。非一時的記憶媒体の例は、リードオンリメモリ（ＲＯＭ）、ランダムアクセスメモリ（ＲＡＭ）、レジスタ、キャッシュメモリ、半導体メモリ装置、内蔵ハードディスクおよび取外可能ディスク装置などの磁気媒体、光磁気媒体、ならびにＣＤ－ＲＯＭディスクおよびデジタル多用途ディスク（ＤＶＤ）などの光学媒体などを含む。

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】