《数字摄影测量学》课堂实验报告(遥感)

河南理工大学测绘学院

《数字摄影测量学》教学实验报告

（专业必修课）

┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄

实验成绩：

评语:

指导老师签名：

2012年月日

实习报告一：点特征提取

一、实验的目的与要求

通过实习掌握点特征提取的原理、方法和过程，并能用计算机编程实现。

二、实验设备

1、计算机一台

2、立体数字图像

3、VC 开发平台

三、实验步骤

1.编写Forstner点特征提取算子的软件，并进行调试。

2.进行BMP数字图像的点特征提取。

3.提取结果分析，评价特征点的提取的情况。

代码如下：double WINAPI Det(double N[], int n)

{

if(n==4)

return (N[0]*N[3]-N[1]*N[2]);

else

return 0;

}

//计算矩阵N之迹（即对角线元素之和）的值

//参数：大小为n的数组、数组宽度

//返回值：矩阵N之迹

double WINAPI Tr(double N[], int w)

{

double sum=0.0;

for(int i=0;i

sum += *(N+i*w+i);

return sum;

}

//判断该点是否为初选点

//参数：大小为n的数组、数组维数、阈值

//返回值：是，返回真；否，返回假

BOOL WINAPI Mid(int dg[], int n,int T)

{

BOOL flag=FALSE;

int i=0,num=0;

for(;i

if(*(dg+i)>T)

num++;

if(num>1)

flag=TRUE;

return flag;

}

///////////////////////////////////////////////////////////////////////////

//特征提取的函数（仅对灰度图象）

//参数：块图象数据指针或数组、块图象宽、块图象高、局部非最大窗口宽、局部非最大窗口高//返回值：标记特征点位的二值化块图象指针

void WINAPI ForstnerCha(BYTE image[],BYTE* newimage, int width, int height, int win_w, int win_h) {

static int T=10;//提取初选点的阈值

static double Tq=0.40;//阈值

int i,j,l,m,n;//计数器

int dg[4];//上下左右四个方向的差分绝对值

// BYTE* newimage=new BYTE[width*height];//为保存特征点开辟的数组

BYTE t,t1,t2,t3,t4;//替换值变量

//提取初选点

for(i=0;i

for(j=0;j

{

if(i==0 || j==0 || i==height-1 || j==width-1)

*(newimage+i*width+j)=0;//窗口边缘点置零

else

{

t=*(image+i*width+j);

t1=*(image+(i+1)*width+j);

t2=*(image+i*width+(j+1));

t3=*(image+(i-1)*width+j);

t4=*(image+i*width+(j-1));

dg[0]=abs(t-t1);//向下的差分绝对值

dg[1]=abs(t-t2);//向右的差分绝对值

dg[2]=abs(t-t3);//向上的差分绝对值

dg[3]=abs(t-t4);//向左的差分绝对值

//若四个方向的差分绝对值有任意两个大于阈值，为真，否则为假

if( Mid(dg,4,T) )

*(newimage+i*width+j)=255;//若为真

else

*(newimage+i*width+j)=0;//若为假

}

}

float* w=new float[width*height];//存放各点权值的数组

int k=3/2;

double N[4];//窗口的灰度协方差矩阵

double gu,gv;//像素的Robert`s梯度

float q;//兴趣值

//提取极值点

for(i=0;i

for(j=0;j

{

*(w+i*width+j)=0.0;//权值数组各元素置零

if(*(newimage+i*width+j)!=0)

{

for(l=0;l<4;l++)

N[l]=0;//灰度协方差矩阵置零

for(m=i-k;m<=i+k-1;m++)

for(n=j-k;n<=j+k-1;n++)

{

//该像素的Robert`s梯度

gu=*(image+(m+1)*width+(n+1))-*(image+m*width+n);

gv=*(image+m*width+(n+1))-*(image+(m+1)*width+n);

//计算灰度协方差矩阵

N[0]+=pow(gu,2.0);

N[1]+=gu*gv;

N[2]=N[1];

N[3]+=pow(gv,2.0);

}

q=(float)(4.0*Det(N,4)/pow(Tr(N,2),2.0));//计算兴趣值

if(q>Tq)

*(w+i*width+j)=(float)( Det(N,4)/Tr(N,2) );//若q大于阈值，计算其权值else

*(newimage+i*width+j)=0;

}

}

int w_h,w_w;

int x,y;

float dmin;

int h_l=win_h,w_l=win_w;

int h_k=h_l/2,w_k=w_l/2;

//局部非最大提取特征点

for(i=0;i

for(j=0;j

{

//控制每次局部非最大的窗口大小

if(j+w_l >= width)

w_w=width;

else

w_w=j+w_l;

if(i+h_l >= height)

w_h=height;

else

w_h=i+h_l;

x=i; y=j;

dmin=*(w+x*width+y);

float temp;

//选取权值最大点为特征点，其余点置零

for(m=i;m

for(n=j;n

{

if(*(newimage+m*width+n)!=0)

*(newimage+m*width+n)=0;

temp=*(w+m*width+n);

if( temp > dmin )

{

dmin=temp;

x=m;

y=n;

}

}

if(*(w+x*width+y)!=0)

*(newimage+x*width+y)=255;

}

delete []w;

// return newimage;//返回记录的特征点

}

四、实验结果与体会

结果处理如上图所示。

通过实习掌握点特征提取的原理、方法和过程，并能读懂别人编的代码，在本次的实习之中我觉得自己懂的东西确实的很少，还得加强努力，相信在不久的将来我也能独立的写出属于自己的代码。

实习报告二：影像匹配

一、实验目的

掌握相关系数影像匹配的原理和基本过程，并能编程实现。

二、实验设备

1、计算机一台

2、立体数字图像

3、VC 开发平台

三、实验步骤

1．使用VC++编写基于相关系数进行影像匹配的程序并进行调试。

2.进行影像匹配的运算。

3. 分析影像匹配结果，和其他函数进行对比，得出匹配质量较好，可以消除灰度的拉伸变化影响，在最小二乘法中可以判断其相关系数，根据其大小进行平差是否。

编写代码如下：

CPoint CRS_Shuguang View::lufei(unsigned char *imgBuf, int width, int height, unsigned char *blockBuf, int blockW, int blockH)

{

int lineByte=(w idth+3)/4*4;

Cxishu shu;

shu.m_xishu=1.0;

shu.DoModal();

UpdateData(true);

//标准图像每行字节数

int blockLineByte=(blockW+3)/4*4;

//循环变量

int i, j, i0, j0;

//中间变量，每搜索一个位置记录该位置上的相关性测度，本函数用对应像素差的//绝对值并求和作为相关性测度，所以sum值最小的位置为最佳匹配位置

//用min记录具有最小sum的值, 及具有最小sum的位置

int minX, minY;

//中间变量

double mum=0,sum=0,t1, t2;

double sum1=0;

double mum1=0;

double scarce=0,row,mm,max;

for(i0=0;i0

{

for(j0=0;j0

{

t2=*(blockBuf+i0*blockLineByte+j0);

sum+=t2;

}

}

sum/=(blockH*blockW);

for(i0=0;i0

{

for(j0=0;j0

{

t2=*(blockBuf+i0*blockLineByte+j0);

sum1=sum1+(t2-sum)*(t2-sum);

}

}

for(i=0;i

for(j=0;j

{

mum=0; mum1=0;scarce=0;

//搜索匹配窗口内的像素，(i0,j0)的步长也可以视情况和需要来定

for(i0=0;i0

{

for(j0=0;j0

{

t1=*(imgBuf+(i+i0)*lineByte+j+j0);

//t2=*(blockBuf+i0*blockLineByte+j0);

mum += t1;;

}

}

mum/=(blockH*blockW);

/

for(i0=0;i0

{

for(j0=0;j0

{

/

t1=*(imgBuf+(i+i0)*lineByte+j+j0);

mum1+=(t1-mum)*(t1-mum);

}

}

for(i0=0;i0

{

for(j0=0;j0

{

/

t1=*(imgBuf+(i+i0)*lineByte+j+j0);

t2=*(blockBuf+i0*blockLineByte+j0);

scarce += (t2-sum)*(t1-mum);

}

}

mm= sqrt((sum1*mum1));

row=scarce/ mm;

//第一个搜索区域，此时也是最小sum值的位置if(i==0&&j==0)

{

max=row;

minX=0;

minY=0;

}

else

{//已经不是第一个搜索位置了

if(max

{//出现一个新的最小位置，记录此时的最小值及其位置坐标

max=row;

minX=j;

minY=i;

}

}

}

if(fabs(max)>shu.m_xishu)

{//返回匹配位置

CPoint pos(minX,minY);

return pos;

}

else{ CPoint pos(0,0);

return pos;}

}

四、实验结果与体会

左灯笼

匹配结果，如图所示。（此图的匹配点位于左下角的角点，利用矩形内的灰度变化进行相关

系数的匹配）

通过本次的实习我掌握了相关系数影像匹配的原理和基本过程，我能读懂了别人编辑的代码，在本次的实习中我遇到了很多的困难，通过问老师和同学都解决了。所以我今后还得继续的努力争取写出属于自己的代码。