ncl run through all outputs

This post was from a previous version of the WRF&MPAS-A Support Forum. New replies have been disabled and if you have follow up questions related to this post, then please start a new thread from the forum home page.

N

Nefeli

New member
Dear NCL users,

I am currently new user of ncl. I manage to plot my wind and slp results from one of the output files, but I would like to plot them for each of the outputs. I understood that I have to use the addfiles but I am getting confused on how to loop through all the output files. Please see bellow what I am trying to do:

load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"

begin
;---Open the WRF data file
dir = "/home/c1735271/Dec2013TestRun/d01/"
files = systemfunc("cd "+ dir+" ; ls wrfout*")
a = addfiles(dir + files,"r")
times = wrf_user_list_times(a) ; get times in the file
it = 0 ; only one time step on this file

print("Working on time: " + times(it) )

;---Get the variables we need
slp = wrf_user_getvar(a,"slp",it) ; sea level pressure
wrf_smooth_2d( slp, 3 ) ; smooth slp
u = wrf_user_getvar(a,"ua",it) ; 3D U at mass points
v = wrf_user_getvar(a,"va",it) ; 3D V at mass points

u10 = u(0,:,:) ; Use lowest level at time 0
v10 = v(0,:,:)

utotal=sqrt(u10^2+v10^2)
delete([/u,v,u10,v10/]) ; Don't need these anymore

;---Attach some descriptive attributes

utotal@units = "m/s"
utotal@description = "Wind at 10m"

;---Open x11 window or png file
wks = gsn_open_wks("png","wrf_slp_uv")

;---Set some options common across all plots
cres = True
cres@MainTitle = "SLP_Wind"
cres@TimeLabel = times(it) ; Set Valid time to use on plots

;---Plot options for sea level pressure
res = cres
res@cnLineColor = "NavyBlue"
res@cnLineLabelBackgroundColor = -1 ; transparent
res@cnLineThicknessF = 2.0 ; default is 1.0
contour_slp = wrf_contour(a[0],wks,slp,res)
delete(res)

;---Plot options for wind countourf
res = cres
res@cnFillOn = True
contour_utot = wrf_contour(a[0],wks,utotal,res)
;res@FieldTitle = "Wind" ; overwrite Field Title
delete(res)

;---Overlay plots on map and draw
mpres = True
mpres@mpGeophysicalLineColor = "Black"
mpres@mpNationalLineColor = "Black"
mpres@mpUSStateLineColor = "Black"
mpres@mpGridLineColor = "Black"
mpres@mpLimbLineColor = "Black"
mpres@mpPerimLineColor = "Black"
mpres@mpGeophysicalLineThicknessF = 1.0
mpres@mpGridLineThicknessF = 1.0
mpres@mpLimbLineThicknessF = 1.0
mpres@mpNationalLineThicknessF = 0.6
mpres@mpUSStateLineThicknessF = 0.6
mpres@mpDataResolution = "FinestResolution"
pltres = True

ov = wrf_map_overlays(a[0],wks,(/contour_slp,contour_utot/),pltres,mpres)
end

Any help on the loop would be highly appreciated!

Thank you in advance,

Nefeli
 
Please see bellow an update of what I was trying:

load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"

begin

diri = "/home/c1735271/Dec2013TestRun/d01/"
fili = systemfunc("cd "+ diri+" ; ls wrfout*")
nfili = dimsizes(fili)
a=addfiles(diri+fili,"r")
times = wrf_user_list_times(a)
it=0

do nf=0,nfili-1
;---Get the variables we need
slp = wrf_user_getvar(a,"slp",it) ; sea level pressure
wrf_smooth_2d( slp, 3 ) ; smooth slp
u = wrf_user_getvar(a,"ua",it) ; 3D U at mass points
v = wrf_user_getvar(a,"va",it) ; 3D V at mass points

u10 = u(0,:,:) ; Use lowest level at time 0
v10 = v(0,:,:)

utotal=sqrt(u10^2+v10^2)
delete([/u,v,u10,v10/]) ; Don't need these anymore

;---Attach some descriptive attributes
utotal@units = "m/s"
utotal@description = "Wind at 10m"

;---Open x11 window or png file
wks = gsn_open_wks("png","wrf_slp_uv")

;---Set some options common across all plots
cres = True
cres@MainTitle = "SLP_Wind"
cres@TimeLabel = times(it) ; Set Valid time to use on plots

;---Plot options for sea level pressure
res = cres
res@cnLineColor = "NavyBlue"
res@cnLineLabelBackgroundColor = -1 ; transparent
res@cnLineThicknessF = 2.0 ; default is 1.0
contour_slp = wrf_contour(a[0],wks,slp,res)
delete(res)

;---Plot options for wind vectors
res = cres
res@cnFillOn = True
contour_utot = wrf_contour(a[0],wks,utotal,res)
delete(res)

;---Overlay plots on map and draw
mpres = True
mpres@mpGeophysicalLineColor = "Black"
mpres@mpNationalLineColor = "Black"
mpres@mpUSStateLineColor = "Black"
mpres@mpGridLineColor = "Black"
mpres@mpLimbLineColor = "Black"
mpres@mpPerimLineColor = "Black"
mpres@mpGeophysicalLineThicknessF = 1.0
mpres@mpGridLineThicknessF = 1.0
mpres@mpLimbLineThicknessF = 1.0
mpres@mpNationalLineThicknessF = 0.6
mpres@mpUSStateLineThicknessF = 0.6
mpres@mpDataResolution = "FinestResolution"
pltres= True

ov = wrf_map_overlays(a[0],wks,(/contour_slp,contour_utot/),pltres,mpres)
end do
end
 
Hi,

Your problem here is that you are not actually looping through all the files. You will want to add a do loop in order to make a plot. I have attached the code here that works for my files. I added a do loop the length of the times list. I loop through each file and make a new plot for each time.
Code: 
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"

begin
;---Open the WRF data file
dir = "./"
files = systemfunc("cd "+ dir+" ; ls wrfout_d01_2013*")
a = addfiles(dir + files+".nc","r")
times = wrf_user_list_times(a) ; get times in the file
it = dimsizes(times)

do i=0,it-1
;it = 0 ; only one time step on this file

print("Working on time: " + times(i) )

;---Get the variables we need
slp = wrf_user_getvar(a,"slp",i) ; sea level pressure
wrf_smooth_2d( slp, 3 ) ; smooth slp
u = wrf_user_getvar(a,"ua",i) ; 3D U at mass points
v = wrf_user_getvar(a,"va",i) ; 3D V at mass points

u10 = u(0,:,:) ; Use lowest level at time 0
v10 = v(0,:,:)

utotal=sqrt(u10^2+v10^2)
delete([/u,v,u10,v10/]) ; Don't need these anymore

;---Attach some descriptive attributes

utotal@units = "m/s"
utotal@description = "Wind at 10m"

;---Open x11 window or png file
wks = gsn_open_wks("X11","wrf_slp_uv_"+times(i))

;---Set some options common across all plots
cres = True
cres@MainTitle = "SLP_Wind"
cres@TimeLabel = times(i) ; Set Valid time to use on plots

;---Plot options for sea level pressure
res = cres
res@cnLineColor = "NavyBlue"
res@cnLineLabelBackgroundColor = -1 ; transparent
res@cnLineThicknessF = 2.0 ; default is 1.0
;contour_slp = wrf_contour(a[0],wks,slp,res)
contour_slp = wrf_contour(a[i],wks,slp,res)
delete(res)

;---Plot options for wind countourf
res = cres
res@cnFillOn = True
;contour_utot = wrf_contour(a[0],wks,utotal,res)
contour_utot = wrf_contour(a[i],wks,utotal,res)
;res@FieldTitle = "Wind" ; overwrite Field Title
delete(res)

;---Overlay plots on map and draw
mpres = True
mpres@mpGeophysicalLineColor = "Black"
mpres@mpNationalLineColor = "Black"
mpres@mpUSStateLineColor = "Black"
mpres@mpGridLineColor = "Black"
mpres@mpLimbLineColor = "Black"
mpres@mpPerimLineColor = "Black"
mpres@mpGeophysicalLineThicknessF = 1.0
mpres@mpGridLineThicknessF = 1.0
mpres@mpLimbLineThicknessF = 1.0
mpres@mpNationalLineThicknessF = 0.6
mpres@mpUSStateLineThicknessF = 0.6
mpres@mpDataResolution = "FinestResolution"
pltres = True

;ov = wrf_map_overlays(a[0],wks,(/contour_slp,contour_utot/),pltres,mpres)
ov = wrf_map_overlays(a[i],wks,(/contour_slp,contour_utot/),pltres,mpres)
end do
end
 
You must log in or register to reply here.
Back
Top