[ncl-talk] Trying to plot streamlines using WRF data, streamlines are not in the correct part of the plot

Craig Tierney - NOAA Affiliate craig.tierney at noaa.gov
Thu Apr 9 15:17:20 MDT 2015


Hello,

After some good help from some list members, I believe that I am setting up
the resources for plotting WRF data without the WRF functions.   However, I
still get no streamlines in my plot, or at least not in the correct place.

I am using NCL 6.1.2.  My WRF domain is close to CONUS.  The domain is
LambertConformal and I am setting the mapping properties stored in the
Global section of the NetCDF file to the appropriate variables for the LC
mapping.   To draw the streamlines, I am basing my code off of the example
here:

https://www.ncl.ucar.edu/Applications/Scripts/gsn_stream_4.ncl

I read in the data, massage it, then try and plot it.  When I plot the
speed, u1, or v1 arrays with contour:

        contour = gsn_csm_contour_map(wks,speed,res)

 I get what I would want, a contour of that variable with the bounds of the
plot that of my domain.  When I try the same thing to create a streamline
plot:

       plot = gsn_streamline_scalar_map(wks,u1,v1,speed,res)

I get nothing.  When I remove setting the Corners of the domain by
commenting out the following:

 res at mpRightCornerLatF                 = lat2d(0,0)
 res at mpRightCornerLonF                = lon2d(0,0)
 res at mpLeftCornerLatF                   = lat2d(nlat-1,nlon-1)
 res at mpLeftCornerLonF                  = lon2d(nlat-1,nlon-1)
 res at mpLimitMode                           = "Corners"

I get streamlines, but they are between 0E and 225E and 0N to 90N.  For the
CONUS domain I am using the corners are about (138W,21N) and (298W,47N).
Something is being calculated, but not what I was expecting.

Any ideas are appreciated.

Thanks,
Craig

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"

begin
;---Open WRF output file.
  dir      =
"/scratch2/portfolios/BMC/winds/Craig.Tierney/theia/wrf/windswrf/trunk/conusplus-wrf-6day/201308010000/nowfp-nocu/wrf/"
  filename = "wrfout_d01_2013-08-04_00:00:00"
  a = addfile(dir + filename + ".nc","r")

;---Read terrain height and lat/lon off file.
  it        = 0     ; first time step
  du1       = wrf_user_getvar(a,"U10",it)    ; Terrain elevation
  dv1       = wrf_user_getvar(a,"V10",it)    ; Terrain elevation
  dlon2d=wrf_user_getvar(a,"XLONG",it)
  dlat2d=wrf_user_getvar(a,"XLAT",it)

; Down sample hack
 d=dimsizes(du1)
 nx=d(0)
 ny=d(1)
 s=8
 u1=du1(0:nx-1:s,0:ny-1:s)
 v1=dv1(0:nx-1:s,0:ny-1:s)
 lat2d=dlat2d(0:nx-1:s,0:ny-1:s)
 lon2d=dlon2d(0:nx-1:s,0:ny-1:s)

 u1 at lon2d=lon2d
 u1 at lat2d=lat2d
 v1 at lon2d=lon2d
 v1 at lat2d=lat2d

  speed=sqrt(u1*u1+v1*v1)
  speed at lat2d=lat2d
  speed at lon2d=lon2d

  wks = gsn_open_wks("png","wrf_gsn")

;---Set some basic plot options
  res               = True

  res at gsnMaximize   = True   ; maximize plot in frame

  res at tiMainString  = filename

  ;res at mpProjection  = "CylindricalEquidistant"
  res at mpProjection  = "LambertConformal"

      res at mpLambertParallel1F = a at TRUELAT1         ; two parallels
      res at mpLambertParallel2F = a at TRUELAT2
      res at mpLambertMeridianF  = a at STAND_LON        ; central meridian
      res at mpCenterLonF        = a at CEN_LON
      res at mpCenterLatF       = a at CEN_LAT

    d=dimsizes(u1)
    nlat=d(0)
    nlon=d(1)

;;;; Comment this section out, and we get streamlines, but in the wrong
place
;;;; This is why they don't show up in the zoomed in plot
;    res at mpRightCornerLatF                 = lat2d(0,0)
;    res at mpRightCornerLonF                = lon2d(0,0)
;    res at mpLeftCornerLatF                   = lat2d(nlat-1,nlon-1)
;    res at mpLeftCornerLonF                  = lon2d(nlat-1,nlon-1)
;    res at mpLimitMode                           = "Corners"

;---Additional resources desired
  res at pmTickMarkDisplayMode = "Always"   ; nicer tickmarks

  res at mpDataBaseVersion     = "MediumRes"       ; better and more map
outlines
  res at mpDataSetName         = "Earth..4"
  res at mpOutlineBoundarySets = "AllBoundaries"
  res at mpOutlineOn           = True

  res at lbOrientation         = "Vertical"
  res at tiMainOffsetYF        = -0.03           ; Move the title down

;; This is where we decide what to plot
  ptype="contour"
;  ptype="streamline"

  if (ptype .eq. "contour") then
print("Plotting: "+ptype)
res at cnFillOn      = True
res at cnFillPalette = "OceanLakeLandSnow"
res at cnLinesOn     = False
;contour = gsn_csm_contour_map(wks,speed,res)
;contour = gsn_csm_contour_map(wks,u1,res)
contour = gsn_csm_contour_map(wks,v1,res)
  else
if (ptype .eq. "streamline") then
print("Plotting: "+ptype)
res at pmTickMarkDisplayMode    = "Always"
  res at pmLabelBarDisplayMode    = "Always"
res at pmLabelBarOrthogonalPosF = -0.02
res at pmLabelBarWidthF         = 0.1
res at lbPerimOn                = False

plot = gsn_streamline_scalar_map(wks,u1,v1,speed,res)
else
print("Unable to plot either contour or streamline: "+ptype)
        exit
end if
  end if

end
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150409/e765950a/attachment.html 


More information about the ncl-talk mailing list