make_rockpred.C
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 // Make noextrap predictions for FD rock MC.
3 //
4 // Options for fhc/rhc and numu/nue/both
5 //
6 // Andrew Sutton - ats2mk@virginia.edu
7 // April 2020
8 ////////////////////////////////////////////////////////////////////////////////
9 
10 
11 #pragma once
12 
13 #include "CAFAna/Core/Binning.h"
14 #include "CAFAna/Core/Cut.h"
15 #include "CAFAna/Core/Loaders.h"
16 #include "CAFAna/Core/Utilities.h"
17 
21 #include "CAFAna/Cuts/SpillCuts.h"
22 
24 #include "CAFAna/Vars/Vars.h"
28 
31 
32 #include "TFile.h"
33 
34 #include <iostream>
35 #include <cmath>
36 
37 using namespace ana;
38 
39 
41 {
42  // What are we doing?
43  std::string horn = "";
44  if(options.find("fhc") != std::string::npos)
45  horn = "fhc";
46  else if(options.find("rhc") != std::string::npos)
47  horn = "rhc";
48  else {
49  std::cout << "\n***** No horn set *****" << std::endl;
50  std::abort();
51  }
52  std::cout << "\n***** Set horn to: " + horn + " *****" << std::endl;
53 
54  std::string ana = "";
55  if(options.find("nue") != std::string::npos)
56  ana = "nue";
57  else if(options.find("numu") != std::string::npos)
58  ana = "numu";
59  else if(options.find("both") != std::string::npos)
60  ana = "both";
61  else {
62  std::cout << "\n***** No analysis set *****" << std::endl;
63  std::abort();
64  }
65  std::cout << "\n***** Set ana to: " + ana + " *****" << std::endl;
66 
67  // Defnames and loader
68  std::string defNS = "prod_caf_R19-11-18-prod5reco.o_fd_genie_N1810j0211a_nonswap_" + horn + "_nova_v08_full_v1";
69  std::string defFS = "prod_caf_R19-11-18-prod5reco.o_fd_genie_N1810j0211a_fluxswap_" + horn + "_nova_v08_full_v1";
70 
75 
76  // Doing some setup here
80  int nquant = 4;
81 
82  // Now for some real stuff
83  std::vector<IPrediction*> preds_nue;
84  std::vector<IPrediction*> preds_numu;
85 
86  if (ana == "numu" || ana == "both") {
87  // All quants first
88  NoExtrapGenerator *gen = new NoExtrapGenerator(axes[0], cuts[0], weight);
89  preds_numu.push_back(gen->Generate(loaders).release());
90 
91  // get the quantiles
92  std::string infile_quant = (std::string)std::getenv("NUMUDATA_DIR")+"lib/ana2020/Quantiles/quantiles_" + horn + "_full_numu2020.root";
93  TFile* infile = new TFile(infile_quant.c_str());
94  TH2* FDSpec2D = (TH2*)infile->FindObjectAny("FDSpec2D");
95  std::vector<Cut> quantCuts = QuantileCutsFromTH2(FDSpec2D, axes[0], kHadEFracAxis, nquant);
96  infile->Close();
97  delete infile;
98 
99  for(int q = 0; q < nquant; ++q) {
100  NoExtrapGenerator *qgen = new NoExtrapGenerator(axes[0], cuts[0] && quantCuts[q], weight);
101  preds_numu.push_back(qgen->Generate(loaders).release());
102  }
103  }
104 
105  if (ana == "nue" || ana == "both") {
106  NoExtrapGenerator *gen0 = new NoExtrapGenerator(axes[1], cuts[1], weight);
107  NoExtrapGenerator *gen1 = new NoExtrapGenerator(axes[2], cuts[1], weight);
108  preds_nue.push_back(gen0->Generate(loaders).release());
109  preds_nue.push_back(gen1->Generate(loaders).release());
110  }
111 
112  loaders.Go();
113 
114  TFile* outFile = new TFile(("fdrock_" + horn + "_" + ana + "2020.root").c_str(), "RECREATE");
115 
116  if(ana == "nue" || ana == "both") {
117  preds_nue[0]->SaveTo(outFile, "pred_FDRock_nue");
118  preds_nue[1]->SaveTo(outFile, "pred_FDRock_nue_merge");
119  }
120 
121  if(ana == "numu" || ana == "both") {
122  preds_numu[0]->SaveTo(outFile, "pred_FDRock_numu");
123 
124  for(int q = 1; q < nquant+1; ++q) {
125  preds_numu[q]->SaveTo(outFile, Form("pred_FDRock_numu_EhadFracQuant%d", q));
126  }
127  }
128 
129  outFile->Close();
130 }
131 
Far Detector at Ash River.
Definition: SREnums.h:11
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const int nquant
Definition: getContProf.C:25
const Var kPPFXFluxCVWgt
Definition: PPFXWeights.h:16
const Var weight
Collection of SpectrumLoaders for many configurations.
Definition: Loaders.h:23
const HistAxis kHadEFracAxis("E_{had.} / E_{#nu}", Binning::Simple(200, 0, 1), kHadEFrac)
HistAxis that implements Hadronic Energy fraction binning used by L.Vinton to derive Hadronic Energy ...
Definition: HistAxes.h:30
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
TFile * outFile
Definition: PlotXSec.C:135
string infile
const HistAxis kNue2020AxisMergedPeripheral("NuE Energy / Analysis Bin", kNue2020BinningMergedPeripheral, kNue2020AnaBinMergedPeripheral)
Definition: NueCuts2020.h:200
std::string getenv(std::string const &name)
const HistAxis kNue2020Axis("NuE Energy / Analysis Bin", kNue2020Binning, kNue2020AnaBin)
Use this Axis for Ana2020, official Axis.
Definition: NueCuts2020.h:195
const Cut kNumu2020FD
Definition: NumuCuts2020.h:59
const HistAxis kNumuCCOptimisedAxis("Reconstructed Neutrino Energy (GeV)", kNumuCCEOptimisedBinning, kCCE)
HistAxis that implements optimised numuCCE from L. Vinton. See docdb 16332. This was close to &#39;custC&#39;...
Definition: HistAxes.h:24
OStream cout
Definition: OStream.cxx:6
std::unique_ptr< IPrediction > Generate(Loaders &loaders, const SystShifts &shiftMC=kNoShift) const override
std::vector< Loaders * > loaders
Definition: syst_header.h:386
std::vector< Cut > QuantileCutsFromTH2(TH2 *quantileHist, const HistAxis &independentAxis, const HistAxis &quantileAxis, const unsigned int &numQuantiles, const bool verbose)
: Do the same as the QuantileCuts function but taking in the TH2 instead of making it...
void make_rockpred(std::string options)
Definition: make_rockpred.C:40
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
const Cut kNue2020FDAllSamples
Definition: NueCuts2020.h:84
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
const Var kXSecCVWgt2020
Definition: XsecTunes.h:106
enum BeamMode string