make_nueFDprediction_kinematics_RHC_REW_pTExtrap.C
Go to the documentation of this file.
1 // Analysis
6 
7 // Core
8 #include "CAFAna/Core/Binning.h"
9 #include "CAFAna/Core/Cut.h"
10 #include "CAFAna/Core/HistAxis.h"
12 #include "CAFAna/Core/Spectrum.h"
15 
16 // Cuts
17 #include "CAFAna/Cuts/Cuts.h"
21 #include "CAFAna/Cuts/SpillCuts.h"
24 
25 // Predict
28 
29 // Vars
30 #include "CAFAna/Vars/XsecTunes.h"
31 #include "CAFAna/Vars/TruthVars.h"
34 #include "CAFAna/Vars/Vars.h"
36 
38 
39 #include <iostream>
40 #include <iomanip>
41 
42 using namespace ana;
43 
44 //sample : ALL/CORE/PERI and the corresponding Weight file name.
46 {
47  const std::string outfilename = ("FDprediction_kinematics_RHC_REW_pTExtrap"+sample+".root").c_str();
48 
49  //For grid
50  // std::string wfile_name = ("/pnfs/nova/scratch/users/zvallari/AcceptSyst/"+wfile).c_str();
51  // TFile *weight_file = TFile::Open(pnfs2xrootd(wfile_name).c_str() );
52 
53  //For non-grid
54  std::string wfile_name = ("/nova/ana/users/zvallari/AcceptSyst2020/"+wfile).c_str();
55  TFile *weight_file = TFile::Open(wfile_name.c_str());
56 
57  TH1 *h = (TH1*)weight_file->Get(("trueQ2_weight_"+sample).c_str());
58  int nbins = h->GetNbinsX();
59  double bins[nbins], binslow[nbins], binshigh[nbins];
60  for(int i = 1; i <= nbins; ++i){
61  bins[i-1] = h->GetBinContent(i);
62  binslow[i-1] = h->GetBinLowEdge(i);
63  binshigh[i-1] = h->GetBinLowEdge(i) + h->GetBinWidth(i);
64  }
65 
66 const Var kTrueQ2Weight([nbins, &binslow, &binshigh, &bins](const caf::SRProxy *sr){
67  if(sr->hdr.det == caf::kFARDET) return 1.0;
68  for(int i = 0; i < nbins; ++i){
69  if((kRecoQ2(sr) > binslow[i]) && (kRecoQ2(sr) < binshigh[i])){
70  return bins[i];
71  }
72  }
73  return 1.0;
74  }
75  );
76 
77  TH1 *hptp = (TH1*)weight_file->Get(("PtP_weight_"+sample).c_str());
78  int nbins_ptp = hptp->GetNbinsX();
79  double bins_ptp[nbins_ptp], binslow_ptp[nbins_ptp], binshigh_ptp[nbins_ptp];
80  for(int i = 1; i <= nbins_ptp; ++i){
81  bins_ptp[i-1] = hptp->GetBinContent(i);
82  binslow_ptp[i-1] = hptp->GetBinLowEdge(i);
83  binshigh_ptp[i-1] = hptp->GetBinLowEdge(i) + hptp->GetBinWidth(i);
84  }
85 
86  const Var kPtPWeight([nbins_ptp, &binslow_ptp, &binshigh_ptp, &bins_ptp](const caf::SRProxy *sr){
87  if(sr->hdr.det == caf::kFARDET) return 1.0;
88  for(int i = 0; i < nbins_ptp; ++i){
89  if((kPtP(sr) > binslow_ptp[i]) && (kPtP(sr) < binshigh_ptp[i])){return bins_ptp[i];}
90  }
91  return 1.0;
92  }
93  );
94 
95  TH1 *hcos = (TH1*)weight_file->Get(("CosNumi_weight_"+sample).c_str());
96  int nbins_cos = hptp->GetNbinsX();
97  double bins_cos[nbins_cos], binslow_cos[nbins_cos], binshigh_cos[nbins_cos];
98  for(int i = 1; i <= nbins_cos; ++i){
99  bins_cos[i-1] = hcos->GetBinContent(i);
100  binslow_cos[i-1] = hcos->GetBinLowEdge(i);
101  binshigh_cos[i-1] = hcos->GetBinLowEdge(i) + hcos->GetBinWidth(i);
102  }
103 
104  const Var kCosWeight([nbins_cos, &binslow_cos, &binshigh_cos, &bins_cos](const caf::SRProxy *sr){
105  if(sr->hdr.det == caf::kFARDET) return 1.0;
106  for(int i = 0; i < nbins_cos; ++i){
107  if((kCosNumi(sr) > binslow_cos[i]) && (kCosNumi(sr) < binshigh_cos[i])){return bins_cos[i];}
108  }
109  return 1.0;
110  }
111  );
112 
113  // close weights file
114  weight_file->Close();
115 
116  const unsigned int nQuants = 3;
117  bool isRHC = true;
118  std::vector<Cut> PtQuantCutsND = GetNueQuantCuts2020( isRHC, caf::kNEARDET, nQuants, ana::kExtrapPt);
119  std::vector<Cut> PtQuantCutsFD = GetNueQuantCuts2020( isRHC, caf::kFARDET, nQuants, ana::kExtrapPt);
120 
123 
124  struct GenDef{
125  const IPredictionGenerator *gen;
126  const TString cutname;
127  const TString varname;
128  };
129 
130  std::vector <GenDef> gens;
131  std::vector <IPrediction * > predictions;
132 
134  sample.c_str(), "nueAxis_NoExtrap"});
135 
136  //Loop over samples
137  for ( unsigned int quant = 0; quant < nQuants; ++quant ){
140  kNue2020FDAllSamples && PtQuantCutsFD[quant],
141  kNumu2020ND && PtQuantCutsND[quant],
142  kNoShift,
143  kXSecCVWgt2020*kPPFXFluxCVWgt*kTrueQ2Weight),
144  sample.c_str(), Form("nueAxis_NueSignalExtrap_Q2Weight_Quant%i", quant+1)});
145 
148  kNue2020FDAllSamples && PtQuantCutsFD[quant],
149  kNumu2020ND && PtQuantCutsND[quant],
150  kNoShift,
151  kXSecCVWgt2020*kPPFXFluxCVWgt*kCosWeight),
152  sample.c_str(), Form("nueAxis_NueSignalExtrap_CosWeight_Quant%i", quant+1)});
153 
156  kNue2020FDAllSamples && PtQuantCutsFD[quant],
157  kNumu2020ND && PtQuantCutsND[quant],
158  kNoShift,
159  kXSecCVWgt2020*kPPFXFluxCVWgt*kPtPWeight),
160  sample.c_str(), Form("nueAxis_NueSignalExtrap_PtPWeight_Quant%i", quant+1)});
161  }
162 
163  for(auto & gen:gens){
164  predictions.push_back(gen.gen->Generate(loaders).release());
165  }
166 
167  loaders.Go();
168 
169  TFile* file = new TFile(outfilename.c_str(),"RECREATE");
170  for (int i = 0 ; i < (int)gens.size(); ++i){
171  auto dir = file->GetDirectory(gens[i].cutname);
172  if(!dir) dir = file->mkdir(gens[i].cutname);
173  predictions[i]->SaveTo(dir, gens[i].varname);
174  }
175 
176  file->Close();
177 
178 }
Near Detector underground.
Definition: SREnums.h:10
Far Detector at Ash River.
Definition: SREnums.h:11
std::map< TString, IPredictionGenerator * > gens
Definition: syst_header.h:387
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const Var kPPFXFluxCVWgt
Definition: PPFXWeights.h:16
void make_nueFDprediction_kinematics_RHC_REW_pTExtrap(const std::string sample, const std::string wfile)
caf::Proxy< caf::SRHeader > hdr
Definition: SRProxy.h:2137
Proxy for caf::StandardRecord.
Definition: SRProxy.h:2126
const Cut kNumu2020ND
Definition: NumuCuts2020.h:57
const HistAxis kNumuCCOptimisedAxis2020("Reconstructed Neutrino Energy (GeV)", kNumuCCEOptimisedBinning, kNumuE2020)
Definition: HistAxes.h:25
Generates FD-only predictions (no extrapolation)
void Go()
Call Go() on all the loaders.
Definition: Loaders.cxx:162
void SetSpillCut(const SpillCut &cut)
Definition: Loaders.cxx:121
string outfilename
knobs that need extra care
const int nbins
Definition: cellShifts.C:15
For nominal spectra and reweighting systs (xsec/flux)
Definition: Prod5Loaders.h:101
const HistAxis kNue2020Axis("NuE Energy / Analysis Bin", kNue2020Binning, kNue2020AnaBin)
Use this Axis for Ana2020, official Axis.
Definition: NueCuts2020.h:195
Generates extrapolated Nue signal-only predictions.
caf::StandardRecord * sr
std::vector< Cut > GetNueQuantCuts2020(const bool isRHC, const caf::Det_t det, const unsigned int nquants, const ExtrapVar var)
const SystShifts kNoShift
Definition: SystShifts.cxx:21
const Binning bins
Definition: NumuCC_CPiBin.h:8
TDirectory * dir
Definition: macro.C:5
std::vector< Loaders * > loaders
Definition: syst_header.h:386
TFile * file
Definition: cellShifts.C:17
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
Given loaders and an MC shift, Generate() generates an IPrediction.
const Cut kNue2020FDAllSamples
Definition: NueCuts2020.h:84
const Var kXSecCVWgt2020
Definition: XsecTunes.h:106
caf::Proxy< caf::Det_t > det
Definition: SRProxy.h:232
enum BeamMode string