Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | List of all members
ana::GeniePCASyst Class Reference

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N21-03-03/3FlavorAna/Systs/GeniePCASyst.h"

Inheritance diagram for ana::GeniePCASyst:
ana::ISyst

Public Member Functions

 GeniePCASyst (const std::string &fname, const std::string &shortname, const std::string &latexname)
 
 ~GeniePCASyst ()
 
void Shift (double sigma, caf::SRProxy *sr, double &weight) const override
 Perform the systematic shift. More...
 
TH1D * getHist (int f, int d, int s, int fl, int p)
 
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 Types

enum  { kFHC, kRHC, kNumFluxType }
 
enum  { kND, kFD, kNumDets }
 
enum  { kNueSel, kNumuSel, kNumSels }
 
enum  {
  kSigQE, kSigNonQE, kBkgFlav, kNC,
  kNumFlavors
}
 
enum  ESign { kMinus, kPlus, kNumSigns }
 

Protected Member Functions

void InitializeHistograms () const
 
double GetWeight (const caf::SRProxy *sr, ESign sign) const
 

Protected Attributes

std::string fFileName
 
std::string fHistName
 
TH1D * fHistos [kNumFluxType][kNumDets][kNumSels][kNumFlavors][kNumSigns]
 

Detailed Description

Definition at line 11 of file GeniePCASyst.h.

Member Enumeration Documentation

anonymous enum
protected
Enumerator
kFHC 
kRHC 
kNumFluxType 

Definition at line 27 of file GeniePCASyst.h.

anonymous enum
protected
Enumerator
kND 
kFD 
kNumDets 

Definition at line 28 of file GeniePCASyst.h.

anonymous enum
protected
Enumerator
kNueSel 
kNumuSel 
kNumSels 

Definition at line 29 of file GeniePCASyst.h.

anonymous enum
protected
enum ana::GeniePCASyst::ESign
protected
Enumerator
kMinus 
kPlus 
kNumSigns 

Definition at line 31 of file GeniePCASyst.h.

Constructor & Destructor Documentation

ana::GeniePCASyst::GeniePCASyst ( const std::string fname,
const std::string shortname,
const std::string latexname 
)

Definition at line 111 of file GeniePCASyst.cxx.

114  : ISyst(shortname,latexname),
115  fFileName(fname), fHistName(shortname), fHistos()
116  {
117  }
std::string fHistName
Definition: GeniePCASyst.h:38
TH1D * fHistos[kNumFluxType][kNumDets][kNumSels][kNumFlavors][kNumSigns]
Definition: GeniePCASyst.h:40
std::string fFileName
Definition: GeniePCASyst.h:37
ISyst(const std::string &shortName, const std::string &latexName)
Definition: ISyst.cxx:10
ana::GeniePCASyst::~GeniePCASyst ( )

Definition at line 119 of file GeniePCASyst.cxx.

References fHistos, MECModelEnuComparisons::i, calib::j, kNumDets, kNumFlavors, kNumFluxType, kNumSels, kNumSigns, submit_hadd::l, and m.

120  {
121  for(int i = 0; i < kNumFluxType; ++i){
122  for(int j = 0; j < kNumDets; ++j){
123  for(int k = 0; k < kNumSels; ++k){
124  for(int l = 0; l < kNumFlavors; ++l){
125  for(int m = 0; m < kNumSigns; ++m){
126  delete fHistos[i][j][k][l][m];
127  fHistos[i][j][k][l][m] = 0;
128  }
129  }
130  }
131  }
132  }
133  }
TH1D * fHistos[kNumFluxType][kNumDets][kNumSels][kNumFlavors][kNumSigns]
Definition: GeniePCASyst.h:40
const double j
Definition: BetheBloch.cxx:29

Member Function Documentation

TH1D * ana::GeniePCASyst::getHist ( int  f,
int  d,
int  s,
int  fl,
int  p 
)

Definition at line 278 of file GeniePCASyst.cxx.

References d, MakeMiniprodValidationCuts::f, fHistos, and InitializeHistograms().

279  {
281 
282  return fHistos[f][d][s][fl][p];
283  }
void InitializeHistograms() const
const char * p
Definition: xmltok.h:285
const XML_Char * s
Definition: expat.h:262
TH1D * fHistos[kNumFluxType][kNumDets][kNumSels][kNumFlavors][kNumSigns]
Definition: GeniePCASyst.h:40
Float_t d
Definition: plot.C:236
double ana::GeniePCASyst::GetWeight ( const caf::SRProxy sr,
ESign  sign 
) const
protected

Definition at line 192 of file GeniePCASyst.cxx.

References om::cerr, fillBadChanDBTables::det, caf::Proxy< caf::SRHeader >::det, allTimeWatchdog::endl, energy, fHistName, fHistos, flux, make_syst_table_plots::h, caf::Proxy< caf::StandardRecord >::hdr, InitializeHistograms(), caf::Proxy< caf::SRSpill >::isRHC, kBkgFlav, kFD, kFHC, ana::kIsNC, ana::kIsNue, ana::kIsNumu, ana::kIsQE, kNC, kND, ana::kNue2020FD, ana::kNue2020ND, kNueSel, ana::kNumu2020FD, ana::kNumu2020ND, kNumuSel, kRHC, kSigNonQE, kSigQE, caf::Proxy< caf::StandardRecord >::mc, caf::Proxy< caf::SRTruthBranch >::nnu, caf::Proxy< caf::SRTruthBranch >::nu, canMan::sign(), and caf::Proxy< caf::StandardRecord >::spill.

Referenced by Shift().

194  {
196 
197  if(sr->mc.nnu == 0) return 1;
198 
199  int flux = sr->spill.isRHC ? kRHC : kFHC;
200 
201  int det = 0;
202  switch (sr->hdr.det){
203  case 1 : det = kND; break;
204  case 2 : det = kFD; break;
205  default: std::cerr <<"Wrong detector; ignore" << std::endl; return 1;
206  };
207 
208  // 2018 remnants
209  //const Cut NumuNDSel = kNumuCutND2018;
210  //const Cut NueNDSel = kNue2018NDCVNSsb;
211  //const Cut NumuFDSel = kNumuCutFD2018;
212  //const Cut NueFDSel = kNue2018FDAllSamples;
213 
214  // 2020 version
215  const Cut NumuNDSel = kNumu2020ND;
216  const Cut NueNDSel = kNue2020ND;
217  const Cut NumuFDSel = kNumu2020FD;
218  const Cut NueFDSel = kNue2020FD;
219 
220  int sel = -1;
221  if(det == kND){
222  if(NumuNDSel(sr)) sel = kNumuSel;
223  if(NueNDSel(sr)) sel = kNueSel;
224  }
225  if(det == kFD){
226  if(NumuFDSel(sr)) sel = kNumuSel;
227  if(NueFDSel(sr)) sel = kNueSel;
228  }
229  if(sel == -1) return 1;
230 
231  int flav = -1;
232  if(sel == kNumuSel){
233  if(kIsNumu(sr) && !kIsNC(sr) && kIsQE(sr)) flav = kSigQE;
234  if(kIsNumu(sr) && !kIsNC(sr) && !kIsQE(sr)) flav = kSigNonQE;
235  if(kIsNC(sr)) flav = kNC;
236  }
237  if(sel == kNueSel){
238  if(kIsNue(sr) && !kIsNC(sr) && kIsQE(sr)) flav = kSigQE;
239  if(kIsNue(sr) && !kIsNC(sr) && !kIsQE(sr)) flav = kSigNonQE;
240  if(kIsNumu(sr) && !kIsNC(sr)) flav = kBkgFlav;
241  if(kIsNC(sr)) flav = kNC;
242  }
243  if(flav == -1) return 1;
244 
245 
246  TH1D* h = fHistos[flux][det][sel][flav][sign];
247  if (h == 0){
248  std::cerr << fHistName+": Can't find desired histogram; ignore"<< std:: endl;
249  return 1;
250  }
251 
252  double energy = sr->mc.nu[0].E;
253 
254  if (energy > h->GetXaxis()->GetXmin() &&
255  energy < h->GetXaxis()->GetXmax() ){
256  return h->Interpolate(energy);
257  }
258 
259  if (energy > h->GetXaxis()->GetXmax())
260  return h->GetBinContent(h->GetNbinsX());
261 
262  return 1;
263  }
const Cut kIsQE
Definition: TruthCuts.cxx:104
caf::Proxy< caf::SRSpill > spill
Definition: SRProxy.h:2143
void InitializeHistograms() const
std::string fHistName
Definition: GeniePCASyst.h:38
caf::Proxy< caf::SRHeader > hdr
Definition: SRProxy.h:2137
const Cut kIsNumu([](const caf::SRProxy *sr){return(sr->mc.nnu==1 &&abs(sr->mc.nu[0].pdg)==14);})
Definition: TruthCuts.h:10
caf::Proxy< std::vector< caf::SRNeutrino > > nu
Definition: SRProxy.h:618
const Cut kNumu2020ND
Definition: NumuCuts2020.h:57
OStream cerr
Definition: OStream.cxx:7
caf::Proxy< short int > nnu
Definition: SRProxy.h:617
const Cut kNue2020ND
Definition: NueCuts2020.h:178
Loaders::FluxType flux
TH1D * fHistos[kNumFluxType][kNumDets][kNumSels][kNumFlavors][kNumSigns]
Definition: GeniePCASyst.h:40
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?
Definition: TruthCuts.h:8
double energy
Definition: plottest35.C:25
_Cut< caf::SRProxy > Cut
Representation of a cut (selection) to be applied to a caf::StandardRecord object.
Definition: Cut.h:96
const Cut kIsNue([](const caf::SRProxy *sr){return(sr->mc.nnu==1 &&abs(sr->mc.nu[0].pdg)==12);})
Definition: TruthCuts.h:9
const Cut kNumu2020FD
Definition: NumuCuts2020.h:59
const Cut kNue2020FD
Definition: NueCuts2020.h:65
caf::Proxy< caf::SRTruthBranch > mc
Definition: SRProxy.h:2138
caf::Proxy< bool > isRHC
Definition: SRProxy.h:1376
def sign(x)
Definition: canMan.py:197
caf::Proxy< caf::Det_t > det
Definition: SRProxy.h:232
void ana::GeniePCASyst::InitializeHistograms ( ) const
protected

Definition at line 135 of file GeniePCASyst.cxx.

References om::cerr, fillBadChanDBTables::det, allTimeWatchdog::endl, fFileName, fHistName, fHistos, fin, flux, compareCafs::histName, kBkgFlav, kFD, kFHC, kMinus, kNC, kND, kNueSel, kNumuSel, kPlus, kRHC, kSigNonQE, kSigQE, and canMan::sign().

Referenced by getHist(), and GetWeight().

136  {
137  if(fHistos[0][0][0][0][0]) return;
138 
139  TFile fin(fFileName.c_str(),"READ");
140  if(fin.IsZombie()){
141  std::cerr << "Warning: couldn't open " << fFileName << std::endl;
142  return;
143  }
144 
145  TIter iterHist(gDirectory->GetListOfKeys());
146  TKey* keyHist;
147  int foundHisto=0;
148 
149  while((keyHist = (TKey*)iterHist())) {
150  TString histName = keyHist->GetName(); // Get a histogram name
151  if(histName.Contains(fHistName)){ //is it the desired syst
152  foundHisto++;
153 
154  int flux=0;
155  if(histName.Contains("FHC")) flux = kFHC;
156  if(histName.Contains("RHC")) flux = kRHC;
157 
158  int det = 0;
159  if(histName.Contains("ND")) det = kND;
160  if(histName.Contains("FD")) det = kFD;
161 
162  int sel = 0;
163  if(histName.Contains("Numu")) sel = kNumuSel;
164  if(histName.Contains("Nue" )) sel = kNueSel;
165 
166  int flav = 0;
167  if(histName.Contains("SigQE" )) flav = kSigQE;
168  if(histName.Contains("SigNonQE")) flav = kSigNonQE;
169  if(histName.Contains("Other" )) flav = kBkgFlav;
170  if(histName.Contains("NC" )) flav = kNC;
171 
172  int sign = 0;
173  if(histName.Contains("minus")) sign = kMinus;
174  if(histName.Contains("plus" )) sign = kPlus;
175 
176  //store relevant histograms
177  fHistos[flux][det][sel][flav][sign] = (TH1D*) fin.Get(histName)->Clone();
178  // disassociate it from the file it came from
179  fHistos[flux][det][sel][flav][sign]->SetDirectory(0);
180  }
181  }
182 
183  if (foundHisto==0) {
184  std::cerr << "Genie PC systematic histogram " << fHistName
185  << " is not in file " << fFileName <<"; abort!" << std::endl;
186  abort();
187  }
188 
189  fin.Close();
190  }
TString fin
Definition: Style.C:24
std::string fHistName
Definition: GeniePCASyst.h:38
OStream cerr
Definition: OStream.cxx:7
Loaders::FluxType flux
TH1D * fHistos[kNumFluxType][kNumDets][kNumSels][kNumFlavors][kNumSigns]
Definition: GeniePCASyst.h:40
std::string fFileName
Definition: GeniePCASyst.h:37
def sign(x)
Definition: canMan.py:197
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
void ana::GeniePCASyst::Shift ( double  sigma,
caf::SRProxy sr,
double &  weight 
) const
overridevirtual

Perform the systematic shift.

Override this function if your systematic depends on non-SRNeutrino quantities. If it is SRNeutrino-only, implement the other function, and let this default forward to you when necessary.

Parameters
sigmaNumber of sigma to shift record by
srThe record to inspect and alter
weightScale this weight for reweighting systematics

Reimplemented from ana::ISyst.

Definition at line 265 of file GeniePCASyst.cxx.

References GetWeight(), kMinus, kPlus, std::max(), and canMan::sign().

267  {
268  if(sigma == 0) return;
269 
270  // Need different histogram for plus or minus sigma
271  ESign sign = kPlus;
272  if(sigma<0) sign = kMinus;
273 
274  weight *=1+(GetWeight(sr, sign)-1.)*TMath::Abs(sigma);
275  weight = std::max(0., weight);
276  }
T max(const caf::Proxy< T > &a, T b)
const Var weight
double GetWeight(const caf::SRProxy *sr, ESign sign) const
double sigma(TH1F *hist, double percentile)
def sign(x)
Definition: canMan.py:197
virtual const std::string& ana::ISyst::ShortName ( ) const
inlinefinalvirtualinherited
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 ana::BeamSyst, demo::DemoSyst1, ana::GenericSystComponentScale< T >, ana::GenericSystComponentScale< T >, ana::GenericSystComponentScale< T >, ana::ReinteractionSyst, 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::string ana::GeniePCASyst::fFileName
protected

Definition at line 37 of file GeniePCASyst.h.

Referenced by InitializeHistograms().

std::string ana::GeniePCASyst::fHistName
protected

Definition at line 38 of file GeniePCASyst.h.

Referenced by GetWeight(), and InitializeHistograms().

TH1D* ana::GeniePCASyst::fHistos[kNumFluxType][kNumDets][kNumSels][kNumFlavors][kNumSigns]
mutableprotected

Definition at line 40 of file GeniePCASyst.h.

Referenced by getHist(), GetWeight(), InitializeHistograms(), and ~GeniePCASyst().


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