JP5504426B2 - Object motion estimation apparatus, object motion estimation method, program, and recording medium - Google Patents

Object motion estimation apparatus, object motion estimation method, program, and recording medium Download PDF

Info

Publication number
JP5504426B2
JP5504426B2 JP2011501533A JP2011501533A JP5504426B2 JP 5504426 B2 JP5504426 B2 JP 5504426B2 JP 2011501533 A JP2011501533 A JP 2011501533A JP 2011501533 A JP2011501533 A JP 2011501533A JP 5504426 B2 JP5504426 B2 JP 5504426B2
Authority
JP
Japan
Prior art keywords
optical flow
luminance
gradient
equation
time
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2011501533A
Other languages
Japanese (ja)
Other versions
JPWO2010098160A1 (en
Inventor
昇 瀬部
英沢 延山
修一 榎田
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Kyushu Institute of Technology NUC
Original Assignee
Kyushu Institute of Technology NUC
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Kyushu Institute of Technology NUC filed Critical Kyushu Institute of Technology NUC
Priority to JP2011501533A priority Critical patent/JP5504426B2/en
Publication of JPWO2010098160A1 publication Critical patent/JPWO2010098160A1/en
Application granted granted Critical
Publication of JP5504426B2 publication Critical patent/JP5504426B2/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/20—Analysis of motion
    • G06T7/269—Analysis of motion using gradient-based methods

Description

本発明は、物体運動推定装置、物体運動推定方法、プログラム及び記録媒体に関し、特に、撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定装置等に関する。   The present invention relates to an object motion estimation device, an object motion estimation method, a program, and a recording medium, and more particularly, to an object motion estimation device that estimates an image of an object by performing image analysis processing on a captured image obtained by imaging an object. .

オプティカルフローは、観測者と物体との相対運動によって生じる画像上の見かけの速度分布である。オプティカルフローは、通常、画像上の各点での見かけの速度ベクトルとして表現される。   The optical flow is an apparent velocity distribution on the image generated by the relative motion between the observer and the object. The optical flow is usually expressed as an apparent velocity vector at each point on the image.

従来の勾配法では、動画像処理において、物体の画像の輝度が移動後も保存されることを仮定し、その保存則を表現したオプティカルフロー方程式と呼ばれる式(1)の偏微分方程式をもとにオプティカルフローを推定している（非特許文献１参照）。ただし、Ｉ(x,y,t)は輝度、（ｘ，ｙ）は画像上の座標、ｔは時刻、［ｕ ｖ］Tはオプティカルフローを表す。 In the conventional gradient method, it is assumed that the brightness of the object image is preserved after moving in the moving image processing, and the partial differential equation (1) called the optical flow equation expressing the conservation law is used. The optical flow is estimated (see Non-Patent Document 1). However, I (x, y, t) represents luminance, (x, y) represents coordinates on the image, t represents time, and [u v] T represents optical flow.

式(1)だけでは２変数ｕ、ｖを決定する方程式としては拘束条件が不足する。そこで、注目点の近傍でオプティカルフローがすべて等しいと仮定することにより、式(2)というオプティカルフローに関する線形方程式が得られる（非特許文献５参照）。ただし、添字はＩtx=（∂2／∂t∂x）Ｉのように偏微分を表す。この線形方程式がオプティカルフローの本質的な推定原理を与えている。ここに現れる輝度のヘッセ行列(3)は、画像中の物体の角と密接に関連した量であり、線形方程式の解きやすさと関連した判別式det(H)という量は推定されたオプティカルフローの信頼度としても利用されている。 The equation (1) alone has insufficient constraints as an equation for determining the two variables u and v. Therefore, by assuming that the optical flows are all equal in the vicinity of the point of interest, a linear equation related to the optical flow of Expression (2) is obtained (see Non-Patent Document 5). However, the subscript represents partial differentiation as I tx = (∂ 2 / ∂t∂x) I. This linear equation gives the essential principle of optical flow estimation. The luminance Hessian matrix (3) that appears here is a quantity closely related to the angle of the object in the image, and the quantity determinant det (H) related to the ease of solving the linear equation is the quantity of the estimated optical flow. It is also used as reliability.

B．K．P．Horn、外１名著，“Determining optical flow.”Artificial Intelligence，17:185-203，1981.B. K. P. Horn, 1 other author, “Determining optical flow.” Artificial Intelligence, 17: 185-203, 1981. S．J．Julier、外１名著，“Unscented filtering and nonlinear estimation.”Proceedings of IEEE，92(3):401-422，2004.S. J. Julier, 1 other author, “Unscented filtering and nonlinear estimation.” Proceedings of IEEE, 92 (3): 401-422, 2004. A．J．Krener、外１名著.“Nonlinear observers with linearizable error dynamics.”SIAM Journal of control and optimization，23(2):197-216，1985.A. J. Krener, 1 other author, “Nonlinear observers with linearizable error dynamics.” SIAM Journal of control and optimization, 23 (2): 197-216, 1985. K．Lucas、外１名著，“An interactive image registration techniques with an application to stereo vision.”In Proceedings of 7th IJCAL，pages 674-679，1981.K. Lucas, 1 other author, “An interactive image registration techniques with an application to stereo vision.” In Proceedings of 7th IJCAL, pages 674-679, 1981. A．Verri、外２名著，“Differential techniques for optical flow.”Journal of optical Society of America A，7(5):912-922，1990.A. Verri, 2 other authors, “Differential techniques for optical flow.” Journal of optical Society of America A, 7 (5): 912-922, 1990. Y．Wu、外２名著，“Unscented Kalman filtering for additive noise case：Augmented versus nonaugmented.”IEEE Signal Processing Letters，12(5):357-360，2005.Y. Wu, 2 other authors, “Unscented Kalman filtering for additive noise case: Augmented versus nonaugmented.” IEEE Signal Processing Letters, 12 (5): 357-360, 2005. 松永、外２名著.“カルマンフィルタを用いたオプティカルフローの抽出に関する研究. ”電子情報通信学会技術研究報告．PRMU，パターン認識・メディア理解，97(501):15-20，1998.Matsunaga and 2 other authors. "Research on optical flow extraction using Kalman filter." IEICE technical report. PRMU, Pattern recognition / Media understanding, 97 (501): 15-20, 1998.

オプティカルフローの研究の方向性としては、空間的広がりを利用して空間差分による雑音の影響を低減する研究（Lucasら（非特許文献４参照）など）、ヘッセ行列の正則性が失われてオプティカルフローの推定の信頼性が劣化する画像上の領域に対して、オプティカルフローの滑らかさなどの新たな拘束条件を要請し信頼性の高い部分のオプティカルフローを用いて補間する大域法（非特許文献１など参照）などが挙げられる。   The direction of optical flow research is to reduce the effects of noise due to spatial differences using spatial spread (Lucas et al. (See Non-Patent Document 4), etc.), and the regularity of the Hessian matrix is lost. Global method that interpolates using optical flow of highly reliable part by requesting new constraints such as smoothness of optical flow for regions on image where reliability of flow estimation deteriorates (non-patent literature) 1) and the like.

しかしながら、多くの改善が提案されているにもかかわらず、勾配法の基本原理はHornらの提案（非特許文献１参照）から25年間以上不変のままである。   However, although many improvements have been proposed, the basic principle of the gradient method has remained unchanged for more than 25 years since the Horn et al. Proposal (see Non-Patent Document 1).

従来の勾配法の問題点としては、例えば、従来法では、オプティカルフローの並進性（つまり、空間勾配が０）を仮定していることがある。オプティカルフローの並進性の仮定は、本来未知であるオプティカルフローに対する仮定であり、可能な限り緩めることが重要である。   As a problem of the conventional gradient method, for example, in the conventional method, there is a case where the translational property of the optical flow (that is, the spatial gradient is 0) is assumed. The assumption of translational nature of the optical flow is an assumption for the optical flow that is originally unknown, and it is important to relax as much as possible.

ゆえに、本発明は、従来の勾配法における問題点を克服し、撮像画像における画像解析処理によるオプティカルフローの推定に適した物体動作推定装置等を提供することを目的とする。   Therefore, an object of the present invention is to provide an object motion estimation apparatus and the like suitable for estimating an optical flow by image analysis processing on a captured image, overcoming the problems in the conventional gradient method.

請求項１に係る発明は、撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定装置において、前記撮像画像の注目点における輝度勾配を演算する輝度勾配演算手段と、前記輝度勾配から時間勾配を演算する時間勾配演算手段と、連立方程式を解く方程式演算手段とを有する画像解析処理手段を備え、前記連立方程式は、前記撮像手段と前記物体との相対運動の速度分布であるオプティカルフローについて、前記撮像画像における前記物体の輝度の時間に関する１階全微分がゼロという時不変性が仮定されることにより、前記オプティカルフローを変数とし、前記撮像画像の画素の輝度の空間勾配及びその時間微分を係数とするだけでなく、さらに、前記撮像手段の撮像軸に対して垂直な平面上を前記物体が移動するという並進性が仮定されないことにより、前記オプティカルフローの空間勾配をも変数とし、前記輝度及びその空間勾配の時間に関して２階微分した値をも係数とするものであり、前記画像解析処理手段が有する方程式演算手段が、前記連立方程式を解くことにより、前記オプティカルフローを推定すると同時に、前記推定されるオプティカルフローの空間勾配をも推定するものである。   According to the first aspect of the present invention, in the object motion estimation device that estimates the motion of the object by performing image analysis processing on the captured image obtained by capturing the object by the imaging unit, the brightness gradient that calculates the brightness gradient at the point of interest of the captured image An image analysis processing unit including a calculation unit, a time gradient calculation unit that calculates a time gradient from the luminance gradient, and an equation calculation unit that solves simultaneous equations, wherein the simultaneous equations are relative to the imaging unit and the object. With respect to the optical flow that is a velocity distribution of motion, invariance is assumed when the first-order total derivative with respect to the luminance time of the object in the captured image is zero, so that the optical flow is a variable, and the pixel of the captured image In addition to using the spatial gradient of the brightness and the time derivative thereof as a coefficient, the plane is perpendicular to the imaging axis of the imaging means. Since the translational nature of the body moving is not assumed, the spatial gradient of the optical flow is also used as a variable, and the value obtained by second-order differentiation with respect to the luminance and the time of the spatial gradient is also used as the coefficient. The equation calculation means included in the processing means estimates the optical flow by solving the simultaneous equations, and at the same time, estimates the spatial gradient of the estimated optical flow.

請求項２に係る発明は、請求項１記載の物体運動推定装置であって、前記撮像画像は、二次元の画像であり、少なくとも３画素×３画素の大きさの画像であって、撮像時刻が異なる３枚の画像であり、撮像画像記憶手段に記憶されており、前記輝度勾配演算手段は、前記撮像画像の注目点Ｐ（ｘ0，ｙ0）における輝度Ｉの輝度勾配Ｉx（ｘ0，ｙ0，ｔ）、Ｉy（ｘ0，ｙ0，ｔ）、Ｉxx（ｘ0，ｙ0，ｔ）、Ｉxy（ｘ0，ｙ0，ｔ）及びＩyy（ｘ0，ｙ0，ｔ）を演算して、式(eq1)に示される観測量ｚ（ｔ）を求めるものであり、前記時間勾配演算手段は、前記輝度Ｉ及び前記輝度勾配Ｉx、Ｉy、Ｉxx、Ｉxy及びＩyyから時間勾配Ｉt、Ｉtx、Ｉty、Ｉtxx、Ｉtxy、Ｉtyy、Ｉtt、Ｉttx、Ｉtty、Ｉttxx、Ｉttxy及びＩttyyを演算するものであり、前記方程式演算手段は、線形方程式(eq2)を解くことにより、前記オプティカルフロー［ｕ ｖ］Tを推定すると同時に、その空間勾配［ｕx ｖx ｕy ｖy］Tも推定するものである。 The invention according to claim 2 is the object motion estimation device according to claim 1, wherein the captured image is a two-dimensional image, and is an image having a size of at least 3 pixels × 3 pixels, and an imaging time Are stored in the captured image storage means, and the brightness gradient calculating means is the brightness gradient I x (x of the brightness I at the point of interest P (x 0 , y 0 ) of the captured image. 0, y 0, t), I y (x 0, y 0, t), I xx (x 0, y 0, t), I xy (x 0, y 0, t) and I yy (x 0, y 0 , t) is calculated to obtain the observation amount z (t) shown in the equation (eq1), and the time gradient calculating means is configured to calculate the luminance I and the luminance gradients I x , I y , I Calculate the time gradients I t , I tx , I ty , I txx , I txy , I tyy , I tt , I ttx , I tty , I ttxx , I ttxy and I ttyy from xx , I xy and I yy The equation calculation means estimates the optical flow [u v] T by solving the linear equation (eq2), and also estimates the spatial gradient [u x v x u y v y ] T To do.

請求項３に係る発明は、撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定装置において、推定モデルに対して、状態変数の逐次推定を行う状態推定処理手段を有する画像解析処理手段を備え、前記推定モデルは、前記撮像手段と前記物体との相対運動の速度分布であるオプティカルフローについて、前記撮像画像における前記物体の輝度の時間に関する１階全微分がゼロという時不変性が仮定されることにより、前記オプティカルフローを状態変数に含むだけでなく、前記撮像手段の撮像軸に対して垂直な平面上を前記物体が移動するという並進性が仮定されないことにより、前記オプティカルフローの空間勾配をも状態変数に含み、さらに、前記撮像画像上の注目点の空間近傍の輝度を観測量とする観測方程式を含むものであり、前記状態推定処理手段が、前記状態変数の逐次推定を行うことにより、前記オプティカルフローを推定すると同時に、前記推定されるオプティカルフローの空間勾配をも推定するものである。 Invention is carried out in the object motion estimating device a captured image pickup means has captured the object image analysis processing to estimate the motion of the object, for the estimated model, the sequential estimation of the state variables according to claim 3 An image analysis processing unit having a state estimation processing unit, wherein the estimation model is a first floor related to a time of luminance of the object in the captured image with respect to an optical flow which is a velocity distribution of relative motion between the imaging unit and the object. By assuming invariance when the total derivative is zero, not only the optical flow is included in the state variable, but also the translational property that the object moves on a plane perpendicular to the imaging axis of the imaging means. By not being assumed, the spatial gradient of the optical flow is also included in the state variable, and the luminance near the space of the point of interest on the captured image is also observed. The state estimation processing means estimates the optical flow by performing sequential estimation of the state variable, and simultaneously estimates the spatial gradient of the estimated optical flow. is there.

請求項４に係る発明は、請求項３に記載の物体運動推定装置であって、前記観測方程式は、システム雑音と観測雑音という雑音を加味したものであり、前記状態変数の逐次推定は、推定値の共分散を近似するカルマンフィルタを用いて前記状態変数を逐次推定するものである。 The invention according to claim 4 is the object motion estimation apparatus according to claim 3 , wherein the observation equation takes into account noises of system noise and observation noise, and the sequential estimation of the state variable is an estimation The state variables are successively estimated using a Kalman filter that approximates the covariance of values.

請求項５に係る発明は、撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定方法において、輝度勾配演算手段が、前記撮像画像の注目点における輝度勾配を演算するステップと、時間勾配演算手段が、前記輝度勾配から時間勾配を演算するステップと、方程式演算手段が、前記撮像手段と前記物体との相対運動の速度分布であるオプティカルフローについて、前記撮像画像における前記物体の輝度の時間に関する１階全微分がゼロという時不変性が仮定されることにより、前記オプティカルフローを変数とし、前記撮像画像の画素の輝度の空間勾配及びその時間微分を係数とするだけでなく、さらに、前記撮像手段の撮像軸に対して垂直な平面上を前記物体が移動するという並進性が仮定されないことにより、前記オプティカルフローの空間勾配をも変数とし、前記輝度及びその空間勾配の時間に関して２階微分した値をも係数とする連立方程式を解くステップを含むことにより、前記オプティカルフローを推定すると同時に、前記推定されるオプティカルフローの空間勾配をも推定する。   According to a fifth aspect of the present invention, in the object motion estimation method for estimating the motion of the object by performing image analysis processing on the captured image obtained by capturing the object by the image capturing unit, the brightness gradient calculating unit includes the brightness at the target point of the captured image. A step of calculating a gradient, a step of calculating a time gradient from the luminance gradient, a step of calculating a gradient from the luminance gradient, and an equation calculating unit for an optical flow that is a velocity distribution of relative motion between the imaging unit and the object, By assuming invariance when the first-order total derivative of the luminance of the object in the captured image is zero, the optical flow is a variable, the spatial gradient of the luminance of the pixel of the captured image and its temporal derivative are coefficients. Furthermore, the translational property that the object moves on a plane perpendicular to the imaging axis of the imaging means is not assumed. And a step of solving a simultaneous equation having a spatial gradient of the optical flow as a variable and a second-order differential value with respect to the luminance and the time of the spatial gradient as a coefficient, and simultaneously estimating the optical flow. The spatial gradient of the estimated optical flow is also estimated.

【００１６】
【請求項６】
撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定方法において、状態推定処理手段が、前記撮像手段と前記物体との相対運動の速度分布であるオプティカルフローについて、前記撮像画像における前記物体の輝度の時間に関する１階全微分がゼロという時不変性が仮定されることにより、前記オプティカルフローを状態変数に含むだけでなく、前記撮像手段の撮像軸に対して垂直な平面上を前記物体が移動するという並進性が仮定されないことにより、前記オプティカルフローの空間勾配をも状態変数に含み、さらに、前記撮像画像上の注目点の空間近傍の輝度を観測量とする観測方程式を含む推定モデルに対して、前記状態変数の逐次推定を行うステップを含むことにより、前記オプティカルフローを推定すると同時に、前記推定されるオプティカルフローの空間勾配をも推定する。 [0016]
[Claim 6]
In an object motion estimation method for estimating a motion of the object by performing image analysis processing on a captured image obtained by imaging the object by the imaging unit, the state estimation processing unit is an optical that is a velocity distribution of relative motion between the imaging unit and the object. Assuming that the flow is invariant when the first-order total derivative with respect to the luminance time of the object in the captured image is zero, the optical flow is not only included in the state variable, but is also included in the imaging axis of the imaging unit. Since the translational property that the object moves on a plane perpendicular to the plane is not assumed, the spatial gradient of the optical flow is also included in the state variable, and the luminance in the vicinity of the point of interest on the captured image is observed. A step of sequentially estimating the state variable with respect to an estimation model including an observation equation as a quantity; And simultaneously estimating a row, also to estimate the spatial gradient of the optical flow the estimation.

請求項７に係る発明は、コンピュータを、請求項１から４のいずれかの物体運動推定装置として動作させるためのプログラムである。   The invention according to claim 7 is a program for causing a computer to operate as the object motion estimation device according to any one of claims 1 to 4.

請求項８に係る発明は、請求項７記載のプログラムを記録したコンピュータ読み取り可能な記録媒体である。   The invention according to claim 8 is a computer-readable recording medium on which the program according to claim 7 is recorded.

なお、請求項２において、線形方程式（ｅｑ２）は過制約の線形方程式である。線形方程式（ｅｑ２）の解は、例えば、最小二乗法などの推定法によって求める。   In claim 2, the linear equation (eq2) is an overconstrained linear equation. The solution of the linear equation (eq2) is obtained by an estimation method such as a least square method.

本願の各請求項に係る発明（以下、「本願発明」という。）は、データとして取得された撮像画像に関わる物理的性質を利用した処理を具体的に行うものであり、オプティカルフローの並進性の仮定を必要とせず、また、画像輝度のヘッセ行列の正則性を必ずしも要請せずに、オプティカルフローの推定が可能となる。従来にも、オプティカルフローが並進的でない場合にも、例えば単純拡大や回転など、その運動に関する事前情報が分かっている場合には、オプティカルフローの推定が可能であることは示されている（非特許文献５参照）。しかし、撮像画像の各点の輝度以外に、オプティカルフローに関する事前情報が必要である点には変わりがない。本願発明によれば、このような事前情報を必要とせず、オプティカルフローの推定が可能である。   The invention according to each claim of the present application (hereinafter referred to as “the present invention”) specifically performs processing using physical properties related to a captured image acquired as data, and translates optical flow. Therefore, it is possible to estimate the optical flow without requiring the regularity of the Hessian matrix of the image luminance. Traditionally, it has been shown that even if the optical flow is not translational, it is possible to estimate the optical flow if the prior information about the motion is known, for example, simple enlargement or rotation (non-null) (See Patent Document 5). However, there is no change in that prior information regarding the optical flow is required in addition to the luminance of each point of the captured image. According to the present invention, optical flow can be estimated without requiring such prior information.

さらに、本願発明によれば、オプティカルフローの推定と同時に、その空間勾配をも推定することが可能となる。オプティカルフローの空間勾配はベクトル場の回転や発散の情報を含み、三次元空間における物体の運動復元を行う上でとても重要な量である。本願発明は、空間勾配という重要な量について、事前に仮定を必要とせず、かつ、オプティカルフローと同時に推定可能とするものである。   Furthermore, according to the present invention, it is possible to estimate the spatial gradient simultaneously with the estimation of the optical flow. The spatial gradient of the optical flow contains information on the rotation and divergence of the vector field and is a very important quantity for restoring the motion of the object in the three-dimensional space. The present invention makes it possible to estimate an important quantity called a spatial gradient at the same time as an optical flow without requiring an assumption in advance.

このように、従来の勾配法では、並進的な動き（つまり、カメラ光軸に垂直な平面上を回転もせずに等速に動くこと）を仮定していた。そのため、光軸に沿った動きや回転を含む運動を計測する際は原理的に精度が劣化することとなる。これに対して、本願発明は、同時に空間勾配をも推定（補正）することができ、原理的に精度が向上する。このように、画像の回転や拡大が直接検出でき、画像の回転や拡大のオプティカルフローへの影響が補正されるので、回転や拡大がある場合のオプティカルフローの精度が向上する。さらに、回転や拡大の影響を用いることで、オプティカルフローが計算可能な画像が増えることとなる。   Thus, in the conventional gradient method, translational movement (that is, moving at a constant speed without rotating on a plane perpendicular to the optical axis of the camera) is assumed. Therefore, when measuring a motion including movement and rotation along the optical axis, the accuracy deteriorates in principle. On the other hand, the present invention can estimate (correct) the spatial gradient at the same time, and the accuracy is improved in principle. As described above, the rotation and enlargement of the image can be directly detected, and the influence of the rotation and enlargement of the image on the optical flow is corrected. Therefore, the accuracy of the optical flow in the case of the rotation and enlargement is improved. Furthermore, by using the effects of rotation and enlargement, the number of images for which optical flow can be calculated increases.

また、従来の勾配法の問題点として、雑音対策に重きを置くあまり、推定原理の本質が見えにくくなっていることがある。例えば、不足する方程式の数を補うために空間的な広がりを利用する方法もあるが、これは、考えている領域内でオプティカルフローが均一であるという仮定を置く。しかし、この仮定は十分小さな領域を考える場合に成り立つが、領域を小さくすると式(1)の方程式の独立性が失われる。領域の大きさによる仮定の成立性と方程式の独立性が相反する。時空間で離散化されたデータを前提とする方法であり、本質的な推定可能性を考えるには不向きである。   Another problem with the conventional gradient method is that the essence of the estimation principle is difficult to see due to the emphasis on noise countermeasures. For example, some methods use spatial spread to compensate for the number of missing equations, which puts the assumption that the optical flow is uniform within the region under consideration. However, this assumption holds when considering a sufficiently small area, but if the area is reduced, the independence of the equation (1) is lost. There is a contradiction between the validity of assumptions based on the size of the region and the independence of equations. This method is based on data that is discretized in space-time, and is not suitable for considering the essential estimability.

本願発明によれば、あくまで注目点の近傍だけの情報、言い換えるならば注目点の輝度及びその時空間勾配の情報のみからオプティカルフローを推定することが可能となる。   According to the present invention, it is possible to estimate an optical flow only from information only in the vicinity of the attention point, in other words, only information on the luminance of the attention point and its spatiotemporal gradient.

さらに、従来の勾配法には、実装上の問題点がある。従来の勾配法では、輝度情報の時空間微分の陽な計算が必要であり、画像情報からえるためには差分演算が不可避であった。この差分演算は外乱に弱く、オプティカルフローの推定の信頼性に大きく影響する。また、これらの差分情報から線形連立方程式を解くためにヘッセ行列の正則性などが要請される。ここでも逆行列計算による推定の劣化が問題となる。   Furthermore, the conventional gradient method has a problem in implementation. The conventional gradient method requires explicit calculation of the spatio-temporal differentiation of luminance information, and the difference calculation is inevitable in order to obtain from the image information. This difference calculation is vulnerable to disturbances and greatly affects the reliability of optical flow estimation. Further, in order to solve the linear simultaneous equations from these difference information, regularity of the Hessian matrix is required. Again, degradation of estimation by inverse matrix calculation is a problem.

本願の請求項４及び５に係る発明によれば、時空間微分を用いず、かつ、陽な逆行列計算を用いることなく推定を行うことが可能となる。このように、例えばカルマンフィルタを用いる等による雑音対策により、逆行列の陽な計算を回避することができ、結果としてオプティカルフローの計算精度が向上することとなる。   According to the inventions according to claims 4 and 5 of the present application, it is possible to perform estimation without using space-time differentiation and without using a positive inverse matrix calculation. Thus, for example, by taking measures against noise by using a Kalman filter or the like, explicit calculation of the inverse matrix can be avoided, and as a result, the calculation accuracy of the optical flow is improved.

本願発明の実施例である物体運動推定装置１の概略ブロック図である。It is a schematic block diagram of the object motion estimation apparatus 1 which is an Example of this invention. 図１の物体運動推定装置１による画像処理の一例を示すフロー図である。It is a flowchart which shows an example of the image process by the object motion estimation apparatus 1 of FIG. 実施例２において、瞬時推定の結果を示す図である。In Example 2, it is a figure which shows the result of instantaneous estimation. 実施例２において、カルマンフィルタを用いた推定の結果を示す図である。In Example 2, it is a figure which shows the result of the estimation using a Kalman filter.

以下、図面を参照して、本願発明の実施例について説明する。ただし、本願発明は、記載した実施例に限定されるものではない。   Embodiments of the present invention will be described below with reference to the drawings. However, the present invention is not limited to the described embodiments.

図１は、本願発明の実施例である物体運動推定装置１の概略ブロック図である。物体運動推定装置１は、例えばドライブレコーダである。   FIG. 1 is a schematic block diagram of an object motion estimation apparatus 1 that is an embodiment of the present invention. The object motion estimation apparatus 1 is, for example, a drive recorder.

物体運動推定装置１は、画像を撮像する撮像部３と、撮像部３により撮像された撮像画像を記憶する撮像画像記憶部５と、撮像画像上の各点について、時間又は／及び輝度の量子化を行い計測画像とする量子化部７と、計測画像を記憶する計測画像記憶部９と、計測画像に平滑化を行い平滑化画像とする平滑化部１１と、平滑化画像を記憶する平滑化画像記憶部１３と、量子化画像に基づいてオプティカルフローを推定する画像解析処理部１５と、推定されたオプティカルフローを記憶するオプティカルフロー記憶部１７を備える。   The object motion estimation device 1 includes an imaging unit 3 that captures an image, a captured image storage unit 5 that stores a captured image captured by the imaging unit 3, and a quantum of time or / and luminance for each point on the captured image. Quantization unit 7 for converting to a measurement image, measurement image storage unit 9 for storing the measurement image, smoothing unit 11 for smoothing the measurement image to obtain a smoothed image, and smoothing for storing the smoothed image A quantized image storage unit 13, an image analysis processing unit 15 for estimating an optical flow based on the quantized image, and an optical flow storage unit 17 for storing the estimated optical flow.

撮像部３は、例えば車載カメラである。車載カメラは、基本的に進行方向に向けられるので光軸に沿った動きが避けられない。また、ドライブレコーダのような広角画像では、画像の歪みによる拡大が避けられない。そのため、撮像画像には、大きな誤差が含まれることとなる。   The imaging unit 3 is an in-vehicle camera, for example. Since the in-vehicle camera is basically directed in the traveling direction, movement along the optical axis is inevitable. Further, in a wide-angle image such as a drive recorder, enlargement due to image distortion is inevitable. Therefore, a large error is included in the captured image.

量子化部７は、例えば、整数の時刻ｔで、撮像画像の（ｘ，ｙ）のｘとｙが整数の点について、小数点以下を切り捨てて量子化を行い、これを計測画像とするものである。量子化部７の処理により、量子化誤差が加わる。   For example, the quantization unit 7 performs quantization by rounding off the fractional part of a point where x and y of (x, y) of the captured image are integers at an integer time t, and this is used as a measurement image. is there. A quantization error is added by the processing of the quantization unit 7.

平滑化部１１は、計測画像の各時刻において、例えばガウスフィルタを用いて平滑化を行い、平滑化画像とする。   The smoothing unit 11 performs smoothing using, for example, a Gaussian filter at each time of the measurement image to obtain a smoothed image.

画像解析処理部１３は、平滑化画像の輝度情報を用いてオプティカルフローの推定を行う。撮像画像の注目点における輝度勾配を演算する輝度勾配演算部２１と、輝度勾配から時間勾配を演算する時間勾配演算部２３と、時間勾配を係数とし、オプティカルフロー及びその空間勾配を変数とする連立方程式を解く方程式演算部２５と、撮像画像の注目点の空間近傍の輝度情報を観測量としてフィルタを用いてオプティカルフロー及びその空間勾配を逐次推定する状態推定処理部２７を有する。以下では、画像解析処理部１３の動作を中心に、具体的に説明する。   The image analysis processing unit 13 performs optical flow estimation using luminance information of the smoothed image. A luminance gradient computing unit 21 that computes a luminance gradient at a point of interest in a captured image, a time gradient computing unit 23 that computes a time gradient from the luminance gradient, and a combination of an optical flow and its spatial gradient as variables, using the time gradient as a coefficient. An equation calculation unit 25 that solves an equation, and a state estimation processing unit 27 that sequentially estimates an optical flow and its spatial gradient using a filter with luminance information in the vicinity of the point of interest of the captured image as an observation amount. Below, it demonstrates concretely centering around operation | movement of the image analysis process part 13. FIG.

まず、オプティカルフロー方程式の一般化について説明する。ここで、時刻ｔにおける画像上の位置(x,y)における輝度をＩ(x,y,t)とし、Ｉ(x,y,t)は必要とされる回数だけ微分可能であると仮定する。一方、オプティカルフローを［ｕ(x,y,t) ｖ(x,y,t)］Tで表し、こちらも必要とされる回数だけ微分可能であると仮定する。なお、本実施例では、偏微分をＩt(=∂I／∂t)のように添字にて表すものとする。First, generalization of the optical flow equation will be described. Here, it is assumed that the luminance at the position (x, y) on the image at time t is I (x, y, t), and I (x, y, t) can be differentiated as many times as necessary. . On the other hand, the optical flow is represented by [u (x, y, t) v (x, y, t)] T , which is also assumed to be differentiable as many times as necessary. In this embodiment, the partial differentiation is represented by a subscript such as I t (= ∂I / ∂t).

輝度不変の拘束式(1)は、式(4)のように、輝度Ｉ(x,y,t)の時空間における［ｕ ｖ １］T方向への勾配が０であると解釈することもできる。The luminance invariant constraint equation (1) may be interpreted as the gradient in the [u v 1] T direction in the space-time of the luminance I (x, y, t) is zero as in the equation (4). it can.

時空間のある点Ｐ(x0,y0,t0)において、式(4)が成立しているとする。さらに、点Ｐの近傍点(x0+Δx,y0+Δy,t0+Δt)においてもオプティカルフローに沿った移動後も輝度が変化しないという条件は、近傍点における方向微分が０という条件として表される。式(5)は、近傍点における輝度、オプティカルフロー及びその空間微分を、微小量Δx，Δy，Δtで近似した条件である。ここで微小量Δx，Δy，Δtの２次の項を無視すると、任意の近傍で式(5)が成立するためには微小量Δx，Δy，Δtの１次の項の係数がそれぞれ０でなければならない。これより、近傍における各点でオプティカルフローに沿った移動後の輝度不変の条件は、式(6a)〜(6c)の３式が成立することとなる。It is assumed that the equation (4) is established at a point P (x 0 , y 0 , t 0 ) in time and space. Further, the condition that the luminance does not change after movement along the optical flow also in the vicinity point (x 0 + Δx, y 0 + Δy, t 0 + Δt) of the point P is that the directional derivative at the vicinity point is 0. Represented as: Equation (5) is a condition that approximates the luminance, optical flow, and spatial differentiation thereof at neighboring points with minute amounts Δx, Δy, and Δt. If the second-order terms of the minute amounts Δx, Δy, and Δt are ignored here, the coefficient of the first-order term of the minute amounts Δx, Δy, and Δt is 0 in order to satisfy the equation (5) in an arbitrary neighborhood. There must be. As a result, the three conditions of Equations (6a) to (6c) are satisfied as the condition of the luminance invariance after movement along the optical flow at each point in the vicinity.

式(6)の一般化されたオプティカルフロー恒等式は、画像上の注目点の近傍において輝度不変であるための条件式をテーラー展開し、その１次の項までで打ち切った条件である。そのように解釈すると、従来の輝度不変の拘束式(4)は、テーラー展開の０次の項に相当することが分かる。また、式(6)の条件は従来のオプティカルフロー方程式(4)をx，y，tについて偏微分することによっても導出ができ、テーラー展開のより高次の項に対応する輝度不変条件も同様に容易に導出することができる。式(7a)〜(7f)は、本実施例で用いる２次の項に対応する条件である。   The generalized optical flow identity of equation (6) is a condition in which the conditional expression for invariance of luminance in the vicinity of the point of interest on the image is Taylor-expanded and truncated to the first order term. When interpreted in that way, it can be seen that the conventional luminance invariant constraint equation (4) corresponds to the 0th-order term of the Taylor expansion. The condition of equation (6) can also be derived by partial differentiation of the conventional optical flow equation (4) with respect to x, y, and t, as well as the luminance invariant conditions corresponding to higher order terms in the Taylor expansion. Can be easily derived. Expressions (7a) to (7f) are conditions corresponding to the second-order terms used in this embodiment.

続いて、オプティカルフローの推定モデルについて説明する。本実施例では、オプティカルフロー方程式の拘束条件を一階非線形連立微分方程式として記述し、動的システムとして扱う。そして、動的システムの状態の可観測性の概念を使って、オプティカルフロー及びその空間勾配を同時に推定する推定原理及び推定可能条件を明らかにする。   Next, an optical flow estimation model will be described. In this embodiment, the constraint condition of the optical flow equation is described as a first-order nonlinear simultaneous differential equation and handled as a dynamic system. Then, using the concept of observability of the state of the dynamic system, the estimation principle and estimable conditions for simultaneously estimating the optical flow and its spatial gradient are clarified.

動的システムとして扱う場合、テーラー展開の次数で制約方程式を選択するのではなく、推定に用いる空間勾配の次数、推定するオプティカルフローの時空間勾配によってシステムを記述するための制約式を選択すべきである。そこで、式(8)の輝度の２次までの空間勾配ｉを用いて、式(9)のオプティカルフローとその空間勾配ｍを推定する。   When treating as a dynamic system, instead of selecting a constraint equation based on the order of Taylor expansion, a constraint equation for describing the system should be selected based on the order of the spatial gradient used for estimation and the spatiotemporal gradient of the optical flow to be estimated. It is. Therefore, the optical flow of Equation (9) and its spatial gradient m are estimated using the spatial gradient i up to the second order of the luminance of Equation (8).

推定に用いない輝度の３次以上の空間勾配と、推定を行わないオプティカルフローの２次以上の空間勾配は、０と仮定する。さらに、オプティカルフローは時不変であると仮定する。すなわち、ｕt＝０，ｖt＝０とする。この仮定は、オプティカルフローの空間勾配も時不変であることを意味する。It is assumed that the spatial gradient of the third or higher order of luminance that is not used for estimation and the second or higher spatial gradient of the optical flow that is not estimated are zero. Further assume that the optical flow is time invariant. That is, u t = 0 and v t = 0. This assumption means that the spatial gradient of the optical flow is also time-invariant.

