遥感影像分类方法试验报告

1、实验报告目录 TOC o 1-5 h z HYPERLINK l bookmark2 o Current Document 1实验目的 4 HYPERLINK l bookmark4 o Current Document 2实验数据 4 HYPERLINK l bookmark6 o Current Document 3实验内容 4 HYPERLINK l bookmark8 o Current Document 4实验步骤 5 HYPERLINK l bookmark10 o Current Document 对人口矢量数据（shapefile ）进行投影转换 5Census.shp 文件投

2、影坐标的检查 5将投影坐标转换为 WGS_1984_UTM_Zone_16N 6 HYPERLINK l bookmark32 o Current Document 对遥感影像进行几何精校正（以经过投影变换的人口矢量数据为基准）.6Census.shp 在 ENVI 软件的加载 6 HYPERLINK l bookmark36 o Current Document 对遥感影像进行几何精校正（以矢量数据为基准） 7 HYPERLINK l bookmark46 o Current Document 用矢量图层对遥感影像进行裁剪 10 HYPERLINK l bookmark57 o Curren

3、t Document 将Pan波段和多光谱波段进行融合，并对融合效果进行定性和定量评价.11两种融合方法的原理 11进行 Gram-Schmidt Spectral Sharpening 融合 11融合效果进行定性评价 14融合效果进行定量评价（软件提供的计算方法） 15融合效果进行定量评价（Matlab编程at算） 16遥感影像融合定量分析代码 20 HYPERLINK l bookmark119 o Current Document 生成住房密度栅格影像 23两表的连接 23计算房屋密度 24 HYPERLINK l bookmark130 o Current Document 直接栅格化

4、 25IDW 插值 25对房屋密度图进行重分类 26 HYPERLINK l bookmark142 o Current Document 将住房密度栅格影像作为额外的通道与ETM多光谱波段进行叠加 26 HYPERLINK l bookmark146 o Current Document 监督分类（融合方法为 HSV波段为5, 4, 3） 27打开Google Earth 影像作为监督分类的参照 27建立兴趣区 29 HYPERLINK l bookmark174 o Current Document 训练样区的选择 30 HYPERLINK l bookmark176 o Current

5、Document 训练样区的评价 31执行监督分类 33 HYPERLINK l bookmark202 o Current Document 分类后处理 35 HYPERLINK l bookmark222 o Current Document 评价结果分析 37分类结果面积统计 38分类结果 41 HYPERLINK l bookmark248 o Current Document 分类结果评价与分析 41 HYPERLINK l bookmark250 o Current Document 未加入房屋密度图层的分类结果评价与分析 41加入IDW插值房屋密度图层的分类结果评价与分析 42加

6、入直接栅格化房屋密度图层的分类结果评价与分析 43加入重分类后IDW插值房屋密度图层的分类结果评价与分析 44 HYPERLINK l bookmark268 o Current Document 从总精度与Kappa系数对分类结果进行评价 45 HYPERLINK l bookmark292 o Current Document 分类结果总体评价 46与其他训练样区的分类精度和Kappa系数的计算 48 HYPERLINK l bookmark308 o Current Document 决策树分类 49决策树分类原理 49数据预处理 49指数的计算 51执行决策树 54不同参数设置的对比

7、575实验体会 60实验中存在的问题 60软件平台使用 63实验总结 631实验目的掌握ArcGIS10和ENVI4.7对遥感图像处理的基本操作与原理熟悉几何精校正的方法，掌握 ENVI软件对遥感影像进行几何精校正掌握全色波段与多光谱波段的融合方法和原理，学会对融合效果进行定性定量分析熟悉掌握ArcGIS的栅格化方法和IDW插值方法熟悉监督分类的方法和基本原理，掌握ENVI软件中进行监督分类了解监督分类后评价过程，对分类结果进行精度评价和分析掌握Erdas的空间建模方法以及原理了解RuleGen算法，掌握决策树分类方法2实验数据带属性数据的shapefile : Census.shp 带有陆地

8、面积字段的矢量图层: GoogleEarth_ 原始拼接：GE1005211134.jpg 研究区域的多光谱波段数据：Stack_b1-6162-7.img 研究区域的全色波段数据：b8.img监督分类参照影像：Google Earth3实验内容对人口矢量数据（shapefile ）进行投影转换： WGS_1984_UTM_Zone_16N对遥感影像进行几何精校正（以经过投影变换的人口矢量数据为基准）：（1）对多光谱波段（30 m空间分辨率）进行几何精校正（小于 0.25个像元）；（2）对Pan波段（15 m空间分辨率）进行几何精校正（小于 0.25个像元）；将Pan波段和多光谱波段进行融合（

