Trajan demo

From OpenDrift 2.0, analysis and plotting of results from OpenDrift simulations will be handled by a new, standalone package: Trajan https://github.com/OpenDrift/trajan

This example creates a test dataset, and demonstrates its anlysis using Trajan

import os
from datetime import datetime, timedelta
import matplotlib.pyplot as plt
import xarray as xr
import trajan as ta
from opendrift.readers import reader_netCDF_CF_generic
from opendrift.models.openoil import OpenOil

Create test dataset with OpenDrift

o = OpenOil(loglevel=20)

# Add forcing
reader_arome = reader_netCDF_CF_generic.Reader(o.test_data_folder() +
    '16Nov2015_NorKyst_z_surface/arome_subset_16Nov2015.nc')
reader_norkyst = reader_netCDF_CF_generic.Reader(o.test_data_folder() +
    '16Nov2015_NorKyst_z_surface/norkyst800_subset_16Nov2015.nc')
o.add_reader([reader_norkyst, reader_arome])

# Seeding some particles
o.seed_elements(lon=4.4, lat=60.1, radius=1000, number=1000,
                time=reader_arome.start_time)

# Running model
o.run(end_time=reader_norkyst.end_time, outfile='openoil.nc')
12:55:13 INFO    opendrift.models.basemodel:515: OpenDriftSimulation initialised (version 1.11.13 / v1.11.13-99-gd2132d3)
12:55:13 INFO    opendrift.readers.reader_netCDF_CF_generic:102: Opening dataset: /root/project/tests/test_data/16Nov2015_NorKyst_z_surface/arome_subset_16Nov2015.nc
/opt/conda/envs/opendrift/lib/python3.11/site-packages/pyproj/crs/crs.py:1286: UserWarning: You will likely lose important projection information when converting to a PROJ string from another format. See: https://proj.org/faq.html#what-is-the-best-format-for-describing-coordinate-reference-systems
  proj = self._crs.to_proj4(version=version)
12:55:13 INFO    opendrift.readers.reader_netCDF_CF_generic:325: Detected dimensions: {'time': 'time', 'x': 'x', 'y': 'y'}
12:55:13 INFO    opendrift.readers.reader_netCDF_CF_generic:102: Opening dataset: /root/project/tests/test_data/16Nov2015_NorKyst_z_surface/norkyst800_subset_16Nov2015.nc
/opt/conda/envs/opendrift/lib/python3.11/site-packages/pyproj/crs/crs.py:1286: UserWarning: You will likely lose important projection information when converting to a PROJ string from another format. See: https://proj.org/faq.html#what-is-the-best-format-for-describing-coordinate-reference-systems
  proj = self._crs.to_proj4(version=version)
12:55:14 INFO    opendrift.readers.reader_netCDF_CF_generic:325: Detected dimensions: {'x': 'X', 'y': 'Y', 'z': 'depth', 'time': 'time'}
12:55:14 INFO    opendrift.models.openoil.openoil:1710: Oil type not specified, using default: GENERIC BUNKER C
12:55:14 INFO    opendrift.models.openoil.adios.dirjs:90: Querying ADIOS database for oil: GENERIC BUNKER C
12:55:14 INFO    opendrift.models.openoil.openoil:1719: Using density 971.1 and viscosity 0.0005020658058702914 of oiltype GENERIC BUNKER C
12:55:14 INFO    opendrift.models.basemodel.environment:218: Adding a dynamical landmask with max. priority based on assumed maximum speed of 1.3 m/s. Adding a customised landmask may be faster...
12:55:17 INFO    opendrift.models.basemodel.environment:245: Fallback values will be used for the following variables which have no readers:
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_surface_height: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    upward_sea_water_velocity: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_surface_wave_significant_height: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_surface_wave_stokes_drift_x_velocity: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_surface_wave_stokes_drift_y_velocity: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_surface_wave_period_at_variance_spectral_density_maximum: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_surface_wave_mean_period_from_variance_spectral_density_second_frequency_moment: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_ice_area_fraction: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_ice_x_velocity: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_ice_y_velocity: 0.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_water_temperature: 10.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_water_salinity: 34.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    sea_floor_depth_below_sea_level: 10000.000000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    ocean_vertical_diffusivity: 0.020000
12:55:17 INFO    opendrift.models.basemodel.environment:248:    ocean_mixed_layer_thickness: 50.000000
12:55:18 INFO    opendrift.models.basemodel:936: Using existing reader for land_binary_mask
12:55:18 INFO    opendrift.models.basemodel:947: All points are in ocean
12:55:18 INFO    opendrift.models.openoil.openoil:687: Oil-water surface tension is 0.035282 Nm
12:55:18 INFO    opendrift.models.openoil.openoil:700: Max water fraction not available for GENERIC BUNKER C, using default
12:55:18 INFO    opendrift.models.basemodel:2038: 2015-11-16 00:00:00 - step 1 of 66 - 1000 active elements (0 deactivated)
12:55:18 INFO    opendrift.models.basemodel:2038: 2015-11-16 01:00:00 - step 2 of 66 - 1000 active elements (0 deactivated)
12:55:18 INFO    opendrift.models.basemodel:2038: 2015-11-16 02:00:00 - step 3 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 03:00:00 - step 4 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 04:00:00 - step 5 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 05:00:00 - step 6 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 06:00:00 - step 7 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 07:00:00 - step 8 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 08:00:00 - step 9 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 09:00:00 - step 10 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 10:00:00 - step 11 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 11:00:00 - step 12 of 66 - 1000 active elements (0 deactivated)
12:55:19 INFO    opendrift.models.basemodel:2038: 2015-11-16 12:00:00 - step 13 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 13:00:00 - step 14 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 14:00:00 - step 15 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 15:00:00 - step 16 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 16:00:00 - step 17 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 17:00:00 - step 18 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 18:00:00 - step 19 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 19:00:00 - step 20 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 20:00:00 - step 21 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 21:00:00 - step 22 of 66 - 1000 active elements (0 deactivated)
12:55:20 INFO    opendrift.models.basemodel:2038: 2015-11-16 22:00:00 - step 23 of 66 - 1000 active elements (0 deactivated)
12:55:21 INFO    opendrift.models.basemodel:2038: 2015-11-16 23:00:00 - step 24 of 66 - 1000 active elements (0 deactivated)
12:55:21 INFO    opendrift.models.basemodel:2038: 2015-11-17 00:00:00 - step 25 of 66 - 1000 active elements (0 deactivated)
12:55:21 INFO    opendrift.models.basemodel:2038: 2015-11-17 01:00:00 - step 26 of 66 - 1000 active elements (0 deactivated)
12:55:21 WARNING opendrift.readers.basereader.structured:324: Data block from /root/project/tests/test_data/16Nov2015_NorKyst_z_surface/norkyst800_subset_16Nov2015.nc not large enough to cover element positions within timestep. Buffer size (8) must be increased. See `Variables.set_buffer_size`.
12:55:21 INFO    opendrift.models.basemodel:2038: 2015-11-17 02:00:00 - step 27 of 66 - 1000 active elements (0 deactivated)
12:55:21 WARNING opendrift.readers.basereader.structured:324: Data block from /root/project/tests/test_data/16Nov2015_NorKyst_z_surface/norkyst800_subset_16Nov2015.nc not large enough to cover element positions within timestep. Buffer size (8) must be increased. See `Variables.set_buffer_size`.
/root/project/opendrift/readers/interpolation/interpolators.py:17: RuntimeWarning: overflow encountered in cast
  data[mask] = np.finfo(np.float64).min
12:55:21 WARNING opendrift.readers.basereader.structured:324: Data block from /root/project/tests/test_data/16Nov2015_NorKyst_z_surface/arome_subset_16Nov2015.nc not large enough to cover element positions within timestep. Buffer size (4) must be increased. See `Variables.set_buffer_size`.
12:55:21 INFO    opendrift.models.basemodel:2038: 2015-11-17 03:00:00 - step 28 of 66 - 1000 active elements (0 deactivated)
12:55:21 WARNING opendrift.readers.basereader.structured:324: Data block from /root/project/tests/test_data/16Nov2015_NorKyst_z_surface/norkyst800_subset_16Nov2015.nc not large enough to cover element positions within timestep. Buffer size (8) must be increased. See `Variables.set_buffer_size`.
12:55:24 INFO    opendrift.models.basemodel:2038: 2015-11-17 04:00:00 - step 29 of 66 - 998 active elements (2 deactivated)
12:55:24 INFO    opendrift.models.basemodel:2038: 2015-11-17 05:00:00 - step 30 of 66 - 991 active elements (9 deactivated)
12:55:24 INFO    opendrift.models.basemodel:2038: 2015-11-17 06:00:00 - step 31 of 66 - 980 active elements (20 deactivated)
12:55:24 INFO    opendrift.models.basemodel:2038: 2015-11-17 07:00:00 - step 32 of 66 - 963 active elements (37 deactivated)
12:55:24 INFO    opendrift.models.basemodel:2038: 2015-11-17 08:00:00 - step 33 of 66 - 920 active elements (80 deactivated)
12:55:24 INFO    opendrift.models.basemodel:2038: 2015-11-17 09:00:00 - step 34 of 66 - 871 active elements (129 deactivated)
12:55:24 INFO    opendrift.models.basemodel:2038: 2015-11-17 10:00:00 - step 35 of 66 - 820 active elements (180 deactivated)
12:55:25 INFO    opendrift.models.basemodel:2038: 2015-11-17 11:00:00 - step 36 of 66 - 784 active elements (216 deactivated)
12:55:25 INFO    opendrift.models.basemodel:2038: 2015-11-17 12:00:00 - step 37 of 66 - 754 active elements (246 deactivated)
12:55:25 INFO    opendrift.models.basemodel:2038: 2015-11-17 13:00:00 - step 38 of 66 - 704 active elements (296 deactivated)
12:55:25 INFO    opendrift.models.basemodel:2038: 2015-11-17 14:00:00 - step 39 of 66 - 568 active elements (432 deactivated)
12:55:25 INFO    opendrift.models.basemodel:2038: 2015-11-17 15:00:00 - step 40 of 66 - 412 active elements (588 deactivated)
12:55:25 INFO    opendrift.models.basemodel:2038: 2015-11-17 16:00:00 - step 41 of 66 - 325 active elements (675 deactivated)
12:55:25 INFO    opendrift.models.basemodel:2038: 2015-11-17 17:00:00 - step 42 of 66 - 285 active elements (715 deactivated)
12:55:25 INFO    opendrift.models.basemodel:2038: 2015-11-17 18:00:00 - step 43 of 66 - 251 active elements (749 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-17 19:00:00 - step 44 of 66 - 235 active elements (765 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-17 20:00:00 - step 45 of 66 - 221 active elements (779 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-17 21:00:00 - step 46 of 66 - 204 active elements (796 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-17 22:00:00 - step 47 of 66 - 198 active elements (802 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-17 23:00:00 - step 48 of 66 - 189 active elements (811 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 00:00:00 - step 49 of 66 - 181 active elements (819 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 01:00:00 - step 50 of 66 - 178 active elements (822 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 02:00:00 - step 51 of 66 - 178 active elements (822 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 03:00:00 - step 52 of 66 - 174 active elements (826 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 04:00:00 - step 53 of 66 - 172 active elements (828 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 05:00:00 - step 54 of 66 - 171 active elements (829 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 06:00:00 - step 55 of 66 - 171 active elements (829 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 07:00:00 - step 56 of 66 - 171 active elements (829 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 08:00:00 - step 57 of 66 - 168 active elements (832 deactivated)
12:55:26 INFO    opendrift.models.basemodel:2038: 2015-11-18 09:00:00 - step 58 of 66 - 166 active elements (834 deactivated)
12:55:27 INFO    opendrift.models.basemodel:2038: 2015-11-18 10:00:00 - step 59 of 66 - 162 active elements (838 deactivated)
12:55:27 INFO    opendrift.models.basemodel:2038: 2015-11-18 11:00:00 - step 60 of 66 - 160 active elements (840 deactivated)
12:55:27 INFO    opendrift.models.basemodel:2038: 2015-11-18 12:00:00 - step 61 of 66 - 159 active elements (841 deactivated)
12:55:27 INFO    opendrift.models.basemodel:2038: 2015-11-18 13:00:00 - step 62 of 66 - 157 active elements (843 deactivated)
12:55:27 INFO    opendrift.models.basemodel:2038: 2015-11-18 14:00:00 - step 63 of 66 - 156 active elements (844 deactivated)
12:55:27 INFO    opendrift.models.basemodel:2038: 2015-11-18 15:00:00 - step 64 of 66 - 156 active elements (844 deactivated)
12:55:27 INFO    opendrift.models.basemodel:2038: 2015-11-18 16:00:00 - step 65 of 66 - 155 active elements (845 deactivated)
12:55:27 INFO    opendrift.models.basemodel:2038: 2015-11-18 17:00:00 - step 66 of 66 - 154 active elements (846 deactivated)
12:55:27 INFO    opendrift.export.io_netcdf:121: Wrote 67 steps to file openoil.nc
/opt/conda/envs/opendrift/lib/python3.11/site-packages/numpy/ma/core.py:467: RuntimeWarning: invalid value encountered in cast
  fill_value = np.array(fill_value, copy=False, dtype=ndtype)

Demonstrating analysis and visualisation of the output dataset, independently of OpenDrift code

if not os.path.exists('openoil.nc'):
    raise ValueError('Please run create_test_dataset.py first')

Importing a trajectory dataset from a simulation with OpenDrift decode_coords is needed so that lon and lat are not interpreted as coordinate variables

d = xr.open_dataset('openoil.nc', decode_coords=False)
# Requirement that status>=0 is needed since non-valid points are not masked in OpenDrift output
d = d.where(d.status>=0)  # only active particles

Displaying a basic plot of trajectories

d.traj.plot(land='mask')
ax = plt.gca()
ax.set_title('Adding custom title')
plt.show()
Adding custom title

Demonstrating how the Xarray Dataset can be modified, allowing for more flexibility than can be provided through the plotting method of OpenDrift

Extracting only the first 100 elements, and every 4th output time steps:

dsub = d.isel(trajectory=range(0, 100), time=range(0, len(d.time), 4))
dsub.traj.plot(land='h')
plt.show()
example trajan

With several plots on the same figure

d.traj.plot(color='red', alpha=0.01, land='mask')  # Plotting individual trajectories in red
dmean = d.mean('trajectory', skipna=True) # Overlaying a "mean" trajectory in black
dmean.traj.plot(color='k', linewidth=5)
# Showing the a sub-period of the mean trajectory in yellow
dmean.sel(time=slice('2015-11-17', '2015-11-17 12')).traj.plot(color='yellow', linewidth=5)
plt.tight_layout()
plt.show()
example trajan

Total running time of the script: (0 minutes 49.555 seconds)

Gallery generated by Sphinx-Gallery