これらの仮定のもと、画像上のある点Ｐ(x0,y0)における輝度、オプティカルフロー及びその空間勾配を用いて、状態ｑ(t)を式(10)と選ぶ。そして、これら変数の時間微分を記述する拘束式として式(4)，(6a)，(6b)，(7a)，(7b)，(7c)を選ぶことで、オプティカルフローの推定モデルが式(11a)〜(11d)のように得られる。ただし、Ｅm×n，Ｏm×nは、それぞれｍ×ｎの単位行列，零行列を表す。なお、関数の引数ｘ0，ｙ0，ｔは適宜省略する。Under these assumptions, the state q (t) is selected as Equation (10) using the luminance, the optical flow and the spatial gradient at a certain point P (x 0 , y 0 ) on the image. Then, by choosing the equations (4), (6a), (6b), (7a), (7b), (7c) as constraint equations that describe the time derivatives of these variables, the optical flow estimation model becomes the equation ( 11a) to (11d). Here, E m × n and O m × n represent an m × n unit matrix and a zero matrix, respectively. The function arguments x 0 , y 0 , and t are omitted as appropriate.

システムの可観測性とは、有限時間区間の出力の観測によって状態を一意に決定できるシステムの性質である。具体的には、式(12)の非線形微分方程式で表される動的システムを考える。ただし、状態ｑ(t)∈Ｒn，出力z(t)∈Ｒmである。このとき、非特許文献３を参照して、可観測性について定義し、補題１が成立する。System observability is the property of a system that can uniquely determine the state by observing the output in a finite time interval. Specifically, consider a dynamic system represented by the nonlinear differential equation of equation (12). However, the state q (t) εR n and the output z (t) εR m . At this time, observability is defined with reference to Non-Patent Document 3, and Lemma 1 is established.

定義：システム(12)が与えられたとする。任意の初期時刻ｔ0、初期状態ｑ0＝ｑ(t0)に対して有限の時刻ｔ1が存在し、観測出力ｚ(t)（t0≦t≦t1）と初期状態がｑ0の近傍にあるという知識から初期状態ｑ(t0)が一意に決定できるとき、システム(12)は状態ｑ0において可観測であるという。Definition: Given system (12). An arbitrary initial time t 0 , a finite time t 1 exists for the initial state q 0 = q (t 0 ), the observed output z (t) (t 0 ≦ t ≦ t 1 ), and the initial state is q 0. The system (12) is said to be observable in state q 0 when the initial state q (t 0 ) can be uniquely determined from the knowledge that it is in the vicinity of.

補題１：以下の命題(i)、(ii)及び(iii)は等価である。
(i) 与えられたシステム(12)は状態ｑ0において可観測である。
(ii) ｑ(t)＝ｑ0において（実ベクトルの意味で）、式(13)である。
(iii) ｑ(t)＝ｑ0において、式(14)である。
ただし、Lie微分は式(15)と、Lie微分の微分形式は式(16)と与えられる。
具体的な状態ｑは、式(17)の非線形連立代数方程式をｑについて解くことによって決定することができる。Lemma 1: The following propositions (i), (ii) and (iii) are equivalent.
(i) A given system (12) is observable in state q 0 .
(ii) When q (t) = q 0 (in the sense of a real vector), equation (13) is obtained.
(iii) Formula (14) when q (t) = q 0 .
However, Lie differentiation is given by equation (15), and the differential form of Lie differentiation is given by equation (16).
The specific state q can be determined by solving the nonlinear simultaneous algebraic equation (17) for q.

式(11)で定義される推定モデルの場合は、出力ｚ(t)から式(10)の状態ｑ(t)が一意に決定できることである。ここで、ｑ(t)はオプティカルフロー［ｕ ｖ］Tを含むので、可観測であればオプティカルフローが一意に決定できることを意味する。オプティカルフローの推定モデル(11)の非線形動的システムとしての可観測性は、Lie微分を用いて調べることができる。式(11d)よりLie微分の微分形式dhは式(18)である。微分形式dhが上記の構造を持つので、可観測性を調べるためには微分形式dLk fh（k=1，...，11)の第7〜12列の列フルランク性を調べればよい。以下では、行列Ａの第7列以降に対応する部分行列を（Ａ）2で表す。In the case of the estimation model defined by Expression (11), the state q (t) of Expression (10) can be uniquely determined from the output z (t). Here, since q (t) includes the optical flow [u v] T , it means that the optical flow can be uniquely determined if it is observable. The observability of the optical flow estimation model (11) as a nonlinear dynamic system can be examined using Lie derivative. From equation (11d), the differential form dh of the Lie derivative is equation (18). Since the differential form dh has the above structure, in order to investigate the observability, the column full rank property of the 7th to 12th columns of the differential form dL k f h (k = 1, ..., 11) is examined. Good. Hereinafter, a partial matrix corresponding to the seventh column and subsequent columns of the matrix A is represented by (A) 2 .

ここで、式(11)のシステムに対する式(16)を計算すると、式(19)となる。ただし、Ｎ、Ｊ、Ｊtは、それぞれ、式(20)、(21)、(22)である。これを用いて実際に微分形式のランクを計算すると、式(23)であるものの式(24)となり、一般の状態ｑ(t)においてシステムは可観測となり、オプティカルフローとその空間勾配ｍが同時に推定可能である。Here, when equation (16) is calculated for the system of equation (11), equation (19) is obtained. However, N, J, and Jt are formulas (20), (21), and (22), respectively. When the rank of the differential form is actually calculated by using this, the equation (23) becomes the equation (24), and the system becomes observable in the general state q (t), and the optical flow and its spatial gradient m are simultaneously measured. It can be estimated.

ただし、特別な状態（例えばｑ(t)＝Ｏ12×1）における可観測性は必ずしも保証されないので、以下、可観測性が成立するために輝度、オプティカルフローが満たすべき条件について考察を与える。式(19)の右辺一つ目の行列が常に正則であることより、式(25)が成り立ち、［ＪT Ｊt T］Tが列フルランクとなるために状態ｑ(t)が満たすべき条件を明らかにすればよい。However, since observability in a special state (for example, q (t) = O 12 × 1 ) is not necessarily guaranteed, consideration will be given to the conditions that the luminance and the optical flow must satisfy in order to establish observability. Since the first matrix on the right side of equation (19) is always regular, equation (25) holds, and [J T J t T ] T has column full rank, so the state q (t) must be satisfied. The conditions should be clarified.

定理１：式(24)が満たされ式(11)のシステムが可観測となるための十分条件は、以下の(i)、(ii)の２つの条件のいずれかが成立することである。
(i) 輝度のヘッセ行列Ｈが式(26)を満たし、さらに式(27a)、(27b)、(27c)のいずれかの付帯条件を満たすこと（ただし、αは任意の定数）。
(ii) 輝度の１次の空間勾配とその時間微分が式(28)を満たすこと。
なお, 定理１の(i)の３つの付帯条件が満たされない場合でも、オプティカルフロー［ｕ ｖ］Tの推定は可能である。Theorem 1: A sufficient condition for satisfying equation (24) and observability of the system of equation (11) is that either of the following two conditions (i) or (ii) is satisfied.
(i) The luminance Hessian matrix H satisfies the equation (26) and further satisfies the incidental condition of any one of the equations (27a), (27b), and (27c) (where α is an arbitrary constant).
(ii) The first-order spatial gradient of brightness and its time derivative satisfy equation (28).
Even if the three incidental conditions of (i) of Theorem 1 are not satisfied, the optical flow [u v] T can be estimated.

定理１(i)の式(26)の条件は、画像の輝度のみから推定可能性が決まる条件で、従来から知られている輝度のヘッセ行列の正則性と等価である。一方、(ii)の式(28)の条件は、画像の輝度勾配の時間変化によって推定可能性が決まるものであり、これまで知られていない新しい条件である。この条件(28)が十分条件として加わることによって、従来法よりも推定可能条件が緩和され、推定可能な状況、つまり輝度とオプティカルフローの組み合わせが増える。   The condition of Equation (26) in Theorem 1 (i) is a condition in which the estimation possibility is determined only from the luminance of the image, and is equivalent to the regularity of the Hessian matrix of luminance known conventionally. On the other hand, the condition of equation (28) in (ii) is a new condition that has not been known so far, because the possibility of estimation is determined by the temporal change of the luminance gradient of the image. When this condition (28) is added as a sufficient condition, the estimable condition is relaxed compared to the conventional method, and the estimable situation, that is, the combination of luminance and optical flow increases.

さらに、この式(28)の条件の特筆すべき点として、輝度のヘッセ行列Ｈが零行列であったとしても推定可能な場合がある。すなわち、式(29)は、式(28)の行列を具体的に書き下したものである。この行列は、Ｈ＝Ｏであったとしても、式(30)の場合を除いて一般にランクは２となり、オプティカルフローの推定が行える。式(30)は単純な拡大縮小の場合のオプティカルフロー（特に、並進はα＝０）を表しているので、単純拡大縮小以外の場合は、式(28)の条件によってもオプティカルフローが推定できることを示している。実際、後に数値実験により示すように、画像全体でＨ＝Ｏであるような画像でも、ある点を中心として一定角速度で回転するような場合には、この推定モデルを用いるとオプティカルフローの推定が可能である。この事実は、従来の推定法（非特許文献５参照）ではＨ≠Ｏであってもdet（H）＝０となるだけで推定ができなくなるのとは大きく異なり、この推定モデルの頑健性を示しているといえる。   Further, as a special point of the condition of the equation (28), there are cases where it is possible to estimate even if the luminance Hessian matrix H is a zero matrix. In other words, Equation (29) is a concretely written down matrix of Equation (28). Even if this matrix is H = O, the rank is generally 2 except for the case of Expression (30), and the optical flow can be estimated. Since Equation (30) represents the optical flow in the case of simple scaling (especially translation is α = 0), the optical flow can be estimated by the condition of Equation (28) in cases other than simple scaling. Is shown. In fact, as will be shown later by numerical experiments, even if an image with H = O in the whole image is rotated at a constant angular velocity around a certain point, this estimation model can be used to estimate the optical flow. Is possible. This fact is greatly different from the conventional estimation method (see Non-Patent Document 5), in which even if H ≠ O, estimation cannot be performed only by det (H) = 0, the robustness of this estimation model is reduced. It can be said that it shows.

さらに、輝度の時空間勾配からオプティカルフローとその空間勾配を求める具体的な計算式を示す。微分形式が列フルランクになるために必要な微分が２回であったことから、瞬時的に推定を行うためには、観測量ｚ(t)の時間に関する２階微分までの情報が必要となる。自明に決定できるｚ(t)を除き、（ｄ／ｄｔ）z(t)，（ｄ2／ｄｔ2）z(t)について式(17)を具体的に書き下すと式(31)となる。なお、式(31)は、請求項の(eq2)に対応する。ここで、式(21)，(22)で定義される行列Ｊ，Ｊtの要素はすべて輝度の時空間勾配であり、観測量ｚ(t)から求めることが可能である。したがって、式(31)の方程式は、推定すべき変数に関して線形方程式となっている。この過制約の線形方程式を、例えば最小二乗法などによって解くことによって、オプティカルフローだけでなく、その空間勾配の瞬時推定も可能である。Furthermore, the specific calculation formula which calculates | requires an optical flow and its spatial gradient from the spatiotemporal gradient of a brightness | luminance is shown. Since the differentiation required for the differential form to be full rank is twice, in order to estimate instantaneously, information up to the second derivative with respect to the time of the observed quantity z (t) is required. Become. Excluding z (t), which can be determined in a self-evident manner, equation (17) is specifically written as (d / dt) z (t), (d 2 / dt 2 ) z (t), and equation (31) is obtained. . Equation (31) corresponds to (eq2) in the claims. Here, equation (21), can be determined from a defined by the matrix J in (22), the elements of J t is the spatial gradients at all brightness, observed quantity z (t). Therefore, the equation (31) is a linear equation with respect to the variable to be estimated. By solving this over-constrained linear equation by, for example, the least square method, it is possible to instantaneously estimate not only the optical flow but also its spatial gradient.

