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

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N20-12-01/CAFAna/XSec/GenieMultiverseSyst.h"

Inheritance diagram for ana::GenieMultiverseNormalizedSpectra:
ana::GenieMultiverseSpectra

Public Types

enum  BandOptions {
  kNoBand, kBandFromMean, kBandFromNominal, kBandFromMPV,
  kBandFromMedian
}
 

Public Member Functions

 GenieMultiverseNormalizedSpectra (unsigned int nuniverses, SpectrumLoaderBase &loader, const HistAxis &axis, const Cut &cut, const SystShifts &shift, const Var &wei=kUnweighted, std::string config_pathname="knob_config.txt")
 
 GenieMultiverseNormalizedSpectra (unsigned int nuniverses, SpectrumLoaderBase &loader, const HistAxis &axis, const Cut &cut, const SystShifts &shift, const Var &wei, std::vector< const ISyst * > systs)
 
 GenieMultiverseNormalizedSpectra (unsigned int nuniverses, std::vector< std::unique_ptr< ana::Spectrum >> &spectra, bool fromfile=false)
 
const SpectrumUpperSigma (BandOptions opt=kBandFromNominal)
 
const SpectrumLowerSigma (BandOptions opt=kBandFromNominal)
 
void SaveTo (TDirectory *dir, const std::string &name)
 Save the Spectrum structure of all universe to a TFile directory. More...
 
TH1 * ToTH1 () const
 
TH1 * ToTH1 (double pot) const
 
TGraphAsymmErrors * ToAreaNormalizedTH1 (double pot, int col=-1, int errCol=-1, float headroom=1.3, bool newaxis=true) const
 
TH2 * ToTH2 () const
 
TH2 * ToTH2 (double pot) const
 
const SpectrumNominal () const
 return the nominal universe More...
 
std::vector< const Spectrum * > AllUniverses () const
 return all universes More...
 
const SpectrumUpperRMS () const
 return RMS extremes More...
 
const SpectrumLowerRMS () const
 
const SpectrumMean () const
 
const SpectrumUpperSigma (BandOptions opt=kBandFromNominal) const
 
const SpectrumLowerSigma (BandOptions opt=kBandFromNominal) const
 
void SaveTo (TDirectory *dir, const std::string &name) const
 Save the Spectrum structure of all universe to a TFile directory. More...
 
unsigned int GetNUniverses ()
 Get the number of universes in an object of this class. More...
 
const auto & GetShiftTable () const
 

Static Public Member Functions

static std::unique_ptr< GenieMultiverseNormalizedSpectraLoadFrom (TDirectory *dir, const std::string &name)
 Load the Spectrum structure of all universe from a TFile directory. More...
 
static std::pair< TH2 *, TH2 * > GetSigmaHistograms (std::vector< TH2 * >, BandOptions opt=kBandFromNominal)
 Calculate error band given 2D histograms. More...
 

Protected Member Functions

void BuildKnobConfigTable (std::string config_pathname)
 
SpectrumCreateUniverse (SpectrumLoaderBase &loader, const HistAxis &axis, const Cut &cut, std::vector< const ISyst * > systs, const Var &wei=kUnweighted)
 
SpectrumCreateUniverse (SpectrumLoaderBase &loader, const NuTruthHistAxis &axis, const NuTruthCut &cut, const SystShifts &shift=kNoShift, const NuTruthVar &wei=kNuTruthUnweighted)
 
SpectrumCreateUniverse (SpectrumLoaderBase &loader, const NuTruthHistAxis &axis, const NuTruthCut &cut, std::vector< const ISyst * > systs, const NuTruthVar &wei=kNuTruthUnweighted)
 

Protected Attributes

std::vector< Spectrum * > fSpectra
 
int fUniverseSeed = 1001
 
BandOptions fCurBandOpt = kNoBand
 
knob_sampling_mode_set fKnobConfigTable
 

Private Member Functions

SpectrumCreateUniverse (SpectrumLoaderBase &loader, const HistAxis &axis, const Cut &cut, const SystShifts &shift=kNoShift, const Var &wei=kUnweighted)
 
void NormalizeSpectra ()
 Normalize spectra. It has to be done AFTER SpectrumLoader::Go() is called. More...
 

Private Attributes

double fNominalCount
 
bool fSpectraNormalized = false
 A flag signaling whether the spectra are normalized already. More...
 

Detailed Description

Class for shape-only systematic effect.

Definition at line 289 of file GenieMultiverseSyst.h.

Member Enumeration Documentation

Constructor & Destructor Documentation

ana::GenieMultiverseNormalizedSpectra::GenieMultiverseNormalizedSpectra ( unsigned int  nuniverses,
SpectrumLoaderBase loader,
const HistAxis axis,
const Cut cut,
const SystShifts shift,
const Var wei = kUnweighted,
std::string  config_pathname = "knob_config.txt" 
)
ana::GenieMultiverseNormalizedSpectra::GenieMultiverseNormalizedSpectra ( unsigned int  nuniverses,
SpectrumLoaderBase loader,
const HistAxis axis,
const Cut cut,
const SystShifts shift,
const Var wei,
std::vector< const ISyst * >  systs 
)

Definition at line 807 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::CreateUniverse(), and MECModelEnuComparisons::i.

815  {
816  // universe with index 0 is the nominal universe
817  fSpectra.push_back(new Spectrum(loader, axis, cut, shift, wei));
818  for(unsigned int i = 0; i < nuniverses-1; i++)
820  }
std::vector< Spectrum * > fSpectra
std::vector< SystGroupDef > systs
Definition: syst_header.h:385
Spectrum * CreateUniverse(SpectrumLoaderBase &loader, const HistAxis &axis, const Cut &cut, const SystShifts &shift=kNoShift, const Var &wei=kUnweighted)
loader
Definition: demo0.py:10
std::vector< float > Spectrum
Definition: Constants.h:610
const Cut cut
Definition: exporter_fd.C:30
const int nuniverses
ana::GenieMultiverseNormalizedSpectra::GenieMultiverseNormalizedSpectra ( unsigned int  nuniverses,
std::vector< std::unique_ptr< ana::Spectrum >> &  spectra,
bool  fromfile = false 
)
inline

Member Function Documentation

std::vector<const Spectrum*> ana::GenieMultiverseSpectra::AllUniverses ( ) const
inlineinherited

return all universes

Definition at line 104 of file GenieMultiverseSyst.h.

References dir, ana::LoadFrom(), MECModelEnuComparisons::opt, ana::SaveTo(), and string.

Referenced by make_muonid_opt(), and make_vertex_optimiz().

104 {return {fSpectra.begin(), fSpectra.end()};}
std::vector< Spectrum * > fSpectra
void ana::GenieMultiverseSpectra::BuildKnobConfigTable ( std::string  config_pathname)
protectedinherited

function for building knob configuration table

Definition at line 422 of file GenieMultiverseSyst.cxx.

References om::cerr, allTimeWatchdog::endl, exit(), MakeMiniprodValidationCuts::f, ana::GenieMultiverseSpectra::fKnobConfigTable, and allTimeWatchdog::index.

Referenced by ana::GenieMultiverseSpectra::GenieMultiverseSpectra().

423  {
424  // Create an empty property tree object
425  //~ using boost::property_tree::ptree;
426  //~ ptree pt;
427 
428  // check file existence
429  ifstream f(fpn.c_str());
430  if(!f.good()){
431  cerr << "ERROR: GENIE knob configuration file \033[1;31m" << fpn << "\033[0m is not found." << endl;
432  exit(0);
433  }
434 
435  // read configuration from the json file
436  // This line is where the boost bug emerges.
437  //~ boost::property_tree::json_parser::read_json(fpn, pt);
438 
439  string knobname;
440  int sample_mode;
441  int index = 0;
442  while(f >> knobname >> sample_mode)
443  fKnobConfigTable.insert(knob_sampling_mode(index++,knobname,sample_mode));
444 
445  // print out key values
446  //~ for(auto& nameIt: fKnobConfigTable.get<0>()) cout << nameIt.id << " " << nameIt.name << " " << nameIt.mode << endl;
447  }
OStream cerr
Definition: OStream.cxx:7
knob_sampling_mode_set fKnobConfigTable
exit(0)
Spectrum * ana::GenieMultiverseSpectra::CreateUniverse ( SpectrumLoaderBase loader,
const HistAxis axis,
const Cut cut,
std::vector< const ISyst * >  systs,
const Var wei = kUnweighted 
)
protectedinherited

Function to create one universe

Definition at line 486 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::fShiftTable, ana::GenieMultiverseSpectra::fUniverseSeed, and central_limit::rand.

491  {
492  map<const ISyst*, double> knob_weight_table;
493 
494  for (const auto & syst: systs)
495  {
496  TRandom3 rand(fUniverseSeed++);
497  knob_weight_table[syst] = rand.Gaus();
498  }
499 
500  auto spec = new Spectrum(loader, axis, cut, knob_weight_table, wei);
501  fShiftTable.push_back(std::move(knob_weight_table));
502  return spec;
503  }
std::vector< SystGroupDef > systs
Definition: syst_header.h:385
std::vector< std::map< const ISyst *, double > > fShiftTable
loader
Definition: demo0.py:10
std::vector< float > Spectrum
Definition: Constants.h:610
const Cut cut
Definition: exporter_fd.C:30
Spectrum * ana::GenieMultiverseSpectra::CreateUniverse ( SpectrumLoaderBase loader,
const NuTruthHistAxis axis,
const NuTruthCut cut,
const SystShifts shift = kNoShift,
const NuTruthVar wei = kNuTruthUnweighted 
)
protectedinherited

Function to create one universe for NuTruth variables.

Definition at line 506 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::fKnobConfigTable, ana::GenieMultiverseSpectra::fUniverseSeed, ana::getAllXsecNuTruthSysts_2017(), ana::GetGenieKnobSyst(), central_limit::rand, and systs.

