diff --git a/pycbc/inference/models/gated_gaussian_noise.py b/pycbc/inference/models/gated_gaussian_noise.py index 8852ecc8583..33e7ecc27e1 100644 --- a/pycbc/inference/models/gated_gaussian_noise.py +++ b/pycbc/inference/models/gated_gaussian_noise.py @@ -46,9 +46,6 @@ def __init__(self, variable_params, data, low_frequency_cutoff, psds=None, static_params=None, highpass_waveforms=False, **kwargs): # we'll want the time-domain data, so store that self._td_data = {} - # cache the current projection for debugging - self.current_proj = {} - self.current_nproj = {} # cache the overwhitened data self._overwhitened_data = {} # cache the current gated data @@ -383,7 +380,6 @@ def _lognl(self): self._det_lognls.clear() # get the times of the gates gate_times = self.get_gate_times() - self.current_nproj.clear() for det, invpsd in self._invpsds.items(): norm = self.det_lognorm(det) gatestartdelay, dgatedelay = gate_times[det] @@ -395,7 +391,6 @@ def _lognl(self): gated_dt = data.gate(gatestartdelay + dgatedelay/2, window=dgatedelay/2, copy=True, invpsd=invpsd, method='paint') - self.current_nproj[det] = (gated_dt.proj, gated_dt.projslc) # convert to the frequency series gated_d = gated_dt.to_frequencyseries() # overwhiten @@ -546,7 +541,6 @@ def _loglikelihood(self): # get the times of the gates gate_times = self.get_gate_times() logl = 0. - self.current_proj.clear() for det, h in wfs.items(): invpsd = self._invpsds[det] norm = self.det_lognorm(det) @@ -562,7 +556,6 @@ def _loglikelihood(self): gated_res = res.gate(gatestartdelay + dgatedelay/2, window=dgatedelay/2, copy=True, invpsd=invpsd, method='paint') - self.current_proj[det] = (gated_res.proj, gated_res.projslc) gated_rtilde = gated_res.to_frequencyseries() # overwhiten gated_rtilde *= invpsd diff --git a/pycbc/inference/models/marginalized_gaussian_noise.py b/pycbc/inference/models/marginalized_gaussian_noise.py index 6db6ade80e4..0a1342d0432 100644 --- a/pycbc/inference/models/marginalized_gaussian_noise.py +++ b/pycbc/inference/models/marginalized_gaussian_noise.py @@ -229,7 +229,7 @@ def __init__(self, variable_params, generator_class=generator.FDomainDetFrameTwoPolNoRespGenerator, gates=self.gates, **kwargs['static_params']) else: - # create a waveform generator for each ifo respestively + # create a waveform generator for each ifo respectively self.waveform_generator = {} for det in self.data: self.waveform_generator[det] = create_waveform_generator( @@ -386,7 +386,7 @@ def __init__(self, variable_params, data, low_frequency_cutoff, psds=None, generator_class=generator.FDomainDetFrameTwoPolGenerator, gates=self.gates, **kwargs['static_params']) else: - # create a waveform generator for each ifo respestively + # create a waveform generator for each ifo respectively self.waveform_generator = {} for det in self.data: self.waveform_generator[det] = create_waveform_generator( diff --git a/pycbc/strain/gate.py b/pycbc/strain/gate.py index e294beeae87..6c2a70f59c3 100644 --- a/pycbc/strain/gate.py +++ b/pycbc/strain/gate.py @@ -168,6 +168,7 @@ def gate_and_paint(data, lindex, rindex, invpsd, copy=True): # Copy the data and zero inside the hole if copy: data = data.copy() + # Here's ambiguity about when gate end time exactly is, rindex-1 or rindex? data[lindex:rindex] = 0 # get the over-whitened gated data @@ -178,6 +179,4 @@ def gate_and_paint(data, lindex, rindex, invpsd, copy=True): proj = linalg.solve_toeplitz(tdfilter[:(rindex - lindex)], owhgated_data[lindex:rindex]) data[lindex:rindex] -= proj - data.projslc = (lindex, rindex) - data.proj = proj return data diff --git a/pycbc/types/timeseries.py b/pycbc/types/timeseries.py index 910922065ad..3bf94065473 100644 --- a/pycbc/types/timeseries.py +++ b/pycbc/types/timeseries.py @@ -595,15 +595,31 @@ def gate(self, time, window=0.25, method='taper', copy=True, # Uses the hole-filling method of # https://arxiv.org/pdf/1908.05644.pdf from pycbc.strain.gate import gate_and_paint + from pycbc.waveform.utils import apply_fd_time_shift if invpsd is None: # These are some bare minimum settings, normally you # should probably provide a psd invpsd = 1. / self.filter_psd(self.duration/32, self.delta_f, 0) lindex = int((time - window - self.start_time) / self.delta_t) - rindex = lindex + int(2 * window / self.delta_t) + rindex = int((time + window - self.start_time) / self.delta_t) lindex = lindex if lindex >= 0 else 0 rindex = rindex if rindex <= len(self) else len(self) - return gate_and_paint(data, lindex, rindex, invpsd, copy=False) + rindex_time = float(self.start_time + rindex * self.delta_t) + offset = rindex_time - (time + window) + if offset == 0: + return gate_and_paint(data, lindex, rindex, invpsd, copy=False) + else: + # time shift such that gate end time lands on a specific data sample + fdata = data.to_frequencyseries() + fdata = apply_fd_time_shift(fdata, offset + fdata.epoch, copy=False) + # gate and paint in time domain + data = fdata.to_timeseries() + data = gate_and_paint(data, lindex, rindex, invpsd, copy=False) + # shift back to the original time + fdata = data.to_frequencyseries() + fdata = apply_fd_time_shift(fdata, -offset + fdata.epoch, copy=False) + tdata = fdata.to_timeseries() + return tdata elif method == 'hard': tslice = data.time_slice(time - window, time + window) tslice[:] = 0