PredictionCombinePeriods.cxx
Go to the documentation of this file.
2 
5 
6 #include "OscLib/IOscCalc.h"
7 
8 #include "TDirectory.h"
9 #include "TH2.h"
10 #include "TObjString.h"
11 #include "TVectorD.h"
12 
13 namespace ana
14 {
15  REGISTER_LOADFROM("PredictionCombinePeriods", IPrediction, PredictionCombinePeriods);
16 
17  //----------------------------------------------------------------------
18  PredictionCombinePeriods::PredictionCombinePeriods(const std::vector<std::pair<const IPrediction*, double>>& preds)
19  : fPreds(preds)
20  {
21  assert(!fPreds.empty());
22  }
23 
24  //----------------------------------------------------------------------
26  {
27  // It's pretty normal for the elements of fPred to have come from the stack
28  // for(auto it: fPreds) delete it.first;
29  }
30 
31  //----------------------------------------------------------------------
33  {
34  return PredictComponent(calc,
37  Sign::kBoth);
38  }
39 
40  //----------------------------------------------------------------------
42  const SystShifts& syst) const
43  {
44  return PredictComponentSyst(calc, syst,
47  Sign::kBoth);
48  }
49 
50  //----------------------------------------------------------------------
53  Flavors::Flavors_t flav,
55  Sign::Sign_t sign) const
56  {
58  flav, curr, sign);
59  }
60 
61  //----------------------------------------------------------------------
64  const SystShifts& syst,
65  Flavors::Flavors_t flav,
67  Sign::Sign_t sign) const
68  {
69  Eigen::ArrayXd asum;
70 
71  std::vector<std::string> labels;
72  std::vector<Binning> bins;
73 
74  double potsum = 0;
75  double livesum = 0;
76  for(auto it: fPreds){
77  const IPrediction* pred = it.first;
78  const double targetPOT = it.second;
79  const Spectrum s = pred->PredictComponentSyst(calc, syst,
80  flav, curr, sign);
81 
82  Eigen::ArrayXd h = s.GetEigen(targetPOT);
83  potsum += targetPOT;
84  // No. FD MC is going to have zero (unknown) livetime anyway. If someone
85  // needs this livetime accounting, they can come and think hard about it
86  // and then reinstate it.
87  // livesum += s.Livetime()*targetPOT/s.POT();
88 
89  if(asum.size() > 0){
90  asum += h;
91  }
92  else{
93  asum = h;
94  labels = s.GetLabels();
95  bins = s.GetBinnings();
96  }
97  }
98 
99  return Spectrum(std::move(asum),
100  HistAxis(labels, bins),
101  potsum, livesum);
102  }
103 
104  //----------------------------------------------------------------------
106  {
107  // TODO this is a total copy-paste job from the above.
108 
109  Eigen::MatrixXd asum;
110 
111  std::vector<std::string> labels;
112  std::vector<Binning> bins;
113 
114  double potsum = 0;
115  double livesum = 0;
116  for(auto it: fPreds){
117  const IPrediction* pred = it.first;
118  const double targetPOT = it.second;
119  const OscillatableSpectrum s = pred->ComponentCC(from, to);
120 
121  Eigen::MatrixXd h = s.GetEigen(targetPOT);
122  potsum += targetPOT;
123  // No. FD MC is going to have zero (unknown) livetime anyway. If someone
124  // needs this livetime accounting, they can come and think hard about it
125  // and then reinstate it.
126  // livesum += s.Livetime()*targetPOT/s.POT();
127 
128  if(asum.size() > 0){
129  asum += h;
130  }
131  else{
132  asum = h;
133  labels = s.GetLabels();
134  bins = s.GetBinnings();
135  }
136  }
137 
138  return OscillatableSpectrum(std::move(asum),
139  HistAxis(labels, bins),
140  potsum, livesum);
141  }
142 
143  //----------------------------------------------------------------------
144  //nc
146  {
149  }
151  {
154  }
156  {
159  }
160 
161  //end nc
162 
163  //----------------------------------------------------------------------
164  void PredictionCombinePeriods::SaveTo(TDirectory* dir, const std::string& name) const
165  {
166  TDirectory* tmp = gDirectory;
167 
168  dir = dir->mkdir(name.c_str()); // switch to subdir
169  dir->cd();
170 
171  TObjString("PredictionCombinePeriods").Write("type");
172 
173  int idx = 0;
174  for(auto it: fPreds)
175  it.first->SaveTo(dir, TString::Format("pred%d", idx++).Data());
176 
177  TVectorD pots(fPreds.size());
178  for(unsigned int i = 0; i < fPreds.size(); ++i)
179  pots[i] = fPreds[i].second;
180 
181  pots.Write("pots");
182 
183  dir->Write();
184  delete dir;
185 
186  tmp->cd();
187  }
188 
189  //----------------------------------------------------------------------
190  std::unique_ptr<PredictionCombinePeriods> PredictionCombinePeriods::
191  LoadFrom(TDirectory* dir, const std::string& name)
192  {
193  dir = dir->GetDirectory(name.c_str()); // switch to subdir
194  assert(dir);
195 
196  TObjString* tag = (TObjString*)dir->Get("type");
197  assert(tag);
198  assert(tag->GetString() == "PredictionCombinePeriods");
199  delete tag;
200 
201  TVectorD* pots = (TVectorD*)dir->Get("pots");
202  assert(pots);
203 
204  std::vector<std::pair<const IPrediction*, double>> preds;
205 
206  for(int i = 0; i < pots->GetNrows(); ++i){
207  const IPrediction* pred = ana::LoadFrom<IPrediction>(dir, TString::Format("pred%d", i).Data()).release();
208  assert(pred);
209  preds.emplace_back(pred, (*pots)[i]);
210  }
211 
212  delete dir;
213 
214  return std::unique_ptr<PredictionCombinePeriods>(new PredictionCombinePeriods(preds));
215  }
216 }
const std::vector< Binning > & GetBinnings() const
Definition: Spectrum.h:256
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
_HistAxis< Var > HistAxis
Definition: HistAxis.h:103
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
General interface to oscillation calculators.
Definition: StanTypedefs.h:22
virtual Spectrum PredictComponentSyst(osc::IOscCalc *calc, const SystShifts &syst, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const
Definition: IPrediction.cxx:79
const Eigen::ArrayXd & GetEigen() const
NB these don&#39;t have POT scaling. For expert high performance ops only!
Definition: Spectrum.h:188
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:40
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:18
osc::NoOscillations noosc
Definition: NueSystFuncs.h:44
std::vector< float > Spectrum
Definition: Constants.h:570
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
REGISTER_LOADFROM("BENDecomp", IDecomp, BENDecomp)
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:255
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