Python netCDF4.num2date() Examples

def set_time(self, time):
        '''
        Converts time data into a pandas date object.

        Parameters
        ----------
        time: netcdf
            Contains time information.

        Returns
        -------
        pandas.DatetimeIndex
        '''
        times = num2date(time[:].squeeze(), time.units,
                         only_use_cftime_datetimes=False,
                         only_use_python_datetimes=True)
        self.time = pd.DatetimeIndex(pd.Series(times), tz=self.location.tz) 

def num2date(times, units, calendar=None):
    """Convert an array of integers into datetime objects.

    This function optimizes the num2date function of python-netCDF4 if the
    standard calendar is used.

    Args:
        times: An array of integers representing timestamps.
        units: A string with the format "{unit} since {epoch}",
            e.g. "seconds since 1970-01-01T00:00:00".
        calendar: (optional) Standard is gregorian. If others are used,
            netCDF4.num2date will be called.

    Returns:
        Either an array of numpy.datetime64 objects (if standard gregorian
        calendar is used), otherwise an array of python datetime objects.
    """
    try:
        unit, epoch = units.split(" since ")
    except ValueError:
        raise InvalidUnitString("Could not convert to datetimes!")

    if calendar is None:
        calendar = "gregorian"
    else:
        calendar = calendar.lower()

    if calendar != "gregorian":
        return netCDF4.num2date(times, units, calendar).astype(
            "M8[%s]" % unit_mapper[unit])

    # Numpy uses the epoch 1970-01-01 natively.
    converted_data = times.astype("M8[%s]" % unit_mapper[unit])

    # numpy.datetime64 cannot read certain time formats while pandas can.
    epoch = pd.Timestamp(epoch).to_datetime64()

    # Maybe there is another epoch used?
    if epoch != np.datetime64("1970-01-01"):
        converted_data -= np.datetime64("1970-01-01") - epoch
    return converted_data 

def _times(dataset, variable):
        """Find times related to variable in dataset"""
        time_dimension = variable.dimensions[0]
        coordinates = variable.coordinates.split()
        for c in coordinates:
            if c.startswith("time"):
                try:
                    var = dataset.variables[c]
                    return netCDF4.num2date(var[:], units=var.units)
                except KeyError:
                    pass
        for v, var in dataset.variables.items():
            if len(var.dimensions) != 1:
                continue
            if v.startswith("time"):
                d = var.dimensions[0]
                if d == time_dimension:
                    return netCDF4.num2date(var[:], units=var.units) 

def datespan_file(self, file):
        """
        return the just the day that this argo file covers
        """
        nc = seapy.netcdf(file)
        try:
            d = netCDF4.num2date(nc.variables['JULD'][0],
                                 nc.variables['JULD'].units)
            st = datetime.datetime(*d.timetuple()[:3])
            en = datetime.datetime(*d.timetuple()[:3] + (23, 59, 59))
        except:
            st = en = None
            pass
        finally:
            nc.close()
            return st, en 

def _valid_times(dataset, variable):
        """Search dataset for time axis"""
        var = dataset.variables[variable]
        for d in var.dimensions:
            if d.startswith('time'):
                if d in dataset.variables:
                    tvar = dataset.variables[d]
                    return np.array(
                        netCDF4.num2date(tvar[:], units=tvar.units),
                        dtype='datetime64[s]')
        coords = var.coordinates.split()
        for c in coords:
            if c.startswith('time'):
                tvar = dataset.variables[c]
                return np.array(
                    netCDF4.num2date(tvar[:], units=tvar.units),
                    dtype='datetime64[s]') 

def _map_ini_times_to_paths(self, paths):
        """
        .. note:: Potentially expensive I/O operation
        """
        mapping = defaultdict(list)
        for path in paths:
            initial_time = _initial_time(path)
            try:
                with netCDF4.Dataset(path) as dataset:
                    if initial_time is None:
                        var = dataset.variables["forecast_reference_time"]
                        initial_time = netCDF4.num2date(var[:], units=var.units)
            except (FileNotFoundError, KeyError) as ex:
                pass
            mapping[self.key(initial_time)].append(path)
        return mapping 

def date2num(dates, units, calendar=None):
    """Convert an array of integer into datetime objects.

    This function optimizes the date2num function of python-netCDF4 if the
    standard calendar is used.

    Args:
        dates: Either an array of numpy.datetime64 objects (if standard
            gregorian calendar is used), otherwise an array of python
            datetime objects.
        units: A string with the format "{unit} since {epoch}",
            e.g. "seconds since 1970-01-01T00:00:00".
        calendar: (optional) Standard is gregorian. If others are used,
            netCDF4.num2date will be called.

    Returns:
        An array of integers.
    """
    if calendar is None:
        calendar = "gregorian"
    else:
        calendar = calendar.lower()

    if calendar != "gregorian":
        return netCDF4.date2num(dates, units, calendar)

    try:
        unit, epoch = units.split(" since ")
    except ValueError:
        raise InvalidUnitString("Could not convert to numeric values!")

    converted_data = \
        dates.astype("M8[%s]" % unit_mapper[unit]).astype("int")

    # numpy.datetime64 cannot read certain time formats while pandas can.
    epoch = pd.Timestamp(epoch).to_datetime64()

    if epoch != np.datetime64("1970-01-01"):
        converted_data -= np.datetime64("1970-01-01") - epoch
    return converted_data 

def findLastYearInNCTime(ncTimeVariable):

    # last datetime
    last_datetime = nc.num2date(ncTimeVariable[len(ncTimeVariable) - 1],\
                                ncTimeVariable.units,\
                                ncTimeVariable.calendar) 
    
    return last_datetime.year 

def findFirstYearInNCTime(ncTimeVariable):

    # first datetime
    first_datetime = nc.num2date(ncTimeVariable[0],\
                                ncTimeVariable.units,\
                                ncTimeVariable.calendar) 
    
    return first_datetime.year 

def calendar(date_in, calendar_type='proleptic_gregorian'):
    """ Get date or number of steps from input.

    Get date from input string using one of the available formats or get time step number from input number or string.
    Used to get the date from CalendarDayStart (input) in the settings xml

    :param date_in: string containing a date in one of the available formats or time step number as number or string
    :param calendar_type:
    :rtype: datetime object or float number
    :returns: date as datetime or time step number as float
    :raises ValueError: stop if input is not a step number AND it is in wrong date format
    """
    try:
        # try reading step number from number or string
        return float(date_in)
    except ValueError:
        # try reading a date in one of available formats
        try:
            _t_units = "hours since 1970-01-01 00:00:00"  # units used for date type conversion (datetime.datetime -> calendar-specific if needed)
            date = parse_time_string(date_in, dayfirst=True)[0]  # datetime.datetime type
            step = date2num(date, _t_units, calendar_type)  # float type
            return num2date(step, _t_units, calendar_type)  # calendar-dependent type from netCDF4.netcdftime._netcdftime module
        except:
            # if cannot read input then stop
            msg = "Wrong step or date format in XML settings file\n Input {}".format(date_in)
            raise LisfloodError(msg) 

def convert_file(self, file, title="OSTIA SST Obs"):
        """
        Load an OSTIA file and convert into an obs structure
        """
        # Load OSTIA Data
        nc = seapy.netcdf(file)
        lon = nc.variables["lon"][:]
        lat = nc.variables["lat"][:]
        dat = np.ma.masked_outside(np.squeeze(
            nc.variables["analysed_sst"][:]) - 273.15,
            self.temp_limits[0], self.temp_limits[1])
        err = np.ma.masked_outside(np.squeeze(
            nc.variables["analysis_error"][:]), 0.01, 2.0)
        dat[err.mask] = np.ma.masked
        time = seapy.roms.num2date(
            nc, "time", records=[0], epoch=self.epoch)[0]
        nc.close()
        if self.grid.east():
            lon[lon < 0] += 360
        lon, lat = np.meshgrid(lon, lat)
        good = dat.nonzero()
        lat = lat[good]
        lon = lon[good]
        data = [seapy.roms.obs.raw_data("TEMP", "SST_OSTIA", dat.compressed(),
                                        err[good], self.temp_error)]
        # Grid it
        return seapy.roms.obs.gridder(self.grid, time, lon, lat, None,
                                      data, self.dt, title) 

def convert_file(self, file, title="NAVO SST Obs"):
        """
        Load a NAVO map file and convert into an obs structure
        """
        import re
        import sys

        nc = seapy.netcdf(file)
        lon = nc.variables["lon"][:]
        lat = nc.variables["lat"][:]
        dat = np.ma.masked_outside(np.squeeze(nc.variables["analysed_sst"][:]) - 273.15,
                                   self.temp_limits[0], self.temp_limits[1])
        err = np.ma.array(np.squeeze(
            nc.variables["analysis_error"][:]), mask=dat.mask)

        # this is an analyzed product and provides errors as a function
        # of space and time directly the temperature is the bulk
        # temperature (ie at around 4m depth, below the e-folding depths of
        # sunlight in the ocean so the product does not have a diuranl cycle
        # (ie you don;t have to worry about hourly variations)
        time = seapy.roms.num2date(
            nc, "time", records=[0], epoch=self.epoch)[0]
        nc.close()

        # here we set the depth to be 4 m below the surface
        if self.grid.east():
            lon[lon < 0] += 360
        lon, lat = np.meshgrid(lon, lat)
        good = dat.nonzero()
        lat = lat[good]
        lon = lon[good]
        data = [seapy.roms.obs.raw_data("TEMP", self.provenance, dat.compressed(),
                                        err[good], self.temp_error)]
        # Grid it
        obs = seapy.roms.obs.gridder(self.grid, time, lon, lat, None,
                                     data, self.dt, depth_adjust=True, title=title)
        obs.z *= 0
        obs.depth = -self.depth * np.ones(len(obs.depth))
        return obs 

def convert_file(self, file, title="REMSS SST Obs"):
        """
        Load an REMSS file and convert into an obs structure
        """
        # Load REMSS Data
        nc = seapy.netcdf(file)
        lon = nc.variables["lon"][:]
        lat = nc.variables["lat"][:]
        dat = np.ma.masked_outside(np.squeeze(
            nc.variables["sea_surface_temperature"][:]) - 273.15,
            self.temp_limits[0], self.temp_limits[1])
        err = np.ma.masked_outside(np.squeeze(
            nc.variables["sses_standard_deviation"][:]), 0.01, 2.0)
        dat[err.mask] = np.ma.masked

        # Check the data flags
        if self.check_qc_flags:
            flags = np.ma.masked_not_equal(
                np.squeeze(nc.variables["quality_level"][:]), 5)
            dat[flags.mask] = np.ma.masked
        else:
            dat = np.ma.masked_where(
                np.squeeze(nc.variables["quality_level"][:]).data == 1, dat)

        # Grab the observation time
        time = seapy.roms.num2date(nc, "time", records=[0])[0] - self.epoch
        dtime = nc.variables["sst_dtime"][:]
        time = np.squeeze((time.total_seconds() + dtime) * seapy.secs2day)
        nc.close()
        if self.grid.east():
            lon[lon < 0] += 360
        good = dat.nonzero()
        data = [seapy.roms.obs.raw_data("TEMP", self.provenance,
                                        dat.compressed(),
                                        err[good], self.temp_error)]
        # Grid it
        return seapy.roms.obs.gridder(self.grid, time[good], lon[good], lat[good],
                                      None, data, self.dt, title) 

