Generate ShakeMap GMFs from `shake_result.hdf` instead of `grid.xml` and `uncertainty.xml`

[raoanirudh]

Motivation

The engine currently obtains ShakeMap ground motions primarily from grid.xml and uncertainty.xml. The resulting data contain the conditioned median and conditional total standard deviation at each grid point.

The Silva&Horspool implementation then constructs a new covariance matrix by applying spatial and cross-IMT correlation models to those total standard deviations. In particular, the spatial covariance is constructed from shakemap['std'] in hazardlib/shakemap/gmfs.py.

This does not distinguish between within-event and between-event residual components. Consequently, a spatial correlation model calibrated for within-event residuals (Jayaram and Baker, 2009) is effectively applied to the total uncertainty, while the spatially common between-event contribution is not sampled separately. A hardcoded cross-IMT correlation model is also applied to the total residuals.

This was a reasonable approximation when only total ShakeMap uncertainty was available through uncertainty.xml, but modern ShakeMap products provide the information needed for uncertainty component-aware, ShakeMap-consistent realizations.

As opposed to grid.xml and uncertainty.xml, the shake_result.hdf product contains several additional realization bookkeeping for each IMT, including:

• C
• Sigma_HH_YD
• add_uncertainty
• sta_per_ix
• sta_phi
• station coordinates

These arrays were extended to all IMTs in shakemap-modules 1.1.4 and are described in the source as the information used to reconstruct realizations. See:

• ShakeMap output documentation
• ShakeMap realization-data change
• Current HDF storage implementation
• Current MVN conditioning implementation

Proposed change

Switch to using shake_result.hdf as the preferred source for ShakeMap-based GMF calculations and reconstruct the posterior covariance represented by the information available in the HDF rather than imposing a new correlation structure on the total std grids from uncertainty.xml.

The workflow should:

1. Prefer download/shake_result.hdf when retrieving a ShakeMap by USGS ID.
2. Read the conditioned mean and uncertainty components for the requested IMTs.
3. Reconstruct the station-conditioned within-event covariance from the recorded correlation model and station bookkeeping.
4. Reconstruct the conditional between-event covariance from the saved C and Sigma_HH_YD factors.
5. Generate GMFs around the conditioned ShakeMap mean using the common correlation/sampling machinery used by regular scenario and event-based calculations.
6. Use the same correlation model used by ShakeMap for conditioning, normally Loth and Baker (2013), for generating the gmf realizations.
7. Preserve the existing XML and NPY inputs as a legacy compatibility path. When correlated realizations are generated from marginal XML uncertainties, document and warn that this is an approximation rather than exact posterior sampling.

Existing oq run, oq engine --run, and oq shakemap2gmfs workflows should continue to function as before.

Backward compatibility

• Current XML and NPY inputs should remain supported.
• Existing uncorrelated or median-only or legacy calculations should retain their current results.
• Scenarios starting from ShakeMaps for new/recent events should be based on shake_result.hdf
• Older or incomplete products can continue using the legacy xml files or fail clearly when exact posterior sampling is explicitly required.

Relationship to the correlation refactor

This work depends on the common covariance/sampling interfaces and implementation of the Loth and Baker (2013) spatial-cross-IMT implementation proposed in #11230

[raoanirudh]

For each individual IMT, the shake_result.hdf file contains enough information to reconstruct the station-conditioned within-event covariance, given that we know the Loth and Baker (2013) spatial-cross-correlation model is used by ShakeMap:

$$W_{YY\mid D}=W_{YY}-W_{YD}W_{DD}^{-1}W_{DY}$$

It also contains the C and Sigma_HH_YD factors needed to reconstruct the conditional between-event covariance.

However, ShakeMap stores separate arrays for each output IMT and does not provide the posterior cross-IMT covariance blocks. Therefore, exact joint multi-IMT posterior sampling does not seem to be possible from the current ShakeMap outputs.

A reasonable approximation could be to:

1. generate jointly correlated unconditional within-event fields using Loth and Baker (2013)
2. apply each IMT’s own ShakeMap conditioning operator to its component of that shared realization
3. add the conditional median supplied by ShakeMap

This would preserve the reconstructed posterior spatial covariance of every individual IMT, ideally also remain positive semidefinite by construction, and reduce to the ordinary Loth and Baker model when no stations are used for conditioning. The resulting cross-IMT posterior dependence would be an approximation, since ShakeMap does not export the joint posterior.