Adding plot script for ECC

git-svn-id: https://svn.ligo.caltech.edu/svn/bayeswave/trunk@85 c56465c9-8126-4a4f-9d7d-ac845eff4865

postprocess/plot_ecc.py

+import matplotlib
+matplotlib.use('Agg')
+import numpy as np
+import acor
+import matplotlib.pyplot as plt
+from glob import glob
+from pylal import bayespputils
+
+# -- Undoe wonky plt parameters
+plt.rcParams['figure.figsize'] = 8, 6
+
+filename = 'chains/signal_1_extchain.dat.0'
+data = np.loadtxt(filename, unpack=True)
+ecc = data[5]
+ra = data[2]
+sindec = data[3]
+print ecc
+
+# -- Plot elipticity chain
+plt.figure()
+plt.plot(ecc, '.')
+plt.title('ECC chain')
+plt.savefig('ecc.png')
+plt.close()
+
+# -- Histogram the chain
+plt.figure()
+plt.hist(ecc, bins=30)
+plt.title('Histogram of ecc, ecc=0 is linear')
+plt.savefig('ecc_hist.png')
+plt.close()
+
+# -- Thin the chain
+acl, mean, sigma = acor.acor(ecc)
+print "The ACL of ECC is {0}".format(acl)
+raacl, mean, sigma = acor.acor(ra)
+print "The ACL of RA is {0}".format(raacl)
+
+ecc = ecc[::raacl]
+ra = ra[::raacl]
+sindec = sindec[::raacl]
+
+# -- Plot the thinned chain
+plt.figure()
+plt.plot(ecc,'.')
+plt.title('Thin Ellipticity chain')
+plt.savefig('thin_ecc.png')
+plt.close()
+
+# -- Plot ACF of thin chain
+plt.figure()
+acf = bayespputils.autocorrelation(ecc)
+plt.plot(acf, '.')
+plt.title('Auto-correlation function of thinned ECC')
+plt.savefig('ecc_acf.png')
+plt.close()
+
+# -- scatter ecc and ra
+plt.figure()
+plt.plot(ra, ecc, '.')
+plt.xlabel('RA')
+plt.ylabel('ECC')
+plt.savefig('ra_ecc.png')
+plt.close()
+
+# -- scatter ecc and sindec
+plt.figure()
+plt.plot(sindec, ecc, '.')
+plt.xlabel('sin(DEC)')
+plt.ylabel('ECC')
+plt.savefig('dec_ecc.png')
+plt.close()
+
+# -- scatter ra and sin dec
+plt.figure()
+plt.scatter(ra, sindec, c=np.abs(ecc))
+plt.xlabel('ra')
+plt.ylabel('sin(dec)')
+plt.title('Color shows abs(ecc)')
+cb = plt.colorbar()
+cb.set_label('abs(ecc)')
+plt.savefig('ecc_radec.png')
+plt.close()
+
+
+# -- make web page
+eccweb = open('ecc.html', 'w')
+figlist = glob('*ecc*.png')
+for fig in figlist:
+
+    eccweb.write('<img src="{0}" width=600>'.format(fig))
+
+eccweb.close()
+print "Thank you, come again"
+    

postprocess/s6/plot_evidence.py

@@ -17,10 +17,11 @@ import sky_area.sky_area_clustering as sky
 import skyview
 
 # --------------
-# Plot Evidence
+# Parameters
 # --------------
 # -- Set skymap resolution
 NSIDE = 32
+mdcfile = '/home/jkanner/baysewave/svn/trunk/burstinj/s6/BurstMDC-ELPTC_S6-Log.txt'
 
 # -- Close any open figures
 plt.close('all')
@@ -29,7 +30,7 @@ plt.close('all')
 # Read S6 MDC log file
 # -----------------
 print "Reading MDC log (this may take a few minutes) ..."
-mdc = ft.Mdc()
+mdc = ft.Mdc(filename=mdcfile)
 print "I read the mdc log!"