def get_reftime(nc, epoch=default_epoch):
    """
    Given a ROMS netCDF4 file, return the reference time for the file. This
    is the timebase of the record dimension in the format:
    "<units> since <reftime>"

    Parameters
    ----------
    nc : netCDF4 dataset
        Input ROMS file
    epoch_str : string, optional
        If lacking units, use this string as the units

    Returns
    -------
    timebase : datetime
        datetime of the origin for the file
    time : string
        name of variable used to generate the base (None if default)
    """
    try:
        tvar = get_timevar(nc)
        calendar, _ = _get_calendar(nc.variables[tvar])

        return netCDF4.num2date(0, nc.variables[tvar].units,
                                calendar=calendar), tvar
    except AttributeError:
        return epoch, None 

def opendap(fetch):
    """Decorator for fetching data from Opendap servers."""
    @wraps(fetch)
    def wrapper(*args, **kwargs):
        url, varname, bbox, dt = fetch(*args, **kwargs)
        ds = open_url(url)
        for var in ds.keys():
            if var.lower().startswith("lon") or var.lower() == "x":
                lonvar = var
            if var.lower().startswith("lat") or var.lower() == "y":
                latvar = var
            if var.lower().startswith("time") or var.lower() == "t":
                timevar = var
        lat = ds[latvar][:].data
        lon = ds[lonvar][:].data
        lon[lon > 180] -= 360
        res = abs(lat[0]-lat[1])  # assume rectangular grid
        i1, i2, j1, j2 = datasets.spatialSubset(np.sort(lat)[::-1], np.sort(lon), res, bbox)
        t = ds[timevar]
        tt = netcdf4.num2date(t[:].data, units=t.units)
        ti = [tj for tj in range(len(tt)) if resetDatetime(tt[tj]) >= dt[0] and resetDatetime(tt[tj]) <= dt[1]]
        if len(ti) > 0:
            lati = np.argsort(lat)[::-1][i1:i2]
            loni = np.argsort(lon)[j1:j2]
            if len(ds[varname].data[0].shape) > 3:
                data = ds[varname].data[0][ti[0]:ti[-1]+1, 0, lati[0]:lati[-1]+1, loni[0]:loni[-1]+1]
            else:
                data = ds[varname].data[0][ti[0]:ti[-1]+1, 0, lati[0]:lati[-1]+1, loni[0]:loni[-1]+1]
            dt = tt[ti]
        else:
            data = None
            dt = None
        lat = np.sort(lat)[::-1][i1:i2]
        lon = np.sort(lon)[j1:j2]
        return data, lat, lon, dt
    return wrapper 

def netcdf(fetch):
    """Decorator for fetching NetCDF files (local or from Opendap servers)."""
    @wraps(fetch)
    def wrapper(*args, **kwargs):
        url, varname, bbox, dt = fetch(*args, **kwargs)
        ds = netcdf4.Dataset(url)
        for var in ds.variables:
            if var.lower().startswith("lon") or var.lower() == "x":
                lonvar = var
            if var.lower().startswith("lat") or var.lower() == "y":
                latvar = var
            if var.lower().startswith("time") or var.lower() == "t":
                timevar = var
        lat = ds.variables[latvar][:]
        lon = ds.variables[lonvar][:]
        lon[lon > 180] -= 360
        res = abs(lat[0]-lat[1])  # assume rectangular grid
        i1, i2, j1, j2 = datasets.spatialSubset(np.sort(lat)[::-1], np.sort(lon), res, bbox)
        t = ds.variables[timevar]
        tt = netcdf4.num2date(t[:], units=t.units)
        ti = [tj for tj in range(len(tt)) if resetDatetime(tt[tj]) >= dt[0] and resetDatetime(tt[tj]) <= dt[1]]
        if len(ti) > 0:
            lati = np.argsort(lat)[::-1][i1:i2]
            loni = np.argsort(lon)[j1:j2]
            if len(ds.variables[varname].shape) > 3:
                data = ds.variables[varname][ti, 0, lati, loni]
            else:
                data = ds.variables[varname][ti, lati, loni]
            dt = tt[ti]
        else:
            data = None
            dt = None
        lat = np.sort(lat)[::-1][i1:i2]
        lon = np.sort(lon)[j1:j2]
        return data, lat, lon, dt
    return wrapper 

