get_numu_data_histogram.C
Go to the documentation of this file.
2 #include "CAFAna/Core/Spectrum.h"
7 #include "CAFAna/Vars/HistAxes.h"
9 
10 using namespace ana;
11 
12 #include "TFile.h"
13 #include <iostream>
14 
15 void get_numu_data_histogram(bool isFHC = true)
16 {
17  /// File for making data histograms for fit
18  /// Use it for extracting OLD data hitstograms from decafs made in 2019 release
19 
20  bool getfiles=0;
21  std::string horn_tag = "fhc";
22  if (!isFHC) horn_tag = "rhc";
23 
24  if (getfiles){
25 
26  std::string fnameData = "defname: prod4_sumrestricteddecaf_fd_numi_fhc_full_goodruns_numu2019";
27  //if (!isFHC) fnameData = "defname: prod4_sumrestricteddecaf_fd_numi_rhc_Ana2018_goodruns_numu2019";
28  if (!isFHC) fnameData = "defname: prod4_sumrestricteddecaf_fd_numi_rhc_full_goodruns_numu2019";
29 
30  SpectrumLoader loaderData(fnameData);
31  loaderData.SetSpillCut(kStandardSpillCuts);
32 
33  std::string infile_quant = "/cvmfs/nova.opensciencegrid.org/externals/numudata/v00.00/NULL/ana2018/Quantiles/quantiles__"+horn_tag+"_full__numu2018.root";
34  TFile* infile = new TFile(infile_quant.c_str());
35  TH2* FDSpec2D = (TH2*)infile->FindObjectAny("FDSpec2D");
36  std::vector<Cut> HadEFracQuantCuts = QuantileCutsFromTH2(FDSpec2D, kNumuCCOptimisedAxis, kHadEFracAxis, 4);
37  infile->Close();
38  delete infile;
39 
40  std::vector<Spectrum*> data;
41  data.push_back(new Spectrum(loaderData, kNumuCCOptimisedAxis, kNumuCutFD2018&&kInBeamSpill, kNoShift));
42  for(auto thisCut : HadEFracQuantCuts){
43  data.push_back(new Spectrum(loaderData, kNumuCCOptimisedAxis, kNumuCutFD2018 && kInBeamSpill && thisCut));
44  }
45 
46  loaderData.Go();
47 
48  TFile* fout = new TFile(("numu_data_2019_"+horn_tag+".root").c_str(), "RECREATE");
49  for(int quant=0; quant<5; quant++){
50  data[quant]-> SaveTo(&fout->, Form("data_q%d", quant));
51  }
52  delete fout;
53  }
54 // in order to get spectra from *.root file do:
55  if (!getfiles){
56 
57  TFile * file = new TFile (("numu_data_2019_"+horn_tag+"_test.root").c_str(), "READ");
58  std::vector<Spectrum*> data;
59  for(int quant=0; quant<5; quant++){
60  auto data = Spectrum::LoadFrom(file, Form("data_q%d", quant));
61  TH1* hist = data->ToTH1(data->POT());
62  cout<<" data: "<<hist->Integral()<<" POT "<<data->POT()<<endl;
63  }
64  }
65 
66 }
Cuts and Vars for the 2020 FD DiF Study.
Definition: vars.h:6
void SaveTo(const osc::IOscCalc &x, TDirectory *dir, const std::string &name)
void get_numu_data_histogram(bool isFHC=true)
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
void SetSpillCut(const SpillCut &cut)
const XML_Char const XML_Char * data
Definition: expat.h:268
string infile
static std::unique_ptr< Spectrum > LoadFrom(TDirectory *dir, const std::string &name)
Definition: Spectrum.cxx:535
virtual void Go() override
Load all the registered spectra.
std::vector< float > Spectrum
Definition: Constants.h:610
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
static bool isFHC
const SystShifts kNoShift
Definition: SystShifts.cxx:21
OStream cout
Definition: OStream.cxx:6
Collaborates with Spectrum and OscillatableSpectrum to fill spectra from CAF files.
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...
TFile * file
Definition: cellShifts.C:17
const Cut kInBeamSpill([](const caf::SRProxy *sr){if(sr->spill.run > util::kLastBadTimingRun) return kInBeamSpill_main(sr);else return kInBeamSpill_main(sr)||kInBeamSpill_shifted(sr);}, [](const caf::SRSpillProxy *spill){if(spill->run > util::kLastBadTimingRun) return kInBeamSpill_main.Livetime(spill);else return kInBeamSpill_main.Livetime(spill)+kInBeamSpill_shifted.Livetime(spill);}, [](const caf::SRSpillProxy *spill) -> double{return spill->spillpot;})
Does the event fall inside the window we call the beam spill?
Definition: TimingCuts.h:8
const Cut kNumuCutFD2018
Definition: NumuCuts2018.h:39
const SpillCut kStandardSpillCuts
Apply this unless you&#39;re doing something special.
Definition: SpillCuts.h:49
enum BeamMode string