Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Friends | List of all members
ana::Spectrum Class Reference

Representation of a spectrum in any variable, with associated POT. More...

#include "/cvmfs/nova.opensciencegrid.org/externals/cafanacore/v01.11/src/CAFAna/Core/Spectrum.h"

Inheritance diagram for ana::Spectrum:
ana::_preview::NewSpectrum

Public Types

enum  ESparse { kDense, kSparse }
 

Public Member Functions

template<class T >
 Spectrum (SpectrumLoaderBase &loader, const _HistAxis< _Var< T >> &axis, const _Cut< T > &cut, const SystShifts &shift=kNoShift, const _Var< T > &wei=Unweighted< T >(), ESparse sparse=kDense)
 One constructor to rule them all. More...
 
template<class T >
 Spectrum (const std::string &label, const Binning &bins, SpectrumLoaderBase &loader, const _Var< T > &var, const _Cut< T > &cut, const SystShifts &shift=kNoShift, const _Var< T > &wei=Unweighted< T >(), ESparse sparse=kDense)
 
template<class T >
 Spectrum (SpectrumLoaderBase &loader, const _HistAxis< _MultiVar< T >> &axis, const _Cut< T > &cut, const SystShifts &shift=kNoShift, const _Var< T > &wei=Unweighted< T >())
 The only MultiVar variant available. More...
 
 Spectrum (Eigen::ArrayXd &&h, const LabelsAndBins &axis, double pot, double livetime)
 Makes a spectrum from an eigen array. More...
 
 Spectrum (Eigen::ArrayXstan &&h, const LabelsAndBins &axis, double pot, double livetime)
 Makes a spectrum from an eigen array of stan vars. More...
 
template<class T >
 Spectrum (SpectrumLoaderBase &loader, const _HistAxis< _Var< T >> &xAxis, const _HistAxis< _Var< T >> &yAxis, const _Cut< T > &cut, const SystShifts &shift=kNoShift, const _Var< T > &wei=Unweighted< T >(), ESparse sparse=kDense)
 2D Spectrum taking 2 HistAxis More...
 
template<class T >
 Spectrum (const std::string &xLabel, const std::string &yLabel, SpectrumLoaderBase &loader, const Binning &binsx, const _Var< T > &varx, const Binning &binsy, const _Var< T > &vary, const _Cut< T > &cut, const SystShifts &shift=kNoShift, const _Var< T > &wei=Unweighted< T >(), ESparse sparse=kDense)
 2D Spectrum of two Vars More...
 
template<class T >
 Spectrum (SpectrumLoaderBase &loader, const _HistAxis< _Var< T >> &xAxis, const _HistAxis< _Var< T >> &yAxis, const _HistAxis< _Var< T >> &zAxis, const _Cut< T > &cut, const SystShifts &shift=kNoShift, const _Var< T > &wei=Unweighted< T >(), ESparse sparse=kDense)
 3D Spectrum taking 3 HistAxis More...
 
template<class T >
 Spectrum (const std::string &xLabel, const std::string &yLabel, const std::string &zLabel, SpectrumLoaderBase &loader, const Binning &binsx, const _Var< T > &varx, const Binning &binsy, const _Var< T > &vary, const Binning &binsz, const _Var< T > &varz, const _Cut< T > &cut, const SystShifts &shift=kNoShift, const _Var< T > &wei=Unweighted< T >(), ESparse sparse=kDense)
 3D Spectrum of three Vars More...
 
 Spectrum (Hist &&hist, const LabelsAndBins &axis, double pot, double livetime)
 Expert constructor for ReweightableSpectrum et al. More...
 
virtual ~Spectrum ()
 
 Spectrum (const Spectrum &rhs)
 
 Spectrum (Spectrum &&rhs)
 
Spectrumoperator= (const Spectrum &rhs)
 
Spectrumoperator= (Spectrum &&rhs)
 
void Fill (double x, double w=1)
 
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. More...
 
TH1D * ToTH1 (double exposure, EExposureType expotype, EBinType bintype=kBinContent) const
 Histogram made from this Spectrum, scaled to some exposure. More...
 
