CN104915968B

CN104915968B - A kind of optical flow velocity method of estimation based on navier stokes equations

Info

Publication number: CN104915968B
Application number: CN201510254923.4A
Authority: CN
Inventors: 齐越; 左庆
Original assignee: Beihang University
Current assignee: Beihang University
Priority date: 2015-05-19
Filing date: 2015-05-19
Publication date: 2017-07-21
Anticipated expiration: 2035-05-19
Also published as: CN104915968A

Abstract

The invention discloses a kind of optical flow velocity method of estimation based on navier stokes equations.This method constructs new energy-optimised function, and navier stokes equations are deformed and are incorporated among energy equation, and whole energy equation includes the constraint of four aspects：Image consistency, velocity field is smooth, velocity field couple in time and velocity field Incoercibility.New energy-optimised equation proposed by the invention can carry out optical flow velocity estimation simultaneously to whole sequence, but need the variable solved very many, in order to stably solve energy-optimised function, we used a variety of optimization means, optimize including multiresolution from coarse to fine, alternative optimization, multiplier alternating direction optimizes, a variety of optimization methods such as point of safes optimization.

Description

A kind of optical flow velocity method of estimation based on Navier Stokes equation

Technical field

The invention belongs to computer vision and field of Computer Graphics, specifically a kind of estimation image or three-dimensional The method of flow field velocity in volume data, this method can be used for video tracking, in terms of air motion estimation and fluid analysis.

Background technology

With developing rapidly for computer technology, it is true in computer animation that people increasingly focus on natural phenomena Property, in field of Computer Graphics, such as cigarette, fire, the simulation of water or other fluid are essential.One class method is completely by physics Equation is simulated, but modulation parameter is very cumbersome, and another kind of is that the stream in true environment is gathered using equipment Volume data, is then analyzed and is reused to it.The present invention is a link of Equations of The Second Kind method, can be for collecting Image or volume data carry out velocity field estimation.

The motion of fluid can be retouched with famous Navier Stokes equation (Navier-Stokes Equation) State：

Wherein u represents velocity field, and f represents external force, and p represents pressure, and t represents the time.And the visible oBject in fluid is (as contaminated Material, or cigarette particle etc., use I is represented) can be described with equation below：

The optical flow analysis method of early stage only takes into account the uniformity matching of formula (3), i.e. image, this method (and after Come improved method on this basis) velocity field can be substantially obtained, but this velocity field does not utilize Na Wei-Stokes side Any constraint of journey, it is accurate not enough for fluid.Formula (2) is incorporated into optical flow velocity estimation by later Gregson Among model so that the flow field of acquisition meets incompressible characteristic.This method using two consecutive frames for being estimated For velocity field, good result can be obtained, however, be that there is very strong coupling for the motion in time of fluid, This coupling is described by formula (1).

The content of the invention

To overcome disadvantages mentioned above, it is an object of the invention to provide a kind of optical flow velocity field method of estimation of global optimum, This method to whole sequence can meet the characteristics of motion of fluid between optical flow velocity estimation, adjacent two frame simultaneously, and whole Sequence also meets the characteristics of motion of fluid.

In order to achieve the above object, the present invention proposes a kind of optical flow velocity estimation based on Navier Stokes equation Method.A brand-new energy optimizing model based on Na Wei-Stokes is first proposed, has then used and has made meticulously very much The optimisation strategy of work carries out velocity field solution.Implementation method is：First, the golden word of multiresolution is constructed to each frame of fluid sequence Tower, the pyramidal number of plies is specified by user, but the speed that standard is lowest resolution image does not exceed a length in pixels； Then, each level of optimization from coarse to fine；In a level every time, optimize each frame according to vertical sequence alternate Velocity field, the iterations of alternative optimization is specified by user；In each frame speed of alternative optimization, by this frame speed place Corresponding energy model splits into three parts：Image consistency is smooth with velocity field, the time coupling of velocity field, velocity field Incoercibility, with multiplier alternating direction optimization method to this three parts alternative optimization, iterations is specified by user.

Brief description of the drawings

Fig. 1 shows the broad flow diagram of the optical flow velocity method of estimation of the invention based on Navier Stokes equation；

Fig. 2 shows the algorithm of multiplier alternating direction optimization method (ADMM) of the present invention.

Embodiment

As shown in figure 1, the optical flow velocity method of estimation based on Navier Stokes equation of the present invention is using following step Suddenly：

(1) multi-resolution pyramid is constructed to each frame of image sequence.The accuracy of arithmetic speed and result is taken into account, I 0.5 zoom factor has been used to construct pyramid, and the pyramidal number of plies is specified by user according to actual conditions, is specified Standard be to allow length of the optical flow velocity field in the image of lowest resolution not exceed a pixel.

(2) it is successively secondary from lowest resolution to highest resolution to optimize.Optimized first lowest resolution, in computing most Start, the initial value of all velocity fields is set to 0, and by optimizing equation below come for each frame speed field computation initial value：

E(u_i)=E_data(u_i)+α²E_sm(u_i)

E_dataFor image consistency energy, specific formula is：

Locus is x, and integral domain is Ω, I_iFor the i-th two field picture.

E_smFor velocity field smoothed energy, specific formula isα is weight coefficient, wherein, N For sequence length, i.e. number of image frames, u_iFor the i-th frame speed, alpha, gamma is weight coefficient.

Although this initial value is also differed farther out from optimal solution, it supports follow-up optimization process enough.

(3) each frame speed is alternately optimized.When optimizing the i-th frame speed, the velocity field of other frames is accordingly to be regarded as Known quantity, is consequently formed the energy theorem specifically designed for the i-th frame speed：

(4) formula (4) is deformed, is divided into three parts and is optimized with multiplier alternating direction optimization method.First It is that, into unconfined formula, divergence constraint switchs to energy term by the constraints conversion of formula (4) Then three parts are respectively：

F₁(u_i)=E_data(u_i)+α²E_sm(u_i),

F₃(u_i)=η E_div(u_i).

Wherein η is coefficient.Multiplier alternating direction optimization method can operate (proximal operator) table with near-end Show, the algorithm seen in Fig. 2.

Near-end handling function is defined as

(a) optimization method of Part I isHere b X3-y1 in corresponding diagram 2 in algorithm, correspondence after formula expansion；

Euler-Lagrange equation is first obtained with the calculus of variations, is then deployed using point of safes alternative manner and first order Taylor Solved.

(b) optimization method of Part II isHere X3-y2 in c respective figures 2 in algorithm, correspondence after formula expansion；

Using the method for point of safes iteration, its Euler-Lagrange equation is during+1 iteration of kth：

During i=1

During i=N-1

1<i<During N-1

Computation sequence is Finally optimal solution is：

During i=1,

During i=N-1,

1<i<During N-1,

(c) optimization method of Part III isHere In d corresponding diagrams 2 in algorithm(it is now formula (4) optimization when coefficient η tends to infinity Target), its optimal solution is equivalent to u_i=P (d).

Claims

1. a kind of optical flow velocity method of estimation based on Navier Stokes equation, it is characterised in that comprise the following steps：

(1) multi-resolution pyramid is constructed for each two field picture, optimizes the optical flow velocity field of each level from coarse to finely；

(2) each frame speed of whole sequence is alternately optimized, when optimizing a certain frame speed, by the speed of other frames Regard constant as in field；

(3) velocity field currently to be optimized, its corresponding energy theorem can be turned into three parts, be optimized using multiplier alternating direction Method (ADMM) is optimized；Part I energy is image consistency energy and velocity field smoothed energy, Part II energy For the coupling energy of velocity field in time, Part III energy is the Incoercibility energy of velocity field；

(4) using point of safes alternative manner optimization Part I energy；

(5) using point of safes alternative manner optimization Part II energy；

(6) carry out velocity projections to optimize Part III energy using Poisson's equation；

Using following light stream energy optimization model, it is specially：

\begin{matrix} E (u_{1}, u_{2}, ..., u_{N - 1}) = Σ_{i = 1}^{N - 1} E_{d a t a} (u_{i}) + α^{2} (Σ_{i}^{N - 1} E_{s m} (u_{i}) + γ Σ_{i = 1}^{N - 2} E_{t e m} (u_{i}, u_{i + 1})) \\ \begin{matrix} s . t . & &ForAll; i &Element; [0, N - 1], & &dtri; \cdot u_{i} = 0 \end{matrix} \end{matrix} - - - (1)

Wherein, N is sequence length, i.e. number of image frames, u_iFor the i-th frame speed, alpha, gamma is weight coefficient；

E_dataFor image consistency energy, specific formula is：

E_{d a t a} (u_{i}) = {&Integral;}_{Ω} | | I_{i} (x) - I_{i + 1} (x + u_{i}) | |_{2}^{2} d Ω,

I_iFor the i-th two field picture, locus is x, and integral domain is Ω,

E_smFor velocity field smoothed energy, specific formula is

E_temEnergy is coupled for the time, specific formula is,

E_{t e m} (u_{i}, u_{i + 1}) = {&Integral;}_{Ω} | | P (u_{i + 1} (x + u_{i}) - f_{i} (x)) - u_{i} (x) | |_{2}^{2} d Ω - - - (2)

E_{t e m} (u_{i - 1}, u_{i}) = {&Integral;}_{Ω} | | P (u_{i - 1} (x - u_{i}) + f_{i} (x)) - u_{i} (x) | |_{2}^{2} d Ω - - - (3)

Speed is thrown into incompressible subspace by P () expressions, the mode i.e. u projected using Poisson_out=P (u_in) be calculated asFormer frame, a later frame velocity field are represented respectively Velocity field after advection, f_iRepresent the external force of the i-th frame.