;Load vital libraries for NCL-WRF Functions. 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/wrf/WRFUserARW.ncl" ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;Begin Program and extract data to plot with. begin ;a=addfile(filename,"r") a=addfile("/tera9/brianjs/model/WRF_3.5/Thompson_MYJ/A_5_24_2007_Thompson_MYJ/wrfout_d01_2007-05-23_12:00:00.nc","r") ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;Specify retrieval data for all times and initiate TIME LOOP. times = wrf_user_getvar(a,"times",-1) ntimes = dimsizes(times) times_in_file = a->Times ; The specific pressure levels that we want the data interpolated to. pressure_levels = (/925,900,875,850/) ; pressure levels to plot - in mb nlevels = dimsizes(pressure_levels) ; number of height levels do it = 0,ntimes-1 ; TIME LOOP do level = 0,nlevels-1 ; LOOP OVER LEVELS title_file = chartostring(times_in_file(it,0:3)) +"_"+ chartostring(times_in_file(it,5:6)) +"_"+ chartostring(times_in_file(it,8:9)) +"_"+ chartostring(times_in_file(it,11:12)) height_file = pressure_levels(level) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;Establish Plots (type of plot, title, color map for contours, indications for manual changes). wks = gsn_open_wks("png","Upper_Air_Convergence_"+title_file+"_"+height_file+"mb") res = True res@MainTitle = "WRF Simulated Winds ~C~at Constant Pressure Field" pltres = True mpres = True mpres@mpGeophysicalLineColor = "Black" mpres@mpNationalLineColor = "Black" mpres@mpUSStateLineColor = "Black" mpres@mpGridAndLimbOn = False mpres@mpGridLineColor = "Black" mpres@mpLimbLineColor = "Black" mpres@mpPerimLineColor = "Black" mpres@mpGeophysicalLineThicknessF = 2.0 mpres@mpGridLineThicknessF = 2.0 mpres@mpLimbLineThicknessF = 2.0 mpres@mpNationalLineThicknessF = 2.0 mpres@mpUSStateLineThicknessF = 2.0 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; print("Working on time: " + times(it) ) res@TimeLabel = times(it) ; Set Valid time to use on plots ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; First get the variables we will need ua = wrf_user_getvar(a,"ua",-1) va = wrf_user_getvar(a,"va",-1) p = wrf_user_getvar(a, "pressure",it) ; pressure is our vertical coordinate z = wrf_user_getvar(a, "z",it) ; grid point height all_lat = wrf_user_getvar(a,"lat",-1) all_lon = wrf_user_getvar(a,"lon",-1) lat = all_lat(0,:,0) lon = all_lon(0,0,:) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;Calculate Needed Variables for wind plotting pressure = pressure_levels(level) ;Interpolate/Calculate wind barbs u_plane = wrf_user_intrp3d(ua(it,:,:,:),p,"h",pressure,0.,False) v_plane = wrf_user_intrp3d(va(it,:,:,:),p,"h",pressure,0.,False) u_plane = u_plane*1.94386 ; kts v_plane = v_plane*1.94386 ; kts u_plane@units = "kts" v_plane@units = "kts" ;Calculate wind magnitude for contours. wnd = u_plane wnd = sqrt((u_plane^2)+(v_plane^2)) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;Calculate Convergence conv = u_plane conv = uv2dv_cfd(u_plane,v_plane,lat,lon,2) conv_scaled = conv/1e-5 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Add windbarbs uv_res = res uv_res@vcRefAnnoOn = False ; turns off the ref vector uv_res@vcRefLengthF = 0.035 ; set length of ref vector uv_res@vcGlyphStyle = "WindBarb" ; turn on wind barbs uv_res@gsnSpreadColorStart = 48;-1 ; don't use added gray uv_res@gsnSpreadColorEnd = -2 ; don't use added gray uv_res@vcGlyphStyle = "WindBarb" ; choose wind barbs uv_res@vcWindBarbLineThicknessF = 3;3;2. ; wind barb thickness ; Plotting options for Wind Vectors opts = uv_res opts@FieldTitle = "Winds" ; overwrite Field Title opts@NumVectors = 13 ; wind barb density vector = wrf_vector(a,wks,u_plane,v_plane,opts) delete(opts) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;draw contours of convergence cnres = True cnres@gsnDraw = False cnres@gsnFrame = False cnres@cnLevelSelectionMode = "ExplicitLevels" cnres@cnLevels = (/-5,-2.5,2.5,5/) cnres@cnLineThicknessF = 2.4 cnres@cnLineLabelInterval = 1. cnres@cnLineLabelsOn = True cnres@cnInfoLabelOn = False cnres@cnLineColor = "blue" cnres@FieldTitle = "Divergence (+) and Convergence (-)" cnres@UnitLabel = "10^-5 s^-1" wrf_smooth_2d (conv_scaled, 230) ;Applies necessary smoothing contour = wrf_contour(a,wks,conv_scaled,cnres) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Overlay contours on top of wind barbs and plot data plot = wrf_map_overlays(a,wks,(/contour,vector/),pltres,mpres) end do ; END OF LEVEL LOOP ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; end do ; END OF TIME LOOP end