GEE|实现K-T变换并导出指定波段

基础知识

缨帽变换特征是通过缨帽变换(Tasseled Cap Trasform,TCT,又称为 K-T变换)得到的特征分量,是一种特殊的主成分分析方法(PCA)。但与 PCA 所不同的是,缨帽变换具有固定的变换矩阵。缨帽变换根据这个固定的变换矩阵将原始影像投影综合变换到具有物理意义的亮度(Brightness)、绿度(Greenness)和湿度(Wetness)特征向量三维特征空间,充分反映了裸土岩石、植被覆盖度和水分信息。这个变换的过程达到了减少特征维数、增强影像信息的效果。其变换公式如下:

式中,Y 表示缨帽变换后栅格影像;C 表示变换矩阵系数; X 表示变换前影像。

缨帽变换只需要确定变换矩阵系数。变换矩阵系数与影像的传感器相关。不同的传感器影像具有不同的变换矩阵系数。

代码实现

代码使用之前需要先将矢量边界导入:

实现缨帽变换的代码如下:

//去云函数
function maskS2clouds(image) {
  var qa = image.select('QA60');
 
  // Bits 10 and 11 are clouds and cirrus, respectively.
  var cloudBitMask = 1 << 10;
  var cirrusBitMask = 1 << 11;
 
  // Both flags should be set to zero, indicating clear conditions.
  var mask = qa.bitwiseAnd(cloudBitMask).eq(0)
      .and(qa.bitwiseAnd(cirrusBitMask).eq(0));
 
  return image.updateMask(mask).divide(10000);
}
var coefficients = ee.Array([  
  [0.0356, 0.0822, 0.1360, 0.2611, 0.2964, 0.3338, 0.3877, 0.3895, 0.0949, 0.3882, 0.1366, 0.4750],  
  [-0.0635, -0.1128, -0.1680, -0.3480, -0.3303, 0.0852, 0.3302, 0.3165, 0.0467, -0.4578,  -0.4064, 0.3625],  
  [0.0649, 0.1363, 0.2802, 0.3072, 0.5288, 0.1379, -0.0001, -0.0807, -0.0302,-0.4064, -0.5602,-0.1389]
]);  
//影像转换为array  
var image =  ee.ImageCollection('COPERNICUS/S2_SR')
                  .filterDate('2021-01-01','2021-12-31')
                  .filterBounds(aoi)
                  .filter(ee.Filter.lt('CLOUDY_PIXEL_PERCENTAGE',10))
                  .map(maskS2clouds)
                  .select(['B1','B2','B3','B4','B5','B6','B7','B8','B9','B11','B12','B8A']);  
 
              
var arrayImage1D = image.mean().toArray();  
//print("arrayImage1D",arrayImage1D)
Map.addLayer(arrayImage1D,{},"arrayImage1D")
 
 
var arrayImage2D = arrayImage1D.toArray(1);  
//print("arrayImage2D",arrayImage2D)
Map.addLayer(arrayImage2D,{},"arrayImage2D")
 
 
//相乘变换  
var componentsImage = ee.Image(coefficients)   
                        .matrixMultiply(arrayImage2D)   
                        .arrayProject([0])   
                        .arrayFlatten([[  
                          'brightness', 'greenness', 'wetness' 
                        ]]);  

var vizParams = {  
  "bands": ['brightness', 'greenness', 'wetness'],  
  "min": -0.1,   
  "max": [0.5, 0.1, 0.1]  
};  
Map.centerObject(image.mean(), 7);  
Map.centerObject(aoi)
Map.addLayer(componentsImage, vizParams, "componentsImage"); 
print("componentsImage",componentsImage)

Export.image.toDrive({
      image:componentsImage.select(["greenness"]),
      description:'GREEN',
      scale:30,
      maxPixels: 1e13,
      region:aoi,
      crs: "EPSG:32649",
      fileFormat: 'GeoTIFF',
      formatOptions: {
        cloudOptimized: true
      }
  });

运行结果

posted @ 2022-08-18 16:29  Weltㅤ  阅读(524)  评论(0编辑  收藏  举报