Namespaces | Functions | Variables
wrong_sign_uncertainty.C File Reference
#include "CAFAna/Core/Spectrum.h"
#include "CAFAna/Core/SpectrumLoader.h"
#include "CAFAna/Cuts/SpillCuts.h"
#include "CAFAna/Vars/XsecTunes.h"
#include "CAFAna/Vars/PPFXWeights.h"
#include "CAFAna/XSec/Multiverse.h"
#include "CAFAna/XSec/SystematicDef.h"
#include "CAFAna/XSec/GenieMultiverseSyst.h"
#include "CAFAna/Systs/XSecSystLists.h"
#include "CAFAna/Analysis/Exposures.h"
#include "CAFAna/Analysis/Style.h"
#include "CAFAna/Systs/BeamSysts.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TLegend.h"
#include "TLegendEntry.h"
#include "TPaletteAxis.h"
#include "THStack.h"
#include "NDAna/nuebarcc_inc/NuebarCCIncCuts.h"
#include "NDAna/nuebarcc_inc/NuebarCCIncVars.h"
#include "NDAna/nuebarcc_inc/NuebarCCIncBins.h"
#include "NDAna/nuebarcc_inc/NuebarCCIncExtra.h"
#include "NDAna/muonid/NDXSecMuonPID.h"

Go to the source code of this file.

Namespaces

 nuebarccinc
 

Functions

void load_libs_muonid ()
 
TH2 * ToTH2 (const TH1 *h, const HistAxis &xaxis, const HistAxis &yaxis)
 
TH1 * MostConservative (TH1 *up, TH1 *down, TH1 *nominal)
 
void PlotMultiverse (Multiverse &mv, TH1 *nominal, std::string title, std::string name, TH1 *error=0)
 
void Plot1D (TH1 *h1, std::string h1_label, std::map< std::string, TH1 * > systs, std::vector< std::string > to_plot, std::string name, double xmax=-1)
 
void Plot2D (TH1 *total, std::map< std::string, TH1 * > systs, std::vector< std::string > to_plot, std::string basename)
 
void abs (TH1 *hist)
 
void MovePalette (TH1 *hist)
 
void wrong_sign_uncertainty (bool make_plots=true, std::string input_file_name="/nova/ana/users/ddoyle/NuebarCCInc/WrongSign/wrong_sign_uncertainty.root", std::string plot_dump="/nova/ana/users/ddoyle/NuebarCCInc/WrongSign/plots")
 

Variables

std::vector< std::string > all_systs = {"ppfx", "genie", "calib_shape", "calib_shift", "ll", "cherenkov", "beam_focusing", "stat"}
 
std::vector< std::string > beam_systs = {"ppfx", "beam_focusing"}
 
std::map< std::string, nuebarccinc::Style > syst_styles
 
std::map< std::string, double > xmax
 
std::map< string, const ana::BeamSyst * > beamSystematics
 

Function Documentation

void abs ( TH1 *  hist)

Definition at line 716 of file wrong_sign_uncertainty.C.

References std::abs(), febshutoff_auto::val, submit_syst::x, submit_syst::y, and test::z.

Referenced by fuzz::FuzzyKVertex::AddProng(), AddUpDwUncerts(), AddUpDwUncerts2D(), mcchk::ShowerAnaCheck::analyze(), remid::ReMIdTrain::analyze(), slicer::S4DParamCalc::analyze(), beamlinereco::ToFSingleCounterAnalysis::analyze(), mcchk::CosmicAna::analyze(), febstat::EventFEBStatus::analyze(), murem::MuonRemoveAna::analyze(), mcchk::NeutrinoAna::analyze(), novaddt::SNAna::analyze(), zcl::MMCheater::analyze(), fuzz::FuzzyKValidate::analyze(), bpfit::BPFTmvaTrainer::analyze(), ncs::Xeff::analyze(), slid::LIDTraining::analyze(), ncs::NCAna::analyze(), comi::Leana::AnaNHit(), lem::MakePID::AvgAndFit(), slid::ParticleIDAlg::CalcCellPlaneTransverseDedx(), ana::CalcPrimE(), cheat::BackTracker::CalcTotalEscapingEnergy(), Calculate2DChi2(), trk::KalmanTrackMerge::CanJoinTracks(), osc::_OscCalcDMP< T >::ChannelCacheIdx(), osc::OscCalcPMNSOptEigen::ChannelCacheIdx(), trk::KalmanTrackAna::CheckRecoTracks(), checkCalibrationCSV.NovaCalibTest::checkRunCoverage(), genie::NievesQELCCPXSec::CNCTCLimUcalc(), om::PlotViewer::ContentProjectionTH2F(), gibuu::GiBUURegen::CopyGenieEvent(), genie::NievesQELCCPXSec::deltaLindhard(), novaddt::HoughTrackMaker::do_like_zukai(), plot_xsec_1d::do_single_shift_1bin(), plot_xsec_1d::do_up_down_1bin(), evd::HitTower::Draw(), ana::SystMakerShift::DrawBkg(), ana::SystMaker::DrawBkg(), ana::SystMakerShift::DrawShift(), ana::SystMakerShift::DrawSig(), ana::SystMaker::DrawSig(), evgen::CosmicsGen::Enhance(), beamlineutil::BeamlineMagneticField::FieldFromCurrent(), calib::StopperThreshold::FillHist(), slidt::FillHists(), std_candles.StdCandlesPlotter::FillMichelEPlots(), caf::FillNumuInfo(), lem::FillPotential(), calib::StopperThreshold::FillTree(), calib::EnergyStandardCandles::FillTruthInfo(), murem::MuonRemove::FillTruthInfo(), comi::CanaFilter::filter(), util::TimeFilter::filter(), mono::MonopoleTrack::filter(), novaddt::MichelETrigger::filter(), trk::KalmanTrack::FilterTracker(), caf::FindAndAddMichels(), generic.tools.HistogramTools::findMaxAsymmetry(), common_tools::findMaxAsymmetry(), generic.tools.HistogramTools::findMaxSkew(), common_tools::findMaxSkew(), bsf::BremShowerFilter::findShowerByTruth(), trk::KalmanTrack::FindTracks(), nuesand::FillNueSandbox::FracAngChanges(), lem::FracChargeMatched(), GetAbsolute(), GetAbsolute2D(), GetBinByEdge(), nuesand::FillNueSandbox::GetECF(), gibuu::GiBUURegen::GetEvent(), cvn::GetInteractionType(), cvn::GetInteractionTypeFromSlice(), GetMaxUpDown(), GetMaxUpDown2D(), numusand::NumuSandFxs::getMissingE(), cvn::GetMultiplicityMap(), getPlots(), getPlotsFromNtupleFD(), getPlotsFromNtupleND(), gibuu::GiBUURegen::GetRecordList(), ana::PredictionWriter::GetShiftName(), GetShifts2017(), ana::GetShifts2018(), ana::GetShifts2020(), ana::FitDelta24InPiUnitsSterile::GetValue(), ana::MINERvA_MEC_Wgt_Var::GetWeight(), mono::Cluster::good_track(), novaddt::smt::Trigger::good_track(), DecTreeLoader::HandleRecord(), ana::hasBPFE(), genie::BaryonResonanceDecayer::HasEvolvedBRs(), ana::MultiHistContainer::Hist(), lem::DistanceMap::InvDist(), genie::BaryonResonanceDecayer::IsDelta(), ana::IsElectronPng(), ana::IsMuonPng(), ana::IsPionPng(), ana::IsProtonPng(), lem::EventSummary::IsSig(), caf::IsTrueDiF(), calib::RockMuonStopperSelection::IsTrulyContainedStopper(), calib::RockMuonStopperSelection::IsTrulyInFiducial(), calib::ThroughgoingSelection::IsTrulyUncontained(), skim::EvaluatorNumuCCpi::KeepEventSlices(), skim::EvaluatorNumuCCpi::KeepSlice(), kinematics(), ana::kMinPngNhit(), ana::kNue2017BasicPartFunc(), load_libs_muonid(), readFlux::Loop(), mono::Cluster::make(), make_dst_cosrejbdttrain(), make_fc_mass_and_oct_nersc_2018(), make_fc_mass_and_oct_nersc_2019(), make_mass_and_oct_fc_input_2019(), cvn::MakeParticleLabel(), airshower::AirSlicer::MakeTrackSlices(), ana::FluxDecomp::MakeWeightsNumuFromKaon(), ana::BENDecomp::MakeWeightsNumuFromKaon(), MakingProfileHadCC(), hough::MultiHough2P::Map(), dt::ViewMerger::MatchScore(), lem::FindLEMMatches::MatchToVertex(), lem::LEM::MatchToVertex(), genie::geometry::ROOTGeomAnalyzer::MaxPathLengthsFluxMethod(), evd::SimulationDrawer::MCTruthTrajectoriesAnyD(), ana::MuonECat(), novaddt::DBSlicer::NeighborScore(), fnex::ShifterAndWeighter::NueAcceptBkg2018Weight(), fnex::ShifterAndWeighter::NueAcceptSignalKin2018Weight(), fnex::ShifterAndWeighter::NueExtrapBkg2017Weight(), fnex::ShifterAndWeighter::NueExtrapSig2017Weight(), NuMu2020_TrimCAFs(), ana::GetBPFMuonID::operator()(), ana::CCFlavSel::operator()(), ana::GetBestBPFTrack::operator()(), ana::NCFlavSel::operator()(), oscW(), osc::OscCalcDumb::P(), osc::OscCalcPMNS_CPT::P(), osc::OscCalcPMNS_NSI::P(), osc::_OscCalcPMNS< T >::P(), osc::_OscCalcPMNSOpt< T >::P(), osc::OscCalcSterile::P(), osc::OscCalcGeneral::P(), osc::OscCalcPMNSOptEigen::P(), cvn::ParticlesSliceClassify(), PlotSystErrorBand(), ppfx_make_systs(), ppfx_smooth_weights_save(), Print(), plot_xsec_2d::process_ckv_calibshape(), plot_xsec_2d::process_focusing(), plot_xsec_2d::process_light_calibration_mue(), demo.ScoreSD::ProcessHits(), ana::SystMakerShift::ProcessShift(), ana::Prod3AbsCalibLoaders::Prod3AbsCalibLoaders(), ana::Prod3LightLevelLoaders::Prod3LightLevelLoaders(), ana::Prod4AbsCalibLoaders::Prod4AbsCalibLoaders(), ana::Prod4LightLevelLoaders::Prod4LightLevelLoaders(), ana::Prod5AbsCalibLoaders::Prod5AbsCalibLoaders(), ana::Prod5LightLevelLoaders::Prod5LightLevelLoaders(), cheat::ProngCheater::produce(), cheat::ShowerCheater::produce(), slid::Recluster::produce(), numusand::FillSandbox::produce(), calib::FiberCalibration::produce(), lem::AlignLibToVtx::produce(), zcl::SPCluster::produce(), evgen::CosmicPionGen::produce(), beamlinereco::CFDHitReco::produce(), evgen::CORSIKAGen::produce(), dif::DiFShowerFinder::produce(), cheat::MCCheater::produce(), g4rwgt::MakeG4WeightTable::produce(), murem::MuonRemove::produce(), caf::CAFMaker::produce(), evgen::CosmicsGen::ProjectMuonsToDetectorBigBox(), ana::ProngTruthCut(), PlotUtilityNew::pull_NormData(), murem::MuonRemove::RemoveByEfficiency(), murem::MuonRemove::RemoveByProngCVN(), murem::MuonRemove::RemoveByTruth(), resolution_event_tree(), RockFilter_D(), novaddt::smt::Trigger::run_algorithm(), evgen::SingleGen::Sample(), screen_to_activecell(), simb::MCTruth::SetNeutrino(), ana::NueExtrapSystBkg2017::Shift(), ana::NueAcceptSystBkg2018FHC::Shift(), ana::NueExtrapSystSignalKin2017::Shift(), ana::NuWROSyst::Shift(), ana::NueAcceptSystSignalKin2018FHC::Shift(), ana::NueAcceptSystSignalKin2020FHC::Shift(), ana::NumuSummedSmallGENIESyst::Shift(), ana::NueAcceptSystBkg2018RHC::Shift(), ana::BeamNueScaleSyst::Shift(), ana::NueAcceptSystSignalKin2018RHC::Shift(), ana::NumuCCScaleSyst::Shift(), ana::NueAcceptSystSignalKin2020RHC::Shift(), ana::NumuTauContaminationSyst::Shift(), signFunc(), ncs::NCAna::SimpleOscProb(), lem::Library::SimpleSurvivalProb(), lem::SimpleSurvivalProb(), murem::MuonRemove::SliceSelection(), murem::MRCCAna::SliceSelection(), cheat::BackTracker::SlicesToMCTruthsTable(), lem::LEMSummarizer::SliceToLEMInput(), osc::_PMNS< T >::SortEigenvalues(), cmf::ShifterAndWeighter::SumSmallXSec2017Weight(), fnex::ShifterAndWeighter::SumSmallXSecNumu2017Weight(), test_flux(), test_fluxhadr_prod_weights_Flux(), test_genieweights(), caf_tools::tightBeamQualityCut(), ana::TopHatPriorInFitVar(), mono::Track3D::Track3D(), cheat::BackTracker::TrackIDToMCTruth(), cheat::BackTracker::TrackIDToMotherParticle(), cheat::BackTracker::TrackIDToParticle(), train_dectree_caf(), TrimCAFs(), trimmubarid(), trimvar(), nuebarccinc::true_electron_png(), slid::NuEEnergyAlg::VertexEnergy(), HandyFuncs::VerticalRange(), ana::VisibleHadE(), and wrong_sign_uncertainty().

717 {
718  for(auto x = 1; x <= hist->GetNbinsX(); x++) {
719  for(auto y = 1; y <= hist->GetNbinsY(); y++) {
720  for(auto z = 1; z <= hist->GetNbinsZ(); z++) {
721  auto val = hist->GetBinContent(x, y, z);
722  hist->SetBinContent(x, y, z, std::abs(val));
723  }
724  }
725  }
726 }
std::abs
float abs(float number)
Definition: d0nt_math.hpp:39
submit_syst.y
y
Definition: submit_syst.py:249
submit_syst.x
x
Definition: submit_syst.py:249
test.z
z
Definition: test.py:28
febshutoff_auto.val
val
Definition: febshutoff_auto.py:277
void load_libs_muonid ( )

Definition at line 32 of file wrong_sign_uncertainty.C.

References abs(), exit(), make_syst_table_plots::h, h1, analysePickle::hist, MostConservative(), MovePalette(), Plot1D(), Plot2D(), PlotMultiverse(), runNovaSAM::ret, string, systs, plotROC::title, produceInfoGainMatrix::total, ToTH2(), xaxis, xmax, and yaxis.

Referenced by wrong_sign_uncertainty().

33 {
34  const int ret = gSystem->Load("libNDAnamuonid");
35  if(ret != 0) exit(1);
36 }
exit
exit(0)
TH1 * MostConservative ( TH1 *  up,
TH1 *  down,
TH1 *  nominal 
)

Definition at line 821 of file wrong_sign_uncertainty.C.

References std::abs(), std::max(), runNovaSAM::ret, ana::UniqueName(), submit_syst::x, submit_syst::y, and test::z.

Referenced by load_libs_muonid(), and wrong_sign_uncertainty().

822 {
823  auto ret = (TH1*) up->Clone(UniqueName().c_str());
824  auto NX = up->GetNbinsX();
825  auto NY = up->GetNbinsY();
826  auto NZ = up->GetNbinsZ();
827  for(auto x = 1; x <= NX; x++) {
828  for(auto y = 1; y <= NY; y++) {
829  for(auto z = 1; z <= NZ; z++) {
830  auto nom = nominal->GetBinContent(x, y, z);
831  auto shift_up = std::abs(nom - up->GetBinContent(x, y, z));
832  auto shift_down = std::abs(nom - down->GetBinContent(x, y, z));
833 
834  ret->SetBinContent(x, y, z, std::max(shift_up, shift_down) + nom);
835  }
836  }
837  }
838  return ret;
839 }
std::max
T max(const caf::Proxy< T > &a, T b)
Definition: BasicTypesProxy.h:414
std::abs
float abs(float number)
Definition: d0nt_math.hpp:39
submit_syst.y
y
Definition: submit_syst.py:249
submit_syst.x
x
Definition: submit_syst.py:249
test.z
z
Definition: test.py:28
ana::UniqueName
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:29
void MovePalette ( TH1 *  hist)

Definition at line 841 of file wrong_sign_uncertainty.C.

Referenced by load_libs_muonid(), Plot2D(), and wrong_sign_uncertainty().

842 {
843  hist->Draw("colz");
844  gPad->Update();
845  auto palette = (TPaletteAxis*) hist->GetListOfFunctions()->FindObject("palette");
846  palette->SetX1NDC(0.81);
847  palette->SetX2NDC(0.85);
848 }
void Plot1D ( TH1 *  h1,
std::string  h1_label,
std::map< std::string, TH1 * >  systs,
std::vector< std::string >  to_plot,
std::string  name,
double  xmax = -1 
)

Definition at line 790 of file wrong_sign_uncertainty.C.

References plot_validation_datamc::c, MECModelEnuComparisons::leg, std::max(), maxy, and syst_styles.

Referenced by load_libs_muonid(), and wrong_sign_uncertainty().

790  {
791  TCanvas * c = new TCanvas();
792  TLegend * leg = new TLegend();
793 
794  // find max y
795  double maxy = -999;
796  for(auto syst : to_plot)
797  maxy = std::max(maxy, systs.at(syst)->GetMaximum());
798 
799  maxy = std::max(maxy, h1->GetMaximum());
800 
801 
802  if(xmax > 0)
803  h1->GetXaxis()->SetRangeUser(0, xmax);
804 
805  h1->GetYaxis()->SetRangeUser(0, maxy * 1.2);
806  h1->SetLineColor(kBlack);
807  h1->Draw("hist");
808  leg->AddEntry(h1, h1_label.c_str(), "l");
809 
810  for(auto syst : to_plot ) {
811  if(syst == "nominal") continue;
812  systs.at(syst)->SetLineColor(syst_styles.at(syst).color);
813  systs.at(syst)->SetLineStyle(syst_styles.at(syst).line_style);
814  leg->AddEntry(systs.at(syst), syst.c_str(), "l");
815  systs.at(syst)->Draw("hist same");
816  }
817  h1->Draw("hist same");
818  leg->Draw();
819  c->Print(name.c_str());
820 }
std::max
T max(const caf::Proxy< T > &a, T b)
Definition: BasicTypesProxy.h:414
name
const XML_Char * name
Definition: expat.h:151
xmax
std::map< std::string, double > xmax
Definition: wrong_sign_uncertainty.C:79
systs
std::vector< SystGroupDef > systs
Definition: syst_header.h:385
maxy
double maxy
Definition: plot_datapredictions.C:15
plot_validation_datamc.c
c
Definition: plot_validation_datamc.py:198
syst_styles
std::map< std::string, nuebarccinc::Style > syst_styles
Definition: wrong_sign_uncertainty.C:50
h1
TH1F * h1
Definition: nslices_potspill.C:4
void Plot2D ( TH1 *  total,
std::map< std::string, TH1 * >  systs,
std::vector< std::string >  to_plot,
std::string  basename 
)

Definition at line 760 of file wrong_sign_uncertainty.C.

References plot_validation_datamc::c, nuebarccinc::costhetabins, nuebarccinc::eelecbins, MovePalette(), ana::ToTH2Helper(), and ana::UniqueName().

Referenced by estimate_energy(), GetHistsFD(), GetHistsND(), load_libs_muonid(), and wrong_sign_uncertainty().

761 {
762 
763  TCanvas * c = new TCanvas();
764  c->SetRightMargin(0.2);
765  std::map<std::string, TH1*> hists2d;
766  hists2d["total"] = ToTH2Helper(std::unique_ptr<TH1>((TH1*) total->Clone(UniqueName().c_str())),
767  nuebarccinc::costhetabins,
768  nuebarccinc::eelecbins);
769  for(auto syst : to_plot) {
770  hists2d[syst] = ToTH2Helper(std::unique_ptr<TH1>((TH1*) systs.at(syst)->Clone(UniqueName().c_str())),
771  nuebarccinc::costhetabins,
772  nuebarccinc::eelecbins);
773  }
774 
775  for(auto hist_it = hists2d.begin(); hist_it != hists2d.end(); hist_it++) {
776  hist_it->second->SetTitle("");
777  hist_it->second->GetXaxis()->SetTitle("True cos#theta_{e}");
778  hist_it->second->GetYaxis()->SetTitle("True Electron Energy (GeV)");
779  hist_it->second->GetXaxis()->SetRangeUser(0.85, 1);
780  hist_it->second->GetYaxis()->SetRangeUser(0, 6);
781  hist_it->second->SetMaximum(1);
782  hist_it->second->SetMarkerColor(kWhite);
783  MovePalette(hist_it->second);
784  hist_it->second->Draw("colz text");
785  c->Print((basename + "_" + hist_it->first + ".pdf").c_str());
786  }
787 
788 }
systs
std::vector< SystGroupDef > systs
Definition: syst_header.h:385
nuebarccinc::eelecbins
const ana::Binning eelecbins
Definition: NuebarCCIncBins.h:19
nuebarccinc::costhetabins
const ana::Binning costhetabins
Definition: NuebarCCIncBins.h:23
MovePalette
void MovePalette(TH1 *hist)
Definition: wrong_sign_uncertainty.C:841
plot_validation_datamc.c
c
Definition: plot_validation_datamc.py:198
ana::ToTH2Helper
TH2 * ToTH2Helper(std::unique_ptr< TH1 > h1, const Binning &binsx, const Binning &binsy, ana::EBinType bintype)
Helper for ana::ToTH2.
Definition: UtilsExt.cxx:135
ana::UniqueName
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:29
void PlotMultiverse ( Multiverse &  mv,
TH1 *  nominal,
std::string  title,
std::string  name,
TH1 *  error = 0 
)

Definition at line 728 of file wrong_sign_uncertainty.C.

References std::abs(), plot_validation_datamc::c, Draw(), make_mec_shifts_plots::GetMaximum(), ana::Multiverse::GetUniversesHist(), hists, kBlue, std::max(), maxy, submit_hadd::u, and submit_syst::x.

Referenced by load_libs_muonid(), and wrong_sign_uncertainty().

729 {
730  auto hists = mv.GetUniversesHist(kCyan);
731 
732  double maxy = -1000;
733  for(auto ihist = 0u; ihist < hists.size(); ihist++)
734  maxy = std::max(maxy, hists[ihist]->GetMaximum());
735 
736  TCanvas * c = new TCanvas();
737  nominal->GetYaxis()->SetRangeUser(0, 1.2 * maxy);
738  nominal->SetLineColor(kBlue);
739  nominal->SetTitle(title.c_str());
740 
741  nominal->Draw("hist");
742 
743  for(auto ihist = 0u; ihist < hists.size(); ihist++)
744  hists[ihist]->Draw("hist same");
745 
746 
747  if(error) {
748  for(auto x = 1; x <= nominal->GetNbinsX(); x++) {
749  nominal->SetBinError(x, std::abs(error->GetBinContent(x) - nominal->GetBinContent(x)));
750  }
751  nominal->Draw("hist same e1");
752  }
753  else
754  nominal->Draw("hist same");
755 
756 
757  c->Print(name.c_str());
758 }
std::max
T max(const caf::Proxy< T > &a, T b)
Definition: BasicTypesProxy.h:414
name
const XML_Char * name
Definition: expat.h:151
Draw
tree Draw("slc.nhit")
maxy
double maxy
Definition: plot_datapredictions.C:15
hists
TString hists[nhists]
Definition: bdt_com.C:3
std::abs
float abs(float number)
Definition: d0nt_math.hpp:39
submit_syst.x
x
Definition: submit_syst.py:249
submit_hadd.u
u
Definition: submit_hadd.py:83
plot_validation_datamc.c
c
Definition: plot_validation_datamc.py:198
ana::Multiverse::GetUniversesHist
const std::vector< TH1 * > GetUniversesHist()
Definition: Multiverse.cxx:120
make_mec_shifts_plots.GetMaximum
def GetMaximum(hist)
Definition: make_mec_shifts_plots.py:20
kBlue
enum BeamMode kBlue
TH2* ToTH2 ( const TH1 *  h,
const HistAxis &  xaxis,
const HistAxis &  yaxis 
)

Referenced by ana::GenieMultiverseSpectra::GenieMultiverseSpectra(), and load_libs_muonid().

void wrong_sign_uncertainty ( bool  make_plots = true,
std::string  input_file_name = "/nova/ana/users/ddoyle/NuebarCCInc/WrongSign/wrong_sign_uncertainty.root",
std::string  plot_dump = "/nova/ana/users/ddoyle/NuebarCCInc/WrongSign/plots" 
)

Definition at line 97 of file wrong_sign_uncertainty.C.

References abs(), all_systs, file_size_ana::axes, allInOneTrainingPlots::axis, beam_systs, beamSystematics, plot_validation_datamc::c, nuebarccinc::costhetabins, Cut(), update_sam_good_runs_metadata::cuts, den, ana::Multiverse::Divide(), nuebarccinc::eelecbins, genie::utils::style::Format(), ana::getAllXsecSysts_2020(), ana::GetkPPFXFluxUnivWgt(), ana::Multiverse::GetMinusOneSigmaShift(), ana::Multiverse::GetPlusOneSigmaShift(), ana::GenieMultiverseParameters::GetSystShifts(), hists, MECModelEnuComparisons::i, ana::kAna2020RHCPOT, kGreen, nuebarccinc::kMuonIDProd4Cut, ana::kNoShift, nuebarccinc::kNuebarCC, nuebarccinc::kNueCC, ana::kPPFXFluxCVWgt, nuebarccinc::kPreselectionLoose, ana::kSelection, ana::kStandardSpillCuts, nuebarccinc::kTrueElectronEVsCosStandardAxis, nuebarccinc::kTrueNeutrinoEnergyStandardAxis, nuebarccinc::kTrueQ2StandardAxis, ana::kXSecCVWgt2020, MECModelEnuComparisons::leg, load_libs_muonid(), loaders, ana::LoadFrom(), make_plots(), MostConservative(), MovePalette(), getGoodRuns4SAM::n, num, output, Plot1D(), Plot2D(), PlotMultiverse(), pot, pots, cet::pow(), nuebarccinc::PROD5_MC_RHC_CALIB_DOWN, nuebarccinc::PROD5_MC_RHC_CALIB_SHAPE, nuebarccinc::PROD5_MC_RHC_CALIB_UP, nuebarccinc::PROD5_MC_RHC_CHERENKOV, nuebarccinc::PROD5_MC_RHC_LIGHT_DOWN, nuebarccinc::PROD5_MC_RHC_LIGHT_UP, nuebarccinc::PROD5_MC_RHC_NOMINAL, std::sqrt(), string, ana::ToTH2Helper(), ana::UniqueName(), febshutoff_auto::val, submit_syst::x, xmax, submit_syst::y, and test::z.

100 {
101  load_libs_muonid();
102  std::map<std::string, const HistAxis*> axes;
103  axes["electron_kin"] = new HistAxis(nuebarccinc::kTrueElectronEVsCosStandardAxis);
104  axes["enu" ] = new HistAxis(nuebarccinc::kTrueNeutrinoEnergyStandardAxis);
105  axes["q2" ] = new HistAxis(nuebarccinc::kTrueQ2StandardAxis );
106 
107  std::map<std::string, const Cut*> cuts;
108  cuts["nue" ] = new Cut(nuebarccinc::kNueCC );
109  cuts["nuebar" ] = new Cut(nuebarccinc::kNuebarCC);
110  cuts["nunubar"] = new Cut(nuebarccinc::kNuebarCC || nuebarccinc::kNueCC);
111 
112  std::map<std::string, SpectrumLoader*> loaders;
113  loaders["nominal" ] = 0;
114  loaders["calib_shift_up"] = 0;
115  loaders["calib_shift_dw"] = 0;
116  loaders["ll_up"] = 0;
117  loaders["ll_dw"] = 0;
118  loaders["cherenkov"] = 0;
119  loaders["calib_shape"] = 0;
120  if(!make_plots) {
121  loaders["nominal" ] = new SpectrumLoader(nuebarccinc::PROD5_MC_RHC_NOMINAL);
122  loaders["calib_shift_up"] = new SpectrumLoader(nuebarccinc::PROD5_MC_RHC_CALIB_UP);
123  loaders["calib_shift_dw"] = new SpectrumLoader(nuebarccinc::PROD5_MC_RHC_CALIB_DOWN);
124  loaders["calib_shape"] = new SpectrumLoader(nuebarccinc::PROD5_MC_RHC_CALIB_SHAPE);
125  loaders["ll_up"] = new SpectrumLoader(nuebarccinc::PROD5_MC_RHC_LIGHT_UP);
126  loaders["ll_dw"] = new SpectrumLoader(nuebarccinc::PROD5_MC_RHC_LIGHT_DOWN);
127  loaders["cherenkov"] = new SpectrumLoader(nuebarccinc::PROD5_MC_RHC_CHERENKOV);
128  const Var kCVWeight = kPPFXFluxCVWgt * kXSecCVWgt2020;
129 
130  const Cut kSelection =
131  nuebarccinc::kPreselectionLoose &&
132  nuebarccinc::kMuonIDProd4Cut;
133 
134  std::map<std::string, std::map<std::string, std::map<std::string, Spectrum*> > > spectra;
135  for(auto load_it = loaders.begin(); load_it != loaders.end(); load_it++) {
136  load_it->second->SetSpillCut(kStandardSpillCuts);
137  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
138  for(auto cut_it = cuts.begin(); cut_it != cuts.end(); cut_it++) {
139  spectra[load_it->first][cut_it->first][ax_it->first] = new Spectrum(*load_it->second,
140  *ax_it->second,
141  kSelection && *cut_it->second,
142  kNoShift,
143  kCVWeight);
144  }
145  }
146  }
147 
148  // beam systs
149  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
150  for(auto cut_it = cuts.begin(); cut_it != cuts.end(); cut_it++) {
151  for(auto beam_it = beamSystematics.begin(); beam_it != beamSystematics.end(); beam_it++) {
152  auto name_up = beam_it->first + "_up";
153  auto name_dw = beam_it->first + "_dw";
154 
155  spectra[name_up][cut_it->first][ax_it->first] = new Spectrum(*loaders["nominal"],
156  *ax_it->second,
157  kSelection && *cut_it->second,
158  ana::SystShifts(beam_it->second, +1),
159  kCVWeight);
160 
161  spectra[name_dw][cut_it->first][ax_it->first] = new Spectrum(*loaders["nominal"],
162  *ax_it->second,
163  kSelection && *cut_it->second,
164  ana::SystShifts(beam_it->second, -1),
165  kCVWeight);
166  }
167  }
168  }
169 
170  const std::vector<Var> ppfx = GetkPPFXFluxUnivWgt();
171  GenieMultiverseParameters genie_multiverse(100, getAllXsecSysts_2020());
172  std::vector<SystShifts> genie_shifts = genie_multiverse.GetSystShifts();
173 
174  std::map<std::string, std::map<std::string, Multiverse *> > mv_ppfx;
175  std::map<std::string, std::map<std::string, Multiverse *> > mv_genie;
176  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
177  for(auto cut_it = cuts.begin(); cut_it != cuts.end(); cut_it++) {
178  mv_ppfx[cut_it->first][ax_it->first] = new Multiverse(*loaders["nominal"],
179  *ax_it->second,
180  kSelection && *cut_it->second,
181  kNoShift,
182  ppfx,
183  kXSecCVWgt2020);
184  mv_genie[cut_it->first][ax_it->first] = new Multiverse(*loaders["nominal"],
185  *ax_it->second,
186  kSelection && *cut_it->second,
187  genie_shifts,
188  kPPFXFluxCVWgt);
189  }
190  }
191 
192  for(auto load_it = loaders.begin(); load_it != loaders.end(); load_it++)
193  load_it->second->Go();
194 
195 
196  TFile* output = new TFile("wrong_sign_uncertainty.root", "recreate");
197  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
198  for(auto cut_it = cuts.begin(); cut_it != cuts.end(); cut_it++) {
199  mv_ppfx .at(cut_it->first).at(ax_it->first)->SaveTo(output, ("ppfx_" + cut_it->first + "_" + ax_it->first).c_str());
200  mv_genie.at(cut_it->first).at(ax_it->first)->SaveTo(output, ("genie_" + cut_it->first + "_" + ax_it->first).c_str());
201  for(auto load_it = loaders.begin(); load_it != loaders.end(); load_it++) {
202  spectra.at(load_it->first).at(cut_it->first).at(ax_it->first)->SaveTo(output, TString::Format("%s_%s_%s",
203  load_it->first.c_str(),
204  cut_it ->first.c_str(),
205  ax_it ->first.c_str()).Data());
206  }
207  for(auto beam_it = beamSystematics.begin(); beam_it != beamSystematics.end(); beam_it++) {
208  spectra.at(beam_it->first + "_up").at(cut_it->first).at(ax_it->first)->SaveTo(output, TString::Format("%s_%s_%s",
209  cut_it ->first.c_str(),
210  ax_it ->first.c_str()).Data());
211  spectra.at(beam_it->first + "_dw").at(cut_it->first).at(ax_it->first)->SaveTo(output, TString::Format("%s_%s_%s",
212  (beam_it->first + "_dw").c_str(),
213  cut_it ->first.c_str(),
214  ax_it ->first.c_str()).Data());
215  }
216  }
217  }
218 
219 
220 
221  output->Close();
222  }
223  else {
224  TFile * input = TFile::Open(input_file_name.c_str());
225 
226  double MCPOT = ((TH1*) input->Get("nominal_nue_q2/pot"))->GetBinContent(1);
227  std::map<std::string, double> pots;
228  std::map<std::string, std::map<std::string, Multiverse *> > mv_ppfx;
229  std::map<std::string, std::map<std::string, Multiverse *> > mv_genie;
230  std::map<std::string, std::map<std::string, std::map<std::string, TH1*> > > hists;
231  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
232  for(auto cut_it = cuts.begin(); cut_it != cuts.end(); cut_it++) {
233  mv_ppfx [ax_it->first][cut_it->first] = Multiverse::LoadFrom(input, ("ppfx_" + cut_it->first + "_" + ax_it->first).c_str())).release(;
234  mv_genie[ax_it->first][cut_it->first] = Multiverse::LoadFrom(input, ("genie_" + cut_it->first + "_" + ax_it->first).c_str())).release(;
235 
236  mv_ppfx .at(ax_it->first).at(cut_it->first)->Scale(kAna2020RHCPOT / MCPOT);
237  mv_genie.at(ax_it->first).at(cut_it->first)->Scale(kAna2020RHCPOT / MCPOT);
238 
239  for(auto load_it = loaders.begin(); load_it != loaders.end(); load_it++) {
240  auto spec = Spectrum::LoadFrom(input, TString::Format("%s_%s_%s",
241  load_it->first.c_str(),
242  cut_it ->first.c_str(),
243  ax_it ->first.c_str()).Data());
244  pots[load_it->first] = spec->POT();
245  auto pot = kAna2020RHCPOT;//spec->POT();
246  hists[ax_it->first][cut_it->first][load_it->first] = spec->ToTH1(pot);
247  }
248 
249  for(auto beam_it = beamSystematics.begin(); beam_it != beamSystematics.end(); beam_it++) {
250  auto spec_up = Spectrum::LoadFrom(input, TString::Format("%s_%s_%s",
251  (beam_it->first + "_up").c_str(),
252  cut_it ->first .c_str(),
253  ax_it ->first .c_str()).Data());
254  auto spec_dw = Spectrum::LoadFrom(input, TString::Format("%s_%s_%s",
255  (beam_it->first + "_dw").c_str(),
256  cut_it ->first .c_str(),
257  ax_it ->first .c_str()).Data());
258 
259  hists[ax_it->first][cut_it->first][beam_it->first + "_up"] = spec_up->ToTH1(kAna2020RHCPOT);
260  hists[ax_it->first][cut_it->first][beam_it->first + "_dw"] = spec_dw->ToTH1(kAna2020RHCPOT);
261  }
262  }
263  }
264 
265  // make fractions
266  std::map<std::string, std::map<std::string, TH1*> > nuebar_fractions;
267  std::map<std::string, Multiverse*> mv_ppfx_nuebar_fractions;
268  std::map<std::string, Multiverse*> mv_genie_nuebar_fractions;
269  std::map<std::string, Multiverse*> mv_ppfx_nue_nuebar_ratios;
270  std::map<std::string, Multiverse*> mv_genie_nue_nuebar_ratios;
271  std::map<std::string, std::map<std::string, TH1*> > nue_nuebar_ratios;
272  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
273  mv_ppfx_nuebar_fractions [ax_it->first] = new Multiverse(*mv_ppfx .at(ax_it->first).at("nuebar"));
274  mv_genie_nuebar_fractions[ax_it->first] = new Multiverse(*mv_genie.at(ax_it->first).at("nuebar"));
275 
276  *mv_ppfx_nuebar_fractions.at(ax_it->first) /= *mv_ppfx .at(ax_it->first).at("nunubar");
277  *mv_genie_nuebar_fractions.at(ax_it->first) /= *mv_genie.at(ax_it->first).at("nunubar");
278 
279  mv_ppfx_nue_nuebar_ratios[ax_it->first] = new Multiverse(*mv_ppfx.at(ax_it->first).at("nue"));
280  *mv_ppfx_nue_nuebar_ratios.at(ax_it->first) /= *mv_ppfx.at(ax_it->first).at("nuebar");
281 
282  mv_genie_nue_nuebar_ratios[ax_it->first] = new Multiverse(*mv_genie.at(ax_it->first).at("nue"));
283  *mv_genie_nue_nuebar_ratios.at(ax_it->first) /= *mv_genie.at(ax_it->first).at("nuebar");
284 
285  for(auto loader_it = loaders.begin(); loader_it != loaders.end(); loader_it++) {
286  nuebar_fractions[ax_it->first][loader_it->first] = (TH1*) hists.at(ax_it->first).at("nuebar").at(loader_it->first)->Clone(UniqueName().c_str());
287  nuebar_fractions.at(ax_it->first).at(loader_it->first)->Divide(hists.at(ax_it->first).at("nunubar").at(loader_it->first));
288 
289  nue_nuebar_ratios[ax_it->first][loader_it->first] = (TH1*) hists.at(ax_it->first).at("nue").at(loader_it->first)->Clone(UniqueName().c_str());
290  nue_nuebar_ratios.at(ax_it->first).at(loader_it->first)->Divide(hists.at(ax_it->first).at("nuebar").at(loader_it->first));
291  }
292 
293  // two sided file systs
294  for(std::string syst : {"calib_shift", "ll"}) {
295  TH1 * tmp_up = (TH1*) hists.at(ax_it->first).at("nuebar").at(syst + "_up")->Clone(UniqueName().c_str());
296  TH1 * tmp_dw = (TH1*) hists.at(ax_it->first).at("nuebar").at(syst + "_dw")->Clone(UniqueName().c_str());
297  tmp_up->Divide(hists.at(ax_it->first).at("nunubar").at(syst + "_up"));
298  tmp_dw->Divide(hists.at(ax_it->first).at("nunubar").at(syst + "_dw"));
299 
300  nuebar_fractions[ax_it->first][syst] = MostConservative(tmp_up, tmp_dw, nuebar_fractions.at(ax_it->first).at("nominal"));
301 
302  delete tmp_up;
303  delete tmp_dw;
304 
305  tmp_up = (TH1*) hists.at(ax_it->first).at("nue").at(syst + "_up")->Clone(UniqueName().c_str());
306  tmp_dw = (TH1*) hists.at(ax_it->first).at("nue").at(syst + "_dw")->Clone(UniqueName().c_str());
307  tmp_up->Divide(hists.at(ax_it->first).at("nuebar").at(syst + "_up"));
308  tmp_dw->Divide(hists.at(ax_it->first).at("nuebar").at(syst + "_dw"));
309 
310  nue_nuebar_ratios[ax_it->first][syst] = MostConservative(tmp_up, tmp_dw, nue_nuebar_ratios.at(ax_it->first).at("nominal"));
311  }
312 
313  // treat all beam systs as two sided even though connor said there are some that are one sided.
314  for(auto beam_it = beamSystematics.begin(); beam_it != beamSystematics.end(); beam_it++) {
315  TH1 * tmp_up = (TH1*) hists.at(ax_it->first).at("nuebar").at(beam_it->first + "_up")->Clone(UniqueName().c_str());
316  TH1 * tmp_dw = (TH1*) hists.at(ax_it->first).at("nuebar").at(beam_it->first + "_dw")->Clone(UniqueName().c_str());
317  tmp_up->Divide(hists.at(ax_it->first).at("nunubar").at(beam_it->first + "_up"));
318  tmp_dw->Divide(hists.at(ax_it->first).at("nunubar").at(beam_it->first + "_dw"));
319 
320  nuebar_fractions[ax_it->first][beam_it->first] = MostConservative(tmp_up, tmp_dw, nuebar_fractions.at(ax_it->first).at("nominal"));
321 
322  delete tmp_up;
323  delete tmp_dw;
324 
325  tmp_up = (TH1*) hists.at(ax_it->first).at("nue") .at(beam_it->first + "_up")->Clone(UniqueName().c_str());
326  tmp_dw = (TH1*) hists.at(ax_it->first).at("nue") .at(beam_it->first + "_dw")->Clone(UniqueName().c_str());
327  tmp_up->Divide(hists.at(ax_it->first).at("nuebar").at(beam_it->first + "_up"));
328  tmp_dw->Divide(hists.at(ax_it->first).at("nuebar").at(beam_it->first + "_dw"));
329 
330  nue_nuebar_ratios[ax_it->first][beam_it->first] = MostConservative(tmp_up, tmp_dw, nue_nuebar_ratios.at(ax_it->first).at("nominal"));
331 
332  }
333 
334  // cherenkov
335  nuebar_fractions.at(ax_it->first)["cherenkov"] = (TH1*) hists.at(ax_it->first).at("nuebar").at("cherenkov")->Clone(UniqueName().c_str());
336  nuebar_fractions.at(ax_it->first).at("cherenkov")->Divide(hists.at(ax_it->first).at("nunubar").at("cherenkov"));
337 
338  nue_nuebar_ratios.at(ax_it->first)["cherenkov"] = (TH1*) hists.at(ax_it->first).at("nue").at("cherenkov")->Clone(UniqueName().c_str());
339  nue_nuebar_ratios.at(ax_it->first).at("cherenkov")->Divide(hists.at(ax_it->first).at("nuebar").at("cherenkov"));
340 
341  // calib shape
342  nuebar_fractions.at(ax_it->first)["calib_shape"] = (TH1*) hists.at(ax_it->first).at("nuebar").at("calib_shape")->Clone(UniqueName().c_str());
343  nuebar_fractions.at(ax_it->first).at("calib_shape")->Divide(hists.at(ax_it->first).at("nunubar").at("calib_shape"));
344 
345  nue_nuebar_ratios.at(ax_it->first)["calib_shape"] = (TH1*) hists.at(ax_it->first).at("nue").at("calib_shape")->Clone(UniqueName().c_str());
346  nue_nuebar_ratios.at(ax_it->first).at("calib_shape")->Divide(hists.at(ax_it->first).at("nuebar").at("calib_shape"));
347 
348 
349  // multiverses
350  TH1 * ppfx_up = mv_ppfx_nuebar_fractions.at(ax_it->first)->GetPlusOneSigmaShift (nuebar_fractions.at(ax_it->first).at("nominal"));
351  TH1 * ppfx_dw = mv_ppfx_nuebar_fractions.at(ax_it->first)->GetMinusOneSigmaShift(nuebar_fractions.at(ax_it->first).at("nominal"));
352 
353  TH1 * genie_up = mv_genie_nuebar_fractions.at(ax_it->first)->GetPlusOneSigmaShift (nuebar_fractions.at(ax_it->first).at("nominal"));
354  TH1 * genie_dw = mv_genie_nuebar_fractions.at(ax_it->first)->GetMinusOneSigmaShift(nuebar_fractions.at(ax_it->first).at("nominal"));
355 
356  nuebar_fractions[ax_it->first]["ppfx"] = MostConservative(ppfx_up , ppfx_dw , nuebar_fractions.at(ax_it->first).at("nominal"));
357  nuebar_fractions[ax_it->first]["genie"] = MostConservative(genie_up, genie_dw, nuebar_fractions.at(ax_it->first).at("nominal"));
358 
359 
360  delete ppfx_up;
361  delete ppfx_dw;
362  delete genie_up;
363  delete genie_dw;
364 
365  ppfx_up = mv_ppfx_nue_nuebar_ratios.at(ax_it->first)->GetPlusOneSigmaShift (nue_nuebar_ratios.at(ax_it->first).at("nominal"));
366  ppfx_dw = mv_ppfx_nue_nuebar_ratios.at(ax_it->first)->GetMinusOneSigmaShift(nue_nuebar_ratios.at(ax_it->first).at("nominal"));
367  nue_nuebar_ratios.at(ax_it->first)["ppfx"] = MostConservative(ppfx_up, ppfx_dw, nue_nuebar_ratios.at(ax_it->first).at("nominal"));
368 
369  genie_up = mv_genie_nue_nuebar_ratios.at(ax_it->first)->GetPlusOneSigmaShift (nue_nuebar_ratios.at(ax_it->first).at("nominal"));
370  genie_dw = mv_genie_nue_nuebar_ratios.at(ax_it->first)->GetMinusOneSigmaShift(nue_nuebar_ratios.at(ax_it->first).at("nominal"));
371  nue_nuebar_ratios.at(ax_it->first)["genie"] = MostConservative(genie_up, genie_dw, nue_nuebar_ratios.at(ax_it->first).at("nominal"));
372  }
373 
374  for(auto it1 = nuebar_fractions.begin(); it1 != nuebar_fractions.end(); it1++)
375  for(auto it2 = nuebar_fractions.at(it1->first).begin(); it2 != nuebar_fractions.at(it1->first).end(); it2++)
376  it2->second->GetYaxis()->SetTitle("#bar{#nu_{e}} Fraction");
377 
378  // uncertainties
379  std::map<std::string, std::map<std::string, TH1*> > nuebar_fractions_uncert_abs;
380  std::map<std::string, std::map<std::string, TH1*> > nuebar_fractions_uncert_frac;
381  std::map<std::string, std::map<std::string, TH1*> > nue_nuebar_ratios_uncert_abs;
382  std::map<std::string, std::map<std::string, TH1*> > nue_nuebar_ratios_uncert_frac;
383  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
384  for(std::string syst : {"ppfx", "genie", "calib_shift", "calib_shape", "ll", "cherenkov"} ) {
385  nuebar_fractions_uncert_abs[ax_it->first][syst] = (TH1*) nuebar_fractions.at(ax_it->first).at(syst)->Clone(UniqueName().c_str());
386  nuebar_fractions_uncert_abs[ax_it->first][syst]->Add(nuebar_fractions.at(ax_it->first).at("nominal"), -1);
387  abs(nuebar_fractions_uncert_abs[ax_it->first][syst]);
388 
389  nuebar_fractions_uncert_frac[ax_it->first][syst] = (TH1*) nuebar_fractions_uncert_abs.at(ax_it->first).at(syst)->Clone(UniqueName().c_str());
390  nuebar_fractions_uncert_frac.at(ax_it->first).at(syst)->Divide(nuebar_fractions.at(ax_it->first).at("nominal"));
391 
392 
393  // ratios
394  nue_nuebar_ratios_uncert_abs[ax_it->first][syst] = (TH1*) nue_nuebar_ratios.at(ax_it->first).at(syst)->Clone(UniqueName().c_str());
395  nue_nuebar_ratios_uncert_abs[ax_it->first][syst]->Add(nue_nuebar_ratios.at(ax_it->first).at("nominal"), -1);
396  abs(nue_nuebar_ratios_uncert_abs[ax_it->first][syst]);
397 
398  nue_nuebar_ratios_uncert_frac[ax_it->first][syst] = (TH1*) nue_nuebar_ratios_uncert_abs.at(ax_it->first).at(syst)->Clone(UniqueName().c_str());
399  nue_nuebar_ratios_uncert_frac.at(ax_it->first).at(syst)->Divide(nue_nuebar_ratios.at(ax_it->first).at("nominal"));
400  }
401 
402  // beam uncertainties
403  for(auto beam_it = beamSystematics.begin(); beam_it != beamSystematics.end(); beam_it++) {
404  std::string syst = beam_it->first;
405  nuebar_fractions_uncert_abs[ax_it->first][syst] = (TH1*) nuebar_fractions.at(ax_it->first).at(syst)->Clone(UniqueName().c_str());
406  nuebar_fractions_uncert_abs[ax_it->first][syst]->Add(nuebar_fractions.at(ax_it->first).at("nominal"), -1);
407  abs(nuebar_fractions_uncert_abs[ax_it->first][syst]);
408 
409  nuebar_fractions_uncert_frac[ax_it->first][syst] = (TH1*) nuebar_fractions_uncert_abs.at(ax_it->first).at(syst)->Clone(UniqueName().c_str());
410  nuebar_fractions_uncert_frac.at(ax_it->first).at(syst)->Divide(nuebar_fractions.at(ax_it->first).at("nominal"));
411 
412  // ratios
413  nue_nuebar_ratios_uncert_abs[ax_it->first][syst] = (TH1*) nue_nuebar_ratios.at(ax_it->first).at(syst)->Clone(UniqueName().c_str());
414  nue_nuebar_ratios_uncert_abs[ax_it->first][syst]->Add(nue_nuebar_ratios.at(ax_it->first).at("nominal"), -1);
415  abs(nue_nuebar_ratios_uncert_abs[ax_it->first][syst]);
416 
417  nue_nuebar_ratios_uncert_frac[ax_it->first][syst] = (TH1*) nue_nuebar_ratios_uncert_abs.at(ax_it->first).at(syst)->Clone(UniqueName().c_str());
418  nue_nuebar_ratios_uncert_frac.at(ax_it->first).at(syst)->Divide(nue_nuebar_ratios.at(ax_it->first).at("nominal"));
419  }
420 
421  // combine beam uncertainties into one
422  nuebar_fractions_uncert_frac.at(ax_it->first)["beam_focusing"] = (TH1*) nuebar_fractions_uncert_frac.at(ax_it->first).at("2kA")->Clone(UniqueName().c_str());
423  nuebar_fractions_uncert_frac.at(ax_it->first).at("beam_focusing")->Reset();
424 
425  nue_nuebar_ratios_uncert_frac.at(ax_it->first)["beam_focusing"] = (TH1*) nue_nuebar_ratios_uncert_frac.at(ax_it->first).at("2kA")->Clone(UniqueName().c_str());
426  nue_nuebar_ratios_uncert_frac.at(ax_it->first).at("beam_focusing")->Reset();
427 
428  auto NX = nuebar_fractions_uncert_frac.at(ax_it->first).at("beam_focusing")->GetNbinsX();
429  auto NY = nuebar_fractions_uncert_frac.at(ax_it->first).at("beam_focusing")->GetNbinsY();
430  auto NZ = nuebar_fractions_uncert_frac.at(ax_it->first).at("beam_focusing")->GetNbinsZ();
431  for(auto x = 1; x <= NX; x++) {
432  for(auto y = 1; y <= NY; y++) {
433  for(auto z = 1; z <= NZ; z++) {
434  double rat = 0;
435  double val = 0;
436  for(auto beam_it = beamSystematics.begin(); beam_it != beamSystematics.end(); beam_it++) {
437  rat += std::pow(nue_nuebar_ratios_uncert_frac.at(ax_it->first).at(beam_it->first)->GetBinContent(x, y, z), 2);
438  val += std::pow(nuebar_fractions_uncert_frac.at(ax_it->first).at(beam_it->first)->GetBinContent(x, y, z), 2);
439  }
440  nue_nuebar_ratios_uncert_frac.at(ax_it->first).at("beam_focusing")->SetBinContent(x, y, z, std::sqrt(rat));
441  nuebar_fractions_uncert_frac.at(ax_it->first).at("beam_focusing")->SetBinContent(x, y, z, std::sqrt(val));
442  }
443  }
444  }
445 
446  // MC statistics
447  nuebar_fractions_uncert_abs .at(ax_it->first)["stat"] = (TH1*) nuebar_fractions.at(ax_it->first).at("nominal")->Clone(UniqueName().c_str());
448  nuebar_fractions_uncert_frac.at(ax_it->first)["stat"] = (TH1*) nuebar_fractions.at(ax_it->first).at("nominal")->Clone(UniqueName().c_str());
449 
450  NX = nuebar_fractions.at(ax_it->first).at("nominal")->GetNbinsX();
451  NY = nuebar_fractions.at(ax_it->first).at("nominal")->GetNbinsY();
452  NZ = nuebar_fractions.at(ax_it->first).at("nominal")->GetNbinsZ();
453  for(auto ix = 1; ix <= NX; ix++) {
454  for(auto iy = 1; iy <= NY; iy++) {
455  for(auto iz = 1; iz <= NZ; iz++) {
456  // rescale back to full MC statistics
457 
458  auto num = hists.at(ax_it->first).at("nuebar") .at("nominal")->GetBinContent(ix, iy, iz) * MCPOT / kAna2020RHCPOT;
459  auto den = hists.at(ax_it->first).at("nunubar").at("nominal")->GetBinContent(ix, iy, iz) * MCPOT / kAna2020RHCPOT;
460  auto err = std::sqrt(1/num + 1/den);
461  auto nom = num / den;
462  nuebar_fractions_uncert_abs .at(ax_it->first).at("stat")->SetBinContent(ix, iy, iz, nom * err);
463  nuebar_fractions_uncert_frac.at(ax_it->first).at("stat")->SetBinContent(ix, iy, iz, err );
464 
465  }
466  }
467  }
468  }
469 
470  // sum in quadrature
471  std::map<std::string, TH1*> nsel_uncert;
472  std::map<std::string, TH1*> mc_frac_uncert;
473  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
474  nsel_uncert [ax_it->first] = (TH1*) hists.at(ax_it->first).at("nuebar").at("nominal")->Clone(UniqueName().c_str());
475  mc_frac_uncert[ax_it->first] = (TH1*) nuebar_fractions.at(ax_it->first).at("nominal")->Clone(UniqueName().c_str());
476  nsel_uncert .at(ax_it->first)->Reset();
477  mc_frac_uncert.at(ax_it->first)->Reset();
478  mc_frac_uncert.at(ax_it->first)->SetTitle("#bar{#nu}_{e} Fraction Uncertainty");
479  nsel_uncert .at(ax_it->first)->SetTitle("Selected #bar{#nu}_{e} Uncertainty");
480 
481  if(ax_it->first == "electron_kin")
482  mc_frac_uncert.at(ax_it->first)->GetZaxis()->SetTitle("Fractional Uncertainty");
483  else
484  mc_frac_uncert.at(ax_it->first)->GetYaxis()->SetTitle("Fractional Uncertainty");
485 
486  auto NX = nuebar_fractions.at(ax_it->first).at("nominal")->GetNbinsX();
487  auto NY = nuebar_fractions.at(ax_it->first).at("nominal")->GetNbinsY();
488  auto NZ = nuebar_fractions.at(ax_it->first).at("nominal")->GetNbinsZ();
489 
490  for(auto ix = 1; ix <= NX; ix++) {
491  for(auto iy = 1; iy <= NY; iy++) {
492  for(auto iz = 1; iz <= NZ; iz++) {
493  double val = 0;
494  for(std::string syst : {"ppfx", "genie", "calib_shape", "calib_shift", "ll", "cherenkov", "beam_focusing", "stat"}) {
495  val += std::pow(nuebar_fractions_uncert_frac.at(ax_it->first).at(syst)->GetBinContent(ix, iy, iz), 2);
496  }// syst
497  mc_frac_uncert.at(ax_it->first)->SetBinContent(ix, iy, iz, std::sqrt(val));
498 
499  auto n = hists.at(ax_it->first).at("nunubar").at("nominal")->GetBinContent(ix, iy, iz);
500  nsel_uncert.at(ax_it->first)->SetBinContent(ix, iy, iz, std::sqrt(1/n + val));
501 
502  } // iz
503  } // iy
504  } // ix
505  } // ax_it
506 
507  mc_frac_uncert.at("enu")->GetXaxis()->SetRangeUser(0, 6);
508  mc_frac_uncert.at("q2" )->GetXaxis()->SetRangeUser(0, 2);
509  nue_nuebar_ratios.at("enu").at("nominal")->GetXaxis()->SetRangeUser(0, 6);
510  nue_nuebar_ratios.at("q2").at("nominal") ->GetXaxis()->SetRangeUser(0, 2);
511 
512  Plot1D(mc_frac_uncert.at("enu"), "Total", nuebar_fractions_uncert_frac.at("enu"), all_systs, plot_dump + "/enu_fractional_uncertainties_all.pdf");
513  Plot1D(mc_frac_uncert.at("q2") , "Total", nuebar_fractions_uncert_frac.at("q2") , all_systs, plot_dump + "/q2_fractional_uncertainties_all.pdf" );
514  Plot1D(mc_frac_uncert.at("enu"), "Total", nuebar_fractions_uncert_frac.at("enu"), beam_systs, plot_dump + "/enu_fractional_uncertainties_beam.pdf");
515  Plot1D(mc_frac_uncert.at("q2") , "Total", nuebar_fractions_uncert_frac.at("q2") , beam_systs, plot_dump + "/q2_fractional_uncertainties_beam.pdf" );
516  for(auto ax_it = axes.begin(); ax_it != axes.end(); ax_it++) {
517  // Plot1D(nuebar_fractions .at(ax_it->first).at("nominal"), "Nominal", nuebar_fractions.at(ax_it->first) , all_systs, plot_dump + "/" + ax_it->first + "_nuebar_fraction_all.pdf" , xmax.at(ax_it->first));
518  // Plot1D(nue_nuebar_ratios.at(ax_it->first).at("nominal"), "Nominal", nue_nuebar_ratios.at(ax_it->first), all_systs, plot_dump + "/" + ax_it->first + "_nue_nuebar_ratios_all.pdf", xmax.at(ax_it->first));
519  //
520  // Plot1D(nuebar_fractions .at(ax_it->first).at("nominal"), "Nominal", nuebar_fractions.at(ax_it->first) , beam_systs, plot_dump + "/" + ax_it->first + "_nuebar_fraction_beam.pdf" , xmax.at(ax_it->first));
521  // Plot1D(nue_nuebar_ratios.at(ax_it->first).at("nominal"), "Nominal", nue_nuebar_ratios.at(ax_it->first), beam_systs, plot_dump + "/" + ax_it->first + "_nue_nuebar_ratios_beam.pdf", xmax.at(ax_it->first));
522 
523  // Plot1D(hists.at(ax_it->first).at("nuebar" ).at("nominal"), "Nominal", hists.at(ax_it->first).at("nuebar" ), plot_dump + "/" + ax_it->first + "_nuebar_events.pdf" , xmax.at(ax_it->first));
524  // Plot1D(hists.at(ax_it->first).at("nue" ).at("nominal"), "Nominal", hists.at(ax_it->first).at("nue" ), plot_dump + "/" + ax_it->first + "_nue_events.pdf" , xmax.at(ax_it->first));
525  // Plot1D(hists.at(ax_it->first).at("nunubar").at("nominal"), "Nominal", hists.at(ax_it->first).at("nunubar"), plot_dump + "/" + ax_it->first + "_nunubar_events.pdf", xmax.at(ax_it->first));
526 
527  PlotMultiverse(*mv_ppfx .at(ax_it->first).at("nuebar" ), hists.at(ax_it->first).at("nuebar" ).at("nominal"), "Selected #bar{#nu_{e}} PPFX" , plot_dump + "/" + ax_it->first + "_nuebar_events_ppfx.pdf" );
528  PlotMultiverse(*mv_genie.at(ax_it->first).at("nuebar" ), hists.at(ax_it->first).at("nuebar" ).at("nominal"), "Selected #bar{#nu_{e}} GENIE", plot_dump + "/" + ax_it->first + "_nuebar_events_genie.pdf");
529  PlotMultiverse(*mv_ppfx .at(ax_it->first).at("nue" ), hists.at(ax_it->first).at("nue" ).at("nominal"), "Selected #nu_{e} PPFX" , plot_dump + "/" + ax_it->first + "_nue_events_ppfx.pdf" );
530  PlotMultiverse(*mv_genie.at(ax_it->first).at("nue" ), hists.at(ax_it->first).at("nue" ).at("nominal"), "Selected #nu_{e} GENIE", plot_dump + "/" + ax_it->first + "_nue_events_genie.pdf");
531  PlotMultiverse(*mv_ppfx .at(ax_it->first).at("nunubar"), hists.at(ax_it->first).at("nunubar").at("nominal"), "Selected #bar{#nu_{e}} + #nu_{e} PPFX" , plot_dump + "/" + ax_it->first + "_nunubar_events_ppfx.pdf" );
532  PlotMultiverse(*mv_genie.at(ax_it->first).at("nunubar"), hists.at(ax_it->first).at("nunubar").at("nominal"), "Selected #bar{#nu_{e}} + #nu_{e} GENIE", plot_dump + "/" + ax_it->first + "_nunubar_events_genie.pdf");
533 
534  PlotMultiverse(*mv_ppfx_nuebar_fractions .at(ax_it->first), nuebar_fractions.at(ax_it->first).at("nominal"), "#bar{#nu_{e}} Fraction PPFX" , plot_dump + "/" + ax_it->first + "_nuebar_fraction_ppfx.pdf",
535  nuebar_fractions.at(ax_it->first).at("ppfx"));
536 
537  PlotMultiverse(*mv_genie_nuebar_fractions.at(ax_it->first), nuebar_fractions.at(ax_it->first).at("nominal"), "#bar{#nu_{e}} Fraction GENIE", plot_dump + "/" + ax_it->first + "_nuebar_fraction_genie.pdf",
538  nuebar_fractions.at(ax_it->first).at("genie"));
539 
540  PlotMultiverse(*mv_ppfx_nue_nuebar_ratios.at(ax_it->first), nue_nuebar_ratios.at(ax_it->first).at("nominal"), "#nu_{e} / #bar{#nu}_{e} PPFX", plot_dump + "/" + ax_it->first + "_nue_nuebar_ratio_ppfx.pdf",
541  nue_nuebar_ratios.at(ax_it->first).at("ppfx"));
542 
543  PlotMultiverse(*mv_genie_nue_nuebar_ratios.at(ax_it->first), nue_nuebar_ratios.at(ax_it->first).at("nominal"), "#nu_{e} / #bar{#nu}_{e} GENIE", plot_dump + "/" + ax_it->first + "_nue_nuebar_ratio_genie.pdf",
544  nue_nuebar_ratios.at(ax_it->first).at("genie"));
545 
546 
547  }
548  gStyle->SetPaintTextFormat(".2f");
549  Plot2D(mc_frac_uncert.at("electron_kin"), nuebar_fractions_uncert_frac.at("electron_kin"), all_systs, plot_dump + "/electron_kin_fractional_uncertainty");
550 
551  TCanvas * c = new TCanvas();
552  nsel_uncert.at("q2")->GetYaxis()->SetTitle("Fractional Uncertainty");
553  nsel_uncert.at("q2")->GetXaxis()->SetRangeUser(0, xmax.at("q2"));
554  nsel_uncert.at("q2")->Draw("hist");
555  c->Print(TString::Format("%s/%s_nsel_uncert.pdf", plot_dump.c_str(), "q2").Data());
556 
557  nsel_uncert.at("enu")->GetYaxis()->SetTitle("Fractional Uncertainty");
558  nsel_uncert.at("enu")->GetXaxis()->SetRangeUser(0, xmax.at("enu"));
559  nsel_uncert.at("enu")->Draw("hist");
560  c->Print(TString::Format("%s/%s_nsel_uncert.pdf", plot_dump.c_str(), "enu").Data());
561 
562 
563  c->SetRightMargin(0.2);
564  auto nsel_uncert_electron_kin_2d = ToTH2Helper(std::unique_ptr<TH1>((TH1*) nsel_uncert.at("electron_kin")->Clone(UniqueName().c_str())),
565  nuebarccinc::costhetabins,
566  nuebarccinc::eelecbins);
567  nsel_uncert_electron_kin_2d->GetZaxis()->SetTitle("Fractional Uncertainty");
568 
569  nsel_uncert_electron_kin_2d->GetXaxis()->SetRangeUser(0.85, 1);
570  nsel_uncert_electron_kin_2d->GetYaxis()->SetRangeUser(0, 6);
571  nsel_uncert_electron_kin_2d->GetZaxis()->SetRangeUser(0, 1);
572  nsel_uncert_electron_kin_2d->SetMarkerColor(kWhite);
573  nsel_uncert_electron_kin_2d->GetXaxis()->SetTitle("True cos#theta_{e}");
574  nsel_uncert_electron_kin_2d->GetYaxis()->SetTitle("True Electron Energy (GeV)");
575  MovePalette(nsel_uncert_electron_kin_2d);
576  nsel_uncert_electron_kin_2d->Draw("colz text");
577 
578  c->Print(TString::Format("%s/%s_nsel_uncert.pdf", plot_dump.c_str(), "electron_kin").Data());
579 
580  auto nnunubar_2d = ToTH2Helper(std::unique_ptr<TH1>((TH1*) hists.at("electron_kin").at("nunubar").at("nominal")->Clone(UniqueName().c_str())),
581  nuebarccinc::costhetabins,
582  nuebarccinc::eelecbins);
583 
584  nnunubar_2d->SetTitle("Selected N_{#bar{#nu}_{e}} + N_{#nu_{e}}");
585  nnunubar_2d->GetXaxis()->SetRangeUser(0.85, 1);
586  nnunubar_2d->GetYaxis()->SetRangeUser(0, 6);
587  nnunubar_2d->SetMarkerColor(kWhite);
588  nnunubar_2d->GetXaxis()->SetTitle("True cos#theta_{e}");
589  nnunubar_2d->GetYaxis()->SetTitle("True Electron Energy (GeV)");
590  MovePalette(nnunubar_2d);
591  nnunubar_2d->Draw("colz text");
592 
593  c->Print(TString::Format("%s/%s_nsel_nunubar.pdf", plot_dump.c_str(), "electron_kin").Data());
594 
595  auto nom_nuebar_frac_2d = ToTH2Helper(std::unique_ptr<TH1>((TH1*) nuebar_fractions.at("electron_kin").at("nominal")->Clone(UniqueName().c_str())),
596  nuebarccinc::costhetabins,
597  nuebarccinc::eelecbins);
598 
599  nom_nuebar_frac_2d->SetTitle("s_{#bar{#nu}_{e}}");
600  nom_nuebar_frac_2d->GetXaxis()->SetRangeUser(0.85, 1);
601  nom_nuebar_frac_2d->GetYaxis()->SetRangeUser(0, 6);
602  nom_nuebar_frac_2d->GetZaxis()->SetRangeUser(0, 1);
603  nom_nuebar_frac_2d->SetMarkerColor(kWhite);
604  nom_nuebar_frac_2d->GetXaxis()->SetTitle("True cos#theta_{e}");
605  nom_nuebar_frac_2d->GetYaxis()->SetTitle("True Electron Energy (GeV)");
606  MovePalette(nom_nuebar_frac_2d);
607  nom_nuebar_frac_2d->Draw("colz text");
608 
609  c->Print(TString::Format("%s/%s_nuebar_fraction_2d.pdf", plot_dump.c_str(), "electron_kin").Data());
610 
611 
612 
613  // calculate required MC POT that would make MC stat. uncert
614  // negligible -- ~1%
615 
616  std::map<std::string, TH1*> req_pot_1d;
617  std::map<std::string, TH1*> req_pot_2d;
618  std::map<std::string, TH1*> req_pot_wrt_prod5;
619  auto target_uncert = 0.01;
620  for(auto loader_it = loaders.begin(); loader_it != loaders.end(); loader_it++) {
621  req_pot_1d[loader_it->first] = (TH1*) nuebar_fractions.at("electron_kin").at(loader_it->first)->Clone(UniqueName().c_str());
622  req_pot_1d.at(loader_it->first)->Reset();
623 
624  for(auto i = 1; i <= req_pot_1d.at(loader_it->first)->GetNbinsX(); i++) {
625  double s = nuebar_fractions.at("electron_kin").at(loader_it->first)->GetBinContent(i);
626  double n = hists.at("electron_kin").at("nunubar").at(loader_it->first)->GetBinContent(i);
627  double pot = pots.at(loader_it->first);
628 
629  if(s != 0. && s != 1)
630  req_pot_1d.at(loader_it->first)->SetBinContent(i, (1+s)/s * std::pow(target_uncert, -2) * pot / n);
631  }
632 
633 
634  req_pot_2d[loader_it->first] = ToTH2Helper(std::unique_ptr<TH1>((TH1*) req_pot_1d.at(loader_it->first)->Clone(UniqueName().c_str())),
635  nuebarccinc::costhetabins,
636  nuebarccinc::eelecbins);
637  req_pot_wrt_prod5[loader_it->first] = (TH1*) req_pot_2d.at(loader_it->first)->Clone(UniqueName().c_str());
638  req_pot_wrt_prod5.at(loader_it->first)->Scale(1./pots.at(loader_it->first));
639 
640 
641  gStyle->SetPaintTextFormat("1.1g");
642  gStyle->SetPalette(kDeepSea);
643  c->SetLogz();
644  req_pot_2d.at(loader_it->first)->SetTitle(TString::Format("Required MC POT (%.1f%% Stat. Uncert)", target_uncert * 100));
645  req_pot_2d.at(loader_it->first)->GetXaxis()->SetRangeUser(0.85, 1);
646  req_pot_2d.at(loader_it->first)->GetYaxis()->SetRangeUser(0, 6);
647  req_pot_2d.at(loader_it->first)->SetMarkerColor(kWhite);
648  req_pot_2d.at(loader_it->first)->GetXaxis()->SetTitle("True cos#theta_{e}");
649  req_pot_2d.at(loader_it->first)->GetYaxis()->SetTitle("True Electron Energy (GeV)");
650  req_pot_2d.at(loader_it->first)->GetZaxis()->SetRangeUser(1e15, 1e25);
651  MovePalette(req_pot_2d.at(loader_it->first));
652  req_pot_2d.at(loader_it->first)->Draw("colz text");
653 
654 
655  c->Print(TString::Format("%s/req_pot_%s.pdf", plot_dump.c_str(), loader_it->first.c_str()).Data());
656 
657  req_pot_wrt_prod5.at(loader_it->first)->SetTitle("Required MC POT / Prod5 POT");
658  req_pot_wrt_prod5.at(loader_it->first)->GetXaxis()->SetRangeUser(0.85, 1);
659  req_pot_wrt_prod5.at(loader_it->first)->GetYaxis()->SetRangeUser(0, 6);
660  req_pot_wrt_prod5.at(loader_it->first)->SetMaximum(500);
661  req_pot_wrt_prod5.at(loader_it->first)->SetMarkerColor(kWhite);
662  req_pot_wrt_prod5.at(loader_it->first)->GetXaxis()->SetTitle("True cos#theta_{e}");
663  req_pot_wrt_prod5.at(loader_it->first)->GetYaxis()->SetTitle("True Electron Energy (GeV)");
664  MovePalette(req_pot_wrt_prod5.at(loader_it->first));
665  req_pot_wrt_prod5.at(loader_it->first)->Draw("colz text");
666 
667 
668  gStyle->SetPaintTextFormat(".1f");
669  c->SetLogz(false);
670  c->Print(TString::Format("%s/req_pot_wrt_prod5%s.pdf", plot_dump.c_str(), loader_it->first.c_str()).Data());
671  }
672 
673  delete c;
674  c = new TCanvas();
675  TLegend * leg = 0;
676  auto kNuebarColor = kGreen;
677  auto kNueColor = kMagenta;
678  // 1D stacked nnunubar plots
679  for(auto loader_it = loaders.begin(); loader_it != loaders.end(); loader_it++) {
680  for(std::string axis : {"enu", "q2"}) {
681  if(leg)
682  delete leg;
683  leg = new TLegend(0.6, 0.89, 0.8, 0.7);
684 
685  THStack * stack = new THStack(UniqueName().c_str(), "");
686  hists.at(axis).at("nue") .at(loader_it->first)->SetFillColor(kNueColor);
687  hists.at(axis).at("nuebar").at(loader_it->first)->SetFillColor(kNuebarColor);
688  hists.at(axis).at("nue") .at(loader_it->first)->SetLineColor(kBlack);
689  hists.at(axis).at("nuebar").at(loader_it->first)->SetLineColor(kBlack);
690  hists.at(axis).at("nue") .at(loader_it->first)->SetLineStyle(1);
691  hists.at(axis).at("nuebar").at(loader_it->first)->SetLineStyle(1);
692 
693  stack->Add(hists.at(axis).at("nue" ).at(loader_it->first));
694  stack->Add(hists.at(axis).at("nuebar").at(loader_it->first));
695  leg->AddEntry(hists.at(axis).at("nue" ).at(loader_it->first), "#nu_{e}" , "f");
696  leg->AddEntry(hists.at(axis).at("nuebar").at(loader_it->first), "#bar{#nu}_{e}", "f");
697 
698 
699  stack->Draw("hist");
700  stack->GetXaxis()->SetRangeUser(0, xmax.at(axis));
701  stack->GetXaxis()->SetTitle(hists.at(axis).at("nue").at(loader_it->first)->GetXaxis()->GetTitle());
702  stack->GetYaxis()->SetTitle("Events");
703  leg->Draw("same");
704  c->Print(TString::Format("%s/%s_nunubar_stacked_%s.pdf",
705  plot_dump.c_str(),
706  axis.c_str(),
707  loader_it->first.c_str()).Data());
708 
709  }
710  }
711 
712 
713  }
714 }
output
ofstream output
Definition: plot_prediction_systs.C:32
ana::HistAxis
_HistAxis< Var > HistAxis
Definition: HistAxis.h:103
xmax
std::map< std::string, double > xmax
Definition: wrong_sign_uncertainty.C:79
ana::Multiverse::GetPlusOneSigmaShift
Spectrum * GetPlusOneSigmaShift(const Spectrum *)
Definition: Multiverse.cxx:443
nuebarccinc::eelecbins
const ana::Binning eelecbins
Definition: NuebarCCIncBins.h:19
ana::kPPFXFluxCVWgt
const Var kPPFXFluxCVWgt
Definition: PPFXWeights.h:16
nuebarccinc::kNuebarCC
const ana::Cut kNuebarCC
Definition: NuebarCCIncCuts.h:87
pots
int pots
Definition: common_convert.C:43
ana::SystShifts
Simple record of shifts applied to systematic parameters.
Definition: SystShifts.h:20
beam_systs
std::vector< std::string > beam_systs
Definition: wrong_sign_uncertainty.C:48
den
Float_t den
Definition: plot.C:36
nuebarccinc::costhetabins
const ana::Binning costhetabins
Definition: NuebarCCIncBins.h:23
std::sqrt
T sqrt(T number)
Definition: d0nt_math.hpp:156
cet::pow
constexpr T pow(T x)
Definition: pow.h:75
Plot2D
void Plot2D(TH1 *total, std::map< std::string, TH1 * > systs, std::vector< std::string > to_plot, std::string basename)
Definition: wrong_sign_uncertainty.C:760
nuebarccinc::PROD5_MC_RHC_LIGHT_UP
const std::string PROD5_MC_RHC_LIGHT_UP
Definition: NuebarCCIncExtra.h:22
abs
void abs(TH1 *hist)
Definition: wrong_sign_uncertainty.C:716
nuebarccinc::PROD5_MC_RHC_CALIB_DOWN
const std::string PROD5_MC_RHC_CALIB_DOWN
Definition: NuebarCCIncExtra.h:20
hists
TString hists[nhists]
Definition: bdt_com.C:3
nuebarccinc::kTrueElectronEVsCosStandardAxis
const ana::HistAxis kTrueElectronEVsCosStandardAxis
Definition: NuebarCCIncVars.h:175
nuebarccinc::PROD5_MC_RHC_LIGHT_DOWN
const std::string PROD5_MC_RHC_LIGHT_DOWN
Definition: NuebarCCIncExtra.h:23
ana::LoadFrom
std::unique_ptr< T > LoadFrom(TDirectory *dir, const std::string &label)
Definition: LoadFromFile.h:17
ana::Multiverse::GetMinusOneSigmaShift
Spectrum * GetMinusOneSigmaShift(const Spectrum *)
Definition: Multiverse.cxx:450
ana::GetkPPFXFluxUnivWgt
std::vector< Var > GetkPPFXFluxUnivWgt()
Definition: PPFXWeights.cxx:53
submit_syst.y
y
Definition: submit_syst.py:249
s
const XML_Char * s
Definition: expat.h:262
ana::getAllXsecSysts_2020
std::vector< const ISyst * > getAllXsecSysts_2020()
Definition: XSecSystLists.cxx:235
Cut
void Cut(double x)
Definition: plot_outliers.C:1
submit_syst.x
x
Definition: submit_syst.py:249
MovePalette
void MovePalette(TH1 *hist)
Definition: wrong_sign_uncertainty.C:841
load_libs_muonid
void load_libs_muonid()
Definition: wrong_sign_uncertainty.C:32
Plot1D
void Plot1D(TH1 *h1, std::string h1_label, std::map< std::string, TH1 * > systs, std::vector< std::string > to_plot, std::string name, double xmax=-1)
Definition: wrong_sign_uncertainty.C:790
nuebarccinc::PROD5_MC_RHC_NOMINAL
const std::string PROD5_MC_RHC_NOMINAL
Definition: NuebarCCIncExtra.h:18
nuebarccinc::PROD5_MC_RHC_CALIB_UP
const std::string PROD5_MC_RHC_CALIB_UP
Definition: NuebarCCIncExtra.h:19
nuebarccinc::PROD5_MC_RHC_CALIB_SHAPE
const std::string PROD5_MC_RHC_CALIB_SHAPE
Definition: NuebarCCIncExtra.h:21
plot_validation_datamc.c
c
Definition: plot_validation_datamc.py:198
pot
#define pot
Definition: template_basic.C:49
cmf::Spectrum
std::vector< float > Spectrum
Definition: Constants.h:610
test.z
z
Definition: test.py:28
nuebarccinc::kPreselectionLoose
const ana::Cut kPreselectionLoose
Definition: NuebarCCIncCuts.h:57
ana::kAna2020RHCPOT
const double kAna2020RHCPOT
Definition: Exposures.h:235
nuebarccinc::PROD5_MC_RHC_CHERENKOV
const std::string PROD5_MC_RHC_CHERENKOV
Definition: NuebarCCIncExtra.h:24
ana::kNoShift
const SystShifts kNoShift
Definition: SystShifts.cxx:22
nuebarccinc::kTrueQ2StandardAxis
const ana::HistAxis kTrueQ2StandardAxis
Definition: NuebarCCIncVars.h:203
ana::kSelection
const Cut kSelection
Definition: SINCpi0_Cuts.h:328
ana::SpectrumLoader
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
Definition: SpectrumLoader.h:21
num
int num
Definition: f2_nu.C:119
loaders
std::vector< Loaders * > loaders
Definition: syst_header.h:386
ana::ToTH2Helper
TH2 * ToTH2Helper(std::unique_ptr< TH1 > h1, const Binning &binsx, const Binning &binsy, ana::EBinType bintype)
Helper for ana::ToTH2.
Definition: UtilsExt.cxx:135
nuebarccinc::kMuonIDProd4Cut
const ana::Cut kMuonIDProd4Cut
Definition: datamc.C:28
file_size_ana.axes
axes
Definition: file_size_ana.py:282
ana::kStandardSpillCuts
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
update_sam_good_runs_metadata.cuts
string cuts
Definition: update_sam_good_runs_metadata.py:10
genie::utils::style::Format
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
nuebarccinc::kNueCC
const ana::Cut kNueCC
Definition: NuebarCCIncCuts.h:86
PlotMultiverse
void PlotMultiverse(Multiverse &mv, TH1 *nominal, std::string title, std::string name, TH1 *error=0)
Definition: wrong_sign_uncertainty.C:728
ana::_Cut< caf::SRProxy >
ana::_Var< caf::SRProxy >
nuebarccinc::kTrueNeutrinoEnergyStandardAxis
const ana::HistAxis kTrueNeutrinoEnergyStandardAxis
Definition: NuebarCCIncVars.h:196
MostConservative
TH1 * MostConservative(TH1 *up, TH1 *down, TH1 *nominal)
Definition: wrong_sign_uncertainty.C:821
kGreen
enum BeamMode kGreen
make_plots
void make_plots(TFile *f, TH2 *hVsRun, TGraph *g, int run, std::string suffix, double xpos, TH2 *&cut)
Definition: plot_cosmic_eff.C:38
ana::UniqueName
std::string UniqueName()
Return a different string each time, for creating histograms.
Definition: Utilities.cxx:29
febshutoff_auto.val
val
Definition: febshutoff_auto.py:277
ana::Multiverse::Divide
void Divide(Multiverse &, bool=false)
Definition: Multiverse.cxx:139
all_systs
std::vector< std::string > all_systs
Definition: wrong_sign_uncertainty.C:47
beamSystematics
std::map< string, const ana::BeamSyst * > beamSystematics
Definition: wrong_sign_uncertainty.C:83
ana::kXSecCVWgt2020
const Var kXSecCVWgt2020
Definition: XsecTunes.h:106
string
enum BeamMode string

Variable Documentation

std::vector<std::string> all_systs = {"ppfx", "genie", "calib_shape", "calib_shift", "ll", "cherenkov", "beam_focusing", "stat"}

Definition at line 47 of file wrong_sign_uncertainty.C.

Referenced by make_cm_pullterm_pdf(), MakeCovarSim(), and wrong_sign_uncertainty().

std::vector<std::string> beam_systs = {"ppfx", "beam_focusing"}

Definition at line 48 of file wrong_sign_uncertainty.C.

Referenced by wrong_sign_uncertainty().

std::map<string, const ana::BeamSyst*> beamSystematics
Initial value:
= {
{"2kA", &ana::kBeamHornCurrent},
{"02mmBeamSpotSize", &ana::kBeamSpotSize},
{"1mmBeamShiftX", &ana::kBeamPosX},
{"1mmBeamShiftY", &ana::kBeamPosY},
{"3mmHorn1X", &ana::kBeamH1PosX},
{"3mmHorn1Y", &ana::kBeamH1PosY},
{"3mmHorn2X", &ana::kBeamH2PosX},
{"3mmHorn2Y", &ana::kBeamH2PosY},
{"7mmTargetZ", &ana::kBeamTarget},
{"MagneticFieldinDecayPipe", &ana::kBeamMagField},
{"1mmHornWater", &ana::kBeamGeomWater}
}
ana::kBeamTarget
const BeamSyst kBeamTarget((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"7mmTargetZ","+/- 7mm Target z Position")
Target z position shift +/-7mm.
Definition: BeamSysts.h:113
ana::kBeamHornCurrent
const BeamSyst kBeamHornCurrent((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"2kA","2kAHornCurrent","+/- 2kA Horn Current")
Horn Current +/-2kA.
Definition: BeamSysts.h:97
ana::kBeamGeomWater
const BeamSyst kBeamGeomWater((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"1mmHornWater","+/- 1mm water on Horn 1")
Water layer on horn 1: +/- 1mm.
Definition: BeamSysts.h:128
ana::kBeamH1PosX
const BeamSyst kBeamH1PosX((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"3mmHorn1X","+/- 3mm Horn 1 X Position")
Horn 1 and 2 position +-3mm in X and Y separately.
Definition: BeamSysts.h:107
ana::kBeamPosY
const BeamSyst kBeamPosY((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"1mmBeamShiftY","Beam Position Y")
Definition: BeamSysts.h:104
ana::kBeamH1PosY
const BeamSyst kBeamH1PosY((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"3mmHorn1Y","+/- 3mm Horn 1 Y Position")
Definition: BeamSysts.h:108
ana::kBeamSpotSize
const BeamSyst kBeamSpotSize((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"0.2mmBeamSpotSize","0p2mmBeamSpotSize"," 1.3 +/- 0.2 mm Spot Size")
Beam Spot Size 1.3 +/- 0.2 mm both XY.
Definition: BeamSysts.h:100
ana::kBeamH2PosX
const BeamSyst kBeamH2PosX((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"3mmHorn2X","+/- 3mm Horn 2 X Position")
Definition: BeamSysts.h:109
ana::kBeamPosX
const BeamSyst kBeamPosX((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"1mmBeamShiftX","Beam Position X")
Beam position on target +-1 mm, X/Y separately.
Definition: BeamSysts.h:103
ana::kBeamH2PosY
const BeamSyst kBeamH2PosY((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"3mmHorn2Y","+/- 3mm Horn 2 Y Position")
Definition: BeamSysts.h:110
ana::kBeamMagField
const BeamSyst kBeamMagField((FindCAFAnaDir()+"/data/beam/TABeamSyst_2018Dec11.root"),"Magnetic Field in Decay Pipe","MagnFieldDecayPipe","Magnetic Field in Decay Pipe")
Constant magnetic field in decay pipe.
Definition: BeamSysts.h:116

Definition at line 83 of file wrong_sign_uncertainty.C.

Referenced by wrong_sign_uncertainty().

std::map<std::string, nuebarccinc::Style> syst_styles
Initial value:
= {{"genie" , {"", kCyan , 1}},
{"ppfx" , {"", kMagenta , 1}},
{"ll" , {"", kOrange , 1}},
{"calib_shift" , {"", kGreen , 1}},
{"cherenkov", {"", kBlue+1 , 1}},
{"ll_up" , {"", kOrange , 9}},
{"ll_dw" , {"", kOrange , 2}},
{"calib_shift_up" , {"", kGreen , 9}},
{"calib_shift_dw" , {"", kGreen , 2}},
{"calib_shape" , {"", kRed-4 , 1}},
{"stat" , {"", kBlack , 2}},
{"2kA" , {"", kCyan, 1}},
{"02mmBeamSpotSize", {"", kCyan + 1, 1}},
{"1mmBeamShiftX", {"", kMagenta + 2, 1}},
{"1mmBeamShiftY", {"", kMagenta+2, 2}},
{"3mmHorn1X", {"", kViolet, 1}},
{"3mmHorn1Y", {"", kViolet, 2}},
{"3mmHorn2X", {"", kAzure, 1}},
{"3mmHorn2Y", {"", kAzure, 2}},
{"7mmTargetZ", {"", kBlue-6, 1}},
{"MagneticFieldinDecayPipe", {"", kGreen+2, 1}},
{"1mmHornWater", {"", kYellow-1, 1}},
{"beam_focusing", {"", kViolet+1}}
}
kOrange
enum BeamMode kOrange
kRed
enum BeamMode kRed
kViolet
enum BeamMode kViolet
kGreen
enum BeamMode kGreen
kBlue
enum BeamMode kBlue

Definition at line 50 of file wrong_sign_uncertainty.C.

Referenced by Plot1D().

std::map<std::string, double> xmax
Initial value:
= {{"electron_kin", -1.},
{"enu" , 6.},
{"q2" , 2.}}

Definition at line 79 of file wrong_sign_uncertainty.C.

Referenced by genie::Spline::Add(), calib::HitEfficiency::analyze(), showere::ShowerEnergyAna::analyze(), evd::TZProjPad::AutoZoom(), evd::TZProjPad::AutoZoomBox(), evd::TZProjPad::AutoZoomHelper(), evd::TZProjPad::AutoZoomSlice(), evd::TZProjPad::AutoZoomTruth(), ana::BayesianSurface::BayesianSurface(), novaddt::DAQHeaderAna::beginJob(), comi::NumiFilteringAna::beginJob(), genie::BLI2DNonUnifGrid::BLI2DNonUnifGrid(), genie::BLI2DUnifGrid::BLI2DUnifGrid(), ana::BayesianSurface::BuildHist(), genie::Spline::BuildSpline(), CalcFiberLoopCorr(), genie::SKKinematicsGenerator::CalculateKin_AtharSingleKaon(), genie::COHKinematicsGenerator::CalculateKin_ReinSehgal(), genie::DFRKinematicsGenerator::ComputeMaxXSec(), genie::DISKinematicsGenerator::ComputeMaxXSec(), genie::SKKinematicsGenerator::ComputeMaxXSec(), genie::DMDISKinematicsGenerator::ComputeMaxXSec(), trk::KalmanTrackAna::Containment(), DCSPlotMaker(), demo_CPT(), genie::utils::kinematics::DISImportanceSamplingEnvelope(), genie::Spline::Divide(), Draw1DProfile(), DrawFrame(), DrawMedianAndUniverses(), DrawSliceCanvas(), ana::CovMxSurface::DrawSurfacePoint(), DrawSurfacePoint(), ana::ExtendGraphToTop(), showere::ShowerEnergyFilterMC::filter(), calib::AttenFit::FitQuality(), calib::AttenuationFit::FitQuality(), evd::XZYZProjectionsView::FixAspectRatio(), genie::Spline::GetAsTGraph(), novaddt::MoonShadow::GetHisto(), moonshadowana::MoonShadowAna::GetHisto(), moonshadowana::MoonShadowAnaHough::GetHisto(), GoodDataSelDate(), GoodDataSelDuration(), GoodDataSelEmptySpillFrac(), GoodDataSelGoodDB(), GoodDataSelLiveTime(), GoodDataSelMipRate(), GoodDataSelNumSlices(), GoodDataSelRuns(), GoodDataSelTimingPeakEnd(), GoodDataSelTimingPeakStart(), GoodDataSelTrkFrac(), genie::BLI2DUnifGrid::Init(), genie::BLI2DNonUnifGrid::Init(), remid::ReMIdValidate::IsContained(), remid::RecoMuon::IsTrackContained(), jointsensitivity(), load_libs_muonid(), ana::LogLikelihood(), main(), trk::CosmicTrackAna::Make1DRecoHistograms(), trk::CosmicTrackAna::Make2DRecoHistograms(), MakeBkgEMPlot(), ana::MakeHistCanvasReady_Quant(), MakeLegend(), genie::Spline::Multiply(), PID2DAxis(), PlotEfficiency(), plotTrainTest(), genie::Spline::Print(), PrintPlot(), numusand::FillSandbox::produce(), request_edarea_size(), resolutionplotter(), resolutionplotterbias(), sensitivity2020(), evd::TZProjPad::ShowFull(), evd::TZProjPad::ShowPartial(), SliceLegend(), ana::SurfaceKrige::SurfaceKrige(), test_mc_rejection_method_lin(), test_mc_rejection_method_log(), TrimGraph(), wrong_sign_uncertainty(), and genie::BLI2DGrid::ZMax().

Generated on Thu Dec 3 2020 05:15:02 for NOvA by   doxygen 1.8.11