def _downloadVariable(varname, dbname, dt, bbox):
    """Download specific variable from the MERRA Reanalysis dataset."""
    # FIXME: Grid is not rectangular, but 0.5 x 0.625 degrees
    log = logging.getLogger(__name__)
    res = 0.5
    try:
        url = "http://goldsmr4.sci.gsfc.nasa.gov:80/dods/M2T1NXSLV"
        ds = netcdf.Dataset(url)
        lat = ds.variables["lat"][:]
        lon = ds.variables["lon"][:]
        lon[lon > 180] -= 360.0
        i1, i2, j1, j2 = datasets.spatialSubset(np.sort(lat)[::-1], np.sort(lon), res, bbox)
        data = np.zeros((i2-i1, j2-j1))
        lati = np.argsort(lat)[::-1][i1:i2]
        loni = np.argsort(lon)[j1:j2]
        t = ds.variables["time"]
        tt = netcdf.num2date(t[:], units=t.units)
        ti = np.where(tt == dt)[0][0]
        if varname == "tmax":
            hdata = ds.variables["t2m"][ti:ti+24, lati, loni]
            data = np.amax(hdata, axis=0) - 273.15
        elif varname == "tmin":
            hdata = ds.variables["t2m"][ti:ti+24, lati, loni]
            data = np.amin(hdata, axis=0) - 273.15
        elif varname in ["wind"]:
            hdata = np.sqrt(ds.variables["u10m"][ti:ti+24, lati, loni]**2 + ds.variables["v10m"][ti:ti+24, lati, loni]**2)
            data = np.mean(hdata, axis=0)
        lat = np.sort(lat)[::-1][i1:i2]
        lon = np.sort(lon)[j1:j2]
        filename = dbio.writeGeotif(lat, lon, res, data)
        table = "{0}.merra".format(varname)
        dbio.ingest(dbname, filename, dt, table)
        log.info("Imported {0} in {1}".format(tt[ti].strftime("%Y-%m-%d"), table))
        os.remove(filename)
    except:
        log.warning("Cannot import MERRA dataset for {0}!".format(dt.strftime("%Y-%m-%d"))) 

def findLastYearInNCTime(ncTimeVariable):

    # last datetime
    last_datetime = nc.num2date(ncTimeVariable[len(ncTimeVariable) - 1],\
                                ncTimeVariable.units,\
                                ncTimeVariable.calendar) 
    
    return last_datetime.year 

def __read_var(self, file, var):
        ds = Dataset(file, 'r')
        self.nx = len(ds.dimensions[self.config.xdim])
        self.ny = len(ds.dimensions[self.config.ydim])
        self.nt = len(ds.dimensions[self.config.tdim])

        self.x = ds.variables[self.config.xname][:]
        self.y = ds.variables[self.config.yname][:]

        # Sort out the dimensions.
        if self.config.clip:
            alldims = {}
            for key, val in list(ds.dimensions.items()):
                alldims[key] = (0, len(val))
            vardims = ds.variables[var].dimensions

            for clipname in self.config.clip:
                clipdims = self.config.clip[clipname]
                common = set(alldims.keys()).intersection([clipname])
                for k in common:
                    alldims[k] = clipdims
            dims = [alldims[d] for d in vardims]


        self.data = np.flipud(np.squeeze(ds.variables[var][
                dims[0][0]:dims[0][1],
                dims[1][0]:dims[1][1],
                dims[2][0]:dims[2][1],
                dims[3][0]:dims[3][1]
                ]))

        self.time = ds.variables[self.config.tname][:]
        self.Times = []
        for t in self.time:
            self.Times.append(num2date(
                t,
                'seconds since {}'.format(self.config.basedate),
                calendar=self.config.calendar
                ))

        ds.close() 

