replaced tabs by spaces

modules/video/src/tvl1flow.cpp

@@ -111,7 +111,7 @@ protected:
     int medianFiltering;
 
 private:
-    void procOneScale(const Mat_<float>& I0, const Mat_<float>& I1, Mat_<float>& u1, Mat_<float>& u2);
+   void procOneScale(const Mat_<float>& I0, const Mat_<float>& I1, Mat_<float>& u1, Mat_<float>& u2, Mat_<float>& u3);
 
     bool procOneScale_ocl(const UMat& I0, const UMat& I1, UMat& u1, UMat& u2);
 
@@ -139,6 +139,7 @@ private:
 
         Mat_<float> v1_buf;
       Mat_<float> v2_buf;
+      Mat_<float> v3_buf;
 
         Mat_<float> p11_buf;
         Mat_<float> p12_buf;
@@ -350,7 +351,7 @@ OpticalFlowDual_TVL1::OpticalFlowDual_TVL1()
     nscales        = 5;
     warps          = 5;
     epsilon        = 0.01;
-	gamma		   = 1.;
+   gamma         = 0.01;
     innerIterations = 30;
     outerIterations = 10;
     useInitialFlow = false;
@@ -406,6 +407,7 @@ void OpticalFlowDual_TVL1::calc(InputArray _I0, InputArray _I1, InputOutputArray
 
     dm.v1_buf.create(I0.size());
    dm.v2_buf.create(I0.size());
+   dm.v3_buf.create(I0.size());
 
     dm.p11_buf.create(I0.size());
     dm.p12_buf.create(I0.size());
@@ -462,7 +464,7 @@ void OpticalFlowDual_TVL1::calc(InputArray _I0, InputArray _I1, InputOutputArray
     for (int s = nscales - 1; s >= 0; --s)
     {
         // compute the optical flow at the current scale
-        procOneScale(dm.I0s[s], dm.I1s[s], dm.u1s[s], dm.u2s[s]);
+        procOneScale(dm.I0s[s], dm.I1s[s], dm.u1s[s], dm.u2s[s], dm.u3s[s]);
 
         // if this was the last scale, finish now
         if (s == 0)
@@ -1027,36 +1029,45 @@ void estimateV(const Mat_<float>& I1wx, const Mat_<float>& I1wy, const Mat_<floa
 ////////////////////////////////////////////////////////////
 // estimateU
 
-float estimateU(const Mat_<float>& v1, const Mat_<float>& v2, const Mat_<float>& div_p1, const Mat_<float>& div_p2, Mat_<float>& u1, Mat_<float>& u2, float theta)
+float estimateU(const Mat_<float>& v1, const Mat_<float>& v2, const Mat_<float>& v3, 
+            const Mat_<float>& div_p1, const Mat_<float>& div_p2, const Mat_<float>& div_p3, 
+            Mat_<float>& u1, Mat_<float>& u2, Mat_<float>& u3,
+            float theta, float gamma)
 {
     CV_DbgAssert( v2.size() == v1.size() );
+   CV_DbgAssert( v3.size() == v1.size() );
     CV_DbgAssert( div_p1.size() == v1.size() );
     CV_DbgAssert( div_p2.size() == v1.size() );
+    CV_DbgAssert( div_p3.size() == v1.size() );
     CV_DbgAssert( u1.size() == v1.size() );
     CV_DbgAssert( u2.size() == v1.size() );
+    CV_DbgAssert( u3.size() == v1.size() );
 
     float error = 0.0f;
     for (int y = 0; y < v1.rows; ++y)
     {
         const float* v1Row = v1[y];
       const float* v2Row = v2[y];
+      const float* v3Row = v3[y];
         const float* divP1Row = div_p1[y];
       const float* divP2Row = div_p2[y];
+      const float* divP3Row = div_p3[y];
 
         float* u1Row = u1[y];
       float* u2Row = u2[y];
+      float* u3Row = u3[y];
 
         for (int x = 0; x < v1.cols; ++x)
         {
             const float u1k = u1Row[x];
          const float u2k = u2Row[x];
+         const float u3k = u3Row[x];
 
             u1Row[x] = v1Row[x] + theta * divP1Row[x];
          u2Row[x] = v2Row[x] + theta * divP2Row[x];
+         u3Row[x] = v3Row[x] + theta * divP3Row[x];
          
-			//u3
-
-            error += (u1Row[x] - u1k) * (u1Row[x] - u1k) + (u2Row[x] - u2k) * (u2Row[x] - u2k);
+         error += (u1Row[x] - u1k) * (u1Row[x] - u1k) + (u2Row[x] - u2k) * (u2Row[x] - u2k) + (u3Row[x] - u3k) * (u3Row[x] - u3k);
         }
     }
 
@@ -1074,10 +1085,14 @@ struct EstimateDualVariablesBody : ParallelLoopBody
     Mat_<float> u1y;
     Mat_<float> u2x;
    Mat_<float> u2y;
+   Mat_<float> u3x;
+   Mat_<float> u3y;
     mutable Mat_<float> p11;
     mutable Mat_<float> p12;
     mutable Mat_<float> p21;
    mutable Mat_<float> p22;
+   mutable Mat_<float> p31;
+   mutable Mat_<float> p32;
     float taut;
 };
 
@@ -1089,38 +1104,55 @@ void EstimateDualVariablesBody::operator() (const Range& range) const
         const float* u1yRow = u1y[y];
         const float* u2xRow = u2x[y];
       const float* u2yRow = u2y[y];
+      const float* u3xRow = u3x[y];
+      const float* u3yRow = u3y[y];
 
         float* p11Row = p11[y];
         float* p12Row = p12[y];
         float* p21Row = p21[y];
       float* p22Row = p22[y];
+      float* p31Row = p31[y];
+      float* p32Row = p32[y];
 
         for (int x = 0; x < u1x.cols; ++x)
         {
             const float g1 = static_cast<float>(hypot(u1xRow[x], u1yRow[x]));
          const float g2 = static_cast<float>(hypot(u2xRow[x], u2yRow[x]));
+         const float g3 = static_cast<float>(hypot(u3xRow[x], u3yRow[x]));
 
             const float ng1  = 1.0f + taut * g1;
          const float ng2 =  1.0f + taut * g2;
+         const float ng3 =  1.0f + taut * g3;
 
             p11Row[x] = (p11Row[x] + taut * u1xRow[x]) / ng1;
             p12Row[x] = (p12Row[x] + taut * u1yRow[x]) / ng1;
             p21Row[x] = (p21Row[x] + taut * u2xRow[x]) / ng2;
          p22Row[x] = (p22Row[x] + taut * u2yRow[x]) / ng2;
+         p31Row[x] = (p31Row[x] + taut * u3xRow[x]) / ng3;
+         p32Row[x] = (p32Row[x] + taut * u3yRow[x]) / ng3;
         }
     }
 }
 
-void estimateDualVariables(const Mat_<float>& u1x, const Mat_<float>& u1y, const Mat_<float>& u2x, const Mat_<float>& u2y,
-                           Mat_<float>& p11, Mat_<float>& p12, Mat_<float>& p21, Mat_<float>& p22, float taut)
+void estimateDualVariables(const Mat_<float>& u1x, const Mat_<float>& u1y, 
+                     const Mat_<float>& u2x, const Mat_<float>& u2y,
+                     const Mat_<float>& u3x, const Mat_<float>& u3y,
+                           Mat_<float>& p11, Mat_<float>& p12, 
+                     Mat_<float>& p21, Mat_<float>& p22, 
+                     Mat_<float>& p31, Mat_<float>& p32,
+                     float taut)
 {
     CV_DbgAssert( u1y.size() == u1x.size() );
     CV_DbgAssert( u2x.size() == u1x.size() );
+    CV_DbgAssert( u3x.size() == u1x.size() );
     CV_DbgAssert( u2y.size() == u1x.size() );
+    CV_DbgAssert( u3y.size() == u1x.size() );
     CV_DbgAssert( p11.size() == u1x.size() );
     CV_DbgAssert( p12.size() == u1x.size() );
     CV_DbgAssert( p21.size() == u1x.size() );
     CV_DbgAssert( p22.size() == u1x.size() );
+    CV_DbgAssert( p31.size() == u1x.size() );
+    CV_DbgAssert( p32.size() == u1x.size() );
 
     EstimateDualVariablesBody body;
 
@@ -1128,10 +1160,14 @@ void estimateDualVariables(const Mat_<float>& u1x, const Mat_<float>& u1y, const
     body.u1y = u1y;
     body.u2x = u2x;
    body.u2y = u2y;
+   body.u3x = u3x;
+   body.u3y = u3y;
     body.p11 = p11;
     body.p12 = p12;
     body.p21 = p21;
    body.p22 = p22;
+   body.p31 = p31;
+   body.p32 = p32;
     body.taut = taut;
 
     parallel_for_(Range(0, u1x.rows), body);
@@ -1225,7 +1261,7 @@ bool OpticalFlowDual_TVL1::procOneScale_ocl(const UMat& I0, const UMat& I1, UMat
     return true;
 }
 
-void OpticalFlowDual_TVL1::procOneScale(const Mat_<float>& I0, const Mat_<float>& I1, Mat_<float>& u1, Mat_<float>& u2)
+void OpticalFlowDual_TVL1::procOneScale(const Mat_<float>& I0, const Mat_<float>& I1, Mat_<float>& u1, Mat_<float>& u2, Mat_<float>& u3)
 {
     const float scaledEpsilon = static_cast<float>(epsilon * epsilon * I0.size().area());
 
@@ -1250,6 +1286,7 @@ void OpticalFlowDual_TVL1::procOneScale(const Mat_<float>& I0, const Mat_<float>
 
     Mat_<float> v1 = dm.v1_buf(Rect(0, 0, I0.cols, I0.rows));
    Mat_<float> v2 = dm.v2_buf(Rect(0, 0, I0.cols, I0.rows));
+   Mat_<float> v3 = dm.v3_buf(Rect(0, 0, I0.cols, I0.rows));
 
     Mat_<float> p11 = dm.p11_buf(Rect(0, 0, I0.cols, I0.rows));
     Mat_<float> p12 = dm.p12_buf(Rect(0, 0, I0.cols, I0.rows));
@@ -1266,7 +1303,7 @@ void OpticalFlowDual_TVL1::procOneScale(const Mat_<float>& I0, const Mat_<float>
 
     Mat_<float> div_p1 = dm.div_p1_buf(Rect(0, 0, I0.cols, I0.rows));
    Mat_<float> div_p2 = dm.div_p2_buf(Rect(0, 0, I0.cols, I0.rows));
-	Mat_<float> div_p3 = dm.div_p2_buf(Rect(0, 0, I0.cols, I0.rows));
+   Mat_<float> div_p3 = dm.div_p3_buf(Rect(0, 0, I0.cols, I0.rows));
 
     Mat_<float> u1x = dm.u1x_buf(Rect(0, 0, I0.cols, I0.rows));
     Mat_<float> u1y = dm.u1y_buf(Rect(0, 0, I0.cols, I0.rows));
@@ -1298,21 +1335,23 @@ void OpticalFlowDual_TVL1::procOneScale(const Mat_<float>& I0, const Mat_<float>
             for (int n_inner = 0; error > scaledEpsilon && n_inner < innerIterations; ++n_inner)
             {
                 // estimate the values of the variable (v1, v2) (thresholding operator TH)
-                estimateV(I1wx, I1wy, u1, u2, grad, rho_c, v1, v2, l_t);
+                estimateV(I1wx, I1wy, u1, u2, u3, grad, rho_c, v1, v2, v3, l_t, gamma);
 
-                // compute the divergence of the dual variable (p1, p2)
+                // compute the divergence of the dual variable (p1, p2, p3)
                 divergence(p11, p12, div_p1);
             divergence(p21, p22, div_p2);
+            divergence(p31, p32, div_p3);
 
                 // estimate the values of the optical flow (u1, u2)
-                error = estimateU(v1, v2, div_p1, div_p2, u1, u2, static_cast<float>(theta));
+                error = estimateU(v1, v2, v3, div_p1, div_p2, div_p3, u1, u2, u3, static_cast<float>(theta), gamma);
 
                 // compute the gradient of the optical flow (Du1, Du2)
                 forwardGradient(u1, u1x, u1y);
             forwardGradient(u2, u2x, u2y);
+            forwardGradient(u3, u3x, u3y);
 
-                // estimate the values of the dual variable (p1, p2)
-                estimateDualVariables(u1x, u1y, u2x, u2y, p11, p12, p21, p22, taut);
+                // estimate the values of the dual variable (p1, p2, p3)
+                estimateDualVariables(u1x, u1y, u2x, u2y, u3x, u3y, p11, p12, p21, p22, p31, p32, taut);
             }
         }
     }
@@ -1402,6 +1441,8 @@ CV_INIT_ALGORITHM(OpticalFlowDual_TVL1, "DenseOpticalFlow.DualTVL1",
                                        "inner iterations (between outlier filtering) used in the numerical scheme");
                   obj.info()->addParam(obj, "outerIterations", obj.outerIterations, false, 0, 0,
                                        "outer iterations (number of inner loops) used in the numerical scheme");
+                  obj.info()->addParam(obj, "gamma", obj.gamma, false, 0, 0,
+                                       "coefficient for additional Ali term");
                   obj.info()->addParam(obj, "useInitialFlow", obj.useInitialFlow))
 
 } // namespace