- 积分
- 40
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2018-11-5
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
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"
begin
; ;***********************************************
; ; get variable names from grib file
; ;***********************************************
; grib_in = addfile("./fnl_20190524_06_00.grib2","r")
; names = getfilevarnames(grib_in); extract all variable names
; ;***********************************************
; ; create output netcdf file
; ;***********************************************
; system("rm fnl_20190524_06_00.nc") ; remove any pre-existing file
; ncdf_out = addfile("fnl_20190524_06_00.nc" ,"c") ; create output netCDF file
; ;***********************************************
; ; loop through variables and output each to netcdf
; ;***********************************************
; do i = 0, dimsizes(names)-1
; ncdf_out->$names(i)$ = grib_in->$names(i)$
; end do
; ;***********************************************
a=15
b=60
c=75
d=125
;---Open file and read data
dir="/home/*/fnl/newtqxst/"
fils = systemfunc ("ls " + dir + "*.nc" )
f=addfiles(fils, "r")
ListSetType(f, "join")
print(f)
; filename = "fnl_20190524_06_00.nc"
; f = addfile(filename,"r")
; lon = f[:]->lon_0(:,{c:d})
; printVarSummary(lon)
; lat = f[:]->lat_0(:,{a:b})
lv=f[:]->lv_ISBL0
print(lv)
; h= a->HGT_P0_L100_GLL0(:,{15:60},{75:125})
; u = a->UGRD_P0_L100_GLL0(:,{15:60},{75:125})
; v = a->VGRD_P0_L100_GLL0(:,{15:60},{75:125})
; T = a->TMP_P0_L100_GLL0(:,{15:60},{75:125})-273.15
; P = a->PRES_P0_L1_GLL0({15:60},{75:125})/100
h= f[:]->HGT_P0_L100_GLL0(:,:,{a:b},{c:d})
u = f[:]->UGRD_P0_L100_GLL0(:,:,{a:b},{c:d})
v = f[:]->VGRD_P0_L100_GLL0(:,:,{a:b},{c:d})
; T = f->TMP_P0_L100_GLL0(:,{a:b},{c:d})-273.15
P = f[:]->PRES_P0_L1_GLL0(:,{a:b},{c:d})/100
;---Start the graphics
printVarSummary(v)
do i=0,2
v!2="lat"
v!3="lon"
lat= (/f[:]->lat_0(i,a:b)/)
lat@units="degrees_north"
lon=(/f[:]->lon_0(i,c:d)/)
lon@units="degrees_east"
v&lat=lat
v&lon=lon
u!2="lat"
u!3="lon"
lat= (/f[:]->lat_0(i,a:b)/)
lat@units="degrees_north"
lon=(/f[:]->lon_0(i,c:d)/)
lon@units="degrees_east"
u&lat=lat
u&lon=lon
P!1="lat"
P!2="lon"
lat= (/f[:]->lat_0(i,a:b)/)
lat@units="degrees_north"
lon=(/f[:]->lon_0(i,c:d)/)
lon@units="degrees_east"
P&lat=lat
P&lon=lon
timess=(/"NCEP GFS 1~S~o~N~*1~S~o~N~ VALID TIME: 2019/05/22/14hrs","NCEP GFS 1~S~o~N~*1~S~o~N~ VALID TIME: 2019/05/23/14hrs","NCEP GFS 1~S~o~N~*1~S~o~N~ VALID TIME: 2019/05/24/14hrs"/)
plot1 = new (3 , "graphic")
plot2 = new (3 , "graphic")
LZ_id = new (3 , "graphic")
txid1 = new (3 , "graphic")
; do j = 30,30
wks = gsn_open_wks("x11","072800850")
gsn_define_colormap(wks,"MPL_seismic")
;---Set some resources
res = True
res@cnLineLabelFont= "times-bold"
res@cnLineLabelFontHeightF=0.011
res@gsnDraw = False
res@gsnFrame = False
; res@tfDoNDCOverlay = True ;very important for overlay
res@gsnMaximize = True ; maximize plot in frame
res@cnFillOn =False ; turn on color fill
res@gsnAddCyclic = False
res@cnLinesOn = True
res@mpOutlineOn = True
res@mpDataBaseVersion = "MediumRes"
res@mpDataSetName = "Earth..4"
res@mpOutlineSpecifiers = "China:states"
res@mpOutlineBoundarySets = "AllBoundaries"
; res@mpFillOn =True
res@mpLandFillColor = "antiquewhite"
res@mpOceanFillColor = "cadetblue1"
res@mpMinLatF = a ; min lat
res@mpMaxLatF = b ; max lat
res@mpMinLonF = c ; min lon
res@mpMaxLonF = d ; max lon
res@mpGeophysicalLineThicknessF= 2.0
res@mpNationalLineThicknessF = 2.0
res@gsnLeftString = lv(0,30)/100+"hPa Wind&Height&T"
res@gsnRightString = "2015/07/28/0(UTC)"
res@gsnRightStringFontHeightF = .02
res@gsnLeftStringFontHeightF = .02
res@lbOrientation = "vertical"
res@lbTitleOn = True ; turn on title
res@lbTitleString = "Celsius degree"
res@lbTitleFontHeightF = .015 ; make title smaller
res@gsnStringFont = "times-bold"
res@pmTickMarkDisplayMode = "Always"
res@tmXBTickSpacingF = 40
res@tiXAxisString = timess(i)
res@tiXAxisFontHeightF = .015
res@tiXAxisFont = "times-bold"
res@tmXBLabelFont = "times-bold"
res@tmYLLabelFont = "times-bold"
; res@gsnDraw = False
; res@gsnFrame = False
res@gsnAddCyclic = False
; res@cnFillOn = False
; res@cnLineThicknessF = 1
; res@cnLineColor = "red"
; res@cnInfoLabelOn=False
; res@cnFillPalette="MPL_hot"
res@gsnLeftString=""
res@gsnRightString=""
; res@cnLevelSelectionMode = "ExplicitLevels" ;绘制等值线
; res@cnLevels = (/0,10,40/)
; res@cnFillColors = (/0,40,80,110,128/)
; plot0 = gsn_csm_contour_map(wks,T(i,:,:),res)
; res1=True
; res1=res
res@gsnDraw = False
res@gsnFrame = False
res@mpFillOn =True
; res@tfDoNDCOverlay = True ;very important for overlay
res@gsnAddCyclic = False
res@cnFillOn = False
res@cnLinesOn = True
res@cnLineThicknessF = 3
; res@cnLineLabelBackgroundColor=-1
res@cnSmoothingOn = True
res@cnLineColor = "black"
res@cnInfoLabelOn=False
res@gsnLeftString=""
res@gsnRightString=""
res@cnLineLabelPlacementMode="constant"
; res1@gsnMaximize=True
res2 = True ; plot mods desired
res2@gsnDraw = False
res2@gsnFrame = False
res2@gsnAddCyclic = False ;为了避免全球图中格林威治子午线的间隙
;res2@vcRefMagnitudeF = 10. ; make vectors larger
;res2@vcRefLengthF = 0.050 ; ref vec length
res2@vcGlyphStyle = "WindBarb" ; select wind barbs
res2@vcMinDistanceF = 0.025 ; thin out windbarbs
res2@vcRefMagnitudeF = 4.
res2@vcRefLengthF = 0.020
res2@vcWindBarbScaleFactorF = 2.5
res2@vcWindBarbLineThicknessF =2.
res2@gsnLeftString=""
res2@gsnRightString=""
res2@vcRefAnnoString2On = False
res2@lbLabelBarOn = False
res2@vcRefAnnoOn = False
plot1(i) = gsn_csm_contour_map(wks,P(i,::9,::9),res)
plot2(i)=gsn_csm_vector(wks,u(i,30,:,:),v(i,30,:,:),res2)
overlay(plot1(i),plot2(i))
resP = True ; modify the panel plot
resP@gsnMaximize = True ; maximize panel area
resP@lbLabelFontHeightF =0.035
resP@lbBoxEndCapStyle ="TriangleBothEnds"
resP@lbTitleFontHeightF=0.013
resP@gsnPanelLabelBar = True ; add common colorbar
resP@lbLabelFontHeightF=0.015
resP@lbTitleFontHeightF=0.021
resP@lbLabelFont = "times-bold"
resP@gsnPanelFigureStrings=(/"(a)","(b)","(c)"/)
resP@amJust="TopLeft"
resP@gsnPanelFigureStringsFontHeightF=0.017
resP@gsnPanelFigureStringsFont="times-bold"
; resP@gsnPanelFigureStringsPerimOn=False
; resP@gsnPanelFigureStringsBackgroundFillColor=-1
; resP@lbTitleString="Divergence"
resP@lbTitleFont ="times-bold"
resP@pmLabelBarOrthogonalPosF=-0.03
lat0 = (/36.072559/)
lon0 = (/103.82007/)
tstrs = (/"LANZHOU "/)
mkres = True
txres = True
mkres@gsMarkerIndex = 16 ; filled dot
mkres@gsMarkerColor = "green"
txres@txFontHeightF = 0.017
txres@txJust = "CenterRight"
txres@txFontColor = "green"
txres@txFont ="times-bold"
; mkid1 = gsn_add_polymarker(wks,plot0,lon0,lat0,mkres)
txid1(i) = gsn_add_text(wks,plot1(i),tstrs,lon0,lat0,txres)
LZ_shp_name="/home/*/wrf/WPSV3CP3/lanzhou/lanzhou.shp"
lnres = True
lnres@gsLineColor = "green"
lnres@gsLineThicknessF = 1.5
LZ_id(i)=gsn_add_shapefile_polylines(wks,plot1(i),LZ_shp_name,lnres)
; overlay(plot0,plot1)
end do
draw(plot1)
frame(wks)
gsn_panel(wks,plot1,(/1,3/),resP)
end
|
|