diff --git a/examples/grid_glitch_search.py b/examples/grid_glitch_search.py new file mode 100644 index 0000000000000000000000000000000000000000..1d6fc7587872f4d3fc025ced3571e39531b27af6 --- /dev/null +++ b/examples/grid_glitch_search.py @@ -0,0 +1,43 @@ +from pyfstat import GridGlitchSearch, Writer +import numpy as np + +phi = 0 +F0 = 2.4343 +F1 = -1e-10 +F2 = 0 +delta_phi = 0 +delta_F0 = 1e-7 +delta_F1 = 0 +duration = 100*86400 +dtglitch = 100*43200 +Alpha = 5e-3 +Delta = 6e-2 +tstart = 700000000 +tend = tstart+duration +tglitch = tstart + dtglitch +tref = tglitch + +glitch_data = Writer(label='glitch', outdir='data', delta_phi=delta_phi, + delta_F0=delta_F0, tref=tref, tstart=tstart, + delta_F1=delta_F1, phi=phi, F0=F0, F1=F1, F2=F2, + duration=duration, dtglitch=dtglitch, Alpha=Alpha, + Delta=Delta) +glitch_data.make_data() + +F0s = [F0] +F1s = [F1] +F2s = [F2] +m = 1e-3 +dF = np.sqrt(12 * m) / (np.pi * duration) +delta_F0s = [-1e-6*F0, 1e-6*F0, dF] +delta_F1s = [delta_F1] +dT = duration / 10. +tglitchs = [tstart+dT, tend-dT, duration/100.] +Alphas = [Alpha] +Deltas = [Delta] +grid = GridGlitchSearch('grid', 'data', glitch_data.label, glitch_data.outdir, + F0s, F1s, F2s, delta_F0s, delta_F1s, tglitchs, Alphas, + Deltas, tref, tstart, tend) +grid.run() +grid.plot_2D('delta_F0', 'tglitch') +print grid.get_max_twoF() diff --git a/pyfstat.py b/pyfstat.py index f0f69acf4069a6beb70e1ca74b1aaa7e009965ec..639f0d1ca21df962f087461cb5e6a5e0c1699c4e 100755 --- a/pyfstat.py +++ b/pyfstat.py @@ -169,7 +169,9 @@ class ComputeFstat(object): self.sft_filepath = self.sftdir+'/*_'+self.sftlabel+"*sft" SFTCatalog = lalpulsar.SFTdataFind(self.sft_filepath, constraints) names = list(set([d.header.name for d in SFTCatalog.data])) - logging.info('Loaded data from detectors {}'.format(names)) + logging.info( + 'Loaded data from detectors {} matching pattern {}'.format( + names, self.sft_filepath)) logging.info('Initialising ephems') ephems = lalpulsar.InitBarycenter(self.earth_ephem, self.sun_ephem) @@ -185,6 +187,9 @@ class ComputeFstat(object): fB = fA + (numBins-1)*SFTCatalog.data[0].header.deltaF self.minCoverFreq = fA + 0.5 self.maxCoverFreq = fB - 0.5 + logging.info('Min/max cover freqs not provided, using ' + '{} and {}, est. from SFTs'.format( + self.minCoverFreq, self.maxCoverFreq)) self.FstatInput = lalpulsar.CreateFstatInput(SFTCatalog, self.minCoverFreq, @@ -347,9 +352,9 @@ class MCMCGlitchSearch(BaseSearchClass): @initializer def __init__(self, label, outdir, sftlabel, sftdir, theta, tref, tstart, tend, nsteps=[100, 100, 100], nwalkers=100, ntemps=1, - nglitch=0, minCoverFreq=29, maxCoverFreq=31, scatter_val=1e-4, - betas=None, detector=None, dtglitchmin=20*86400, - earth_ephem=None, sun_ephem=None): + nglitch=0, minCoverFreq=None, maxCoverFreq=None, + scatter_val=1e-4, betas=None, detector=None, + dtglitchmin=20*86400, earth_ephem=None, sun_ephem=None): """ Parameters label, outdir: str @@ -904,7 +909,7 @@ class GridGlitchSearch(BaseSearchClass): def __init__(self, label, outdir, sftlabel=None, sftdir=None, F0s=[0], F1s=[0], F2s=[0], delta_F0s=[0], delta_F1s=[0], tglitchs=None, Alphas=[0], Deltas=[0], tref=None, tstart=None, tend=None, - minCoverFreq=29, maxCoverFreq=31, write_after=1000, + minCoverFreq=None, maxCoverFreq=None, write_after=1000, earth_ephem=None, sun_ephem=None): """ Parameters @@ -946,8 +951,8 @@ class GridGlitchSearch(BaseSearchClass): if os.path.isdir(outdir) is False: os.mkdir(outdir) self.out_file = '{}/{}_gridFS.txt'.format(self.outdir, self.label) - self.keys = ['F0', 'F1', 'F2', 'delta_F0', 'delta_F1', 'tglitch', - 'Alpha', 'Delta'] + self.keys = ['F0', 'F1', 'F2', 'Alpha', 'Delta', 'delta_F0', + 'delta_F1', 'tglitch'] def get_array_from_tuple(self, x): if len(x) == 1: @@ -995,7 +1000,7 @@ class GridGlitchSearch(BaseSearchClass): counter = 0 data = [] for vals in self.input_data: - FS = self.search.compute_glitch_fstat(*vals) + FS = self.search.compute_glitch_fstat_single(*vals) data.append(list(vals) + [FS]) if counter > self.write_after: