karlmsmith
commented
6 years ago Comment by @AnsleyManke on 24 May 2007 22:54 UTC This has been fixed; see the benchmark script
err601_irregular_axis.jnl
Closed karlmsmith closed 6 years ago
karlmsmith
commented
6 years ago Comment by @AnsleyManke on 24 May 2007 22:54 UTC This has been fixed; see the benchmark script
err601_irregular_axis.jnl
karlmsmith
commented
6 years ago Comment by @AnsleyManke on 26 Jun 2007 22:34 UTC Further comments, and a further fix -
The code in cd_get_1_axis.F implements TM_FPEQ_EPS as described in previous comments, butI am now changing it so we DO NOT recompute epsilon for each delta value. Instead compute it once, for the largest coordinate.
epsilon = epsilon_23 * 2.*(ABS(line_mem(ipte)) / delta )
DO 360 i=ipt1+2,ipte
delta = line_mem(i) - line_mem(i-1)
IF (.NOT. TM_FPEQ_EPS(epsilon, first_delta, delta) ) GOTO 380360 CONTINUE
In addition, use TM_FPEQ to check whether a longitude axis is 360 degrees or less long. Previously, we had a simple IF (axwwlen .LE. 360.D0) where axwwlen is the result of a call to double precision function TM_WW_AXLEN(iaxis) which computes the entire axis span to the cell edges. Allowing for some computational fuzzyness means that the axis, for instance, in "http://dods.jpl.nasa.gov/dods-bin/nph-hdf/pub/sea_surface_temperature/GHRSST/data/L4/GLOB/UKMO/OSTIA/2007/170/20070619-UKMO-L4HRfnd-GLOB-v01-fv02-OSTIA.nc.bz2" is detected as a modulo longitude axis.
The above fix also means that the axes of this GHRSST data are found to be regularly-spaced.
karlmsmith
commented
6 years ago Comment by @AnsleyManke on 27 Aug 2008 17:27 UTC Further problems with detecting axes as regular or not. See the user-list thread "Error in plotting time series data with 6.1?" http://www.pmel.noaa.gov/maillists/tmap/ferret_users/fu_2008/msg00440.html
The code is still wrong. If the coordinates themselves are large, the relative-error computation detects even large delta-coordinate variations as regularly spaced. Here are coordinates from an axis with t0=1-jan-1970
17-APR-2007 11:22:51 / 68912: 1176808960.
17-APR-2007 11:25:04 / 68913: 1176809088.
17-APR-2007 11:27:17 / 68914: 1176809216.
27-APR-2007 23:39:56 / 68915: 1177717248.
27-APR-2007 23:42:07 / 68916: 1177717376.
27-APR-2007 23:44:08 / 68917: 1177717504.
27-APR-2007 23:46:20 / 68918: 1177717632.
27-APR-2007 23:48:36 / 68919: 1177717760.first-delta is 128. delta between 17-APR-2007 11:27:17 and 27-APR-2007 23:39:56 is 6809216. eps is 2.2
The routine tm_fpeq_eps had this test, where a = first_delta and b is the delta we're checking: TM_FPEQ_EPS = ABS(a-b) .LE. ABS(b)*eps
Well, if b is much larger than a, this will be computed as TRUE. Instead shouldn't it be
TM_FPEQ_EPS = (ABS(a-b) .LE. ABS(b)*eps) .AND.
. (ABS(a-b) .LE. ABS(a)*eps)In addition, if the file has the attribute point_spacing = "irregular" we should take it at its word and treat the axis as irregular.
karlmsmith
commented
6 years ago Comment by @AnsleyManke on 18 Mar 2009 21:49 UTC Make a different test in double-precision arithmetic if the coordinates come in as double precision, doing the tests in double precision. If the coordinates come in as single precision, test in single precision.
The netcdf file examples above in the report are also correctly detected as irregular.
Also tested this by defining an axis with origin at 1-jan-0001, in seconds
yes? define axis/t=1-jan-2008:"1-jan-2008:00:02:20":1/units=seconds/t0="1-jan-0001" tax
yes? let cc = L[gt=tax]
yes? save/file=cc.nc/clobber ccUse ncdump to make a cdl file and edit the cdl file, removing one coordinate value in the time axis and one value of the variable cc. Run ncgen to create a new netcdf file, ccirreg.nc. Previous Ferret executables incorrectly detected the time axis as regular, with the fix it's detected as irregular.
karlmsmith
commented
6 years ago Attachment from @AnsleyManke on 23 Jan 2007 20:16 UTC example netCDF file to demonstrate the behavior z_wide_variation.nc.zip
Reported by @AnsleyManke on 23 Jan 2007 20:15 UTC This file, a subset of one sent us by Ehouarn Millour, "Ferret bug report (about vertical coordinate)" on 22-Jan-2007, has this ncdump output for its z axis
----- output of ncdump -v altitude temp_pb.nc ----------------- altitude = 1400, 1100, 800, 565, 400, 280, 200, 140, 100, 70, 50, 35, 25, 18, 13, 9, 6, 3.632, 1.804, 0.8957, 0.4448, 0.2209, 0.1097, 0.05447, 0.02705, 0.01343, 0.00667, 0.003312, 0.001645, 0.0008168, 0.0004056, 0.0002014, 0.0001, 4.967e-05, 2.466e-05, 1.225e-05, 6.082e-06, 3.02e-06, 1.5e-06, 7.448e-07, 3.699e-07, 1.837e-07, 9.12e-08, 4.529e-08, 2.249e-08, 1.117e-08, 5.546e-09, 2.754e-09, 1.368e-09, 6.792e-10 ;
yes? use z_wide_variation.nc NOTE: Units on axis "ALTITUDE" are not recognized: Pa NOTE: They will not be convertible: yes? show grid/z atmos GRID GOR1 name axis # pts start end LONGITUDE47_47 LONGITUDE 1mr 78.75E 78.75E LATITUDE14_14 LATITUDE 1 r 41.249S 41.249S ALTITUDE Z (Pa) 50 r 6.7916E-10 1400 TIME TIME 1mr 01-JAN 04:00 01-JAN 04:00
In Ferret v6.0 the axis is detected as being actually regular! and a regular axis substituted for it.
In an attempt to better capture what "machine precision" in this version, we tried a new test for regular axis:
OLD TEST
Check if axis is actually regularly spaced (though not flagged as "even")
. SNGL(line_mem(i)-line_mem(i-1))) ) GOTO 380 360 CONTINUE (save axis as regular) 380 (irregular axis processing)
NEW TEST ( see the thread Re: [Fwd: "GHRSST-compliant" netCDF?] in my mail directory fer_netcdf )
Check if axis is actually regularly spaced (though not flagged as "even") epsilon_23 = 2.**(-23)
. SNGL(line_mem(i)-line_mem(i-1))) ) GOTO 380 360 CONTINUE
What happens with the example axis above, is that the epsilon we compute is huge, because the max coord divided by the first delta is 1400 / (1.368e-09 - 6.792e-10). All of this computed after the coordinates are reversed. Then of course TM_FPEQ_EPS says, the difference between the current delta and the first delta is smaller than epsilon*delta, so the axis must be equally-spaced....
Safer though more computing if we compute the relative epsilon for each coordinate
360 CONTINUE
Migrated-From: http://dunkel.pmel.noaa.gov/trac/ferret/ticket/1483