ORB Test Hanson

之前介绍了ORB，一种具备旋转不变形的局部特征描述子。OpenCV2.3中提供了实现，但是缺少使用例程。下面是一个简单的样例程序。

随便拍了两张图片作为测试图像。

下面上下两图分别为模板图像和查询图像：


提取左右图特征：

?
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
Mat img1 = imread(image_filename1, 0);
Mat img2 = imread(image_filename2, 0);
//GaussianBlur(img1, img1, Size(5, 5), 0);
//GaussianBlur(img2, img2, Size(5, 5), 0);
  
<A class="st_tag internal_tag" title="标签 orb 下的日志" href="http://www.cvchina.info/tag/orb/" rel=tag>ORB</A> orb1(3000, <A class="st_tag internal_tag" title="标签 orb 下的日志" href="http://www.cvchina.info/tag/orb/" rel=tag>ORB</A>::CommonParams(1.2, 8));
<A class="st_tag internal_tag" title="标签 orb 下的日志" href="http://www.cvchina.info/tag/orb/" rel=tag>ORB</A> orb2(100, <A class="st_tag internal_tag" title="标签 orb 下的日志" href="http://www.cvchina.info/tag/orb/" rel=tag>ORB</A>::CommonParams(1.2, 1));
  
vector keys1, keys2;
Mat descriptors1, descriptors2;
  
<A class="st_tag internal_tag" title="标签 orb 下的日志" href="http://www.cvchina.info/tag/orb/" rel=tag>orb</A>1(img1, Mat(), keys1, descriptors1, false);
printf("tem feat num: %d\n", keys1.size());
  
int64 st, et;
st = cvGetTickCount();
<A class="st_tag internal_tag" title="标签 orb 下的日志" href="http://www.cvchina.info/tag/orb/" rel=tag>orb</A>2(img2, Mat(), keys2, descriptors2, false);
et = cvGetTickCount();
printf("<A class="st_tag internal_tag" title="标签 orb 下的日志" href="http://www.cvchina.info/tag/orb/" rel=tag>orb</A>2 extraction time: %f\n", (et-st)/(double)cvGetTickFrequency()/1000.);
printf("query feat num: %d\n", keys2.size());

注：模板图像在多尺度提取特征，查询图像只在提取原始尺度上的特征。

做穷举式的最近邻检索：

?
1
2
3
4
5
6
7
8
9
10
11
// find matches
vector matches;
  
st = cvGetTickCount();
//for(int i = 0; i < 10; i++){
naive_nn_search2(keys1, descriptors1, keys2, descriptors2, matches);
//}
et = cvGetTickCount();
  
printf("match time: %f\n", (et-st)/(double)cvGetTickFrequency()/1000.);
printf("matchs num: %d\n", matches.size());

hamming距离测算通过查找表实现：

?
1
2
3
4
5
6
7
8
unsigned int hamdist2(unsigned char* a, unsigned char* b, size_t size)
{
HammingLUT lut;
  
unsigned int result;
result = lut((a), (b), size);
return result;
}

绘图：

?
1
2
3
4
5
6
Mat showImg;
drawMatches(img2, keys2, img1, keys1, matches, showImg, CV_RGB(0, 255, 0), CV_RGB(0, 0, 255));
string winName = "Matches";
namedWindow( winName, WINDOW_AUTOSIZE );
imshow( winName, showImg );
waitKey();

估计单应矩阵，计算重投影误差：

?
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
Mat homo;
  
st = cvGetTickCount();
homo = findHomography(pt1, pt2, Mat(), CV_RANSAC, 5);
et = cvGetTickCount();
printf("ransac time: %f\n", (et-st)/(double)cvGetTickFrequency()/1000.);
  
printf("homo\n"
"%f %f %f\n"
"%f %f %f\n"
"%f %f %f\n",
homo.at(0,0), homo.at(0,1), homo.at(0,2),
homo.at(1,0), homo.at(1,1), homo.at(1,2),
homo.at(2,0),homo.at(2,1),homo.at(2,2));
  
vector
reproj;
reproj.resize(pt1.size());
  
perspectiveTransform(pt1, reproj, homo);
  
Mat diff;
diff = Mat(reproj) - Mat(pt2);
  
int inlier = 0;
double err_sum = 0;
for(int i = 0; i < diff.rows; i++){
float* ptr = diff.ptr(i);
float err = ptr[0]*ptr[0] + ptr[1]*ptr[1];
if(err < 25.f){
inlier++;
err_sum += sqrt(err);
}
}
printf("inlier num: %d\n", inlier);
printf("ratio %f\n", inlier / (float)(diff.rows));
printf("mean reprojection error: %f\n", err_sum / inlier);

结果分析：

tem feat num: 743
orb2 extraction time: 1.672435
query feat num: 100
match time: 3.698276
matchs num: 8
ransac time: 143.570586
homo
0.974942 0.410833 4.426035
-0.182418 0.828115 52.742661
0.001191 0.000144 1.000000
inlier num: 8
ratio 1.000000
mean reprojection error: 0.976777

可见最近邻检索是系统的瓶颈，（进行了743*100次hamming距离（32bytes）计算。）一个简单的优化如下，分段计算hamming距离，先计算前16byte的hamming距离，如超过某一阈值，则直接认为非候选，如小于某阈值，则继续进行后一半16bytes的距离计算。（粗略估计可以减少30%+的最近邻查询时间）。更复杂的办法是使用LSH，此处按下不提，有空再续。

完整代码如下：

?
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
#include "opencv2/objdetect/objdetect.hpp"
#include "opencv2/features2d/features2d.hpp"
#include "opencv2/highgui/highgui.hpp"
#include "opencv2/calib3d/calib3d.hpp"
#include "opencv2/imgproc/imgproc_c.h"
#include "opencv2/imgproc/imgproc.hpp"
  
#include
#include
#include
  
using namespace std;
using namespace cv;
  
char* image_filename1 = "apple_vinegar_0.png";
char* image_filename2 = "apple_vinegar_2.png";
  
unsigned int hamdist(unsigned int x, unsigned int y)
{
unsigned int dist = 0, val = x ^ y;
  
// Count the number of set bits
while(val)
{
++dist;
val &= val - 1;
}
  
return dist;
}
  
unsigned int hamdist2(unsigned char* a, unsigned char* b, size_t size)
{
HammingLUT lut;
  
unsigned int result;
result = lut((a), (b), size);
return result;
}
  
void naive_nn_search(vector& keys1, Mat& descp1,
vector& keys2, Mat& descp2,
vector& matches)
{
for( int i = 0; i < (int)keys2.size(); i++){
unsigned int min_dist = INT_MAX;
int min_idx = -1;
unsigned char* query_feat = descp2.ptr(i);
for( int j = 0; j < (int)keys1.size(); j++){
unsigned char* train_feat = descp1.ptr(j);
unsigned int dist =  hamdist2(query_feat, train_feat, 32);
  
if(dist < min_dist){
min_dist = dist;
min_idx = j;
}
}
  
//if(min_dist <= (unsigned int)(second_dist * 0.8)){
if(min_dist <= 50){
matches.push_back(DMatch(i, min_idx, 0, (float)min_dist));
}
}
}
  
void naive_nn_search2(vector& keys1, Mat& descp1,
vector& keys2, Mat& descp2,
vector& matches)
{
for( int i = 0; i < (int)keys2.size(); i++){
unsigned int min_dist = INT_MAX;
unsigned int sec_dist = INT_MAX;
int min_idx = -1, sec_idx = -1;
unsigned char* query_feat = descp2.ptr(i);
for( int j = 0; j < (int)keys1.size(); j++){
unsigned char* train_feat = descp1.ptr(j);
unsigned int dist =  hamdist2(query_feat, train_feat, 32);
  
if(dist < min_dist){
sec_dist = min_dist;
sec_idx = min_idx;
min_dist = dist;
min_idx = j;
}else if(dist < sec_dist){
sec_dist = dist;
sec_idx = j;
}
}
  
if(min_dist <= (unsigned int)(sec_dist * 0.8) && min_dist <=50){
//if(min_dist <= 50){
matches.push_back(DMatch(i, min_idx, 0, (float)min_dist));
}
}
}
  
int main(int argc, char* argv[])
{
Mat img1 = imread(image_filename1, 0);
Mat img2 = imread(image_filename2, 0);
//GaussianBlur(img1, img1, Size(5, 5), 0);
//GaussianBlur(img2, img2, Size(5, 5), 0);
  
ORB orb1(3000, ORB::CommonParams(1.2, 8));
ORB orb2(100, ORB::CommonParams(1.2, 1));
  
vector keys1, keys2;
Mat descriptors1, descriptors2;
  
orb1(img1, Mat(), keys1, descriptors1, false);
printf("tem feat num: %d\n", keys1.size());
  
int64 st, et;
st = cvGetTickCount();
orb2(img2, Mat(), keys2, descriptors2, false);
et = cvGetTickCount();
printf("orb2 extraction time: %f\n", (et-st)/(double)cvGetTickFrequency()/1000.);
printf("query feat num: %d\n", keys2.size());
  
// find matches
vector matches;
  
st = cvGetTickCount();
//for(int i = 0; i < 10; i++){
naive_nn_search2(keys1, descriptors1, keys2, descriptors2, matches);
//}
et = cvGetTickCount();
  
printf("match time: %f\n", (et-st)/(double)cvGetTickFrequency()/1000.);
printf("matchs num: %d\n", matches.size());
  
Mat showImg;
drawMatches(img2, keys2, img1, keys1, matches, showImg, CV_RGB(0, 255, 0), CV_RGB(0, 0, 255));
string winName = "Matches";
namedWindow( winName, WINDOW_AUTOSIZE );
imshow( winName, showImg );
waitKey();
  
vector
pt1;
vector
pt2;
  
for(int i = 0; i < (int)matches.size(); i++){
pt1.push_back(Point2f(keys2[matches[i].queryIdx].pt.x, keys2[matches[i].queryIdx].pt.y));
  
pt2.push_back(Point2f(keys1[matches[i].trainIdx].pt.x, keys1[matches[i].trainIdx].pt.y));
}
  
Mat homo;
  
st = cvGetTickCount();
homo = findHomography(pt1, pt2, Mat(), CV_RANSAC, 5);
et = cvGetTickCount();
printf("ransac time: %f\n", (et-st)/(double)cvGetTickFrequency()/1000.);
  
printf("homo\n"
"%f %f %f\n"
"%f %f %f\n"
"%f %f %f\n",
homo.at(0,0), homo.at(0,1), homo.at(0,2),
homo.at(1,0), homo.at(1,1), homo.at(1,2),
homo.at(2,0),homo.at(2,1),homo.at(2,2));
  
vector
reproj;
reproj.resize(pt1.size());
  
perspectiveTransform(pt1, reproj, homo);
  
Mat diff;
diff = Mat(reproj) - Mat(pt2);
  
int inlier = 0;
double err_sum = 0;
for(int i = 0; i < diff.rows; i++){
float* ptr = diff.ptr(i);
float err = ptr[0]*ptr[0] + ptr[1]*ptr[1];
if(err < 25.f){
inlier++;
err_sum += sqrt(err);
}
}
printf("inlier num: %d\n", inlier);
printf("ratio %f\n", inlier / (float)(diff.rows));
printf("mean reprojection error: %f\n", err_sum / inlier);
  
return 0;
}
 
分类: 练习 标签: hamming, homography, lsh, opencv2.3, orb, 单应矩阵
  1. boil
    2011年9月26日00:38 | #1
    回复 | 引用
     

    sf，不错，学习了～

     
  2. 壮壮
    2011年9月26日07:37 | #2
    回复 | 引用
     

    试了下，能用，不过希望改下源码显示插件，太难懂了，有的地方被过滤掉了，还要自己猜！！！

     
  3. cvchina
    2011年9月26日08:58 | #3
    回复 | 引用
     

    壮壮 :

    试了下，能用，不过希望改下源码显示插件，太难懂了，有的地方被过滤掉了，还要自己猜！！！

    确实。谢谢建议。

     
  4. 老杨
    2011年9月26日10:42 | #4
    回复 | 引用
     

    你可以学习一下我的博客怎么显示代码的么？

     
  5. gxf
    2011年9月26日10:48 | #5
    回复 | 引用
     

    ORB是用什么做检测子的呢？

     
  6. cvchina
    2011年9月26日11:43 | #6
    回复 | 引用
     

    gxf :

    ORB是用什么做检测子的呢？

    FAST

     
  7. jsj063
    2011年9月26日17:32 | #7
    回复 | 引用
     

    我一直条不同啊，能否把您的工程发给俺一份啊？？gx_cheng@163.com@壮壮

     
  8. janet
    2011年9月27日05:47 | #8
    回复 | 引用
     

    新人注册, 感谢LZ分享基于最新ORB的实验结果,学习了!

     
  9. 阿棠
    2011年9月27日10:27 | #9
    回复 | 引用
     

    效果不错，不过实时视频，初步估计还是不行。

     
  10. pandawang
    2011年9月27日10:35 | #10
    回复 | 引用
     

    显示homo矩阵元素那里一直调不通。。。什么原因呢

     
  11. pandawang
    2011年9月27日10:36 | #11
    回复 | 引用
     

    homo.at(0,0), 这个总是报错

     
  12. pandawang
    2011年9月28日11:41 | #12
    回复 | 引用
     

    已经好了，homo.at(0,0)才可以。还有就是虽然orb检测时间快了些，但是match的时间增加了，因此实时还是无法做到

     
  13. 阿七
    2011年9月28日12:16 | #13
    回复 | 引用
     

    一开始就遇到error C2955: ‘std::vector’ : use of class template requires template argument list 唔 不明白

     
  14. cvchina
    2011年9月28日13:05 | #14
    回复 | 引用
     

    @pandawang
    match可以优化的。

     
  15. Anonymous
    2011年9月28日19:05 | #15
    回复 | 引用
     

    ……难道模板类的尖括号都被当做HTML语言的标签解析掉了 – -|||

     
  16. vinjn
    2011年9月28日22:31 | #16
    回复 | 引用
     

    我整理了下代码，建了个VS2008的工程，放在github上
    https://github.com/vinjn/CreativeCoding/tree/master/OpenCVExperiments/OrbTest

     
  17. leo
    2011年9月29日21:47 | #17
    回复 | 引用
     

    @vinjn

    thx

     
  18. pandawang
    2011年9月30日16:19 | #18
    回复 | 引用
     

    @cvchina
    即使优化了，感觉做实时也是有困难的。而且现在还碰到一个问题：同样的两张图片用opencv2.1的find_obj.cpp跑，速度比用opencv2.3.0跑的快。而且快的不止一倍。不知道你们遇到过没有？很奇怪

     
  19. merrylifeng
    2011年11月19日17:11 | #19
    回复 | 引用
     

    为何我的运行时间那么慢，和你的结果分析的时间相差一个数量级，很奇怪，希望能够给予解释？非常感谢。。。

     
  20. chuck
    2011年11月19日20:32 | #20
    回复 | 引用
     

    你的代码是VS2008的吗？我运行不了！我是VS2008@vinjn

     
  21. chuck
    2011年11月19日20:54 | #21
    回复 | 引用
     

    在么，能不能把你调通过的代码发我看看？ORB的shmily0923@126.com。上面那个给的VS2008的我反正没打开。@merrylifeng

     
  22. chuck
    2011年11月19日21:30 | #22
    回复 | 引用
     

    阿七 :一开始就遇到error C2955: ‘std::vector’ : use of class template requires template argument list 唔 不明白

    解决了吗？

     
  23. cvchina
    2011年11月20日11:13 | #23
    回复 | 引用
     

    你结果是啥？我怎么分析。。。@merrylifeng

     
  24. vinjn
    2011年11月21日23:58 | #24
    回复 | 引用
     

    我是vs2008
    @chuck

     
  25. captain
    2011年11月28日21:20 | #25
    回复 | 引用
     

    看了ORB原文，使用了您的ORB例子。请教您几个问题：1）对于左右图像的feature数目，您是怎么定的；2）金字塔的层数，您是如何选取呢。收获良多，谢谢！

     
  26. cvchina
    2011年11月29日17:58 | #26
    回复 | 引用
     

    这几个参数都没有做什么特别的考虑。
    唯一考虑的是将模板图像的feature数目提高，层数增多，同时降低索引图像的feature数目，层数。目的是降低匹配速度。

     
  27. rookiepig
    2011年12月7日09:49 | #27
    回复 | 引用
     

    void operator()(const Mat& img, const Mat& mask,
    vector& keypoints,
    cv::Mat& descriptors,
    bool useProvidedKeypoints=false) const;

    有没有人尝试过 使用自己提供的keypoint？

     
  28. 若众
    2011年12月13日21:19 | #28
    回复 | 引用
     

    @rookiepig
    这个应该是给出keypoint后计算相应的descriptor，没有本质区别

【推广】 免费学中医，健康全家人
原文地址：https://www.cnblogs.com/scnucs/p/2294193.html