Functions | Variables
fiducial_accounting.C File Reference
#include <bitset>
#include "CAFAna/Cuts/SpillCuts.h"
#include "CAFAna/Core/SpectrumLoader.h"
#include "CAFAna/Core/SpectrumLoaderBase.h"
#include "CAFAna/Core/HistAxis.h"
#include "CAFAna/Core/Loaders.h"
#include "CAFAna/Vars/Vars.h"
#include "CAFAna/Core/Binning.h"
#include "CAFAna/Core/Cut.h"
#include "CAFAna/Cuts/Cuts.h"
#include "3FlavorAna/Cuts/NumuCuts.h"
#include "CAFAna/Cuts/NueCutsSecondAna.h"
#include "3FlavorAna/Cuts/NueCuts2017.h"
#include "CAFAna/Decomp/ProportionalDecomp.h"
#include "3FlavorAna/Decomp/NueDecomp.h"
#include "CAFAna/Decomp/NCDecomp.h"
#include "3FlavorAna/Decomp/NumuDecomp.h"
#include "CAFAna/Extrap/ModularExtrap.h"
#include "CAFAna/Prediction/PredictionExtrap.h"
#include "CAFAna/Prediction/PredictionNoExtrap.h"
#include "TFile.h"
#include "OscLib/OscCalcDumb.h"
#include "CAFAna/Core/Spectrum.h"
#include "TF1.h"
#include "TH1.h"
#include "TH2.h"
#include "TEfficiency.h"
#include "CAFAna/Experiment/SingleSampleExperiment.h"
#include "TCanvas.h"
#include "CAFAna/Fit/FrequentistSurface.h"
#include "CAFAna/Vars/FitVars.h"
#include "CAFAna/Core/Utilities.h"
#include "Utilities/func/MathUtil.h"
#include "CAFAna/Cuts/AnalysisMasks.h"

Go to the source code of this file.

Functions

void GetExtrap (PredictionNoExtrap *all, PredictionNoExtrap *sel, osc::IOscCalc *calc, TString title)
 
void fiducial_accounting ()
 

Variables

const Cut kBoxCutIdeal ([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return false;if((planesA.second-planesA.first+1)/64< 14) return false;if(sr->mc.nu[0].vtx.x< -100) return false;if(sr->mc.nu[0].vtx.x > 100) return false;if(sr->mc.nu[0].vtx.y< -100) return false;if(sr->mc.nu[0].vtx.y > 100) return false;if(sr->mc.nu[0].vtx.z< 1000) return false;if(sr->mc.nu[0].vtx.z > 3500) return false;return true;})
 
Cut kBoxCut ([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;if(sr->mc.nu[0].vtx.x< -762) return false;if(sr->mc.nu[0].vtx.x > 763) return false;if(sr->mc.nu[0].vtx.y< -762) return false;if(sr->mc.nu[0].vtx.y > 763) return false;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return false;if((planesA.second-planesA.first+1)/64< 4) return false;int first=planesA.first/64;int last=(planesA.second+1)/64;if(sr->mc.nu[0].vtx.z< ((first)*426.0)) return false;if(sr->mc.nu[0].vtx.z >((last)*426.0)) return false;return true;})
 
const Var kBlock ([](const caf::SRProxy *sr){int size=-1;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return size;size=(planesA.second-planesA.first+1)/64;return size;})
 

Function Documentation

void fiducial_accounting ( )

Definition at line 127 of file fiducial_accounting.C.

References allInOneTrainingPlots::axis, calc, file, GetExtrap(), ana::Loaders::Go(), ana::kBeam, kBlock, ana::kBoxCut, ana::kBoxCutIdeal, caf::kFARDET, ana::Loaders::kFluxSwap, ana::Loaders::kMC, ana::kNue2017CorePresel, kSelectionCut, ana::kStandardSpillCuts, ana::Loaders::SetLoaderPath(), ana::Loaders::SetSpillCut(), and ana::Binning::Simple().

128 {
129 
130  TH1::SetDefaultSumw2();
131  TH2::SetDefaultSumw2();
132 
134 
135  // SpectrumLoader loaderFDSwap(
136  // "prod_caf_R17-03-01-prod3reco.l_fd_genie_fluxswap_fhc_nova_v08_full_v1" );
137  //
138  // NullLoader loaderFDNonSwap;
139  Loaders loader_fdmc;
140  loader_fdmc.SetLoaderPath("prod_caf_R17-03-01-prod3reco.l_fd_genie_fluxswap_fhc_nova_v08_full_v1",
141  caf::kFARDET, Loaders::kMC, kBeam, Loaders::kFluxSwap);
142 
143  const HistAxis axis(
144  "diblocks", Binning::Simple(15,-0.5,14.5), kBlock );
145 
146  // make prediction for all true neutrinos in fiducial volume
147  PredictionNoExtrap* predAll = new PredictionNoExtrap(
148  loader_fdmc,
149  axis, kBoxCut
150  );
151  //prediction for all true neutrinos in fiducial volume that pass preselection
152  PredictionNoExtrap* predSel = new PredictionNoExtrap(
153  loader_fdmc,
154  axis, kSelectionCut&&kBoxCut
155  );
156  // make prediction for all true neutrinos in small box at center of detector
157  PredictionNoExtrap* predIdealAll = new PredictionNoExtrap(
158  loader_fdmc,
160  );
161  //prediction for all true neutrinos in fiducial volume that pass preselection in small box at center of detector
162  PredictionNoExtrap* predIdealSel = new PredictionNoExtrap(
163  loader_fdmc,
164  axis, kSelectionCut&&kBoxCutIdeal
165  );
166 
167  loader_fdmc.SetSpillCut(kStandardSpillCuts);
168 
169 
170 
171  // Load from files
172  // loaderFDNonSwap.Go();
173  // loaderFDSwap.Go();
174  loader_fdmc.Go();
175 
176  // Open output file
177  TFile file("fd_fiducial_mass_efficiency.root","RECREATE");
178  // Save things in output file
179 
180 
181 
183 
184  //caclulate preselection efficieny as a function of diblock
185  GetExtrap(predAll,predSel, &calc,"FidMass");
186  //caclulate preselection efficieny for small window in ideal detector
187  GetExtrap(predIdealAll,predIdealSel, &calc,"SmallWindow");
188 
189  file.Close();
190 
191 }
Far Detector at Ash River.
Definition: SREnums.h:11
Collection of SpectrumLoaders for many configurations.
Definition: Loaders.h:23
const Color_t kMC
osc::OscCalcDumb calc
const Cut kSelectionCut
void Go()
Call Go() on all the loaders.
Definition: Loaders.cxx:162
void SetSpillCut(const SpillCut &cut)
Definition: Loaders.cxx:121
const Cut kBoxCutIdeal([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return false;if((planesA.second-planesA.first+1)/64< 14) return false;if(sr->mc.nu[0].vtx.x< -100) return false;if(sr->mc.nu[0].vtx.x > 100) return false;if(sr->mc.nu[0].vtx.y< -100) return false;if(sr->mc.nu[0].vtx.y > 100) return false;if(sr->mc.nu[0].vtx.z< 1000) return false;if(sr->mc.nu[0].vtx.z > 3500) return false;return true;})
Simple oscillation probability calculator that has no solar term or mass hierarchy or delta so it&#39;s s...
Definition: OscCalcDumb.h:16
Cut kBoxCut([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;if(sr->mc.nu[0].vtx.x< -762) return false;if(sr->mc.nu[0].vtx.x > 763) return false;if(sr->mc.nu[0].vtx.y< -762) return false;if(sr->mc.nu[0].vtx.y > 763) return false;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return false;if((planesA.second-planesA.first+1)/64< 4) return false;int first=planesA.first/64;int last=(planesA.second+1)/64;if(sr->mc.nu[0].vtx.z< ((first)*426.0)) return false;if(sr->mc.nu[0].vtx.z >((last)*426.0)) return false;return true;})
void GetExtrap(PredictionNoExtrap *all, PredictionNoExtrap *sel, osc::IOscCalc *calc, TString title)
const Cut kNue2017CorePresel
Full preselection for 2017 "core" sample with Cut-based CosRej.
Definition: NueCuts2017.h:105
TFile * file
Definition: cellShifts.C:17
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
Prediction that just uses FD MC, with no extrapolation.
const Var kBlock([](const caf::SRProxy *sr){int size=-1;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return size;size=(planesA.second-planesA.first+1)/64;return size;})
void SetLoaderPath(const std::string &path, caf::Det_t det, DataMC datamc, DataSource src=kBeam, SwappingConfig swap=kNonSwap)
Configure loader via wildcard path.
Definition: Loaders.cxx:25
void GetExtrap ( PredictionNoExtrap all,
PredictionNoExtrap sel,
osc::IOscCalc calc,
TString  title 
)

Definition at line 96 of file fiducial_accounting.C.

References MECModelEnuComparisons::i, calib::j, ana::Sign::kBoth, ana::Current::kCC, ana::Flavors::kNuMuToNuE, kRed, ana::Spectrum::POT(), ana::PredictionExtrap::Predict(), ana::PredictionExtrap::PredictComponent(), and ana::Spectrum::ToTH1().

Referenced by fiducial_accounting(), make_extrap_figure_hists(), and xsec_extrap_plots().

97 {
98  TH1::SetDefaultSumw2();
99  TH2::SetDefaultSumw2();
100  double kPOTFD = all->Predict(calc).POT();
101  TH1* hSig = all->PredictComponent(calc,
103  Current::kCC,
104  Sign::kBoth).ToTH1(kPOTFD );
105  hSig->SetLineColor(kRed-10);
106 
107  TH1* hSigSel = sel->PredictComponent(calc,
109  Current::kCC,
110  Sign::kBoth).ToTH1(kPOTFD );
111 
112  TH1D* hN = new TH1D("hN_"+title,"Nue CC selection efficiency;diblocks;efficiency",11,3.5,14.5);
113  TH1D* hD = new TH1D("hD_"+title,"Nue CC selection efficiency;diblocks;efficiency",11,3.5,14.5);
114  for(unsigned int i=5;i<16;++i){
115  for(unsigned int j=0;j<hSigSel->GetBinContent(i);++j){
116  hN->Fill(i-1);
117  }
118  for(unsigned int j=0;j<hSig->GetBinContent(i);++j){
119  hD->Fill(i-1);
120  }
121  }
122  hN->Write();
123  hD->Write();
124 
125 }
Spectrum PredictComponent(osc::IOscCalc *calc, Flavors::Flavors_t flav, Current::Current_t curr, Sign::Sign_t sign) const override
enum BeamMode kRed
(&#39; appearance&#39;)
Definition: IPrediction.h:18
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
Charged-current interactions.
Definition: IPrediction.h:39
Spectrum Predict(osc::IOscCalc *calc) const override
const double j
Definition: BetheBloch.cxx:29
double POT() const
Definition: Spectrum.h:227

Variable Documentation

const Var kBlock([](const caf::SRProxy *sr){int size=-1;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return size;size=(planesA.second-planesA.first+1)/64;return size;})

Referenced by fiducial_accounting().

Cut kBoxCut([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;if(sr->mc.nu[0].vtx.x< -762) return false;if(sr->mc.nu[0].vtx.x > 763) return false;if(sr->mc.nu[0].vtx.y< -762) return false;if(sr->mc.nu[0].vtx.y > 763) return false;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return false;if((planesA.second-planesA.first+1)/64< 4) return false;int first=planesA.first/64;int last=(planesA.second+1)/64;if(sr->mc.nu[0].vtx.z< ((first)*426.0)) return false;if(sr->mc.nu[0].vtx.z >((last)*426.0)) return false;return true;})
const Cut kBoxCutIdeal([](const caf::SRProxy *sr){if(sr->mc.nnu==0) return false;std::bitset< 14 > binary(sr->hdr.dibmask);std::pair< int, int > planesA=calcFirstLastLivePlane(sr->slc.firstplane, binary);std::pair< int, int > planesB=calcFirstLastLivePlane(sr->slc.lastplane, binary);if((planesA.first!=planesB.first)||(planesA.second!=planesB.second)) return false;if((planesA.second-planesA.first+1)/64< 14) return false;if(sr->mc.nu[0].vtx.x< -100) return false;if(sr->mc.nu[0].vtx.x > 100) return false;if(sr->mc.nu[0].vtx.y< -100) return false;if(sr->mc.nu[0].vtx.y > 100) return false;if(sr->mc.nu[0].vtx.z< 1000) return false;if(sr->mc.nu[0].vtx.z > 3500) return false;return true;})