Spectrum.cxx
Go to the documentation of this file.
1 #include "CAFAna/Core/Spectrum.h"
2 
3 #include "CAFAna/Core/Ratio.h"
4 
5 #include "TDirectory.h"
6 #include "TH2.h"
7 #include "TH3.h"
8 #include "THnSparse.h"
9 #include "TObjString.h"
10 #include "TRandom3.h"
11 
12 #include <cassert>
13 #include <iostream>
14 
15 namespace
16 {
17  namespace util{template<class T> T sqr(const T& x){return x*x;}}
18 }
19 
20 namespace ana
21 {
22  //----------------------------------------------------------------------
24  : fHist(Hist::Uninitialized()), fPOT(0), fLivetime(0), fAxis(axis)
25  {
26  const Binning bins1D = fAxis.GetBins1D();
27 
28  if(sparse == kSparse){
29  fHist = Hist::ZeroSparse(bins1D.NBins());
30  }
31  else{
32  fHist = Hist::Zero(bins1D.NBins());
33  }
34  }
35 
36  //----------------------------------------------------------------------
37  Spectrum::Spectrum(Eigen::ArrayXd&& h,
38  const LabelsAndBins& axis,
39  double pot, double livetime)
40  : fHist(Hist::Adopt(std::move(h))), fPOT(pot), fLivetime(livetime), fAxis(axis)
41  {
42  }
43 
44  //----------------------------------------------------------------------
46  const LabelsAndBins& axis,
47  double pot, double livetime)
48  : fHist(Hist::AdoptStan(std::move(h))), fPOT(pot), fLivetime(livetime), fAxis(axis)
49  {
50  }
51 
52  //----------------------------------------------------------------------
54  fHist(rhs.fHist),
55  fPOT(rhs.fPOT),
56  fLivetime(rhs.fLivetime),
57  fAxis(rhs.fAxis)
58  {
59  assert(rhs.fReferences.empty()); // Copying with pending loads is unexpected
60  }
61 
62  //----------------------------------------------------------------------
64  {
65  if(this == &rhs) return *this;
66 
67  fHist = rhs.fHist;
68  fPOT = rhs.fPOT;
69  fLivetime = rhs.fLivetime;
70  fAxis = rhs.fAxis;
71 
72  assert(fReferences.empty()); // Copying with pending loads is unexpected
73 
74  return *this;
75  }
76 
77  //----------------------------------------------------------------------
79  {
80  // Unregister self from anything that might still want to fill us
81  for(Spectrum** ref: fReferences) *ref = 0;
82  }
83 
84  //----------------------------------------------------------------------
85  TH1D* Spectrum::ToTH1(double exposure,
86  EExposureType expotype,
87  EBinType bintype) const
88  {
89  // Could have a file temporarily open
90  DontAddDirectory guard;
91 
92  TH1D* ret = fHist.ToTH1(fAxis.GetBins1D());
93 
94  ret->GetXaxis()->SetTitle(fAxis.GetLabel1D().c_str());
95  ret->GetYaxis()->SetTitle("Events");
96 
97  if(expotype == kPOT){
98  const double pot = exposure;
99  if(fPOT){
100  ret->Scale(pot/fPOT);
101  }
102  else{
103  // Allow zero POT if there are also zero events
104  if(ret->Integral() > 0){
105  std::cout << "Error: Spectrum with " << ret->Integral()
106  << " entries has zero POT, no way to scale to "
107  << exposure << " POT.";
108  if(fLivetime > 0){
109  std::cout << " Spectrum has " << fLivetime << " seconds livetime. "
110  << "Did you mean to pass kLivetime to ToTH1()?";
111  }
112  std::cout << std::endl;
113  abort();
114  }
115  }
116  }
117  if(expotype == kLivetime){
118  const double livetime = exposure;
119  if(fLivetime){
120  ret->Scale(livetime/fLivetime);
121  }
122  else{
123  // Allow zero exposure if there are also zero events
124  if(ret->Integral() > 0){
125  std::cout << "Error: Spectrum with " << ret->Integral()
126  << " entries has zero livetime, no way to scale to "
127  << livetime << " seconds.";
128  if(fPOT > 0){
129  std::cout << " Spectrum has " << fPOT << " POT. "
130  << "Did you mean to pass kPOT to ToTH1()?";
131  }
132  std::cout << std::endl;
133  abort();
134  }
135  }
136  }
137 
138  if(bintype == kBinDensity) ret->Scale(1, "width");
139 
140  // Allow GetMean() and friends to work even if this histogram never had any
141  // explicit Fill() calls made.
142  if(ret->GetEntries() == 0) ret->SetEntries(1);
143 
144  return ret;
145  }
146 
147  //----------------------------------------------------------------------
148  TH1D* Spectrum::ToTH1(double exposure, Color_t col, Style_t style,
149  EExposureType expotype,
150  EBinType bintype) const
151  {
152  // Could have a file temporarily open
153  DontAddDirectory guard;
154 
155  TH1D* ret = ToTH1(exposure, expotype, bintype);
156 
157  ret->SetLineColor(col);
158  ret->SetMarkerColor(col);
159  ret->SetLineStyle(style);
160 
161  return ret;
162  }
163 
164  //----------------------------------------------------------------------
165  TH2* Spectrum::ToTH2(double exposure, EExposureType expotype, EBinType bintype) const
166  {
167  if(NDimensions() != 2){
168  std::cout << "Error: This Spectrum does not appear to be 2D." << std::endl;
169  abort();
170  }
171 
172  TH2* ret = ana::ToTH2(*this, exposure, expotype, GetBinnings()[0], GetBinnings()[1]);
173 
174  ret->GetXaxis()->SetTitle(GetLabels()[0].c_str());
175  ret->GetYaxis()->SetTitle(GetLabels()[1].c_str());
176 
177  if(bintype == kBinDensity) ret->Scale(1, "width");
178 
179  // Allow GetMean() and friends to work even if this histogram never had any
180  // explicit Fill() calls made.
181  if(ret->GetEntries() == 0) ret->SetEntries(1);
182 
183  return ret;
184  }
185 
186  //----------------------------------------------------------------------
187  TH3* Spectrum::ToTH3(double exposure, EExposureType expotype, EBinType bintype) const
188  {
189  if(NDimensions() != 3){
190  std::cout << "Error: This Spectrum does not appear to be 3D." << std::endl;
191  abort();
192  }
193 
194  TH3* ret = ana::ToTH3(*this, exposure, expotype,
195  GetBinnings()[0], GetBinnings()[1], GetBinnings()[2]);
196 
197  ret->GetXaxis()->SetTitle(GetLabels()[0].c_str());
198  ret->GetYaxis()->SetTitle(GetLabels()[1].c_str());
199  ret->GetZaxis()->SetTitle(GetLabels()[2].c_str());
200 
201  if(bintype == kBinDensity) ret->Scale(1, "width");
202 
203  // Allow GetMean() and friends to work even if this histogram never had any
204  // explicit Fill() calls made.
205  if(ret->GetEntries() == 0) ret->SetEntries(1);
206 
207  return ret;
208  }
209 
210  //----------------------------------------------------------------------
211  TH1* Spectrum::ToTHX(double exposure, EExposureType expotype, EBinType bintype) const
212  {
213  if(NDimensions() == 2) return ToTH2(exposure, expotype, bintype);
214  if(NDimensions() == 3) return ToTH3(exposure, expotype, bintype);
215  return ToTH1(exposure, expotype, bintype);
216  }
217 
218  //----------------------------------------------------------------------
219  Eigen::ArrayXd Spectrum::GetEigen(double exposure, EExposureType expotype) const
220  {
221  if(expotype == kPOT)
222  return (exposure/fPOT) * fHist.GetEigen();
223  else
224  return (exposure/fLivetime) * fHist.GetEigen();
225  }
226 
227  //----------------------------------------------------------------------
228  Eigen::ArrayXstan Spectrum::GetEigenStan(double exposure, EExposureType expotype) const
229  {
230  if(expotype == kPOT)
231  return (exposure/fPOT) * fHist.GetEigenStan();
232  else
233  return (exposure/fLivetime) * fHist.GetEigenStan();
234  }
235 
236  //----------------------------------------------------------------------
237  void Spectrum::Scale(double c)
238  {
239  fHist.Scale(c);
240  }
241 
242  //----------------------------------------------------------------------
244  {
245  fHist.Scale(c);
246  }
247 
248  //----------------------------------------------------------------------
249  double Spectrum::Integral(double exposure, double* err,
250  EExposureType expotype) const
251  {
252  const double ratio = (expotype == kPOT) ? exposure/fPOT : exposure/fLivetime;
253 
254  if(err){
255  *err = 0;
256 
257  for(int i = 0; i < fHist.GetNbinsX()+2; ++i){
258  *err += util::sqr(fHist.GetBinError(i));
259  }
260  *err = sqrt(*err) * ratio;
261  }
262 
263  return fHist.Integral() * ratio;
264  }
265 
266  //----------------------------------------------------------------------
267  double Spectrum::Mean() const
268  {
269  const Binning bins = fAxis.GetBins1D();
270 
271  if(fHist.GetBinContent(0) != 0){
272  std::cout << "Spectrum::Mean(): Warning ignoring underflow bin content " << fHist.GetBinContent(0) << std::endl;
273  }
274 
275  if(fHist.GetBinContent(bins.NBins()+1) != 0){
276  std::cout << "Spectrum::Mean(): Warning ignoring overflow bin content " << fHist.GetBinContent(bins.NBins()+1) << std::endl;
277  }
278 
279  double mean = 0;
280  double W = 0;
281  for(int i = 1; i <= bins.NBins(); ++i){
282  const double w = fHist.GetBinContent(i);
283  W += w;
284  const double x0 = bins.Edges()[i-1];
285  const double x1 = bins.Edges()[i];
286  mean += w * (x0+x1)/2;
287  }
288 
289  return mean/W;
290  }
291 
292  //----------------------------------------------------------------------
293  void Spectrum::Fill(double x, double w)
294  {
295  fHist.Fill(fAxis.GetBins1D(), x, w);
296  // TODO Pull binning out of Hist entirely and just update an index?
297  }
298 
299  //----------------------------------------------------------------------
300  Spectrum Spectrum::MockData(double pot, int seed) const
301  {
302  Spectrum ret = FakeData(pot);
303 
304  TRandom3 rnd(seed); // zero seeds randomly
305 
306  for(int i = 0; i < ret.fHist.GetNbinsX()+2; ++i){
307  ret.fHist.SetBinContent(i, rnd.Poisson(ret.fHist.GetBinContent(i)));
308  }
309 
310  // Drop old errors, which are based on the MC statistics, and create new
311  // ones that are based on the prediction for the data
312  ret.fHist.ResetErrors();
313 
314  return ret;
315  }
316 
317  //----------------------------------------------------------------------
319  {
320  Spectrum ret = *this;
321  if(fPOT > 0) ret.fHist.Scale(pot/fPOT);
322  ret.fPOT = pot;
323 
324  // Drop old errors, which are based on the MC statistics, and create new
325  // ones that are based on the prediction for the data
326  ret.fHist.ResetErrors();
327 
328  return ret;
329  }
330 
331  //----------------------------------------------------------------------
333  {
334  Spectrum ret = *this;
335 
336  if(fPOT > 0) ret.fHist.Scale(pot/fPOT);
337  ret.fPOT = pot;
338 
339  // overwrite livetime for fake data
340  // PlusEqualsHelper will take over the rest
341  ret.fLivetime = livetime;
342 
343  ret.fHist.ResetErrors();
344 
345  return ret;
346  }
347 
348  //----------------------------------------------------------------------
350  {
351  return AsimovData(pot);
352  }
353 
354  //----------------------------------------------------------------------
355  Spectrum Spectrum::FakeData(double pot, double livetime) const
356  {
357  return AsimovData(pot, livetime);
358  }
359 
360  //----------------------------------------------------------------------
362  {
363  fHist.Reset();
364  }
365 
366  //----------------------------------------------------------------------
368  {
369  fReferences.erase(ref);
370  }
371 
372  //----------------------------------------------------------------------
374  {
375  fReferences.insert(ref);
376  }
377 
378  //----------------------------------------------------------------------
380  {
381  // In this case it would be OK to have no POT/livetime
382  if(rhs.fHist.Initialized() && rhs.fHist.Integral() == 0) return *this;
383 
384  if((!fPOT && !fLivetime) || (!rhs.fPOT && !rhs.fLivetime)){
385  std::cout << "Error: can't sum Spectrum with no POT or livetime: "
386  << fPOT << " " << rhs.fPOT << " " << fLivetime << " " << rhs.fLivetime
387  << std::endl;
388  abort();
389  }
390 
391  if(!fLivetime && !rhs.fPOT){
392  std::cout << "Error: can't sum Spectrum with POT ("
393  << fPOT << ") but no livetime and Spectrum with livetime ("
394  << rhs.fLivetime << " sec) but no POT." << std::endl;
395  abort();
396  }
397 
398  if(!fPOT && !rhs.fLivetime){
399  std::cout << "Error: can't sum Spectrum with livetime ("
400  << fLivetime << " sec) but no POT and Spectrum with POT ("
401  << rhs.fPOT << ") but no livetime." << std::endl;
402  abort();
403  }
404 
405  // And now there are still a bunch of good cases to consider
406 
407  if(fPOT && rhs.fPOT){
408  // Scale by POT when possible
409  fHist.Add(rhs.fHist, sign*fPOT/rhs.fPOT);
410 
411  if(fLivetime && rhs.fLivetime){
412  // If POT/livetime ratios match, keep regular lifetime, otherwise zero
413  // it out.
414  if(AlmostEqual(fLivetime*rhs.fPOT, rhs.fLivetime*fPOT))
415  fLivetime = 0;
416  }
417  if(!fLivetime && rhs.fLivetime){
418  // If the RHS has a livetime and we don't, copy it in (suitably scaled)
419  fLivetime = rhs.fLivetime * fPOT/rhs.fPOT;
420  }
421  // Otherwise, keep our own livetime (if any)
422 
423  return *this;
424  }
425 
426  if(fLivetime && rhs.fLivetime){
427  // Scale by livetime, the only thing in common
428  fHist.Add(rhs.fHist, sign*fLivetime/rhs.fLivetime);
429 
430  if(!fPOT && rhs.fPOT){
431  // If the RHS has a POT and we don't, copy it in (suitably scaled)
432  fPOT = rhs.fPOT * fLivetime/rhs.fLivetime;
433  }
434  // Otherwise, keep our own POT (if any)
435 
436  return *this;
437  }
438 
439  // That should have been all the cases. I definitely want to know what
440  // happened if it wasn't.
441  std::cout << "Spectrum::operator+=(). How did we get here? "
442  << fPOT << " " << fLivetime << " "
443  << rhs.fPOT << " " << rhs.fLivetime << std::endl;
444  abort();
445  }
446 
447  //----------------------------------------------------------------------
449  {
450  return PlusEqualsHelper(rhs, +1);
451  }
452 
453  //----------------------------------------------------------------------
455  {
456  Spectrum ret = *this;
457  ret += rhs;
458  return ret;
459  }
460 
461  //----------------------------------------------------------------------
463  {
464  return PlusEqualsHelper(rhs, -1);
465  }
466 
467  //----------------------------------------------------------------------
469  {
470  Spectrum ret = *this;
471  ret -= rhs;
472  return ret;
473  }
474 
475  //----------------------------------------------------------------------
477  {
478  fHist.Multiply(rhs.fHist);
479  return *this;
480  }
481 
482  //----------------------------------------------------------------------
484  {
485  Spectrum ret = *this;
486  ret *= rhs;
487  return ret;
488  }
489 
490  //----------------------------------------------------------------------
492  {
493  fHist.Divide(rhs.fHist);
494  return *this;
495  }
496 
497  //----------------------------------------------------------------------
499  {
500  Spectrum ret = *this;
501  ret /= rhs;
502  return ret;
503  }
504 
505  //----------------------------------------------------------------------
506  void Spectrum::SaveTo(TDirectory* dir, const std::string& name) const
507  {
508  TDirectory* tmp = gDirectory;
509 
510  dir = dir->mkdir(name.c_str()); // switch to subdir
511  dir->cd();
512 
513  TObjString("Spectrum").Write("type");
514 
516  TH1D hPot("", "", 1, 0, 1);
517  hPot.Fill(.5, fPOT);
518  hPot.Write("pot");
519  TH1D hLivetime("", "", 1, 0, 1);
520  hLivetime.Fill(.5, fLivetime);
521  hLivetime.Write("livetime");
522 
523  for(unsigned int i = 0; i < NDimensions(); ++i){
524  TObjString(GetLabels()[i].c_str()).Write(TString::Format("label%d", i).Data());
525  GetBinnings()[i].SaveTo(dir, TString::Format("bins%d", i).Data());
526  }
527 
528  dir->Write();
529  delete dir;
530 
531  tmp->cd();
532  }
533 
534  //----------------------------------------------------------------------
535  std::unique_ptr<Spectrum> Spectrum::LoadFrom(TDirectory* dir, const std::string& name)
536  {
537  dir = dir->GetDirectory(name.c_str()); // switch to subdir
538  assert(dir);
539 
540  DontAddDirectory guard;
541 
542  TObjString* tag = (TObjString*)dir->Get("type");
543  assert(tag);
544  assert(tag->GetString() == "Spectrum");
545  delete tag;
546 
547  TH1* hPot = (TH1*)dir->Get("pot");
548  assert(hPot);
549  TH1* hLivetime = (TH1*)dir->Get("livetime");
550  assert(hLivetime);
551 
552  std::vector<std::string> labels;
553  std::vector<Binning> bins;
554  for(int i = 0; ; ++i){
555  const std::string subname = TString::Format("bins%d", i).Data();
556  TDirectory* subdir = dir->GetDirectory(subname.c_str());
557  if(!subdir) break;
558  delete subdir;
559  bins.push_back(*Binning::LoadFrom(dir, subname));
560  TObjString* label = (TObjString*)dir->Get(TString::Format("label%d", i));
561  labels.push_back(label ? label->GetString().Data() : "");
562  delete label;
563  }
564 
565  std::unique_ptr<Spectrum> ret = std::make_unique<Spectrum>(Eigen::ArrayXd(), LabelsAndBins(labels, bins), hPot->GetBinContent(1), hLivetime->GetBinContent(1));
566  ret->fHist = Hist::FromDirectory(dir);
567 
568  delete hPot;
569  delete hLivetime;
570 
571  delete dir;
572 
573  return ret;
574  }
575 }
static std::unique_ptr< Binning > LoadFrom(TDirectory *dir, const std::string &name)
Definition: Binning.cxx:230
const std::vector< Binning > & GetBinnings() const
Definition: Spectrum.h:264
const XML_Char * name
Definition: expat.h:151
void SetBinContent(int i, double x)
Definition: Hist.cxx:363
void Write(const Binning &bins) const
Definition: Hist.cxx:583
Filter events based on their run/event numbers.
const Eigen::ArrayXd & GetEigen() const
Definition: Hist.h:53
TH2 * ToTH2(double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Spectrum must be 2D to obtain TH2.
Definition: Spectrum.cxx:165
Spectrum operator*(const Ratio &rhs) const
Definition: Spectrum.cxx:483
void Add(const Hist &rhs, double scale=1)
Definition: Hist.cxx:455
Represent the binning of a Spectrum&#39;s x-axis.
Definition: Binning.h:16
void AddLoader(Spectrum **)
Definition: Spectrum.cxx:373
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
LabelsAndBins fAxis
Definition: Spectrum.h:291
double Integral() const
Definition: Hist.cxx:252
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:148
const SpillVar fPOT([](const caf::SRSpillProxy *spill){return spill->spillpot;})
subdir
Definition: cvnie.py:7
void ResetErrors()
Definition: Hist.cxx:342
Spectrum & PlusEqualsHelper(const Spectrum &rhs, int sign)
Helper for operator+= and operator-=.
Definition: Spectrum.cxx:379
Float_t x1[n_points_granero]
Definition: compare.C:5
T sqrt(T number)
Definition: d0nt_math.hpp:156
Spectrum & operator+=(const Spectrum &rhs)
Definition: Spectrum.cxx:448
void RemoveLoader(Spectrum **)
Definition: Spectrum.cxx:367
Definition: Hist.h:29
std::set< Spectrum ** > fReferences
Things that point at this Spectrum. Maintained by SpectrumLoader.
Definition: Spectrum.h:289
static Hist FromDirectory(TDirectory *dir)
Definition: Hist.cxx:158
const Eigen::ArrayXd & GetEigen() const
NB these don&#39;t have POT scaling. For expert high performance ops only!
Definition: Spectrum.h:189
int GetNbinsX() const
Definition: Hist.cxx:228
Divide bin contents by bin widths.
Definition: UtilsExt.h:31
TH1D * ToTH1(const Binning &bins) const
Definition: Hist.cxx:207
void Clear()
Definition: Spectrum.cxx:361
double fPOT
Definition: Spectrum.h:285
Float_t tmp
Definition: plot.C:36
void Fill(double x, double w=1)
Definition: Spectrum.cxx:293
def ratio(spec1, spec2, nbins, pot, binrange=(0, 1))
double Integral(double exposure, double *err=0, EExposureType expotype=kPOT) const
Return total number of events scaled to pot.
Definition: Spectrum.cxx:249
Spectrum & operator-=(const Spectrum &rhs)
Definition: Spectrum.cxx:462
T sqr(T x)
More efficient square function than pow(x,2)
Definition: MathUtil.h:23
void Scale(double c)
Multiply this spectrum by a constant c.
Definition: Spectrum.cxx:237
Spectrum MockData(double pot, int seed=0) const
Mock data is FakeData with Poisson fluctuations applied.
Definition: Spectrum.cxx:300
TH1 * ToTHX(double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Definition: Spectrum.cxx:211
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
const char * label
static Hist ZeroSparse(int nbins)
Definition: Hist.cxx:51
EBinType
Definition: UtilsExt.h:28
bool Initialized() const
Definition: Hist.h:39
static std::unique_ptr< Spectrum > LoadFrom(TDirectory *dir, const std::string &name)
Definition: Spectrum.cxx:535
const std::string & GetLabel1D() const
unsigned int seed
Definition: runWimpSim.h:102
bool AlmostEqual(double a, double b, double eps)
Definition: UtilsExt.cxx:40
Spectrum & operator/=(const Ratio &rhs)
Definition: Spectrum.cxx:491
Spectrum FakeData(double pot) const
Synonymous with AsimovData(). Retained for compatibility.
Definition: Spectrum.cxx:349
Spectrum operator/(const Ratio &rhs) const
Definition: Spectrum.cxx:498
const Binning & GetBins1D() const
Appropriate binning and labelling for that 1D Var.
Int_t col[ntarg]
Definition: Style.C:29
EExposureType
For use as an argument to Spectrum::ToTH1.
Definition: UtilsExt.h:35
Spectrum AsimovData(double pot) const
Asimov data is a MC spectrum scaled to the POT expected in the data.
Definition: Spectrum.cxx:318
#define pot
Spectrum & operator=(const Spectrum &rhs)
Definition: Spectrum.cxx:63
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
double fLivetime
Definition: Spectrum.h:286
void Scale(double s)
Definition: Hist.cxx:296
void Reset()
Definition: Hist.cxx:382
void SaveTo(TDirectory *dir, const std::string &name) const
Definition: Spectrum.cxx:506
unsigned int NDimensions() const
Definition: Spectrum.h:262
std::vector< float > Spectrum
Definition: Constants.h:728
Represent the ratio between two spectra.
Definition: Ratio.h:8
const std::vector< double > & Edges() const
Definition: Binning.h:34
OStream cout
Definition: OStream.cxx:6
void Fill(const Binning &bins, double x, double w)
Definition: Hist.cxx:265
const Binning bins
Definition: NumuCC_CPiBin.h:8
Spectrum operator+(const Spectrum &rhs) const
Definition: Spectrum.cxx:454
TH3 * ToTH3(double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Spectrum must be 3D to obtain TH3.
Definition: Spectrum.cxx:187
void Divide(const Hist &rhs)
Definition: Hist.cxx:531
Eigen::Array< stan::math::var, Eigen::Dynamic, 1 > ArrayXstan
Definition: StanUtils.h:7
TDirectory * dir
Definition: macro.C:5
Hist fHist
Definition: Ratio.h:56
double livetime
Definition: saveFDMCHists.C:21
int NBins() const
Definition: Binning.h:29
double GetBinContent(int i) const
Definition: Hist.cxx:349
static Hist Zero(int nbins)
Definition: Hist.cxx:41
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
double T
Definition: Xdiff_gwt.C:5
virtual ~Spectrum()
Definition: Spectrum.cxx:78
const std::vector< std::string > & GetLabels() const
Definition: Spectrum.h:263
Prevent histograms being added to the current directory.
Definition: UtilsExt.h:46
double GetBinError(int i) const
Definition: Hist.cxx:241
Float_t w
Definition: plot.C:20
#define W(x)
Spectrum operator-(const Spectrum &rhs) const
Definition: Spectrum.cxx:468
TH3 * ToTH3(const Spectrum &s, double exposure, ana::EExposureType expotype, const Binning &binsx, const Binning &binsy, const Binning &binsz, ana::EBinType bintype)
Same as ToTH2, but with 3 dimensions.
Definition: UtilsExt.cxx:162
TH2 * ToTH2(const Spectrum &s, double exposure, ana::EExposureType expotype, const Binning &binsx, const Binning &binsy, ana::EBinType bintype)
For use with Var2D.
Definition: UtilsExt.cxx:115
def sign(x)
Definition: canMan.py:197
Spectrum & operator*=(const Ratio &rhs)
Definition: Spectrum.cxx:476
const Eigen::ArrayXstan & GetEigenStan() const
Definition: Spectrum.h:190
void Multiply(const Hist &rhs)
Definition: Hist.cxx:479
const Eigen::ArrayXstan & GetEigenStan() const
Definition: Hist.h:54
double Mean() const
Return mean of 1D histogram.
Definition: Spectrum.cxx:267
enum BeamMode string