Python读取NC数据画风场和急流区

jl2587t

Python读取NC数据画风场和急流区Python读取NC数据画风场和急流区
作者今天不想说话，注意点都在代码注释里了，各位看官凑活看哈

1. #知识点   1.多图共用一个色标
   
2. #         2.四张一样的图做循环，不用一张一张画
   
3. 
   
4. #共用色标的原理，其实就是多加一个子图，把色标放子图里，然后调整一下位置就好了
   
5. 
   
6. 
   
7. #循环的原理就是建立四个子图，同样沟通步骤循环四遍
   
8. #因为四张图是一摸一样的，只有数据不同，所以通过循环4次直接产生四张图
   
9. 
   
10. 
    
11. 
    
12. 
    
13. #处理数据的包
    
14. import xarray as xr
    
15. import numpy as np
    
16. #画图的包
    
17. import matplotlib.pyplot as plt
    
18. #地图的包
    
19. import cartopy.crs as ccrs
    
20. import my_class
    
21. 
    
22. 
    
23. 
    
24. 
    
25. 
    
26. #打开nc文件
    
27. #这些都是只需要执行一次的东西
    
28. ds = xr.open_dataset('download_201309.nc')
    
29. lat = ds.latitude
    
30. lon = ds.longitude
    
31. time = ds.time
    
32. u = ds['u']
    
33. v = ds['v']
    
34. lat_range = lat[(lat>=22) & (lat<=44)]
    
35. lon_range = lon[(lon>=102) & (lon<=135)]
    
36. fig = plt.figure(figsize=(18,9))
    
37. 
    
38. #设置4个子图，并且放到列表里面，方便循环的时候用
    
39. plt.subplots_adjust(left=0.07, right=0.90, top=0.95, bottom=0.05,wspace=0.0009,hspace=0.1)
    
40. ax_a = fig.add_subplot(2,2,1,projection = ccrs.PlateCarree())
    
41. ax_b = fig.add_subplot(2,2,2,projection = ccrs.PlateCarree())
    
42. ax_c = fig.add_subplot(2,2,3,projection = ccrs.PlateCarree())
    
43. ax_d = fig.add_subplot(2,2,4,projection = ccrs.PlateCarree())
    
44. 
    
45. ax_list = [ax_a,ax_b,ax_c,ax_d]
    
46. abcd = ['(a)','(b)','(c)','(d)']#图右上角的abcd
    
47. i = 0#循环变量
    
48. 
    
49. 
    
50. 
    
51. 
    
52. #开始循环，以时间为循环的变量，i作为辅助
    
53. for time_i in time:
    
54.     #取数据
    
55.     u_region = u.sel(longitude=lon_range, latitude=lat_range,time = time_i)
    
56.     v_region = v.sel(longitude=lon_range, latitude=lat_range,time = time_i)
    
57.     #算风速
    
58.     wind_speed = np.sqrt(u_region * u_region + v_region * v_region)
    
59.     #清洗风速，小于24的都置零
    
60.     e = (wind_speed > 23)
    
61.     wind_speed_region = np.where(e,wind_speed,0)
    
62. 
    
63. 
    
64.     #i的作用来了，利用i得到相应的子图
    
65.     ax_list[i].set_xticks(np.arange(102,136,3),crs=ccrs.PlateCarree())
    
66.     ax_list[i].set_yticks(np.arange(22,45,2),crs=ccrs.PlateCarree())
    
67.     xticks_str = ['102', '105', '108', '111', '114', '117', '120','123', '126', '129', '132', '135°E']
    
68.     ax_list[i].set_xticklabels(xticks_str,fontsize = 11)
    
69.     yticks_str = ['22   ','24    ','26    ','28    ','30    ','32    ','34    ','36    ','38    ','40    ','42    ','44°N']
    
70.     ax_list[i].set_yticklabels(yticks_str,fontsize = 11)
    
71.     my_class.readshapefile('bou2_4l.shp',linewidth=1,ax=ax_list[i])
    
72. 
    
73.     #填色
    
74.     #色标可以自己选择风格，也可以自定义，在之前的代码里面有
    
75.     cf_rh = ax_list[i].contourf(lon_range,lat_range ,wind_speed_region,cmap='hot_r',extend='both')
    
76. 
    
77.     #风杆
    
78.     ax_list[i].barbs(lon_range[::4],lat_range[::4],u_region[::4,::4],v_region[::4,::4],barbcolor=['k'],linewidth=0.5,length=5,barb_increments=dict(half=2, full=4, flag=20))
    
79. 
    
80.     #写abcd
    
81.     ax_list[i].text(0.05, 0.95, abcd[i],transform=ax_list[i].transAxes, fontsize=11)
    
82. 
    
83.     #一次循环结束 i+1
    
84.     i += 1
    
85. 
    
86. 
    
87. 
    
88. 
    
89. #色条额外添加一个子图，【左侧距离，下侧距离，宽度，高度】
    
90. position=fig.add_axes([0.90,0.3,0.01,0.4])
    
91. cb=fig.colorbar(cf_rh,cax=position,orientation='vertical')
    
92. font = {'family' : 'serif',#字体
    
93.         'color'  : 'k',#颜色
    
94.         'weight' : 'normal',#自己百度吧
    
95.         'size'   : 11,#字体大小
    
96.         }
    
97. cb.set_label('    m/s',rotation=360,fontdict=font)
    
98. 
    
99. 
    
100. 
     
101. 
     
102. plt.show()

复制代码

原帖位于：https://blog.csdn.net/weixin_423 ... 1001.2014.3001.5502

sms

学习学习{:eb502:}

紫藤萝

my_class在哪里呢，谢谢

luckyxiaogaga

{:eb502:}代码简洁，效率高，感谢楼主分享

bahzhong

my_class可以共享吗，谢谢！

sjysjxyf66

感谢分享！厉害厉害！

ironchen

感谢分享！!!

ゼヅ欢チニ

多谢代码，学习了