Add empirical calculations to intensity measures - #127
Conversation
![gemini-code-assist[bot]](https://avatars.githubusercontent.com/in/956858?s=60&v=4){kind=small}
There was a problem hiding this comment.
Code Review
This pull request introduces the capability to calculate empirical intensity measures from ground motion models using OpenQuake wrappers. It updates dependencies, default parameters, and schemas, and adds helper methods to average multi-fault parameters (rakes, magnitudes) weighted by fault moment. The im_calc.py script is refactored to calculate source-to-site distances, site parameters, and empirical IMs, outputting the results as an xr.DataTree in NetCDF format. The review feedback highlights critical issues where running the script with empirical=False or with broadband data lacking vs30 coordinates will result in NameError or AttributeError exceptions.
Important
The consumer version of Gemini Code Assist on GitHub is being sunset. Starting June 18, 2026, new organization installations will be blocked, and all code review activity will officially cease on July 17, 2026.
For more details on the timeline and next steps, please review the Help Documentation.
There was a problem hiding this comment.
Pull request overview
Refactors intensity-measure (IM) calculation output to an xarray.DataTree netCDF structure and adds optional empirical (GMM-based) IM calculations driven by per-realisation configuration (tectonic type + model list), carrying through per-station parameters like Vs30 and derived basin depths.
Changes:
- Refactor
im-calcoutput from a singlexarray.Datasetto a structuredxarray.DataTree, with per-IM datasets and shared station/source metadata. - Add empirical IM calculation via
oq_wrapperfor supported IMs/tectonic types and store results under{im}/empirical/{model}. - Add
vs30to broadband waveform outputs and introduce realisation/schema/default support for empirical configuration.
Reviewed changes
Copilot reviewed 6 out of 7 changed files in this pull request and generated 4 comments.
Show a summary per file
| File | Description |
|---|---|
| workflow/scripts/im_calc.py | Major refactor to DataTree output; adds empirical GMM evaluation and metadata/unit annotation. |
| workflow/scripts/bb_sim.py | Adds vs30 coordinate to broadband netCDF output for downstream site parameter usage. |
| workflow/schemas.py | Adds EMPIRICAL_PARAMETERS schema for validating empirical config in realisations/defaults. |
| workflow/realisations.py | Adds moment-weighted averaging helpers and EmpiricalParameters realisation configuration. |
| workflow/default_parameters/root/defaults.yaml | Introduces default empirical configuration (active_shallow + NSHM2022). |
| uv.lock | Lockfile change related to dependency resolution (cffi/pycparser marker). |
| pyproject.toml | Adds netCDF4 dependency to enforce safe import order vs OpenQuake/HDF5 stack. |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
|
|
||
| Ground motion models describe a rupture with a single magnitude, rake, | ||
| dip and depth. Multi-fault realisations are collapsed into these by | ||
| averaging each fault's contribution, weighted by its moment. |
There was a problem hiding this comment.
Based on what we have done before for empirical results for multi fault events we averaged based on the area of the fault not the moment. There was a few discussions regarding what to use and we landed on using area, however could be open to reconsider this adjustment if there has been a new development for this.
There was a problem hiding this comment.
Ah didn't see this, I will change the calculation to use that averaging.
| The site parameters, with basin depths estimated using the Chiou | ||
| and Youngs (2008) relations. | ||
| """ | ||
| z1pt0 = chiou_young_08_calc_z1p0(vs30) # ty: ignore[invalid-argument-type] |
There was a problem hiding this comment.
Do we want to estimate for all of the sites? Wondering if we have any "real sites" in this where we have a better measured / estimate value of z1.0 / 2.5 etc to use those instead from the site database.
There was a problem hiding this comment.
Yeah this is a stand-in at the moment, I intend to replace this with actual site database values at a later date.
5508803 to
6c4b53c
Compare
![github-actions[bot]](https://avatars.githubusercontent.com/in/15368?s=60&v=4){kind=small}
| @@ -14,12 +14,13 @@ dependencies = [ | |||
| "im-calculation @ git+https://github.com/ucgmsim/IM_calculation@no_parallel", | |||
There was a problem hiding this comment.
| "im-calculation @ git+https://github.com/ucgmsim/IM_calculation@no_parallel", | |
| "im-calculation @ git+https://github.com/ucgmsim/IM_calculation.git@no_parallel", |
| @@ -14,12 +14,13 @@ dependencies = [ | |||
| "im-calculation @ git+https://github.com/ucgmsim/IM_calculation@no_parallel", | |||
There was a problem hiding this comment.
| "im-calculation @ git+https://github.com/ucgmsim/IM_calculation@no_parallel", | |
| "im-calculation @ git+https://github.com/ucgmsim/IM_calculation.git@no_parallel", |
6c4b53c to
74f7c58
Compare
6e72aee to
ed8468f
Compare
786f419 to
d4a7ec7
Compare
d4a7ec7 to
e3d3997
Compare
e3d3997 to
e55d9fc
Compare
`im-calc` now also evaluates the ground motion models named in the
realisation, writing them alongside the simulated measures at
`{im}/empirical/{model}` in the same data tree. Having both in one file
is the point: comparing a simulation against a GMM should not require
joining two files on station and period.
A multi-fault realisation has to be reduced to the single rupture a GMM
describes -- one magnitude, rake, dip, ztor, zbot -- which
`calculate_source_parameters` does by moment-weighted averaging, with
rakes and dips averaged as unit vectors rather than as angles. Site
parameters come from vs30 via the Chiou & Young z1pt0/z2pt5
estimators, and the distances are the ones already computed for the
simulated measures, so both sides see identical geometry.
Combinations a model does not support are skipped with a warning rather
than failing the run -- the NSHM2022 logic tree configures pSA but not
PGA or PGV, and that is not an error.
Two deliberate awkwardnesses:
`netCDF4` is imported before anything that pulls in OpenQuake. OpenQuake
brings h5py, which is linked against a different HDF5 build than netCDF4
is, and whichever loads second cannot open our waveform files. The
import order is load bearing, not stylistic.
`EmpiricalParameters` types its fields as `Any` and its schema validates
plain strings rather than `oq_wrapper.constants` enum members. Importing
that module pulls in OpenQuake, which is expensive and must be
precompiled; the strings are resolved to enum members inside the
calculation instead, so `workflow.realisations` stays cheap to import.
vs30 is loaded eagerly rather than left as a dask array -- it is one
float per station, and letting it stay lazy propagates dask into every
station coordinate.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
e55d9fc to
30ced21
Compare
Adds empirical calculations to intensity measure outputs where empirical models support a given tectonic type and intensity measure.
To support this structure, I have also refactored im calc so that it finally outputs in xarray data tree format instead of dataset format. The upshot of doing this is that we can uncouple the components so that each intensity measure carries only components it actually computes. Especially for models like NSHM2022 that only support pSA it's kind of silly to have every intensity measure carry all-NaN EAS and empirical model components.