\ ; Example script to produce dbz plots for a WRF real-data run, ; with the ARW coordinate dynamics option. ; November 2008 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/wrf/WRFUserARW.ncl" load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRF_contributed.ncl" begin ; ; The WRF ARW input file. ; This needs to have a ".nc" appended, so just do it. a = addfile("/scratch/summit/nilu4970/WRFTEST/WRFV3/test/em_real/wrfout_d03_2015-06-13_12:00:00.nc","r") ; a = addfile("/Users/nicholasluchetti/Build_WRF/WRFV3/test/em_quarter_ss/wrfout_d01_simulation_1.nc","r") ; b = addfile("/Users/nicholasluchetti/Build_WRF/WRFV3/test/em_quarter_ss/wrfout_d01_control.nc","r") ; c = addfile("/Users/nicholasluchetti/Build_WRF/WRFV3/test/em_quarter_ss/wrfout_d01_simulation_4.nc","r") ; c = addfile("/Users/nicholasluchetti/Build_WRF/WRFV3/test/em_quarter_ss/wrfout_d01_simulation_3.nc","r") ; c = addfile("/Users/nicholasluchetti/Build_WRF/WRFV3/test/em_quarter_ss/wrfout_d01_simulation_2.nc","r") ; We generate plots, but what kind do we prefer? ; type = "x11" type = "pdf" ; type = "ps" ; type = "ncgm" ; wks = gsn_open_wks(type,"plt_ts_simulation_2_240.pdf") ; wks = gsn_open_wks(type,"plt_ts_simulation_3_100.pdf") wks = gsn_open_wks(type,"plt_ts_061315_terrain.pdf") ; gsn_define_colormap(wks,"WhViBlGrYeOrReWh") ; Overwrite the standard color map ; Set some basic resources res = True res@MainTitle = "XPIA 06/13/15" llres = True ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;s Time = wrf_times_c(a->Times, 0) print(Time) restick = True restick@ttmFormat = "%H:%M:%S" time_axis_labels(Time,res,restick) ; Set desired LocX, LocY, or equivalent lat lon llres = True llres@returnInt = True locij = wrf_user_ll_to_ij(a, -105.0285,40.0003, llres) locij = locij -1 locX = locij(0) locY = locij(1) ; Point 1 ; locX = 85 ; locY = 120 ; Point 2 ; locX = 90 ; locY = 31 ;Horizontal velocity time-series opts_Vel = res opts_Vel@tiXAxisString = "Time (hr)" opts_Vel@tiYAxisString = "10m Wind Speed (m/s)" opts_Vel@tiMainString = "10m Wind Speed (m/s) at BAO" opts_Vel@xyLineThicknesses = (/3.0,3.0,3.0/) opts_Vel@xyLineColors = (/"red","blue","green"/) opts_Vel@pmLegendDisplayMode = "Always" ; turn on legend opts_Vel@pmLegendWidthF = 0.07 ; Change width and opts_Vel@pmLegendHeightF = 0.07 ; height of legend. opts_Vel@lgLabelFontHeightF = .012 ; change font height opts_Vel@pmLegendOrthogonalPosF = -0.99 ; move up slightly opts_Vel@lgBoxMinorExtentF = 0.6 opts_Vel@xyExplicitLegendLabels = (/"10m Wind Speed"/) opts_Vel@pmLegendParallelPosF = 0.23 opts_Vel@tmXTOn = False opts_Vel@tmYLFormat = "f" u = wrf_user_getvar(a,"U10",-1) ; ua is u averaged to mass points v = wrf_user_getvar(a,"V10",-1) Vel = sqrt(u^2 + v^2) printVarSummary(Vel) Vel_point_a = Vel(72:,locY,locX) ; u_b = wrf_user_getvar(b,"ua",-1) ; ua is u averaged to mass points ; v_b = wrf_user_getvar(b,"va",-1) ; Vel_b = sqrt(u_b^2 + v_b^2) ; Vel_point_b = Vel_b(:,0,locY,locX) ; u_c = wrf_user_getvar(c,"ua",-1) ; ua is u averaged to mass points ; v_c = wrf_user_getvar(c,"va",-1) ; Vel_c = sqrt(u_c^2 + v_c^2) ; Vel_point_c = Vel_c(:,0,locY,locX) ;Vel_final = (/Vel_point_a,Vel_point_b, Vel_point_c/) Vel_final = (/Vel_point_a/) vel_plot =gsn_csm_xy(wks,Time(72:),Vel_final,opts_Vel) ; Vertical Velocity time series plot ; nlvmax = 29 ; Z_h = wrf_user_getvar(a,"height",-1) ; Zter = wrf_user_getvar(f,"ter" ,-1) ; Zcat = Z(:,0:nlvmax-1,:,:) ;do kk = 0, nlvmax-1 ; Zcat(:,kk,:,:) = Zcat(:,kk,:,:) - Zter ;end do ; print(Zcat(0,:,29,29)) height_levels = (/ 825./) nlevels = dimsizes(height_levels) print(nlevels) height_levels_2 = (/ 820./) height_levels_3 = (/ 815./) height_levels_4 = (/ 750./) w = wrf_user_getvar(a,"wa",-1) printVarSummary(w) z = wrf_user_getvar(a, "z",-1) p = wrf_user_getvar(a, "pressure",-1) w_plane = wrf_user_intrp3d( w,p,"h",height_levels,0.,False) w_plane = w_plane printVarSummary(w_plane) opts_w = res opts_w@tiXAxisString = "Time (hr)" opts_w@tiYAxisString = "Vertical Velocity (m/s)" opts_w@tiMainString = "Vertical Velocity (m/s) at BAO" opts_w@xyLineThicknesses = (/ 3.0, 3.0, 3.0, 3.0/) opts_w@pmLegendDisplayMode = "Always" ; turn on legend opts_w@pmLegendWidthF = 0.07 ; Change width and opts_w@pmLegendHeightF = 0.07 ; height of legend. opts_w@lgLabelFontHeightF = .012 ; change font height opts_w@pmLegendOrthogonalPosF = -1.20 ; move up slightly opts_w@lgBoxMinorExtentF = 0.6 opts_w@xyExplicitLegendLabels = (/"40m","100m","300m","1000m"/) opts_w@pmLegendParallelPosF = 0.23 opts_w@xyLineColors = (/"red","blue","green","orange"/) w_point = w_plane(72:,locY,locX) printVarSummary(w_point) w_plane_2 = wrf_user_intrp3d( w,p,"h",height_levels_2,0.,False) w_point_2 = w_plane_2(72:,locY,locX) w_plane_3 = wrf_user_intrp3d( w,p,"h",height_levels_3,0.,False) w_point_3 = w_plane_3(72:,locY,locX) w_plane_4 = wrf_user_intrp3d( w,p,"h",height_levels_4,0.,False) w_point_4 = w_plane_4(72:,locY,locX) ;Plotting multiple simulations ; w_b = wrf_user_getvar(b,"wa",-1) ; z_b = wrf_user_getvar(b, "z",-1) ; p_b = wrf_user_getvar(b, "pressure",-1) ; w_plane_b = wrf_user_intrp3d(w_b,z_b,"h",height_levels,0.,False) ; w_plane_b = w_plane_b ; w_point_b = w_plane_b(:,locY,locX) ;w_c = wrf_user_getvar(c,"wa",-1) ;z_c = wrf_user_getvar(c, "z",-1) ;p_c = wrf_user_getvar(c, "pressure",-1) ;w_plane_c = wrf_user_intrp3d(w_c,z_c,"h",height_levels,0.,False) ;w_plane_c = w_plane_c ;w_point_c = w_plane_c(:,locY,locX) ;w_final = (/w_point,w_point_b,w_point_c/) w_final = (/w_point,w_point_2,w_point_3,w_point_4/) w_plot = gsn_csm_xy(wks,Time(72:),w_final,opts_w) ;Velocity plots at different heights ua = wrf_user_getvar(a,"ua",-1) va = wrf_user_getvar(a,"va",-1) u_plane = wrf_user_intrp3d( ua,p,"h",height_levels,0.,False) v_plane = wrf_user_intrp3d( va,p,"h",height_levels,0.,False) Vel2_1 = sqrt(u_plane^2 + v_plane^2) vel2_point_1 = Vel2_1(72:,locY,locX) opts_Vel2 = res opts_Vel2@tiXAxisString = "Time (hr)" opts_Vel2@tiYAxisString = "Wind Speed (m/s)" opts_Vel2@tiMainString = "Wind Speed (m/s) at BAO" opts_Vel2@xyLineThicknesses = (/3.0,3.0,3.0,3.0/) opts_Vel2@xyLineColors = (/"red","blue","green","orange"/) opts_Vel2@pmLegendDisplayMode = "Always" ; turn on legend opts_Vel2@pmLegendWidthF = 0.07 ; Change width and opts_Vel2@pmLegendHeightF = 0.07 ; height of legend. opts_Vel2@lgLabelFontHeightF = .012 ; change font height opts_Vel2@pmLegendOrthogonalPosF = -0.99 ; move up slightly opts_Vel2@lgBoxMinorExtentF = 0.6 opts_Vel2@xyExplicitLegendLabels = (/"40m","100m","300m","1000m"/) opts_Vel2@pmLegendParallelPosF = 0.23 u_plane_2 = wrf_user_intrp3d(ua,p,"h",height_levels_2,0.,False) v_plane_2 = wrf_user_intrp3d(va,p,"h",height_levels_2,0.,False) Vel2_2 = sqrt(u_plane_2^2 + v_plane_2^2) vel2_point_2 = Vel2_2(72:,locY,locX) u_plane_3 = wrf_user_intrp3d(ua,p,"h",height_levels_3,0.,False) v_plane_3 = wrf_user_intrp3d(va,p,"h",height_levels_3,0.,False) Vel2_3 = sqrt(u_plane_3^2 + v_plane_3^2) vel2_point_3 = Vel2_3(72:,locY,locX) u_plane_4 = wrf_user_intrp3d(ua,p,"h",height_levels_4,0.,False) v_plane_4 = wrf_user_intrp3d(va,p,"h",height_levels_4,0.,False) Vel2_4 = sqrt(u_plane_4^2 + v_plane_4^2) vel2_point_4 = Vel2_4(72:,locY,locX) vel2_final = (/vel2_point_1,vel2_point_2,vel2_point_3,vel2_point_4/) vel2_plot = gsn_csm_xy(wks,Time(72:),vel2_final,opts_Vel2) ;Temperature time series plot th = wrf_user_getvar(a,"th",-1) ; get temperature (C) tc2 = th-273.16 ; T2 in C tf2 = 1.8*tc2+32. ; Turn temperature into Fahrenheit tf2@description = "Temperature" tf2@units = "F" th_plane = wrf_user_intrp3d(tc2,p,"h",height_levels,0.,False) th_point = th_plane(72:,locY,locX) opts_th = res opts_th@tiMainString = "Temperature at BAO" ;res@xyMarkLineMode = "MarkLines" opts_th@tiXAxisString = "Time (hr)" opts_th@tiYAxisString = "Temperature (C)" opts_th@xyLineThicknesses = (/ 3.0, 3.0, 3.0, 3.0/) opts_th@pmLegendDisplayMode = "Always" ; turn on legend opts_th@pmLegendWidthF = 0.07 ; Change width and opts_th@pmLegendHeightF = 0.07 ; height of legend. opts_th@lgLabelFontHeightF = .012 ; change font height opts_th@pmLegendOrthogonalPosF = -0.45 ; move up slightly opts_th@lgBoxMinorExtentF = 0.6 opts_th@xyExplicitLegendLabels = (/"40m","100m","300m","1000m"/) opts_th@pmLegendParallelPosF = 0.23 opts_th@xyLineColors = (/"red","blue","green","orange"/) th_plane_2 = wrf_user_intrp3d(tc2,p,"h",height_levels_2,0.,False) th_point_2 = th_plane_2(72:,locY,locX) th_plane_3 = wrf_user_intrp3d(tc2,p,"h",height_levels_3,0.,False) th_point_3 = th_plane_3(72:,locY,locX) th_plane_4 = wrf_user_intrp3d(tc2,p,"h",height_levels_4,0.,False) th_point_4 = th_plane_4(72:,locY,locX) ;Use below for different simulations ; th_b = wrf_user_getvar(b,"th",-1) ; get temperature (C) ; tc2_b = th_b-273.16 ; T2 in C ; tf2_b = 1.8*tc2_b+32. ; Turn temperature into Fahrenheit ; tf2_b@description = "Temperature" ; tf2_b@units = "F" ; th_plane_b = wrf_user_intrp3d(tf2_b,z_b,"h",height_levels,0.,False) ; th_point_b = th_plane_b(:,locY,locX) ; th_c = wrf_user_getvar(c,"th",-1) ; get temperature (C) ; tc2_c = th_c-273.16 ; T2 in C ; tf2_c = 1.8*tc2_c+32. ; Turn temperature into Fahrenheit ; tf2_c@description = "Temperature" ; tf2_c@units = "F" ;th_plane_c = wrf_user_intrp3d(tf2_c,z_c,"h",height_levels,0.,False) ;th_point_c = th_plane_c(:,locY,locX) ;th_final = (/th_point,th_point_b,th_point_c/) th_final = (/th_point,th_point_2,th_point_3,th_point_4/) th_plot = gsn_csm_xy(wks,Time(72:),th_final,opts_th) ; Relative Humidity time series rh=wrf_user_getvar(a,"rh",-1) rh_plane = wrf_user_intrp3d(rh,p,"h",height_levels,0.,False) rh@description = "Relative Humidity (%)" rh_point = rh_plane(72:,locY,locX) rh_plane_2 = wrf_user_intrp3d(rh,p,"h",height_levels_2,0.,False) rh_point_2 = rh_plane_2(72:,locY,locX) rh_plane_3 = wrf_user_intrp3d(rh,p,"h",height_levels_3,0.,False) rh_point_3 = rh_plane_3(72:,locY,locX) rh_plane_4 = wrf_user_intrp3d(rh,p,"h",height_levels_4,0.,False) rh_point_4 = rh_plane_4(72:,locY,locX) opts_rh = res opts_rh@tiMainString = "Relative Humidity (%) at BAO" ;res@xyMarkLineMode = "MarkLines" opts_rh@tiXAxisString = "Time (hr)" opts_rh@tiYAxisString = "Relative Humidity (%)" opts_rh@xyLineThicknesses = (/ 3.0, 3.0, 3.0, 3.0/) opts_rh@pmLegendDisplayMode = "Always" ; turn on legend opts_rh@pmLegendWidthF = 0.071 ; Change width and opts_rh@pmLegendHeightF = 0.07 ; height of legend. opts_rh@lgLabelFontHeightF = .012 ; change font height opts_rh@pmLegendOrthogonalPosF = -1.05 ; move up slightly opts_rh@lgBoxMinorExtentF = 0.6 opts_rh@xyExplicitLegendLabels = (/"40m","100m","300m","1000m"/) opts_rh@pmLegendParallelPosF = 0.23 opts_rh@xyLineColors = (/"red","blue","green","orange"/) ;for multiple simulations ;rh_b=wrf_user_getvar(b,"rh",-1) ;rh_b@description = "Relative Humidity" ;rh_point_b = rh_b(:,0,locY,locX) ;rh_c=wrf_user_getvar(c,"rh",-1) ;rh_c@description = "Relative Humidity" ;rh_point_c = rh_c(:,0,locY,locX) ;rh_final = (/rh_point,rh_point_b,rh_point_c/) rh_final = (/rh_point,rh_point_2,rh_point_3,rh_point_4/) rh_plot = gsn_csm_xy(wks,Time(72:),rh_final,opts_rh) ;Wind Direction Time Series ua = wrf_user_getvar(a,"ua",-1) va = wrf_user_getvar(a,"va",-1) u_plane = wrf_user_intrp3d( ua,p,"h",height_levels,0.,False) v_plane = wrf_user_intrp3d( va,p,"h",height_levels,0.,False) wdir = wind_direction(u_plane,v_plane,0) opts_dir = res opts_dir@tiMainString = "Wind Direction at BAO" ;res@xyMarkLineMode = "MarkLines" opts_dir@tiXAxisString = "Time (hr)" opts_dir@tiYAxisString = "Wind Direction (Degrees)" opts_dir@xyLineThicknesses = (/ 3.0,3.0,3.0,3.0/) opts_dir@pmLegendDisplayMode = "Always" ; turn on legend opts_dir@pmLegendWidthF = 0.07 ; Change width and opts_dir@pmLegendHeightF = 0.07 ; height of legend. opts_dir@lgLabelFontHeightF = .012 ; change font height opts_dir@pmLegendOrthogonalPosF = -0.90 ; move up slightly opts_dir@lgBoxMinorExtentF = 0.6 opts_dir@xyExplicitLegendLabels = (/"40m","100m","300m","1000m"/) opts_dir@pmLegendParallelPosF = 0.23 opts_dir@xyLineColors = (/"red","blue","green","orange"/) dir_point = wdir(72:,locY,locX) u_plane_2 = wrf_user_intrp3d( ua,p,"h",height_levels_2,0.,False) v_plane_2 = wrf_user_intrp3d( va,p,"h",height_levels_2,0.,False) wdir2 = wind_direction(u_plane_2,v_plane_2,0) dir_point_2 = wdir2(72:,locY,locX) u_plane_3 = wrf_user_intrp3d( ua,p,"h",height_levels_3,0.,False) v_plane_3 = wrf_user_intrp3d( va,p,"h",height_levels_3,0.,False) wdir3 = wind_direction(u_plane_3,v_plane_3,0) dir_point_3 = wdir3(72:,locY,locX) u_plane_4 = wrf_user_intrp3d( ua,p,"h",height_levels_4,0.,False) v_plane_4 = wrf_user_intrp3d( va,p,"h",height_levels_4,0.,False) wdir4 = wind_direction(u_plane_4,v_plane_4,0) dir_point_4 = wdir4(72:,locY,locX) ;Plotting for multiple simulations ;wdir_b = wind_direction(u_b,v_b,0) ;dir_point_b = wdir_b(:,0,locY,locX) ;wdir_c = wind_direction(u_c,v_c,0) ;dir_point_c = wdir_c(:,0,locY,locX) ;dir_final = (/dir_point,dir_point_b,dir_point_c/) dir_final = (/dir_point,dir_point_2,dir_point_3,dir_point_4/) dir_plot = gsn_csm_xy(wks,Time(72:),dir_final,opts_dir) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; end