Py学习  »  Python

Python计算整层水汽通量散度

happy科研 • 2 年前 • 1234 次点击  

本文作者信息、相关参考及致谢见文末。


1 计算水汽通量散度

import matplotlib.pyplot as pltimport cartopy.crs as ccrsfrom cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatterfrom pylab import *import numpy as npimport xarray as xrimport metpy.calc as mpcalcimport metpy.constants as constantsimport cmapsfrom matplotlib import rcParamsconfig = {"font.family":'Times New Roman',"font.size":16,"mathtext.fontset":'stix'}rcParams.update(config)u=xr.open_dataset('F:/Rpython/lp36/data/water/uwnd.2018.nc')v=xr.open_dataset('F:/Rpython/lp36/data/water/vwnd.2018.nc')sh=xr.open_dataset('F:/Rpython/lp36/data/water/shum.2018.nc')lat=sh.latlon=sh.lonlev=sh.level[0:8]dx, dy = mpcalc.lat_lon_grid_deltas(lon, lat)div_qv_gc=np.zeros((3,8,73,144))qv_u_gc=np.zeros((3,8,73,144))qv_v_gc=np.zeros((3,8,73,144))# qv_u_gc1=np.zeros((8,73,144))# qv_v_gc1=np.zeros((8,73,144))# for j in pd.date_range(data1.loc[i]['start_day'],data1.loc[i]['end_day'],freq='D'):  fu=u['uwnd'].loc['2018-1-18':'2018-1-20',1000:300,:,:]fv=v['vwnd'].loc['2018-1-18':'2018-1-20',1000:300,:,:]# fu1=u['uwnd'].loc['2018-1-18':'2018-1-20',1000:300,:,:].mean(axis=0)# fv1=v['vwnd'].loc['2018-1-18':'2018-1-20',1000:300,:,:].mean(axis=0)fq=sh['shum'].loc['2018-1-18':'2018-1-20',1000:300,:,:]#kg/kg转换成g/kg# fq1=sh['shum'].loc['2018-1-18':'2018-1-20',1000:300,:,:].mean(axis=0)qv_u_gc[:,:,:,:]=fu*fqqv_v_gc[:,:,:,:]=fv*fq# 计算每次过程的每层水汽通量散度for j in range(3):    for k in range(lev.shape[0]):        div_qv_gc[j,k,:,:] = mpcalc.divergence(u = qv_u_gc[j,k,:,:]/9.8,v=qv_v_gc[j,k,:,:]/9.8,dx = dx ,dy = dy)#计算整层水汽通量散度div_total_gc=np.zeros((3,73,144))qv_u_gc_total=np.zeros((3,73,144))qv_v_gc_total=np.zeros((3,73,144))for l in range(3):    div_total_gc[l,:,:]= np.trapz(div_qv_gc[l,::-1],lev[::-1],axis=0)    qv_u_gc_total[l,:,:]=np.trapz(qv_u_gc[l,::-1],lev[::-1],axis=0)    qv_v_gc_total[l,:,:]=np.trapz(qv_v_gc[l,::-1],lev[::-1],axis=0)div_mean=div_total_gc.mean(axis=0)*10000000qu_mean= qv_u_gc_total.mean(axis=0)/9.8*100qv_mean= qv_v_gc_total.mean(axis=0)/9.8*100   div_ano1=xr.DataArray(div_mean , coords=[lat,lon], dims=['lat','lon'],name='div') div_ano1.to_netcdf('F:/Rpython/lp36/data/water/water2/vaper_div_2018.nc')     #计算整层散度通量qu_ano=xr.DataArray(qu_mean, coords=[lat,lon], dims=['lat','lon'],name='qu') qv_ano=xr.DataArray(qv_mean, coords=[lat,lon], dims=['lat','lon'],name='qv')qu_ano.to_netcdf('F:/Rpython/lp36/data/water/water2/qu_2018.nc')  qv_ano.to_netcdf('F:/Rpython/lp36/data/water/water2/qv_2018.nc')

2 作图

2.1 contourf着色作图

cf=ax.contourf(cycle_LON,cycle_LAT,cycle_div_mean,levels=levs,cmap=cmaps.NCV_bright,extend='both')


2.2 streamplot流线作图

sp=ax.streamplot(lon2,lat,qv_u_gc_total,qv_v_gc_total,color='black')


2.3 quiver作图

q=ax.quiver(LON[::4],LAT[::4],qv_u_gc


    
_total[1:73,:][::4],qv_v_gc_total[1:73,:][::4])


本文作者:南信大 cl

作者邮箱:1343293715@qq.com

特别致谢 南信大 摸鱼咯 博士。

参考链接:https://www.jianshu.com/u/9293eb1f7254



Python社区是高质量的Python/Django开发社区
本文地址:http://www.python88.com/topic/121736
 
1234 次点击