(57)【特許請求の範囲】
【請求項１】
解析対象となる動画像データを、前記動画像データの３次元信号と、周波数及び初期位相を用いた位相と振幅とによって表される正弦波モデル信号との差の二乗和が最小値になるような前記周波数，前記振幅及び前記初期位相を、非周期信号のフーリエ変換式のパラメータとして求め、時空間スペクトルを抽出し、抽出された前記時空間スペクトル分布を、平面をなすスペクトル群毎に切り分けして１つ又は複数の平面群を推定し、推定された前記平面群の傾きを求めることで、前記動画像内の１つ又は複数の移動物体の動きを推定する動き推定方法であって、
前記時空間スペクトルに対して最小二乗法を用いたクラスタリングを行い、１つ又は複数の平面群を推定するに際し、前記時空間スペクトルが所定の平面に属するものとした場合の存在確率変数を１又は小数に設定し、最小二乗法を用いたクラスタリングを行うことを特徴とする動き推定方法。
【請求項２】
解析対象となる動画像データを入力する動画像入力手段と、入力された前記動画像データの３次元信号と、周波数及び初期位相を用いた位相と振幅とによって表される正弦波モデル信号との差の二乗和が最小値になるような前記周波数，前記振幅及び前記初期位相を、非周期信号のフーリエ変換式のパラメータとして求め、時空間スペクトルを抽出する時空間スペクトル抽出手段と、抽出された前記時空間スペクトル分布を、平面をなすスペクトル群毎に切り分けして１つ又は複数の平面群を推定する平面推定手段と、推定された前記平面群の傾きを求める傾き算出手段とが備わっており、前記動画像内の１つ又は複数の移動物体の動きを推定する動き推定装置であって、
前記時空間スペクトルに対して最小二乗法を用いたクラスタリングを行い、１つ又は複数の平面群を推定するに際し、前記時空間スペクトルが所定の平面に属するものとした場合の存在確率変数を１又は小数に設定し、最小二乗法を用いたクラスタリングを行う構成であることを特徴とする動き推定装置。
【請求項３】
解析対象となる動画像データを入力する動画像入力手段と、入力された前記動画像データの３次元信号と、周波数及び初期位相を用いた位相と振幅とによって表される正弦波モデル信号との差の二乗和が最小値になるような前記周波数，前記振幅及び前記初期位相を、非周期信号のフーリエ変換式のパラメータとして求め、時空間スペクトルを抽出する時空間スペクトル抽出手段、前記時空間スペクトル分布を、平面をなすスペクトル群毎に切り分けして１つ又は複数の平面群を推定する平面推定手段、及び、推定された前記平面群の傾きを求める傾き算出手段としてコンピュータを機能させ、前記動画像内の１つ又は複数の移動物体の動きを推定するコンピュータ実行可能な動き推定プログラムであって、
前記時空間スペクトルに対して最小二乗法を用いたクラスタリングを行い、１つ又は複数の平面群を推定するに際し、前記時空間スペクトルが所定の平面に属するものとした場合の存在確率変数を１又は小数に設定し、最小二乗法を用いたクラスタリングを行うことを特徴とする動き推定プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、動画像内の移動物体の動き解析や動画像の圧縮符号化の分野に関し、特に、動画像内の１つ又は複数の移動物体の動きを推定する動き推定方法、動き推定装置、及び動き推定プログラムに関する。
【背景技術】
【０００２】
動画像の圧縮符号化において、参照フレームからの移動物体の動き量を予測する動き予測技術は、動画像データの記述量を大幅に削減する上で極めて重要な技術である。
【０００３】
動き予測に広く使われている方法としては、ＭＰＥＧ−４  ＡＶＣ／Ｈ．２６４として規格化されている方法が知られている。ＭＰＥＧ−４  ＡＶＣ／Ｈ．２６４の動画像符号化技術で用いられる動き予測技術においては、動画像から得られるデータを時間領域又は空間領域において解析することによって動き予測を行っている。すなわち、この動き予測技術は、２次元の画素ブロックに対して整数精度の直交変換を行う解析手法であり、また、その解析は、動画像データを構成する前後のフレームを利用し、これら前後のフレームの画素情報から予測した予測誤差画像に対して行っている。
【０００４】
しかしながら、このような動き予測は、画素ブロック単位でのパターン認識に基づくものであることから、異なる方向へ移動する複数の物体や、複雑に形状を変えながら移動する物体については対応することができず、さらに別の手法を必要とする。
【０００５】
このような別の手法としては、動画像を、高さ、幅、時間からなる３次元信号として捉えて３次元解析し、得られた時空間スペクトルを利用して動画像内の移動物体の動きを予測する手法がある。
【０００６】
一般的に、移動物体を含む動画像を３次元解析して得られた時空間スペクトルのうち、高エネルギーを有するスペクトルに注目すると、これらのスペクトルは、３次元周波数領域空間において平面状に分布する平面群となる。この平面群の構成は、動画像内の物体の動きによって決まり、平面群の傾きを求めることによって動き予測が可能となる。したがって、複数のそれぞれ異なる動きをする移動物体を含む動画像の場合には、移動物体の数と平面の数とが等しくなるため、各平面の傾きをそれぞれ求めることにより、個々の物体の動きを予測することが可能となる。
【０００７】
しかしながら、周波数解析手法として通常用いられる高速フーリエ変換（Fast
Fourier Transform；ＦＦＴ）や離散フーリエ変換（Discrete Fourier Transform；ＤＦＴ）に代表される、等間隔な周波数分解幅を有する解析手法においては、分解能がフレーム長に依存し、動画像の解析において十分な分解能が得られないことから、分解能不足に起因して解析できない動きに対応するためには、さらに別の手法を組み合わせる必要がある。
【０００８】
例えば、非特許文献１には、ＦＦＴを用いて動画像を３次元解析して時空間スペクトルを得た後、時空間スペクトルで表される平面群の傾きを決定する方法として直線群フィルタを用いる手法が開示されている。
【０００９】
また、非特許文献２に記載された技術においては、平面推定を前提としたＦＦＴアルゴリズムを利用して時空間周波数スペクトルを得ることによって動画像を３次元解析し、平面推定時にファジー推論を使用してスペクトルデータの誤差を補正することによって解析精度の向上を図っている。
【００１０】
一方、特許文献１には、本願発明者の一部が考案した周波数解析手法として、非周期信号の解析手法であるNon-Harmonic Analysis（ＮＨＡ）が開示されている。このＮＨＡは、解析対象信号と、周波数ｆ’及び初期位相φ’を用いた位相と振幅Ａ’とによって表される正弦波モデル信号との差の二乗和が最小値になるような周波数ｆ’、振幅Ａ’、及び初期位相φ’を、非周期信号のフーリエ変換式のパラメータとして算出するものである。
【先行技術文献】
【特許文献】
【００１１】
【特許文献１】国際公開第２００９／０３８０５６号
【非特許文献】
【００１２】
【非特許文献１】Akira KOJIMA and Jun-ichiHISHIGAMI, “Motion Detection using 3D-FFT Spectrum”, ITEC’92, 1992年
【非特許文献２】C. E. Erdem, G. Karabulut, E.Yanmaz and E. Anarim, “Motion Estimation in thefrequency domain using fuzzy c-planes clustering”, vol.10, pp. 1873-1879, 2001年12月
【発明の概要】
【発明が解決しようとする課題】
【００１３】
しかしながら、非特許文献１に記載された技術においては、得られる時空間スペクトルの精度が向上するわけではなく、移動物体の予測精度にも限界がある。また、非特許文献２に記載された技術においては、ファジー推論のために複雑な処理を必要とし、計算コストが高くなるという問題があった。
【００１４】
本発明は、このような実情に鑑みてなされたものであり、高精度に且つ少ない演算量で、動画像を３次元解析して得られた時空間スペクトルを利用して動画像内の移動物体の動きを推定することができる動き推定方法、動き推定装置、及び動き推定プログラムを提供することを目的とする。
【課題を解決するための手段】
【００１５】
本発明の動き推定方法は、解析対象となる動画像データを、前記動画像データの３次元信号と、周波数及び初期位相を用いた位相と振幅とによって表される正弦波モデル信号との差の二乗和が最小値になるような前記周波数，前記振幅及び前記初期位相を、非周期信号のフーリエ変換式のパラメータとして求め、時空間スペクトルを抽出し、抽出された前記時空間スペクトル分布を、平面をなすスペクトル群毎に切り分けして１つ又は複数の平面群を推定し、推定された前記平面群の傾きを求めることで、前記動画像内の１つ又は複数の移動物体の動きを推定する動き推定方法であって、
前記時空間スペクトルに対して最小二乗法を用いたクラスタリングを行い、１つ又は複数の平面群を推定するに際し、前記時空間スペクトルが所定の平面に属するものとした場合の存在確率変数を１又は小数に設定し、最小二乗法を用いたクラスタリングを行うことを特徴とする。
【００１６】
本発明の動き推定装置は、解析対象となる動画像データを入力する動画像入力手段と、入力された前記動画像データの３次元信号と、周波数及び初期位相を用いた位相と振幅とによって表される正弦波モデル信号との差の二乗和が最小値になるような前記周波数，前記振幅及び前記初期位相を、非周期信号のフーリエ変換式のパラメータとして求め、時空間スペクトルを抽出する時空間スペクトル抽出手段と、抽出された前記時空間スペクトル分布を、平面をなすスペクトル群毎に切り分けして１つ又は複数の平面群を推定する平面推定手段と、推定された前記平面群の傾きを求める傾き算出手段とが備わっており、前記動画像内の１つ又は複数の移動物体の動きを推定する動き推定装置であって、
前記時空間スペクトルに対して最小二乗法を用いたクラスタリングを行い、１つ又は複数の平面群を推定するに際し、前記時空間スペクトルが所定の平面に属するものとした場合の存在確率変数を１又は小数に設定し、最小二乗法を用いたクラスタリングを行う構成であることを特徴とする。
【００１７】
本発明の動き推定プログラムは、解析対象となる動画像データを入力する動画像入力手段と、入力された前記動画像データの３次元信号と、周波数及び初期位相を用いた位相と振幅とによって表される正弦波モデル信号との差の二乗和が最小値になるような前記周波数，前記振幅及び前記初期位相を、非周期信号のフーリエ変換式のパラメータとして求め、時空間スペクトルを抽出する時空間スペクトル抽出手段、前記時空間スペクトル分布を、平面をなすスペクトル群毎に切り分けして１つ又は複数の平面群を推定する平面推定手段、及び、推定された前記平面群の傾きを求める傾き算出手段としてコンピュータを機能させ、前記動画像内の１つ又は複数の移動物体の動きを推定するコンピュータ実行可能な動き推定プログラムであって、
前記時空間スペクトルに対して最小二乗法を用いたクラスタリングを行い、１つ又は複数の平面群を推定するに際し、前記時空間スペクトルが所定の平面に属するものとした場合の存在確率変数を１又は小数に設定し、最小二乗法を用いたクラスタリングを行うことを特徴とする動き推定プログラム。
【００１８】
このような本発明にかかる動き推定方法、動き推定装置、及び動き推定プログラムが実装された装置においては、周波数分解能が解析窓長に依存しない周波数解析手法を動画像の３次元信号解析に用いることにより、従来の周波数解析手法を用いた場合に比べて高精度に時空間スペクトルを抽出することができ、動画像内の移動物体の動きが周波数領域空間に精度良く反映される。そのため、本発明にかかる動き推定方法、動き推定装置、及び動き推定プログラムが実装された装置においては、得られた時空間スペクトルのうち、高エネルギーを有するスペクトルの理論平面からの誤差を極めて少なくすることができる。
【発明の効果】
【００１９】
本発明においては、周波数解析によって得られた時空間スペクトル分布を高精度に求めることができるため、動きの平面を推定する工程では複雑な処理を必要とせず、高精度に且つ少ない演算量で動画像内の移動物体の動きを推定することができる。
【図面の簡単な説明】
【００２０】
【図１】本発明の実施の形態として示す動き推定装置の構成を示すブロック図である。
【図２】本周波数解析手法とＤＦＴとＧＨＡとの違いを説明するための図であり、各手法の誤差を求めた結果を示す図である。
【図３】本発明の実施の形態として示す動き推定装置において、動画像内の１つ又は複数の移動物体の動きを推定する際の一連の処理を示すフローチャートである。
【図４】第１の検証にて用いた動画像を示す図である。
【図５】第１の検証において物体の移動速度Ｖｘを１．０ピクセル／フレームとした場合に、３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図６】第１の検証において物体の移動速度Ｖｘを３．３ピクセル／フレームとした場合に、３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図７】第１の検証において物体の移動速度Ｖｘを３．３ピクセル／フレームとした場合に、本周波数解析手法を用いて求めた時空間スペクトル分布を示す図である。
【図８】第２の検証にて用いた動画像を示す図である。
【図９】第２の検証において３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図１０】第２の検証において本周波数解析手法を用いて求めた時空間スペクトル分布を示す図である。
【図１１】第３の検証にて用いた動画像を示す図である。
【図１２】第３の検証において３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図１３】第３の検証において本周波数解析手法を用いて求めた時空間スペクトル分布を示す図である。
【図１４】図１３の時空間スペクトル分布の視点を変えて平面を描画した様子を示す図である。
【図１５】類似する複数の動きを含む動画像について３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図１６】類似する複数の動きを含む動画像について本周波数解析手法を用いて求めた時空間スペクトル分布を示す図である。
【図１７】図１５の場合とは異なる複数の動きを含む動画像について３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図１８】図１６の場合とは異なる複数の動きを含む動画像について３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図１９】第４の検証にて用いた動画像を示す図である。
【図２０】第４の検証において本周波数解析手法を用いて求めた１００本の時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【図２１】第４の検証において本周波数解析手法を用いて求めた１２５本の時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【図２２】第４の検証において本周波数解析手法を用いて求めた１６０本の時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【図２３】第４の検証において本周波数解析手法を用いて求めた１８０本の時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【図２４】第４の検証において本周波数解析手法を用いて求めた２００本の時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【図２５】第４の検証において本周波数解析手法を用いて求めた２２５本の時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【図２６】第４の検証において本周波数解析手法を用いて求めた２５０本の時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【図２７】第５の検証にて用いた動画像を示す図である。
【図２８】第５の検証において３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図２９】第５の検証において本周波数解析手法を用いて求めた時空間スペクトル分布を示す図である。
【図３０】第６の検証にて用いた動画像を示す図である。
【図３１】第６の検証において３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図３２】第６の検証において本周波数解析手法を用いて求めた時空間スペクトル分布を示す図である。
【図３３】第７の検証にて用いた動画像を示す図である。
【図３４】第７の検証において３次元ＦＦＴを用いて求めた時空間スペクトル分布を示す図である。
【図３５】第７の検証において本周波数解析手法を用いて求めた時空間スペクトル分布を示す図である。
【図３６】図３４に示す時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【図３７】図３５に示す時空間スペクトルに基づいて推定した平面を描画した様子を示す図である。
【発明を実施するための形態】
【００２１】
以下、本発明を適用した具体的な実施の形態について図面を参照しながら詳細に説明する。
【００２２】
この実施の形態は、動画像内の１つ又は複数の移動物体の動きを推定する動き推定装置である。特に、この動き推定装置は、非線形方程式を解くことでフーリエ係数を推定することによって周波数分解能が解析窓長に依存しない新たな周波数解析手法を利用して動画像の３次元信号解析を行うものである。
【００２３】
［動き推定装置の構成］
動き推定装置は、例えばコンピュータ等から構成され、図１に示すように、各部を統括的に制御するＣＰＵ（Central Processing Unit）１１と、各種プログラムを含む各種情報を格納する読み取り専用のＲＯＭ（Read Only Memory）１２と、ワークエリアとして機能するＲＡＭ（Random
Access Memory）１３と、各種情報を読み出し及び／又は書き込み可能に記憶する記憶部１４と、ユーザインターフェースとしての図示しない所定の操作デバイスを介した入力操作の処理及び制御を行う入力操作制御部１５と、各種情報を表示する表示部１６とを備える。
【００２４】
ＣＰＵ１１は、記憶部１４等に格納されている各種アプリケーションプログラムをはじめとする各種プログラムを実行し、各部を統括的に制御する。
【００２５】
ＲＯＭ１２は、各種プログラムをはじめとする各種情報を格納している。このＲＯＭ１２に格納されている情報は、ＣＰＵ１１の制御のもとに読み出される。
【００２６】
ＲＡＭ１３は、ＣＰＵ１１が各種プログラムを実行する際のワークエリアとして機能し、ＣＰＵ１１の制御のもとに、各種情報を一時記憶するとともに、記憶している各種情報を読み出す。
【００２７】
記憶部１４は、本発明にかかる動き推定プログラム等のアプリケーションプログラムの他、解析対象となる動画像データをはじめとする各種情報を記憶する。この記憶部１４としては、例えば、ハードディスクや不揮発性メモリ等を用いることができる。また、記憶部１４には、本体に対して着脱可能とされるフレキシブルディスクやメモリカード等の記憶媒体に対して、各種情報の読み出し及び／又は書き込みを行うドライブ装置も含まれる。この記憶部１４に記憶されている各種情報は、ＣＰＵ１１の制御のもとに読み出される。
【００２８】
入力操作制御部１５は、例えば、キーボード、マウス、キーパッド、赤外線リモートコントローラ、スティックキー、又はプッシュボタンといった、ユーザインターフェースとしての図示しない所定の操作デバイスを介した入力操作を受け付け、操作内容を示す制御信号をＣＰＵ１１に対して供給する。
【００２９】
表示部１６は、例えば、液晶ディスプレイ（Liquid Crystal
Display；ＬＣＤ）、プラズマ・ディスプレイ・パネル（Plasma Display Panel；ＰＤＰ）、有機エレクトロルミネッセンス（Organic ElectroLuminescent）ディスプレイ、又はＣＲＴ（Cathode Ray Tube）といった、各種表示デバイスであり、ＣＰＵ１１の制御のもとに各種情報を表示する。例えば、表示部１６は、ＣＰＵ１１によって動き推定プログラムが起動されると、その画面を表示し、入力された解析対象としての動画像データや動き推定結果等を表示する。
【００３０】
このような各部を備える動き推定装置は、ＣＰＵ１１の制御のもとに、動き推定プログラムを実行すると、ＣＰＵ１１の制御のもとに、入力された動画像データの周波数解析を行うことによって時空間スペクトルを抽出し、得られた時空間スペクトルについて最小二乗法を用いたクラスタリングを行うことにより、時空間スペクトル分布を、平面をなすスペクトル群毎に切り分けして１つ又は複数の平面群を推定する。なお、周波数解析の対象となる信号、すなわち、解析対象となる動画像データは、図示しない動画像入力部を介してＣＰＵ１１に入力される。例えば、動き推定装置は、テレビ映像を録画して得られた動画像データの圧縮符号化を行う場合には、当該動き推定装置を含む録画機とテレビとを接続する所定のインターフェースを介して解析対象としての動画像データを入力する。すなわち、動画像入力部は、解析対象としての動画像データをＣＰＵ１１に入力させる機能を有する部位である。なお、動画像入力部は、アナログ信号を入力した場合には、Ａ／Ｄ変換を行ってディジタル信号に変換する機能をあわせ持つことはいうまでもない。このとき、動画像入力部は、必要に応じてアンチエイリアシングフィルタを含むＡ／Ｄ変換器としてもよい。動き推定装置は、ＣＰＵ１１の制御のもとに、このようにして入力された解析対象としての動画像データの周波数解析を行うことによって動き推定を行い、圧縮符号化したデータ等を、図示しない出力部を介して記憶部１４に記憶させたり、その他の機器に出力したりする。
【００３１】
［周波数解析アルゴリズム］
まず、動き推定装置における一連の動き推定アルゴリズムの説明に先立って、動き推定を行う際に利用する周波数解析アルゴリズムについて詳述する。なお、動画像データは３次元信号であることはいうまでもないが、ここでは、説明の便宜上、１次元の解析対象信号に対する周波数解析アルゴリズムについて説明するものとする。
【００３２】
動き推定装置に適用する周波数解析手法（以下、本周波数解析手法という。）においては、次式（１）に示す非周期信号のフーリエ変換式の周波数パラメータを求める問題を非線形方程式の最適解を求める問題に置き換えている。
【００３３】
【数１】
0005467346.tif 000002
【００３４】
具体的には、本周波数解析手法においては、次式（２）に示すように、解析対象信号ｘ（ｎ）と正弦波モデル信号との差の二乗和で表される非線形方程式の最適解として、この非線形方程式の右辺が最小値になるような周波数ｆ’、振幅Ａ’、及び初期位相φ’を求める。なお、次式（２）において、Ｌはフレーム長（解析窓長）であり、ｆｓはサンプリング周波数［Ｈｚ］である。本周波数解析手法においては、このような最小二乗法によって非線形方程式の最適解を求める問題に帰着させることにより、解析窓の影響やエイリアシングの影響がなくなり、解析窓長が、１周期未満であってもよく、周期の整数倍でなくてもよく、さらには、不等間隔であってもよい等、柔軟な周波数解析処理を実現することが可能となる。
【００３５】
【数２】
0005467346.tif 000003
【００３６】
さて、上式（２）に示す非線形方程式の最適解を実際に求めるにあたっては、以下のような方法をとることができる。
【００３７】
本周波数解析手法においては、振幅Ａ’、周波数ｆ’、及び初期位相φ’のそれぞれについて適切な初期値を求め、これら初期値から非線形方程式の解法を用いて最適解に収束させる。この非線形問題では、上式（２）をコスト関数とする最小化問題とする。なお、適切な初期値は、離散フーリエ変換（Discrete Fourier Transform；ＤＦＴ）やウェーブレット変換等の任意の周波数変換を行ったり、フィルタリングを行うことによっておおよその見当をつけたりする等、既存の任意の方法を適用して求めることができる。
【００３８】
まず、本周波数解析手法においては、上式（２）における正弦波モデル信号の位相を構成する周波数パラメータｆ’，φ’について、いわゆる最急降下法を適用し、周波数パラメータｆｍ’，φｍ’を次式（３）及び次式（４）によって求める。
【００３９】
【数３】
0005467346.tif 000004
【数４】
0005467346.tif 000005
【００４０】
なお、上式（３）及び上式（４）においては、次式（５）と略している。また、μｍは、いわゆる減速法に基づく重み係数であり、各漸化式によって求められるコスト関数を単調減少数列にするために、適時０〜１の値をとる。
【００４１】
【数５】
0005467346.tif 000006
【００４２】
周波数パラメータｆｍ’，φｍ’を求めることができれば、上式（２）における正弦波モデル信号の係数としての周波数パラメータＡ’を一意に求めることができるため、本周波数解析手法においては、次式（６）によって周波数パラメータＡｍ’を収束させる。
【００４３】
【数６】
0005467346.tif 000007
【００４４】
本周波数解析手法においては、これら一連の計算を反復して行うことにより、振幅Ａ’、周波数ｆ’、及び初期位相φ’を高精度に収束させることができる。特に、本周波数解析手法においては、上式（２）における正弦波モデル信号の位相を構成する周波数パラメータｆ’，φ’と、係数としての周波数パラメータＡ’とを別個に求めることにより、計算を簡便に行うことができる。
【００４５】
しかしながら、最急降下法は、比較的広い範囲から収束するものの、１回の反復では精度が低く、収束するまでに時間を要する。
【００４６】
そこで、本周波数解析手法においては、最急降下法を適用して周波数パラメータｆｍ’，φｍ’をある程度まで収束させた後、さらに、いわゆるニュートン法を適用して高精度に収束させるのが望ましい。具体的には、本周波数解析手法においては、ニュートン法として、次式（７）及び次式（８）に示す漸化式によって周波数パラメータｆｍ’，φｍ’を求める。
【００４７】
【数７】
0005467346.tif 000008
【数８】
0005467346.tif 000009
【００４８】
ただし、上式（７）及び上式（８）において、Ｊは次式（９）とし、次式（１０）と略している。また、νｍもμｍと同様に減速法に基づく重み係数であり、適時０〜１の値をとる。
【００４９】
【数９】
0005467346.tif 000010
【数１０】
0005467346.tif 000011
【００５０】
本周波数解析手法においては、上式（７）及び上式（８）によって周波数パラメータｆｍ’，φｍ’を求めた後、最急降下法と同様に、上式（６）によって周波数パラメータＡｍ’を収束させ、この一連の計算をさらに反復して行う。
【００５１】
このように、本周波数解析手法においては、最急降下法とニュートン法とを組み合わせたハイブリッド型の解法を用いることにより、高速に且つ高精度に周波数パラメータＡ’，ｆ’，φ’を推定することができる。
【００５２】
また、本周波数解析手法においては、解析対象信号ｘ（ｎ）が複合正弦波の場合であっても、逐次減算処理することにより、近似的にスペクトルパラメータを導出することができる。ここで、解析対象信号ｘ（ｎ）が複数の正弦波の和であり、次式（１１）のように表されているとする。
【００５３】
【数１１】
0005467346.tif 000012
【００５４】
パーセヴァル（Parseval）の定理より、解析対象信号ｘ（ｎ）の周波数ｆｋと正弦波モデル信号の周波数パラメータｆ’とが全く一致しない場合、すなわち、次式（１２）である場合には、上式（２）に示す非線形方程式は次式（１３）となる。また、周波数パラメータｆ’，φ’の組が、周波数ｆｋ及び初期位相φｋの組のいずれかに一致する場合には、上式（２）に示す非線形方程式は次式（１４）となる。さらに、振幅Ａｊが周波数パラメータＡ’とも一致した場合には、解析対象信号から推定スペクトルに関する周波数成分を完全に消去することができる。そのため、最適解を求める問題は、周波数に対して独立であり、解析対象信号から順次個別に推定すれば、複数の正弦波で表される信号にも応用することができる。
【００５５】
【数１２】
0005467346.tif 000013
【数１３】
0005467346.tif 000014
【数１４】
0005467346.tif 000015
【００５６】
すなわち、本周波数解析手法においては、解析対象信号ｘ（ｎ）が複合正弦波の場合であっても、逐次残差信号に対して同様に処理を行い、複数の正弦波を抽出することができる。
【００５７】
音声信号や音響信号等の信号を複合正弦波によって表現するためには、これまで多くのスペクトル数（正弦波の数）が必要であったが、本周波数解析手法においては、そのような信号であっても僅かなスペクトル数で誤差なく表現することができる。すなわち、信号をより少ないスペクトル数で表現可能であることは、情報圧縮の用途に有効であることを示している。
【００５８】
［本周波数解析手法の有効性］
以下、本周波数解析手法の有効性について具体的に説明する。
【００５９】
本周波数解析手法は、非線形方程式の最適解を求めることにより、正弦波モデル信号の周波数ｆ’、振幅Ａ’、及び初期位相φ’を高速に且つ高精度に求めることができる。具体的な精度を立証するために、本願発明者は、ＤＦＴと、ＤＦＴの発展型のうち最も解析精度が高いといわれているＧＨＡ（Generalized Harmonic Analysis）とを比較対象として精度の検証を行った。
【００６０】
なお、ＤＦＴやＧＨＡは、１つの解析窓長に見かけ上複数の窓長を持たせていることから、周波数分解能が解析窓長に依存するが、その分解周波数が有限長であり、解析対象信号の周波数が分解周波数以外の周波数となった場合には解析することができず、解析対象信号が正確に解析できる周波数と異なる場合には、最も近い分解周波数の他に、その周辺に小さなスペクトルの周波数（側帯波成分）が現れ、複数の周波数が出現してしまう。
【００６１】
このような現象が本周波数解析手法においても生じるか否かについて、すなわち、本周波数解析手法の周波数分解能を検証するために、解析窓長を１秒（１０２４サンプル）とした１次元の非常に短い単一正弦波を解析し、各手法によって正弦波を１本抽出して元の信号との二乗誤差を調べた。その結果を図２に示す。
【００６２】
図２に示すように、ＤＦＴにおいては、基本周波数の整数倍以外の周波数における解析精度の悪化がみられた。また、ＧＨＡにおいては、１Ｈｚ以上の周波数ではＤＦＴと比べて２〜５桁程度の精度向上がみられた。これに対して、本周波数解析手法においては、１Ｈｚ以上の周波数ではＤＦＴと比べて１０桁以上、ＧＨＡと比べて５桁以上の精度向上がみられた。すなわち、本周波数解析手法は、既存の周波数解析手法と比べて１０万〜１００億倍以上の精度向上がみられた。特に、１Ｈｚ以下の周波数を正確に推定することができるということは、解析窓長を超えた長い周期信号であっても解析可能であることを示している。
【００６３】
このように、本周波数解析手法は、最も解析精度が高いといわれているＧＨＡと比べても驚くべき高精度に解析を行うことができるものである。動き推定装置は、このような周波数解析手法を利用して、図３に示すような一連の処理を行うことにより、入力された動画像内の１つ又は複数の移動物体の動きを推定する。
【００６４】
［動画像からの領域切り出し］
まず、動き推定装置は、図３に示すように、ステップＳ１において、ＣＰＵ１１の制御のもとに、図示しない動画像入力部を介して入力されてＲＡＭ１３等のメモリに記憶されているオリジナルの動画像データから移動物体が含まれる領域を切り出す。ここで、切り出した領域の横方向のピクセル数、縦方向のピクセル数、フレーム数を、それぞれ、Ｎｘ，Ｎｙ，Ｎｚとし、各要素をｉ，ｊ，ｋとする。なお、ｉ，ｊ，ｋは、画像を再現するためのインデックスであり、実数である。また、切り出した３次元信号からなる動画像データを、ξｏｂｊ（ｉ，ｊ，ｋ）とする。この切り出された動画像データは、ＲＡＭ１３等に記憶される。
【００６５】
［時空間スペクトルの抽出］
続いて、動き推定装置は、ステップＳ２において、ＣＰＵ１１の制御のもとに、次式（１５）に示すように、上述した周波数解析手法を利用して３次元信号からなる動画像データξｏｂｊ（ｉ，ｊ，ｋ）からＬ本の時空間スペクトルを抽出する。次式（１５）において、ｆｘｓ，ｆｙｓ，ｆｚｓは、それぞれ、動画像の横軸方向、縦軸方向、時間軸方向のサンプリング周波数［Ｈｚ］であり、Ａ’，ｆｘ’，ｆｙ’，ｆｚ’，φ’は、それぞれ、実際に抽出するスペクトルの振幅、各軸に対応する周波数、初期位相である。動き推定装置は、３次元信号である動画像データを次式（１５）に示すような正弦波モデル関数を用いて表し、実際の信号と、次式（１５）に示す正弦波モデル信号との差が最小となるように、次式（１６）を用いてパラメータを変化させ、各周波数を求める。
【００６６】
【数１５】
0005467346.tif 000016
【数１６】
0005467346.tif 000017
【００６７】
そして、動画像データは、次式（１７）に示すように、時空間スペクトルの合成によって表現することができるため、動き推定装置は、ＣＰＵ１１の制御のもとに、動画像データからＬ本の時空間スペクトルを抽出する。このとき、ｌ番目の時空間スペクトルを次式（１８）に示すように表す。
【００６８】
【数１７】
0005467346.tif 000018
【数１８】
0005467346.tif 000019
【００６９】
［時空間スペクトルのクラスタリング］
そして、動き推定装置は、ステップＳ３において、ＣＰＵ１１の制御のもとに、時空間スペクトルのクラスタリングを行う。ここでは、Ｉ個の平面が存在すると仮定し、ｌ番目の時空間スペクトルがｉ番目の平面に属するものとした場合の存在確率変数をＰｉ（ｌ）とする。なお、ここでは、確率変数Ｐは、０又は１の値をとるものとし、当初は０で初期化されているものとする。
【００７０】
まず、動き推定装置は、ＣＰＵ１１の制御のもとに、時空間スペクトルＳｌから２本の時空間スペクトルＳｕ，Ｓｖを選択する。なお、ｕ，ｖは、次式（１９）に示すように集合Ｌに含まれ、且つ、互いに異なる値である。
【００７１】
【数１９】
0005467346.tif 000020
【００７２】
続いて、動き推定装置は、ＣＰＵ１１の制御のもとに、次式（２０）に示すように、選択した時空間スペクトルＳｕ，Ｓｖと原点とによって構成される平面の法線ベクトルｎｕ，ｖを求める。
【００７３】
【数２０】
0005467346.tif 000021
【００７４】
続いて、動き推定装置は、ＣＰＵ１１の制御のもとに、次式（２１）に示すように、求めた法線ベクトルｎｕ，ｖを用いて平面の方程式を仮定する
【００７５】
【数２１】
0005467346.tif 000022
【００７６】
続いて、動き推定装置は、ＣＰＵ１１の制御のもとに、次式（２２）に示すように、上式（２１）によって求めた平面と、時空間スペクトルＳｕ，Ｓｖを除く時空間スペクトルＳｍとの距離Ｄｕ，ｖ（ｍ）をそれぞれ求める。なお、ｍは、集合Ｌに含まれ、且つ、ｕ，ｖとは異なる値である。
【００７７】
【数２２】
0005467346.tif 000023
【００７８】
そして、動き推定装置は、ＣＰＵ１１の制御のもとに、次式（２３）に示すように、求めた距離Ｄｕ，ｖ（ｍ）が所定の閾値ε以下となる時空間スペクトルを探索し、その確率変数Ｐｉ（ｍ）を１とする。
【００７９】
【数２３】
0005467346.tif 000024
【００８０】
なお、ここでは、確率変数を１とすることにより、平面に完全に属するものとして演算を行っているが、動き推定装置においては、小数の確率変数を設定し、平面に近い点とみなす閾値を設け、この閾値を用いた演算を行うようにしてもよく、この閾値を次回の演算の際の算出パラメータの１つとしてもよい。
【００８１】
そして、動き推定装置は、ＣＰＵ１１の制御のもとに、上式（２３）を満たす時空間スペクトルＳｍに対して最小二乗法を適用し、次式（２４）に示す平面の方程式を計算し、推定した平面の傾きから移動物体の横方向及び縦方向の速度Ｖｘ，Ｖｙを推定する。なお、次式（４）におけるｄは、ｆｚ軸の切片である。
【００８２】
【数２４】
0005467346.tif 000025
【００８３】
［時空間スペクトルからの平面推定］
動き推定装置は、ステップＳ４において、ＣＰＵ１１の制御のもとに、以上のようにして時空間スペクトルから平面を推定し、平面ｉに属していない時空間スペクトルが存在する場合には、上式（２３）を満たしていない時空間スペクトルを用いて、ステップＳ３における時空間スペクトルのクラスタリングを繰り返す。
【００８４】
動き推定装置は、このような一連の処理を行うことにより、入力された動画像内の１つ又は複数の移動物体の動きを推定することができる。特に、動き推定装置は、時空間スペクトルの抽出に上述した本周波数解析手法を利用することにより、理論平面からの誤差が少ない平面の推定を行うことができる。そのため、動き推定装置は、直線群フィルタやファジー推論を用いたクラスタリング法ではなく、平面推定の簡易的な手法である最小二乗法を用いたクラスタリング法を適用することが可能となり、高精度に且つ少ない演算量で動き推定を行うことができる。
【００８５】
［動き推定装置の有効性］
以下、動き推定装置の有効性について具体的に説明する。
［第１の検証］
まず、図４に示すように、横方向のピクセル数、縦方向のピクセル数、フレーム数を、それぞれ、Ｎｘ，Ｎｙ，Ｎｚとし、初期フレームの左端にある黒い縦棒状の物体が、時間の進行にともない右方に移動するような動画像を考える。移動速度Ｖｘは、連続するフレーム間を物体が移動するピクセル数として定義し、フレームの右方又は上方に物体が移動する場合を正とする。速度がフレームあたりのピクセル数であり、Ｎｚ＞（Ｎｘ／Ｖｘ）である場合には、物体は、ｎ番目のフレームまで存在する。縦棒状の物体を用いていることから、縦方向の周波数ｆｙはゼロである。したがって、動画像は、周波数ｆｘ，ｆｚによって特徴付けられることになる。本願発明者は、３次元の平面ｆｘ−ｆｚ上のスペクトル構造を求め、振幅及び周波数ｆｘ，ｆｚの変化を観察した。この場合、時空間スペクトルは、平面ｆｘ−ｆｚ上で直線状に分布し、次式（２５）に示すように定義される。なお、次式（２５）において、ｄは、ｆｚ軸の切片であり、速度Ｖｘは、理論移動の一次方程式における周波数ｆｘの係数である。この直線は、理論直線と称される。
【００８６】
【数２５】
0005467346.tif 000026
【００８７】
ここで、３次元ＦＦＴによる周波数解析を実行した場合、移動物体の速度及び解析窓は、スペクトルの分布が理論直線の形態をとるかどうかに影響する。また、３次元ＦＦＴの分解能も解析窓に依存する。Ｎｚ＝Ｎｘ／Ｖｘが一定の場合には、移動物体は、動画像を構成する全てのフレームにおいて存在し、動画像データは略周期的なものとなる。また、分解能及びフレーム数が整数であることから、速度もまた整数となる。したがって、この場合には、３次元ＦＦＴにおいて高エネルギーを有するスペクトルは理論直線上にあることになる。しかしながら、Ｎｚ＜Ｎｘ／Ｖｘ又はＮｚ＞Ｎｘ／Ｖｘである場合には、解析データは、横方向ｘ及び時間軸方向ｚについて値がゼロであるデータと等しくなる。このような場合には、解析窓内でデータが周期的でないことから、スペクトルは、理論直線上の複数個所で隆起したものとなる。
【００８８】
実際に、図４に示した動画像として、サイズを（Ｎｘ×Ｎｙ×Ｎｚ）＝（６４×６４×６４）とし、物体の大きさを６４ピクセル×８ピクセルとした場合において、物体の移動速度Ｖｘを１．０ピクセル／フレームとした場合と、３．３ピクセル／フレームとした場合とについて時空間スペクトルの分布を求めた。なお、スペクトルの振幅は、全てのスペクトルが０ｄＢよりも大きくなるように最大値を４０ｄＢに設定した。また、初期位相を１に、周波数ｆｘ，ｆｙ，ｆｚをそれぞれ６４Ｈｚに正規化した。
【００８９】
図５に、物体の移動速度Ｖｘを１．０ピクセル／フレームとした場合における３次元ＦＦＴの平面ｆｘ−ｆｚ上の時空間スペクトルを示す。図５における破線は、理論直線であり、この直線を正確に求めることができれば正確な移動速度を得ることができる。図５においては、３次元ＦＦＴの分解能が１であり、物体の移動速度に一致していることから、高エネルギーを有するスペクトルが直線を形成している。しかしながら、物体の移動速度Ｖｘが、３次元ＦＦＴの分解能に等しくない３．３ピクセル／フレームである場合には、図６に示すように、スペクトルが理論直線上のみに存在するのではなくなり、その周囲にも存在するものとなる。したがって、正確な速度を得るために、理論直線上にないスペクトルのエネルギーを考慮することが重要である。
【００９０】
これに対して、物体の移動速度Ｖｘを３．３ピクセル／フレームとした場合において、上述した本周波数解析手法を用いて求めた時空間スペクトルは、図７に示すように、理論直線上に存在し、その周囲には存在しないものとなった。なお、図示しないが、物体の移動速度Ｖｘが１．０ピクセル／フレームである場合にも、同様に、時空間スペクトルが理論直線上に存在する分布が得られた。
【００９１】
これらの結果を定量的に評価するために、本願発明者は、３次元ＦＦＴ及び本周波数解析手法のそれぞれについての時空間スペクトルと理論直線との分散を求めた。分散σ２は、次式（２６）によって求めた。なお、次式（２６）において、ｆｘ’，ｆｚ’は、それぞれ、３次元ＦＦＴ又は本周波数解析手法のスペクトルパラメータであり、Ｌはスペクトルの本数である。
【００９２】
【数２６】
0005467346.tif 000027
【００９３】
この分散σ２が小さいほど、理論直線に対する誤差が小さくなる。３次元ＦＦＴの場合には、分散σ２は５６０．０となり、本周波数解析手法の場合には、０．０２６８と極めて小さい値となった。本願発明者は、振幅が約４０ｄＢから６０ｄＢのスペクトルについて、３次元ＦＦＴの場合には分散σ２が約５７０．０に収束する一方で、本周波数解析手法の場合には分散σ２が約０．０３に収束することを確認している。
【００９４】
したがって、本周波数解析手法によるスペクトルは、ほぼ理論直線上に存在することがわかる。移動物体の周期が解析窓の長さと等しい場合には、３次元ＦＦＴによるスペクトルは直線を形成するが、逆に、そうでない場合には、スペクトルは直線の形態をとらない。本周波数解析手法を適用した動き推定装置は、解析窓の大きさにかかわらず常に正確なスペクトルを得ることができ、高精度に動き推定を行うことができることがわかる。
【００９５】
［第２の検証］
第１の検証においては、物体が１方向にのみ移動する場合を示したが、２方向に移動する場合の例として、図８に示すように、初期フレームの左下にある黒い正方形状の物体が、時間の進行にともない縦横方向に均一な速度で右上に移動するような動画像を考えた。このような物体の移動は、図４に示したような物体の移動において、フレームの下方から上方に移動する横棒状物体の重ね合わせによって表現される。横方向の移動速度をＶｘとし、縦方向の移動速度をＶｙとした場合、３次元ＦＦＴにおいて高エネルギーを有する時空間スペクトルは、３次元周波数領域空間における平面上の共役スペクトルを除いてはこの平面上にあることになる。この平面は、次式（２７）のような理論平面として定義される。したがって、２つの速度Ｖｘ，Ｖｙは、理論平面方程式における周波数ｆｘ，ｆｚの係数を計算することによって求められる。
【００９６】
【数２７】
0005467346.tif 000028
【００９７】
ここで、ＮｚがＮｘ／Ｖｘ，Ｎｙ／Ｖｙの最小公倍数と等しい場合には、動画像データは略周期的なものとなり、３次元ＦＦＴによるスペクトルの大部分は、理論平面の近くに存在することになる。しかしながら、ＮｚがＮｘ／Ｖｘ，Ｎｙ／Ｖｙの最小公倍数と等しくない場合には、時空間スペクトルは、理論平面に限らずどこにでも存在すると考えられる。したがって、３次元ＦＦＴによる時空間スペクトルは、解析窓の影響を受けると推測される。本願発明者は、上述した本周波数解析手法を用い、本周波数解析手法の分解能をテストするために速度を整数に設定することにより、この問題の解決を試みた。
【００９８】
実際に、図８に示した動画像として、サイズを（Ｎｘ×Ｎｙ×Ｎｚ）＝（６４×６４×６４）とし、物体の大きさを８ピクセル×８ピクセルとした場合において、物体の移動速度をＶｘ＝Ｖｙ＝３．３ピクセル／フレームとした場合について時空間スペクトルの分布を求めた。なお、３次元周波数領域におけるスペクトル分布の観察を明確にするために、スペクトルの振幅は、全てのスペクトルが０ｄＢよりも大きくなるように最大値を２０ｄＢに設定した。
【００９９】
図９及び図１０に、それぞれ、３次元周波数領域において３次元ＦＦＴ及び本周波数解析手法を用いて求めた時空間スペクトルを示す。図９及び図１０において、３つの破線によって形成される平面が理論平面であり、上式（２７）におけるパラメータＶｘ，Ｖｙ、すなわち、物体の移動速度は３．３である。これら図９及び図１０は、求めた時空間スペクトルと理論平面との関係を示している。
【０１００】
３次元ＦＦＴによる時空間スペクトルの一部は、理論平面上に存在するが、大部分のスペクトルは、理論平面の近くに存在している。一方、本周波数解析手法による時空間スペクトルは、その大部分が理論平面上に存在している。分散σ２は、次式（２８）によって求めることができる。
【０１０１】
【数２８】
0005467346.tif 000029
【０１０２】
３次元ＦＦＴの場合には、分散σ２は１１．６６となり、本周波数解析手法の場合には、０．０１と極めて小さい値となった。また、特に図示しないが、物体の大きさが３２ピクセル×３２ピクセルと大きくなった場合には、分散σ２は、対数的に増加する結果が得られた。
【０１０３】
このように、本周波数解析手法を適用した動き推定装置は、２方向の速度を有する物体の移動についての正確な時空間スペクトルを得るためにも有益であることがわかる。
【０１０４】
［第３の検証］
つぎに、同時に移動する２つの物体を含む動画像の例として、図１１に示すように、初期フレームの左下及び右下にある黒い正方形状の２つの物体ａ，ｂが、それぞれ、時間の進行にともない、分解能に一致していない異なる速度で右上及び左上に移動するような動画像を考えた。なお、動画像のサイズは、（Ｎｘ×Ｎｙ×Ｎｚ）＝（６４×６４×６４）であり、物体の大きさは、８ピクセル×８ピクセルである。
【０１０５】
図１２及び図１３に、それぞれ、３次元周波数領域において３次元ＦＦＴ及び本周波数解析手法を用いて求めた時空間スペクトルを示す。
【０１０６】
図１２から、３次元ＦＦＴによる時空間スペクトルは、分解能が不十分であり、スペクトルを正確に抽出することができないため、分解能に一致していない速度で移動している物体を含む動画像の場合には、スペクトルのエネルギーが理想とする場所の周囲に漏れてしまい、平面が積層して厚みを有する層のようなスペクトル分布となることがわかる。このような層状のスペクトル分布が得られた場合には、これらのスペクトルを分離して平面を推定することが困難となる。すなわち、ＤＦＴやＦＦＴに代表される等間隔の周波数分解能を有する解析手法においては、その等間隔の周波数に一致する周波数成分であれば正確なスペクトルを得ることができるものの、そうでない場合には正確なスペクトルを得ることができず、理論平面の周囲の周波数に群れたスペクトルが出現し、層状のスペクトル分布となってしまう。
【０１０７】
一方、本周波数解析手法による時空間スペクトルは、図１３に示すように、２つの物体に対応する２つの平面状に分布したものとなり、ばらつきもほとんどみられない。このような２つの平面は、上述したように、最小二乗法を用いて容易に推定することができる。図１４に、２つの平面がよく観察できるように視点を変えて平面を描画した様子を示している。この図から、２つの平面が高精度に推定されていることがわかる。動き推定装置においては、推定した平面のｆｘ軸及びｆｙ軸の傾きをそれぞれ求めることにより、物体の横方向及び縦方向の速度を推定することができる。
【０１０８】
なお、上述した３次元ＦＦＴによる時空間スペクトルにおける層状のスペクトル分布は、特に類似する複数の動きを含む動画像の場合に顕著にみられる。本願発明者は、このような類似する複数の動きを含む動画像における本周波数解析手法の有効性を確認する実験も行っている。
【０１０９】
動画像のサイズを（Ｎｘ×Ｎｙ×Ｎｚ）＝（６４×６４×６４）とし、物体の大きさを８ピクセル×８ピクセルとした場合において、先に図８に示したように、初期フレームにおいて左下にある２つの黒い正方形状の物体が、時間の進行にともない縦横方向に均一な速度で右上に移動するような動画像を対象とした。２つの物体をａ，ｂとすると、物体ａの移動速度は、Ｖｘ＝１．０ピクセル／フレーム、Ｖｙ＝１．７ピクセル／フレームであり、物体ｂの移動速度は、Ｖｘ＝１．０ピクセル／フレーム、Ｖｙ＝２．３ピクセル／フレームである。すなわち、これら２つの物体ａ，ｂは、同じ場所を始点とし、時間の進行とともに、右上方向に移動しながら徐々に高さ方向に離間していく動きをするものである。なお、スペクトルの振幅は、全てのスペクトルが０ｄＢよりも大きくなるように最大値を２８ｄＢに設定した。
【０１１０】
３次元ＦＦＴ及び本周波数解析手法の双方とも２００本のスペクトルを用いて分布を求めたところ、それぞれ、図１５及び図１６に示すような結果が得られた。そして、この結果から推定した平面の傾き、すなわち、物体ａ，ｂのそれぞれの移動速度を求めると、３次元ＦＦＴの場合には、
物体ａについて、Ｖｘ＝１．０、Ｖｙ＝１．８５４５
物体ｂについて、Ｖｘ＝１．０、Ｖｙ＝１．５４８４
が得られたのに対して、本周波数解析手法の場合には、
物体ａについて、Ｖｘ＝０．９９９９４、Ｖｙ＝１．７２１１
物体ｂについて、Ｖｘ＝０．９８７１２、Ｖｙ＝２．３１０８
が得られた。すなわち、３次元ＦＦＴを用いた場合には、縦方向の移動速度Ｖｙの精度が悪く、本周波数解析手法を用いた場合には極めて高精度に動きを推定していることがわかる。
【０１１１】
同様に、物体ａの移動速度を、Ｖｘ＝１．０ピクセル／フレーム、Ｖｙ＝２．３ピクセル／フレームに変化させ、物体ｂの移動速度を、Ｖｘ＝１．７ピクセル／フレーム、Ｖｙ＝３．２ピクセル／フレームに変化させて同様にスペクトル分布を求めた。ここでも、スペクトルの振幅は、全てのスペクトルが０ｄＢよりも大きくなるように最大値を２８ｄＢに設定した。
【０１１２】
３次元ＦＦＴの場合には３００本のスペクトル、本周波数解析手法の場合には２００本のスペクトルを用いて分布を求めたところ、それぞれ、図１７及び図１８に示すような結果が得られた。そして、この結果から推定した平面の傾き、すなわち、物体ａ，ｂのそれぞれの移動速度を求めると、３次元ＦＦＴの場合には、
物体ａについて、Ｖｘ＝１．００３４、Ｖｙ＝２．１５２５
物体ｂについて、Ｖｘ＝１．８０５３、Ｖｙ＝３．５６９３
が得られたのに対して、本周波数解析手法の場合には、
物体ａについて、Ｖｘ＝１．０１０６、Ｖｙ＝２．３１１６
物体ｂについて、Ｖｘ＝１．７５０２、Ｖｙ＝３．２３０４
が得られ、この場合にも、３次元ＦＦＴを用いた場合に比べ、本周波数解析手法を用いた場合の方が良好な結果が得られた。
【０１１３】
このように、本周波数解析手法を適用した動き推定装置は、類似する複数の動きを含む動画像における動き推定を行う場合にも極めて有効である。
【０１１４】
［第４の検証］
つぎに、同時に移動する３つの物体を含む動画像の例として、図１９に示すように、初期フレームの右上、左下及び右下にある黒い正方形状の３つの物体ａ，ｂ，ｃが、それぞれ、時間の進行にともない、異なる速度で左下、右上及び左上に移動するような動画像を考えた。なお、動画像のサイズは、（Ｎｘ×Ｎｙ×Ｎｚ）＝（６４×６４×６４）であり、物体の大きさは、８ピクセル×８ピクセルである。また、物体ａの移動速度は、Ｖｘ＝−１．０ピクセル／フレーム、Ｖｙ＝−１．０ピクセル／フレームであり、物体ｂの移動速度は、Ｖｘ＝１．３ピクセル／フレーム、Ｖｙ＝２．７ピクセル／フレームであり、物体ｃの移動速度は、Ｖｘ＝−３．３ピクセル／フレーム、Ｖｙ＝２．５ピクセル／フレームである。このような動画像について、本周波数解析手法を用いて時空間スペクトルを求め、使用するスペクトルの本数を変えながら平面を推定した。
【０１１５】
まず、１００本のスペクトルを用いた場合には、図２０に示すような結果が得られ、推定した３つの平面の傾き、すなわち、物体ａ，ｂ，ｃのそれぞれの移動速度は、
物体ａについて、Ｖｘ＝−１．００３８、Ｖｙ＝−１．００１６
物体ｂについて、Ｖｘ＝１．２３５３、Ｖｙ＝２．６９３３
が得られたが、物体ｃの移動速度については、スペクトルの本数が少なく求めることができなかった。
【０１１６】
また、１２５本のスペクトルを用いた場合には、図２１に示すような結果が得られ、物体ａ，ｂ，ｃのそれぞれの移動速度は、
物体ａについて、Ｖｘ＝−１．００３８、Ｖｙ＝−１．００２４
物体ｂについて、Ｖｘ＝１．２５９３、Ｖｙ＝２．７１６６
物体ｃについて、Ｖｘ＝−３．２３３２、Ｖｙ＝２．４５２
が得られた。
【０１１７】
さらに、１６０本のスペクトルを用いた場合には、図２２に示すような結果が得られ、物体ａ，ｂ，ｃのそれぞれの移動速度は、
物体ａについて、Ｖｘ＝−１．００５３、Ｖｙ＝−１．００１５
物体ｂについて、Ｖｘ＝１．２９４８、Ｖｙ＝２．７０１２
物体ｃについて、Ｖｘ＝−３．２８９７、Ｖｙ＝２．４５９４
が得られた。
【０１１８】
さらにまた、１８０本のスペクトルを用いた場合には、図２３に示すような結果が得られ、物体ａ，ｂ，ｃのそれぞれの移動速度は、
物体ａについて、Ｖｘ＝−１．００３、Ｖｙ＝−１．００１３
物体ｂについて、Ｖｘ＝１．３０８５、Ｖｙ＝２．６９７４
物体ｃについて、Ｖｘ＝−３．２９０６、Ｖｙ＝２．４５６１
が得られた。
【０１１９】
また、２００本のスペクトルを用いた場合には、図２４に示すような結果が得られ、物体ａ，ｂ，ｃのそれぞれの移動速度は、
物体ａについて、Ｖｘ＝−１．００２４、Ｖｙ＝−１．０００８
物体ｂについて、Ｖｘ＝１．３１０４、Ｖｙ＝２．７００７
物体ｃについて、Ｖｘ＝−３．２８７４、Ｖｙ＝２．４５９５
が得られた。
【０１２０】
さらに、２２５本のスペクトルを用いた場合には、図２５に示すような結果が得られ、物体ａ，ｂ，ｃのそれぞれの移動速度は、
物体ａについて、Ｖｘ＝−１．００２３、Ｖｙ＝−１．００１１
物体ｂについて、Ｖｘ＝１．３０９３、Ｖｙ＝２．６９９４
物体ｃについて、Ｖｘ＝−３．２９９８、Ｖｙ＝２．４４９４
が得られた。
【０１２１】
そして、２５０本のスペクトルを用いた場合には、図２６に示すような結果が得られ、物体ａ，ｂ，ｃのそれぞれの移動速度は、
物体ａについて、Ｖｘ＝−１．００１９、Ｖｙ＝−１．００１
物体ｂについて、Ｖｘ＝１．３０９１、Ｖｙ＝２．６８３
物体ｃについて、Ｖｘ＝−３．２７５６、Ｖｙ＝２．４７４７
が得られた。
【０１２２】
すなわち、本周波数解析手法を適用した動き推定装置の精度は、時空間スペクトルの本数に依存し、図１９に示した動画像の場合には、約１８０本程度の時空間スペクトルを用いれば、十分な精度を実現できることがわかる。なお、この十分な精度を実現可能な時空間スペクトルの本数は、移動物体の速度や移動方向等、動画像の内容によって異なる。換言すれば、動き推定装置においては、動画像に応じて適切な時空間スペクトルを用いることにより、演算量を削減しながらも高精度に動き推定を行うことが可能である。
【０１２３】
［第５の検証］
つぎに、移動しない背景がある動画像の前景としての物体が移動する場合の例として、図２７に示すように、初期フレームの左下にある黒い正方形状の物体が時間の進行にともない、右上に移動するような動画像を考えた。なお、動画像のサイズは、（Ｎｘ×Ｎｙ×Ｎｚ）＝（６４×６４×６４）であり、物体の大きさは、８ピクセル×８ピクセルであり、物体の移動速度は、Ｖｘ＝Ｖｙ＝３．３ピクセル／フレームである。また、背景は、所定の繰り返しパターンから構成されている。具体的には、背景は、振幅Ａ＝１２８（０〜２５５）、位相φ＝０、周波数ｆｘ＝８、ｆｙ＝０、ｆｚ＝１のパターンと、振幅Ａ＝１２８（０〜２５５）、位相φ＝０、周波数ｆｘ＝０、ｆｙ＝８、ｆｚ＝１のパターンとから構成されている。このような動画像について、３次元ＦＦＴ及び本周波数解析手法のそれぞれを用いて時空間スペクトルを求め、平面を推定した。なお、スペクトルの振幅は、全てのスペクトルが０ｄＢよりも大きくなるように、３次元ＦＦＴの場合には最大値を２５．０７ｄＢ、本周波数解析手法の場合には２２．７５ｄＢに設定した。
【０１２４】
まず、３次元ＦＦＴによる１０１本の時空間スペクトルの分布を求めると、図２８に示すように、背景のパターンの影響を受けて断続的な分布が得られた。そして、この結果から、上述した閾値εを変化させながら物体の移動速度を求めると、
閾値ε＝０．１の場合、Ｖｘ＝３．４６９７、Ｖｙ＝３．４６９７
閾値ε＝０．２の場合、Ｖｘ＝３．３８２４、Ｖｙ＝３．３８２４
閾値ε＝０．３〜０．５の場合、Ｖｘ＝３．２８５７、Ｖｙ＝３．２８５７
といったように、閾値εにかかわらず誤差が大きい結果が得られた。
【０１２５】
一方、本周波数解析手法による１００本の時空間スペクトルの分布を求めると、図２９に示すように、背景のパターンにかかわらず連続的な分布が得られた。そして、この結果から、上述した閾値εを変化させながら物体の移動速度を求めると、
閾値ε＝０．１の場合、Ｖｘ＝３．３０９３、Ｖｙ＝３．３０３３
閾値ε＝０．２の場合、Ｖｘ＝３．３０６、Ｖｙ＝３．３０５３
閾値ε＝０．３〜０．５の場合、Ｖｘ＝３．３０３７、Ｖｙ＝３．３０３１
といったように、誤差が極めて小さい結果が得られた。
【０１２６】
このように、本周波数解析手法を適用した動き推定装置は、背景の影響を受けにくく、高精度に動き推定を行うことが可能である。
【０１２７】
［第６の検証］
つぎに、カラーの動画像において前景が移動する場合の例として、図３０に示すように、人物が左方から右方に移動するような動画像を考えた。なお、動画像のサイズは、（Ｎｘ×Ｎｙ×Ｎｚ）＝（６４×６４×６４）であり、人物の移動速度は、Ｖｘ＝１．６４１ピクセル／フレームである。また、動画像は、所定の温度分布からなる環境において撮像した熱画像から構成されている。このような動画像について、３次元ＦＦＴ及び本周波数解析手法のそれぞれを用いて時空間スペクトルを求め、平面を推定した。なお、スペクトルの振幅は、全てのスペクトルが０ｄＢよりも大きくなるように、３次元ＦＦＴの場合には最大値を１９．４ｄＢ、本周波数解析手法の場合には２２．５８ｄＢに設定した。
【０１２８】
まず、３次元ＦＦＴによる１０１本の時空間スペクトルの分布を求めると、図３１に示すような結果が得られ、この結果から、上述した閾値εを変化させながら人物の移動速度を求めると、
閾値ε＝０．４の場合、Ｖｘ＝１．４８３８
閾値ε＝０．５の場合、Ｖｘ＝１．５７３２
といったように、閾値εにかかわらず誤差が大きい結果が得られた。
【０１２９】
一方、本周波数解析手法による１００本の時空間スペクトルの分布を求めると、図３２に示すような結果が得られ、この結果から、上述した閾値εを変化させながら人物の移動速度を求めると、
閾値ε＝０．４の場合、Ｖｘ＝１．６４５２
閾値ε＝０．５の場合、Ｖｘ＝１．７０２７
といったように、誤差が極めて小さい結果が得られた。
【０１３０】
このように、本周波数解析手法を適用した動き推定装置は、カラー動画像の場合にも高精度に動き推定を行うことが可能である。特に、動き推定装置においては、動画像に応じて適切な閾値εを用いることにより、動き推定の精度を向上させることができる。
【０１３１】
［第７の検証］
つぎに、自然動画像に適用した場合の例として、図３３に示すように、いわゆる「Coast
Guard」という標準動画像を用いた。この動画像は、前景であるボートはほとんど静止した状態であるが、カメラワークによって背景が左方に移動するように観察されるものである。この実験では、１３０フレーム目から１６１フレーム目までの３２フレーム間におけるボートの移動量を推定した。
【０１３２】
まず、３次元ＦＦＴによる時空間スペクトルの分布を求めると、図３４に示すように、周波数ｆｘ，ｆｙの低周波領域に誤差と思われるばらつきが顕著にみられた。一方、本周波数解析手法による１００本の時空間スペクトルの分布を求めると、図３５に示すように、ばらつきがほとんどみられない結果が得られた。そして、これらの結果から、最小二乗法を用いたクラスタリング法によって平面を推定すると、３次元ＦＦＴ及び本周波数解析手法の場合についてそれぞれ図３６及び図３７に示す結果が得られた。求めた平面の傾きに基づいて速度を求めると、次表１に示す結果が得られた。
【０１３３】
【表１】
0005467346.tif 000030
【０１３４】
３次元ＦＦＴの場合には、背景については一応推定ができていたが、その推定速度結果に誤差があり、また、ボートについては推定することができず、誤った速度を推定していることがわかる。一方、本周波数解析手法の場合には、誤差が小さいことが確認された。
【０１３５】
このように、本周波数解析手法を適用した動き推定装置は、自然動画像の場合にも高精度に動き推定を行うことが可能である。
【０１３６】
［動き推定装置の効果］
以上説明したように、本周波数解析手法を適用した動き推定装置は、動画像符号化において動き補償や動き予測、補間といった処理を行う場合に、演算量を大幅に削減し且つ高精度に処理を行うことができる、したがって、この動き推定装置は、動画像内に含まれる物体の動きを推定・予測したり、フレーム補間等を行ったりする際に利用して好適である。
【０１３７】
また、この動き推定装置は、動画像符号化に応用することにより、従来に比べて大幅な高圧縮化を図ることができる。さらに、この動き推定装置は、大規模な変更を加えることなく、処理の質的・量的変化に対応することができる、スケーラブルな符号化技術の開発に利用可能である。これにより、この動き推定装置は、昨今増加傾向にあるネットワークトラフィックの低減に貢献することができる。
【０１３８】
なお、本発明は、上述した実施の形態に限定されるものではない。例えば、上述した実施の形態では、動き推定装置によってソフトウェアによる周波数解析を行うものとして説明したが、本発明は、本周波数解析手法を含む動き推定処理のアルゴリズムを実装したＤＳＰ（Digital Signal Processor）等、積和演算を行うことが可能であればハードウェアによっても実現することができる。
【０１３９】
このように、本発明は、その趣旨を逸脱しない範囲で適宜変更が可能であることはいうまでもない。
【符号の説明】
【０１４０】
１１  ＣＰＵ
１２  ＲＯＭ
１３  ＲＡＭ
１４  記憶部
１５  入力操作制御部
１６  表示部
【図１】
0005467346.tif 000031
【図２】
0005467346.tif 000032
【図３】
0005467346.tif 000033
【図４】
0005467346.tif 000034
【図５】
0005467346.tif 000035
【図６】
0005467346.tif 000036
【図７】
0005467346.tif 000037
【図８】
0005467346.tif 000038
【図９】
0005467346.tif 000039
【図１０】
0005467346.tif 000040
【図１１】
0005467346.tif 000041
【図１２】
0005467346.tif 000042
【図１３】
0005467346.tif 000043
【図１４】
0005467346.tif 000044
【図１５】
0005467346.tif 000045
【図１６】
0005467346.tif 000046
【図１７】
0005467346.tif 000047
【図１８】
0005467346.tif 000048
【図１９】
0005467346.tif 000049
【図２０】
0005467346.tif 000050
【図２１】
0005467346.tif 000051
【図２２】
0005467346.tif 000052
【図２３】
0005467346.tif 000053
【図２４】
0005467346.tif 000054
【図２５】
0005467346.tif 000055
【図２６】
0005467346.tif 000056
【図２７】
0005467346.tif 000057
【図２８】
0005467346.tif 000058
【図２９】
0005467346.tif 000059
【図３０】
0005467346.tif 000060
【図３１】
0005467346.tif 000061
【図３２】
0005467346.tif 000062
【図３３】
0005467346.tif 000063
【図３４】
0005467346.tif 000064
【図３５】
0005467346.tif 000065
【図３６】
0005467346.tif 000066
【図３７】
0005467346.tif 000067