9、自选至少一种融合算法），并对融合效果进行定性和定量评价；生成住房密度栅格影像：（1）直接栅格化；IDW 插值；将住房密度栅格影像作为额外的通道(或波段)与ETM侈光谱波段进行叠加；进行监督分类和分类后处理( Post-Classification , Expert Rules )利用ERDA函件的空间建模(Spatial Modeler )进行水体信息(MNDW旨数)和植 被信息(NDVI指数)的提取；利用“自动阈值决策树分类算法”进 Marion County的土地利用/覆盖分类信息提取(使用的数据：原始各波段 +MNDWI+NDVI+ISODATA,或其他有益的波段组合) 探讨“自动阈值决

10、策树分类算法”中的各个参数意义及如何设置更合理对分类结果进行评价与分析对分类结果进行精度评价和分析；4实验步骤对人口矢量数据(shapefile )进行投影转换Census.shp文件投影坐标的检查根据实验要求，人口矢量数据(shapefile )进行投影坐标应为： WGS_1984_UTM_Zone_16N在ArcGIS软件的图层右击 Properties ,在Layer Properties 的Source下查看投影信息，如图1。得到Census.shp的投影坐标为：GCS_North_American_1983 ,与实验要求不符合，需进行投影转换。Layer PropertieiI 餐图

11、1将投影坐标转换为 WGS_1984_UTM_Zone_16N在Arcgis软件的工具箱中的Define Projection 工具，设置输入数据为：Census.shp ,坐标系统为GCS_WGS_198麻工具箱中的工具，设置输入数据：Census.shp ,导入遥感影像的投影坐标系，即GCS_WGS_1984如图2）。图2对遥感影像进行几何精校正（以经过投影变换的人口矢量数据为基准）Census.shp 在 ENVI 软件的加载在ENVI软件中，File - Open Vector File ,选择 Census.shp ,设置好参数，生成evf文件（如图3）。图30K :对遥感影像进行几

12、何精校正（以矢量数据为基准）在ENVI分别打开遥感影像img和矢量文件vef,选择Map - Registration - Select GCPs :image to map,以矢量数据为基准，设置好投影，如图 4。RegistTiartJDnRjgorDUE OrthorEdilicErbonOrrfiorectificafrienMtwickirgSelect GO% Image ：d nugeSelect Image m MapGEareiEnerKt fram Input Geomdry4 fWarp from GCP Image to linageWarp From GCPs； Im

13、age to MapAutcmaLc fieg?trjtion: Image tD bn呼在ENVI中，在zoom窗口下采集控制点，这次实验采集的控制点数为13个,控制点的主要定位在道路与道路之间的交叉点,如图5,其控制点的 RMS Error为0.246390,如图6,如图7为20个控制点的采集情况。图5r Im啊母 to Map GCP U5t1三5File Optioris如Inap r ii4c x工附解甘Yl七dim 尸飞加打丫 Error XErrcr 丫郦I 1#14-5BZ066. 04S12652,依丽.各L379.5TT40.Z5671379. 331ED 0137-0.

14、23日4. 2391*58365. 104409712. 3E 1416. W1479.001416. 98Sf 147S. 110D. 1130. USS1*34-551698. 5731526. 7E 1395. 23L8&4.4I1395 DTT2t J8B4 13E-D 2tB0P 26640,3300570131.304.4003. St 1009.071S21.931009.20211522. ISieO 132*T0.23160. 2669*5+曲L,E4制旺30. 514g 5M03 3rm】丽口 TTOsQ 40T30.0602QTIBQ566633.604.3936T& 2

15、692.932013.078ft3.0&7720I2.7&710-15774). 30290.341S*7+50654315543-94329. 5 1556 79l瓯兆1556 &S1H9D S9S2-D 1519Q 3】健Q. 3457559711.77W0HI5J5.7EG61.9S1765.97662.01261163. 120,OS2fi0. 154120.16BS5T1675.494432566, IQW. S3l迎571060 56S336Z 516= Q 0553P 05350, OTTO#10*57&4J60.11435055. 2JC1220. ?!1220. 503i 1

16、66. 6ME CL 空轴0. 0150. 26*15值通554388401292.21139.031?92 21103B 975C0 0719-0. 05500. 0905#12+S&SBi工的401256.769.0717M. 527&8.T9201759.前氓-0. 2T8fl0. 13050.33E5* LI* 1口红UQ. | 西加七.K&de Lisi选择校正参数输出结果，在Ground Control Points Selection窗口选择 Option - WarpFile，如图8;数学模型为Polynomial ,设定参数为2,从采样方法为最临近法，如图 9。File Op

