PredictionCombinePeriods.cxx
Go to the documentation of this file.
2 
4 
5 #include "OscLib/IOscCalc.h"
6 
7 #include "TDirectory.h"
8 #include "TH2.h"
9 #include "TObjString.h"
10 #include "TVectorD.h"
11 
12 namespace ana
13 {
14  //----------------------------------------------------------------------
15  PredictionCombinePeriods::PredictionCombinePeriods(const std::vector<std::pair<const IPrediction*, double>>& preds)
16  : fPreds(preds)
17  {
18  assert(!fPreds.empty());
19  }
20 
21  //----------------------------------------------------------------------
23  {
24  // It's pretty normal for the elements of fPred to have come from the stack
25  // for(auto it: fPreds) delete it.first;
26  }
27 
28  //----------------------------------------------------------------------
30  {
31  return PredictComponent(calc,
34  Sign::kBoth);
35  }
36 
37  //----------------------------------------------------------------------
39  const SystShifts& syst) const
40  {
41  return PredictComponentSyst(calc, syst,
44  Sign::kBoth);
45  }
46 
47  //----------------------------------------------------------------------
50  Flavors::Flavors_t flav,
52  Sign::Sign_t sign) const
53  {
55  flav, curr, sign);
56  }
57 
58  //----------------------------------------------------------------------
61  const SystShifts& syst,
62  Flavors::Flavors_t flav,
64  Sign::Sign_t sign) const
65  {
66  Eigen::ArrayXd asum;
67 
68  std::vector<std::string> labels;
69  std::vector<Binning> bins;
70 
71  double potsum = 0;
72  double livesum = 0;
73  for(auto it: fPreds){
74  const IPrediction* pred = it.first;
75  const double targetPOT = it.second;
76  const Spectrum s = pred->PredictComponentSyst(calc, syst,
77  flav, curr, sign);
78 
79  Eigen::ArrayXd h = s.GetEigen(targetPOT);
80  potsum += targetPOT;
81  // No. FD MC is going to have zero (unknown) livetime anyway. If someone
82  // needs this livetime accounting, they can come and think hard about it
83  // and then reinstate it.
84  // livesum += s.Livetime()*targetPOT/s.POT();
85 
86  if(asum.size() > 0){
87  asum += h;
88  }
89  else{
90  asum = h;
91  labels = s.GetLabels();
92  bins = s.GetBinnings();
93  }
94  }
95 
96  return Spectrum(std::move(asum),
97  HistAxis(labels, bins),
98  potsum, livesum);
99  }
100 
101  //----------------------------------------------------------------------
103  {
104  // TODO this is a total copy-paste job from the above.
105 
106  Eigen::MatrixXd asum;
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  Eigen::MatrixXd h = s.GetEigen(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(asum.size() > 0){
126  asum += h;
127  }
128  else{
129  asum = h;
130  labels = s.GetLabels();
131  bins = s.GetBinnings();
132  }
133  }
134 
135  return OscillatableSpectrum(std::move(asum),
136  HistAxis(labels, bins),
137  potsum, livesum);
138  }
139 
140  //----------------------------------------------------------------------
141  //nc
143  {
146  }
148  {
151  }
153  {
156  }
157 
158  //end nc
159 
160  //----------------------------------------------------------------------
161  void PredictionCombinePeriods::SaveTo(TDirectory* dir, const std::string& name) const
162  {
163  TDirectory* tmp = gDirectory;
164 
165  dir = dir->mkdir(name.c_str()); // switch to subdir
166  dir->cd();
167 
168  TObjString("PredictionCombinePeriods").Write("type");
169 
170  int idx = 0;
171  for(auto it: fPreds)
172  it.first->SaveTo(dir, TString::Format("pred%d", idx++).Data());
173 
174  TVectorD pots(fPreds.size());
175  for(unsigned int i = 0; i < fPreds.size(); ++i)
176  pots[i] = fPreds[i].second;
177 
178  pots.Write("pots");
179 
180  dir->Write();
181  delete dir;
182 
183  tmp->cd();
184  }
185 
186  //----------------------------------------------------------------------
187  std::unique_ptr<PredictionCombinePeriods> PredictionCombinePeriods::
188  LoadFrom(TDirectory* dir, const std::string& name)
189  {
190  dir = dir->GetDirectory(name.c_str()); // switch to subdir
191  assert(dir);
192 
193  TObjString* tag = (TObjString*)dir->Get("type");
194  assert(tag);
195  assert(tag->GetString() == "PredictionCombinePeriods");
196  delete tag;
197 
198  TVectorD* pots = (TVectorD*)dir->Get("pots");
199  assert(pots);
200 
201  std::vector<std::pair<const IPrediction*, double>> preds;
202 
203  for(int i = 0; i < pots->GetNrows(); ++i){
204  const IPrediction* pred = ana::LoadFrom<IPrediction>(dir, TString::Format("pred%d", i).Data()).release();
205  assert(pred);
206  preds.emplace_back(pred, (*pots)[i]);
207  }
208 
209  delete dir;
210 
211  return std::unique_ptr<PredictionCombinePeriods>(new PredictionCombinePeriods(preds));
212  }
213 }
const std::vector< Binning > & GetBinnings() const
Definition: Spectrum.h:268
const XML_Char * name
Definition: expat.h:151
Pass neutrinos through unchanged.
Definition: IOscCalc.h:40
OscillatableSpectrum ComponentCC(int from, int to) const override
Spectrum PredictSyst(osc::IOscCalc *calc, const SystShifts &syst) const override
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
set< int >::iterator it
Antineutrinos-only.
Definition: IPrediction.h:50
static std::unique_ptr< PredictionCombinePeriods > LoadFrom(TDirectory *dir, const std::string &name)
int pots
Simple record of shifts applied to systematic parameters.
Definition: SystShifts.h:20
virtual Spectrum PredictComponentSyst(osc::IOscCalc *calc, const SystShifts &syst, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const
const Eigen::ArrayXd & GetEigen() const
NB these don&#39;t have POT scaling. For expert high performance ops only!
Definition: Spectrum.h:200
PredictionCombinePeriods(const std::vector< std::pair< const IPrediction *, double >> &preds)
Spectrum ComponentNCTotal() const override
Float_t tmp
Definition: plot.C:36
static SystShifts Nominal()
Definition: SystShifts.h:34
osc::OscCalcDumb calc
std::vector< std::pair< const IPrediction *, double > > fPreds
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:33
const std::vector< std::string > & GetLabels() const
const XML_Char * s
Definition: expat.h:262
Spectrum PredictComponentSyst(osc::IOscCalc *calc, const SystShifts &syst, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const override
Interactions of both types.
Definition: IPrediction.h:42
virtual OscillatableSpectrum ComponentCC(int from, int to) const
const std::vector< Binning > & GetBinnings() const
Spectrum ComponentNCAnti() const override
std::unique_ptr< IPrediction > LoadFrom< IPrediction >(TDirectory *dir, const std::string &label)
Definition: IPrediction.cxx:28
osc::NoOscillations noosc
Definition: NueSystFuncs.h:44
std::vector< float > Spectrum
Definition: Constants.h:527
Spectrum PredictComponent(osc::IOscCalc *calc, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const override
Neutrinos-only.
Definition: IPrediction.h:49
const Binning bins
Definition: NumuCC_CPiBin.h:8
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
TDirectory * dir
Definition: macro.C:5
GenericHistAxis< Var > HistAxis
Definition: HistAxis.h:111
Neutral-current interactions.
Definition: IPrediction.h:40
assert(nhit_max >=nhit_nbins)
Eigen::MatrixXd GetEigen(double pot) const
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
const std::vector< std::string > & GetLabels() const
Definition: Spectrum.h:267
All neutrinos, any flavor.
Definition: IPrediction.h:26
Spectrum with true energy information, allowing it to be oscillated
void SaveTo(TDirectory *dir, const std::string &name) const override
def sign(x)
Definition: canMan.py:197
Spectrum Predict(osc::IOscCalc *calc) const override