Creating the initial conditions

Creating the initial conditions#

[1]:

from roms_tools import Grid, InitialConditions

We start by creating a grid that contains Iceland and has a horizontal resolution of 2 km and 100 vertical layers.

[2]:

grid = Grid(
    nx=500,
    ny=500,
    size_x=1000,
    size_y=1000,
    center_lon=-20,
    center_lat=65,
    rot=0,
    N=100,
    topography_source={
        "name": "SRTM15",
        "path": "/anvil/projects/x-ees250129/Datasets/SRTM15/SRTM15_V2.6.nc",
    },
)

Next, we specify the time that we want to make the initial conditions for.

[3]:

from datetime import datetime

[4]:

ini_time = datetime(2012, 1, 2, 12, 0, 0) # noon on January 2, 2012

Physical initial conditions from GLORYS#

In this section, we use GLORYS data to create our physical initial conditions, i.e., temperature, salinity, sea surface height, and velocities. (We will learn how to add biogeochemical initial conditions further down in the notebook.) Our downloaded GLORYS data sits at the following location.

[5]:

path = "/anvil/projects/x-ees250129/Datasets/GLORYS/NA/2012/mercatorglorys12v1_gl12_mean_20120102.nc"

Note

It would also be okay to provide a filename that contains data for more than just the day of interest. ROMS-Tools will pick out the correct day (and complain if the day of interest is not in the provided filename.) Or we can even use wildcards, such as path="/global/cfs/projectdirs/m4746/Datasets/GLORYS/NA/2012/*.nc". Note, however, that ROMS-Tools will operate more efficiently when the filename is as specific as possible.

You can also download your own GLORYS data, which must span the desired ROMS domain and ini_time. For more details, please refer to this page.

We can now create the InitialConditions object.

[6]:

%%time

initial_conditions = InitialConditions(
    grid=grid,
    ini_time=ini_time,
    source={"name": "GLORYS", "path": path},
    model_reference_date=datetime(2000, 1, 1), # this is the default
    use_dask=True,
)

CPU times: user 32 s, sys: 900 ms, total: 32.9 s
Wall time: 8.05 s

Note

In the cell above, we created our initial conditions with use_dask = True. This enables Dask, a Python library designated to facilitate scalable, out-of-memory data processing by distributing computations across multiple threads or processes. Here you can learn more about using Dask with ROMS-Tools.

The initial conditions variables are held in an xarray.Dataset that is accessible via the .ds property.

[7]:

initial_conditions.ds

[7]:

<xarray.Dataset> Size: 508MB
Dimensions:     (ocean_time: 1, s_rho: 100, eta_rho: 502, xi_rho: 502,
                 xi_u: 501, eta_v: 501, s_w: 101)
Coordinates:
  * ocean_time  (ocean_time) float64 8B 3.788e+08
    abs_time    (ocean_time) datetime64[ns] 8B 2012-01-02T12:00:00
