VS2010+GDAL+openCV的遥感图像K均值算法的实现和植被指数的计算 联系客服

发布时间 : 星期二 文章VS2010+GDAL+openCV的遥感图像K均值算法的实现和植被指数的计算更新完毕开始阅读83e0880b0640be1e650e52ea551810a6f524c8ba

像作为样本存储在该样本矩阵中,每一个波段为一列,行数为影像的长宽之积。同时定义整型样本标识矩阵,矩阵大小为样本矩阵的一列。构造如下函数:

::imgKMeans(CvMat* sample,IplImage *pResult,int nClusters) Sample为样本矩阵,pResult,clusters为输出的聚类结果,nClusters为聚类数。采用openCV工具的cvKmeans2()函数可以实现K均值聚类:

cvKMeans2(sample,nClusters,clusters ,cvTermCriteria(CV_TERMCRIT_ITER ,20,1.0));

cvTermCriteria(CV_TERMCRIT_ITER ,20,1.0)为迭代终止规则,表示迭代次数与迭代精度的控制。

分类类数可以控制在[2,20]之间,分类太少没有意义,分类数太多往往将同类地物分开,效果不好。不同分类数的K均值聚类结果如图4。

2.3 植被指数的计算

本次实习中实现了四中植被指数的计算过程,包括归一化差分植被指数(NDNI),比值植被指数(RVI),土壤调节植被指数(SAVI)和差值环境植被指数(DVI)。各植被指数的计算公式如上文所述。在进行计算前,已经将各个波段的数据一次存在一个指针数组之中,计算哪几个波段就将哪几个波段的数据读出即可。

2.3.1 计算NDVI指数

依次进行矩阵元素级减法,加法和除法运算: cvSub(image[3],image[2],res1);//减法 cvAdd(image[3],image[2],res2);//加法 cvDiv(res1,res2,res);//除法

image[3]表示图像的近红外波段,image[2]表示图像的红光波段。得出的res值的大小即为NDVI指数。为了直观的显示NDVI指数的效果,结合NDVI不同值所表示的植被覆盖情况:

A. -1<=NDVI<0 表示地面覆盖为云、水、雪等,对可见光高反射; B. NDVI=0

表示有岩石或裸土等,NIR和R近似相等;

C. 0

则可以对大于0值的NDVI对应的区域按照NDVI值的大小在绿色通道上赋不同大小的像素值,得到NDVI图像特征:所有绿色区域为植被覆盖区域,灰色区域没有植被覆盖;绿色越亮的区域植被覆盖越茂密,绿色越暗的区域植被覆盖度越低。分类结果见图5。

5

2.3.2 计算RVI

RVI的计算方法为图像的近红外波段与红光波段的比值。RVI值的范围仍然在原图的像素范围之内。计算比较简单:

cvDiv(image[3],image[2],res);//除法

由于植被的RVI通常大于2,即RVI>=2为植被覆盖区域,设为绿色,同NDVI,绿色越亮的区域植被覆盖越茂密,绿色越暗的区域植被覆盖度越低;其他区域为灰度。分类结果见图5。 2.3.3 计算SAVI

SAVI=((NIR-R)/(NIR+R+L))(1+L)。L是随着植被密度变化的参数,取值范围从0到1,当植被覆盖度很高时L值为0,很低时L值为1。很明显L=0时,SAVI=NDVI。但对某幅图像进行植被指数计算时并不知道其植被覆盖度如何,故采用L=0.5:

cvSub(image[3],image[2],res1);//相减 cvAdd(image[3],image[2],res2);//相加 CvScalar nplus;nplus.val[0]=0.5f; cvAddS(res2,nplus,res2);//加0.5 cvDiv(res1,res2,res);//除法 …

double s = cvGetReal2D(res,i,j); s = s*1.5f;//每个像素值乘以1.5 …

其计算过程类似于NDVI,故显示方法同NDVI指数。 2.3.4计算DVI

DVI=NIR-R。即近红外波段减去红光波段的值: cvSub(image[3],image[2],res);//相减

相减之后得出的值就是DVI指数,为了直观显示和可视化分析,也将DVI指数拉伸到[0,255]并显示出来:

CvScalar ps;

double s = cvGetReal2D(res,i,j);

s=(s-mm[0])/(mm[1]-mm[0])*255.0f;//进行最大最小值拉伸 ps.val[0] = s; …

6

显示结果如图5。

3 实验结果与分析

3.1图像的读取与显示

下图是打开16位4波段img格式的影像后的显示结果,左图是单波段显示波段1的结果,右图是3波段(1,4,3)组合显示:

图1 波段1(左图),波段1,4,3(右图)

读取其他格式如8位3波段*.bmp格式图像:

图2单波段(左图),三波段(右图)

7

3.2K均值聚类

由分类结果可知,分类结果过多会导致相同地物被分为多类,所以合理选择分类数很必要。且随着分类图像波段数的增多以及图像长宽的变大都会使分类时间变长。

图3 进行K均值聚类的4波段图像(显示:1,2,3)

图4K均值算法依次分4,6,10,15类的结果

8