def initial_time_netcdf4(self, path):
        with netCDF4.Dataset(path) as dataset:
            try:
                var = dataset.variables["forecast_reference_time"]
                result = netCDF4.num2date(var[:], units=var.units)
            except KeyError:
                result = None
        return result 

def netcdf4_strategy(path):
        with netCDF4.Dataset(path) as dataset:
            var = dataset.variables["forecast_reference_time" ]
            values = netCDF4.num2date(var[:], units=var.units)
        return values 

def read_times(path):
    """Read time axis from a file"""
    with netCDF4.Dataset(path) as dataset:
        var = dataset.variables["time"]
        times = netCDF4.num2date(var[:], units=var.units)
    return np.array([forest.util.to_datetime(t) for t in times], dtype=object) 

def _load_cube(self, path, variable, lon0, lat0, pressure=None):
        """ Constrain data loading to points required """
        cube = iris.load_cube(path, variable)
        # reference longitude axis by "axis='X'" and latitude axis as axis='Y',
        # to accommodate various types of coordinate system.
        # e.g. 'grid_longitude'. See iris.utils.guess_coord_axis.
        if cube.coord(axis='X').points[-1] > 180.0:
            # get circular longitude values
            lon0 = iris.analysis.cartography.wrap_lons(np.asarray(lon0), 0, 360)
        # Construct constraint
        coord_values={cube.coord(axis='X').standard_name: lon0,
                      cube.coord(axis='Y').standard_name: lat0,                     
                      }
        if pressure is not None and 'pressure' in [coord.name() for coord in cube.coords()]:
            ptol = 0.01 * pressure
            coord_values['pressure'] = (
                lambda cell: (pressure - ptol) < cell < (pressure + ptol)
            )
        cube = cube.extract(iris.Constraint(coord_values=coord_values))
        assert cube is not None
        # Get validity times and data values
        # list the validity times as datetime objects
        time_coord = cube.coord('time')
        times = time_coord.units.num2date(time_coord.points).tolist()
        values = cube.data
        return {
            "x": times,
            "y": values} 

def __init__(self, paths):
        self.paths = paths
        self.table = defaultdict(list)
        for path in paths:
            time = _initial_time(path)
            if time is None:
                try:
                    with netCDF4.Dataset(path) as dataset:
                        var = dataset.variables["forecast_reference_time"]
                        time = netCDF4.num2date(var[:], units=var.units)
                except KeyError:
                    continue
            self.table[self.key(time)].append(path) 

def test_locator_times(tmpdir):
    path = str(tmpdir / "eida50_20190417.nc")
    with netCDF4.Dataset(path, "w") as dataset:
        _eida50(dataset, TIMES, LONS, LATS)
    result = eida50.Locator.load_time_axis(path)
    with netCDF4.Dataset(path) as dataset:
        var = dataset.variables["time"]
        expect = netCDF4.num2date(var[:], units=var.units)
    np.testing.assert_array_equal(expect, result) 

def findFirstYearInNCTime(ncTimeVariable):

    # first datetime
    first_datetime = nc.num2date(ncTimeVariable[0],\
                                ncTimeVariable.units,\
                                ncTimeVariable.calendar) 
    
    return first_datetime.year 

