[ncl-talk] climo_3.ncl + loop + conventional subscripts

Dennis Shea shea at ucar.edu
Sun Oct 27 15:01:26 MDT 2019


The following is untested. However, it gives you the idea.
It is your responsibility to look up documentation and look at the printed
output.



 fils = systemfunc ("ls precip.*.nc") ; file paths
 a    = addfiles (fils, "r")
 printVarSummary (a)

;************************************************
; Read the file
;************************************************

  ymdStrt = 19790101                            ; climatology start
year/month/day
  ymdLast = 20151231                             ;      "      last
 year/month/day
  time    = a[:]->time
  yyyymmdd= cd_calendar(time, -2)       ; note -2
print(yyyymmdd)
  ntStrt = ind(yyyymmdd.eq.ymdStrt)                ; start time index
  ntLast = ind(yyyymmdd.eq.ymdLast)              ; last  time index
print("ntStrt="+ntStrt+"  ntLast="+ntLast)
   prc   = a[:]->precip(ntStrt:ntLast,:,:)               ; read precip from
all files

   printVarSummary (prc)
   printMinMax (prc,0)
  print("=====")
                   ; you decide what you need
   prcMonth = *calculate_monthly_values(*
<http://www.ncl.ucar.edu/Document/Functions/Contributed/calculate_monthly_values.shtml>prc,
"avg", 0, False)   ; monthly mean
 ;;prcMonth = *calculate_monthly_values(*
<http://www.ncl.ucar.edu/Document/Functions/Contributed/calculate_monthly_values.shtml>prc,
"sum", 0, False)  ; monthly total

   printVarSummary (prcMonth)
   printMinMax (prcMonth,0)
   print("=====")
   opt = True                          ; examine data distribute
[information only]
   opt at PrintStat = True
   stat_prc = stat_dispersion(prcMonth, opt )

   prcClm   = *clmMonTLL*
<http://www.ncl.ucar.edu/Document/Functions/Contributed/clmMonTLL.shtml>(prcMonth)
; Compute monthly climatology
   printVarSummary (prcClm)
   printMinMax (prcClm,0)

====
Add the plot code.
Look carefully at: http://www.ncl.ucar.edu/Applications/HiResPrc.shtml

EG: The following may be inappropriate for your dataset but this is the
idea.

 res at cnLevelSelectionMode = "ExplicitLevels"
 res at cnLevels             = (/ 0.01, 0.02, 0.04, 0.08, 0.16, \
                               0.32, 0.64, 0.96/)
 res at cnFillColors         = (/"white","cyan", "green","yellow",\
                              "darkorange","red","magenta","purple",\
                              "black"/)


On Sun, Oct 27, 2019 at 2:48 PM Dennis Shea <shea at ucar.edu> wrote:

