几个要点

1、原始数据是站点数据，故要用插值法转成格点数据

2、涉及到指定地图范围的色斑图，必须用到地图白化，借鉴原文：哔哩哔哩 搜索摸鱼气象&       Pthon 白化绘图基础，这里就不粘贴了，涉嫌广告。

3、绘制等值线用函数Axes.contour(*args, data=None, **kwargs)

     等值线填图用函数Axes.contourf(*args, data=None, **kwargs)

先贴一张成品图：强降水次数分布图

代码涉及到的数据有一个降水的excel文档，和一个本市的分县shp文档

# -*- coding: utf-8 -*-import numpy as npimport pandas as pdfrom scipy.interpolate import Rbf # 径向基函数 ： 将站点信息插到格点上 用于绘制等值线import matplotlib.pyplot as pltimport matplotlib.colors as colorsimport matplotlib as mplimport cartopy.crs as ccrs # 默认投影from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatterimport geopandas as gpdfrom matplotlib.path import Pathfrom cartopy.mpl.patch import geos_to_pathplt.rcParams['font.sans-serif'] = ['SimHei'] # 用来正常显示中文plt.rcParams['axes.unicode_minus'] = False # 用来正常显示负号# 数据准备df = pd.read_excel(r'C:\……\rain.xlsx') # 读取Excellon = df['Lon']lat = df['Lat']count = df['PRE_1h_Count']# 绘图准备olon = np.linspace(114.5, 115.75, 100) # 经纬坐标，0.05°分辨率 114°到117°0.05分辨率是有60个格点olat = np.linspace(29.5, 30.5, 80) # 纬度坐标，0.05°分辨率olon, olat = np.meshgrid(olon, olat) # 生成坐标网格 meshgrid网格化func = Rbf(lon, lat, count, function='linear') # 插值函数 调用Rbf插值函数中的 cubic 插值法linearcount_new = func(olon, olat) # 插值# count_new[count_new < 0] = 0 # 如果插值有可能涉及到负数，那么将其设置为0# 画布及绘图声明fig = plt.figure(figsize=(12, 8), facecolor='#666666', edgecolor='Blue', frameon=False) # 画布proj = ccrs.PlateCarree()ax = fig.add_subplot(111, projection=proj) # 绘图区# 色彩定制clevs = [10, 20, 30, 40, 50] # 自定义颜色列表cdict = ['#87CEFA', '#1E90FF', '#4169E1','#000080'] # 自定义颜色列表 '#A9F090','#40B73F','#63B7FF','#0000FE','#FF00FC','#850042'my_cmap = colors.ListedColormap(cdict) # 自定义颜色映射 color-mapnorm = mpl.colors.BoundaryNorm(clevs, my_cmap.N) # 基于离散区间生成颜色映射索引# 绘制等值线、等值线填色cf = ax.contourf(olon, olat, count_new, clevs, transform=proj, cmap=my_cmap, norm=norm)# c1 = ax.contour(olon, olat, count_new, levels=[35], colors='purple', linewidths=3, linestyles='--', transform=proj)position = fig.add_axes([0.77, 0.15, 0.03, 0.2]) # 位置【左，下，宽。高】plt.colorbar(cf, cax=position) # 颜色参照表position.set_yticklabels((10, 20, 30, 40, 50))ax.xaxis.set_major_formatter(LongitudeFormatter(zero_direction_label=True))ax.yaxis.set_major_formatter(LatitudeFormatter())ax.set_xticks(np.arange(114.4, 115.8, 0.2), crs=ccrs.PlateCarree())ax.set_yticks(np.arange(29.4, 30.4, 0.2), crs=ccrs.PlateCarree())# ax.gridlines() # 显示背景线# 边界白化shp = gpd.read_file(r'C:\……\420200.shp',encoding='utf-8') # utf 8 编码 或gbk编码 a = shp['geometry'][shp['name'].isin(['**区','**区','**区','**区','**县','**市'])]# 生成裁剪路径 – 关键操作1path_clip = Path.make_compound_path(*geos_to_path(a.to_list()))# 将裁剪路径应用到图层 – 关键操作2[collection.set_clip_path(path_clip, transform=ax.transData) for collection in cf.collections]# 绘制多边形边缘线ax.add_geometries(a, crs=ccrs.PlateCarree(), facecolor='none', edgecolor='black')plt.show()

 之前一直研究一个maskout.py文件，总是不能正常白化，干脆去学习了下白化的原理，其实也不难，就是程序最后几排的代码，关键要对自己的shp文件有一个了解，才能得出正确的边界，并且白化用的python包，安装的有时候会不成功，需要多试几次。