http://blog.sina.com.cn/s/blog_916b71bb0100wbdh.html 算法的中文简介
http://cvlab.epfl.ch/research/detect/brief/ 算法的homepage
这个算法比较简单,
1.利用已有的算法(比如说fast-9算法)找到角点,在作者的示例程序中使用surf算法来搜寻角点。
2.在找到的角点的周围规定一个S*S的patch,对patch所在的S*S区域进行高斯模糊。
3.论文中介绍了5种算法对patch中的找到nd个(x,y)个像素对(其中x and y代表一个像素点的位置)进行取样。其中第二种算法,作者认为是最好的,并且在给中的代码中仅给出了这种算法。
翻译:(X,Y)均服从均值为0,方差为1/25*s*s的正太分布,实验是在各向同性的高斯分布下抽样的。我们发现:在方差为1/25*s*s参数下可以取得最好的识别率。
4.抽样完毕后根据:
计算出每个tao,然后根据:
计算出特征描述子。特征描述子是一个binary序列。
5.在进行匹配的时候使用hamming距离计算相似度。
论文代码分析:
TestSampler< TestLocT >::sample(tests_, DESC_LEN, INIT_PATCH_SZ, 1);
下面这个函数产生高斯取样。
double randNormal(const double std, const double mu)
{
double x1, x2, w, y1;
static double y2;
static bool use_last = false;
if (use_last) /* use value from previous call */
{
y1 = y2;
use_last = false;
}
else
{
do {
x1 = 2.0 * cvRandReal(&LocalStorage::getInst().rng_obj) - 1.0;
x2 = 2.0 * cvRandReal(&LocalStorage::getInst().rng_obj) - 1.0;
w = x1 * x1 + x2 * x2;
} while ( w >= 1.0 );
w = sqrt( (-2.0 * log( w ) ) / w );
y1 = x1 * w;
y2 = x2 * w;
use_last = true;
}
return mu + y1 * std;
}
} // namespace utils
下面的代码产生brief描述子。把描述子放在res中。
template< typename KPT_T >
void BRIEF::getBRIEF(const std::vector< KPT_T >& kpts,
std::vector< std::bitset<DESC_LEN> >& res,
const IplImage* precomp_int_img)
{
assert(kpts[0].image->nChannels == 1);
// Make sure all points are from the same image
const IplImage* img = kpts[0].image;
assert(img);
#if DEBUG
for (int i=1; i<(int)kpts.size(); ++i)
assert(kpts[i].image == img);
#endif
// Pre-compute integral image
// #if USE_INTEGRAL_IMAGE
// if (!precomp_int_img && cur_int_img_.src != img) // Update needed for our integral image?
// {
// if (!cur_int_img_.src)
// cur_int_img_.data = new float[img->width*img->height];
// else if (cur_int_img_.src->width*cur_int_img_.src->height < img->width*img->height)
// {
// delete [] cur_int_img_.data;
// cur_int_img_.data = new float[img->width*img->height];
// }
//
// cur_int_img_.src = img;
// utils::computeIntegralImage(img, cur_int_img_.data);
// }
// #endif
// Recompute bounding box, if needed
if (act_patch_sz_ <= 0)
recomputeActPatchSz();
// Compute all descriptors
std::bitset<DESC_LEN> some_desc;
res.reserve(kpts.size());
for (int i=0; i<(int)kpts.size(); ++i)
{
const KPT_T &kp = kpts[i];
#if DEBUG
if (!checkAllTestInsideImage(kp)) {
some_desc.reset();
res.push_back(some_desc);
continue;
}
#endif
// Rotate and scale tests, if information in keypoints valid
// Note: Mind the INIT_PATCH_SZ and BOX_SZ params or this will segfault
#if SOS_SCALE != SOS_DO_NOT_INTERPRET && SOS_ORI != SOS_DO_NOT_INTERPRET
rotateAndScaleTests(interpretOri(kp.ori), interpretScale(kp.scale), true);
#define KERNEL_SZ act_kernel_sz_
#elif SOS_SCALE != SOS_DO_NOT_INTERPRET
scaleTests(interpretScale(kp.scale), true);
#define KERNEL_SZ act_kernel_sz_
#elif SOS_ORI != SOS_DO_NOT_INTERPRET
rotateTests(interpretOri(kp.ori), true);
#define KERNEL_SZ INIT_KERNEL_SZ
#else
#define KERNEL_SZ INIT_KERNEL_SZ
#endif
// Smooth patch, compute descriptor
#if USE_INTEGRAL_IMAGE
assert(!precomp_int_img || precomp_int_img->depth == IPL_DEPTH_64F);
const float* ii_data = cur_int_img_.data;
const int ii_step = cur_int_img_.src->width;
#if SMOOTH_ENTIRE_PATCH // smooth entire patch, then apply tests
#if DEBUG
if (!checkSmoothAreaInsideIntImage(kp)) {
some_desc.reset();
res.push_back(some_desc);
continue;
}
#endif
// Compute UNNORMALIZED smoothed patch around (kp.x, kp.y), ok since only
// relative values are of importance
// Using "area = B + D - A - C" where arrangement is B A
// 对patch区域进行高斯模糊。 C D
float *p_smooth = smooth_;
for (int v=0; v<act_patch_sz_; ++v)
{
// cix: smoothing window central pixel offset
int cix = (kp.y-act_patch_sz_/2+v)*ii_step + kp.x-act_patch_sz_/2;
const float *p_a = &ii_data[cix - KERNEL_SZ/2*ii_step + KERNEL_SZ/2];
const float *p_b = &ii_data[cix - KERNEL_SZ/2*ii_step - KERNEL_SZ/2];
const float *p_c = &ii_data[cix + KERNEL_SZ/2*ii_step - KERNEL_SZ/2];
const float *p_d = &ii_data[cix + KERNEL_SZ/2*ii_step + KERNEL_SZ/2];
for (int u=0; u<act_patch_sz_; ++u, ++p_smooth, ++p_a, ++p_b, ++p_c, ++p_d)
*p_smooth = ((*p_b + *p_d - *p_a - *p_c));
}
// Ready for tests
res.push_back(some_desc);
applyTests(smooth_, res[res.size()-1], tests_);
#else // smooth around test locations only
#if DEBUG
if (!checkSmoothAreaForTestInsideIntImage(kp)) {
some_desc.reset();
res.push_back(some_desc);
continue;
}
#endif
res.push_back(some_desc);
std::bitset<DESC_LEN> &bits = res[res.size()-1];
const TestLocT *pxy = tests_;//获得高斯采样的的数组,这样可以获得高的速度,一次计算后,以后查表。
for (int j=0; j<DESC_LEN; ++j, pxy+=4)
{
// cix: central pixel offset
const int cix_1 = (kp.y + pxy[1])*ii_step + kp.x+pxy[0];
const int cix_2 = (kp.y + pxy[3])*ii_step + kp.x+pxy[2];
// Compute tests按照公式(1)计算。
bits[j] = ( ii_data[cix_1 - KERNEL_SZ/2*ii_step - KERNEL_SZ/2] +
ii_data[cix_1 + KERNEL_SZ/2*ii_step + KERNEL_SZ/2] -
ii_data[cix_1 - KERNEL_SZ/2*ii_step + KERNEL_SZ/2] -
ii_data[cix_1 + KERNEL_SZ/2*ii_step - KERNEL_SZ/2] )
>
( ii_data[cix_2 - KERNEL_SZ/2*ii_step - KERNEL_SZ/2] +
ii_data[cix_2 + KERNEL_SZ/2*ii_step + KERNEL_SZ/2] -
ii_data[cix_2 - KERNEL_SZ/2*ii_step + KERNEL_SZ/2] -
ii_data[cix_2 + KERNEL_SZ/2*ii_step - KERNEL_SZ/2] );
}
#endif
#else // not using integral image but cvSmooth instead
#ifdef DEBUG
if (!checkPatchInsideImage(kp)) {
some_desc.reset();
res.push_back(some_desc);
continue;
}
#endif
// Need to copy patch data before smoothing
assert(img->depth == (int)IPL_DEPTH_8U || img->depth == (int)IPL_DEPTH_8S); // or memcpy will fail
const int step = img->widthStep;
uchar* start = (uchar*)&img->imageData[(kp.y-act_patch_sz_/2)*step + kp.x-act_patch_sz_/2];
for (int j=0; j<act_patch_sz_; ++j, start += step)
memcpy(patch_8u_->imageData + j*act_patch_sz_, start, act_patch_sz_);
cvSmooth(patch_8u_, patch_8u_2_, CV_BLUR, KERNEL_SZ, KERNEL_SZ);
// Ready for tests
res.push_back(some_desc);
applyTests((uchar*)patch_8u_2_->imageData, res[res.size()-1], tests_);
#endif
}
#undef KERNEL_SZ
}
下面的代码对计算出的描述子进行match,以下代码本人没有看懂。
template< typename KPT_T, int DESC_LEN >
struct BRIEFMatcher
{
// Matching with L/R check
void matchLeftRight(const std::vector< std::bitset< DESC_LEN > > feat_left,
const std::vector< std::bitset< DESC_LEN > > feat_right,
const std::vector< KPT_T > kpts_left,
const std::vector< KPT_T > kpts_right,
std::insert_iterator< std::vector< KPT_T > > lmatch_ins,
std::insert_iterator< std::vector< KPT_T > > rmatch_ins)
{
int *dist = new int[MAX(feat_left.size(), feat_right.size())];
int *left_nn = new int[feat_left.size()];
int *right_nn = new int[feat_right.size()];
// Matching with left-right check
int min_ix;
for (int i=0; i<(int)feat_left.size(); ++i)
{
for (int k=0; k<(int)feat_right.size(); ++k)
dist[k] = (int)(feat_left[i] ^ feat_right[k]).count();
utils::minVect(dist, (int)feat_right.size(), &min_ix);
left_nn[i] = min_ix;
}
for (int i=0; i<(int)feat_right.size(); ++i)
{
for (int k=0; k<(int)feat_left.size(); ++k)
dist[k] = (int)(feat_right[i] ^ feat_left[k]).count();
utils::minVect(dist, (int)feat_left.size(), &min_ix);
right_nn[i] = min_ix;
}
for (int i=0; i<(int)feat_left.size(); ++i)
if (right_nn[left_nn[i]] == i)
{
*lmatch_ins = kpts_left[i];
*rmatch_ins = kpts_right[left_nn[i]];
++lmatch_ins; ++rmatch_ins;
}
delete [] dist;
delete [] left_nn;
delete [] right_nn;
}
};