|
|
|
|
@ -44,7 +44,7 @@ |
|
|
|
|
//M*/ |
|
|
|
|
|
|
|
|
|
__kernel void centeredGradientKernel(__global const float* src, int src_col, int src_row, int src_step, |
|
|
|
|
__global float* dx, __global float* dy, int dx_step) |
|
|
|
|
__global float* dx, __global float* dy, int dx_step) |
|
|
|
|
{ |
|
|
|
|
int x = get_global_id(0); |
|
|
|
|
int y = get_global_id(1); |
|
|
|
|
@ -53,13 +53,6 @@ __global float* dx, __global float* dy, int dx_step) |
|
|
|
|
{ |
|
|
|
|
int src_x1 = (x + 1) < (src_col -1)? (x + 1) : (src_col - 1); |
|
|
|
|
int src_x2 = (x - 1) > 0 ? (x -1) : 0; |
|
|
|
|
|
|
|
|
|
//if(src[y * src_step + src_x1] == src[y * src_step+ src_x2]) |
|
|
|
|
//{ |
|
|
|
|
// printf("y = %d\n", y); |
|
|
|
|
// printf("src_x1 = %d\n", src_x1); |
|
|
|
|
// printf("src_x2 = %d\n", src_x2); |
|
|
|
|
//} |
|
|
|
|
dx[y * dx_step+ x] = 0.5f * (src[y * src_step + src_x1] - src[y * src_step+ src_x2]); |
|
|
|
|
|
|
|
|
|
int src_y1 = (y+1) < (src_row - 1) ? (y + 1) : (src_row - 1); |
|
|
|
|
@ -97,24 +90,24 @@ __kernel void warpBackwardKernel(__global const float* I0, int I0_step, int I0_c |
|
|
|
|
int u2_offset_x, |
|
|
|
|
int u2_offset_y) |
|
|
|
|
{ |
|
|
|
|
const int x = get_global_id(0); |
|
|
|
|
const int y = get_global_id(1); |
|
|
|
|
int x = get_global_id(0); |
|
|
|
|
int y = get_global_id(1); |
|
|
|
|
|
|
|
|
|
if(x < I0_col&&y < I0_row) |
|
|
|
|
{ |
|
|
|
|
//const float u1Val = u1(y, x); |
|
|
|
|
const float u1Val = u1[(y + u1_offset_y) * u1_step + x + u1_offset_x]; |
|
|
|
|
//const float u2Val = u2(y, x); |
|
|
|
|
const float u2Val = u2[(y + u2_offset_y) * u2_step + x + u2_offset_x]; |
|
|
|
|
//float u1Val = u1(y, x); |
|
|
|
|
float u1Val = u1[(y + u1_offset_y) * u1_step + x + u1_offset_x]; |
|
|
|
|
//float u2Val = u2(y, x); |
|
|
|
|
float u2Val = u2[(y + u2_offset_y) * u2_step + x + u2_offset_x]; |
|
|
|
|
|
|
|
|
|
const float wx = x + u1Val; |
|
|
|
|
const float wy = y + u2Val; |
|
|
|
|
float wx = x + u1Val; |
|
|
|
|
float wy = y + u2Val; |
|
|
|
|
|
|
|
|
|
const int xmin = ceil(wx - 2.0f); |
|
|
|
|
const int xmax = floor(wx + 2.0f); |
|
|
|
|
int xmin = ceil(wx - 2.0f); |
|
|
|
|
int xmax = floor(wx + 2.0f); |
|
|
|
|
|
|
|
|
|
const int ymin = ceil(wy - 2.0f); |
|
|
|
|
const int ymax = floor(wy + 2.0f); |
|
|
|
|
int ymin = ceil(wy - 2.0f); |
|
|
|
|
int ymax = floor(wy + 2.0f); |
|
|
|
|
|
|
|
|
|
float sum = 0.0f; |
|
|
|
|
float sumx = 0.0f; |
|
|
|
|
@ -126,7 +119,7 @@ __kernel void warpBackwardKernel(__global const float* I0, int I0_step, int I0_c |
|
|
|
|
{ |
|
|
|
|
for (int cx = xmin; cx <= xmax; ++cx) |
|
|
|
|
{ |
|
|
|
|
const float w = bicubicCoeff(wx - cx) * bicubicCoeff(wy - cy); |
|
|
|
|
float w = bicubicCoeff(wx - cx) * bicubicCoeff(wy - cy); |
|
|
|
|
|
|
|
|
|
//sum += w * tex2D(tex_I1 , cx, cy); |
|
|
|
|
int2 cood = (int2)(cx, cy); |
|
|
|
|
@ -140,30 +133,30 @@ __kernel void warpBackwardKernel(__global const float* I0, int I0_step, int I0_c |
|
|
|
|
} |
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
const float coeff = 1.0f / wsum; |
|
|
|
|
float coeff = 1.0f / wsum; |
|
|
|
|
|
|
|
|
|
const float I1wVal = sum * coeff; |
|
|
|
|
const float I1wxVal = sumx * coeff; |
|
|
|
|
const float I1wyVal = sumy * coeff; |
|
|
|
|
float I1wVal = sum * coeff; |
|
|
|
|
float I1wxVal = sumx * coeff; |
|
|
|
|
float I1wyVal = sumy * coeff; |
|
|
|
|
|
|
|
|
|
I1w[y * I1w_step + x] = I1wVal; |
|
|
|
|
I1wx[y * I1w_step + x] = I1wxVal; |
|
|
|
|
I1wy[y * I1w_step + x] = I1wyVal; |
|
|
|
|
|
|
|
|
|
const float Ix2 = I1wxVal * I1wxVal; |
|
|
|
|
const float Iy2 = I1wyVal * I1wyVal; |
|
|
|
|
float Ix2 = I1wxVal * I1wxVal; |
|
|
|
|
float Iy2 = I1wyVal * I1wyVal; |
|
|
|
|
|
|
|
|
|
// store the |Grad(I1)|^2 |
|
|
|
|
grad[y * I1w_step + x] = Ix2 + Iy2; |
|
|
|
|
|
|
|
|
|
// compute the constant part of the rho function |
|
|
|
|
const float I0Val = I0[y * I0_step + x]; |
|
|
|
|
float I0Val = I0[y * I0_step + x]; |
|
|
|
|
rho[y * I1w_step + x] = I1wVal - I1wxVal * u1Val - I1wyVal * u2Val - I0Val; |
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
static float readImage(__global const float *image, const int x, const int y, const int rows, const int cols, const int elemCntPerRow) |
|
|
|
|
static float readImage(__global float *image, int x, int y, int rows, int cols, int elemCntPerRow) |
|
|
|
|
{ |
|
|
|
|
int i0 = clamp(x, 0, cols - 1); |
|
|
|
|
int j0 = clamp(y, 0, rows - 1); |
|
|
|
|
@ -185,24 +178,24 @@ __kernel void warpBackwardKernelNoImage2d(__global const float* I0, int I0_step, |
|
|
|
|
int I1_step, |
|
|
|
|
int I1x_step) |
|
|
|
|
{ |
|
|
|
|
const int x = get_global_id(0); |
|
|
|
|
const int y = get_global_id(1); |
|
|
|
|
int x = get_global_id(0); |
|
|
|
|
int y = get_global_id(1); |
|
|
|
|
|
|
|
|
|
if(x < I0_col&&y < I0_row) |
|
|
|
|
{ |
|
|
|
|
//const float u1Val = u1(y, x); |
|
|
|
|
const float u1Val = u1[y * u1_step + x]; |
|
|
|
|
//const float u2Val = u2(y, x); |
|
|
|
|
const float u2Val = u2[y * u2_step + x]; |
|
|
|
|
//float u1Val = u1(y, x); |
|
|
|
|
float u1Val = u1[y * u1_step + x]; |
|
|
|
|
//float u2Val = u2(y, x); |
|
|
|
|
float u2Val = u2[y * u2_step + x]; |
|
|
|
|
|
|
|
|
|
const float wx = x + u1Val; |
|
|
|
|
const float wy = y + u2Val; |
|
|
|
|
float wx = x + u1Val; |
|
|
|
|
float wy = y + u2Val; |
|
|
|
|
|
|
|
|
|
const int xmin = ceil(wx - 2.0f); |
|
|
|
|
const int xmax = floor(wx + 2.0f); |
|
|
|
|
int xmin = ceil(wx - 2.0f); |
|
|
|
|
int xmax = floor(wx + 2.0f); |
|
|
|
|
|
|
|
|
|
const int ymin = ceil(wy - 2.0f); |
|
|
|
|
const int ymax = floor(wy + 2.0f); |
|
|
|
|
int ymin = ceil(wy - 2.0f); |
|
|
|
|
int ymax = floor(wy + 2.0f); |
|
|
|
|
|
|
|
|
|
float sum = 0.0f; |
|
|
|
|
float sumx = 0.0f; |
|
|
|
|
@ -213,7 +206,7 @@ __kernel void warpBackwardKernelNoImage2d(__global const float* I0, int I0_step, |
|
|
|
|
{ |
|
|
|
|
for (int cx = xmin; cx <= xmax; ++cx) |
|
|
|
|
{ |
|
|
|
|
const float w = bicubicCoeff(wx - cx) * bicubicCoeff(wy - cy); |
|
|
|
|
float w = bicubicCoeff(wx - cx) * bicubicCoeff(wy - cy); |
|
|
|
|
|
|
|
|
|
int2 cood = (int2)(cx, cy); |
|
|
|
|
sum += w * readImage(tex_I1, cood.x, cood.y, I0_col, I0_row, I1_step); |
|
|
|
|
@ -223,24 +216,24 @@ __kernel void warpBackwardKernelNoImage2d(__global const float* I0, int I0_step, |
|
|
|
|
} |
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
const float coeff = 1.0f / wsum; |
|
|
|
|
float coeff = 1.0f / wsum; |
|
|
|
|
|
|
|
|
|
const float I1wVal = sum * coeff; |
|
|
|
|
const float I1wxVal = sumx * coeff; |
|
|
|
|
const float I1wyVal = sumy * coeff; |
|
|
|
|
float I1wVal = sum * coeff; |
|
|
|
|
float I1wxVal = sumx * coeff; |
|
|
|
|
float I1wyVal = sumy * coeff; |
|
|
|
|
|
|
|
|
|
I1w[y * I1w_step + x] = I1wVal; |
|
|
|
|
I1wx[y * I1w_step + x] = I1wxVal; |
|
|
|
|
I1wy[y * I1w_step + x] = I1wyVal; |
|
|
|
|
|
|
|
|
|
const float Ix2 = I1wxVal * I1wxVal; |
|
|
|
|
const float Iy2 = I1wyVal * I1wyVal; |
|
|
|
|
float Ix2 = I1wxVal * I1wxVal; |
|
|
|
|
float Iy2 = I1wyVal * I1wyVal; |
|
|
|
|
|
|
|
|
|
// store the |Grad(I1)|^2 |
|
|
|
|
grad[y * I1w_step + x] = Ix2 + Iy2; |
|
|
|
|
|
|
|
|
|
// compute the constant part of the rho function |
|
|
|
|
const float I0Val = I0[y * I0_step + x]; |
|
|
|
|
float I0Val = I0[y * I0_step + x]; |
|
|
|
|
rho[y * I1w_step + x] = I1wVal - I1wxVal * u1Val - I1wyVal * u2Val - I0Val; |
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
@ -253,38 +246,35 @@ __kernel void estimateDualVariablesKernel(__global const float* u1, int u1_col, |
|
|
|
|
__global float* p12, |
|
|
|
|
__global float* p21, |
|
|
|
|
__global float* p22, |
|
|
|
|
const float taut, |
|
|
|
|
float taut, |
|
|
|
|
int u2_step, |
|
|
|
|
int u1_offset_x, |
|
|
|
|
int u1_offset_y, |
|
|
|
|
int u2_offset_x, |
|
|
|
|
int u2_offset_y) |
|
|
|
|
{ |
|
|
|
|
|
|
|
|
|
//const int x = blockIdx.x * blockDim.x + threadIdx.x; |
|
|
|
|
//const int y = blockIdx.y * blockDim.y + threadIdx.y; |
|
|
|
|
const int x = get_global_id(0); |
|
|
|
|
const int y = get_global_id(1); |
|
|
|
|
int x = get_global_id(0); |
|
|
|
|
int y = get_global_id(1); |
|
|
|
|
|
|
|
|
|
if(x < u1_col && y < u1_row) |
|
|
|
|
{ |
|
|
|
|
int src_x1 = (x + 1) < (u1_col - 1) ? (x + 1) : (u1_col - 1); |
|
|
|
|
const float u1x = u1[(y + u1_offset_y) * u1_step + src_x1 + u1_offset_x] - u1[(y + u1_offset_y) * u1_step + x + u1_offset_x]; |
|
|
|
|
float u1x = u1[(y + u1_offset_y) * u1_step + src_x1 + u1_offset_x] - u1[(y + u1_offset_y) * u1_step + x + u1_offset_x]; |
|
|
|
|
|
|
|
|
|
int src_y1 = (y + 1) < (u1_row - 1) ? (y + 1) : (u1_row - 1); |
|
|
|
|
const float u1y = u1[(src_y1 + u1_offset_y) * u1_step + x + u1_offset_x] - u1[(y + u1_offset_y) * u1_step + x + u1_offset_x]; |
|
|
|
|
float u1y = u1[(src_y1 + u1_offset_y) * u1_step + x + u1_offset_x] - u1[(y + u1_offset_y) * u1_step + x + u1_offset_x]; |
|
|
|
|
|
|
|
|
|
int src_x2 = (x + 1) < (u1_col - 1) ? (x + 1) : (u1_col - 1); |
|
|
|
|
const float u2x = u2[(y + u2_offset_y) * u2_step + src_x2 + u2_offset_x] - u2[(y + u2_offset_y) * u2_step + x + u2_offset_x]; |
|
|
|
|
float u2x = u2[(y + u2_offset_y) * u2_step + src_x2 + u2_offset_x] - u2[(y + u2_offset_y) * u2_step + x + u2_offset_x]; |
|
|
|
|
|
|
|
|
|
int src_y2 = (y + 1) < (u1_row - 1) ? (y + 1) : (u1_row - 1); |
|
|
|
|
const float u2y = u2[(src_y2 + u2_offset_y) * u2_step + x + u2_offset_x] - u2[(y + u2_offset_y) * u2_step + x + u2_offset_x]; |
|
|
|
|
float u2y = u2[(src_y2 + u2_offset_y) * u2_step + x + u2_offset_x] - u2[(y + u2_offset_y) * u2_step + x + u2_offset_x]; |
|
|
|
|
|
|
|
|
|
const float g1 = hypot(u1x, u1y); |
|
|
|
|
const float g2 = hypot(u2x, u2y); |
|
|
|
|
float g1 = hypot(u1x, u1y); |
|
|
|
|
float g2 = hypot(u2x, u2y); |
|
|
|
|
|
|
|
|
|
const float ng1 = 1.0f + taut * g1; |
|
|
|
|
const float ng2 = 1.0f + taut * g2; |
|
|
|
|
float ng1 = 1.0f + taut * g1; |
|
|
|
|
float ng2 = 1.0f + taut * g2; |
|
|
|
|
|
|
|
|
|
p11[y * p11_step + x] = (p11[y * p11_step + x] + taut * u1x) / ng1; |
|
|
|
|
p12[y * p11_step + x] = (p12[y * p11_step + x] + taut * u1y) / ng1; |
|
|
|
|
@ -299,8 +289,8 @@ static float divergence(__global const float* v1, __global const float* v2, int |
|
|
|
|
|
|
|
|
|
if (x > 0 && y > 0) |
|
|
|
|
{ |
|
|
|
|
const float v1x = v1[y * v1_step + x] - v1[y * v1_step + x - 1]; |
|
|
|
|
const float v2y = v2[y * v2_step + x] - v2[(y - 1) * v2_step + x]; |
|
|
|
|
float v1x = v1[y * v1_step + x] - v1[y * v1_step + x - 1]; |
|
|
|
|
float v2y = v2[y * v2_step + x] - v2[(y - 1) * v2_step + x]; |
|
|
|
|
return v1x + v2y; |
|
|
|
|
} |
|
|
|
|
else |
|
|
|
|
@ -328,30 +318,25 @@ __kernel void estimateUKernel(__global const float* I1wx, int I1wx_col, int I1wx |
|
|
|
|
__global const float* p22, /*int p22_step,*/ |
|
|
|
|
__global float* u1, int u1_step, |
|
|
|
|
__global float* u2, |
|
|
|
|
__global float* error, const float l_t, const float theta, int u2_step, |
|
|
|
|
__global float* error, float l_t, float theta, int u2_step, |
|
|
|
|
int u1_offset_x, |
|
|
|
|
int u1_offset_y, |
|
|
|
|
int u2_offset_x, |
|
|
|
|
int u2_offset_y, |
|
|
|
|
char calc_error) |
|
|
|
|
{ |
|
|
|
|
|
|
|
|
|
//const int x = blockIdx.x * blockDim.x + threadIdx.x; |
|
|
|
|
//const int y = blockIdx.y * blockDim.y + threadIdx.y; |
|
|
|
|
|
|
|
|
|
int x = get_global_id(0); |
|
|
|
|
int y = get_global_id(1); |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
if(x < I1wx_col && y < I1wx_row) |
|
|
|
|
{ |
|
|
|
|
const float I1wxVal = I1wx[y * I1wx_step + x]; |
|
|
|
|
const float I1wyVal = I1wy[y * I1wx_step + x]; |
|
|
|
|
const float gradVal = grad[y * I1wx_step + x]; |
|
|
|
|
const float u1OldVal = u1[(y + u1_offset_y) * u1_step + x + u1_offset_x]; |
|
|
|
|
const float u2OldVal = u2[(y + u2_offset_y) * u2_step + x + u2_offset_x]; |
|
|
|
|
float I1wxVal = I1wx[y * I1wx_step + x]; |
|
|
|
|
float I1wyVal = I1wy[y * I1wx_step + x]; |
|
|
|
|
float gradVal = grad[y * I1wx_step + x]; |
|
|
|
|
float u1OldVal = u1[(y + u1_offset_y) * u1_step + x + u1_offset_x]; |
|
|
|
|
float u2OldVal = u2[(y + u2_offset_y) * u2_step + x + u2_offset_x]; |
|
|
|
|
|
|
|
|
|
const float rho = rho_c[y * I1wx_step + x] + (I1wxVal * u1OldVal + I1wyVal * u2OldVal); |
|
|
|
|
float rho = rho_c[y * I1wx_step + x] + (I1wxVal * u1OldVal + I1wyVal * u2OldVal); |
|
|
|
|
|
|
|
|
|
// estimate the values of the variable (v1, v2) (thresholding operator TH) |
|
|
|
|
|
|
|
|
|
@ -370,31 +355,31 @@ __kernel void estimateUKernel(__global const float* I1wx, int I1wx_col, int I1wx |
|
|
|
|
} |
|
|
|
|
else if (gradVal > 1.192092896e-07f) |
|
|
|
|
{ |
|
|
|
|
const float fi = -rho / gradVal; |
|
|
|
|
float fi = -rho / gradVal; |
|
|
|
|
d1 = fi * I1wxVal; |
|
|
|
|
d2 = fi * I1wyVal; |
|
|
|
|
} |
|
|
|
|
|
|
|
|
|
const float v1 = u1OldVal + d1; |
|
|
|
|
const float v2 = u2OldVal + d2; |
|
|
|
|
float v1 = u1OldVal + d1; |
|
|
|
|
float v2 = u2OldVal + d2; |
|
|
|
|
|
|
|
|
|
// compute the divergence of the dual variable (p1, p2) |
|
|
|
|
|
|
|
|
|
const float div_p1 = divergence(p11, p12, y, x, I1wx_step, I1wx_step); |
|
|
|
|
const float div_p2 = divergence(p21, p22, y, x, I1wx_step, I1wx_step); |
|
|
|
|
float div_p1 = divergence(p11, p12, y, x, I1wx_step, I1wx_step); |
|
|
|
|
float div_p2 = divergence(p21, p22, y, x, I1wx_step, I1wx_step); |
|
|
|
|
|
|
|
|
|
// estimate the values of the optical flow (u1, u2) |
|
|
|
|
|
|
|
|
|
const float u1NewVal = v1 + theta * div_p1; |
|
|
|
|
const float u2NewVal = v2 + theta * div_p2; |
|
|
|
|
float u1NewVal = v1 + theta * div_p1; |
|
|
|
|
float u2NewVal = v2 + theta * div_p2; |
|
|
|
|
|
|
|
|
|
u1[(y + u1_offset_y) * u1_step + x + u1_offset_x] = u1NewVal; |
|
|
|
|
u2[(y + u2_offset_y) * u2_step + x + u2_offset_x] = u2NewVal; |
|
|
|
|
|
|
|
|
|
if(calc_error) |
|
|
|
|
{ |
|
|
|
|
const float n1 = (u1OldVal - u1NewVal) * (u1OldVal - u1NewVal); |
|
|
|
|
const float n2 = (u2OldVal - u2NewVal) * (u2OldVal - u2NewVal); |
|
|
|
|
float n1 = (u1OldVal - u1NewVal) * (u1OldVal - u1NewVal); |
|
|
|
|
float n2 = (u2OldVal - u2NewVal) * (u2OldVal - u2NewVal); |
|
|
|
|
error[y * I1wx_step + x] = n1 + n2; |
|
|
|
|
} |
|
|
|
|
} |
|
|
|
|
|
Roman Donchenko
11 years ago
committed by
OpenCV Buildbot