GEE(python版API)提取采样点附近5*5共25个像元的反射率值方法1(建立缓冲区buffer版)

建立在另一篇《python版GEE批量提取采样点最临近时间内landsat中心像元的反射率》的基础上来求的,由于上一篇已经求出来影像ID以及采样点中心像元最中心的坐标经纬度ImageCoords,根据像元之间的距离建立缓冲区圆,将圆变为矩形再变为数组,提取周围25个像元的反射率。

如图,5*5=25个像元(landsat一个像元分辨率为30m),知采样点中心像元最中心O点的坐标经纬度ImageCoords,建立一个缓冲区圆使其能够缓冲到附近的25个像元(不能超过25个像元要求建立的缓冲区圆半径小于OA距离,即小于75m;不能少于25个像元要求建立的缓冲区圆半径大于OB,即大于45倍根2;所以缓冲区半径要求大于等于64,小于75)

开始来写代码

1、导入必要的库

import ee
import geemap
import tqdm    #显示进度条
import pandas as pd
import numpy as np
geemap.set_proxy(****)#端口号
ee.Initialize()

2*去云和转换计算,可要可不要,我写这段的时候没有要,因为上一篇笔记已经去云了,先记下来万一以后用得着

def maskL457sr(image):    #定义一个影像去云,mask掉L457质量不好的函数
    qaMask = image.select('QA_PIXEL').bitwiseAnd(int('11111', 2)).eq(0)    #以QA波段为标准的mask,int('11111', 2)为31
    saturationMask = image.select('QA_RADSAT').eq(0)        #用于检查 QA_RADSAT 波段的每个像素值是否等于 0,0的部分为饱和,mask掉
    # Apply the scaling factors to the appropriate bands.
    opticalBands = image.select('SR_B.').multiply(0.0000275).add(-0.2)     #将常见的光谱波段进行转换计算
    thermalBand = image.select('ST_B6').multiply(0.00341802).add(149.0)    #热红外波段转换计算
    # Replace the original bands with the scaled ones and apply the masks.
    return image.addBands(opticalBands, None, True) \
                .addBands(thermalBand, None, True) \
                .updateMask(qaMask) \
                .updateMask(saturationMask)     #将换算过的SR-B。    热红外波段添加到影像中,更新mask后的QA_PIXEL,QA_RADSAT,返回影像

3、建立缓冲区,后面需要3*3,7*7也可以依次类推改一改参数就好了

df = pd.read_excel(r'C:\Users\带有影像ID和像元中心点经纬度的表格.xlsx')
bands = ['SR_B1', 'SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B7']    #将这几个波段放进字典bands中
datas = []    #准备一个表格存储数据
for index, row in tqdm.tqdm(df.iterrows()):   #tqdm.tqdm()创建进度条对象,df.iterrows()逐行遍历表格,返回每一行索引index和当前行数据row
    try:
        # 获取image影像
        image = ee.Image(row['ImageId'])
        # 定义选择反射率波段
        # 对反射率值进行缩放(官方提供)
        image = image.select(bands).multiply(0.0000275).add(-0.2)
        # 获取点坐标位置
        point = ee.Geometry.Point([row['centrol_lon'], row['centrol_lat']])
        # 以该点为中心生成5*5的矩阵区域
        area_of_interest = point.buffer(64).bounds()  #这个值不一定是64,64到74每一个数都行,一般来说,前提是这里提供的经纬度必须得是像元最中心
        # 根据5*5矩形区域将该范围影像的值转为矩阵     图像对象转换为 NumPy 数组
        im_data = geemap.ee_to_numpy(image,region=area_of_interest, scale=30)
        # 验证一遍
        if im_data.shape[0]*im_data.shape[1] == 25:
            datas.append({"number":row['number'],"coordinates":(row['longitude'], row['latitude']), "matrix":im_data})
    except Exception as e:
        print(f"{row['number']} 点矩阵大小完全不符合要求")
    # print(datas)

4、提取反射率到表格里去

# 创建一个空的 DataFrame 用于存储所有数据
all_data = []

for data in datas:
    x, y = data['coordinates']
    matrix = data['matrix']
    number = data['number']
    # print(f"Number: {number}, Coordinates: ({x}, {y}), Matrix: {matrix}")
    
    # 添加当前点的数据
    for i in range(5):  # 遍历每一行
        for j in range(5):  # 遍历每一列
            row = [number,x, y, i + 1, j + 1]  # 添加点坐标和像素位置,从1开始
            row.extend(matrix[i, j, :])  # 添加该像素的所有波段值,matrix为三位数组,第三围表示波段
            
            all_data.append(row)
    
    # 添加空行
    all_data.append([''] * (5 + len(bands)))

# 创建 DataFrame

columns = ['number','Longitude', 'Latitude', 'Row', 'Col'] + bands
df = pd.DataFrame(all_data, columns=columns)

# 写入 Excel 文件
output_file = r'C:\Users\需要保存的表格.xlsx'
df.to_excel(output_file, index=False)

print(f"Data has been written to {output_file}")

读取原始数据的表格是这样的前面的经纬度是采样点经纬度,后面的经纬度是开头提到的上一篇笔记提的最中心的经纬度,也可以两篇笔记的代码结合一下

结果表格是这样的,有行列信息,能清楚指导25个像元位置

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值