22 bool nd = opt.Contains(
"neardet");
26 if (nd) loaders.
SetND(
true);
33 std::vector<std::tuple<SpectrumLoaderBase*, const Cut, std::string>>
components;
41 components.push_back(std::make_tuple(nonswap,
kIsBeamNue &&
kIsAntiNu,
"cc_nuebar_to_nuebar" ));
42 components.push_back(std::make_tuple(nonswap,
kIsNumuCC && !
kIsAntiNu,
"cc_numu_to_numu" ));
43 components.push_back(std::make_tuple(nonswap,
kIsNumuCC &&
kIsAntiNu,
"cc_numubar_to_numubar"));
44 components.push_back(std::make_tuple(nonswap,
kIsNC,
"nc" ));
53 components.push_back(std::make_tuple(swap,
kIsSig && !
kIsAntiNu,
"cc_numu_to_nue" ));
54 components.push_back(std::make_tuple(swap,
kIsSig &&
kIsAntiNu,
"cc_numubar_to_nuebar" ));
56 components.push_back(std::make_tuple(swap,
kIsNumuApp &&
kIsAntiNu,
"cc_nuebar_to_numubar" ));
65 std::vector<Var> kPPFXFluxUnivWgts;
66 for (
size_t i = 0;
i < 100; ++
i) {
69 if (sr->
mc.
nnu != 1)
return 1.f;
70 if (sr->
mc.
nu[0].rwgt.ppfx.vuniv[
i] <= 0)
return 1.f;
71 return (
float)sr->
mc.
nu[0].rwgt.ppfx.vuniv[
i];
73 kPPFXFluxUnivWgts.push_back(tempPPFXWgt);
88 std::vector<Spectrum*> nominal(components.size());
89 std::vector<std::vector<Spectrum*> > ppfx(components.size(),
90 std::vector<Spectrum*>(kPPFXFluxUnivWgts.size(),
nullptr));
91 std::vector<std::vector<Spectrum*> > genie_up(components.size(),
92 std::vector<Spectrum*>(systs.size(),
nullptr));
93 std::vector<std::vector<Spectrum*> > genie_down(components.size(),
94 std::vector<Spectrum*>(systs.size(),
nullptr));
98 for (
size_t c = 0;
c < components.size(); ++
c) {
102 nominal[
c] =
new Spectrum(*std::get<0>(components[
c]),
103 axis, sel && std::get<1>(components[c]),
kNoShift, kReweight);
105 for (
size_t u = 0;
u < 100; ++
u)
106 ppfx[c][
u] =
new Spectrum(*std::get<0>(components[c]),
109 for (
size_t s = 0;
s < systs.size(); ++
s) {
113 genie_up[
c][
s] =
new Spectrum(*std::get<0>(components[c]),
114 axis, sel && std::get<1>(components[c]), shift_up, kReweight);
117 genie_down[
c][
s] =
new Spectrum(*std::get<0>(components[c]),
118 axis, sel && std::get<1>(components[c]), shift_down, kReweight);
128 TFile*
f = TFile::Open(
"SystRatios.root",
"recreate");
130 TDirectory*
dir = f->mkdir(
"nominal");
131 for (
size_t c = 0;
c < components.size(); ++
c)
132 nominal[
c]->
SaveTo(dir, std::get<2>(components[
c]).c_str());
134 dir = f->mkdir(
"ppfx");
135 for (
size_t u = 0;
u < 100; ++
u) {
136 TDirectory* udir = dir->mkdir(Form(
"u%zu",
u+1));
137 for (
size_t c = 0; c < components.size(); ++
c)
138 ppfx[c][
u]->
SaveTo(udir, std::get<2>(components[c]).c_str());
141 dir = f->mkdir(
"genie");
142 for (
size_t s = 0;
s < systs.size(); ++
s) {
143 TDirectory* sdir = dir->mkdir(systs[
s]->ShortName().c_str());
144 TDirectory* updir = sdir->mkdir(
"up");
145 TDirectory* downdir = sdir->mkdir(
"down");
146 for (
size_t c = 0; c < components.size(); ++
c) {
147 genie_up[
c][
s]->SaveTo(updir, std::get<2>(components[c]).c_str());
148 genie_down[
c][
s]->SaveTo(downdir, std::get<2>(components[c]).c_str());
Near Detector underground.
const Cut kIsNumuCC(CCFlavSel(14, 14))
Select CC .
Far Detector at Ash River.
Cuts and Vars for the 2020 FD DiF Study.
std::vector< const ISyst * > getAllXsecSysts_2018()
Get master XSec syst list for 2018 analyses.
void SaveTo(const osc::IOscCalc &x, TDirectory *dir, const std::string &name)
Simple record of shifts applied to systematic parameters.
caf::Proxy< caf::SRHeader > hdr
Proxy for caf::StandardRecord.
caf::Proxy< std::vector< caf::SRNeutrino > > nu
const Cut kIsAntiNu([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;assert(sr->mc.nnu==1);return sr->mc.nu[0].pdg< 0;})
Is this truly an antineutrino?
const Cut kIsBeamNue(CCFlavSel(12, 12))
Select CC .
const Cut kIsNumuApp(CCFlavSel(14, 12))
Select CC .
caf::Proxy< short int > nnu
void Go()
Call Go() on all the loaders.
void SetSpillCut(const SpillCut &cut)
const Cut kIsTauFromE(CCFlavSel(16, 12))
Select CC .
const Cut kIsTauFromMu(CCFlavSel(16, 14))
Select CC .
const Cut kIsNC([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;assert(sr->mc.nnu==1);return!sr->mc.nu[0].iscc;})
Is this a Neutral Current event?
SpectrumLoaderBase & GetLoader(caf::Det_t det, DataMC datamc, DataSource src=kBeam, SwappingConfig swap=kNonSwap)
Retrieve a specific loader.
void MakeSystRatios(TString opt)
const Cut kIsSig(CCFlavSel(12, 14))
Select CC .
const Cut kNus18ND
Full Nus18 ND analysis selection.
std::vector< float > Spectrum
const SystShifts kNoShift
Base class for the various types of spectrum loader.
caf::Proxy< caf::SRTruthBranch > mc
const SpillCut kStandardSpillCuts
Apply this unless you're doing something special.
For nominal spectra and reweighting systs (xsec/flux)
void SetShift(const ISyst *syst, double shift, bool force=false)