Skip to content

Commit 273e335

Browse files
committed
RF: check fp error in offset for offset scaling
When trying to scale integers only with inter, we need to make sure the inter as cast to the scaling type does in fact keep the integers in the range of the output type. We also want to make sure that inter generates numbers with the smallest possible absolute magnitude to avoid floating point error in the resulting numbers.
1 parent dd2b7f7 commit 273e335

File tree

1 file changed

+19
-3
lines changed

1 file changed

+19
-3
lines changed

nibabel/arraywriters.py

Lines changed: 19 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -464,9 +464,25 @@ def _iu2iu(self):
464464
t_min, t_max = np.iinfo(out_dtype).min, np.iinfo(out_dtype).max
465465
type_range = as_int(t_max) - as_int(t_min)
466466
mn2mx = mx - mn
467-
if mn2mx <= type_range: # offset enough?
468-
self.inter = mn - t_min
469-
return
467+
if mn2mx <= type_range: # might offset be enough?
468+
if t_min == 0: # uint output - take min to 0
469+
# decrease offset with floor_exact, meaning mn >= t_min after
470+
# subtraction. But we may have pushed the data over t_max,
471+
# which we check below
472+
inter = floor_exact(mn - t_min, self.scaler_dtype)
473+
else: # int output - take midpoint to 0
474+
# ceil below increases inter, pushing scale up to 0.5 towards
475+
# -inf, because ints have abs min == abs max + 1
476+
midpoint = mn + as_int(np.ceil(mn2mx / 2.0))
477+
# Floor exact decreases inter, so pulling scaled values more
478+
# positive. This may make mx - inter > t_max
479+
inter = floor_exact(midpoint, self.scaler_dtype)
480+
# Need to check still in range after floor_exact-ing
481+
int_inter = as_int(inter)
482+
assert mn - int_inter >= t_min
483+
if mx - int_inter <= t_max:
484+
self.inter = inter
485+
return
470486
# Try slope options (sign flip) and then range scaling
471487
super(SlopeInterArrayWriter, self)._iu2iu()
472488

0 commit comments

Comments
 (0)