Skip to content

Add dingo_generate_multibanded_domain script#369

Open
annalena-k wants to merge 3 commits into
mainfrom
generate-mfd
Open

Add dingo_generate_multibanded_domain script#369
annalena-k wants to merge 3 commits into
mainfrom
generate-mfd

Conversation

@annalena-k

Copy link
Copy Markdown
Collaborator

Summary
This PR adds a new CLI script to automatically generate a MultibandedFrequencyDomain (MFD) waveform dataset settings file from an existing UniformFrequencyDomain (UFD) settings file.

Previously, determining the MFD band nodes required manual tuning in a notebook. The new script automates this by:

  1. Generating waveforms once under an extreme prior (minimum chirp mass, geocent_time = 0.12 s) to stress-test the decimation.
  2. Searching over the whitened waveform difference threshold using a two-phase strategy (geometric walk from a central value, followed by bisection) until the desired median mismatch is reached.
  3. Saving the resulting MFD settings as a new YAML file next to the input, with _ufd replaced by _mfd in the filename.

Waveforms are generated only once and reused across all search iterations, keeping the computational cost minimal.

Changes

  • dingo/gw/dataset/generate_multibanded_domain.py — new script (replaces manual notebook workflow)
  • dingo/gw/dataset/evaluate_multibanded_domain.py — refactored to use shared utilities
  • dingo/gw/dataset/_multibanded_domain_utils.py — new internal module shared by both scripts, providing build_extreme_prior and rint_mismatch_stats
  • tests/gw/test_multibanded_domain_generation.py — unit tests for all new functions
  • pyproject.toml — registered dingo_generate_multibanded_domain entry point

Usage
Command:

dingo_generate_multibanded_domain --settings_file settings_wfd_ufd_32s.yaml --num_samples 1000 --target_median_mismatch 0.0005 --num_processes 0

Example output:

dingo_generate_multibanded_domain --settings_file settings_wfd_ufd_32s.yaml --num_samples 1000 --target_median_mismatch 0.0005 --num_processes 0
Prior (extreme settings for stress-testing multibanding):
  mass_1: Constraint(minimum=1.0, maximum=500.0, name='mass_1', latex_label='$m_1$', unit=None)
  mass_2: Constraint(minimum=1.0, maximum=500.0, name='mass_2', latex_label='$m_2$', unit=None)
  chirp_mass: DeltaFunction(peak=10.0, name=None, latex_label=None, unit=None)
  mass_ratio: bilby.gw.prior.UniformInComponentsMassRatio(minimum=0.05, maximum=1.0, name='mass_ratio', latex_label='$q$', unit=None, boundary=None, equal_mass=False)
  phase: Uniform(minimum=0.0, maximum=6.283185307179586, name='phase', latex_label='$\\phi$', unit=None, boundary='periodic')
  a_1: Uniform(minimum=0.0, maximum=0.99, name='a_1', latex_label='$a_1$', unit=None, boundary=None)
  a_2: Uniform(minimum=0.0, maximum=0.99, name='a_2', latex_label='$a_2$', unit=None, boundary=None)
  tilt_1: Sine(minimum=0.0, maximum=3.141592653589793, name='tilt_1', latex_label='$\\theta_1$', unit=None, boundary=None)
  tilt_2: Sine(minimum=0.0, maximum=3.141592653589793, name='tilt_2', latex_label='$\\theta_2$', unit=None, boundary=None)
  phi_12: Uniform(minimum=0.0, maximum=6.283185307179586, name='phi_12', latex_label='$\\Delta\\phi$', unit=None, boundary='periodic')
  phi_jl: Uniform(minimum=0.0, maximum=6.283185307179586, name='phi_jl', latex_label='$\\phi_{JL}$', unit=None, boundary='periodic')
  theta_jn: Sine(minimum=0.0, maximum=3.141592653589793, name='theta_jn', latex_label='$\\theta_{JN}$', unit=None, boundary=None)
  luminosity_distance: DeltaFunction(peak=100.0, name=None, latex_label=None, unit=None)
  geocent_time: DeltaFunction(peak=0.12, name=None, latex_label=None, unit=None)

Generating 1000 waveforms at 2x resolution...
Setting spin_conversion_phase = 0.0. Using this value for the phase parameter for conversion to cartesian spins.
Generating dataset of size 1000
Computing waveform differences per decimation factor...

Searching for threshold targeting median mismatch ≤ 5.00e-04:
          Step     Threshold   Median mismatch   Bands   Compression
          init     5.000e-03        5.1820e-06       7         27.7x
          up 1     1.000e-02        1.9095e-05       7         30.3x
          up 2     4.000e-02        2.4247e-04       7         36.7x
          up 3     3.200e-01        7.7831e-03       7         50.9x

  Refining bracket [4.000e-02, 3.200e-01]:
      bisect 1     1.131e-01        1.3485e-03       7         42.9x
      bisect 2     6.727e-02        6.0093e-04       7         39.7x
      bisect 3     5.187e-02        3.6988e-04       7         38.0x
      bisect 4     5.907e-02        4.7102e-04       7         38.9x
      bisect 5     6.304e-02        5.3825e-04       7         39.3x
      bisect 6     6.102e-02        4.8410e-04       7         39.0x
      bisect 7     6.202e-02        5.1435e-04       7         39.2x
      bisect 8     6.152e-02        5.0439e-04       7         39.0x
      bisect 9     6.127e-02        4.8410e-04       7         39.0x
     bisect 10     6.140e-02        4.9848e-04       7         39.0x
     bisect 11     6.146e-02        4.9848e-04       7         39.0x
  Nodes unchanged; bisection converged.

Final MFD:
  Nodes          : [20.0, 30.3125, 40.8125, 52.4375, 64.9375, 79.9375, 97.9375, 2047.9375]
  delta_f_initial: 0.03125
  Num bands      : 7
  Compression    : 39.00x

Mismatches between UFD waveforms and MFD waveforms interpolated to UFD.
This is a conservative estimate of the MFD performance when training networks.
num_samples = 1000
  Mean mismatch = 0.0006994324356612301
  Standard deviation = 0.0006637966397081304
  Max mismatch = 0.006980803438541017
  Median mismatch = 0.0004984823001107919
  Percentiles:
    99    -> 0.002959556774394364
    99.9  -> 0.006171660205806462
    99.99 -> 0.006888877547291416

Saved MFD settings to: /path/to/settings_wfd_mfd_32s.yaml

annalena-k and others added 3 commits March 18, 2026 14:55
…n based on a uniform frequency domain waveform settings file.
Add optional token_size to the multibanded-domain generator so it can produce
domains for the Dingo-T1 transformer, whose tokenizer requires each band to span
a whole number of token_size-bin tokens.

- get_band_nodes_for_adaptive_decimation: new min_mfd_bins_per_band argument;
  bands advance by dec_factor * min_mfd_bins_per_band and are asserted to be
  integer multiples of it. Default 1 reproduces the previous NPE banding exactly.
- Thread token_size through _build_mfd_for_threshold and
  generate_multibanded_domain_settings, exposed as --token_size (use 16 for T1).
- Honour the waveform_generator new_interface flag instead of hardcoding the new
  interface, matching generate_dataset.py (standard WaveformGenerator is the
  default used by IMRPhenomXPHM / Dingo-T1).
- MultibandedFrequencyDomain.frequency_mask now returns a boolean array so it can
  be used directly for index selection.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
build_extreme_prior pinned chirp_mass to a scalar, which bilby turns into a
DeltaFunction. Sampling such a prior together with the mass_1/mass_2 Constraint
priors (which forces constrained sampling) raised "non-broadcastable output
operand with shape ()" in PriorDict.sample, so the generator could not run on any
mass-constrained prior (independent of the chirp-mass value). Pin chirp_mass with
a negligibly narrow Uniform at the minimum instead; it samples cleanly while
keeping every draw at the minimum chirp mass.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
The 1x ASD is obtained by taking every other sample of the 2x ASD, so only
``ufd_2x`` is needed. A spectral artefact (bump) in the frequency range 477–483 Hz
is suppressed by replacing those values with the ASD value at 477 Hz, preventing it
from artificially inflating the waveform difference in that band.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess this feature is specific to the ASD aLIGO_ZERO_DET_high_P_asd.txt that is used by default, but would not be the case for other ASDs which will in general have other features at different frequencies (e.g., violin modes) and _get_asd_file can in principle return a user-supplied ASD. Would it be possible to check that the procedure is robust when using some different ASD, for example, something from a real event which has more features than a synthetic ASD?


CLI usage::

python -m dingo.gw.dataset.generate_multibanded_domain \\

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This could also mention the entry point dingo_generate_multibanded_domain?

"from a waveform dataset settings file with a UniformFrequencyDomain.",
)
parser.add_argument(
"--settings_file",

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The arguments here use underscores, but the CLI usage at the top of the file uses hyphens, this should be made consistent.

Overall most of the DINGO code uses hyphens, but generate_dataset.pyuses underscores, so maybe it's best to keep this file consistent with that one, and update the CLI usage at the top?

)

print("Computing waveform differences per decimation factor...")
diffs, freqs = compute_waveform_difference_per_decimation_factor(

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could this have the option of passing difference_over_full_window? Otherwise that option isn't really used.

frequencies_per_decimation_factor,
diffs_per_decimation_factor,
):
f_max = f[np.argmax(diff < threshold)]

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe this could be something like

below_threshold = np.where(diff < threshold)[0]
if len(below_threshold) == 0:
    break
f_max = f[below_threshold[0]]

If the threshold is very small and the decimation factor error is never below the threshold, this would otherwise return f_max = f[0] and give the maximum decimation, while we want the opposite.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants