Public Member Functions | Static Public Member Functions | Public Attributes | Private Member Functions | Private Attributes | Static Private Attributes | List of all members
NeutrinoFluxReweight::ThinTargetMC Class Reference

A class to manage the MC value for thin target. More...

#include <ThinTargetMC.h>

Public Member Functions

double getMCval_pC_X (double incP, double xf, double pt, int pdgcode)
 MC value for this HP production. More...
 
double getMCxs_pC_piK (int genid, double inc_mom)
 Get the MC roduction cross-section pC->pi, K: More...
 
double getMCxs_pC_nucleon (int genid, int pdg, double inc_mom)
 Get the MC roduction cross-section pC->n, p: More...
 

Static Public Member Functions

static ThinTargetMCgetInstance ()
 

Public Attributes

std::vector< std::vector< TH2F * > > hTTScl
 Vector of the scaling histograms: More...
 
std::vector< TH1F * > hTTScl_n
 Vector of the scaling histograms for neutrons: More...
 

Private Member Functions

 ThinTargetMC ()
 

Private Attributes

TFile * fpC_x [8]
 
std::vector< TH2D * > vpC_x [8]
 
std::vector< TH1D * > vpC_n
 
std::vector< std::stringspart_prod
 
std::vector< int > mom_inc
 
std::vector< std::stringspart_qe_corr
 
TFile * fqe_corr [2]
 
std::vector< TH2D * > vqe_corr_p
 
std::vector< TH1D * > vqe_corr_n
 
std::vector< TFile * > fTTscale
 

Static Private Attributes

static ThinTargetMCinstance = 0
 

Detailed Description

A class to manage the MC value for thin target.

Definition at line 14 of file ThinTargetMC.h.

Constructor & Destructor Documentation

NeutrinoFluxReweight::ThinTargetMC::ThinTargetMC ( )
private

Definition at line 11 of file ThinTargetMC.cpp.

11  {
12 
13  const char* ppfxDir = getenv("PPFX_DIR");
14  char dirData[400];
15  sprintf(dirData,"%s/data",ppfxDir);
16  spart_prod.push_back("pip");
17  spart_prod.push_back("pim");
18  spart_prod.push_back("kap");
19  spart_prod.push_back("kam");
20  spart_prod.push_back("prt");
21  spart_prod.push_back("klong");
22  spart_prod.push_back("kshort");
23  spart_prod.push_back("neu"); const int idx_neu = 7;
24  mom_inc.push_back(12);mom_inc.push_back(20);
25  mom_inc.push_back(31);mom_inc.push_back(40);
26  mom_inc.push_back(50);mom_inc.push_back(60);
27  mom_inc.push_back(70);mom_inc.push_back(80);
28  mom_inc.push_back(90);mom_inc.push_back(100);
29  mom_inc.push_back(110);mom_inc.push_back(120);
30 
31  int part_size = spart_prod.size();
32  int mom_size = mom_inc.size();
33 
34  for(int i=0;i<part_size;i++){
35  if(i!=idx_neu){
36  fpC_x[i] = new TFile(Form("%s/MC/FTFP/invxs_%s_FTFP_BERT.root",dirData,spart_prod[i].c_str()),"read");
37  for(int j=0;j<mom_size;j++)vpC_x[i].push_back((TH2D*)fpC_x[i]->Get(Form("xF_pT_%dGeV",mom_inc[j])));
38  }
39  else if(i==idx_neu){
40  fpC_x[i] = new TFile(Form("%s/MC/FTFP/yield_%s_FTFP_BERT.root",dirData,spart_prod[i].c_str()),"read");
41  for(int j=0;j<mom_size;j++)vpC_n.push_back((TH1D*)fpC_x[i]->Get(Form("dndxf_%dGeV",mom_inc[j])));
42  }
43  }
44  //QEL fraction:
45  spart_qe_corr.push_back("prt");
46  spart_qe_corr.push_back("neu");
47  int qe_size = spart_qe_corr.size();
48  for(int i=0;i<qe_size;i++){
49  if(i<qe_size-1){
50  fqe_corr[i] = new TFile(Form("%s/MC/FTFP/invxs_qe_corr_%s.root",dirData,spart_qe_corr[i].c_str()),"read");
51  for(size_t j=0;j<mom_inc.size();j++)vqe_corr_p.push_back((TH2D*)fqe_corr[i]->Get(Form("frac_prod_xF_pT_%dGeV",mom_inc[j])));
52  }
53  else if(i==(qe_size-1)){
54  fqe_corr[i] = new TFile(Form("%s/MC/FTFP/yield_qe_corr_%s.root",dirData,spart_qe_corr[i].c_str()),"read");
55  for(int j=0;j<mom_size;j++)vqe_corr_n.push_back((TH1D*)fqe_corr[i]->Get(Form("frac_prod_xf_%dGeV",mom_inc[j])));
56  }
57  }
58 
59  //Scaling:
60  const int Nscl = 5;
61  const char* cscl_parts[Nscl] = {"pip","pim","kap","kam","prt"};
62  const int Nmomscl = 11;
63  const int mom_scale[Nmomscl] = {12,20,31,40,50,60,70,80,100,120,158};
64  for(Int_t ipart=0;ipart<Nscl;ipart++){
65  ThinTargetMC::fTTscale.push_back(new TFile(Form("%s/SCALING/%s_scaling.root",dirData,cscl_parts[ipart])));
66  }
67  ThinTargetMC::fTTscale.push_back(new TFile(Form("%s/SCALING/%s_scaling.root",dirData,"neu")));
68  //Loading scale histograms:
69  for(Int_t ipart=0;ipart<Nscl;ipart++){
70  std::vector<TH2F*> tmp_scl;
71  for(int imom=0;imom<Nmomscl;imom++){
72  tmp_scl.push_back((TH2F*)ThinTargetMC::fTTscale[ipart]->Get(Form("xF_pT_%dGeV",mom_scale[imom])));
73  }
74  hTTScl.push_back(tmp_scl);
75  tmp_scl.clear();
76  }
77  //neutron scaling
78  for(int imom=0;imom<Nmomscl;imom++){
79  hTTScl_n.push_back((TH1F*)ThinTargetMC::fTTscale[Nscl]->Get(Form("xF_%dGeV",mom_scale[imom])));
80  }
81 
82  }
std::vector< std::string > spart_qe_corr
Definition: ThinTargetMC.h:45
std::vector< TH2D * > vqe_corr_p
Definition: ThinTargetMC.h:47
std::vector< TH1F * > hTTScl_n
Vector of the scaling histograms for neutrons:
Definition: ThinTargetMC.h:33
std::vector< TH1D * > vpC_n
Definition: ThinTargetMC.h:39
std::vector< std::vector< TH2F * > > hTTScl
Vector of the scaling histograms:
Definition: ThinTargetMC.h:31
std::vector< TH2D * > vpC_x[8]
Definition: ThinTargetMC.h:38
std::vector< TFile * > fTTscale
Definition: ThinTargetMC.h:50
fInnerVessel push_back(Point(-578.400000, 0.000000, 0.000000))
std::string getenv(std::string const &name)
Definition: getenv.cc:15
std::vector< std::string > spart_prod
Definition: ThinTargetMC.h:41
std::vector< TH1D * > vqe_corr_n
Definition: ThinTargetMC.h:48

Member Function Documentation

ThinTargetMC * NeutrinoFluxReweight::ThinTargetMC::getInstance ( )
static

Definition at line 209 of file ThinTargetMC.cpp.

209  {
210  if (instance == 0) instance = new ThinTargetMC;
211  return instance;
212  }
static ThinTargetMC * instance
Definition: ThinTargetMC.h:52
double NeutrinoFluxReweight::ThinTargetMC::getMCval_pC_X ( double  incP,
double  xf,
double  pt,
int  pdgcode 
)

MC value for this HP production.

Definition at line 84 of file ThinTargetMC.cpp.

84  {
85 
86  //check:
87  if(incP<12)return -1;
88  if(pdgcode!=211 && pdgcode!=-211 && pdgcode!=321 && pdgcode!=-321 && pdgcode!=2212 && pdgcode!=2112 && pdgcode!=130 && pdgcode!=310)return -1;
89  //idx:
90  int idx_part = -1;
91  int idx_qe_corr = -1;
92  int idx_lowp = -1;
93  int idx_hip = -1;
94  for(size_t i=0;i<mom_inc.size()-1;i++){
95  if(incP>=double(mom_inc[i]) && incP<=double(mom_inc[i+1])){
96  idx_lowp=i;
97  idx_hip =i+1;
98  }
99  }
100  if(idx_lowp==-1){
101  idx_lowp = mom_inc.size()-2;
102  idx_hip = mom_inc.size()-1;
103  }
104  if(pdgcode == 211)idx_part=0;
105  if(pdgcode == -211)idx_part=1;
106  if(pdgcode == 321)idx_part=2;
107  if(pdgcode == -321)idx_part=3;
108  if(pdgcode == 2212){
109  idx_part = 4;
110  idx_qe_corr = 0;
111  }
112  if(pdgcode == 2112){
113  idx_part = 7;
114  idx_qe_corr = 1;
115  }
116  if(pdgcode == 130)idx_part=5;
117  if(pdgcode == 310)idx_part=6;
118 
119  double mcval = 0.0;
120  double qe_corr = 1.0;
121  if(idx_part<7){
122  int binp = vpC_x[idx_part][idx_lowp]->FindBin(xf,pt);
123  double mclow = vpC_x[idx_part][idx_lowp]->GetBinContent(binp);
124  double mchi = vpC_x[idx_part][idx_hip]->GetBinContent(binp);
125  mcval = mclow + (incP-double(mom_inc[idx_lowp]))*(mchi-mclow)/(double(mom_inc[idx_hip])-double(mom_inc[idx_lowp]));
126  if(idx_qe_corr==0){
127  int binqe = vqe_corr_p[idx_lowp]->FindBin(xf,pt);
128  double qelow = vqe_corr_p[idx_lowp]->GetBinContent(binqe);
129  double qehi = vqe_corr_p[idx_hip]->GetBinContent(binqe);
130  qe_corr = qelow + (incP-double(mom_inc[idx_lowp]))*(qehi-qelow)/(double(mom_inc[idx_hip])-double(mom_inc[idx_lowp]));
131  }
132  }
133  else if(idx_part==7){
134  int binp = vpC_n[idx_lowp]->FindBin(xf);
135  double mclow = vpC_n[idx_lowp]->GetBinContent(binp);
136  double mchi = vpC_n[idx_hip]->GetBinContent(binp);
137  mcval = mclow + (incP-double(mom_inc[idx_lowp]))*(mchi-mclow)/(double(mom_inc[idx_hip])-double(mom_inc[idx_lowp]));
138  int binqe = vqe_corr_n[idx_lowp]->FindBin(xf);
139  double qelow = vqe_corr_n[idx_lowp]->GetBinContent(binqe);
140  double qehi = vqe_corr_n[idx_hip]->GetBinContent(binqe);
141  qe_corr = qelow + (incP-double(mom_inc[idx_lowp]))*(qehi-qelow)/(double(mom_inc[idx_hip])-double(mom_inc[idx_lowp]));
142  }
143 
144  mcval /=qe_corr;
145 
146  //check:
147  if(mcval<1.e-12 || mcval!=mcval)return -1;
148 
149  return mcval;
150 
151  }
std::vector< TH2D * > vqe_corr_p
Definition: ThinTargetMC.h:47
std::vector< TH1D * > vpC_n
Definition: ThinTargetMC.h:39
std::vector< TH2D * > vpC_x[8]
Definition: ThinTargetMC.h:38
const double e
std::vector< TH1D * > vqe_corr_n
Definition: ThinTargetMC.h:48
double NeutrinoFluxReweight::ThinTargetMC::getMCxs_pC_nucleon ( int  genid,
int  pdg,
double  inc_mom 
)

Get the MC roduction cross-section pC->n, p:

Definition at line 182 of file ThinTargetMC.cpp.

182  {
183  double xx[13] ={12,20,31,40,50,60,70,80,90,100,110,120,158};
184  double yy[13] ={153386793./197812683.,160302538./197811564.,164508480./197831250.,166391359./197784915.,
185  167860919./197822312.,168882647./197807739.,169681805./197803099.,170311264./197811098.,
186  170860912./197822002.,171309291./197834756.,171651963./197811822.,171991260./197823012.,
187  172902228./197804669.};
188 
189  int idx_lowp = -1;
190  int idx_hip = -1;
191  for(int i=0;i<12;i++){
192  if(inc_mom>=double(xx[i]) && inc_mom<double(xx[i+1])){
193  idx_lowp=i;
194  idx_hip =i+1;
195  }
196  }
197  double frac_low = yy[idx_lowp];
198  double frac_hi = yy[idx_hip];
199  double frac_m = frac_low + (inc_mom-double(xx[idx_lowp]))*(frac_hi-frac_low)/(double(xx[idx_hip])-double(xx[idx_lowp]));
200 
201  if(genid==0 && pdg==2212)return frac_m*243.2435;
202  if(genid>0 && pdg==2212)return frac_m;
203  if(genid==0 && pdg==2112)return frac_m;
204  if(genid>0 && pdg==2112)return frac_m/243.2435;
205  return 1.0;
206 
207  }
double NeutrinoFluxReweight::ThinTargetMC::getMCxs_pC_piK ( int  genid,
double  inc_mom 
)

Get the MC roduction cross-section pC->pi, K:

Definition at line 154 of file ThinTargetMC.cpp.

154  {
155  double xx[13] ={12,20,31,40,50,60,70,80,90,100,110,120,158};
156  double yy[13] ={153386793./197812683.,160302538./197811564.,164508480./197831250.,166391359./197784915.,
157  167860919./197822312.,168882647./197807739.,169681805./197803099.,170311264./197811098.,
158  170860912./197822002.,171309291./197834756.,171651963./197811822.,171991260./197823012.,
159  172902228./197804669.};
160 
161  int idx_lowp = -1;
162  int idx_hip = -1;
163  for(int i=0;i<12;i++){
164  if(inc_mom>=double(xx[i]) && inc_mom<double(xx[i+1])){
165  idx_lowp=i;
166  idx_hip =i+1;
167  }
168  }
169  double frac_low = yy[idx_lowp];
170  double frac_hi = yy[idx_hip];
171  double frac_m = frac_low + (inc_mom-double(xx[idx_lowp]))*(frac_hi-frac_low)/(double(xx[idx_hip])-double(xx[idx_lowp]));
172 
173  if(genid==0)return frac_m*243.2435;
174  else if(genid>0)return frac_m;
175  else{
176  std::cout<<"Something is wrong with gen "<<std::endl;
177  return 1.0;
178  }
179 
180  }
QTextStream & endl(QTextStream &s)

Member Data Documentation

TFile* NeutrinoFluxReweight::ThinTargetMC::fpC_x[8]
private

Definition at line 37 of file ThinTargetMC.h.

TFile* NeutrinoFluxReweight::ThinTargetMC::fqe_corr[2]
private

Definition at line 46 of file ThinTargetMC.h.

std::vector<TFile*> NeutrinoFluxReweight::ThinTargetMC::fTTscale
private

Definition at line 50 of file ThinTargetMC.h.

std::vector< std::vector<TH2F*> > NeutrinoFluxReweight::ThinTargetMC::hTTScl

Vector of the scaling histograms:

Definition at line 31 of file ThinTargetMC.h.

std::vector<TH1F*> NeutrinoFluxReweight::ThinTargetMC::hTTScl_n

Vector of the scaling histograms for neutrons:

Definition at line 33 of file ThinTargetMC.h.

ThinTargetMC * NeutrinoFluxReweight::ThinTargetMC::instance = 0
staticprivate

Definition at line 52 of file ThinTargetMC.h.

std::vector<int> NeutrinoFluxReweight::ThinTargetMC::mom_inc
private

Definition at line 42 of file ThinTargetMC.h.

std::vector<std::string> NeutrinoFluxReweight::ThinTargetMC::spart_prod
private

Definition at line 41 of file ThinTargetMC.h.

std::vector<std::string> NeutrinoFluxReweight::ThinTargetMC::spart_qe_corr
private

Definition at line 45 of file ThinTargetMC.h.

std::vector<TH1D*> NeutrinoFluxReweight::ThinTargetMC::vpC_n
private

Definition at line 39 of file ThinTargetMC.h.

std::vector<TH2D*> NeutrinoFluxReweight::ThinTargetMC::vpC_x[8]
private

Definition at line 38 of file ThinTargetMC.h.

std::vector<TH1D*> NeutrinoFluxReweight::ThinTargetMC::vqe_corr_n
private

Definition at line 48 of file ThinTargetMC.h.

std::vector<TH2D*> NeutrinoFluxReweight::ThinTargetMC::vqe_corr_p
private

Definition at line 47 of file ThinTargetMC.h.


The documentation for this class was generated from the following files: