prongcvn_optimization.C
Go to the documentation of this file.
1 #pragma once
2 
4 
5 // analysis cuts/vars
10 
11 // for data pot
13 
14 // Multiverse includes
16 #include "CAFAna/Vars/XsecTunes.h"
20 
21 using namespace ana;
22 
24 {
25 
26  std::map<TString, SpectrumLoader *> loaders;
27  loaders["kNominal"] = new SpectrumLoader(NominalMC_entire_RHC_prod4);
28  loaders["kCalibNeg"] = new SpectrumLoader(DWNCalibFiles_full_RHC_prod4);
29  loaders["kCalibPos"] = new SpectrumLoader(UPCalibFiles_full_RHC_prod4);
30  loaders["kCalibShape"] = new SpectrumLoader(ShapeCalibFiles_full_RHC_prod4);
31  loaders["kLightUpCalibDown"] = new SpectrumLoader(UPLightFiles_full_RHC_prod4);
32  loaders["kLightDownCalibUp"] = new SpectrumLoader(DWNLightFiles_full_RHC_prod4);
33  loaders["kCherenkov"] = new SpectrumLoader(CherenkovFiles_full_RHC_prod4);
34 
35  // define a hist axis to optimize
36  const HistAxis * opti_axis = new HistAxis("Max ProngCVN Electron Score",
37  Binning::Simple(50, 0, 1),
38  kvProngCVN);
39  // Define selection and signal cuts
40  const Cut * kSignal = new Cut(CutFromNuTruthCut(kIsNuebarCCST && kTrueDetectorST));
41 
42  const Cut * kSelection = new Cut(kcDQ &&
43  kcFrontPlanes &&
44  kcNHits &&
45  kTrueDetector &&
47  kcNueCCIncFiducialLoose &&
48  kcMuonIDCut);
49 
50  // central value weight to be applied to all samples
52 
53  // nominal definition
54  SystematicDef * sNominal =
55  new SystematicDef("sNominal",
56  loaders["kNominal"],
57  opti_axis,
58  kCVWeight);
59 
60  // systematic definitions
61  SystematicDef * sLightlevelUpDown =
62  new SystematicDef("sLightlevelUpDown",
63  loaders["kLightUpCalibDown"],
64  loaders["kLightDownCalibUp"],
65  opti_axis,
66  kCVWeight);
67  SystematicDef * sCalibPosNeg =
68  new SystematicDef("sCalibPosNeg",
69  loaders["kCalibPos"],
70  loaders["kCalibNeg"],
71  opti_axis,
72  kCVWeight);
73  SystematicDef * sCalibShape =
74  new SystematicDef("sCalibShape",
75  loaders["kCalibShape"],
76  opti_axis,
77  kCVWeight);
78  SystematicDef * sCherenkov =
79  new SystematicDef("sCherenkov",
80  loaders["kCherenkov"],
81  opti_axis,
82  kCVWeight);
83 
84  // Setup the genie multiverses
85  // We take flux uncertainty to be constant, so no need for PPFX universe
86  GenieMultiverseParameters genie_multiverse(100, getAllXsecSysts_2018());
87  std::vector<SystShifts> genie_shifts = genie_multiverse.GetSystShifts();
88 
89  SystematicDef * sGenieMultiverse =
90  new SystematicDef("sGenieMultiverse",
91  loaders["kNominal"],
92  opti_axis,
93  kCVWeight);
94 
95  // create the optimizer.
96  CutOptimization * optimize = new CutOptimization(opti_axis,
97  kSignal,
98  kSelection);
99  optimize->DefineNominal(sNominal);
100  optimize->DefineSystematic(sLightlevelUpDown);
101  optimize->DefineSystematic(sCalibPosNeg);
102  optimize->DefineSystematic(sCalibShape);
103  optimize->DefineSystematic(sCherenkov);
104  optimize->DefineMultiverseSystematic(sGenieMultiverse, genie_shifts);
105 
106  // optionally set spill cuts for all loaders
107  optimize->SetSpillCuts(kStandardSpillCuts);
108 
109  // Let the loaders go
110  optimize->Go();
111 
112  // save to file
113  TFile * output = new TFile("prongcvn_optimization.root", "recreate");
114  optimize->SaveTo(output, "prongcvn_optimize");
115  output->Close();
116 
117 }
const Cut kcNueCCIncContainment([](const caf::SRProxy *sr){if(sr->vtx.elastic[0].fuzzyk.nshwlid< 1) return false;TVector3 start=sr->vtx.elastic[0].fuzzyk.png[0].shwlid.start;TVector3 stop=sr->vtx.elastic[0].fuzzyk.png[0].shwlid.stop;for(uint ix=0;ix< sr->vtx.elastic[0].fuzzyk.nshwlid;ix++){TVector3 start_muoncatcher=sr->vtx.elastic[0].fuzzyk.png[ix].shwlid.start;TVector3 stop_muoncatcher=sr->vtx.elastic[0].fuzzyk.png[ix].shwlid.stop;if(std::max(start_muoncatcher.Z(), stop_muoncatcher.Z()) > 1250.0) return false;}if(sr->sel.nuecosrej.distallpngtop< 50) return false;if(sr->sel.nuecosrej.distallpngbottom< 30) return false;if(sr->sel.nuecosrej.distallpngeast< 50) return false;if(sr->sel.nuecosrej.distallpngwest< 30) return false;if(sr->sel.nuecosrej.distallpngfront< 150) return false;return true;})
ofstream output
const std::string NominalMC_entire_RHC_prod4
const Cut kcFrontPlanes
Definition: NueCCIncCuts.h:117
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const Cut kcNHits
Definition: NueCCIncCuts.h:118
std::vector< const ISyst * > getAllXsecSysts_2018()
Get master XSec syst list for 2018 analyses.
const Var kPPFXFluxCVWgt
Definition: PPFXWeights.h:16
const NuTruthCut kTrueDetectorST([](const caf::SRNeutrinoProxy *sr){return(sr->vtx.X()< tvtxmax.X()&& sr->vtx.X() > tvtxmin.X()&& sr->vtx.Y() > tvtxmin.Y()&& sr->vtx.Y()< tvtxmax.Y()&& sr->vtx.Z() > tvtxmin.Z()&& sr->vtx.Z()< tvtxmax.Z());})
Definition: NumuCCIncCuts.h:32
const ana::Var kCVWeight
void prongcvn_optimization()
const std::string DWNLightFiles_full_RHC_prod4
void DefineNominal(SystematicDef *nominal)
const Cut kcMuonIDCut
Definition: NueCCIncCuts.h:123
const std::string ShapeCalibFiles_full_RHC_prod4
const Cut kcDQ
Definition: NueCCIncCuts.h:114
const std::string UPCalibFiles_full_RHC_prod4
void SaveTo(TDirectory *dir, const std::string &name) const
void DefineMultiverseSystematic(SystematicDef *syst, std::vector< Var > mv_weights)
_Var< caf::SRProxy > Var
Definition: Var.h:7
const Var kvProngCVN([](const caf::SRProxy *sr){if(sr->vtx.nelastic< 1) return-1000.f;if(sr->vtx.elastic[0].fuzzyk.npng< 1) return-1000.f;int nProngs=sr->vtx.elastic[0].fuzzyk.npng;float maxProngCVNe=-5.;for(int iProng=0;iProng< nProngs;iProng++){if(sr->vtx.elastic[0].fuzzyk.png[iProng].cvnpart.electronid > maxProngCVNe){maxProngCVNe=sr->vtx.elastic[0].fuzzyk.png[iProng].cvnpart.electronid;}}return maxProngCVNe;})
GenericSystematicDef< Cut, Var > SystematicDef
void SetSpillCuts(const SpillCut &)
const std::string UPLightFiles_full_RHC_prod4
const Cut kSignal
Definition: SINCpi0_Cuts.h:383
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
_Cut< caf::SRProxy, caf::SRSpillProxy > Cut
Definition: Cut.h:9
void DefineSystematic(SystematicDef *syst)
const Var kXSecCVWgt2018
Definition: XsecTunes.h:48
std::vector< Loaders * > loaders
Definition: syst_header.h:386
std::vector< SystShifts > GetSystShifts()
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
const std::string DWNCalibFiles_full_RHC_prod4
const std::string CherenkovFiles_full_RHC_prod4
const Cut kTrueDetector
static Binning Simple(int n, double lo, double hi, const std::vector< std::string > &labels={})
Definition: Binning.cxx:107
_HistAxis< Var > HistAxis
Definition: HistAxis.h:10
Cut CutFromNuTruthCut(const NuTruthCut &stc)
Definition: Cut.cxx:7