Re: NCL not plotting first 200 m in the vertical cross section

Guiting,

The problem in Basit's case, if I remember correctly, is that all the data were missing in that location, so the plot was being drawn correctly.

If you think NCL is not plotting your data correctly, then you need to provide us with your script and any necessary data files. Dennis Shea has already asked for these offline.

Thanks,

--Mary

On Mar 14, 2011, at 2:00 AM, Song Guiting (Dr) wrote:

> Dear Mary,
>
> I also encounter the same problem. My NCL version is 6.00 beta.
>
> Best regards,
> Guiting
>
>
> From: ncl-talk-bounces@ucar.edu [mailto:ncl-talk-bounces@ucar.edu] On Behalf Of Mary Haley
> Sent: Monday, February 14, 2011 11:18 PM
> To: Basit Khan
> Cc: ncl-talk@ucar.edu
> Subject: Re: NCL not plotting first 200 m in the vertical cross section
>
> This might be due to a bug we fixed awhile back. What version of NCL (ncl -V) are you running?
> If you can provide me with the file, then I can try your script and see if the bug is fixed.
>
> --Mary
>
>
>
> On Feb 13, 2011, at 3:15 PM, Basit Khan wrote:
>
>
> Hi,
> I am trying to plot vertical cross section (height x latitude) from a WRF output file. The problem is that NCL not plotting anything in the first 200 m while I have at least 11 vertical levels in the first 200 m asl. I am using terrain following coordinates system and most of the cross section area is either sea or terrain less than 50 m high asl. I will greatly appreciate if somebody could help me fix this problem. The ncl code is copied below and the output of the ncl code is attached.
>
> Regards,
>
> Basit
> Centre for Atmospheric Research
> Department of Geography,
> University of Canterbury, Private Bag 4800,
> Christchurch 8140, New Zealand.
> Tel: +643 364 2987 ext. 3088
> email: basit.khan@canterbury.ac.nz
>
> -------------------------------------------------------------------------------
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
>
> begin
> ;
> a = addfile("/hpc/bluefern/WRF/basit/basit-archive/b-haze/bh30/wrfout_d03_2006-05-30_12:00:00.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_CrossSection4")
>
>
> ; Set some basic resources
> res = True
> res@MainTitle = "REAL-TIME WRF"
> res@Footer = False
>
> pltres = True
>
> ;;;;;;;;;;;;;;;; HIGHRES RANGS/GSHHS COASTLINE PLOT RESOURCES ;;;;;;;;;;;;;;;;;;;;;;;;
> mpres =True
> mpres@gsnMaximize = True
> mpres@mpOutlineOn = True
> mpres@mpFillOn = False
> mpres@mpProjection = "Mercator"
> mpres@mpLimitMode = "LatLon"
> mpres@mpMinLonF = 174.4061
> mpres@mpMaxLonF = 175.2502
> mpres@mpMinLatF = -37.27655
> mpres@mpMaxLatF = -36.60195
> mpres@mpOutlineOn = True
> mpres@mpGridLineColor = -1
> mpres@gsnDraw = False ; don't draw the plots
> mpres@gsnFrame = False ; don't advance the frame
> mpres@tmYROn = False
> mpres@tmXTOn = False
> mpres@mpGeophysicalLineColor= "Red"
> mpres@mpLimbLineColor = "Blue"
> mpres@mpGeophysicalLineThicknessF = 1.5
> mpres@mpDataBaseVersion = "HighRes" ; Use the high-res database
> mpres@pmTickMarkDisplayMode = "Always" ; Turn on map tickmarks.
> ; mpres@tiMainString = "Using the RANGS-GSHHS coastal database"
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> FirstTime = True
> FirstTimeMap = True
> times = wrf_user_list_times(a) ; get times in the file
> ntimes = dimsizes(times) ; number of times in the file
>
> mdims = getfilevardimsizes(a,"P") ; get some dimension sizes for the file; print(mdims)=97x50x150x150
> nd = dimsizes(mdims) ; incase of 'P' the print(nd)=4 cause P has four dims
>
> xlat = wrf_user_getvar(a, "XLAT",0)
> xlon = wrf_user_getvar(a, "XLONG",0)
> ter = wrf_user_getvar(a, "HGT",0)
>
> ;---------------------------------------------------------------
>
> do it = 58,74,1 ; TIME LOOP
>
> 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
> z = wrf_user_getvar(a, "height",it) ; grid point height
>
> if ( FirstTime ) then ; get height info for labels
> zmin = 0
> zmax = 1. ; We are only interested in the first 1000 m
> nz = floattoint(zmax + 10)
> end if
>
>
> ;---------------------------------------------------------------
>
> plane = new(2,float)
> plane = (/ mdims(nd-1)/2, mdims(nd-2)/2 /) ; pivot point is center of domain (x,y)
> opts = False
>
>
> angle = 0.
> X_plane = wrf_user_intrp2d(xlat,plane,angle,opts)
> X_desc = "latitude"
>
> rh_plane = wrf_user_intrp3d(rh,z,"v",plane,angle,opts) ; "v" is for vertical and "h" is for horinzontal interporlation
> tc_plane = wrf_user_intrp3d(tc,z,"v",plane,angle,opts) ; Array of vertical coordninates to interpolate
>
>
> ; Find the index where 6km is - only need to do this once
> if ( FirstTime ) then
>
> zz = wrf_user_intrp3d(z,z,"v",plane,angle,opts) ; printVarSummary(zz) = dimsizesV x Dimsize H = 100 x 149
>
> b = ind(zz(:,0) .gt. zmax*1000. ) ; (zmax *1000 --> conver km into metres
>
> 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 + 10)
> ; nx = (xmax-xmin)/2
> ;---------------------------------------------------------------
> ; print("the xmin value is "+xmin)
> ; print("the xmax value is "+xmax)
>
> ; Options for XY Plots
> opts_xy = res
> 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@tmXBValues = fspan(0,xspan,nx) ; Create tick marks
> opts_xy@tmXBLabels = sprintf("%6.2f",fspan(xmin,xmax,nx)) ; Create labels
> opts_xy@tmXBLabelFontHeightF = 0.015
> opts_xy@tmYLMode = "Explicit"
> opts_xy@tmYLValues = fspan(0,zspan,nz) ; Create tick marks
> opts_xy@tmYLLabels = sprintf("%3.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_rh = opts_xy
> opts_rh@ContourParameters = (/ 10., 90., 10. /)
> opts_rh@pmLabelBarOrthogonalPosF = -0.1
> opts_rh@cnFillOn = True
> opts_rh@cnFillColors = (/"White","White","White", \
> "White","Chartreuse","Green", \
> "Green3","Green4", \
> "ForestGreen","PaleGreen4"/)
>
> ; 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 = (/ 1. /)
>
>
> ; Get the contour info for the rh and temp
> contour_tc = wrf_contour(a,wks,tc_plane(0:zmax_pos,:),opts_tc)
> contour_rh = wrf_contour(a,wks,rh_plane(0:zmax_pos,:),opts_rh)
>
> ;---------------------------------------------------------------
>
> ; MAKE PLOTS
>
> if (FirstTimeMap) then
> lat_plane = wrf_user_intrp2d(xlat,plane,angle,opts)
> lon_plane = wrf_user_intrp2d(xlon,plane,angle,opts)
> mpres = True
> pltres = True
> pltres@FramePlot = False
> optsM = res
> optsM@NoHeaderFooter = True
> optsM@cnFillOn = True
> optsM@lbTitleOn = False
> contour = wrf_contour(a,wks,ter,optsM)
> plot = wrf_map_overlays(a,wks,(/contour/),pltres,mpres)
> lnres = True
> lnres@gsLineThicknessF = 3.0
> lnres@gsLineColor = "Red"
>
> do ii = 0,dimsX(0)-2
> gsn_polyline(wks,plot,(/lon_plane(ii),lon_plane(ii+1)/),(/lat_plane(ii),lat_plane(ii+1)/),lnres)
> end do
>
> frame(wks)
> delete(lon_plane)
> delete(lat_plane)
> pltres@FramePlot = True
> end if
>
> plot = wrf_overlays(a,wks,(/contour_rh,contour_tc/),pltres) ; plot x-section
>
> ; Delete options and fields, so we don't have carry over
> delete(opts_xy)
> delete(opts_tc)
> delete(opts_rh)
> delete(tc_plane)
> delete(rh_plane)
> delete(X_plane)
>
> end do ; make next cross section
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> FirstTimeMap = False
> ; end do ; END OF TIME LOOP
>
> end
>
> This email may be confidential and subject to legal privilege, it may
> not reflect the views of the University of Canterbury, and it is not
> guaranteed to be virus free. If you are not an intended recipient,
> please notify the sender immediately and erase all copies of the message
> and any attachments.
>
> Please refer to http://www.canterbury.ac.nz/emaildisclaimer for more
> information.
>
> <vcs_ncl.eps>_______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
> CONFIDENTIALITY: This email is intended solely for the person(s) named and may be confidential and/or privileged. If you are not the intended recipient, please delete it, notify us and do not copy, use, or disclose its content. Thank you.
>
> Towards A Sustainable Earth: Print Only When Necessary
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Wed Mar 16 08:38:31 2011