Reading gating files and applying the documented gates

As a pre-processing step, the PyCBC search uses an inverted-Tukey window to miitgate the effect of loud, non-Gaussian features in the data. This is further described in:

Usman et al 2016 Class. Quantum Grav. 33 215004 - https://dx.doi.org/10.1088/0264-9381/33/21/215004

A subset of these times are the times listed in the txt files

  • H1-O3_GATES_1238166018-31197600.txt
  • L1-O3_GATES_1238166018-31197600.txt

These times in these files were chosen based on auxiliary monitors of overflows in the digital-to-analog converters used to control the positions of the test masses. The gated times (i.e. the time period where the data is zeroed) are time segments where these monitors recorded an overflow were. The "central time" and "suggested half-width of zero time" were chosen to fully cover these time seconds. The final gating parameter, the "suggested taper time" was chosen to be 0.5 to balance the cost of impacting more data with the window function versus introducing additional artifacts into the data.

The syntax of the files themselves is

{central time} {suggested half-width of zero time} {suggested taper time}

with each row containing the parameters of a single gate.

This notebook provides an example of how to read in and apply one of the suggested gates.

In [1]:
import numpy
import pycbc.frame
import matplotlib.pyplot as plt
%matplotlib inline
In [2]:
# read in the list of gate parameters
gate_times, gate_widths, gate_tapers = numpy.loadtxt('H1-O3_GATES_1238166018-31197600.txt',unpack=True)

# pick parameters for a specific gate
gate_num = 29
gate_time = gate_times[gate_num]
gate_width = gate_widths[gate_num]
gate_taper = gate_tapers[gate_num]

print(gate_time, gate_width, gate_taper)

1238434615.41 0.34375 0.5
In [3]:
# download the required data file

!wget https://www.gw-openscience.org/archive/data/O3a_4KHZ_R1/1238368256/H-H1_GWOSC_O3a_4KHZ_R1-1238433792-4096.gwf

--2021-10-19 15:27:09--  https://www.gw-openscience.org/archive/data/O3a_4KHZ_R1/1238368256/H-H1_GWOSC_O3a_4KHZ_R1-1238433792-4096.gwf
Resolving www.gw-openscience.org (www.gw-openscience.org)... 131.215.113.73
Connecting to www.gw-openscience.org (www.gw-openscience.org)|131.215.113.73|:443... connected.
HTTP request sent, awaiting response... 301 Moved Permanently
Location: /archive/data/O3a_4KHZ_R1/1238368256/H-H1_GWOSC_O3a_4KHZ_R1-1238433792-4096.gwf/ [following]
--2021-10-19 15:27:09--  https://www.gw-openscience.org/archive/data/O3a_4KHZ_R1/1238368256/H-H1_GWOSC_O3a_4KHZ_R1-1238433792-4096.gwf/
Reusing existing connection to www.gw-openscience.org:443.
HTTP request sent, awaiting response... 200 OK
Length: 129377846 (123M) [application/octet-stream]
Saving to: ‘H-H1_GWOSC_O3a_4KHZ_R1-1238433792-4096.gwf.1’

H-H1_GWOSC_O3a_4KHZ 100%[===================>] 123.38M   122MB/s    in 1.0s    

2021-10-19 15:27:10 (122 MB/s) - ‘H-H1_GWOSC_O3a_4KHZ_R1-1238433792-4096.gwf.1’ saved [129377846/129377846]
In [4]:
# read the strain data

data = pycbc.frame.read_frame('H-H1_GWOSC_O3a_4KHZ_R1-1238433792-4096.gwf', 'H1:GWOSC-4KHZ_R1_STRAIN', 
                              int(gate_time)-128, int(gate_time)+128)
In [5]:
# plot the original data

plot = plt.figure(figsize=(12,4))
ax = plot.gca()
ax.plot(data.sample_times-gate_time, data, label='Raw data')
ax.legend()
ax.set_xlim(-4, 4)
ax.set_xlabel('Time [seconds] from center of gate')
ax.set_ylabel('[strain]')
ax.grid(True)
plot.show()
In [6]:
# plot a spectrogram of the original data using the q-transform

whitened_data = data.whiten(4,2).time_slice(gate_time-4, gate_time+4)
times, freqs, qvals = whitened_data.qtransform(qrange=[20,20],frange=[10,500])

plot = plt.figure(figsize=(12,4))
ax = plot.gca()
im = ax.pcolormesh(times - gate_time, freqs, qvals, cmap='viridis', shading='auto', vmin=0,vmax=25)

ax.set_yscale('log')
ax.set_xlabel('Time [seconds] from center of gate')
ax.set_ylim(20, 500)
ax.set_ylabel('Frequency [Hz]')
ax.grid(True, axis='both', which='both', alpha=0.75)
plot.colorbar(im, label='Normalized energy')
plot.show()
In [7]:
# gate the data

gated_data = data.gate(time=gate_time, window=gate_width, taper_width=gate_taper, method='taper')
In [8]:
# plot the gated data

plot = plt.figure(figsize=(12,4))
ax = plot.gca()
ax.plot(gated_data.sample_times-gate_time, gated_data, label='Raw data, gated')
ax.legend()
ax.set_xlim(-4, 4)
ax.set_xlabel('Time [seconds] from center of gate')
ax.set_ylabel('[strain]')
ax.grid(True)
plot.show()
In [9]:
# plot a spectrogram of the gated data using the q-transform

whitened_gated_data = gated_data.whiten(4,2).time_slice(gate_time-4, gate_time+4)
times, freqs, qvals = whitened_gated_data.qtransform(qrange=[20,20],frange=[10,500])

plot = plt.figure(figsize=(12,4))
ax = plot.gca()
im = ax.pcolormesh(times - gate_time, freqs, qvals, cmap='viridis', shading='auto', vmin=0,vmax=25)

ax.set_xlabel('Time [seconds] from center of gate')
ax.set_yscale('log')
ax.set_ylim(20, 500)
ax.set_ylabel('Frequency [Hz]')
ax.grid(True, axis='both', which='both', alpha=0.75)
plot.colorbar(im, label='Normalized energy')
plot.show()
In [10]:
# also plot a spectrogram of the gated data using the q-transform in gwpy

from gwpy.timeseries import TimeSeries

gated_data_gwpy = TimeSeries.from_pycbc(gated_data)

qscan = gated_data_gwpy.crop(gate_time-32,gate_time+32)\
                  .q_transform(outseg=(gate_time-4, gate_time+4), qrange=(20,20))

plot = qscan.plot(figsize=(12,4))
ax = plot.gca()
ax.set_xscale('seconds')
ax.set_yscale('log')
ax.set_ylim(20, 500)
ax.set_ylabel('Frequency [Hz]')
ax.grid(True, axis='y', which='both')
ax.colorbar(cmap='viridis', label='Normalized energy',vmin=0,vmax=25)
plot.show()