また、実装の際に必要とされる画像データは、輝度の２次の空間勾配を計算する必要があるので、各時刻において３画素×３画素以上必要である。さらに、輝度の空間勾配の時間に関する２次微分まで必要であるので、時系列方向に３枚以上の画像が必要である。つまり、式(25)に従ってオプティカルフローとその空間勾配を推定するために必要とされる画像サイズは、少なくとも３画素(ｘ方向)×３画素(ｙ方向)×３(ｔ方向)である。   Moreover, since it is necessary to calculate the secondary spatial gradient of the luminance, the image data required for mounting requires 3 pixels × 3 pixels or more at each time. Furthermore, since it is necessary up to the second derivative with respect to the time of the spatial gradient of luminance, three or more images are required in the time series direction. In other words, the image size required to estimate the optical flow and its spatial gradient according to Equation (25) is at least 3 pixels (x direction) × 3 pixels (y direction) × 3 (t direction).

図２を参照して、実際の計算手法の一例を説明する。時刻ｔ、座標（ｘ，ｙ）における輝度をＩ(x,y,t)とし、時刻ｔ0における注目点（ｘ0，ｙ0）におけるオプティカルフロー［ｕ ｖ］Tとその空間勾配を求める。An example of an actual calculation method will be described with reference to FIG. The luminance at time t and coordinates (x, y) is I (x, y, t), and the optical flow [u v] T and its spatial gradient at the point of interest (x 0 , y 0 ) at time t 0 are obtained.

まず、撮像部３が、３画素以上×３画素以上の画像を、時間を異にして３枚以上撮影して撮像画像とし（ステップＳＴ１）、量子化部７が撮像画像上の各点について、時間又は／及び輝度の量子化を行って計測画像とし（ステップＳＴ２）、平滑化部１１が計測画像に平滑化を行い平滑化画像とする（ステップＳＴ３）。ここで、具体的に説明するために、用いるデータは、I(x,y,t)のうち、ｘがｘ0−１、ｘ0、ｘ0＋１、ｙがｙ0−１、ｙ0、ｙ0＋１、ｔがｔ0−２、ｔ0−１、ｔ0を組み合わせた２７画素の輝度であるとする。First, the imaging unit 3 captures three or more images of three pixels or more × 3 pixels or more at different times to form a captured image (step ST1), and the quantization unit 7 Time or / and luminance are quantized to obtain a measurement image (step ST2), and the smoothing unit 11 smoothes the measurement image to obtain a smoothed image (step ST3). Here, in order to explain specifically, the data to be used is I (x, y, t), where x is x 0 -1, x 0 , x 0 +1, y is y 0 -1, y 0 , It is assumed that y 0 +1, t is the luminance of 27 pixels combining t 0 -2, t 0 -1, and t 0 .

まず、輝度勾配演算部２１は、輝度勾配を求める（ステップＳＴ４）。計算で求めたデータは、輝度自体も得られたデータから一番確からしい値に修正するため、観測データと区別するためＧと表す。これが、式(32)である。なお、式(32a)、(32b)、(32c)、(32d)、(32e)、(32f)についてはｔ＝ｔ0−２、ｔ0−１、ｔ0である。First, the luminance gradient calculation unit 21 obtains a luminance gradient (step ST4). The data obtained by calculation is expressed as G in order to distinguish it from observation data in order to correct the luminance itself from the obtained data to the most probable value. This is Equation (32). Note that t = t 0 -2, t 0 -1, and t 0 for the expressions (32a), (32b), (32c), (32d), (32e), and (32f).

次に、時間勾配演算部２３は、このデータから時間勾配を求める（ステップＳＴ５）。これが、式(33)である。なお、式(33a)、(33b)、(33c)、(33g)、(33h)、(33i)についてはｔ＝ｔ0−１、ｔ0である。そして、方程式演算部２５は、オプティカルフローとその空間勾配を式（34）の連立方程式を解くことにより求める（ステップＳＴ６）。画像処理装置１３は、求めたオプティカルフロー及びその空間勾配を、オプティカルフロー記憶部１５に保存する（ステップＳＴ７）。Next, the time gradient calculation unit 23 obtains a time gradient from this data (step ST5). This is Equation (33). It should be noted that t = t 0 -1, t 0 for the expressions (33a), (33b), (33c), (33g), (33h), and (33i). Then, the equation calculation unit 25 obtains the optical flow and its spatial gradient by solving the simultaneous equations of the equation (34) (step ST6). The image processing device 13 stores the obtained optical flow and its spatial gradient in the optical flow storage unit 15 (step ST7).

なお、代表的な仮定としては、時不変性の仮定以外に、並進性の仮定（すなわち、式(35)が成立するとの仮定）がある。提案したものは、時不変性は仮定し、並進性は仮定しないものである。そこで、これらの各仮定をした場合としない場合の組み合わせについて、提案のものと比較する。   As typical assumptions, there is a translational assumption (that is, an assumption that Equation (35) holds) in addition to the time-invariant assumption. The proposed one assumes time invariance and does not assume translation. Therefore, the combination with and without each of these assumptions is compared with the proposed one.

まず、時不変性の仮定だけでなく、並進性も仮定した場合に対して提案したものと比較した場合、提案したものは、定理１において説明したように、推定可能条件が緩和され、推定可能な状況、つまり輝度とオプティカルフローの組み合わせが増える。   First, when compared with the proposed case for not only the assumption of time invariance but also for the translational nature, the proposed one has a relaxed estimable condition and can be estimated as explained in Theorem 1. Situation, that is, the combination of luminance and optical flow increases.

次に、時不変性の仮定をせず、並進性の仮定をした場合について検討する。この推定モデルは、比較的時間変動の大きなオプティカルフローの推定に適したものといえる。   Next, let us consider the case of translational assumptions without the assumption of time invariance. This estimation model can be said to be suitable for estimation of optical flows with relatively large temporal fluctuations.

画像上のある点Ｐ(x0,y0)における輝度、オプティカルフロー及びその空間勾配を用いて、状態ｑ(t)を式(36)として、オプティカルフローの推定モデルを式(37)のようにおく。オプティカルフローの推定モデル(36)の非線形システムとしての可観測性を調べると、推定モデル(37)が可観測となるための必要十分条件は式(38)となる。また、微分形式が列フルランクになるために必要な微分が２回であったことから、瞬時的に推定を行うためには、観測量ｚ(t)の時間に関する２階微分までの情報が必要である。これまでの推定モデルと同様に推定するための代数方程式を求めると式(39)となる。この（過制約の）線形方程式を解くことによって（誤差や雑音の影響を無視すれば）瞬時推定を行うことが可能である。オプティカルフロー［ｕ ｖ］Tだけであれば、上式の第１〜３行目の線形方程式を解くことでも推定が行えるが、第４行目の制約式によって［ｕt ｖt］Tの情報を［ｕ ｖ］Tの推定に反映させることができる。一方、この推定モデルの可観測性条件(38)は、従来の可観測性条件よりも厳しい条件である。これは時不変性の仮定を設けなかったことの代償である。さらに、一般に推定モデルの状態数の増加は推定精度の劣化、計算時間の増加に結びつくので、この推定モデルが従来の推定モデルより有効であるかどうかについては注意が必要である。Using the brightness at the point P (x 0 , y 0 ) on the image, the optical flow, and its spatial gradient, the state q (t) is expressed as equation (36), and the optical flow estimation model is expressed as equation (37) Put it in. When the observability of the optical flow estimation model (36) as a nonlinear system is examined, the necessary and sufficient condition for the estimation model (37) to be observable is Equation (38). In addition, since the differentiation required for the differential form to be a full rank column was twice, in order to estimate instantaneously, the information up to the second derivative with respect to the time of the observable z (t) is necessary. When an algebraic equation for estimation is obtained in the same manner as the previous estimation models, Equation (39) is obtained. It is possible to perform instantaneous estimation (by ignoring the influence of errors and noise) by solving this (over-constrained) linear equation. If only the optical flow [u v] T is used, estimation can also be performed by solving the linear equations in the first to third rows of the above equation, but information on [u t v t ] T by the constraint equation in the fourth row. Can be reflected in the estimation of [u v] T. On the other hand, the observability condition (38) of this estimation model is stricter than the conventional observability condition. This is the price of not making the time invariance assumption. Furthermore, since an increase in the number of states of an estimation model generally leads to a decrease in estimation accuracy and an increase in calculation time, it is necessary to pay attention to whether this estimation model is more effective than a conventional estimation model.

なお、時不変性の仮定と並進性の仮定のいずれか一方のみを仮定した発明は、次のように特定される。すなわち、物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定装置において、前記物体が撮像された少なくとも３枚の撮像画像に対し、時間に関して２階微分した情報を用いてオプティカルフローという情報を画像解析処理して推定する画像解析処理手段を備えた、物体運動装置、である。なお、時不変性と並進性を共に仮定したものは１階微分した情報を用い、時不変性と並進性を共に仮定しないものは３階微分した情報を用いる。   The invention that assumes only one of the time-invariant assumption and the translational assumption is specified as follows. That is, in an object motion estimation apparatus that estimates the motion of the object by performing image analysis processing on the captured image obtained by capturing the object, information obtained by performing second-order differentiation with respect to time with respect to at least three captured images obtained by capturing the object. It is an object motion apparatus provided with an image analysis processing means that uses an image analysis process to estimate information called optical flow. Note that information obtained by first-order differentiation is used for both assuming time invariance and translation, and information obtained by third-order differentiation is used for not assuming both time invariance and translation.

次に、時不変性も並進性も仮定しない、最も一般的な推定モデルでは、オプティカルフローが時不変かつ単純拡大・縮小という基本的な状況で推定モデルが輝度Ｉによらず不可観測になる。他の推定モデルでは、輝度に対して可観測であるための制約が課されていた。これに対して、この推定モデルは、輝度に関わらずオプティカルフローの時空間勾配に対して条件が課される点で大きく異なる。実質的にこの推定モデルでは、オプティカルフローの推定ができないことを意味している。   Next, in the most general estimation model that assumes neither time-invariance nor translation, the estimation model becomes unobservable regardless of the luminance I in a basic situation where the optical flow is time-invariant and simple enlargement / reduction. Other estimation models impose restrictions on observability of brightness. On the other hand, this estimation model is greatly different in that a condition is imposed on the spatiotemporal gradient of the optical flow regardless of the luminance. In effect, this estimation model means that optical flow cannot be estimated.

このように、時不変性を仮定しない推定モデルは、時不変性を仮定しないモデルに比べ可観測となるための条件が厳しく、実画像への適用には注意が必要である。   As described above, an estimation model that does not assume time invariance has stricter conditions for observability than a model that does not assume time invariance, and attention must be paid to application to an actual image.

続いて、カルマンフィルタによる推定について説明する。瞬時推定では、輝度の時空間微分の陽な計算が必要であり、画像情報からこれらの情報を得るためには差分演算が不可避であった。この差分演算は外乱に弱く、オプティカルフローの推定の信頼性に大きく影響する。   Next, estimation by the Kalman filter will be described. The instantaneous estimation requires explicit calculation of the spatiotemporal derivative of luminance, and the difference calculation is inevitable in order to obtain such information from the image information. This difference calculation is vulnerable to disturbances and greatly affects the reliability of optical flow estimation.

一方、カルマンフィルタは、雑音が加わる動的システムの出力の時系列データから逐次的にシステムの状態を推定する。カルマンフィルタによる推定は、時間解像度は犠牲にするものの、時間微分の陽な計算は不要となり、式(31)の線形方程式を解く際の逆行列計算も不要となり、外乱に強く精度の高い推定が可能となる。   On the other hand, the Kalman filter sequentially estimates the state of the system from the time series data of the output of the dynamic system to which noise is added. Although estimation by the Kalman filter sacrifices time resolution, explicit calculation of time differentiation is not required, and inverse matrix calculation is not required when solving the linear equation of Equation (31). It becomes.

オプティカルフローの推定にカルマンフィルタを用い推定精度の向上を試みることは既に行われている（非特許文献７など参照）。しかし、これらの手法は、あくまでも瞬時的な方法に基づいて求まったオプティカルフローの推定値に、時間方向の滑らかさを要請するためにカルマンフィルタを用いているだけであり、式(31)の線形方程式を解くという意味で本質的には瞬時的な推定方法となんら変わらない。   Attempts have been made to improve estimation accuracy using a Kalman filter for optical flow estimation (see Non-Patent Document 7, etc.). However, these methods only use the Kalman filter to request the smoothness in the time direction to the optical flow estimation value obtained based on the instantaneous method, and the linear equation of Equation (31) is used. Is essentially the same as the instantaneous estimation method.