17、tions | HelpWarp Displayed B渺tLWarp File,Warp Delayed Band (as Image to Map.岫Reverse脉训抑p 7i.二kt Degree (RSI Only)=Auto Predict/ R#gistrjti&n Param*t看如图10、11为几何校正前后，矢量图层与遥感影像吻合度的对比，可以明显看出，经过几何校正后的遥感影像与矢量图层吻合程度有明显的改善，有部分水体边界不吻合， 这主要是由于水体会随时间而改变；而街区与道路吻合程度良好。图10图1113如图12为对多光谱波段（30 m空间分辨率）进行几何精校正后的遥感影像；

18、如图 为对Pan波段（15 m空间分辨率）进行几何精校正后的遥感影像。图12图13用矢量图层对遥感影像进行裁剪在ENVI中打开Census.shp文件，将该图层转换成ROI,如图14。*1 IectfK4.% J= 回PMa* 前 O0ww H 而Open VkDck FH*, fCraate New Layer.力*1 kJ牲r4 to Te口值整.RffsbEr? Lasers From Templffir.Export A。证泮r 2 ROtf.EjqKrt A6rt Ldr CD 珈）芦网.Calculate Bufer Zone,Cwd图14通过ROI进行裁剪遥感影像，选才B Bas

19、ic Tools -Subset Data via ROIs ;选择转换好的ROI进行裁剪，如图15。图1510同样，对pan波段的遥感影像进行裁剪，得到遥感影像如图16。图164.3 将Pan波段和多光谱波段进行融合，并对融合效果进行定性和定量评价两种融合方法的原理Gram-Schmidt可以对具有高分辨率的高光谱数据进行锐化。从低分辨率的波谱波段中模拟出一个全色波段。对该全色波段和波谱波段进行Gram-Schmidt变换，其中模拟的全色波段被作为第一个波段。用Gram-Schmidt变换后的第一个波段替换高空间分辨率的全色波段。应用Gram-Schmidt反变换构成pan锐化后的波谱波段。

20、Gram-Schmidt Spectral Sharpening方法进行图像增强能够比较好的保留原多光谱图像的光谱信息，使遥感影像的融合保留多光谱影像的增强效果。用PC可以对具有高空间分辨率的光谱图像进行锐化。先对多光谱数据进行主成分变换。用高分辨率波段替换第一主成分波段，在此之前，高分辨率波段已被缩放匹配到第一主成分波段，从而避免波谱信息失真。进行主成分逆变换。函数自动地用最近邻、 双线性或三次卷积技术将多光谱数据重 采样到高分辨率像元尺寸。进行 Gram-Schmidt Spectral Sharpening 融合在 ENVI 软件中，选择 Transform - Image Sharpe

21、ning - Gram-Schmidt SpectralSharpening ,在 Select Low Spatial Resolution Multi Band Input File 中选择 Stack_b1-6162-7_CJ.img 多光谱波段，Spatial Subset 为 Full Scene , Spectral Sunset 为8个波段，如图17。11图17在 Select High Spatial Resolution Pan Input File 段，如图18。窗口中选择 b8_CJ.img全色波Q1 Sdect High Spatigl RewluiiiDn Pamj.

22、 t5u.icl input Bu-ii：Z g) 5Uck_bl-6I.E?-7_CJ:期工刖隼+ffM-p 口用MT 口初IKrsi jfQ.iy*r口AOIKniztffwpiliiTtr口KOIEiEiztffrpCLiytr、口hOIKi-siitfftrpili-ytr口KOIAn in：力树pLkytr:口MIKifiit4tMeCL中衬,口101Kaljlt4tHiId44r；+ 曜 Hifii 3日固WJJ JU I J；ai3 X FkTP 1133 EBl Sail e 僧自心布CBmd1 :LT l02L032_CK2CffiMB06 JLG. TU) 3E5k_H6L

23、fl2T. 阴除T】 tnd Inf ctbjiI L-anm i ： =TDins 2225 x 22 I I LEQSLe Byti 1.903,350 IrLiii fil Typt Ekfl StuxdrdTunlaa Image Sharpening - PC Spectral Sharpening在 Select Low Spatial Resolution Multi Band Input File 中选择 Stack_b1-6162-7_CJ.img 多光谱波段，Spatial Subset 为 Full Scene , Spectral Sunset 为 8 个波段，如图 2

