使用matplotlib畫色斑圖

幾個要點

1、原始數據是站點數據，故要用插值法轉成格點數據

2、涉及到指定地圖范圍的色斑圖，必須用到地圖白化，借鑒原文：嗶哩嗶哩 搜索摸魚氣象&? ? ? ?Pthon 白化繪圖基礎，這里就不粘貼了，涉嫌廣告。

3、繪制等值線用函數Axes.contour(*args,?data=None,?**kwargs)

? ? ?等值線填圖用函數Axes.contourf(*args,?data=None,?**kwargs)

先貼一張成品圖：強降水次數分布圖

代碼涉及到的數據有一個降水的excel文檔，和一個本市的分縣shp文檔

# -*- coding: utf-8 -*- import numpy as np import pandas as pd from scipy.interpolate import Rbf # 徑向基函數 ： 將站點信息插到格點上 用于繪制等值線import matplotlib.pyplot as plt import matplotlib.colors as colors import matplotlib as mplimport cartopy.crs as ccrs # 默認投影 from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter import geopandas as gpd from matplotlib.path import Path from 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') # 讀取Excel lon = 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 插值法linear count_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-map norm = 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(['**區','**區','**區','**區','**縣','**市'])] # 生成裁剪路徑 – 關鍵操作1 path_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包，安裝的有時候會不成功，需要多試幾次。