脚本说明-version2

一、q-cross-ori.ncl

该脚本用于绘制比湿（g*kg-1）随时间的垂直剖面分布。根据定义的时间范围、空间范围，读取era5数据绘图。其中时间范围自定义，空间范围可以是该时间范围下暴雨点的分布（匹配晨杰的暴雨文件）、或台风位置（匹配台风路径文件）、或自定义。

        1.函数库挂载

;============（区域平均或单点）比湿随时间的垂直剖面分布=================
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/contrib/time_axis_labels.ncl"
load "/Users/yuxiaoyu/progm/time-cross/get-data-fuc.ncl"

最后一行挂载的get-data-fuc.ncl文件，路径更改为当前电脑该文件所在的路径。

        2.数据路径

;============================= file path ==========================================
;era5数据路径
cd = "/Users/yuxiaoyu/progm/data/";era5 data path
;图片输出路径
out_dir = "/Users/yuxiaoyu/progm/time-cross/"
;暴雨数据
file2 = cd + "3hPre100(2).csv";100mm以上
;台风数据
fiTY = cd + "201909.txt";201909利奇马

（建议把era5数据和晨杰计算的暴雨数据放在一个路径下，为cd路径）

cd = 当前电脑下era5数据路径

out_dir = 绘制的图文件输出路径

file2 = 晨杰计算的暴雨数据

fiTY = 台风路径数据，台风路径数据格式：

        3.数据自定义说明

时间范围：（era5数据、暴雨数据、台风数据）

;============================== data range=============================================
;定义需要的时间范围
initial_time = 2019080808                                                          ;start time 北京时
gap = ispan(0,3*24,6)                                                              ;时间间隔，总共3天，逐6h
;-----era5数据level层选取，-1为全选---
; level = (/1000,925,850,800,700,600,500,400,300,200/)
; level = -1
level = (/1000,950,900,850,800,700,600,500,400,300,200/)                           ;era5数据level层选取，-1为全选inter = 4   

initial_time：起始时间，自定义，格式yyyymmddhh，北京时。

gap：选取时间段，自定义，格式1：数组：ispan(st,ed,inter)；格式2：数值，读取的时间为initial_time+gap时刻。

level：垂直层。-1为era5原始垂直层全选。若不按原始数据垂直层则自定义。

inter：图中底部横坐标时间间隔（h）。

;rainstorm 数据
st_pp = get_data_rs(file2,initial_time,gap,"rain","3h")          

站点数据：用get_data_rs函数读取，get_data_rs(站点数据路径，起始时间，时间间隔范围，数据类型，3h或1h）。其中数据类型：rain，station，latlon，time。

此处用的是晨杰计算的3h暴雨数据，若用1h暴雨，最后的”3h“改为”1h"。

数据格式：st_pp:[num] 1维数据，定义的时间范围内出现3h100mm以上暴雨站点(num)的雨量

st_id、st_time同st_pp，st_loc：[2]×[num]：[0]lat，[1]lon。

台风数据：用get_data_ty函数读取，得到ty：[6]×[num]：[0]lat,[1]lon,[2]dir,[3]windspeed [4]u and [5] v with unit:m/s

定义时间范围内所有台风路径信息。

空间范围：（era5数据）

①自定义

area = -1

area：格式1(/lat0,lon0,lat1,lon1/)，表示选取（lat0-lat1，lon0-lon1）范围内数据；

格式2(/lat0,lat1,-1/)，表示选取（lat0-lat1）范围内所有经度数据；(/-1,lon0,lon1/)，表示选取（lon0-lon1）范围内所有纬度数据；

格式3(/lat,lon/)，表示选取（lat,lon）单点数据；

格式4 area= -1，表示用era5原始数据范围。

②根据暴雨站点经纬度选取era5数据范围

A.若选取暴雨站点范围：

利用前面已经读取的定义时间内暴雨站点经纬度信息：st_loc

lat0 = min(st_loc(0,:))
lat1 = max(st_loc(0,:))
lon0 = min(st_loc(1,:))
lon1 = max(st_loc(1,:))
;范围
area = (/lat0,lon0,lat1,lon1/);

B.若选取暴雨站点单点：

station_id = "k3004"
indid = str_match_ind_ic(st_id,station_id)
if(any(ismissing(indid)))thenprint("Error: There is no '"+ station_id+"' station in rainstorm data!!")exit
elseif(dimsizes(indid).gt.1)thenlatlon = st_loc(:,indid(0))elselatlon = st_loc(:,indid)end if
end if
print(latlon)
;单站
area = (/latlon(0),latlon(1)/)    ;

③根据台风路径的经纬度选取era5数据范围

A.若选取台风该时段内路径的范围：

lat0 = min(ty(0,:))
lat1 = max(ty(0,:))
lon0 = min(ty(1,:))
lon1 = max(ty(1,:))
area = (/lat0,lon0,lat1,lon1/);

B.若选取某一时刻台风位置单点：

time_pick = "1008"
indid = str_match_ind_ic(ty_time,time_pick)
if(any(ismissing(indid)))thenprint("Error: There is no "+ time_pick+" time in ty data!!")exit
elseif(dimsizes(indid).gt.1)thenlat_ty = ty(0,indid(0))lon_ty = ty(1,indid(0))elselat_ty = ty(0,indid)lon_ty = ty(1,indid)end if
end if
;;单点
area = (/lat_ty,lon_ty/)    ;

每种定义area方式不共存，选取一种方式则需要隐掉其余方式，最后有且仅有一个area。

若用范围，后面data calculate 模块中用区域平均比湿，隐掉单点；若用单点，则用单点比湿，并隐掉另一种。

;区域平均比湿
q_area = wgt_areaave_Wrap(q(:,:,{lat0:lat1},{lon0:lon1}),1,1,0)   ;区域平均
q_area_d = q_area(lev|:,time|:)                       ;绘图格式;单点比湿
; q_point = q(:,:,{latlon(0)},{latlon(1)})              ;单点
; q_point_d = q_point(lev|:,time|:)    

        4.读取era5数据

u = get_data(cd,"u",initial_time,gap,level,area)

用get_data函数读取，get_data(era5数据路径，数据类型，起始时间，时间间隔，垂直层次，空间范围)

起始时间可以与前面一致也可自定义，注意格式为北京时，时间间隔同理。

其中数据类型：

“t2m","u10","v10","cape","convec_pp","div","hgt","large_pp","mean_tn_lwr","mslp","pv","rh","sst","q" ,"s_lhf","surface_net_sr","surface_net_sr_clear","surface_p","s_shf","t","total_pp","u","v","w","vort"

对应era5数据库里的数据。

        5.绘图模块

wks =   gsn_open_wks("png",out_dir+"q-cross")

图片文件名”q-cross“，可自定义。

绘图resource基本不用更改，需注意若level层次有所更改，不是-1且不是原始脚本，则该行resource需更改：

res@tmYLLabels = (/"1000"," ","900"," ","850"," ","700","600","500","400","300","200"/)

该行表示绘图中纵坐标垂直层的标注，大小应与level一致，不标注的层次用” “代替。

二、tlogp-ori.ncl

该脚本可选取站点或自定义位置绘制T-LOGP图。

        1.函数库挂载        

与q-cross-ori脚本一样，更改get-data-fuc文件路径。        

        2.数据自定义

和q-cross-ori脚本一致。

        3.era5数据读取

u = get_data(cd,"u",initial_time,gap2,level,area)

t-logp图为单一时刻图，故但gap2非数组，必须为数值，initial_time可与前面一致也可自行定义。

其余无需变动。

三、出图结果

q-cross运行结果应出如下图：

横坐标：8日08时-11日08时，每隔6小时

t-lop应得到如下：

 站点k3004于2019年8月8日14时的压力-对数图。蓝线：露点廓线，黑线：温度廓线，红虚线：CAPE廓线