diff options
| -rwxr-xr-x | plot_sens.py | 18 | ||||
| -rw-r--r-- | submitter/mcmc_dag.py | 4 | ||||
| -rw-r--r-- | submitter/mcmc_submit.sub | 1 | ||||
| -rw-r--r-- | submitter/sens_dag.py | 18 | ||||
| -rw-r--r-- | submitter/sens_submit.sub | 4 | ||||
| -rw-r--r-- | utils/fr.py | 2 | ||||
| -rw-r--r-- | utils/gf.py | 17 | ||||
| -rw-r--r-- | utils/plot.py | 14 |
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(): |
