爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 10619|回复: 24

利用hysplit分析后向水汽轨迹

[复制链接]

新浪微博达人勋

发表于 2021-3-8 14:42:17 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
请教一下各位前辈,老师,和同学们,我参照王老师的帖子成功把ERA5的nc数据转换为驱动hysplit的arl格式,由于自己是要研究水汽轨迹,就新增了6小时降水;比湿,

轨迹

轨迹

密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2021-3-8 14:44:39 | 显示全部楼层
请问一下,这些轨迹就是水汽轨迹吗?
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2021-3-8 14:45:49 | 显示全部楼层
所以你的问题呢?最近我也在用hysplit,我是批量运行的hysplit里的ecm2arl.exe实现ERA5转arl的
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2021-3-8 14:48:40 | 显示全部楼层
但为什么图层属性里面没有比湿的变化
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2021-3-8 15:03:36 | 显示全部楼层
龍_彳亍_天_下 发表于 2021-3-8 14:45
所以你的问题呢?最近我也在用hysplit,我是批量运行的hysplit里的ecm2arl.exe实现ERA5转arl的

我的问题是,我画出来的这些轨迹代表的是水汽的轨迹吗?(ERA5数据里包含比湿和降水)但是单点每条轨迹的图层属性,只有经纬度和高度,没有属性变化,模式输出的文件里也没有比湿
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2021-3-8 15:43:39 | 显示全部楼层
这是我的输出文件
微信图片_20210308154309.png
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2021-3-12 20:30:54 | 显示全部楼层
还得找到每条轨迹对应的起始点水汽分析吧   PS:数据格式怎么转呢
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2021-3-16 10:19:58 | 显示全部楼层
淡水河 发表于 2021-3-12 20:30
还得找到每条轨迹对应的起始点水汽分析吧   PS:数据格式怎么转呢

王老师的脚本:
这是王老师的代码,我只把文件名改成我自己的
# Convert ERA5 GRIB data to ARL data

#---- Set data folder
datadir = r'E:/ERA5/down'
#---- Set output data file
outfn = os.path.join(datadir, 'test_era5_grib_diff1.arl')
#---- Read a GRIB data file
infn3d = addfile('{}/201901_press.grib'.format(datadir))
infn2d = addfile('{}/201901sur.grib  '.format(datadir))

print 'GRIB data file has been opened...'
#---- Set output ARL data file
arlf = addfile(outfn, 'c', dtype='arl')

#---- Set variable and level list
#---- Variable names in ERA5 data file
gvar2d = ['Surface_pressure_surface','2_metre_temperature_surface','10_metre_U_wind_component_surface',\
    '10_metre_V_wind_component_surface','Boundary_layer_height_surface','Convective_available_potential_energy_surface',\
    'Instantaneous_eastward_turbulent_surface_stress_surface','Instantaneous_northward_turbulent_surface_stress_surface']

gvar3d = ['Geopotential_isobaric','Temperature_isobaric','Vertical_velocity_isobaric',\
    'U_component_of_wind_isobaric','V_component_of_wind_isobaric','Relative_humidity_isobaric']

#---- Corresponding variable names in ARL data file
avar2d = ['PRSS','T02M','U10M','V10M','PBLH','CAPE','UMOF','VMOF']
avar3d = ['HGTS','TEMP','WWND','UWND','VWND','RELH']
#--- Add DIFF fields - difference between the original data and the packed data
avar3d_diff = list(avar3d)
avar3d_diff.append('DIFW')

#---- Set parameters of ARL data file
gv = infn3d['Geopotential_isobaric']
nx = gv.dimlen(gv.ndim - 1)
ny = gv.dimlen(gv.ndim - 2)
levels = gv.dimvalue(gv.ndim - 3)[::-1]
nz = len(levels)
arlf.setlevels(levels)
arlf.set2dvar(avar2d)
for l in levels:
    arlf.set3dvar(avar3d_diff)
arlf.setx(gv.dimvalue(gv.ndim - 1))
arlf.sety(gv.dimvalue(gv.ndim - 2))

#---- Write ARL data file
tNum = infn3d.timenum()
fhour = 0
for t in range(0, tNum):
    print 'Time index: ' + str(t)
    atime = infn3d.gettime(t)
    print atime.strftime('%Y-%m-%d %H:00')
    dhead = arlf.getdatahead(infn3d.proj, 'RSMC', 2, fhour)
    #Pre-write index record without checksum - will be over-write latter
    arlf.writeindexrec(atime, dhead)
    #Checksum list
    ksumlist = []
    # Write 2d variables
    ksums = []
    for avname,gvname in zip(avar2d, gvar2d):
        gdata = infn2d[gvname][t,:,:]
        if avname == 'PRSS':
            gdata = gdata * 0.01
        ksum = arlf.writedatarec(atime, 0, avname, fhour, 99, gdata)
        ksums.append(ksum)
    ksumlist.append(ksums)
    # Write 3d variables
    for lidx in range(0, nz):
        ksums = []
        llidx = nz - lidx - 1
        print(lidx,llidx)
        for avname,gvname in zip(avar3d, gvar3d):
            gdata = infn3d[gvname][t,llidx,:,:]
            if avname == 'WWND':
                gdata = gdata * 0.01
                difw = arlf.diff_origin_pack(gdata)
            elif avname == 'SPHU':
                gdata = gdata * 1000.
            elif avname == 'HGTS':
                gdata = gdata / 9.80665
            ksum = arlf.writedatarec(atime, lidx + 1, avname, fhour, 99, gdata)
            ksums.append(ksum)
        ksum = arlf.writedatarec(atime, lidx + 1, 'DIFW', fhour, 99, difw)
        ksums.append(ksum)
        ksumlist.append(ksums)
    #Re-write index record with checksum
    arlf.writeindexrec(atime, dhead, ksumlist)
    fhour += 1
arlf.close()
print 'Finished!'
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2021-9-7 10:38:35 | 显示全部楼层
shizhaoyang 发表于 2021-3-8 15:43
这是我的输出文件

请问您解决了这个问题吗?我的输出也没有看到比湿,但是在文献中有计算水汽源地贡献率的公式,是利用比湿计算的。请问您现在会了吗?可以帮忙解答一下吗,万分感谢!!
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2021-9-14 21:18:58 | 显示全部楼层
Sharon_Liu 发表于 2021-9-7 10:38
请问您解决了这个问题吗?我的输出也没有看到比湿,但是在文献中有计算水汽源地贡献率的公式,是利用比湿 ...

用的meteoinfo,将MeteoInfo_3.1.5\MeteoInfo\plugins\TrajStat\working中tm_sphu = 1,应该就可以得到比湿了
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表