    const bool withinImageBounds = xOut < output_width; // check whether x exceed the boundary
    const bool withinGridBounds = blockIdx.x*blockDim.y + threadIdx.x / 2 < output_width; // check whether
 
+   if(!withinImageBounds || !withinGridBounds) return;
+
    const int width = inputImages_width;
    const int height = inputImages_height;
 
-   float yf,xf;
-
-   __shared__ float gridData[32];
-   __shared__ float gridData_r[32];
-   __shared__ float gridData_b[32];
-   __shared__ float gridData_br[32];
-
-   if (threadIdx.y==0 && withinGridBounds)
-   {
-     # if __CUDA_ARCH__>=200
-       printf("%d \n", grids_strideWidth);
-     #endif
-      gridData[threadIdx.x] = grids_data[b*grids_strideBatch + yOut*grids_strideHeight + xOut*grids_strideWidth + threadIdx.x];
+   int xOut_l = xOut > 0 ? (xOut - 1) : xOut;
+   int yOut_t = yOut > 0 ? (yOut - 1) : yOut;
+   int xOut_r = xOut < (output_width - 1) ? (xOut + 1) : xOut;
+   int yOut_b = yOut < (output_height - 1) ? (yOut + 1) : yOut;
+
+   float yf = grids_data[b*grids_strideBatch + yOuy*grids_strideHeight + xOut*grids_strideWidth];
+   float xf = grids_data[b*grids_strideBatch + yOut*grids_strideHeight + xOut*grids_strideWidth + 1];
+
+   // get coordinates for four corners of target response map in the source response map
+   float yf_tl = grids_data[b*grids_strideBatch + yOut_t*grids_strideHeight + xOut_l*grids_strideWidth];
+   float xf_tl = grids_data[b*grids_strideBatch + yOut_t*grids_strideHeight + xOut_l*grids_strideWidth + 1];
+
+   float yf_tr = grids_data[b*grids_strideBatch + yOut_t*grids_strideHeight + xOut_r*grids_strideWidth];
+   float xf_tr = grids_data[b*grids_strideBatch + yOut_t*grids_strideHeight + xOut_r*grids_strideWidth + 1];
+
+   float yf_bl = grids_data[b*grids_strideBatch + yOut_b*grids_strideHeight + xOut_l*grids_strideWidth];
+   float xf_bl = grids_data[b*grids_strideBatch + yOut_b*grids_strideHeight + xOut_l*grids_strideWidth + 1];
+
+   float yf_br = grids_data[b*grids_strideBatch + yOut_b*grids_strideHeight + xOut_r*grids_strideWidth];
+   float xf_br = grids_data[b*grids_strideBatch + yOut_b*grids_strideHeight + xOut_r*grids_strideWidth + 1];
+
+   // compute the bottom left coordinates of source response map
+   float xs_tl = min(min(min(xf_tl, xf_tr), xf_bl), xf_br);
+   float ys_tl = min(min(min(yf_tl, yf_tr), yf_bl), yf_br);
+
+   float xs_br = max(max(max(xf_tl, xf_tr), xf_bl), xf_br);
+   float ys_br = max(max(max(yf_tl, yf_tr), yf_bl), yf_br);
+
+   int yi_t, xi_l;
+   float yWeight_t, xWeight_l;
+
+   int yi_b, xi_r;
+   float yWeight_b, xWeight_r;
+
+   // compute the top left, top right, bottom left, bottom right corner for current grid
+   // compute the nearest bottom right coordiate in source map of top left grid coordiate
+   getBottomRight(xs_tl, inputImages_width, xi_l, xWeight_l);
+   getBottomRight(ys_tl, inputImages_height, yi_t, yWeight_t);
+   // compute the nearest top left coordiate in source map of bottom right grid coordiate
+   getTopLeft(xs_br, inputImages_width, xi_r, xWeight_r);
+   getTopLeft(ys_br, inputImages_height, yi_b, yWeight_b);
+
+   // compute average mask value in region [topleft, bottom right]
+   float weights[32];
+   float m = 0;
+   float weight_sum = 0;
+   int id_point = 0;
+   for (int y = yi_t; y <= yi_b; ++y) {
+     if (!between(y, 0, height-1)) continue;
+     for (int x = xi_l; x <= xi_r; +=x) {
+       if (!between(x, 0, width-1)) continue;
+       float weight = __expf((x - xf) * (x - xf) + (y - yf) * (y - yf));
+       weight_sum += weight;
+       weights[id_point] = weight;
+       int address = masks_strideBatch * b + masks_strideHeight * y + masks_strideWidth * x;
+       m += weight * masks_data[address];
+       ++id_point;
+     }
    }
-   __syncthreads(); // synchronize to ensure all threads begin to get data after storing
-
-   if(!withinImageBounds) return;
-   yf = gridData[threadIdx.y*2];
-   xf = gridData[threadIdx.y*2+1];
-
-   int yInTopLeft, xInTopLeft;
-   float yWeightTopLeft, xWeightTopLeft;
-
-   // compute the top left, top right, bottom left, bottom right corner for current grid
-   getTopLeft(xf, inputImages_width, xInTopLeft, xWeightTopLeft);
-   getTopLeft(yf, inputImages_height, yInTopLeft, yWeightTopLeft);
-
-   xWeightTopLeft = 0.5;
-   yWeightTopLeft = 0.5;
-
-   // compute the address (location) for output
-   const int outAddress = output_strideBatch * b + output_strideHeight * yOut + output_strideWidth * xOut;
-
-
-   const int inTopLeftAddress = inputImages_strideBatch * b + inputImages_strideHeight * yInTopLeft + inputImages_strideWidth * xInTopLeft;
-   const int inTopRightAddress = inTopLeftAddress + inputImages_strideWidth;
-   const int inBottomLeftAddress = inTopLeftAddress + inputImages_strideHeight;
-   const int inBottomRightAddress = inBottomLeftAddress + inputImages_strideWidth;
-
-   const int inTopLeftMaskAddress = masks_strideBatch * b + masks_strideHeight * yInTopLeft + masks_strideWidth * xInTopLeft;
-   const int inTopRightMaskAddress = inTopLeftMaskAddress + masks_strideWidth;
-   const int inBottomLeftMaskAddress = inTopLeftMaskAddress + masks_strideHeight;
-   const int inBottomRightMaskAddress = inBottomLeftMaskAddress + masks_strideWidth;
+   m /= weight_sum;
 
    float v=0;
    float inTopLeft=0;

                                         float* canvas_data, int canvas_strideBatch, int canvas_strideYX, int canvas_strideHeight, int canvas_strideWidth,
                                         float* gradCanvas_data, int gradCanvas_strideBatch, int gradCanvas_strideYX, int gradCanvas_strideHeight, int gradCanvas_strideWidth,
                                         float* gradOutput_data, int gradOutput_strideBatch, int gradOutput_strideChannels, int gradOutput_strideHeight, int gradOutput_strideWidth,
                                         int inputImages_channels, int inputImages_height, int inputImages_width, int gradOutput_width)
{
   // each (32,16) block 16 output pixels (for coalescing the grid read)
   // x,y = coordinates

                                         float* masks_data, int masks_strideBatch, int masks_strideYX, int masks_strideHeight, int masks_strideWidth,
                                         float* canvas_data, int canvas_strideBatch, int canvas_strideYX, int canvas_strideHeight, int canvas_strideWidth,
                                         float* output_data, int output_strideBatch, int output_strideChannels, int output_strideHeight, int output_strideWidth,
                                         int inputImages_channels, int inputImages_height, int inputImages_width, int output_height, int output_width)
{
   // each (32,16) block 16 output pixels (for coalescing the grid read)
   // x,y = coordinates (xOut = blockIdx.x*16+blockDim.y+threadIdx.y)
