Evaluate_BDTMLP_Algs_Spectra.C
Go to the documentation of this file.
1 // A macro to test the results that I got from training the Cos Rej BDT.
2 // something like TMVA.
3 
4 #include "CAFAna/Core/Binning.h"
5 #include "CAFAna/Core/Spectrum.h"
6 #include "CAFAna/Core/SpectrumLoader.h"
7 //#include "CAFAna/Cuts/Cuts.h"
8 #include "CAFAna/Cuts/QuantileCuts.h"
9 //#include "CAFAna/Cuts/NumuCuts2018.h"
10 //#include "CAFAna/Cuts/SpillCuts.h"
11 
12 #include "CAFAna/Vars/BPFVars.h"
13 #include "CAFAna/Vars/BPFEnergyVars.h"
14 #include "CAFAna/Vars/CosmicRejBDTVars.h"
15 #include "CAFAna/Vars/GenieWeights.h"
16 #include "CAFAna/Vars/HistAxes.h"
17 #include "CAFAna/Vars/PPFXWeights.h"
18 #include "CAFAna/Vars/Vars.h"
19 #include "CAFAna/Vars/NumuVars.h"
20 
21 #include "CosRej/training/Ana20/Evaluate_header.h"
22 
23 #include "Utilities/rootlogon.C"
24 
25 #include "TTree.h"
26 #include "TString.h"
27 
28 using namespace ana;
29 
30 // -----------------------------------------------------------
31 // -----------------------------------------------------------
32 // --- Now for my actual module.
33 // -----------------------------------------------------------
34 // -----------------------------------------------------------
35 void Evaluate_BDTMLP_Algs_Spectra( std::string sIdentifier, bool isFHC, bool isBPF, bool Quant, bool isEval = true ) {
36  // -----------------------------------------------------------
37  // --- First off, lets set some global things and make sure we know what we're running over...
38  // -----------------------------------------------------------
39  std::string sFHC_RHC = (isFHC == true ? "fhc" :"rhc" );
40  std::string sisEval = (isEval == true ? "Test" :"Train" );
41  std::string sisBPF = (isBPF == true ? "BPF" :"Kal" );
42  std::string sisData = "Cosmics";
43 
44  if (!isEval) {
45  std::cout << "\n\nYou're running over the training datasets....Are you sure you want to do this?"
46  << "\nThey should be reserved for training, not evaluating!"
47  << std::endl;
48  }
49 
50  std::string OutName = "BDTSpectra-"+sisData+"-"+sFHC_RHC+"-"+sisBPF+"-"+sIdentifier+"-"+sisEval+".root";
51 
52  // And write them out!
53  std::cout << "\n\nLots of input parameters, so I'm going to write them out now....;"
54  << "\n\t isFHC " << (isFHC == true ? "true --> Run in FHC " :"false --> Run in RHC" ) << " ==> sFHC_RHC is " << sFHC_RHC
55  << "\n\t isBPF " << (isBPF == true ? "true --> BPF training" :"false --> Kal training" ) << " ==> sisBPF is " << sisBPF
56  << "\n\t Only making spectra, therefore running over data."
57  << "\n\t Identifier is " << sIdentifier
58  << "\n\t isEval " << (isEval == true ? "true --> Evaluation dataset":"false --> Training dataset" ) << " ==> sisEval is " << sisEval
59  << "\n\t Quant " << (Quant == true ? "true --> Make quantile plot":"false --> Just All Quant" )
60  << "\n\n OutName is therefore " << OutName
61  << "\n\n" << std::endl;
62 
63  // -----------------------------------------------------------
64  // --- Now lets get our Quartile cuts
65  // -----------------------------------------------------------
66  /*
67  std::string fdspecfile = (std::string)std::getenv("NUMUDATA_DIR")+"ana2018/Quantiles/quantiles__"+sFHC_RHC+"_full__numu2018.root";
68  //std::string fdspecfile = "/pnfs/nova/persistent/analysis/numu/Ana2018/official/Quantiles/quantiles__"+sFHC_RHC+"_full__numu2018.root";
69  TFile* inFile = TFile::Open( fdspecfile.c_str() );
70  TH2 *FDSpec2D = (TH2*)inFile->FindObjectAny( "FDSpec2D" ); //TString("FD_full_") + TString(sFHC_RHC) );
71  // How many quantiles?
72  const int NHadEFracQuantiles = 4; // defines how many divisions of hadEFrac are used
73  // Make my cuts.
74  //std::vector<Cut> HadEFracQuantCuts = QuantileAndPIDCutsFromTH2(FDSpec2D, kNumuCCOptimisedAxis, kHadEFracAxis, NHadEFracQuantiles, isFHC);
75  std::vector<Cut> HadEFracQuantCuts = QuantileCutsFromTH2(FDSpec2D, kNumuCCOptimisedAxis, kHadEFracAxis, NHadEFracQuantiles);
76  */
77 
78  // -----------------------------------------------------------
79  // --- Lets figure out which definition I want to run over....
80  // -----------------------------------------------------------
81  std::string MyDef = "karlwarb-prod4reco-cosmicrej-fd-cosmics-loose-"+sFHC_RHC+"-"+sIdentifier+"_"+sisEval;
82  //std::string MyDef = "karlwarb-prod4reco-cosmicrej-fd-cosmics-"+sFHC_RHC+"-"+sIdentifier+"_"+sisEval;
83 
84  // --- And declare my loader.
85  SpectrumLoader loader( MyDef );
86  loader.SetSpillCut(kStandardSpillCuts);
87 
88  // -----------------------------------------------------------
89  // --- What cuts do I want to use?
90  // -----------------------------------------------------------
91  Cut BaseCut = kNumuQuality && kNumuContainFD2017 && kNumuPID2018 && kInCosmicTimingWindow;
92  std::vector<Cut> Cuts;
93  std::vector<std::string> CutNames;
94  std::cout << "Going into my function." << std::endl;
95 
96  MakeCutVec( sIdentifier, isFHC, isBPF, Quant, BaseCut, Cuts, CutNames );
97  unsigned int NumCuts = Cuts.size();
98 
99  std::cout << "Out of my function. CutNames has size " << NumCuts << std::endl;
100  for (unsigned int cc=0; cc<NumCuts; ++cc) { std::cout << "\t Name["<<cc<<"] is " << CutNames[cc] << std::endl; }
101 
102  // -----------------------------------------------------------
103  // --- Any custom Hist Axes?
104  // -----------------------------------------------------------
105  const unsigned int NumBins = 100;
106  const Binning NumHitBins = Binning::Simple(NumBins, 0 , 1000 );
107  const Binning LengthBins = Binning::Simple(NumBins, 0 , 20 );
108  const Binning RatioBins = Binning::Simple(NumBins, -1.01, 1.01 );
109  const Binning XRangeBins = Binning::Simple(NumBins, -8 , 8 );
110  const Binning YRangeBins = Binning::Simple(NumBins, -8 , 8 );
111  const Binning ZRangeBins = Binning::Simple(NumBins, -0.1 , 60.1 );
112 
113  // -----------------------------------------------------------
114  // --- Define my spectra.
115  // -----------------------------------------------------------
116  // Energy quanitities
117  //Spectrum *sTrueNuE [NumCuts];
118  Spectrum *sRecoNuE [NumCuts];
119  //Spectrum *sNuMuPtOP [NumCuts];
120  /*
121  // --- My New Variables.
122  Spectrum *sCosRejPID;
123  Spectrum *sFHCMLPPID;
124  Spectrum *sRHCMLPPID;
125  Spectrum *sFHCBDTPID;
126  Spectrum *sRHCBDTPID;
127  // -- Some 2D spectra vs Reco E and True E
128  // Reco
129  Spectrum *sRecoE_CosRejPID;
130  Spectrum *sRecoE_FHCMLPPID;
131  Spectrum *sRecoE_RHCMLPPID;
132  Spectrum *sRecoE_FHCBDTPID;
133  Spectrum *sRecoE_RHCBDTPID;
134  // True
135  Spectrum *sTrueE_CosRejPID;
136  Spectrum *sTrueE_FHCMLPPID;
137  Spectrum *sTrueE_RHCMLPPID;
138  Spectrum *sTrueE_FHCBDTPID;
139  Spectrum *sTrueE_RHCBDTPID;
140  */
141  // -----------------------------------------------------------
142  // --- Create the spectra classes
143  // -----------------------------------------------------------
144  for (unsigned int cc=0; cc<NumCuts; ++cc) {
145  Cut ThisCut = Cuts[cc];
146  // Energy quanitities
147  //sTrueNuE [cc] = new Spectrum( "True Neutrino E", kNumuCCEOptimisedBinning, loader, kTrueE , ThisCut );
148  sRecoNuE [cc] = new Spectrum( "Reco Neutrino E", kNumuCCEOptimisedBinning, loader, kCCE , ThisCut );
149  //sNuMuPtOP [cc] = new Spectrum( "Muon Pt Over P" , RatioBins , loader, kNumuMuonPtP, ThisCut );
150  }
151  /*
152  // --- My New Variables.
153  sCosRejPID = new Spectrum( "Ana18 CosRej" , RatioBins, loader, kNumuContPID , BaseCut );
154  sFHCMLPPID = new Spectrum( "FHC MLP CosRej", RatioBins, loader, kFHCMLPCosRej, BaseCut );
155  sRHCMLPPID = new Spectrum( "RHC MLP CosRej", RatioBins, loader, kRHCMLPCosRej, BaseCut );
156  sFHCBDTPID = new Spectrum( "FHC BDT CosRej", RatioBins, loader, kFHCBDTCosRej, BaseCut );
157  sRHCBDTPID = new Spectrum( "RHC BDT CosRej", RatioBins, loader, kRHCBDTCosRej, BaseCut );
158  // -- Some 2D spectra vs Reco E and True E
159  // Reco
160  sRecoE_CosRejPID = new Spectrum( "RecoE vs Ana18 CosRej" , loader, kNumuCCEOptimisedBinning, kCCE , RatioBins, kNumuContPID , BaseCut );
161  sRecoE_FHCMLPPID = new Spectrum( "RecoE vs FHC MLP CosRej", loader, kNumuCCEOptimisedBinning, kCCE , RatioBins, kFHCMLPCosRej, BaseCut );
162  sRecoE_RHCMLPPID = new Spectrum( "RecoE vs RHC MLP CosRej", loader, kNumuCCEOptimisedBinning, kCCE , RatioBins, kRHCMLPCosRej, BaseCut );
163  sRecoE_FHCBDTPID = new Spectrum( "RecoE vs FHC BDT CosRej", loader, kNumuCCEOptimisedBinning, kCCE , RatioBins, kFHCBDTCosRej, BaseCut );
164  sRecoE_RHCBDTPID = new Spectrum( "RecoE vs RHC BDT CosRej", loader, kNumuCCEOptimisedBinning, kCCE , RatioBins, kRHCBDTCosRej, BaseCut );
165  // True
166  sTrueE_CosRejPID = new Spectrum( "TrueE vs Ana18 CosRej" , loader, kNumuCCEOptimisedBinning, kTrueE, RatioBins, kNumuContPID , BaseCut );
167  sTrueE_FHCMLPPID = new Spectrum( "TrueE vs FHC MLP CosRej", loader, kNumuCCEOptimisedBinning, kTrueE, RatioBins, kFHCMLPCosRej, BaseCut );
168  sTrueE_RHCMLPPID = new Spectrum( "TrueE vs RHC MLP CosRej", loader, kNumuCCEOptimisedBinning, kTrueE, RatioBins, kRHCMLPCosRej, BaseCut );
169  sTrueE_FHCBDTPID = new Spectrum( "TrueE vs FHC BDT CosRej", loader, kNumuCCEOptimisedBinning, kTrueE, RatioBins, kFHCBDTCosRej, BaseCut );
170  sTrueE_RHCBDTPID = new Spectrum( "TrueE vs RHC BDT CosRej", loader, kNumuCCEOptimisedBinning, kTrueE, RatioBins, kRHCBDTCosRej, BaseCut );
171  */
172  // --- Do it!
173  loader.Go();
174 
175  // -----------------------------------------------------------
176  // --- Finally, write everything to my output file.
177  // -----------------------------------------------------------
178  TFile *outFile = new TFile( OutName.c_str(), "RECREATE" );
179  outFile->cd();
180  for (unsigned int cc=0; cc<NumCuts; ++cc) {
181  // Energy quanitities
182  //sTrueNuE [cc] -> SaveTo( outFile->mkdir( TString("TrueNuE" )+TString(CutNames[cc]) ) );
183  sRecoNuE [cc] -> SaveTo( outFile->mkdir( TString("RecoNuE_" )+TString(CutNames[cc]) ) );
184  //sNuMuPtOP [cc] -> SaveTo( outFile->mkdir( TString("MuPtOverP")+TString(CutNames[cc]) ) );
185  }
186  /*
187  // --- My New Variables.
188  sCosRejPID -> SaveTo( outFile->mkdir( TString("Ana18CosRej" ) ) );
189  sFHCMLPPID -> SaveTo( outFile->mkdir( TString("FHCMLPCosRej") ) );
190  sRHCMLPPID -> SaveTo( outFile->mkdir( TString("RHCMLPCosRej") ) );
191  sFHCBDTPID -> SaveTo( outFile->mkdir( TString("FHCBDTCosRej") ) );
192  sRHCBDTPID -> SaveTo( outFile->mkdir( TString("RHCBDTCosRej") ) );
193  // -- Some 2D spectra vs Reco E and True E
194  // Reco
195  sRecoE_CosRejPID -> SaveTo( outFile->mkdir( TString("RecoE_vs_Ana18CosRej" ) ) );
196  sRecoE_FHCMLPPID -> SaveTo( outFile->mkdir( TString("RecoE_vs_FHCMLPCosRej") ) );
197  sRecoE_RHCMLPPID -> SaveTo( outFile->mkdir( TString("RecoE_vs_RHCMLPCosRej") ) );
198  sRecoE_FHCBDTPID -> SaveTo( outFile->mkdir( TString("RecoE_vs_FHCBDTCosRej") ) );
199  sRecoE_RHCBDTPID -> SaveTo( outFile->mkdir( TString("RecoE_vs_RHCBDTCosRej") ) );
200  // True
201  sTrueE_CosRejPID -> SaveTo( outFile->mkdir( TString("TrueE_vs_Ana18CosRej" ) ) );
202  sTrueE_FHCMLPPID -> SaveTo( outFile->mkdir( TString("TrueE_vs_FHCMLPCosRej") ) );
203  sTrueE_RHCMLPPID -> SaveTo( outFile->mkdir( TString("TrueE_vs_RHCMLPCosRej") ) );
204  sTrueE_FHCBDTPID -> SaveTo( outFile->mkdir( TString("TrueE_vs_FHCBDTCosRej") ) );
205  sTrueE_RHCBDTPID -> SaveTo( outFile->mkdir( TString("TrueE_vs_RHCBDTCosRej") ) );
206  */
207  // --- And now close and we'll be done!
208  outFile->Close();
209 
210 } // --- DONE! ---
ana::ZRangeBins
const Binning ZRangeBins
Definition: ComparisonHints.h:36
ana::Binning
Represent the binning of a Spectrum&#39;s x-axis.
Definition: Binning.h:16
ana
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
ana::SaveTo
void SaveTo(const osc::IOscCalc &x, TDirectory *dir, const std::string &name)
Definition: LoadFromFile.cxx:105
ana::YRangeBins
const Binning YRangeBins
Definition: ComparisonHints.h:35
ana::SpectrumLoaderBase::SetSpillCut
void SetSpillCut(const SpillCut &cut)
Definition: SpectrumLoaderBase.cxx:64
ana::Spectrum
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:40
Evaluate_BDTMLP_Algs_Spectra
void Evaluate_BDTMLP_Algs_Spectra(std::string sIdentifier, bool isFHC, bool isBPF, bool Quant, bool isEval=true)
Definition: Evaluate_BDTMLP_Algs_Spectra.C:35
ana::kInCosmicTimingWindow
const Cut kInCosmicTimingWindow
Is the event far from the start and ends of the spill ? For FD cosmic selection.
Definition: TimingCuts.cxx:165
outFile
TFile * outFile
Definition: PlotXSec.C:135
ana::kNumuContainFD2017
const Cut kNumuContainFD2017
Definition: NumuCuts2017.h:21
ana::MakeCutVec
void MakeCutVec(std::string sIdentifier, bool isFHC, bool isBPF, bool Quants, Cut &kBaseCut, std::vector< Cut > &Cuts, std::vector< std::string > &CutNames)
Definition: Evaluate_header.h:20
Evaluate_header.h
ana::kCCE
const Var kCCE
Definition: NumuVars.h:21
ana::SpectrumLoader::Go
virtual void Go() override
Load all the registered spectra.
Definition: SpectrumLoader.cxx:66
demo0.loader
loader
Definition: demo0.py:10
cmf::Spectrum
std::vector< float > Spectrum
Definition: Constants.h:610
isFHC
static bool isFHC
Definition: draw_plots_util.h:46
om::cout
OStream cout
Definition: OStream.cxx:6
ana::NumHitBins
const Binning NumHitBins
Definition: ComparisonHints.h:22
ana::SpectrumLoader
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
Definition: SpectrumLoader.h:21
ana::LengthBins
const Binning LengthBins
Definition: ComparisonHints.h:37
cc
void cc()
Definition: test_ana.C:28
ana::kNumuPID2018
const Cut kNumuPID2018([](const caf::SRProxy *sr){std::cout<< "ERROR::kNumuPID2018, cutting on both cvnProd3Train and cvn2017."<< " Neither branch exists anymore. Returning False."<< std::endl;abort();return false;})
Definition: NumuCuts2018.h:22
ana::kNumuCCEOptimisedBinning
const Binning kNumuCCEOptimisedBinning
Optimised binning for numuCCE from L. Vinton. See docdb 16332. This was close to &#39;custC&#39; in that talk...
Definition: Binnings.cxx:28
isBPF
bool isBPF
Definition: Plot_Evals_BDTMLP_Algs.C:40
ana::NumBins
const unsigned int NumBins
Definition: ComparisonHints.h:21
ana::kStandardSpillCuts
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
ana::_Cut< caf::SRProxy >
ana::kNumuQuality
const Cut kNumuQuality
Definition: NumuCuts.h:18
ana::Binning::Simple
static Binning Simple(int n, double lo, double hi, const std::vector< std::string > &labels={})
Definition: Binning.cxx:107
ana::RatioBins
const Binning RatioBins
Definition: ComparisonHints.h:23
CutNames
std::vector< std::string > CutNames
Definition: MakeCutFlow.C:49
string
enum BeamMode string
Generated on Wed Dec 2 2020 05:11:44 for NOvA by   doxygen 1.8.11