Commit 6cbeb585 authored by ABHIRUP GHOSH's avatar ABHIRUP GHOSH Committed by John Douglas Veitch
Browse files

updates on the IMR consistency codes

lalinference/python/imrtgr_imr_consistency_test.py:
1. Pass the range of epsilon, sigma
2. Contour labels

lalinference/python/lalinference/imrtgr/result.html
Webpage cleaning up: remove the middle panel that is unused.
parent cea710f4
......@@ -21,6 +21,8 @@ import pickle, gzip
import sys
from lalinference import git_version
from scipy.stats import gaussian_kde #rahul: for KDE implementation
from matplotlib import rc
import matplotlib
matplotlib.rc('text.latex', preamble = '\usepackage{txfonts}')
......@@ -102,8 +104,12 @@ if __name__ == '__main__':
parser.add_option("-w", "--waveform", dest="waveform", help="waveform used for recovery")
parser.add_option("-d", "--debug-plots", dest="debug_plots", help="debug plots")
parser.add_option("--N_bins", type="int", dest="N_bins", default=201, help="number of bins (default=201)")
parser.add_option("--dMfbyMf_lim", type="float", dest="dMfbyMf_lim", default=1., help="absolute value of limit for range of dMfbyMf_vec, defined as [-dMfbyMf_lim, +dMfbyMf_lim]")
parser.add_option("--dchifbychif_lim", type="float", dest="dchifbychif_lim", default=1., help="absolute value of limit for range of dchifbychif_vec, defined as [-dchifbychif_lim, +dchifbychif_lim]")
parser.add_option("--use_KDE", type="int", dest="MfafKDE", help="use KDE or not after getting samples of Mf, af")
(options, args) = parser.parse_args()
MfafKDE = options.MfafKDE
insp_post = options.insp_post
ring_post = options.ring_post
imr_post = options.imr_post
......@@ -126,7 +132,8 @@ if __name__ == '__main__':
print('Recovery approximant not provided. To have it displayed on the results page, please pass command line option --waveform.')
N_bins = int(options.N_bins) # Number of grid points along either axis (dMfbyMf, dchifbychif) for computation of the posteriors
dMfbyMf_lim = float(options.dMfbyMf_lim)
dchifbychif_lim = float(options.dchifbychif_lim)
lalinference_datadir = os.getenv('LALINFERENCE_DATADIR')
# create output directory and copy the script
......@@ -257,19 +264,47 @@ if __name__ == '__main__':
Mf_intp = (Mf_bins[:-1] + Mf_bins[1:])/2.
chif_intp = (chif_bins[:-1] + chif_bins[1:])/2.
# compute the 2D posterior distributions for the inspiral, ringodwn and IMR analyses
P_Mfchif_i, Mf_bins, chif_bins = np.histogram2d(Mf_i, chif_i, bins=(Mf_bins, chif_bins), normed=True)
P_Mfchif_r, Mf_bins, chif_bins = np.histogram2d(Mf_r, chif_r, bins=(Mf_bins, chif_bins), normed=True)
P_Mfchif_imr, Mf_bins, chif_bins = np.histogram2d(Mf_imr, chif_imr, bins=(Mf_bins, chif_bins), normed=True)
# transpose to go from (X,Y) indexing returned by np.histogram2d() to array (i,j) indexing for further
# computations. From now onwards, different rows (i) correspond to different values of Mf and different
# columns (j) correspond to different values of chif
P_Mfchif_i = P_Mfchif_i.T
P_Mfchif_r = P_Mfchif_r.T
P_Mfchif_imr = P_Mfchif_imr.T
print 'useKDE=',MfafKDE
if MfafKDE==1:
print 'replacing lal P(Mfaf) with its KDE pdf'
M_i,C_i=np.meshgrid(Mf_intp,chif_intp)
joint_data=np.vstack([Mf_i,chif_i]);kernel=gaussian_kde(joint_data)
f_i = lambda x,y:kernel.evaluate([x,y])
print "for inspiral kernel",kernel.integrate_box([-Mf_lim,-chif_lim],[Mf_lim,chif_lim])
P_Mfchif_i = np.vectorize(f_i)(M_i,C_i)/kernel.integrate_box([-Mf_lim,-chif_lim],[Mf_lim,chif_lim])
joint_data=np.vstack([Mf_r,chif_r]);kernel=gaussian_kde(joint_data)#;M_i,C_i=np.meshgrid(Mf_bins,chif_bins)
f_r = lambda x,y:kernel.evaluate([x,y])
print "for post-inspiral kernel",kernel.integrate_box([-Mf_lim,-chif_lim],[Mf_lim,chif_lim])
P_Mfchif_r = np.vectorize(f_r)(M_i,C_i)/kernel.integrate_box([-Mf_lim,-chif_lim],[Mf_lim,chif_lim])
joint_data=np.vstack([Mf_imr,chif_imr]);kernel=gaussian_kde(joint_data)#;M_i,C_i=np.meshgrid(Mf_bins,chif_bins)
f_imr = lambda x,y:kernel.evaluate([x,y])
print "for imr kernel",kernel.integrate_box([-Mf_lim,-chif_lim],[Mf_lim,chif_lim])
P_Mfchif_imr = np.vectorize(f_imr)(M_i,C_i)/kernel.integrate_box([-Mf_lim,-chif_lim],[Mf_lim,chif_lim])
#rahul: end KDE of Mf,af
elif MfafKDE==0:
print 'using default samples, NOKDE'
# compute the 2D posterior distributions for the inspiral, ringodwn and IMR analyses
P_Mfchif_i, Mf_bins, chif_bins = np.histogram2d(Mf_i, chif_i, bins=(Mf_bins, chif_bins), normed=True)
P_Mfchif_r, Mf_bins, chif_bins = np.histogram2d(Mf_r, chif_r, bins=(Mf_bins, chif_bins), normed=True)
P_Mfchif_imr, Mf_bins, chif_bins = np.histogram2d(Mf_imr, chif_imr, bins=(Mf_bins, chif_bins), normed=True)
# transpose to go from (X,Y) indexing returned by np.histogram2d() to array (i,j) indexing for further
# computations. From now onwards, different rows (i) correspond to different values of Mf and different
# columns (j) correspond to different values of chif
#rahul:Transpose is forbidden if using KDE
P_Mfchif_i = P_Mfchif_i.T
P_Mfchif_r = P_Mfchif_r.T
P_Mfchif_imr = P_Mfchif_imr.T
print 'computed P_Mfchif without using KDE'
###############################################################################################
......@@ -313,8 +348,8 @@ if __name__ == '__main__':
P_Mfchif_r_interp_object = scipy.interpolate.interp2d(Mf_intp, chif_intp, P_Mfchif_r, fill_value=0., bounds_error=False)
# defining limits of delta_Mf/Mf and delta_chif/chif.
dMfbyMf_vec = np.linspace(-1.0, 1.0, N_bins)
dchifbychif_vec = np.linspace(-1.0, 1.0, N_bins)
dMfbyMf_vec = np.linspace(-dMfbyMf_lim, dMfbyMf_lim, N_bins)
dchifbychif_vec = np.linspace(-dchifbychif_lim, dchifbychif_lim, N_bins)
# compute the P(dMf/Mf, dchif/chif) by evaluating the integral
diff_dMfbyMf = np.mean(np.diff(dMfbyMf_vec))
......@@ -642,23 +677,27 @@ if __name__ == '__main__':
plt.ylabel('$\chi_f$')
plt.grid()
CSi.levels = np.asarray(CSi.levels)
CSr.levels = np.asarray(CSr.levels)
CSimr.levels = np.asarray(CSimr.levels)
strs_i = [ 'inspiral', 'inspiral' ]
strs_r = [ 'ringdown', 'ringdown' ]
strs_imr = [ 'IMR', 'IMR' ]
fmt_i = {}
fmt_r = {}
fmt_imr = {}
#for l,s in zip(CSi.levels, strs_i):
#fmt_i[l] = s
#for l,s in zip(CSr.levels, strs_r):
#fmt_r[l] = s
#for l,s in zip(CSimr.levels, strs_imr):
#fmt_imr[l] = s
for l,s in zip(CSi.levels, strs_i):
fmt_i[l] = s
for l,s in zip(CSr.levels, strs_r):
fmt_r[l] = s
for l,s in zip(CSimr.levels, strs_imr):
fmt_imr[l] = s
## Label every other level using strings
#plt.clabel(CSi,CSi.levels[::2],inline=True,fmt=fmt_i,fontsize=14, use_clabeltext=True)
#plt.clabel(CSr,CSr.levels[::2],inline=True,fmt=fmt_r,fontsize=14, use_clabeltext=True)
#plt.clabel(CSimr,CSimr.levels[::2],inline=True,fmt=fmt_imr,fontsize=10)
plt.clabel(CSi,CSi.levels[::2],inline=True,fmt=fmt_i,fontsize=14, use_clabeltext=True)
plt.clabel(CSr,CSr.levels[::2],inline=True,fmt=fmt_r,fontsize=14, use_clabeltext=True)
plt.clabel(CSimr,CSimr.levels[::2],inline=True,fmt=fmt_imr,fontsize=10)
plt.savefig('%s/img/IMR_overlap.png'%(out_dir), dpi=300)
plt.savefig('%s/img/IMR_overlap_thumb.png'%(out_dir), dpi=72)
......@@ -706,8 +745,8 @@ if __name__ == '__main__':
plt.plot(0, 0, 'k+', ms=12, mew=2)
plt.xlabel('$\Delta M_f/M_f$')
plt.ylabel('$\Delta \chi_f/\chi_f$')
plt.xlim([-1.,1.])
plt.ylim([-1.,1.])
plt.xlim([-dMfbyMf_lim,dMfbyMf_lim])
plt.ylim([-dchifbychif_lim,dchifbychif_lim])
plt.grid()
plt.subplot2grid((3,3), (1,2), rowspan=2)
......
......@@ -22,16 +22,12 @@
<div class="report">
<h1 id = "summary" class="heading">Summary</h1>
<table style="width:70%">
<table style="width:50%">
<tr>
<td>
<img src="img/IMR_overlap_thumb.png" alt="(Mf,chif plot)" class="graph-image">
<img src="img/IMR_overlap.png" alt="(Mf,chif plot)" class="graph-image-big">
</td>
<td>
<img src="img/dMfdchif_thumb.png" alt="(delta_Mf, delta_chif plot)" class="graph-image">
<img src="img/dMfdchif.png" alt="(Mf,chif plot)" class="graph-image-big">
</td>
<td>
<img src="img/dMfbyMfdchifbychif_thumb.png" alt="(delta_Mf/Mf, delta_chif/chif plot)" class="graph-image">
<img src="img/dMfbyMfdchifbychif.png" alt="(Mf,chif plot)" class="graph-image-big">
......@@ -39,7 +35,6 @@
</tr>
<tr>
<td class="posteriors-graph-caption graph-caption">The plot shows the 68% and 95% confidence regions on the posteriors on the mass (M<sub>f</sub>) and dimensionless spin (&chi;<sub>f</sub>) of the final black hole estimated from the inspiral and ringdown parts of the signal, as well as the full IMR signal.</td>
<td class="posteriors-graph-caption graph-caption">Posteriors on the parameters (&Delta;M<sub>f</sub>, &Delta;&chi;<sub>f</sub>) describing the difference in the estimates of the mass and the spin of the final black hole estimated from the inspiral and ringdown parts of the signal. The GR prediction is marked by the + sign. </td>
<td class="posteriors-graph-caption graph-caption">Posteriors on the fractional parameters (&Delta;M<sub>f</sub>/M<sub>f</sub>, &Delta;&chi;<sub>f</sub>/&chi;<sub>f</sub>). The GR prediction is marked by the + sign.</td>
</tr>
</table>
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment