CN103700117A - Robust optical flow field estimating method based on TV-L1 variation model - Google Patents

Abstract

The invention discloses a robust optical flow field estimating method based on a TV-L1 variation model. The robust optical flow field estimating method comprises the following steps: firstly, performing structural texture resolution on an input image, and establishing an optical flow calculating model based on the TV-L1; secondly, establishing an image pyramid, calculating optical flow on the lowest image resolution layer by a discretized alternating iteration method, further calculating with a calculated value as an initial value of a next higher resolution layer till the highest resolution layer, namely the original image resolution, and accelerating the algorithm by using a GPU (graphic processing unit) so as to improve the instantaneity of the algorithm; finally, calculating the error of the algorithm by using an optical flow error evaluating function. In the robust optical flow field estimating method, the input image is processed by a structural texture resolving method and a texture image is applied to optical flow calculation, so that influence of an image shadow caused by illumination variation on the calculation is avoided; by the robust optical flow field estimating method based on the TV-L1 variation model, segmenting smoothness of the image is kept and the optical flow calculating precision and optical flow calculating speed are improved.

Description

A kind of based on TV-L 1the robust Optical flow estimation method of Variation Model

Technical field

The invention belongs to computer vision field, relate to a kind of based on TV-L ¹the robust Optical flow estimation method of (single order data item) Variation Model.

Background technology

Light stream refers to the instantaneous velocity of space motion object pixel motion on observation imaging surface, and the information that it is comprising moving object, so can utilize it to understand the motion conditions of object.Light stream comprises following three key elements: the one, and the motion that light stream is produced, namely velocity field; The 2nd, can carry information and there is the supporting body of optical characteristics, such as the picture element with gray scale; The 3rd, object is projected to the plane of delineation from scene, make object can be observed projection.Optical flow computation is one of important research field of computer vision and image processing, in fields such as space flight, military affairs, medical science, industry, all has a wide range of applications.

Calculating optical flow field is at present mainly the variational method, its starting point is to change obtaining of optical flow field into a global energy functional minimization problem, and the foundation of energy functional model is the key of this variational method, its mathematical model is mainly comprised of data item and level and smooth two parts, data item mainly comprises various normal logically false establishing on duty, such as gray scale conservation hypothesis (light stream fundamental equation), gradient conservation hypothesis etc., the constraint condition consisting of these conservation hypothesis is in optical flow computation process, to determine the principal element of motion model.And level and smooth Xiang Ze is the various level and smooth and sectionally smooth strategy that has reflected optical flow field, and can make this model obtain unique solution.

The people such as Horn and Schunk applies to the variational method in optical flow computation at first, creationary grey scale pixel value and two-dimension speed field is linked together, and introduces the smoothness constraint equation of light stream, obtains classical optical flow computation method.The method is mainly utilized gray scale conservation hypothesis, and very short in the time interval of adjacent two two field pictures, variation of image grayscale is very little, supposes that gray-scale value is constant, be I(x+dx, y+dy, t+dt)=I(x, y, t), this formula is carried out to single order Taylor expansion, can derive gradient constraint equation, that is:

$\frac{\∂I}{\∂x}\mathrm{dx}+\frac{\∂I}{\∂y}\mathrm{dy}+\frac{\∂I}{\∂t}\mathrm{dt}=0$

Order $u=\frac{\mathrm{dx}}{\mathrm{dt}},v=\frac{\mathrm{dy}}{\mathrm{dt}}$ Represent that respectively this point is at x, the light stream component in y direction, $I_x=\frac{\∂I}{\∂x},I_y=\frac{\∂I}{\∂y},I_t=\frac{\∂I}{\∂t}$ Distinguish representative image gray scale with respect to x, y, the local derviation of t, above formula can be written as following basic optical flow constraint equation:

I _xu+I _yv+I _t＝0

Optical flow constraint equation contains two unknown quantity u and v, therefore cannot determine unique solution.For this reason, Horn & Schunck introduces an overall smoothness constraint to light stream to above formula and solves light stream, thereby obtains following energy functional:

$E=\∫\∫[{(I_{x}u+I_{y}v+I_t)}^2+{\mathrm{\α}}^2({|\▿u|}^2+{|\▿u|}^2)]\mathrm{dxdy}$

Above formula is basic optical flow computation method.But the method can not keep sectionally smooth well, the spill point problem of deal with data item steadily, and by this energy functional of iterative, want to obtain more stable result, must, through thousands of up to a hundred iteration, reduce computing velocity.

Application number is that 201310174158.6 patent has proposed a kind of quick the Computation of Optical Flow based on error Distributed-tier grid, solving on the real time problems of optical flow computation, utilize multiple grid method to solve energy model, but the basic optical flow computation model of the method utilization, just on derivation algorithm, improved computing velocity, for model, do not improved.Want fundamentally to put forward the high-precision optical flow computation of simultaneously accelerating, must on original model, modify, to obtain the better computation model of effect.

Summary of the invention

For the precision and the real time problems that exist in basic optical flow computation method, the present invention proposes a kind of based on TV-L ¹the robust Optical flow estimation method of Variation Model, introduces single order data item (being L1 norm), can keep well sectionally smooth, and can accelerate the computing velocity of light stream.

Provide the principle of the Computation of Optical Flow of the present invention below.

It is as follows that variational method based on gray scale conservation solves the energy functional of light stream:

${\∫}_{\mathrm{\Ω}}\{\mathrm{\λ\φ}(I_0\left(x\right)-I_1(x+u\left(x\right)))+\mathrm{\ψ}(u,\▿u,...)\}\mathrm{dx}$

Wherein, u (x) is the light stream (x direction and the light stream of y direction) on two-dimensional directional, φ (I ₀(x)-I ₁(x+u (x))) be the data penalty term of image,

for level and smooth.Parameter lambda is the weight coefficient between data item and level and smooth item.If select φ (x)=x ²with

this energy functional is the basic optical flow computation method of Horn & Schunck.

Select

this energy functional is based on TV-L ¹the light stream energy functional of (the total variation method of single order data item), as follows:

$E={\∫}_{\mathrm{\Ω}}\left\{\mathrm{\λ}\right|I_0\left(x\right)-I_1(x+u\left(x\right))|+|\▿u\left|\right\}\mathrm{dx}$

Above-mentioned energy functional formula seems simply, but solves but very difficultly, and main cause is that data item and level and smooth are not continuously differentiable.In order to address this problem, can be by

with functional expression that below can be micro-, replace, with

(ε is a very little constant, and while preventing differentiate, denominator is the zero calculation overflow that causes) just can utilize this energy functional of Numerical Methods Solve like this.But the method has been introduced the margin of error, therefore can affect the precision of optical flow computation.Adopting primal dual algorithm to replace this energy functional of iterative can effectively avoid introducing the margin of error and cause solving out of true.

A kind of based on TV-L ¹the robust Optical flow estimation method of Variation Model, its technical characterictic mainly comprises following steps:

Step 1, input image sequence.

Step 2, carries out structural texture decomposition to image.

In actual applications, the factor that affects optical flow computation precision comprise moving object large change in displacement, weak texture region, block and intensity of illumination variation etc.Wherein, the impact that intensity of illumination changes is particularly evident, this important hypothesis of gray scale conservation of widespread use is no longer set up, thereby be difficult to obtain optical flow field accurately.

Picture breakdown technology is absorbed in and from image, is extracted useful, interested information.The method of using structural texture to decompose solves intensity of illumination and changes, and the problem such as consequent shade.Its theoretical foundation is, image can be regarded as by structure division (mainly comprising geological information in image, as striped, edge etc.) and texture part (mainly comprising the small scale detailed information that some have cyclophysis or oscillating characteristic) and forms.Image is after structure-texture decomposes, and intensity of illumination changes the negative effect producing and only appears in structural drawing.Input quantity using texture image as optical flow computation process, can avoid the impact of illumination variation on result of calculation.

ROF (Rudin, Osher, the Fatemi) denoising model of utilization based on total variation carried out structural texture decomposition.For gray level image I (x), its structure division I _ssolving model be:

$\underset{I_S}{\min }{\∫}_{\mathrm{\Ω}}\left\{\right|\▿I_S|+\frac{1}{2\mathrm{\θ}}{(I_S-I)}^2\}\mathrm{dx}$

In formula, θ is a very little constant, only has the I of working as in optimizing process _swhile approaching with I, just can make energy functional obtain minimum value, I is original image gray-scale value.

Utilize primal dual algorithm to minimize this energy functional.Introduce I _sdual variable p _i(i=1,2), adopt dual variable p=(p ₁, p ₂) ^titerative solution equation:

I _S＝I+θdiv?p

Iterative formula is:

${\tilde{p}}^{n+1}=p^n+\frac{\mathrm{\τ}}{\mathrm{\θ}}(\▿(I+\mathrm{\θ}\mathrm{div}{p}^n))$

$p^{n+1}=\frac{{\tilde{p}}^{n+1}}{\max \{1,|{\tilde{p}}^{n+1}\left|\right\}}$

Wherein, p ⁰=0, τ≤1/4.

The texture component I of image _t(x) equal the poor of original graph and structure component, that is: I _t(x)=I (x)-I _s(x).

By the texture maps obtaining after decomposing in order to follow-up optical flow computation.Practice shows, the method can reduce the impact of illumination variation on optical flow computation, improves the precision solving, and operation is also more efficient.

Step 3, sets up based on TV-L ¹the energy functional model of Variation Model.

Based on TV-L ¹the energy functional of model is as follows:

$E={\∫}_{\mathrm{\Ω}}\left\{\mathrm{\λ}\right|I_1(x+u\left(x\right))-I_0\left(x\right)|+|\▿u\left|\right\}\mathrm{dx}$

By image I ₁at x+u ₀near carry out linearization, to I ₁(x+u (x)) single order Taylor expansion:

$I_1(x+u\left(x\right))=I_1(x+u_0)+(u-u_0)\▿I_1(x+u_0)$

Fixing u ₀and utilize linear-apporximation to replace I ₁(x+u (x)), TV-L ¹energy functional is write as following form:

$E={\&Integral;}_{\mathrm{\&Omega;}}\left\{\mathrm{\&lambda;}\right|u\&dtri;{\mathrm{<figure-callout\; id="1"\; label="I"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">I</figure-callout>}}_1+I_1(x+u_0)-u_0\&dtri;I_1-I_0|+|\&dtri;u\left|\right\}\mathrm{dx}$

With ρ (u), represent

introduce auxiliary variable v, by TV-L ¹energy functional is write as following convex function form:

$E_{\mathrm{\&theta;}}={\&Integral;}_{\mathrm{\&Omega;}}\left\{\right|\&dtri;u|+\frac{1}{2\mathrm{\&theta;}}{(u-v)}^2+\mathrm{\&lambda;}|\mathrm{\&rho;}\left(v\right)\left|\right\}\mathrm{dx}$

Wherein, θ is a very little constant, in the process of iteration, only has when u and v approach, can make above formula energy functional obtain minimum value.

This energy functional is write as to the mathematical model of concrete two dimensional form:

$E_{\mathrm{\&theta;}}={\&Integral;}_{\mathrm{\&Omega;}}\left\{\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\right|\&dtri;u_i|+\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\frac{1}{2\mathrm{\&theta;}}{(u_i-v_i)}^2+\mathrm{\&lambda;}|\mathrm{\&rho;}\left(v\right)\left|\right\}\mathrm{dx}$

Wherein, u ₁and u ₂represent respectively the light stream of x direction and the light stream of y direction.Utilize the alternately alternative manner of primal dual algorithm to optimize above-mentioned model, can try to achieve light stream.

Step 4, utilizes alternately alternative manner to solve energy model.

(1) for i(i=1,2), fixing v ₁and v ₂, solve u ₁and u ₂, optimize with drag:

$\underset{u_i}{\min }{\&Integral;}_{\mathrm{\&Omega;}}\left\{\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\right|\&dtri;u_i|+\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\frac{1}{2\mathrm{\&theta;}}{(u_i-v_i)}^2\}\mathrm{dx}$

This model is the image denoising model based on ROF, solves the following formula of this model utilization:

u _i＝v _i+θdiv?p _i

U _ithe dual variable of (i=1,2) is p _i(i=1,2), the iterative formula that solves dual variable is as follows:

${\tilde{p}}^{n+1}=p^n+\frac{\mathrm{\&tau;}}{\mathrm{\&theta;}}(\&dtri;(I+\mathrm{\&theta;}\mathrm{div}{p}^n))$

$p^{n+1}=\frac{{\tilde{p}}^{n+1}}{\max \{1,|{\tilde{p}}^{n+1}\left|\right\}}$

Wherein, p ⁰=0, τ≤1/4.

(2) for i(i=1,2), fixing u ₁and u ₂, solve v ₁and v ₂, optimize with drag:

$\underset{v}{\min }{\&Integral;}_{\mathrm{\&Omega;}}\{\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\frac{1}{2\mathrm{\&theta;}}{(u_i-v_i)}^2+\mathrm{\&lambda;}|\mathrm{\&rho;}\left(v\right)\left|\right\}\mathrm{dx}$

Solve the method that this model need to use pointwise to solve, its solution procedure is used following threshold method:

$v=u+\left\{\begin{array}{ccc}\mathrm{\&lambda;\&theta;}\&dtri;{\mathrm{<figure-callout\; id="1"\; label="I"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">I</figure-callout>}}_1& \mathrm{if}& \mathrm{\&rho;}\left(u\right)<-\mathrm{\&lambda;\&theta;}{|\&dtri;{\mathrm{<figure-callout\; id="1"\; label="I"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">I</figure-callout>}}_1|}^2\\ -\mathrm{\&lambda;\&theta;}\&dtri;{\mathrm{<figure-callout\; id="1"\; label="I"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">I</figure-callout>}}_1& \mathrm{if}& \mathrm{\&rho;}\left(u\right)>\mathrm{\&lambda;\&theta;}{|\&dtri;{\mathrm{<figure-callout\; id="1"\; label="I"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">I</figure-callout>}}_1|}^2\\ -\mathrm{\&rho;}\left(u\right)\&dtri;{\mathrm{<figure-callout\; id="1"\; label="I"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">I</figure-callout>}}_1/{|\&dtri;{\mathrm{<figure-callout\; id="1"\; label="I"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">I</figure-callout>}}_1|}^2& \mathrm{if}& \left|\mathrm{\&rho;}\left(u\right)\right|\&le;\mathrm{\&lambda;\&theta;}{|\&dtri;I_1|}^2\end{array}\right.$

In the specific implementation process of program solution, due to pixel light stream value gradient

very difficult with solving of dual variable div p, take the method for value solving of discretize to solve energy model, method is as follows:

Utilize consequent difference method to solve div p, utilize preceding paragraph difference method to calculate the gradient of light stream.

For one wide be that in N, the high view data scope for M, pixel (i, j) is located gradient preceding paragraph difference discrete be expressed as:

${(\&dtri;u)}_{i,\mathrm{<figure-callout\; id="1"\; label="j"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">j</figure-callout>}}^1=\left\{\begin{array}{ccc}{u}_{i+1,j}-u_{i,j},& \mathrm{if}& i<N\\ 0,& \mathrm{if}& i=N\end{array}\right.$

${(\&dtri;u)}_{i,\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">j</figure-callout>}}^2=\left\{\begin{array}{ccc}{u}_{i,j+1}-u_{i,j},& \mathrm{if}& i<M\\ 0,& \mathrm{if}& i=M\end{array}\right.$

The consequent difference discrete of divergence operator is expressed as:

${\left(\mathrm{div}p\right)}_{i,j}=\left\{\begin{array}{ccc}{p}_{i,j}^1-p_{i-1,\mathrm{<figure-callout\; id="1"\; label="j"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">j</figure-callout>}}^1,& \mathrm{if}& 1<i<N\\ p_{i,\mathrm{<figure-callout\; id="1"\; label="j"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">j</figure-callout>}}^1,& \mathrm{if}& i=1\\ -p_{i-1,\mathrm{<figure-callout\; id="1"\; label="j"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">j</figure-callout>}}^1,& \mathrm{if}& i=N\end{array}\right.+\left\{\begin{array}{ccc}{p}_{i,j}^2-p_{i,j-1}^2,& \mathrm{if}& 1<j<M\\ p_{i,\mathrm{<figure-callout\; id="2"\; label="j"\; filenames="HDA0000419724490000021.png"\; state="\{\{state\}\}">j</figure-callout>}}^2,& \mathrm{if}& j=1\\ -p_{i,j-1}^2,& \mathrm{if}& j=M\end{array}\right.$

Numerical evaluation is carried out (for example image-region of rectangle) in regular grid, utilizes computing power and parallel processing capability that GPU is powerful, and the solution procedure of energy functional is accelerated, and makes optical flow computation have good real-time.

In linearization procedure, because the single order Taylor expansion of data item is only applicable to the projection coordinate of little displacement, the solution procedure of light stream is restricted.In order to overcome this difficulty, the solution procedure of energy model has been taked by slightly to smart pyramid algorith: generate the pyramid diagram picture that a series of thicknesses are different, on the thickest tomographic image yardstick, utilize TV-L ¹variation Model is asked for light stream, and the initial solution using this solution on the thinner image of lower one deck, repeats this step, until the thinnest resolution tomographic image, i.e. original image.

Compared with prior art, the present invention has the following advantages:

Utilize structural texture decomposition method to process input picture, texture maps is applied in optical flow computation, avoided the shade that in image, illumination variation produces on calculating the impact causing.Utilization is based on TV-L ¹the Optical flow estimation method of Variation Model, can well keep the sectionally smooth of image, makes result of calculation more accurate, and can improve the computing velocity of light stream.

Accompanying drawing explanation

Fig. 1 is based on TV-L ¹the robust optical flow computation method flow diagram of Variation Model;

Fig. 2 solves TV-L for replacing alternative manner ¹energy functional process flow diagram;

Fig. 3 is two groups of image sequences and light stream color diagram and polar plot in Middlebury java standard library.

Embodiment

Below in conjunction with the drawings and specific embodiments, the present invention will be further described.

Hardware configuration of the present invention is for being equipped with the PC of i3-3220CPU and GT630GPU, and running environment is Windows7 operating system and Visual Studio2008 software platform.

The technical scheme that the present invention takes is: first the image of input is carried out to structural texture decomposition, the texture part of asking for applies in light stream estimation below, sets up based on TV-L ¹optical flow computation model, then set up 4～5 tomographic image pyramids, on minimum image resolution ratio layer, by the method for the alternately iteration after discretize, calculate light stream, respectively in the hope of the worthwhile initial value of doing next floor height layers of resolution continue to calculate, until highest resolution layer (being original image resolution), utilize GPU to accelerate this algorithm to improve algorithm real-time, finally utilize light stream error assessment function to calculate the error of this algorithm.Fig. 1 is method flow of the present invention, specifically comprises following step:

Step 1, input two continuous frames view data.

Step 2, carries out structural texture decomposition to image.

Step 3, sets up based on TV-L ¹the energy functional model of Variation Model.

Step 4, takes alternately alternative manner to solve energy model, and process flow diagram as shown in Figure 2.

Provide an application example of the present invention below.

Select that two groups of test patterns propose the present invention based on TV-L ¹the robust Optical flow estimation method of Variation Model is verified, as shown in Figure 3, two groups of images are all from the test pattern sequence (a1 in the Middlebury java standard library extensively adopting in the world, b1 and a2, b2 is respectively the two continuous frames of RubberWhale and two groups of image sequences of Hydrangea), the light stream color diagram drawing is respectively as shown in c1, c2, and light stream vector figure is respectively as shown in d1, d2.

In order to compare with prior art, adopt respectively traditional optical flow computation model improve one's methods (method one) and the present invention proposes based on TV-L ¹the robust optical flow computation method (method two) of Variation Model is tested, and experiment still adopts a1, b1 in Fig. 3 and a2, two groups of image sequences of b2.

Adopt the estimation of error mode of average angle error AAE (Average Angular Error) and average end point error AEPE (Average Endpoint Error) to evaluate the precision that two kinds of distinct methods calculate.

The computing formula of AAE is as follows:

$\mathrm{AAE}=\frac{1}{N_0}\underset{i=1}{\overset{N_0}{\mathrm{\&Sigma;}}}{\mathrm{\&phi;}}_e\left(i\right)$

${\mathrm{\&phi;}}_e\left(i\right)=\arccos \left[\frac{u_{1i}^c{u}_{1i}^e+u_{2i}^c{u}_{2i}^e+k^2}{\sqrt{{\left(u_{1i}^c\right)}^2+{\left(u_{2i}^c\right)}^2+k^2}\sqrt{{\left(u_{1i}^e\right)}^2+{\left(u_{2i}^e\right)}^2+k^2}}\right]$

Wherein, N ₀the total pixel number that represents a two field picture,

the standard light flow vector that represents i pixel,

the light stream vector of i the pixel that expression calculates, k represents the frame number of being separated by.AAE has reflected the degree of the whole departure standard light flow vector field, light stream vector field of calculating.

The computing formula of AEPE is:

$\mathrm{AEPE}=\frac{1}{N_0}\underset{i=1}{\overset{N_0}{\mathrm{\&Sigma;}}}{\mathrm{\&phi;}}_{\mathrm{ep}}\left(i\right)$

Wherein:

${\mathrm{\&phi;}}_{\mathrm{ep}}\left(i\right)=\sqrt{{(u_{1i}^c-u_{1i}^e)}^2+{(u_{2i}^c-u_{2i}^e)}^2}$

AEPE is in order to weigh the vector length of optical flow field of calculating and the error between the vector length in standard light flow field.

Table 1 has provided AAE, the AEPE of two kinds of methods and has calculated the time used.As shown in Table 1, the present invention is based on TV-L ¹aAE, the AEPE of the robust Optical flow estimation method of Variation Model and computing time be improving one's methods lower than traditional optical flow computation model all, compared with prior art, can not only improve the precision of optical flow algorithm, can also improve the speed of optical flow computation, there is good real-time.

The contrast of table 1 the present invention and the prior art error of calculation and speed

Claims (4)

1. one kind based on TV-L ¹the robust Optical flow estimation method of Variation Model, is characterized in that comprising following steps:

Step 1, input image sequence;

Step 2, and image is carried out to structural texture decomposition;

ROF (Rudin, Osher, the Fatemi) denoising model of utilization based on total variation carried out structural texture decomposition; For gray level image I (x), its structure division I _ssolving model be:

$\underset{I_S}{\min }{\&Integral;}_{\mathrm{\&Omega;}}\left\{\right|\&dtri;I_S|+\frac{1}{2\mathrm{\&theta;}}{(I_S-I)}^2\}\mathrm{dx}$

In formula, θ is a very little constant, only has the I of working as in optimizing process _swhile approaching with I, just can make energy functional obtain minimum value, I is original image gray-scale value;

Utilize primal dual algorithm to minimize this energy functional; Introduce I _sdual variable p _i(i=1,2), adopt dual variable p=(p ₁, p ₂) ^titerative solution equation:

I _S＝I+θdiv?p

Iterative formula is:

${\tilde{p}}^{n+1}=p^n+\frac{\mathrm{\&tau;}}{\mathrm{\&theta;}}(\&dtri;(I+\mathrm{\&theta;}\mathrm{div}{p}^n))$

$p^{n+1}=\frac{{\tilde{p}}^{n+1}}{\max \{1,|{\tilde{p}}^{n+1}\left|\right\}}$

Wherein, p ⁰=0, τ≤1/4;

The texture component I of image _t(x) equal the poor of original graph and structure component, that is: I _t(x)=I (x)-I _s(x);

Step 3, sets up based on TV-L ¹the energy functional model of Variation Model;

Based on TV-L ¹the energy functional of model is as follows:

E＝∫ _Ω{λ|I ₁(x+u(x))-I ₀(x)|+|▽u|}dx

By image I ₁at x+u ₀near carry out linearization, to I ₁(x+u (x)) single order Taylor expansion:

I ₁(x+u(x))＝I ₁(x+u ₀)+(u-u ₀)▽I ₁(x+u ₀)

Fixing u ₀and utilize linear-apporximation to replace I ₁(x+u (x)), TV-L ¹energy functional is write as following form:

E＝∫ _Ω{λ|u▽I ₁+I ₁(x+u ₀)-u ₀▽I ₁-I ₀|+|▽u|}dx

With ρ (u), represent I ₁(x+u ₀)+(u-u ₀) ▽ I ₁-I ₀, introduce auxiliary variable v, by TV-L ¹energy functional is write as following convex function form:

$E={\&Integral;}_{\mathrm{\&Omega;}}\left\{\mathrm{\&lambda;}\right|I_1(x+u\left(x\right))-I_0\left(x\right)|+|\&dtri;u\left|\right\}\mathrm{dx}$

Wherein, θ is a very little constant, in iterative process, only have when u and v approaching in, can make above formula energy functional obtain minimum value;

This energy functional is write as to the mathematical model of concrete two dimensional form:

$E_{\mathrm{\&theta;}}={\&Integral;}_{\mathrm{\&Omega;}}\left\{\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\right|\&dtri;u_i|+\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\frac{1}{2\mathrm{\&theta;}}{(u_i-v_i)}^2+\mathrm{\&lambda;}|\mathrm{\&rho;}\left(v\right)\left|\right\}\mathrm{dx}$

Wherein, u ₁and u ₂represent respectively the light stream of x direction and the light stream of y direction; Utilize the alternately alternative manner of primal dual algorithm to optimize above-mentioned model, can try to achieve light stream;

Step 4, utilizes alternately alternative manner to solve energy model;

(1) for i(i=1,2), fixing v ₁and v ₂, solve u ₁and u ₂, optimize with drag:

$\underset{u_i}{\min }{\&Integral;}_{\mathrm{\&Omega;}}\left\{\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\right|\&dtri;u_i|+\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\frac{1}{2\mathrm{\&theta;}}{(u_i-v_i)}^2\}\mathrm{dx}$

This model is the image denoising model based on ROF, solves the following formula of this model utilization:

u _i＝v _i+θdiv?p _i

U _ithe dual variable of (i=1,2) is p _i(i=1,2), the iterative formula that solves dual variable is as follows:

${\tilde{p}}^{n+1}=p^n+\frac{\mathrm{\&tau;}}{\mathrm{\&theta;}}(\&dtri;(I+\mathrm{\&theta;}\mathrm{div}{p}^n))$

$p^{n+1}=\frac{{\tilde{p}}^{n+1}}{\max \{1,|{\tilde{p}}^{n+1}\left|\right\}}$

Wherein, p ⁰=0, τ≤1/4;

(2) for i(i=1,2), fixing u ₁and u ₂, solve v ₁and v ₂, optimize with drag:

$\underset{v}{\min }{\&Integral;}_{\mathrm{\&Omega;}}\{\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}\frac{1}{2\mathrm{\&theta;}}{(u_i-v_i)}^2+\mathrm{\&lambda;}|\mathrm{\&rho;}\left(v\right)\left|\right\}\mathrm{dx}$

Solve the method that this model need to use pointwise to solve, its solution procedure is used following threshold method:

$v=u+\left\{\begin{array}{cc}\mathrm{\&lambda;\&theta;}\&dtri;I_1& \mathrm{if\&rho;}\left(u\right)<-\mathrm{\&lambda;\&theta;}{|\&dtri;I_1|}^2\\ -\mathrm{\&lambda;\&theta;}\&dtri;I_1& \mathrm{if\&rho;}\left(u\right)>\mathrm{\&lambda;\&theta;}{|\&dtri;I_1|}^2\\ -\mathrm{\&rho;}\left(u\right)\&dtri;I_1/{|\&dtri;I_1|}^2& \mathrm{if}\left|\mathrm{\&rho;}\left(u\right)\right|\&le;\mathrm{\&lambda;\&theta;}{|\&dtri;I_1|}^2\end{array}\right..$

2. according to claim 1 a kind of based on TV-L ¹the robust Optical flow estimation method of Variation Model, is characterized in that, in order to solve pixel light stream value gradient (▽ u) _i,jsolve difficult problem with dual variable div p, take the method for value solving of discretize to solve energy model, method is as follows:

Utilize consequent difference method to solve div p, utilize preceding paragraph difference method to calculate the gradient of light stream;

For one wide be that in N, the high view data scope for M, pixel (i, j) is located gradient

preceding paragraph difference discrete be expressed as:

${(\&dtri;u)}_{i,j}^1=\left\{\begin{array}{cc}{u}_{i+1,j}-u_{i,j},& \mathrm{ifi}<N\\ 0,& \mathrm{ifi}=N\end{array}\right.$

${(\&dtri;u)}_{i,j}^2=\left\{\begin{array}{cc}{u}_{i,j+1}-u_{i,j},& \mathrm{ifi}<M\\ 0,& \mathrm{ifi}=M\end{array}\right.$

The consequent difference discrete of divergence operator is expressed as:

${\left(\mathrm{div}p\right)}_{i,j}=\left\{\begin{array}{cc}{p}_{i,j}^1-p_{i-1,j}^1,& \mathrm{if}1<i<N\\ p_{i,j}^1,& \mathrm{ifi}=1\\ -p_{i-1,j}^1,& \mathrm{ifi}=N\end{array}\right.+\left\{\begin{array}{cc}{p}_{i,j}^2-p_{i,j-1}^2,& \mathrm{if}1<j<M\\ p_{i,j}^2,& \mathrm{ifj}=1\\ -p_{i,j-1}^2,& \mathrm{ifj}=M\end{array}\right..$

3. according to claim 1 a kind of based on TV-L ¹the robust Optical flow estimation method of Variation Model, it is characterized in that, numerical evaluation is carried out in regular grid, utilize computing power that GPU is powerful and and parallel processing capability, solution procedure to energy functional is accelerated, and makes optical flow computation have good real-time.

4. according to claim 1 a kind of based on TV-L ¹the robust Optical flow estimation method of Variation Model, it is characterized in that, for solving single order Taylor expansion due to data item, be only applicable to the problem that the projection coordinate of little displacement is restricted light stream solution procedure, the solution procedure of energy model has been taked by slightly to smart pyramid algorith: generate the pyramid diagram picture that a series of thicknesses are different, on the thickest tomographic image yardstick, utilize TV-L ¹variation Model is asked for light stream, and the initial solution using this solution on the thinner image of lower one deck, repeats this step, until the thinnest resolution tomographic image, i.e. original image.

Images