【ProPlot库】ProPlot3兰伯特投影-可添加刻度(三)
虽然 cartopy 下的 Plate Carrée 投影使用方便，但中纬度下使用 Lambert 投影能更好的呈现真实的地图。用一个正圆锥切于或割于球面，将地球面投影到圆锥面上，然后沿一母线展开成平面。下图是使用proplot绘制的最终效果：
在 proplot 中，可在以下链接找到相关投影名称表，其中兰伯特投影简称 'lcc'。
https://proplot.readthedocs.io/en/latest/api/proplot.constructor.Proj.html#proj-table
下面直观感受一下两种投影的不同：
Import proplot as plot
fig,axs = plot.subplots(ncols=2,width=5,height=3,projection=['cyl','lcc'])
fig.format(reso='hi',coast=True,metalinewidth=2,coastlinewidth=0.5)
axs[0].format(title='Equidistant Cylindrical (Plate Carrée)')
axs[1].format(title='Lambert Conformal Conic')
也可以使用 plot.Proj 来具体定义投影的参数和中央经线位置。同时可在format 里设置 labels=True (这是直接利用cartopy的 gridlines 方法)快速绘制坐标。
import proplot as plot
proj = plot.Proj('lcc', lon_0=105)
fig,axs = plot.subplots(ncols=2,width=5,height=3,projection=['cyl',proj])
fig.format(reso='hi',coast=True,metalinewidth=2,coastlinewidth=0.5,
lonlim=(80,130),latlim=(15,55),labels=True)
axs[0].format(title='Equidistant Cylindrical (Plate Carrée)')
axs[1].format(title='Lambert Conformal Conic')
也可以直接使用 cartopy 提供的 cartopy.mpl.geoaxes.GeoAxes.gridlines 方法添加 grid：
import proplot as plot
proj = plot.Proj('lcc', lon_0=105)
fig,axs = plot.subplots(ncols=1,width=5,height=3,projection=proj)
axs.format(reso='hi',coast=True,metalinewidth=2,coastlinewidth=0.5,
lonlim=(80,130),latlim=(15,55),title='Lambert Conformal Conic',titlesize=15,titlepad=10)
gl=axs[0].gridlines(
xlocs=np.arange(-180, 180 + 1, 10), ylocs=np.arange(-90, 90 + 1, 10),
draw_labels=True, x_inline=False, y_inline=False,
linewidth=0.5, linestyle='--', color='gray')
gl.top_labels = False
gl.right_labels = False
gl.rotate_labels =0
如果只想用 proplot 的 format 函数，也可以达到相同的效果
import proplot as plot
import numpy as np
proj = plot.Proj('lcc', lon_0=105)
fig,axs = plot.subplots(ncols=1,width=5,height=3,projection=proj)
axs.format(reso='hi',coast=True,metalinewidth=2,coastlinewidth=0.5,
lonlim=(80,130),latlim=(15,55),title='Lambert Conformal Conic',titlesize=15,titlepad=10,grid=True,
latlocator=np.arange(-180, 180 + 1, 10),labels=True,lonlabels=True,gridlabelsize=10,gridlabelpad=10,
gridlinestyle='--', gridcolor='grey',gridlinewidth=0.5)
#https://proplot.readthedocs.io/en/latest/api/proplot.axes.GeoAxes.html?highlight=geo#proplot.axes.GeoAxes
在上述例子中，gridlabelsize=10 , gridlabelpad=10 分别控制经纬度标签的大小和离axis的距离；gridlinestyle='--', gridcolor='grey' , gridlinewidth=0.5 来改变 gridline 的风格。追求快速出图的朋友可以采用这种方案绘制底图。
添加自定义刻度
虽然 cartopy 提供了 gridlines 的方法，但经度的标签出现在了 y 轴上（80°E），并且和matplotlib 的 matplotlib.axes.Axes.set_xticks 等并不相融，难以给出 tickmarker ；同时，ncl 可以为顶部添加 ticks label。网上的解决方案并不多，自己尝试了很多办法。下面提供一种可以给 cartopy 兰伯特投影添加 major、minor ticks 的方法(略微有些复杂)。
该方法部分参照：
https://zhajiman.github.io/post/cartopy_lambert/
https://gist.github.com/ajdawson/dd536f786741e987ae4e代码如下：
def add_alt(self,**kwargs):
locator = self._make_inset_locator([0, 0, 1, 1], self.transAxes)
ax = plot.CartesianAxes(self.figure, locator(self, None).bounds, **kwargs)
ax.set_axes_locator(locator)
cx=self.add_child_axes(ax) # to facilitate tight layout
cx.grid(False)
limsx=self.get_xlim()
limsy=self.get_ylim()
cx.set_xlim(limsx)
cx.set_ylim(limsy)
return cx
def find_side(ls, side):
"""
Given a shapely LineString which is assumed to be rectangular, return the
line corresponding to a given side of the rectangle.
"""
minx, miny, maxx, maxy = ls.bounds
points = {'left': [(minx, miny), (minx, maxy)],
'right': [(maxx, miny), (maxx, maxy)],
'bottom': [(minx, miny), (maxx, miny)],
'top': [(minx, maxy), (maxx, maxy)],}
return sgeom.LineString(points[side])
def lambert_xticks(ax,ax2, ticks,po='bottom',formatter = LongitudeFormatter()):
"""Draw ticks on the bottom x-axis of a Lambert Conformal projection."""
""" ax用于计算交点，ax2用于呈现ticks """
te = lambda xy: xy[0]
lc = lambda t, n, b: np.vstack((np.zeros(n) + t, np.linspace(b[2]+0.1*(b[2]-b[3]), b[3]-0.1*(b[2]-b[3]), n))).T
xticks, xticklabels = _lambert_ticks(ax, ticks, po, lc, te)
#ax.xaxis.tick_bottom()
ax2.set_xticks(xticks)
#print(xticks)
#formatter = LongitudeFormatter()
ax2.grid(False)
ax2.set_xticklabels([formatter(xtick) for xtick in xticklabels])
def lambert_yticks(ax,ax2, ticks,po='left',formatter=LatitudeFormatter()):
"""Draw ricks on the left y-axis of a Lamber Conformal projection."""
te = lambda xy: xy[1]
lc = lambda t, n, b: np.vstack((np.linspace(b[0], b[1], n), np.zeros(n) + t)).T
yticks, yticklabels = _lambert_ticks(ax, ticks, po, lc, te)
#ax.yaxis.tick_left()
ax2.set_yticks(yticks)
#formatter = LatitudeFormatter()
ax2.set_yticklabels([formatter(ytick) for ytick in yticklabels])
ax2.grid(False)
def _lambert_ticks(ax, ticks, tick_location, line_constructor, tick_extractor):
"""Get the tick locations and labels for an axis of a Lambert Conformal projection."""
outline_patch = sgeom.LineString(ax.spines['geo'].get_path().vertices.tolist())
axis = find_side(outline_patch, tick_location)
n_steps = 30
extent = ax.get_extent(ccrs.PlateCarree())
_ticks = []
for t in ticks:
xy = line_constructor(t, n_steps, extent)
proj_xyz = ax.projection.transform_points(ccrs.Geodetic(), xy[:, 0], xy[:, 1])
xyt = proj_xyz[..., :2]
ls = sgeom.LineString(xyt.tolist())
locs = axis.intersection(ls)
if not locs:
tick = [None]
else:
tick = tick_extractor(locs.xy)
#print(tick[0])
_ticks.append(tick[0])
# Remove ticks that aren't visible:
ticklabels = copy(ticks)
while True:
try:
index = _ticks.index(None)
except ValueError:
break
_ticks.pop(index)
ticklabels.pop(index)
return _ticks, ticklabels
def lambert_minorxticks(ax,ax2, ticks,po='bottom',formatter = LongitudeFormatter()):
"""Draw ticks on the bottom x-axis of a Lambert Conformal projection."""
te = lambda xy: xy[0]
lc = lambda t, n, b: np.vstack((np.zeros(n) + t, np.linspace(b[2]+0.1*(b[2]-b[3]), b[3]-0.1*(b[2]-b[3]), n))).T
xticks, xticklabels = _lambert_ticks(ax, ticks, po, lc, te)
ax2.grid(False)
ax2.xaxis.set_minor_locator(mticker.FixedLocator(xticks))
def lambert_minoryticks(ax,ax2, ticks,po='left',formatter=LatitudeFormatter()):
"""Draw ricks on the left y-axis of a Lamber Conformal projection."""
te = lambda xy: xy[1]
lc = lambda t, n, b: np.vstack((np.linspace(b[0], b[1], n), np.zeros(n) + t)).T
yticks, yticklabels = _lambert_ticks(ax, ticks, po, lc, te)
ax2.yaxis.set_minor_locator(mticker.FixedLocator(yticks))
ax2.grid(False)
#以下是画图的代码:
proj = plot.Proj('lcc', lon_0=105)
proj2 = plot.Proj('lcc', lon_0=115)
fig, axs = plot.subplots(ncols=2,refwidth=5,refheight=5.5,projection=proj,outerpad=4,share=False)
xticks = np.arange(65, 185+1, 20).tolist()
yticks = np.arange(0, 50+1, 10).tolist()
xmiticks = np.arange(65, 185+1, 4 ).tolist()
ymiticks = np.arange(0, 50+1, 2.5).tolist()
xticks2 = np.arange(65, 185, 10)
yticks2 = np.arange(0, 50, 10)
axs.format(grid=True,gridalpha=0.5,gridcolor='grey',lonlim=(80,130),latlim=(15,55),reso='hi',labels=False,coast=True,metalinewidth=2,
lonlocator=xticks,latlocator=yticks
,gridlinestyle='--',gridlinewidth = 2)
fig.canvas.draw() #在proplot中，如果没有插入shape地图，则要加上这句，原因未知（matplotlib不用）
ax=axs[0]
ax.format(title='without top and right ticks',titlesize=20,titlepad=10)
lambert_xticks(ax,ax, xticks)
lambert_yticks(ax,ax,yticks)
lambert_minorxticks(ax,ax,xmiticks)
lambert_minoryticks(ax,ax,ymiticks)
ax.tick_params(which='major', width=1.5, length=8,bottom=True,left=True,right=False,top=False , )
ax.tick_params(labelsize=15, pad=5)
ax.tick_params(which='minor', width=0.8,length=4, color='k',bottom=True,left=True,right=False,top=False ,)
#添加南海子图
ins = ax.inset_axes([0.825,0.001,0.2,0.2],projection =proj2) #proplot中inset_axes可以自己选择子图的projection，非常便捷
ins.format(grid=True,coast=True,lonlim=(105,125),latlim=(0,25),metalinewidth=1.5,reso='hi',coastlinewidth=0.2,
lonlocator=xticks2, latlocator=yticks2,gridlinestyle='--',gridlinewidth = 1)
ax=axs[1]
lambert_xticks(ax,ax, xticks)
lambert_yticks(ax,ax,yticks)
lambert_minorxticks(ax,ax,xmiticks)
lambert_minoryticks(ax,ax,ymiticks)
#添加secondary axis （CartesianAxes projection） 在matplotlib中可直接用ax.secondary_xaxis 和 ax.secondary_yaxis
#也可以直接ax2 = ax.inset_axes([0.,0.00,1,1],projection =proj)
ax2=add_alt(ax)
lambert_xticks(ax,ax2, xticks,'top') #由ax判断交点，在ax2上绘制
lambert_yticks(ax,ax2, yticks,'right')
lambert_minorxticks(ax,ax2,xmiticks,'top')
lambert_minoryticks(ax,ax2,ymiticks,'right')
ax.tick_params(which='major', width=1.5, length=8,bottom=True,left=True,right=False,top=False , )
ax.tick_params(labelsize=15, pad=5)
ax.tick_params(which='minor', width=0.8,length=4, color='k',bottom=True,left=True,right=False,top=False ,)
ax2.tick_params(which='major', width=1.5, length=8,bottom=False,left=False,right=True,top=True,color='k' ,labelleft=False,labelbottom=False
,labelright=True, labeltop=True)
ax2.tick_params(labelsize=15, pad=5)
ax2.tick_params(which='minor',width=0.8, length=4,bottom=False,left=False,right=True,top=True,color='k')
#添加南海子图
ins = ax.inset_axes([0.825,0.001,0.2,0.2],projection =proj2)
ins.format(grid=True,coast=True,lonlim=(105,125),latlim=(0,25),metalinewidth=1.5,reso='hi',coastlinewidth=0.2,
lonlocator=xticks2, latlocator=yticks2,gridlinestyle='--',gridlinewidth = 1)
利用 Shapley 判断网格线与一个轴线相交的位置，进而找到刻度的位置。传入set_xticks 和 set_yticks 中，再利用ax.set_xticklabels 添加刻度。在proplot中,目前我只找到了用ax.tick_params添加geoaxes 的tickmarker 的例子，在format 里给地图设置刻度标签的方法可能还未实现。
目前网上较多的案例是添加左侧和底部的majortick，而由于 top 和right 轴的坐标范围发生变化，于是模仿Secondary Axis 方法，叠加新的“图层”来为顶部和右侧添加处于不同位置的 ticks，add_alt(self,**kwargs) 是专门给 proplot 改编的，如果想用matplotlib 绘制，便可以直接采用axes.Axes.secondary_xaxis 和axes.Axes.secondary_yaxis ，为geoAxes 添加一个笛卡尔坐标系的轴，再使用 lambert_xticks 函数添加标签。
在实际画图中也可以只选择采用左图的方案。另外，proplot 的inset_axes 方法可以添加任意投影(proj参数)的小子图，比matplotlib 的inset_axes 方便很多。
再添加如下代码，完成等值线和 colorbar 的绘制；当然，也可也不画顶部和右侧的 tick，或是 ticklabel 和 tickmarker。
m=ax.contourf(T[0],cmap=cmaps.BlueWhiteOrangeRed,lw=0.2,levels=12,vmax=316,vmin=264)
cb=ax.colorbar(m,loc='b',format='%.0f',pad=0.0, shrink=1 ,drawedges=True,width=0.2,linewidth=1,
ticklabelsize=12,label='')
m2=ins.contourf(T[0],cmap=cmaps.BlueWhiteOrangeRed,lw=0.2,levels=12,vmax=316,vmin=260)
ax.format(grid=False,ltitle='2m Temperature',titlesize=15,titlepad=15,rtitle='(K)',suptitle='Your Second Plot' ,suptitlesize=20,suptitlepad=15)
本文参与 腾讯云自媒体同步曝光计划，分享自微信公众号。
原始发表：2021-12-20，如有侵权请联系 cloudcommunity@tencent.com 删除
评论
登录后参与评论
推荐阅读
腾讯云开发者
