Python读取.grd高程数据

import xarray as xr
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from geopy.distance import geodesic

def extracted_mete_etopo(lon_i, lat_i, lon_arr, lat_arr, aim_mete_arr, distance):
    # 选取一个更小的经纬度框
    # 根据所选的研究区域的纬度范围，可以选取±1°的范围内
    
    lon_index = np.argmin(np.abs(lon_i - lon_arr))
    # lon_index[0][0]
    # lon_index[-1][0]
    lat_index = np.argmin(np.abs(lat_i - lat_arr))
    # lat_index[0][0]
    # lat_index[-1][0]
    region_dat = aim_mete_arr[lat_index, lon_index]
    
    return region_dat

file = 'ETOPO1_Bed_g_gmt4.grd'
data = xr.open_dataset(file)

lon_arr = np.array(data.x)
lat_arr = np.array(data.y)
z = np.array(data.z)

# % 读取文件和数据
loc_file = 'Ross_ADPE_track.csv'
loc_df = pd.read_csv(loc_file)

distance = 10
# % 下面开始循环提取值
ext_mete_list = []
for i in range(len(loc_df)):
    print(i)
    lon_i = loc_df.iloc[i,4]
    lat_i = loc_df.iloc[i,5]
    mete_arr = z
    ext_mete = extracted_mete_etopo(lon_i, lat_i, lon_arr, lat_arr, mete_arr, distance)
    ext_mete_list.append(ext_mete)

# % 导出数据
out_put = pd.DataFrame({'z':ext_mete_list})
out_put.to_csv('z_'+str(distance)+'km.csv')