[ncl-talk] Temp advection, wrf_contour fails
Ryan Connelly
rconne01 at gmail.com
Fri Apr 22 13:22:47 MDT 2016
Hi all,
It seems like the grad_latlon_cfd function will work once I convert from
curvilinear grid to rectilinear, which seems straightforward enough with
NCL's ESMF regridding function.
Only problem now is... contributed.ncl_640 isn't where it's supposed to
be. I need to download it so that I can use the grad_latlon_cfd function
in V6.4.0, as described here:
http://ncl.ucar.edu/Applications/gradients.shtml, but the link to it is a
404 page: http://www.cgd.ucar.edu/~shea/contributed.ncl_640
Thanks,
Ryan
On Thu, Apr 21, 2016 at 12:24 PM, Dennis Shea <shea at ucar.edu> wrote:
> I think you would have to ask wrf-help at ucar.edu
>
> I don't think this is an NCL-issue per-se. NCL is providing an interface
> to fortran subroutines provided by the WRF project.
>
> You state that the "latitude off by ~1 deg" ? Not sure how you get this.
>
> I'm not a WRF user/expert so I can not be of much further help.
>
> ==
> Yes ... people have calculated gradients on curvilinear grids but they
> use ESMF regridding to interpolate to a rectilinear grid; perform the
> gradient calculations and then reinterpolate back to the original
> curvilinear grid.
>
> Not gradients related per se but it give you an idea of how to proceed ...
>
> https://www.ncl.ucar.edu/Applications/ESMF.shtml
> Example 5
> Example 29
>
> Example 30 is for a NARR grid (curvilinear). It illustrates that the
> errors of back and forth interpolation are small.: original values are
> 200-300 .... differeines are 0.05. That's pretty good.
>
> D
>
>
>
>
>
> On Thu, Apr 21, 2016 at 9:55 AM, Ryan Connelly <rconne01 at gmail.com> wrote:
>
>> Yeah that makes sense why wrf_user_intrp3d would fail. Your second
>> suggestion got me lat/lon minima and maxima that are close, but not close
>> enough, to the corners of my domain (latitude off by ~1 deg).
>>
>> I'm starting to think this isn't physically possible in NCL. The fact
>> that WRF is a curvilinear grid means that I couldn't even manually code up
>> conditionals and loops to calculate the gradient components, right? Since
>> advancing from i=2 to i=4, e.g., isn't moving due west-east.
>>
>> Has no one ever computed advections of WRF variables in NCL before?
>>
>>
>>
>> On Wed, Apr 20, 2016 at 6:30 PM, Dennis Shea <shea at ucar.edu> wrote:
>>
>>> Likely, I was not clear. The following file is needed.
>>>
>>> wrfout_d01_2016-02-04_12:00:00
>>>
>>> ===
>>> You should look at:
>>> lon_plane = wrf_user_intrp3d(lon,p,"h",pressure,0.,False)
>>> lat_plane = wrf_user_intrp3d(lat,p,"h",pressure,0.,False)
>>>
>>> printVarSummary(lon_plane)
>>> printMinMax(lon_plane,0)
>>>
>>> printVarSummary(lat_plane)
>>> printMinMax(lat_plane,0)
>>>
>>> Really, since lat/lon have no vertical extent .... how can you use
>>> wrf_user_intrp3d?
>>>
>>> https://www.ncl.ucar.edu/Document/Functions/WRF_arw/wrf_user_intrp3d.shtml
>>>
>>> The 1st arguments is:
>>>
>>> *var3d*
>>>
>>> Data on model levels that will be interpolated. The rightmost dimensions
>>> of this array is *nz* x *ny* x *nx*.
>>>
>>> I am not a WRF user ... perhaps someone more knowledgeable can provide
>>> an answer?
>>>
>>> Likely, you could 'trick' lat/lon to multiple dimensions via
>>>
>>> lat3d = conform(p, lat, (/1,2/))
>>> lon3d = conform(p, lon, (/1,2/))
>>>
>>> lon_plane = wrf_user_intrp3d(lon3d,p,"h",pressure,0.,False)
>>> lat_plane = wrf_user_intrp3d(lat3d,p,"h",pressure,0.,False)
>>>
>>>
>>>
>>>
>>> ===
>>> re: "realistic values" for dTdx instead of inf's:
>>>
>>> (0) min=-15.1425 max=18.5633
>>> (0) min=-24.1597 max=15.1911
>>>
>>> ===
>>> The actual gradients will likely *not* be these nice numbers.
>>>
>>> The (say): (T(n+1)-T(n-1))/(LAT(n+1)-LAT(n-1))
>>>
>>> (20-10)/(15-12.5)
>>>
>>> would yield 10/2.5 = 4 [degC/degLat] or degF in your case
>>>
>>> however, most gradients are (say) degC/m_or_km
>>>
>>> so that 2.5 degrees of lat would have to be converted to m or km.
>>>
>>> ===
>>> You are aware that the WRF grid is curvilinear ... not rectilinear?
>>>
>>>
>>>
>>> On Wed, Apr 20, 2016 at 11:49 AM, Ryan Connelly <rconne01 at gmail.com>
>>> wrote:
>>>
>>>> Hi Dennis,
>>>>
>>>> My file is now in ../incoming/wrf_temp_adv_rconne01.ncl.
>>>>
>>>> It doesn't like dimension reduction because the resulting argument
>>>> isn't the right size:
>>>>
>>>> fatal:center_finite_diff_n: r must either be a scalar, a 1D array the
>>>> same length as the dim-th dimemsion of q, or the same size as q
>>>>
>>>>
>>>> It also throws a fatal on the declaration of lon_plane and/or lat_plane:
>>>>
>>>> fatal:Subscript out of range, error in subscript #0
>>>> fatal:An error occurred reading dims
>>>> fatal:["Execute.c":8128]:Execute: Error occurred at or near line 224 in
>>>> file $NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl
>>>>
>>>> fatal:["Execute.c":8128]:Execute: Error occurred at or near line 97 in
>>>> file wrf_temp_adv.ncl
>>>>
>>>>
>>>> My best success so far has been with your first suggestion, merely
>>>> changing the last arg of the dTDx declaration from 0 to 1, which gave me
>>>> realistic values for dTdx instead of inf's:
>>>>
>>>> (0) min=-15.1425 max=18.5633
>>>> (0) min=-24.1597 max=15.1911
>>>>
>>>> But that ultimately still hangs up, though interestingly enough, it now
>>>> gets one step farther. Before, wrf_contour was where it stopped. Now, it
>>>> gets past wrf_contour, but fails on plot=.
>>>>
>>>> Thanks,
>>>> Ryan
>>>>
>>>> On Wed, Apr 20, 2016 at 1:19 PM, Ryan Connelly <rconne01 at gmail.com>
>>>> wrote:
>>>>
>>>>> Hi Mary,
>>>>>
>>>>> That shows that dTdx is the problem. (dTdx is first, dTdy second:)
>>>>>
>>>>> (0) min=-inf max=inf
>>>>> (0) min=-24.1597 max=15.1911
>>>>>
>>>>> I'll work on Dennis's suggestion now.
>>>>>
>>>>> Ryan
>>>>>
>>>>> On Wed, Apr 20, 2016 at 11:59 AM, Dennis Shea <shea at ucar.edu> wrote:
>>>>>
>>>>>> NOTE: I have not looked at this in detail.
>>>>>>
>>>>>> [1]
>>>>>> You could not use 'advect_variable' because it requires "[t]he array *must
>>>>>> be global and ordered south to north*." Also, this assumes a
>>>>>> rectilinear grid (eg: fixed or gaussian)
>>>>>>
>>>>>> Your data are not global :-(
>>>>>> Yourhorizontal grid is not rectilinear.
>>>>>>
>>>>>> Note: I could change ''advect_variable" which uses spherical
>>>>>> harmonics (hence, the global grid requirement) to derive the horizontal and
>>>>>> latitudinal gradients tusing the 6.4.0
>>>>>>
>>>>>> https://www.ncl.ucar.edu/Document/Functions/Contributed/grad_latlon_cfd.shtml
>>>>>> BUT that would take a day or two.
>>>>>>
>>>>>> [2]
>>>>>> center_finite_diff "the issue" .. yes and no
>>>>>>
>>>>>> The error occurs in center_finite_diff but the issue is, I speculate,
>>>>>> the arguments. Specifically, the use of lat and lon
>>>>>>
>>>>>> ====
>>>>>> lon = wrf_user_getvar(a,"lon",it)
>>>>>> lat = wrf_user_getvar(a,"lat",it)
>>>>>>
>>>>>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>>>>>> ; Interpolate to isobaric planes
>>>>>>
>>>>>> pressure = 850.
>>>>>>
>>>>>> tc_plane = wrf_user_intrp3d(tc,p,"h",pressure,0.,False)
>>>>>> z_plane = wrf_user_intrp3d(z_dec,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)
>>>>>>
>>>>>>
>>>>>> ; Define variables needed to calculate temp advection
>>>>>>
>>>>>> dTdx = center_finite_diff_n(tc_plane,lon,False,0,0)
>>>>>> dTdy = center_finite_diff_n(tc_plane,lat,False,0,0)
>>>>>> ===
>>>>>>
>>>>>> Maybe,
>>>>>>
>>>>>>
>>>>>> dTdx = center_finite_diff_n(tc_plane,lon,False,0,1) ; <==
>>>>>> dTdy = center_finite_diff_n(tc_plane,lat,False,0,0)
>>>>>>
>>>>>> Or [281] x [353]=>nlat,mlon
>>>>>>
>>>>>> dTdx =
>>>>>> center_finite_diff_n(tc_plane,lon(nlat/2,:),False,0,0)
>>>>>> dTdy = center_finite_diff_n(tc_plane,lat(:,mlon/2),False,0,0)
>>>>>>
>>>>>> Note that NCL's 'dimension reduction' (elimination of degenerate
>>>>>> dimensions), lon(nlat/2,:) and lat(:,mlon/2) are one dimensional.
>>>>>> Use printVarSummary( lon(nlat/2,:) ) to verify
>>>>>>
>>>>>> [3]
>>>>>>
>>>>>> [2] is a thumb in the dike approach. I think it would be best (if
>>>>>> possible) to calculate
>>>>>>
>>>>>> lon_plane =
>>>>>> lat_plane =
>>>>>>
>>>>>> dTdx = center_finite_diff_n(tc_plane,lon_plane,False,0,0)
>>>>>> dTdy = center_finite_diff_n(tc_plane,lat_plane,False,0,0)
>>>>>>
>>>>>> [4]
>>>>>> Maybe you could send a sample filw (one only) to
>>>>>>
>>>>>> ftp ftp.cgd.ucar.edu
>>>>>> anonymous
>>>>>> your_email
>>>>>> cd incoming
>>>>>> put .....
>>>>>> quit
>>>>>>
>>>>>> then let ncl_talk know the name of the file.
>>>>>>
>>>>>> On Tue, Apr 19, 2016 at 4:21 PM, Mary Haley <haley at ucar.edu> wrote:
>>>>>>
>>>>>>> Hi Ryan,
>>>>>>>
>>>>>>> That makes sense why you can't upgrade.
>>>>>>>
>>>>>>> As for the "inf", this is definitely an issue. To further determine
>>>>>>> where the problem is, I suggest doing a printMinMax on the individual
>>>>>>> variables in the calculation:
>>>>>>>
>>>>>>> printMinMax(dTdx,0)
>>>>>>> printMinMax(dTdy,0)
>>>>>>>
>>>>>>> etc.
>>>>>>>
>>>>>>> You may want to do this on every variable that you are getting with
>>>>>>> "wrf_user_getvar" to make sure those all look okay.
>>>>>>>
>>>>>>> If it looks like center_finite_diff is the issue, then revisit the
>>>>>>> documentation to make sure all the input parameters are being input
>>>>>>> correctly:
>>>>>>>
>>>>>>>
>>>>>>> http://www.ncl.ucar.edu/Document/Functions/Built-in/center_finite_diff.shtml
>>>>>>>
>>>>>>> --Mary
>>>>>>>
>>>>>>> On Tue, Apr 19, 2016 at 3:53 PM, Ryan Connelly <rconne01 at gmail.com>
>>>>>>> wrote:
>>>>>>>
>>>>>>>> Hi Mary,
>>>>>>>>
>>>>>>>> I use NCL that's installed on a shared server where WRF is
>>>>>>>> configured. I don't suspect I have permissions to remove an old version.
>>>>>>>> I recently updated to 6.3.0 on my local machine, but then I have to copy
>>>>>>>> the wrfout over, which is a bit time-consuming since I'm running a triple
>>>>>>>> nest down to 1.333 km! :O
>>>>>>>>
>>>>>>>> These extra commands that you suggested certainly give a clue that
>>>>>>>> I'm not calculating things right:
>>>>>>>>
>>>>>>>> Variable: temp_adv
>>>>>>>> Type: float
>>>>>>>> Total Size: 396772 bytes
>>>>>>>> 99193 values
>>>>>>>> Number of Dimensions: 2
>>>>>>>> Dimensions and sizes: [281] x [353]
>>>>>>>> Coordinates:
>>>>>>>> Number Of Attributes: 1
>>>>>>>> _FillValue : 9.96921e+36
>>>>>>>> (0) min=-inf max=inf
>>>>>>>>
>>>>>>>>
>>>>>>>> 281x353 is my grid size, so that's good, but the infinities
>>>>>>>> probably are not...
>>>>>>>>
>>>>>>>> Thanks,
>>>>>>>> Ryan
>>>>>>>>
>>>>>>>> On Tue, Apr 19, 2016 at 10:38 AM, Mary Haley <haley at ucar.edu>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> Hi Ryan,
>>>>>>>>>
>>>>>>>>> Is there a reason why you don't want to upgrade from V6.1.2? The
>>>>>>>>> version you have is over 3 years old.
>>>>>>>>>
>>>>>>>>> I can't be sure of why your temp_adv variable is not plotting,
>>>>>>>>> without actually being able to run the script.
>>>>>>>>>
>>>>>>>>> However, rather than printing the whole temp_adv array, what does
>>>>>>>>> the following report:
>>>>>>>>>
>>>>>>>>> printVarSummary(temp_adv)
>>>>>>>>> printMinMax(temp_adv,0)
>>>>>>>>>
>>>>>>>>> Sometimes this will give you a clue where the problem is, for
>>>>>>>>> example, your min/max values are off scale, or the size of the array is not
>>>>>>>>> what you were expecting.
>>>>>>>>>
>>>>>>>>> You can also try plotting temp_adv with a basic contour plot call:
>>>>>>>>>
>>>>>>>>> plot = gsn_csm_contour(wks,temp_adv,False)
>>>>>>>>>
>>>>>>>>> just to make sure that the problem isn't with wrf_contour itself.
>>>>>>>>>
>>>>>>>>> --Mary
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> On Mon, Apr 18, 2016 at 12:49 PM, Ryan Connelly <
>>>>>>>>> rconne01 at gmail.com> wrote:
>>>>>>>>>
>>>>>>>>>> Hi,
>>>>>>>>>>
>>>>>>>>>> Running 6.1.2, so can't use advect_variable without upgrading.
>>>>>>>>>> So instead I have...
>>>>>>>>>>
>>>>>>>>>> ;do it = 12,30,3 ; TIME LOOP
>>>>>>>>>> do it = 12,13,1
>>>>>>>>>>
>>>>>>>>>> print("Working on time: " + times(it) )
>>>>>>>>>> res at TimeLabel = times(it) ; Set Valid time to use on plots
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>>>>>>>>>> ; First get the variables we will need
>>>>>>>>>>
>>>>>>>>>> tc = wrf_user_getvar(a,"tc",it) ; 3D tc
>>>>>>>>>> td = wrf_user_getvar(a,"td",it) ; 3D td
>>>>>>>>>> u = wrf_user_getvar(a,"ua",it) ; 3D U at mass points
>>>>>>>>>> v = wrf_user_getvar(a,"va",it) ; 3D V at mass points
>>>>>>>>>> p = wrf_user_getvar(a, "pressure",it) ; pressure is our
>>>>>>>>>> vertical coordinate
>>>>>>>>>> z = wrf_user_getvar(a,"z",it) ; Full model height in
>>>>>>>>>> meters
>>>>>>>>>> z_dec = z/10. ; Height in decameters
>>>>>>>>>> lon = wrf_user_getvar(a,"lon",it)
>>>>>>>>>> lat = wrf_user_getvar(a,"lat",it)
>>>>>>>>>>
>>>>>>>>>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>>>>>>>>>> ; Interpolate to isobaric planes
>>>>>>>>>>
>>>>>>>>>> pressure = 850.
>>>>>>>>>>
>>>>>>>>>> tc_plane = wrf_user_intrp3d(tc,p,"h",pressure,0.,False)
>>>>>>>>>> z_plane = wrf_user_intrp3d(z_dec,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)
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> ; Define variables needed to calculate temp advection
>>>>>>>>>>
>>>>>>>>>> dTdx = center_finite_diff_n(tc_plane,lon,False,0,0)
>>>>>>>>>> dTdy = center_finite_diff_n(tc_plane,lat,False,0,0)
>>>>>>>>>>
>>>>>>>>>> temp_adv = u_plane*dTdx + v_plane*dTdy
>>>>>>>>>>
>>>>>>>>>> print(temp_adv)
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>>>>>>>>>>
>>>>>>>>>> ; Plotting options for Tc
>>>>>>>>>> opts = res
>>>>>>>>>> print("Defined ops")
>>>>>>>>>> opts at cnFillOn = True
>>>>>>>>>> print("Turned fill on")
>>>>>>>>>> opts at cnLabelMasking = True
>>>>>>>>>> opts at cnInfoLabelOn = False
>>>>>>>>>> opts at cnLineLabelPerimOn = False
>>>>>>>>>> opts at ContourParameters = (/ -60., 60., 3./)
>>>>>>>>>> print("Accepted contour parameters")
>>>>>>>>>> ;opts at gsnSpreadColorEnd = -3 ; End third from the last
>>>>>>>>>> color in color map
>>>>>>>>>> contour_tc = wrf_contour(a,wks,temp_adv,opts) ; <- Breaks
>>>>>>>>>> right here
>>>>>>>>>> print("wrf_contour called")
>>>>>>>>>> delete(opts)
>>>>>>>>>> print("opts deleted")
>>>>>>>>>>
>>>>>>>>>> print("Got past plotting options for Tc")
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> ; MAKE PLOTS
>>>>>>>>>> ;plot =
>>>>>>>>>> wrf_map_overlays(a,wks,(/contour_tc,contour_td,contour_z/),pltres,mpres)
>>>>>>>>>> plot = wrf_map_overlays(a,wks,(/contour_tc/),pltres,mpres)
>>>>>>>>>> ;plot =
>>>>>>>>>> wrf_map_overlays(a,wks,(/contour_td,vector/),pltres,mpres)
>>>>>>>>>>
>>>>>>>>>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>>>>>>>>>>
>>>>>>>>>> end do ; END OF TIME LOOP
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> When it prints out the full grid of temp_adv, or random indexes
>>>>>>>>>> of the array, I get reasonable values and no reason to suspect the grid is
>>>>>>>>>> undefined anywhere. Yet as you can see from my print statements, the call
>>>>>>>>>> to wrf_contour fails. It does not throw an error; it just hangs up. I
>>>>>>>>>> checked with just tc_plane instead, and that completed fine, so it has to
>>>>>>>>>> be an issue with the temp_adv grid I calculated.
>>>>>>>>>>
>>>>>>>>>> Thanks,
>>>>>>>>>> Ryan
>>>>>>>>>>
>>>>>>>>>> --
>>>>>>>>>> Ryan Connelly
>>>>>>>>>> M.S. Student in Atmospheric Sciences, Stony Brook University
>>>>>>>>>> B.S. in Meteorology with Minors in Mathematics and GIS, Valparaiso
>>>>>>>>>> University
>>>>>>>>>> rconne01 at gmail.com
>>>>>>>>>> ryan.connelly at stonybrook.edu
>>>>>>>>>>
>>>>>>>>>> _______________________________________________
>>>>>>>>>> ncl-talk mailing list
>>>>>>>>>> ncl-talk at ucar.edu
>>>>>>>>>> List instructions, subscriber options, unsubscribe:
>>>>>>>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> --
>>>>>>>> Ryan Connelly
>>>>>>>> M.S. Student in Atmospheric Sciences, Stony Brook University
>>>>>>>> B.S. in Meteorology with Minors in Mathematics and GIS, Valparaiso
>>>>>>>> University
>>>>>>>> rconne01 at gmail.com
>>>>>>>> ryan.connelly at stonybrook.edu
>>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> _______________________________________________
>>>>>>> ncl-talk mailing list
>>>>>>> ncl-talk at ucar.edu
>>>>>>> List instructions, subscriber options, unsubscribe:
>>>>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>>>>
>>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> Ryan Connelly
>>>>> M.S. Student in Atmospheric Sciences, Stony Brook University
>>>>> B.S. in Meteorology with Minors in Mathematics and GIS, Valparaiso
>>>>> University
>>>>> rconne01 at gmail.com
>>>>> ryan.connelly at stonybrook.edu
>>>>>
>>>>
>>>>
>>>>
>>>> --
>>>> Ryan Connelly
>>>> M.S. Student in Atmospheric Sciences, Stony Brook University
>>>> B.S. in Meteorology with Minors in Mathematics and GIS, Valparaiso
>>>> University
>>>> rconne01 at gmail.com
>>>> ryan.connelly at stonybrook.edu
>>>>
>>>
>>>
>>
>>
>> --
>> Ryan Connelly
>> M.S. Student in Atmospheric Sciences, Stony Brook University
>> B.S. in Meteorology with Minors in Mathematics and GIS, Valparaiso
>> University
>> rconne01 at gmail.com
>> ryan.connelly at stonybrook.edu
>>
>
>
--
Ryan Connelly
M.S. Student in Atmospheric Sciences, Stony Brook University
B.S. in Meteorology with Minors in Mathematics and GIS, Valparaiso
University
rconne01 at gmail.com
ryan.connelly at stonybrook.edu
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160422/283cfaeb/attachment.html
More information about the ncl-talk
mailing list