Classes | Public Types | Public Member Functions | Protected Member Functions | Static Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
ana::NuWROSyst Class Reference

Reweight events to match NuWRO. More...

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N20-12-01/CAFAna/Systs/NuWROSysts.h"

Inheritance diagram for ana::NuWROSyst:
ana::ISyst

Classes

struct  Key_t
 

Public Types

enum  EReweightVar {
  EReweightVar::kEnu, EReweightVar::kX, EReweightVar::kY, EReweightVar::kQ2,
  EReweightVar::kMode
}
 The variable to reweight as a function of. More...
 

Public Member Functions

 NuWROSyst (EReweightVar rwVar)
 
 ~NuWROSyst ()
 
void Shift (double sigma, caf::SRProxy *sr, double &weight) const override
 +1 means match nuWRO, -1 moves in the opposite direction More...
 
 NuWROSyst (const NuWROSyst &)=delete
 
NuWROSystoperator= (const NuWROSyst &)=delete
 
virtual const std::stringShortName () const final
 The name printed out to the screen. More...
 
virtual const std::stringLatexName () const final
 The name used on plots (ROOT's TLatex syntax) More...
 
virtual void TruthShift (double sigma, caf::SRNeutrinoProxy *nu, double &weight) const
 
virtual bool IsGenieReweight () const
 GENIE reweights can only provide +/-1,2sigma. More...
 

Protected Member Functions

std::string ShortNameHelper (EReweightVar v) const
 
std::string LatexNameHelper (EReweightVar v) const
 

Static Protected Member Functions

static void InitializeHistograms ()
 

Protected Attributes

EReweightVar fReweightVar
 

Static Protected Attributes

static std::map< Key_t, TH1 * > fgGenieHists
 
static std::map< Key_t, TH1 * > fgNuWroHists
 
static int fgNInstances = 0
 

Detailed Description

Reweight events to match NuWRO.

Definition at line 12 of file NuWROSysts.h.

Member Enumeration Documentation

The variable to reweight as a function of.

Enumerator
kEnu 
kX 
kY 
kQ2 
kMode 

Definition at line 16 of file NuWROSysts.h.

16  {
17  kEnu,
18  kX,
19  kY,
20  kQ2,
21  kMode
22  };
const Var kMode([](const caf::SRProxy *sr){return(sr->mc.nnu==0)?-1:int(sr->mc.nu[0].mode);})
Neutrino interaction mode.
Definition: Vars.h:99
True neutrino energy (GeV)
const Var kY([](const caf::SRProxy *sr){float tmp=0.f;if(sr->mc.nu.empty()) return tmp;tmp=sr->mc.nu[0].y;return tmp;})
True vertex x position (cm)

Constructor & Destructor Documentation

ana::NuWROSyst::NuWROSyst ( EReweightVar  rwVar)

Definition at line 28 of file NuWROSysts.cxx.

References fgNInstances.

29  : ISyst(ShortNameHelper(rwVar), LatexNameHelper(rwVar)),
30  fReweightVar(rwVar)
31  {
32  ++fgNInstances;
33  }
static int fgNInstances
Definition: NuWROSysts.h:57
EReweightVar fReweightVar
Definition: NuWROSysts.h:38
std::string ShortNameHelper(EReweightVar v) const
Definition: NuWROSysts.cxx:127
ISyst(const std::string &shortName, const std::string &latexName)
Definition: ISyst.cxx:10
std::string LatexNameHelper(EReweightVar v) const
Definition: NuWROSysts.cxx:140
ana::NuWROSyst::~NuWROSyst ( )

Definition at line 36 of file NuWROSysts.cxx.

References fgNInstances.

37  {
38  --fgNInstances;
39  // Crashes at shutdown. Maybe because load_libs means everything winds up
40  // loaded twice? Leak isn't that big, just leave it.
41  /*
42  if(fgNInstances == 0){
43  // Only clean up when there are no instances left
44  for(auto it: fgGenieHists) delete it.second;
45  for(auto it: fgNuWroHists) delete it.second;
46  // In case a new instance is created
47  fgGenieHists.clear();
48  fgNuWroHists.clear();
49  }
50  if(fgNInstances < 0){
51  std::cerr << "NuWROSysts::fgNInstances is " << fgNInstances
52  << ". This should never happen" << std::endl;
53  }
54  */
55  }
static int fgNInstances
Definition: NuWROSysts.h:57
ana::NuWROSyst::NuWROSyst ( const NuWROSyst )
delete

Member Function Documentation

void ana::NuWROSyst::InitializeHistograms ( )
staticprotected

Definition at line 58 of file NuWROSysts.cxx.

References ana::assert(), cc(), ana::Progress::Done(), MakeMiniprodValidationCuts::f, fgGenieHists, fgNuWroHists, genie::utils::style::Format(), hists, MECModelEnuComparisons::i, caf::kDIS, kEnu, kMode, kQ2, caf::kQE, caf::kRes, kX, kY, submit_nova_art::mode, make_root_from_grid_output::pdg, cacheDefinitionData::prog, genie::utils::kinematics::Q2(), ana::Progress::SetProgress(), make_root_from_grid_output::tr, submit_syst::x, and submit_syst::y.

Referenced by Shift().

59  {
60  // Already called
61  if(!fgGenieHists.empty()) return;
62 
63  Progress prog("Loading NuWRO histograms");
64 
65  for(int pdg: {12, 14}){
66 
67  TFile f(TString::Format("$NOVA_ANA/caf/genie_vs_nuwro_%s.root", (pdg == 12) ? "nue" : "numu"));
68 
69  assert(!f.IsZombie());
70 
71  DontAddDirectory guard;
72 
73  for(bool genie: {false, true}){
74  TTree* tr = (TTree*)f.Get(genie ? "genie" : "nuwro");
75  assert(tr);
76  std::map<Key_t, TH1*>& hists = genie ? fgGenieHists : fgNuWroHists;
77 
78  for(int mode = 0; mode < 4; ++mode){
79  for(bool cc: {false, true}){
80  hists[{pdg, mode, cc, EReweightVar::kEnu }] = new TH1F("", "", 100, 0, 10);
81  hists[{pdg, mode, cc, EReweightVar::kX }] = new TH1F("", "", 100, 0, 5);
82  hists[{pdg, mode, cc, EReweightVar::kY }] = new TH1F("", "", 100, 0, 1);
83  hists[{pdg, mode, cc, EReweightVar::kQ2 }] = new TH1F("", "", 100, 0, 5);
84  hists[{pdg, mode, cc, EReweightVar::kMode}] = new TH1F("", "", 1, 0, 1);
85  }
86  }
87 
88  bool qel, res, dis, cc;
89  tr->SetBranchAddress("qel", &qel);
90  tr->SetBranchAddress("res", &res);
91  tr->SetBranchAddress("dis", &dis);
92  tr->SetBranchAddress("cc", &cc);
93 
94  double Ev, x, y, Q2;
95  tr->SetBranchAddress("Ev", &Ev);
96  tr->SetBranchAddress("x", &x);
97  tr->SetBranchAddress("y", &y);
98  tr->SetBranchAddress("Q2", &Q2);
99 
100  const unsigned int I = tr->GetEntries();
101  for(unsigned int i = 0; i < I; ++i){
102  tr->GetEntry(i);
103 
104  int mode = -1;
105  if(qel) mode = caf::kQE;
106  if(res) mode = caf::kRes;
107  if(dis) mode = caf::kDIS;
108  if(mode == -1) continue; // coherent?
109 
110  hists[{pdg, mode, cc, EReweightVar::kEnu }]->Fill(Ev);
111  hists[{pdg, mode, cc, EReweightVar::kX }]->Fill(x );
112  hists[{pdg, mode, cc, EReweightVar::kY }]->Fill(y );
113  hists[{pdg, mode, cc, EReweightVar::kQ2 }]->Fill(Q2);
114  hists[{pdg, mode, cc, EReweightVar::kMode}]->Fill(.5);
115 
116  if(i%10000 == 0) prog.SetProgress(i/(4.*I) +
117  (genie ? .25 : 0) +
118  ((pdg == 14) ? .5 : 0));
119  } // end for i
120  } // end for genie
121  } // end for pdg
122 
123  prog.Done();
124  }
THE MAIN GENIE PROJECT NAMESPACE
Definition: GeneratorBase.h:8
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:991
static std::map< Key_t, TH1 * > fgNuWroHists
Definition: NuWROSysts.h:56
TString hists[nhists]
Definition: bdt_com.C:3
void cc()
Definition: test_ana.C:28
assert(nhit_max >=nhit_nbins)
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
static std::map< Key_t, TH1 * > fgGenieHists
Definition: NuWROSysts.h:56
virtual bool ana::ISyst::IsGenieReweight ( ) const
inlinevirtualinherited

GENIE reweights can only provide +/-1,2sigma.

Reimplemented in ana::SummedSyst.

Definition at line 56 of file ISyst.h.

56 {return false;}
virtual const std::string& ana::ISyst::LatexName ( ) const
inlinefinalvirtualinherited

The name used on plots (ROOT's TLatex syntax)

Definition at line 30 of file ISyst.h.

References ana::ISyst::fLatexName, ana::ISyst::Shift(), sigma(), sr, and ana::weight.

Referenced by ana::PredictionInterp::DebugPlotColz(), GetGENIEShiftLabels(), ana::NuISyst::SaveTo(), SystsGENIEAna(), and WriteSystName().

30 {return fLatexName;}
std::string fLatexName
Definition: ISyst.h:60
std::string ana::NuWROSyst::LatexNameHelper ( EReweightVar  v) const
protected

Definition at line 140 of file NuWROSysts.cxx.

References ana::assert(), kEnu, kMode, kQ2, kX, and kY.

141  {
142  switch(v){
143  case EReweightVar::kEnu : return "NuWRO reweight on neutrino energy";
144  case EReweightVar::kX : return "NuWRO reweight on hadronic x";
145  case EReweightVar::kY : return "NuWRO reweight on hadronic y";
146  case EReweightVar::kQ2 : return "NuWRO reweight on Q^{2}";
147  case EReweightVar::kMode: return "NuWRO reweight on interaction mode";
148  }
149  assert(0 && "Unknown reweight");
150  }
assert(nhit_max >=nhit_nbins)
NuWROSyst& ana::NuWROSyst::operator= ( const NuWROSyst )
delete
void ana::NuWROSyst::Shift ( double  sigma,
caf::SRProxy sr,
double &  weight 
) const
overridevirtual

+1 means match nuWRO, -1 moves in the opposite direction

Reimplemented from ana::ISyst.

Definition at line 153 of file NuWROSysts.cxx.

References abs(), ana::assert(), caf::Proxy< caf::SRNeutrino >::E, fgGenieHists, fgNuWroHists, fReweightVar, GetXaxis(), InitializeHistograms(), caf::Proxy< caf::SRNeutrino >::iscc, caf::kDIS, kEnu, findDuplicateFiles::key, kMode, kQ2, caf::kQE, kX, kY, caf::Proxy< caf::StandardRecord >::mc, caf::Proxy< caf::SRNeutrino >::mode, caf::Proxy< caf::SRTruthBranch >::nnu, caf::Proxy< caf::SRTruthBranch >::nu, caf::Proxy< caf::SRNeutrino >::pdg, caf::Proxy< caf::SRNeutrino >::q2, caf::Proxy< caf::SRNeutrino >::x, and caf::Proxy< caf::SRNeutrino >::y.

155  {
157 
158  if(sr->mc.nnu == 0) return;
159  //Leave weight unaltered
160  if(sigma == 0) return;
161 
162  const caf::SRNeutrinoProxy& nu = sr->mc.nu[0];
163 
164  if(nu.mode < caf::kQE || nu.mode > caf::kDIS) return; // probably COH
165 
166  if(abs(nu.pdg) != 12 && abs(nu.pdg) != 14) return; // don't reweight taus
167 
168  const Key_t key = {abs(nu.pdg), nu.mode, nu.iscc, fReweightVar};
169 
170  TH1* hgen = fgGenieHists[key];
171  TH1* hnuw = fgNuWroHists[key];
172 
173  assert(hgen && hnuw);
174 
175  double nuvar = -1;
176  if(fReweightVar == EReweightVar::kEnu ) nuvar = nu.E;
177  if(fReweightVar == EReweightVar::kX ) nuvar = nu.x;
178  if(fReweightVar == EReweightVar::kY ) nuvar = nu.y;
179  if(fReweightVar == EReweightVar::kQ2 ) nuvar = nu.q2;
180  if(fReweightVar == EReweightVar::kMode) nuvar = .5;
181  assert(nuvar != -1);
182 
183  // Off our axis range. Don't reweight
184  if(nuvar < hgen->GetXaxis()->GetXmin()) return;
185  if(nuvar > hgen->GetXaxis()->GetXmax()) return;
186 
187  const double wGenie = hgen->GetBinContent(hgen->FindBin(nuvar));
188  const double wNuWro = hnuw->GetBinContent(hnuw->FindBin(nuvar));
189 
190  // Either division by zero or 0/0. Don't reweight
191  if(wGenie == 0) return;
192 
193  // Do the actual reweight
194  weight *= 1 + (wNuWro/wGenie - 1)*sigma;
195 
196  if(weight < 0) weight = 0;
197  }
caf::Proxy< float > x
Definition: SRProxy.h:575
static void InitializeHistograms()
Definition: NuWROSysts.cxx:58
const Var weight
Proxy for caf::SRNeutrino.
Definition: SRProxy.h:510
static std::map< Key_t, TH1 * > fgNuWroHists
Definition: NuWROSysts.h:56
caf::Proxy< std::vector< caf::SRNeutrino > > nu
Definition: SRProxy.h:618
void abs(TH1 *hist)
caf::Proxy< int > mode
Definition: SRProxy.h:543
caf::Proxy< short int > nnu
Definition: SRProxy.h:617
caf::Proxy< float > E
Definition: SRProxy.h:520
EReweightVar fReweightVar
Definition: NuWROSysts.h:38
correl_xv GetXaxis() -> SetDecimals()
caf::Proxy< bool > iscc
Definition: SRProxy.h:538
caf::Proxy< float > y
Definition: SRProxy.h:577
double sigma(TH1F *hist, double percentile)
caf::Proxy< caf::SRTruthBranch > mc
Definition: SRProxy.h:2138
assert(nhit_max >=nhit_nbins)
caf::Proxy< float > q2
Definition: SRProxy.h:557
caf::Proxy< short int > pdg
Definition: SRProxy.h:552
static std::map< Key_t, TH1 * > fgGenieHists
Definition: NuWROSysts.h:56
virtual const std::string& ana::ISyst::ShortName ( ) const
inlinefinalvirtualinherited
std::string ana::NuWROSyst::ShortNameHelper ( EReweightVar  v) const
protected

Definition at line 127 of file NuWROSysts.cxx.

References ana::assert(), kEnu, kMode, kQ2, kX, and kY.

128  {
129  switch(v){
130  case EReweightVar::kEnu : return "nuwroEnu";
131  case EReweightVar::kX : return "nuwrox";
132  case EReweightVar::kY : return "nuwroy";
133  case EReweightVar::kQ2 : return "nuwroq2";
134  case EReweightVar::kMode: return "nuwroqmode";
135  }
136  assert(0 && "Unknown reweight");
137  }
assert(nhit_max >=nhit_nbins)
virtual void ana::ISyst::TruthShift ( double  sigma,
caf::SRNeutrinoProxy nu,
double &  weight 
) const
inlinevirtualinherited

For systematics that deal only with the neutrino truth and not any reconstruction/PID details. Systematics defined this way will work on nuTree-derived spectra too (e.g. denominators of efficiencies).

Reimplemented in demo::DemoSyst1, ana::BeamSyst, ana::GenericSystComponentScale< T >, ana::GenericSystComponentScale< T >, ana::GenericSystComponentScale< T >, and ana::NOvARwgtSyst.

Definition at line 46 of file ISyst.h.

Referenced by ana::ISyst::Shift().

49  {
50  // Implement this function if your systematic depends only
51  // SRNeutrino. Left blank by default, since systematics using other
52  // information can do nothing sensible to the nuTree.
53  }

Member Data Documentation

std::map< NuWROSyst::Key_t, TH1 * > ana::NuWROSyst::fgGenieHists
staticprotected

Definition at line 56 of file NuWROSysts.h.

Referenced by InitializeHistograms(), and Shift().

int ana::NuWROSyst::fgNInstances = 0
staticprotected

Definition at line 57 of file NuWROSysts.h.

Referenced by NuWROSyst(), and ~NuWROSyst().

std::map< NuWROSyst::Key_t, TH1 * > ana::NuWROSyst::fgNuWroHists
staticprotected

Definition at line 56 of file NuWROSysts.h.

Referenced by InitializeHistograms(), and Shift().

EReweightVar ana::NuWROSyst::fReweightVar
protected

Definition at line 38 of file NuWROSysts.h.

Referenced by Shift().


The documentation for this class was generated from the following files: