# Scalefactors and intercepts¶

SPM Analyze and nifti1 images have scalefactors. nifti1 images also have intercepts. If `A` is an array in memory, and `S` is the array that will be written to disk, then:

```R = (A - intercept) / scalefactor
```

and `R == S` if `R` is already the data dtype we need to write.

If we load the image from disk, we exactly recover `S` (and `R`). To get something approximating `A` (say `Aprime`) we apply the intercept and scalefactor:

```Aprime = (S * scalefactor) + intercept
```

In a perfect world `A` would be exactly the same as `Aprime`. However `scalefactor` and `intercept` are floating point values. With floating point, if `r = (a - b) / c; p = (r * c) + b` it is not necessarily true that `p == a`. For example:

```>>> import numpy as np
>>> a = 10
>>> b = np.e
>>> c = np.pi
>>> r = (a - b) / c
>>> p = (r * c) + b
>>> p == a
False
```

So there will be some error in this reconstruction, even when `R` is the same type as `S`.

More common is the situation where `R` is a different type from `S`. If `R` is of type `r_dtype`, `S` is of type `s_dtype` and `cast_function(R, dtype)` is some function that casts `R` to the desired type `dtype`, then:

```R = (A - intercept) / scalefactor
S = cast_function(R, s_dtype)
R_prime = cast_function(S, r_dtype)
A_prime = (R_prime * scalefactor) + intercept
```

The type of `R` will depend on what numpy did for upcasting ```A, intercept, scalefactor```.

In order that `cast_function(S, r_dtype)` can best reverse ```cast_function(R, s_dtype)```, the second needs to know the type of `R`, which is not stored. The type of `R` depends on the types of `A` and of `intercept, scalefactor`. We don’t know the type of `A` because it is not stored.

`R` is likely to be a floating point type because of the application of scalefactor and intercept. If `(intercept, scalefactor)` are not the identity (0, 1), then we can ensure that `R` is at minimum the type of the ```intercept, scalefactor``` by making these be at least 1D arrays, so that floating point types will upcast in `R = (A - intercept) / scalefactor`.

The cast of `R` to `S` and back to `R_prime` can lose resolution if the types of `R` and `S` have different resolution.

Our job is to select:

• scalefactor

• intercept

• `cast_function`

such that we minimize some measure of difference between `A` and `A_prime`.