ENVI去云补丁Haze_Tool使用说明

Haze Tool User Guide 使用方式............................................................................................................................. 2 云层厚度检测..................................................................................................................... 2 HOT13 ........................................................................................................................ 2 操作步骤............................................................................................................. 3 HOT123 ...................................................................................................................... 4 操作步骤............................................................................................................. 4 案例............................................................................................................................. 5 云层厚度完善..................................................................................................................... 8 Haze perfection TM .................................................................................................... 9 操作步骤............................................................................................................. 9 Haze perfection QB................................................................................................... 11 Maskandinterpolation ....................................................................................... 12 案例........................................................................................................... 13 Thresholdandinterpolation ................................................................................ 13 fill_sink ............................................................................................................. 14 flatten_peak....................................................................................................... 15 操作步骤................................................................................................... 15 云层去除........................................................................................................................... 17 Dark Substract........................................................................................................... 17 操作步骤........................................................................................................... 20 Hist match ................................................................................................................. 22 操作步骤........................................................................................................... 22 Cloud Point ............................................................................................................... 22 操作步骤........................................................................................................... 23 案例................................................................................................................... 28 鸣谢................................................................................................................................... 29 引文................................................................................................................................... 29 使用方式 该模块是在 ENVI4.4 中二次开发而成。使用时只需将文件置于 ***\ITT\IDL64\products\envi44\save_add文件夹内，然后运行 ENVI，Basic Tool菜单 中会出现一个 Haze tool按钮。绝对不要修改文件名，否则无效。 Haze tool主要由三部分组成：云层厚度检测（haze detection），云层厚度完善（haze perfection），云层去除（haze removal）。这三部分前后连贯，每一步的结果都会影响 到最终的除云效果，而且在参数的选择上主观性较大，这也是这个模块的缺点即不 够自动化。 有问题联系作者请加 QQ：27126797 或者 e-mail联系：whogamble@qq.com 除云案例请见作者博客：idlandrs.bokee.com 云层厚度检测 一副影像中往往云层厚薄不一，因此不同厚度云覆盖下的地需要恢复的强度 不一；当然，完全遮挡地面的厚云不在考虑范围内。在该模块中云层厚度检测有两 个指数：HOT13和 HOT123。两者都是相对厚度检测，而不是绝对的光学厚度检测， HOT123是作者在 HOT13的基础上的改进。 HOT13 HOT13由加拿大遥感研究中心的 Zhang ying提出的(Zhang et al. 2002) （原称 HOT，13是作者加上去的以示与 HOT123的区别）。根据地物在蓝色（TM1）和红 色（TM3）波段的高度相关性，在特征空间里绝大部分像素分布在晴空线上。云的 存在会使得云下地物的光谱偏离这条晴空线，云越厚，偏离越大。HOT13等于偏移 距离。 晴空线 : b1sinΦ-b3cosΦ-a=0 (1) HOT= b1sinΦ-b3cosΦ-a (2) Φ 是晴空线的倾角, a晴空线的截距, b1，b3 分别是蓝色(TM1)和红色(TM3)波段。 通过人工选择无云区域的 1、3波段回归，得到晴空线。 操作步骤 1、首先打开需要处理的影像，建立无云区域的 ROI作为对照，并保存（之后其 他操作仍然要用到）. 2、Basic Tool->haze tool->haze detection->HOT13 3、弹出对话框选择需要处理的影像，并且在 Spectral subset中选择蓝色和 红色两个波段。 4、选择一个之前就建立并且已经打开的无云区域 ROI，点击 OK即在内存中生 成一个 HOT13 图像。为了减少内存占用，生成的 HOT13 是放大十倍的 int 格式。 HOT123 在很多情况下, 蓝色（TM1）和红色（TM3）高度相关，相关系数大于 0.9。但 是，当地物更加复杂的时候往往就不成立。比如除了植被，同时存在大量的土壤和 水体，水体非常浑浊等等。在 Quickbird 等高分辨率影像中，除柏油和水泥以外的 其他彩色地物也会降低HOT13的准确性。因此作者对其进行了改进，提出了HOT123, 即利用可见光波段（TM1\2\3）提取云层厚度。 HOT123=k1*b1+k2*b2+k3*b3-b K1\k2\k3\b 这 4 个参数的值使得 |Mean_cloud-Mean_clear|/SD_clear 最小。 Mean_cloud和Mean_clear分别是有云区域（选择的区域的云越厚越好）和无云区域 的 HOT123平均值，SD_clear是无云区域的 HOT123差。 HOT123 要满足有云区域和无云区域的厚度差值尽可能大，而又要使无云区域 （背景）的方差尽可能小。 操作步骤 主要步骤和 HOT13基本一致，不过波段选择不是蓝色和红色两个波段，而是蓝 色、绿色、红色三个波段，即 TM图像上的第一第二第三波段。最后一步的不同之 处是，HOT123需要选择两个 ROI：无云区域 ROI和厚云区域 ROI（其中的厚云是 指将地物完全遮挡住的云层，但是如果图像中没有此类厚云，也可以选择相对最厚 的薄云来代替）。无云区域的选择要尽量涵盖各种土地利用类型，厚云区域无此要求。 为了减少内存占用，生成的 HOT123是放大十倍的 int格式。 案例 左边是 431组合的假彩色图，中间是 HOT123，右边是 HOT13。根据作者处理 大量影像后的经验，推荐使用 HOT123。当没有厚云存在，云非常薄的情况下，有 的时候可以选择 HOT13。 云层厚度完善 虽然 HOT13\HOT123尽量地突出云层信息，抑制背景信息，但是仍然有不少无 云区域的云层厚度值偏离零很大，主要发生在水、浪花、土壤、雪、建筑区等地物 上。以 HOT13举例说明，在如下的 13波段特征空间里（横坐标是第一波段，纵坐 标是第三波段），无云区域用灰色表示，厚云区域用黑色表示；无云区域形状如同飞 机。一翼代表偏小的云层厚度检测，一翼代表偏大的云层厚度检测。云层厚度的完 善实际上就是把无云区域的这辆飞机的两翼去掉。 因此，当地物很复杂的时候就需要对 HOT值进行一些额外处理，以修正这些偏 差。当然，云层厚度完善中的几种方法都是非必须的，可以根据影像的特征酌情选 用，也可不用。 Haze perfection QB适用于几乎所有光学遥感影像，如 Quickbird等 4波段影像 （当然也适用于 TM 等更多波段影像），这类影像会有几个可见光波段被云污染。 Haze perfection TM 只适用于 TM等含有中红外波段的影像，中红外波段往往很少 受薄云影响。不管选择哪一个，处理后的结果都需要重新将无云区域 ROI内的 HOT 平均值归零；因此，必须确保在进行云层厚度完善操作之前内存中有无云区域 ROI 的存在，否则最后将不会进行归零操作。 我们将以以下案例举例说明。左边是有云 TM影像，中间是无云的 TM影像（城 市区域），右边是有云的 QUICKBIRD 影像。可见有很多特别暗的背景的云层厚度 是负值，也有很多特别亮而不是云的背景。 Haze perfection TM 这种方法的灵感来自 Liang shunlin的聚类匹配方法(Liang et al. 2001; Liang et al. 2002)。作者稍作修改用于 HOT值的完善上。该方法以利用 TM4\5\7波段（忽略这 些波段的云污染）进行非监督分类的结果作为输入参数，通过计算每类地物的平均 HOT值并将其减去以达到归零的目的。 操作步骤 1、检查 4\5\7波段是否受云的影响，如果影响可以忽略，则利用 TM4\5\7波段 进行非监督分类。类的数量看地物复杂程度随意决定，推荐 20-50类。 2、选择需要处理的 HOT影像（必须有无云区域的 ROI）。 3、选择对应的非监督分类图 4、选择无云区域 ROI 5、决定是否要对某些方差特别大的类型进行插值操作。如果不需要，则选择 Cancel；如果需要，则在相对应的类型前面打勾，然后点击 OK。结果会在 内存中生成。 以下是案例的结果，这步操作的效果取决于地物类型的光谱特征。当存在某几种地 物类型，其 4\5\7波段相似，而 1\2\3可见光波段相差很大的时候，这一步操作就会 对这几种地物类型失效。 Haze perfection QB 点击 Haze Perfection QB，会弹出对话框要求选择要处理的 HOT图像。然后， 弹出一个多项选择的对话框，选择合适处理方法。可以选择一种，也可以同时选择 多种，作者不推荐同时选择多种方法处理。主要有四种方法可供选择： maskandinterpolation、thresholdandinterpolation、fill_sink、flatten_peak。其中 fill_sink 和 flatten_peak非常吃内存，是 haze tool最大能处理的影像的限制步骤。 Maskandinterpolation 选择一个手动勾勒的需要掩膜并插值的异常区域 ROI进行插值。该操作是人工 方法，操作简单，当异常较少的时候比其他 haze perfection方法要快。 运行完毕后，将跳出对话框要求选择无云区域 ROI进行归零操作以纠正改变后 的 HOT带来的偏差。 案例 如下，左边有个多边形区域被手工勾勒并且插值，右边是原始图像。 Thresholdandinterpolation 通过设定一个阈值范围[minimum, maximum]，将范围之外的 HOT值掩模并进行 插值。插值方式采用周围像素取平均的方式。这种方法较简单，因为 HOT 值偏小 的地物受到云的影响依然可以有较大的 HOT 值。因此，这种方法现在作者很少采 用，仅仅作为一种异常值去除的方法（去除厚云，和一些极端负值）。 运算需要输入 4个参数：minimum(最小值),maximum（最大值）,kernel size for interpolation (插值窗口大小)，kernel size for smooth (平滑窗口大小)。前三个参数用 于插值，最后一个参数用于结果图像的平滑滤波。默认情况下，minimum=整幅影像 最小值，maximum=整幅影像最大值，kernel size for interpolation=7, 不进行平滑滤 波。 运行完毕后，将跳出对话框要求选择无云区域 ROI进行归零操作以纠正改变后 的 HOT带来的偏差。 fill_sink 该方法来自于水文领域中对地形图的处理(Planchon and Darboux 2002)：填洼算 法。我们把 HOT图像看作地形图，HOT值偏小的地物往往会在地形图上形成一个 一个的洼地。洼地虽然比周围的值要小，一个洼地可以在海拔很低的地方，也可以 在一座高山上。这也是作者很少利用之前的 maskandinterpolation的原因。 运行完毕后，将跳出对话框要求选择无云区域 ROI进行归零操作以纠正改变后 的 HOT带来的偏差。 案例 如下图，左边是经过 fill_sink处理的 HOT，右边是未经过处理的 HOT。很明显， fill_sink对暗斑去除很有效。 flatten_peak 除了 HOT值偏小的地物，还有 HOT值偏大的地物。一个像素的 HOT值很大， 有可能是该像素被云覆盖，也有可能是本身 HOT 值偏大。作者利用云层厚度具有 渐变的特征而地物边界骤变的特征将这两者区别开来。 操作步骤 1、预处理。在弹出对话框中选择 morph reconstruction。并在接下去弹出的参数 对话框中输入一个合适整数 n, 进行 n次形态学侵蚀操作。n 的大小要保证绝大多数 偏大的地物能够被侵蚀完全。过大，则徒耗计算量。 预处理步骤在内存中生成两个图像 max_change_of_successive_erosion和 total_change_after_reconstruction。这两个图像作为输入图像用于下一步真正的 flatten_peak操作。 2、在预处理完后再重新 flatten_peak，不过这回在弹出的对话框中选择 choose previous results 。 并 在 接 下 来 的 对 话 框 中 依 次 选 择 预 处 理 中 生 成 的 max_change_of_successive_erosion和 total_change_after_reconstruction。 3、对 max_change_of_successive_erosion和 total_change_after_reconstruction分别 设定最大阈值 maxchange_limit 和 totalchange_limit。将同时大于这两阈值的区域进 行掩模并插值。Kernel size for mask interpolation和 kernel size for smooth和之前提到 的意义相同。 4、运行完毕后，将跳出对话框要求选择无云区域 ROI进行归零操作以纠正改变 后的 HOT带来的偏差。 当然，haze perfection操作的目的尽可能地突出云的信息，抑制背景信息。这里提供 的几种方法的最终结果并不一定完美，有待将来开发更多更好的算法。 云层去除 完成了云层厚度检测和云层厚度完善之后即可进行云层去除。这里提供了三种 方法：Dark Substract (暗物质法), Hist match（直方图匹配法）, Cloud Point （云点 法）。Cloud Point 是在暗物质法的基础上的改进。以上的几种方法都是按照云层厚 度对图像进行分层除云处理。 Dark Substract Dark substract即 HOT的提出这 Zhang ying采用的方法(Zhang et al. 2002)，主要 是对 HOT图像分层分割，每个像素各个波段的 DN值减去该像素所在层各个波段的 最小值。 然后考虑到最小值不够稳定，因此作者采用 Dal Moro推荐的固定百分点方法代 替最小值作为下限(Dal Moro and Halounova 2007)。 作者不推荐使用这种除云方法，因为云层不仅能够增强较暗地物的亮度，也能降低 较亮地物的亮度，从而缩小云下地物的对比度。云的模糊效果不仅仅是加性的，还有乘 性的效果。而且云的厚度应该是渐变的，人为将其按照厚度分割成一层层会导致结果的 分块效应。 操作步骤 1、确保无云区域 ROI（必须）以及需要处理的区域 ROI（可选）已经在内存中。 然后选择需要除云的影像，并且在 Spectral subset 中选择需要除云的波段。作者推 荐处理可见光波段和近红外波段。 2、选择经过云层厚度检测和云层厚度完善（可选）而得到的 HOT图像。 3、选择无云区域 ROI。 4、选择需要进行除云处理的区域，如果需要对整幅影像进行处理则选择 Cancel。 5、对图像云层进行切割分别处理。[minimum,maximum]范围之外将不被处理， 默认情况下，minimum=最小 HOT，maximum=最大 HOT。Increment表示切割云层 的厚度，默认情况下 increment=1，太大则有明显的分块效应，太小则每层的像素数 量太少而结果不够稳定。Lowerandupperlimit是作为下限的固定百分点，默认情况下 等于 2%，0%则表示最小值作为下限。 6、处理结果自动输出到内存，并附带一张各层下限的对应波段的最小值曲线。 曲线波动越小效果越好。如下图，说明上一步定义的 minimum 不合适，minumum 至少要大于-4。 Hist match Hist match是分层直方图匹配方法，即将 HOT图像分层分割，然后以无云区域 ROI作为参照对每一层进行直方图匹配操作。这种方法在地物分布较均匀的情况下 效果较好。当然，结果也会由于分层的原因导致分块效应。 操作步骤 与 Dark substract的操作步骤雷同。 Cloud Point Cloud point是 dark substract的改进，相同之处是仍然要对 HOT图像进行分层分 割，不同之处是每层同时提取上限和下限两个值。通过上限和下限分别回归得到两 条回归直线的交点即云点（代表最厚的云），如下图（蓝色波段为例）。任何一个像 素点都将通过云点进行点投影到 HOT=0的竖直直线上实现除云处理。 a是虚拟的云点，是各层的上限回归直线和下限回归直线的焦点。b 点是除云前 的像素，c是除云后的像素，c点的 HOT值等于 0. 操作步骤 前两步与 Dark substract的操作步骤雷同。 3、选择有云区域 ROI用于计算云点，ROI的选择直接决定了最终的除云质量。 该 ROI的建立的目的在于得到稳定的云点，因此要求上限点和下限点稳定而且回归 直线的回归系数尽可能的大。所以手工建立的有云区域要求：1、不同厚度的云覆盖， 该区域内的 HOT值范围尽可能大，使有更多的下限和上线点用于回归。2、区域内 的地物尽量分布均匀，如此可以减少地物造成的影响。如果选择 Cancel，表明利用 整幅影像计算云点，不仅计算量大而且效果不如人工选择更合适的区域。 4、对图像云层进行切割分别处理。[minimum,maximum]范围之外将不被处理， 默认情况下，minimum=最小 HOT，maximum=最大 HOT。Increment表示切割云层 的厚度，默认情况下 increment=1，。Lowerandupperlimit是作为下限的固定百分点， 默认情况下等于 2%，那么上限的固定百分点是 98%。0%则表示最小值作为下限， 最大值作为上限。 确定后，会输出各层的结果图： 各层的像素个数 各层的下限值 各层的上限值 各层的上限值减去下限值 6、根据上一步的结果选择用于计算云点的 HOT范围。1、要求范围内各层的像 素不能太少（根据各层的像素个数图），2、要求各层的上限值减去下限值随 着 HOT的增加逐渐稳定地减少（根据各层的上限值减去下限值图）。如上面 的例子，选择范围最小值不能是-3以下，最大值不能超过 15。 以下案例说明了第 3步，有云区域 ROI选择的重要性。下图的有云区域的 ROI 选择比较合适。 下图的有云区域 ROI选择不合适。造成的错误是：蓝色、绿色、红色波段随着 HOT的增加对比度增加，这个显然与现实不符。 案例 Cloud point的结果 Dark substract的结果 鸣谢 南京师范大学蒋圣同学在模块测试中的无私帮助。 引文 Dal Moro, G., & Halounova, L. (2007). Haze removal for high-resolution satellite data: A case study. International Journal of Remote Sensing, 28, 2187-2205 Liang, S.L., Fang, H.L., & Chen, M.Z. (2001). Atmospheric correction of landsat ETM+ land surface imagery - Part I: Methods. Ieee Transactions on Geoscience and Remote Sensing, 39, 2490-2498 Liang, S.L., Fang, H.L., Morisette, J.T., Chen, M.Z., Shuey, C.J., Walthall, C.L., & Daughtry, C.S.T. (2002). Atmospheric correction of landsat ETM plus land surface imagery - Part II: Validation and applications. Ieee Transactions on Geoscience and Remote Sensing, 40, 2736-2746 Planchon, O., & Darboux, F. (2002). A fast, simple and versatile algorithm to fill the depressions of digital elevation models. Catena, 46, 159-176 Zhang, Y., Guindon, B., & Cihlar, J. (2002). An image transform to characterize and compensate for spatial variations in thin cloud contamination of Landsat images. Remote Sensing of Environment, 82, 173-187