Fixing NetCDF Values That Look Wrong: scale_factor, _FillValue and Units

On this page

Problem statement

The file opens, the map has the right shape, and the numbers are wrong: sea temperatures of 3000, a global mean of −infinity, rainfall that turns negative, a Kelvin field that other software reads as entirely missing. Almost every case comes from the same small set of attributes — scale_factor and add_offset, _FillValue and missing_value, valid_min and valid_max, units — being applied by one reader and not another, or being carried into a file where they no longer fit the data.

Measured on NOAA OISST v2.1 sea surface temperature, NOAA CPC precipitation and the NCEP/NCAR Reanalysis 1:

Reading OISST through rasterio returned the stored integers: a maximum of 3513 and a mean of 617.16, including −999 land cells. Decoded, the mean was 14.254 °C.
CPC precipitation read without masking had a mean of −infinity, because missing cells hold −9.97 × 10³⁶.
A 30-day accumulation written with OISST's inherited 16-bit encoding came back capped at 327.6; its true maximum was 1053.9, and 294,074 cells turned negative.
A Kelvin field written the same way read correctly in xarray but was entirely masked in netCDF4-python, which applied the old valid_min and valid_max of −300 and 4500.

Quick answer

Find the symptom, then check the attribute behind it:

symptom                                      likely cause                           fix
values 100x too large, integers              scale_factor not applied               read with xarray, or apply scales
huge negative or -inf means                  missing_value / _FillValue not masked  mask before arithmetic
values ~273 too high                         units are Kelvin                       check the units attribute
written values capped, wrapped or negative   inherited int16 packing                clear .encoding before writing
file reads as all missing elsewhere          stale valid_min / valid_max            remove or update the attributes
values rounded to whole numbers              least_significant_digit in encoding    clear .encoding before writing

Then compare stored and decoded values directly:

import xarray as xr

raw = xr.open_dataset("oisst-avhrr-v02r01.20240715.nc", mask_and_scale=False)["sst"]
decoded = xr.open_dataset("oisst-avhrr-v02r01.20240715.nc")["sst"]
print(raw.dtype, float(raw.max()), "->", decoded.dtype, float(decoded.max()), decoded.encoding.get("scale_factor"))

int16 3513.0 -> float32 35.130001068115234 0.01

Six symptoms, four attributes: most wrong values are a decoding step applied in one place and not another.

Compare stored and decoded, then check units

1. Compare stored and decoded values

Open the variable twice, with mask_and_scale=False and with the default, and compare minimum, maximum, mean and data type (Example 1). For one day of OISST the stored values were 16-bit integers from −999 to 3513; decoded, they were 32-bit floats up to 35.13 with land as NaN. If your numbers match the stored column, a decoding step is missing.

2. Check which reader applied the decoding

xarray and netCDF4-python apply the scale factor and mask fill values by default; netCDF4 returned a masked array with a mean of 14.254 and 345,650 masked cells. GDAL and rasterio do not: reading netcdf:file.nc:sst gave int16 values, nodata = -999, scales = (0.01,), and a raw mean of 617.16. The metadata was there; applying it is the caller's job (Example 3). rioxarray's open_rasterio on the same subdataset failed outright with a CoordinateValidationError about the zlev coordinate.

3. Mask `missing_value` as well as `_FillValue`

CPC precipitation marks missing cells with missing_value = -9.96921e+36. Read without masking, one day's mean was −infinity and its minimum −9.97 × 10³⁶; decoded, the mean was 0.995 mm. xarray masks both attributes. Code that masks only _FillValue, or only NaN, lets the sentinel into every sum.

4. Check the units before the values

The NCEP 6-hourly surface file stores temperature in degK: its 2024 mean was 278.92. The monthly file stores degC. Nothing converts between them. A difference of about 273 is almost always a units problem, and a factor of 86,400 is almost always a rate per second against a rate per day.

5. Clear the encoding when the meaning of the values changes

The in-place row is the dangerous one, because it looks like ordinary arithmetic.

xarray keeps a variable's on-disk encoding through some operations and drops it through others. On OISST, isel, sel, copy, load, in-place += and assignment to .values kept dtype=int16, scale_factor=0.01; arithmetic that creates a new array, where, fillna, astype, clip, mean and rolling dropped it. Anything written with the kept encoding is packed into 16 bits with a 0.01 step, which holds −327.68 to 327.67. A 30-day accumulation with values up to 1053.9 came back with a maximum of 327.6 and 294,074 negative cells. Set da.encoding = {} whenever the values are no longer the quantity the encoding was designed for (Example 2).

6. Update or remove stale valid ranges

OISST's valid_min and valid_max, −300 and 4500, describe the stored integers of a Celsius field. Adding 273.15 in place and writing kept both attributes and the packing. xarray read the result back with a mean of 287.40; netCDF4-python, which applies valid ranges, masked all 1,036,800 cells, and rasterio saw nothing but −999. When a variable's values change, fix its attributes as well as its encoding.

7. Watch for precision loss on re-save

Re-saving the NCEP monthly file with its default encoding produced a 4.9 MB file instead of 29.6 MB: the source encoding carried least_significant_digit = 0, and every value was rounded to a whole degree — at the worst cell 22.50 became 22.00, a change of 0.5. The 2024 global mean barely moved, from 14.8730 to 14.8731, which is why the damage goes unnoticed. The same re-saved file then refused a second write with ValueError: Variable 'air' has conflicting _FillValue (nan) and missing_value.

8. Check GeoTIFF exports for scale tags

rio.to_raster on a decoded OISST field wrote int16 values with a 0.01 scale tag, because the encoding travelled with the data. rioxarray's default read of that GeoTIFF returned integers up to 3513 with the scale in attrs; with mask_and_scale=True it returned floats up to 35.13. See converting NetCDF to GeoTIFF for writing physical values instead.

The packing that suited a daily temperature silently capped a 30-day accumulation.

Packing, fill values and units across four NOAA files

The data used on this page

Four NOAA files: one day of OISST v2.1 sea surface temperature (15 July 2024) from NCEI; the 2024 daily CPC global precipitation grid; and the NCEP/NCAR Reanalysis 1 surface air temperature, both the 2024 six-hourly file and the long-term monthly means, from NOAA PSL. PSL drops long transfers, so the setup resumes and retries.

import os
import subprocess

PSL = ["https://downloads.psl.noaa.gov/Datasets/", "https://psl.noaa.gov/thredds/fileServer/Datasets/"]
FILES = {
    "oisst-avhrr-v02r01.20240715.nc": ["https://www.ncei.noaa.gov/data/sea-surface-temperature-optimum-interpolation/"
                                       "v2.1/access/avhrr/202407/oisst-avhrr-v02r01.20240715.nc"],
    "cpc_precip.2024.nc": [p + "cpc_global_precip/precip.2024.nc" for p in PSL],
    "air.sig995.2024.nc": [p + "ncep.reanalysis/surface/air.sig995.2024.nc" for p in PSL],
    "air.mon.mean.nc": [p + "ncep.reanalysis/Monthlies/surface/air.mon.mean.nc" for p in PSL],
}
for name, urls in FILES.items():
    for url in urls * 5:
        if os.path.exists(name):
            break
        if subprocess.run(["curl", "-sfL", "-C", "-", "--retry", "3", "-o", name + ".part", url]).returncode == 0:
            os.replace(name + ".part", name)
    if not os.path.exists(name):
        raise RuntimeError(f"could not download {name}")
print({name: f"{os.path.getsize(name) / 1e6:.1f} MB" for name in FILES})

{'oisst-avhrr-v02r01.20240715.nc': '1.5 MB', 'cpc_precip.2024.nc': '66.1 MB', 'air.sig995.2024.nc': '22.2 MB', 'air.mon.mean.nc': '29.6 MB'}

Example 1 — stored against decoded values for any variable

import numpy as np
import xarray as xr


def diagnose(path, variable):
    raw = xr.open_dataset(path, mask_and_scale=False)[variable]
    dec = xr.open_dataset(path)[variable]
    enc = dec.encoding
    print(f"{path} {variable}: stored {raw.dtype}, decoded {dec.dtype}, units {dec.attrs.get('units')!r}")
    print(f"  scale_factor {enc.get('scale_factor')}, add_offset {enc.get('add_offset')}, "
          f"_FillValue {enc.get('_FillValue')}, missing_value {enc.get('missing_value')}")
    with np.errstate(over="ignore", invalid="ignore"):
        print(f"  stored:  min {float(raw.min()):.6g}, max {float(raw.max()):.6g}, mean {float(raw.mean()):.6g}")
    print(f"  decoded: min {float(dec.min()):.6g}, max {float(dec.max()):.6g}, mean {float(dec.mean()):.6g}, "
          f"missing {float(dec.isnull().mean()):.1%}")


diagnose("oisst-avhrr-v02r01.20240715.nc", "sst")
diagnose("cpc_precip.2024.nc", "precip")
diagnose("air.sig995.2024.nc", "air")

oisst-avhrr-v02r01.20240715.nc sst: stored int16, decoded float32, units 'Celsius'
  scale_factor 0.009999999776482582, add_offset 0.0, _FillValue -999, missing_value None
  stored:  min -999, max 3513, mean 617.162
  decoded: min -1.8, max 35.13, mean 14.2542, missing 33.3%
cpc_precip.2024.nc precip: stored float32, decoded float32, units 'mm'
  scale_factor None, add_offset None, _FillValue None, missing_value -9.969209968386869e+36
  stored:  min -9.96921e+36, max 795.932, mean -inf
  decoded: min 0, max 795.932, mean 1.42416, missing 64.1%
air.sig995.2024.nc air: stored float32, decoded float32, units 'degK'
  scale_factor None, add_offset None, _FillValue None, missing_value -9.969209968386869e+36
  stored:  min 190.6, max 324.3, mean 278.922
  decoded: min 190.6, max 324.3, mean 278.922, missing 0.0%

The CPC figures cover the whole year, 366 daily grids; the NCEP file declares a missing_value that no cell uses.

Example 2 — write a derived field without inherited packing

def check_roundtrip(da, path):
    da.to_netcdf(path)
    back = xr.open_dataset(path)[da.name]
    print(f"{path}: written as {da.encoding.get('dtype', da.dtype)}, read back max {float(back.max()):.1f}, "
          f"min {float(back.min()):.1f}, largest error {float(abs(back - da).max()):.1f}")


sst = xr.open_dataset("oisst-avhrr-v02r01.20240715.nc")["sst"]
accumulation = sst.copy(deep=True)
accumulation.values[:] = accumulation.values * 30      # new quantity, old int16 x 0.01 encoding
check_roundtrip(accumulation, "accumulation_inherited.nc")

accumulation.encoding = {}
accumulation.attrs = {"long_name": "30-day SST accumulation", "units": "degC day"}
check_roundtrip(accumulation, "accumulation_clean.nc")

accumulation_inherited.nc: written as int16, read back max 327.6, min -327.6, largest error 1310.7
accumulation_clean.nc: written as float32, read back max 1053.9, min -54.0, largest error 0.0

Replacing the attributes matters as much as clearing the encoding: the old valid_min and valid_max would mask most of the new values in readers that apply them.

Example 3 — apply the scale when reading through rasterio

import rasterio

with rasterio.open("netcdf:oisst-avhrr-v02r01.20240715.nc:sst") as src:
    band = src.read(1, masked=True)                       # masks cells equal to nodata (-999)
    values = band * src.scales[0] + src.offsets[0]
    print(f"stored max {int(band.max())}, nodata {src.nodata}, masked {int(band.mask.sum()):,}; "
          f"scaled max {float(values.max()):.2f}, scaled mean {float(values.mean()):.3f}")

stored max 3513, nodata -999.0, masked 345,650; scaled max 35.13, scaled mean 14.254

Example 4 — the other claims on this page, checked

import warnings

import netCDF4
import rioxarray

warnings.filterwarnings("ignore")
path = "oisst-avhrr-v02r01.20240715.nc"

with netCDF4.Dataset(path) as nc:                      # step 2
    v = nc["sst"][:]
    print(f"netCDF4: mean {float(v.mean()):.3f}, masked {int(v.mask.sum()):,}")
with rasterio.open(f"netcdf:{path}:sst") as src:
    print(f"rasterio: dtype {src.dtypes[0]}, nodata {src.nodata}, scales {src.scales}, "
          f"raw mean {float(src.read(1).mean()):.2f}")
try:
    rioxarray.open_rasterio(f"netcdf:{path}:sst")
    print("rioxarray: opened")
