, or view it on GitHub. Notebook version corresponds with v3-support.
TimeStamps and calendars#
[1]:
import warnings
from glob import glob
import numpy as np
import parcels
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 = parcels.download_example_dataset("WOA_data")
tempfield = parcels.Field.from_netcdf(
glob(f"{example_dataset_folder}/woa18_decav_*_04.nc"),
"t_an",
{"lon": "lon", "lat": "lat", "time": "time"},
)
C:\Users\asche\Desktop\po-code\parcels_dev\parcels\parcels\field.py:502: FileWarning: File C:\Users\asche\AppData\Local\parcels\parcels\Cache\WOA_data\woa18_decav_t01_04.nc could not be decoded properly by xarray (version 2024.6.0). It will be opened with no decoding. Filling values might be wrongly parsed.
with _grid_fb_class(
C:\Users\asche\Desktop\po-code\parcels_dev\parcels\parcels\field.py:330: FileWarning: File C:\Users\asche\AppData\Local\parcels\parcels\Cache\WOA_data\woa18_decav_t01_04.nc could not be decoded properly by xarray (version 2024.6.0). It will be opened with no decoding. Filling values might be wrongly parsed.
with _grid_fb_class(
---------------------------------------------------------------------------
KeyError Traceback (most recent call last)
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\coding\times.py:322, in decode_cf_datetime(num_dates, units, calendar, use_cftime)
321 try:
--> 322 dates = _decode_datetime_with_pandas(flat_num_dates, units, calendar)
323 except (KeyError, OutOfBoundsDatetime, OutOfBoundsTimedelta, OverflowError):
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\coding\times.py:256, in _decode_datetime_with_pandas(flat_num_dates, units, calendar)
255 time_units, ref_date = _unpack_netcdf_time_units(units)
--> 256 time_units = _netcdf_to_numpy_timeunit(time_units)
257 try:
258 # TODO: the strict enforcement of nanosecond precision Timestamps can be
259 # relaxed when addressing GitHub issue #7493.
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\coding\times.py:118, in _netcdf_to_numpy_timeunit(units)
117 units = f"{units}s"
--> 118 return {
119 "nanoseconds": "ns",
120 "microseconds": "us",
121 "milliseconds": "ms",
122 "seconds": "s",
123 "minutes": "m",
124 "hours": "h",
125 "days": "D",
126 }[units]
KeyError: 'months'
During handling of the above exception, another exception occurred:
ValueError Traceback (most recent call last)
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\coding\times.py:216, in _decode_cf_datetime_dtype(data, units, calendar, use_cftime)
215 try:
--> 216 result = decode_cf_datetime(example_value, units, calendar, use_cftime)
217 except Exception:
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\coding\times.py:324, in decode_cf_datetime(num_dates, units, calendar, use_cftime)
323 except (KeyError, OutOfBoundsDatetime, OutOfBoundsTimedelta, OverflowError):
--> 324 dates = _decode_datetime_with_cftime(
325 flat_num_dates.astype(float), units, calendar
326 )
328 if (
329 dates[np.nanargmin(num_dates)].year < 1678
330 or dates[np.nanargmax(num_dates)].year >= 2262
331 ):
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\coding\times.py:240, in _decode_datetime_with_cftime(num_dates, units, calendar)
238 if num_dates.size > 0:
239 return np.asarray(
--> 240 cftime.num2date(num_dates, units, calendar, only_use_cftime_datetimes=True)
241 )
242 else:
File src\\cftime\\_cftime.pyx:587, in cftime._cftime.num2date()
File src\\cftime\\_cftime.pyx:101, 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 c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\conventions.py:440, in decode_cf_variables(variables, attributes, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
439 try:
--> 440 new_vars[k] = decode_cf_variable(
441 k,
442 v,
443 concat_characters=concat_characters,
444 mask_and_scale=mask_and_scale,
445 decode_times=decode_times,
446 stack_char_dim=stack_char_dim,
447 use_cftime=use_cftime,
448 decode_timedelta=decode_timedelta,
449 )
450 except Exception as e:
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\conventions.py:291, in decode_cf_variable(name, var, concat_characters, mask_and_scale, decode_times, decode_endianness, stack_char_dim, use_cftime, decode_timedelta)
290 if decode_times:
--> 291 var = times.CFDatetimeCoder(use_cftime=use_cftime).decode(var, name=name)
293 if decode_endianness and not var.dtype.isnative:
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\coding\times.py:987, in CFDatetimeCoder.decode(self, variable, name)
986 calendar = pop_to(attrs, encoding, "calendar")
--> 987 dtype = _decode_cf_datetime_dtype(data, units, calendar, self.use_cftime)
988 transform = partial(
989 decode_cf_datetime,
990 units=units,
991 calendar=calendar,
992 use_cftime=self.use_cftime,
993 )
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\coding\times.py:226, in _decode_cf_datetime_dtype(data, units, calendar, use_cftime)
221 msg = (
222 f"unable to decode time units {units!r} with {calendar_msg!r}. Try "
223 "opening your dataset with decode_times=False or installing cftime "
224 "if it is not installed."
225 )
--> 226 raise ValueError(msg)
227 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.
The above exception was the direct cause of the following exception:
ValueError Traceback (most recent call last)
File ~\Desktop\po-code\parcels_dev\parcels\parcels\tools\converters.py:281, in convert_xarray_time_units(ds, time)
280 try:
--> 281 da2 = xr.decode_cf(da2)
282 except ValueError:
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\conventions.py:581, in decode_cf(obj, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
579 raise TypeError("can only decode Dataset or DataStore objects")
--> 581 vars, attrs, coord_names = decode_cf_variables(
582 vars,
583 attrs,
584 concat_characters,
585 mask_and_scale,
586 decode_times,
587 decode_coords,
588 drop_variables=drop_variables,
589 use_cftime=use_cftime,
590 decode_timedelta=decode_timedelta,
591 )
592 ds = Dataset(vars, attrs=attrs)
File c:\Users\asche\miniconda3\envs\parcels_dev\Lib\site-packages\xarray\conventions.py:451, in decode_cf_variables(variables, attributes, concat_characters, mask_and_scale, decode_times, decode_coords, drop_variables, use_cftime, decode_timedelta)
450 except Exception as e:
--> 451 raise type(e)(f"Failed to decode variable {k!r}: {e}") from e
452 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)
Cell In[2], line 2
1 example_dataset_folder = parcels.download_example_dataset("WOA_data")
----> 2 tempfield = parcels.Field.from_netcdf(
3 glob(f"{example_dataset_folder}/woa18_decav_*_04.nc"),
4 "t_an",
5 {"lon": "lon", "lat": "lat", "time": "time"},
6 )
File ~\Desktop\po-code\parcels_dev\parcels\parcels\field.py:540, in Field.from_netcdf(cls, filenames, variable, dimensions, indices, grid, mesh, timestamps, allow_time_extrapolation, time_periodic, deferred_load, **kwargs)
535 raise RuntimeError("Multiple files given but no time dimension specified")
537 if grid is None:
538 # Concatenate time variable to determine overall dimension
539 # across multiple files
--> 540 time, time_origin, timeslices, dataFiles = cls.collect_timeslices(
541 timestamps, data_filenames, _grid_fb_class, dimensions, indices, netcdf_engine, netcdf_decodewarning
542 )
543 grid = Grid.create_grid(lon, lat, depth, time, time_origin=time_origin, mesh=mesh)
544 grid.timeslices = timeslices
File ~\Desktop\po-code\parcels_dev\parcels\parcels\field.py:333, in Field.collect_timeslices(timestamps, data_filenames, _grid_fb_class, dimensions, indices, netcdf_engine, netcdf_decodewarning)
329 for fname in data_filenames:
330 with _grid_fb_class(
331 fname, dimensions, indices, netcdf_engine=netcdf_engine, netcdf_decodewarning=netcdf_decodewarning
332 ) as filebuffer:
--> 333 ftime = filebuffer.time
334 timeslices.append(ftime)
335 dataFiles.append([fname] * len(ftime))
File ~\Desktop\po-code\parcels_dev\parcels\parcels\fieldfilebuffer.py:221, in NetcdfFileBuffer.time(self)
219 @property
220 def time(self):
--> 221 return self.time_access()
File ~\Desktop\po-code\parcels_dev\parcels\parcels\fieldfilebuffer.py:231, in NetcdfFileBuffer.time_access(self)
228 return np.array([None])
230 time_da = self.dataset[self.dimensions["time"]]
--> 231 convert_xarray_time_units(time_da, self.dimensions["time"])
232 time = (
233 np.array([time_da[self.dimensions["time"]].data])
234 if len(time_da.shape) == 0
235 else np.array(time_da[self.dimensions["time"]])
236 )
237 if isinstance(time[0], datetime.datetime):
File ~\Desktop\po-code\parcels_dev\parcels\parcels\tools\converters.py:283, in convert_xarray_time_units(ds, time)
281 da2 = xr.decode_cf(da2)
282 except ValueError:
--> 283 raise RuntimeError(
284 "Xarray could not convert the calendar. If you're using from_netcdf, "
285 "try using the timestamps keyword in the construction of your Field. "
286 "See also the tutorial at https://docs.oceanparcels.org/en/latest/examples/tutorial_timestamps.html"
287 )
288 ds[time] = da2[time]
RuntimeError: Xarray could not convert the calendar. If you're 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(f"2001-{m:02d}-15") for m in range(1, 13)]), axis=1
)
And then we can add the timestamps as an extra argument
[4]:
with warnings.catch_warnings():
warnings.simplefilter("ignore", parcels.FileWarning)
tempfield = parcels.Field.from_netcdf(
glob(f"{example_dataset_folder}/woa18_decav_*_04.nc"),
"t_an",
{"lon": "lon", "lat": "lat", "time": "time"},
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 warnings.catch_warnings() with warnings.simplefilter("ignore", parcels.FileWarning) to wrap the FieldSet.from_nemo() call above. This is to silence an expected warning because the time dimension in the coordinates.nc file can’t be decoded by xarray.