>
> [0]
> As noted on ncl-talk many time the first rule of data processing is *Look
> at your data!"
> The temperature data were *monthly means* [K].
> The precipitation data are* daily* totals [mm/day]. How do I know? I
> looked at the README file at the URL. You should do that.
>
> [1]
> a =* addfiles*
> <http://www.ncl.ucar.edu/Document/Functions/Built-in/addfiles.shtml>
> (fils, "r")
>
> You have the following which is appropriate for an *addfile*
> <http://www.ncl.ucar.edu/Document/Functions/Built-in/addfile.shtml>
> reference
>        time   = a->time
>
> However, you are using *addfiles*
> <http://www.ncl.ucar.edu/Document/Functions/Built-in/addfiles.shtml> It
> should be:
>        time   = a*[:]*->time   ; *[:]* syntax is appropriate for an
> *addfiles* reference
>
> [2]
> Did you read the documentation for *cd_calendar *
> <http://www.ncl.ucar.edu/Document/Functions/Built-in/cd_calendar.shtml>as
> I suggested? The whole point of using *cd_calendar* is to return
> 'human-readable' start and end times from which the actual start/stop time
> index values could be determined.. Your script uses the actual time values
> [ 692496 , 1016808 ] mixed with human-readable yyyymm.  Further, you used
>
>     yyyymm = cd_calendar(time, *-1*)         ; -1 returns YYYYMM
>
> The precipitation files contain [as noted] daily values. You must use a
> human readable time like YYYYMMDD [eg: 19790101]. The us
>
>     yyyymmdd = cd_calendar(time, *-2*)         ; -2 returns YYYYMMDD
>
> This converts the actual time values [eg:      692496] to 'human-readable'
> yyyymmdd [year=>yyyy, month=>mm, day=dd>* 692496=>197901*01] values.
>
>
> ===
> Also, your script will need to be changed to accommodate the precipitation
> variable plot values
> There are many examples of how to do this.
>
> eg: *http://www.ncl.ucar.edu/Applications/HiResPrc.shtml*
> <http://www.ncl.ucar.edu/Applications/HiResPrc.shtml>
> ===
> In a previous email, you mentioned that you converted degrees Kelvin to C.
> It is best to self-document the  variable with its new units.
>
>    tmp = tmp-273.15
>    tmp at units = "degC"  ; this will be automatically plotted
> ===
> Good Luck
>
>
>
>
>
>
> On Sat, Oct 26, 2019 at 8:19 AM zoe jacobs <zoejacobs1990 at gmail.com>
> wrote:
>
>> Dear all,
>> I would like to have climatology precipitation like temperature
>> climatology discussed through previous emails. However, there are
>> precipitation files for each year from 1979- 2015 (
>> ftp://ftp.cdc.noaa.gov/Datasets/cpc_global_precip/). In the first step,
>> I tried to merge all those files as one, and then use Dennis script to plot
>> a map. I faced some errors and need your advice please. the script I am
>> running has been attached.Many thanks in advance,
>> Best regardsn Sat, Oct 26, 2019 at 3:20 PM zoe jacobs <
>> zoejacobs1990 at gmail.com> wrote:
>>
>>> I just noticed that I could use res at gsnRightString = " deg C" .
>>> Once again, many thanks for your guide.
>>> Best wishes,
>>>
>>>
>>> On Fri, Oct 25, 2019 at 11:07 PM Dennis Shea <shea at ucar.edu> wrote:
>>>
>>>> I forgot to mention that plots can be done in 'portrait' [default] or
>>>> 'landscape' mode.
>>>>
>>>> You should experiment to see which is appropriate for your needs.
>>>>
>>>> resP at gsnPaperOrientation = "landscape"        ; "portrait" is default
>>>>
>>>>   do nmo=0,11                                   ; loop over the months
>>>>      res at gsnLeftString     = months(nmo)
>>>>      plot(nmo) = gsn_csm_contour_map(wks,tmpClm(nmo,:,:), res)  ;
>>>> create plot
>>>>   end do
>>>>
>>>>   if (resP at gsnPaperOrientation .eq."landscape") then
>>>>       gsn_panel(wks,plot,(/3,4/),resP)
>>>> ; 3 rows x 4 columns
>>>>   else
>>>>       gsn_panel(wks,plot,(/4,3/),resP)            ; portrait         ;
>>>> 4 rows x 3 columns
>>>>   end if
>>>>
>>>> On Fri, Oct 25, 2019 at 12:15 PM Dennis Shea <shea at ucar.edu> wrote:
>>>>
>>>>> Hello,
>>>>>
>>>>> We do not usually do this but this is offline from ncl-talk.
>>>>>
>>>>> [1]
>>>>> Conventional subscripts start at 1 for Fortran and Matlab.  Hence, for
>>>>> 12 months [nmos=12]
>>>>>      Fortran: do nmo=1,nmos
>>>>>      Matlab: for 1:nmos   or 1:1:nmos  [I am not a Matlab user so this
>>>>> is a guess.]
>>>>> NCL 'conventional' subscripts start at 0 [like C/C++/IDL/Python]
>>>>>      NCL: do nmo=0,ntim-1
>>>>>
>>>>> [2] Dave illustrated 'coordinate subscripting' which uses the *{*...
>>>>> *}* syntax.
>>>>>
>>>>> The script I am attaching uses conventional subscripting. Just like
>>>>> Matlab/Fortran except it uses a range fro 0 to 11 [12 elements].
>>>>>
>>>>> [3] Please carefully read the documentation for *clmMonTLL*
>>>>> <http://www.ncl.ucar.edu/Document/Functions/Contributed/clmMonTLL.shtml>,
>>>>> *cd_calendar*
>>>>> <http://www.ncl.ucar.edu/Document/Functions/Built-in/cd_calendar.shtml>
>>>>> ,* ind*
>>>>> <http://www.ncl.ucar.edu/Document/Functions/Built-in/ind.shtml>
>>>>>
>>>>> [4] Look at the output from '*printVarSummary*
>>>>> <http://www.ncl.ucar.edu/Document/Functions/Built-in/printVarSummary.shtml>'.
>>>>> USE *printVarSummary *frequently.
>>>>>      Note the dimensions and added attributes:
>>>>>      time_op_ncl : Climatology: 29 years
>>>>>      info : function clmMonTLL: contributed.ncl
>>>>>
>>>>> [5] NCL offer a large number of color tables [palettes]:
>>>>>
>>>>> *http://www.ncl.ucar.edu/Document/Graphics/color_table_gallery.shtml*
>>>>> <http://www.ncl.ucar.edu/Document/Graphics/color_table_gallery.shtml>
>>>>>      Using a color map that emphasizes features can be very useful.
>>>>>
>>>>> [6] You must invest the time to learn any new language.
>>>>> Karin Meier-Fliesher and Michael Bottinger [DKRZ] wrote a wonderful
>>>>> tutorial.
>>>>>
>>>>> *http://www.ncl.ucar.edu/Document/Manuals/NCL_User_Guide/*
>>>>> <http://www.ncl.ucar.edu/Document/Manuals/NCL_User_Guide/>
>>>>>
>>>>> I suggest you read it.
>>>>>
>>>>> Good luck
>>>>> ==============
>>>>> *%>* ncl  zoe_jacobs.hgcn_cams.ncl
>>>>>
>>>>> will produce a png file.
>>>>>
>>>>>
>>>>>
>>>>> On Fri, Oct 25, 2019 at 11:16 AM Dave Allured - NOAA Affiliate via
>>>>> ncl-talk <ncl-talk at ucar.edu> wrote:
>>>>>
>>>>>> 1.  In your example, please notice that the month loop has a stride
>>>>>> of 3.  Therefore it is making plots for only four months:  January, April,
>>>>>> July, and October.
>>>>>>
>>>>>> If you want all 12 months, remove the stride.  Then make sure the
>>>>>> array "plot" is dimensioned 12.  I think you can interchangeably use either
>>>>>> a 1-D or 2-D plot array (i.e. 3*4), your choice.
>>>>>>
>>>>>> The secondary index "i" is used to write the four plots into
>>>>>> positions 0, 1, 2, 3 in the "plot" graphics array.  If you are making 12
>>>>>> plots, you do not really need to use a secondary index.
>>>>>>
>>>>>> Panel plotting uses a graphics array containing multiple plots.
>>>>>> Please study basic examples and documentation for making panel plots.
>>>>>>
>>>>>> 2.  That data file has full coordinates.  Therefore you can use
>>>>>> either conventional or coordinate subscripting, or mixed, your choice.
>>>>>>
>>>>>> I suggest using one of the date functions with coordinate
>>>>>> subscripting, to index the time subset that you want.  Something like this:
>>>>>>
>>>>>>     time_units = f->air&time at units
>>>>>>     time1 = cd_inv_calendar (year1, 1, 1, 0, 0, 0, time_units, 0)
>>>>>>     time2 = cd_inv_calendar (year2, 12, 31, 23, 0, 0, time_units, 0)
>>>>>>     air_subset = f->air({time1:time2},:,:)
>>>>>>
>>>>>> Use printVarSummary to ensure that the subset has the dimension sizes
>>>>>> that you expect.  Then you can proceed to compute the climatology.
>>>>>>
>>>>>>
>>>>>> On Fri, Oct 25, 2019 at 10:14 AM zoe jacobs via ncl-talk <
>>>>>> ncl-talk at ucar.edu> wrote:
>>>>>>
>>>>>>> Hi all,
>>>>>>> regarding climo_3.ncl (
>>>>>>> https://www.ncl.ucar.edu/Applications/Scripts/climo_3.ncl) ,
>>>>>>> I have 2 questions:
>>>>>>> 1. I need to plot all months on one panel (say 3*4), and cannot
>>>>>>> understand the logic behind the below loop, which used in the climo_3.ncl
>>>>>>> script :
>>>>>>>
>>>>>>> i = -1                                        ; Climatologies
>>>>>>>   do nmo=0,11,3                                 ; loop over the months
>>>>>>>      i = i+1
>>>>>>>      res at gsnCenterString   = months(nmo)+":"+time(0)/100 +"-"+ time(ntim-1)/100
>>>>>>>      plot(i) = gsn_csm_contour_map(wks,prcClm(nmo,:,:), res)  ; create plot
>>>>>>>   end do
>>>>>>>
>>>>>>> How  does it  work??/
>>>>>>>
>>>>>>> 2. I would like to show climotology temperature from 1987- 2015. Data which I am using is air.mon.mean.nc  ( https://www.esrl.noaa.gov/psd/data/gridded/data.ghcncams.html)  and it is using conventional subscripts. I am not familiar with conventional subscripts. So how can I convert coordinate subscripting to conventional subscripts?!!
>>>>>>>
>>>>>>> Please kindly advice me .
>>>>>>>
>>>>>>> Many thanks in advance,
>>>>>>>
>>>>>>> Best regards,
>>>>>>>
>>>>>>> _______________________________________________
>>>>>> ncl-talk mailing list
>>>>>> ncl-talk at ucar.edu
>>>>>> List instructions, subscriber options, unsubscribe:
>>>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>>
>>>>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20191027/f91ff9f0/attachment.html>


More information about the ncl-talk mailing list