except Exception as exc:
    print(f"rioxarray: {type(exc).__name__}")

cpc = xr.open_dataset("cpc_precip.2024.nc", mask_and_scale=False)["precip"].isel(time=0)   # step 3
with np.errstate(over="ignore", invalid="ignore"):
    print(f"CPC one day raw: min {float(cpc.min()):.3g}, mean {float(cpc.mean())}; "
          f"decoded mean {float(xr.open_dataset('cpc_precip.2024.nc')['precip'].isel(time=0).mean()):.3f}")

sst = xr.open_dataset(path)["sst"]                      # step 5: which operations keep the encoding
ops = {"isel": lambda d: d.isel(lat=slice(0, 10)), "copy": lambda d: d.copy(), "load": lambda d: d.load(),
       "d + 1": lambda d: d + 1, "where": lambda d: d.where(d > 0), "fillna": lambda d: d.fillna(0),
       "astype": lambda d: d.astype("float64"), "mean": lambda d: d.mean("time")}
print("encoding kept: " + ", ".join(f"{k} {'yes' if f(sst).encoding.get('scale_factor') else 'no'}"
                                    for k, f in ops.items()))

kelvin = xr.open_dataset(path)["sst"]                   # step 6: valid ranges that no longer fit
kelvin.values[:] += 273.15
kelvin.to_netcdf("sst_kelvin.nc")
back = xr.open_dataset("sst_kelvin.nc")["sst"]
with netCDF4.Dataset("sst_kelvin.nc") as nc:
    masked = int(nc["sst"][:].mask.sum())
print(f"Kelvin file: xarray mean {float(back.mean()):.2f}; netCDF4 masked {masked:,} of {back.size:,} "
      f"(valid_min {back.encoding.get('valid_min', kelvin.attrs.get('valid_min'))}, "
      f"valid_max {back.encoding.get('valid_max', kelvin.attrs.get('valid_max'))})")

mon = xr.open_dataset("air.mon.mean.nc")                # step 7: re-saving with the source encoding
mon.to_netcdf("air_resaved.nc")
again = xr.open_dataset("air_resaved.nc")
print(f"re-save: {os.path.getsize('air.mon.mean.nc') / 1e6:.1f} MB -> {os.path.getsize('air_resaved.nc') / 1e6:.1f} MB; "
      f"least_significant_digit {mon.air.encoding.get('least_significant_digit')}; "
      f"largest change {float(abs(again.air - mon.air).max()):.2f}")
y24 = lambda d: float(d.air.sel(time="2024").weighted(np.cos(np.deg2rad(d.lat))).mean())   # area-weighted
print(f"2024 global mean {y24(mon):.4f} -> {y24(again):.4f}")
diff = abs(again.air - mon.air)
where = np.unravel_index(int(diff.argmax()), diff.shape)
print(f"for example {float(mon.air.values[where]):.2f} became {float(again.air.values[where]):.2f}")
inherited = xr.open_dataset("accumulation_inherited.nc")["sst"]
print(f"inherited-encoding accumulation: {int((inherited < 0).sum()):,} cells read back negative")
try:
    again.to_netcdf("air_resaved_twice.nc")
except ValueError as exc:
    print("second write:", str(exc).splitlines()[0])

sst.squeeze().rio.set_spatial_dims(x_dim="lon", y_dim="lat").rio.write_crs(4326).rio.to_raster("sst.tif")  # step 8
t = rioxarray.open_rasterio("sst.tif")
t2 = rioxarray.open_rasterio("sst.tif", mask_and_scale=True)
print(f"GeoTIFF: dtype {t.dtype}, max {float(t.max()):.0f}, scale in attrs {t.attrs.get('scale_factor')}; "
      f"mask_and_scale max {float(t2.max()):.2f}")

netCDF4: mean 14.254, masked 345,650
rasterio: dtype int16, nodata -999.0, scales (0.009999999776482582,), raw mean 617.16
rioxarray: CoordinateValidationError
CPC one day raw: min -9.97e+36, mean -inf; decoded mean 0.995
encoding kept: isel yes, copy yes, load yes, d + 1 no, where no, fillna no, astype no, mean no
Kelvin file: xarray mean 287.40; netCDF4 masked 1,036,800 of 1,036,800 (valid_min -300, valid_max 4500)
re-save: 29.6 MB -> 4.9 MB; least_significant_digit 0; largest change 0.50
2024 global mean 14.8730 -> 14.8731
for example 22.50 became 22.00
inherited-encoding accumulation: 294,074 cells read back negative
second write: Variable 'air' has conflicting _FillValue (nan) and missing_value (-9.969209968386869e+36). Cannot encode data.
GeoTIFF: dtype int16, max 3513, scale in attrs 0.009999999776482582; mask_and_scale max 35.13