511  {
512  const vector<const ISyst*> systs = getAllXsecNuTruthSysts_2017();
513 
514  map<const ISyst*, double> knob_weight_table;
515  for(auto& idIt: fKnobConfigTable.get<0>())
516  {
517  TRandom3 rand(fUniverseSeed++);
518  rwgt::ReweightLabel_t knoblabel = (rwgt::ReweightLabel_t)idIt.id;
519  double nsigma = 0;
520  if(idIt.mode == 0) continue;
521  else if(idIt.mode == 1) nsigma = rand.Gaus();
522  else if(idIt.mode == 2) nsigma = rand.Integer(2);
523  // For tags later than 2017/08/14, when the underlying framework
524  // underwent an overhaul, this conditional check is needed.
525 
526  // Note that in GetGenieKnobSyst function in XSecSysts.cxx,
527  // rwgt::kReweightNull is excluded. Therefore a call with argument
528  // rwgt::kReweightNull will lead to seg fault. Skip it!
529  if(knoblabel == rwgt::kReweightNull) continue;
530  if(find(systs.begin(), systs.end(), GetGenieKnobSyst(knoblabel)) != systs.end())
531  knob_weight_table[GetGenieKnobSyst(knoblabel)] = nsigma;
532  }
533 
534  return new Spectrum(loader, axis, cut, SystShifts(knob_weight_table), wei);
535  }
std::vector< SystGroupDef > systs
Definition: syst_header.h:385
knob_sampling_mode_set fKnobConfigTable
loader
Definition: demo0.py:10
std::vector< float > Spectrum
Definition: Constants.h:610
const Cut cut
Definition: exporter_fd.C:30
std::vector< const ISyst * > getAllXsecNuTruthSysts_2017()
Get master XSec syst list for 2017 analyses (NuTruthSyst variant)
const NOvARwgtSyst * GetGenieKnobSyst(rwgt::ReweightLabel_t knobIdx, std::string altName, std::string altLabel)
Convenience function to get a GENIE knob syst. (Allows using the GENIE knob name & description as the...
Definition: XSecSysts.cxx:119
Spectrum * ana::GenieMultiverseSpectra::CreateUniverse ( SpectrumLoaderBase loader,
const NuTruthHistAxis axis,
const NuTruthCut cut,
std::vector< const ISyst * >  systs,
const NuTruthVar wei = kNuTruthUnweighted 
)
protectedinherited

Function to create one universe for NuTruth variables.

Definition at line 538 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::fUniverseSeed, ana::GetGenieKnobSyst(), ana::kReweightLabels(), and central_limit::rand.

543  {
544  map<const ISyst*, double> knob_weight_table;
545  for(auto knoblabel : ana::kReweightLabels)
546  {
547  TRandom3 rand(fUniverseSeed++);
548  // Note that in GetGenieKnobSyst function in XSecSysts.cxx,
549  // rwgt::kReweightNull is excluded. Therefore a call with argument
550  // rwgt::kReweightNull will lead to seg fault. Skip it!
551  if(knoblabel == rwgt::kReweightNull) continue;
552  if(find(systs.begin(), systs.end(), GetGenieKnobSyst(knoblabel)) != systs.end())
553  knob_weight_table[GetGenieKnobSyst(knoblabel)] = rand.Gaus();
554  }
555 
556  return new Spectrum(loader, axis, cut, SystShifts(knob_weight_table), wei);
557  }
std::vector< SystGroupDef > systs
Definition: syst_header.h:385
const std::list< rwgt::ReweightLabel_t > kReweightLabels({rwgt::fReweightMaNCEL, rwgt::fReweightEtaNCEL, rwgt::fReweightNormCCQE, rwgt::fReweightNormCCQEenu, rwgt::fReweightMaCCQEshape, rwgt::fReweightMaCCQE, rwgt::fReweightVecCCQEshape, rwgt::fReweightNormCCRES, rwgt::fReweightMaCCRESshape, rwgt::fReweightMvCCRESshape, rwgt::fReweightMaCCRES, rwgt::fReweightMvCCRES, rwgt::fReweightNormNCRES, rwgt::fReweightMaNCRESshape, rwgt::fReweightMvNCRESshape, rwgt::fReweightMaNCRES, rwgt::fReweightMvNCRES, rwgt::fReweightMaCOHpi, rwgt::fReweightR0COHpi, rwgt::fReweightRvpCC1pi, rwgt::fReweightRvpCC2pi, rwgt::fReweightRvpNC1pi, rwgt::fReweightRvpNC2pi, rwgt::fReweightRvnCC1pi, rwgt::fReweightRvnCC2pi, rwgt::fReweightRvnNC1pi, rwgt::fReweightRvnNC2pi, rwgt::fReweightRvbarpCC1pi, rwgt::fReweightRvbarpCC2pi, rwgt::fReweightRvbarpNC1pi, rwgt::fReweightRvbarpNC2pi, rwgt::fReweightRvbarnCC1pi, rwgt::fReweightRvbarnCC2pi, rwgt::fReweightRvbarnNC1pi, rwgt::fReweightRvbarnNC2pi, rwgt::fReweightAhtBY, rwgt::fReweightBhtBY, rwgt::fReweightCV1uBY, rwgt::fReweightCV2uBY, rwgt::fReweightAhtBYshape, rwgt::fReweightBhtBYshape, rwgt::fReweightCV1uBYshape, rwgt::fReweightCV2uBYshape, rwgt::fReweightNormDISCC, rwgt::fReweightRnubarnuCC, rwgt::fReweightDISNuclMod, rwgt::fReweightNC, rwgt::fReweightAGKY_xF1pi, rwgt::fReweightAGKY_pT1pi, rwgt::fReweightFormZone, rwgt::fReweightMFP_pi, rwgt::fReweightMFP_N, rwgt::fReweightFrCEx_pi, rwgt::fReweightFrInel_pi, rwgt::fReweightFrAbs_pi, rwgt::fReweightFrPiProd_pi, rwgt::fReweightFrCEx_N, rwgt::fReweightFrInel_N, rwgt::fReweightFrAbs_N, rwgt::fReweightFrPiProd_N, rwgt::fReweightCCQEPauliSupViaKF, rwgt::fReweightCCQEMomDistroFGtoSF, rwgt::fReweightBR1gamma, rwgt::fReweightBR1eta, rwgt::fReweightTheta_Delta2Npi, rwgt::fReweightZNormCCQE, rwgt::fReweightZExpA1CCQE, rwgt::fReweightZExpA2CCQE, rwgt::fReweightZExpA3CCQE, rwgt::fReweightZExpA4CCQE, rwgt::fReweightAxFFCCQEshape})
loader
Definition: demo0.py:10
std::vector< float > Spectrum
Definition: Constants.h:610
const Cut cut
Definition: exporter_fd.C:30
const NOvARwgtSyst * GetGenieKnobSyst(rwgt::ReweightLabel_t knobIdx, std::string altName, std::string altLabel)
Convenience function to get a GENIE knob syst. (Allows using the GENIE knob name & description as the...
Definition: XSecSysts.cxx:119
Spectrum * ana::GenieMultiverseNormalizedSpectra::CreateUniverse ( SpectrumLoaderBase loader,
const HistAxis axis,
const Cut cut,
const SystShifts shift = kNoShift,
const Var wei = kUnweighted 
)
private

Definition at line 822 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseParameters::fKnobConfigTable, ana::GenieMultiverseParameters::fUniverseSeed, ana::getAllXsecNuTruthSysts_2017(), ana::GetGenieKnobSyst(), central_limit::rand, and systs.

828  {
829  map<const ISyst*, double> knob_weight_table;
830  for(auto& idIt: fKnobConfigTable.get<0>())
831  {
832  TRandom3 rand(fUniverseSeed++);
833  rwgt::ReweightLabel_t knoblabel = (rwgt::ReweightLabel_t)idIt.id;
834  double nsigma = 0;
835  if(idIt.mode == 0) continue;
836  else if(idIt.mode == 1) nsigma = rand.Gaus();
837  else if(idIt.mode == 2) nsigma = rand.Integer(2);
838  // For tags later than 2017/08/14, when the underlying framework
839  // underwent an overhaul, this conditional check is needed.
840  vector<const ISyst*> systs = getAllXsecNuTruthSysts_2017();
841  if(find(systs.begin(), systs.end(), GetGenieKnobSyst(knoblabel)) != systs.end())
842  knob_weight_table[GetGenieKnobSyst(knoblabel)] = nsigma;
843  }
844  Spectrum* normSpec = new Spectrum(loader, axis, cut, knob_weight_table, wei);
845 
846  return normSpec;
847  }
std::vector< SystGroupDef > systs
Definition: syst_header.h:385
knob_sampling_mode_set fKnobConfigTable
loader
Definition: demo0.py:10
std::vector< float > Spectrum
Definition: Constants.h:610
const Cut cut
Definition: exporter_fd.C:30
std::vector< const ISyst * > getAllXsecNuTruthSysts_2017()
Get master XSec syst list for 2017 analyses (NuTruthSyst variant)
const NOvARwgtSyst * GetGenieKnobSyst(rwgt::ReweightLabel_t knobIdx, std::string altName, std::string altLabel)
Convenience function to get a GENIE knob syst. (Allows using the GENIE knob name & description as the...
Definition: XSecSysts.cxx:119
unsigned int ana::GenieMultiverseSpectra::GetNUniverses ( )
inlineinherited

Get the number of universes in an object of this class.

Definition at line 119 of file GenieMultiverseSyst.h.

References MECModelEnuComparisons::opt.

119 {return fNUniverses;};
unsigned int fNUniverses
number of universes to generate
const auto& ana::GenieMultiverseSpectra::GetShiftTable ( ) const
inlineinherited

Definition at line 124 of file GenieMultiverseSyst.h.

124 { return fShiftTable; };
std::vector< std::map< const ISyst *, double > > fShiftTable
pair< TH2 *, TH2 * > ana::GenieMultiverseSpectra::GetSigmaHistograms ( std::vector< TH2 * >  huniv,
BandOptions  opt = kBandFromNominal 
)
staticinherited

Calculate error band given 2D histograms.

Definition at line 185 of file GenieMultiverseSyst.cxx.

References a, ana::GenieMultiverseSpectra::BinSigma(), plot_validation_datamc::Clone(), MECModelEnuComparisons::i, ana::GenieMultiverseSpectra::kBandFromMean, ana::GenieMultiverseSpectra::kBandFromMedian, ana::GenieMultiverseSpectra::kBandFromMPV, ana::GenieMultiverseSpectra::kBandFromNominal, ana::GenieMultiverseSpectra::kNoBand, cet::sqlite::max(), extractScale::mean, min(), and nbinsx.

186  {
187  // create up and down histograms
188  TH2* hup = (TH2*)huniv[0]->Clone("hup");
189  TH2* hdown = (TH2*)huniv[0]->Clone("hdown");
190 
191  // loop through all bins
192  int nbinsx = hup->GetNbinsX();
193  int nbinsy = hup->GetNbinsY();
194  for(int ix = 0; ix <= nbinsx; ix++){
195  for(int iy = 0; iy <= nbinsy; iy++)
196  {
197  // collect counts from all universes
198  vector<float> counts_in_one_bin;
199  for(unsigned int iu = 0; iu < huniv.size(); iu++)
200  counts_in_one_bin.push_back(huniv[iu]->GetBinContent(ix, iy));
201 
202  // determine from where the error band is drawn
203  double pivot = 0.;
204  double mean = TMath::Mean(counts_in_one_bin.begin(), counts_in_one_bin.end());
205  if(opt == kBandFromMean) pivot = mean;
206  else if(opt == kBandFromNominal || opt == kNoBand) pivot = huniv[0]->GetBinContent(ix, iy);
207  else if(opt == kBandFromMedian){
208  float* a = &counts_in_one_bin[0];
209  pivot = TMath::Median(counts_in_one_bin.size(), a);
210  }
211  else if(opt == kBandFromMPV){
212  double a[counts_in_one_bin.size()];
213  for(unsigned int i = 0; i < counts_in_one_bin.size(); i++) a[i] = counts_in_one_bin[i];
214  double min = *min_element(counts_in_one_bin.begin(), counts_in_one_bin.end());
215  double max = *max_element(counts_in_one_bin.begin(), counts_in_one_bin.end());
216  TKDE kde(counts_in_one_bin.size(), a, min, max, "kUnbinned");
217  pivot = kde.GetFunction(counts_in_one_bin.size(), min, max)->GetMaximumX();
218  }
219 
220  // draw band
221  hup->SetBinContent(ix, iy, BinSigma(counts_in_one_bin, 1, pivot));
222  hdown->SetBinContent(ix, iy, BinSigma(counts_in_one_bin, -1, pivot));
223  }
224  }
225 
226  return pair<TH2*, TH2*>(hup, hdown);
227  }
static float BinSigma(std::vector< float > events, float nsigmas, float mean=-1)
const double a
static float min(const float a, const float b, const float c)
Definition: absgeo.cxx:45
Int_t nbinsx
Definition: plot.C:23
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
std::unique_ptr< GenieMultiverseNormalizedSpectra > ana::GenieMultiverseNormalizedSpectra::LoadFrom ( TDirectory *  dir,
const std::string name 
)
static

Load the Spectrum structure of all universe from a TFile directory.

Definition at line 849 of file GenieMultiverseSyst.cxx.

References ana::assert(), dir, MECModelEnuComparisons::i, and art::to_string().

850  {
851  dir = dir->GetDirectory(name.c_str()); // switch to subdir
852  assert(dir);
853 
854  TObjString* universes = (TObjString*)dir->Get("nUniverses");
855  unsigned int nUniverses = universes->GetString().Atoi();
856 
857  std::vector<std::unique_ptr<ana::Spectrum>> spectra;
858  for(unsigned int i = 0; i < nUniverses; i++){
859  spectra.push_back(ana::LoadFrom<Spectrum>(dir, "universe_"+std::to_string(i)));
860  }
861 
862  delete dir;
863 
864  return std::make_unique<GenieMultiverseNormalizedSpectra>(nUniverses, spectra, true);
865  }
const XML_Char * name
Definition: expat.h:151
TDirectory * dir
Definition: macro.C:5
assert(nhit_max >=nhit_nbins)
std::string to_string(ModuleType mt)
Definition: ModuleType.h:32
const Spectrum * ana::GenieMultiverseSpectra::LowerRMS ( ) const
inherited
const Spectrum * ana::GenieMultiverseSpectra::LowerSigma ( BandOptions  opt = kBandFromNominal) const
inherited
const Spectrum * ana::GenieMultiverseNormalizedSpectra::LowerSigma ( BandOptions  opt = kBandFromNominal)

Definition at line 867 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::LowerSigma().

868  {
871  }
bool fSpectraNormalized
A flag signaling whether the spectra are normalized already.
const Spectrum * LowerSigma(BandOptions opt=kBandFromNominal) const
void NormalizeSpectra()
Normalize spectra. It has to be done AFTER SpectrumLoader::Go() is called.
const Spectrum * ana::GenieMultiverseSpectra::Mean ( ) const
inherited
const Spectrum * ana::GenieMultiverseSpectra::Nominal ( ) const
inherited

return the nominal universe

Definition at line 244 of file GenieMultiverseSyst.cxx.

References om::cerr, e, allTimeWatchdog::endl, ana::GenieMultiverseSpectra::fSpectra, and stan::math::out_of_range().

Referenced by make_muonid_opt(), make_vertex_optimiz(), MakePlots(), ana::GenieMultiverseSpectra::ToAreaNormalizedTH1(), ana::GenieMultiverseSpectra::ToTH1(), ana::GenieMultiverseSpectra::ToTH2(), xsec_tot_uncert_optimization(), and xsec_uncertainty_per_bin().

245  {
246  Spectrum* nominal_u = 0;
247  try {
248  nominal_u = fSpectra.at(0);
249  } catch(const std::out_of_range& e) {
250  cerr << "Nominal universe is accessed before existence." << endl;
251  }
252 
253  return nominal_u;
254  }
std::vector< Spectrum * > fSpectra
OStream cerr
Definition: OStream.cxx:7
std::vector< float > Spectrum
Definition: Constants.h:610
void out_of_range(const char *function, int max, int index, const char *msg1="", const char *msg2="")
Float_t e
Definition: plot.C:35
void ana::GenieMultiverseNormalizedSpectra::NormalizeSpectra ( )
private

Normalize spectra. It has to be done AFTER SpectrumLoader::Go() is called.

Definition at line 873 of file GenieMultiverseSyst.cxx.

References MECModelEnuComparisons::i, ana::Spectrum::Integral(), POT, ana::Spectrum::POT(), scale, and ana::Spectrum::Scale().

874  {
875  fNominalCount = fSpectra[0]->Integral(fSpectra[0]->POT());
876  for(unsigned int i = 1; i < fSpectra.size(); i++)
877  {
878  Spectrum* normSpec = fSpectra[i];
879  // Scale this spectrum to the area of the nominal spectrum.
880  double scale = fNominalCount/normSpec->Integral(normSpec->POT());
881  normSpec->Scale(scale);
882  }
883  fSpectraNormalized = true;
884  }
std::vector< Spectrum * > fSpectra
bool fSpectraNormalized
A flag signaling whether the spectra are normalized already.
Double_t scale
Definition: plot.C:25
std::vector< float > Spectrum
Definition: Constants.h:610
std::vector< double > POT
void ana::GenieMultiverseSpectra::SaveTo ( TDirectory *  dir,
const std::string name 
) const
inherited

Save the Spectrum structure of all universe to a TFile directory.

Definition at line 655 of file GenieMultiverseSyst.cxx.

References confusionMatrixTree::count, dir, ana::GenieMultiverseSpectra::fSpectra, it, tmp, and art::to_string().

Referenced by SaveTo().

655  {
656  TDirectory *tmp = gDirectory;
657 
658  dir = dir->mkdir(name.c_str()); // switch to subdir
659  dir->cd();
660 
661  //~ TVectorD universes(1);
662  //~ universes[0] = {(Double_t)fSpectra.size()};
663  //~ universes.Write("nUniverses");
664  // Could not get the above code to store a resonable unsigned int value,
665  // which leads to failure when executing hadd_cafana.
666  // Among the hadd_cafana data types, TObjString types are also not added.
667  // Gonna try this.
668  TObjString universes(Form("%lu",fSpectra.size()));
669  universes.Write("nUniverses");
670 
671  int count = 0;
672  for(auto& it: fSpectra){
673  it->SaveTo(dir, "universe_"+std::to_string(count));
674  count++;
675  }// end loop over spectra
676 
677  dir->Write();
678  delete dir;
679 
680  tmp->cd();
681  }// end of SaveTo
std::vector< Spectrum * > fSpectra
const XML_Char * name
Definition: expat.h:151
set< int >::iterator it
Float_t tmp
Definition: plot.C:36
TDirectory * dir
Definition: macro.C:5
std::string to_string(ModuleType mt)
Definition: ModuleType.h:32
void ana::GenieMultiverseNormalizedSpectra::SaveTo ( TDirectory *  dir,
const std::string name 
)

Save the Spectrum structure of all universe to a TFile directory.

Definition at line 886 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::SaveTo().

887  {
890  }
const XML_Char * name
Definition: expat.h:151
bool fSpectraNormalized
A flag signaling whether the spectra are normalized already.
void SaveTo(TDirectory *dir, const std::string &name) const
Save the Spectrum structure of all universe to a TFile directory.
TDirectory * dir
Definition: macro.C:5
void NormalizeSpectra()
Normalize spectra. It has to be done AFTER SpectrumLoader::Go() is called.
TGraphAsymmErrors * ana::GenieMultiverseSpectra::ToAreaNormalizedTH1 ( double  pot,
int  col = -1,
int  errCol = -1,
float  headroom = 1.3,
bool  newaxis = true 
) const
inherited

Definition at line 256 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::fLowerRMS, ana::GenieMultiverseSpectra::fUpperRMS, MECModelEnuComparisons::g, hi(), kRed, lo(), ana::GenieMultiverseSpectra::Nominal(), std::sqrt(), std::swap(), ana::Spectrum::ToTH1(), w, and submit_syst::y.

258  {
259  if(col == -1){
260  col = kRed;
261  errCol = kRed-10;
262  }
263  else if(errCol == -1) errCol = col-7; // hopefully a lighter version
264 
265  TH1* nom = Nominal()->ToTH1(pot);
266 
267  std::vector<TH1*> ups, dns;
268  ups.push_back(fUpperRMS->ToTH1(pot));
269  dns.push_back(fLowerRMS->ToTH1(pot));
270 
271  nom->SetLineColor(col);
272  nom->GetXaxis()->CenterTitle();
273  nom->GetYaxis()->CenterTitle();
274  //~ if(newaxis) nom->Draw("hist ]["); // Set the axes up
275 
276  double yMax = nom->GetBinContent(nom->GetMaximumBin());
277 
278  TGraphAsymmErrors* g = new TGraphAsymmErrors;
279 
280  for(int binIdx = 0; binIdx < nom->GetNbinsX()+2; ++binIdx){
281  const double y = nom->GetBinContent(binIdx);
282  g->SetPoint(binIdx, nom->GetXaxis()->GetBinCenter(binIdx), y);
283 
284  const double w = nom->GetXaxis()->GetBinWidth(binIdx);
285 
286  double errUp = 0;
287  double errDn = 0;
288 
289  for(unsigned int systIdx = 0; systIdx < ups.size(); ++systIdx){
290  double hi = ups[systIdx]->GetBinContent(binIdx)-y;
291  double lo = y-dns[systIdx]->GetBinContent(binIdx);
292 
293  if(lo <= 0 && hi <= 0) std::swap(lo, hi);
294 
295  errUp += hi*hi;
296  errDn += lo*lo;
297 
298  // TODO: what happens if they're both high or both low?
299  } // end for systIdx
300 
301 
302  g->SetPointError(binIdx, w/2, w/2, sqrt(errDn), sqrt(errUp));
303  } // end for i
304 
305  g->SetFillColor(errCol);
306  //~ g->Draw("e2 same");
307  g->GetYaxis()->SetRangeUser(0, headroom*yMax);
308  nom->GetYaxis()->SetRangeUser(0, headroom*yMax);
309 
310  //~ nom->Draw("hist ][ same");
311  return g;
312  }
TSpline3 lo("lo", xlo, ylo, 12,"0")
enum BeamMode kRed
Spectrum * fLowerRMS
spectra for lower rms
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
T sqrt(T number)
Definition: d0nt_math.hpp:156
void swap(art::HLTGlobalStatus &lhs, art::HLTGlobalStatus &rhs)
TSpline3 hi("hi", xhi, yhi, 18,"0")
Int_t col[ntarg]
Definition: Style.C:29
#define pot
Spectrum * fUpperRMS
spectra for upper rms
Float_t w
Definition: plot.C:20
const Spectrum * Nominal() const
return the nominal universe
TH1 * ana::GenieMultiverseSpectra::ToTH1 ( ) const
inherited

return the nominal histogram with the RMS of all other universes as error band

Definition at line 314 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::Nominal(), and POT.

315  {
316  return ToTH1(Nominal()->POT());
317  }
std::vector< double > POT
const Spectrum * Nominal() const
return the nominal universe
TH1 * ana::GenieMultiverseSpectra::ToTH1 ( double  pot) const
inherited

Definition at line 325 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::fSpectra, MECModelEnuComparisons::i, it, kRed, nbinsx, ana::GenieMultiverseSpectra::Nominal(), push_back(), and ana::Spectrum::ToTH1().

326  {
327  TH1* hNominal = Nominal()->ToTH1(pot);
328  vector<TH1*> hUniverses;
329  map< int, vector<float> > counts_in_bins;
330  for(auto& it: fSpectra)
331  {
332  TH1* curhist = it->ToTH1(pot);
333  hUniverses.push_back(curhist);
334  int nbinsx = curhist->GetNbinsX();
335  for(int i = 1; i <= nbinsx; i++)
336  counts_in_bins[i].push_back(curhist->GetBinContent(i));
337  }
338 
339  for(auto& it: counts_in_bins)
340  {
341  hNominal->SetBinError(it.first, TMath::RMS(it.second.begin(), it.second.end()));
342  //~ cout << TMath::RMS(it.second.begin(), it.second.end()) << endl;
343  }
344  hNominal->SetLineColor(kRed);
345  hNominal->SetMarkerColor(kRed);
346 
347  return hNominal;
348  }
std::vector< Spectrum * > fSpectra
enum BeamMode kRed
set< int >::iterator it
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
base_types push_back(int_type())
#define pot
Int_t nbinsx
Definition: plot.C:23
const Spectrum * Nominal() const
return the nominal universe
TH2 * ana::GenieMultiverseSpectra::ToTH2 ( ) const
inherited

Definition at line 319 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::Nominal(), and POT.

320  {
321  return ToTH2(Nominal()->POT());
322  }
std::vector< double > POT
const Spectrum * Nominal() const
return the nominal universe
TH2 * ana::GenieMultiverseSpectra::ToTH2 ( double  pot) const
inherited

Definition at line 351 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::fSpectra, MECModelEnuComparisons::i, it, calib::j, make_pair(), nbinsx, ana::GenieMultiverseSpectra::Nominal(), and ana::Spectrum::ToTH2().

352  {
353  TH2* hNominal = Nominal()->ToTH2(pot);
354  vector<TH2*> hUniverses;
355  map< std::pair<int,int>, vector<float> > counts_in_bins;
356  for(auto& it: fSpectra)
357  {
358  TH2* curhist = it->ToTH2(pot);
359  hUniverses.push_back(curhist);
360  int nbinsx = curhist->GetNbinsX();
361  int nbinsy = curhist->GetNbinsY();
362  for(int i = 1; i <= nbinsx; i++)
363  {
364  for(int j = 1; j<= nbinsy; j++)
365  {
366  counts_in_bins[std::make_pair(i,j)].push_back(curhist->GetBinContent(i,j));
367  }
368  }
369  }
370 
371  for(auto& it: counts_in_bins)
372  {
373  hNominal->SetBinError(it.first.first,it.first.second, TMath::RMS(it.second.begin(), it.second.end()));
374  //~ cout << TMath::RMS(it.second.begin(), it.second.end()) << endl;
375  }
376 
377  return hNominal;
378  }
std::vector< Spectrum * > fSpectra
TH2 * ToTH2(double exposure, EExposureType expotype=kPOT, EBinType bintype=kBinContent) const
Spectrum must be 2D to obtain TH2.
Definition: Spectrum.cxx:165
set< int >::iterator it
std::pair< Spectrum *, CheatDecomp * > make_pair(SpectrumLoaderBase &loader_data, SpectrumLoaderBase &loader_mc, HistAxis *axis, Cut *cut, const SystShifts &shift, const Var &wei)
Definition: DataMCLoad.C:336
#define pot
const double j
Definition: BetheBloch.cxx:29
Int_t nbinsx
Definition: plot.C:23
const Spectrum * Nominal() const
return the nominal universe
const Spectrum * ana::GenieMultiverseSpectra::UpperRMS ( ) const
inherited

return RMS extremes

Definition at line 381 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::fBoundariesFound, ana::GenieMultiverseSpectra::FindBandBoundaries(), and ana::GenieMultiverseSpectra::fUpperRMS.

382  {
384  return fUpperRMS;
385  }
Spectrum * fUpperRMS
spectra for upper rms
void FindBandBoundaries(BandOptions opt=kBandFromNominal) const
const Spectrum * ana::GenieMultiverseSpectra::UpperSigma ( BandOptions  opt = kBandFromNominal) const
inherited
const Spectrum * ana::GenieMultiverseNormalizedSpectra::UpperSigma ( BandOptions  opt = kBandFromNominal)

Definition at line 892 of file GenieMultiverseSyst.cxx.

References ana::GenieMultiverseSpectra::UpperSigma().

893  {
896  }
bool fSpectraNormalized
A flag signaling whether the spectra are normalized already.
const Spectrum * UpperSigma(BandOptions opt=kBandFromNominal) const
void NormalizeSpectra()
Normalize spectra. It has to be done AFTER SpectrumLoader::Go() is called.

Member Data Documentation

BandOptions ana::GenieMultiverseSpectra::fCurBandOpt = kNoBand
mutableprotectedinherited
knob_sampling_mode_set ana::GenieMultiverseSpectra::fKnobConfigTable
protectedinherited
double ana::GenieMultiverseNormalizedSpectra::fNominalCount
private

Definition at line 323 of file GenieMultiverseSyst.h.

std::vector<Spectrum*> ana::GenieMultiverseSpectra::fSpectra
protectedinherited
bool ana::GenieMultiverseNormalizedSpectra::fSpectraNormalized = false
private

A flag signaling whether the spectra are normalized already.

Definition at line 326 of file GenieMultiverseSyst.h.

int ana::GenieMultiverseSpectra::fUniverseSeed = 1001
protectedinherited

Random seed sequence. Need to have each universe with the same physical parameters. Start from 1001, giving way to the PPFX multi-universe.

Definition at line 132 of file GenieMultiverseSyst.h.

Referenced by ana::GenieMultiverseSpectra::CreateUniverse().


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