; Example script to produce plots for a WRF real-data run, ; with the ARW coordinate dynamics option. ; Plot data on a cross section ; This script will plot data at a set angle through a specified point ; This script adds lon/lat info along X-axis load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl" load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl" load "$NCARG_ROOT/lib/ncarg/nclscripts/contrib/cd_string.ncl" external GET_IJ "./get_ij.so" ; external file required to find center location of cross section undef("set_res1_list") function set_res1_list() begin res = True res@cnFillOn = True res@cnLevelSelectionMode = "ExplicitLevels" ; use explicit levels res@cnLevels = ispan(-9,57,3) ; set the contour levels cmap = read_colormap_file("NCV_rainbow2") cmap := cmap(50:, :) cmap(0,:) = 1. ; Make first color white res@cnFillPalette = cmap ; set color map ; res@pmLabelBarOrthogonalPosF = -0.1 res@lbOrientation = "vertical" ; vertical label bar ; opts_fft@lbLabelBarOn = False res@gsnDraw = False ; don't draw res@gsnFrame = False ; don't advance frame res@cnLinesOn = False ; this gives a warning (once) that it's not a valid setting, but it is. ; I think it happens because I reuse the opts for different variables. res@cnLineLabelsOn = False ; opts_fft@lbLabelBarOn = False res@cnInfoLabelOn = False return(res) end undef("set_res2_list") function set_res2_list() begin res = True res@cnFillOn = True res@cnLinesOn = False res@cnLineLabelsOn = False res@lbOrientation = "vertical" ; vertical label bar ; res@lbLabelBarOn = False res@gsnDraw = False ; don't draw res@gsnFrame = False ; don't advance frame res@cnLevelSelectionMode = "ExplicitLevels" ; use explicit levels cmap = read_colormap_file("NCV_rainbow2") cmap(0,:) = 1. ; Make first color white res@cnFillPalette = cmap ; res@cnFillPalette="NCV_rainbow2" res@cnLevels = (/-2.0,-1.0,-0.75,-0.5,-0.25,0,0.25,0.5,0.75,1.0,1.25,1.50,1.75,2.0,\ 2.5,3.0,3.5,4.0,4.5,5.0,7.5,10.0/) ; res@cnFillColors = (/97,90,80,70,63,55,108,115,122,129,157,161,169,176,183,190,197,204,211, \ res@cnFillColors = (/0,90,80,70,63,55,108,115,122,129,157,161,169,176,183,190,197,204,211, \ 218,225,232,239,247,252/) return(res) end undef("set_res3_list") function set_res3_list() begin res = True res@cnFillOn = True res@cnLinesOn = False res@cnLineLabelsOn = False res@lbOrientation = "vertical" ; vertical label bar ; res@lbLabelBarOn = False res@gsnDraw = False ; don't draw res@gsnFrame = False ; don't advance frame ; res@cnFillPalette = "MPL_gist_ncar" cmap = read_colormap_file("MPL_gist_ncar") cmap(0,:) = 1. ; Make first color white res@cnFillPalette = cmap res@cnLevelSelectionMode = "ExplicitLevels" ; use explicit levels res@cnLevels = (/ 0.70,0.71,0.72,0.73,0.74,0.75,0.76,0.77,0.78,0.79,0.80,0.81,\ 0.82,0.83,0.84,0.85,0.86,0.87,0.88,\ 0.89,0.90,0.91,0.92,0.93,0.94,0.95,0.96,\ 0.97,0.98,0.99,1.0/) res@cnFillColors = (/ 0,17,21,25,29,33,37,41,\ 45,49,53,57,61,65,69,73,75,77,\ 81 ,83 ,85,\ 87,89,91,93,95,126,122,117,115,113,111,110/) res@cnLinesOn = False ;--- These are required for gsn_csm_contour_map plotting ;res = wrf_map_resources(f,res) res@tfDoNDCOverlay = True res@gsnAddCyclic = False return(res) end undef("set_res4_list") function set_res4_list() begin res = True res@cnFillOn = True res@cnLinesOn = False res@cnLineLabelsOn = False res@lbOrientation = "vertical" ; vertical label bar ; res@lbLabelBarOn = False res@gsnDraw = False ; don't draw res@gsnFrame = False ; don't advance frame res@cnLevelSelectionMode = "ExplicitLevels" ; use explicit levels cmap = read_colormap_file("NCV_rainbow2") cmap(0,:) = 1. ; Make first color white res@cnFillPalette = cmap ; res@cnLevels = (/ 0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.2,1.8/) res@cnLevels = (/-5,-2,-1.0,-0.5,-0.25,0,0.01,0.02,0.06,0.1,0.2,0.3,0.4,0.5,0.6,\ 0.7,0.8,0.9,1.0,2.0,3.0,5.0/) res@cnFillColors = (/0,90,80,70,63,55,108,115,122,129,157,161,169,176,183,190,197,204,211, \ 218,225,232,239,247,251,252,253/) return(res) end begin plot_rad1 = new ( 4, graphic ) plot_rad2 = new ( 4, graphic ) plot_rad3 = new ( 4, graphic ) ; ; The WRF ARW input file. ; This needs to have a ".nc" appended, so just do it. a = addfile("wrfout_d02_2013-02-08_15:30:00.nc","r") ; my additional code -- set bounds xlat = wrf_user_getvar(a, "XLAT",0) xlon = wrf_user_getvar(a, "XLONG",0) ter = wrf_user_getvar(a, "HGT",0) dims2d=dimsizes(xlat) ; lats = (/ 39.00, 42.50 /) ; lons = (/ -75.0, -69.87 /) ; lats = (/ 39.00, 42.50 /) ; lons = (/ -76.0, -70.0/) lats = (/ 39.50, 42.00 /) lons = (/ -75.0, -71.0/) loc = wrf_user_ll_to_ij(a, lons, lats, True) ; loc(0,;) is west-east (x) ; loc(1,:) is south-north (y) ; subtract one since we want to use it as an index in NCL x_start = loc(0,0) - 1 x_end = loc(0,1) - 1 y_start = loc(1,0) - 1 y_end = loc(1,1) - 1 print("x_start = " + x_start) print("x_end = " + x_end) print("y_start = " + y_start) print("y_end = " + y_end) ; latitude longitude of Upton, NY xlat_zoom = xlat(y_start:y_end,x_start:x_end) xlon_zoom = xlon(y_start:y_end,x_start:x_end) ter_zoom = ter(y_start:y_end,x_start:x_end) dims_zoom = dimsizes(xlat_zoom) lat_beg = 40.8682 ; lon_beg = 72.8792 ; new point further west lon_beg = 73.25 i_loc = 1 j_loc = 1 GET_IJ::get_ij(xlat_zoom,xlon_zoom,lat_beg,lon_beg,i_loc,j_loc,dims_zoom(0),dims_zoom(1)) i_point = i_loc-1 j_point = j_loc-1 print("i_point = " + i_point) print("j_point = " + j_point) print(xlat(j_point,i_point)) print(xlon(j_point,i_point)) ; We generate plots, but what kind do we prefer? ; type = "x11" type = "png" ; wtype = "pdf" type@wkWidth = 2500 ; produces a better looking PNG type@wkHeight = 2500 ; type = "pdf" ; type = "ps" ; type = "ncgm" wks = gsn_open_wks(type,"plt_radar_temp_cross_sections_90_deg") ; wks = gsn_open_wks(type,"plt_radar_temp_cross_sections_00_deg") ; Set some basic resources ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; FirstTime = True FirstTimeMap = True times = wrf_user_getvar(a,"times",-1) ; get times in the file ntimes = dimsizes(times) ; number of times in the file yr = str_get_cols(times, 0, 3) year = stringtoint(yr) mnth = str_get_cols(times, 5, 6) month = stringtoint(mnth) dy = str_get_cols(times, 8, 9) day = stringtoint(dy) hr = str_get_cols(times, 11, 12) hour = stringtoint(hr) mn = str_get_cols(times, 14, 15) minute = stringtoint(mn) print("hour = " + hour) ; not needed with zoom ; mdims = getfilevardimsizes(a,"P") ; get some dimension sizes for the file ; nd = dimsizes(mdims) opts_rad = set_res1_list() opts_zdr = set_res2_list() opts_crs = set_res3_list() opts_kdp = set_res4_list() ; ;--------------------------------------------------------------- ; do it = 0,ntimes-1,2 ; TIME LOOP do it = 15,23,1 ; do it = 16,22,1 ; do it = 16,16,1 it_val = it new_time_units = "hours since 1800-01-01 00:00" new_time_units = new_time_units new_time_array = cd_inv_calendar( year(it_val), month(it_val), day(it_val), \ hour(it_val), minute(it_val), 0, new_time_units,0) print("new_time_array = " + new_time_array) ; format = "" ;; defaults to "%H%M IST %d %c %Y" format = "%H%M UTC %d %c %Y" ; format = "%d %c %Y" new_time_array = new_time_array time_var=cd_string(new_time_array,format) print("time_var = " + time_var) print("Working on time: " + times(it) ) ; res@TimeLabel = times(it) ; Set Valid time to use on plots tc = wrf_user_getvar(a,"tc",it) ; T in C ; rh = wrf_user_getvar(a,"rh",it) ; relative humidity fft = wrf_user_getvar(a,"fft_zh",it) ; total reflectivity ff1 = wrf_user_getvar(a,"ff1_zh",it) ; liquid ff2 = wrf_user_getvar(a,"ff2_zh",it) ; plates ff3 = wrf_user_getvar(a,"ff3_zh",it) ; columns ff4 = wrf_user_getvar(a,"ff4_zh",it) ; dendrites ff5 = wrf_user_getvar(a,"ff5_zh",it) ; snow ff6 = wrf_user_getvar(a,"ff6_zh",it) ; graupel ff7 = wrf_user_getvar(a,"ff7_zh",it) ; hail z = wrf_user_getvar(a, "z",it) ; grid point height zdr = wrf_user_getvar(a,"fft_zdr",it) crs = wrf_user_getvar(a,"fft_crs",it) kdp = wrf_user_getvar(a,"fft_kdp",it) tc_zoom = tc(:,y_start:y_end,x_start:x_end) fft_zoom = fft(:,y_start:y_end,x_start:x_end) ff1_zoom = ff1(:,y_start:y_end,x_start:x_end) ff2_zoom = ff2(:,y_start:y_end,x_start:x_end) ff3_zoom = ff3(:,y_start:y_end,x_start:x_end) ff4_zoom = ff4(:,y_start:y_end,x_start:x_end) ff5_zoom = ff5(:,y_start:y_end,x_start:x_end) ff6_zoom = ff6(:,y_start:y_end,x_start:x_end) ff7_zoom = ff7(:,y_start:y_end,x_start:x_end) z_zoom = z(:,y_start:y_end,x_start:x_end) zdr_zoom = zdr(:,y_start:y_end,x_start:x_end) crs_zoom = crs(:,y_start:y_end,x_start:x_end) kdp_zoom = kdp(:,y_start:y_end,x_start:x_end) if ( FirstTime ) then ; get height info for labels zmin = 0. zmax = 6. ; We are only interested in the first 6km nz = floattoint(zmax + 1) end if ;--------------------------------------------------------------- ; do ip = 1, 3 ; we are doing 3 plots do ip = 2, 2 ; we are doing 3 plots ; do ip = 1, 1 ; we are doing 3 plots ; all with the pivot point (plane) in the center of the domain ; at angles 0, 45 and 90 ; ; | ; angle=0 is | ; | ; plane = new(2,float) ; note this part of the code assumes that the arrays have been created (x,y) ; plane = (/ mdims(nd-1)/2, mdims(nd-2)/2 /) ; pivot point is center of domain (x,y) plane = (/ i_point, j_point /) ; pivot point is Upton opts = False if(ip .eq. 1) then angle = 90. X_plane = wrf_user_intrp2d(xlon_zoom,plane,angle,opts) X_desc = "longitude" end if if(ip .eq. 2) then angle = 0. X_plane = wrf_user_intrp2d(xlat_zoom,plane,angle,opts) X_desc = "latitude" end if if(ip .eq. 3) then ; angle = 45. angle = 135. X_plane = wrf_user_intrp2d(xlon_zoom,plane,angle,opts) X_desc = "longitude" end if fft_plane = wrf_user_intrp3d(fft_zoom,z_zoom,"v",plane,angle,opts) ff1_plane = wrf_user_intrp3d(ff1_zoom,z_zoom,"v",plane,angle,opts) ff2_plane = wrf_user_intrp3d(ff2_zoom,z_zoom,"v",plane,angle,opts) ff3_plane = wrf_user_intrp3d(ff3_zoom,z_zoom,"v",plane,angle,opts) ff4_plane = wrf_user_intrp3d(ff4_zoom,z_zoom,"v",plane,angle,opts) ff5_plane = wrf_user_intrp3d(ff5_zoom,z_zoom,"v",plane,angle,opts) ff6_plane = wrf_user_intrp3d(ff6_zoom,z_zoom,"v",plane,angle,opts) ff7_plane = wrf_user_intrp3d(ff7_zoom,z_zoom,"v",plane,angle,opts) tc_plane = wrf_user_intrp3d(tc_zoom,z_zoom,"v",plane,angle,opts) zdr_plane = wrf_user_intrp3d(zdr_zoom,z_zoom,"v",plane,angle,opts) crs_plane = wrf_user_intrp3d(crs_zoom,z_zoom,"v",plane,angle,opts) kdp_plane = wrf_user_intrp3d(kdp_zoom,z_zoom,"v",plane,angle,opts) printVarSummary(fft_plane) ; Find the index where 6km is - only need to do this once if ( FirstTime ) then zz = wrf_user_intrp3d(z_zoom,z,"v",plane,angle,opts) b = ind(zz(:,0) .gt. zmax*1000. ) zmax_pos = b(0) - 1 if ( abs(zz(zmax_pos,0)-zmax*1000.) .lt. abs(zz(zmax_pos+1,0)-zmax*1000.) ) then zspan = b(0) - 1 else zspan = b(0) end if delete(zz) delete(b) FirstTime = False end if ; X-axis lables dimsX = dimsizes(X_plane) xmin = X_plane(0) xmax = X_plane(dimsX(0)-1) xspan = dimsX(0)-1 nx = floattoint( (xmax-xmin)/2 + 3) ; ; xmin_new = nx-xmin ; xmax_new = nx+xmin ; print("xmin_new = " + xmin_new) ; print("xmax_new = " + xmax_new) ;--------------------------------------------------------------- ; Options for XY Plots opts_xy = True opts_xy@tiXAxisString = X_desc opts_xy@tiYAxisString = "Height (km)" opts_xy@cnMissingValPerimOn = True opts_xy@cnMissingValFillColor = 0 opts_xy@cnMissingValFillPattern = 11 opts_xy@tmXTOn = False opts_xy@tmYROn = False opts_xy@tmXBMode = "Explicit" opts_xy@tmXBValues = fspan(0,xspan,nx) ; Create tick marks opts_xy@tmXBLabels = sprintf("%.1f",fspan(xmin,xmax,nx)) ; Create labels print("opts_xy@tmXBLabels" + opts_xy@tmXBLabels) opts_xy@tmXBLabelFontHeightF = 0.015 opts_xy@tmYLMode = "Explicit" opts_xy@tmYLValues = fspan(0,zspan,nz) ; Create tick marks opts_xy@tmYLLabels = sprintf("%.1f",fspan(zmin,zmax,nz)) ; Create labels opts_xy@tiXAxisFontHeightF = 0.020 opts_xy@tiYAxisFontHeightF = 0.020 opts_xy@tmXBMajorLengthF = 0.02 opts_xy@tmYLMajorLengthF = 0.02 opts_xy@tmYLLabelFontHeightF = 0.015 ; opts_xy@PlotOrientation = tc_plane@Orientation ; Plotting options for RH opts_fft = opts_xy ; opts_fft = True opts_fft@cnFillOn = True ; opts_fft@ContourParameters = (/ -10., 60., 5. /) opts_fft@cnLevelSelectionMode = "ExplicitLevels" ; use explicit levels opts_fft@cnLevels = ispan(-9,57,3) ; set the contour levels cmap := read_colormap_file("NCV_rainbow2") cmap := cmap(50:, :) cmap(0,:) = 1. ; Make first color white opts_fft@cnFillPalette = cmap ; set color map opts_fft@pmLabelBarOrthogonalPosF = -0.1 ; opts_fft@lbLabelBarOn = False opts_fft@gsnDraw = False ; don't draw opts_fft@gsnFrame = False ; don't advance frame opts_fft@cnLinesOn = False ; this gives a warning (once) that it's not a valid setting, but it is. ; I think it happens because I reuse the opts for different variables. opts_fft@cnLineLabelsOn = False opts_fft@lbLabelBarOn = False opts_fft@cnInfoLabelOn = False opts_ff1 = opts_fft opts_ff2 = opts_fft opts_ff3 = opts_fft opts_ff4 = opts_fft opts_ff5 = opts_fft opts_ff6 = opts_fft opts_ff7 = opts_fft ; this gives an not valid at this time, but it is valid. Strange. opts_fft@gsnLeftString = "Total Radar Reflectivity" opts_fft@gsnRightString = "dBZ" opts_ff1@gsnLeftString = "Radar Reflectivity of Liquid Water" opts_ff2@gsnLeftString = "Radar Reflectivity of Plates" opts_ff3@gsnLeftString = "Radar Reflectivity of Columns" opts_ff4@gsnLeftString = "Radar Reflectivity of Dendrites" opts_ff5@gsnLeftString = "Radar Reflectivity of Snow" opts_ff6@gsnLeftString = "Radar Reflectivity of Graupel" opts_ff7@gsnLeftString = "Radar Reflectivity of Hail" ;;; opts_xy opts_rad@tiXAxisString = X_desc opts_rad@tiYAxisString = "Height (km)" opts_rad@cnMissingValPerimOn = True opts_rad@cnMissingValFillColor = 0 opts_rad@cnMissingValFillPattern = 11 opts_rad@tmXTOn = False opts_rad@tmYROn = False opts_rad@tmXBMode = "Explicit" opts_rad@tmXBValues = fspan(0,xspan,nx) ; Create tick marks opts_rad@tmXBLabels = sprintf("%.1f",fspan(xmin,xmax,nx)) ; Create labels print("opts_rad@tmXBLabels" + opts_rad@tmXBLabels) opts_rad@tmXBLabelFontHeightF = 0.015 opts_rad@tmYLMode = "Explicit" opts_rad@tmYLValues = fspan(0,zspan,nz) ; Create tick marks opts_rad@tmYLLabels = sprintf("%.1f",fspan(zmin,zmax,nz)) ; Create labels opts_rad@tiXAxisFontHeightF = 0.020 opts_rad@tiYAxisFontHeightF = 0.020 opts_rad@tmXBMajorLengthF = 0.02 opts_rad@tmYLMajorLengthF = 0.02 opts_rad@tmYLLabelFontHeightF = 0.015 opts_zdr@tiXAxisString = X_desc opts_zdr@tiYAxisString = "Height (km)" opts_zdr@cnMissingValPerimOn = True opts_zdr@cnMissingValFillColor = 0 opts_zdr@cnMissingValFillPattern = 11 opts_zdr@tmXTOn = False opts_zdr@tmYROn = False opts_zdr@tmXBMode = "Explicit" opts_zdr@tmXBValues = fspan(0,xspan,nx) ; Create tick marks opts_zdr@tmXBLabels = sprintf("%.1f",fspan(xmin,xmax,nx)) ; Create labels print("opts_zdr@tmXBLabels" + opts_zdr@tmXBLabels) opts_zdr@tmXBLabelFontHeightF = 0.015 opts_zdr@tmYLMode = "Explicit" opts_zdr@tmYLValues = fspan(0,zspan,nz) ; Create tick marks opts_zdr@tmYLLabels = sprintf("%.1f",fspan(zmin,zmax,nz)) ; Create labels opts_zdr@tiXAxisFontHeightF = 0.020 opts_zdr@tiYAxisFontHeightF = 0.020 opts_zdr@tmXBMajorLengthF = 0.02 opts_zdr@tmYLMajorLengthF = 0.02 opts_zdr@tmYLLabelFontHeightF = 0.015 ; opts_zdr@PlotOrientation = tc_plane@Orientation opts_crs@tiXAxisString = X_desc opts_crs@tiYAxisString = "Height (km)" opts_crs@cnMissingValPerimOn = True opts_crs@cnMissingValFillColor = 0 opts_crs@cnMissingValFillPattern = 11 opts_crs@tmXTOn = False opts_crs@tmYROn = False opts_crs@tmXBMode = "Explicit" opts_crs@tmXBValues = fspan(0,xspan,nx) ; Create tick marks opts_crs@tmXBLabels = sprintf("%.1f",fspan(xmin,xmax,nx)) ; Create labels print("opts_crs@tmXBLabels" + opts_crs@tmXBLabels) opts_crs@tmXBLabelFontHeightF = 0.015 opts_crs@tmYLMode = "Explicit" opts_crs@tmYLValues = fspan(0,zspan,nz) ; Create tick marks opts_crs@tmYLLabels = sprintf("%.1f",fspan(zmin,zmax,nz)) ; Create labels opts_crs@tiXAxisFontHeightF = 0.020 opts_crs@tiYAxisFontHeightF = 0.020 opts_crs@tmXBMajorLengthF = 0.02 opts_crs@tmYLMajorLengthF = 0.02 opts_crs@tmYLLabelFontHeightF = 0.015 ; opts_crs@PlotOrientation = tc_plane@Orientation opts_kdp@tiXAxisString = X_desc opts_kdp@tiYAxisString = "Height (km)" opts_kdp@cnMissingValPerimOn = True opts_kdp@cnMissingValFillColor = 0 opts_kdp@cnMissingValFillPattern = 11 opts_kdp@tmXTOn = False opts_kdp@tmYROn = False opts_kdp@tmXBMode = "Explicit" opts_kdp@tmXBValues = fspan(0,xspan,nx) ; Create tick marks opts_kdp@tmXBLabels = sprintf("%.1f",fspan(xmin,xmax,nx)) ; Create labels print("opts_kdp@tmXBLabels" + opts_kdp@tmXBLabels) opts_kdp@tmXBLabelFontHeightF = 0.015 opts_kdp@tmYLMode = "Explicit" opts_kdp@tmYLValues = fspan(0,zspan,nz) ; Create tick marks opts_kdp@tmYLLabels = sprintf("%.1f",fspan(zmin,zmax,nz)) ; Create labels opts_kdp@tiXAxisFontHeightF = 0.020 opts_kdp@tiYAxisFontHeightF = 0.020 opts_kdp@tmXBMajorLengthF = 0.02 opts_kdp@tmYLMajorLengthF = 0.02 opts_kdp@tmYLLabelFontHeightF = 0.015 ; opts_kdp@PlotOrientation = tc_plane@Orientation ; Plotting options for Temperature opts_tc = opts_xy opts_tc@cnInfoLabelZone = 1 opts_tc@cnInfoLabelSide = "Top" opts_tc@cnInfoLabelPerimOn = True opts_tc@cnInfoLabelOrthogonalPosF = -0.00005 ; opts_tc@ContourParameters = (/ 5. /) opts_tc@cnLevelSelectionMode = "ExplicitLevels" ; use explicit levels opts_tc@cnLevels = ispan(-16,8,2) ; set the contour levels opts_tc@gsnDraw = False ; don't draw opts_tc@gsnFrame = False ; don't advance frame opts_tc@cnInfoLabelOn = False opts_tc@gsnContourZeroLineThicknessF = 3.5 ; Get the contour info for the rh and temp plot_tc1 = gsn_contour(wks,tc_plane(0:zmax_pos,:),opts_tc) plot_tc2 = gsn_contour(wks,tc_plane(0:zmax_pos,:),opts_tc) plot_tc3 = gsn_contour(wks,tc_plane(0:zmax_pos,:),opts_tc) plot_tc4 = gsn_contour(wks,tc_plane(0:zmax_pos,:),opts_tc) plot_rad1(0) = gsn_csm_contour(wks,fft_plane(0:zmax_pos,:),opts_fft) plot_rad2(0) = gsn_csm_contour(wks,ff1_plane(0:zmax_pos,:),opts_ff1) plot_rad3(0) = gsn_csm_contour(wks,fft_plane(0:zmax_pos,:),opts_rad) plot_rad1(1) = gsn_csm_contour(wks,ff2_plane(0:zmax_pos,:),opts_ff2) plot_rad1(2) = gsn_csm_contour(wks,ff3_plane(0:zmax_pos,:),opts_ff3) plot_rad2(1) = gsn_csm_contour(wks,ff5_plane(0:zmax_pos,:),opts_ff5) plot_rad2(2)= gsn_csm_contour(wks,ff6_plane(0:zmax_pos,:),opts_ff6) plot_rad1(3) = gsn_csm_contour(wks,ff4_plane(0:zmax_pos,:),opts_ff4) plot_rad2(3) = gsn_csm_contour(wks,ff7_plane(0:zmax_pos,:),opts_ff7) plot_rad3(1) = gsn_csm_contour(wks,zdr_plane(0:zmax_pos,:),opts_zdr) plot_rad3(2) = gsn_csm_contour(wks,crs_plane(0:zmax_pos,:),opts_crs) plot_rad3(3) = gsn_csm_contour(wks,kdp_plane(0:zmax_pos,:),opts_kdp) ;--------------------------------------------------------------- ; MAKE PLOTS ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; FirstTimeMap = False ; Panel the WRF plots. Type = "QC + QR" pnlres = True ; pnlres@txString = t@description + " (" + t@units + ")" ; pnlres@txString = Type + " at " + radar_deg + " deg elevation (g/cm^3)" ; pnlres@txString = "Vertical Cross Sections" pnlres@txString = time_var pnlres@gsnPanelYWhiteSpacePercent = 5 ; Add white space b/w plots. pnlres@gsnPanelLabelBar = True ; Turn on common labelbar pnlres@lbLabelAutoStride = True ; Spacing of lbar labels. pnlres@lbBoxMinorExtentF = 0.13 pnlres@lbLabelFontHeightF = .015 pnlres@lbPerimOn = False pnlres@lbLabelFont = "Helvetica" ; label font pnlres@lbTitleOn = True ; turn on title pnlres@lbTitleString = "dBZ" ; title string ; pnlres@lbTitlePosition = "Top" ; title position pnlres@lbTitlePosition = "Bottom" ; title position pnlres@lbTitleFontHeightF= .015 ; make title smaller pnlres@lbTitleDirection = "Across" ; title direction pnlres@pmLabelBarWidthF = 0.5 pnlres@pmLabelBarHeightF = 0.15 pnlres@lbTopMarginF = 0.001 pnlres@pmLabelBarOrthogonalPosF = 0.02 ; position wrt plot pnlres@lbTitleOffsetF = -0.4 overlay(plot_rad1(0),plot_tc1) ; NhlRemoveOverlay(plot_rad1(0),plot_tc,False) overlay(plot_rad1(1),plot_tc2) ; NhlRemoveOverlay(plot_rad1(1),plot_tc,False) overlay(plot_rad1(2),plot_tc3) ; NhlRemoveOverlay(plot_rad1(2),plot_tc,False) overlay(plot_rad1(3),plot_tc4) gsn_panel(wks,(/plot_rad1/),(/2,2/),pnlres) NhlRemoveOverlay(plot_rad1(0),plot_tc1,False) NhlRemoveOverlay(plot_rad1(1),plot_tc2,False) NhlRemoveOverlay(plot_rad1(2),plot_tc3,False) NhlRemoveOverlay(plot_rad1(3),plot_tc4,False) overlay(plot_rad2(0),plot_tc1) overlay(plot_rad2(1),plot_tc2) overlay(plot_rad2(2),plot_tc3) overlay(plot_rad2(3),plot_tc4) gsn_panel(wks,(/plot_rad2/),(/2,2/),pnlres) NhlRemoveOverlay(plot_rad2(0),plot_tc1,False) NhlRemoveOverlay(plot_rad2(1),plot_tc2,False) NhlRemoveOverlay(plot_rad2(2),plot_tc3,False) NhlRemoveOverlay(plot_rad2(3),plot_tc4,False) overlay(plot_rad3(0),plot_tc1) overlay(plot_rad3(1),plot_tc2) overlay(plot_rad3(2),plot_tc3) overlay(plot_rad3(3),plot_tc4) pnlres@gsnPanelLabelBar = False ; Turn on common labelbar gsn_panel(wks,(/plot_rad3/),(/2,2/),pnlres) delete(opts_xy) delete(opts_tc) delete(opts_fft) delete(tc_plane) delete(fft_plane) delete(ff1_plane) delete(ff2_plane) delete(ff3_plane) delete(ff4_plane) delete(ff5_plane) delete(ff6_plane) delete(ff7_plane) delete(X_plane) delete(zdr_plane) delete(crs_plane) delete(kdp_plane) end do ; make next cross section ; Delete options and fields, so we don't have carry over end do ; END OF TIME LOOP end