CAFReweightExample.C
Go to the documentation of this file.
1 //Example macro for demonstrating the implementation of caf based reweighting
2 //First run CAFReweight/loadlibs.C
3 //Currently I don't think that this will run under the cafe. It doesn't appear
4 //that the cafe loads the correct libraries that are needed for reweighting. That
5 //can certainly be done if there is a desire to do it. So if reweighting is
6 //needed for a macro run as a bare root using the above loadlibs.C macro.
7 
8 // C/C++ includes
9 #include <math.h>
10 #include <vector>
11 #include <fstream>
12 
13 #include "TTree.h"
14 #include "TH1D.h"
15 #include "TFile.h"
16 #include "TLorentzVector.h"
17 #include "TMath.h"
18 #include "TVector.h"
19 #include "TCanvas.h"
20 
22 #include "StandardRecord/SRGenerator.h"
25 
26 
27 
28 void CAFReweightExample(char *nIn, char *nOut) {
29 
30  cafrwgt::CAFReweight* fGrwgt[3];
31  TH1D *wgt_dist[3];
32  for(int i = 0; i < 3; i++) {
33 
34  double sigma = (double)(i+1);
35  fGrwgt[i] = new cafrwgt::CAFReweight();
36  fGrwgt[i]->AddReweightValue(rwgt::fReweightMaCCQE, sigma);
37  fGrwgt[i]->AddReweightValue(rwgt::fReweightMaCCRES, sigma);
38  fGrwgt[i]->AddReweightValue(rwgt::fReweightMaNCRES, sigma);
39  fGrwgt[i]->AddReweightValue(rwgt::fReweightRvpCC1pi, sigma);
40  fGrwgt[i]->AddReweightValue(rwgt::fReweightRvnCC1pi, sigma);
41  fGrwgt[i]->AddReweightValue(rwgt::fReweightRvpCC2pi, sigma);
42  fGrwgt[i]->AddReweightValue(rwgt::fReweightRvnCC2pi, sigma);
43  fGrwgt[i]->AddReweightValue(rwgt::fReweightRvpNC1pi, sigma);
44  fGrwgt[i]->AddReweightValue(rwgt::fReweightRvnNC1pi, sigma);
45  fGrwgt[i]->AddReweightValue(rwgt::fReweightRvpNC2pi, sigma);
46  fGrwgt[i]->AddReweightValue(rwgt::fReweightRvnNC2pi, sigma);
47 
48  fGrwgt[i]->Configure();
49 
50  }
51  wgt_dist[0] = new TH1D("h1", "Weight 3#sigma", 200, -2, 2);
52  wgt_dist[1] = new TH1D("h2", "Weight 2#sigma", 200, -2, 2);
53  wgt_dist[2] = new TH1D("h3", "Weight 3#sigma", 200, -2, 2);
54 
55  TFile *fIn = new TFile(nIn, "read");
56  TTree *cafT = (TTree*)fIn->Get("recTree");
57 
59  cafT->SetBranchAddress("rec", &sr);
60 
61  int nEntries = cafT->GetEntries();
62 
63  for(int iEntry = 0; iEntry < nEntries; iEntry++) {
64  cafT->GetEntry(iEntry);
65  std::cout << "Loop over tree entries" << std::endl;
66 
67  if(iEntry % 50000 == 0) {
68  std::cout << iEntry << " of " << nEntries << std::endl;
69  }
70 
71 
72  double wgt[3];
73  for(int i = 0; i < 3; i++) {
74  wgt[i] = fGrwgt[i]->CalcWeight(sr->mc.gen);
75  if(sr->mc.gen.size() > 0) {
76  std::cout << "Found a Neutrino Interaction" << std::endl;
77  wgt_dist[i]->Fill(wgt[i]);
78  }
79  }
80 
81  }
82 
83  TFile *fOut = new TFile(nOut, "recreate");
84  fOut->cd();
85  for(int i = 0; i < 3; i++) {
86  wgt_dist[i]->Write();
87  }
88  fOut->Close();
89 
90 
91 }
void AddReweightValue(rwgt::ReweightLabel_t rLabel, double value, unsigned int i=0)
Definition: CAFReweight.cxx:66
double CalcWeight(const std::vector< caf::SRGenerator > &gen, unsigned int i=0)
caf::StandardRecord * sr
const ana::Var wgt
double sigma(TH1F *hist, double percentile)
Wrapper for reweight neutrino interactions with GENIE.
OStream cout
Definition: OStream.cxx:6
The StandardRecord is the primary top-level object in the Common Analysis File trees.
SRTruthBranch mc
Truth branch for MC: energy, flavor, etc.
void Configure(int i=-1)
Definition: CAFReweight.cxx:74
void CAFReweightExample(char *nIn, char *nOut)