Python画南极地图-借花献佛

对方的凤飞飞

看了Masterpiece的帖子http://bbs.06climate.com/forum.php?mod=viewthread&tid=95503,利用Python画高程图,之前我利用GMT画了一个南极高程地图(左图),所以想要利用python画南极地图,但是想要使用自己喜欢的任意,利用的是http://bbs.06climate.com/forum.php?mod=viewthread&tid=43521的colormap,我基本是在Masterpiece的帖子上进行修改,具体请看代码.但是有一个问题是m.drawmeridians设置为什么不能显示出-90和90,目前表示不会设置(右图),gmt这点就非常灵活.

1. import numpy as np
   
2. import matplotlib.pyplot as plt
   
3. from mpl_toolkits.basemap import Basemap
   
4. import netCDF4 as nc
   
5. import matplotlib as mpl
   
6. import cmaps
   
7. 
   
8. topo_file='ETOPO2v2g_f4.nc'
   
9. data=nc.Dataset(topo_file)
   
10. topo=data.variables['z'][:,:]
    
11. lat=data.variables['y'][:]
    
12. lon=data.variables['x'][:]
    
13. 
    
14. fig,ax = plt.subplots(figsize=(14,9)) #建立绘图平台
    
15. m = Basemap(projection='spstere',boundinglat=-60,lon_0=180,resolution='l')
    
16. m.drawcoastlines(linewidth=0.2, color='gray',zorder=3)
    
17. 
    
18. parallels = np.arange(-90,-50,10) #纬线
    
19. m.drawparallels(parallels,labels=[True,True,False,True])#第三个不能改,第一个,第二个和第四个无所谓
    
20. meridians = np.arange(0,360,45) #经线
    
21. m.drawmeridians(meridians,labels=[False,False,True,True])
    
22. lons, lats = np.meshgrid(lon,lat) #经纬度2维化
    
23. x, y = m(lons, lats) #投影映射
    
24. cmap=cmaps.WhiteBlueGreenYellowRed
    
25. cmap.set_under([255/255,255/255,255/255]) #低于0的填色为白色
    
26. lev=np.arange(0,4500,200)
    
27. cf=m.contourf(x,y,topo,levels=lev,cmap=cmap,extend='both')
    
28. cb=plt.colorbar(cf, ax=ax,shrink=0.52,aspect=30,pad=0.05,orientation='horizontal') #色标
    
29. cb.ax.tick_params(labelsize=10,pad=2,direction='in') #色标tick
    
30. 
    
31. plt.savefig('global_etopo.png',dpi=600,bbox_inches='tight') #存图
    
32. plt.show()

复制代码

同时也欢迎大神指教,怎么添加-90和90的坐标.如何去掉外面的方框呢?

Masterpiece

第一个问题，南/北半球投影有一个参数round可以设置和图1一样的圆盘图，就是程序中的第15行加入：
m = Basemap(projection='spstere',boundinglat=-60,lon_0=180,resolution='l',round=True)

第二个问题，因为原程序是没有第一个问题中的设置圆盘，所以axis是原始的方框，所以第21行：
m.drawmeridians(meridians,labels=[False,False,True,True])
之中的,labels=[False,False,True,True]还是会按照方框axis的位置去标注这个经度labels
并且这个labels后边的列表位置分别是[左，右，上，下]，所以你可以看到左右的两边是没有labels的，而左右两边正好是90W和90E的位置，于是这两个就空着了。解决办法也简单，把labels全部设为True即可。

如果建立投影的时候就已经弄了round=True的话，这时候在"左右"或"上下"的位置设置True即显示全部的经度刻度标签。

对方的凤飞飞

Masterpiece 发表于 2020-6-20 22:33
第一个问题，南/北半球投影有一个参数round可以设置和图1一样的圆盘图，就是程序中的第15行加入：
m = Bas ...

高,真的高.

不用回头

楼主，请问怎么把纬度值也标上去呀？