From 2625c7b68325438254c638082f5b68a68dc185f5 Mon Sep 17 00:00:00 2001 From: Wei Changfeng <2681968849@qq.com> Date: Tue, 22 Oct 2019 13:21:02 +0100 Subject: [PATCH] Upload New File --- bbh_test_strain.ipynb | 193 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 193 insertions(+) create mode 100644 bbh_test_strain.ipynb diff --git a/bbh_test_strain.ipynb b/bbh_test_strain.ipynb new file mode 100644 index 0000000..1534cc7 --- /dev/null +++ b/bbh_test_strain.ipynb @@ -0,0 +1,193 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "#bilby_bbh_test_onewaveform\n", + "#Version 1.01\n", + "#more strain data from different interferometers\n", + "import numpy as np\n", + "import bilby\n", + "import matplotlib.pyplot as plt\n", + "from gwpy.timeseries import TimeSeries" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "outdir = 'outdir_bbh_test_strain'\n", + "label = 'bbh_test_strain'" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "sampling_frequency = 16384.\n", + "dt = 1.0/sampling_frequency" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "#Load data.\n", + "data = np.loadtxt('./mdc-master/injection-h_m1_L0.9_l2m2_r300.dat')\n", + "N = len(data)\n", + "duration = N*dt\n", + "time = np.arange(0,duration,dt)\n", + "\n", + "\n", + "#The strain data Daniel gave me only include one column. This is to add time column to data.\n", + "f = open('injection-h_m1_L0.9_l2m2_r300.txt','w')\n", + "for i in range(0,N):\n", + " f.write(\"%.20e %.20e\\n\" % (time[i],data[i])) \n", + " \n", + "f.close()\n", + "plt.plot(time,data)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "#Translate \"txt\" file to \"gwf\" file, because BILBY only reads \"gwf\" file.\n", + "\n", + "\n", + "#Read new data which has 2 columns as a timeseries.\n", + "gwffile = TimeSeries.read('injection-h_m1_L0.9_l2m2_r300.txt')\n", + "gwffile.t0 = .0\n", + "\n", + "#Give channel a name.\n", + "gwffile.name='H1:GWOSC-4KHZ_R1_STRAIN'\n", + "\n", + "gwffile.write('injection-h_m1_L0.9_l2m2_r300.gwf')\n", + "\n", + "#This is to check if \"gwf\" file is succeed to produce.\n", + "#data_2=TimeSeries.read('injection-h_m1_L0.9_l2m2_r300.gwf','H1:GWOSC-4KHZ_R1_STRAIN')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# inject the signal into H1 interferometer.\n", + "#ifos = bilby.gw.detector.get_empty_interferometer('H1')\n", + "ifos = bilby.gw.detector.InterferometerList(['H1'])\n", + "\n", + "ifos[0].set_strain_data_from_frame_file(frame_file='injection-h_m1_L0.9_l2m2_r300.gwf',\n", + " sampling_frequency=sampling_frequency, duration=duration, \n", + " start_time=.0,channel='H1:GWOSC-4KHZ_R1_STRAIN')\n", + "#ifos[0].set_strain_data_from_power_spectral_densities(sampling_frequency=sampling_frequency, duration=duration, start_time=-0.5)\n", + "print(ifos[0])\n", + "#Usually for 1 event, we have strain data from 2 or more interferometers, so the following code is needed:\n", + "#for i in range(len(ifos)):\n", + "# ifos[i].set_strain_data_from_frame_file(frame_file='new_injection-h_m16_L0.18_l2m2_r300.gwf',sampling_frequency=sampling_frequency, duration=duration, start_time=.0,channel='')\n", + "# #ifos[i].strain_data_from_gwpy_time_series('injection-h_m16_L0.18_l2m2_r300.gwf')\n", + "# ifos[i].set_strain_data_from_power_spectral_densities(\n", + "# sampling_frequency=sampling_frequency, duration=duration,\n", + "# start_time=-0.5)\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Generate waveform\n", + "waveform_arguments = dict(waveform_approximant='IMRPhenomPv2',\n", + " reference_frequency=50., minimum_frequency=10.)\n", + "\n", + "\n", + "waveform_generator = bilby.gw.WaveformGenerator(\n", + " duration=duration, sampling_frequency=sampling_frequency,\n", + " frequency_domain_source_model=bilby.gw.source.lal_binary_black_hole,\n", + " waveform_arguments=waveform_arguments)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "#Use default BBH priors.\n", + "priors = bilby.core.prior.PriorDict(filename='bilby_gw_prior_files_binary_black_holes.prior')\n", + "\n", + "\n", + "priors['geocent_time'] = bilby.prior.Uniform(\n", + " minimum=.0-0.1,\n", + " maximum=.0+0.1,\n", + " name='geocent_time', latex_label='$t_c$', unit='$s$')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "likelihood = bilby.gw.GravitationalWaveTransient(\n", + " interferometers=ifos, waveform_generator=waveform_generator, priors=priors,\n", + " distance_marginalization=True, phase_marginalization=True, time_marginalization=True)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "result = bilby.run_sampler(\n", + " likelihood=likelihood, priors=priors, sampler='dynesty', npoints=1000,\n", + " injection_parameters=None, outdir=outdir, label=label)\n", + "\n", + "# Make a corner plot.\n", + "result.plot_corner()\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.6.8" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} -- GitLab