24、0。图20在 Select High Spatial Resolution Pan Input File窗口中选择 b8_CJ.img全色波段，如图21。Select SpatialInput Fi eStLtrt Input Bud压5第+ 印d * 国 W,CJBOI Rai;三下作IBW据畀段 2225 x 220fl n 1 解布工电 I班3】35c bflts.Fil* Typt5t 100 ;说明了融合后影像的亮度信息提高了。但PC和GS两者的Value通过图上很难分辨出区别。iao50号M_ 200320SampleHorizorrtol Profile250图 26 X Pro

25、fileMrii3 PpR加国h 口HnlpVertical ProfileG Fl Vertical Prodi1#2 件也H IWik25C28D Q Q D- 5 fl 5II53 Vertal ProflfeFd/ EdH PptK?rvw Pforl %nftKri Hplp#3 Vertical PraflleView Stitisties Fik$um BdfbdtGrrwjte Rsndmra Sjmpe3u a EgMiMu潜Ep D*Ka Lajrer StackingCctft Da |BSQ Blk 日晒Slrrtch DataEuikiki如而.tahs&csChan

26、g：e Detection TmlD/ “小28。图 27 Y Profile在 ENVI 4.7 软件中，选择 Basic Tools -Stayistics -Compute Statistics15图28如表1为原始影像、PC融合后影像和GS融合后影像的部分统计参数。晶:稣康信强白感比米图29蹒雁融理驰晡均值方差band281,3476921.59512原始影像band379.2623232.45205band498.1691620. 65547band281,3970921.85553PC融合ban.d379.6408333, 08145band4IT20121.71396band28

27、1.3107819.67364暨融合bandST9. 2934927.6280920. 623%band498.182OE表1从图29、30,两张折线图上我们可以直观的了解原始影像、PC融合后影像和 GS融合后影像在亮度信息上的关系。在均值上，PC融合影像和GS融合影像均有所增加，但是从幅度上来看，PC融合影像的均值增加的幅度较大，较GS融合影像更有利于目视判读。在方差上，与原始影像相比较 PC融合影像的方差有所增加，而GS融合影像的方差大幅下降，不利于目 视判读。因此，在亮度指标上最好的为PC融合，其次为 GS融合。融合效果进行定量评价 (Matlab编程计算)此次实验从融合影像的亮度信息，

28、清晰程度，光谱保持程度，信息丰富程度等多角度进行评价分析，相比传统的单一定量评价全面，能够减少评价的随机性， 使得定量评价更加科学全面。主要通过 4方面进行统计分析：亮度信息，针对融合后影像亮度信息进行评价， 主要包括均值和标准差等指标。清晰度，评价融合影像的清晰程度，主要包括平均梯度和 空间频率等指标。光谱信息，评价融合后影像的光谱变形情况 ，包括扭曲程度，相关系数等 指标。信息量，评价融合后信息量的保持或增加程度，主要包括信息嫡,交叉嫡和互信息等指标。以下公式中，影像函数为 Z(X, Y),影像的行数和列数分别为 M和N,影像的大小则 为 MX N。164.3.6.1 亮度信息均值(z)影

29、像均值是像素的灰度平均值 ，对人眼反映为平均亮度。如果均值适中，则影像效果良好 其定义为：M N、 z(Xi, yj) i 1 j 1 Z M N标准差(b)标准差反映了影像灰度相对于灰度平均值的离散情况，在某种程度上，标准差也可用来评价影像反差的大小。若标准差大，则影像灰度级分布分散，影像的反差大，可以看出更多的信息，其公式为： TOC o 1-5 h z HYPERLINK l bookmark71 o Current Document M N2工工(Z(Xi - yj ) -之)二i ,M N方差(D)方差反映了灰度相对于灰度均值的离散情况，方差越大，则灰度级分布越分散。此时, 图像中所

30、有灰度级出现的概率越趋于相等，从而包含的信息量趋于越大。方差计算公式为:n2D = Pi(DNi -DN) i 14.3.6.2 清晰度平均才度(G)平均梯度可敏感地反映影像对微小细节反差表达的能力，可用来评价影像的清晰程度, 还可同时反映出影像中微小细节反差和纹理变化特征，其计算公式为：jZ(xi,yj) 2 JZ(x,yj) 21(M -1)(N -1)(-)(-)二 Xi二 yi2空间频率(SF)空间频率反映了一幅影像空间的总体活跃程度，空间频率越大，说明融合效果越好。 它包括空间行频率 RF和空间列频率 CF组成，其公式为：17RF11M N i4% Z(x,yj).Z(x,yj,)2

31、j 41 M NCF=：M N：yj)-Z(Qyj)2总体的空间频率值取 RF和CF的均方根，即:222SF = . RF2 CF2光谱信息扭曲程度(D)扭曲程度D是指融合影像F像素灰度平均值与源影像 A像素灰度平均值之差，也可以说是融合影像与源影像的差值影像的灰度平均值，它的表达式为：M |F(x,yi) - A(Xi, yi)j 1它反映融合影像和源多光谱影像在光谱信息上的差异大小和光谱特性变化的平均程度。相关系数(p)融合影像与源影像的相关系数能反映两幅影像光谱特征的相似程度，其定义如下:M N|F(Xi, yj) - f A(Xi, yj) - a:二 y jj回N.必 |F(Xi,y

32、j)-f A(x, yj) -a其中，f和a分别为融合影像与源影像的均值。通过比较融合前后的影像相关系数可以看出多光谱影像的光谱信息的改变程度o信息量信息I嫡(E)影像的信息嫡值是衡量影像信息丰富程度的一个重要指标。对于一幅单独的影像，其灰度分布为p=p0,R,，R,pL_J,Pi为灰度值等于i的像素数与影像总像素数之比。根据Shannon信息论白定理，一幅影像的信息嫡为：L 4E 二 pi log 2 pii =0融合影像的嫡越大，融合影像所含的信息越丰富，融合质量越好。交叉嫡(C)18交叉嫡可以用来测定两幅影像灰度分布的信息差异。设源影像和融合影像的灰度分布分别为：p=p(), Pi,，P

