保姆级教程:在GEE里用Landsat 5数据一键计算亮度、绿度、湿度(附完整代码)
2026/6/9 10:58:47 网站建设 项目流程

零基础玩转GEE:Landsat 5缨帽变换实战指南

当你第一次听说"缨帽变换"这个专业术语时,是不是感觉像在听天书?别担心,今天我们就用最接地气的方式,带你从零开始在Google Earth Engine(GEE)平台上实现Landsat 5数据的缨帽变换分析。不需要复杂的数学基础,不需要理解矩阵运算原理,只需要跟着步骤操作,你就能轻松获取亮度、绿度和湿度三个关键指标的可视化结果。

1. 准备工作与环境搭建

在开始之前,确保你已经拥有一个可用的GEE账号。访问 Google Earth Engine官网 注册即可免费使用。登录后,点击"Code Editor"进入编程界面,这里就是我们今天的主战场。

GEE平台提供了丰富的遥感数据集,其中Landsat系列是最常用的地表观测数据之一。我们今天使用的Landsat 5数据已经过大气校正(TOA),包含以下波段:

波段名称波长范围(nm)主要用途
B1450-520蓝光波段
B2520-600绿光波段
B3630-690红光波段
B4760-900近红外波段
B51550-1750短波红外1
B72080-2350短波红外2

提示:在GEE中,Landsat 5的波段编号与Landsat 8/9有所不同,使用时需特别注意对应关系。

2. 数据加载与预处理

现在让我们开始加载数据。在代码编辑器中,首先清除默认的示例代码,然后输入以下内容:

// 定义研究区域(以北京为例) var roi = ee.Geometry.Point([116.404, 39.915]); // 加载Landsat 5 TOA数据 var image = ee.Image("LANDSAT/LT05/C01/T1_TOA/LT05_044034_20081011") .select(['B1', 'B2', 'B3', 'B4', 'B5', 'B7']); // 设置地图中心并显示原始影像 Map.centerObject(roi, 10); Map.addLayer(image, {bands: ['B4', 'B3', 'B2'], min: 0, max: 0.3}, '原始影像');

这段代码做了三件事:

  1. 定义了研究区域坐标(北京天安门)
  2. 加载了2008年10月11日的Landsat 5数据
  3. 以假彩色方式显示影像(近红外-红-绿波段组合)

运行后,你应该能在右侧地图窗口看到北京市区的卫星影像。如果想让影像更清晰,可以调整minmax参数的值。

3. 缨帽变换核心实现

缨帽变换的关键在于那个神秘的变换矩阵。对于Landsat 5数据,科学家们已经通过大量实验确定了最佳系数,我们直接拿来用即可:

// 定义Landsat 5缨帽变换系数矩阵 var coefficients = ee.Array([ [0.3037, 0.2793, 0.4743, 0.5585, 0.5082, 0.1863], // 亮度系数 [-0.2848, -0.2435, -0.5436, 0.7243, 0.0840, -0.1800], // 绿度系数 [0.1509, 0.1973, 0.3279, 0.3406, -0.7112, -0.4572] // 湿度系数 ]);

接下来是重头戏——矩阵运算。别被吓到,GEE已经帮我们封装好了相关方法:

// 将影像转换为数组形式 var arrayImage1D = image.toArray(); var arrayImage2D = arrayImage1D.toArray(1); // 执行矩阵乘法运算 var componentsImage = ee.Image(coefficients) .matrixMultiply(arrayImage2D) .arrayProject([0]) .arrayFlatten([['brightness', 'greenness', 'wetness']]);

这段代码完成了缨帽变换的核心计算:

  1. toArray()将影像转换为数学上的数组形式
  2. matrixMultiply()执行矩阵乘法运算
  3. arrayFlatten()将结果重新转换为影像格式

4. 结果可视化与导出

现在我们已经得到了亮度、绿度和湿度三个分量,如何直观地查看它们呢?

// 设置可视化参数 var vizParams = { bands: ['brightness', 'greenness', 'wetness'], min: -0.1, max: [0.5, 0.1, 0.1] }; // 显示缨帽变换结果 Map.addLayer(componentsImage, vizParams, '缨帽变换结果');

为了让结果更清晰,我们可以单独查看每个分量:

// 单独显示各分量 Map.addLayer(componentsImage.select('brightness'), {min: -0.1, max: 0.5, palette: ['black', 'white']}, '亮度分量'); Map.addLayer(componentsImage.select('greenness'), {min: -0.1, max: 0.1, palette: ['brown', 'yellow', 'green']}, '绿度分量'); Map.addLayer(componentsImage.select('wetness'), {min: -0.1, max: 0.1, palette: ['red', 'blue']}, '湿度分量');

如果需要将结果导出到Google Drive,可以使用以下代码:

// 导出结果到Google Drive Export.image.toDrive({ image: componentsImage, description: 'TCT_Result', folder: 'GEE_Exports', fileNamePrefix: 'L5_TCT', region: roi.buffer(10000).bounds(), scale: 30, maxPixels: 1e13 });

5. 实际应用与参数调整

缨帽变换的三个分量各有其独特的应用价值:

  • 亮度分量:反映地表反射率,可用于识别裸土、城市建筑等
  • 绿度分量:与植被覆盖度高度相关,是植被监测的重要指标
  • 湿度分量:反映地表水分状况,可用于干旱监测、湿地识别等

如果你想分析不同时期的土地变化,只需修改影像ID即可。例如,分析1990年的北京:

// 加载1990年影像 var image1990 = ee.Image("LANDSAT/LT05/C01/T1_TOA/LT05_044034_19900915") .select(['B1', 'B2', 'B3', 'B4', 'B5', 'B7']); // 应用相同的缨帽变换 var array1990 = image1990.toArray().toArray(1); var components1990 = ee.Image(coefficients) .matrixMultiply(array1990) .arrayProject([0]) .arrayFlatten([['brightness', 'greenness', 'wetness']]); // 比较2008年与1990年的绿度变化 var greennessChange = componentsImage.select('greenness') .subtract(components1990.select('greenness')); Map.addLayer(greennessChange, {min: -0.05, max: 0.05, palette: ['red', 'white', 'green']}, '绿度变化(2008-1990)');

6. 常见问题排查

在实际操作中,你可能会遇到以下问题:

  1. 影像加载失败

    • 检查影像ID是否正确
    • 确认该影像在研究区域内有数据覆盖
  2. 计算结果异常

    • 确认波段选择顺序与系数矩阵匹配
    • 检查影像是否经过大气校正(TOA或SR)
  3. 可视化效果不佳

    • 调整minmax参数值
    • 尝试不同的颜色方案(palette)
  4. 导出任务未启动

    • 检查Google Drive存储空间
    • 确认导出区域不过大

注意:GEE对免费用户有一定的计算资源限制,处理大范围或长时间序列数据时,建议先在小范围测试。

7. 完整代码清单

为了方便使用,这里提供完整的可运行代码:

// 定义研究区域 var roi = ee.Geometry.Point([116.404, 39.915]); // 加载Landsat 5数据 var image = ee.Image("LANDSAT/LT05/C01/T1_TOA/LT05_044034_20081011") .select(['B1', 'B2', 'B3', 'B4', 'B5', 'B7']); // 显示原始影像 Map.centerObject(roi, 10); Map.addLayer(image, {bands: ['B4', 'B3', 'B2'], min: 0, max: 0.3}, '原始影像'); // 定义缨帽变换系数 var coefficients = ee.Array([ [0.3037, 0.2793, 0.4743, 0.5585, 0.5082, 0.1863], [-0.2848, -0.2435, -0.5436, 0.7243, 0.0840, -0.1800], [0.1509, 0.1973, 0.3279, 0.3406, -0.7112, -0.4572] ]); // 执行缨帽变换 var arrayImage1D = image.toArray(); var arrayImage2D = arrayImage1D.toArray(1); 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.addLayer(componentsImage, vizParams, '缨帽变换结果'); // 导出设置 Export.image.toDrive({ image: componentsImage, description: 'TCT_Result', folder: 'GEE_Exports', fileNamePrefix: 'L5_TCT', region: roi.buffer(10000).bounds(), scale: 30, maxPixels: 1e13 });

在实际项目中,我发现最容易出错的地方是波段选择顺序与系数矩阵的对应关系。有一次我误将SWIR1和SWIR2波段顺序颠倒,导致湿度分量结果完全异常。因此建议在正式分析前,先用小范围测试确认各分量结果是否符合预期。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询