Potentially serious ncgen bug

  • Stubaan

    Stubaan - 2014-04-09

    Hi folks

    I have a series of nc files containing monthly average rates of nitrogen deposition, calculated over decades. My workflow works fine for all files covering decades from 1850-1859 through to and including 1960-1969. However, the process breaks down for all files covering decades from 1970-1979 to 2100-2109.

    I have found that ncgen appears to exclude the exponent when it is e-10 when working on cdl files that I created using ncks --cdl.

    The original data files can be found here: https://sourceforge.net/projects/ipcc3hourly/files/NCO%20ncgen/

    Following is my workflow, which I encourage you to reproduce:

    1) Average the 12 monthly values to obtain an annual average rate, per decade

    ncra ndep_ocn_1960-1969_1.9x2.5_FD.nc ndep_ocn_1960-1969_1.9x2.5_FD_ncra.nc

    2) Sum the forms of deposition to obtain a total nitrogen deposition

    ncap2 -s 'Ndep=NHx_wet_deposition+NHx_dry_deposition+NOy_wet_deposition+NOy_dry_deposition' ndep_ocn_1960-1969_1.9x2.5_FD_ncra.nc ndep_ocn_1960-1969_1.9x2.5_FD_ncra_ncap2.nc

    At this juncture, for the case of 1970-1979, the values are still correct:

    lon(45)=110 lat(47)=-2.84210526316 time(1) Ndep(6669)=2.72749e-11 kg(N)/m2/s
    !!-> lon(46)=112.5 lat(47)=-2.84210526316 time(1) Ndep(6670)=1.19372e-10 kg(N)/m2/s
    !!-> lon(47)=115 lat(47)=-2.84210526316 time(1) Ndep(6671)=1.56785e-10 kg(N)/m2/s
    lon(48)=117.5 lat(47)=-2.84210526316 time(1) Ndep(6672)=7.39291e-12 kg(N)/m2/s

    3) The original files have a time dimension, but no time variable or associated values. One goal is to catenate all these decadal averages into a ~150-year file, so I therefore manually input a time variable and specify the appropriate value (units: "months since 01-01-1850")

    ncks --cdl -v Ndep ndep_ocn_1960-1969_1.9x2.5_FD_ncra_ncap2.nc > ndep_ocn_1960-1969_1.9x2.5_FD_ncra_ncap2.cdl

    3.1) Insert into the .cdl file under "variable:"

    double time(time) ;
    time:long_name = "time" ;
    time:units = "months since 01-01-1850" ;

    3.2) Insert into the .cdl file under "data:"

    time = 1380; // 1380 (months since 01-01-1850 for 1960-69; 1500 for 1970-79

    4) Now, regenerating the nc file from the edited cdl file:

    ncgen -k netCDF-4 -b -o ndep_ocn_1960-1969.nc ndep_ocn_1960-1969_1.9x2.5_FD_ncra_ncap2.cdl

    At this juncture, for the case of 1970-1979, the values have been corrupted:

    lon(45)=110 lat(47)=-2.84210526316 time(1)=1500 Ndep(6669)=2.72749e-11 kg(N)/m2/s
    !!-> lon(46)=112.5 lat(47)=-2.84210526316 time(1)=1500 Ndep(6670)=0.119372 kg(N)/m2/s
    !!-> lon(47)=115 lat(47)=-2.84210526316 time(1)=1500 Ndep(6671)=0.156785 kg(N)/m2/s
    lon(48)=117.5 lat(47)=-2.84210526316 time(1)=1500 Ndep(6672)=7.39291e-12 kg(N)/m2/s

    It appears that 1970-1970 is the first file to contain data points at e-10 (previously they range from e-11 to e-15), and that all files from this point on contain data at the e-10 order of magnitude. The sites where these values occur are highly localized around the lat47,lon46-47 region.

    The values themselves shouldn't be an issue though, and it seems that ncgen is somehow interpreting e-10 as e-1, replacing 1.19372e-10 with 0.119372 for example.

    The work around is to
    1) create a cdl of the dimensions only (ncks -c --cdl), which I can then manually edit for time variables and values
    2) create a nc file with ncgen (has no variables)
    3) append the variable from my intermediate file (*_ncra_ncap2.nc) to the new nc file.

    Still, one has to find the glitch before even knowing the workaround is necessary.

  • Charlie Zender

    Charlie Zender - 2014-04-09

    Thank you for reporting this. It's probably either a ncks --cdl or a ncgen bug. It would be most helpful if you could provide a small ncks-generated CDL file that ncgen appears to convert incorrectly. Just hyperslab the *#@! out of the
    file you have until all that remains are a few values that give the weird behavior. That's what we need to fix it. Thanks, cz

  • Stubaan

    Stubaan - 2014-04-11

    Perform the following steps to obtain a small version that manifests the bug discussed:

    1) Download file ndep_ocn_1970-1979_1.9x2.5_FD_ncra_ncap2.nc from https://sourceforge.net/projects/ipcc3hourly/files/NCO%20ncgen/

    2) ncks -O -F --cdl -v Ndep -d lat,47 -d lon,40,50,1 ndep_ocn_1970-1979_1.9x2.5_FD_ncra_ncap2.nc > ndep_ocn_1970-1979_1.9x2.5_FD_ncra_ncap2_bug.cdl

    3) ncgen -O -k netCDF-4 -b -o ndep_ocn_1970-1979_bug.nc ndep_ocn_1970-1979_1.9x2.5_FD_ncra_ncap2_bug.cdl

    If you look at the resulting values you will see that values with e-10 in step 2) are converted to decimals with equivalent e-1 by step 3).

  • Charlie Zender

    Charlie Zender - 2014-04-12

    Thank you for reporting this. Apparently this is an ncks --cdl bug, not an ncgen bug. The error occurs in step 2 not step 3. ncks hyperslabs correctly but for some reason ncks --cdl truncates the e-10 exponent (ncks without --cdl does fine) when printing. A subtle yet serious bug in the formatting of floating point CDL output. Will patch and post ASAP.

  • Charlie Zender

    Charlie Zender - 2014-04-12

    Thanks to your testcase I have found the bug, which is in the routine that trims trailing zeros from printed numbers. Will have patch within a few days...

  • Charlie Zender

    Charlie Zender - 2014-04-13

    The fix is now in the daily snapshot. From the draft 4.4.4 ANNOUNCE:
    A. Fix bug that caused CDL and XML printing to truncate trailing zeros
    of floating point exponents. Before this bugfix, ncks --cdl would
    print 1.100e-10 as 1.100e-1 instead of 1.1e-10. All CDL and XML
    printing prior to 4.4.4 mis-prints floating point numbers with
    exponents that end in zero (i.e., +/- 10, 20, 30 ... 200, 210...).
    Users of CDL and XML mode are urged to upgrade. This bug changed
    affected values by some multiple of 10 orders of magnitude.
    Hence it is likely that users affected by this bug would notice it
    when carefully performing any statistics, yet possible they would
    not for fields with large dynamic ranges. Thanks to Stubaan for report.

  • Stubaan

    Stubaan - 2014-04-14

    Excellent work! Thank you.


Log in to post a comment.