GEE学习笔记 九十二：Sentinel-2 最新去云方法总结

    关于影像去云实际上是非常复杂的算法来实现的，但是我们不是专业研究如何识别云的，所以我们采用被人已经研究好的成果来实现我们的需求就可以。在GEE如何上实现Sentinel-2去云，之前也专门提到过，这篇文章就是对目前最新的去云方法做一个总结介绍。

---

 

下面使用例子的原始影像截图如下：

 

 

第一种方法：使用QA波段去云

    这是我们最常用的方法，具体原理就是利用QA60波段标记实现去云，具体代码如下：

var s2 = ee.ImageCollection("COPERNICUS/S2"),

    point = /* color: #98ff00 */ee.Geometry.Point([116.20553071765003, 39.404020061278715]);



/**

*

* remove cloud by QA bands

*

* */

function rmCloudByQA(image) {

  var qa = image.select('QA60');

  var cloudBitMask = 1 << 10;

  var cirrusBitMask = 1 << 11;

  var mask = qa.bitwiseAnd(cloudBitMask).eq(0)

               .and(qa.bitwiseAnd(cirrusBitMask).eq(0));

  return image.updateMask(mask);

}



function main() {

  var startDate = "2019-4-4";

  var endDate = "2019-4-8";

  Map.centerObject(point, 8);

  var s2Imgs = s2.filterDate(startDate, endDate)

                 .filterBounds(point);

  Map.addLayer(s2Imgs.first(), {min:0, max:3000, bands:["B4", "B3", "B2"]}, "raw", false);

  s2Imgs = s2Imgs.map(rmCloudByQA);

  Map.addLayer(s2Imgs.first(), {min:0, max:3000, bands:["B4", "B3", "B2"]}, "cloud");

}



main();

 

去云的效果如下图：

 

代码分析：

    这段代码中去云的方法是：rmCloudByQA()其他的代码也都是基础代码，我就不在过多解释。

 

第二种方法：计算云量分数去云

    这种方法比较有意思，是一个国外的大牛利用波段组合以及NDSI指数计算得到当前像素是不是云，具体代码如下：

var s2 = ee.ImageCollection("COPERNICUS/S2"),

    point = /* color: #98ff00 */ee.Geometry.Point([116.20553071765003, 39.404020061278715]);



var _cloudScore = function(img) {

  var rescale = function(img, exp, thresholds) {

    return img.expression(exp, {img: img})

        .subtract(thresholds[0]).divide(thresholds[1] - thresholds[0]);

  };

  var score = ee.Image.constant(1.0);

  score = score.min(rescale(img, 'img.blue', [0.1, 0.3]));

  score = score.min(rescale(img, 'img.red + img.green + img.blue', [0.2, 0.8]));

  score = score.min(rescale(img, 'img.nir + img.swir1 + img.swir2', [0.3, 0.8]));

  var ndsi = img.normalizedDifference(['green', 'swir1']);

  return score.min(rescale(ndsi, 'img', [0.8, 0.6]));

};



/**

*

* remove cloud by Score bands

*

* */

function rmCloudByScore(image, thread) {

  var preBands = ["B2","B3","B4","B8","B11","B12"];

  var newBands = ['blue','green','red','nir','swir1','swir2'];

  var score = _cloudScore(image.select(preBands, newBands));

  score = score.multiply(100).byte().rename('cloud');

  return image.addBands(score)

              .updateMask(score.lte(thread));

}



/**

*

* scale image

*

* */

function scaleImage(image) {

  var time_start = image.get("system:time_start");

  image = image.divide(10000);

  image = image.set("system:time_start", time_start);

  return image;

}



function main() {

  var startDate = "2019-4-4";

  var endDate = "2019-4-8";

  Map.centerObject(point, 8);

  var s2Imgs = s2.filterDate(startDate, endDate)

                 .filterBounds(point)

                 .map(scaleImage)

                 .map(function(image) {

                    return rmCloudByScore(image, 30);

                 });

  Map.addLayer(s2Imgs.first(), {min:0, max:0.3, bands:["B4", "B3", "B2"]}, "cloud");

}



main();

 

运行代码效果：

代码分析：

    观察这段代码，有一点需要注意的是：调用去云方法rmCloudByScore(image, threshold)需要传入两个参数，一个是影像另外一个是识别云的阈值（也就是大于多大的值就认为是云，范围是0-100）。还有就是调用这个方法之间首先需要将影像的波段做缩放和将波段名称重新命名为指定的名称，这是因为这个算法内部主要逻辑方法_cloudScore(xxx)它内部计算时候取的都是特殊值，具体参考上述代码。

第三种方法：使用s2cloudless算法去云

    这是目前最新的一种去云方法，s2cloudless算法是一种机器学习算法，Sentinel-2影像通过它可以计算得到云掩膜文件。这个算法是由SENTINEL Hub开源的一种算法，地址是：

https://github.com/sentinel-hub/sentinel2-cloud-detector

    在GEE上我们不需要自己去用这个算法来计算云掩膜，GEE官方目前正在利用这个算法来生成一个新的数据集，这个数据集就是Sentinel-2云掩膜数据集（“COPERNICUS/S2_CLOUD_PROBABILITY”）。需要注意的是由于目前这个数据集依然在生产中，所以不一定所有的地方和所有时间段都存在，我们在使用的时候一定要注意。

 

下面以具体的例子说明如何使用这个数据集：

var s2 = ee.ImageCollection("COPERNICUS/S2"),

    s2_cloud = ee.ImageCollection("COPERNICUS/S2_CLOUD_PROBABILITY"),

    point = /* color: #98ff00 */ee.Geometry.Point([116.20553071765003, 39.404020061278715]);



/**

*

* remove cloud by probability bands

*

* */

function rmCloudByProbability(image, thread) {

  var prob = image.select("probability");

  return image.updateMask(prob.lte(thread));

}



/**

*

* scale image

*

* */

function scaleImage(image) {

  var time_start = image.get("system:time_start");

  image = image.divide(10000);

  image = image.set("system:time_start", time_start);

  return image;

}



function getMergeImages(primary, secondary) {

  var join = ee.Join.inner();

  var filter = ee.Filter.equals({

    leftField: 'system:index',

    rightField: 'system:index'

  });

  var joinCol = join.apply(primary, secondary, filter);

  joinCol = joinCol.map(function(image) {

    var img1 = ee.Image(image.get("primary"));

    var img2 = ee.Image(image.get("secondary"));

    return img1.addBands(img2);

  });

  return ee.ImageCollection(joinCol);

}



function main() {

  var startDate = "2019-4-4";

  var endDate = "2019-4-8";

  Map.centerObject(point, 8);

  var s2Imgs1 = s2.filterDate(startDate, endDate)

                  .filterBounds(point)

                  .map(scaleImage);

  var s2Imgs2 = s2_cloud.filterDate(startDate, endDate)

                        .filterBounds(point);

  var s2Imgs = getMergeImages(s2Imgs1, s2Imgs2);

  s2Imgs = s2Imgs.map(function(image) {

    return rmCloudByProbability(image, 50);

  });

  Map.addLayer(s2Imgs.first(), {min:0, max:0.3, bands:["B4", "B3", "B2"]}, "cloud");

}



main();

 

运行效果：

 

代码分析：

    从上述代码可以看到，由于这是两个数据集，所以我们首先需要将筛选后的两个数据（s2Imgs1和s2Imgs2）利用join合并为一个数据集，具体方法参考getMergeImages(xxx)。合并后的影像中每一幅都包含一个叫做”probability“的波段，这个波段计算的就是云掩膜的信息（值范围是0-100）。程序中我们利用rmCloudByProbability(image, threshold)方法传入影像和云量阈值，就可以实现去云操作。

 

其实还有一种方法就是利用GEE封装好的模块算法来实现去云，具体的方法是：

但是这个方法有一定的缺陷，大家想要用的话可以自行研究一下。

 

    总结一下，这里面介绍了目前常用的哨兵2的去云方法，也展示了对于同一景影像不同去云方法最终实现的效果，大家可以通过对比来选择使用具体的方法来做，也可以同时结合多种算法来实现去云。需要说明的一点就是去云操作本身就没有完全完美的算法，每一种算法都有自己的缺陷和不足。以上面第三种最新的算法来说，它在大部分区域上计算云掩膜都比较不错，但是在某些特定的水体上却有问题，会把无云的地方计算的云量分数非常高。

 

 

---

大家如果有问题需要交流或者有项目需要合作，可以微信联系我，加微信好友请留言加上“GEE”。

知乎专栏：https://zhuanlan.zhihu.com/c_123993183

CSDN：https://blog.csdn.net/shi_weihappy