def convert_file(self, file, title="AVISO SLA Track Obs"):
        """
        Load an AVISO file and convert into an obs structure
        """
        # Load AVISO Data
        nc = seapy.netcdf(file)
        lon = nc.variables["longitude"][:]
        lat = nc.variables["latitude"][:]
        slaname = 'SLA' if 'SLA' in nc.variables.keys() else 'sla_filtered'
        dat = nc.variables[slaname][:]
        time = seapy.roms.num2date(nc, "time", epoch=self.epoch)
        nc.close()

        # make them into vectors
        lat = lat.ravel()
        lon = lon.ravel()
        dat = dat.ravel()
        err = np.ones(dat.shape) * _aviso_sla_errors.get(self.provenance, 0.1)

        if not self.grid.east():
            lon[lon > 180] -= 360

        good = dat.nonzero()
        data = [seapy.roms.obs.raw_data("ZETA", self.provenance,
                                        dat[good], err[good], err[0])]
        # Grid it
        obs = seapy.roms.obs.gridder(self.grid, time, lon[good], lat[good], None,
                                     data, self.dt, title)

        # Apply the model mean ssh to the sla data
        if self.ssh_mean is not None and obs is not None:
            m, p = seapy.oasurf(self.grid.I, self.grid.J, self.ssh_mean,
                                obs.x, obs.y, nx=1, ny=1, weight=7)
            obs.value += m

        # Duplicate the observations before and after as per the repeat
        # time unless it is zero
        if self.repeat and obs:
            prior = obs.copy()
            after = obs.copy()
            prior.time -= self.repeat / 24
            after.time += self.repeat / 24
            obs.add(prior)
            obs.add(after)

        return obs 

def checknetcdf(name, start, end):
    """ Check available time steps in netCDF input file
    
    Check available timesteps in netCDF file. Get first and last available timestep in netCDF file and compare with
    first and last computation timestep of the model.
    It can use sub-daily steps.
    
    :param name: string containing path and name of netCDF file
    :param start: initial date or step number of model simulation
    :param end: final date or step of model simulation  
    :return: none
    :raises Exception: stop if netCDF maps do not cover simulation time period
    """
    settings = LisSettings.instance()
    binding = settings.binding
    filename = name + ".nc"
    nf1 = iterOpenNetcdf(filename, "Netcdf map stacks: \n", "r")

    # read information from netCDF file
    t_steps = nf1.variables['time'][:]    # get values for timesteps ([  0.,  24.,  48.,  72.,  96.])
    t_unit = nf1.variables['time'].units  # get unit (u'hours since 2015-01-01 06:00:00')
    t_cal = get_calendar_type(nf1)
    if t_cal != binding['calendar_type']:
        warnings.warn(calendar_inconsistency_warning(filename, t_cal, binding['calendar_type']))

    # get date of first available timestep in netcdf file
    date_first_step_in_ncdf = num2date(t_steps[0], units=t_unit, calendar=t_cal)
    # get date of last available timestep in netcdf file
    date_last_step_in_ncdf = num2date(t_steps[-1], units=t_unit, calendar=t_cal)

    nf1.close()
    #calendar date start (CalendarDayStart)
    begin = calendar(binding['CalendarDayStart'], binding['calendar_type'])
    DtSec = float(binding['DtSec'])
    DtDay = DtSec / 86400.
    # Time step, expressed as fraction of day (same as self.var.DtSec and self.var.DtDay)

    date_first_sim_step = calendar(start, binding['calendar_type'])
    if type(date_first_sim_step) is not datetime.datetime:
        date_first_sim_step = begin + datetime.timedelta(days=(date_first_sim_step - 1) * DtDay)
    if (date_first_sim_step < date_first_step_in_ncdf):
        msg = "First simulation time step is before first time step in netCDF input data file \n" \
              "File name: "+ filename +"\n" \
              "netCDF start date: "+ date_first_step_in_ncdf.strftime('%d/%m/%Y %H:%M') +"\n" \
              "simulation start date: "+ date_first_sim_step.strftime('%d/%m/%Y %H:%M')
        raise LisfloodError(msg)

    date_last_sim_step = calendar(end, binding['calendar_type'])
    if type(date_last_sim_step) is not datetime.datetime:
        date_last_sim_step = begin + datetime.timedelta(days=(date_last_sim_step - 1) * DtDay)
    if (date_last_sim_step > date_last_step_in_ncdf):
        msg = "Last simulation time step is after last time step in netCDF input data file \n" \
              "File name: " + filename +"\n" \
              "netCDF last date: " + date_last_step_in_ncdf.strftime('%d/%m/%Y %H:%M') +"\n" \
              "simulation last date: " + date_last_sim_step.strftime('%d/%m/%Y %H:%M')
        raise LisfloodError(msg)
    return 

def assemble_matches(filename):
    """
    Read a DOMS netCDF file and return a list of matches.
    
    Parameters
    ----------
    filename : str
        The DOMS netCDF file name.
    
    Returns
    -------
    matches : list
        List of matches. Each list element is a dictionary.
        For match m, netCDF group GROUP (SatelliteData or InsituData), and
        group variable VARIABLE:
        matches[m][GROUP]['matchID']: MatchedRecords dimension ID for the match
        matches[m][GROUP]['GROUPID']: GROUP dim dimension ID for the record
        matches[m][GROUP][VARIABLE]: variable value 
    """
    
    try:
        # Open the netCDF file
        with Dataset(filename, 'r') as doms_nc:
            # Check that the number of groups is consistent w/ the MatchedGroups
            # dimension
            assert len(doms_nc.groups) == doms_nc.dimensions['MatchedGroups'].size,\
                ("Number of groups isn't the same as MatchedGroups dimension.")
            
            matches = []
            matched_records = doms_nc.dimensions['MatchedRecords'].size
            
            # Loop through the match IDs to assemble matches
            for match in range(0, matched_records):
                match_dict = OrderedDict()
                # Grab the data from each platform (group) in the match
                for group_num, group in enumerate(doms_nc.groups):
                    match_dict[group] = OrderedDict()
                    match_dict[group]['matchID'] = match
                    ID = doms_nc.variables['matchIDs'][match][group_num]
                    match_dict[group][group + 'ID'] = ID
                    for var in doms_nc.groups[group].variables.keys():
                        match_dict[group][var] = doms_nc.groups[group][var][ID]
                    
                    # Create a UTC datetime field from timestamp
                    dt = num2date(match_dict[group]['time'],
                                  doms_nc.groups[group]['time'].units)
                    match_dict[group]['datetime'] = dt
                LOGGER.info(match_dict)
                matches.append(match_dict)
            
            return matches
    except (OSError, IOError) as err:
        LOGGER.exception("Error reading netCDF file " + filename)
        raise err 

def convert_file(self, file, title="REMSS SST Obs"):
        """
        Load an REMSS file and convert into an obs structure
        """
        # Load REMSS Data
        nc = seapy.netcdf(file)
        lon = nc.variables["lon"][:]
        lat = nc.variables["lat"][:]
        dat = np.ma.masked_outside(np.squeeze(
            nc.variables["sea_surface_temperature"][:]) - 273.15,
            self.temp_limits[0], self.temp_limits[1])
        err = np.ma.masked_outside(np.squeeze(
            nc.variables["SSES_standard_deviation_error"][:]), 0.01, 2.0)
        dat[err.mask] = np.ma.masked

        # Check the data flags
        flags = np.ma.masked_not_equal(
            np.squeeze(nc.variables["rejection_flag"][:]), 0)
        dat[flags.mask] = np.ma.masked
        err[flags.mask] = np.ma.masked

        # Grab the observation time
        time = seapy.roms.num2date(nc, "time", epoch=self.epoch)
        sst_time = nc.variables["sst_dtime"][:] * seapy.secs2day
        for n, i in enumerate(time):
            sst_time[n, :, :] += i
        sst_time[dat.mask] = np.ma.masked

        # Set up the coordinate
        lon, lat = np.meshgrid(lon, lat)
        lon = np.ma.masked_where(dat.mask, seapy.adddim(lon, len(time)))
        lat = np.ma.masked_where(dat.mask, seapy.adddim(lat, len(time)))

        nc.close()

        if self.grid.east():
            lon[lon < 0] += 360
        data = [seapy.roms.obs.raw_data("TEMP", self.provenance,
                                        dat.compressed(),
                                        err.compressed(), self.temp_error)]
        # Grid it
        return seapy.roms.obs.gridder(self.grid, sst_time.compressed(),
                                      lon.compressed(), lat.compressed, None,
                                      data, self.dt, title)