そこで、提案した推定モデル(11)にカルマンフィルタを適用し、オプティカルフローを推定する方法について説明する。   Therefore, a method for estimating an optical flow by applying a Kalman filter to the proposed estimation model (11) will be described.

まず、空間勾配を用いない推定モデルについて説明する。カルマンフィルタの特徴は、時間微分を用いる必要がないことである。一方、提案した推定モデルの観測量は、輝度の空間微分を含んでいる。そのため、これらの推定モデルをそのまま用いては、空間微分の計算は不可欠である。   First, an estimation model that does not use a spatial gradient will be described. A feature of the Kalman filter is that it is not necessary to use time differentiation. On the other hand, the observation amount of the proposed estimation model includes a spatial differential of luminance. Therefore, using these estimation models as they are, the calculation of spatial differentiation is indispensable.

そこで、推定モデルの観測方程式を修正し、カルマンフィルタによる逐次推定において空間微分の計算を不要となる推定モデルを提案する。具体的には、観測量として注目点Ｐ(x0,y0)の空間勾配の代わりに、注目点Ｐ(x0,y0)の空間近傍の輝度Ｉを用いることを考え、新たな観測量〜ｚ(t)として、式(40)という、注目点Ｐ(x0,y0)の３×３近傍の輝度を考える。Therefore, we propose an estimation model that corrects the observation equation of the estimation model and eliminates the need to calculate spatial differentials in successive estimation using the Kalman filter. Specifically, instead of the spatial gradient of the target point P (x 0, y 0) as the observation quantity, consider the use of brightness I of the spatial vicinity of the target point P (x 0, y 0), a new observation As a quantity to z (t), a luminance in the vicinity of 3 × 3 of the point of interest P (x 0 , y 0 ), which is expressed by Equation (40), is considered.

近傍点の輝度Ｉ(x0+Δx,y0+Δy)（Δx,Δy＞0）は式(41)と近似できるので、〜ｚ(t)は式(42)のように近似できる。ただし、ｎは状態ｑの次数（又は要素数）である。ここで、式(42b)の行列〜Ｈのランクが常に６であり、式(43)によって〜ｚ(t)からｚ(t)が復元できる。このことから、システムの可観測性が維持されることが示せる。このことより、空間勾配の代わりに近傍点の輝度情報を観測量とすることでカルマンフィルタが構成できる。Luminance I in the vicinity of point (x0 + Δx, y0 + Δy ) (Δx, Δy> 0) so can be approximated as equation (41), ~ z (t) can be approximated as equation (42). However, n is the order (or the number of elements) of the state q. Here, the rank of the matrix ~ H in the expression (42b) is always 6, and z (t) can be restored from ~ z (t) by the expression (43). This shows that the observability of the system is maintained. Thus, a Kalman filter can be configured by using the luminance information of neighboring points as an observation amount instead of the spatial gradient.

なお、ここでは簡単のため３×３近傍を考えたが、一般のｋ×ｌ近傍（ただし、３×３以上）や矩形以外の近傍に対しても式(42b)の〜Ｈは簡単に計算でき、そのランクが６であればカルマンフィルタが構成できる。Note that although 3 × 3 neighborhoods are considered here for simplicity, ~ H in Equation (42b) is easily calculated for general k × l neighborhoods (however, 3 × 3 or more) and neighborhoods other than rectangles. If the rank is 6, a Kalman filter can be configured.

次に、図１の画像解析処理部１３の状態推定処理部２７による、カルマンフィルタを用いた推定処理について説明する。雑音を付加した式(44)のモデルを考える。ただし、ｆ(q)，〜Ｈは式(11c)，(42b)で与えられ、システム雑音ｕ(t)，観測雑音ｗ(t)は、それぞれ共分散行列がＱ，Ｒのガウス雑音と仮定する。この枠組みにおいて、カルマンフィルタによる逐次推定を行えばオプティカルフローの推定が行える。Next, estimation processing using the Kalman filter by the state estimation processing unit 27 of the image analysis processing unit 13 in FIG. 1 will be described. Consider the model of equation (44) with added noise. However, f (q), ~ H are given by equations (11c) and (42b), and system noise u (t) and observation noise w (t) are assumed to be Gaussian noises with covariance matrices of Q and R, respectively. To do. In this framework, optical flow can be estimated by performing successive estimation using a Kalman filter.

オプティカルフローの推定モデルは非線形システムであり、その非線形性も強い。そこで、非線形性を線形近似する拡張カルマンフィルタではなく、推定値の共分散を近似する無香料カルマンフィルタ（unscented Kalman filter）（非特許文献２参照）を用いて状態推定を行い、オプティカルフローを推定することを提案する。無香料カルマンフィルタの特徴として、システム外乱を状態に含めて推定を行うことにより、推定精度を向上させることができること（非特許文献６参照）、時間方向の離散化において、ルンゲ・クッタ法などの高次の積分公式を用いることが可能で、これによっても推定精度を向上させることができることが挙げられる。これらは推定に必要な計算時間とのトレードオフになるが、システム外乱を状態に含めたモデルを採用し、数値積分公式は非線形性を考慮して２次以上の積分公式を用いるべきである。   The estimation model of optical flow is a nonlinear system, and its nonlinearity is strong. Therefore, state estimation is performed by using an unscented Kalman filter (see Non-Patent Document 2) that approximates the covariance of the estimated value instead of an extended Kalman filter that linearly approximates nonlinearity, and optical flow is estimated. Propose. As a feature of the unscented Kalman filter, the estimation accuracy can be improved by performing the estimation including the system disturbance (see Non-Patent Document 6), and in the discretization in the time direction, such as the Runge-Kutta method It is possible to use the following integration formula, which can also improve the estimation accuracy. These are trade-offs with the calculation time required for estimation, but a model including system disturbance is adopted, and a numerical integration formula should use a quadratic or higher order integration formula in consideration of nonlinearity.

以上をまとめると、本実施例では、輝度不変の仮定と局所的な輝度情報からのみによるオプティカルフロー推定可能性を動的システムの可観測性の立場から示した。オプティカルフローに対する並進性の仮定を緩めるために、輝度不変の拘束式の拡張を行った。導出した拘束式は、画像上の注目点の近傍において輝度不変であるための条件式をテーラー展開し、その各項の係数が０となるための条件である。そのように解釈すると、従来の輝度不変の拘束式はテーラー展開の０次の項に相当することが分かる。   In summary, in this embodiment, the possibility of optical flow estimation based only on the assumption of luminance invariance and local luminance information is shown from the standpoint of observability of a dynamic system. In order to relax the translational assumption for optical flow, we extended the invariant constraint with brightness. The derived constraint formula is a condition for the Taylor expansion of a conditional formula for invariance of luminance in the vicinity of the point of interest on the image, and for the coefficient of each term to be zero. When interpreted in that way, it can be seen that the conventional luminance invariant constraint equation corresponds to the 0th-order term of the Taylor expansion.

さらに、本実施例では、オプティカルフローに対する並進性の仮定を緩めるだけでなく、オプティカルフローとその空間勾配を同時に推定するという立場から、動的システムとしての推定モデルを提案した。そして、動的システムの可観測性の立場から、推定可能性に対して考察を行い、輝度のヘッセ行列が零行列であってもオプティカルフローとその空間勾配の推定が行える場合があることを示した。   Furthermore, in this embodiment, an estimation model as a dynamic system was proposed from the standpoint of not only relaxing the assumption of translational property for optical flow but also estimating optical flow and its spatial gradient simultaneously. And from the standpoint of observability of dynamic systems, we consider the possibility of estimation and show that even if the luminance Hessian is zero, the optical flow and its spatial gradient may be estimated. It was.

また、ヘッセ行列などの逆行列計算や、時空間勾配の陽な計算を排除するために、推定モデルの出力方程式を修正し、カルマンフィルタを用いた数値的に安定な推定法を提案した。   In addition, in order to eliminate inverse matrix calculations such as the Hessian matrix and explicit calculation of the spatiotemporal gradient, the output equation of the estimation model was modified and a numerically stable estimation method using a Kalman filter was proposed.

なお、カルマンフィルタは一例であり、他の逐次推定手法を用いてもよい。   Note that the Kalman filter is an example, and other sequential estimation methods may be used.

また、例えば、テーラー展開のより高次の項まで用いた条件式に基づく推定モデルを導出したり、オプティカルフローの時変性へ対応したりするようにしてもよい。   In addition, for example, an estimation model based on a conditional expression that uses higher-order terms of Taylor expansion may be derived, or optical flow time-varying may be dealt with.

本実施例では、具体的な例を通して、提案した推定モデルの有効性を示す。   In this embodiment, the effectiveness of the proposed estimation model is shown through a specific example.

元画像の輝度Ｉ(x,y,t)を式(45)のように与える。この元画像は、輝度の空間勾配が平坦な（ヘッセ行列が零行列の）画像を、点(0,0)を中心に一定角速度ω（≠０）で回転させたものである。この元画像の座標(x,y)におけるオプティカルフロー及びその空間勾配は、時不変性の仮定の元では一意に決まり、式(46)である。一方、輝度の時空間勾配は式(47)であり、輝度のヘッセ行列が零行列となるので、従来の局所推定法ではオプティカルフローの推定は行えない。   The luminance I (x, y, t) of the original image is given as in the formula (45). This original image is obtained by rotating an image with a flat luminance spatial gradient (the Hessian matrix is zero) at a constant angular velocity ω (≠ 0) around the point (0,0). The optical flow and its spatial gradient at the coordinates (x, y) of the original image are uniquely determined under the assumption of time invariance, and are given by Equation (46). On the other hand, the spatio-temporal gradient of the luminance is Equation (47), and the luminance Hessian matrix is a zero matrix, so that the optical flow cannot be estimated by the conventional local estimation method.

次に、理論計算について説明する。ここでは、式(45)の元画像の輝度及びその空間勾配から、オプティカルフローが理論的に計算可能かどうかを検討する。   Next, theoretical calculation will be described. Here, it is examined whether the optical flow can be theoretically calculated from the luminance of the original image in Equation (45) and its spatial gradient.

まず、Hornら（1981）の大域法を適用する。Hornら（1981）では、ある領域Ｄについて式(48)を最小化する。ただし、γ＞０は平滑化の強さを調整するパラメータである。領域Ｄ＝{(x,y)|-a≦x≦a，-b≦y≦b}（a,b＞0），時刻は一般性を失わずｔ＝０を考える。式(47a)、(47b)より、ｔ＝０で上記の評価関数は式(49)となり、ｕ，γの選び方に関わらずｖx＝ｖy＝０、つまり、ｖ＝const.が解となる。このことは、ｖ（x,y,0）が一意に決定できないこと、決定したとしても真値(46)とは異なることを意味する。また、ｕ(x,y,0)についても、領域Ｄが仮定のように矩形の場合はｕx＝０となるようなｕのうち式(50)を最小化するものを選べばよい。変分法を用いて解くと式(51)が最小解となる。このｕ(x,y,0)は滑らかさの重みγがγ→０で真値に一致するものの、γ＞０である限りはやはり真値とは一致しない。以上より、この数値例ではHornら（1981）の大域法を用いてもオプティカルフローが一意に定められないこと、定められたとしても真値とは一致しないことがわかる。First, the global method of Horn et al. (1981) is applied. Horn et al. (1981) minimize Equation (48) for a region D. However, γ> 0 is a parameter for adjusting the strength of smoothing. Region D = {(x, y) | -a≤x≤a, -b≤y≤b} (a, b> 0), and time is considered not to lose generality t = 0. From the equations (47a) and (47b), the above evaluation function becomes the equation (49) when t = 0, and v x = v y = 0, that is, v = const. Become. This means that v (x, y, 0) cannot be uniquely determined, and even if it is determined, it is different from the true value (46). As for u (x, y, 0), if the region D is rectangular as assumed, u may be selected so that u x = 0, which minimizes the expression (50). When solving using the variational method, Equation (51) becomes the minimum solution. Although u (x, y, 0) matches the true value when the smoothness weight γ is γ → 0, it does not match the true value as long as γ> 0. From the above, it can be seen that in this numerical example, the optical flow cannot be uniquely determined even if the global method of Horn et al. (1981) is used, and even if it is determined, it does not match the true value.

