0107] 利用以上DEM数据中标记阴影区域可以排除水体单元中的阴影。因为这些阴影区 域往往是陡峭的斜坡，因此排除大范围湖水等静态水体存在的可能，而对于流动的水体如 瀑布等，其斜坡顶部非阴影区域一定存在河流等水体。因此，在本发明中完全被标记的阴影 区域包含在内的水体单元将被视为阴影而排除，这样就大大降低了水体提取过程中错检的 概率。至此完成本发明目标水体信息的提取。
[0108] 本发明未详细阐述部分属于本领域公知技术。
【主权项】
1. 一种遥感卫星多光谱影像的水体信息自动提取方法，其特征在于：包括水体指数提 取阶段、水体指数阈值分割阶段和山体阴影剔除阶段， 所述水体指数提取阶段步骤如下： (1) 选取一景无云的国产卫星多光谱遥感影像，根据该遥感影像的绿光波段的灰度值 计算表观反射率值P ，根据近红外波段的灰度值计算表观反射率Pnir，计算过程如下： 先按照如下公式将影像的灰度值转换为表观辐亮度La: La=GainXDN+Bias， 式中Gain为增益，DN为灰度值，Bias为偏移量； 然后，将表观辐射亮度La转换为表观反射率Pa，公式如下：式中d是日-地距离订正因子，Es是大气外太阳光谱辐照度，0s是太阳天顶角，这些 参数可从影像头文件和中国资源卫星应用中心网站获得； (2) 根据步骤（1)中计算出的遥感影像的绿光波段的表观反射率P 和近红外波段 的表观反射率Pnir获取归一化水体指数影像NDWI，公式如下：该归一化水体指数影像包括水体与非水体目标，水体与非水体目标一般在水体指数影 像直方图上呈双峰分布；由于归一化水体指数影像NDWI是范围为[_1，1]的离散浮点数，对 归一化水体指数影像NDWI乘以100并四舍五入，得到水体指数影像，该水体指数影像像元 值是范围为[-100,100]的整数，至此完成水体指数提取阶段； 所述水体指数阈值分割阶段步骤如下： (3) 将步骤⑵中水体指数影像中像元的值从小到大依次排列，得到水体指数影像中 像元的值的集合； (4) 设定水体指数影像的阈值T。，该T。为步骤（3)水体指数影像中像元的值的集合中 的第i个值，小于阈值T。的值为非水体目标，大于阈值T。的值为水体目标；i= 1…P，P为 大于1的正整数；水体目标的像元数占水体指数影像的所有像元数的比例为所有水体 目标的像元的平均值为非水体目标的像元数占水体指数影像的所有像元数的比例为 ?:，所有非水体目标的像元的平均值为y:，水体指数影像的平均像元值为y，水体与非水 体目标间的方差记为g，i从1到N，依次计算水体指数影像的平均像元值y，公式如下 \i =Q0X U0+Q1X U1 然后计算水体与非水体目标间的方差记为g，公式如下： g= OQXCyo-yy+OiXCyfy)2 记录每组阈值T。对应的水体与非水体目标间的方差g，将最大的水体与非水体目标间 的方差g对应的阈值T。，作为全局阈值T; (5) 将水体指数影像中大于等于步骤（4)的全局阈值T的像元设置为水体目标，将水体 指数影像中小于步骤（4)的全局阈值T的像元设置为非水体目标；去除孤立的单个水体目 标，将所有相邻接的水体目标设置为一个水体单元； 对每个水体单元的边界各向外每次扩展一个像元，直到扩展的像元数之和与水体单元 像元数的差最小，此时扩展的水体单元作为疑似水体单元； (6) 将疑似水体单元中像元的值从小到大依次排列，得到疑似水体单元NDWI中像元的 值的集合； (7) 设定疑似水体单元的阈值Tj,该Tj为步骤（6)疑似水体单元中像元的值的集合中 的第j个值，小于阈值Tj的值为非水体目标，大于阈值Tj的值为水体目标；j= 1…L，L为 大于1的正整数；水体目标的像元数占疑似水体单元的所有像元数的比例为《j。，所有水体 目标的像元的平均值为yi。，非水体目标的像元数占疑似水体单元的所有像元数的比例为 所有非水体目标的像元的平均值为為V，疑似水体单元的平均像元值为h，水体与非 水体目标间的方差记为gj，j从1到L，依次计算疑似水体单元的平均像元值Iij，公式如下 y j= w j0X y j0+o nX y n 然后计算水体与非水体目标间的方差记为gj，公式如下： Sj= w j〇_V-j)2+ ? ji X ( y J1-v-j)2 记录每组阈值Tj与水体与非水体目标间的方差gj，将最大的水体与非水体目标间的方 差gj对应的阈值Tj，作为局部阈值Th; (8) 将疑似水体单元中大于等于步骤（7)的局部阈值1；的像元设置为新的疑似水体单 元，再对该新的疑似水体单元的边界向外每次扩展一个像元，直到扩展的像元数之和与新 的疑似水体单元像元数的差最小，此时扩展的新的水体单元再次作为疑似水体单元； (9) 将步骤（8)得到的疑似水体单元，代替步骤(6)中的疑似水体单元，重复进行步骤 (6)至步骤（8)，直至步骤（8)得到的新的疑似水体单元与扩展之前的疑似水体单元像元数 之差小于扩展之前的疑似水体单元像元数的1/100为止，此时，设此新的疑似水体单元为 最终疑似水体单元，该最终疑似水体单元可能包含山体阴影； 所述山体阴影剔除阶段步骤如下： (10) 将步骤（9)中新的疑似水体单元在国产卫星多光谱遥感影像相同地理范围的DEM 数据中标记出来，形成标记的疑似水体信息，DEM为国产卫星多光谱遥感影像范围内的地面 高程数据； (11) 从步骤（1)中的一景无云的国产卫星多光谱遥感影像的元数据中获取太阳方位 角L，即太阳光线在地平面上的投影与当地子午线的夹角； (12) 计算DEM数据中每一个像元点对给定太阳方位角L方向的坡度 (13) 从步骤（1)中的一景无云的国产卫星多光谱遥感影像的元数据中获取太阳高度 角0，将DEM数据中每一个像元的太阳方位角L方向的坡度^与太阳高度角0进行比较， 如果(6)0，则DEM数据中的此像元以及太阳光线经过此像元在地平面上的投影为阴影区 ±或，将DEM数据中所有阴影区域标记为疑似阴影信息； (14) 在DEM数据中，对比步骤（10)所有标记的疑似水体信息与步骤（13)的标记的疑 似阴影信息，如果标记的疑似水体信息被完全包含在标记的疑似阴影信息里，则视此标记 的疑似水体信息为阴影信息，则从步骤（9)中最终疑似水体单元中去除对应的阴影信息， 最终完成水体信息提取。2.根据权利要求1所述的一种遥感卫星多光谱影像的水体信息自动提取方法，其特征 在于：所述步骤（4)中在国产遥感卫星影像中归一化水体指数影像NDWI的水体与非水体的 初始值阈值一般为0. 1，因此此处设置T。的初始值为10。3.根据权利要求1所述的一种遥感卫星多光谱影像的水体信息自动提取方法，其特征 在于：所述步骤（12)中计算DEM数据中每一个像元点对给定太阳方位角L方向的坡度0M 公式如下：式中巾SL与北方向顺时针方向夹角，范围取[0°，360° ] 是反映单位长度上高 程的升降值，即某点梯度方向的坡度，P表达式为：式中，坡度P范围取[-90°，90° ]，fxSX轴方向的地形曲线函数，fy*Y轴方向的 地形曲线函数，使用三阶反距离平方权差分的计算表达式表示，公式为：式中，1为DEM-个3X3像元窗口中某像元点的高程，i为此3X3像元窗口各像元点 的编号，i= 1，2,…9,r为DEM像元的分辨率。
【专利摘要】一种遥感卫星多光谱影像的水体信息自动提取方法如下：(1)对国产遥感卫星多光谱影像的表观反射率数据进行归一化水体指数提取，形成归一化水体指数影像NDWI；(2)对归一化水体指数影像进行基于最大类间方差算法的全局阈值分割初步提取水体信息；(3)对初步的水体信息进行局部缓冲区扩张和阈值分割，获取稳定的水体单元结果集合；(4)根据高程数据获取坡度和阴影信息并对水体单元结果集合进行阴影移除，最终完成水体信息提取。本发明采用计算机自动解译的方法，利用从全局分割到局部迭代分割的思想，对归一化水体指数流程化地进行水体提取，提取精度高，工程容易。
【IPC分类】G06F19/00, G01S7/48, G01S17/89
【公开号】CN105046087
【申请号】CN201510472553
【发明人】吕争, 傅俏燕, 王小燕, 隋正伟, 高青山
【申请人】中国资源卫星应用中心
【公开日】2015年11月11日
【申请日】2015年8月4日