make_RHC_WrongSign_Nue.C
Go to the documentation of this file.
1 // ----------------------------------------------------------------
2 // Make the Predictions at ND for nue, including data
3 // See docDB-29919 for tech note on WS study
4 // ----------------------------------------------------------------
5 
6 #include "CAFAna/Core/Binning.h"
7 #include "CAFAna/Core/Spectrum.h"
8 #include "CAFAna/Core/SpectrumLoader.h"
9 #include "CAFAna/Cuts/Cuts.h"
10 #include "CAFAna/Cuts/TruthCuts.h"
11 #include "CAFAna/Cuts/SpillCuts.h"
12 #include "3FlavorAna/Cuts/NueCuts2018.h"
13 #include "3FlavorAna/Cuts/WrongSignCuts2018.h"
14 #include "3FlavorAna/Cuts/NumuCuts2018.h"
15 #include "3FlavorAna/Cuts/QuantileCuts.h"
16 #include "CAFAna/Vars/Vars.h"
17 #include "3FlavorAna/Vars/NumuVars.h"
18 #include "3FlavorAna/Vars/NueVars.h"
19 #include "CAFAna/Vars/GenieWeights.h"
20 #include "CAFAna/Vars/PPFXWeights.h"
21 #include "3FlavorAna/Vars/WrongSignBDTVars.h"
22 
23 #include "CAFAna/Prediction/PredictionNoExtrap.h"
24 #include "CAFAna/Prediction/PredictionNoOsc.h"
25 #include "OscLib/OscCalcPMNSOpt.h"
26 
27 #include "TCanvas.h"
28 #include "TH2.h"
29 #include "TLegend.h"
30 #include "TLatex.h"
31 #include "TFile.h"
32 #include "TProfile.h"
33 #include "TString.h"
34 
35 #include "Utilities/rootlogon.C"
36 
37 using namespace ana;
38 
39 void make_RHC_WrongSign_Nue()
40 {
41  // ND MC data set, Prod 4 RHC
42  // nue
43  std::string fname_nueRHC = "prod_sumdecaf_R17-11-14-prod4reco.e_nd_genie_nonswap_rhc_nova_v08_full_v1_nue2018";
44  SpectrumLoader ld_nueRHC(fname_nueRHC);
45  ld_nueRHC.SetSpillCut(kStandardSpillCuts);
46  // numu
47  std::string fname_numuRHC = "prod_sumdecaf_R17-11-14-prod4reco.e_nd_genie_nonswap_rhc_nova_v08_full_v1_numu2018";
48  SpectrumLoader ld_numuRHC(fname_numuRHC);
49  ld_numuRHC.SetSpillCut(kStandardSpillCuts);
50 
51  // ND real data set, all periods concats (I think this is periods 4, 6)
52  std::string rhcDataNue = "prod_sumdecaf_R17-09-05-prod4recopreview.f_nd_numi_rhc_full_v1_addShortSimpleCVN_goodruns_nue2018";
53  SpectrumLoader ld_dataRHCNue(rhcDataNue);
54  ld_dataRHCNue.SetSpillCut(kStandardSpillCuts);
55 
56  std::string rhcDataNumu = "prod_sumdecaf_R17-09-05-prod4recopreview.f_nd_numi_rhc_full_v1_addShortSimpleCVN_goodruns_numu2018";
57  SpectrumLoader ld_dataRHCNumu(rhcDataNumu);
58  ld_dataRHCNumu.SetSpillCut(kStandardSpillCuts);
59 
60  // Standard weights
61  Var wei = kXSecCVWgt2018 * kPPFXFluxCVWgt;
62 
63  // Extra Cuts and Vars for analysis
64  // -------------------------------------------------------------------------------------------------
65  const Cut kCutHighCVNSSe = kCVNSSe >= 0.98; // Ana 2018 Hi CVN bin (docDB-27317)
66  const Cut kCutLowCVNSSe = kCVNSSe >= 0.89; // Ana 2018 Lo CVN bin has this as low edge
67 
68  const Cut kSelectRSNumuND0 = kNumuCutND2018 && kCVNFSProtonScore2018<0.49;
69  const Cut kSelectRSNumuND1 = kNumuCutND2018 && !kCVNProngProtonCtNumu;
70  const Cut kSelectRSNumuND2 = kNumuCutND2018 && kAntiNumuBDTCVN>0.11;
71 
72  const Cut kSelectWSNumuND0 = kNumuCutND2018 && kCVNFSProtonScore2018>=0.49;
73  const Cut kSelectWSNumuND1 = kNumuCutND2018 && kCVNProngProtonCtNumu;
74  const Cut kSelectWSNumuND2 = kNumuCutND2018 && kAntiNumuBDTCVN<=0.11;
75 
76  // Nue-style Binnings for BDT
77  const Var kNueBDT3Bin([](const caf::SRProxy *sr)
78  {
79  if( kAntiNueBDTCVN(sr)>0. ) return float(0.5);
80  else return float(1.5);
81 
82  return float(-1.5);
83  });
84 
85  const Var kNueBDT3AnaBinning([&](const caf::SRProxy *sr)
86  {
87  int selBin = kNueBDT3Bin(sr);
88 
89  float nuE = kNueEnergy2018(sr);
90  int nuEBin = nuE/0.5;
91  assert(nuEBin <= 8 && "An event with nuE > 4.5 should never happen");
92 
93  int anaBin = 9*selBin + nuEBin;
94  return anaBin;
95  });
96 
97  const HistAxis kAxisFSProton("NuE Energy / Analysis Bin",Binning::Simple(18,0,18),kNueFSProtonAnaBinning);
98  const HistAxis kAxisCVNProng("NuE Energy / Analysis Bin",Binning::Simple(18,0,18),kNueCVNProngAnaBinning);
99  const HistAxis kAxisBDT3("NuE Energy / Analysis Bin",Binning::Simple(18,0,18),kNueBDT3AnaBinning);
100  // -------------------------------------------------------------------------------------------------
101 
102  // Predictions for nue RHC
103  PredictionNoOsc pNueFSProtonHiCVN(ld_nueRHC,kAxisFSProton,kNue2018NDCVNSsb && kCutHighCVNSSe,kNoShift,wei);
104  PredictionNoOsc pNueCVNProngHiCVN(ld_nueRHC,kAxisCVNProng,kNue2018NDCVNSsb && kCutHighCVNSSe,kNoShift,wei);
105  PredictionNoOsc pNueBDT3HiCVN(ld_nueRHC,kAxisBDT3,kNue2018NDCVNSsb && kCutHighCVNSSe,kNoShift,wei);
106  // Splitting into NC from nu and antinu using spectra for this analysis
107  Spectrum sNueFSProtonHiCVNNC(ld_nueRHC,
108  kAxisFSProton,
109  kNue2018NDCVNSsb && kCutHighCVNSSe && kIsNC && !kIsAntiNu,
110  kNoShift,
111  wei);
112  Spectrum sNueFSProtonHiCVNNCAnti(ld_nueRHC,
113  kAxisFSProton,
114  kNue2018NDCVNSsb && kCutHighCVNSSe && kIsNC && kIsAntiNu,
115  kNoShift,
116  wei);
117  Spectrum sNueCVNProngHiCVNNC(ld_nueRHC,
118  kAxisCVNProng,
119  kNue2018NDCVNSsb && kCutHighCVNSSe && kIsNC && !kIsAntiNu,
120  kNoShift,
121  wei);
122  Spectrum sNueCVNProngHiCVNNCAnti(ld_nueRHC,
123  kAxisCVNProng,
124  kNue2018NDCVNSsb && kCutHighCVNSSe && kIsNC && kIsAntiNu,
125  kNoShift,
126  wei);
127  Spectrum sNueBDT3HiCVNNC(ld_nueRHC,
128  kAxisBDT3,
129  kNue2018NDCVNSsb && kCutHighCVNSSe && kIsNC && !kIsAntiNu,
130  kNoShift,
131  wei);
132  Spectrum sNueBDT3HiCVNNCAnti(ld_nueRHC,
133  kAxisBDT3,
134  kNue2018NDCVNSsb && kCutHighCVNSSe && kIsNC && kIsAntiNu,
135  kNoShift,
136  wei);
137 
138  PredictionNoOsc pNueFSProtonLoCVN(ld_nueRHC,kAxisFSProton,kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe,kNoShift,wei);
139  PredictionNoOsc pNueCVNProngLoCVN(ld_nueRHC,kAxisCVNProng,kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe,kNoShift,wei);
140  PredictionNoOsc pNueBDT3LoCVN(ld_nueRHC,kAxisBDT3,kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe,kNoShift,wei);
141  // Splitting into NC from nu and antinu using spectra for this analysis
142  Spectrum sNueFSProtonLoCVNNC(ld_nueRHC,
143  kAxisFSProton,
144  kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe && kIsNC && !kIsAntiNu,
145  kNoShift,
146  wei);
147  Spectrum sNueFSProtonLoCVNNCAnti(ld_nueRHC,
148  kAxisFSProton,
149  kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe && kIsNC && kIsAntiNu,
150  kNoShift,
151  wei);
152  Spectrum sNueCVNProngLoCVNNC(ld_nueRHC,
153  kAxisCVNProng,
154  kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe && kIsNC && !kIsAntiNu,
155  kNoShift,
156  wei);
157  Spectrum sNueCVNProngLoCVNNCAnti(ld_nueRHC,
158  kAxisCVNProng,
159  kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe && kIsNC && kIsAntiNu,
160  kNoShift,
161  wei);
162  Spectrum sNueBDT3LoCVNNC(ld_nueRHC,
163  kAxisBDT3,
164  kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe && kIsNC && !kIsAntiNu,
165  kNoShift,
166  wei);
167  Spectrum sNueBDT3LoCVNNCAnti(ld_nueRHC,
168  kAxisBDT3,
169  kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe && kIsNC && kIsAntiNu,
170  kNoShift,
171  wei);
172 
173  // ------------------------------------------------------------------------------------------------
174  // Nue DATA
175  Spectrum sNueData_FSProtonHiCVN(ld_dataRHCNue,kAxisFSProton,kNue2018NDCVNSsb && kCutHighCVNSSe,kNoShift,kUnweighted);
176  Spectrum sNueData_CVNProngHiCVN(ld_dataRHCNue,kAxisCVNProng,kNue2018NDCVNSsb && kCutHighCVNSSe,kNoShift,kUnweighted);
177  Spectrum sNueData_BDT3HiCVN(ld_dataRHCNue,kAxisBDT3,kNue2018NDCVNSsb && kCutHighCVNSSe,kNoShift,kUnweighted);
178  Spectrum sNueData_FSProtonLoCVN(ld_dataRHCNue,kAxisFSProton,kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe,kNoShift,kUnweighted);
179  Spectrum sNueData_CVNProngLoCVN(ld_dataRHCNue,kAxisCVNProng,kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe,kNoShift,kUnweighted);
180  Spectrum sNueData_BDT3LoCVN(ld_dataRHCNue,kAxisBDT3,kNue2018NDCVNSsb && kCutLowCVNSSe && !kCutHighCVNSSe,kNoShift,kUnweighted);
181  // ------------------------------------------------------------------------------------------------
182 
183  // Press play
184  ld_nueRHC.Go();
185  ld_dataRHCNue.Go();
186 
187  // Save the predictions
188  TFile *outFile = new TFile("./PredictionsAndDataWrongSignNue.root","RECREATE");
189  outFile->cd();
190  pNueFSProtonHiCVN.SaveTo(outFile, "pNueFSProtonHiCVN");
191  pNueCVNProngHiCVN.SaveTo(outFile, "pNueCVNProngHiCVN");
192  pNueBDT3HiCVN.SaveTo(outFile, "pNueBDT3HiCVN");
193  pNueFSProtonLoCVN.SaveTo(outFile, "pNueFSProtonLoCVN");
194  pNueCVNProngLoCVN.SaveTo(outFile, "pNueCVNProngLoCVN");
195  pNueBDT3LoCVN.SaveTo(outFile, "pNueBDT3LoCVN");
196 
197  sNueFSProtonHiCVNNC.SaveTo(outFile, "sNueFSProtonHiCVNNC");
198  sNueFSProtonHiCVNNCAnti.SaveTo(outFile, "sNueFSProtonHiCVNNCAnti");
199  sNueCVNProngHiCVNNC.SaveTo(outFile, "sNueCVNProngHiCVNNC");
200  sNueCVNProngHiCVNNCAnti.SaveTo(outFile, "sNueCVNProngHiCVNNCAnti");
201  sNueBDT3HiCVNNC.SaveTo(outFile, "sNueBDT3HiCVNNC");
202  sNueBDT3HiCVNNCAnti.SaveTo(outFile, "sNueBDT3HiCVNNCAnti");
203  sNueFSProtonLoCVNNC.SaveTo(outFile, "sNueFSProtonLoCVNNC");
204  sNueFSProtonLoCVNNCAnti.SaveTo(outFile, "sNueFSProtonLoCVNNCAnti");
205  sNueCVNProngLoCVNNC.SaveTo(outFile, "sNueCVNProngLoCVNNC");
206  sNueCVNProngLoCVNNCAnti.SaveTo(outFile, "sNueCVNProngLoCVNNCAnti");
207  sNueBDT3LoCVNNC.SaveTo(outFile, "sNueBDT3LoCVNNC");
208  sNueBDT3LoCVNNCAnti.SaveTo(outFile, "sNueBDT3LoCVNNCAnti");
209 
210  // DATA
211  sNueData_FSProtonHiCVN.SaveTo(outFile, "sNueData_FSProtonHiCVN");
212  sNueData_CVNProngHiCVN.SaveTo(outFile, "sNueData_CVNProngHiCVN");
213  sNueData_BDT3HiCVN.SaveTo(outFile, "sNueData_BDT3HiCVN");
214  sNueData_FSProtonLoCVN.SaveTo(outFile, "sNueData_FSProtonLoCVN");
215  sNueData_CVNProngLoCVN.SaveTo(outFile, "sNueData_CVNProngLoCVN");
216  sNueData_BDT3LoCVN.SaveTo(outFile, "sNueData_BDT3LoCVN");
217 
218  outFile->Close();
219 }
ana
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
ana::kPPFXFluxCVWgt
const Var kPPFXFluxCVWgt
Definition: PPFXWeights.h:16
ana::kNumuCutND2018
const Cut kNumuCutND2018
Definition: NumuCuts2018.h:41
caf::Proxy< caf::StandardRecord >
Proxy for caf::StandardRecord.
Definition: SRProxy.h:2126
ana::kAntiNueBDTCVN
const Var kAntiNueBDTCVN
Definition: WrongSignBDTVars.h:43
ana::SpectrumLoaderBase::SetSpillCut
void SetSpillCut(const SpillCut &cut)
Definition: SpectrumLoaderBase.cxx:64
check_time_usage.float
float
Definition: check_time_usage.py:6
ana::kNueEnergy2018
const Var kNueEnergy2018([](const caf::SRProxy *sr){if(!sr->vtx.elastic.IsValid) return-1.0;if(sr->vtx.elastic.fuzzyk.npng< 1) return-1.0;if(kIsRHC(sr)) return kNueEnergyRHC(sr);else return kNueEnergyFHC(sr);})
Definition: NueEnergy2018.h:25
ana::Energy
const Cut Energy
Definition: NCPi0SingleProngCuts.h:313
sr
caf::StandardRecord * sr
Definition: FillTrainingTrees.h:18
ana::kCVNFSProtonScore2018
const Var kCVNFSProtonScore2018([](const caf::SRProxy *sr){return CVNFinalStateScore2018(sr, 2212);})
Proton score from CVN Final State labels, using Prod3Train CVN.
Definition: CVNFinalStates.h:22
ana::kCVNSSe
const Var kCVNSSe([](const caf::SRProxy *sr){throw std::runtime_error("kCVNSSe is no longer available. Fix your macro so you don't use it.");return-5.;})
2018 nue PID
Definition: Vars.h:52
ana::kCVNProngProtonCtNumu
const Cut kCVNProngProtonCtNumu([](const caf::SRProxy *sr){ int count=0;if(sr->vtx.elastic.fuzzyk.npng==1) return false;for(unsigned int i=0;i< sr->vtx.elastic.fuzzyk.npng;++i){if(sr->vtx.elastic.fuzzyk.png[i].cvnpart.pdgmax!=2212) continue;if(util::pythag(sr->vtx.elastic.vtx.X()-sr->vtx.elastic.fuzzyk.png[i].start.X(), sr->vtx.elastic.vtx.Y()-sr->vtx.elastic.fuzzyk.png[i].start.Y(), sr->vtx.elastic.vtx.Z()-sr->vtx.elastic.fuzzyk.png[i].start.Z()) > 20) continue;if(sr->vtx.elastic.fuzzyk.png[i].cvnpart.maxval< 0.72) continue;count++;}return count >0;})
Definition: WrongSignCuts2018.h:17
make_RHC_WrongSign_Nue
void make_RHC_WrongSign_Nue()
Definition: make_RHC_WrongSign_Nue.C:39
ana::SpectrumLoader
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
Definition: SpectrumLoader.h:21
ana::kAntiNumuBDTCVN
const Var kAntiNumuBDTCVN
Definition: WrongSignBDTVars.h:55
ana::kXSecCVWgt2018
const Var kXSecCVWgt2018
Definition: XsecTunes.h:48
ana::assert
assert(nhit_max >=nhit_nbins)
ana::kStandardSpillCuts
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
ana::_Cut< caf::SRProxy >
ana::_Var< caf::SRProxy >
string
enum BeamMode string
Generated on Mon Jan 18 2021 14:12:49 for NOvA by   doxygen 1.8.11