Re: WRF data fro streamlines

From: Adam Phillips <asphilli_at_nyahnyahspammersnyahnyah>
Date: Tue, 23 Oct 2007 15:28:57 -0600

Erik,

Back 4 emails ago you printed this output from printVarSummary(u_plane)
at my request:

Variable: u_plane
Type: float
Total Size: 60604 bytes
             15151 values
Number of Dimensions: 2
Dimensions and sizes: [109] x [139]
Coordinates:
Number Of Attributes: 6
   lon2D : <ARRAY of 15260 elements>
   lat2D : <ARRAY of 15260 elements>
   description : u Velocity
   units : m/s
   _FillValue : -999999
   PlotLevelID : 925 hPa

Notice that the attached lon2D/lat2D arrays have 15260 elements, but
your u_plane grid has 15151 elements (109x139). So you have a mismatch
between your 2D coordinates and the arrays that you are plotting. I took
that to be the reason why you were getting strange plotting results.
Note that the difference betyween 15260 and 15151 is 109, the # of
latitude points in your u_plane array. Thus, I suspected that your
u_plane and v_plane were on staggered grids.

You still have this issue..

Cut your script down to this:
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
a=addfile("/Volumes/Data_and_Models/Model-Output/Athena/WRF-SOP3-Athena_Phys-3-1-1-1.nc","r")
u = wrf_user_getvar(a,"ua",it) ; u averaged to mass points
v = wrf_user_getvar(a,"va",it) ; v averaged to mass points
p = wrf_user_getvar(a, "pressure",it) ; pressure
pressure = 700.
u_plane = wrf_user_intrp3d( u,p,"h",pressure,0.,False)
v_plane = wrf_user_intrp3d( v,p,"h",pressure,0.,False)
u_plane_at_lat2D = a->XLAT_U(0,:,:) ; do you need XLAT_U, or XLAT?
u_plane_at_lon2D = a->XLONG_U(0,:,:) ; do you need XLONG_U, or XLONG?
v_plane_at_lat2D = a->XLAT_V(0,:,:) ; do you need XLAT_V, or XLAT?
v_plane_at_lon2D = a->XLONG_V ; do you need XLONG_V, or XLONG?

printVarSummary(u)
printVarSummary(u_plane)
printVarSummary(a->XLAT_U)
printVarSummary(a->XLAT)

wks = gsn_open_wks("ps","test")
mpres = True ; Create map background
mpres_at_gsnMaximize = True
mpres_at_tiMainString = "Streamlines" ; some tit
mpres_at_mpGeophysicalLineThicknessF = 3.0
mpres_at_gsnAddCyclic = False
plot = gsn_csm_streamline_map(wks,u_plane,v_plane,mpres)
mpres_at_mpMinLatF = -30.
mpres_at_mpMaxLatF = 20.
mpres_at_mpMinLonF = -40.
mpres_at_mpMaxLonF = 40.
plot2 = gsn_csm_streamline_map(wks,u_plane,v_plane,mpres)
end

(I am not sure why you are specifying XLAT_U/XLONG_U when your winds are
not on staggered grids, you may need to replace with XLAT/XLONG.)

Send back to ncl-talk all outputted information, including _all_ error
messages.. Adam

Erik Noble wrote:
> You are correct that WRF, by default, outputs the u/v variables on
> staggered grids.
>
> From the printVarSummary's, you can see that u and v are the same.
> You give the suggestion that I should not have to use the example from
> the 7th example at http://www.ncl.ucar.edu/Applications/wrflc.shtml.
> You suggested that I use this method (below).
>
> So what should use to plot streamlines? I am at square one again.
>
> Two days ago I posted to NCL-talk the same observation that I am
> including below (with my code) noting that NCL incorrectly plots
> regional data directly from WRF using the regular map functions.
> Although the WRF team provided contributed WRF-NCL functions, you have
> trouble plotting anything outside of what is given. (streamlines for
> example).
> Am I the only person that is having this trouble?
>
> I want to plot streamlines. The streamlines are spread all over the
> globe, when they are only supposed to represent a region (here, West
> Africa).
>
> If I use the "gsn_csm_streamline" function (streamline plot draw without
> a map) because the U and V data is strictly just that, "regional
> data".
> When I try to use the gsn_csm_streamline_map function along with a
> specific a subregion, NCL incorrectly assumes the data to represent
> the globe and zooms in on the regional data.
>
> Pictures and previous notes are included.
>
> Any suggestions? Thank you for the help.
> Sincerely,
> Erik Noble
>
>
> ---------- Forwarded message ----------
> From: Adam Phillips <asphilli_at_cgd.ucar.edu <mailto:asphilli_at_cgd.ucar.edu>>
> Date: Oct 22, 2007 4:27 PM
> Subject: Re: Help with regional streamline data drawn to a
> corresponding map
> To: Erik Noble <enoble_at_giss.nasa.gov <mailto:enoble_at_giss.nasa.gov>>
>
>
> Hi Erik,
>
> The u and v winds that come out of the WRF model are not on the same
> grid. You can see this yourself by doing a printVarSummary(u_plane)
> followed by a printVarSummary(v_plane).. Thus, the problems that you're
> having with NCL not plotting anything for you. The solution to this is
> to apply the method outlined in example 7 of the WRF lambert conformal page:
>
> http://www.ncl.ucar.edu/Applications/wrflc.shtml
>
> u = 0.5*(U(nt,kl,:,0:mlonU-2)+U(nt,kl,:,1:mlonU-1))
> v = 0.5*(V(nt,kl,0:nlatV-2,:)+V(nt,kl,1:nlatV-1,:))
>
> u and v would then be on the same grid, and you should be good to go...
> Adam
>
>
> From: Adam Phillips <asphilli_at_cgd.ucar.edu <mailto:asphilli_at_cgd.ucar.edu>>
> Date: Oct 19, 2007 12:15 PM
> Subject: Re: Help with regional streamline data drawn to a
> corresponding map
> To: Erik Noble <enoble_at_giss.nasa.gov <mailto:enoble_at_giss.nasa.gov>>
> Cc: ncl-talk_at_ucar.edu <mailto:ncl-talk_at_ucar.edu>
>
>
> Hi Erik,
>
> If you want a map underneath your streamlines, then
> gsn_csm_streamline_map should be the function that you use. I believe
> the WRF model uses 2D latitudes/longitudes, and thus in NCL you need to set:
> u_plane_at_lat2D = a->XLAT(0,:,:)
> u_plane_at_lon2D = a->XLONG(0,:,:)
> and the same for v_plane...
> Good luck,
> Adam
>
>> *From: *Erik Noble <enoble_at_giss.nasa.gov <mailto:enoble_at_giss.nasa.gov>>
>> *Date: *October 19, 2007 8:36:57 PM EDT
>> *To: *Adam Phillips <asphilli_at_cgd.ucar.edu <mailto:asphilli_at_cgd.ucar.edu>>
>> *Subject: **Re: [ncl-talk] Help with regional streamline data drawn to
>> a* *corresponding ** **map*
>
> Hi Adam. Thank you for this help again. I seem to always have trouble
> when it comes to using map functions that are not in the NCL WRFWRW
> package with WRF data.
>
> Following your suggestions, I get the following map result (pic13) ,
> which is incorrect. The streamlines are spread all over the globe, when
> they are only supposed to represent a region (here, West Africa).
> If i try to provide a subregion as well using something like:
> lonW = -35
> lonE = 35
> latN = 35
> latS = -25
> mpres_at_mpMinLonF = lonW ; select a subregion
> mpres_at_mpMaxLonF = lonE
> mpres_at_mpMinLatF = latS
> mpres_at_mpMaxLatF = latN
>
> I get the following result (pic12) which is even more incorrect, as it
> is just a zoomed-in version.
> How else can I get the data to correctly display over its particular region?
> Sincerely,
> Erik
>
> On Oct 19, 2007, at 12:15 PM, Adam Phillips wrote:
>
>
>
> Hi Erik,
>
> If you want a map underneath your streamlines, then
> gsn_csm_streamline_map should be the function that you use. I believe
> the WRF model uses 2D latitudes/longitudes, and thus in NCL you need to set:
> u_plane_at_lat2D = a->XLAT(0,:,:)
> u_plane_at_lon2D = a->XLONG(0,:,:)
> and the same for v_plane...
> Good luck,
> Adam
>
> Erik Noble wrote:
> Could I please have some help with plotting the regional U and V data
> onto a WRF map?
> Below is NCL code that takes the U and V data from WRF model output
> and invokes the NCL streamline function.
> I use the "gsn_csm_streamline" function (streamline plot draw without
> a map) because the U and V data is strictly just that, "regional
> data".
> When I try to use the gsn_csm_streamline_map_ce function along with a
> specific a subregion, NCL incorrectly assumes the data to represent
> the globe and zooms in on the regional data.
> Is there a way to plot the streamlines first and then, using the lat
> and lons from the data itself, display an outline of the map/continent
> beneath?
> Attached is a picture of the resulting 925 mb streamlines over West
> Africa, where you can already see a feint outline of Africa below.
> any suggestions?
> Thank you,
> Erik
> ; Example script to produce stremaline plots for a WRF real-data run,
> ; with the ARW coordinate dynamics option.
> 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
> ;
> ; The WRF ARW input file.
> ; This needs to have a ".nc" appended, so just do it.
> a =
> addfile("/Volumes/Data_and_Models/Model-Output/Athena/WRF-SOP3-Athena_Phys-3-1-1-1.nc","r")
> ; We generate plots, but what kind do we prefer?
> type = "x11"
> ; type = "pdf"
> ; type = "ps"
> ; type = "ncgm"
> wks = gsn_open_wks(type,"WRF_SOP3_Run-2_Sep10-13_3_1_1_Streamlines")
> ; Set some Basic Plot options
> ARWres = True
> ARWres_at_MainTitle = "REAL-TIME WRF"
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ;What times and how many time steps are in the data set?
> times = wrf_user_list_times(a) ; get 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 = (/925., 850., 700., 500., 300./) ; pressure levels
> to plot
> nlevels = dimsizes(pressure_levels) ; number of pressure
> levels
> do it = 324,324 ;;;;;Time step you want!!!!
> print("Working on time: " + times(it) )
> ARWres_at_TimeLabel = times(it) ; Set Valid time to use on plots
> mpres = True ; Create map background
> mpres_at_gsnMaximize = True
> mpres = ARWres
> mpres_at_tiMainString = "Streamlines" ; some tit
> ;mpres_at_FieldTitle = "Circulation" ; overwrite Field Title
> mpres_at_mpOutlineOn = True
> mpres_at_mpGeophysicalLineColor = "Black"
> mpres_at_mpGeophysicalLineThicknessF = "3.0"
> mpres_at_mpGridLineColor = "Black"
> mpres_at_mpLimbLineColor = "Black"
> mpres_at_mpNationalLineColor = "Black"
> mpres_at_mpPerimLineColor = "Black"
> mpres_at_mpUSStateLineColor = "Black"
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ;Variables
> p = wrf_user_getvar(a, "pressure",it) ; pressure is our vertical
> coordinate
> u = wrf_user_getvar(a,"ua",it) ; u averaged to mass points
> v = wrf_user_getvar(a,"va",it) ; v averaged to mass points
> z = wrf_user_getvar(a, "z",it) ; grid point height
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> do level = 0,nlevels-1 ; LOOP OVER LEVELS
> pressure = pressure_levels(level)
> z_plane = wrf_user_intrp3d( z,p,"h",pressure,0.,False)
> ; wrf_smooth_2d( z_plane, 1 )
> u_plane = wrf_user_intrp3d( u,p,"h",pressure,0.,False)
> v_plane = wrf_user_intrp3d( v,p,"h",pressure,0.,False)
> ; u_plane = u_plane*1.94386 ; kts
> ; v_plane = v_plane*1.94386 ; kts
> ; u_plane_at_units = "kts"
> ; v_plane_at_units = "kts"
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ;Plot
> if ( pressure .eq. 925 ) then
> plot = gsn_csm_streamline(wks,u_plane,v_plane,mpres)
> end if
> end do
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> end do ; END OF TIME LOOP
> end
> ------------------------------------------------------------------------
> ------------------------------------------------------------------------
> _______________________________________________
>
> On Oct 23, 2007, at 12:32 PM, Adam Phillips wrote:
>
>> Erik,
>>
>> In future emails if you can tell us which line of your code is
>> specified in the error messages you receive that would be helpful,
>> otherwise we are stuck counting... Also, using printVarSummary to
>> print out summaries of each of the variables from the the offending
>> line is also a good debugging method...
>>
>> Anyway, the error message you are receiving is due to your nlatV =
>> dimV(2) line.. It should actually be nlatV = dimV(1)...
>>
>> But more importantly, if the order of your script is indeed this:
>>
>> u = wrf_user_getvar(a,"ua",it) ; u averaged to mass points
>> v = wrf_user_getvar(a,"va",it) ; v averaged to mass points
>> z = wrf_user_getvar(a, "z",it) ; grid point height
>> printVarSummary(u)
>> printVarSummary(v)
>>
>> and the output from the printVarSummary's are this:
>>
>> >> Variable: u
>> >> Type: float
>> >> Total Size: 1636308 bytes
>> >> 409077 values
>> >> Number of Dimensions: 3
>> >> Dimensions and sizes: [27] x [109] x [139]
>> >> Coordinates:
>> >> Number Of Attributes: 2
>> >> description : u Velocity
>> >> units : m/s
>> >>
>> >>
>> >> Variable: v
>> >> Type: float
>> >> Total Size: 1636308 bytes
>> >> 409077 values
>> >> Number of Dimensions: 3
>> >> Dimensions and sizes: [27] x [109] x [139]
>> >> Coordinates:
>> >> Number Of Attributes: 2
>> >> description : v Velocity
>> >> units : m/s
>>
>> then u/v are likely on the same grid and you do not have to do the
>> U = 0.5*(u(level,:,0:mlonU-2)+u(level,:,1:mlonU-1))
>> V = 0.5*(v(level,0:nlatV-2,:)+v( level,1:nlatV-1,:))
>>
>> coding. I would be surprised if this was the case because WRF by
>> default outputs the u/v variables on staggered grids.. Adam
>>
>> Dennis Shea wrote:
>>> I f you look at the Example 7 carefully you will not that
>>> on the right hand sude the variable ar capitalized [U,V]
>>> while in the left they are small [u,v].
>>> Given that you are already using small u,v ... I would merely
>>> change the left hand side and plot the U,V
>>> U = 0.5*(u(level,:,0:mlonU-2)+u(level,:,1:mlonU-1))
>>> V = 0.5*(v(level,0:nlatV-2,:)+v( level,1:nlatV-1,:))
>>> copy_VarAtts (u,U)
>>> copy_varAtts (v,V)
>>> lat2d = a->XLAT(0,:,:)
>>> lon2d = a->XLONG(0,:,)
>>> U_at_lat2d = lat2d
>>> U_at_lon2d = lon2d
>>> V_at_lat2d = lat2d
>>> V_at_lon2d = lon2d
>>> Erik Noble wrote:
>>>> Hi. I rewrote the portion of the my code, following the the method
>>>> outlined in example 7 of the WRF lambert conformal page:
>>>>
>>>> http://www.ncl.ucar.edu/Applications/wrflc.shtml as follows:
>>>>
>>>> do it = 324,324 ;;;;;Time step you want!!!!
>>>> ;Variables
>>>> p = wrf_user_getvar(a, "pressure",it) ; pressure is our
>>>> vertical coordinate
>>>> u = wrf_user_getvar(a,"ua",it) ; u averaged to mass points
>>>> v = wrf_user_getvar(a,"va",it) ; v averaged to mass points
>>>> z = wrf_user_getvar(a, "z",it) ; grid point height
>>>> dimU = dimsizes(u) ; demo "dimsizes"
>>>> mlonU = dimU(2) ; number of "west_east_stag"
>>>> longitudes
>>>> dimV = dimsizes(v)
>>>> nlatV = dimV(2)
>>>> printVarSummary(u)
>>>> printVarSummary(v)
>>>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>>>> do level = 0,nlevels-1 ; LOOP OVER LEVELS
>>>>
>>>> pressure = pressure_levels(level)
>>>> u = 0.5*(u(level,:,0:mlonU-2)+u(level,:,1:mlonU-1))
>>>> v = 0.5*(v(level,0:nlatV-2,:)+v( level,1:nlatV-1,:))
>>>> And this is the error that I am getting;
>>>>
>>>> noble:/Volumes/Data_and_Models/ncl_scripts enoble$ ncl
>>>> wrf_Streamlines.ncl
>>>> Copyright (C) 1995-2007 - All Rights Reserved
>>>> University Corporation for Atmospheric Research
>>>> NCAR Command Language Version 4.3.1
>>>> The use of this software is governed by a License Agreement.
>>>> See http://www.ncl.ucar.edu/ for more details.
>>>> (0) Working on time: 2006-09-10_12:00:00
>>>>
>>>>
>>>> Variable: u
>>>> Type: float
>>>> Total Size: 1636308 bytes
>>>> 409077 values
>>>> Number of Dimensions: 3
>>>> Dimensions and sizes: [27] x [109] x [139]
>>>> Coordinates:
>>>> Number Of Attributes: 2
>>>> description : u Velocity
>>>> units : m/s
>>>>
>>>>
>>>> Variable: v
>>>> Type: float
>>>> Total Size: 1636308 bytes
>>>> 409077 values
>>>> Number of Dimensions: 3
>>>> Dimensions and sizes: [27] x [109] x [139]
>>>> Coordinates:
>>>> Number Of Attributes: 2
>>>> description : v Velocity
>>>> units : m/s
>>>> fatal:Number of dimensions on right hand side do not match number of
>>>> dimension
>>>> n left hand side
>>>> fatal:Execute: Error occurred at or near line 70 in file
>>>> wrf_Streamlines.ncl
>>>>
>>>> noble:/Volumes/Data_and_Models/ncl_scripts enoble$
>>>>
>>>> How do I fix this?
>>>> -Erik
>>>> _______________________________________________
>>>> ncl-talk mailing list
>>>> ncl-talk_at_ucar.edu <mailto:ncl-talk_at_ucar.edu>
>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>> _______________________________________________
>>> ncl-talk mailing list
>>> ncl-talk_at_ucar.edu <mailto:ncl-talk_at_ucar.edu>
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>> --
>> --------------------------------------------------------------
>> Adam Phillips asphilli_at_ucar.edu <mailto:asphilli_at_ucar.edu>
>> National Center for Atmospheric Research tel: (303) 497-1726
>> ESSL/CGD/CAS fax: (303) 497-1333
>> P.O. Box 3000
>> Boulder, CO 80307-3000 http://www.cgd.ucar.edu/cas/asphilli
>>
>

-- 
--------------------------------------------------------------
Adam Phillips			             asphilli_at_ucar.edu
National Center for Atmospheric Research   tel: (303) 497-1726
ESSL/CGD/CAS                               fax: (303) 497-1333
P.O. Box 3000				
Boulder, CO 80307-3000	  http://www.cgd.ucar.edu/cas/asphilli
_______________________________________________
ncl-talk mailing list
ncl-talk_at_ucar.edu
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Tue Oct 23 2007 - 15:28:57 MDT

This archive was generated by hypermail 2.2.0 : Wed Oct 24 2007 - 09:14:00 MDT