📄 image.cpp
字号:
// Image.cpp: implementation of the CImg class.
//
//////////////////////////////////////////////////////////////////////
#include "stdafx.h"
#include "disparity 200712.h"
#include "Img.h"
#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif
//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////
CImg::CImg()
{
kp = new CKeypoint[N_KP];
n_kp = 0;
}
CImg::~CImg()
{
}
void CImg::CreateImg(char *filename)
{
printf("读取图像...");
image_c = cvLoadImage( filename, 1); //读彩色图像
image = cvLoadImage( filename, 0); //强制转换为单色
width = image->width; //保存图像的尺寸
height = image->height;
image32f = cvCreateImage ( cvSize(width,height), IPL_DEPTH_32F, 1);
cvConvert( image, image32f); //转换为32位浮点图像
//构造高斯图像组Gau[*]
printf("ok\n");
for (int i=0; i<LEVEL; i++)
{
printf("构造高斯图像组...%d%%\r",(int)(i*100.0/LEVEL));
Gau[i] = cvCreateImage ( cvSize(width,height), IPL_DEPTH_32F, 1);
cvSmooth( image32f, Gau[i], CV_GAUSSIAN, 0, 0, pow(RATE,i));
}
//计算高斯图像的差分Dog[*]
printf("构造高斯图像组...ok \n");
printf("构造高斯图像的差分...");
for (i=0; i<LEVEL-1; i++)
{
Dog[i] = cvCreateImage ( cvSize(width,height), IPL_DEPTH_32F, 1);
cvSub( Gau[i+1], Gau[i], Dog[i] );
}
printf("ok\n");
}
void CImg::FindKeypoint(void)
{
int margin; //图像边缘留出的空白宽度
double TH;
for(int i=1; i<LEVEL-2; i++)//高斯模糊等级
{
//显示提示
printf("查找特征点...%d%%\r",(int)(i*100.0/(LEVEL-1)) );
margin = int(6*pow(RATE,i)); //搜索特征点时要在图像边缘留出一定宽度，以防止出界
TH = 0; //域值，但似乎作用不大，所以设为0
for (int x=margin; x<width-margin; x++)
{
for (int y=margin; y<height-margin; y++)
{
//查找在 3x3x3 临域中的极值
if ( (ELEM(Dog[i],x,y) - ELEM(Dog[i],x+1,y) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x-1,y) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x+1,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x-1,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x+1,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x-1,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x+1,y) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x-1,y) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x+1,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x-1,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x+1,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x-1,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x,y) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x+1,y) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x-1,y) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x+1,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x-1,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x+1,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x-1,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x,y+1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x,y-1) > TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x,y)> TH )
|| (ELEM(Dog[i],x,y) - ELEM(Dog[i],x+1,y) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x-1,y) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x+1,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x-1,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x+1,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x-1,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i],x,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x+1,y) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x-1,y) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x+1,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x-1,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x+1,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x-1,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i+1],x,y) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x+1,y) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x-1,y) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x+1,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x-1,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x+1,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x-1,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x,y+1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x,y-1) < -TH
&& ELEM(Dog[i],x,y) - ELEM(Dog[i-1],x,y)< -TH ) )
//下面这一段是把图像的像素当作正六边形（蜂巢形）来处理，但效果不理想
/*if(
(ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i],x+1,y) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i],x-1,y) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i],x+0.5,y+0.5*sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i],x+0.5,y-0.5*sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i],x-0.5,y+0.5*sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i],x-0.5,y-0.5*sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i+1],x,y+1/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i+1],x+0.5,y-0.5/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i+1],x+0.5,y-0.5/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i-1],x,y+1/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i-1],x+0.5,y-0.5/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)>ELEM_IN(Dog[i-1],x+0.5,y-0.5/sqrt(3)) ) ||
(ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i],x+1,y) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i],x-1,y) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i],x+0.5,y+0.5*sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i],x+0.5,y-0.5*sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i],x-0.5,y+0.5*sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i],x-0.5,y-0.5*sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i+1],x,y+1/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i+1],x+0.5,y-0.5/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i+1],x+0.5,y-0.5/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i-1],x,y+1/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i-1],x+0.5,y-0.5/sqrt(3)) &&
ELEM_IN(Dog[i],x,y)<ELEM_IN(Dog[i-1],x+0.5,y-0.5/sqrt(3)) ) )*/
{
if ( fabs(ELEM(Dog[i],x,y)) > 3 )//这一行用来排除特征不明显的点
{
//以下用来排除曲率太小的点（边缘点）
double Dxx,Dyy,Dxy,Tr,Det,r,threshold;
double rate = pow(RATE,i);
Dxx = (ELEM_IN(Dog[i],x+rate,y) + ELEM_IN(Dog[i],x-rate,y) -2*ELEM_IN(Dog[i],x,y));
Dyy = (ELEM_IN(Dog[i],x,y+rate) + ELEM_IN(Dog[i],x,y-rate) -2*ELEM_IN(Dog[i],x,y));
Dxy = (( ELEM_IN(Dog[i],x+rate,y+rate) + ELEM_IN(Dog[i],x-rate,y-rate) - ELEM_IN(Dog[i],x-rate,y+rate) - ELEM_IN(Dog[i],x+rate,y-rate)))/4;
Tr = Dxx + Dyy;
Det = Dxx*Dyy - Dxy*Dxy;
r=5;
threshold = (r+1)*(r+1)/r;
if ( fabs(Tr*Tr/Det) < threshold )//(Tr*Tr/Det)表示的是特征点的曲率，其绝对值大于一定值才保留该点
{
//最后剩下的特征点,保存下来
kp[n_kp].level = i; //保存高斯模糊的等级
kp[n_kp].x = x;
kp[n_kp].y = y; //保存特征点的坐标
if(n_kp<N_KP-1)
{
n_kp++; //特征点个数加一（防止出界N_KP）
}
}
}
}
}
}
}
//显示提示
printf("特征点查找...ok \n共找到%d个特征点\n", n_kp);
}
void CImg::DrawKeypoint(void)
{
for (int n=0; n<n_kp; n++)
{
int x=kp[n].x;
int y=kp[n].y;
float length = (float)pow(RATE,kp[n].level); //方向矢量的长度
//if(kp[n].level==1)
{
//在彩色图像上用白色“+”表示出特征点的位置
cvLine( image_c, cvPoint(x-1,y), cvPoint(x+1,y), cvScalar(255,255,255) );
cvLine( image_c, cvPoint(x,y-1), cvPoint(x,y+1), cvScalar(255,255,255) );
//画出特征点的方向矢量
cvLine( image_c, cvPoint(x,y), cvPoint( x+(int)(cos(kp[n].direction*PI/18)*length), y+(int)(sin(kp[n].direction*PI/18)*length)), cvScalar(255,255,255) );
}
}
}
void CImg::ShowImg(char *title)
{
cvNamedWindow(title, 1 );
cvShowImage( title, image_c ); //显示彩色图像
}
//计算特征点的主方向和副方向
void CImg::Direction4Kp(void)
{
float x,y;
int level;
float rate;
float temp; //临时变量
float m,q; //m表示模，q表示角度
int n_kp_t = n_kp; //由Dog的极值所得到的特征点的数量（因为有可能增加）
printf("计算特征点的方向...");
for (int n=0; n<n_kp_t; n++)
{
level = kp[n].level;
rate = 2*(float)pow(RATE,level);
//构造方向直方图
for (x=-2; x<=2; x+=.5)
{
for (y=-2; y<=2; y+=.5)
{
if( x*x+y*y < 4)
{
m = (float)sqrt(pow(ELEM_IN(Gau[level],kp[n].x+(x+1)*rate,kp[n].y+y*rate)-ELEM_IN(Gau[level],kp[n].x+(x-1)*rate,kp[n].y+y*rate),2)+pow(ELEM_IN(Gau[level],kp[n].x+x*rate,kp[n].y+(y+1)*rate)-ELEM_IN(Gau[level],kp[n].x+x*rate,kp[n].y+(y-1)*rate),2));
q = (float)atan2(ELEM_IN(Gau[level],kp[n].x+x*rate,kp[n].y+(y+1)*rate)-ELEM_IN(Gau[level],kp[n].x+x*rate,kp[n].y+(y-1)*rate), ELEM_IN(Gau[level],kp[n].x+(x+1)*rate,kp[n].y+y*rate)-ELEM_IN(Gau[level],kp[n].x+(x-1)*rate,kp[n].y+y*rate));
kp[n].orient[(int)floor((q+2*PI)/(PI/18))%36] += m*(float)exp(-(x*x+y*y)/2);
}
}
}
//选出最大方向
temp = kp[n].orient[0];
for (int d=1; d<36; d++)
{
if (kp[n].orient[d]>temp)
{
temp = kp[n].orient[d];
kp[n].direction = d;
}
}
//*/选出副方向
for (d=0; d<36; d++)
{
if ( d!=kp[n].direction && kp[n].orient[d]>0.9*kp[n].orient[kp[n].direction])
{
kp[n_kp].level = kp[n].level;
kp[n_kp].x = kp[n].x;
kp[n_kp].y = kp[n].y;
kp[n_kp].direction = d;
if(n_kp<N_KP-1)
{
n_kp++;
}
}
}//*///end of 选出副方向
}//end of for
printf("ok\n");
printf("特征点的数量增加到%d个\n",n_kp);
}
void CImg::ConstructDescripter(void)
{
double region_kp[12][12]; //特征点的 12x12 的临域
double direction; //主方向(0-2*PI)
double x,y; //临域中每一点在原图中的对应坐标
double rate;
printf("提取特征...");
for (int n=0; n<n_kp; n++)
{
direction = kp[n].direction*PI/18; //主方向(0-2*PI)
rate = 1.5*pow(RATE,kp[n].level); //特征点所在的等级的比例
//为特征点构造一个12x12的临域
for (int i=0; i<12; i++)
{
for (int j=0; j<12; j++)
{
//一些几何变换
x = (i-5.5*rate)*cos(direction)+(j-5.5*rate)*sin(direction) + kp[n].x;
y = -(i-5.5*rate)*sin(direction)+(j-5.5*rate)*cos(direction) + kp[n].y;
if ( x<=0 || x>=width-1 || y<=0 || y>=height-1)
{
//如果落在图像边界之外则给一个平均值128
region_kp[i][j] = 128;
}
else
{
//如果在图像边界之内则采用插值
region_kp[i][j] = ELEM_IN(Gau[kp[n].level], x, y);
}
}
}//end of 构造临域
//计算特征点的128维特征
for (int x=0; x<4; x++)
{
for (int y=0; y<4; y++)
{
int i = x*3+1;
int j = y*3+1;
//注释掉的一段：文献上的方法，其实效果差不多
/* kp[n].feature[0+(y*4+x)*8] = region_kp[i][j] - region_kp[i+1][j];
kp[n].feature[1+(y*4+x)*8] = region_kp[i][j] - region_kp[i][j+1];
kp[n].feature[2+(y*4+x)*8] = region_kp[i][j] - region_kp[i-1][j];
kp[n].feature[3+(y*4+x)*8] = region_kp[i][j] - region_kp[i][j-1];
kp[n].feature[4+(y*4+x)*8] = region_kp[i][j] - region_kp[i+1][j+1];
kp[n].feature[5+(y*4+x)*8] = region_kp[i][j] - region_kp[i+1][j-1];
kp[n].feature[6+(y*4+x)*8] = region_kp[i][j] - region_kp[i-1][j+1];
kp[n].feature[7+(y*4+x)*8] = region_kp[i][j] - region_kp[i-1][j-1];
}
}*///end of 计算特征
kp[n].feature[0+(y*4+x)*4] = region_kp[i][j] - region_kp[i+1][j];
kp[n].feature[1+(y*4+x)*4] = region_kp[i][j] - region_kp[i][j+1];
kp[n].feature[2+(y*4+x)*4] = region_kp[i][j] - region_kp[i-1][j];
kp[n].feature[3+(y*4+x)*4] = region_kp[i][j] - region_kp[i][j-1];
}
}
{
kp[n].feature[64] = region_kp[0][6] - region_kp[11][6];
kp[n].feature[65] = region_kp[6][0] - region_kp[6][11];
kp[n].feature[66] = region_kp[0][3] - region_kp[6][3];
kp[n].feature[67] = region_kp[3][0] - region_kp[3][6];
kp[n].feature[68] = region_kp[0][9] - region_kp[6][9];
kp[n].feature[69] = region_kp[9][0] - region_kp[9][6];
kp[n].feature[70] = region_kp[3][6] - region_kp[3][11];
kp[n].feature[71] = region_kp[6][3] - region_kp[11][3];
kp[n].feature[72] = region_kp[9][6] - region_kp[9][11];
kp[n].feature[73] = region_kp[6][9] - region_kp[11][9];
}
}//end of for
// printf("ok\n");
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -