一个巨简单绘图掩膜海洋的方法

MissLH88

之前一直是mask array掩膜海洋上数据的，今天突然发现，
只要使用地图投影包catropy自带的海洋填色功能
然后使用zorder把这个填色的图层调增到最外面即可。
核心代码如下：

import cartopy.feature as cfeature

ax.add_feature(cfeature.OCEAN,facecolor='w',linewidth=0,zorder=100)

附赠我完整的绘图代码：
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
from netCDF4 import Dataset
import netCDF4 as nc
from wrf import getvar,disable_xarray,ALL_TIMES
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import numpy.ma as ma
from cartopy.mpl.gridliner import LATITUDE_FORMATTER, LONGITUDE_FORMATTER
from matplotlib.pyplot import MultipleLocator
import nclcmaps


plt.rcParams['font.sans-serif']=['Times New Roman']
matplotlib.rcParams.update({'font.size': 8})

    ref_lat = 32
    ref_lon = 98
    true_lat1 = 30
    true_lat2 = 60
    false_easting = (180-1)/2*30000
    false_northing = (99.5-1)/2*30000
   
      
    proj_lambert = ccrs.LambertConformal(
                    central_longitude=ref_lon,
                    central_latitude=ref_lat,
                    standard_parallels=(true_lat1,true_lat2),
                    cutoff=-30,
                    false_easting=false_easting,
                    false_northing=false_northing,
                )
    ## 创建坐标系
    fig = plt.figure(figsize=(4.5,2.8), dpi=200)  # 创建页面
    ax = fig.add_axes([0.1,0.15,0.85,0.88], projection=proj_lambert)  
    ax.set_extent([0, false_easting*2, 0, false_northing*2], crs=proj_lambert)
    ax.coastlines(linewidth=0.8,resolution='50m',zorder=101)
    ax.add_feature(cfeature.LAKES.with_scale('50m'),linewidth=0.5,edgecolor='black',facecolor='none')
    ax.add_feature(cfeature.OCEAN,facecolor='w',linewidth=0,zorder=100)
        
    # 添加shp文件
    import cartopy.io.shapereader as shpreader
    filepath ='D:\LHwork\shpfile\中国GIS地图\国家基础地理数据/hyd1_4l.shp'
    crs = ccrs.PlateCarree()
    reader = shpreader.Reader(filepath)
    geoms = reader.geometries()
    ax.add_geometries(geoms, crs, lw=0.5, fc='none')
    reader.close()
   
    ## --设置网格属性, 不画默认的标签
    gl=ax.gridlines(draw_labels=True,linestyle=":",linewidth=0.6 ,x_inline=False, y_inline=False,color='k', zorder=103)
   
    ## 关闭上面和右边的经纬度显示
    gl.top_labels=False #关闭上部经纬标签                                 
    # gl.bottom_labels = False
    # gl.left_labels = False
    import matplotlib.ticker as mticker
    gl.right_labels=False
    gl.xformatter = LONGITUDE_FORMATTER  #使横坐标转化为经纬度格式            
    gl.yformatter = LATITUDE_FORMATTER        
   
                                    
    gl.xlocator=mticker.FixedLocator(np.arange(60,160,10))                              
    gl.ylocator=mticker.FixedLocator(np.arange(10,60,10))
    gl.rotate_labels = False
    gl.xlabel_style={'size':10,}
                      # 'rotation':'horizontal','rotation_mode':'default',
                      # 'horizontalalignment':'center','verticalalignment':'center'}#修改经纬度字体大小                             
    gl.ylabel_style={'size':10,}#'rotation':'horizontal'}
    ax.spines['geo'].set_linewidth(0.6)#调节边框粗细
   
    plt.plot([110,122], [28,28], 'r--',linewidth=1,transform=ccrs.PlateCarree(),zorder=105)
    plt.plot([110,122], [32,32], 'r--',linewidth=1,transform=ccrs.PlateCarree(),zorder=105)
    plt.plot([110,110], [28,32], 'r--',linewidth=1,transform=ccrs.PlateCarree(),zorder=105)
    plt.plot([122,122], [28,32], 'r--',linewidth=1,transform=ccrs.PlateCarree(),zorder=105)
   
    # SM
    levels = (np.arange(0,0.625,0.025))
    v = np.around(np.arange(0,0.7,0.1),1)
    # SM_diff
    levels = (np.arange(-0.050,0.06,0.01))
    v = np.around(np.arange(-0.050,0.06,0.01),2)   
   
    data = ma.masked_array(data_pert-data_ctrl,mask = b)
   
    cs = ax.contourf(lon,lat,data,transform=ccrs.PlateCarree(),
                      levels=levels,cmap= nclcmaps.cmap('MPL_BrBG'),#)#MPL_BuGn#CBR_wet,#'terrain', 'OceanLakeLandSnow'
                       extend='both')
    # over_t = (224/255,224/255,224/255)
    # under_t = (0, 71/255 , 71/255 )
    # under_t = (1,1,1)
    # cs.cmap.set_over(over_t) #'fuchsia'
    # cs.cmap.set_under(under_t) #21
    # cs.changed()
   
   
    cax = fig.add_axes([0.1,0.07,0.8,0.06])
    bar = plt.colorbar(cs,cax = cax,orientation='horizontal',
                        ticks=v,drawedges=False)
   
    bar.ax.tick_params(labelcolor='k',labelsize=8,pad=0.5)
    bar.outline.set_linewidth(0.3)
    bar.ax.set_xticklabels(v,fontsize=8,fontname="Times New Roman")
   

    fig.text(0.90,0.028,"($m^3/m^3$)",fontsize=6) # bsf
    tail = 'ctrl'
    # tail = 'rand'+str(num)
    tail = 'diff'+str(num)
    plt.savefig(year+'_d'+str(day)+'_h'+str(hour)+'_z'+str(level)+tail+".png",dpi=800)

灭火器

没图没真相

MissLH88

D:\LHwork\wrfout\code\1988_d0_h0_z0ctrl.png

MissLH88

MissLH88 发表于 2023-11-28 12:37

之前画类似的图，海洋的地方都没有问题，就是这个时刻的数据就算把海洋都mask了，绘图的时候还是会有海洋。所以就尝试直接用cartopy的海洋填色试试看，居然可以。但是就是要注意那些需要在海洋外面的东西，比如经纬度的网格，还有我加红标注的区域这些绘制的zorder需要大于海洋填色的zorder，不然会被挡住。

好运人

置敦，插眼