[ncl-talk] WRF streamlines in ncl
Zilore Mumba
zmumba at gmail.com
Fri Oct 7 02:01:27 MDT 2022
Rick,
Thank you once again. wrf_gsn_8.ncl does exactly what I was trying to do. I
do not know how I missed it. I spent the whole of Thursday 6th Oct. trying
to make my code work.
best regards.
On Fri, Oct 7, 2022 at 7:29 AM Zilore Mumba <zmumba at gmail.com> wrote:
> Thanks Rick, I will try that example.
>
> On Thu, Oct 6, 2022 at 10:37 PM Rick Brownrigg <brownrig at ucar.edu> wrote:
>
>> Hi Zilore,
>>
>> There's a WRF examples page that's a bit obscurely referenced from the
>> main examples pages:
>>
>> http://ncl.ucar.edu/Applications/wrfgsn.shtml#ex8
>>
>> Example #8 in particular shows how to set up the coordinates. I am not an
>> expert at this, but hopefully there's enough there to get you going.
>>
>> Rick
>>
>>
>> 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
>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
>>>
>>> begin
>>> ;
>>> ; Read in data.
>>> ;;;;;;;;;;;;;;
>>> ;
>>> f =
>>> addfile("/home/zmumba/DA/OUTPUT/2022092800/noda/wrfout_d01_2022-09-28_00:00:
>>> 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 =
>>> ("$HOME/DA/SCRIPTS/06_Utility_Files/Shapefiles/Lesotho/lso_admbnda_adm1_FAO_MLGCA_2019.shp")
>>>
>>> ;Set shapefile resources
>>> ;***********************
>>> shpres = True
>>> shpres at gsLineThicknessF = 1 ; increase line
>>> thickness
>>> shpres at gsLineColor = "Black" ; line
>>> colorgsLineThicknessF
>>> shpres at gsLineDashPattern = 1
>>>
>>> shp= gsn_add_shapefile_polylines(wks,plot,shp_filename,shpres)
>>>
>>> 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/20221007/6ed038a4/attachment.htm>
More information about the ncl-talk
mailing list