# [ncl-talk] WRF streamlines in ncl

Dennis Shea shea at ucar.edu
Thu Oct 6 16:06:18 MDT 2022

```1] NCL calculates streamlines ONLY for global rectilinear grids using
spherical harmonics.
*    https://www.ncl.ucar.edu/Applications/stream.shtml*
<https://www.ncl.ucar.edu/Applications/stream.shtml>

NCL does not provide any functions for creating streamlines for WRF
winds, which, I believe, are on a curvilinear staggered grid.

I am not aware of WRF streamlines being plotted in *ANY* language
(Matlab, Python, IDL, GrADS, ferret, ...).

[2] Perhaps, a WRF post-processing group could help.

[3] Generally, speaking, the following is not correct:

lat2d    = f->XLAT(0,:,:)
lon2d    = f->XLONG(0,:,:)

lat      = lat2d(:,0)           ; create classic 1D coordinate arrays
lon      = lon2d(0,:)

[4] **IF IT WAS CORTRECT**, you would have to do

lat at units= "degrees_north"
lon at units= "degrees_east"
lat!0    = "lat"
lon!0    = "lon"
lat&lat  =  lat
lon&lon  =  lon

[5] Hypothetically, if the grid was  {u,v}_plane  was 2-dimensional
u_plane!0   = "lat"
u_plane!1   = "lon"
u_plane&lat =  lat
u_plane&lon =  lon
printVarSummary(uplane)

same for v_plane

[5]

plot = gsn_csm_streamline_map(wks,u_plane,v_plane,res)

On Thu, Oct 6, 2022 at 2:04 PM Zilore Mumba via ncl-talk <
ncl-talk at mailman.ucar.edu> wrote:

>
> I am wondering if there is anybody drawing streamlines with WRF output. My
> search has yielded nothing and my attempt has been disastrous.
> I am having two problems:
> 1. how do I attach lon/lat units to the u, v components
> (check_for_y_lat_coord: Warning: Data either does not contain
> a valid latitude coordinate array or doesn't contain one at all... errors).
> See code pasted below.
> 2. The streamlines are drawn, but in the wrong place over a global map.
> If there is anyone who is able to draw streamlines, I will appreciate to
> learn from them.
>
> my doe so far, below
>
> begin
> ;
> ;;;;;;;;;;;;;;
> ;
>   f =
> 00.nc","r")
>
>    ; What times and how many time steps are in the data set?
>    ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>    times = wrf_user_getvar(f,"times",-1)  ; get all times in the file
>    ntimes = dimsizes(times)         ; number of times in the file
>
>    ; The specific pressure levels that we want the data interpolated to.
>    ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>    pressure_levels = (/ 850., 700., 500., 300./)   ; pressure levels to
> plot
>    nlevels         = dimsizes(pressure_levels)     ; number of pressure
> levels
>
>
>     tc = wrf_user_getvar(f,"tc",-1)        ; T in C
>     u  = wrf_user_getvar(f,"ua",-1)        ; u averaged to mass points
>     v  = wrf_user_getvar(f,"va",-1)        ; v averaged to mass points
>     p  = wrf_user_getvar(f, "pressure",-1) ; pressure is our vertical
> coordinate
>     z  = wrf_user_getvar(f, "z",-1)        ; grid point height
>     rh = wrf_user_getvar(f,"rh",-1)        ; relative humidity
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
>      ;lat2d    = wrf_user_getvar(f,"XLAT",0)
>      ;lon2d    = wrf_user_getvar(f,"XLONG",0)
>
>      lat2d    = f->XLAT(0,:,:)
>      lon2d    = f->XLONG(0,:,:)
>
>      lat      = lat2d(:,0)           ; create classic 1D coordinate arrays
>      lon      = lon2d(0,:)
>
>       ;u&lat at units = "degrees_north"
>       ;u&lon at units = "degrees_east"
>       ;v&lat at units = "degrees_north"
>       ;v&lon at units = "degrees_east"
>
>      ;lat at units= "degrees_north"
>      ;lon at units= "degrees_east"
>      ;lat!0    = "lat"
>      ;lon!0    = "lon"
>      ;lat&lat  =  lat
>      ;lon&lon  =  lon
>
>
>       ;pressure = pressure_levels(level)
>  ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>       pressure = pressure_levels(0)
>
>       tc_plane = wrf_user_intrp3d(tc,p,"h",pressure,0.,False)
>       z_plane  = wrf_user_intrp3d( z,p,"h",pressure,0.,False)
>       rh_plane = wrf_user_intrp3d(rh,p,"h",pressure,0.,False)
>       u_plane  = wrf_user_intrp3d( u,p,"h",pressure,0.,False)
>       v_plane  = wrf_user_intrp3d( v,p,"h",pressure,0.,False)
>
>       spd     = (u_plane*u_plane + v_plane*v_plane)^(0.5) ; m/sec
>       spd at description = "Wind Speed"
>       spd at units = "m/s"
>       u_plane = u_plane*1.94386     ; kts
>       v_plane = v_plane*1.94386     ; kts
>       u_plane at units = "kts"
>       v_plane at units = "kts"
>
>      wks  = gsn_open_wks("x11","gsn_csm_streamline_map")
>
>      res                  = True
>      res at gsnAddCyclic     = False
>      res at gsnMaximize      = True
>      res at gsnDraw          = False         ; don't draw yet
>      res at gsnFrame         = False         ; don't advance yet
>
>      res at tiMainString    = "WRF Streamlines"           ; draw a title
>      res at gsnLeftString   = ""
>      res at gsnRightString  = ""
>      res at gsnCenterString = " "
>
>      res at mpProjection  = "Mercator"
>
>     ;---Zoom in on map and plot again
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>      ;res at cnLineThicknessF = 2            ; make lines thicker, Default: 1
>      ;res at mpMinLatF      = min(lat)   ; zoom in on lat/lon area
>      ;res at mpMaxLatF      = max(lat)
>      ;res at mpMinLonF      = min(lon)
>      ;res at mpMaxLonF      = max(lon)
>
>      res at mpDataBaseVersion     = "MediumRes"       ; better and more map
> outlines
>      res at mpDataSetName         = "Earth..4"
>      res at mpOutlineBoundarySets = "AllBoundaries"
>      res at mpOutlineOn           = True
>      res at mpNationalLineThicknessF    = 2.5          ; turn on country
> boundaries
>
>      res at lbOrientation         = "Vertical"
>      res at tiMainOffsetYF        = -0.03           ; Move the title down
>
>      plot = gsn_csm_streamline_map(wks,u_plane,v_plane,res)
>
>     ;Define shapefiles
>     ;*****************
>      shp_filename =
>
>     ;Set shapefile resources
>     ;***********************
>      shpres                    =  True
>      shpres at gsLineThicknessF   =  1             ; increase line thickness
>      shpres at gsLineColor        =  "Black"       ; line
> colorgsLineThicknessF
>      shpres at gsLineDashPattern  =  1
>
>
>      draw(plot)
>      frame(wks)
> end
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at mailman.ucar.edu
> List instructions, subscriber options, unsubscribe:
> https://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20221006/648989d8/attachment.htm>
```