[ncl-talk] Help in saving values in looping

David Brown dbrown at ucar.edu
Wed Jun 7 14:24:02 MDT 2017


Hi Andreas,
There are a couple of suggestions on the doc page for the where
function about how to avoid the division by 0 issue:

Caveat: NCL evaluates the entire expression inside the where function,
which will cause problems if it's an expression NCL would normally
fail on. For example, if you try to execute:

  yinv = where(y.ne.0, 1./y, 0.)

and y contains values equal to 0, you will get a fatal error:

  fatal:divide: Division by 0, Can't continue

There are two approaches to use where.

yinv = 1. / where(y.ne.0,y,y at _FillValue)

Set y to a _FillValue value where it is equal to 0, and then can do the divide:

  y    = where(y.ne.0,y,y at _FillValue)
  yinv = 1. / y

If you further want to set y back to where it originally contained zeros:

  y = where(ismissing(y),0,y)


Hope this helps.
 -dave

On Wed, Jun 7, 2017 at 1:16 PM, Andreas Chrysanthou <eeac at leeds.ac.uk> wrote:
> Following up on the loop I created, when I tried to add some criteria based
> on the slope of the linear trend and the regression coefficients in order to
> save them to a new array I got the error:
>
> fatal:divide: Division by 0, Can't continue
> fatal:Div: operator failed, can't continue
>
> I want to store the calculated value only when it meets one of the following
> criteria.
> Any ideas on how to improve the syntax? Or ways to bypass the problem?
>
> I’ve read on the NCL website in the lazy evaluation section that: "If the
> left operand is an array, you cannot depend on a left-hand side test to
> guard against possible error conditions resulting from evaluating an
> expression on the right side. For example, you cannot avoid division by 0
> and the consequences of an error result…”
>
> The snippet of the code is:
>
> finarr_nies = new((/22,22,1/),typeof(value_test))
> finarr_nies!0 = "startyear"
> finarr_nies&startyear = ispan(1980,2001,1)
> finarr_nies!1 = "endyear"
> finarr_nies&endyear = ispan(1989,2010,1)
>
> do i=1980,2001,1
>    do j=1989,2010,1
>         ntimes3      =  ispan(i,j,1)
>         start_date   = i
>         end_date    = j
>         test_rc_nies = regline_stats(ntimes3,value_test_nies({i:j}))
>  if (j-i .ge. 9) then
>             finarr_nies({i},{j},0) = where(test_rc_nies at b(1) .gt. 0 .and.
> test_rc_nies at b95(0) .gt. 0 , finarr_nies({i},{j},0), value_test at _FillValue )
>         else
>             finarr_nies({i},{j},0) = where(test_rc_nies at b(1) .lt. 0 .and.
> test_rc_nies at b95(1) .lt. 0 , finarr_nies({i},{j},0), value_test at _FillValue )
>         end if
>         delete(ntimes3)
>     end do
> end do
>
>
> Cheers,
> Andreas
>
> On 6 Jun 2017, at 13:22, Andreas Chrysanthou <eeac at leeds.ac.uk> wrote:
>
> I really appreciate your help Adam and Dennis.
>
> I could manage to extract the one value I wanted from the slope and the
> regression coefficients, pass it on to the new array outside of the loop and
> improved the code by writing the loop with variables instead of constants.
>
> Thanks so much.
>
>
>
> On 6 Jun 2017, at 01:29, Dennis Shea <shea at ucar.edu> wrote:
>
> I am not sure what is meant by: "There was one not insignificant thing I
> didn’t mention"
>
> ---
> As you note:
>
> "the regline_stats function for the b calculates  both the b(0) +b(1*)x of
> the linear regression, while the b95 calculates both the 2.5% and 97.5%
> regression coefficient confidence intervals. So when it tries to save to
> that cell, there are two values for the @b and @b95."
>
> Change:
>         trend_slope=test_rc at b*10  ; prints slope of the linear trend for the
> selected range of years
>         reg_coeff=test_rc at b95         ; prints 2.5% and 97.5% regression
> coefficient confidence intervals
>
> To extract only one value
>         trend_slope=test_rc at b(1) *10  ; prints slope of the linear trend for
> the selected range of years
>         reg_coeff=test_rc at b95(1)          ; prints 2.5% and 97.5% regression
> coefficient confidence intervals
>
> ---
>
> "4 different columns (start_date, end_date,trend_slope,reg_coeff)
>
> ---
> The following create 4 (well 3) different columns. I may give you an idea
> how to proceed.
> Also, good programming practice is to use variables and not constants.
>
> nyrSpan   = 10
> nyrCrit   = 10
> nyrStrt1  = 1980
> nyrStrt2  = 2001
> nyrLast1  = nyrStrt1+nyrSpan-1
> nyrLast2  = 2010
>
> NTIM      = 1000    ; arbitrary max possible # elements
> strtYear  = new( NTIM, "integer")
> lastYear  = new( NTIM, "integer")
> trend     = new( NTIM, "float")
> rc        = new( NTIM, "float")
>
> nt        = -1
> do nyrStrt=nyrStrt1,nyrStrt2
>    do nyrLast=nyrLast1,nyrLast2
>       nyrs = nyrLast-nyrStrt+1
>       if (nyrs.ge.nyrCrit) then
>           test_rc =
> regline_stats(tofloat(ispan(1,nyrSpan,1)),random_normal(-1,1,nyrSpan))
>           nt = nt+1
>           strtYear(nt) = nyrStrt
>           lastYear(nt) = nyrLast
>           spanYear     = nyrLast-nyrStrt+1
>           rc(nt)       = (/ test_rc /)   ; (/.../) no meta data ; same as
> b(1)
>          ;trend(nt)    = ?????
>           print(strtYear(nt)+"  "+lastYear(nt)+"  "+spanYear+"  "+rc(nt))
>        end if
>    end do
> end do
>
>  print("=================================================")
>  ntim = nt+1    ; # of elements
>  print("ntim="+ntim)
>
>
>
> On Mon, Jun 5, 2017 at 11:10 AM, Andreas Chrysanthou <eeac at leeds.ac.uk>
> wrote:
>>
>> There was one not insignificant thing I didn’t mention.
>>
>> the regline_stats function for the b calculates  both the b(0) +b(1*)x of
>> the linear regression, while the b95 calculates both the 2.5% and 97.5%
>> regression coefficient confidence intervals.
>> So when it tries to save to that cell, there are two values for the @b and
>> @b95.
>>
>> That’s why I got the error of : "fatal:Dimension sizes on right hand side
>> of assignment do not match dimension sizes of left hand side"
>>
>> I only want the slope of the liner trend (b(1)) and my criterion if I need
>> to save these values is that the confidence intervals should not be negative
>> (that would probably need to be included in the if loop).
>> But I don’t know how to keep only one value of those pair of calculated
>> stats.
>>
>> Cheers,
>> Andreas
>>
>>
>> On 5 Jun 2017, at 17:17, Adam Phillips <asphilli at ucar.edu> wrote:
>>
>> Hi Andreas,
>> As there are singular values for reg_coef and trend_slope for each start
>> and end time, I think it would be best to simply create a 3D array
>> dimensioned (start_date, end_date,2), with the last dimension (0) =
>> trend_slope and (1) representing reg_coeff. If you must have a 4-dimensional
>> array then you should create an array that is dimensioned 22 x 22 x 2 x 2,
>> but (:,:,1,0) and (:,:,0,1) will be missing. In the example below I show how
>> to create both types of arrays:
>>
>> finarr = new((/22,22,2/),typeof(value_test))
>> finarr!0 = "startyear"
>> finarr&startyear = ispan(1980,2001,1)
>> finarr!1 = "endyear"
>> finarr&endyear = ispan(1989,2010,1)
>> finarrZ = new((/22,22,2,2/),typeof(value_test))
>> finarrZ!0 = "startyear"
>> finarrZ&startyear = ispan(1980,2001,1)
>> finarrZ!1 = "endyear"
>> finarrZ&endyear = ispan(1989,2010,1)
>>
>> do i=1980,2001,1
>>    do j=1989,2010,1
>>         ntimes=ispan(i,j,1)
>>         if (j-i .ge. 9 ) then
>>            test_rc=regline_stats(ntimes,value_test({i:j}))
>>            print(test_rc)
>>            start_date=i
>>            end_date=j
>>            finarr({i},{j},0) = test_rc at b*10  ; prints slope of the linear
>> trend for the selected range of years
>>            finarr({i},{j},1) = test_rc at b95         ; prints 2.5% and 97.5%
>> regression coefficient confidence intervals
>>            finarrZ({i},{j},0,0) = test_rc at b*10
>>            finarrZ({i},{j},1,1) =  test_rc at b95
>>        else
>>            print(“this age range has not been selected, since its less
>> than 9 years")
>>        end if
>>        delete(ntimes)
>>    end do
>> end do
>>
>> Hope that helps. If you have any further questions please respond to the
>> ncl-talk email list.
>> Adam
>>
>>
>> On Mon, Jun 5, 2017 at 9:48 AM, Andreas Chrysanthou <eeac at leeds.ac.uk>
>> wrote:
>>>
>>> Hi NCL users,
>>>
>>>
>>> I created a loop for calculating a linear trend sensitivity for different
>>> start and end dates.
>>>
>>> My timeseries span over 1980-2010 and I want to calculate the linear
>>> trend of the timeseries for when the start - end date are more than >=
>>> 9years.
>>>
>>> (1980-1989, 1980-1990, 1980-1991,..., 1980-2010
>>>
>>> 1981-1990, 1981-1991, ..., 1981-2010
>>>
>>> ...
>>>
>>> ...
>>>
>>> ...
>>>
>>> 2000-2009, 2001-2010
>>>
>>> 2001-2010)
>>>
>>>
>>>
>>> I’ve created the loops that calculate those values but I need to save
>>> them in an array as 4 different columns (start_date,
>>> end_date,trend_slope,reg_coeff)
>>> My aim is to plot those as a a contour box plot with x axis the start
>>> date and y axis the end date, and the calculated (i’m gonna filter out those
>>> values with some criteria) trend values based on a colormap.
>>>
>>> Can you help me in order to save the values to that new array so I can
>>> plot it?
>>>
>>> A snippet of the code follows:
>>>
>>> do i=1980,2001,1
>>>    do j=1989,2010,1
>>>         ntimes=ispan(i,j,1)
>>> if (j-i .ge. 9 ) then
>>>         test_rc=regline_stats(ntimes,value_test({i:j}))
>>>         print(test_rc)
>>>         start_date=i
>>>         end_date=j
>>>         trend_slope=test_rc at b*10  ; prints slope of the linear trend for
>>> the selected range of years
>>>         reg_coeff=test_rc at b95         ; prints 2.5% and 97.5% regression
>>> coefficient confidence intervals
>>> else print(“this age range has not been selected, since its less than 9
>>> years")
>>>     end if
>>>       delete(ntimes)
>>> end do
>>> end do
>>>
>>> Cheers,
>>> Andreas
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>> Andreas
>>>
>>>
>>>
>>> _______________________________________________
>>> ncl-talk mailing list
>>> ncl-talk at ucar.edu
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>
>>
>>
>>
>> --
>> Adam Phillips
>> Associate Scientist,  Climate and Global Dynamics Laboratory, NCAR
>> www.cgd.ucar.edu/staff/asphilli/   303-497-1726
>>
>>
>>
>> _______________________________________________
>> ncl-talk mailing list
>> ncl-talk at ucar.edu
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>
>
>
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>


More information about the ncl-talk mailing list