Hi NCL Community,
I have a working cross section of PM data, wind, and temperature. However,
I'd like to cut off the graph at 3km in height and below. When I do this by
putting in a value at opts@trYMaxF, the Yaxis labels get messed up, and I
think it's inaccurately labeling the y-axis because the terrain height and
other contours don't match as the once did. I'd also like the make the
tickmarks more prevalent on the y-axis. Attached is my code and the plots,
one with a cut off yaxis at opts@trYMaxF=24 (? Was experimenting with
values) and other with the entire WRF column. Any help is appreciated,
thanks!
Michael Pirhalla
Graduate Research Assistant
University of Alaska Fairbanks
;--------------------------------------------------------------------------
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/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/contrib/ut_string.ncl"
;load "color_map.ncl"
;---------------------------------------------------------------------------
begin
;What type of plot do we want?
type = "x11"
;type = "ps"
;type = "pdf"
wks = gsn_open_wks(type,"WRF_PM_Cross_Section")
gsn_define_colormap(wks,"BlAqGrYeOrReVi200")
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
; READ WRF OUTPUT FILE AND INCLUDE ALL CASES
; all_files = systemfunc ("ls
/center/d/molders/WITH/wrfout_d01_200*")+".nc"
; print(all_files)
; a = addfiles(all_files,"r")
; print(a)
; ListSetType (a, "cat") ; Combine all files
; idim=dimsizes(all_files)
a =
addfile("/center/d/molders/WITH/wrfout_d01_2008-07-19_00:00:00"+".nc","r")
;---------------------------------------------------------------------------
FirstTime = True
times = wrf_user_getvar(a,"times",-1) ; get times in the file
ntimes = dimsizes(times) ; number of times in the file
print("ntimes"+ntimes)
;---------------------------------------------------------------------------
-
; Set some basic resources
res = True
res@MainTitle = "REAL-TIME WRF" ; main title of plot
res@Footer = True ; footer information
pltres = True
;---------------------------------------------------------------------------
-- ; HOURLY LOOP OF PARTICULATE MATTER do it = 0,ntimes-1 ; loop over all hours ;do it=11,11 ; choose specific hour print("Working on time: " + times(it) ) res@TimeLabel = times(it) ; Set Valid time to use on plots ; First get the variables we will need tc = wrf_user_getvar(a,"tc",it) ; T in C z = wrf_user_getvar(a,"z",it) ; grid point height hgt = wrf_user_getvar(a,"HGT",it) if(FirstTime) then ; get height for labels zmin = 0. zmax = max(z) / 1000 ; make to km nz = floattoint(zmax/2 + 1) FirstTime = False end if ;z_values=fspan(zmin,zmax,100) ;print(z_values) ;exit ; Get PM10 data from WRF/Chem pm = a->PM10(it,:,:,:) pm@description = "PM~B~10~N~ Concentration" pm@units = "~F33~m~F21~g/m~S~3~N~" ; Get the wind data uu=wrf_user_getvar(a,"ua",it) ; u component of wind on mass points vv=wrf_user_getvar(a,"va",it) ; v component of wind on mass points ;uu = a->U10(it,:,:) ; u wind component at 10m ;vv = a->V10(it,:,:) ; v wind component at 10m ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;maxlat=59.1125 ;maxlon=-137.228 ;minlat=58.3728 ;minlon=-135.739 ;lats =(/maxlat,minlat/) ;lons = (/maxlon,minlon/) ; loc = wrf_user_latlon_to_ij(a,lats,lons) ; ys = loc(0,0) ; Starting y (j) location ; ye = loc(1,0) ; Ending y (j) location ; xs = loc(0,1) ; Starting x (i) location ; xe = loc(1,1) ; Ending x (i) location do ip = 3,3 ; we are doing 3 plots, specifying different start and end points opts = True ; setting start and end times/points plane = new(4,float) if(ip .eq. 1) then plane = (/51,42, 52,62/) ; start x;y & end x;y point (X-sec top left diagonal) ;The values should be ordered as (xstart,xend,ystart,yend). ; plane = (/xs,ys, xe,ye/) end if if(ip .eq. 2) then plane = (/48,49, 50,66/) ; start x;y & end x;y point (X-sec top right diagonal) end if if(ip .eq. 3) then plane = (/43,55, 57,57/) end if ; Interpolate WRF Model data vertically angle = 0. ; Establish starting and end points in domain ; plane = new(4,float) ; plane = (/ xs,ys, xe,ye /) ; start x;y & end x;y point pm_plane = wrf_user_intrp3d(pm,z,"v",plane,angle,opts) ; "v" vertical cross section tc_plane = wrf_user_intrp3d(tc,z,"v",plane,angle,opts) u_plane = wrf_user_intrp3d(uu,z,"v",plane,angle,opts) v_plane = wrf_user_intrp3d(vv,z,"v",plane,angle,opts) u_plane=u_plane*1.94386 ; convert to kts v_plane=v_plane*1.94386 u_plane@units = "kts" v_plane@units = "kts" xcoord = ispan(0,dimsizes(pm_plane(0,:)), 1) ;takes the size of the plane (in this case 12 grid points) dim = dimsizes(pm_plane) ; find the data span for use in labels zspan = dim(0) ; Plotting Options for PM opts_xy = res opts_xy@tiYAxisString = "Altitude (km)" opts_xy@tiXAxisString = "Cross Section Distance (km)" opts_xy@tiXAxisFontHeightF = 0.015 opts_xy@tiYAxisOn = True opts_xy@tmXBMode = "Explicit" ;turn manual X axis labels on opts_xy@tmXBValues = xcoord ;labels the start to end grid points as stated above opts_xy@tmXBLabels = xcoord*7 ;take values and mult by 7 since grid cells are 7km width opts_xy@tmXBLabelFontHeightF = 0.01 opts_xy@tmXTOn = False ;turns off top plot tickmarks opts_xy@tmYROn = False ;turns off right tickmarks opts_xy@gsnMaximize = True ;Plot maximized in workstation ; yxis options for Height/PM opts_xy@tmYLMode = "Explicit" opts_xy@tmYLValues = fspan(0,zspan,nz) opts_xy@tmYLLabels = sprintf("%.1f",fspan(zmin,zmax,nz)) ; opts_xy@trYMaxF = 24 ; not sure how this works? ; Specifies a minimum Y coordinate value that defines the lower bound of the Y-Axis coordinate system. ; PM10 Plotting Commands opts_xy@cnFillOn = True ;Fill in the PM Contours opts_xy@cnFillMode = "AreaFill" ; Smooth them opts_xy@cnLinesOn = False ; turn off PM contour lines opts_xy@cnLineLabelsOn = False ; turn off PM contour labels ;options for size of plot dimensions opts_xy@vpXF = 0.105 ; left edge of the view objects bounding box opts_xy@vpWidthF=0.75 ; width of the viewing box opts_xy@vpYF = 0.84355 ; top edge location opts_xy@vpHeightF = 0.51 ; height of the viewing box ;--Temperature lotting Options opts_tc = res opts_tc@cnInfoLabelOn = False ; turn off contour label opts_tc@ContourParameters = (/ 2. /) ; label temp every 5 degrees opts_tc@cnLineLabelInterval = 1 ; every contour labeled opts_tc@cnLineLabelDensityF = .5 ; only one label per contour line ; opts_tc@gsnMaximize = True ; opts_tc@tiXAxisOn = False ;options for size of plot dimensions NECESSARY opts_tc@vpXF = 0.105 opts_tc@vpWidthF=.75 opts_tc@vpYF = 0.84355 opts_tc@vpHeightF = 0.51 ;Same Y Plotting options to restrict data to certain level opts_tc@tmYLMode = "Explicit" opts_tc@tmYLValues = fspan(0,zspan,nz) opts_tc@tmYLLabels = sprintf("%.1f",fspan(zmin,zmax,nz)) ; opts_tc@trYMaxF = 24 ; cut off plot at this level? ;---Wind Plotting Options opts_w=res opts_w@FieldTitle = "Wind Speed" ; Overwrite Field Title opts_w@vcMinDistanceF = 0.05 ; thin out windbarbs opts_w@NumVectors = 20 ; Wind barb density opts_w@gsnMaximize = True ;options for size of plot dimensions NECESSARY opts_w@vpXF = 0.105 opts_w@vpWidthF=.75 opts_w@vpYF = 0.84355 opts_w@vpHeightF = 0.51 ;Same Y Plotting options to restrict data to certain level opts_w@tmYLMode = "Explicit" opts_w@tmYLValues = fspan(0,zspan,35) opts_w@tmYLLabels = sprintf("%.1f",fspan(zmin,zmax,35)) ; opts_w@trYMaxF = 24 ; cut off plot at this level? ; Contour information for PM, Temp, Wind contour_pm = wrf_contour(a,wks,pm_plane,opts_xy) contour_tc = wrf_contour(a,wks,tc_plane,opts_tc) vector=wrf_vector(a,wks,u_plane,v_plane,opts_w) ; RESOURCES FOR FINAL OVERLAY PLOT pltres@gsnMaximize = True pltres@vpXF = 0.105 pltres@vpWidthF=0.75 pltres@vpYF = 0.84355 pltres@vpHeightF = 0.51 ; Plot the graph plot = wrf_overlays(a,wks,(/contour_pm,contour_tc,vector/),pltres) ; delete option fields to avoid carry over delete(opts_tc) delete(opts_xy) delete(tc_plane) delete(pm_plane) end do ; make the next cross section ;----------------------------------------- end do ; end of time loop ;---------------------------------------- end
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
This archive was generated by hypermail 2.1.8 : Wed Dec 04 2013 - 20:42:38 MST