Why tools disagree about the same file

Why the same file gives different numbers in different tools

The NetCDF format stores integers and attributes; the rule that scale_factor multiplies the integers is a convention, not part of the storage. xarray, netCDF4-python and CDO apply it. GDAL exposes it as band metadata and leaves applying it to you. Any tool that reads the raw array sees stored values, and a correct file can therefore produce wrong numbers.

Why packing is dangerous on output

Packing is designed for one quantity with a known range: OISST's int16 with a 0.01 step covers about ±327 °C, generous for sea temperature. The encoding does not know that the variable now holds an accumulation, an anomaly in hundredths or a Kelvin value. Out-of-range values are written as whatever the integer conversion produces, and the file gives no warning on reading.

Why valid ranges bite later

valid_min and valid_max are ignored by xarray on reading, so a file with stale limits looks fine in xarray. Readers that follow the CF recommendation to mask values outside the valid range — netCDF4-python by default — throw away everything outside the old limits. The error appears for whoever uses the file next.

Why rounding hides in averages

Rounding every cell to a whole degree adds errors of up to half a degree that are nearly symmetric, so area means and trends change in the fourth decimal place. Maps, extremes, gradients and anything computed from differences between neighbours change far more.

One file, four readers, three different sets of numbers.

Offsets, unsigned packing and odd fill values

add_offset is usually 0 in these files; when it is not, forgetting it shifts every value by a constant.
Unsigned packing uses _Unsigned = "true" on signed integer variables; readers that ignore it get negative numbers.
Several fill values can coexist: _FillValue, missing_value and values outside valid_range may all mean "no data".
Float fill values such as 1e20 survive into means as huge numbers when masking is switched off.
Accumulations in GRIB and model output may be totals since the forecast start rather than per step.
Units strings vary: degC, Celsius, degree_Celsius and K all appear; normalise before comparing.
Encoding from a source file is reused by to_netcdf unless replaced; pass encoding= or clear it explicitly.

Internal links

CF conventions explained: how a NetCDF file says what its numbers mean — the attributes involved
NetCDF and gridded data explained: dimensions, variables and attributes — stored versus decoded values
How to open a NetCDF file in Python with xarray — decoding switches
How to convert NetCDF to GeoTIFF in Python — packing in exports
Rasterio returns the wrong values: NoData, scaling and dtype fixes — the GDAL side
Fixing a NetCDF map that is shifted, flipped or split at 180° — when the values are right but misplaced
How to read GRIB weather forecast data in Python — units and accumulations in GRIB
The raster data model explained: bands, dtype, NoData and the transform — scale and offset in rasters

FAQ

Why are my NetCDF values 100 times too large?

The values are packed integers and the reader did not apply scale_factor. OISST stores 0.01 °C steps, so a stored 3513 means 35.13 °C; xarray applies the factor, rasterio does not.

Why is the mean of my NetCDF variable negative infinity?

Missing cells hold a sentinel such as −9.97 × 10³⁶ that was not masked. Read with xarray's default decoding, or mask missing_value and _FillValue before any arithmetic.

Why did my written NetCDF values get capped or turn negative?

The variable kept its source file's int16 packing. A 30-day accumulation written that way was capped at 327.6; set .encoding = {} before writing derived quantities.

Why does another program see my NetCDF file as all missing?

Stale valid_min and valid_max attributes. A Kelvin field written with Celsius-era limits was fully masked by netCDF4-python while xarray read it correctly.

Why were my values rounded after saving with xarray?

The source encoding included least_significant_digit. Re-saving the NCEP file rounded every value to a whole degree; clear the encoding or pass a new one.

Does GDAL apply scale_factor when reading NetCDF?

No. It reports the scale and offset as band metadata; multiply by src.scales and add src.offsets yourself, or read with mask_and_scale=True in rioxarray.