# 【OpenCV】SIFT原理与源码分析：方向赋值

### 《SIFT原理与源码分析》系列文章索引：http://blog.csdn.net/xiaowei_cqu/article/details/8069548

// 计算梯度直方图
float omax = calcOrientationHist(gauss_pyr[o*(nOctaveLayers+3) + layer],
Point(c1, r1),
SIFT_ORI_SIG_FCTR * scl_octv,
hist, n);

# 源码

// Computes a gradient orientation histogram at a specified pixel
// 计算特定点的梯度方向直方图
static float calcOrientationHist( const Mat& img, Point pt, int radius,
float sigma, float* hist, int n )
{
//len：2r+1也就是以r为半径的圆（正方形）像素个数

float expf_scale = -1.f/(2.f * sigma * sigma);
AutoBuffer<float> buf(len*4 + n+4);
float *X = buf, *Y = X + len, *Mag = X, *Ori = Y + len, *W = Ori + len;
float* temphist = W + len + 2;

for( i = 0; i < n; i++ )
temphist[i] = 0.f;

// 图像梯度直方图统计的像素范围
for( i = -radius, k = 0; i <= radius; i++ )
{
int y = pt.y + i;
if( y <= 0 || y >= img.rows - 1 )
continue;
{
int x = pt.x + j;
if( x <= 0 || x >= img.cols - 1 )
continue;

float dx = (float)(img.at<short>(y, x+1) - img.at<short>(y, x-1));
float dy = (float)(img.at<short>(y-1, x) - img.at<short>(y+1, x));

X[k] = dx; Y[k] = dy; W[k] = (i*i + j*j)*expf_scale;
k++;
}
}

len = k;

// compute gradient values, orientations and the weights over the pixel neighborhood
exp(W, W, len);
fastAtan2(Y, X, Ori, len, true);
magnitude(X, Y, Mag, len);

// 计算直方图的每个bin
for( k = 0; k < len; k++ )
{
int bin = cvRound((n/360.f)*Ori[k]);
if( bin >= n )
bin -= n;
if( bin < 0 )
bin += n;
temphist[bin] += W[k]*Mag[k];
}

// smooth the histogram
// 高斯平滑
temphist[-1] = temphist[n-1];
temphist[-2] = temphist[n-2];
temphist[n] = temphist[0];
temphist[n+1] = temphist[1];
for( i = 0; i < n; i++ )
{
hist[i] = (temphist[i-2] + temphist[i+2])*(1.f/16.f) +
(temphist[i-1] + temphist[i+1])*(4.f/16.f) +
temphist[i]*(6.f/16.f);
}

// 得到主方向
float maxval = hist[0];
for( i = 1; i < n; i++ )
maxval = std::max(maxval, hist[i]);

return maxval;
}