-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathqflux_separate_vector_calc.ncl
95 lines (85 loc) · 3.09 KB
/
qflux_separate_vector_calc.ncl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
; Program
; This NCL script uses vertical integration of qflux to calcaluate its divergent and
; rotational wind components, and then writes them to netCDF file
; History:
; i2015/1/29 MeteoBoy4 First Release
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"
begin
;-------Open file to read
dir="/run/media/MeteoBoy4/Data/MData/ERA-Interim"
system("rm -f "+dir+"/1985-2014_ALEV/Monthly/Calculation/Qflux_500hpa/vint_qflux_sepa_vector.nc")
efile=addfile(dir+"/1985-2014_ALEV/Monthly/Calculation/Qflux_500hpa/qu.nc","r")
nfile=addfile(dir+"/1985-2014_ALEV/Monthly/Calculation/Qflux_500hpa/qv.nc","r")
qu=efile->qu(:,0,::-1,:) ;Pay attentions to the dimensions !
qv=nfile->qv(:,0,::-1,:)
dims=dimsizes(qu)
ntime=dims(0)
nlat=dims(1)
nlon=dims(2)
latitude=qu&latitude
longitude=qu&longitude
time=qu&time
printVarSummary(qu)
print("Ready to calculate!")
;--------Calculate the divergent and rotational wind components
div=uv2dvF_Wrap(qu,qv)
print("Divergence Done!")
uvd=dv2uvF_Wrap(div)
print("Divergent Components Done!")
uvr=uvd
uvr(0,:,:,:)=qu-uvd(0,:,:,:)
uvr(1,:,:,:)=qv-uvd(1,:,:,:)
print("Rotational Components Derived!")
sfvp=uvd
sfvp=uv2sfvpF(qu,qv)
print("Functions Done!")
delete(qu)
delete(qv)
;-------Write the resulats into netCDF file
system("free")
fout=addfile(dir+"/1985-2014_ALEV/Monthly/Calculation/Qflux_500hpa/vint_qflux_sepa_vector.nc","c")
fileAtt=True
fileAtt@title="The divergent and rotational wind components of vertical integration of qflux on 500hpa"
fileAtt@Conventions="None"
fileAtt@creation_date=systemfunc("date")
setfileoption(fout,"DefineMode",True)
fileattdef(fout,fileAtt)
dimNames=(/"time","latitude","longitude"/)
dimSizes=(/ntime,nlat,nlon/)
dimUnlim=(/False,False,False/)
filedimdef(fout,dimNames,dimSizes,dimUnlim)
filevardef(fout,"time",typeof(time),"time")
filevardef(fout,"latitude",typeof(latitude),"latitude")
filevardef(fout,"longitude",typeof(longitude),"longitude")
filevardef(fout,"diu",typeof(uvd),getvardims(uvd(0,:,:,:)))
filevardef(fout,"div",typeof(uvd),getvardims(uvd(1,:,:,:)))
filevardef(fout,"rou",typeof(uvr),getvardims(uvr(0,:,:,:)))
filevardef(fout,"rov",typeof(uvr),getvardims(uvr(1,:,:,:)))
filevardef(fout,"sf",typeof(sfvp),getvardims(sfvp(0,:,:,:)))
filevardef(fout,"vp",typeof(sfvp),getvardims(sfvp(1,:,:,:)))
filevarattdef(fout,"time",time)
filevarattdef(fout,"latitude",latitude)
filevarattdef(fout,"longitude",longitude)
filevarattdef(fout,"diu",uvd(0,:,:,:))
filevarattdef(fout,"div",uvd(1,:,:,:))
filevarattdef(fout,"rou",uvr(0,:,:,:))
filevarattdef(fout,"rov",uvr(1,:,:,:))
filevarattdef(fout,"sf",sfvp(0,:,:,:))
filevarattdef(fout,"vp",sfvp(1,:,:,:))
fout->time=(/time/)
fout->latitude=(/latitude/)
fout->longitude=(/longitude/)
fout->diu=(/uvd(0,:,:,:)/)
fout->div=(/uvd(1,:,:,:)/)
fout->rou=(/uvr(0,:,:,:)/)
fout->rov=(/uvr(1,:,:,:)/)
fout->sf=(/sfvp(0,:,:,:)/)
fout->vp=(/sfvp(1,:,:,:)/)
print("Written to netCDF format!")
delete(div)
delete(uvd)
delete(uvr)
delete(sfvp)
end