三門峽建設銀行網(wǎng)站緬甸今日新聞
問題
當繪圖時,往往并不需要繪制整塊區(qū)域,而是想聚焦于局部地區(qū),此時我們需要繪制扇形圖。
在cartopy中,只提供普通正方形的框架,如果我們需要其他,邊界,需要自己去繪制,最常見的是使用set_boundary繪制邊界,set_extend規(guī)定邊界,如:
latmin = 10
latmax = 60
lonmin = 70
lonmax = 140
lats = np.linspace(latmax, latmin, latmax - latmin + 1)
lons = np.linspace(lonmin, lonmax, lonmax - lonmin + 1)
vertices = [(lon, latmin) for lon in range(lonmin, lonmax + 1, 1)] + \[(lon, latmax) for lon in range(lonmax, lonmin - 1, -1)]
boundary = mpath.Path(vertices)
ax.set_boundary(boundary, transform=ccrs.PlateCarree())
ax.set_extent(extent)
需要注意的是,在使用set_boundary裁剪邊界后,需要使用set_extend規(guī)定邊界,否則會出現(xiàn)如下情況:
但是,這種方法有著自己的bug與缺陷:比如
- 需要自己添加網(wǎng)格和經(jīng)緯度標簽
- set_boundary和 set_extend一起使用時,繪制到邊界消失,如:cartopy二次曲線外觀保持
- 添加坐標時可能會有A LinearRing must have at least 3 coordinate tuples的報錯。
解決方式
實際上,這個問題原因還是由于投影轉換的問題,在set_extend時,繪制的上下邊界仍然是方形、未被正確投影的邊界,與我們的set_boundary存在沖突,最根本的原因還是在于cartopy對于投影計算的一些缺陷。
這類問題的解決在:GeoAxes.set_extent on non-cylindrical projections得到了充分的討論與解決,主要是使用了這幾行代碼:
lon1, lon2, lat1, lat2 = [-20, 20, 50, 80]
rect = mpath.Path([[lon1, lat1], [lon2, lat1],[lon2, lat2], [lon1, lat2], [lon1, lat1]]).interpolated(50)
proj_to_data = ccrs.PlateCarree()._as_mpl_transform(ax) - ax.transData
rect_in_target = proj_to_data.transform_path(rect)
將我們繪制的邊界進行了轉換。
我自己繪制到圖形示例如下:
import matplotlib.pyplot as plt
import matplotlib.path as mpath
import cartopy.crs as ccrs
import cartopy.mpl.ticker as ctk
import cartopy.feature as cfeature
import cmapscmap = cmaps.BlueDarkRed18
leftlon, rightlon, lowerlat, upperlat = [0,120,65,85]rect = mpath.Path([[leftlon, lowerlat], [rightlon, lowerlat],[rightlon, upperlat], [leftlon, upperlat], [leftlon, lowerlat]]).interpolated(50)proj=ccrs.NearsidePerspective(central_longitude=(leftlon+rightlon)*0.5,central_latitude=(lowerlat+upperlat)*0.5)
fig= plt.figure(figsize=(8,8),dpi=300)#設置畫布大小
ax= fig.add_axes([0.2,0.3,0.5,0.5],projection =proj)
#ax.set_extent(extent, ccrs.PlateCarree())#這樣截出來是方形的
ax.add_feature(cfeature.COASTLINE.with_scale('110m'))proj_to_data = ccrs.PlateCarree()._as_mpl_transform(ax) - ax.transData
rect_in_target = proj_to_data.transform_path(rect)
ax.set_boundary(rect_in_target)
ax.set_xlim(rect_in_target.vertices[:,0].min(), rect_in_target.vertices[:,0].max())
ax.set_ylim(rect_in_target.vertices[:,1].min(), rect_in_target.vertices[:,1].max())gl=ax.gridlines(draw_labels=True, x_inline=False, y_inline=False, linestyle='dashed')
gl.top_labels=False
gl.right_labels=False
gl.rotate_labels=False
gl.xlocator=ctk.LongitudeLocator(4)
gl.ylocator=ctk.LatitudeLocator(6)
gl.xformatter=ctk.LongitudeFormatter(zero_direction_label=True)
gl.yformatter=ctk.LatitudeFormatter()
ax.set_title('AMJ-pc1 & SON-SIC',loc='center',fontsize =12)c1 = ax.contourf(lon1,lat1, s, zorder=0,levels =np.arange(-0.09,0.12,0.03),extend = 'both',cmap=cmap,transform=ccrs.PlateCarree())
c1b =ax.contourf(lon1,lat1, p,[0,0.1 ,1], zorder=1,hatches=['...', None],colors="none",transform=ccrs.PlateCarree())position=fig.add_axes([0.2, 0.2, 0.5, 0.02])
fig.colorbar(c1,cax=position,orientation='horizontal')#,format='%.2f',)
plt.show()
另外我們需要指出的是:**該方法不適用于極地投影,即NorthPolarStereo,由于NorthPolarStereo本身投影特性只需一個參數(shù),本身并不適合。
**
極地局部繪制(不推薦)
我們繪制極地投影時,同樣也是使用set_boundary繪制圓形邊界,那么當我們想要繪制扇形時,可以通過只繪制一部分的圓形,通過調整繪制圓的參數(shù)來局部繪制。
詳細例子可見:極地局部圖繪制
在我進行實驗后發(fā)現(xiàn)它存在一些缺陷:
1、難以準確截出自己需要的區(qū)域
2、圖形位置不好確定。
3、同樣不好添加網(wǎng)格標簽。
該方法本質是:繪制整個圓形邊界,而只顯示部分區(qū)域,在繪圖時會給人不過完整的感覺。
直接繪制如下圖,需要自行裁剪或調整繪制:
請根據(jù)喜好自行選擇,