Fisher and DALI Forecasts

gwforge_fisher forecasts how well a network will measure the parameters of a source. It works with any waveform bilby.gw.WaveformGenerator can produce.

At order = 1 this is the Fisher matrix of Vallisneri (2008),

$$\Gamma_{ab} = \langle \partial_a h \mid \partial_b h\rangle, \qquad \Sigma = \Gamma^{-1},$$

with the noise-weighted inner product \(\langle a \mid b\rangle = 4,\mathrm{Re}\int \tilde a^* \tilde b / S_n(f),\mathrm{d}f\) summed over the detectors of the network.

At order = 2 it adds the doublet-DALI correction of Sellentin, Quartin & Amendola (2014), as applied to gravitational waves by GWDALI:

$$\log L = -\tfrac{1}{2} F_{ab}\Delta^a\Delta^b - \tfrac{1}{2} G_{abc}\Delta^a\Delta^b\Delta^c - \tfrac{1}{8} H_{abcd}\Delta^a\Delta^b\Delta^c\Delta^d,$$

with \(F_{ab} = \langle h_{,a}|h_{,b}\rangle\), \(G_{abc} = \langle h_{,ab}|h_{,c}\rangle\) and \(H_{abcd} = \langle h_{,ab}|h_{,cd}\rangle\).

Note

order = 3 (triplet-DALI) is not implemented and raises an error. Grouping the expansion by highest derivative, the terms above already exhaust everything expressible with first and second derivatives; the next group consists entirely of third derivatives of the waveform.

Quick start

gwforge_fisher --config-file fisher.ini
[Injection_Parameters]
; either write the source out here ...
chirp_mass = 28.0
symmetric_mass_ratio = 0.24
chi_1 = 0.3
chi_2 = -0.2
luminosity_distance = 1500.0
theta_jn = 0.7
psi = 0.9
phase = 1.3
geocent_time = 1893024018.0
ra = 1.375
dec = -1.2108

; ... or point at a gwforge_population file instead
; injection-file = bbh.h5
; start-index = 0
; end-index = 100

[Waveform_Generator]
waveform-approximant = SEOBNRv5HM_ROM
waveform-minimum-frequency = 20
frequency-domain-source-model = lal_binary_black_hole

[IFOS]
detectors = ['CE20', 'CE40', 'ET']
sampling-frequency = 4096

[Fisher]
order = 1

[Output]
output-file = fisher.pkl
plot-corner = True

Important

Point injection-file at a catalogue and the redundant parametrisations are stripped for you. A gwforge_population file carries mass_1, mass_2, chirp_mass, symmetric_mass_ratio, mass_ratio, total_mass and the _source variants, all at once. bilby’s convert_to_lal_binary_black_hole_parameters prefers the component masses and never looks at chirp_mass, so a Fisher over (chirp_mass, symmetric_mass_ratio) would displace keys the waveform ignores: the mass derivatives came back exactly zero, the Fisher was singular in the mass block, and nothing raised — \(\sigma(\mathcal{M})\) was wrong by six orders of magnitude. GWForge.fisher.parameters.strip_shadowed_parameters() now removes whatever duplicates a varied parameter, and logs what it dropped. The same shadowing applies to spins, where a_1 wins over chi_1.

If you build a FisherMatrix yourself, this is handled in its constructor; if you assemble derivatives directly, call the helper first.

frequency-domain-source-model accepts lal_binary_black_hole, lal_binary_neutron_star, gwsignal_binary_black_hole, lal_eccentric_binary_black_hole_no_spins and EFPE_binary_black_hole now. You can plugin any waveform generator pretty easily.

Three named parameter sets are available, and picking the right one matters: DEFAULT_PARAMETERS (aligned spin), PRECESSING_PARAMETERS and ECCENTRIC_PRECESSING_PARAMETERS. Using the aligned-spin set with a precessing waveform is a silent mistake — bilby’s conversion turns chi_1/chi_2 into a_1/a_2 with both tilts zero, so the waveform generates happily and you forecast an aligned-spin binary.

SEOBNRv5HM/SEOBNRv5PHM need pyseobnr (pip install -e .[waveforms]), and pyEFPEHM is installed from git. pyEFPEHM ignores minimum_frequency: set f22_start in the waveform arguments instead, and set it below the analysis band, because its eccentric harmonics radiate below the quadrupole start frequency. Similarly, for TEOBResumS install teobresums.

The segment length is chosen per source from GWForge.conversion.get_safe_signal_durations(), the same estimate gwforge_optimal_snr uses, so a forecast and an SNR see the same signal.

Step sizes (Claudio speaking)

Finite-difference steps are calibrated per source so that one step changes the waveform by a fixed fraction, \(\lVert\Delta h\rVert/\lVert h\rVert \approx 10^{-2}\).

This matters more than it sounds. LAL places each higher mode’s onset at a mass-dependent frequency, so differencing across that edge contributes a spurious term growing like \(1/\text{step}\). With a step an order of magnitude too small, \(\sigma(\eta)\) for IMRPhenomXHM is wrong by 150%; with the calibration on, rescaling the starting guess by a factor of thirty moves every \(\sigma\) by less than half a percent.

Measuring the waveform’s own numerical noise

A finite difference cannot see beneath the noise of the function it differences. GWForge measures that noise directly, before computing any derivative, using the ECNoise procedure of Moré & Wild (2011): sample a smooth scalar functional of the waveform at closely spaced points and build the finite difference table. A smooth function’s \(k\)-th differences shrink like \(h^k\); noise does not shrink at all.

The functional is the noise-weighted overlap \(\langle h(x)\mid h_0\rangle\) — smooth in the parameter (unlike \(\lVert\Delta h\rVert\), which has a kink at zero) and PSD-weighted exactly as the Fisher matrix is. It is kept complex, because the real part is stationary at the fiducial point for a phase-like parameter.

WaveformDerivatives.noise_ratios reports the noise divided by the change one step produces — the approximate relative error of that derivative — and it is stored in the results pickle so the caveat travels with the forecast. Measured worst case per model:

model

worst relative error

IMRPhenomXHM

4e-6 — machine precision

SEOBNRv5HM

3e-3 on masses and spins; 5e-14 on theta_jn and phase

SEOBNRv5PHM

2e-4 on the precession angles

pyEFPEHM

2e-4 on chirp_mass

The SEOBNRv5HM row is the one that shows the estimator is measuring something real: masses and spins pass through the adaptive EOB ODE solve and pick up scatter, while theta_jn and phase are applied analytically outside it and come out ten orders of magnitude cleaner.

The estimate is used only to raise a step that is too small to clear the noise. For a closed-form waveform it never binds.

Note

Noise is not the dominant limit in practice. Scanning IMRPhenomXHM’s chi_1 step from 3e-4 to 1e-2 moves sigma(chi_1) by 0.03%; the same scan for SEOBNRv5HM moves it by 70%, with only a narrow flat region around 3e-4 to 1e-3. That is ten times larger than the measured noise, so it is truncation — the derivative genuinely varying with the step — not scatter. A Fisher forecast with an EOB waveform is correspondingly less well determined than one with a closed-form model, and no step-selection strategy changes that.

Override the calibration if you need to:

[Fisher]
calibrate-steps = False
step-sizes = {'chirp_mass': 1e-4, 'chi_1': 1e-3}
target-change = 1e-2

Response corrections

[Fisher]
earth-rotation = True
finite-size = True

Both default to True, matching how gwforge_inject projects the same signal. Turning both off recovers bilby’s static long-wavelength response.

Priors

[Priors]
prior-file = my.prior          ; optional at order 1, REQUIRED at order 2
enforce-physicality = False
  • Order 1 — optional. The reported covariance is always the pure-likelihood \(\Gamma^{-1}\). A Gaussian prior additionally contributes \(1/\sigma^2\) to the diagonal; a uniform prior carries no Fisher information in its interior.

  • Order 2 — required. The DALI posterior is a truncated Taylor expansion, which is only meaningful inside bounds; without them the cubic and quartic terms produce spurious modes that are artefacts of the truncation.

enforce-physicality additionally clips priors to physical ranges (\(q \le 1\), \(|\chi| < 1\), angles in range).

Warning

A prior narrower than the Fisher \(\sigma\) silently sets the posterior width rather than bounding it. GWForge warns when a prior reaches fewer than three sigma from the fiducial value.

Output

The results pickle is a list of one record per source, each holding the fiducial parameters, the network and per-detector Fisher matrices, the covariance, the sigmas, the optimal SNRs, the condition number and inversion error, the calibrated steps — plus doublet, quartic and posterior samples at order 2.

Order 2 posteriors are non-Gaussian, so they have no covariance matrix and no closed-form marginals. A corner plot is a set of marginals, so GWForge samples log L_DALI with emcee through bilby.run_sampler. All the waveform evaluations happen while building the tensors; evaluating the DALI likelihood afterwards is pure tensor contraction.

Two waveforms, one binary

validation/compare_waveforms.py forecasts the same aligned-spin binary with IMRPhenomXHM and SEOBNRv5HM and overlays the two error ellipses. Both models describe the source equally well, so the difference between them is a lower bound on the waveform systematic in any Fisher forecast.

python validation/compare_waveforms.py --outdir waveform_comparison

For a 28 M☉ chirp mass source at network SNR 329 on CE40+CE20+ET:

quantity

agreement

network SNR

0.07%

luminosity_distance, ra, dec, theta_jn, psi

0.6 – 1.6%

chirp_mass

35%

geocent_time, phase

factor 1.8 – 2.4

chi_1, chi_2

factor 4 – 5

The SNRs agree almost exactly and the distance/orientation block overlays closely; the mass–spin–time–phase block, which is where the two models actually differ, does not. This is a genuine waveform systematic, not a numerical artefact — the measured noise in both models is three orders of magnitude too small to explain it.