ENVI 5.3 + 哨兵2号（Sentinel-2）L2A提取健康水体和不健康水体

大学是学GIS的，最近因为需要而搞起了遥感影像，不过GIS和RS不分家嘛！多了解学习没有坏处！针对学习过程中遇到的问题和获得的收获，特发文一篇以作记录，积累经验。

【1.数据下载】

首先需要下载数据，参考的文章使用的是哨兵2号卫星的遥感影像，所以我也使用相同的卫星影像。在欧空局（ESA）网站上（Sentinel-2 - Missions - Sentinel Online - Sentinel Online）可以下载到哨兵系列卫星的数据，还是比较方便的，从【Data Access】-【Access to sentinel data】可以进入影像下载页面。

Level-1C 是经过正射校正和几何精校正的大气表观反射率产品，并没有进行大气校正。

Level-2A 主要包含经过大气校正和辐射定标的大气底层反射率数据（Bottom-of-Atmosphere corrected reflectance）

下图为哨兵2号的波段介绍。

下载数据完成后，影像数据保存在[GRANULE]-[IMG_DATA]文件夹下

【2.数据处理】

在网上查了下，ENVI 5.3可以直接打开哨兵系列的影像（jp2格式），后续需要在envi里另存为envi支持的.dat格式，再用layer stacking工具进行波段合成就可以了。

为了方便，我选择将10m分辨率里的波段2，3，4，8一次性合成，后续进行波段运算的时候就比较方便。

[文件]-[另存为]-[另存为...（ENVI,NITF,TIFF）]，选择想要转换的文件（我这里是用10m分辨率下的TCI影像，即True Colour Image为例，大家做的时候选择B02,B03,B04,B08即可），命名输出文件即可。

注意，图中NDBI应该改为NDBWI，当时没注意写错了。 

 将文件都转为.dat格式后，进行波段合成，该工具在envi工具箱的栅格数据管理文件夹下，叫波段组合，英文版自行对照...

将那四个波段import，会自动读取信息，命名输出文件即可。

 波段融合完成。

【3.提取水体公式】

根据文章：基于python实现高分二号遥感影像水体提取与水质反演（黑臭水体与水体富营养化）_一只盒饭的博客-CSDN博客_遥感影像水体提取

得知可以通过光谱影像中的红绿波段来识别正常水体和黑臭水体 ，但是通过波段计算后发现，虽然最后结果能区分开水体和建筑和裸地，但是不能很好的区分开水体和植被，所以还需要计算一个NDVI指数来将植被剔除。

NDVI 的计算公式非常简单，为近红外（NIR）与红光（Red）之差与之和的比值，即：

NDVI = (NIR - Red) / (NIR + Red)，对应到哨兵2号即为第八波段（NIR）和第四波段（Red）

此公式的结果介于 -1 和 +1 之间，如果红色通道的反射率较低，而近红外通道的反射率较高，则将产生较高的 NDVI 值；反之亦然。

【4.水体提取】

公式已经列出来了，现在就需要确定一个提取的阈值，来保证最后结果大部分都为水体，没有其他地物干扰。这一个阈值需要针对实际来确定，结合他人文章和自我判断，暂时将：

0<NDBWI<0.05

范围内的水体识别为黑臭水体，不正常水体。将：

NDBWI≥0.05

识别为正常水体。

注意，因为还没将水体和植被区分开，所以还需要添加一个条件：

NDVI≤0

来将植被和水体区分开。

最后的波段计算公式为：

不健康水体：

(((float(b8)-float(b4))/(float(b8)+float(b4)) le 0) and  ((float(b3)-float(b4))/(float(b3)+float(b4)) lt 0.05) and  ((float(b3)-float(b4))/(float(b3)+float(b4)) gt 0.00))*1.0

健康水体：

(((float(b8)-float(b4))/(float(b8)+float(b4)) le 0) and  ((float(b3)-float(b4))/(float(b3)+float(b4)) ge 0.05))*1.0

附：

 结果：

不健康水体

 健康水体

【5.补充】

1.按照网上查阅的资料加上自己的理解，曲线的数值是大气底层反射率（不知道能不能理解为地物反射率），但是在envi里被默认放大了10000倍。

2023年4月10日更新

        在某乎上不经意间看到了一篇文章，里面提到了Landsat Collection 1 Level 2影像的像元值（pixel value）需要乘以0.0001（scale factor）才能转化为地表反射率的值。我已经忘记了做这次实验用的是不是Landsat C1L2数据，但是我认为就是这个原因。

附链接：Landsat Collection 1 和Collection 2数据集（徐丹丹，安德帅） - 知乎

 2.对最终结果新建密度分割，想要输出密度分割结果为shp文件时，envi报错。

具体原因是啥不清楚，全英文路径也不行，最后是按照贴吧老哥的办法，新建了一个window账户，在该账户内操作才成功，比较麻烦，但是居然成功了，逆天操作。

在ArcMap里面打开：

附贴吧老哥链接：求助，辐射定标和导出shp文件都有这个问题【envi吧】_百度贴吧

 3.由于本人非专业人士，若其中有错误或漏洞或本人理解有误的地方，还请各位不吝赐教，大家一起交流学习！

4.这位大佬写的哨兵系列的介绍很详细：2.Sentinel(哨兵系列)详细介绍以及下载方式 - 简书

5.另一位大佬的文章：哨兵系列卫星（哨兵1号-12号）介绍与下载教程 - 知乎 

6.补充一个问题：就是识别出的水体还是有一部分零零散散的噪声点，这些是误差的点，目前还不知道怎么去除，ENVI里的图像平滑好像用不了。

20230915

对于噪声点，一个比较简单的办法是在ArcGIS里面筛选出面积小于一定量的要素，删除即可，只要保守一点误差应该不会很大