pint.residuals.WidebandDMResiduals
- class pint.residuals.WidebandDMResiduals(toas=None, model=None, residual_type='toa', unit=Unit('s'), subtract_mean=True, use_weighted_mean=True, track_mode=None, use_abs_phase=True)[source]
Bases:
Residuals
Residuals for independent DM measurement (i.e. Wideband TOAs).
This class manages the DM residuals from data that includes direct DM measurements associated with the TOAs.
pint.residuals.WidebandTOAResiduals
combines one of these objects with apint.residuals.Residuals
object.The values of interest are probably best accessed using the
.resids
property, and the uncertainty using the.get_data_error()
.- Variables:
dm_data (
astropy.units.Quantity
) – The DM data extracted from the TOAs.dm_error (
astropy.units.Quantity
) – The DM uncertainties extracted from the TOAs.
- Parameters:
toas (
pint.toa.TOAs
) – TOAs. They should include DM measurement data.model (
pint.models.timing_model.TimingModel
) – The timing model.
Methods
Return the weighted chi-squared for the model and toas.
calc_phase_mean
([weighted])Calculate mean phase of residuals, optionally weighted
calc_phase_resids
([subtract_mean, ...])Compute timing model residuals in pulse phase.
calc_resids
()calc_time_mean
([calctype, weighted])Calculate mean time of residuals, optionally weighted
calc_time_resids
([calctype, subtract_mean, ...])Compute timing model residuals in time (seconds).
Compute whitened timing residuals (dimensionless).
d_Ndiag_d_param
(param)Derivative of the white noise covariance matrix diagonal elements w.r.t.
d_lnlikelihood_d_ECORR
(param)d_lnlikelihood_d_Ndiag
()d_lnlikelihood_d_param
(param)ecorr_average
([use_noise_model])Uses the ECORR noise model time-binning to compute "epoch-averaged" residuals.
get_PSR_freq
([calctype])Return pulsar rotational frequency in Hz.
get_data_error
([scaled])Get data errors.
Get the independent measured DM data from TOA flags.
Compute the log-likelihood for the model and TOAs.
Compute weighted RMS of the residuals in time.
update
()Recalculate everything in residuals class after changing model or TOAs
update_model
(new_model, **kwargs)Up date DM models from a new PINT timing model
Attributes
Compute chi-squared as needed and cache the result.
Reduced chi-squared.
Return number of degrees of freedom for the DM model.
Return the weighted reduced chi-squared for the model and toas.
Residuals in time units.
Get pure value of the residuals use the given base unit.
- property resids
Residuals in time units.
- property resids_value
Get pure value of the residuals use the given base unit.
- property dof
Return number of degrees of freedom for the DM model.
- get_data_error(scaled=True)[source]
Get data errors.
- Parameters:
scaled (bool, optional) – If errors get scaled by the noise model.
- calc_chi2()[source]
Return the weighted chi-squared for the model and toas.
If the errors on the TOAs are independent this is a straightforward calculation, but if the noise model introduces correlated errors then obtaining a meaningful chi-squared value requires a Cholesky decomposition. This is carried out using the :method:`~pint.residuals.Residuals._calc_gls_chi2` helper function.
The return value here is available as self.chi2, which will not redo the computation unless necessary.
The chi-squared value calculated here is suitable for use in downhill minimization algorithms and Bayesian approaches.
Handling of problematic results - degenerate conditions explored by a minimizer for example - may need to be checked to confirm that they correctly return infinity.
If lognorm=True is given, the log-normalization-factor of the likelihood function will also be returned.
- Parameters:
lognorm (bool) – If True, return the the log-normalization-factor of the likelihood function along with the chi2 value.
- Returns:
chi2 if lognorm is False
(chi2, log_norm) if lognorm is True
- get_dm_data()[source]
Get the independent measured DM data from TOA flags.
This is to extract DM and uncertainty data from its representation in the flags on TOAs.
FIXME: there should be a
set_dm_data
function.- Returns:
valid_dm (numpy.ndarray) – Independent measured DM data from TOA line. It only returns the DM values that is present in the TOA flags.
valid_error (numpy.ndarray) – The error associated with DM values in the TOAs.
valid_index (list) – The TOA with DM data index.
- update_model(new_model, **kwargs)[source]
Up date DM models from a new PINT timing model
- Parameters:
new_model (pint.timing_model.TimingModel) –
- calc_phase_mean(weighted=True)
Calculate mean phase of residuals, optionally weighted
- Parameters:
weighted (bool, optional) –
- Return type:
- calc_phase_resids(subtract_mean=None, use_weighted_mean=None, use_abs_phase=None)
Compute timing model residuals in pulse phase.
if
subtract_mean
oruse_weighted_mean
is None, will use the values set for the object itself- Parameters:
subtract_mean (bool or None, optional) – Subtract the mean of the residuals. This is ignored if the PhaseOffset component is present in the model. Default is to use the class attribute.
use_weighted_mean (bool or None, optional) – Whether to use weighted mean for mean subtraction. Default is to use the class attribute.
use_abs_phase (bool or None, optional) – Whether to use absolute phase (w.r.t. the TZR TOA). Default is to use the class attribute.
- Return type:
- calc_time_mean(calctype='taylor', weighted=True)
Calculate mean time of residuals, optionally weighted
- Parameters:
calctype (str, optional) – Calculation time for phase to time conversion. See
pint.residuals.Residuals.calc_time_resids()
for details.weighted (bool, optional) –
- Return type:
- calc_time_resids(calctype='taylor', subtract_mean=None, use_weighted_mean=None, use_abs_phase=None)
Compute timing model residuals in time (seconds).
Converts from phase residuals to time residuals using several possible ways to calculate the frequency.
If
subtract_mean
oruse_weighted_mean
is None, will use the values set for the object itself- Parameters:
calctype ({'taylor', 'modelF0', 'numerical'}) – Type of calculation. If calctype == “modelF0”, then simply the
F0
parameter from the model. If calctype == “numerical”, then try a numerical derivative If calctype == “taylor”, evaluate the frequency with a Taylor seriessubtract_mean (bool or None, optional) – Subtract the mean of the residuals. This is ignored if the PhaseOffset component is present in the model. Default is to use the class attribute.
use_weighted_mean (bool or None, optional) – Whether to use weighted mean for mean subtraction. Default is to use the class attribute.
use_abs_phase (bool or None, optional) – Whether to use absolute phase (w.r.t. the TZR TOA). Default is to use the class attribute.
- Returns:
residuals
- Return type:
- calc_whitened_resids()
Compute whitened timing residuals (dimensionless).
Whitened residuals are computed by subtracting the correlated noise realization from the time residuals and normalizes the result using scaled TOA uncertainties.
This requires the noise_resids attribute to be set. This is usually available in a post-fit residuals object.
Example usage:
>>> ftr = Fitter.auto(toas, model) >>> ftr.fit_toas() >>> res = ftr.resids >>> white_res = res.calc_whitened_resids()
- property chi2
Compute chi-squared as needed and cache the result.
- property chi2_reduced
Reduced chi-squared.
- d_Ndiag_d_param(param)
Derivative of the white noise covariance matrix diagonal elements w.r.t. a white noise parameter (EFAC or EQUAD).
- ecorr_average(use_noise_model=True)
Uses the ECORR noise model time-binning to compute “epoch-averaged” residuals.
Requires ECORR be used in the timing model. If
use_noise_model
is true, the noise model terms (EFAC, EQUAD, ECORR) will be applied to the TOA uncertainties, otherwise only the raw uncertainties will be used.Returns a dictionary with the following entries:
mjds Average MJD for each segment
freqs Average topocentric frequency for each segment
time_resids Average residual for each segment, time units
noise_resids Dictionary of per-noise-component average residual
errors Uncertainty on averaged residuals
indices List of lists giving the indices of TOAs in the original TOA table for each segment
- get_PSR_freq(calctype='modelF0')
Return pulsar rotational frequency in Hz.
- Parameters:
calctype ({'modelF0', 'numerical', 'taylor'}) – Type of calculation. If calctype == “modelF0”, then simply the
F0
parameter from the model. If calctype == “numerical”, then try a numerical derivative If calctype == “taylor”, evaluate the frequency with a Taylor series- Returns:
freq – Either the single
F0
in the model or the spin frequency at the moment of each TOA.- Return type:
- lnlikelihood()
Compute the log-likelihood for the model and TOAs.
- property reduced_chi2
Return the weighted reduced chi-squared for the model and toas.
- update()
Recalculate everything in residuals class after changing model or TOAs