python对landsat8数据进行辐射标定和大气校正

import numpy as np
import rasterio
from rasterio.plot import show
from rasterio.warp import reproject, Resampling
 
# 设定输入和输出的影像路径
input_tif_path = 'path_to_input_landsat8_tif.tif'
output_tif_path = 'path_to_output_calibrated_tif.tif'
 
# 设置大气校正参数
atmospheric_correction_coefficients = [
    -0.000003892, 0.000024572, 0.000002235, -0.000000093, 0.000006173,
    -0.000000236, -0.000029017, 0.000011665, -0.000000641, 0.000014978,
    0.000000039, -0.000000019, 0.000000194, -0.000000005, 0.000000014
]
 
# 读取影像数据
with rasterio.open(input_tif_path) as src:
    data = src.read()
    profile = src.profile
    profile.update(dtype=rasterio.float32)
 
# 进行大气校正
# 假设影像波段顺序为 [B2, B3, B4, B5, B6, B7]
bands_to_correct = [data[i] for i in [2, 3, 4, 5, 6, 7]]
corrected_bands = [band - np.polyval(atmospheric_correction_coefficients, band) for band in bands_to_correct]
 
# 将校正后的波段写入新的影像
profile.update(count=len(corrected_bands))
with rasterio.open(output_tif_path, 'w', **profile) as dst:
    dst.write(corrected_bands)
 
# 输出结果可视化
with rasterio.open(output_tif_path) as src:
    show(src.read(), transform=src.transform)

这段代码使用了rasterio库来读取和写入GeoTiff格式的影像数据，并使用numpy来执行大气校正的多项式计算。代码中的atmospheric_correction_coefficients是假设的大气校正多项式的系数，实际应用中需要根据实际的卫星传感器参数进行调整。代码示例中的波段索引[2, 3, 4, 5, 6, 7]假设是Landsat 8数据的标准波段顺序，实际应用中需要根据具体数据进行调整。最后，代码使用rasterio.plot.show函数将校正后的影像可视化。