一方、提案した推定モデルの瞬時推定式(31)は、式(47c)に注意すると式(52)となる。ここで、両辺とも０となる自明な方程式は除くとともに、方程式の順番を入れ替えている。上式は３組の２次元線形連立代数方程式となり、式(47a)，(47b)を代入することで式(53)の３つの式と等価になる。いずれの式もｔによらず常に一意に解が定まり、式(46)の理論値と一致する。この事実は、大域法に対する提案した推定モデルの優位性を示している。   On the other hand, the instantaneous estimation equation (31) of the proposed estimation model becomes equation (52) when attention is paid to equation (47c). Here, the obvious equations that are 0 on both sides are excluded, and the order of the equations is changed. The above equation becomes three sets of two-dimensional linear simultaneous algebraic equations, and is equivalent to the three equations of equation (53) by substituting equations (47a) and (47b). In any equation, the solution is always uniquely determined regardless of t and coincides with the theoretical value of equation (46). This fact shows the superiority of the proposed estimation model over the global method.

次に、数値計算について説明する。ここでは、式(45)の元画像のパラメータをβ＝0.2、ω＝π／100、Ｉ0＝128とし、時空間及び輝度の量子化を行った計測画像に対して提案手法を適用し、その有効性を示す。Next, numerical calculation will be described. Here, the proposed method is applied to a measurement image in which the spatio-temporal and luminance are quantized with the parameters of the original image of Equation (45) being β = 0.2, ω = π / 100, and I 0 = 128, The effectiveness is shown.

具体的には、図１の量子化部７が、式(45)で与えられる元画像を撮像画像として、ｔが整数の時刻において(x,y)が整数の点で小数点以下を切り捨てて量子化を行い、これを計測画像とする。そして、図１の平滑化部１１が、計測画像の各時刻において、σ＝３のガウスフィルタを用いて平滑化を行う。画像解析処理部１５は、この平滑化画像の輝度情報を用いてオプティカルフローの推定を行う。   Specifically, the quantization unit 7 in FIG. 1 uses the original image given by Equation (45) as a captured image, and at the time when t is an integer, (x, y) is rounded down to the integer and the quantum is cut off. This is used as a measurement image. Then, the smoothing unit 11 in FIG. 1 performs smoothing using a Gaussian filter with σ = 3 at each time of the measurement image. The image analysis processing unit 15 estimates the optical flow using the luminance information of the smoothed image.

まず、提案した推定モデルを用いて瞬時推定を行った結果を示す。図３は、座標(60,80)におけるオプティカルフロー［ｕ ｖ］Tの推定値を400フレーム（２回転）分示している。量子化誤差の影響で推定値に大きな誤差が加わっており、推定がうまくできているとは言い難い。これは元画像が非常に特徴の少ない画像であるため、量子化誤差程度でも逆行列計算に大きく影響するためと考えられる。First, the results of instantaneous estimation using the proposed estimation model are shown. FIG. 3 shows an estimated value of the optical flow [u v] T at coordinates (60, 80) for 400 frames (two rotations). A large error is added to the estimated value due to the influence of the quantization error, and it is difficult to say that the estimation is successful. This is presumably because the original image is an image with very few features, and even the degree of quantization error greatly affects the inverse matrix calculation.

次に、量子化誤差に対応するため、カルマンフィルタを用いて推定を行った結果を示す。推定は提案した推定モデルである式(11)に、式(42)による出力方程式の修正を施したモデルを用いる。このモデルに対して雑音が印可された式(44)を考え、無香料カルマンフィルタを用いて推定を行う。ただし、式(44)における雑音の共分散行列はそれぞれ式(54)、式(55)とし、カルマンフィルタで用いる状態の推定値の初期値はすべて０としている。また、オプティカルフローのモデルの数値積分は古典的ルンゲ・クッタの４次の公式を用い、かつ無香料カルマンフィルタのモデルとしては、状態にシステム雑音を含め拡張した推定モデル（非特許文献６参照）で実験を行った。   Next, the results of estimation using a Kalman filter in order to deal with quantization errors will be shown. For the estimation, a model obtained by correcting the output equation according to Equation (42) to Equation (11), which is the proposed estimation model, is used. Considering Equation (44) in which noise is applied to this model, estimation is performed using an unscented Kalman filter. However, the noise covariance matrices in Equation (44) are Equation (54) and Equation (55), respectively, and the initial values of the estimated values used in the Kalman filter are all zero. The numerical integration of the optical flow model uses the classic Runge-Kutta fourth-order formula, and the unscented Kalman filter model is an estimated model that includes system noise in the state (see Non-Patent Document 6). The experiment was conducted.

図４は、座標(60, 80)におけるオプティカルフロー［ｕ ｖ］Tとその空間勾配ｕx、ｕy、ｖx、ｖyの推定値を400フレーム（２回転）分示している。初期値が大きくずれているため収束まで時間がかかっているものの、精度良く推定できていることが分かる。また、オプティカルフローの空間勾配も同時に推定できていることが確認できる。これより、提案した推定モデル及びカルマンフィルタを用いた推定法の有効性を示すことができた。FIG. 4 shows the estimated values of the optical flow [u v] T at the coordinates (60, 80) and the spatial gradients u x , u y , v x , and v y for 400 frames (two rotations). It can be seen that although the initial value is greatly deviated, it takes a long time to converge, but it can be estimated accurately. It can also be confirmed that the spatial gradient of the optical flow can be estimated at the same time. Thus, the effectiveness of the proposed estimation model and the estimation method using the Kalman filter could be demonstrated.

オプティカルフローは、画像による物体の運動計測に用いることができる汎用性の高い技術である。具体的には、例えば、自動車の車載カメラによる移動体の運動計測がある。車載カメラは、基本的に進行方向に向けられるので光軸に沿った動きが避けられない。また、ドライブレコーダのような広角画像では、画像の歪みによる拡大が避けられない。本願発明によれば、このような場合にも高精度でオプティカルフローが計算できる。さらに、移動体の一般的な検出も可能であり、例えば、車載カメラによる飛び出しや追い越し車両の警告、建物における侵入者検出、運動するカメラの座標補正なども可能である。   Optical flow is a highly versatile technique that can be used to measure the movement of an object using images. Specifically, for example, there is a motion measurement of a moving body by an in-vehicle camera of an automobile. Since the in-vehicle camera is basically directed in the traveling direction, movement along the optical axis is inevitable. Further, in a wide-angle image such as a drive recorder, enlargement due to image distortion is inevitable. According to the present invention, the optical flow can be calculated with high accuracy even in such a case. Furthermore, a general detection of a moving body is also possible. For example, it is possible to use a vehicle-mounted camera to jump out, warn an overtaking vehicle, detect an intruder in a building, or correct the coordinates of a moving camera.

また、既に組みあがってしまった回転機械の回転速度を、回転軸に何かを設置することなく測定することもできる。従来の角速度センサは回転軸に取り付けるものが主流である。本願発明によれば、必ずしも回転中心を含まない回転体の時系列画像から回転角速度及び回転の中心までの距離を計測することが可能になる。   It is also possible to measure the rotational speed of a rotating machine that has already been assembled without installing anything on the rotating shaft. Conventional angular velocity sensors are mainly attached to a rotating shaft. According to the present invention, it is possible to measure the rotational angular velocity and the distance from the rotation center to the time-series image of the rotating body that does not necessarily include the rotation center.

さらに、同様に物体の運動の並進性が仮定できない、流体計測（特に渦、乱流）などにも有用である。   Furthermore, it is also useful for fluid measurement (especially vortex and turbulent flow) where the translational motion of the object cannot be assumed.

１ 物体運動推定装置、１３ 画像解析処理部、２１ 輝度勾配演算部、２３ 時間勾配演算部、２５ 方程式演算部、２７ 状態推定処理部   DESCRIPTION OF SYMBOLS 1 Object motion estimation apparatus, 13 Image analysis processing part, 21 Luminance gradient calculating part, 23 Time gradient calculating part, 25 Equation calculating part, 27 State estimation processing part

Claims (8)

撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定装置において、
前記撮像画像の注目点における輝度勾配を演算する輝度勾配演算手段と、
前記輝度勾配から時間勾配を演算する時間勾配演算手段と、
連立方程式を解く方程式演算手段とを有する
画像解析処理手段を備え、
前記連立方程式は、
前記撮像手段と前記物体との相対運動の速度分布であるオプティカルフローについて、前記撮像画像における前記物体の輝度の時間に関する１階全微分がゼロという時不変性が仮定されることにより、前記オプティカルフローを変数とし、前記撮像画像の画素の輝度の空間勾配及びその時間微分を係数とするだけでなく、さらに、
前記撮像手段の撮像軸に対して垂直な平面上を前記物体が移動するという並進性が仮定されないことにより、前記オプティカルフローの空間勾配をも変数とし、前記輝度及びその空間勾配の時間に関して２階微分した値をも係数とするものであり、
前記画像解析処理手段が有する方程式演算手段が、前記連立方程式を解くことにより、前記オプティカルフローを推定すると同時に、前記推定されるオプティカルフローの空間勾配をも推定する物体運動推定装置。 In an object motion estimation apparatus for estimating motion of an object by performing image analysis processing on a captured image obtained by imaging the object by an imaging means,
A luminance gradient calculating means for calculating a luminance gradient at a point of interest of the captured image;
Time gradient calculating means for calculating a time gradient from the luminance gradient;
Image analysis processing means having equation calculation means for solving simultaneous equations,
The simultaneous equations are
The optical flow, which is the velocity distribution of relative motion between the imaging means and the object, is assumed to be invariant when the first-order total derivative with respect to the time of luminance of the object in the captured image is assumed to be zero. As a variable, the spatial gradient of the luminance of the pixel of the captured image and its time derivative as a coefficient,
Since the translational property that the object moves on a plane perpendicular to the imaging axis of the imaging means is not assumed, the spatial gradient of the optical flow is also a variable, and the second order with respect to the luminance and the time of the spatial gradient. The differentiated value is also used as a coefficient.
An object motion estimation apparatus in which an equation calculation unit included in the image analysis processing unit estimates the optical flow and solves the spatial gradient of the estimated optical flow by solving the simultaneous equations. 前記撮像画像は、二次元の画像であり、少なくとも３画素×３画素の大きさの画像であって、撮像時刻が異なる３枚の画像であり、撮像画像記憶手段に記憶されており、
前記輝度勾配演算手段は、前記撮像画像の注目点Ｐ（ｘ0，ｙ0）における輝度Ｉの輝度勾配Ｉx（ｘ0，ｙ0，ｔ）、Ｉy（ｘ0，ｙ0，ｔ）、Ｉxx（ｘ0，ｙ0，ｔ）、Ｉxy（ｘ0，ｙ0，ｔ）及びＩyy（ｘ0，ｙ0，ｔ）を演算して、式(eq1)に示される観測量ｚ（ｔ）を求めるものであり、
前記時間勾配演算手段は、前記輝度Ｉ及び前記輝度勾配Ｉx、Ｉy、Ｉxx、Ｉxy及びＩyyから時間勾配Ｉt、Ｉtx、Ｉty、Ｉtxx、Ｉtxy、Ｉtyy、Ｉtt、Ｉttx、Ｉtty、Ｉttxx、Ｉttxy及びＩttyyを演算するものであり、
前記方程式演算手段は、線形方程式(eq2)を解くことにより、前記オプティカルフロー［ｕ ｖ］Tを推定すると同時に、その空間勾配［ｕx ｖx ｕy ｖy］Tも推定するものである、請求項１記載の物体運動推定装置。
The captured image is a two-dimensional image, is an image having a size of at least 3 pixels × 3 pixels, and is three images with different imaging times, and is stored in the captured image storage means,
The luminance gradient calculation means is a luminance gradient I x (x 0 , y 0 , t), I y (x 0 , y 0 , t) of the luminance I at the point of interest P (x 0 , y 0 ) of the captured image. , I xx (x 0 , y 0 , t), I xy (x 0 , y 0 , t) and I yy (x 0 , y 0 , t), and the observed quantity shown in the equation (eq1) z (t) is obtained,
The time gradient calculating means, said luminance I and the intensity gradient I x, I y, I xx, I xy and I yy from time gradient I t, I tx, I ty , I txx, I txy, I tyy, I tt , I ttx , I tty , I ttxx , I ttxy, and I ttyy ,
The equation calculation means estimates the optical flow [u v] T by solving the linear equation (eq2), and simultaneously estimates the spatial gradient [u x v x u y v y ] T. The object motion estimation apparatus according to claim 1.
撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定装置において、
推定モデルに対して、状態変数の逐次推定を行う状態推定処理手段
を有する画像解析処理手段を備え、
前記推定モデルは、
前記撮像手段と前記物体との相対運動の速度分布であるオプティカルフローについて、前記撮像画像における前記物体の輝度の時間に関する１階全微分がゼロという時不変性が仮定されることにより、前記オプティカルフローを状態変数に含むだけでなく、
前記撮像手段の撮像軸に対して垂直な平面上を前記物体が移動するという並進性が仮定されないことにより、前記オプティカルフローの空間勾配をも状態変数に含み、
さらに、前記撮像画像上の注目点の空間近傍の輝度を観測量とする観測方程式を含むものであり、
前記状態推定処理手段が、前記状態変数の逐次推定を行うことにより、前記オプティカルフローを推定すると同時に、前記推定されるオプティカルフローの空間勾配をも推定する物体運動推定装置。 In an object motion estimation apparatus for estimating motion of an object by performing image analysis processing on a captured image obtained by imaging the object by an imaging means,
For the estimated model, with an image analysis processing means having a state estimation processing unit that sequentially estimates the state variables,
The estimation model is
The optical flow, which is the velocity distribution of relative motion between the imaging means and the object, is assumed to be invariant when the first-order total derivative with respect to the time of luminance of the object in the captured image is assumed to be zero. As a state variable,
By not assuming the translational property that the object moves on a plane perpendicular to the imaging axis of the imaging means, the spatial gradient of the optical flow is also included in the state variable,
Furthermore, it includes an observation equation in which the luminance in the vicinity of the point of interest on the captured image is an observation amount,
An object motion estimation apparatus in which the state estimation processing means estimates the optical flow by sequentially estimating the state variables and simultaneously estimates a spatial gradient of the estimated optical flow. 前記観測方程式は、システム雑音と観測雑音という雑音を加味したものであり、
前記状態変数の逐次推定は、推定値の共分散を近似するカルマンフィルタを用いて前記状態変数を逐次推定するものである、
請求項３に記載の物体運動推定装置。 The observation equation takes into account system noise and observation noise,
The sequential estimation of the state variable is to sequentially estimate the state variable using a Kalman filter that approximates the covariance of the estimated value.
The object motion estimation apparatus according to claim 3 . 撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定方法において、
輝度勾配演算手段が、前記撮像画像の注目点における輝度勾配を演算するステップと、
時間勾配演算手段が、前記輝度勾配から時間勾配を演算するステップと、
方程式演算手段が、前記撮像手段と前記物体との相対運動の速度分布であるオプティカルフローについて、前記撮像画像における前記物体の輝度の時間に関する１階全微分がゼロという時不変性が仮定されることにより、前記オプティカルフローを変数とし、前記撮像画像の画素の輝度の空間勾配及びその時間微分を係数とするだけでなく、さらに、前記撮像手段の撮像軸に対して垂直な平面上を前記物体が移動するという並進性が仮定されないことにより、前記オプティカルフローの空間勾配をも変数とし、前記輝度及びその空間勾配の時間に関して２階微分した値をも係数とする連立方程式を解くステップ
を含むことにより、前記オプティカルフローを推定すると同時に、前記推定されるオプティカルフローの空間勾配をも推定する物体運動推定方法。 In the object motion estimation method for estimating the motion of the object by performing image analysis processing on a captured image obtained by imaging the object by the imaging unit,
A luminance gradient calculating means calculating a luminance gradient at a point of interest in the captured image;
A time gradient calculating means for calculating a time gradient from the luminance gradient;
For the optical flow that is the velocity distribution of the relative motion between the imaging means and the object, the equation calculation means is assumed to be invariant when the first-order total derivative with respect to the luminance time of the object in the captured image is zero. By using the optical flow as a variable, not only the spatial gradient of the brightness of the pixel of the captured image and its time derivative as a coefficient, but also the object on a plane perpendicular to the imaging axis of the imaging means By not assuming the translational nature of movement, including the step of solving simultaneous equations having the spatial gradient of the optical flow as a variable and the second derivative of the luminance and the time of the spatial gradient as a coefficient An object that estimates the optical flow and at the same time estimates the spatial gradient of the estimated optical flow Motion estimation method. 撮像手段が物体を撮像した撮像画像を画像解析処理して前記物体の運動を推定する物体運動推定方法において、
状態推定処理手段が、前記撮像手段と前記物体との相対運動の速度分布であるオプティカルフローについて、前記撮像画像における前記物体の輝度の時間に関する１階全微分がゼロという時不変性が仮定されることにより、前記オプティカルフローを状態変数に含むだけでなく、前記撮像手段の撮像軸に対して垂直な平面上を前記物体が移動するという並進性が仮定されないことにより、前記オプティカルフローの空間勾配をも状態変数に含み、さらに、前記撮像画像上の注目点の空間近傍の輝度を観測量とする観測方程式を含む推定モデルに対して、前記状態変数の逐次推定を行うステップ
を含むことにより、前記オプティカルフローを推定すると同時に、前記推定されるオプティカルフローの空間勾配をも推定する物体運動推定方法。 In the object motion estimation method for estimating the motion of the object by performing image analysis processing on a captured image obtained by imaging the object by the imaging unit,
For the optical flow that is the velocity distribution of the relative motion between the imaging means and the object, the state estimation processing means is assumed to be invariant when the first-order total derivative with respect to the time of luminance of the object in the captured image is zero. Thus, not only the optical flow is included in the state variable, but also the translational property that the object moves on a plane perpendicular to the imaging axis of the imaging means is not assumed, so that the spatial gradient of the optical flow is reduced. Further including a step of performing sequential estimation of the state variable with respect to an estimation model including an observation equation having an observation amount of luminance near the space of the point of interest on the captured image. An object motion estimation method for estimating an optical flow and simultaneously estimating a spatial gradient of the estimated optical flow. コンピュータを、請求項１から４のいずれかの物体運動推定装置として動作させるためのプログラム。   A program for operating a computer as the object motion estimation apparatus according to claim 1. 請求項７記載のプログラムを記録したコンピュータ読み取り可能な記録媒体。   A computer-readable recording medium on which the program according to claim 7 is recorded.
JP2011501533A 2009-02-24 2010-01-25 Object motion estimation apparatus, object motion estimation method, program, and recording medium Active JP5504426B2 (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2011501533A JP5504426B2 (en) 2009-02-24 2010-01-25 Object motion estimation apparatus, object motion estimation method, program, and recording medium

Applications Claiming Priority (4)

Application Number Priority Date Filing Date Title
JP2009041360 2009-02-24
JP2009041360 2009-02-24
PCT/JP2010/050924 WO2010098160A1 (en) 2009-02-24 2010-01-25 Object motion inferring device, object motion inferring method, program, and recording medium
JP2011501533A JP5504426B2 (en) 2009-02-24 2010-01-25 Object motion estimation apparatus, object motion estimation method, program, and recording medium

Publications (2)

Publication Number Publication Date
JPWO2010098160A1 JPWO2010098160A1 (en) 2012-08-30
JP5504426B2 true JP5504426B2 (en) 2014-05-28

Family

ID=42665370

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2011501533A Active JP5504426B2 (en) 2009-02-24 2010-01-25 Object motion estimation apparatus, object motion estimation method, program, and recording medium

Country Status (3)

Country Link
US (1) US8805023B2 (en)
JP (1) JP5504426B2 (en)
WO (1) WO2010098160A1 (en)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP5986689B2 (en) * 2012-11-06 2016-09-06 アルカテル−ルーセント System and method for processing visual information for event detection
CN110415276B (en) * 2019-07-30 2022-04-05 北京字节跳动网络技术有限公司 Motion information calculation method and device and electronic equipment
CN113870308B (en) * 2021-09-07 2024-03-22 广西科技大学 Weak and small target detection method for spiral gradient optimization estimation

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2000242797A (en) * 1999-02-18 2000-09-08 Toyota Motor Corp Method for detecting movement of image and object detecting device
JP2008203538A (en) * 2007-02-20 2008-09-04 National Univ Corp Shizuoka Univ Image display system

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5991459A (en) * 1992-01-22 1999-11-23 Eastman Kodak Company Method of modifying a time-varying image sequence by estimation of velocity vector fields
US6671324B2 (en) * 2001-04-16 2003-12-30 Mitsubishi Electric Research Laboratories, Inc. Estimating total average distortion in a video with variable frameskip
US7558402B2 (en) * 2003-03-07 2009-07-07 Siemens Medical Solutions Usa, Inc. System and method for tracking a global shape of an object in motion
US7689035B2 (en) * 2005-06-17 2010-03-30 The Regents Of The University Of California Methods for identifying, separating and editing reflection components in multi-channel images and videos
US7711147B2 (en) * 2006-07-28 2010-05-04 Honda Motor Co., Ltd. Time-to-contact estimation device and method for estimating time to contact
US8264614B2 (en) * 2008-01-17 2012-09-11 Sharp Laboratories Of America, Inc. Systems and methods for video processing based on motion-aligned spatio-temporal steering kernel regression
JP5547739B2 (en) * 2008-10-15 2014-07-16 イノベイティブ テクノロジー ディストリビューターズ エルエルシー Digital processing method and system for optical flow determination

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2000242797A (en) * 1999-02-18 2000-09-08 Toyota Motor Corp Method for detecting movement of image and object detecting device
JP2008203538A (en) * 2007-02-20 2008-09-04 National Univ Corp Shizuoka Univ Image display system

Also Published As

Publication number Publication date
WO2010098160A1 (en) 2010-09-02
US20110299739A1 (en) 2011-12-08
US8805023B2 (en) 2014-08-12
JPWO2010098160A1 (en) 2012-08-30

Similar Documents

Publication Publication Date Title
Jana et al. Computer vision‐based real‐time cable tension estimation in Dubrovnik cable‐stayed bridge using moving handheld video camera
JP4814669B2 (en) 3D coordinate acquisition device
Kearney et al. Optical flow estimation: An error analysis of gradient-based methods with local optimization
US10217271B2 (en) Method and apparatus for compensating for movement of head-mounted display
EP2128821B1 (en) Moving object recognizing apparatus
US20100073366A1 (en) Model generation apparatus and method
US20180075614A1 (en) Method of Depth Estimation Using a Camera and Inertial Sensor
Besnard et al. A space–time approach in digital image correlation: Movie-DIC
JP5504426B2 (en) Object motion estimation apparatus, object motion estimation method, program, and recording medium
KR20160020172A (en) Method and Apparatus for detecting lane of road
Shi et al. Deep learning based virtual point tracking for real-time target-less dynamic displacement measurement in railway applications
CN103875020A (en) Resolving homography decomposition ambiguity based on orientation sensors
CN109074480A (en) Method, computing device, driver assistance system and the motor vehicles of rolling shutter effect in the image of environmental area for detecting motor vehicles
Hu et al. Real-time video stabilization for fast-moving vehicle cameras
Wang et al. Completely non-contact modal testing of full-scale bridge in challenging conditions using vision sensing systems
JP4734568B2 (en) Method and apparatus for determining moving object measurement point on image
Singh et al. Methods to estimate subpixel level small motion from video of vibrating cutting tools
Shiba et al. A fast geometric regularizer to mitigate event collapse in the contrast maximization framework
JP2011064639A (en) Distance measuring device and distance measuring method
KR101791166B1 (en) Apparatus and Method for Estimation of Spatial information of an object
US20060020562A1 (en) Apparatus and method for estimating optical flow
WO2019097576A1 (en) Measurement system, correction processing device, correction processing method and computer-readable recording medium
CN113447111B (en) Visual vibration amplification method, detection method and system based on morphological component analysis
JP3985610B2 (en) Vehicle traveling path recognition device
Raizada et al. Comparative analysis of image processing schemes to register motion from video of vibrating cutting tools

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20130109

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20131010

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20140217

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20140219

R150 Certificate of patent or registration of utility model

Ref document number: 5504426

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150

R250 Receipt of annual fees

Free format text: JAPANESE INTERMEDIATE CODE: R250

R250 Receipt of annual fees

Free format text: JAPANESE INTERMEDIATE CODE: R250

R250 Receipt of annual fees

Free format text: JAPANESE INTERMEDIATE CODE: R250