params.upright = 0;

            float descriptor_dir;

            if (params->upright == 0) {

                icvResizeHaarPattern( dx_s, dx_t, NX, 4, grad_wav_size, sum->cols );

                icvResizeHaarPattern( dy_s, dy_t, NY, 4, grad_wav_size, sum->cols );

                for( kk = 0, nangle = 0; kk < nOriSamples; kk++ )

                {

                    const int* ptr;

                    float vx, vy;

                    x = cvRound( center.x + apt[kk].x*s - (float)(grad_wav_size-1)/2 );

                    y = cvRound( center.y + apt[kk].y*s - (float)(grad_wav_size-1)/2 );

                    if( (unsigned)y >= (unsigned)(sum->rows - grad_wav_size) ||

                        (unsigned)x >= (unsigned)(sum->cols - grad_wav_size) )

                        continue;

                    ptr = sum_ptr + x + y*sum_cols;

                    vx = icvCalcHaarPattern( ptr, dx_t, 2 );

                    vy = icvCalcHaarPattern( ptr, dy_t, 2 );

                    X[nangle] = vx*aptw[kk]; Y[nangle] = vy*aptw[kk];

                    nangle++;

                }

                if ( nangle == 0 )

                {

                    /* No gradient could be sampled because the keypoint is too

                    * near too one or more of the sides of the image. As we

                    * therefore cannot find a dominant direction, we skip this

                    * keypoint and mark it for later deletion from the sequence. */

                    kp->size = -1;

                }

                matX.cols = matY.cols = _angle.cols = nangle;

                cvCartToPolar( &matX, &matY, 0, &_angle, 1 );

                float bestx = 0, besty = 0, descriptor_mod = 0;

                for( i = 0; i < 360; i += ORI_SEARCH_INC )

                    float sumx = 0, sumy = 0, temp_mod;

                    for( j = 0; j < nangle; j++ )

                        int d = std::abs(cvRound(angle[j]) - i);

                        if( d < ORI_WIN/2 || d > 360-ORI_WIN/2 )

                        {

                            sumx += X[j];

                            sumy += Y[j];

                        }

                    temp_mod = sumx*sumx + sumy*sumy;

                    if( temp_mod > descriptor_mod )

                    {

                        descriptor_mod = temp_mod;

                        bestx = sumx;

                        besty = sumy;

                    }

                descriptor_dir = cvFastArctan( besty, bestx );

                kp->dir = descriptor_dir;

            } else {

                descriptor_dir = 0;

                kp->dir = 0;

            if (params->upright == 0) {

                float sin_dir = sin(descriptor_dir);

                float cos_dir = cos(descriptor_dir) ;

                /* Subpixel interpolation version (slower). Subpixel not required since

                the pixels will all get averaged when we scale down to 20 pixels */

                /*  

                float w[] = { cos_dir, sin_dir, center.x,

                -sin_dir, cos_dir , center.y };

                CvMat W = cvMat(2, 3, CV_32F, w);

                cvGetQuadrangleSubPix( img, &win, &W );

                */

                /* Nearest neighbour version (faster) */

                float win_offset = -(float)(win_size-1)/2;

                float start_x = center.x + win_offset*cos_dir + win_offset*sin_dir;

                float start_y = center.y - win_offset*sin_dir + win_offset*cos_dir;

                uchar* WIN = win.data.ptr;

                for( i = 0; i < win_size; i++, start_x += sin_dir, start_y += cos_dir )

                    float pixel_x = start_x;

                    float pixel_y = start_y;

                    for( j = 0; j < win_size; j++, pixel_x += cos_dir, pixel_y -= sin_dir )

                    {

                        int x = std::min(std::max(cvRound(pixel_x), 0), img->cols-1);

                        int y = std::min(std::max(cvRound(pixel_y), 0), img->rows-1);

                        WIN[i*win_size + j] = img->data.ptr[y*img->step + x];

                    }

            } else {

                /* extract rect - slightly optimized version of the code above

                   TODO: find faster code, as this is simply an extract rect operation, 

                         e.g. by using cvGetSubRect, problem is the border processing */

                float win_offset = -(float)(win_size-1)/2;

                int start_x = cvRound(center.x + win_offset);

                int start_y = cvRound(center.y + win_offset);

                uchar* WIN = win.data.ptr;

                for( i=0; i<win_size; i++, start_y++ )

                {

                    int pixel_x = start_x;

                    for( j=0; j<win_size; j++, pixel_x++ )

                    {

                        x = MAX( pixel_x, 0 );

                        y = MAX( start_y, 0 );

                        x = MIN( x, img->cols-1 );

                        y = MIN( y, img->rows-1 );

                        WIN[i*win_size + j] = img->data.ptr[y*img->step+x];

                    }

                }               

	    cv::SURFInvoker invoker(&params, keypoints, descriptors, img, sum);

	    invoker(cv::BlockedRange(0, N));

    upright = 0;

SURF::SURF(double _threshold, int _nOctaves, int _nOctaveLayers, bool _extended, bool _upright)

    upright = _upright;