TimeStamps and calendars#
[1]:
from glob import glob
import numpy as np
from parcels import Field, download_example_dataset
Some NetCDF files, such as for example those from the World Ocean Atlas, have time calendars that can’t be parsed by xarray. These result in a ValueError: unable to decode time units, for example when the calendar is in ‘months since’ a particular date.
In these cases, a workaround in Parcels is to use the timestamps argument in Field (or FieldSet) creation. Here, we show how this works for example temperature data from the World Ocean Atlas in the Pacific Ocean
The following cell raises an error, since the calendar of the World Ocean Atlas data is in “months since 1955-01-01 00:00:00”
[2]:
example_dataset_folder = download_example_dataset("WOA_data")
tempfield = Field.from_netcdf(
glob(f"{example_dataset_folder}/woa18_decav_*_04.nc"),
"t_an",
{"lon": "lon", "lat": "lat", "time": "time"},
)
WARNING: File /Users/erik/Library/Caches/parcels/WOA_data/woa18_decav_t10_04.nc could not be decoded properly by xarray (version 2022.11.0).
It will be opened with no decoding. Filling values might be wrongly parsed.
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/coding/times.py:275, in decode_cf_datetime(num_dates, units, calendar, use_cftime)
274 try:
--> 275 dates = _decode_datetime_with_pandas(flat_num_dates, units, calendar)
276 except (KeyError, OutOfBoundsDatetime, OutOfBoundsTimedelta, OverflowError):
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/coding/times.py:216, in _decode_datetime_with_pandas(flat_num_dates, units, calendar)
215 delta, ref_date = _unpack_netcdf_time_units(units)
--> 216 delta = _netcdf_to_numpy_timeunit(delta)
217 try:
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/coding/times.py:110, in _netcdf_to_numpy_timeunit(units)
109 units = f"{units}s"
--> 110 return {
111 "nanoseconds": "ns",
112 "microseconds": "us",
113 "milliseconds": "ms",
114 "seconds": "s",
115 "minutes": "m",
116 "hours": "h",
117 "days": "D",
118 }[units]
KeyError: 'months'
During handling of the above exception, another exception occurred:
ValueError Traceback (most recent call last)
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/coding/times.py:180, in _decode_cf_datetime_dtype(data, units, calendar, use_cftime)
179 try:
--> 180 result = decode_cf_datetime(example_value, units, calendar, use_cftime)
181 except Exception:
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/coding/times.py:277, in decode_cf_datetime(num_dates, units, calendar, use_cftime)
276 except (KeyError, OutOfBoundsDatetime, OutOfBoundsTimedelta, OverflowError):
--> 277 dates = _decode_datetime_with_cftime(
278 flat_num_dates.astype(float), units, calendar
279 )
281 if (
282 dates[np.nanargmin(num_dates)].year < 1678
283 or dates[np.nanargmax(num_dates)].year >= 2262
284 ):
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/coding/times.py:202, in _decode_datetime_with_cftime(num_dates, units, calendar)
200 if num_dates.size > 0:
201 return np.asarray(
--> 202 cftime.num2date(num_dates, units, calendar, only_use_cftime_datetimes=True)
203 )
204 else:
File src/cftime/_cftime.pyx:580, in cftime._cftime.num2date()
File src/cftime/_cftime.pyx:98, in cftime._cftime._dateparse()
ValueError: 'months since' units only allowed for '360_day' calendar
During handling of the above exception, another exception occurred:
ValueError Traceback (most recent call last)
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/conventions.py:523, in decode_cf_variables(variables, attributes, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
522 try:
--> 523 new_vars[k] = decode_cf_variable(
524 k,
525 v,
526 concat_characters=concat_characters,
527 mask_and_scale=mask_and_scale,
528 decode_times=decode_times,
529 stack_char_dim=stack_char_dim,
530 use_cftime=use_cftime,
531 decode_timedelta=decode_timedelta,
532 )
533 except Exception as e:
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/conventions.py:369, in decode_cf_variable(name, var, concat_characters, mask_and_scale, decode_times, decode_endianness, stack_char_dim, use_cftime, decode_timedelta)
368 if decode_times:
--> 369 var = times.CFDatetimeCoder(use_cftime=use_cftime).decode(var, name=name)
371 dimensions, data, attributes, encoding = variables.unpack_for_decoding(var)
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/coding/times.py:699, in CFDatetimeCoder.decode(self, variable, name)
698 calendar = pop_to(attrs, encoding, "calendar")
--> 699 dtype = _decode_cf_datetime_dtype(data, units, calendar, self.use_cftime)
700 transform = partial(
701 decode_cf_datetime,
702 units=units,
703 calendar=calendar,
704 use_cftime=self.use_cftime,
705 )
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/coding/times.py:190, in _decode_cf_datetime_dtype(data, units, calendar, use_cftime)
185 msg = (
186 f"unable to decode time units {units!r} with {calendar_msg!r}. Try "
187 "opening your dataset with decode_times=False or installing cftime "
188 "if it is not installed."
189 )
--> 190 raise ValueError(msg)
191 else:
ValueError: unable to decode time units 'months since 1955-01-01 00:00:00' with 'the default calendar'. Try opening your dataset with decode_times=False or installing cftime if it is not installed.
During handling of the above exception, another exception occurred:
ValueError Traceback (most recent call last)
File ~/Codes/ParcelsCode/parcels/tools/converters.py:266, in convert_xarray_time_units(ds, time)
265 try:
--> 266 da2 = xr.decode_cf(da2)
267 except ValueError:
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/conventions.py:659, in decode_cf(obj, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
657 raise TypeError("can only decode Dataset or DataStore objects")
--> 659 vars, attrs, coord_names = decode_cf_variables(
660 vars,
661 attrs,
662 concat_characters,
663 mask_and_scale,
664 decode_times,
665 decode_coords,
666 drop_variables=drop_variables,
667 use_cftime=use_cftime,
668 decode_timedelta=decode_timedelta,
669 )
670 ds = Dataset(vars, attrs=attrs)
File ~/opt/anaconda3/envs/py3_parcels/lib/python3.8/site-packages/xarray/conventions.py:534, in decode_cf_variables(variables, attributes, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
533 except Exception as e:
--> 534 raise type(e)(f"Failed to decode variable {k!r}: {e}")
535 if decode_coords in [True, "coordinates", "all"]:
ValueError: Failed to decode variable 'time': unable to decode time units 'months since 1955-01-01 00:00:00' with 'the default calendar'. Try opening your dataset with decode_times=False or installing cftime if it is not installed.
During handling of the above exception, another exception occurred:
RuntimeError Traceback (most recent call last)
/Users/erik/Codes/ParcelsCode/docs/examples/tutorial_timestamps.ipynb Cell 5 in <cell line: 2>()
<a href='vscode-notebook-cell:/Users/erik/Codes/ParcelsCode/docs/examples/tutorial_timestamps.ipynb#W4sZmlsZQ%3D%3D?line=0'>1</a> example_dataset_folder = download_example_dataset("WOA_data")
----> <a href='vscode-notebook-cell:/Users/erik/Codes/ParcelsCode/docs/examples/tutorial_timestamps.ipynb#W4sZmlsZQ%3D%3D?line=1'>2</a> tempfield = Field.from_netcdf(
<a href='vscode-notebook-cell:/Users/erik/Codes/ParcelsCode/docs/examples/tutorial_timestamps.ipynb#W4sZmlsZQ%3D%3D?line=2'>3</a> glob(f"{example_dataset_folder}/woa18_decav_*_04.nc"),
<a href='vscode-notebook-cell:/Users/erik/Codes/ParcelsCode/docs/examples/tutorial_timestamps.ipynb#W4sZmlsZQ%3D%3D?line=3'>4</a> "t_an",
<a href='vscode-notebook-cell:/Users/erik/Codes/ParcelsCode/docs/examples/tutorial_timestamps.ipynb#W4sZmlsZQ%3D%3D?line=4'>5</a> {"lon": "lon", "lat": "lat", "time": "time"},
<a href='vscode-notebook-cell:/Users/erik/Codes/ParcelsCode/docs/examples/tutorial_timestamps.ipynb#W4sZmlsZQ%3D%3D?line=5'>6</a> )
File ~/Codes/ParcelsCode/parcels/field.py:443, in Field.from_netcdf(cls, filenames, variable, dimensions, indices, grid, mesh, timestamps, allow_time_extrapolation, time_periodic, deferred_load, **kwargs)
438 raise RuntimeError('Multiple files given but no time dimension specified')
440 if grid is None:
441 # Concatenate time variable to determine overall dimension
442 # across multiple files
--> 443 time, time_origin, timeslices, dataFiles = cls.collect_timeslices(timestamps, data_filenames,
444 _grid_fb_class, dimensions,
445 indices, netcdf_engine, netcdf_decodewarning)
446 grid = Grid.create_grid(lon, lat, depth, time, time_origin=time_origin, mesh=mesh)
447 grid.timeslices = timeslices
File ~/Codes/ParcelsCode/parcels/field.py:276, in Field.collect_timeslices(timestamps, data_filenames, _grid_fb_class, dimensions, indices, netcdf_engine, netcdf_decodewarning)
273 for fname in data_filenames:
274 with _grid_fb_class(fname, dimensions, indices, netcdf_engine=netcdf_engine,
275 netcdf_decodewarning=netcdf_decodewarning) as filebuffer:
--> 276 ftime = filebuffer.time
277 timeslices.append(ftime)
278 dataFiles.append([fname] * len(ftime))
File ~/Codes/ParcelsCode/parcels/fieldfilebuffer.py:214, in NetcdfFileBuffer.time(self)
212 @property
213 def time(self):
--> 214 return self.time_access()
File ~/Codes/ParcelsCode/parcels/fieldfilebuffer.py:224, in NetcdfFileBuffer.time_access(self)
221 return np.array([None])
223 time_da = self.dataset[self.dimensions['time']]
--> 224 convert_xarray_time_units(time_da, self.dimensions['time'])
225 time = np.array([time_da[self.dimensions['time']].data]) if len(time_da.shape) == 0 else np.array(time_da[self.dimensions['time']])
226 if isinstance(time[0], datetime.datetime):
File ~/Codes/ParcelsCode/parcels/tools/converters.py:268, in convert_xarray_time_units(ds, time)
266 da2 = xr.decode_cf(da2)
267 except ValueError:
--> 268 raise RuntimeError('Xarray could not convert the calendar. If you''re using from_netcdf, '
269 'try using the timestamps keyword in the construction of your Field. '
270 'See also the tutorial at https://docs.oceanparcels.org/en/latest/'
271 'examples/tutorial_timestamps.html')
272 ds[time] = da2[time]
RuntimeError: Xarray could not convert the calendar. If youre using from_netcdf, try using the timestamps keyword in the construction of your Field. See also the tutorial at https://docs.oceanparcels.org/en/latest/examples/tutorial_timestamps.html
However, we can create our own numpy array of timestamps associated with each of the 12 snapshots in the netcdf file
[3]:
timestamps = np.expand_dims(
np.array([np.datetime64("2001-%.2d-15" % m) for m in range(1, 13)]), axis=1
)
And then we can add the timestamps as an extra argument
[4]:
tempfield = Field.from_netcdf(
glob(f"{example_dataset_folder}/woa18_decav_*_04.nc"),
"t_an",
{"lon": "lon", "lat": "lat", "time": "time"},
netcdf_decodewarning=False,
timestamps=timestamps,
)
Note, by the way, that adding the time_periodic=True argument to Field.from_netcdf() will also mean that the climatology can be cycled for multiple years.
Furthermore, note that we used netcdf_decodewarning=False in the FieldSet.from_netcdf() call above. This is to silence an expected warning because the time dimension in the coordinates.nc file can’t be decoded by xarray.