|
|
|
@ -209,10 +209,10 @@ namespace tvl1flow |
|
|
|
|
|
|
|
|
|
__global__ void estimateUKernel(const PtrStepSzf I1wx, const PtrStepf I1wy, |
|
|
|
|
const PtrStepf grad, const PtrStepf rho_c, |
|
|
|
|
const PtrStepf p11, const PtrStepf p12, |
|
|
|
|
const PtrStepf p21, const PtrStepf p22, |
|
|
|
|
const PtrStepf p31, const PtrStepf p32, |
|
|
|
|
PtrStepf u1, PtrStepf u2, PtrStepf u3, PtrStepf error, |
|
|
|
|
const PtrStepf p11, const PtrStepf p12, |
|
|
|
|
const PtrStepf p21, const PtrStepf p22, |
|
|
|
|
const PtrStepf p31, const PtrStepf p32, |
|
|
|
|
PtrStepf u1, PtrStepf u2, PtrStepf u3, PtrStepf error, |
|
|
|
|
const float l_t, const float theta, const float gamma, const bool calcError) |
|
|
|
|
{ |
|
|
|
|
const int x = blockIdx.x * blockDim.x + threadIdx.x; |
|
|
|
@ -225,8 +225,8 @@ namespace tvl1flow |
|
|
|
|
const float I1wyVal = I1wy(y, x); |
|
|
|
|
const float gradVal = grad(y, x); |
|
|
|
|
const float u1OldVal = u1(y, x); |
|
|
|
|
const float u2OldVal = u2(y, x); |
|
|
|
|
const float u3OldVal = u3(y, x); |
|
|
|
|
const float u2OldVal = u2(y, x); |
|
|
|
|
const float u3OldVal = u3(y, x); |
|
|
|
|
|
|
|
|
|
const float rho = rho_c(y, x) + (I1wxVal * u1OldVal + I1wyVal * u2OldVal + gamma * u3OldVal); |
|
|
|
|
|
|
|
|
@ -234,61 +234,61 @@ namespace tvl1flow |
|
|
|
|
|
|
|
|
|
float d1 = 0.0f; |
|
|
|
|
float d2 = 0.0f; |
|
|
|
|
float d3 = 0.0f; |
|
|
|
|
float d3 = 0.0f; |
|
|
|
|
|
|
|
|
|
if (rho < -l_t * gradVal) |
|
|
|
|
{ |
|
|
|
|
d1 = l_t * I1wxVal; |
|
|
|
|
d2 = l_t * I1wyVal; |
|
|
|
|
d3 = l_t * gamma; |
|
|
|
|
d3 = l_t * gamma; |
|
|
|
|
} |
|
|
|
|
else if (rho > l_t * gradVal) |
|
|
|
|
{ |
|
|
|
|
d1 = -l_t * I1wxVal; |
|
|
|
|
d2 = -l_t * I1wyVal; |
|
|
|
|
d3 = -l_t * gamma; |
|
|
|
|
d3 = -l_t * gamma; |
|
|
|
|
} |
|
|
|
|
else if (gradVal > numeric_limits<float>::epsilon()) |
|
|
|
|
{ |
|
|
|
|
const float fi = -rho / gradVal; |
|
|
|
|
d1 = fi * I1wxVal; |
|
|
|
|
d2 = fi * I1wyVal; |
|
|
|
|
d3 = fi * gamma; |
|
|
|
|
d3 = fi * gamma; |
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
const float v1 = u1OldVal + d1; |
|
|
|
|
const float v2 = u2OldVal + d2; |
|
|
|
|
const float v3 = u3OldVal + d3; |
|
|
|
|
const float v2 = u2OldVal + d2; |
|
|
|
|
const float v3 = u3OldVal + d3; |
|
|
|
|
|
|
|
|
|
// compute the divergence of the dual variable (p1, p2) |
|
|
|
|
|
|
|
|
|
const float div_p1 = divergence(p11, p12, y, x); |
|
|
|
|
const float div_p2 = divergence(p21, p22, y, x); |
|
|
|
|
const float div_p3 = divergence(p31, p32, y, x); |
|
|
|
|
const float div_p2 = divergence(p21, p22, y, x); |
|
|
|
|
const float div_p3 = divergence(p31, p32, y, x); |
|
|
|
|
|
|
|
|
|
// estimate the values of the optical flow (u1, u2) |
|
|
|
|
|
|
|
|
|
const float u1NewVal = v1 + theta * div_p1; |
|
|
|
|
const float u2NewVal = v2 + theta * div_p2; |
|
|
|
|
const float u3NewVal = v3 + theta * div_p3; |
|
|
|
|
const float u2NewVal = v2 + theta * div_p2; |
|
|
|
|
const float u3NewVal = v3 + theta * div_p3; |
|
|
|
|
|
|
|
|
|
u1(y, x) = u1NewVal; |
|
|
|
|
u2(y, x) = u2NewVal; |
|
|
|
|
u3(y, x) = u3NewVal; |
|
|
|
|
u2(y, x) = u2NewVal; |
|
|
|
|
u3(y, x) = u3NewVal; |
|
|
|
|
|
|
|
|
|
if (calcError) |
|
|
|
|
{ |
|
|
|
|
const float n1 = (u1OldVal - u1NewVal) * (u1OldVal - u1NewVal); |
|
|
|
|
const float n2 = (u2OldVal - u2NewVal) * (u2OldVal - u2NewVal); |
|
|
|
|
const float n3 = 0;// (u3OldVal - u3NewVal) * (u3OldVal - u3NewVal); |
|
|
|
|
const float n2 = (u2OldVal - u2NewVal) * (u2OldVal - u2NewVal); |
|
|
|
|
const float n3 = 0;// (u3OldVal - u3NewVal) * (u3OldVal - u3NewVal); |
|
|
|
|
error(y, x) = n1 + n2 + n3; |
|
|
|
|
} |
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
void estimateU(PtrStepSzf I1wx, PtrStepSzf I1wy, |
|
|
|
|
PtrStepSzf grad, PtrStepSzf rho_c, |
|
|
|
|
PtrStepSzf p11, PtrStepSzf p12, PtrStepSzf p21, PtrStepSzf p22, PtrStepSzf p31, PtrStepSzf p32, |
|
|
|
|
PtrStepSzf u1, PtrStepSzf u2, PtrStepSzf u3, PtrStepSzf error, |
|
|
|
|
PtrStepSzf p11, PtrStepSzf p12, PtrStepSzf p21, PtrStepSzf p22, PtrStepSzf p31, PtrStepSzf p32, |
|
|
|
|
PtrStepSzf u1, PtrStepSzf u2, PtrStepSzf u3, PtrStepSzf error, |
|
|
|
|
float l_t, float theta, float gamma, bool calcError) |
|
|
|
|
{ |
|
|
|
|
const dim3 block(32, 8); |
|
|
|
@ -306,8 +306,8 @@ namespace tvl1flow |
|
|
|
|
|
|
|
|
|
namespace tvl1flow |
|
|
|
|
{ |
|
|
|
|
__global__ void estimateDualVariablesKernel(const PtrStepSzf u1, const PtrStepf u2, const PtrStepSzf u3, |
|
|
|
|
PtrStepf p11, PtrStepf p12, PtrStepf p21, PtrStepf p22, PtrStepf p31, PtrStepf p32, const float taut) |
|
|
|
|
__global__ void estimateDualVariablesKernel(const PtrStepSzf u1, const PtrStepf u2, const PtrStepSzf u3, |
|
|
|
|
PtrStepf p11, PtrStepf p12, PtrStepf p21, PtrStepf p22, PtrStepf p31, PtrStepf p32, const float taut) |
|
|
|
|
{ |
|
|
|
|
const int x = blockIdx.x * blockDim.x + threadIdx.x; |
|
|
|
|
const int y = blockIdx.y * blockDim.y + threadIdx.y; |
|
|
|
@ -321,26 +321,26 @@ namespace tvl1flow |
|
|
|
|
const float u2x = u2(y, ::min(x + 1, u1.cols - 1)) - u2(y, x); |
|
|
|
|
const float u2y = u2(::min(y + 1, u1.rows - 1), x) - u2(y, x); |
|
|
|
|
|
|
|
|
|
const float u3x = u3(y, ::min(x + 1, u1.cols - 1)) - u3(y, x); |
|
|
|
|
const float u3y = u3(::min(y + 1, u1.rows - 1), x) - u3(y, x); |
|
|
|
|
const float u3x = u3(y, ::min(x + 1, u1.cols - 1)) - u3(y, x); |
|
|
|
|
const float u3y = u3(::min(y + 1, u1.rows - 1), x) - u3(y, x); |
|
|
|
|
|
|
|
|
|
const float g1 = ::hypotf(u1x, u1y); |
|
|
|
|
const float g2 = ::hypotf(u2x, u2y); |
|
|
|
|
const float g3 = ::hypotf(u3x, u3y); |
|
|
|
|
const float g2 = ::hypotf(u2x, u2y); |
|
|
|
|
const float g3 = ::hypotf(u3x, u3y); |
|
|
|
|
|
|
|
|
|
const float ng1 = 1.0f + taut * g1; |
|
|
|
|
const float ng2 = 1.0f + taut * g2; |
|
|
|
|
const float ng3 = 1.0f + taut * g3; |
|
|
|
|
const float ng2 = 1.0f + taut * g2; |
|
|
|
|
const float ng3 = 1.0f + taut * g3; |
|
|
|
|
|
|
|
|
|
p11(y, x) = (p11(y, x) + taut * u1x) / ng1; |
|
|
|
|
p12(y, x) = (p12(y, x) + taut * u1y) / ng1; |
|
|
|
|
p21(y, x) = (p21(y, x) + taut * u2x) / ng2; |
|
|
|
|
p22(y, x) = (p22(y, x) + taut * u2y) / ng2; |
|
|
|
|
p31(y, x) = (p31(y, x) + taut * u3x) / ng3; |
|
|
|
|
p32(y, x) = (p32(y, x) + taut * u3y) / ng3; |
|
|
|
|
p22(y, x) = (p22(y, x) + taut * u2y) / ng2; |
|
|
|
|
p31(y, x) = (p31(y, x) + taut * u3x) / ng3; |
|
|
|
|
p32(y, x) = (p32(y, x) + taut * u3y) / ng3; |
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
void estimateDualVariables(PtrStepSzf u1, PtrStepSzf u2, PtrStepSzf u3, PtrStepSzf p11, PtrStepSzf p12, PtrStepSzf p21, PtrStepSzf p22, PtrStepSzf p31, PtrStepSzf p32, float taut) |
|
|
|
|
void estimateDualVariables(PtrStepSzf u1, PtrStepSzf u2, PtrStepSzf u3, PtrStepSzf p11, PtrStepSzf p12, PtrStepSzf p21, PtrStepSzf p22, PtrStepSzf p31, PtrStepSzf p32, float taut) |
|
|
|
|
{ |
|
|
|
|
const dim3 block(32, 8); |
|
|
|
|
const dim3 grid(divUp(u1.cols, block.x), divUp(u1.rows, block.y)); |
|
|
|
|