[ncl-talk] Temp advection, wrf_contour fails

Dennis Shea shea at ucar.edu
Wed Apr 20 16:30:47 MDT 2016


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


More information about the ncl-talk mailing list