python根据点经纬度来提取遥感影像中对应的像元值（单张&多张影像）

问题：现在我有一些点坐标数据和一幅带有地理信息的遥感影像，如何提取这些点对应的遥感影像的光谱信息呢？

废话不多说，直接上码：

先导入一些包

import rasterio
import numpy as np
import os
import pandas as pd
from pyproj import Transformer

先定义函数

# 提取单张影像的像素值的函数
def extract_pixel_values_single_image(image_path, points_df):
    with rasterio.open(image_path) as src:
        transform = src.transform
        band1 = src.read(1)
        vals = []
        for _, row in points_df.iterrows():
            col, row = ~transform * (row['Longitude_degree'], row['Latitude_degree'])
            col, row = int(col), int(row)
            vals.append(band1[row, col])
    return vals

单张影像提取

如果是单张影像，直接利用这个函数即可，函数中band1 = src.read(1)是读取的第一个波段，参数image_path是遥感影像的路径，points_df是对应的点数据，这里Longitude_degree是经度所在列名，Latitude_degree是纬度所在列名

如：

多张影像提取

如果是多张影像提取，可以利用循环读取文件夹，结合单张影像提取函数来读取：

for folder_name1 in ['fn1','fn2']:
    for folder_name2 in ['fn11', 'fn22', 'fn33', 'fn44', 'fn55', 'fn66']:
        for continent in ['c1','c2','c3','c4','c5']:
            for year in range(2000, 2020):  # 要处理的年份
                csv_path =r'H:\1.Datasets\7.LandCover\Global\\'+folder_name1+'\\data_'+folder_name2+'\\'+continent+'\\data_'+folder_name2+'_'+continent+'_All_' + str(year) + '_' + str(year + 1) + '_LCdata.csv'
                df = pd.read_csv(csv_path)

                # 数据类型及其对应文件夹路径
                data_folders = {
                    '2m_air_temperature_mean':r"H:\1.Datasets\global_temperature_2m_mean\Global_Annual_temperature_2m_mean_" + str(year) + '.tif',
                    '2m_air_temperature_max':r"H:\1.Datasets\global_temperature_2m_max\Global_Annual_temperature_2m_max_" + str(year) + '.tif',
                    '2m_air_temperature_min':r"H:\1.Datasets\global_temperature_2m_min\Global_Annual_temperature_2m_min_" + str(year) + '.tif',
                    'precipitation_sum':r"H:\1.Datasets\global_total_precipitation\Global_Annual_Total_Precipitation_" + str(year) + '.tif',
                    'precipitation_max':r"H:\1.Datasets\global_max_precipitation\Global_Annual_max_Precipitation_" + str(year) + '.tif',
                    'precipitation_min':r"H:\1.Datasets\global_min_precipitation\Global_Annual_min_Precipitation_" + str(year) + '.tif'
                 }

                # 遍历每种数据类型及其文件夹
                for data_type, image_path in data_folders.items():
                    if os.path.exists(image_path):
                        pixel_values = extract_pixel_values_single_image(image_path, df)
                        df[data_type] = pixel_values
                    else:
                        print(f"File not found: {image_path}")

                output_fn = str(year) + '_' + str(year + 1) + '_'+ folder_name2 +'.csv'
                df.to_csv(r'H:\1.Datasets\7.LandCover\Global合并\\'+folder_name1+'\\data_'+folder_name2+'\\'+continent+'\\' + output_fn, index=False)
                print(year,folder_name1,folder_name2,continent, '完成')

前面的4个for循环主要是为了构建文件名，以及构建一个包含不同特征影像的data_folders，便于输入下面的像元提取循环。

因为我有好几层文件夹的数据都要处理，而且每个文件夹内时间范围是2000-2019，20个文件，所以文件名的构建就较为复杂。各位在用的时候可以根据自己情况进行更改。