33、i，Pl/d q =q0,q1,，qi,，q(_,则交叉嫡定义为C = Pi log 艮 i=0qi交叉嫡是评价两幅影像差别的重要指标，它直接反映了两幅影像对应像素的差异。 对融合影像前的源影像和融合结果影像求交叉嫡 ,即可得到融合影像与源影像的差异。差异越小， 则该融合方法从源影像提取的信息量越多。互信息(MI)互信息是信息论中的一个重要概念，它可作为两个变量之间相关性的量度，或一个变量包含另一个变量的信息量的量度。F与A, B的互信息分别表示为 MI FA和MI FB :LMI FA =、k=0L10g27L JMI fb : 1 k =0L 1pFB (k, I)PFB(k,i)log2

34、 -口pF(k)pB(i)式中：PA, Pb和PF分别是a, B和C的概率密度；pFA(k,i)和PFB(k,i)分别代表两组影像的联合概率密度。互信息的值越大，表示融合影像从源影像中获取的信息越丰富，融合效果越好。通过定量评价分析对 2种融合方法进行评价。评价内容包括亮度信息、清晰度、光谱信息和信息量等4个方面。试验证明这些评价参数能够很好地反映遥感影像融合效果。融合结果参数评价指标见表所示。亮度信息清晰度均值力左标准差平均梯度空间频率b873.0791488.295122.09743.46829.0288Stack90.48971488.73338.58417.645620.0627PC9

35、0.67711464.02738.26266.107315.398GS90.53151474.09938.39406.134415.4558光谱信息信息量扭曲程度相关系数信息嫡交叉嫡互信息b85.5133Stack6.8639PC0.18750.87336.87130.23680.0409GS0.04190.88586.90860.77940.0405191）亮度信息由表第一列亮度信息可以看出2种融合方法的均值有所降低，方差和标准差有所增加；2）清晰度由表第二列可以看出,3种融合方法的平均梯度和空间频率都小于原始的多光谱影像，说明融合后的结果影像减少了细节纹理信息，使影像的清晰度降低了。3）光

36、谱信息依据表中的第三列可以看出 GS融合扭曲程度最小，而相关系数最大，因此GSI虫合的光谱 信息保持得最好。4）信息量通过信息嫡可以看出 GS融合所得的结果信息量最大；从交叉嫡看到，GS融合和PC融合的结果交叉嫡依次减小，说明融合后的结果影像与源影像对应像素差异也依次减小，即GS融合保留原始图像信息量最多 ，PC融合保留最少；从互信息也可以得到同样规律。而通过定 性评价分析不能得出图像携带的信息量，它仅仅是通过比较分析图像的亮度信息进行比较评价的。4.3.7遥感影像融合定量分析代码Dfusion=imread(C:UsersAdministratorDesktopRS1GS.tif);Dlow

37、=imread(C:UsersAdministratorDesktopRS1Stack_b1-6162-7_CJ.tif);Dhigh=imread(C:UsersAdministratorDesktopRS1b8_CJ.tif);rh,ch=size(Dhigh);均值%扭曲程度%标准差方差%求差方差信息嫡%相关系数Dlowh=imresize(Dlow,rh,ch,bicubic);mean=mean2(Dlow);%meanf=mean2(Dfusion);DIF=abs(meanf-mean);Std=std2(Dfusion);Std2=std2(Dlowh);Ds=StdA2;%D

38、l=Std2A2;D3=abs(Ds-Dl);p=imhist(Dfusion(:),8)/numel(Dfusion(:);r=entropy(Dfusion(:);%c=corr2(Dfusion(:),Dlowh(:);h1=diff（Dfusion）;%求影像差分，反映图像清晰度h=mean2(h1);g1=diff(Dlow);g=mean2(g1);G=h-g;%求融合后影像与原影像差分差值20s=size(size(Dlow); if s(2)=3;% 判断是灰度图还是 RGB f1=rgb2gray(Dlow);f2=rgb2gray(Dfusion); elsef1=Dlow

39、;f2=Dfusion; end G1=double(f1); G2=double(f2); m1,n1=size(G1); m2,n2=size(G2); m2=m1; n2=n1;X1=zeros(1,256); X2=zeros(1,256); result=0;%计两图各灰度级像素for i=1:m1for j=1:n1X1(G1(i,j)+1)=X1(G1(i,j)+1)+1;X2(G2(i,j)+1)=X2(G2(i,j)+1)+1; end end for k=1:256P1(k)=X1(k)/(m1*n1); P2(k)=X2(k)/(m1*n1); if(P1(k)=0)&(

40、P2(k)=0) result=P1(k)*10g2(P1(k)/P2(k)+result;end end f=results=size(size(Dlow);%if s(2)=3;% 判断是灰度图还是 RGB a=rgb2gray(Dlow);b=rgb2gray(Dfusion); end a=double(Dlow); b=double(Dfusion); Ma,Na = size(a); Mb,Nb = size(b); M=min(Ma,Mb);交叉嫡互信息21N=min(Na,Nb); % 初始化直方图数组hab = zeros(256,256);ha = zeros(1,256)

41、;hb = zeros(1,256); % 归一化if max(max(a)=min(min(a)a = (a-min(min(a)/(max(max(a)-min(min(a);elsea = zeros(M,N);endif max(max(b)-min(min(b)b = (b-min(min(b)/(max(max(b)-min(min(b);elseb = zeros(M,N);enda = double(int16(a*255)+1;b = double(int16(b*255)+1;嗾计直方图for i=1:Mfor j=1:Nindexx = a(i,j);indexy = b

42、(i,j);hab(indexx,indexy) = hab(indexx,indexy)+1;%联合直方图ha(indexx) = ha(indexx)+1;%a图直方图hb(indexy) = hb(indexy)+1;%b图直方图endend%十算联合信息嫡hsum = sum(sum(hab);index = find(hab=0);p = hab/hsum;Hab = sum(sum(-p(index).*10g(p(index);刈算 a 图信息嫡 hsum = sum(sum(ha);index = find(ha=0); p = ha/hsum;Ha = sum(sum(-p(

43、index).*log(p(index);%十算b图信息嫡hsum = sum(sum(hb);index = find(hb=0);p = hb/hsum;Hb = sum(sum(-p(index).*log(p(index);%十算a和b的互信息mi = Ha+Hb-Hab22%十算a和b的归一化互信息 mil = hab/(Ha+Hb);img = double(img);%平均梯度% Get the size of imgr,c,b = size(img);dx = 1;dy = 1;for k = 1 : bband = img(:,:,k);dzdx,dzdy = gradien

44、t(band,dx,dy);s = sqrt(dzdx ,A 2 + dzdy .A2) ./ 2);g(k) = sum(sum(s) / (r - 1) * (c - 1);endoutval = mean2(g);m,n=size(Dfusion);腕间频率s=size(size(Dfusion);if s(1)=3;%判断是灰度图还是RGBimg=rgb2gray(Dfusion);endimg =double(img);rf=0.00;cf=0.00;for i=1:mfor j=2:nrf=rf+(img(i,j)-img(i,j-1)*(img(i,j)-img(i,j-1);e

45、ndendRF=rf/(m*n);for i=2:mfor j=1:ncf=cf+(img(i,j)-img(i-1,j)*(img(i,j)-img(i-1,j);endendCF=cf/(m*n);MSF=(RF+CF) 1/2);4.4 生成住房密度栅格影像两表的连接在ArcGIS软件中，选择 Open Attribute Table中查看属性数据，其中，没有该图层的面积字段。因此需要进行有陆地面积的表格进行连接，如图 31。23| 1QL4| 1603751Joins arid Relates卜Related Tablesn -;Jain Remov&-1 lihi Create Gr

46、aph.Add Table to Layout二Relate.Remove Relate(s) ： l图31在Join Data数据中选择两表连接的属性分别为STFID和BLKIDFP00 ,如图32。图324.4.2计算房屋密度创建房屋密度字段 hou_den,并进行字段计算，公式为：hou_den=HSE_UNITS*1000000/ALAND00其中，HSE_UNITS为各街区房屋单元数，ALAND00为各街陆地面积，单位平方米。乘 以1000000是为将房屋密度的单位转换为个/平方千米，如图33。三 Sort Ascending Sort Descending Advanced Sor

47、ting.- Summarize,Statistics.Pwwa e-b sc-pt OMRunrMdx：d_20 向_uMM_u6的00 .BiKn尹Rj 九 3g 舞炯 7 tHitatetaOCLUiXUCi：_LB4R ?_LK4fc*a .HTC匚口 工20曲_ ij9W7_ubWkK . Lftj M IW)7. U1C1N4-d ：_ LBCT ? _L4fc*a .rLHCET dJBKL UMHPAWodM .A1AMR 4b3QMaIW7_hrthiKfaOCI.AMATTF Time 罩忖6b便0 541mgi BattField Calculater.,.Calcula

48、te GeometfyMI ISIwwCMttbkMh Cvnaia pfn.hou d rrom wrs由胃比FvletadateiJiToCAO甲TkjCuIImeJ.4I 十 i-TnUcvfk 口 四 丹IToOASE丁ToOe-o-d a-ta ba s e-*i 4.-I oKMLTu Ftu , liLrrA. dTr tea R hi fa.I,Hr|.,b IFF ts-k R jv1I-F-e-atuine to R.arte-rFloat to Raster a*, Pol mt to RLiaste-r图344.4.4 IDW 插值在ArcGIS工具箱中选择工具 Fea

49、ture to Point ,不勾选Inside ,将选择多边形的中心点作为生成的点，如图 35。图35在ArcGIS软件的空间分析模块中，打开IDW插值方法的窗口，对房屋密度的点图进行IDW个 Kriging:-f：. Nlatural NeighkorSplineSpline with Barriers插值，选择字段为 hou_den,输出栅格大小15m,如图36。自Spatial Analyst Tools Conditional) Densityj 由Distance 由“Extrartion1300三个级别。IDW插值房屋密度图和直接栅格化房屋密度图的重分类结果分 别如图38、图39

50、。ResJDW | D-4OORes_Ra0-400X3W图38图39将住房密度栅格影像作为额外的通道与ETM多光谱波段进行叠加在ArcGIS 10软件中，分别将两种方法生成的住房密度栅格影像数据导出，格式为dat ,如图40。26I I争 copyX Hmm息 Op*n Attribute Table.Wfii and餐 Zgnrt T Uer甲 Zoom To Raster IRfsobtonMTbk 加卜 RjngeDatoEdft FtihjnH%)gave As Layer Kk.，/ /rte Ui胖t Package.-廿 PrcjMttei.7 Enport Dsi%. Vie

51、w rem 口吊pion-Epabai HewMweEitpon： Rmter Dag - idExieffiflUhArrUTtfiO MtEf Dewet Mgn(|* HwIf Ofltwr： OOthiOulsuiRjerUgftrttaiNrqrccRQQS.K SpOLims.NdOLam k54wreSmdiFWCflrfi JKMipreHed =t 更dart QbH, tsp, ngM. baHnrr)fwerh-RM国历闺(.ru4.go, WfLM二舐工 VK2、E:酷三不出/出富3赤lb#幅料极需klw 由 t5zQiidCYCL1M0：图40在Arcgis 10 软

52、件中，分别把 ras、Re_ras、idw和Re_idw导出成dat文件。在ENVI 4.7 软件中，分别打开影像和房屋密度栅格图，选择 Basic Tools - Layer Stacking ,如图 41。图41在 Layer Stacking Parameters 窗口中，点击 Import File 将 idw.dat 和 HSV文件导入,设置保存路径和投影坐标（如图42）。图42监督分类（融合方法为HSV波段为5, 4, 3）打开Google Earth影像作为监督分类的参照27使用Envi的Google Earth Bridge将遥感影像和街区矢量图层导入Google Earth

53、,以便在Google Earth寻找参考影像时更快定位。在 ENVI 4.7 软件中，选择 Spectral - SPEAR Tools - Google Earth Bridge43。SPEAR Tools Anaimaly Delcdi&nTarget OeH沁n WizardChange DetectionGoogle Edrlh BridgeSpectral Uhrarie&*RegistrationSpKtral 勾Mrsk Iridependeril Cm时站m鬲图43在Google Earth Bridge 窗口中，选择添加的遥感影像大三下数据处理HSV,如图44,点击NEXT进

54、入下一步,设置遥感影像显示的参数，选择3波段进行彩色显示，如图 45,点击NEXT进入下一步，同样添加矢量图层，如图作业5、4、46,保存路径为大三下作业 数据处理 栅格数据 Google Earth Bridge.kmlFiU1ap trwHUjer.玲 MTlwinhil FsrirnDitpciL Flitct 部m. age，StrrUk Ri#： tia E * M 强血?由*ktl Jtfri 5.Ctl-ir 1 5t 3 1 9Frdiri i 1|.图45图44图46284.6.2 建立兴趣区在ENVI 4.7软件中，分别打开 HSV.img和重分类后的直接栅格化房屋密度图，

55、并为两个窗口建立其同步，如图 47。Link Di印Geographic Link.Geographic LinkZ Profile (Spectrum).ROJThLInterarfive Stretching.-Save mage图47Select 1the Di splays far theGeographic UAkDisplay *1 皿炉 Di知1纾戋On _t打开兴趣区模版，选择 Overlay - Region of Interest,打开 ROI Tool ,如图48。L (R:HSV Sharp R ：ROI Mask (Wap ( OvcrLy Enhance Tools

56、 Wind Annotation.Classification.Contour Lines. Density Sli” Grid Lines.Region of ntErest.H1 Vectors.*i rcm图48在ROI Tool窗口中，我们可以对 ROI Name和Color进行编辑，双击 ROI Name的区域，可以对训练样区进行名字的修改，右击Color可以选择训练样区的颜色，如图 49, windows选才Z Zoom窗口，即只能在 Zoom中选择训练样区；在 #1 Zoom窗口中，我们鼠标左键选择训练样区，右击确定，第二次右击即保存该训练样区，如图50。46/to J Dele

57、te Part |:上山 GrciwROIUIBlueYellowCyanMagentaMaroon图4929图50在ROI Tool窗口中，点击 New Region建立其新的 ROI样区，Goto可以跳转各个训练 样区，就可以对其进行修改，在实验中，我们将分成11类，进行监督分类，如： water,ransportation, forest, commercial, grass, industrial, crop land, fallow, residential-H, residential-M, residential-L，共 11 类，如图 51。图514.6.3训练样区的选择训练样

58、区的选择依据，如图 52, (a) Commercial , (b) industial , ( c) transportation ,(d) high-density residential, ( e ) medium-density residential, (f ) low-density residential30图524.6.4训练样区的评价4.6.4.1分离性的定量分析在 ROI Tool 窗口中，选择 Options - Compute ROI Separability,计算样本的可分离性。如图 53,表示各个样本类型之间的可分离性，用 Jeffries-Matusita, Tr

59、ansformedDivergence参数表示，这两个参数的值在02.0之间，大于1.9说明样本之间可分离性好;小于1.8大于1.4属于合格样本；小于1.4需要重新选择样本；小于 1,考虑将两类样本合成一类样本。图5331R.OI,parabiilHy Calculationvl YvIIqhI J Ml p 3 Ml truisportVhstifcj -lEfi pointsIcLlwS 337 points 44Gli 4 JltfjHilriZ 2TBK.V二id.卓Tdi5由0.243 373 points T3.ys Erieeii 377 p4ins carap_luid Eur

60、pLol 3912 poLnts fall opt FiurrlaS 4006 points f pr#5t i酎曾ua飘我paHl3Atty BLut L825 paintsUTJitr of 11ws iselectied,电力1 All工l电粕口p皿 BDI fil*.CIau All 工te n-D Vi?Uflliier nput ROh而 | ROI.Type Option5 Hdp Saw ROIs.,.Sdlfict Iripul ROIsRfrsore ROIs.Subtet Data via ROkExport ROIs to EVF.a, Export ROIs to S