[ncl-talk] Suspected bug in calculate_daily_values

[Dave Macleod]

Hi,

Long time user, first time emailer ... thanks for all the hard work on NCL! FYI, I believe I have found a bug in calculate_daily_values.

The argument "nDim" is meant to identify the time dimension, but in fact the function only works correctly when nDim=0 (in the case of 2,3,4-dimensional data), and nDim=1 (in the case of 5D data).

I have copied a MWE below to demonstrate. But I have traced the problem to point of the statistical operation in the function (e.g. line 18402 in /contributed.ncl, NCL v 6.4.0):

> xStat(:,nt,:,:,:)  = dim_avg_n( x(:,it,:,:,:), nDim)

Here "nDim" specifies the dimension on which to carry out the averaging and "it" is the index corresponding to day "nt" . Since "it" is hard-coded on dimension #1, any value of nDim not equal to 1 will not produce the correct result.

Similarly for the 2,3,4D case,
> xStat(nt,:,:,:)  = dim_avg_n( x(it,:,:,:), nDim)

...requires nDim=0.

A quick fix would be to modify the 5D case such that it require nDim=0, then add a error condition that nDim=0. Or something more clever to preserve the functionality of nDim! (I think there is also an issue where the xStat is defined according to dimx).

NB, from the code I think the same error may affect calculate_monthly_values and possibly calculate_segment_values, though I haven't tested this.

All the best,
Dave MacLeod
-------------------------------------
Postdoctoral Research Scientist
Atmospheric Oceanic and Planetary Physics, University of Oxford


MWE demonstrating the problem described above:

load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
begin

; 2D, 3D, 4D data: only works if nDim=0
x4D=random_uniform(-1.0,1.0,(/4,4,4,4/))        ; dummy data
nDim=0                                          ; choose a dimension to be 'time' - CHANGE FROM A VALUE DIFFERENT FROM ZERO TO BREAK
x4D!nDim="time"                         ; add metadata; 4 points, over 2 days
x4D&time=fspan(0.25,1.75,4)
x4D&time at units="days since 2017-09-12"
printVarSummary(x4D)
x4DDay = calculate_daily_values (x4D, "avg",nDim, False)
printVarSummary(x4DDay)


; 5D data: only works if nDim=1
x5D=random_uniform(-1.0,1.0,(/4,4,4,4,4/))      ; dummy data
nDim=1                                          ; choose a dimension to be 'time' - CHANGE TO A VALUE DIFFERENT FROM ONE TO BREAK
x5D!nDim="time"                         ; add metadata; 4 points, over 2 days
x5D&time=fspan(0.25,1.75,4)
x5D&time at units="days since 2017-09-12"
printVarSummary(x5D)
x5DDay = calculate_daily_values (x5D, "avg",nDim, False)
printVarSummary(x5DDay)
end
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20170912/66b844cc/attachment.html>