makeRealDataFluxes.C
Go to the documentation of this file.
1 // makeFakeData.C
2 
3 #include "CAFAna/Core/Cut.h"
4 #include "CAFAna/Vars/Vars.h"
5 #include "CAFAna/Core/Binning.h"
6 #include "CAFAna/Cuts/Cuts.h"
8 #include "CAFAna/Core/Spectrum.h"
12 #include "CAFAna/Systs/Systs.h"
15 
18 
24 
25 #include "TFile.h"
26 #include "TH1.h"
27 #include "TH2.h"
28 
29 #include <iostream>
30 #include <cmath>
31 #include <sstream>
32 #include <string>
33 #include <fstream>
34 #include <iomanip>
35 #include <map>
36 
37 using namespace ana;
38 
39 std::map<std::string, ana::HistAxis> unfoldingVars{
40  {"data_eavail", kRecoMuKEVsCosVsEavailStandardAxis},
41  {"data_enu", kRecoEStandardAxis},
42  {"data_q2", kRecoQ2StandardAxis}
43 };
44 
45 std::map<std::string, ana::Cut> recoCuts{
46  {"data_eavail", kAllNumuCCCuts},
47  {"data_enu", kAllNumuCC1DCuts},
48  {"data_q2", kAllNumuCC1DCuts}
49 };
50 
52 const ana::Var ppfxXSecWeights = ana::VarFromNuTruthVar(ppfxXSecWeightsST, 1);
53 const vector<pair<const string, const ana::Var>> weights{
54  {"unwgt", kUnweighted},
55  {"ppfxwgt", kPPFXFluxCVWgt},
56  {"ppfxxsecwgt", ppfxXSecWeights}
57 };
58 
59 void makeRealDataFluxes(std::string outputfile = "fRealDataFluxes.root"){
60 
61  std::string data_def = "def_snapshot prod_caf_R17-11-14-prod4reco.remid-hotfix.b_nd_numi_fhc_full_v1_addShortSimpleCVN_goodruns";
62  std::string nominal_def = "def_snapshot nd_fhc_remid-hotfix_caf_minus_muonid_training_minus_fakedata";
63 
64  SpectrumLoader * dataloader = new SpectrumLoader(data_def);
65  SpectrumLoader * nomloader = new SpectrumLoader(nominal_def);
66 
67  std::map<std::string, std::map<std::string, Spectrum*> > dataSpectrumsReco;
68  ////////////// Fake data /////////////////////////////////////////////////
69  for(std::pair<std::string, ana::HistAxis> unfoldingVar : unfoldingVars){
70  for(std::pair<const string, const ana::Var> weight : weights){
71  ana::HistAxis& recoAxis = unfoldingVar.second;
72  Spectrum * dataReco = new Spectrum(*dataloader, recoAxis, recoCuts.at(unfoldingVar.first), kNoShift, weight.second);
73  dataSpectrumsReco[unfoldingVar.first][weight.first] = dataReco;
74  }
75  }
76 
77  ////////////// Fluxes ////////////////////////////////////////////////////
78  const int wanted_pdg = 14;
79  const Binning FluxBinning = Binning::Simple(480, 0.0, 120.0);
80 
81  // Vertex default to the definitions in NumuCCIncCuts.h
82  TVector3 flux_vtxmin = vtxmin;
83  TVector3 flux_vtxmax = vtxmax;
84 
85  Spectrum * fluxData = DeriveFlux(*dataloader, FluxBinning, wanted_pdg, &flux_vtxmin, &flux_vtxmax, kNoShift, kPPFXFluxCVWgtST);
86  Spectrum * fluxCV = DeriveFlux(*nomloader, FluxBinning, wanted_pdg, &flux_vtxmin, &flux_vtxmax, kNoShift, kPPFXFluxCVWgtST);
87 
88  ////////////// Go() //////////////////////////////////////////////////////
89  dataloader->Go();
90  nomloader->Go();
91 
92  ////////////// Saving ////////////////////////////////////////////////////
93  TFile * fOut = new TFile(outputfile.c_str(), "RECREATE");
94  TDirectory * dir;
95 
96  for (std::pair<std::string, std::map<std::string, Spectrum*> > dataReco : dataSpectrumsReco){
97  // Data Spectrums
98  dir = fOut->mkdir(dataReco.first.c_str());
99  dir = dir->mkdir("reco");
100  dir->cd();
101 
102  std::map<std::string, Spectrum*>& weightSpecs = dataReco.second;
103  for(std::pair<std::string, Spectrum*> weightSpec : weightSpecs)
104  weightSpec.second->SaveTo(dir->mkdir(weightSpec.first.c_str()));
105  }
106 
107  // Flux
108  dir = fOut->mkdir("fluxes");
109  dir->cd();
110  fluxData->SaveTo(dir->mkdir("Data"));
111  fluxCV->SaveTo(dir->mkdir("CV"));
112 }
const HistAxis kRecoQ2StandardAxis("Reco Q2 (GeV)", q2bins, kRecoq2)
Represent the binning of a Spectrum&#39;s x-axis.
Definition: Binning.h:16
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
const Var kPPFXFluxCVWgt
Definition: PPFXWeights.h:16
const Var weight
const TVector3 vtxmin(-130,-176, 225)
void makeRealDataFluxes(std::string outputfile="fRealDataFluxes.root")
std::vector< double > Spectrum
Definition: Constants.h:746
Representation of a spectrum in any variable, with associated POT.
Definition: Spectrum.h:39
std::map< std::string, ana::HistAxis > unfoldingVars
const vector< pair< const string, const ana::Var > > weights
virtual void Go() override
Load all the registered spectra.
void SaveTo(TDirectory *dir, const std::string &name) const
Definition: Spectrum.cxx:534
const NuTruthVar kXSecCVWgt2018_smallerDISScale_NT
Definition: XsecTunes.h:47
const SystShifts kNoShift
Definition: SystShifts.cxx:22
const Cut kAllNumuCC1DCuts
Var VarFromNuTruthVar(const NuTruthVar &stv, double _default)
Definition: Var.cxx:8
const ana::NuTruthVar ppfxXSecWeightsST
Spectrum * DeriveFlux(SpectrumLoaderBase &loader, const Binning &bins, int pdg, const TVector3 *min, const TVector3 *max, const SystShifts &shift, const NuTruthVar weight)
Definition: Flux.cxx:58
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
TDirectory * dir
Definition: macro.C:5
const HistAxis kRecoMuKEVsCosVsEavailStandardAxis("Reco. T_{#mu} vs cos #{theta} vs Available Energy (GeV)", angvsmukevseavailbins, kRecoMuKEVsCosVsEavail)
std::map< std::string, ana::Cut > recoCuts
const TVector3 vtxmax(160, 160, 1000)
Template for Vars applied to any type of object.
Definition: FwdDeclare.h:12
const HistAxis kRecoEStandardAxis("Reconstructed Neutrino Energy (GeV)", enubins, kRecoE)
static Binning Simple(int n, double lo, double hi, const std::vector< std::string > &labels={})
Definition: Binning.cxx:107
const Var kUnweighted
Definition: Var.h:14
const NuTruthVar kPPFXFluxCVWgtST([](const caf::SRNeutrinoProxy *nu){ if(nu->rwgt.ppfx.cv!=nu->rwgt.ppfx.cv){return 1.f;}if(nu->rwgt.ppfx.cv >90){return 1.f;}return float(nu->rwgt.ppfx.cv);})
weight events with the flux PPFX Central value correction.
Definition: PPFXWeights.h:12
const Cut kAllNumuCCCuts
const ana::Var ppfxXSecWeights
enum BeamMode string