Dimensions without coordinates: s_rho, eta_rho, xi_rho, xi_u, eta_v, s_w
Data variables:
    temp        (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    salt        (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    u           (ocean_time, s_rho, eta_rho, xi_u) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    v           (ocean_time, s_rho, eta_v, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    zeta        (ocean_time, eta_rho, xi_rho) float32 1MB -0.7428 ... -0.8718
    ubar        (ocean_time, eta_rho, xi_u) float32 1MB dask.array<chunksize=(1, 409, 409), meta=np.ndarray>
    vbar        (ocean_time, eta_v, xi_rho) float32 1MB dask.array<chunksize=(1, 409, 409), meta=np.ndarray>
    w           (ocean_time, s_w, eta_rho, xi_rho) float32 102MB 0.0 0.0 ... 0.0
    Cs_r        (s_rho) float32 400B -0.992 -0.9753 ... -8.89e-05 -9.874e-06
    Cs_w        (s_w) float32 404B -1.0 -0.9837 -0.9667 ... -3.95e-05 0.0
Attributes:
    title:                                ROMS initial conditions file create...
    roms_tools_version:                   10000.dev361+g2e6c47bd8
    ini_time:                             2012-01-02 12:00:00
    model_reference_date:                 2000-01-01 00:00:00
    adjust_depth_for_sea_surface_height:  False
    source:                               GLORYS
    theta_s:                              5.0
    theta_b:                              2.0
    hc:                                   300.0

Let’s make some plots! As an example, let’s have a look at the temperature field temp. It is three-dimensional with horizontal dimensions eta_rho and xi_rho, and vertical dimension s_rho.

We first want to plot different layers of the temperature field, i.e., slice along the vertical dimension s.

[8]:

initial_conditions.plot("temp", s=-1)  # plot uppermost layer

[########################################] | 100% Completed | 1.01 sms
_images/initial_conditions_18_1.png 
[9]:

initial_conditions.plot("temp", s=0, depth_contours=True)  # plot bottom layer
_images/initial_conditions_19_0.png

Next, we slice our domain along one of the horizontal dimensions and look at temperature along these sections.

[10]:

initial_conditions.plot("temp", xi=0, layer_contours=True)
_images/initial_conditions_21_0.png

Note that even though we have a total of 100 layers, layer_contours = True will create a plot with a maximum of 10 contours to ensure plot clarity.

[11]:

initial_conditions.plot("temp", eta=50, layer_contours=False)
_images/initial_conditions_23_0.png

We can also plot a depth profile at a certain spatial location.

[12]:

initial_conditions.plot("temp", eta=0, xi=0)
_images/initial_conditions_25_0.png

Finally, we can look at a transect in a certain layer and at a fixed eta_rho (similarly xi_rho).

[13]:

initial_conditions.plot("temp", eta=0, s=-1)
_images/initial_conditions_27_0.png

Plotting velocity fields works similarly.

[14]:

initial_conditions.plot("u", s=-1)  # plot uppermost layer

[########################################] | 100% Completed | 1.51 sms
_images/initial_conditions_29_1.png 
[15]:

initial_conditions.plot("u", eta=0)
_images/initial_conditions_30_0.png 
[16]:

initial_conditions.plot("u", xi=0)
_images/initial_conditions_31_0.png 
[17]:

initial_conditions.plot("v", xi=0)

[########################################] | 100% Completed | 1.51 sms
_images/initial_conditions_32_1.png

Adding Biogeochemical (BGC) Initial Conditions#

To prepare a ROMS simulation with MARBL biogeochemistry (BGC), we need to generate both physical and biogeochemical initial conditions. Since ROMS requires a single initial conditions file, we create both sets of conditions together.

We use the following datasets:

  • Physical Initial Conditions:
    Derived from GLORYS data, including temperature, salinity, sea surface height, and velocities.
  • Biogeochemical (BGC) Initial Conditions:
    A unified BGC climatology with 1° horizontal resolution, that combines multiple observationally and model based sources:

    • Nutrient concentrations (nitrate, phosphate, silicate) and dissolved oxygen from the 2018 World Ocean Atlas

    • Dissolved iron and nitrous oxide from in-situ measurements

    • Dissolved inorganic carbon (DIC) and alkalinity from the global GLODAPv2 product

    • Remaining nutrients (ammonium, nitrite, and organic nitrogen) and dissolved organic matter from CESM model simulations

The BGC climatology can be accessed at the following location:

[18]:

unified_bgc_path = "/anvil/projects/x-ees250129/Datasets/UNIFIED/BGCdataset_v2_filled.nc"

The initial conditions are created as above, but now with additional information about the bgc_source.

[19]:

%%time

initial_conditions_with_unified_bgc = InitialConditions(
    grid=grid,
    ini_time=ini_time,
    source={"name": "GLORYS", "path": path},
    bgc_source={
        "name": "UNIFIED",
        "path": unified_bgc_path,
        "climatology": True,
    },  # bgc_source is optional
    model_reference_date=datetime(2000, 1, 1), # this is the default
    use_dask=True,
)

/anvil/projects/x-ees250129/x-smaticka/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1716: UserWarning: rename 'lon' to 'longitude' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
/anvil/projects/x-ees250129/x-smaticka/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1716: UserWarning: rename 'lat' to 'latitude' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
/anvil/projects/x-ees250129/x-smaticka/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1716: UserWarning: rename 'dep' to 'depth' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
2026-07-09 14:09:58 - WARNING - Optional variables missing (but not critical): ['Lig', 'DIC_ALT_CO2', 'Alk_ALT_CO2']

CPU times: user 26.1 s, sys: 1.56 s, total: 27.7 s
Wall time: 10.5 s

[20]:

initial_conditions_with_unified_bgc.ds

[20]:

<xarray.Dataset> Size: 4GB
Dimensions:      (ocean_time: 1, s_rho: 100, eta_rho: 502, xi_rho: 502,
                  xi_u: 501, eta_v: 501, s_w: 101)
Coordinates:
    abs_time     (ocean_time) datetime64[ns] 8B 2012-01-02T12:00:00
  * ocean_time   (ocean_time) float64 8B 3.788e+08
Dimensions without coordinates: s_rho, eta_rho, xi_rho, xi_u, eta_v, s_w
Data variables: (12/42)
    temp         (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    salt         (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    u            (ocean_time, s_rho, eta_rho, xi_u) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    v            (ocean_time, s_rho, eta_v, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    zeta         (ocean_time, eta_rho, xi_rho) float32 1MB dask.array<chunksize=(1, 502, 502), meta=np.ndarray>
    ubar         (ocean_time, eta_rho, xi_u) float32 1MB dask.array<chunksize=(1, 409, 409), meta=np.ndarray>
    ...           ...
    Lig          (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    DIC_ALT_CO2  (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    ALK_ALT_CO2  (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    w            (ocean_time, s_w, eta_rho, xi_rho) float32 102MB 0.0 ... 0.0
    Cs_r         (s_rho) float32 400B -0.992 -0.9753 ... -8.89e-05 -9.874e-06
    Cs_w         (s_w) float32 404B -1.0 -0.9837 -0.9667 ... -3.95e-05 0.0
Attributes:
    title:                                ROMS initial conditions file create...
    roms_tools_version:                   4.0.0a2.dev29+gdcb35eb01
    ini_time:                             2012-01-02 12:00:00
    model_reference_date:                 2000-01-01 00:00:00
    adjust_depth_for_sea_surface_height:  False
    source:                               GLORYS
    bgc_source:                           UNIFIED
    theta_s:                              5.0
    theta_b:                              2.0
    hc:                                   300.0
[21]:

initial_conditions_with_unified_bgc.plot("PO4", xi=0, layer_contours=True)

[########################################] | 100% Completed | 1.11 sms
_images/initial_conditions_38_1.png 
[22]:

initial_conditions_with_unified_bgc.plot("ALK", eta=0)
_images/initial_conditions_39_0.png

Interpolate BGC variables by density or mixed-layer depth instead of plain depth#

By default (bgc_interpolation_method="depth"), BGC tracers are interpolated onto the model’s vertical grid in depth. Because biogeochemical fields tend to follow water masses rather than fixed depths, two alternatives are available:

  • bgc_interpolation_method="density" interpolates the tracers in potential-density (isopycnal) space. The density coordinate is built from the BGC dataset’s own temperature/salinity for the source profile and from the model’s temperature/salinity for the target.

  • bgc_interpolation_method="density_mld" finds the mixed layer depth (MLD) in the source and target density fields, scales the source mixed layer so its MLD matches the target’s, and then interpolates 1:1 in depth below the MLD. This keeps the mixed layers aligned and preserves the absolute depth of sub-mixed-layer features, while avoiding the surface degeneracy of pure density space.

Both density methods require the BGC source to provide temperature/salinity (e.g. the unified dataset’s temp_WOA/salt_WOA); otherwise interpolation falls back to depth space.

[23]:

%%time

initial_conditions_with_unified_bgc_mld = InitialConditions(
    grid=grid,
    ini_time=ini_time,
    source={"name": "GLORYS", "path": path},
    bgc_source={
        "name": "UNIFIED",
        "path": unified_bgc_path,
        "climatology": True,
    },  # bgc_source is optional
    model_reference_date=datetime(2000, 1, 1), # this is the default
    use_dask=True,
    bgc_interpolation_method="density_mld",
)

/home/x-kthyng/packages/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1539: UserWarning: rename 'lon' to 'longitude' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
/home/x-kthyng/packages/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1539: UserWarning: rename 'lat' to 'latitude' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
/home/x-kthyng/packages/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1539: UserWarning: rename 'dep' to 'depth' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
2026-06-26 19:39:37 - WARNING - Optional variables missing (but not critical): ['Lig', 'DIC_ALT_CO2', 'Alk_ALT_CO2']

CPU times: user 3min 28s, sys: 2.65 s, total: 3min 30s
Wall time: 11.4 s

[24]:

initial_conditions_with_unified_bgc_mld.plot("ALK", eta=0)
_images/initial_conditions_42_0.png

Initializing from ROMS Restart Files#

When running a nested ROMS simulation, it is necessary to initialize the child simulation using a restart file from a previously run parent simulation. For a detailed guide on preparing nested simulations, see Nested Simulations.

Typically, you will have already created both the parent and child grid objects, as described in the nesting guide linked above. To initialize the child simulation, we first create the corresponding grid objects for both simulations.

In this example, the parent grid covers the northeast Pacific at a horizontal resolution of 3km and with 60 vertical layers.

[25]:

parent_grid = Grid(filename="/anvil/projects/x-ees250129/Datasets/ROMSOutput/eastpac3km/eastpac_3km_grd.nc")

[26]:

parent_grid.plot()
_images/initial_conditions_45_0.png

We will use the following restart file, produced by the parent grid simulation for December 31, 2023, at 08:25:20, to initialize the child grid.

[27]:

restart_date = datetime(2023, 12, 31, 8, 25, 20)
restart_file = "/anvil/projects/x-ees250129/Datasets/ROMSOutput/eastpac3km/eastpac3km_rst.20231231082520.nc"

The child grid covers the Salish Sea at a horizontal resolution of 600m and 100 vertical layers. Here, we create a standard Grid that has not been adapted to the parent grid yet.

[28]:

child_grid = Grid(
    nx = 960,
    ny = 640,
    size_x = 576,
    size_y = 384,
    center_lon=-125.0,
    center_lat=49.2,
    rot=-30.0,
    N=100,
    hc=300.0,
)

[29]:

child_grid.plot()
_images/initial_conditions_50_0.png

Before we make the initial conditions for the child simulation, we need to smooth the child grid boundaries into the parent grid. See the nesting documentation.

[30]:

from roms_tools import align_grids
child_grid = align_grids(parent_grid=parent_grid, child_grid=child_grid, verbose=True)

2026-06-26 19:39:55 - INFO - No `boundaries` provided. Using mask-based defaults: {'south': True, 'east': False, 'north': False, 'west': True}
2026-06-26 19:43:00 - INFO - === Modifying the child mask ===
2026-06-26 19:43:33 - INFO - Total time: 33.032 seconds
2026-06-26 19:43:33 - INFO - ================================================================================================
2026-06-26 19:43:33 - INFO - === Modifying the child topography ===
2026-06-26 19:44:06 - INFO - Total time: 32.955 seconds
2026-06-26 19:44:06 - INFO - ================================================================================================

We now create the initial conditions for the child simulations as follows.

[31]:

%%time

initial_conditions_from_roms = InitialConditions(
    grid=child_grid,
    ini_time=restart_date,
    source={"name": "ROMS", "grid": parent_grid, "path": restart_file},
    use_dask=True,
    bgc_source={
        "name": "ROMS",
        "grid": parent_grid,
        "path": restart_file,
    },
)

CPU times: user 14min 12s, sys: 3.18 s, total: 14min 15s
Wall time: 46.8 s

Note

In the cell above, we created the initial conditions using both source and bgc_source with the same information (the ROMS parent grid and restart file). If you do not provide bgc_source, ROMS-Tools will only read and prepare the physical variables from the restart file.

[32]:

initial_conditions_from_roms.ds

[32]:

<xarray.Dataset> Size: 9GB
Dimensions:      (ocean_time: 1, eta_rho: 642, xi_rho: 962, s_rho: 100,
                  xi_u: 961, eta_v: 641, s_w: 101)
Coordinates:
  * ocean_time   (ocean_time) float64 8B 7.573e+08
    abs_time     (ocean_time) datetime64[ns] 8B 2023-12-31T08:25:20
Dimensions without coordinates: eta_rho, xi_rho, s_rho, xi_u, eta_v, s_w
Data variables: (12/42)
    zeta         (ocean_time, eta_rho, xi_rho) float32 2MB dask.array<chunksize=(1, 642, 962), meta=np.ndarray>
    temp         (ocean_time, s_rho, eta_rho, xi_rho) float32 247MB dask.array<chunksize=(1, 100, 642, 962), meta=np.ndarray>
    salt         (ocean_time, s_rho, eta_rho, xi_rho) float32 247MB dask.array<chunksize=(1, 100, 642, 962), meta=np.ndarray>
    u            (ocean_time, s_rho, eta_rho, xi_u) float32 247MB dask.array<chunksize=(1, 100, 642, 961), meta=np.ndarray>
    v            (ocean_time, s_rho, eta_v, xi_rho) float32 247MB dask.array<chunksize=(1, 100, 641, 962), meta=np.ndarray>
    ubar         (ocean_time, eta_rho, xi_u) float32 2MB dask.array<chunksize=(1, 642, 961), meta=np.ndarray>
    ...           ...
    diazFe       (ocean_time, s_rho, eta_rho, xi_rho) float32 247MB dask.array<chunksize=(1, 100, 642, 962), meta=np.ndarray>
    spCaCO3      (ocean_time, s_rho, eta_rho, xi_rho) float32 247MB dask.array<chunksize=(1, 100, 642, 962), meta=np.ndarray>
    zooC         (ocean_time, s_rho, eta_rho, xi_rho) float32 247MB dask.array<chunksize=(1, 100, 642, 962), meta=np.ndarray>
    w            (ocean_time, s_w, eta_rho, xi_rho) float32 250MB 0.0 ... 0.0
    Cs_r         (s_rho) float32 400B -0.992 -0.9753 ... -8.89e-05 -9.874e-06
    Cs_w         (s_w) float32 404B -1.0 -0.9837 -0.9667 ... -3.95e-05 0.0
Attributes:
    title:                                ROMS initial conditions file create...
    roms_tools_version:                   10000.dev361+g2e6c47bd8
    ini_time:                             2023-12-31 08:25:20
    model_reference_date:                 2000-01-01 00:00:00
    adjust_depth_for_sea_surface_height:  True
    source:                               ROMS
    bgc_source:                           ROMS
    theta_s:                              5.0
    theta_b:                              2.0
    hc:                                   300.0

Let’s make some plots as before.

[33]:

initial_conditions_from_roms.plot("temp", s=-1)

[########################################] | 100% Completed | 8.47 sms
_images/initial_conditions_58_1.png 
[34]:

initial_conditions_from_roms.plot("ALK", eta=0)
_images/initial_conditions_59_0.png

Strict vs. flexible initial time#

The daily GLORYS dataset provides one time stamp per day, which in our locally available dataset is always at 12:00 UTC (noon). Similarly, the restart files typically contain 1–2 time stamps that may be at irregular times (as seen above).

But what happens if we request an initial time that does not exactly match one of these available time stamps?

[35]:

midnight_ini_time = datetime(2012, 1, 2, 0, 0, 0) # midnight on (i.e., the start of) January 2, 2012

[36]:

initial_conditions_strict_midnight = InitialConditions(
    grid=grid,
    ini_time=midnight_ini_time,
    source={"name": "GLORYS", "path": path},
    use_dask=True,
)

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[40], line 1
----> 1 initial_conditions_strict_midnight = InitialConditions(
      2     grid=grid,
      3     ini_time=midnight_ini_time,
      4     source={"name": "GLORYS", "path": path},

File <string>:15, in __create_fn__.<locals>.__init__(self, grid, ini_time, source, bgc_source, model_reference_date, allow_flex_time, use_dask, chunks, initial_slice_bounds, bypass_validation, bgc_interpolation_method, use_density_interpolation)

File ~/packages/roms-tools/roms_tools/setup/initial_conditions.py:219, in InitialConditions.__post_init__(self)
    216 self._input_checks()
    218 processed_fields = {}
--> 219 processed_fields = self._process_data(processed_fields, type="physics")
    221 if self.bgc_source is not None:
    222     processed_fields = self._process_data(processed_fields, type="bgc")

File ~/packages/roms-tools/roms_tools/setup/initial_conditions.py:247, in InitialConditions._process_data(self, processed_fields, type)
    244 def _process_data(self, processed_fields, type="physics"):
    245     target_coords = get_target_coords(self.grid)
--> 247     data = self._get_data(forcing_type=type)
    248     data.choose_subdomain(
    249         target_coords,
    250         unchunk_lateral_dims=True,
    251     )
    252     # Enforce double precision to ensure reproducibility

File ~/packages/roms-tools/roms_tools/setup/initial_conditions.py:655, in InitialConditions._get_data(self, forcing_type)
    652 else:
    653     self.adjust_depth_for_sea_surface_height = False
--> 655     data = data_type(
    656         filename=source_dict["path"],  # type: ignore
    657         start_time=self.ini_time,
    658         climatology=source_dict["climatology"],  # type: ignore
    659         allow_flex_time=self.allow_flex_time,
    660         use_dask=self.use_dask,
    661         chunks=self.chunks,
    662         initial_slice_bounds=self.initial_slice_bounds,
    663     )
    665 return data

File <string>:19, in __create_fn__.<locals>.__init__(self, filename, start_time, end_time, dim_names, var_names, opt_var_names, climatology, has_encoded_times, needs_lateral_fill, use_dask, chunks, initial_slice_bounds, read_zarr, allow_flex_time, apply_post_processing, ds_loader_fn)

File ~/packages/roms-tools/roms_tools/datasets/lat_lon_datasets.py:189, in LatLonDataset.__post_init__(self)
    187 if "time" in self.dim_names and self.start_time is not None:
    188     ds = self.add_time_info(ds)
--> 189     ds = self.select_relevant_times(ds)
    191     if self.dim_names["time"] != "time":
    192         ds = ds.rename({self.dim_names["time"]: "time"})

File ~/packages/roms-tools/roms_tools/datasets/lat_lon_datasets.py:353, in LatLonDataset.select_relevant_times(self, ds)
    350 if self.start_time is None:
    351     raise ValueError("select_relevant_times called but start_time is None.")
--> 353 ds = select_relevant_times(
    354     ds=ds,
    355     time_dim=time_dim,
    356     time_coord=time_dim,
    357     start_time=self.start_time,
    358     end_time=self.end_time,
    359     climatology=self.climatology,
    360     allow_flex_time=self.allow_flex_time,
    361 )
    363 return ds

File ~/packages/roms-tools/roms_tools/datasets/utils.py:294, in select_relevant_times(ds, time_dim, time_coord, start_time, end_time, climatology, allow_flex_time)
    290     ds = ds.assign_coords({time_dim: convert_cftime_to_datetime(ds[time_coord])})
    292 if not end_time:
    293     # Assume we are looking for exactly one time record for initial conditions
--> 294     return _select_initial_time(
    295         ds, time_dim, time_coord, start_time, climatology, allow_flex_time
    296     )
    298 if climatology:
    299     return ds

File ~/packages/roms-tools/roms_tools/datasets/utils.py:397, in _select_initial_time(ds, time_dim, time_coord, ini_time, climatology, allow_flex_time)
    394 else:
    395     # Strict match required
    396     if not (ds[time_coord].values == np.datetime64(ini_time)).any():
--> 397         raise ValueError(
    398             f"No exact match found for initial time {ini_time}. Consider setting allow_flex_time to True."
    399         )
    401     ds = ds.sel({time_coord: np.datetime64(ini_time)})
    403 if time_dim not in ds.dims:

ValueError: No exact match found for initial time 2012-01-02 00:00:00. Consider setting allow_flex_time to True.

If you request an initial time that does not exactly match one of these entries, ROMS-Tools will by default raise an error because it cannot find an exact match. This behavior is controlled by the allow_flex_time flag:

  • False (default): requires an exact match to ini_time. A ValueError is raised if no match exists (see above example).

  • True: allows a +24-hour search window starting from ini_time and selects the closest available entry within that window. If no match is found, a ValueError is still raised.

[37]:

initial_conditions_flex_midnight = InitialConditions(
    grid=grid,
    ini_time=midnight_ini_time,
    allow_flex_time=True,
    source={"name": "GLORYS", "path": path},
    use_dask=True,
)

2026-06-26 19:46:48 - WARNING - Selected time entry closest to the specified start_time in +24 hour range: ['2012-01-02T12:00:00.000000000']

The logging output indicates which time stamp was actually selected. To verify, the following output confirms that ROMS-Tools adjusted the requested midnight_ini_time to the nearest available timestamp, which in this case was noon, 12 hours later.

[38]:

initial_conditions_flex_midnight.ds.abs_time

[38]:

<xarray.DataArray 'abs_time' (ocean_time: 1)> Size: 8B
array(['2012-01-02T12:00:00.000000000'], dtype='datetime64[ns]')
Coordinates:
  * ocean_time  (ocean_time) float64 8B 3.788e+08
    abs_time    (ocean_time) datetime64[ns] 8B 2012-01-02T12:00:00
Attributes:
    long_name:  absolute time

The same behavior applies when creating InitialConditions from a ROMS restart file. Instead of specifying the exact time December 31, 2023, at 08:25:20, we simply provide December 31, 2023, at 00:00 (midnight), together with allow_flex_time = True.

[39]:

%%time

initial_conditions_from_roms_flex_time = InitialConditions(
    grid=child_grid,
    ini_time=datetime(2023, 12, 31),
    allow_flex_time=True,
    source={"name": "ROMS", "grid": parent_grid, "path": restart_file},
    use_dask=True,
)

2026-06-26 19:46:52 - WARNING - Selected time entry closest to the specified start_time in +24 hour range: 2023-12-31T08:22:50.000000000

CPU times: user 1min 23s, sys: 1.51 s, total: 1min 24s
Wall time: 16.2 s

[40]:

initial_conditions_from_roms_flex_time.ds.abs_time

[40]:

<xarray.DataArray 'abs_time' (ocean_time: 1)> Size: 8B
array(['2023-12-31T08:22:50.000000000'], dtype='datetime64[ns]')
Coordinates:
  * ocean_time  (ocean_time) float64 8B 7.573e+08
    abs_time    (ocean_time) datetime64[ns] 8B 2023-12-31T08:22:50
Attributes:
    long_name:  absolute time

Saving as NetCDF or YAML file#

We can now save the dataset as a NetCDF file.

[41]:

filepath = "/anvil/scratch/x-smaticka/initial_conditions/my_initial_conditions.nc"

[42]:

%time initial_conditions_with_unified_bgc.save(filepath)

2026-07-09 14:10:50 - INFO - Writing the following NetCDF files:
/anvil/scratch/x-smaticka/initial_conditions/my_initial_conditions.nc

[########################################] | 100% Completed | 26.88 ss
CPU times: user 3min 8s, sys: 13 s, total: 3min 21s
Wall time: 27.5 s

[42]:

[PosixPath('/anvil/scratch/x-smaticka/initial_conditions/my_initial_conditions.nc')]

We can also export the parameters of our InitialConditions object to a YAML file.

[43]:

yaml_filepath = "/anvil/scratch/x-smaticka/initial_conditions/my_initial_conditions.yaml"

[44]:

initial_conditions_with_unified_bgc.to_yaml(yaml_filepath)

[45]:

# Open and read the YAML file
with open(yaml_filepath, "r") as file:
    file_contents = file.read()

# Print the contents
print(file_contents)

---
roms_tools_version: 4.0.0a2.dev29+gdcb35eb01
---
Grid:
  nx: 500
  ny: 500
  size_x: 1000
  size_y: 1000
  center_lon: -20
  center_lat: 65
  rot: 0
  N: 100
  theta_s: 5.0
  theta_b: 2.0
  hc: 300.0
  topography_source:
    name: SRTM15
    path: /anvil/projects/x-ees250129/Datasets/SRTM15/SRTM15_V2.6.nc
  mask_shapefile: null
  close_narrow_channels: false
  hmin: 5.0
  filename: null
InitialConditions:
  ini_time: '2012-01-02T12:00:00'
  source:
    name: GLORYS
    path: /anvil/projects/x-ees250129/Datasets/GLORYS/NA/2012/mercatorglorys12v1_gl12_mean_20120102.nc
    climatology: false
  bgc_source:
    name: UNIFIED
    path: /anvil/projects/x-ees250129/Datasets/UNIFIED/BGCdataset_v2_filled.nc
    climatology: true
  model_reference_date: '2000-01-01T00:00:00'
  allow_flex_time: false
  chunks: null
  initial_slice_bounds: null
  bypass_validation: false
  bgc_interpolation_method: depth

Creating initial conditions from an existing YAML file#

[46]:

%time the_same_initial_conditions_with_unified_bgc = InitialConditions.from_yaml(yaml_filepath, use_dask=True)

/anvil/projects/x-ees250129/x-smaticka/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1716: UserWarning: rename 'lon' to 'longitude' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
/anvil/projects/x-ees250129/x-smaticka/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1716: UserWarning: rename 'lat' to 'latitude' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
/anvil/projects/x-ees250129/x-smaticka/roms-tools/roms_tools/datasets/lat_lon_datasets.py:1716: UserWarning: rename 'dep' to 'depth' does not create an index anymore. Try using swap_dims instead or use set_index after rename to create an indexed coordinate.
  ds = ds.rename({"lon": "longitude", "lat": "latitude", "dep": "depth"})
2026-07-09 14:11:49 - WARNING - Optional variables missing (but not critical): ['Lig', 'DIC_ALT_CO2', 'Alk_ALT_CO2']

CPU times: user 2min 52s, sys: 1.6 s, total: 2min 54s
Wall time: 16.7 s

[47]:

the_same_initial_conditions_with_unified_bgc.ds

[47]:

<xarray.Dataset> Size: 4GB
Dimensions:      (ocean_time: 1, s_rho: 100, eta_rho: 502, xi_rho: 502,
                  xi_u: 501, eta_v: 501, s_w: 101)
Coordinates:
    abs_time     (ocean_time) datetime64[ns] 8B 2012-01-02T12:00:00
  * ocean_time   (ocean_time) float64 8B 3.788e+08
Dimensions without coordinates: s_rho, eta_rho, xi_rho, xi_u, eta_v, s_w
Data variables: (12/42)
    temp         (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    salt         (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    u            (ocean_time, s_rho, eta_rho, xi_u) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    v            (ocean_time, s_rho, eta_v, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    zeta         (ocean_time, eta_rho, xi_rho) float32 1MB dask.array<chunksize=(1, 502, 502), meta=np.ndarray>
    ubar         (ocean_time, eta_rho, xi_u) float32 1MB dask.array<chunksize=(1, 409, 409), meta=np.ndarray>
    ...           ...
    Lig          (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    DIC_ALT_CO2  (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    ALK_ALT_CO2  (ocean_time, s_rho, eta_rho, xi_rho) float32 101MB dask.array<chunksize=(1, 100, 409, 409), meta=np.ndarray>
    w            (ocean_time, s_w, eta_rho, xi_rho) float32 102MB 0.0 ... 0.0
    Cs_r         (s_rho) float32 400B -0.992 -0.9753 ... -8.89e-05 -9.874e-06
    Cs_w         (s_w) float32 404B -1.0 -0.9837 -0.9667 ... -3.95e-05 0.0
Attributes:
    title:                                ROMS initial conditions file create...
    roms_tools_version:                   4.0.0a2.dev29+gdcb35eb01
    ini_time:                             2012-01-02 12:00:00
    model_reference_date:                 2000-01-01 00:00:00
    adjust_depth_for_sea_surface_height:  False
    source:                               GLORYS
    bgc_source:                           UNIFIED
    theta_s:                              5.0
    theta_b:                              2.0
    hc:                                   300.0
[ ]: