{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Install the PyCBC gravitational-wave analysis toolkit\n* https://github.com/gwastro/pycbc\n\nPyCBC is a python package developed by the community of GW astronomers to help analyze gravitational-wave data, detect signals, and even estimate the parameters of a source binary. It's meant to be accessible and welcomes contributions.","metadata":{}},{"cell_type":"code","source":"!pip install pycbc","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2021-07-04T13:10:55.288668Z","iopub.execute_input":"2021-07-04T13:10:55.289181Z","iopub.status.idle":"2021-07-04T13:11:21.893426Z","shell.execute_reply.started":"2021-07-04T13:10:55.289076Z","shell.execute_reply":"2021-07-04T13:11:21.892352Z"},"collapsed":true,"jupyter":{"outputs_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Get labels","metadata":{}},{"cell_type":"code","source":"import pandas as pd\ndf = pd.read_csv('../input/g2net-gravitational-wave-detection/training_labels.csv')\n\ndef train_filepath_from_id(id):\n    return f'../input/g2net-gravitational-wave-detection/train/{id[0]}/{id[1]}/{id[2]}/{id}.npy'\n\ndf['filepath'] = df.id.apply(train_filepath_from_id)\nprint(df.shape)","metadata":{"execution":{"iopub.status.busy":"2021-07-04T13:11:21.895156Z","iopub.execute_input":"2021-07-04T13:11:21.895448Z","iopub.status.idle":"2021-07-04T13:11:22.795258Z","shell.execute_reply.started":"2021-07-04T13:11:21.895417Z","shell.execute_reply":"2021-07-04T13:11:22.794029Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Plot the q-transform of the data\n\nA constant q-transform (https://en.wikipedia.org/wiki/Constant-Q_transform) is a common method to visualize gravitational-wave data. Almost all time-frequency figures in the gravitaitonal-wave literature use this method. ","metadata":{}},{"cell_type":"code","source":"import pylab, glob, os\nimport numpy as np\nimport pycbc.types\nfrom time import time\n\ndef plot_cqf(filepath, label):    \n    data = np.load(filepath)\n    \n    fig, axes = pylab.subplots(4, 1, figsize=[9,8], dpi=100)\n    pylab.title(f\"{label} / {os.path.basename(filepath)}\")\n    powers = list()\n    for i in range(3):\n        ts = pycbc.types.TimeSeries(data[i, :], epoch=0, delta_t=1.0/2048) \n        \n        # whiten the data (i.e. normalize the noise power at different frequencies)\n        ts = ts.whiten(0.125, 0.125, remove_corrupted=True)\n        \n        # calculate the qtransform\n        t, freq, power = ts.qtransform(.002, logfsteps=100, qrange=(10, 10), frange=(20, 512))\n        powers.append(power)\n        \n    #standardized_powers = [(p -np.mean(p))/np.std(p) for p in powers]\n    mean_of_all_sites = sum(powers)/3\n    powers.append(mean_of_all_sites)\n    \n    for i, power in enumerate(powers):\n        # convert the data to a TimeSeries instance\n        #print(time()-t0)\n        pylab.sca(axes[i])\n        pylab.pcolormesh(t, freq, power, vmax=15, vmin=0)\n        #pylab.xlim(.25, 1.75)\n        pylab.yscale('log')\n        #breakpoint()\n    pylab.show()\n    print()\n\nfor _, row in df[df.target == 1].sample(n=1).iterrows():\n    plot_cqf(row['filepath'], row['target'])\n    \n#for _, row in df[df.target == 0].sample(n=10).iterrows():\n#    plot_cqf(row['filepath'], row['target'])","metadata":{"execution":{"iopub.status.busy":"2021-07-04T13:11:45.623929Z","iopub.execute_input":"2021-07-04T13:11:45.624318Z","iopub.status.idle":"2021-07-04T13:11:51.245731Z","shell.execute_reply.started":"2021-07-04T13:11:45.624289Z","shell.execute_reply":"2021-07-04T13:11:51.244562Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from pathlib import Path\nfrom tqdm.notebook import tqdm\n\ndef compute_and_save_qtransform(filepath: str, output_directory: str):    \n    output_path = Path(output_directory) /  '/'.join(filepath.split('/')[-5:])\n    output_path.parent.mkdir(parents=True, exist_ok=True)\n    data = np.load(filepath)\n    \n    powers = list()\n    for i in range(3):\n        ts = pycbc.types.TimeSeries(data[i, :], epoch=0, delta_t=1.0/2048) \n        ts = ts.whiten(0.125, 0.125, remove_corrupted=True)\n        t, freq, power = ts.qtransform(.002, logfsteps=100, qrange=(10, 10), frange=(20, 512))\n        powers.append(power)\n    all_powers = np.array(powers)\n    np.save(output_path, all_powers)\n    \nfrom joblib import Parallel, delayed\nParallel(n_jobs=2)(delayed(compute_and_save_qtransform)(filepath, '/kaggle/working/') \n                   for filepath in tqdm(df['filepath'].values[:50]))","metadata":{"execution":{"iopub.status.busy":"2021-07-04T13:27:00.901900Z","iopub.execute_input":"2021-07-04T13:27:00.902292Z","iopub.status.idle":"2021-07-04T13:27:09.022627Z","shell.execute_reply.started":"2021-07-04T13:27:00.902249Z","shell.execute_reply":"2021-07-04T13:27:09.021848Z"},"trusted":true},"execution_count":null,"outputs":[]}]}