How to use OceanVal in Python

Validating a simulation involves three steps: register the observational datasets, matchup the simulation output against them, and generate statistics and an HTML report. Always create a new directory before running OceanVal for a new simulation, and run all commands from within it. Once you have reports for multiple simulations, you can also compare them side by side. To validate without writing any Python, use the OceanVal app instead.

Quickstart

The shortest path from a model output directory to a validation report is three function calls. Work from a fresh, empty directory: OceanVal writes its matchup files and report there.

import oceanval

# 1. Register: compare the model variable "thetao" with WOA23 temperature
oceanval.add_gridded_comparison(
    name="temperature",
    model_variable="thetao",
    recipe={"temperature": "woa23"},
    start=2005,
    end=2014,
    climatology=True,
)

# 2. Match: pair model output with the observations
oceanval.matchup(
    sim_dir="/path/to/model/output",
    start=2005,
    end=2014,
    cores=4,
)

# 3. Report: compute statistics and build the HTML report
oceanval.validate()

Replace thetao with the temperature variable name used in your model's NetCDF files. matchup scans sim_dir, reports the file pattern it has identified, and asks you to confirm before matching; the matched data is written to an oceanval_matchups directory. The report is written below oceanval_report and opens in your browser when the build completes. Each step is described in full below.

If no matchups are produced, check the model directory structure, variable names, time resolution, units, and climatology setting. The most common issue with monthly model output and in-situ observations is using daily matching precision; see the time resolution guidance below, or browse the Q&A.

Step 1: Register observational datasets

Register the observational datasets you want to validate against by specifying the location of the data files and any necessary metadata, using oceanval.add_point_comparison and oceanval.add_gridded_comparison.

Setting up gridded observational data

Register a gridded (NetCDF) dataset by pointing OceanVal at a directory of observation files and saying whether the data is a climatology:

oceanval.add_gridded_comparison(
    name="oxygen",
    source="CMEMS",
    source_info="Gridded observations from the Copernicus Marine Environment Monitoring Service",
    short_name="oxygen concentration",
    model_variable="oxygen",
    obs_variable="O2_concentration",
    obs_path="/path/to/obs_data/",
    climatology=False,
)

Required parameters:

• name required

  A name for the dataset, e.g. "temperature".

• source required

  The source of the observational data, e.g. "CMEMS".

• model_variable required

  The name of the model variable to compare against the observations.

• obs_path required

  Path to the directory containing the observational data files.

• obs_variable

  Name of the variable in the observational files. If omitted, OceanVal assumes a single variable is present.

• climatology required

  Whether the observational data is a climatology.

Using built-in recipes for gridded data

Built-in recipe definitions make it easy to register standard observational datasets — especially useful for a standard climatology such as WOA23 or GLODAP, without manually specifying all the metadata.

oceanval.add_gridded_comparison(
    name="temperature",
    source="WOA23",
    model_variable="temp",
    recipe={"temperature": "woa23"},
    start=2005, end=2014,
    climatology=True,
)

This uses the built-in metadata and file locations for the WOA23 temperature climatology. recipe also works for salinity, oxygen, nitrate, phosphate, silicate, chlorophyll and pH — see all recipes. The underlying helper can also be used directly:

recipe = oceanval.parsers.find_recipe({"temperature": "woa23"}, start=2005, end=2014)
print(recipe["source"], recipe["obs_variable"])

Browse the full recipe catalogue to see every supported observational dataset, with the exact call to use for each.

Setting up point (in-situ) observational data

Register an in-situ dataset by pointing OceanVal at a CSV (or directory of CSVs) and naming the model variable to compare it against:

oceanval.add_point_comparison(
    name="nitrate",
    source="ICES",
    source_info="In-situ observations from the International Council for the Exploration of the Sea",
    short_name="nitrate concentration",
    model_variable="temp",
    obs_path="/path/to/obs_data/",
)

Required parameters:

• name required

  A name for the dataset, e.g. "temperature" — used internally to keep track of things. Letters and numbers only.

• source required

  The source of the observational data, e.g. "NOAA".

• model_variable required

  The name of the model variable to compare against the observations.

• obs_path required

  Path to a file or directory containing the observational data. If a directory, ensure it only contains files relevant to this variable — OceanVal recursively uses all NetCDF files inside it.

Instead of obs_path, you can use a built-in point recipe, which downloads the observations for you when oceanval.matchup runs. Ten variables from the ICES Oceanographic database are available this way — see the point data recipes table for the full list and units. The observations are downloaded for the years being matched and the lon_lim/lat_lim area given to matchup (or the model's own extent, if these are not given). The name, source and labels all come from the recipe:

oceanval.add_point_comparison(
    model_variable="votemper",
    recipe={"temperature": "ices"},
    vertical=True,
)

See point data recipes for details.

Validating a variable against more than one dataset

Call oceanval.add_gridded_comparison or oceanval.add_point_comparison once per dataset, with the same name and a different source. Each source is matched up and reported separately:

oceanval.add_gridded_comparison(
    name="temperature",
    model_variable="votemper",
    recipe={"temperature": "cobe2"},
)
oceanval.add_gridded_comparison(
    name="temperature",
    model_variable="votemper",
    recipe={"temperature": "nsbc"},
)

Every call for a variable must use the same model_variable and labels. Registering a source that is already registered for that variable replaces it, so re-running a script or notebook cell is safe. Where a variable has more than one source, the summary page labels it with the source, e.g. "Temperature (COBE2)".

Step 2: Matchup model output with observations

Once your datasets are registered, pair the model output against them with oceanval.matchup:

oceanval.matchup(
    sim_dir="/path/to/simulation/output/",
    start=2000, end=2010,
    cores=4,
    lon_lim=[-80, 0], lat_lim=[20, 60],
    thickness="cell_thickness",
)

Required parameters:

• sim_dir required

  Path to the directory containing the model simulation output files.

• start / end required

  First / last year of the simulation to use for validation.

• cores

  Number of CPU cores to use for parallel processing.

• lon_lim / lat_lim

  Longitude / latitude limits for the validation region, e.g. [-180, 180].

• thickness vertical only

  "z_level" or the exact name of the variable holding the cell thickness, e.g. "e3t" (a variable that only has it in its name, such as O3_c_e3t, is not used) — required for vertical validation. The thickness can have a time dimension or not (a static thickness such as e3t_0), and its vertical coordinate needs a positive attribute ("down" for NEMO). A thickness file must hold only the thickness variable. With "z_level", the model’s levels must be depths (positive down); output on heights (positive up) needs the cell thickness instead.

Note: for monthly-resolution simulations, set point_time_res=["year", "month"] when matching up in-situ observations — otherwise day-of-year will almost never match and few matchups will be found. matchup shows the time resolution of the files it uses for each variable (monthly, or 1d, 5d and so on), and if a point dataset is matched by day against output coarser than daily, it warns you and offers to change point_time_res for all such datasets or for each one before anything is matched up.

Summing simulation output: to compare observations against the sum of multiple model variables, set something like "var1+var2+var3" as model_variable. In the create_recipes window you can tick the variables in the pop-out under the model variable box instead of typing the sum; the variables selected are listed with their long names at the top of the pop-out.

Step 3: Build the validation report

Once matched, compute statistics and generate plots with oceanval.validate, run from the same directory as the matchup step:

oceanval.validate()

This generates and opens an HTML page you can view in a web browser.

Validating along a transect

To see how well the model reproduces the observations along a line, give validate a transect. It must run north–south, with the same longitude at both ends, or east–west, with the same latitude at both ends:

oceanval.validate(transect={"start": [-30, 0], "end": [-30, 65]})

Each gridded matchup’s page then has a section with:

• a map of the transect over the land;

• the monthly climatology of the surface model and observations along the transect, and their difference, by latitude (or longitude) and month — or, if the matchup has a single time step, the values along the line;

• if the matchup was made through the water column (vertical=True), a section of the annual mean, with depth down and latitude (or longitude) across. The matchup’s depths can be unevenly spaced, so they are first interpolated onto 30 evenly spaced depths between the shallowest and the deepest.

The values are extracted with nctoolkit’s to_transect, bilinearly interpolated onto about one point per grid cell of the matchup, so it needs nctoolkit 1.3.6 or later. The model and observations share a viridis colour scale, and the difference is blue–white–red. In the oceanval window, the last page asks whether you want a transect, and does not let you carry on until it is a straight north–south or east–west line.

Choosing the depth bins

Point matchups made through the water column (vertical=True) are summarised by depth, on each point dataset’s page and in the summary. By default the depth ranges are 0–10, 10–30, 30–60, 60–100, 100–150, 150–300, 300–600, 600–1000 and >1000 m. Give validate your own as a list of [min, max] pairs, in metres:

oceanval.validate(depth_bins=[[0, 20], [20, 200], [200, None]])

A depth is in a bin if it is deeper than its min and no deeper than its max (the shallowest bin also takes its min). The deepest bin can have a max of None, for everything below its min. Bins can leave gaps, and observations in a gap are left out, but they cannot overlap. The same bins are used for every point dataset and in the summary. In the oceanval window, the last page asks for them when a point dataset is matched up with Vertical ticked, starting from the defaults: each bin can be changed or removed with its ×, and Add a bin adds another.

Validating some of the matchups

By default the report validates every matchup in oceanval_matchups. To validate only some, give validate their variables, each one for all its sources, or "variable/source" for one source:

oceanval.validate(variables=["temperature", "nitrate/WOA23"])

The names are those of the matchups (the folders and files in oceanval_matchups), and are not case sensitive. Each entry must match at least one matchup, or validate says which matchups there are. The summary covers only the matchups validated, and the others stay in oceanval_matchups for a later report. The oceanval window lists them all, ticked, for you to untick any to leave out.

Matchups from earlier runs

matchup never removes matchups from out_dir: a run that matches up other variables or sources leaves the earlier run’s matchups there, and by default validate validates them too. matchup records which matchups each run made (or kept, with overwrite=False), so previous=False validates only those of the last run, and a list adds the earlier matchups it picks:

oceanval.validate(previous=False)
oceanval.validate(previous=["chlorophyll/OCCCI"])

variables then picks among the matchups validated. The interim report only ever has the matchups of the run building it. When the output directory has matchups already, the OceanVal app lists them, by variable, source and whether they are point or gridded, and asks What do you want to do with existing matchups? Each has two options to tick, one of them: Use (as each is to start with) keeps it as it is, without making it again, and validates it; Replace/Ignore leaves it unused, so it is made again if this run matches it up and left out of the report if not. Ticking an option at the top of its column chooses it for them all. Nothing is removed from the directory. With no report to build, the second option is just Replace, and those used are kept and used as they are. Either way, you must tick to confirm what is done with them before carrying on.

Using your own regions

The report summarises model performance in each region: a map of the regions, the seasonal cycle of the spatial average, and a table of bias, RMSD and correlations. Instead of the built-in "global" or "nwes" regions, you can supply your own netCDF file:

oceanval.validate(subregions="my_regions.nc")

The file should be set up as follows:

• One variable per region, each a 2D map on a longitude/latitude grid. It is regridded to the model grid using nearest neighbour, so it does not need to match the model grid.

• Grid cells in the region are 1, and all other cells are missing values. A 0 is also accepted: it is set to missing and ignored in the regional summaries. Any other value is rejected.

• Each region is named in the report by the variable's long_name attribute, or by the variable name if it has no long_name.

A copy of the file is saved with the results, so oceanval.compare can also compare regions across simulations. When comparing, use the same regions file for every simulation. OceanVal cannot tell two different custom files apart.

Word and PDF versions of the report

The report is HTML by default. If you need something to edit, circulate or drop into a deliverable, ask for a Word version as well:

oceanval.validate(word=True)

This writes oceanval_report.docx — the whole report as a single document — plus one .docx per page, alongside the HTML files. The report itself gains download buttons for them, so you can pick up the Word file straight from the page you are reading.

The Word version is built to be worked on rather than just read:

• Equations are real, editable Word equations rather than pictures, so you can lift a formula straight into a paper or report.

• Figures and tables are centred, the page is A4, and every page carries a page number.

• In the combined document each page becomes a numbered chapter, so figures and tables read as 1.1, 1.2, 2.1 and so on, and cross-references in the text follow the same numbering. The single-page files number their own figures from 1, matching the HTML page they came from.

Word output needs pandoc installed (conda install -c conda-forge pandoc). If it is missing, OceanVal warns and still builds the HTML report.

PDF works the same way, via pdf=True, and the two can be combined. Both are off by default because they add time to the build and most people only need the HTML.

Looking at the results while the matchups run

Matching up a long simulation can take hours. oceanval.matchup can build an interim validation report as it goes, so that you can look at each matchup’s results as soon as it is made:

oceanval.matchup(
    sim_dir="/path/to/simulation/output/",
    start=2000, end=2010,
    live_validation={"lon_lim": [-20, 10], "lat_lim": [40, 65]},
)
oceanval.validate(lon_lim=[-20, 10], lat_lim=[40, 65], pdf=True)

Each matchup’s page is added to the interim report as soon as the matchup is made, and the summary is run again each time. live_validation takes True, for validate’s defaults, or a dict of validate’s report options: out_dir, lon_lim, lat_lim, subregions, fixed_scale, concise, transect, depth_bins and variables (the matchups it leaves out are still made, but have no page). The interim report is HTML only, so PDF, Word and zip versions are left to validate.

It is built in oceanval_interim_report, in out_dir: open oceanval_interim_report/oceanval_report/_build/html/notebooks/summary.html, and reload it to see the latest. Each page says how many of the matchups it has so far. matchup waits for the interim report to be finished before it returns, and the full report is then built by validate, as usual. The OceanVal app does all of this for you.

Comparing validation outputs from multiple simulations

Once you have validation reports for different simulations, compare them with oceanval.compare, which summarises the differences between them:

oceanval.compare(
    model_dict={
        "model_a": "/path/to/model_a",
        "model_b": "/path/to/model_b",
    },
    view=True,
    ask=True,
)

model_dict maps a short name for each model to its validation output directory. The HTML comparison report is written to oceanval_comparison/compare/_build/html/notebooks/comparison_seasonal.html.

Comparison reports take word=True and pdf=True as well, producing the same Word and PDF versions as a validation report.

Analyzing a simulation on its own

oceanval.analyze works out climatologies, phenology and time series from a simulation’s output alone, with no observations, and builds an HTML report of them the way validate builds the validation report:

oceanval.analyze(
    sim_dir="/path/to/simulation",
    variables={"Chlorophyll": "P1_Chl+P2_Chl+P3_Chl+P4_Chl", "Temperature": "votemper"},
    climatology={
        "Temperature": {"start": 2000, "end": 2009, "surface": True, "water_column": "mean"},
        "Chlorophyll": {"start": 2000, "end": 2009, "surface": True, "water_column": "integral"},
    },
    phenology={"Chlorophyll": {"start": 2000, "end": 2009}},
    long_term={"Temperature": {"start": 2000, "end": 2009, "surface": True}},
    seasonal={"Temperature": {"start": 2000, "end": 2009, "surface": True}},
    map_changes={
        "Temperature": {
            "historical": {"start": 1990, "end": 1999},
            "future": [{"start": 2040, "end": 2049}, {"start": 2090, "end": 2099}],
            "plot": "absolute_changes", "change": "percentage", "surface": True,
        },
    },
    transects={
        "Temperature": [
            {"start": [-30, 0], "end": [-30, 65],
             "climatology": {"start": 2000, "end": 2009},
             "vertical": {"start": 2000, "end": 2009},
             "trends": {"start": 2000, "end": 2009, "how": "absolute"}},
            {"start": [-20, 50], "end": [0, 50],
             "trends": {"start": 2000, "end": 2009, "how": "percentage_change"}},
        ],
    },
    stations={"L4": [-4.217, 50.25]},
    subregions="nwes",
    n_dirs_down=2,
    thickness="e3t",
    cores=4,  # there is no default: say how many
)

• Climatology: the mean of each month of the year over the years from start to end, and the annual mean, both mapped, at the surface, at the bottom ("bottom": True: the deepest value that is not missing, from nctoolkit’s bottom method with choice="value") and/or through the water column, as its vertical mean (weighted by the cell thickness) or its depth integral (a concentration per m³ becomes an amount per m²). The water column needs thickness: the thickness variable in the output, a file with it, or "z_level".
• Phenology: the day of the year of the seasonal peak at the surface in each grid cell, from the variable’s daily climatology (its mean on each day of the year), so it needs daily output. Days are counted in the model’s own calendar.
• Long-term time series (long_term): the mean of each year from start to end (the mean of its months’ means), over the model’s domain, over each region of subregions ("nwes", "global" or your own regions .nc file, as in validate), and at each of stations, the nearest model cell with values, within two grid cells. lon_lim and lat_lim cut a box out of the domain for the domain and region means. A year with output missing is the mean of the months there are, and the report marks it.
• Seasonal cycles (seasonal): the mean seasonal cycle: the mean of each month of the year from start to end, over the domain, each region and each station. Each year’s months are not drawn, but are saved in the results. It needs monthly output, or output more often.
• Map temporal changes (map_changes): maps of how a variable changes between a historical period and one or two future periods. For each variable, historical and future (a list of one or two) are dicts of start and end years, and the layers are surface, bottom and/or water_column, as for climatology. Each period is mapped as the mean of its years’ annual means. plot is what is drawn: "absolute" (the historical and future maps), "absolute_changes" (those, and the maps of the changes) or "changes" (the historical map, and the changes), and change is "absolute" (future minus historical) or "percentage" (that as a percentage of the historical value, missing where it is 0). The periods share one colour scale, and the changes another, symmetric about zero. All the maps of a layer are drawn in one figure, each titled by the years it is of. Each period must be years the output has, but the years between them need not be.
• Transects (transects): the model’s values along straight lines that run north-south (the same longitude at both ends) or east-west (the same latitude at both ends), as many as you like, for each variable on its own. Each variable has a list of transects, each a dict of start and end ([lon, lat]) and at least one of climatology (the monthly climatology at the surface), trends (at the surface, a heat map of the position along the transect by year, with how of "absolute", the annual means, "absolute_change", each annual mean minus the first, or "percentage_change", that as a percentage of the first), vertical (the annual mean through the water column, in a section; it needs thickness, which places each cell at its depth, or "z_level") and vertical_trends (the same how, as heat maps at some of the section’s depths), each a dict of its own start and end years; trends and vertical_trends can also have a window, the number of years of a rolling mean of the annual means, worked out with nctoolkit’s rolling_mean after the annual means (default 1: no rolling mean; with a window of N each year is the last of N years, and the first rolling mean is what the changes are from). Each month of each year is averaged first, a year’s mean is the mean of its months’ means, and years without all their months are left out of the trends. The changes need output in at least 2 years (plus the window less one); the app does not offer them for a variable with only one year. The report maps each transect and draws each option. It needs nctoolkit 1.3.6 or later.
• cores has no default: you must say how many cores OceanVal may use. as_missing is a value, or [minimum, maximum], to treat as missing, as in matchup. All spatial means use nctoolkit’s spatial_mean. If nctoolkit warns that the model’s grid has no cell bounds, so that the cell areas cannot be worked out, the mean is worked out again after regridding the grid to a regular lon-lat grid with to_latlon, by nearest neighbour and by bilinear interpolation (the bilinear values are used, and the report says how far the two differ), and the report says that the raw grid did not have enough information to calculate grid cell areas. The means of the regions in subregions are not worked out that way: their files are lon-lat grids already, so the raw model is regridded (bilinear) straight to the region file’s grid, and each region’s mean is taken there. The domain as a whole is called “Full Domain”: it is the mean over the whole model grid, so a region of that name in your regions file (as in the nwes file) is left out rather than shown as a second one. In the time series figures, with more than six regions, the full domain has a panel of its own, rather than a grey line in each region’s.
• figure_code=True gives each figure in the report a Copy code to reproduce figure button at its top right. It copies the Python code that draws the figure from the analysis’s saved results, with its colour scale, labels and colours written out, so you can change the figure and run the code in the conda environment OceanVal ran in. Default False.

Everything is saved in oceanval_analysis/ in out_dir, replacing any analysis already there: the fields as netCDF in results/, with analysis.json saying what is in them, the report in report/, and report.html, a link to its first page.

The report is built a variable at a time: each variable’s page is built, and the overview run again, as soon as the variable is worked out, while the next one is, so the report can be read (from report/_build/html/notebooks/000_overview.html) long before the analysis ends. Until then each page has a note saying how many of the variables are in it. Once the last is in, the note is taken off and report.html is made, with nothing left to build; if the analysis stops with an error, the pages of the variables already worked out are kept.

In the OceanVal app, Analyze simulations on the first page does the same with no Python. The years given must be ones the output has: analyze() refuses a start before the output’s first year or an end after its last, and the app says so beside the box.

Comparing simulations side by side

oceanval.compare_simulations does everything analyze does for two to five simulations, and builds one HTML report with them side by side: the maps are a panel for each simulation on one colour scale, and the time series have a colour for each simulation:

oceanval.compare_simulations(
    simulations=[
        {"name": "Control", "sim_dir": "/path/to/control", "n_dirs_down": 2},
        {"name": "New mixing", "sim_dir": "/path/to/new_mixing", "n_dirs_down": 1},
    ],
    # the model variable holding each variable in each simulation
    variables={"Temperature": {"Control": "votemper", "New mixing": "thetao"}},
    climatology={"Temperature": {"start": 2000, "end": 2009, "surface": True}},
    long_term={"Temperature": {"start": 2000, "end": 2009, "surface": True}},
    unify={"baseline": "Control"},
    cores=4,
)

• simulations are two to five dicts, each with a name, a sim_dir and n_dirs_down, and, if they differ between simulations, exclude, require, time_shift, thickness and as_missing. variables gives the model variable holding each variable in each simulation; a variable must be in the same units in all of them.
• climatology, phenology, long_term, seasonal, map_changes, transects, stations, subregions, lon_lim, lat_lim, thickness, as_missing, time_shift, figure_code and cores are as in analyze, and the same for every simulation. The years must be years every simulation has.
• unify says what to do if the simulations are on different grids. None leaves each on its own. {"lon": [-20, 10], "lat": [40, 65], "res": 0.1} puts all of them on a lon-lat grid of that box and resolution, with nctoolkit’s to_latlon. {"baseline": "Control"} regrids the others to the grid of that simulation, with nctoolkit’s regrid. Either way, a cell with no value in any one simulation, judged from its first time step, is made missing in all of them, so every simulation is compared over exactly the same cells.

Everything is saved in oceanval_simulation_comparison/ in out_dir: results/ has a folder for each simulation, laid out as an analysis is, and comparison.json, which says what is in them; report/ is the report’s book, and report.html a link to its first page. As for analyze, the report is built a variable at a time: each variable is worked out in every simulation, then its page is built while the next is. The app’s Compare simulations box does the same, with no Python.