plot tweaks, bug fixes

5 files changed, 85 insertions, 65 deletions

diff --git a/fr.py b/fr.py
index da561e4..ef0b775 100755
--- a/fr.py
+++ b/fr.py
@@ -207,12 +207,38 @@ def main():
         )
         mcmc_utils.save_chains(samples, outfile)
 
+    raw = np.load(outfile+'.npy')
+    raw[:,4] *= 1E23
+    raw[:,5] *= 1E21
+    ranges = list(llh_paramset.ranges)
+    ranges[4] = [x*1E23 for x in ranges[4]]
+    ranges[5] = [x*1E21 for x in ranges[5]]
+
+    labels = [
+        r'${\rm sin}^2\theta_{12}$',
+        r'${\rm cos}^4\theta_{13}$',
+        r'${\rm sin}^2\theta_{23}$',
+        r'$\delta$',
+        r'$\Delta m_{21}^2\left[10^{-5}\,{\rm eV}^2\right]$',
+        r'$\Delta m_{31}^2\left[10^{-3}\,{\rm eV}^2\right]$',
+        r'$N_{\rm conv}$',
+        r'$N_{\rm prompt}$',
+        r'$N_{\rm muon}$',
+        r'$N_{\rm astro}$',
+        r'$\gamma_{\rm astro}$',
+        r'${\rm log}_{10}\left[\Lambda^{-1}_{('+ \
+        r'{0}'.format(args.dimension)+r')}'+ \
+        misc_utils.get_units(args.dimension)+r'\right]$'
+    ]
+
     plot_utils.chainer_plot(
-        infile       = outfile+'.npy',
+        infile       = raw,
         outfile      = outfile[:5]+outfile[5:].replace('data', 'plots'),
         outformat    = ['pdf'],
         args         = args,
-        llh_paramset = llh_paramset
+        llh_paramset = llh_paramset,
+        labels       = labels,
+        ranges       = ranges
     )
     print "DONE!"
 
diff --git a/submitter/fr_dag.py b/submitter/fr_dag.py
index c1cf11a..b93ad86 100644
--- a/submitter/fr_dag.py
+++ b/submitter/fr_dag.py
@@ -4,8 +4,8 @@ import os
 import numpy as np
 
 scenarios = [
-    ((1, 0, 0), 'OUT'),
-    ((0, 1, 0), 'OET')
+    ((0, 1, 0), 'OET'),
+    ((1, 0, 0), 'OUT')
 ]
 
 dims = [
diff --git a/submitter/sens_dag.py b/submitter/sens_dag.py
index 829c0b3..5d3d6a3 100644
--- a/submitter/sens_dag.py
+++ b/submitter/sens_dag.py
@@ -3,36 +3,27 @@
 import os
 import numpy as np
 
-# sources = [
-#     (1, 2, 0),
-#     (1, 0, 0),
-#     (0, 1, 0),
-# ]
-
-# XLIMS = (0., 0.3)
-XLIMS = (0.89, 1.0)
-
 x_segments = 100
-x_array = np.linspace(0, 1, x_segments)
-sources = []
-for x in x_array:
-    if x >= XLIMS[0] and x <= XLIMS[1]:
-        sources.append([x, 1-x, 0])
+scenarios = [
+    [(0.00, 0.30), 'OET'],
+    [(0.89, 1.00), 'OUT']
+]
 
 dims = [
     6
 ]
 
-textures = [
-    # 'OET'
-    'OUT'
-]
+for i, (xlims, tex) in enumerate(scenarios):
+    x_array = np.linspace(0, 1, x_segments)
+    sources = []
+    for x in x_array:
+        if x >= xlims[0] and x <= xlims[1]:
+            sources.append([x, 1-x, 0])
+    scenarios[i][0] = tuple(sources)
 
 datadir = '/data/user/smandalia/flavour_ratio/data/sensitivity'
 
-# prefix = ''
-# prefix = '_OET'
-prefix = '_OUT'
+prefix = ''
 
 golemfitsourcepath = os.environ['GOLEMSOURCEPATH'] + '/GolemFit'
 condor_script = golemfitsourcepath + '/scripts/flavour_ratio/submitter/sens_submit.sub'
@@ -80,10 +71,10 @@ with open(dagfile, 'w') as f:
     for dim in dims:
         print 'dims', dim
         of_d = datadir + '/DIM{0}/{1}'.format(dim, prefix)
-        for src in sources:
-            print 'source flavour', src
-            for tex in textures:
-                print 'texture', tex
+        for sources, tex in scenarios:
+            print 'texture', tex
+            for src in sources:
+                print 'source flavour', src
                 for r in xrange(GLOBAL_PARAMS['segments']):
                     print 'run', r
                     f.write('JOB\tjob{0}\t{1}\n'.format(job_number, condor_script))
diff --git a/submitter/sens_submit.sub b/submitter/sens_submit.sub
index 6980e8f..7d543d8 100644
--- a/submitter/sens_submit.sub
+++ b/submitter/sens_submit.sub
@@ -1,5 +1,5 @@
 Executable = /data/user/smandalia/GolemTools/sources/GolemFit/scripts/flavour_ratio/sens.py
-Arguments = "--ast $(ast) --data $(data) --dimension $(dimension) --no-bsm $(no_bsm) --datadir $(datadir) --seed $(seed) --source-ratio $(sr0) $(sr1) $(sr2) --threads $(threads) --binning $(binning) --texture $(texture) --segments $(segments) --eval-segment $(eval_segment) --stat-method $(stat_method) --mn-live-points $(mn_live_points) --mn-tolerance $(mn_tolerance) --mn_efficiency $(mn_efficiency) --mn-output $(mn_output) --run-mn $(run_mn)"
+Arguments = "--ast $(ast) --data $(data) --dimension $(dimension) --no-bsm $(no_bsm) --datadir $(datadir) --seed $(seed) --source-ratio $(sr0) $(sr1) $(sr2) --threads $(threads) --binning $(binning) --texture $(texture) --segments $(segments) --eval-segment $(eval_segment) --stat-method $(stat_method) --mn-live-points $(mn_live_points) --mn-tolerance $(mn_tolerance) --mn-efficiency $(mn_efficiency) --mn-output $(mn_output) --run-mn $(run_mn)"
 
 # All logs will go to a single file 
 log    = /scratch/smandalia/flavour_ratio/submitter/logs/job_$(Cluster).log
diff --git a/utils/plot.py b/utils/plot.py
index 01dd2bb..175d0fb 100644
--- a/utils/plot.py
+++ b/utils/plot.py
@@ -78,18 +78,27 @@ elif os.environ.get('GOLEMSOURCEPATH') is not None:
 if 'submitter' in socket.gethostname():
     rc('text', usetex=False)
 
+mpl.rcParams['text.latex.preamble'] = [
+    r'\usepackage{xcolor}',
+    r'\usepackage{amsmath}',
+    r'\usepackage{amssymb}']
+mpl.rcParams['text.latex.unicode'] = True
+
 
 def gen_figtext(args):
     """Generate the figure text."""
     t = r'$'
-    t += r'{\rm Source\:flavour\:ratio}'+r'\:=\:({0})'.format(
+    if args.data is DataType.REAL:
+        t += r'\textbf{IceCube\:Preliminary}' + '$\n$'
+    elif args.data in [DataType.ASIMOV, DataType.REALISATION]:
+        t += r'{\rm\bf IceCube\:Simulation}' + '$\n$'
+        t += r'\rm{Injected\:flavour\:ratio}'+r'\:=\:({0})_\oplus'.format(
+            solve_ratio(args.injected_ratio).replace('_', ':')
+        ) + '$\n$'
+    t += r'{\rm Source\:flavour\:ratio}'+r'\:=\:({0})_S'.format(
         solve_ratio(args.source_ratio).replace('_', ':')
     )
-    if args.data in [DataType.ASIMOV, DataType.REALISATION]:
-        t += '\n' + r'\rm{Injected\:flavour\:ratio}'+r'\:=\:({0})'.format(
-            solve_ratio(args.injected_ratio).replace('_', ':')
-        )
-    t += '$\n' + r'${\rm Dimension}'+r' = {0}$'.format(args.dimension)
+    t += '$\n$' + r'{\rm Dimension}'+r' = {0}$'.format(args.dimension)
     return t
 
 
@@ -382,7 +391,6 @@ def flavour_contour(frs, nbins, coverage, ax=None, smoothing=0.4,
     else:
         return ev_polygon
 
-
 def plot_Tchain(Tchain, axes_labels, ranges):
     """Plot the Tchain using getdist."""
     Tsample = mcsamples.MCSamples(
@@ -390,21 +398,22 @@ def plot_Tchain(Tchain, axes_labels, ranges):
     )
 
     Tsample.updateSettings({'contours': [0.90, 0.99]})
-    Tsample.num_bins_2D=500
-    Tsample.fine_bins_2D=500
-    Tsample.smooth_scale_2D=0.03
+    Tsample.num_bins_2D=10
+    Tsample.fine_bins_2D=50
+    Tsample.smooth_scale_2D=0.05
 
     g = plots.getSubplotPlotter()
     g.settings.num_plot_contours = 2
-    g.settings.axes_fontsize = 10
+    g.settings.axes_fontsize = 9
     g.settings.figure_legend_frame = False
     g.triangle_plot(
-        [Tsample], filled=True,
+        [Tsample], filled=True#, contour_colors=['green', 'lightgreen']
     )
     return g
 
 
-def chainer_plot(infile, outfile, outformat, args, llh_paramset, fig_text=None):
+def chainer_plot(infile, outfile, outformat, args, llh_paramset, fig_text=None,
+                 labels=None, ranges=None):
     """Make the triangle plot."""
     if hasattr(args, 'plot_elements'):
         if not args.plot_angles and not args.plot_elements:
@@ -423,25 +432,26 @@ def chainer_plot(infile, outfile, outformat, args, llh_paramset, fig_text=None):
     if fig_text is None:
         fig_text = gen_figtext(args)
 
-    axes_labels = llh_paramset.labels
-    ranges = llh_paramset.ranges
+    if labels is None: axes_labels = llh_paramset.labels
+    else: axes_labels = labels
+    if ranges is None: ranges = llh_paramset.ranges
 
     if args.plot_angles:
         print "Making triangle plots"
         Tchain = raw
         g = plot_Tchain(Tchain, axes_labels, ranges)
 
-        mpl.pyplot.figtext(0.5, 0.7, fig_text, fontsize=15)
+        mpl.pyplot.figtext(0.6, 0.7, fig_text, fontsize=20)
 
-        for i_ax_1, ax_1 in enumerate(g.subplots):
-            for i_ax_2, ax_2 in enumerate(ax_1):
-                if i_ax_1 == i_ax_2:
-                    itv = interval(Tchain[:,i_ax_1], percentile=90.)
-                    for l in itv:
-                        ax_2.axvline(l, color='gray', ls='--')
-                        ax_2.set_title(r'${0:.2f}_{{{1:.2f}}}^{{+{2:.2f}}}$'.format(
-                            itv[1], itv[0]-itv[1], itv[2]-itv[1]
-                        ), fontsize=10)
+        # for i_ax_1, ax_1 in enumerate(g.subplots):
+        #     for i_ax_2, ax_2 in enumerate(ax_1):
+        #         if i_ax_1 == i_ax_2:
+        #             itv = interval(Tchain[:,i_ax_1], percentile=90.)
+        #             for l in itv:
+        #                 ax_2.axvline(l, color='gray', ls='--')
+        #                 ax_2.set_title(r'${0:.2f}_{{{1:.2f}}}^{{+{2:.2f}}}$'.format(
+        #                     itv[1], itv[0]-itv[1], itv[2]-itv[1]
+        #                 ), fontsize=10)
 
         # if not args.fix_mixing:
         #     sc_index = llh_paramset.from_tag(ParamTag.SCALE, index=True)
@@ -451,13 +461,6 @@ def chainer_plot(infile, outfile, outformat, args, llh_paramset, fig_text=None):
         #         '1E{2}'.format(itv[0], itv[2], itv[1])
         #     )
 
-        if args.data is DataType.REAL:
-            plt.text(0.8, 0.9, 'IceCube Preliminary', color='red', fontsize=15,
-                     ha='center', va='center')
-        elif args.data in [DataType.ASIMOV, DataType.REALISATION]:
-            plt.text(0.8, 0.9, 'IceCube Simulation', color='red', fontsize=15,
-                     ha='center', va='center')
-
         for of in outformat:
             print 'Saving', outfile+'_angles.'+of
             g.export(outfile+'_angles.'+of)
@@ -700,7 +703,7 @@ def plot_table_sens(data, outfile, outformat, args):
                 ))
 
     ax.get_xaxis().set_visible(True)
-    ax.set_xlabel(r'${\rm New\:Physics\:Scale}\:[\:{\rm log}_{10} (\Lambda^{-1}\:/\:{\rm GeV}^{-d+4})\: ]$',
+    ax.set_xlabel(r'${\rm New\:Physics\:Scale}\:[\:{\rm log}_{10} (\Lambda^{-1}_{(d)}\:/\:{\rm GeV}^{-d+4})\: ]$',
                  fontsize=19)
     ax.tick_params(axis='x', labelsize=16)
 
@@ -894,8 +897,8 @@ def plot_x(data, outfile, outformat, args, normalise=False):
     if normalise:
         fig.text(
             0.02, 0.5,
-            r'${\rm New\:Physics\:Scale}\:[\:{\rm log}_{10} \left (\Lambda_{' +
-            r'\:{0}'.format(args.dimension)+r'}\:/\:{\rm M}_{\:\rm Planck}^{\:'+
+            r'${\rm New\:Physics\:Scale}\:[\:{\rm log}_{10} \left (\Lambda^{-1}_{(' +
+            r'\:{0}'.format(args.dimension)+r')}\:\cdot\:{\rm M}_{\:\rm Planck}^{\:'+
             r'{0}'.format(args.dimension-4)+ r'}\right )\: ]$', ha='left',
             va='center', rotation='vertical', fontsize=largesize
         )