[ncl-talk] Temp advection, wrf_contour fails

Ryan Connelly rconne01 at gmail.com
Thu Apr 21 09:55:32 MDT 2016


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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160421/a9cf4189/attachment.html 


More information about the ncl-talk mailing list