xx = 0.0;
   for (i=0;i<64;i++)
   {
      lag[0][i] = (float) ((xx)*(xx-1.0)*(xx-2.0) / -6.0);
      lag[1][i] = (float) ((xx+1.0)*(xx-1.0)*(xx-2.0) / 2.0);
      lag[2][i] = (float) ((xx+1.0)*(xx)*(xx-2.0) / -2.0);
      lag[3][i] = (float) ((xx)*(xx+1.0)*(xx-1.0) / 6.0);
      xx+= 1.0/64.0;
   }
   for (i=0;i<IMAGE_HEIGHT;i++)
   {
      for(j=0;j<IMAGE_WIDTH;j++)
      {
         out_y = (float) i + (taskid * IMAGE_SLICE);
         out_x = (float) j;

         /* Calculate coordinates of corresponding input pixel */
         in_x = (cos(RAD_PER_DEGREE * ANGLE) * (out_x - 256.0)) +
                (sin(RAD_PER_DEGREE * ANGLE) * (256.0 - out_y)) +
                256.0;
         in_y = (sin(RAD_PER_DEGREE * ANGLE) * (out_x - 256.0)) +
                (cos(RAD_PER_DEGREE * ANGLE) * (out_y - 256.0)) +
                256.0;
         row = (int) in_y;
         col = (int) in_x;

         /* Perform Lagrange interpolation (a 4 by 4 window surrounding */
         /* the source input pixel is used to determine output pixel value */
         if ((row >= 1) && (col >= 1) &&
             (row < IMAGE_HEIGHT - 3) && (col < IMAGE_WIDTH - 3))
         {
            yfrac = (int)((in_y - (float)(row)) * 64.0);
            xfrac = (int)((in_x - (float)(col)) * 64.0);

            sum0= (lag[0][xfrac] * input_image[row-1][col-1]) +
                  (lag[1][xfrac] * input_image[row-1][col]) +
                  (lag[2][xfrac] * input_image[row-1][col+1]) +
                  (lag[3][xfrac] * input_image[row-1][col+2]);
            sum1= (lag[0][xfrac] * input_image[row][col-1]) +
                  (lag[1][xfrac] * input_image[row][col]) +
                  (lag[2][xfrac] * input_image[row][col+1]) +
                  (lag[3][xfrac] * input_image[row][col+2]);
            sum2= (lag[0][xfrac] * input_image[row+1][col-1]) +
                  (lag[1][xfrac] * input_image[row+1][col]) +
                  (lag[2][xfrac] * input_image[row+1][col+1]) +
                  (lag[3][xfrac] * input_image[row+1][col+2]);
            sum3= (lag[0][xfrac] * input_image[row+2][col-1]) +
                  (lag[1][xfrac] * input_image[row+2][col]) +
                  (lag[2][xfrac] * input_image[row+2][col+1]) +
                  (lag[3][xfrac] * input_image[row+2][col+2]);

            output_image[i][j] =
                    (unsigned char)((lag[0][yfrac] * sum0) +
                                    (lag[1][yfrac] * sum1) +
                                    (lag[2][yfrac] * sum2) +
                                    (lag[3][yfrac] * sum3));
         } else
            output_image[i][j] = 0;
      }
   }