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 2023.9.0). It will be opened with no decoding. Filling values might be wrongly parsed.
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/coding/times.py:319, in decode_cf_datetime(num_dates, units, calendar, use_cftime)
318 try:
--> 319 dates = _decode_datetime_with_pandas(flat_num_dates, units, calendar)
320 except (KeyError, OutOfBoundsDatetime, OutOfBoundsTimedelta, OverflowError):
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/coding/times.py:253, in _decode_datetime_with_pandas(flat_num_dates, units, calendar)
252 time_units, ref_date = _unpack_netcdf_time_units(units)
--> 253 time_units = _netcdf_to_numpy_timeunit(time_units)
254 try:
255 # TODO: the strict enforcement of nanosecond precision Timestamps can be
256 # relaxed when addressing GitHub issue #7493.
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/coding/times.py:115, in _netcdf_to_numpy_timeunit(units)
114 units = f"{units}s"
--> 115 return {
116 "nanoseconds": "ns",
117 "microseconds": "us",
118 "milliseconds": "ms",
119 "seconds": "s",
120 "minutes": "m",
121 "hours": "h",
122 "days": "D",
123 }[units]
KeyError: 'months'
During handling of the above exception, another exception occurred:
ValueError Traceback (most recent call last)
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/coding/times.py:213, in _decode_cf_datetime_dtype(data, units, calendar, use_cftime)
212 try:
--> 213 result = decode_cf_datetime(example_value, units, calendar, use_cftime)
214 except Exception:
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/coding/times.py:321, in decode_cf_datetime(num_dates, units, calendar, use_cftime)
320 except (KeyError, OutOfBoundsDatetime, OutOfBoundsTimedelta, OverflowError):
--> 321 dates = _decode_datetime_with_cftime(
322 flat_num_dates.astype(float), units, calendar
323 )
325 if (
326 dates[np.nanargmin(num_dates)].year < 1678
327 or dates[np.nanargmax(num_dates)].year >= 2262
328 ):
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/coding/times.py:237, in _decode_datetime_with_cftime(num_dates, units, calendar)
235 if num_dates.size > 0:
236 return np.asarray(
--> 237 cftime.num2date(num_dates, units, calendar, only_use_cftime_datetimes=True)
238 )
239 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 ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/conventions.py:428, in decode_cf_variables(variables, attributes, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
427 try:
--> 428 new_vars[k] = decode_cf_variable(
429 k,
430 v,
431 concat_characters=concat_characters,
432 mask_and_scale=mask_and_scale,
433 decode_times=decode_times,
434 stack_char_dim=stack_char_dim,
435 use_cftime=use_cftime,
436 decode_timedelta=decode_timedelta,
437 )
438 except Exception as e:
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/conventions.py:279, in decode_cf_variable(name, var, concat_characters, mask_and_scale, decode_times, decode_endianness, stack_char_dim, use_cftime, decode_timedelta)
278 if decode_times:
--> 279 var = times.CFDatetimeCoder(use_cftime=use_cftime).decode(var, name=name)
281 if decode_endianness and not var.dtype.isnative:
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/coding/times.py:831, in CFDatetimeCoder.decode(self, variable, name)
830 calendar = pop_to(attrs, encoding, "calendar")
--> 831 dtype = _decode_cf_datetime_dtype(data, units, calendar, self.use_cftime)
832 transform = partial(
833 decode_cf_datetime,
834 units=units,
835 calendar=calendar,
836 use_cftime=self.use_cftime,
837 )
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/coding/times.py:223, in _decode_cf_datetime_dtype(data, units, calendar, use_cftime)
218 msg = (
219 f"unable to decode time units {units!r} with {calendar_msg!r}. Try "
220 "opening your dataset with decode_times=False or installing cftime "
221 "if it is not installed."
222 )
--> 223 raise ValueError(msg)
224 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 ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/conventions.py:569, in decode_cf(obj, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
567 raise TypeError("can only decode Dataset or DataStore objects")
--> 569 vars, attrs, coord_names = decode_cf_variables(
570 vars,
571 attrs,
572 concat_characters,
573 mask_and_scale,
574 decode_times,
575 decode_coords,
576 drop_variables=drop_variables,
577 use_cftime=use_cftime,
578 decode_timedelta=decode_timedelta,
579 )
580 ds = Dataset(vars, attrs=attrs)
File ~/miniconda3/envs/parcels/lib/python3.11/site-packages/xarray/conventions.py:439, in decode_cf_variables(variables, attributes, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
438 except Exception as e:
--> 439 raise type(e)(f"Failed to decode variable {k!r}: {e}")
440 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 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:457, in Field.from_netcdf(cls, filenames, variable, dimensions, indices, grid, mesh, timestamps, allow_time_extrapolation, time_periodic, deferred_load, **kwargs)
452 raise RuntimeError('Multiple files given but no time dimension specified')
454 if grid is None:
455 # Concatenate time variable to determine overall dimension
456 # across multiple files
--> 457 time, time_origin, timeslices, dataFiles = cls.collect_timeslices(timestamps, data_filenames,
458 _grid_fb_class, dimensions,
459 indices, netcdf_engine, netcdf_decodewarning)
460 grid = Grid.create_grid(lon, lat, depth, time, time_origin=time_origin, mesh=mesh)
461 grid.timeslices = timeslices
File ~/Codes/ParcelsCode/parcels/field.py:290, in Field.collect_timeslices(timestamps, data_filenames, _grid_fb_class, dimensions, indices, netcdf_engine, netcdf_decodewarning)
287 for fname in data_filenames:
288 with _grid_fb_class(fname, dimensions, indices, netcdf_engine=netcdf_engine,
289 netcdf_decodewarning=netcdf_decodewarning) as filebuffer:
--> 290 ftime = filebuffer.time
291 timeslices.append(ftime)
292 dataFiles.append([fname] * len(ftime))
File ~/Codes/ParcelsCode/parcels/fieldfilebuffer.py:215, in NetcdfFileBuffer.time(self)
213 @property
214 def time(self):
--> 215 return self.time_access()
File ~/Codes/ParcelsCode/parcels/fieldfilebuffer.py:225, in NetcdfFileBuffer.time_access(self)
222 return np.array([None])
224 time_da = self.dataset[self.dimensions['time']]
--> 225 convert_xarray_time_units(time_da, self.dimensions['time'])
226 time = np.array([time_da[self.dimensions['time']].data]) if len(time_da.shape) == 0 else np.array(time_da[self.dimensions['time']])
227 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 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.