Fill in holes in data.

fixholes.py

#!/usr/bin/env python3

import numpy as np
import gdal
import sys

MV = 8500

gdal.UseExceptions()

fn = sys.argv[1]
fd = gdal.Open(fn, gdal.GA_Update)

data = fd.ReadAsArray()

rb = fd.GetRasterBand(1)
nodata = rb.GetNoDataValue()

mask = (data > MV).any(axis=0)
data[:, mask] = nodata

for i in range(fd.RasterCount):
    rb = fd.GetRasterBand(i + 1)
    rb.WriteArray(data[i, :, :])
    rb.SetNoDataValue(nodata)

    gdal.FillNodata(rb, None, 5, 0, [], None)

fd = None