Skip to content
Snippets Groups Projects
Commit bda75d0e authored by Yifan Wang's avatar Yifan Wang
Browse files

upload notebooks

parent d3dbcce1
No related branches found
No related tags found
No related merge requests found
%% Cell type:code id:8cf15b3c tags:
``` python
%matplotlib notebook
```
%% Cell type:code id:d5332cc4 tags:
``` python
import itertools
import numpy
import matplotlib
from matplotlib import patches
from matplotlib import pyplot
from pycbc import waveform
from pycbc.inference import io
from pycbc import conversions
from scripts import scatter_histograms
```
%% Cell type:code id:0dd8e4c7 tags:
``` python
pyplot.style.use('seaborn-colorblind')
```
%% Cell type:code id:29a5dff2 tags:
``` python
color_cycle = [c['color'] for c in matplotlib.rcParams['axes.prop_cycle']]
colors = itertools.cycle(color_cycle)
# cycle around
for _ in range(5):
next(colors)
```
%% Cell type:markdown id:f7551a5f tags:
# PyCBC Inference Read Samples
%% Cell type:code id:29d23980 tags:
``` python
filenames = {'GW190521': '/work/yifan.wang/gw190521_xrerun/t3-localrun/0_700/result.hdf',
'GW191109': '/work/yifan.wang/testingGR-o3b/t2-o3bevents/o3btgr_output/samples_files/H1L1V1-INFERENCE_GW191109_010717-1126259200-400.hdf',
'GW200220':'/work/yifan.wang/testingGR-o3b/t4-mpi-massive/t2-200220-mpv1000/result.hdf'}
samples = {}
for lbl, fn in filenames.items():
fp = io.loadfile(fn, 'r')
s = fp.read_samples(['srcmchirp', 'q', 'parity_mpvinverse'])
samples[lbl] = s
fp.close()
```
%% Cell type:code id:a94f9251 tags:
``` python
filelabels = {'GW190521': 'GW190521',
'GW191109': 'GW191109',
'GW200220': 'GW200220'
}
```
%% Cell type:code id:fc371486 tags:
``` python
samples
```
%% Output
{'GW190521': rec.array([(-27873.26400656, 83.26078987, 436.04949975, 1.43400538),
(-27875.2432163 , 71.57790771, 541.07468378, 1.64651789),
(-27875.26155477, 74.09536189, 440.63653902, 1.73495024), ...,
(-27880.60033757, 63.38915481, 277.02222192, 1.29027691),
(-27869.3461176 , 95.44301611, 381.03490203, 1.50594504),
(-27877.14788405, 79.98899193, 313.76083404, 1.28679349)],
dtype=[('loglikelihood', '<f8'), ('srcmchirp', '<f8'), ('parity_mpvinverse', '<f8'), ('q', '<f8')]),
'GW191109': rec.array([(-17938.0557928 , 48.73074706, 144.75730556, 1.16390601),
(-17942.74637522, 52.43607729, 329.71287488, 1.02073016),
(-17942.07580232, 41.9276795 , 180.19740616, 2.19670947), ...,
(-17940.19991335, 48.96967047, 405.20982651, 1.45601162),
(-17935.2245303 , 50.17113232, 235.30007999, 1.16774801),
(-17938.96972275, 53.97280568, 244.20575908, 1.28908152)],
dtype=[('loglikelihood', '<f8'), ('srcmchirp', '<f8'), ('parity_mpvinverse', '<f8'), ('q', '<f8')]),
'GW200220': rec.array([(-81869.72084721, 64.4399525 , 511.48166829, 2.134752 ),
(-81863.8137886 , 89.29121051, 418.60442017, 3.00415035),
(-81867.14210206, 67.036858 , 239.90584564, 1.45710404), ...,
(-81863.32826024, 59.46329134, 606.33211039, 2.63692637),
(-81864.82296474, 69.8439018 , 467.65624429, 1.19154618),
(-81866.23224365, 104.03905068, 100.25604593, 2.5061787 )],
dtype=[('loglikelihood', '<f8'), ('srcmchirp', '<f8'), ('parity_mpvinverse', '<f8'), ('q', '<f8')])}
%% Cell type:code id:ea533ca3 tags:
``` python
labels = {'srcmchirp': waveform.parameters.srcmchirp.label,
'q': waveform.parameters.q.label,
'parity_mpvinverse':'$M_\mathrm{PV}^{-1} (\mathrm{GeV}^{-1})$'
}
```
%% Cell type:code id:d1fc164a tags:
``` python
x = 'srcmchirp'
y = 'q'
z = 'parity_mpvinverse'
mins = {x: min(s[x].min() for s in samples.values()),
y: min(s[y].min() for s in samples.values()),
z: min(s[z].min() for s in samples.values())
}
maxs = {x: max(s[x].max() for s in samples.values()),
y: min(max(s[y].max() for s in samples.values()),6),
z: max(s[z].max() for s in samples.values())
}
```
%% Cell type:code id:e95d473c tags:
``` python
mins
```
%% Output
{'srcmchirp': 31.170769850622765,
'q': 1.0000229039873965,
'parity_mpvinverse': 0.006938599419005487}
%% Cell type:code id:6ed4aeed tags:
``` python
maxs
```
%% Output
{'srcmchirp': 149.1743322543928,
'q': 6,
'parity_mpvinverse': 999.9916720638518}
%% Cell type:code id:e84b1343 tags:
``` python
# create a pool for evaluating the kombine kde
from pycbc.pool import choose_pool
pool = choose_pool(10)
```
%% Cell type:code id:37a7a798 tags:
``` python
samples
```
%% Output
{'GW190521': rec.array([(-27873.26400656, 83.26078987, 436.04949975, 1.43400538),
(-27875.2432163 , 71.57790771, 541.07468378, 1.64651789),
(-27875.26155477, 74.09536189, 440.63653902, 1.73495024), ...,
(-27880.60033757, 63.38915481, 277.02222192, 1.29027691),
(-27869.3461176 , 95.44301611, 381.03490203, 1.50594504),
(-27877.14788405, 79.98899193, 313.76083404, 1.28679349)],
dtype=[('loglikelihood', '<f8'), ('srcmchirp', '<f8'), ('parity_mpvinverse', '<f8'), ('q', '<f8')]),
'GW191109': rec.array([(-17938.0557928 , 48.73074706, 144.75730556, 1.16390601),
(-17942.74637522, 52.43607729, 329.71287488, 1.02073016),
(-17942.07580232, 41.9276795 , 180.19740616, 2.19670947), ...,
(-17940.19991335, 48.96967047, 405.20982651, 1.45601162),
(-17935.2245303 , 50.17113232, 235.30007999, 1.16774801),
(-17938.96972275, 53.97280568, 244.20575908, 1.28908152)],
dtype=[('loglikelihood', '<f8'), ('srcmchirp', '<f8'), ('parity_mpvinverse', '<f8'), ('q', '<f8')]),
'GW200220': rec.array([(-81869.72084721, 64.4399525 , 511.48166829, 2.134752 ),
(-81863.8137886 , 89.29121051, 418.60442017, 3.00415035),
(-81867.14210206, 67.036858 , 239.90584564, 1.45710404), ...,
(-81863.32826024, 59.46329134, 606.33211039, 2.63692637),
(-81864.82296474, 69.8439018 , 467.65624429, 1.19154618),
(-81866.23224365, 104.03905068, 100.25604593, 2.5061787 )],
dtype=[('loglikelihood', '<f8'), ('srcmchirp', '<f8'), ('parity_mpvinverse', '<f8'), ('q', '<f8')])}
%% Cell type:code id:737146e4 tags:
``` python
filecolors = {}
fig, axd = (None, None)
for run in samples:
c = next(colors)
filecolors[run] = c
s = samples[run]
fig, axd = scatter_histograms.create_multidim_plot(
['parity_mpvinverse', 'srcmchirp'], s, labels=labels,
mins=mins, maxs=maxs,
plot_scatter=False, plot_contours=True,
marginal_percentiles=[], contour_percentiles=[50, 90],
fill_color=None, contour_color=c, hist_color=c, line_color=c,
contour_ls=['solid', 'dashed'], contour_labels=False, fold_masses=True,
use_kombine=True, kdeargs={'pool': pool, 'max_samples': 200},
fig=fig, axis_dict=axd,
)
# add legend
handles = []
lbls = []
for run, c in filecolors.items():
l = filelabels[run]
handles.append(patches.Patch(color=c, label=l))
lbls.append(l)
fig.legend(loc=(0.7, 0.75), handles=handles, labels=lbls)
fig.savefig('3events.png')
```
%% Cell type:code id:95c5bcdf tags:
``` python
color_cycle = ['#1f77b4', '#ff7f0e']
colors = itertools.cycle(color_cycle)
# cycle around
for _ in range(5):
next(colors)
```
%% Cell type:code id:e69aacc6 tags:
``` python
color_cycle = [c['color'] for c in matplotlib.rcParams['axes.prop_cycle']]
colors = itertools.cycle(color_cycle)
# cycle around
#for _ in range(5):
# next(colors)
```
%% Cell type:code id:5902c7f1 tags:
``` python
for i in range(5):
print(i)
```
%% Output
0
1
2
3
4
%% Cell type:code id:4ff63d8a tags:
``` python
color_cycle
```
%% Output
['#0072B2', '#009E73', '#D55E00', '#CC79A7', '#F0E442', '#56B4E9']
%% Cell type:code id:2b6d8f06 tags:
``` python
[c['color'] for c in matplotlib.rcParams['axes.prop_cycle']]
plt.figure()
plt.plot([1,2],color=c)
```
%% Output
['#0072B2', '#009E73', '#D55E00', '#CC79A7', '#F0E442', '#56B4E9']
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
/local/user/yifan.wang/ipykernel_3285978/1511099479.py in <module>
----> 1 [c['color'] for c in matplotlib.rcParams['axes.prop_cycle']]
2 plt.figure()
3 plt.plot([1,2],color=c)
NameError: name 'matplotlib' is not defined
%% Cell type:code id:87080ba5 tags:
``` python
filecolors = {}
fig, axd = (None, None)
for run in samples:
c = next(colors)
filecolors[run] = c
s = samples[run]
fig, axd = scatter_histograms.create_multidim_plot(
['parity_mpvinverse', 'srcmchirp','q'], s, labels=labels,
mins=mins, maxs=maxs,
plot_scatter=False, plot_contours=True,
marginal_percentiles=[], contour_percentiles=[50, 90],
fill_color=None, contour_color=c, hist_color=c, line_color=c,
contour_ls=['solid', 'dashed'], contour_labels=False, fold_masses=True,
use_kombine=True, kdeargs={'pool': pool, 'max_samples': 200},
fig=fig, axis_dict=axd,
)
# add legend
handles = []
lbls = []
for run, c in filecolors.items():
l = filelabels[run]
handles.append(patches.Patch(color=c, label=l))
lbls.append(l)
fig.legend(loc=(0.7, 0.75), handles=handles, labels=lbls)
fig.show()
fig.savefig('./3events-2.png')
```
%% Cell type:code id:132771ea tags:
``` python
```
......
Source diff could not be displayed: it is too large. Options to address this: view the blob.
This diff is collapsed.
......@@ -2,7 +2,7 @@
"cells": [
{
"cell_type": "code",
"execution_count": 36,
"execution_count": 1,
"id": "expressed-speed",
"metadata": {},
"outputs": [],
......@@ -40,7 +40,7 @@
},
{
"cell_type": "code",
"execution_count": 37,
"execution_count": 2,
"id": "settled-motel",
"metadata": {},
"outputs": [],
......
%% Cell type:code id:expressed-speed tags:
``` python
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from pycbc.inference import io, models
from pycbc import waveform
from pycbc.detector import Detector
from pycbc import cosmology
# PLOTTING OPTIONS
fig_width_pt = 3*246.0 # Get this from LaTeX using \showthe\columnwidth
inches_per_pt = 1.0/72.27 # Convert pt to inch
golden_mean = (np.sqrt(5)-1.0)/2.0 # Aesthetic ratio
fig_width = fig_width_pt*inches_per_pt # width in inches
fig_height = fig_width*golden_mean # height in inches
fig_size = [fig_width,fig_height]
params = { 'axes.labelsize': 24,
'font.family': 'serif',
'font.serif': 'Computer Modern Raman',
'font.size': 24,
'legend.fontsize': 20,
'xtick.labelsize': 24,
'ytick.labelsize': 24,
'axes.grid' : True,
'text.usetex': True,
'savefig.dpi' : 100,
'lines.markersize' : 14,
'figure.figsize': fig_size}
mpl.rcParams.update(params)
```
%% Cell type:code id:settled-motel tags:
``` python
def fn_to_waveform(fn):
# read in the raw posterior file and extract the max loglikelihood parameter value
with io.loadfile(fn, 'r') as fp:
data = fp.read_data()
psds = fp.read_psds()
cp = fp.read_config_file()
samples = fp.read_samples(list(fp['samples'].keys()))
idx = samples['loglikelihood'].argmax()
params = {p: samples[p][idx] for p in samples.fieldnames}
redshift = cosmology.redshift(params['distance'])
params['comoving_volume'] = cosmology.cosmological_quantity_from_redshift(redshift, 'comoving_volume')
# set up the model
model = models.read_from_config(cp, data=data, psds=psds)
model.update(**params)
_ = model.loglikelihood
pol = model.current_stats['maxl_polarization']
#compute the max loglikelihood waveform
wfs = model.waveform_generator.generate(polarization=pol,**model.current_params)
whdata = {}
for det in model.data:
d = model.data[det]
d = d.copy()
d *= model.weight[det] / (4*model.psds[det].delta_f)**0.5
whdata[det] = d.to_timeseries()
waveforms = {}
whitened_waveforms = {}
for detname in wfs:
hs = {}
whitend_hs = {}
det = Detector(detname)
fp, fc = det.antenna_pattern(model.current_params['ra'],
model.current_params['dec'],
pol,
model.current_params['tc'])
hp, hc = wfs[detname]
h = fp*hp + fc*hc
h.resize(len(model.weight[detname]))
wh = h.copy()
wh *= model.weight[detname] / (4*model.psds[detname].delta_f)**0.5
whitened_waveforms[detname] = wh.to_timeseries()
return whdata, whitened_waveforms,model
```
%% Cell type:code id:8691a3c5 tags:
``` python
parity_path = '/work/yifan.wang/testingGR-o3b/snr/GW190521/convert.hdf'
parity_whdata, parity_whwaveform, parity_model = fn_to_waveform(parity_path)
gr_path = '/work/yifan.wang/gw190521_xrerun/t1-mpi-gr/mpirun_output/config_files/convert.hdf'
_, gr_whwaveform,gr_model = fn_to_waveform(gr_path)
```
%% Output
WARNING:root:No sampling_params section read from config file
WARNING:root:WARNING: The following args are not being used for waveform generation: redshift, srcmchirp, comoving_volume, delta_tc, parity_mpvinverse, q, baseapprox, trigger_time
WARNING: fval is not bracketed by func(zmin)=3.630907107506232e-13 Mpc3 and func(zmax)=164042102234.23068 Mpc3. This means either there is no solution, or that there is more than one solution between zmin and zmax satisfying fval = func(z). [astropy.cosmology.funcs]
WARNING:astropy:fval is not bracketed by func(zmin)=3.630907107506232e-13 Mpc3 and func(zmax)=164042102234.23068 Mpc3. This means either there is no solution, or that there is more than one solution between zmin and zmax satisfying fval = func(z).
WARNING: fval is not bracketed by func(zmin)=3.630907107506232e-13 Mpc3 and func(zmax)=3763143105730.4414 Mpc3. This means either there is no solution, or that there is more than one solution between zmin and zmax satisfying fval = func(z). [astropy.cosmology.funcs]
WARNING:astropy:fval is not bracketed by func(zmin)=3.630907107506232e-13 Mpc3 and func(zmax)=3763143105730.4414 Mpc3. This means either there is no solution, or that there is more than one solution between zmin and zmax satisfying fval = func(z).
WARNING:root:No sampling_params section read from config file
WARNING:root:WARNING: The following args are not being used for waveform generation: redshift, srcmchirp, comoving_volume, delta_tc, q, trigger_time
%% Cell type:code id:d27e743c tags:
``` python
xlim=(-0.1, 0.1)
fig = plt.figure(figsize=(10.211706102117061, 6.311181454233049*1.6))
fig.subplots_adjust(hspace=.2)
ax = fig.add_subplot(311);
tc = parity_model.static_params['trigger_time']
det = list(gr_whwaveform.keys())
wdata = parity_whdata[det[0]]
ax.plot(wdata.sample_times-tc, wdata, c='k', alpha=0.2, label='whitened data')
wwf = parity_whwaveform[det[0]]
ax.plot(wwf.sample_times-tc, wwf,label='birefringence')
ax.plot(gr_whwaveform[det[0]].sample_times-tc, gr_whwaveform[det[0]],label='GR')
ax.grid(True)
ax.set_xticklabels([])
ax.set_xlim(*xlim)
ax.set_ylim(-100,100)
ax.text(0.054, 70, 'LIGO Hanford', style='italic',fontsize=20,
bbox={'facecolor':'white', 'alpha':0.5, 'pad':5,'edgecolor':'none'})
ax.set_title('GW190521')
################
bx = fig.add_subplot(312);
wdata = parity_whdata[det[1]]
bx.plot(wdata.sample_times-tc, wdata, c='k', alpha=0.2, label='whitened data')
wwf = parity_whwaveform[det[1]]
bx.plot(wwf.sample_times-tc, wwf,label='birefringence')
bx.plot(gr_whwaveform[det[1]].sample_times-tc, gr_whwaveform[det[1]],label='GR')
bx.set_ylabel('Whitened Strain')
bx.set_xlim(*xlim)
bx.set_ylim(-100,100)
bx.text(0.0515, 70, 'LIGO Livingston', style='italic',fontsize=20,
bbox={'facecolor':'white', 'alpha':0.5, 'pad':0,'edgecolor':'none'})
bx.grid(True)
bx.set_xticklabels([])
################
cx = fig.add_subplot(313);
wdata = parity_whdata[det[2]]
cx.plot(wdata.sample_times-tc, wdata, c='k', alpha=0.2, label='whitened data')
wwf = parity_whwaveform[det[2]]
cx.plot(wwf.sample_times-tc, wwf,label='birefringence')
cx.plot(gr_whwaveform[det[2]].sample_times-tc, gr_whwaveform[det[2]],label='GR')
cx.text(0.07, 70, 'Virgo', style='italic',fontsize=20,
bbox={'facecolor':'white', 'alpha':0.5, 'pad':5,'edgecolor':'none'})
cx.set_xlim(*xlim)
cx.set_ylim(-100,100)
cx.set_xlabel('Time - '+str(tc)+' [s]')
cx.legend(loc='lower right', ncol=3,framealpha=0)
fig.savefig('0521_gw190521.pdf',bbox_inches='tight')
```
%% Cell type:code id:a2769cdb tags:
``` python
parity_model.current_params
```
%% Output
{'inclination': 1.3773388232444717,
'maxl_polarization': 5.132211422080623,
'srcmchirp': 59.988695836431965,
'logwt': -28046.517628006437,
'V1_optimal_snrsq': 4.645117051883857,
'parity_mpvinverse': 571.5796601060491,
'distance': 1333.6349776943568,
'spin1_a': 0.9116850950619036,
'dec': 0.5513436858058255,
'delta_tc': 0.0061986454908000005,
'spin1_polar': 2.076365152215487,
'spin2_azimuthal': 6.273223231326342,
'spin2_polar': 0.9148528363767893,
'ra': 3.3490655420980504,
'logprior': -49.260229458075116,
'spin2_a': 0.09566031928959656,
'coa_phase': 2.930192323180816,
'loglr': 122.35747765652349,
'L1_optimal_snrsq': 163.8200365758998,
'loglikelihood': -28001.80399824368,
'logjacobian': 0.0,
'q': 4.542030730356726,
'maxl_loglr': 125.90973288111286,
'spin1_azimuthal': 2.6322431174435468,
'H1_optimal_snrsq': 69.48418153589377,
'comoving_volume': 5018380778.811694,
'approximant': 'mpvnosmall',
'baseapprox': 'IMRPhenomXPHM',
'f_lower': 20.0,
'f_ref': 20.0,
'trigger_time': 1242442967.4473,
'spin1x': -0.6963825154031204,
'spin1y': 0.3889331472174398,
'spin1z': -0.44153381727155666,
'spin2x': 0.07580444315196702,
'spin2y': -0.0007551945954393175,
'spin2z': 0.05834391798748846,
'tc': 1242442967.4534986,
'redshift': 0.2556724782843754,
'mass1': 194.34779119906074,
'mass2': 42.78874422846471}
%% Cell type:markdown id:7476300f tags:
# SNR
%% Cell type:code id:6bd3b807 tags:
``` python
import h5py
from pycbc.conversions import snr_from_loglr
nongr = h5py.File(parity_path,'r')
gr = h5py.File(gr_path,'r')
```
%% Cell type:code id:2b6eb9a3 tags:
``` python
snr_from_loglr(np.max(nongr['samples']['loglr']))
```
%% Output
15.643367774013592
%% Cell type:code id:358ab6ca tags:
``` python
snr_from_loglr(np.max(gr['samples']['loglr']))
```
%% Output
15.312425524512493
%% Cell type:code id:eae483fd tags:
``` python
```
%% Cell type:code id:expressed-speed tags:
``` python
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from pycbc.inference import io, models
from pycbc import waveform
from pycbc.detector import Detector
from pycbc import cosmology
# PLOTTING OPTIONS
fig_width_pt = 3*246.0 # Get this from LaTeX using \showthe\columnwidth
inches_per_pt = 1.0/72.27 # Convert pt to inch
golden_mean = (np.sqrt(5)-1.0)/2.0 # Aesthetic ratio
fig_width = fig_width_pt*inches_per_pt # width in inches
fig_height = fig_width*golden_mean # height in inches
fig_size = [fig_width,fig_height]
params = { 'axes.labelsize': 24,
'font.family': 'serif',
'font.serif': 'Computer Modern Raman',
'font.size': 24,
'legend.fontsize': 20,
'xtick.labelsize': 24,
'ytick.labelsize': 24,
'axes.grid' : True,
'text.usetex': True,
'savefig.dpi' : 100,
'lines.markersize' : 14,
'figure.figsize': fig_size}
mpl.rcParams.update(params)
```
%% Cell type:code id:settled-motel tags:
``` python
def fn_to_waveform(fn):
# read in the raw posterior file and extract the max loglikelihood parameter value
with io.loadfile(fn, 'r') as fp:
data = fp.read_data()
psds = fp.read_psds()
cp = fp.read_config_file()
samples = fp.read_samples(list(fp['samples'].keys()))
idx = samples['loglikelihood'].argmax()
params = {p: samples[p][idx] for p in samples.fieldnames}
redshift = cosmology.redshift(params['distance'])
params['comoving_volume'] = cosmology.cosmological_quantity_from_redshift(redshift, 'comoving_volume')
# set up the model
model = models.read_from_config(cp, data=data, psds=psds)
model.update(**params)
_ = model.loglikelihood
pol = model.current_stats['maxl_polarization']
#compute the max loglikelihood waveform
wfs = model.waveform_generator.generate(polarization=pol,**model.current_params)
whdata = {}
for det in model.data:
d = model.data[det]
d = d.copy()
d *= model.weight[det] / (4*model.psds[det].delta_f)**0.5
whdata[det] = d.to_timeseries()
waveforms = {}
whitened_waveforms = {}
for detname in wfs:
hs = {}
whitend_hs = {}
det = Detector(detname)
fp, fc = det.antenna_pattern(model.current_params['ra'],
model.current_params['dec'],
pol,
model.current_params['tc'])
hp, hc = wfs[detname]
h = fp*hp + fc*hc
h.resize(len(model.weight[detname]))
wh = h.copy()
wh *= model.weight[detname] / (4*model.psds[detname].delta_f)**0.5
whitened_waveforms[detname] = wh.to_timeseries()
return whdata, whitened_waveforms,model
```
%% Cell type:code id:8691a3c5 tags:
``` python
parity_path = '/work/yifan.wang/testingGR-o3b/snr/GW191109/t2run/convert.hdf'
parity_whdata, parity_whwaveform, parity_model = fn_to_waveform(parity_path)
#gr_path = '/work/yifan.wang/testingGR-o3b/t4-mpi-massive/t2-gr-191109-lowlive/convert.hdf'
gr_path = '/work/yifan.wang/4ogc/release_prod/convertsnr/1212_posterior/GW191109_010717-PYCBC-POSTERIOR-IMRPhenomXPHM.hdf'
_, gr_whwaveform,gr_model = fn_to_waveform(gr_path)
```
%% Output
WARNING:root:No sampling_params section read from config file
WARNING:root:WARNING: The following args are not being used for waveform generation: delta_tc, q, trigger_time, redshift, baseapprox, comoving_volume, parity_mpvinverse, srcmchirp
WARNING:root:No sampling_params section read from config file
WARNING:root:WARNING: The following args are not being used for waveform generation: recalib_phase_l1_2, q, recalib_amplitude_l1_7, redshift, recalib_amplitude_l1_4, recalib_amplitude_l1_6, recalib_amplitude_l1_3, recalib_phase_l1_9, recalib_phase_l1_3, recalib_phase_l1_6, recalib_amplitude_l1_0, recalib_phase_h1_8, recalib_amplitude_h1_8, recalib_phase_h1_7, recalib_phase_h1_9, delta_tc, recalib_phase_h1_4, recalib_phase_h1_0, recalib_amplitude_l1_1, recalib_phase_l1_7, recalib_phase_h1_5, recalib_amplitude_h1_3, recalib_amplitude_h1_1, recalib_amplitude_h1_7, recalib_phase_l1_0, recalib_amplitude_l1_2, recalib_amplitude_h1_5, recalib_amplitude_h1_6, recalib_amplitude_l1_8, recalib_phase_l1_4, trigger_time, recalib_phase_h1_2, recalib_phase_h1_3, comoving_volume, srcmchirp, recalib_phase_h1_1, recalib_amplitude_l1_5, recalib_amplitude_h1_9, recalib_phase_l1_5, recalib_amplitude_l1_9, recalib_phase_l1_8, recalib_amplitude_h1_2, recalib_amplitude_h1_4, recalib_amplitude_h1_0, recalib_phase_h1_6, recalib_phase_l1_1
%% Cell type:code id:d27e743c tags:
``` python
xlim=(-0.1, 0.1)
fig = plt.figure(figsize=(10.211706102117061, 6.311181454233049*1.6))
fig.subplots_adjust(hspace=.2)
ax = fig.add_subplot(311);
tc = parity_model.static_params['trigger_time']
det = list(gr_whwaveform.keys())
wdata = parity_whdata[det[0]]
ax.plot(wdata.sample_times-tc, wdata, c='k', alpha=0.2, label='whitened data')
wwf = parity_whwaveform[det[0]]
ax.plot(wwf.sample_times-tc, wwf,label='birefringence')
ax.plot(gr_whwaveform[det[0]].sample_times-tc, gr_whwaveform[det[0]],label='GR')
ax.grid(True)
ax.set_xticklabels([])
ax.set_xlim(*xlim)
ax.set_ylim(-100,100)
ax.text(0.054, 70, 'LIGO Hanford', style='italic',fontsize=20,
bbox={'facecolor':'white', 'alpha':0.5, 'pad':5,'edgecolor':'none'})
ax.set_title('GW191109')
################
bx = fig.add_subplot(312);
wdata = parity_whdata[det[1]]
bx.plot(wdata.sample_times-tc, wdata, c='k', alpha=0.2, label='whitened data')
wwf = parity_whwaveform[det[1]]
bx.plot(wwf.sample_times-tc, wwf,label='birefringence')
bx.plot(gr_whwaveform[det[1]].sample_times-tc, gr_whwaveform[det[1]],label='GR')
bx.set_ylabel('Whitened Strain')
bx.set_xlim(*xlim)
bx.set_ylim(-100,100)
bx.text(0.0515, 70, 'LIGO Livingston', style='italic',fontsize=20,
bbox={'facecolor':'white', 'alpha':0.5, 'pad':0,'edgecolor':'none'})
bx.grid(True)
bx.set_xlim(*xlim)
bx.set_ylim(-100,100)
bx.set_xlabel('Time - '+str(tc)+' [s]')
ax.legend(loc='lower right', ncol=3,framealpha=0,bbox_to_anchor=(1,-0.1))
fig.savefig('0521_gw191109.pdf',bbox_inches='tight')
```
%% Cell type:markdown id:7476300f tags:
# SNR
%% Cell type:code id:6bd3b807 tags:
``` python
import h5py
from pycbc.conversions import snr_from_loglr
nongr = h5py.File(parity_path,'r')
gr = h5py.File(gr_path,'r')
```
%% Cell type:code id:2b6eb9a3 tags:
``` python
snr_from_loglr(np.max(nongr['samples']['maxl_loglr']))
```
%% Output
16.374853000017954
%% Cell type:code id:358ab6ca tags:
``` python
snr_from_loglr(np.max(gr['samples']['maxl_loglr']))
```
%% Output
16.461053196773342
%% Cell type:code id:cdd6d44c tags:
``` python
```
%% Cell type:code id:eae483fd tags:
``` python
nongr = h5py.File('/work/yifan.wang/testingGR-o3b/snr/GW191109/convert.hdf','r')
```
%% Cell type:code id:6f366b66 tags:
``` python
gr = h5py.File('/work/yifan.wang/testingGR-o3b/t4-mpi-massive/t2-gr-191109-lowlive/convert.hdf','r')
```
%% Cell type:code id:0e541d35 tags:
``` python
snr_from_loglr(np.max(nongr['samples']['maxl_loglr']))
```
%% Output
16.340538761683053
%% Cell type:code id:e7326128 tags:
``` python
snr_from_loglr(np.max(gr['samples']['loglr']))
```
%% Output
16.057563122880396
%% Cell type:code id:d33a01f7 tags:
``` python
gr=h5py.File('/work/yifan.wang/testingGR-o3b/t4-mpi-massive/gr-191109/convert.hdf','r')
snr_from_loglr(np.max(gr['samples']['maxl_loglr']))
```
%% Output
16.521160910966348
%% Cell type:code id:b92e9523 tags:
``` python
file=h5py.File('/work/yifan.wang/testingGR-o3b/t4-mpi-massive/highlive2-GW191109-xphm-mpv1000/convert.hdf','r')
snr_from_loglr(np.max(file['samples']['maxl_loglr']))
```
%% Output
16.37840081629306
%% Cell type:code id:92155527 tags:
``` python
file=h5py.File('/work/yifan.wang/testingGR-o3b/t4-mpi-massive/highlive-GW191109-xphm-mpv1000/convert.hdf','r')
snr_from_loglr(np.max(file['samples']['maxl_loglr']))
```
%% Output
16.375657933055365
%% Cell type:code id:16ec631d tags:
``` python
fn = '/work/yifan.wang/testingGR-o3b/t4-mpi-massive/t2-gr-191109-lowlive/result.hdf'
def printsnr(fn):
with io.loadfile(fn, 'r') as fp:
data = fp.read_data()
psds = fp.read_psds()
cp = fp.read_config_file()
samples = fp.read_samples(list(fp['samples'].keys()))
idx = samples['loglikelihood'].argmax()
params = {p: samples[p][idx] for p in samples.fieldnames}
# set up the model
model = models.read_from_config(cp, data=data, psds=psds)
model.update(**params)
_ = model.loglikelihood
print(snr_from_loglr(model.loglr))
```
%% Cell type:code id:5c0a894f tags:
``` python
printsnr(fn)
```
%% Output
WARNING:root:No sampling_params section read from config file
WARNING:root:WARNING: The following args are not being used for waveform generation: delta_tc, q, trigger_time, redshift, comoving_volume, srcmchirp
16.05759603653594
%% Cell type:code id:bd8c92d1 tags:
``` python
printsnr('/work/yifan.wang/testingGR-o3b/t4-mpi-massive/GW191109-xphm-mpv1000/result.hdf')
```
%% Output
WARNING:root:No sampling_params section read from config file
WARNING:root:WARNING: The following args are not being used for waveform generation: delta_tc, q, trigger_time, redshift, baseapprox, comoving_volume, parity_mpvinverse, srcmchirp
16.094266974440536
%% Cell type:code id:ce51561b tags:
``` python
printsnr('/work/yifan.wang/testingGR-o3b/t4-mpi-massive/highlive2-GW191109-xphm-mpv1000/result.hdf')
```
%% Output
WARNING:root:No sampling_params section read from config file
WARNING:root:WARNING: The following args are not being used for waveform generation: delta_tc, q, trigger_time, redshift, baseapprox, comoving_volume, parity_mpvinverse, srcmchirp
16.136603108612462
%% Cell type:code id:aa022dfd tags:
``` python
printsnr('/work/yifan.wang/testingGR-o3b/t4-mpi-massive/highlive-GW191109-xphm-mpv1000/result.hdf')
```
%% Output
WARNING:root:No sampling_params section read from config file
WARNING:root:WARNING: The following args are not being used for waveform generation: delta_tc, q, trigger_time, redshift, baseapprox, comoving_volume, parity_mpvinverse, srcmchirp
16.13759201817148
%% Cell type:code id:ec31e12b tags:
``` python
printsnr('/work/yifan.wang/testingGR-o3b/t4-mpi-massive/highlive2-gr-GW191109/')
```
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment