ana::ReweightableSpectrum Class Reference

Spectrum with the value of a second variable, allowing for reweighting More...

#include <ReweightableSpectrum.h>

Inheritance diagram for ana::ReweightableSpectrum:

Public Member Functions
	ReweightableSpectrum (SpectrumLoaderBase &loader, const HistAxis &recoAxis, const HistAxis &trueAxis, const Cut &cut, const SystShifts &shift=kNoShift, const Var &wei=kUnweighted)
	ReweightableSpectrum (SpectrumLoaderBase &loader, const HistAxis &recoAxis, const HistAxis &trueAxis, const SpillCut &spillcut, const SliceCut &slicecut, const SystShifts &shift=kNoShift, const Var &wei=kUnweighted)
	ReweightableSpectrum (const Var &rwVar, const std::string &xlabel, const std::string &ylabel, double pot, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax)
	ReweightableSpectrum (const Var &rwVar, TH2 *h, const std::vector< std::string > &labels, const std::vector< Binning > &bins, double pot, double livetime)
	ReweightableSpectrum (const Var &rwVar, std::unique_ptr< TH2D > h, const std::vector< std::string > &labels, const std::vector< Binning > &bins, double pot, double livetime)
virtual	~ReweightableSpectrum ()
	ReweightableSpectrum (const ReweightableSpectrum &rhs)
ReweightableSpectrum &	operator= (const ReweightableSpectrum &rhs)
const Var &	ReweightVar () const
	The variable that will be used to fill the y-axis. More...
void	Fill (double x, double y, double w=1)
TH2D *	ToTH2 (double pot) const
Spectrum	UnWeighted () const
Spectrum	WeightingVariable () const
Spectrum	WeightedBy (const TH1 *weights) const
void	ReweightToTrueSpectrum (const Spectrum &target)
	Rescale bins so that WeightingVariable will return target. More...
void	ReweightToRecoSpectrum (const Spectrum &target)
	Recale bins so that Unweighted will return target. More...
void	Clear ()
void	SaveTo (TDirectory *dir) const
unsigned int	NDimensions () const
std::vector< std::string >	GetLabels () const
std::vector< Binning >	GetBinnings () const

Static Public Member Functions
static std::unique_ptr < ReweightableSpectrum >	LoadFrom (TDirectory *dir)

Protected Member Functions
	ReweightableSpectrum (const std::vector< std::string > &labels, const std::vector< Binning > &bins, const Var &rwVar)
	ReweightableSpectrum (const std::string &label, const Binning &bins, const Var &rwVar)
	ReweightableSpectrum (TH2 *h, const std::vector< std::string > &labels, const std::vector< Binning > &bins, double pot, double livetime)
void	RemoveLoader (SpectrumLoaderBase *)
void	AddLoader (SpectrumLoaderBase *)

Protected Attributes
Var	fRWVar
	What goes on the y axis? More...
TH2D *	fHist
double	fPOT
double	fLivetime
std::vector< std::string >	fLabels
std::vector< Binning >	fBins
std::string	fTrueLabel
std::set< SpectrumLoaderBase * >	fLoaderCount
	This count is maintained by SpectrumLoader, as a sanity check. More...

Friends
class	SpectrumLoaderBase
class	SpectrumLoader
class	NullLoader
class	MRCCLoader

Detailed Description

Spectrum with the value of a second variable, allowing for reweighting

Definition at line 14 of file ReweightableSpectrum.h.

Constructor & Destructor Documentation

ana::ReweightableSpectrum::ReweightableSpectrum	(	SpectrumLoaderBase &	loader,
		const HistAxis &	recoAxis,
		const HistAxis &	trueAxis,
		const Cut &	cut,
		const SystShifts &	shift = kNoShift,
		const Var &	wei = kUnweighted
	)

Definition at line 21 of file ReweightableSpectrum.cxx.

    : ReweightableSpectrum(loader, recoAxis, trueAxis, kNoSpillCut, cut, shift, wei)
  {
  }

ana::ReweightableSpectrum::ReweightableSpectrum	(	SpectrumLoaderBase &	loader,
		const HistAxis &	recoAxis,
		const HistAxis &	trueAxis,
		const SpillCut &	spillcut,
		const SliceCut &	slicecut,
		const SystShifts &	shift = kNoShift,
		const Var &	wei = kUnweighted
	)

Definition at line 32 of file ReweightableSpectrum.cxx.

    : ReweightableSpectrum(recoAxis.GetLabels(), recoAxis.GetBinnings(),
                           trueAxis.GetVars()[0])
  {
    assert(trueAxis.NDimensions() == 1);

    fTrueLabel = trueAxis.GetLabels()[0];

    DontAddDirectory guard;

    // Can't use HistCache here because y-axis is not necessarily
    // TrueEnergyBinning. TODO - that should maybe be generalized.

    const std::string name = UniqueName();

    Binning xbins = fBins[0];
    if(fBins.size() > 1){
      int n = 1;
      for(const Binning& b: fBins) n *= b.NBins();
      xbins = Binning::Simple(n, 0, n);
    }

    const Binning ybins = trueAxis.GetBinnings()[0];


    // Ugh combinatorics
    if(xbins.IsSimple() && ybins.IsSimple())
      fHist =  new TH2D(name.c_str(), "",
                        xbins.NBins(), xbins.Min(), xbins.Max(),
                        ybins.NBins(), ybins.Min(), ybins.Max());

    if(xbins.IsSimple() && !ybins.IsSimple())
      fHist =  new TH2D(name.c_str(), "",
                        xbins.NBins(), xbins.Min(), xbins.Max(),
                        ybins.NBins(), &ybins.Edges()[0]);

    if(!xbins.IsSimple() && ybins.IsSimple())
      fHist =  new TH2D(name.c_str(), "",
                        xbins.NBins(), &xbins.Edges()[0],
                        ybins.NBins(), ybins.Min(), ybins.Max());

    if(!xbins.IsSimple() && !ybins.IsSimple())
      fHist =  new TH2D(name.c_str(), "",
                        xbins.NBins(), &xbins.Edges()[0],
                        ybins.NBins(), &ybins.Edges()[0]);

    loader.AddReweightableSpectrum(*this, recoAxis.GetMultiDVar(), spillcut, slicecut, shift, wei);
  }

ana::ReweightableSpectrum::ReweightableSpectrum	(	const Var &	rwVar,
		const std::string &	xlabel,
		const std::string &	ylabel,
		double	pot,
		int	nbinsx,
		double	xmin,
		double	xmax,
		int	nbinsy,
		double	ymin,
		double	ymax
	)

Definition at line 88 of file ReweightableSpectrum.cxx.

    : ReweightableSpectrum(xlabel,
                           Binning::Simple(nbinsx, xmin, xmax),
                           rwVar)
  {
    DontAddDirectory guard;

    fHist = new TH2D(UniqueName().c_str(), "",
                     nbinsx, xmin, xmax, nbinsy, ymin, ymax);

    fTrueLabel = ylabel;

    // Ensure errors get accumulated properly
    fHist->Sumw2();
  }

ana::ReweightableSpectrum::ReweightableSpectrum	(	const Var &	rwVar,
		TH2 *	h,
		const std::vector< std::string > &	labels,
		const std::vector< Binning > &	bins,
		double	pot,
		double	livetime
	)

