### 使用Python处理遥感数据实现地表温度反演 #### 方法概述 为了通过Python进行遥感数据处理并实现地表温度反演，通常采用单窗算法作为基础方法之一。该过程涉及多个步骤和技术细节，包括但不限于获取卫星影像、预处理数据以及应用特定公式来估算地表温度。 #### 所需库介绍 对于此类任务而言，有几个重要的Python库可以提供帮助： - **rasterio**: 用于读取和写入栅格文件； - **numpy**: 提供高效的数组运算支持； - **scipy**: 实现科学计算功能； - **matplotlib/pyplot**: 绘制图表辅助分析结果可视化； - **gdal (osgeo)**: 支持地理空间数据分析与转换操作。 #### 数据准备阶段 针对Landsat系列传感器所采集的数据集，在实际编程前需要注意区分不同版本之间的差异，比如Collection 1 和 Collection 2 的缩放因子有所不同[^3]。因此，在加载原始数据之后应当先确认其具体版本以便正确调整后续计算逻辑中的参数设置。 #### 单窗算法的应用实例 下面给出一段简单的代码片段展示如何利用上述提到的技术栈完成基本的地表温度反演工作流： ```python import numpy as np from osgeo import gdal, ogr import rasterio def calculate_brightness_temperature(dn_values): """将DN值转换成亮度温度""" k1 = 774.8853 # 常量K1_λ 对于 Landsat 8 BAND10 波段 k2 = 1321.0789 # 常量K2_λ 对于 Landsat 8 BAND10 波段 brightness_temp = (k2 / np.log((k1/dn_values)+1)) return brightness_temp def apply_single_window_algorithm(brightness_temperatures, emissivity=0.99): """基于给定的亮度温度和表面发射率估计真实地表温度""" lambda_val = 10.895e-6 # 发射辐射波长 μm 转换为 m单位 c1_star = 1.19104273e-8 * ((lambda_val*1e6)**(-5)) # 辐射常数C1* c2 = 1.4387752 # 第二类普朗克常数 C2 surface_temperature_kelvin = ( brightness_temperatures / (1 + (c1_star/(emissivity*(brightness_temperatures**2))) * np.exp(c2/(lambda_val*brightness_temperatures))) ) return surface_temperature_kelvin - 273.15 # 将 K 温度转摄氏度显示 with rasterio.open('path_to_your_landsat_image.tif') as src: band_data = src.read(1).astype(float) # 应用缩放因子（此处假设为Collection 2） scale_factor = 0.0001 band_radiance = band_data * scale_factor bt = calculate_brightness_temperature(band_radiance) lst = apply_single_window_algorithm(bt) print(f"Estimated Land Surface Temperature:\n{lst}") ``` 此脚本展示了从读取图像到最终获得地表温度图的过程。注意这里简化了一些可能存在的复杂情况，例如云层遮挡校正等问题并未在此详述。