TH2 * ToTH2 (double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
 Spectrum must be 2D to obtain TH2. More...
 
TH3 * ToTH3 (double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
 Spectrum must be 3D to obtain TH3. More...
 
TH1 * ToTHX (double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
 
bool HasStan () const
 
const Eigen::ArrayXd & GetEigen () const
 NB these don't have POT scaling. For expert high performance ops only! More...
 
const Eigen::ArrayXstanGetEigenStan () const
 
Eigen::ArrayXd GetEigen (double exposure, EExposureType expotype=kPOT) const
 
Eigen::ArrayXstan GetEigenStan (double exposure, EExposureType expotype=kPOT) const
 
double Integral (double exposure, double *err=0, EExposureType expotype=kPOT) const
 Return total number of events scaled to pot. More...
 
double Mean () const
 Return mean of 1D histogram. More...
 
Spectrum MockData (double pot, int seed=0) const
 Mock data is FakeData with Poisson fluctuations applied. More...
 
Spectrum AsimovData (double pot) const
 Asimov data is a MC spectrum scaled to the POT expected in the data. More...
 
Spectrum AsimovData (double pot, double livetime) const
 
Spectrum FakeData (double pot) const
 Synonymous with AsimovData(). Retained for compatibility. More...
 
Spectrum FakeData (double pot, double livetime) const
 Synonymous with AsimovData(). Retained for compatibility. More...
 
double POT () const
 
double Livetime () const
 Seconds. For informational purposes only. No calculations use this. More...
 
void OverridePOT (double newpot)
 DO NOT USE UNLESS YOU ARE 110% CERTAIN THERE ISN'T A BETTER WAY! More...
 
void OverrideLivetime (double newlive)
 DO NOT USE UNLESS YOU ARE 110% CERTAIN THERE ISN'T A BETTER WAY! More...
 
void Clear ()
 
void Scale (double c)
 Multiply this spectrum by a constant c. More...
 
void Scale (const stan::math::var &v)
 
Spectrumoperator+= (const Spectrum &rhs)
 
Spectrum operator+ (const Spectrum &rhs) const
 
Spectrumoperator-= (const Spectrum &rhs)
 
Spectrum operator- (const Spectrum &rhs) const
 
Spectrumoperator*= (const Ratio &rhs)
 
Spectrum operator* (const Ratio &rhs) const
 
Spectrumoperator/= (const Ratio &rhs)
 
Spectrum operator/ (const Ratio &rhs) const
 
void SaveTo (TDirectory *dir, const std::string &name) const
 
unsigned int NDimensions () const
 
const std::vector< std::string > & GetLabels () const
 
const std::vector< Binning > & GetBinnings () const
 

Static Public Member Functions

static Spectrum Uninitialized ()
 
static std::unique_ptr< SpectrumLoadFrom (TDirectory *dir, const std::string &name)
 

Protected Member Functions

 Spectrum ()
 Constructor for Uninitialized() More...
 
 Spectrum (const LabelsAndBins &axis, ESparse sparse=kDense)
 Helper for constructors. More...
 
void RemoveLoader (Spectrum **)
 
void AddLoader (Spectrum **)
 
SpectrumPlusEqualsHelper (const Spectrum &rhs, int sign)
 Helper for operator+= and operator-=. More...
 

Protected Attributes

Hist fHist
 
double fPOT
 
double fLivetime
 
std::set< Spectrum ** > fReferences
 Things that point at this Spectrum. Maintained by SpectrumLoader. More...
 
LabelsAndBins fAxis
 

Friends

class SpectrumLoaderBase
 
class SpectrumSink
 
class SpectrumSinkBase< Spectrum >
 
class Ratio
 

Detailed Description

Representation of a spectrum in any variable, with associated POT.

Definition at line 40 of file Spectrum.h.

Member Enumeration Documentation

Enumerator
kDense 
kSparse 

Definition at line 48 of file Spectrum.h.

Constructor & Destructor Documentation

template<class T >
ana::Spectrum::Spectrum ( SpectrumLoaderBase loader,
const _HistAxis< _Var< T >> &  axis,
const _Cut< T > &  cut,
const SystShifts shift = kNoShift,
const _Var< T > &  wei = UnweightedT >(),
ESparse  sparse = kDense 
)

One constructor to rule them all.

template<class T >
ana::Spectrum::Spectrum ( const std::string label,
const Binning bins,
SpectrumLoaderBase loader,
const _Var< T > &  var,
const _Cut< T > &  cut,
const SystShifts shift = kNoShift,
const _Var< T > &  wei = UnweightedT >(),
ESparse  sparse = kDense 
)
template<class T >
ana::Spectrum::Spectrum ( SpectrumLoaderBase loader,
const _HistAxis< _MultiVar< T >> &  axis,
const _Cut< T > &  cut,
const SystShifts shift = kNoShift,
const _Var< T > &  wei = UnweightedT >() 
)

The only MultiVar variant available.

ana::Spectrum::Spectrum ( Eigen::ArrayXd &&  h,
const LabelsAndBins axis,
double  pot,
double  livetime 
)

Makes a spectrum from an eigen array.

Definition at line 37 of file Spectrum.cxx.

40  : fHist(Hist::Adopt(std::move(h))), fPOT(pot), fLivetime(livetime), fAxis(axis)
41  {
42  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
double fPOT
Definition: Spectrum.h:285
#define pot
double fLivetime
Definition: Spectrum.h:286
double livetime
Definition: saveFDMCHists.C:21
static Hist Adopt(Eigen::ArrayXd &&v)
Definition: Hist.cxx:149
ana::Spectrum::Spectrum ( Eigen::ArrayXstan &&  h,
const LabelsAndBins axis,
double  pot,
double  livetime 
)

Makes a spectrum from an eigen array of stan vars.

Definition at line 45 of file Spectrum.cxx.

49  {
50  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
double fPOT
Definition: Spectrum.h:285
#define pot
double fLivetime
Definition: Spectrum.h:286
static Hist AdoptStan(Eigen::ArrayXstan &&v)
Definition: Hist.cxx:140
double livetime
Definition: saveFDMCHists.C:21
template<class T >
ana::Spectrum::Spectrum ( SpectrumLoaderBase loader,
const _HistAxis< _Var< T >> &  xAxis,
const _HistAxis< _Var< T >> &  yAxis,
const _Cut< T > &  cut,
const SystShifts shift = kNoShift,
const _Var< T > &  wei = UnweightedT >(),
ESparse  sparse = kDense 
)

2D Spectrum taking 2 HistAxis

template<class T >
ana::Spectrum::Spectrum ( const std::string xLabel,
const std::string yLabel,
SpectrumLoaderBase loader,
const Binning binsx,
const _Var< T > &  varx,
const Binning binsy,
const _Var< T > &  vary,
const _Cut< T > &  cut,
const SystShifts shift = kNoShift,
const _Var< T > &  wei = UnweightedT >(),
ESparse  sparse = kDense 
)

2D Spectrum of two Vars

template<class T >
ana::Spectrum::Spectrum ( SpectrumLoaderBase loader,
const _HistAxis< _Var< T >> &  xAxis,
const _HistAxis< _Var< T >> &  yAxis,
const _HistAxis< _Var< T >> &  zAxis,
const _Cut< T > &  cut,
const SystShifts shift = kNoShift,
const _Var< T > &  wei = UnweightedT >(),
ESparse  sparse = kDense 
)

3D Spectrum taking 3 HistAxis

template<class T >
ana::Spectrum::Spectrum ( const std::string xLabel,
const std::string yLabel,
const std::string zLabel,
SpectrumLoaderBase loader,
const Binning binsx,
const _Var< T > &  varx,
const Binning binsy,
const _Var< T > &  vary,
const Binning binsz,
const _Var< T > &  varz,
const _Cut< T > &  cut,
const SystShifts shift = kNoShift,
const _Var< T > &  wei = UnweightedT >(),
ESparse  sparse = kDense 
)

3D Spectrum of three Vars

ana::Spectrum::Spectrum ( Hist &&  hist,
const LabelsAndBins axis,
double  pot,
double  livetime 
)
inline

Expert constructor for ReweightableSpectrum et al.

Definition at line 135 of file Spectrum.h.

139  : fHist(std::move(hist)), fPOT(pot), fLivetime(livetime), fAxis(axis)
140  {
141  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
double fPOT
Definition: Spectrum.h:285
#define pot
double fLivetime
Definition: Spectrum.h:286
double livetime
Definition: saveFDMCHists.C:21
ana::Spectrum::~Spectrum ( )
virtual

Definition at line 78 of file Spectrum.cxx.

References fReferences.

Referenced by Uninitialized().

79  {
80  // Unregister self from anything that might still want to fill us
81  for(Spectrum** ref: fReferences) *ref = 0;
82  }
std::set< Spectrum ** > fReferences
Things that point at this Spectrum. Maintained by SpectrumLoader.
Definition: Spectrum.h:289
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
ana::Spectrum::Spectrum ( const Spectrum rhs)

Definition at line 53 of file Spectrum.cxx.

References ana::assert(), and fReferences.

53  :
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  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
double fPOT
Definition: Spectrum.h:285
double fLivetime
Definition: Spectrum.h:286
assert(nhit_max >=nhit_nbins)
ana::Spectrum::Spectrum ( Spectrum &&  rhs)
ana::Spectrum::Spectrum ( )
inlineprotected

Constructor for Uninitialized()

Definition at line 268 of file Spectrum.h.

References AddLoader(), kDense, PlusEqualsHelper(), RemoveLoader(), and canMan::sign().

Referenced by Uninitialized().

270  fPOT(0), fLivetime(0),
271  fAxis(std::vector<std::string>(), std::vector<Binning>())
272  {
273  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
double fPOT
Definition: Spectrum.h:285
double fLivetime
Definition: Spectrum.h:286
static Hist Uninitialized()
Definition: Hist.h:38
ana::Spectrum::Spectrum ( const LabelsAndBins axis,
ESparse  sparse = kDense 
)
protected

Helper for constructors.

Definition at line 23 of file Spectrum.cxx.

References fAxis, fHist, ana::LabelsAndBins::GetBins1D(), kSparse, ana::Binning::NBins(), ana::Hist::Zero(), and ana::Hist::ZeroSparse().

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  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
double fPOT
Definition: Spectrum.h:285
static Hist ZeroSparse(int nbins)
Definition: Hist.cxx:51
const Binning & GetBins1D() const
Appropriate binning and labelling for that 1D Var.
double fLivetime
Definition: Spectrum.h:286
static Hist Zero(int nbins)
Definition: Hist.cxx:41
static Hist Uninitialized()
Definition: Hist.h:38

Member Function Documentation

void ana::Spectrum::AddLoader ( Spectrum **  ref)
protected

Definition at line 373 of file Spectrum.cxx.

References fReferences.

Referenced by Spectrum().

374  {
375  fReferences.insert(ref);
376  }
std::set< Spectrum ** > fReferences
Things that point at this Spectrum. Maintained by SpectrumLoader.
Definition: Spectrum.h:289
Spectrum ana::Spectrum::AsimovData ( double  pot) const

Asimov data is a MC spectrum scaled to the POT expected in the data.

Use for sensitivity plots and testing fit convergence

Definition at line 318 of file Spectrum.cxx.

References fHist, fPOT, pot, ana::Hist::ResetErrors(), runNovaSAM::ret, and ana::Hist::Scale().

Referenced by FakeData(), and GetEigenStan().

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  }
double fPOT
Definition: Spectrum.h:285
#define pot
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
Spectrum ana::Spectrum::AsimovData ( double  pot,
double  livetime 
) const

Use for sensitivity plots when fake cosmic data is needed. Fake cosmic spectra can be added to FakeData by desired livetime.

Definition at line 332 of file Spectrum.cxx.

References fHist, fLivetime, fPOT, livetime, pot, ana::Hist::ResetErrors(), runNovaSAM::ret, and ana::Hist::Scale().

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  }
double fPOT
Definition: Spectrum.h:285
#define pot
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
double livetime
Definition: saveFDMCHists.C:21
void ana::Spectrum::Clear ( void  )
Spectrum ana::Spectrum::FakeData ( double  pot) const
Spectrum ana::Spectrum::FakeData ( double  pot,
double  livetime 
) const

Synonymous with AsimovData(). Retained for compatibility.

Definition at line 355 of file Spectrum.cxx.

References AsimovData().

356  {
357  return AsimovData(pot, livetime);
358  }
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
double livetime
Definition: saveFDMCHists.C:21
void ana::Spectrum::Fill ( double  x,
double  w = 1 
)

Definition at line 293 of file Spectrum.cxx.

References fAxis, fHist, ana::Hist::Fill(), ana::LabelsAndBins::GetBins1D(), w, and submit_syst::x.

Referenced by ana::SpectrumSink::HandleRecord(), ana::MultiVarSpectrumSink::HandleRecord(), ana::ReweightableSpectrumSink::HandleRecord(), and Uninitialized().

294  {
295  fHist.Fill(fAxis.GetBins1D(), x, w);
296  // TODO Pull binning out of Hist entirely and just update an index?
297  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
const Binning & GetBins1D() const
Appropriate binning and labelling for that 1D Var.
void Fill(const Binning &bins, double x, double w)
Definition: Hist.cxx:265
Float_t w
Definition: plot.C:20
const std::vector<Binning>& ana::Spectrum::GetBinnings ( ) const
inline
const Eigen::ArrayXd& ana::Spectrum::GetEigen ( ) const
inline
Eigen::ArrayXd ana::Spectrum::GetEigen ( double  exposure,
EExposureType  expotype = kPOT 
) const

Definition at line 219 of file Spectrum.cxx.

References fHist, fLivetime, fPOT, ana::Hist::GetEigen(), and ana::kPOT.

220  {
221  if(expotype == kPOT)
222  return (exposure/fPOT) * fHist.GetEigen();
223  else
224  return (exposure/fLivetime) * fHist.GetEigen();
225  }
const Eigen::ArrayXd & GetEigen() const
Definition: Hist.h:53
double fPOT
Definition: Spectrum.h:285
double fLivetime
Definition: Spectrum.h:286
const Eigen::ArrayXstan& ana::Spectrum::GetEigenStan ( ) const
inline
Eigen::ArrayXstan ana::Spectrum::GetEigenStan ( double  exposure,
EExposureType  expotype = kPOT 
) const

Definition at line 228 of file Spectrum.cxx.

References fHist, fLivetime, fPOT, ana::Hist::GetEigenStan(), and ana::kPOT.

229  {
230  if(expotype == kPOT)
231  return (exposure/fPOT) * fHist.GetEigenStan();
232  else
233  return (exposure/fLivetime) * fHist.GetEigenStan();
234  }
double fPOT
Definition: Spectrum.h:285
double fLivetime
Definition: Spectrum.h:286
const Eigen::ArrayXstan & GetEigenStan() const
Definition: Hist.h:54
const std::vector<std::string>& ana::Spectrum::GetLabels ( ) const
inline
bool ana::Spectrum::HasStan ( ) const
inline

Definition at line 187 of file Spectrum.h.

References fHist, and ana::Hist::HasStan().

Referenced by ana::SingleSampleExperiment::PredHistIncCosmics(), and ana::PredictionInterp::ShiftSpectrum().

187 {return fHist.HasStan();}
bool HasStan() const
Definition: Hist.h:52
double ana::Spectrum::Integral ( double  exposure,
double *  err = 0,
EExposureType  expotype = kPOT 
) const

Return total number of events scaled to pot.

Parameters
exposurePOT/livetime to scale to
errThe statistical error on this total (optional)
expotypeWhat the first parameter represents

Definition at line 249 of file Spectrum.cxx.

References fHist, fLivetime, fPOT, ana::Hist::GetBinError(), ana::Hist::GetNbinsX(), MECModelEnuComparisons::i, ana::Hist::Integral(), ana::kPOT, PandAna.reco_validation.prod5_pid_validation::ratio(), util::sqr(), and std::sqrt().

Referenced by BicountEllipse_dCP(), ana::CountingExperiment::ChiSq(), ComputeEfficiency(), ana::DataMCAreaNormalizedRatio(), demo_flat(), ana::SpectrumComponents::DrawLegend(), ana::DataMCPair::DrawMCNormSyst(), DataMCPair::DrawMCNormSyst(), FCContour(), ana::FluxMultiverseSyst::FindSigmaBoundaries(), FitParamEffectsAna(), FitSystEffectsAna(), GetEigenStan(), getHists_FNEX(), GetNueCosmics2017(), GetNueCosmics2018(), GetNueCosmics2019(), ana::GetNueCosmics2020(), GetNueCosmicsFuture(), GetNueData2017(), GetNueData2018(), GetNueData2019(), ana::GetNueData2020(), GetNuePrediction2017(), GetNuePrediction2018(), GetNuePrediction2019(), ana::GetNuePrediction2020(), GetNuePredictionFuture(), ana::GetNumuCosmics2020(), ana::GetNumuData2020(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getTemplateShapeOnly1D(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getTemplateShapeOnly3D(), Integral(), joint_fit_2017_make_fc_slice(), joint_fit_2017_make_fc_surf(), make_fc_mass_and_oct_nersc_2018(), make_fc_mass_and_oct_nersc_2019(), make_fc_mh_nersc_2018(), make_fc_mh_nersc_2019(), make_fc_oct_nersc_2018(), make_fc_oct_nersc_2019(), make_fc_slices_nersc_2018(), make_fc_slices_nersc_2018_stats(), make_fc_slices_nersc_2019(), make_fc_surfaces_nersc_2018(), make_fc_surfaces_nersc_2018_stats(), make_fc_surfaces_nersc_2019(), make_nominal_xs(), make_xs(), MakeHBar(), MakeTable(), MCCompPredictionTable(), ana::GenieMultiverseNormalizedSpectra::NormalizeSpectra(), PeripheralCuts(), ana::PlotWithAreaSystErrorBand(), predEventCountWithSystError(), SpectrumParamEffectsAna(), and ana::UnfoldSVD::Truth().

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  }
double Integral() const
Definition: Hist.cxx:252
T sqrt(T number)
Definition: d0nt_math.hpp:156
int GetNbinsX() const
Definition: Hist.cxx:228
double fPOT
Definition: Spectrum.h:285
def ratio(spec1, spec2, nbins, pot, binrange=(0, 1))
T sqr(T x)
More efficient square function than pow(x,2)
Definition: MathUtil.h:23
double fLivetime
Definition: Spectrum.h:286
double GetBinError(int i) const
Definition: Hist.cxx:241
double ana::Spectrum::Livetime ( ) const
inline

Seconds. For informational purposes only. No calculations use this.

Definition at line 230 of file Spectrum.h.

References fLivetime.

Referenced by ApplyOscillations(), ana::FitInAnaBinsBkgdEstimator::Background(), ana::NumuCC2p2hBkgdEstimator::Background(), BlessedPlotsLoad(), CalcRWithSystsNus17(), Cumulative(), ana::TwoSampleDecomp::Decomp(), FCContour(), fd_plot(), FitParamEffectsAna(), FitSystEffectsAna(), getBNBPlots(), GetFluxError(), getHists_FNEX(), ana::Multiverse::GetNSigmaShift(), GetNueCosmics2017(), GetNueCosmics2018(), GetNueCosmics2019(), ana::GetNueCosmics2020(), GetNueCosmicsFuture(), GetNueData2017(), GetNueData2018(), GetNueData2019(), ana::GetNueData2020(), ana::GetNumuCosmics2020(), ana::GetNumuData2020(), ana::PredictionInterp::InitFits(), make_fc_nus_surfs_nersc_2018(), make_fc_nus_surfs_nersc_2019(), make_nus17_fc_surfs(), make_nus_fc_surfs(), make_xs(), make_xs_1D(), MakeExtrapSurface(), MakeHBar(), makeMatrixElementSurface(), MakeNusPrediction(), ana::MichelDecomp::MCToDCMPComp(), ana::PredictionExtendToPeripheral::MergePeripheral(), mrbrem_get_reweighted_spectra(), mrbrem_plots(), PlotNus17Prediction(), PlotNus17PredSystsData(), PlotNus18Sideband(), PlotNusSensAna01(), ana::PlotSpectra(), PlotSysts(), ana::SingleSampleExperiment::PredHistIncCosmics(), ana::PredictionExtendToPeripheral::ReduceHelperNC(), selection_story_plots(), ShiftedCosmics(), ana::PredictionInterp::ShiftSpectrum(), ana::SingleSampleExperiment::SingleSampleExperiment(), timingPeak(), ana::UnfoldSVD::Truth(), and ana::UnfoldTikhonov::Truth().

230 {return fLivetime;}
double fLivetime
Definition: Spectrum.h:286
std::unique_ptr< Spectrum > ana::Spectrum::LoadFrom ( TDirectory *  dir,
const std::string name 
)
static

Definition at line 535 of file Spectrum.cxx.

References ana::assert(), ana::bins, dir, genie::utils::style::Format(), ana::Hist::FromDirectory(), MECModelEnuComparisons::i, label, PandAna.Demos.pi0_spectra::labels, ana::Binning::LoadFrom(), runNovaSAM::ret, string, cvnie::subdir, and getGoodRuns4SAM::tag.

Referenced by Analyse_GetEfficiency(), Analyse_GetEfficiency_UseNEntries(), bdtstudyplotter(), BlessedPlotsAnaByPeriod(), CalcRWithSysts(), calculateComponentsNumu(), calculateWrongSignNue(), calculateWrongSignNumuQ1(), calculateWrongSignNumuQ2(), calculateWrongSignNumuQ3(), calculateWrongSignNumuQ4(), CalculateXSec(), ana::CompareNDDataMC(), ana::CompareNDDataOneMC(), ana::CompareNDDataTwoMC(), CutTableAna(), CVNphotonHist(), CVNphotonSplitHist(), demo2p5b(), DoThePlots(), doUnfolding(), EHadVisMECpairs(), example_plot(), FDDataMC(), FDDataMCSystBandAna(), FidOptHist(), FillFlavorContainers(), FitParamEffectsAna(), FitSystEffectsAna(), FOMCalc(), GeniePredictionToRoot(), get_data_histogram(), get_numu_data_histogram(), GetCosmics(), Plotter::GetDataPlots(), GetHist(), getHists_FNEX(), GetHistVectors(), Plotter::GetMCPlots(), ana::GetNDComponents(), ana::GetNDMCComponents(), ana::GetNumuCosmics2020(), GetNumuCosmicsFuture(), GetNumuData2017(), GetNumuData2018(), GetNumuData2019(), GetSpectToHist(), ana::InteractionSpectra::InteractionSpectra(), Load1DHistFromSpec(), LoadCosmic(), LoadFakeData(), ana::CountingExperiment::LoadFrom(), jw::TrivialPrediction::LoadFrom(), ana::CheatDecomp::LoadFrom(), ana::FluxReweight::LoadFrom(), ana::NCDecomp::LoadFrom(), ana::NueDecomp::LoadFrom(), ana::SingleSampleExperiment::LoadFrom(), ana::NumuDecomp::LoadFrom(), ana::CrossSectionSpectra::LoadFrom(), ana::TrivialPrediction::LoadFrom(), ana::ProportionalDecomp::LoadFrom(), ana::CutOptimization::LoadFrom(), ana::FakeDecomp::LoadFrom(), ana::TrivialExtrap::LoadFrom(), ana::TrivialCrossSectionAnalysis::LoadFrom(), ana::Multiverse::LoadFrom(), ana::TwoSampleDecomp::LoadFrom(), ana::SpectrumComponents::LoadFrom(), ana::FluxDecomp::LoadFrom(), ana::ModularExtrapSterile::LoadFrom(), ana::BENDecomp::LoadFrom(), ana::DataMCPair::LoadFrom(), ana::RecoReweight::LoadFrom(), ana::LoadFromUpDownSpectra(), LoadMaps(), make_dataMC(), make_eff_plots_areaNorm(), make_estimate_energy(), make_muonid_opt(), make_nominal_xs(), Make_NuMuCC_Inc_XS(), Make_NuMuCC_Inc_XS_v2(), make_vertex_optimiz(), MakeCutFlow(), makeEnergyEstimator(), makeFlatWeight(), MakePlots(), makeResolutionPlots(), MakeSurface(), MakeSurfaceBinningStudy(), MakeThePlots(), makeXSecPlots1D(), makeXSecPlots2D(), makeXSecPlots_TemplateFit(), meanWeight_plot(), mec_nux_tester_2020(), mec_tuning_fitter_2020(), MRDiFStudy_FHC_Step2(), MRDiFStudy_RHC_Step2(), multiverse_efficiency_plot(), MuonCatcherComp_ProdPlots(), nc_bkgd_by_interaction_mode(), NDDataMC(), NDDataMCSystBandAna(), nue_data_mc_validation(), nue_decomp_scales(), nue_fd_mc_validation(), numu_data_mc_validation(), OverrideLivetime(), plot_DataMCComp_numu(), plot_estimate_energy(), plot_ND_DataMC(), plot_ND_DataMC_energybinning(), plot_nd_spectra_2018(), plot_recoE_numu(), plot_spectra_2dplots(), plot_uncertainty(), plotHist_SpectrumCVNID(), PlotResolution(), Plotsidebandfittest(), Plotting_Data2DataComp(), Plotting_Data2DataComp_SingleCanvas(), Plotting_DataAndPrediction(), Plotting_DataSpectra_MakeTGraph(), Plotting_OverlayStudies(), ppfx_smooth_weights_make(), preselection_cutflow(), ProduceCompPlots(), ReMIdHist(), resolutionplotter(), ana::ResolutionScan::ResolutionScan(), ShwZOptHist(), Spec2DToHist(), Spec2DtoHist(), SpecToHist(), SpectrumParamEffectsAna(), ana::TrivialPrediction::TrivialPrediction(), twodvtxcontplotter(), Unfold1D(), Unfold3D(), vertexstudyploter(), xsec_tot_uncert_optimization(), and xsec_uncertainty_per_bin().

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  }
static std::unique_ptr< Binning > LoadFrom(TDirectory *dir, const std::string &name)
Definition: Binning.cxx:230
const XML_Char * name
Definition: expat.h:151
subdir
Definition: cvnie.py:7
static Hist FromDirectory(TDirectory *dir)
Definition: Hist.cxx:158
const char * label
const Binning bins
Definition: NumuCC_CPiBin.h:8
TDirectory * dir
Definition: macro.C:5
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
enum BeamMode string
double ana::Spectrum::Mean ( ) const

Return mean of 1D histogram.

Definition at line 267 of file Spectrum.cxx.

References ana::bins, om::cout, ana::Binning::Edges(), allTimeWatchdog::endl, fAxis, fHist, ana::Hist::GetBinContent(), ana::LabelsAndBins::GetBins1D(), MECModelEnuComparisons::i, extractScale::mean, ana::Binning::NBins(), w, W, and x1.

Referenced by ana::DataMCPair::CreateSystTable(), and GetEigenStan().

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  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
Float_t x1[n_points_granero]
Definition: compare.C:5
const Binning & GetBins1D() const
Appropriate binning and labelling for that 1D Var.
OStream cout
Definition: OStream.cxx:6
const Binning bins
Definition: NumuCC_CPiBin.h:8
double GetBinContent(int i) const
Definition: Hist.cxx:349
Float_t w
Definition: plot.C:20
#define W(x)
Spectrum ana::Spectrum::MockData ( double  pot,
int  seed = 0 
) const

Mock data is FakeData with Poisson fluctuations applied.

Use for low-budget MDCs, or just getting a sense of the expected scale of statistical variation. NB seed = 0 is true random

Definition at line 300 of file Spectrum.cxx.

References FakeData(), fHist, ana::Hist::GetBinContent(), ana::Hist::GetNbinsX(), MECModelEnuComparisons::i, ana::Hist::ResetErrors(), runNovaSAM::ret, generate_hists::rnd, and ana::Hist::SetBinContent().

Referenced by demo4(), demo5(), FCTutorial2020(), fd_plot(), fill_col(), GenerateFutureData(), ana2019::fakedata::Get2019Prediction(), GetEigenStan(), GetMockData(), joint_fit_2017_make_fc_slice(), joint_fit_2017_make_fc_surf(), make_fc_mass_and_oct_nersc_2018(), make_fc_mass_and_oct_nersc_2019(), make_fc_mh_nersc_2018(), make_fc_mh_nersc_2019(), make_fc_nus_surfs_nersc_2018(), make_fc_nus_surfs_nersc_2019(), make_fc_oct_nersc_2018(), make_fc_oct_nersc_2019(), make_fc_slices_nersc_2018(), make_fc_slices_nersc_2018_stats(), make_fc_slices_nersc_2019(), make_fc_surfaces_2020_validation(), make_fc_surfaces_nersc_2018(), make_fc_surfaces_nersc_2018_stats(), make_fc_surfaces_nersc_2019(), make_mockdata_syst_contours(), make_nus17_fc_surfs(), make_nus_fc_surfs(), MakeFakeData(), selection_story_plots(), syst_plot_test(), and test_saloaders().

301  {
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  }
unsigned int seed
Definition: runWimpSim.h:102
Spectrum FakeData(double pot) const
Synonymous with AsimovData(). Retained for compatibility.
Definition: Spectrum.cxx:349
#define pot
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
unsigned int ana::Spectrum::NDimensions ( ) const
inline
Spectrum ana::Spectrum::operator* ( const Ratio rhs) const

Definition at line 483 of file Spectrum.cxx.

References runNovaSAM::ret.

Referenced by OverrideLivetime().

484  {
485  Spectrum ret = *this;
486  ret *= rhs;
487  return ret;
488  }
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
Spectrum & ana::Spectrum::operator*= ( const Ratio rhs)

Definition at line 476 of file Spectrum.cxx.

References ana::Ratio::fHist, fHist, and ana::Hist::Multiply().

Referenced by OverrideLivetime().

477  {
478  fHist.Multiply(rhs.fHist);
479  return *this;
480  }
void Multiply(const Hist &rhs)
Definition: Hist.cxx:479
Spectrum ana::Spectrum::operator+ ( const Spectrum rhs) const

Definition at line 454 of file Spectrum.cxx.

References runNovaSAM::ret.

Referenced by OverrideLivetime().

455  {
456  Spectrum ret = *this;
457  ret += rhs;
458  return ret;
459  }
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
Spectrum & ana::Spectrum::operator+= ( const Spectrum rhs)

Definition at line 448 of file Spectrum.cxx.

References PlusEqualsHelper().

Referenced by OverrideLivetime().

449  {
450  return PlusEqualsHelper(rhs, +1);
451  }
Spectrum & PlusEqualsHelper(const Spectrum &rhs, int sign)
Helper for operator+= and operator-=.
Definition: Spectrum.cxx:379
Spectrum ana::Spectrum::operator- ( const Spectrum rhs) const

Definition at line 468 of file Spectrum.cxx.

References runNovaSAM::ret.

Referenced by OverrideLivetime().

469  {
470  Spectrum ret = *this;
471  ret -= rhs;
472  return ret;
473  }
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
Spectrum & ana::Spectrum::operator-= ( const Spectrum rhs)

Definition at line 462 of file Spectrum.cxx.

References PlusEqualsHelper().

Referenced by OverrideLivetime().

463  {
464  return PlusEqualsHelper(rhs, -1);
465  }
Spectrum & PlusEqualsHelper(const Spectrum &rhs, int sign)
Helper for operator+= and operator-=.
Definition: Spectrum.cxx:379
Spectrum ana::Spectrum::operator/ ( const Ratio rhs) const

Definition at line 498 of file Spectrum.cxx.

References runNovaSAM::ret.

Referenced by OverrideLivetime().

499  {
500  Spectrum ret = *this;
501  ret /= rhs;
502  return ret;
503  }
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268
Spectrum & ana::Spectrum::operator/= ( const Ratio rhs)

Definition at line 491 of file Spectrum.cxx.

References ana::Hist::Divide(), ana::Ratio::fHist, and fHist.

Referenced by OverrideLivetime().

492  {
493  fHist.Divide(rhs.fHist);
494  return *this;
495  }
void Divide(const Hist &rhs)
Definition: Hist.cxx:531
Spectrum & ana::Spectrum::operator= ( const Spectrum rhs)

Definition at line 63 of file Spectrum.cxx.

References ana::assert(), fAxis, fHist, fLivetime, fPOT, and fReferences.

Referenced by Uninitialized().

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  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
std::set< Spectrum ** > fReferences
Things that point at this Spectrum. Maintained by SpectrumLoader.
Definition: Spectrum.h:289
double fPOT
Definition: Spectrum.h:285
double fLivetime
Definition: Spectrum.h:286
assert(nhit_max >=nhit_nbins)
Spectrum& ana::Spectrum::operator= ( Spectrum &&  rhs)
void ana::Spectrum::OverrideLivetime ( double  newlive)
inline
void ana::Spectrum::OverridePOT ( double  newpot)
inline
Spectrum & ana::Spectrum::PlusEqualsHelper ( const Spectrum rhs,
int  sign 
)
protected

Helper for operator+= and operator-=.

Definition at line 379 of file Spectrum.cxx.

References ana::Hist::Add(), ana::AlmostEqual(), om::cout, allTimeWatchdog::endl, fHist, fLivetime, fPOT, ana::Hist::Initialized(), and ana::Hist::Integral().

Referenced by operator+=(), operator-=(), and Spectrum().

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  }
void Add(const Hist &rhs, double scale=1)
Definition: Hist.cxx:455
double fPOT
Definition: Spectrum.h:285
bool AlmostEqual(double a, double b, double eps)
Definition: UtilsExt.cxx:40
double fLivetime
Definition: Spectrum.h:286
OStream cout
Definition: OStream.cxx:6
def sign(x)
Definition: canMan.py:197
double ana::Spectrum::POT ( ) const
inline

Definition at line 227 of file Spectrum.h.

References fPOT.

Referenced by ana::nueccinc::NueCCIncCrossSectionAnalysis::AddEnhancedSample(), ana::DataMCPair::AddExposure(), ApplyOscillations(), caf_numu_fd_validation_data(), CalcChi2(), CalcR(), CalcRWithSystsNus17(), calculateComponentsNumu(), calculateWrongSignNue(), calculateWrongSignNumuQ1(), calculateWrongSignNumuQ2(), calculateWrongSignNumuQ3(), calculateWrongSignNumuQ4(), ana::CountingExperiment::ChiSq(), ana::SingleSampleExperiment::ChiSq(), ana::CovMxExperiment::ChiSq(), CompareMissingLeptons(), ComputeEfficiency(), Cumulative(), ana::DataMCAreaNormalizedRatio(), ana::DataMCComparison(), Plotter::DataMCComparison(), ana::DataMCComparisonAreaNormalized(), ana::DataMCComparisonComponents(), ana::TwoSampleDecomp::Decomp(), demo6(), draw_decomp_plots(), ana::DataMCPair::DrawData(), ana::SpectrumComponents::DrawLegend(), ana::DataMCPair::DrawMCComponents(), ana::DataMCPair::DrawMCNormSyst(), DataMCPair::DrawMCNormSyst(), ana::DataMCPair::DrawMCSyst(), drawPlot(), DrawPreSelectionPlots(), ana::TwoSampleDecomp::DrawSigBkgOverlay(), ana::TwoSampleDecomp::DrawTwoSamplesWithRatios(), EHadVisMECpairs(), example_plot(), FCContour(), fd_plot(), fill_col(), ana::FluxMultiverseSyst::FindSigmaBoundaries(), FitParamEffectsAna(), FitSystEffectsAna(), GetExtrap(), getHists_FNEX(), ana::Multiverse::GetNSigmaShift(), GetNueData2017(), GetNueData2018(), GetNueData2019(), ana::GetNueData2020(), ana::GetNumuData2020(), ana::MichelDecomp::GetSum(), ana::MichelDecomp::GetTemplateContent(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getTemplateShapeOnly1D(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getTemplateShapeOnly3D(), ana::PredictionInterp::InitFits(), ana::SingleSampleExperiment::LogLikelihood(), Make2DPlot(), make_muonid_opt(), make_nus17_fc_surfs(), make_nus_fc_surfs(), make_plots(), make_vertex_optimiz(), make_xs(), make_xs_1D(), makeEnergyEstimator(), makeFlatWeight(), MakeHBar(), makeMatrixElementSurface(), MakePlot(), MakePlots(), ana::MakeQuantileHistogram(), MakeSelectionPlots(), ana::FluxDecomp::MakeWeightsNumuFromKaon(), ana::BENDecomp::MakeWeightsNumuFromKaon(), ana::FluxDecomp::MakeWeightsNumuFromPion(), ana::BENDecomp::MakeWeightsNumuFromPion(), ana::MichelDecomp::MCToDCMPComp(), mec_nux_tester_2020(), mec_tuning_fitter_2020(), ana::PredictionExtendToPeripheral::MergePeripheral(), MichelDecompTest(), mrbrem_get_reweighted_spectra(), mrbrem_plots(), mre_blessed(), mre_comp_split(), multiverse_efficiency_plot(), multiverse_macro(), ana::GenieMultiverseNormalizedSpectra::NormalizeSpectra(), nue_decomp_scales(), nue_decomp_scales_for_make_decomp(), ana::FluxDecomp::NueEstimate(), ana::BENDecomp::NueEstimate(), ana::FluxDecomp::NueEstimateFromKa(), ana::FluxDecomp::NueEstimateFromPi(), ana::BENDecomp::NueEstimateFromPi(), ana::CutOptimization::OptimizedSigmaOverSigma(), plot(), plot_2d(), plot_2d_vars(), plot_3NDvsFD(), plot_3NDvsFD_FHC(), plot_3NDvsFD_RHC(), plot_datamc_ND_numu(), plot_datamc_ND_numu_REW(), plot_diff(), plot_nd_data_mc(), plot_ND_numu_NOMvsREW(), plot_NDvsFD_REW(), plot_NDvsFD_weights(), plot_NDvsFD_weights_FHC(), plot_NDvsFD_weights_RHC(), plot_predictions(), plot_recoE_numu(), plot_time(), ana::PlotAllSelectionDecomposition(), ana::PlotAllSignalEstimates(), ana::CrossSectionAnalysis::PlotBackgroundEstimate(), ana::NumuCCIncAnalysis::PlotBackgroundEstimate(), ana::SingleNucAnalysis::PlotBackgroundEstimate2D(), ana::NumuCC2p2hAnalysis::PlotBackgroundEstimate2D(), ana::CrossSectionAnalysis::PlotData(), ana::NumuCCIncAnalysis::PlotData(), ana::SingleNucAnalysis::PlotData2D(), ana::NumuCC2p2hAnalysis::PlotData2D(), ana::CutOptimization::PlotDebug(), PlotNueDecompFourBins(), PlotNus17PredSystsData(), PlotNus18Sideband(), PlotNusSensAna01(), ana::CrossSectionAnalysis::PlotRecoToTrueMatrix(), ana::SingleNucAnalysis::PlotRecoToTrueMatrix2D(), ana::NumuCC2p2hAnalysis::PlotRecoToTrueMatrix2D(), ana::CrossSectionAnalysis::PlotSignalEstimate(), ana::SingleNucAnalysis::PlotSignalEstimate2D(), ana::NumuCC2p2hAnalysis::PlotSignalEstimate2D(), ana::PlotSpectra(), PlotStack(), ana::PlotStack(), PlotSyst(), PlotSysts(), ana::CrossSectionAnalysis::PlotUnfoldedSignal(), ana::SingleNucAnalysis::PlotUnfoldedSignal2D(), ana::NumuCC2p2hAnalysis::PlotUnfoldedSignal2D(), PlotVertices(), ana::SingleSampleExperiment::PredHistIncCosmics(), ana::NDPredictionSterile::Predict(), ana::NDPredictionSterile::PredictComponent(), ana::Ratio::Ratio(), ana::TwoSampleDecomp::RatioCalc(), reco_minus_true_panels(), ana::PredictionExtendToPeripheral::ReduceHelperNC(), ana::CrossSectionAnalysis::Result(), ana::NumuCCIncAnalysis::Result(), ana::NumuCC2p2hAnalysis::Result1DFluxInt(), ana::SingleNucAnalysis::Result2D(), ana::NumuCC2p2hAnalysis::Result2D(), ana::NumuCCIncAnalysis::Result2D(), ana::NumuCCIncAnalysis::ReturnHists(), ana::MichelDecomp::SaveCompPlots(), ana::FluxDecomp::SavePlots(), ana::BENDecomp::SavePlots(), ana::TruthReweight::SavePlots(), ana::RecoReweight::SavePlots(), ana::FluxDecomp::SavePlotsKa(), ana::BENDecomp::SavePlotsKa(), ana::FluxDecomp::SavePlotsPi(), ana::BENDecomp::SavePlotsPi(), saveS(), saveS1(), ana::MichelDecomp::SaveTempPlots(), selection_story_plots(), ShiftedCosmics(), ana::PredictionInterp::ShiftSpectrum(), ana::SimpleFOM(), SpectrumParamEffectsAna(), TableNueNDComponents(), TableNueNDDataMC(), test_beam_errorband(), test_micheldecomp(), test_sam(), test_sam_project(), test_stanfit_systpulls(), ana::Multiverse::ToHist(), ana::CutOptimization::ToHist(), ana::ICrossSectionAnalysis::ToHist(), ana::UnfoldIterative::Truth(), ana::UnfoldSVD::Truth(), and ana::UnfoldTikhonov::Truth().

227 {return fPOT;}
double fPOT
Definition: Spectrum.h:285
void ana::Spectrum::RemoveLoader ( Spectrum **  ref)
protected

Definition at line 367 of file Spectrum.cxx.

References fReferences.

Referenced by Spectrum().

368  {
369  fReferences.erase(ref);
370  }
std::set< Spectrum ** > fReferences
Things that point at this Spectrum. Maintained by SpectrumLoader.
Definition: Spectrum.h:289
void ana::Spectrum::SaveTo ( TDirectory *  dir,
const std::string name 
) const

Definition at line 506 of file Spectrum.cxx.

References dir, fAxis, fHist, fLivetime, genie::utils::style::Format(), fPOT, GetBinnings(), ana::LabelsAndBins::GetBins1D(), GetLabels(), MECModelEnuComparisons::i, NDimensions(), tmp, and ana::Hist::Write().

Referenced by Ana2017_box_opening_macro(), Ana2017_sb_opening_macro(), Ana2018_box_opening_macro(), angle(), bdtstudyspectrums(), BlessedPlotsLoad(), caf_nue_data_mc(), ccpiinc_mc_studies(), containmentstudy(), CosmicPred(), CutTableLoad(), CVNphoton(), CVNphotonSplit(), datamc_ND_numu_kinematics(), datamc_ND_numu_kinematics_FHC(), datamc_ND_numu_kinematics_FHC_pTBins(), datamc_ND_numu_kinematics_FHC_REW(), datamc_ND_numu_kinematics_FHC_REW_pTBins(), datamc_ND_numu_kinematics_REW(), datamc_ND_numu_kinematics_RHC(), datamc_ND_numu_kinematics_RHC_pTBins(), datamc_ND_numu_kinematics_RHC_REW(), datamc_ND_numu_kinematics_RHC_REW_pTBins(), dataprocess_numuccinc(), estimate_energy(), FDDataMCSystBandLoad(), FidOpt(), FidWShwCuts(), FillSpectra(), FitSystEffectsLoad(), get_cosmic_spectra(), get_data_and_cosmic(), get_data_histogram(), get_numi_data_histogram(), getData(), getFlatWeightSpectra(), getSpectra_ForFitting(), getTimePeakTotal(), make_mockdata_syst_contours(), make_pi0_xcheck(), make_pid(), MakeCosmics(), MakeFakeData(), makeFakeDataFluxes(), MakeNus17Prediction(), MakeNus18SidebandPred(), MakeNusPrediction(), makeRealDataFluxes(), makeXSecPlots_TemplateFit(), meanWeight_macro(), mrbrem_get_initial_spectra(), mrbrem_get_reweighted_spectra(), MRDiFStudy_FHC_Step1(), MRDiFStudy_RHC_Step1(), mre_blessed(), mre_comp_split(), ncpi0HistoGrid2(), ND_DataMC(), NDDataMCSystBandLoad(), nue_decomp_scales(), nus17_box_opening(), nus17_fd_cut_tables(), nus17_fd_cut_tables2D(), nus18_box_opening(), nus_ana01_box_opening(), nus_ana01_sideband_box_opening(), OverrideLivetime(), pion_multiverse(), ReMId(), resolutionscript(), resolutionstudy(), saveFDMCHists(), saveSpectraForUnf(), ana::xsec::UnfoldingVariable::SaveSpectrums(), ana::PredictionXSecTuning::SaveTo(), ana::CountingExperiment::SaveTo(), ana::PredictionNoOsc::SaveTo(), ana::FluxReweight::SaveTo(), ana::CheatDecomp::SaveTo(), jw::TrivialPrediction::SaveTo(), ana::FitInAnaBinsBkgdEstimator::SaveTo(), ana::NumuCC2p2hBkgdEstimator::SaveTo(), ana::nueccinc::NueCCIncMRECorrection::SaveTo(), ana::NCDecomp::SaveTo(), ana::NueDecomp::SaveTo(), ana::SingleSampleExperiment::SaveTo(), ana::NumuDecomp::SaveTo(), ana::CrossSectionSpectra::SaveTo(), ana::TrivialPrediction::SaveTo(), ana::PredictionScaleComp::SaveTo(), ana::ProportionalDecomp::SaveTo(), ana::CutOptimization::SaveTo(), ana::TrivialCrossSectionAnalysis::SaveTo(), ana::FakeDecomp::SaveTo(), ana::TrivialExtrap::SaveTo(), ana::CrossSectionAnalysis::SaveTo(), ana::TwoSampleDecomp::SaveTo(), ana::SpectrumComponents::SaveTo(), ana::NumuCCIncAnalysis::SaveTo(), ana::FluxDecomp::SaveTo(), ana::nueccinc::NueCCIncCrossSectionAnalysis::SaveTo(), ana::MichelDecomp::SaveTo(), ana::BENDecomp::SaveTo(), ana::DataMCPair::SaveTo(), ana::RecoReweight::SaveTo(), ana::nueccinc::NueCCIncEnhancedSamples::SaveTo(), mcmc::SaveToFile(), ShwZOpt(), sidebandfittest(), signal_count(), test_beam_errorband(), test_fluxhadr_prod_weights_Flux(), test_genieweights(), test_predictionscalecomp(), ana::TrivialPrediction::TrivialPrediction(), TrueNCSpectrum(), TrueSpectrumFromKaons(), twodstudyvtxcont(), validationscript(), and vertexstudy().

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  }
const std::vector< Binning > & GetBinnings() const
Definition: Spectrum.h:264
const XML_Char * name
Definition: expat.h:151
void Write(const Binning &bins) const
Definition: Hist.cxx:583
LabelsAndBins fAxis
Definition: Spectrum.h:291
double fPOT
Definition: Spectrum.h:285
Float_t tmp
Definition: plot.C:36
const Binning & GetBins1D() const
Appropriate binning and labelling for that 1D Var.
double fLivetime
Definition: Spectrum.h:286
unsigned int NDimensions() const
Definition: Spectrum.h:262
TDirectory * dir
Definition: macro.C:5
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:263
void ana::Spectrum::Scale ( double  c)
void ana::Spectrum::Scale ( const stan::math::var v)

Definition at line 243 of file Spectrum.cxx.

References fHist, and ana::Hist::Scale().

244  {
245  fHist.Scale(c);
246  }
void Scale(double s)
Definition: Hist.cxx:296
TH1D * ana::Spectrum::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.

Parameters
exposurePOT or livetime (seconds)
colHistogram color (default black)
styleHistogram line style (default solid)
expotypeHow to interpret exposure (kPOT (default) or kLivetime)

Definition at line 148 of file Spectrum.cxx.

References runNovaSAM::ret.

Referenced by ana::nueccinc::NueCCIncCrossSectionAnalysis::AddEnhancedSample(), ana::AddErrorInQuadrature(), ana::PredictionSystJoint2018::AddNormSyst(), ana::PredictionSystJointDemo::AddNormSyst(), ana::PredictionSyst3Flavor2020::AddNormSyst(), ana::PredictionSystNue2017::AddNormSysts(), ana::PredictionSystNumu2017::AddNormSysts(), ana::CovMxManager::AddSystematic(), AnalyzeNus18Pred(), ana::NumuCC2p2hBkgdEstimator::Background(), bin_composition_pie_chart(), ana::covmx::CovarianceMatrix::BuildFullCovMx(), caf_numu_fd_validation_data(), caf_numu_fd_validation_MC(), caf_numu_fd_validation_MC_no_tau(), calcAlphaBetaEachBin(), CalcChi2(), ana::nueccinc::NueCCIncMRECorrection::CalcEfficiencyCorrection(), CalcR(), CalcRWithSystsNus17(), calculateComponentsNumu(), calculateWrongSignNue(), calculateWrongSignNumuQ1(), calculateWrongSignNumuQ2(), calculateWrongSignNumuQ3(), calculateWrongSignNumuQ4(), cc(), check_predinterp(), check_predinterp_numu(), ana::OscCovMxExperiment::ChiSq(), ana::LikelihoodCovMxExperiment::ChiSq(), Compare(), CompareBinningSchemes(), CompareDecompDataMC(), CompareMCCompPrediction(), CompareMissingLeptons(), ComparePredictions(), ana::ComparePredictions(), ana::ModularExtrapComponent::ComparisonPlot(), Cumulative(), CVNCuts(), ana::DataMCComparison(), Plotter::DataMCComparison(), ana::DataMCComparisonAreaNormalized(), ana::DataMCComparisonComponents(), ana::PredictionInterp::DebugPlot(), ana::TwoSampleDecomp::Decomp(), ana::MichelDecomp::Decompose(), demo1(), demo3(), demo5(), demo6(), demo_flat(), demo_nueNumuSysts(), DrawBackgrounds(), DrawBins(), ana::DataMCPair::DrawData(), drawPlot(), DrawPulls(), DrawSensitivityDip(), ana::TwoSampleDecomp::DrawSigBkgOverlay(), DrawSurface(), DrawSystShifts(), drawSystsShiftingNDdata(), drawSystsShiftingNDdata_updatedAna(), ana::TwoSampleDecomp::DrawTwoSamplesWithRatios(), efficiency(), efficiencySA(), EHadVisMECpairs(), example_plot(), FCContour(), fd_plot(), fill_col(), FitParamEffectsAna(), FitSystEffectsAna(), GenerateFutureData(), GeniePredictionToRoot(), GetBackgroundStatisticalUncertainty(), GetBackgroundSystematicUncertainty(), ana::GetBFSystBands(), GetBG(), ana::GetBG(), getBNBPlots(), getContProf(), getContProf_Sensitivity(), GetDenominator(), ana::nueccinc::NueCCIncMRECorrection::getEfficiency1D(), GetEfficiencyDenominator(), GetEfficiencySystematicUncertainty(), GetExtrap(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getFlux1D(), GetFluxError(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getFullTemplate(), getHists_FNEX(), GetNC(), ana::GetNC(), ana::GetNDComponents(), ana::GetNDCompsFromDecomp(), ana::GetNDDecompsFromDecomp(), ana::Multiverse::GetNSigmaShift(), GetNueCosmics2017(), GetNueCosmics2018(), GetNueCosmics2019(), GetNueCosmicsFuture(), GetNumuCosmicsFuture(), ana::nueccinc::NueCCIncMRECorrection::getPreselected1D(), GetQuantilePredictionHist(), ana::nueccinc::NueCCIncMRECorrection::getSelected1D(), GetSelectedStatisticalUncertainty(), GetSig(), GetSpectra(), ana::GetSpectrum(), ana::GetSystBands(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getTemplateShapeOnly1D(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getTemplateShapeOnly3D(), getTimePeakPlots(), Integral(), make_extrap_figure(), make_fc_nus_surfs_nersc_2018(), make_muonid_opt(), make_nus17_fc_surfs(), make_nus_fc_surfs(), make_plots(), make_vertex_optimiz(), make_xs_1D(), MakeCovMx(), makeEnergyEstimator(), MakeFakeData(), makeFlatWeight(), makeMatrixElementSurface(), MakePeriCutPlot(), MakePlot(), MakePlots(), MakeResultPlot(), MakeSurface(), MakeSurfaceBinningStudy(), MakeTable(), ana::FluxDecomp::MakeWeightsNumuFromKaon(), ana::BENDecomp::MakeWeightsNumuFromKaon(), makeXSecPlots1D(), makeXSecPlots2D(), MatrixDeterminant(), MichelDecompTest(), mrbrem_get_reweighted_spectra(), mrbrem_plots(), mre_blessed(), mre_comp_split(), multiverse_efficiency_plot(), multiverse_macro(), NDHists(), nue_decomp_scales(), nue_decomp_scales_for_make_decomp(), nue_pid_effs(), nue_pid_effs_miniprod(), nue_pid_effs_paper_numu_neweff(), numu_cut_flow(), PeripheralCuts(), pi0_xcheck(), plot(), plot_3NDvsFD(), plot_3NDvsFD_FHC(), plot_3NDvsFD_RHC(), plot_datamc_ND_numu(), plot_datamc_ND_numu_REW(), plot_datamcpred(), plot_diff(), plot_kinematics_cafana(), plot_nd_data_mc(), plot_ND_numu_NOMvsREW(), plot_NDvsFD_REW(), plot_NDvsFD_weights(), plot_NDvsFD_weights_FHC(), plot_NDvsFD_weights_RHC(), plot_nue_xsec_pred(), plot_nueFD_Signal_REWvsNOM_FHC(), plot_nueFD_Signal_REWvsNOM_pTExtrap_FHC(), plot_nueFD_Signal_REWvsNOM_pTExtrap_RHC(), plot_nueFD_Signal_REWvsNOM_RHC(), plot_pi0_xcheck(), plot_pid(), plot_predictions(), plot_recoE_numu(), plot_time(), ana::PlotAllSelectionDecomposition(), ana::PlotAllSignalEstimates(), ana::CrossSectionAnalysis::PlotBackgroundEstimate(), ana::NumuCCIncAnalysis::PlotBackgroundEstimate(), PlotComp(), ana::CrossSectionAnalysis::PlotData(), ana::NumuCCIncAnalysis::PlotData(), ana::CutOptimization::PlotDebug(), ana::CrossSectionAnalysis::PlotFluxEstimate(), ana::SingleNucAnalysis::PlotFluxEstimate2D(), ana::NumuCC2p2hAnalysis::PlotFluxEstimate2D(), PlotNuePredictionFourBins(), PlotNumuPredData(), PlotNus17PredSystsData(), PlotNus18Sideband(), PlotNusSensAna01(), PlotPionPlots(), PlotPurEff(), ana::PlotPurEff(), plots(), ana::PlotSpectra(), PlotStack(), ana::PlotStack(), PlotSyst(), PlotSysts(), Plotting_DataAndPrediction(), ana::CrossSectionAnalysis::PlotUnfoldedSignal(), ana::PlotWithSystErrorBand(), ana::PlotWithSystErrorBand_Quant(), ana::PlotWithSystErrorBandTwoPreds(), ana::covmx::CovarianceMatrix::Predict(), PredictCosmic(), ana::SpectrumComponents::Purity(), ana::TwoSampleDecomp::RatioCalc(), reco_minus_true_panels(), resolution2018(), ana::CrossSectionAnalysis::Result(), ana::NumuCCIncAnalysis::Result(), ana::NumuCC2p2hAnalysis::Result1DFluxInt(), ana::SingleNucAnalysis::Result2D(), ana::NumuCC2p2hAnalysis::Result2D(), ana::NumuCCIncAnalysis::Result2D(), ana::NumuCCIncAnalysis::ReturnHists(), Save(), ana::MichelDecomp::SaveCompPlots(), SaveDCMPPlots(), ana::FluxReweight::SavePlots(), ana::FluxDecomp::SavePlots(), ana::BENDecomp::SavePlots(), ana::TruthReweight::SavePlots(), ana::RecoReweight::SavePlots(), ana::FluxDecomp::SavePlotsKa(), ana::BENDecomp::SavePlotsKa(), ana::BENDecomp::SavePlotsPi(), saveS(), ScaleCovarianceMatrix(), selection_story_plots(), ShiftedCosmics(), signal_count(), ana::SimpleFOM(), SpectrumParamEffectsAna(), sterile_demo(), syst_plot_test(), systematics_summary_from_pred_interp(), TableNuePredictionFourBins(), template_basic(), test_ana(), test_beam_errorband(), test_fluxhadr_prod_weights_Flux(), test_genie_systs(), test_genieweights(), test_micheldecomp(), test_nue2017Prediction(), test_nue2018_fitter(), test_nueextrapsyst(), test_nuwro(), test_saloaders(), test_sam(), test_sam_project(), test_stanfit_systpulls(), TestPred(), timingPeak(), ana::FluxMultiverseSyst::ToAreaNormalizedTH1(), ana::GenieMultiverseSpectra::ToAreaNormalizedTH1(), ana::Multiverse::ToHist(), ana::CutOptimization::ToHist(), ana::ICrossSectionAnalysis::ToHist(), ana::FluxMultiverseSyst::ToTH1(), ana::GenieMultiverseSpectra::ToTH1(), ana::ToTH2(), ana::ToTH3(), ToTHX(), Toy_analyses(), ana::UnfoldSVD::Truth(), ana::UnfoldTikhonov::Truth(), Unfold1D(), Unfold3D(), UnfoldInOut(), and Uninitialized().

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  }
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
Int_t col[ntarg]
Definition: Style.C:29
TH1D * ana::Spectrum::ToTH1 ( double  exposure,
EExposureType  expotype,
EBinType  bintype = kBinContent 
) const

Histogram made from this Spectrum, scaled to some exposure.

Parameters
exposurePOT or livetime (seconds)
expotypeHow to interpret exposure (kPOT (default) or kLivetime)

Definition at line 85 of file Spectrum.cxx.

References om::cout, allTimeWatchdog::endl, fAxis, fHist, fLivetime, fPOT, ana::LabelsAndBins::GetBins1D(), ana::LabelsAndBins::GetLabel1D(), ana::kBinDensity, ana::kLivetime, ana::kPOT, livetime, pot, runNovaSAM::ret, and ana::Hist::ToTH1().

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  }
LabelsAndBins fAxis
Definition: Spectrum.h:291
Divide bin contents by bin widths.
Definition: UtilsExt.h:31
TH1D * ToTH1(const Binning &bins) const
Definition: Hist.cxx:207
double fPOT
Definition: Spectrum.h:285
const std::string & GetLabel1D() const
const Binning & GetBins1D() const
Appropriate binning and labelling for that 1D Var.
#define pot
double fLivetime
Definition: Spectrum.h:286
OStream cout
Definition: OStream.cxx:6
double livetime
Definition: saveFDMCHists.C:21
TH2 * ana::Spectrum::ToTH2 ( double  exposure,
EExposureType  expotype = kPOT,
EBinType  bintype = kBinContent 
) const

Spectrum must be 2D to obtain TH2.

Definition at line 165 of file Spectrum.cxx.

References om::cout, allTimeWatchdog::endl, GetBinnings(), GetLabels(), ana::kBinDensity, NDimensions(), runNovaSAM::ret, and ana::ToTH2().

Referenced by ana::nueccinc::NueCCIncMRECorrection::CalcEfficiencyCorrection(), demo6(), ana::nueccinc::NueCCIncCrossSectionAnalysis::doUnfolding1D(), ana::nueccinc::NueCCIncCrossSectionAnalysis::doUnfolding2D(), drawLongTerm(), drawTimePlots(), drawVsPOT(), FD_Data_PosComp(), FD_plots(), GetHistsFD(), GetHistsND(), GetSpectra2D(), getTimePeakPlots(), ana::nueccinc::NueCCIncCrossSectionAnalysis::getUnfoldingMatrix(), HadEFit(), Make2DPlot(), make_quantiles_histogram_2020(), makeEnergyEstimator(), MakePlots(), ana::MakeQuantileHistogram(), mec_nux_tester_2020(), mec_tuning(), mec_tuning_fitter_2020(), MichelDecompTest(), MuonFit(), plot_2d(), plot_2d_vars(), plot_pi0_xcheck(), ana::SingleNucAnalysis::PlotBackgroundEstimate2D(), ana::NumuCC2p2hAnalysis::PlotBackgroundEstimate2D(), ana::SingleNucAnalysis::PlotData2D(), ana::NumuCC2p2hAnalysis::PlotData2D(), ana::SingleNucAnalysis::PlotUnfoldedSignal2D(), ana::NumuCC2p2hAnalysis::PlotUnfoldedSignal2D(), PlotVertices(), PositionComparison(), print_tables(), PrintPlot(), Save(), saveS1(), signal_count(), ana::Multiverse::ToHist(), ana::CutOptimization::ToHist(), ana::ICrossSectionAnalysis::ToHist(), ana::GenieMultiverseSpectra::ToTH2(), ToTHX(), and Uninitialized().

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  }
const std::vector< Binning > & GetBinnings() const
Definition: Spectrum.h:264
Divide bin contents by bin widths.
Definition: UtilsExt.h:31
unsigned int NDimensions() const
Definition: Spectrum.h:262
OStream cout
Definition: OStream.cxx:6
const std::vector< std::string > & GetLabels() const
Definition: Spectrum.h:263
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
TH3 * ana::Spectrum::ToTH3 ( double  exposure,
EExposureType  expotype = kPOT,
EBinType  bintype = kBinContent 
) const

Spectrum must be 3D to obtain TH3.

Definition at line 187 of file Spectrum.cxx.

References om::cout, allTimeWatchdog::endl, GetBinnings(), GetLabels(), ana::kBinDensity, NDimensions(), runNovaSAM::ret, and ana::ToTH3().

Referenced by ana::FluxDecomp::MakeWeightsNumuFromPion(), ana::BENDecomp::MakeWeightsNumuFromPion(), ana::FluxDecomp::NueEstimateFromPi(), ana::BENDecomp::NueEstimateFromPi(), Save(), ana::FluxDecomp::SavePlotsPi(), ana::BENDecomp::SavePlotsPi(), ana::Multiverse::ToHist(), ana::CutOptimization::ToHist(), ana::ICrossSectionAnalysis::ToHist(), ToTHX(), and Uninitialized().

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  }
const std::vector< Binning > & GetBinnings() const
Definition: Spectrum.h:264
Divide bin contents by bin widths.
Definition: UtilsExt.h:31
unsigned int NDimensions() const
Definition: Spectrum.h:262
OStream cout
Definition: OStream.cxx:6
const std::vector< std::string > & GetLabels() const
Definition: Spectrum.h:263
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
TH1 * ana::Spectrum::ToTHX ( double  exposure,
EExposureType  expotype = kPOT,
EBinType  bintype = kBinContent 
) const

Definition at line 211 of file Spectrum.cxx.

References NDimensions(), ToTH1(), ToTH2(), and ToTH3().

Referenced by Uninitialized().

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  }
TH2 * ToTH2(double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Spectrum must be 2D to obtain TH2.
Definition: Spectrum.cxx:165
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
unsigned int NDimensions() const
Definition: Spectrum.h:262
TH3 * ToTH3(double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Spectrum must be 3D to obtain TH3.
Definition: Spectrum.cxx:187
static Spectrum ana::Spectrum::Uninitialized ( )
inlinestatic

The only valid thing to do with such a spectrum is to assign something else into it.

Definition at line 145 of file Spectrum.h.

References col, Fill(), ana::kBinContent, ana::kPOT, operator=(), Spectrum(), ToTH1(), ToTH2(), ToTH3(), ToTHX(), w, submit_syst::x, and ~Spectrum().

Referenced by ana::PredictionAddRock::RockComponent().

145 {return Spectrum();}
Spectrum()
Constructor for Uninitialized()
Definition: Spectrum.h:268

Friends And Related Function Documentation

friend class Ratio
friend

Definition at line 46 of file Spectrum.h.

friend class SpectrumLoaderBase
friend

Definition at line 43 of file Spectrum.h.

friend class SpectrumSink
friend

Definition at line 44 of file Spectrum.h.

friend class SpectrumSinkBase< Spectrum >
friend

Definition at line 45 of file Spectrum.h.

Member Data Documentation

LabelsAndBins ana::Spectrum::fAxis
protected

Definition at line 291 of file Spectrum.h.

Referenced by Fill(), GetBinnings(), GetLabels(), Mean(), NDimensions(), operator=(), SaveTo(), Spectrum(), and ToTH1().

Hist ana::Spectrum::fHist
protected
double ana::Spectrum::fLivetime
protected
double ana::Spectrum::fPOT
protected
std::set<Spectrum**> ana::Spectrum::fReferences
protected

Things that point at this Spectrum. Maintained by SpectrumLoader.

Definition at line 289 of file Spectrum.h.

Referenced by AddLoader(), operator=(), RemoveLoader(), Spectrum(), and ~Spectrum().


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