ENVI实验报告

实验目的：基于环境小卫星的草原荒漠化监测1、专题讲述浑善达克地区位于内蒙古草原阴山北麓锡林郭勒盟高原中部，是亚洲草原荒漠化土地东部边缘区的重要组成部分，经纬度在东经1145511638，北纬41464307之间，平均海拔高度在1100米左右。退化去属温带半干旱区、中温带干旱大陆性季风气候，浑善达克退化土地多为沙地，植被稀少，特别是春季地面回暖解冻，地表裸露，多细沙土，狂风起时沙尘弥漫，形成沙尘天气。近年来频频发生在京津地区春季的沙尘天气与该地区生态环境恶化，人地关系严重失调有关。据统计，京津地区沙尘暴70的沙源来自这个区域。根据该区域不同时期的植被覆盖数据可以实现该区域的植被变化监测。本专题详细介绍了利用环境小卫星CCD-1A 图像反演植被覆盖图的完整流程，专题涉及环境小卫星的数据读取、辐射定标、大气校正、植被覆盖反演模型的建立、遥感反演过程、植被覆盖变化监测等内容。除了使用ENVI 主模块功能外，还需要用到FLAASH 大气校正扩展模块、IDL 开发的环境小卫星数据读取补丁、波段运算等功能。2、处理流程介绍根据环境小卫星CCD 数据特点及草原植被变化监测的要求，采用的技术路线为：先对环境小卫星CCD 数据进行数据预处理：数据读取、辐射定标、大气校正、研究区裁剪，建立反演模型，利用波段运算工具，反演出整个浑善达克地区的归一化植被指数、植被覆盖度，根据不同时期的植被覆盖数据，实现草原土地退化的遥感监测。流程说明：1） 图像获取本文的基础数据包括2006年土地覆盖图数据、环境小卫星CCD-1A 数据。环境小卫星CCD-1A 数据可以在环保部卫星环境应用中心免费下载获取，土地覆盖图数据的制作时间是2006年8月份，本专题选择了影像质量良好的2009年8月11日的环境卫星数据作为数据源。环境小卫星CCD 数据下载地址：2） 数据读取和定标网上免费获取的HJ-1A 卫星CCD 的分发格式为Geotiff，每一个波段为一个Geotiff 文件，并提供一个元数据说明（.XML）。可以依次打开每个波段，并用元数据说明中提供的定标参数用Band Math 工具进行手动定标，再用Layer stcaking 功能将4个波段的图像合成为一个多波段的数据文件。本专题中用的是环境小卫星的读取补丁，直接读取CCD 数据，输出结果为一个多波段的ENVI 标准栅格文件，并带有中心波长等信息，可直接输出定标结果（辐射亮度）。补丁下载地址为：/ESRI/viewthread.php?tid=755753） 图像配准获取的2006年的土地利用分类图已经过精确的地理定位，以该景作为进准影像，对环境小卫星数据进行配准，便于后续的变化监测。4） 大气校正对于环境小卫星数据，提供光谱响应函数，可用FLAASH 大气校正模块中进行大气校正，波谱响应函数下载地址：http//n16/n1115/n1522/n2118/index.html。5）反演模型构建及模型应用6）植被变化监测根据2006年8月和2009年8月的植被覆盖数据，采用波段运算，实现草原植被变化的遥感监测。3、详细处理过程本实验的数据存放在“草原植被变化遥感监测”文件夹内。3.1数据预处理第一步：安装环境小卫星数据处理补丁将ENVI_HJ1A1B_Tools.sav 补丁放在homeITTIDL80productsenvi48save_add 目录下。第二步：数据读取和定标主菜单 - File - Open External File -HJ-1A/1B Tools，打开环境小卫星数据读取补丁。在HJ-1A/1B Tools，打开环境小卫星数据读取补丁。在HJ-1A/1B Tools v3.0 面板中，选择CCD，点击Input File 输入“1-环境小卫星数据HJ1A-CCD2-2-64-20090811-L20000154793154793”文件夹中的.xml 文件中，点击Output Path设置数据的输出路径，勾选“Calibration”“Layer Stacking”两个选项，单击Apply 按钮。第三步：工程区裁剪由于整景影像数据范围非常大，本专题工程区只是其中一小部分，在进行大气校正之前，先将浑善达克以及周边区域裁剪出来。（1）打开上一步生成的文件：HJ1A-CCD2-2-64-20090811-L20000154793_Calbrated_LayerStacking.img，（2）主菜单 - File - Save File As - ENVI Standard，弹出New File Builder 面板，（3）在New File Builder 面板中，单击Import File，弹出的Create New Input File面板，（4）在Create New File Input File面板中，选中Select Input File 列表中的裁剪数据，单击Spatial Subset 按钮，（5）在Select Spatial Subset面板中，单击Image，弹出Subset by Image 对话框，（6）在Subset by Image 对话框中，按住鼠标左键拖动图像中的红色矩形框确定裁剪区域，裁剪出包括浑善达克区域的一部分，单击OK，（7）在Select Spatial Subset面板中，可以看到裁剪区域信息，单击OK，（8）在Create New File Input File对话框中，单击OK，（9）在New File Builder，单击Choose 设置输出文件名20090811-Cal-sub.img及路径，单击OK。完成浑善达克及周边区域的裁剪。第四步：图像配准下面以土地利用图作为基准影像对环境小卫星图像进行图像配准。（1）分别打开和显示基准影像“基于环境小卫星的草原植被变化监测浑善达克2006年8月土地利用分类图.img”。（2）主菜单 - Map - Registration - Select GCPs：Image to Image，打开几何校正模块。（3）选择显示2006年土地利用分类图文件的Display为基准影像（Basic Image），显示环境星文件的Display为待校正影像（Warp Image），点击OK进行采集地面控制点。（4）打开Tools - Link - Geographic link，将两个窗口都选择为on，单击确定，找到定位的大致区域后，再Tools - Link - Geographic link，改为off，关闭链接。（5）在两个Display中找到相同区域，在Zoom窗口中，点击左下角第三个按钮，打开定位十字光标，将十字光标到相同点上，点击Ground Control Points Selection上的Add Point按钮，将当前找到的点加入控制点列表。（6）用同样的方法继续寻找其余的点，当选择控制点的数量达到3时，RMS被自动计算。Ground Control Points Selection上的Predict按钮可用，选择Option Auto Predict，打开自动预测可能。这时在Basic Image上面定位点，Warp Image上会自动预测区域。（7）完成控制点的选择，RMS值小于1个像素，点击Ground Control Points Selection上的File Save Coefficients to ASCII，将控制点保存。（8）在Ground Control Points Selection上选择Option - Warp File（as Image Map），选择校正文件（HJ数据文件）。（9）在校正参数面板中，默认投影参数和像元大小与基准影像一致，30米。（10）重采样选择Nearest Neighor，背景值（Background）为0 。（11）Output Image Extent：默认是根据基准图像大小计算，可以做适当的调整。（12）选择输出路径和文件名，单击OK按钮。第五步：大气校正环境小卫星提供了波普响应函数，以文本形式提供，第一列表示波长（nm），后面四列分别表示4个波段对应波长的波普响应值。需要制作波普曲线来描述波普响应函数，用于大气校正。1） 制作波普曲线（1）主菜单Window Start New Plot Window，打开ENVI Plot Window 面板，在波普绘制窗口中，选择File Input Data ASCII，如图所示，自动将第一列作为x轴，后面4列作为y轴，波长单位选择Nanometers，单击OK。（2）如图，在绘制窗口生成了4条曲线，选择Edit Data Parameters，编辑每条线的名称为b1，b2，b3，b4，便于区分。（3）选择File - Save Plot As -Spectral Library，在Output Plots to Spectral Library面板中，单击Select All Items，单击OK。（4）在Output Spectral Library 面板中，有输出曲线相关参数设置，按默认，选择输出路径和文件名，单击OK，将波普曲线保存为波普库文件：环境1A星CCD2光谱响应.sli 。2）FLAASH大气校正（1）准备数据FLAASH对图像文件有以下几个要求：数据是经过定标后的辐射亮度（辐射率）数据，单位是：（W）/（cm*nm*sr）。数据带有中心波长（wavelength）值，如果是高光谱还必须有波段宽度（FWHM），这两个参数都可以通过编辑头文件信息输入（Edit Header）。数据类型：支持四种数据类型：浮点型（floating）、长整型（long integer）、整型（integer）和无符号型（unsigned int）。数据存储类型：ENVI标准栅格格式文件，且是BIP或者BIL。波普范围：4002500nm。本次用的环境小卫星经过以上处理，已经定标为W*m(-2)*sr(-1)*um(-1)单位、浮点型的辐射率数据，有中心波长信息，下面将BSQ格式转化成BIL格式。选择主菜单Basic Tools - Convert Data（BSQ、BIL、BIP），选择已经定标和配准的数据20090811-cal-jz.img，在Convert File Parameters中，Output Interleave选择BIL，选择Convert In Place：yes，单击OK。（2）设置参数进行FLAASH大气校正主菜单Spectral - FLAASH 打开FLAASH 大气校正模块；点击Input Radiance Image，选择BIL格式的环境小卫星数据20090811-cal-jz.img，在Radiance Scale Factors面板中选择Use single scale factor for all bands，由于定标的辐射量数据与FLAASH的辐射亮度的单位相差10倍，所以在此Single scale factor选择默认：10【要求实验报告中解释为什么是10】，单击OK；设置输出文件及路径设置传感器基本信息设置：成像中心点经纬度、传感高度、成像区域平均高度、成像时间设置，这些都可以从数据头文件中读取HJ1A-CCD2-2-64-20090811-L20000154793.XML。大气模型，选择MLS，气溶胶模型，选择Rural，气溶胶反演方法选择None，能见度设置为40km。参数如下图单击Multispectral Setting按钮，在Filer Function File导入之前做好的光谱响应曲线“环境1A星CCD2光谱响应.sli”，单击OK。单击Advanced Setting按钮，在高级设置中，Tile Size默认的是Cash Size的大小，手动改为100Mb，单击OK。设置好后，在大气校正模块面板中，单击Apply。大气校正完成后，检查大气校正的结果，分别加载校正前后的图像，将两幅影像进行地理链接，移动到植被区域（植被的波普曲线比较特殊），在影像上右键，选择Z Profile（Spectrum）打开光谱曲线窗口，显示两幅图像同一位置的光谱曲线图，如下图。第六步：裁剪浑善达克区（1）显示大气校正后的图像20090811-cal-jz-FLAASH.img，在image窗口选择Overlay - Vectors，打开Vector Parameters面板，选择File - Open Vector File，打开hunshandake.evf。（2）在Available Vectors List 面板中选择该矢量文件，点击Load Selected，选择显示图像的Display，单击OK，矢量叠加在影像上。（3）在Available Vectors List 面板中，选择File -Export Layers to ROI，在Select Data File to Associate with new ROI中，选择Convert all records of an EVF layer to one ROI，单击OK，将矢量转化为一个ROI。（4）在图像窗口，选择Overlay - Region of Interest，打开ROI面板，浑善达克ROI显示在图像上，在ROI Tool面板中，选择File - Subset Data via ROIs，在Select Input File to Subset via ROI面板中，选择20090811-cal-jz-FLAASH.img，Mask pixels outside of ROI选择Yes。（5）在Spatial Subset via ROI Parameters面板中选择浑善达克ROI，Mask pixels outside of ROI选择Yes。（6）设置输出路径及文件名20090811-yanjiuqu.img，单击OK。到此，数据预处理工作已经全部完成，下面介绍如何应用ENVI进行植被覆盖度反演和植被变化监测。3.2植被覆盖度反演第一步：计算归一化植被指数（1）选择主菜单 - File - Transform - NDVI计算模块。（2）选择裁剪后的数据：20090811-yanjiuqu.img。（3）在NDVI Calculation Parameters面板中进行相应设置，环境星波段参数与TM数据相似，应用被植被强吸收的红光波段（环境星第3波段）和被植被强反射的近红外波段（环境星第4波段）计算归一化植被指数。数据类型默认为浮点型。（4）选择文件保存名（20090811-yanjiuqu-NDVI）和路径，单击OK执行。第二步：计算植被覆盖度（1）计算植被覆盖度FC采用的是混合像元分解法，将正景影像的地类大致分为水体、植被和建筑，具体的计算公式如下：FC=(NDVI-NDVIsoil)/(NDVIveg-NDVIsoil)其中，NDVI为归一化植被指数，NDVI为完全是裸土或无植被覆盖区域的NDVIsoil值，NDVIveg则代表完全被植被所覆盖的像元的NDVI值，即纯植被像元的NDVI值。取经验值NDVIveg=0.70和NDVIsoil=0.00，且有当某个像元的NDVI值大于0.70时，FC取值为1；当NDVI小于0.00时，FC取值为0。利用ENVI主菜单- Basic Tools - Band Math，在公式输入栏中输入：（b1 gt 0.7）*1+(b1 1t 0.0)*0+(b1 ge 0 and b1 le 0.7)*(b1-0.0)/(0.7-0.0)b1：选择ENVI图像选择文件保存名（20090811-yanjiuqu-VCI）和路径，单击OK执行，得到植被覆盖度图像。3.3两图像植被覆盖区域提取第一步：植被覆盖区提取土地利用分类类型影像DN值耕地1林地2草地3水域4城市、工矿、居民用地5未利用地6沙地61（1）2009年8月植被覆盖区提取利用ENVI主菜单 Basic Tools - Band Math，在公式输入栏中输入（b1 1e 0.3）*0+（b1 gt 0.3）*1b1：选择HJ1A-CCD2-20090811-cal-sub-flaash-ref-hunshandake-VCI图像，0.3为经验值。选择文件保存名（2009年8月植被覆盖区）和路径，单击OK执行，得到2009年8月植被盖度区。（2）2006年8月植被覆盖区提取获取的2006年8月标准的土地利用分类图，分类代码与地物类型如下表所示，耕地、林地、草地，属于植被覆盖区域。ENVI主菜单- Basic Tools - Band Math，在公式输入栏中输入：（bl