{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.7.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":37077,"databundleVersionId":4333111,"sourceType":"competition"}],"dockerImageVersionId":30301,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport os\nimport h5py\nimport matplotlib.pyplot as plt\nimport tqdm\nimport pickle","metadata":{"execution":{"iopub.status.busy":"2024-06-07T07:33:59.013456Z","iopub.execute_input":"2024-06-07T07:33:59.014709Z","iopub.status.idle":"2024-06-07T07:33:59.022648Z","shell.execute_reply.started":"2024-06-07T07:33:59.014633Z","shell.execute_reply":"2024-06-07T07:33:59.020944Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Check the contents of a HDF file","metadata":{}},{"cell_type":"code","source":"with h5py.File('../input/g2net-detecting-continuous-gravitational-waves/train/001121a05.hdf5','r') as f:\n    for k1 in f:\n        print(\"0:\",k1)\n        for k2 in f[k1]:\n            print(\"1:\",f[k1][k2])\n            for k3 in f[k1][k2]:\n                try:\n                    print(\"2:\",f[k1][k2][k3])\n                except TypeError:\n                    pass","metadata":{"execution":{"iopub.status.busy":"2024-06-07T07:33:59.025276Z","iopub.execute_input":"2024-06-07T07:33:59.025787Z","iopub.status.idle":"2024-06-07T07:33:59.177633Z","shell.execute_reply.started":"2024-06-07T07:33:59.025728Z","shell.execute_reply":"2024-06-07T07:33:59.176024Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Extract abstract of train/test data\nI extract\n1. min/max frequency\n2. H1/L1 SFTs sizes \n3. H1/L1 SFTs timestamps_GPS min/max\n4. H1/L1 SFTs timestamps_GPS maximum interval \n","metadata":{}},{"cell_type":"code","source":"df_train = pd.DataFrame(columns = ['freq min','freq max','H1 SFTs h','H1 SFTs w','H1 time min','H1 time max','H1 max inteval','L1 SFTs h','L1 SFTs w','L1 time min','L1 time max','L1 max inteval'])\ndf_test =  pd.DataFrame(columns = ['freq min','freq max','H1 SFTs h','H1 SFTs w','H1 time min','H1 time max','H1 max inteval','L1 SFTs h','L1 SFTs w','L1 time min','L1 time max','L1 max inteval'])\n","metadata":{"execution":{"iopub.status.busy":"2024-06-07T07:33:59.179164Z","iopub.execute_input":"2024-06-07T07:33:59.179527Z","iopub.status.idle":"2024-06-07T07:33:59.191716Z","shell.execute_reply.started":"2024-06-07T07:33:59.179496Z","shell.execute_reply":"2024-06-07T07:33:59.190286Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# For train data","metadata":{}},{"cell_type":"code","source":"dirname = '/kaggle/input/g2net-detecting-continuous-gravitational-waves/train/'\ntrainfiles = sorted(os.listdir(dirname))\nfilenames = trainfiles # + testfiles\ntimestamps_train_H1 = []\ntimestamps_train_L1 = []\nprint(f\"train files = {len(trainfiles)}\")\nfor filename in tqdm.tqdm(filenames):\n    filepath = os.path.join(dirname, filename)\n    file_id = filename.split('.')[0]\n    f = h5py.File(filepath,'r')\n    H1tinterval = np.array(f[file_id]['H1'][\"timestamps_GPS\"])\n    timestamps_train_H1.append(H1tinterval) # save timestamps (not interval)\n    H1tinterval = H1tinterval[1:]-H1tinterval[0:-1]\n    L1tinterval = np.array(f[file_id]['L1'][\"timestamps_GPS\"])\n    timestamps_train_L1.append(L1tinterval) # save timestamps (not interval)\n    L1tinterval = L1tinterval[1:]-L1tinterval[0:-1]\n    df_train.loc[file_id] =  [min(f[file_id]['frequency_Hz']),max(f[file_id]['frequency_Hz']),\n                    f[file_id]['H1'][\"SFTs\"].shape[0],f[file_id]['H1'][\"SFTs\"].shape[1],\n                              min(f[file_id]['H1'][\"timestamps_GPS\"]),max(f[file_id]['H1'][\"timestamps_GPS\"]),max(H1tinterval),\n                    f[file_id]['L1'][\"SFTs\"].shape[0],f[file_id]['L1'][\"SFTs\"].shape[1],\n                              min(f[file_id]['L1'][\"timestamps_GPS\"]),max(f[file_id]['L1'][\"timestamps_GPS\"]),max(H1tinterval)]\n    f.close()\ndf_train.to_csv('train_summary.csv')\nwith open(\"timestamps_train_H1.pickle\", 'wb') as f: pickle.dump(timestamps_train_H1, f)\nwith open(\"timestamps_train_L1.pickle\", 'wb') as f: pickle.dump(timestamps_train_L1, f)","metadata":{"execution":{"iopub.status.busy":"2024-06-07T07:33:59.194608Z","iopub.execute_input":"2024-06-07T07:33:59.195048Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# For test data","metadata":{}},{"cell_type":"code","source":"dirname = '/kaggle/input/g2net-detecting-continuous-gravitational-waves/test/'\ntestfiles = sorted(os.listdir(dirname))\nfilenames = testfiles \ntimestamps_test_H1 = []\ntimestamps_test_L1 = []\nprint(f\"test files = {len(testfiles)}\")\nfor filename in tqdm.tqdm(filenames):\n    filepath = os.path.join(dirname, filename)\n    file_id = filename.split('.')[0]\n    f = h5py.File(filepath,'r')\n    H1tinterval = np.array(f[file_id]['H1'][\"timestamps_GPS\"])\n    timestamps_test_H1.append(H1tinterval) # save timestamps (not interval)\n    H1tinterval = H1tinterval[1:]-H1tinterval[0:-1]\n    L1tinterval = np.array(f[file_id]['L1'][\"timestamps_GPS\"])\n    timestamps_test_L1.append(L1tinterval) # save timestamps (not interval)\n    L1tinterval = L1tinterval[1:]-L1tinterval[0:-1]\n    df_test.loc[file_id] =  [min(f[file_id]['frequency_Hz']),max(f[file_id]['frequency_Hz']),\n                    f[file_id]['H1'][\"SFTs\"].shape[0],f[file_id]['H1'][\"SFTs\"].shape[1],\n                              min(f[file_id]['H1'][\"timestamps_GPS\"]),max(f[file_id]['H1'][\"timestamps_GPS\"]),max(H1tinterval),\n                    f[file_id]['L1'][\"SFTs\"].shape[0],f[file_id]['L1'][\"SFTs\"].shape[1],\n                              min(f[file_id]['L1'][\"timestamps_GPS\"]),max(f[file_id]['L1'][\"timestamps_GPS\"]),max(H1tinterval)]\n    f.close()\ndf_test.to_csv('test_summary.csv')\nwith open(\"timestamps_test_H1.pickle\", 'wb') as f: pickle.dump(timestamps_test_H1, f)\nwith open(\"timestamps_test_L1.pickle\", 'wb') as f: pickle.dump(timestamps_test_L1, f)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Check the extracted data","metadata":{}},{"cell_type":"code","source":"df_train","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_test","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"SFTs h(frequency) size seems to be all 360, so check it.","metadata":{}},{"cell_type":"code","source":"print(df_train['H1 SFTs h'].unique())\nprint(df_train['L1 SFTs h'].unique())\nprint(df_test['H1 SFTs h'].unique())\nprint(df_test['L1 SFTs h'].unique())\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Distribution of SFTs w(time) size","metadata":{}},{"cell_type":"markdown","source":"SFTs w(time) sizes are different from file to file. So, check the distribution.","metadata":{}},{"cell_type":"code","source":"fig,ax = plt.subplots(2,2,figsize=(12,4))\nplt.subplots_adjust(hspace=0.6)\n_ = ax[0,0].hist(df_train['H1 SFTs w'],range(4000,5000,10))\n_ = ax[0,0].set_title('H1 SFTs w-size (Train)')\n_ = ax[0,1].hist(df_train['L1 SFTs w'],range(4000,5000,10))\n_ = ax[0,1].set_title('L1 SFTs w-size (Train)')\n_ = ax[1,0].hist(df_test['H1 SFTs w'],range(4000,5000,10))\n_ = ax[1,0].set_title('H1 SFTs w-size (Test)')\n_ = ax[1,1].hist(df_test['L1 SFTs w'],range(4000,5000,10))\n_ = ax[1,1].set_title('L1 SFTs w-size (Test)')\n\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Distribution of SFTs min/max frequency","metadata":{}},{"cell_type":"code","source":"fig,ax = plt.subplots(2,2,figsize=(12,4))\nplt.subplots_adjust(hspace=0.6)\n_ = ax[0,0].hist(df_train['freq min'] ,100)\n_ = ax[0,0].set_title('Frequency min [Hz] (Train)')\n_ = ax[0,1].hist(df_train['freq max'] ,100)\n_ = ax[0,1].set_title('Frequency max [Hz] (Train)')\n#_ = ax[0,2].hist(df_train['freq max']-df_train['freq min']  ,100)\n#_ = ax[0,2].set_title('Frequency span [Hz] (Train)')\n\n_ = ax[1,0].hist(df_test['freq min'] ,100)\n_ = ax[1,0].set_title('Frequency min [Hz] (Test)')\n_ = ax[1,1].hist(df_test['freq max'] ,100)\n_ = ax[1,1].set_title('Frequency max [Hz] (Test)')\n#_ = ax[1,2].hist(df_test['freq max']-df_test['freq min']  ,100)\n#_ = ax[1,2].set_title('Frequency span [Hz] (Test)')\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Target's frequency dependency","metadata":{}},{"cell_type":"markdown","source":"Now that I found that frequency window is different from data to data, I should see them grouped by with/without GW(gravitational wave).","metadata":{}},{"cell_type":"code","source":"df_label = pd.read_csv('../input/g2net-detecting-continuous-gravitational-waves/train_labels.csv')\n(df_label.id != df_train.index).sum()\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Be careful. There are -1 in labels as written in the Dataset Description, which means 'currently unable to determine the status'.","metadata":{"execution":{"iopub.status.busy":"2022-10-29T19:38:40.03093Z","iopub.execute_input":"2022-10-29T19:38:40.031235Z","iopub.status.idle":"2022-10-29T19:38:40.038893Z","shell.execute_reply.started":"2022-10-29T19:38:40.031206Z","shell.execute_reply":"2022-10-29T19:38:40.037286Z"}}},{"cell_type":"code","source":"df_label.groupby('target').count()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig , ax = plt.subplots(1,2,figsize=(12,4))\n_ = ax[0].hist(df_train.iloc[df_label.target.values==1]['freq min'],100,label='With GW')\n_ = ax[0].hist(df_train.iloc[df_label.target.values==0]['freq min'],100,label='W/O GW')\n_ = ax[0].legend()\n_ = ax[0].set_title('Frequency min [Hz] (Train)')\n_ = ax[1].hist(df_train.iloc[df_label.target.values==0]['freq min'],100,label='W/O GW')\n_ = ax[1].hist(df_train.iloc[df_label.target.values==1]['freq min'],100,label='With GW')\n_ = ax[1].legend()\n_ = ax[1].set_title('Frequency min [Hz] (Train), Swap order')\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We have to notice, there are only 'with GW' data or 'without GW' data for some frequencies. (Again, frequency windows are just 0.2 Hz)\nLet's magnify some frequency.","metadata":{}},{"cell_type":"code","source":"fig , ax = plt.subplots(2,1,figsize=(12,4))\nplt.subplots_adjust(hspace=0.6)\n_ = ax[0].hist(df_train.iloc[df_label.target.values==1]['freq min'],bins = np.linspace(100,120,100),label='With GW')\n_ = ax[0].hist(df_train.iloc[df_label.target.values==0]['freq min'],bins = np.linspace(100,120,100),label='W/O GW')\n_ = ax[0].legend()\n_ = ax[0].set_title('Frequency min [Hz] (Train)')\n_ = ax[1].hist(df_test['freq min'] ,bins = np.linspace(100,120,100),label='we have to predict With/Without GW')\n_ = ax[1].legend()\n_ = ax[1].set_title('Frequency min [Hz] (Test)')\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"With freqency windows being 0.2Hz, SFTs's frequencies do not overwrap. We have to make predictions without corresponding frequency pairs if we use only provided training dataset.","metadata":{}},{"cell_type":"markdown","source":"# Distribution of timestamps_GPS range\nAlso check the spans of timestamps_GPS. These are narrower than SFTs w-size.","metadata":{}},{"cell_type":"code","source":"fig,ax = plt.subplots(2,2,figsize=(12,4))\nplt.subplots_adjust(hspace=0.6)\n_ = ax[0,0].hist(df_train['H1 time max'] - df_train['H1 time min'],bins = np.linspace(10e6,12.5e6,100))\n_ = ax[0,0].set_title('H1 timestamps_GPS span (Train)')\n_ = ax[0,1].hist(df_train['L1 time max'] - df_train['L1 time min'],bins = np.linspace(10e6,12.5e6,100))\n_ = ax[0,1].set_title('L1 timestamps_GPS span (Train)')\n_ = ax[1,0].hist(df_test['H1 time max'] - df_test['H1 time min'],bins = np.linspace(10e6,12.5e6,100))\n_ = ax[1,0].set_title('H1 timestamps_GPS span (Test)')\n_ = ax[1,1].hist(df_test['L1 time max'] - df_test['L1 time min'],bins = np.linspace(10e6,12.5e6,100))\n_ = ax[1,1].set_title('L1 timestamps_GPS span (Test)')\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Gaps in timestamps_GPS\nThere are gaps in timestamps_GPS due to detectors' maintenacne or connection being offline, etc. Let's look at how gaps are.<br>\nFirst, look at SFTs' intervals of three data from train and test data.","metadata":{}},{"cell_type":"code","source":"for dirname in ['/kaggle/input/g2net-detecting-continuous-gravitational-waves/train/', \n                '/kaggle/input/g2net-detecting-continuous-gravitational-waves/test/']:\n    print(dirname)\n    trainfiles = sorted(os.listdir(dirname))\n    filenames = trainfiles # + testfiles\n    for filename in (filenames[:3]):\n        filepath = os.path.join(dirname, filename)\n        file_id = filename.split('.')[0]\n        f = h5py.File(filepath,'r')\n        H1tinterval = np.array(f[file_id]['H1'][\"timestamps_GPS\"])\n        H1tinterval = H1tinterval[1:]-H1tinterval[0:-1]\n        fig,ax = plt.subplots(2,2,figsize=(10,4))\n        plt.subplots_adjust(hspace=0.6)\n        ax[0,0].plot(H1tinterval)\n        ax[0,0].set_title(f\"{file_id} H1 SFTs interval (Time series)\")\n        ax[0,1].hist(H1tinterval,50)\n        ax[0,1].set_title(f\"{file_id} H1 SFTs interval (Distribution)\")\n\n        L1tinterval = np.array(f[file_id]['L1'][\"timestamps_GPS\"])\n        L1tinterval = L1tinterval[1:]-L1tinterval[0:-1]\n        ax[1,0].plot(L1tinterval)\n        ax[1,0].set_title(f\"{file_id} L1 SFTs interval (Time series)\")\n        ax[1,1].hist(L1tinterval,50)\n        ax[1,1].set_title(f\"{file_id} L1 SFTs interval (Distribution)\")\n        plt.show()\n        print(f\"{file_id} H1 first 5 intervals = {H1tinterval[:5]}  max interval = {H1tinterval.max()}\\n\")\n        print(f\"{file_id} L1 first 5 intervals = {L1tinterval[:5]}  max interval = {L1tinterval.max()}\\n\")\n        f.close()\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Typical intervals are 1800, but all data above have intervals greater than 1800. Having gaps are very common and we have to take this into consideration when we generate training data by ourselves with PyFstat.\nFinally, I plot the distribution of maximum intervals for the train and test dataset. Each data usually have gaps 20 to 30 times longer than typical interval, but sometimes have 50 times longer.","metadata":{}},{"cell_type":"markdown","source":"# Distribution of maximum gap","metadata":{}},{"cell_type":"code","source":"fig,ax = plt.subplots(2,2,figsize=(12,4))\nplt.subplots_adjust(hspace=0.6)\n_ = ax[0,0].hist(df_train['H1 max inteval'] ,bins=np.linspace(0,100000,100))\n_ = ax[0,0].set_title('Max interval (Train)')\n_ = ax[0,1].hist(df_train['L1 max inteval'] ,bins=np.linspace(0,100000,100))\n_ = ax[0,1].set_title('Max interval (Train)')\n\n_ = ax[1,0].hist(df_test['H1 max inteval'] ,bins=np.linspace(0,100000,100))\n_ = ax[1,0].set_title('Max interval (Test)')\n_ = ax[1,1].hist(df_test['L1 max inteval'] ,bins=np.linspace(0,100000,100))\n_ = ax[1,1].set_title('Max interval (Test)')\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}