aboutsummaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
-rwxr-xr-xplot_sens.py18
-rw-r--r--submitter/mcmc_dag.py4
-rw-r--r--submitter/mcmc_submit.sub1
-rw-r--r--submitter/sens_dag.py18
-rw-r--r--submitter/sens_submit.sub4
-rw-r--r--utils/fr.py2
-rw-r--r--utils/gf.py17
-rw-r--r--utils/plot.py14
8 files changed, 52 insertions, 26 deletions
diff --git a/plot_sens.py b/plot_sens.py
index 279d66a..8357ffc 100755
--- a/plot_sens.py
+++ b/plot_sens.py
@@ -145,6 +145,10 @@ def parse_args(args=None):
help='Make sensitivity plots'
)
parser.add_argument(
+ '--data', type=misc_utils.parse_bool, default='False',
+ help='Use HESE7 data'
+ )
+ parser.add_argument(
'--plot-statistic', type=misc_utils.parse_bool, default='False',
help='Plot MultiNest evidence or LLH value'
)
@@ -227,6 +231,7 @@ def main():
infile += '/golemfit/'
elif args.likelihood is Likelihood.GAUSSIAN:
infile += '/gaussian/'
+ # infile += '/DIM{0}/fix_ifr/HESESim'.format(dim)
infile += '/DIM{0}/fix_ifr/'.format(dim)
if args.likelihood is Likelihood.GAUSSIAN:
infile += '{0}/'.format(str(args.sigma_ratio).replace('.', '_'))
@@ -273,6 +278,7 @@ def main():
argsc.dimension = dim
_, scale_region = fr_utils.estimate_scale(argsc)
argsc.scale_region = scale_region
+ # infile = base_infile + '/DIM{0}/fix_ifrHESESim/'.format(dim)
infile = base_infile + '/DIM{0}/fix_ifr/'.format(dim)
if args.likelihood is Likelihood.GAUSSIAN:
infile += '{0}/'.format(str(args.sigma_ratio).replace('.', '_'))
@@ -282,8 +288,11 @@ def main():
argsc.source_ratio = src
finfile = infile +'/{0}/{1}/fr_stat'.format(args.stat_method, args.run_method) \
+ misc_utils.gen_identifier(argsc)
- basename = os.path.dirname(finfile) + '/statrun/' + os.path.basename(finfile)
+ basename = os.path.dirname(finfile) + '/statrun/'
baseoutfile = basename[:5]+basename[5:].replace('data', 'plots')
+ if args.data:
+ baseoutfile += '/data/'
+ baseoutfile += os.path.basename(finfile)
if args.run_method is SensitivityCateg.FULL:
outfile = baseoutfile
plot_utils.plot_statistic(
@@ -331,11 +340,12 @@ def main():
if args.plot:
print 'Plotting sensitivities'
+ basename = args.infile[:5]+args.infile[5:].replace('data', 'plots')
if args.likelihood in [Likelihood.GOLEMFIT, Likelihood.GF_FREQ]:
- basename = args.infile + 'golemfit/{1}'.format(args.likelihood, args.stat_method)
+ baseoutfile = basename + 'golemfit/{1}'.format(args.likelihood, args.stat_method)
+ if args.data: baseoutfile += '/data'
elif args.likelihood is Likelihood.GAUSSIAN:
- basename = args.infile + 'gaussian/{1}'.format(args.likelihood, args.stat_method)
- baseoutfile = basename[:5]+basename[5:].replace('data', 'plots')
+ baseoutfile = basename + 'gaussian/{1}'.format(args.likelihood, args.stat_method)
if args.run_method is SensitivityCateg.FULL:
plot_utils.plot_sens_full(
diff --git a/submitter/mcmc_dag.py b/submitter/mcmc_dag.py
index d06e920..7a0cde9 100644
--- a/submitter/mcmc_dag.py
+++ b/submitter/mcmc_dag.py
@@ -34,11 +34,11 @@ GLOBAL_PARAMS.update(dict(
))
# FR
-dimension = [3, 6]
+dimension = [3]
# dimension = [4, 5, 7, 8]
GLOBAL_PARAMS.update(dict(
threads = 1,
- binning = '1e4 1e7 5',
+ binning = '6e4 1e7 5',
no_bsm = 'False',
scale_region = "1E10",
energy_dependance = 'spectral',
diff --git a/submitter/mcmc_submit.sub b/submitter/mcmc_submit.sub
index c565205..998932e 100644
--- a/submitter/mcmc_submit.sub
+++ b/submitter/mcmc_submit.sub
@@ -22,6 +22,7 @@ request_cpus = 1
Universe = vanilla
Notification = never
++AccountingGroup="sanctioned.$ENV(USER)"
# run on both SL5 and 6
# +WantRHEL6 = True
# +WantSLC6 = False
diff --git a/submitter/sens_dag.py b/submitter/sens_dag.py
index d0b4752..e02a671 100644
--- a/submitter/sens_dag.py
+++ b/submitter/sens_dag.py
@@ -29,24 +29,25 @@ GLOBAL_PARAMS.update(dict(
sens_run = 'True',
run_method = 'fixed_angle', # full, fixed_angle, corr_angle
stat_method = 'bayesian',
- sens_bins = 40,
+ sens_bins = 20,
seed = 'None'
))
# MultiNest
GLOBAL_PARAMS.update(dict(
- mn_live_points = 800,
- mn_tolerance = 0.01,
+ mn_live_points = 300,
+ mn_tolerance = 0.1,
mn_output = './mnrun'
))
# FR
-dimension = [3, 6]
+dimension = [6]
+# dimension = [3, 6]
# dimension = [4, 5, 7, 8]
# dimension = [3, 4, 5, 6, 7, 8]
GLOBAL_PARAMS.update(dict(
threads = 1,
- binning = '1e4 1e7 5',
+ binning = '6e4 1e7 5',
no_bsm = 'False',
scale_region = "1E10",
energy_dependance = 'spectral',
@@ -73,7 +74,7 @@ GLOBAL_PARAMS.update(dict(
plot_statistic = 'True'
))
-outfile = 'dagman_FR_SENS_{0}_{1}_{2}.submit'.format(
+outfile = 'dagman_FR_SENS_{0}_{1}_{2}_sim.submit'.format(
GLOBAL_PARAMS['stat_method'], GLOBAL_PARAMS['run_method'], GLOBAL_PARAMS['likelihood']
)
golemfitsourcepath = os.environ['GOLEMSOURCEPATH'] + '/GolemFit'
@@ -95,7 +96,10 @@ with open(outfile, 'w') as f:
)
for frs in fix_sfr_mfr:
print 'frs', frs
- output = outchain_head + '/fix_ifr/'
+ # output = outchain_head + '/fix_ifr/'
+ # output = outchain_head + '/fix_ifr/HESESim'
+ output = outchain_head + '/fix_ifr/sim'
+ # output = outchain_head + '/fix_ifr/data'
if GLOBAL_PARAMS['likelihood'].lower() == 'gaussian':
output += '{0}/'.format(str(GLOBAL_PARAMS['sigma_ratio']).replace('.', '_'))
output += 'fr_stat'
diff --git a/submitter/sens_submit.sub b/submitter/sens_submit.sub
index f3a3e4c..5234578 100644
--- a/submitter/sens_submit.sub
+++ b/submitter/sens_submit.sub
@@ -22,6 +22,7 @@ request_cpus = 1
Universe = vanilla
Notification = never
++AccountingGroup="sanctioned.$ENV(USER)"
# run on both SL5 and 6
# +WantRHEL6 = True
# +WantSLC6 = False
@@ -31,9 +32,10 @@ Notification = never
# +TransferOutput=""
++NATIVE_OS = True
# Requirements = IS_GLIDEIN && HAS_CVMFS_icecube_opensciencegrid_org && (OpSysAndVer =?= "CentOS6" || OpSysAndVer =?= "RedHat6" || OpSysAndVer =?= "SL6")
# Requirements = IS_GLIDEIN
-Requirements = (OpSysMajorVer =?= 6)
+# Requirements = (OpSysMajorVer =?= 6)
# GO!
queue
diff --git a/utils/fr.py b/utils/fr.py
index 5acd6f8..64e2706 100644
--- a/utils/fr.py
+++ b/utils/fr.py
@@ -248,7 +248,7 @@ def fr_argparse(parser):
help='Fold in the spectral index when using GolemFit'
)
parser.add_argument(
- '--binning', default=[1e4, 1e7, 5], type=float, nargs=3,
+ '--binning', default=[6e4, 1e7, 5], type=float, nargs=3,
help='Binning for spectral energy dependance'
)
parser.add_argument(
diff --git a/utils/gf.py b/utils/gf.py
index b651b5a..7ded152 100644
--- a/utils/gf.py
+++ b/utils/gf.py
@@ -55,7 +55,8 @@ def fit_flags(llh_paramset):
def steering_params(args):
steering_categ = args.ast
- params = gf.SteeringParams()
+ # params = gf.SteeringParams(gf.sampleTag.HESE)
+ params = gf.SteeringParams(gf.sampleTag.MagicTau)
params.quiet = False
params.fastmode = True
params.simToLoad= steering_categ.name.lower()
@@ -65,12 +66,15 @@ def steering_params(args):
# For Tianlu
# params.years = [999]
+ params.minFitEnergy = 1.0e5 # GeV
+
return params
def set_up_as(fitter, params):
print 'Injecting the model', params
- asimov_params = gf.FitParameters(gf.sampleTag.HESE)
+ # asimov_params = gf.FitParameters(gf.sampleTag.HESE)
+ asimov_params = gf.FitParameters(gf.sampleTag.MagicTau)
for parm in params:
asimov_params.__setattr__(parm.name, float(parm.value))
fitter.SetupAsimov(asimov_params)
@@ -81,12 +85,14 @@ def setup_fitter(args, asimov_paramset):
sparams = steering_params(args)
npp = gf.NewPhysicsParams()
fitter = gf.GolemFit(datapaths, sparams, npp)
- set_up_as(fitter, asimov_paramset)
+ # comment to use data
+ # set_up_as(fitter, asimov_paramset)
return fitter
def get_llh(fitter, params):
- fitparams = gf.FitParameters(gf.sampleTag.HESE)
+ # fitparams = gf.FitParameters(gf.sampleTag.HESE)
+ fitparams = gf.FitParameters(gf.sampleTag.MagicTau)
for parm in params:
fitparams.__setattr__(parm.name, float(parm.value))
llh = -fitter.EvalLLH(fitparams)
@@ -95,7 +101,8 @@ def get_llh(fitter, params):
def get_llh_freq(fitter, params):
print 'setting to {0}'.format(params)
- fitparams = gf.FitParameters(gf.sampleTag.HESE)
+ # fitparams = gf.FitParameters(gf.sampleTag.HESE)
+ fitparams = gf.FitParameters(gf.sampleTag.MagicTau)
for parm in params:
fitparams.__setattr__(parm.name, float(parm.value))
fitter.SetFitParametersSeed(fitparams)
diff --git a/utils/plot.py b/utils/plot.py
index 8792cbb..fec8066 100644
--- a/utils/plot.py
+++ b/utils/plot.py
@@ -313,8 +313,8 @@ def plot_sens_full(data, outfile, outformat, args):
null = statistic[min_idx]
if args.stat_method is StatCateg.BAYESIAN:
reduced_ev = -(statistic - null)
- # al = scales[reduced_ev > np.log(10**(3/2.))] # Strong degree of belief
- al = scales[reduced_ev > 0.4] # Testing
+ al = scales[reduced_ev > np.log(10**(3/2.))] # Strong degree of belief
+ # al = scales[reduced_ev > 0.4] # Testing
elif args.stat_method is StatCateg.FREQUENTIST:
reduced_ev = -2*(statistic - null)
al = scales[reduced_ev > 2.71] # 90% CL for 1 DOF via Wilks
@@ -360,7 +360,8 @@ def plot_sens_fixed_angle(data, outfile, outformat, args):
print 'Making FIXED_ANGLE sensitivity plot'
colour = {0:'red', 1:'blue', 2:'green', 3:'purple', 4:'orange', 5:'black'}
- xticks = [r'$\mathcal{O}_{12}$', r'$\mathcal{O}_{13}$', r'$\mathcal{O}_{23}$']
+ # xticks = [r'$\mathcal{O}_{12}$', r'$\mathcal{O}_{13}$', r'$\mathcal{O}_{23}$']
+ xticks = [r'$\mathcal{O}_{e\mu}$', r'$\mathcal{O}_{e\tau}$', r'$\mathcal{O}_{\mu\tau}$']
argsc = deepcopy(args)
for idim in xrange(len(data)):
dim = args.dimensions[idim]
@@ -390,6 +391,7 @@ def plot_sens_fixed_angle(data, outfile, outformat, args):
if args.stat_method is StatCateg.BAYESIAN:
reduced_ev = -(statistic - null)
al = scales[reduced_ev > np.log(10**(3/2.))] # Strong degree of belief
+ # al = scales[reduced_ev > np.log(10**(1/2.))]
elif args.stat_method is StatCateg.FREQUENTIST:
reduced_ev = -2*(statistic - null)
al = scales[reduced_ev > 2.71] # 90% CL for 1 DOF via Wilks
@@ -403,7 +405,8 @@ def plot_sens_fixed_angle(data, outfile, outformat, args):
print 'limit = {0}'.format(lim)
label = '[{0}, {1}, {2}]'.format(*misc_utils.solve_ratio(src))
if lim < yranges[0]: yranges[0] = lim
- if lim > yranges[1]: yranges[1] = lim+4
+ # if lim > yranges[1]: yranges[1] = lim+5
+ if lim > yranges[1]: yranges[1] = lim
line = plt.Line2D(
(ian+1-0.1, ian+1+0.1), (lim, lim), lw=3, color=colour[isrc], label=label
)
@@ -417,7 +420,7 @@ def plot_sens_fixed_angle(data, outfile, outformat, args):
)
try:
- yranges = (myround(yranges[0], up=True), myround(yranges[1], down=True))
+ # yranges = (myround(yranges[0], up=True), myround(yranges[1], down=True))
ax.set_ylim(yranges)
except: pass
@@ -433,7 +436,6 @@ def plot_sens_fixed_angle(data, outfile, outformat, args):
print 'Saving plot as {0}'.format(out+'.'+of)
fig.savefig(out+'.'+of, bbox_inches='tight', dpi=150)
-
for ymaj in ax.yaxis.get_majorticklocs():
ax.axhline(y=ymaj, ls=':', color='gray', alpha=0.4, linewidth=1)
for xmaj in ax.xaxis.get_majorticklocs():