Python绘制流线图因x、y非等间距报错求助

已实现用Python绘制WRF数据的相对湿度（RH）和风羽，但添加流线图时失败，报错提示x、y坐标非等间距。尝试以下操作均未解决：

• 仅完成坐标投影转换：x, y = m(LON[Time_Index],LAT[Time_Index])
• 使用np.arange()/np.linspace()时触发错误：IndexError: too many indices for array: array is 1-dimensional, but 2 were indexed
• 使用meshgrid时出现警告：WARNING: x coordinate not montonically increasing，后续contourf又报错：ValueError: operands could not be broadcast together with shapes (78,129) (10062,10062)

数据由外部脚本从WRF文件读取，需解决坐标等间距处理问题以绘制流线图。

WRF输出的是非规则曲线网格，投影后的坐标自然也不是等间距的，而streamplot要求输入严格的等间距网格。解决核心是将原始风场插值到自定义的规则网格上：

1. 在地图范围内生成等间距的经纬度规则网格
2. 用插值方法将原始U、V风场数据映射到规则网格
3. 将规则经纬度转换为Basemap投影坐标
4. 基于插值后的风场和规则投影坐标绘制流线图

# Import libraries      
import numpy as np
import netCDF4
from read_wrf import *
from matplotlib import pyplot as plt
from mpl_toolkits.basemap import Basemap, cm
from datetime import datetime, timedelta
import shapefile as shp
from scipy.interpolate import griddata

# Read in the file using the read_wrf script
filename = './wrfout_d01_2022-09-28_00:00:00'
LAT, LON, RH, U, V = main(filename)

# Get the times from the model and pressure levels
TIMES = get_wrf_var('Times', filename)
PLevels = np.array([1000, 850, 700, 500, 300, 200])

# Calculate wind speed magnitude
WSPD = np.sqrt(U**2 + V**2)

# Loop over times, and pressure levels
for Time_Index, T in enumerate(TIMES):
    for Height_Index, H in enumerate(PLevels):
        # Declare the figure
        fig = plt.figure(figsize=(10,10))
        ax = plt.subplot(111)

        # Define map extent
        lllon, lllat, urlon, urlat = 15.63819, -34.30278, 36.36181, -23.41068
        
        # Set up Basemap instance
        m = Basemap(
            projection='merc', 
            llcrnrlon=lllon, llcrnrlat=lllat, 
            urcrnrlon=urlon, urcrnrlat=urlat, 
            resolution='h'
        )
            
        # Read shapefile
        m.readshapefile(
            '/home/zmumba/DA/SCRIPTS/06_Utility_Files/Shapefiles/Lesotho/lso_admbnda_adm1_FAO_MLGCA_2019', 
            'lso_admbnda_adm1_FAO_MLGCA_2019'
        )

        # Add Coastlines, States, and Country Boundaries
        m.drawcoastlines()
        m.drawmapboundary(fill_color='white')
        m.drawcountries(
            linewidth=1.25, linestyle='solid', color='#000073',
            antialiased=True, ax=ax, zorder=3
        )
           
        # Add Grid Lines
        # Draw parallels
        parallels = np.arange(-35., -25., 5.)
        m.drawparallels(
            parallels, labels=[1,0,0,0], color='0.25',
            linewidth=0.5, fontsize=10
        )
        m.drawlsmask(land_color='coral', ocean_color='aqua', lakes=True)

        # Draw meridians
        meridians = np.arange(15., 35., 5.)
        m.drawmeridians(
            meridians, labels=[0,0,0,1], color='black',
            linewidth=0.5, fontsize=10
        )

        # ---------------------- 处理坐标与插值 ----------------------
        # 原始WRF网格的经纬度（当前时间步）
        lon_2d = LON[Time_Index]
        lat_2d = LAT[Time_Index]
        # 转换为投影坐标（用于原始数据绘图）
        x_orig, y_orig = m(lon_2d, lat_2d)

        # 生成规则经纬度网格（自定义密度，这里用和原始网格相近的点数）
        lon_regular = np.linspace(lllon, urlon, lon_2d.shape[1])
        lat_regular = np.linspace(lllat, urlat, lon_2d.shape[0])
        lon_reg_grid, lat_reg_grid = np.meshgrid(lon_regular, lat_regular)
        # 转换为投影坐标（用于流线图）
        x_reg, y_reg = m(lon_reg_grid, lat_reg_grid)

        # 将原始风场插值到规则网格
        # 先把原始网格展平为点集
        points = np.column_stack((lon_2d.flatten(), lat_2d.flatten()))
        # 插值U和V
        U_interp = griddata(points, U[Time_Index, Height_Index].flatten(), (lon_reg_grid, lat_reg_grid), method='linear')
        V_interp = griddata(points, V[Time_Index, Height_Index].flatten(), (lon_reg_grid, lat_reg_grid), method='linear')

        # ---------------------- 绘制RH和风羽 ----------------------
        clev = np.arange(0, 110, 5)
        # 用原始投影坐标绘制RH填色
        im = m.contourf(x_orig, y_orig, RH[Time_Index, Height_Index], clev, cmap=plt.cm.BrBG)

        # 绘制风羽（采样原始网格）
        m.barbs(
            x_orig[::10, ::10], y_orig[::10, ::10],
            U[Time_Index, Height_Index, ::10, ::10], V[Time_Index, Height_Index, ::10, ::10],
            sizes=dict(emptybarb=0.25, spacing=.1, height=0.3), 
            flip_barb=True, color='r'
        )
        
        # ---------------------- 绘制流线图 ----------------------
        m.streamplot(
            x_reg, y_reg, U_interp, V_interp,
            linewidth=1., color='blue', density=2.5
        )

        # Add Colorbar
        cbar = m.colorbar(im)
        cbar.set_label('Relative Humidity (%)')

        # Save image and display
        vt = Time_Index  # 根据实际时间偏移逻辑调整
        plt.title(f"RH and Wind at {H} hPa, T+{vt}H")
        plt.savefig(f"V{H}_T+{vt}H.png", bbox_inches='tight')
        plt.show()
        plt.close()

• 添加scipy.interpolate.griddata用于风场插值，解决非规则网格到规则网格的转换问题
• 区分原始投影坐标（用于RH填色和风羽）和规则投影坐标（用于流线图），避免维度不匹配错误
• 修正原代码中的缩进问题（循环内代码未正确缩进）
• 补充了标题和色标标签，提升图件可读性

内容的提问来源于stack exchange，提问作者Zilore Mumba