Summary
write_restart_file (src/roguewave/wavewatch3/io.py) does the following when given a Dataset or Spectrum input:
if isinstance(spectra, (Dataset, Spectrum)):
spectra = spectra.variance_density.values
roguewavespectrum.Spectrum.variance_density is not an alias for the full directional spectrum — for a 2D/directional spectrum it directionally-integrates and returns a 1D frequency-only marginal. Confirmed directly: Spectrum._directionally_integrate (_spectrum.py) does (data_array * self.direction_binwidth).sum(NAME_D, skipna=True) — a standard xarray .sum(dim=...) that drops the direction dimension entirely (not reduced to size 1). This is inconsistent with the sibling function in the same file, write_partial_restart_file, which correctly reads .directional_variance_density.values for the same kind of input.
Impact
Every Spectrum/Dataset this module itself produces is built from directional_variance_density, never variance_density:
RestartFile.__getitem__ (restart_file.py) builds its dataset with "directional_variance_density".
RestartFile.interpolate_in_space does the same.
So calling write_restart_file with any such Spectrum collapses the direction axis via directional integration before the shape checks run — confirmed a directional (n_points, n_freq, n_dir) spectrum comes back as (n_points, n_freq). The very next line does shape[2] != parent_restart_file.number_of_frequencies — but shape now only has 2 dimensions, so this raises IndexError: tuple index out of range instead of the intended ValueError, and the restart file is never written.
This isn't a theoretical path: io.py's own clone_restart_file helper calls write_restart_file(source_restart_file[:], target, source_restart_file, True) — source_restart_file[:] returns exactly this kind of directional Spectrum via __getitem__. Also reproduces in esm-python-scripts's rewrite_restart_file.py's same-resolution branch (input_res == output_res), which passes a Spectrum from __getitem__ straight into write_restart_file.
Scope check — confirmed NOT to affect the newly-merged data-assimilation#152 driver. I checked wavespotter/data-assimilation's drivers/compute_ensemble_mean_ww3_restart_files.py on main directly: it never calls write_restart_file at all. It calls write_partial_restart_file instead, passing a raw DataArray (from xr.concat(...).mean(dim="member") on .directional_variance_density values) — which takes write_partial_restart_file's elif isinstance(spectra, DataArray) branch, bypassing the Dataset/Spectrum branch this bug lives in entirely. So the two known trigger points remain clone_restart_file (unused by any test) and esm-python-scripts's rewrite_restart_file.py (unmaintained, deprioritized per sofarocean/roguewave's current direction).
Why this was never caught
No test in tests/restart_files/test_io.py calls write_restart_file with a Spectrum or Dataset argument — coverage there is limited to write_partial_restart_file and reassemble_restart_file_from_parts. clone_restart_file, the one in-repo caller that would hit this (via __getitem__'s Spectrum output), is imported into the test helpers module but never actually called by any test.
Suggested fix
Change spectra.variance_density.values to spectra.directional_variance_density.values in write_restart_file, matching write_partial_restart_file. Add a test that calls write_restart_file with an actual Spectrum/Dataset (e.g. via clone_restart_file) rather than only a raw numpy.ndarray, to catch this class of regression.
Related
Found while investigating sofarocean/roguewave#14 and the esm-python-scripts regrid call site (same file, io.py, same "correct in one sibling function, wrong in the other" pattern as #14's depth-interpolation bug and #12's directions-axis validation bug).
Summary
write_restart_file(src/roguewave/wavewatch3/io.py) does the following when given aDatasetorSpectruminput:roguewavespectrum.Spectrum.variance_densityis not an alias for the full directional spectrum — for a 2D/directional spectrum it directionally-integrates and returns a 1D frequency-only marginal. Confirmed directly:Spectrum._directionally_integrate(_spectrum.py) does(data_array * self.direction_binwidth).sum(NAME_D, skipna=True)— a standard xarray.sum(dim=...)that drops the direction dimension entirely (not reduced to size 1). This is inconsistent with the sibling function in the same file,write_partial_restart_file, which correctly reads.directional_variance_density.valuesfor the same kind of input.Impact
Every
Spectrum/Datasetthis module itself produces is built fromdirectional_variance_density, nevervariance_density:RestartFile.__getitem__(restart_file.py) builds its dataset with"directional_variance_density".RestartFile.interpolate_in_spacedoes the same.So calling
write_restart_filewith any suchSpectrumcollapses the direction axis via directional integration before the shape checks run — confirmed a directional(n_points, n_freq, n_dir)spectrum comes back as(n_points, n_freq). The very next line doesshape[2] != parent_restart_file.number_of_frequencies— butshapenow only has 2 dimensions, so this raisesIndexError: tuple index out of rangeinstead of the intendedValueError, and the restart file is never written.This isn't a theoretical path:
io.py's ownclone_restart_filehelper callswrite_restart_file(source_restart_file[:], target, source_restart_file, True)—source_restart_file[:]returns exactly this kind of directionalSpectrumvia__getitem__. Also reproduces inesm-python-scripts'srewrite_restart_file.py's same-resolution branch (input_res == output_res), which passes aSpectrumfrom__getitem__straight intowrite_restart_file.Scope check — confirmed NOT to affect the newly-merged
data-assimilation#152driver. I checkedwavespotter/data-assimilation'sdrivers/compute_ensemble_mean_ww3_restart_files.pyonmaindirectly: it never callswrite_restart_fileat all. It callswrite_partial_restart_fileinstead, passing a rawDataArray(fromxr.concat(...).mean(dim="member")on.directional_variance_densityvalues) — which takeswrite_partial_restart_file'selif isinstance(spectra, DataArray)branch, bypassing theDataset/Spectrumbranch this bug lives in entirely. So the two known trigger points remainclone_restart_file(unused by any test) andesm-python-scripts'srewrite_restart_file.py(unmaintained, deprioritized persofarocean/roguewave's current direction).Why this was never caught
No test in
tests/restart_files/test_io.pycallswrite_restart_filewith aSpectrumorDatasetargument — coverage there is limited towrite_partial_restart_fileandreassemble_restart_file_from_parts.clone_restart_file, the one in-repo caller that would hit this (via__getitem__'sSpectrumoutput), is imported into the test helpers module but never actually called by any test.Suggested fix
Change
spectra.variance_density.valuestospectra.directional_variance_density.valuesinwrite_restart_file, matchingwrite_partial_restart_file. Add a test that callswrite_restart_filewith an actualSpectrum/Dataset(e.g. viaclone_restart_file) rather than only a rawnumpy.ndarray, to catch this class of regression.Related
Found while investigating
sofarocean/roguewave#14and theesm-python-scriptsregrid call site (same file,io.py, same "correct in one sibling function, wrong in the other" pattern as#14's depth-interpolation bug and#12's directions-axis validation bug).