Definition at line 109 of file ReweightableSpectrum.cxx.

    : ReweightableSpectrum(labels, bins, rwVar)
  {
    fPOT = pot;
    fLivetime = livetime;

    if(!h){
      fHist = 0;
      return;
    }

    DontAddDirectory guard;

    const TString className = h->ClassName();

    if(className == "TH2D"){
      // Shortcut if types match
      fHist = new TH2D(*((TH2D*)h));
    }
    else{
      const TAxis* ax = h->GetXaxis();
      const TAxis* ay = h->GetYaxis();
      // Must have both or neither
      assert(bool(ax->GetXbins()->GetArray()) == bool(ay->GetXbins()->GetArray()));

      if(ax->GetXbins()->GetArray()){
        fHist = new TH2D(UniqueName().c_str(), "",
                         ax->GetNbins(), ax->GetXbins()->GetArray(),
                         ay->GetNbins(), ay->GetXbins()->GetArray());
      }
      else{
        fHist = new TH2D(UniqueName().c_str(), "",
                         ax->GetNbins(), ax->GetXmin(), ax->GetXmax(),
                         ay->GetNbins(), ay->GetXmin(), ay->GetXmax());
      }

      fHist->Add(h);
    }

    fTrueLabel = h->GetYaxis()->GetTitle();
  }

ana::ReweightableSpectrum::ReweightableSpectrum	(	const Var &	rwVar,
		std::unique_ptr< TH2D >	h,
		const std::vector< std::string > &	labels,
		const std::vector< Binning > &	bins,
		double	pot,
		double	livetime
	)

Definition at line 155 of file ReweightableSpectrum.cxx.

    : ReweightableSpectrum(labels, bins, rwVar)
  {
    fHist = h.release();
    fPOT = pot;
    fLivetime = livetime;

    fTrueLabel = fHist->GetYaxis()->GetTitle();
  }

ana::ReweightableSpectrum::~ReweightableSpectrum ( )

virtual

Definition at line 170 of file ReweightableSpectrum.cxx.

  {
    if(fHist && fHist->GetDirectory()){
      static bool once = true;
      if(once){
        once = false;
        std::cerr << "ReweightableSpectrum's fHist is associated with a directory. How did that happen?" << std::endl;
      }
    }

    // Can't use HistCache here because that expects all 2D spectra have the
    // usual true energy binning on their y-axis. In the case that we're
    // actually an OscillatableSpectrum that destructor will indeed take care
    // of invoking HistCache.
    delete fHist;
  }

ana::ReweightableSpectrum::ReweightableSpectrum

(

const ReweightableSpectrum &

rhs

)

Definition at line 188 of file ReweightableSpectrum.cxx.

    : fRWVar(rhs.fRWVar), fLabels(rhs.fLabels), fBins(rhs.fBins)
  {
    DontAddDirectory guard;

    fHist = new TH2D(*rhs.fHist);

    fPOT = rhs.fPOT;
    fLivetime = rhs.fLivetime;

    assert( rhs.fLoaderCount.empty() ); // Copying with pending loads is unexpected
  }

ana::ReweightableSpectrum::ReweightableSpectrum ( const std::vector< std::string > &  labels, const std::vector< Binning > &  bins, const Var &  rwVar  )

inlineprotected

Definition at line 97 of file ReweightableSpectrum.h.

      : fRWVar(rwVar),
        fHist(0), fPOT(0), fLivetime(0),
        fLabels(labels), fBins(bins)
    {
    }

ana::ReweightableSpectrum::ReweightableSpectrum ( const std::string &  label, const Binning &  bins, const Var &  rwVar  )

inlineprotected

Definition at line 106 of file ReweightableSpectrum.h.

      : fRWVar(rwVar),
        fHist(0), fPOT(0), fLivetime(0),
        fLabels(1, label), fBins(1, bins)
    {
    }

ana::ReweightableSpectrum::ReweightableSpectrum ( TH2 *  h, const std::vector< std::string > &  labels, const std::vector< Binning > &  bins, double  pot, double  livetime  )

inlineprotected

Constructor needed by LoadFrom. Since there's no good way to store a Var, ReweightVar will return nonsense for ReweightableSpectrum that are loaded from a file

Definition at line 118 of file ReweightableSpectrum.h.

      : ReweightableSpectrum(kUnweighted, h, labels, bins, pot, livetime)
    {
    }

Member Function Documentation

void ana::ReweightableSpectrum::AddLoader ( SpectrumLoaderBase *  p )

protected

Definition at line 456 of file ReweightableSpectrum.cxx.

  { fLoaderCount.insert(p); }

void ana::ReweightableSpectrum::Clear

(

)

Definition at line 446 of file ReweightableSpectrum.cxx.

  {
    fHist->Reset();
  }

void ana::ReweightableSpectrum::Fill	(	double	x,
		double	y,
		double	w = 1
	)

Definition at line 249 of file ReweightableSpectrum.cxx.

  {
    fHist->Fill(x, y, w);
  }

std::vector<Binning> ana::ReweightableSpectrum::GetBinnings ( ) const

inline

Definition at line 93 of file ReweightableSpectrum.h.

{return fBins;}

std::vector<std::string> ana::ReweightableSpectrum::GetLabels ( ) const

inline

Definition at line 92 of file ReweightableSpectrum.h.

{return fLabels;}

std::unique_ptr< ReweightableSpectrum > ana::ReweightableSpectrum::LoadFrom ( TDirectory *  dir )

static

Definition at line 485 of file ReweightableSpectrum.cxx.

  {
    TObjString* tag = (TObjString*)dir->Get("type");
    assert(tag);
    assert(tag->GetString() == "ReweightableSpectrum");

    TH2D* spect = (TH2D*)dir->Get("hist");
    assert(spect);
    TH1* hPot = (TH1*)dir->Get("pot");
    assert(hPot);
    TH1* hLivetime = (TH1*)dir->Get("livetime");
    assert(hLivetime);

    std::vector<std::string> labels;
    std::vector<Binning> bins;

    for(int i = 0; ; ++i){
      TDirectory* subdir = dir->GetDirectory(TString::Format("bins%d", i));
      if(!subdir) break;
      bins.push_back(*Binning::LoadFrom(subdir));
      TObjString* label = (TObjString*)dir->Get(TString::Format("label%d", i));
      labels.push_back(label ? label->GetString().Data() : "");
    }

    if(bins.empty() && labels.empty()){
      // Must be an old file. Make an attempt at backwards compatibility.
      bins.push_back(Binning::FromTAxis(spect->GetXaxis()));
      labels.push_back(spect->GetXaxis()->GetTitle());
    }

    return std::make_unique<ReweightableSpectrum>(kUnweighted,
                                                  spect,
                                                  labels, bins,
                                                  hPot->GetBinContent(1),
                                                  hLivetime->GetBinContent(1));
  }

unsigned int ana::ReweightableSpectrum::NDimensions ( ) const

inline

Definition at line 91 of file ReweightableSpectrum.h.

{return fLabels.size();}

ReweightableSpectrum & ana::ReweightableSpectrum::operator=

(

const ReweightableSpectrum &

rhs

)

Definition at line 202 of file ReweightableSpectrum.cxx.

  {
    if(this == &rhs) return *this;

    DontAddDirectory guard;

    fRWVar = rhs.fRWVar;
    fLabels = rhs.fLabels;
    fBins = rhs.fBins;

    delete fHist;
    fHist = new TH2D(*rhs.fHist);
    fPOT = rhs.fPOT;
    fLivetime = rhs.fLivetime;

    assert( fLoaderCount.empty() ); // Copying with pending loads is unexpected

    return *this;
  }

void ana::ReweightableSpectrum::RemoveLoader ( SpectrumLoaderBase *  p )

protected

Definition at line 452 of file ReweightableSpectrum.cxx.

  { fLoaderCount.erase(p); }

void ana::ReweightableSpectrum::ReweightToRecoSpectrum

(

const Spectrum &

target

)

Recale bins so that Unweighted will return target.

Definition at line 410 of file ReweightableSpectrum.cxx.

  {
    // This is a big component of what extrapolations do, so it has been
    // optimized for speed

    Ratio corr(target, UnWeighted());
    std::unique_ptr<TH1D> hcorr(corr.ToTH1());

    assert(hcorr->GetNbinsX() == fHist->GetNbinsX());

    const int X = fHist->GetNbinsX();
    const int Y = fHist->GetNbinsY();

    // Direct access to the bins is faster
    double* histArr = fHist->GetArray();
    double* corrArr = hcorr->GetArray();

    int bin = 0;
    for(int y = 0; y < Y+2; ++y){
      for(int x = 0; x < X+2; ++x){
        // Our loops go over the bins in the order they are internally in
        // fHist, and we do overflows, so we keep up exactly. If you get
        // paranoid, reenable this briefly.

        // assert(bin == fHist->GetBin(x, y));

        histArr[bin] *= corrArr[x];
        ++bin;
      }
    }

    TH1D* todel = hcorr.release();
    HistCache::Delete(todel);
  }

void ana::ReweightableSpectrum::ReweightToTrueSpectrum

(

const Spectrum &

target

)

Rescale bins so that WeightingVariable will return target.

Definition at line 373 of file ReweightableSpectrum.cxx.

  {
    // This is a big component of what extrapolations do, so it has been
    // optimized for speed

    Ratio corr(target, WeightingVariable());
    std::unique_ptr<TH1D> hcorr(corr.ToTH1());

    assert(hcorr->GetNbinsX() == fHist->GetNbinsY());

    const int X = fHist->GetNbinsX();
    const int Y = fHist->GetNbinsY();

    // Direct access to the bins is faster
    double* histArr = fHist->GetArray();
    double* corrArr = hcorr->GetArray();

    int bin = 0;
    for(int y = 0; y < Y+2; ++y){
      const double w = corrArr[y];
      for(int x = 0; x < X+2; ++x){
        // Our loops go over the bins in the order they are internally in
        // fHist, and we do overflows, so we keep up exactly. If you get
        // paranoid, reenable this briefly.

        // assert(bin == fHist->GetBin(x, y));

        histArr[bin] *= w;
        ++bin;
      }
    }

    TH1D* todel = hcorr.release();
    HistCache::Delete(todel);
  }

const Var& ana::ReweightableSpectrum::ReweightVar ( ) const

inline

The variable that will be used to fill the y-axis.

By convention, return zero if the information can't be obtained, and this event will be skipped.

Definition at line 64 of file ReweightableSpectrum.h.

{return fRWVar;}

void ana::ReweightableSpectrum::SaveTo

(

TDirectory *

dir

)

const

Function to save a ReweightableSpectrum to file the fRWVar member is not written to file, so when the spectrum is loaded back from file, ReweightVar should not be accessed, but reweighting still works

Definition at line 460 of file ReweightableSpectrum.cxx.

  {
    TDirectory* tmp = gDirectory;
    dir->cd();

    TObjString("ReweightableSpectrum").Write("type");

    fHist->GetYaxis()->SetTitle(fTrueLabel.c_str());
    fHist->Write("hist");
    TH1D hPot("", "", 1, 0, 1);
    hPot.Fill(.5, fPOT);
    hPot.Write("pot");
    TH1D hLivetime("", "", 1, 0, 1);
    hLivetime.Fill(.5, fLivetime);
    hLivetime.Write("livetime");

    for(unsigned int i = 0; i < fBins.size(); ++i){
      TObjString(fLabels[i].c_str()).Write(TString::Format("label%d", i).Data());
      fBins[i].SaveTo(dir->mkdir(TString::Format("bins%d", i)));
    }

    tmp->cd();
  }

TH2D * ana::ReweightableSpectrum::ToTH2

(

double

pot

)

const

Definition at line 223 of file ReweightableSpectrum.cxx.

  {
    // Could have a file temporarily open
    DontAddDirectory guard;

    TH2D* ret = new TH2D(*fHist);
    ret->SetName(UniqueName().c_str());
    if(fPOT){
      ret->Scale(pot/fPOT);
    }
    else{
      // How did it get events with no POT?
      assert(ret->Integral() == 0);
    }


    std::string label;
    for(const std::string& l: fLabels) label += l + " and ";
    label.resize(label.size()-5); // drop the last "and"
    ret->GetXaxis()->SetTitle(label.c_str());
    ret->GetYaxis()->SetTitle(fTrueLabel.c_str());

    return ret;
  }

Spectrum ana::ReweightableSpectrum::UnWeighted

(

)

const

Definition at line 281 of file ReweightableSpectrum.cxx.

  {
    DontAddDirectory guard;

    // Create a suitably-sized space for the result
    std::unique_ptr<TH1D> h(HistCache::New("", fHist->GetXaxis()));

    ProjectionX(fHist, h.get());

    return Spectrum(std::move(h), fLabels, fBins, fPOT, fLivetime);
  }

Spectrum ana::ReweightableSpectrum::WeightedBy

(

const TH1 *

weights

)

const

Definition at line 333 of file ReweightableSpectrum.cxx.

  {
    // This function is in the inner loop of oscillation fits, so some
    // optimization has been done.

    DontAddDirectory guard;

    assert(ws->GetNbinsX() == fHist->GetNbinsY());

    TAxis* ax = fHist->GetXaxis();
    TH1D* hRet = HistCache::New("", ax);

    const int X = fHist->GetNbinsX();
    const int Y = fHist->GetNbinsY();

    // Direct access to the bins is faster
    double* retArr = hRet->GetArray();
    double* histArr = fHist->GetArray();

    int bin = 0;
    for(int y = 0; y < Y+2; ++y){
      const double w = ws->GetBinContent(y);
      for(int x = 0; x < X+2; ++x){
        // Our loops go over the bins in the order they are internally in
        // fHist, and we do overflows, so we keep up exactly. If you get
        // paranoid, reenable this briefly.

        // assert(bin == fHist->GetBin(x, y));

        retArr[x] += histArr[bin]*w;
        ++bin;
      }
    }

    // TODO: can this all be more efficient?
    return Spectrum(std::unique_ptr<TH1D>(hRet), fLabels, fBins, fPOT, fLivetime);
  }

Spectrum ana::ReweightableSpectrum::WeightingVariable

(

)

const

Definition at line 320 of file ReweightableSpectrum.cxx.

  {
    DontAddDirectory guard;

    // Create a suitably-sized space for the result
    std::unique_ptr<TH1D> h(HistCache::New("", fHist->GetYaxis()));

    ProjectionY(fHist, h.get());

    return Spectrum(std::move(h), fLabels, fBins, fPOT, fLivetime);
  }

Friends And Related Function Documentation

friend class MRCCLoader

friend

Definition at line 20 of file ReweightableSpectrum.h.

friend class NullLoader

friend

Definition at line 19 of file ReweightableSpectrum.h.

friend class SpectrumLoader

friend

Definition at line 18 of file ReweightableSpectrum.h.

friend class SpectrumLoaderBase

friend

Definition at line 17 of file ReweightableSpectrum.h.

Member Data Documentation

std::vector<Binning> ana::ReweightableSpectrum::fBins

protected

Definition at line 136 of file ReweightableSpectrum.h.

TH2D* ana::ReweightableSpectrum::fHist

protected

Definition at line 131 of file ReweightableSpectrum.h.

std::vector<std::string> ana::ReweightableSpectrum::fLabels

protected

Definition at line 135 of file ReweightableSpectrum.h.

double ana::ReweightableSpectrum::fLivetime

protected

Definition at line 133 of file ReweightableSpectrum.h.

std::set<SpectrumLoaderBase*> ana::ReweightableSpectrum::fLoaderCount

protected

This count is maintained by SpectrumLoader, as a sanity check.

Definition at line 141 of file ReweightableSpectrum.h.

double ana::ReweightableSpectrum::fPOT

protected

Definition at line 132 of file ReweightableSpectrum.h.

Var ana::ReweightableSpectrum::fRWVar

protected

What goes on the y axis?

Definition at line 129 of file ReweightableSpectrum.h.

std::string ana::ReweightableSpectrum::fTrueLabel

protected

Definition at line 138 of file ReweightableSpectrum.h.

---

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

• ReweightableSpectrum.h
• ReweightableSpectrum.cxx