diff --git a/CHANGELOG-unreleased.md b/CHANGELOG-unreleased.md index d8ff0a4b7..d87f8ac1e 100644 --- a/CHANGELOG-unreleased.md +++ b/CHANGELOG-unreleased.md @@ -45,6 +45,7 @@ the released changes. - Moved the test in `test_pmtransform_units.py` into a function. - Fixed bug in residual calculation when adding or removing phase wraps - Fix #1766 by correcting logic and more clearly naming argument (clkcorr->undo_clkcorr) +- `make_fake_toas_fromtim` now handles BIPM corrections and wideband DMs correctly. - Fix removal of top-level parameter - Minimal fixes to allow usage of numpy 2.0 ### Removed diff --git a/src/pint/simulation.py b/src/pint/simulation.py index 7e3782cbb..86153039a 100644 --- a/src/pint/simulation.py +++ b/src/pint/simulation.py @@ -514,25 +514,14 @@ def make_fake_toas_fromtim( -------- :func:`make_fake_toas` """ - ephem = ( - model.EPHEM.value - if hasattr(model, "EPHEM") and model.EPHEM.value is not None - else None - ) - planets = ( - model.PLANET_SHAPIRO.value - if hasattr(model, "PLANET_SHAPIRO") and model.PLANET_SHAPIRO.value is not None - else False - ) + ts = pint.toa.get_TOAs(timfile, model=model) - input_ts = pint.toa.get_TOAs(timfile, ephem=ephem, planets=planets) - - if input_ts.is_wideband(): - dm_errors = input_ts.get_dm_errors() - ts = update_fake_dms(model, input_ts, dm_errors, add_noise) + if ts.is_wideband(): + dm_errors = ts.get_dm_errors() + ts = update_fake_dms(model, ts, dm_errors, add_noise) return make_fake_toas( - input_ts, + ts, model=model, add_noise=add_noise, add_correlated_noise=add_correlated_noise, diff --git a/tests/test_fake_toas.py b/tests/test_fake_toas.py index 16c95e46a..199389d2e 100644 --- a/tests/test_fake_toas.py +++ b/tests/test_fake_toas.py @@ -358,6 +358,8 @@ def test_fake_from_timfile(planets): ) r_sim = pint.residuals.Residuals(t_sim, m) + assert t.clock_corr_info["bipm_version"] == t_sim.clock_corr_info["bipm_version"] + assert np.isclose( r.calc_time_resids().std(), r_sim.calc_time_resids().std(), rtol=2 )