[ncl-talk] wrf-ncl titles and subtitles
Zilore Mumba
zmumba at gmail.com
Thu Apr 7 13:13:12 MDT 2022
Firstly, my sincere apologies for this very long mail.
Coming from a weather forecasting background, it is preferable to have the
title on the graphic as, e.g.
WRF Surface Pressure, Winds and Temperatures
Valid at 06Z on 2022-04-07, From 00Z on 2022-04-07"
This quickly tells the forecaster 1) what he/she is looking at, 2) the
validity time and date, and 3) the model initiation time and date.
The current titling and subtitling in wrf-ncl scripts is not convenient for
a forecasting setup.
In the script below I have stripped off the ncl hard-coded titles, mainly
with ”opts at FieldTitle = " " and opts at PlotLevelID = " ", and then setting
all units to “”
Not knowing how to extract part of a file name in ncl, I have used perl to
extract yyyy-mm-dd and hh from wrfout_d01_yyy-mm-dd_hh:mm:ss. These are
written to files and are like below
Variable: data1
Type: string
Total Size: 104 bytes
13 values
Number of Dimensions: 1
Dimensions and sizes: [13]
Coordinates:
Number Of Attributes: 1
_FillValue : missing
(0) 2022-03-28
(1) 2022-03-28
(2) 2022-03-28
(3) 2022-03-28
(4) 2022-03-28
(5) 2022-03-28
(6) 2022-03-28
(7) 2022-03-28
(8) 2022-03-28
(9) 2022-03-28
Variable: data2
Type: string
Total Size: 104 bytes
13 values
Number of Dimensions: 1
Dimensions and sizes: [13]
Coordinates:
Number Of Attributes: 1
_FillValue : missing
(0) 00
(1) 01
(2) 02
(3) 03
(4) 04
(5) 05
(6) 06
(7) 07
(8) 08
(9) 09
(10) 10
(11) 11
(12) 12
If anyone can assist
1. how to read data1 and data2 into 1d subscripted arrays
2. how to correctly position the main title, if possible
I know ncl scripts are being migrated to python, so this may not be
interesting.
; Example script to produce plots for a WRF real-data run,
; with the ARW coordinate dynamics option.
; Interpolating to specified pressure levels
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
begin
;
; The WRF ARW input file.
DATADir = "/home/zmumba/DA/OUTPUT/2022022800/noda/"
FILES = systemfunc (" ls -1 " + DATADir + "wrfout_d01* ")
numFILES = dimsizes(FILES)
files_all = addfiles(FILES+".nc","r")
; Perl script takes the files and extracts yyyy-mm-dd and hh
;-----------------------------------------------------------
system("perl /home/zmumba/DA/SCRIPTS/03_Ncl_Scripts/sort_files.pl")
data1 =
asciiread("/home/zmumba/DA/SCRIPTS/03_Ncl_Scripts/datesfile.txt",-1,"string")
data2 =
asciiread("/home/zmumba/DA/SCRIPTS/03_Ncl_Scripts/timesfile.txt",-1,"string")
; Need to read files above into a 1d array to be used in title
;,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,
;Code below does not work
;-------------------------
; Init_date = date1(0)
; Init_time = date2(0)
;,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,
;do i = 0, numFILES-1
; valid_date(i) = data1(i)
; valid_time(i) = data2(i)
; print(valid_time(i))
;end do
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; We generate plots, but what kind do we prefer?
; type = "x11"
; type = "pdf"
; type = "ps"
; type = "ncgm"
type = "png"
; wks = gsn_open_wks(type,"plt_PressureLevel1")
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; What times and how many time steps are in the data set?
; get time information and strip out the day and hour
times_in_file = files_all[:]->Times
dims = dimsizes (times_in_file)
ntimes = dims(0)
times_to_use = new(dims(0),string) ;create empty array
do i=0,dims(0)-1
times_to_use(i) = chartostring(times_in_file(i,8:12))
print(times_to_use(i))
end do
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; Set some Basic Plot options
res = True
;res at NoHeaderFooter = True
res at MainTitle = "WRF Surface Pressure, Winds and Temperatures ~C~ Valid
at Valid_time(i)Z on Valid_date(i) From Init_time Z on Init_date"
res at tiMainOffsetXF = 6.0
res at InitTime = False
res at Footer = False
pltres = True
mpres = True
mpres at mpFillOn = False ; turn off gray fill
mpres at mpOutlineBoundarySets = "National" ; turn on country
boundaries
mpres at mpGeophysicalLineColor = "Navy" ; color of cont. outlines
mpres at mpGeophysicalLineThicknessF = 1.5 ; thickness of outlines
mpres at mpNationalLineColor = "Black"
mpres at mpNationalLineThicknessF = 2.5 ; turn on country
boundaries
mpres at mpDataBaseVersion = "MediumRes" ; choose higher
resolution
mpres at mpDataSetName = "Earth..4" ; choose most recent
boundaries
mpres at mpMaxLatF = -4 ; choose subregion
mpres at mpMinLatF = -21
mpres at mpMaxLonF = 37
mpres at mpMinLonF = 21
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; The specific pressure levels that we want the data interpolated to.
pressure_levels = (/ 850., 700., 500., 300./) ; pressure levels to plot
nlevels = dimsizes(pressure_levels) ; number of pressure
levels
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; VTlab = (/ "T+00", T+00", "T+06", "T+12", "T+18", "T+24", "T+30",
"T+36", "T+42", "T+48"/) ; Validity for file name
do it = 0,ntimes(0)-1,2 ; TIME LOOP
print("Working on time: " + times_to_use(it) )
;res at TimeLabel = times_to_use(it) ; Set Valid time to use on plots
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; First get the variables we will need
tc = wrf_user_getvar(files_all,"tc",it) ; T in C
u = wrf_user_getvar(files_all,"ua",it) ; u averaged to mass
points
v = wrf_user_getvar(files_all,"va",it) ; v averaged to mass
points
p = wrf_user_getvar(files_all, "pressure",it) ; pressure is our
vertical coordinate
z = wrf_user_getvar(files_all, "z",it) ; grid point height
rh = wrf_user_getvar(files_all,"rh",it) ; relative humidity
tc at units = ""
u at units = ""
v at units = ""
p at units = ""
z at units = ""
rh at units = ""
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
do level = 0,nlevels-1 ; LOOP OVER LEVELS
pressure = pressure_levels(level)
tc_plane = wrf_user_intrp3d(tc,p,"h",pressure,0.,False)
z_plane = wrf_user_intrp3d( z,p,"h",pressure,0.,False)
rh_plane = wrf_user_intrp3d(rh,p,"h",pressure,0.,False)
u_plane = wrf_user_intrp3d( u,p,"h",pressure,0.,False)
v_plane = wrf_user_intrp3d( v,p,"h",pressure,0.,False)
spd = (u_plane*u_plane + v_plane*v_plane)^(0.5) ; m/sec
spd at description = "Wind Speed"
;spd at units = "m/s"
spd at units = ""
u_plane = u_plane*1.94386 ; kts
v_plane = v_plane*1.94386 ; kts
;u_plane at units = "kts"
;v_plane at units = "kts"
u_plane at units = ""
v_plane at units = ""
wks = gsn_open_wks(type,"Plots00-" + "P" + pressure_levels(level) +
"_" + times_to_use(it) + "Z")
; Plotting options for T
opts = res
opts at cnLineColor = "Red"
opts at ContourParameters = (/ 5.0 /)
opts at cnInfoLabelOrthogonalPosF = 0.07 ; offset second label
information
opts at gsnContourLineThicknessesScale = 2.0
opts at FieldTitle = "" ; overwrite Field Title
opts at PlotLevelID = ""
opts at cnInfoLabelString = ""
contour_tc = wrf_contour(files_all[it],wks,tc_plane,opts)
delete(opts)
; Plotting options for RH
opts = res
opts at cnFillOn = True
opts at pmLabelBarOrthogonalPosF = -0.1
opts at ContourParameters = (/ 10., 90., 10./)
opts at FieldTitle = "" ; overwrite Field Title
opts at PlotLevelID = ""
opts at cnInfoLabelString = ""
opts at cnFillColors = (/"White","White","White", \
"White","Chartreuse","Green",\
"Green3","Green4", \
"ForestGreen","PaleGreen4"/)
contour_rh = wrf_contour(files_all[it],wks,rh_plane,opts)
delete(opts)
; Plotting options for Wind Speed
opts = res
opts at cnLineColor = "MediumSeaGreen"
opts at ContourParameters = (/ 10. /)
opts at cnInfoLabelOrthogonalPosF = 0.07 ; offset second label
information
opts at FieldTitle = "" ; overwrite Field Title
opts at PlotLevelID = ""
opts at cnInfoLabelString = ""
opts at gsnContourLineThicknessesScale = 3.0
contour_spd = wrf_contour(files_all[it],wks,spd,opts)
delete(opts)
; Plotting options for Wind Vectors
opts = res
;opts at FieldTitle = "Wind" ; overwrite Field Title
opts at FieldTitle = "" ; overwrite Field Title
opts at PlotLevelID = ""
;opts at cnInfoLabelString = ""
opts at NumVectors = 47 ; wind barb density
vector = wrf_vector(files_all[it],wks,u_plane,v_plane,opts)
delete(opts)
; Plotting options for Geopotential Heigh
opts_z = res
opts_z at cnLineColor = "Blue"
opts_z at FieldTitle = "" ; overwrite Field Title
opts_z at PlotLevelID = ""
;opts at cnInfoLabelString = ""
opts_z at gsnContourLineThicknessesScale = 3.0
; MAKE PLOTS
; wks = gsn_open_wks(type,"plt_PressureLevel1"+level+it)
if ( pressure .eq. 850 ) then ; plot temp, rh, height, wind barbs
opts_z at ContourParameters = (/ 20.0 /)
contour_height = wrf_contour(files_all[it],wks,z_plane,opts_z)
plot =
wrf_map_overlays(files_all[it],wks,(/contour_rh,contour_tc,contour_height, \
vector/),pltres,mpres)
end if
if ( pressure .eq. 700 ) then ; plot temp, height, wind barbs
opts_z at ContourParameters = (/ 30.0 /)
contour_height = wrf_contour(files_all[it],wks, z_plane,opts_z)
plot =
wrf_map_overlays(files_all[it],wks,(/contour_tc,contour_height, \
vector/),pltres,mpres)
end if
if ( pressure .eq. 500 ) then ; plot temp, height, wind barbs
opts_z at ContourParameters = (/ 60.0 /)
contour_height = wrf_contour(files_all[it],wks, z_plane,opts_z)
plot =
wrf_map_overlays(files_all[it],wks,(/contour_tc,contour_height, \
vector/),pltres,mpres)
end if
if ( pressure .eq. 300 ) then ; plot windspeed, height, wind barbs
opts_z at ContourParameters = (/ 60.0 /)
contour_height = wrf_contour(files_all[it],wks, z_plane,opts_z)
plot =
wrf_map_overlays(files_all[it],wks,(/contour_spd,contour_height, \
vector/),pltres,mpres)
end if
delete(opts_z)
; delete(wks)
end do ; END OF LEVEL LOOP
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
end do ; END OF TIME LOOP
end
<https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail>
Virus-free.
www.avast.com
<https://www.avast.com/sig-email?utm_medium=email&utm_source=link&utm_campaign=sig-email&utm_content=webmail>
<#DAB4FAD8-2DD7-40BB-A1B8-4E2AA1F9FDF2>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20220407/a721238e/attachment.html>
More information about the ncl-talk
mailing list