PredictionCombinePeriods.cxx
Go to the documentation of this file.
2 
5 
7 
8 #include "TDirectory.h"
9 #include "TH2.h"
10 #include "TObjString.h"
11 #include "TVectorD.h"
12 
13 namespace ana
14 {
15  //----------------------------------------------------------------------
16  PredictionCombinePeriods::PredictionCombinePeriods(const std::vector<std::pair<const IPrediction*, double>>& preds)
17  : fPreds(preds)
18  {
19  assert(!fPreds.empty());
20  }
21 
22  //----------------------------------------------------------------------
24  {
25  // It's pretty normal for the elements of fPred to have come from the stack
26  // for(auto it: fPreds) delete it.first;
27  }
28 
29  //----------------------------------------------------------------------
31  {
32  return PredictComponent(calc,
35  Sign::kBoth);
36  }
37 
38  //----------------------------------------------------------------------
40  const SystShifts& syst) const
41  {
42  return PredictComponentSyst(calc, syst,
45  Sign::kBoth);
46  }
47 
48  //----------------------------------------------------------------------
51  Flavors::Flavors_t flav,
53  Sign::Sign_t sign) const
54  {
56  flav, curr, sign);
57  }
58 
59  //----------------------------------------------------------------------
62  const SystShifts& syst,
63  Flavors::Flavors_t flav,
65  Sign::Sign_t sign) const
66  {
67  TH1D* hsum = 0;
68 
69  std::vector<std::string> labels;
70  std::vector<Binning> bins;
71 
72  double potsum = 0;
73  double livesum = 0;
74  for(auto it: fPreds){
75  const IPrediction* pred = it.first;
76  const double targetPOT = it.second;
77  const Spectrum s = pred->PredictComponentSyst(calc, syst,
78  flav, curr, sign);
79 
80  TH1D* h = s.ToTH1(targetPOT);
81  potsum += targetPOT;
82  // No. FD MC is going to have zero (unknown) livetime anyway. If someone
83  // needs this livetime accounting, they can come and think hard about it
84  // and then reinstate it.
85  // livesum += s.Livetime()*targetPOT/s.POT();
86 
87  if(hsum){
88  hsum->Add(h);
90  }
91  else{
92  hsum = h;
93  labels = s.GetLabels();
94  bins = s.GetBinnings();
95  }
96  }
97 
98  return Spectrum(std::unique_ptr<TH1D>(hsum), labels, bins, potsum, livesum);
99  }
100 
101  //----------------------------------------------------------------------
103  {
104  // TODO this is a total copy-paste job from the above.
105 
106  std::unique_ptr<TH2D> hsum(nullptr);
107 
108  std::vector<std::string> labels;
109  std::vector<Binning> bins;
110 
111  double potsum = 0;
112  double livesum = 0;
113  for(auto it: fPreds){
114  const IPrediction* pred = it.first;
115  const double targetPOT = it.second;
116  const OscillatableSpectrum s = pred->ComponentCC(from, to);
117 
118  TH2D* h = s.ToTH2(targetPOT);
119  potsum += targetPOT;
120  // No. FD MC is going to have zero (unknown) livetime anyway. If someone
121  // needs this livetime accounting, they can come and think hard about it
122  // and then reinstate it.
123  // livesum += s.Livetime()*targetPOT/s.POT();
124 
125  if(hsum){
126  hsum->Add(h);
128  }
129  else{
130  hsum.reset(h);
131  labels = s.GetLabels();
132  bins = s.GetBinnings();
133  }
134  }
135 
136  return OscillatableSpectrum(hsum.get(), labels, bins, potsum, livesum);
137  }
138 
139  //----------------------------------------------------------------------
140  //nc
142  {
145  }
147  {
150  }
152  {
155  }
156 
157  //end nc
158 
159  //----------------------------------------------------------------------
160  void PredictionCombinePeriods::SaveTo(TDirectory* dir) const
161  {
162  TDirectory* tmp = gDirectory;
163 
164  dir->cd();
165 
166  TObjString("PredictionCombinePeriods").Write("type");
167 
168  int idx = 0;
169  for(auto it: fPreds)
170  it.first->SaveTo(dir->mkdir(TString::Format("pred%d", idx++)));
171 
172  TVectorD pots(fPreds.size());
173  for(unsigned int i = 0; i < fPreds.size(); ++i)
174  pots[i] = fPreds[i].second;
175 
176  pots.Write("pots");
177 
178  tmp->cd();
179  }
180 
181  //----------------------------------------------------------------------
182  std::unique_ptr<PredictionCombinePeriods> PredictionCombinePeriods::
183  LoadFrom(TDirectory* dir)
184  {
185  TObjString* tag = (TObjString*)dir->Get("type");
186  assert(tag);
187  assert(tag->GetString() == "PredictionCombinePeriods");
188 
189  TVectorD* pots = (TVectorD*)dir->Get("pots");
190  assert(pots);
191 
192  std::vector<std::pair<const IPrediction*, double>> preds;
193 
194  for(int i = 0; i < pots->GetNrows(); ++i){
195  const IPrediction* pred = ana::LoadFrom<IPrediction>(dir->GetDirectory(TString::Format("pred%d", i))).release();
196  assert(pred);
197  preds.emplace_back(pred, (*pots)[i]);
198  }
199 
200  return std::unique_ptr<PredictionCombinePeriods>(new PredictionCombinePeriods(preds));
201  }
202 }
Pass neutrinos through unchanged.
OscillatableSpectrum ComponentCC(int from, int to) const override
osc::OscCalculatorDumb calc
Oscillation analysis framework, runs over CAF files outside of ART.
set< int >::iterator it
Antineutrinos-only.
Definition: IPrediction.h:50
General interface to oscillation calculators.
Definition: FwdDeclare.h:15
TH1D * ToTH1(double exposure, Color_t col=kBlack, Style_t style=kSolid, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Histogram made from this Spectrum, scaled to some exposure.
Definition: Spectrum.cxx:553
int pots
Simple record of shifts applied to systematic parameters.
Definition: SystShifts.h:20
PredictionCombinePeriods(const std::vector< std::pair< const IPrediction *, double >> &preds)
Spectrum ComponentNCTotal() const override
static std::unique_ptr< PredictionCombinePeriods > LoadFrom(TDirectory *dir)
Spectrum PredictComponent(osc::IOscCalculator *calc, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const override
Spectrum Predict(osc::IOscCalculator *calc) const override
Float_t tmp
Definition: plot.C:36
static SystShifts Nominal()
Definition: SystShifts.h:34
std::vector< std::pair< const IPrediction *, double > > fPreds
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:37
virtual Spectrum PredictComponentSyst(osc::IOscCalculator *calc, const SystShifts &syst, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const
const XML_Char * s
Definition: expat.h:262
Interactions of both types.
Definition: IPrediction.h:42
std::vector< Binning > GetBinnings() const
std::vector< std::string > GetLabels() const
static void Delete(TH1D *&h)
Definition: HistCache.cxx:215
Spectrum ComponentNCAnti() const override
osc::NoOscillations noosc
Definition: NueSystFuncs.h:44
Neutrinos-only.
Definition: IPrediction.h:49
def sign(x)
Definition: canMan.py:204
Spectrum PredictComponentSyst(osc::IOscCalculator *calc, const SystShifts &syst, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const override
const Binning bins
Definition: NumuCC_CPiBin.h:8
std::unique_ptr< IPrediction > LoadFrom< IPrediction >(TDirectory *dir)
Definition: IPrediction.cxx:38
TDirectory * dir
Definition: macro.C:5
virtual OscillatableSpectrum ComponentCC(int from, int to) const
Definition: IPrediction.h:103
string syst
Definition: plotSysts.py:176
Spectrum PredictSyst(osc::IOscCalculator *calc, const SystShifts &syst) const override
Neutral-current interactions.
Definition: IPrediction.h:40
Both neutrinos and antineutrinos.
Definition: IPrediction.h:52
Standard interface to all prediction techniques.
Definition: IPrediction.h:57
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
std::vector< std::string > GetLabels() const
Definition: Spectrum.h:298
All neutrinos, any flavor.
Definition: IPrediction.h:26
Spectrum with true energy information, allowing it to be oscillated
void SaveTo(TDirectory *dir) const override
std::vector< Binning > GetBinnings() const
Definition: Spectrum.h:299
h
Definition: demo3.py:41
TH2D * ToTH2(double pot) const