pyfstat package
These pages document the full API and set of classes provided by PyFstat.
See here for installation instructions and other general information.
Subpackages
- pyfstat.utils package
- Submodules
- pyfstat.utils.atoms module
- pyfstat.utils.cli module
- pyfstat.utils.converting module
- pyfstat.utils.ephemeris module
- pyfstat.utils.formatting module
- pyfstat.utils.gsl module
- pyfstat.utils.importing module
- pyfstat.utils.io module
- pyfstat.utils.predict module
- pyfstat.utils.runlalsuite module
- pyfstat.utils.sft module
- Module contents
Submodules
pyfstat.core module
The core tools used in pyfstat
- class pyfstat.core.BaseSearchClass(*args, **kwargs)[source]
Bases:
object
The base class providing parent methods to other PyFstat classes.
This does not actually have any ‘search’ functionality, which needs to be added by child classes along with full initialization and any other custom methods.
- binary_keys = ['asini', 'period', 'ecc', 'tp', 'argp']
List of extra parameters for sources in binaries.
- default_search_keys = ['F0', 'F1', 'F2', 'Alpha', 'Delta']
Default order of the traditionally supported search parameter names.
FIXME: these are only used as fallbacks for the deprecated style of passing keys one by one; not needed when using the new parameters dictionary.
- tex_labels = {'Alpha': '$\\alpha$', 'Delta': '$\\delta$', 'F0': '$f$', 'F1': '$\\dot{f}$', 'F2': '$\\ddot{f}$', 'F3': '$\\dddot{f}$', 'argp': '$\\mathrm{argp}$', 'asini': '$\\mathrm{asin}\\,i$', 'delta_F0': '$\\delta f$', 'delta_F1': '$\\delta \\dot{f}$', 'ecc': '$\\mathrm{ecc}$', 'lnBtSG': '$\\ln\\mathcal{B}_{\\mathrm{tS/G}}$', 'log10BSGL': '$\\log_{10}\\mathcal{B}_{\\mathrm{SGL}}$', 'maxTwoF': '$\\max\\widetilde{2\\mathcal{F}}$', 'period': '$P$', 'tglitch': '$t_\\mathrm{glitch}$', 'tp': '$t_p$', 'transient_duration': '$\\Delta T$', 'transient_tstart': '$t_\\mathrm{start}$', 'twoF': '$\\widetilde{2\\mathcal{F}}$'}
Formatted labels used for plot annotations.
- unit_dictionary = {'Alpha': 'rad', 'Delta': 'rad', 'F0': 'Hz', 'F1': 'Hz/s', 'F2': 'Hz/s$^2$', 'F3': 'Hz/s$^3$', 'argp': '', 'asini': '', 'delta_F0': 'Hz', 'delta_F1': 'Hz/s', 'ecc': '', 'period': 's', 'tglitch': 's', 'tp': 's', 'transient_duration': 's', 'transient_tstart': 's'}
Units for standard parameters.
- fmt_detstat = '%.9g'
Standard output precision for detection statistics.
- fmt_doppler = '%.16g'
Standard output precision for Doppler (frequency evolution) parameters.
- set_ephemeris_files(earth_ephem=None, sun_ephem=None)[source]
Set the ephemeris files to use for the Earth and Sun.
NOTE: If not given explicit arguments, default values from utils.get_ephemeris_files() are used.
- Parameters:
earth_ephem (str, optional) – Paths of the two files containing positions of Earth and Sun, respectively at evenly spaced times, as passed to CreateFstatInput
Default:None
sun_ephem (str, optional) – Paths of the two files containing positions of Earth and Sun, respectively at evenly spaced times, as passed to CreateFstatInput
Default:None
- pprint_init_params_dict()[source]
Pretty-print a parameters dictionary for output file headers.
- Returns:
pretty_init_parameters – A list of lines to be printed, including opening/closing “{” and “}”, consistent indentation, as well as end-of-line commas, but no comment markers at start of lines.
- Return type:
list
- get_output_file_header()[source]
Constructs a meta-information header for text output files.
This will include PyFstat and LALSuite versioning, information about when/where/how the code was run, and input parameters of the instantiated class.
- Returns:
header – A list of formatted header lines.
- Return type:
list
- read_par(filename=None, label=None, outdir=None, suffix='par', raise_error=True)[source]
Read a key=val file and return a dictionary.
- Parameters:
filename (str or None, optional) – Filename (path) containing rows of key=val data to read in.
Default:None
label (str or None, optional) – If filename is None, form the file to read as outdir/label.suffix.
Default:None
outdir (str or None, optional) – If filename is None, form the file to read as outdir/label.suffix.
Default:None
suffix (str or None, optional) – If filename is None, form the file to read as outdir/label.suffix.
Default:'par'
raise_error (bool, optional) – If True, raise an error for lines which are not comments, but cannot be read.
Default:True
- Returns:
params_dict – A dictionary of the parsed key=val pairs.
- Return type:
dict
- static translate_keys_to_lal(dictionary)[source]
Convert input keys into lalpulsar convention.
In PyFstat’s convention, input keys (search parameter names) are F0, F1, F2, …, while lalpulsar functions prefer to use Freq, f1dot, f2dot, ….
Since lalpulsar keys are only used internally to call lalpulsar routines, this function is provided so the keys can be translated on the fly.
- Parameters:
dictionary (dict) – Dictionary to translate. A copy will be made (and returned) before translation takes place.
- Returns:
translated_dict – Copy of “dictionary” with new keys according to lalpulsar convention.
- Return type:
dict
- class pyfstat.core.ComputeFstat(*args, **kwargs)[source]
Bases:
BaseSearchClass
Base search class providing an interface to lalpulsar.ComputeFstat.
In most cases, users should be using one of the higher-level search classes from the grid_based_searches or mcmc_based_searches modules instead.
See the lalpulsar documentation at https://lscsoft.docs.ligo.org/lalsuite/lalpulsar/group___compute_fstat__h.html and R. Prix, The F-statistic and its implementation in ComputeFstatistic_v2 ( https://dcc.ligo.org/T0900149/public ) for details of the lalpulsar module and the meaning of various technical concepts as embodied by some of the class’s parameters.
Normally this will read in existing data through the sftfilepattern argument, but if that option is None and the necessary alternative arguments are used, it can also generate simulated data (including noise and/or signals) on the fly.
NOTE that the detection statistics that can be computed from an instance of this class depend on the BSGL, BtSG and transientWindowType arguments given at initialisation. See get_fullycoherent_detstat() and get_transient_detstats() for details. To change what you want to compute, you may need to initialise a new instance with different options.
NOTE for GPU users (tCWFstatMapVersion=”pycuda”): This class tries to conveniently deal with GPU context management behind the scenes. A known problematic case is if you try to instantiate it twice from the same session/script. If you then get some messages like RuntimeError: make_default_context() and invalid device context, that is because the GPU is still blocked from the first instance when you try to initiate the second. To avoid this problem, use context management:
with pyfstat.ComputeFstat( [...], tCWFstatMapVersion="pycuda", ) as search: search.get_fullycoherent_detstat([...])
or manually call the search.finalizer_() method where needed.
- Parameters:
tref (int) – GPS seconds of the reference time.
sftfilepattern (str, optional) – Pattern to match SFTs using wildcards (*?) and ranges [0-9]; multiple patterns can be given separated by colons.
Default:None
minStartTime (int, optional) – Only use SFTs with timestamps starting from within this range, following the XLALCWGPSinRange convention: half-open intervals [minStartTime,maxStartTime].
Default:None
maxStartTime (int, optional) – Only use SFTs with timestamps starting from within this range, following the XLALCWGPSinRange convention: half-open intervals [minStartTime,maxStartTime].
Default:None
Tsft (int, optional) – SFT duration in seconds. Only required if sftfilepattern=None and hence simulted data is generated on the fly.
Default:1800
binary (bool, optional) – If true, search over binary parameters.
Default:False
singleFstats (bool, optional) – If true, also compute the single-detector twoF values.
Default:False
BSGL (bool, optional) –
If true, compute the log10BSGL statistic rather than the twoF value. For details, see Keitel et al (PRD 89, 064023, 2014): https://arxiv.org/abs/1311.5738 Note this automatically sets singleFstats=True as well. Tuning parameters are currently hardcoded:
Default:False
Fstar0=15 for coherent searches.
A p-value of 1e-6 and correspondingly recalculated Fstar0 for semicoherent searches.
Uniform per-detector prior line-vs-Gaussian odds.
BtSG (bool, optional) – If true and transientWindowType is not None, compute the transient \(\ln\mathcal{B}_{\mathrm{tS}/\mathrm{G}}\) statistic from Prix, Giampanis & Messenger (PRD 84, 023007, 2011) (tCWFstatMap marginalised over uniform t0, tau priors). rather than the maxTwoF value.
Default:False
transientWindowType (str, optional) – If rect or exp, allow for the Fstat to be computed over a transient range. (none instead of None explicitly calls the transient-window function, but with the full range, for debugging.) (If not None, will also force atoms regardless of computeAtoms option.)
Default:None
t0Band (int, optional) – Search ranges for transient start-time t0 and duration tau. If >0, search t0 in (minStartTime,minStartTime+t0Band) and tau in (tauMin,2*Tsft+tauBand). If =0, only compute the continuous-wave Fstat with t0=minStartTime, tau=maxStartTime-minStartTime.
Default:None
tauBand (int, optional) – Search ranges for transient start-time t0 and duration tau. If >0, search t0 in (minStartTime,minStartTime+t0Band) and tau in (tauMin,2*Tsft+tauBand). If =0, only compute the continuous-wave Fstat with t0=minStartTime, tau=maxStartTime-minStartTime.
Default:None
tauMin (int, optional) – Minimum transient duration to cover, defaults to 2*Tsft.
Default:None
dt0 (int, optional) – Grid resolution in transient start-time, defaults to Tsft.
Default:None
dtau (int, optional) – Grid resolution in transient duration, defaults to Tsft.
Default:None
detectors (str, optional) – Two-character references to the detectors for which to use data. Specify None for no constraint. For multiple detectors, separate by commas.
Default:None
minCoverFreq (float, optional) – The min and max cover frequency passed to lalpulsar.CreateFstatInput. For negative values, these will be used as offsets from the min/max frequency contained in the sftfilepattern. If either is None, the search_ranges argument is used to estimate them. If the automatic estimation fails and you do not have a good idea what to set these two options to, setting both to -0.5 will reproduce the default behaviour of PyFstat <=1.4 and may be a reasonably safe fallback in many cases.
Default:None
maxCoverFreq (float, optional) – The min and max cover frequency passed to lalpulsar.CreateFstatInput. For negative values, these will be used as offsets from the min/max frequency contained in the sftfilepattern. If either is None, the search_ranges argument is used to estimate them. If the automatic estimation fails and you do not have a good idea what to set these two options to, setting both to -0.5 will reproduce the default behaviour of PyFstat <=1.4 and may be a reasonably safe fallback in many cases.
Default:None
search_ranges (dict, optional) – Dictionary of ranges in all search parameters, only used to estimate frequency band passed to lalpulsar.CreateFstatInput, if minCoverFreq, maxCoverFreq are not specified (==`None`). For actually running searches, grids/points will have to be passed separately to the .run() method. The entry for each parameter must be a list of length 1, 2 or 3: [single_value], [min,max] or [min,max,step].
Default:None
injectSources (dict or str, optional) – Either a dictionary of the signal parameters to inject, or a string pointing to a .cff file defining a signal.
Default:None
injectSqrtSX (float or list or str, optional) – Single-sided PSD values for generating fake Gaussian noise on the fly. Single float or str value: use same for all IFOs. List or comma-separated string: must match len(detectors) and/or the data in sftfilepattern. Detectors will be paired to list elements following alphabetical order.
Default:None
randSeed (int or None, optional) – random seed for on-the-fly noise generation using injectSqrtSX. Setting this to 0 or None is equivalent; both will randomise the seed, following the behaviour of XLALAddGaussianNoise(), while any number not equal to 0 will produce a reproducible noise realisation.
Default:None
assumeSqrtSX (float or list or str, optional) – Don’t estimate noise-floors but assume this (stationary) single-sided PSD. Single float or str value: use same for all IFOs. List or comma-separated string: must match len(detectors) and/or the data in sftfilepattern. Detectors will be paired to list elements following alphabetical order. If working with signal-only data, please set assumeSqrtSX=1 .
Default:None
SSBprec (int, optional) – Flag to set the Solar System Barycentring (SSB) calculation in lalpulsar: 0=Newtonian, 1=relativistic, 2=relativistic optimised, 3=DMoff, 4=NO_SPIN
Default:None
RngMedWindow (int, optional) – Running-Median window size for F-statistic noise normalization (number of SFT bins).
Default:None
tCWFstatMapVersion (str, optional) – Choose between implementations of the transient F-statistic functionality: standard lal implementation, pycuda for GPU version, and some others only for devel/debug.
Default:'lal'
cudaDeviceName (str, optional) – GPU name to be matched against drv.Device output, only for tCWFstatMapVersion=pycuda.
Default:None
computeAtoms (bool, optional) – Request calculation of ‘F-statistic atoms’ regardless of transientWindowType.
Default:False
earth_ephem (str, optional) – Earth ephemeris file path. If None, will check standard sources as per utils.get_ephemeris_files().
Default:None
sun_ephem (str, optional) – Sun ephemeris file path. If None, will check standard sources as per utils.get_ephemeris_files().
Default:None
allowedMismatchFromSFTLength (float, optional) – Maximum allowed mismatch from SFTs being too long [Default: what’s hardcoded in XLALFstatMaximumSFTLength]
Default:None
- init_computefstatistic()[source]
Initialization step for the F-stastic computation internals.
This sets up the special input and output structures the lalpulsar module needs, the ephemerides, optional on-the-fly signal injections, and extra options for multi-detector consistency checks and transient searches.
All inputs are taken from the pre-initialized object, so this function does not have additional arguments of its own.
- estimate_min_max_CoverFreq()[source]
Extract spanned spin-range at reference -time from the template bank.
To use this method, self.search_ranges must be a dictionary of lists per search parameter which can be either [single_value], [min,max] or [min,max,step].
- get_fullycoherent_detstat(F0=None, F1=None, F2=None, Alpha=None, Delta=None, asini=None, period=None, ecc=None, tp=None, argp=None, params=None, tstart=None, tend=None)[source]
Computes the detection statistic(s) fully-coherently at a single point.
Currently supported statistics:
twoF (CW)
log10BSGL (CW or transient)
maxTwoF (transient)
lnBtSG (transient)
All computed statistics are stored as attributes, but only one statistic is returned.
As the basic statistic of this class, twoF is always computed and stored as self.twoF as well, and it is the default return value.
If self.singleFstats, additionally the single-detector 2F-stat values are stored in self.twoFX.
If self.BSGL, the log10BSGL statistic for CWs is additionally stored, and it is returned instead of twoF.
If transient parameters are enabled (self.transientWindowType is set), maxTwoF will always be computed and stored, and returned by default. Depending on the self.BSGL and self.BtSG options, either log10BSGL (a transient version of it, superseding the CW version) or lnBtSG will also be computed, stored, and returned instead of maxTwoF. The full transient-F-stat map is also computed here, but stored in self.FstatMap, not returned.
NOTE: the old way of calling this with explicit [F0,F1,F2,Alpha,Delta,…] parameters is DEPRECATED and may be removed in future versions. Currently, this method can be either called with
a complete set of (F0, F1, F2, Alpha, Delta) (plus optional binary parameters),
OR a params dictionary;
and only the latter version will be supported going forward.
- Parameters:
F0 (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
F1 (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
F2 (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
Alpha (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
Delta (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
asini (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
period (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
ecc (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
tp (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
argp (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
params (dict, optional) – A dictionary defining a parameter space point. See get_fullycoherent_twoF() for more information.
Default:None
tstart (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime. This is only passed on to self.get_transient_detstats(), i.e. only used if self.transientWindowType is set.
Default:None
tend (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime. This is only passed on to self.get_transient_detstats(), i.e. only used if self.transientWindowType is set.
Default:None
- Returns:
stat – A single value of the main detection statistic at the input parameter values.
- Return type:
float
- get_fullycoherent_twoF(F0=None, F1=None, F2=None, Alpha=None, Delta=None, asini=None, period=None, ecc=None, tp=None, argp=None, params=None)[source]
Computes the fully-coherent 2F statistic at a single point.
NOTE: This always uses the full data set as defined when initialising the search object. If you want to restrict the range of data used for a single 2F computation, you need to set a self.transientWindowType and then call self.get_fullycoherent_detstat() with tstart and tend options instead of this function.
NOTE the old way of calling this with explicit (F0,F1,F2,Alpha,Delta,…) parameters is DEPRECATED and may be removed in future versions. Currently, this method can be either called with
a complete set of (F0, F1, F2, Alpha, Delta) (plus optional binary parameters),
OR a params dictionary;
and only the latter version will be supported going forward.
- Parameters:
F0 (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
F1 (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
F2 (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
Alpha (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
Delta (float, optional) – DEPRECATED: Parameters at which to compute the statistic.
Default:None
asini (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
period (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
ecc (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
tp (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
argp (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the statistic.
Default:None
params (dict, optional) – A dictionary defining a parameter space point, with [“F0”,”Alpha”,”Delta”] required as a minimum set of keys. Also supported: Fk with 1<k<lalpulsar.PULSAR_MAX_SPINS and binary parameters asini, period, ecc, tp, argp (basically anything that fits within a PulsarDopplerParams struct).
Default:None
- Returns:
twoF – A single value of the fully-coherent 2F statistic at the input parameter values. Also stored as self.twoF.
- Return type:
float
- get_fullycoherent_single_IFO_twoFs()[source]
Computes single-detector F-stats at a single point.
This requires self.get_fullycoherent_twoF() to be run first.
- Returns:
twoFX – A list of the single-detector detection statistics twoF. Also stored as self.twoFX.
- Return type:
list
- get_fullycoherent_log10BSGL()[source]
Computes the line-robust statistic log10BSGL at a single point.
This requires self.get_fullycoherent_twoF() and self.get_fullycoherent_single_IFO_twoFs() to be run first.
- Returns:
log10BSGL – A single value of the detection statistic log10BSGL at the input parameter values. Also stored as self.log10BSGL.
- Return type:
float
- get_transient_detstats(tstart=None, tend=None)[source]
Computes one or more transient detection statistics at a single point.
This requires self.get_fullycoherent_twoF() to be run first.
All computed statistics will be stored as attributes of self, but only one (twoF, log10BSGL or lnBtSG) will be the return value.
The full transient-F-stat map will also be computed here, but stored in self.FstatMap, not returned.
- Parameters:
tstart (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime.
Default:None
tend (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime.
Default:None
- Returns:
detstat – A single value of the main chosen detection statistic (maxTwoF, log10BSGL or lnBtSG) at the input parameter values.
- Return type:
float
- get_transient_maxTwoFstat(tstart=None, tend=None)[source]
Computes the transient maxTwoF statistic at a single point.
This requires self.get_fullycoherent_twoF() to be run first, and is itself now only a backwards compatibility / convenience wrapper around the more general get_transient_detstats().
The full transient-F-stat map will also be computed here, but stored in self.FstatMap, not returned.
- Parameters:
tstart (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime.
Default:None
tend (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime.
Default:None
- Returns:
maxTwoF – A single value of the detection statistic at the input parameter values. Also stored as self.maxTwoF.
- Return type:
float
- get_transient_log10BSGL()[source]
Computes a transient detection statistic log10BSGL at a single point.
This should normally be called through get_transient_detstats(), but if called stand-alone, it requires self.get_transient_maxTwoFstat() to be run first.
The single-detector 2F-stat values used for that computation (at the index of maxTwoF) are saved in self.twoFXatMaxTwoF, not returned.
- Returns:
log10BSGL – A single value of the detection statistic log10BSGL at the input parameter values. Also stored as self.log10BSGL.
- Return type:
float
- calculate_twoF_cumulative(F0=None, F1=None, F2=None, Alpha=None, Delta=None, asini=None, period=None, ecc=None, tp=None, argp=None, params=None, tstart=None, tend=None, transient_tstart=None, transient_duration=None, num_segments=1000)[source]
Calculate the cumulative twoF over subsets of the observation span.
This means that we consider sub-“segments” of the [tstart,tend] interval, each starting at the overall tstart and with increasing durations, and compute the 2F for each of these, which for a true CW signal should increase roughly with duration towards the full value.
- Parameters:
F0 (float, optional) – DEPRECATED: Parameters at which to compute the cumulative twoF.
Default:None
F1 (float, optional) – DEPRECATED: Parameters at which to compute the cumulative twoF.
Default:None
F2 (float, optional) – DEPRECATED: Parameters at which to compute the cumulative twoF.
Default:None
Alpha (float, optional) – DEPRECATED: Parameters at which to compute the cumulative twoF.
Default:None
Delta (float, optional) – DEPRECATED: Parameters at which to compute the cumulative twoF.
Default:None
asini (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the cumulative 2F.
Default:None
period (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the cumulative 2F.
Default:None
ecc (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the cumulative 2F.
Default:None
tp (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the cumulative 2F.
Default:None
argp (float, optional) – DEPRECATED: Optional: Binary parameters at which to compute the cumulative 2F.
Default:None
params (dict, optional) – A dictionary defining a parameter space point. See get_fullycoherent_twoF() for more information.
Default:None
tstart (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime;. If outside those: auto-truncated.
Default:None
tend (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime;. If outside those: auto-truncated.
Default:None
num_segments (int, optional) – Number of segments to split [tstart,tend] into.
Default:1000
transient_tstart (float or None, optional) – These are not actually used by this function, but just included so a parameters dict can be safely passed.
Default:None
transient_duration (float or None, optional) – These are not actually used by this function, but just included so a parameters dict can be safely passed.
Default:None
- Returns:
cumulative_durations (ndarray of shape (num_segments,)) – Offsets of each segment’s tend from the overall tstart.
twoFs (ndarray of shape (num_segments,)) – Values of twoF computed over [[tstart,tstart+duration] for duration in cumulative_durations].
- predict_twoF_cumulative(F0, Alpha, Delta, h0, cosi, psi, tstart=None, tend=None, num_segments=10, **predict_fstat_kwargs)[source]
Calculate expected 2F, with uncertainty, over subsets of the observation span.
This yields the expected behaviour that calculate_twoF_cumulative() can be compared against: 2F for CW signals increases with duration as we take longer and longer subsets of the total observation span.
- Parameters:
F0 (float) – Parameters at which to compute the cumulative predicted twoF.
Alpha (float) – Parameters at which to compute the cumulative predicted twoF.
Delta (float) – Parameters at which to compute the cumulative predicted twoF.
h0 (float) – Parameters at which to compute the cumulative predicted twoF.
cosi (float) – Parameters at which to compute the cumulative predicted twoF.
psi (float) – Parameters at which to compute the cumulative predicted twoF.
tstart (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime. If outside those: auto-truncated.
Default:None
tend (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime. If outside those: auto-truncated.
Default:None
num_segments (int, optional) – Number of segments to split [tstart,tend] into.
Default:10
predict_fstat_kwargs – Other kwargs to be passed to utils.predict_fstat().
- Returns:
tstart (int) – GPS start time of the observation span.
cumulative_durations (ndarray of shape (num_segments,)) – Offsets of each segment’s tend from the overall tstart.
pfs (ndarray of size (num_segments,)) – Predicted 2F for each segment.
pfs_sigma (ndarray of size (num_segments,)) – Standard deviations of predicted 2F.
- plot_twoF_cumulative(CFS_input, PFS_input=None, tstart=None, tend=None, num_segments_CFS=1000, num_segments_PFS=10, custom_ax_kwargs=None, savefig=False, label=None, outdir=None, **PFS_kwargs)[source]
Plot how 2F accumulates over time.
This compares the accumulation on the actual data set (‘CFS’, from self.calculate_twoF_cumulative()) against (optionally) the average expectation (‘PFS’, from self.predict_twoF_cumulative()).
- Parameters:
CFS_input (dict) – Input arguments for self.calculate_twoF_cumulative() (besides [tstart, tend, num_segments]).
PFS_input (dict, optional) – Input arguments for self.predict_twoF_cumulative() (besides [tstart, tend, num_segments]). If None: do not calculate predicted 2F.
Default:None
tstart (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime. If outside those: auto-truncated.
Default:None
tend (int or None, optional) – GPS times to restrict the range of data used. If None: falls back to self.minStartTime and self.maxStartTime. If outside those: auto-truncated.
Default:None
num_segments_(CFS|PFS) (int) – Number of time segments to (compute|predict) twoF.
custom_ax_kwargs (dict, optional) – Optional axis formatting options.
Default:None
savefig (bool, optional) – If true, save the figure in outdir. If false, return an axis object without saving to disk.
Default:False
label (str, optional) – Output filename will be constructed by appending _twoFcumulative.png to this label. (Ignored unless savefig=true.)
Default:None
outdir (str, optional) – Output folder (ignored unless savefig=true).
Default:None
PFS_kwargs (dict) – Other kwargs to be passed to self.predict_twoF_cumulative().
- Returns:
ax – The axes object containing the plot.
- Return type:
matplotlib.axes._subplots_AxesSubplot, optional
- write_atoms_to_file(fnamebase='', comments='%%')[source]
Save F-statistic atoms (time-dependent quantities) for a given parameter-space point.
- Parameters:
fnamebase (str, optional) – Basis for output filename, full name will be {fnamebase}_Fstatatoms_{dopplerName}.dat where dopplerName is a canonical lalpulsar formatting of the ‘Doppler’ parameter space point (frequency-evolution parameters).
Default:''
comments (str, optional) – Comments marker character(s) to be prepended to header lines. Note that the column headers line (last line of the header before the atoms data) is printed by lalpulsar, with %% as comments marker, so (different from most other PyFstat functions) the default here is %% too.
Default:'%%'
- class pyfstat.core.SemiCoherentSearch(