SIFT特征点匹配中KD-tree与Ransac算法的使用
轉自：http://blog.csdn.net/ijuliet/article/details/4471311
Step1:BBF算法，在KD-tree上找KNN。第一步做匹配咯~
1.什么是KD-tree（fromwiki）
K-Dimension tree，實際上是一棵平衡二叉樹。
一般的KD-tree構造過程：
functionkdtree (list of points pointList, int depth)
{
ifpointList is empty
returnnil;
else {
// Select axis based on depth so thataxis cycles through all valid values
varint axis := depth mod k;
// Sort point list and choose medianas pivot element
selectmedian by axis from pointList;
// Create node and construct subtrees
vartree_node node;
node.location:= median;
node.leftChild:= kdtree(points in pointList before median, depth+1);
node.rightChild:= kdtree(points in pointList after median, depth+1);
returnnode;
}
}
2.BBF算法，在KD-tree上找KNN ( K-nearest neighbor)
BBF(BestBin First)算法，借助優先隊列（這里用最小堆）實現。從根開始，在KD-tree上找路子的時候，錯過的點先塞到優先隊列里，自己先一個勁兒掃到leaf；然后再從隊列里取出目前key值最小的（這里是是ki維上的距離最小者），重復上述過程，一個勁兒掃到leaf；直到隊列找空了，或者已經重復了200遍了停止。
Step1:將img2的features建KD-tree; kd_root = kdtree_build( feat2,n2 );。在這里，ki是選取均方差最大的那個維度，kv是各特征點在那個維度上的median值，features是你率領的整個兒子孫子特征大軍，n是你兒子孫子個數。
?
| /** a? node in a k-d tree */ struct kd_node{ ???? int ki; /**<? partition key index */ ???? double kv; /**<? partition key value */ ???? int leaf; /**<? 1 if node is a leaf, 0 otherwise */ ???? struct feature* features; /**< features at this node */ ???? int n; /**<? number of features */ ???? struct kd_node* kd_left; /**< left child */ ???? struct kd_node* kd_right; /**< right child */ }; |
Step2: 將img1的每個feat到KD-tree里找k個最近鄰，這里k=2。
k= kdtree_bbf_knn( kd_root, feat, 2, &nbrs, KDTREE_BBF_MAX_NN_CHKS );
?
| ???? min_pq = minpq_init(); ???? minpq_insert( min_pq, kd_root, 0 ); ???? while( min_pq->n > 0 && t < max_nn_chks ) //隊列里有東西就繼續搜，同時控制在t<200（即200步內） ???? { ???????? expl = (struct kd_node*)minpq_extract_min(? min_pq ); //取出最小的，front & pop ???????? expl = explore_to_leaf( expl, feat,? min_pq ); //從該點開始，explore到leaf，路過的“有意義的點”就塞到最小隊列min_pq中。 ???????? for( i =? 0; i < expl->n; i++ ) // ???????? { ????????????? tree_feat =? &expl->features[i]; ????????????? bbf_data->old_data =? tree_feat->feature_data; ????????????? bbf_data->d =? descr_dist_sq(feat, tree_feat); //兩feat均方差 ????????????? tree_feat->feature_data =? bbf_data; ????????????? n += insert_into_nbr_array(? tree_feat, _nbrs, n, k ); //按從小到大塞到neighbor數組里，到時候取前k個就是 KNN 咯~ n 每次加1或0，表示目前已有的元素個數 ???????? } ???????? t++; ???? } |
對“有意義的點”的解釋：
?
| struct kd_node* explore_to_leaf( struct? kd_node* kd_node, struct feature* feat, ???????????????????????????????????? struct? min_pq* min_pq )//expl, feat, min_pq { ???? struct kd_node* unexpl, * expl = kd_node; ???? double kv; ???? int ki; ???? while( expl && ! expl->leaf ) ???? { ???????? ki = expl->ki; ???????? kv = expl->kv; ???????? if(? feat->descr[ki] <= kv ) { ????????????? unexpl = expl->kd_right; ????????????? expl = expl->kd_left; //走左邊，右邊點將被記下來 ???????? } ???????? else{ ????????????? unexpl = expl->kd_left; ????????????? expl = expl->kd_right; //走右邊，左邊點將被記下來 ???????? } ???????? minpq_insert( min_pq, unexpl, ABS( kv? - feat->descr[ki] ) ) ;//將這些點插入進來，key鍵值為|kv? - feat->descr[ki]| 即第ki維上的差值 ???? } ???? return expl; } |
?????? Step3: 如果k近鄰找到了(k=2)，那么判斷是否能作為有效特征，d0/d1<0.49就算是咯~
?
| ????????????? d0 = descr_dist_sq( feat,? nbrs[0] );//計算兩特征間squared Euclidian distance ????????????? d1 = descr_dist_sq( feat,? nbrs[1] ); ????????????? if( d0? < d1 * NN_SQ_DIST_RATIO_THR )//如果d0/d1小于閾值0.49 ????????????? { ?????????????????? pt1 = cvPoint( cvRound(? feat->x ), cvRound( feat->y ) ); ?????????????????? pt2 = cvPoint( cvRound(? nbrs[0]->x ), cvRound( nbrs[0]->y ) ); ?????????????????? pt2.y += img1->height; ???????? ???????? cvLine(? stacked, pt1, pt2, CV_RGB(255,0,255), 1, 8, 0 );//畫線 ?????????????????? m++;//matches個數 ?????????????????? feat1[i].fwd_match =? nbrs[0]; ????????????? } |
Step2:通過RANSAC算法來消除錯配，什么是RANSAC先？
1.RANSAC(Random Sample Consensus, 隨機抽樣一致)(from wiki)
該算法做什么呢？呵呵，用一堆數據去搞定一個待定模型，這里所謂的搞定就是一反復測試、迭代的過程，找出一個error最小的模型及其對應的同盟軍（consensusset）。用在我們的SIFT特征匹配里，就是說找一個變換矩陣出來，使得盡量多的特征點間都符合這個變換關系。
算法思想：
input:
data - a set of observations
model - a model that can be fitted todata
n - the minimum number of datarequired to fit the model
k - the maximum number of iterationsallowed in the algorithm
t - a threshold value for determiningwhen a datum fits a model
d - the number of close data valuesrequired to assert that a model fits well to data
output:
best_model - model parameters whichbest fit the data (or nil if no good model is found)
best_consensus_set - data point fromwhich this model has been estimated
best_error - the error of this modelrelative to the data
iterations:= 0
best_model:= nil
best_consensus_set:= nil
best_error:= infinity
whileiterations < k //進行K次迭代
maybe_inliers:= n randomly selected values from data
maybe_model:= model parameters fitted to maybe_inliers
consensus_set:= maybe_inliers
forevery point in data not in maybe_inliers
ifpoint fits maybe_model with an error smaller than t //錯誤小于閾值t
addpoint to consensus_set //成為同盟，加入consensus set
if thenumber of elements in consensus_set is > d //同盟軍已經大于d個人，夠了
(thisimplies that we may have found a good model,
nowtest how good it is)
better_model:= model parameters fitted to all points in consensus_set
this_error:= a measure of how well better_model fits these points
ifthis_error < best_error
(wehave found a model which is better than any of the previous ones,
keepit until a better one is found)
best_model:= better_model
best_consensus_set:= consensus_set
best_error:= this_error
incrementiterations
returnbest_model, best_consensus_set, best_error
2.RANSAC去除錯配：
H= ransac_xform( feat1, n1, FEATURE_FWD_MATCH, lsq_homog, 4,0.01,homog_xfer_err, 3.0, NULL, NULL );
?
| ???? nm = get_matched_features( features, n,? mtype, &matched ); ???? /*? initialize random number generator */ ???? rng = gsl_rng_alloc( gsl_rng_mt19937 ); ???? gsl_rng_set( rng, time(NULL) ); ???? in_min = calc_min_inliers( nm, m,? RANSAC_PROB_BAD_SUPP, p_badxform ); //符合這一要求的內點至少得有多少個 ???? p = pow( 1.0 - pow( in_frac, m ), k ); ???? i = 0; ???? while( p > p_badxform )//p>0.01 ???? { ???????? sample = draw_ransac_sample( matched,? nm, m, rng ); ???????? extract_corresp_pts( sample, m,? mtype, &pts, &mpts ); ???????? M = xform_fn( pts, mpts, m ); ???????? if( ! M? ) ????????????? goto? iteration_end; ???????? in = find_consensus( matched, nm,? mtype, M, err_fn, err_tol, &consensus); ???????? if( in? > in_max )? { ????????????? if(? consensus_max ) ?????????????????? free( consensus_max ); ????????????? consensus_max = consensus; ????????????? in_max = in; ????????????? in_frac = (double)in_max? / nm; ???????? } ???????? else ????????????? free( consensus ); ???????? cvReleaseMat( &M ); iteration_end: ???????? release_mem( pts, mpts, sample ); ???????? p = pow( 1.0 - pow( in_frac, m ), ++k? ); ???? } ???? /*? calculate final transform based on best consensus set */ ???? if( in_max >= in_min ) ???? { ???????? extract_corresp_pts( consensus_max,? in_max, mtype, &pts, &mpts ); ???????? M = xform_fn( pts, mpts, in_max ); ???????? in = find_consensus( matched, nm,? mtype, M, err_fn, err_tol, &consensus); ???????? cvReleaseMat( &M ); ???????? release_mem( pts, mpts, consensus_max? ); ???????? extract_corresp_pts( consensus, in,? mtype, &pts, &mpts ); ???????? M = xform_fn( pts, mpts, in );?????? |
思考中的一些問題：
features間的對應關系，記錄在features->fwd_match里（matching feature from forward
imge）。
1.數據是nm個特征點間的對應關系，由它們產生一個3*3變換矩陣（xform_fn= hsq_homog函數，此要>=4對的對應才可能計算出來咯~），此乃模型model。
2.然后開始找同盟軍（find_consensus函數），判斷除了sample的其它對應關系是否滿足這個模型（err_fn= homog_xfer_err函數，<=err_tol就OK~），滿足則留下。
3.一旦大于當前的in_max，那么該模型就升級為目前最牛的模型。（最最原始的RANSAC是按錯誤率最小走的，我們這會兒已經保證了錯誤率在err_tol范圍內，按符合要求的對應數最大走，盡量多的特征能匹配地上）
4.重復以上3步，直到(1-wm)k <=p_badxform (即0.01)，模型就算找定~
5.最后再把模型和同盟軍定一下，齊活兒~
聲明：以上代碼參考Rob Hess的SIFT實現。
?
其它參考文獻：
1、http://www.cnblogs.com/slysky/archive/2011/11/08/2241247.html
2、http://en.wikipedia.org/wiki/Kd_tree
3、http://www.cnblogs.com/tjulxh/archive/2011/12/31/2308921.html
4、http://grunt1223.iteye.com/blog/961063
?
總結
以上是生活随笔為你收集整理的SIFT特征点匹配中KD-tree与Ransac算法的使用的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： Windows7 libsvm库中gri
- 下一篇： 在Ubuntu11.10中安装配置Ope