Classes | Public Member Functions | Private Types | Private Member Functions | Private Attributes | List of all members
nbw::skzpReweight Class Reference

#include <skzpReweight.h>

Classes

struct  LessThan
 
struct  mapkey
 

Public Member Functions

 skzpReweight (std::string fpath="/nova/data/flux/SKZPdata/fluka05ptxf.root", std::string bpath="/nova/data/flux/SKZPdata/IPNDhists.root", int flag=2)
 
 ~skzpReweight ()
 
void SetParams (std::vector< double > fpar, std::vector< double > bpar)
 
void SetFlukParams (std::vector< double > fpar)
 
void SetBeamParams (std::vector< double > bpar)
 
double GetWeight (const simb::MCFlux *mcf, double Enu, int det, int beam)
 
double GetFlukWeight (const simb::MCFlux *mcf)
 
double GetFlukWeight (int ptype, double pT, double pz)
 
double GetBeamWeight (int ntype, double Enu, int det=1, int beam=2)
 

Private Types

typedef std::map< double, double > WeightMap_t
 

Private Member Functions

void FlukConfig ()
 
Conventions::ParticleType_t GeantToEnum (int ptype)
 
std::string PartEnumToString (Conventions::ParticleType_t ptype)
 
void BeamConfig ()
 
std::string GetHname (int inu, int eff, int beam, int det)
 
void FillVector (TH1D *hist, int ntype, int eff, int beam, int det)
 
void FillVector (TH1F *hist, int ntype, int eff, int beam, int det)
 
std::string BeamSysToString (Conventions::BeamSys_t bstype)
 
std::string BeamTypeToString (Conventions::BeamType_t btype)
 
std::string DetTypeToString (Conventions::DetType_t dtype)
 

Private Attributes

std::vector< double > fFPar
 
std::string fFpath
 
std::vector< Conventions::ParticleType_tfPlist
 
std::map< Conventions::ParticleType_t, TH2F * > fPTPZ
 
std::map< Conventions::ParticleType_t, TH2F * > fWeightedPTPZ
 
std::map< Conventions::ParticleType_t, TH2F * > fWeightHist
 
std::map< Conventions::ParticleType_t, double > fMeanPT
 
std::map< Conventions::ParticleType_t, double > fN
 
std::map< Conventions::ParticleType_t, double > fNWeighted
 
std::map< Conventions::ParticleType_t, double > fWeightedMeanPT
 
std::map< Conventions::ParticleType_t, int > fNBinsY
 
std::map< Conventions::ParticleType_t, int > fNBinsX
 
std::vector< double > fBPar
 
std::string fBpath
 
TFile * fBeamSysFile
 
std::map< mapkey, WeightMap_t, LessThanfBeamSysMap
 
int fBflag
 

Detailed Description

Definition at line 23 of file skzpReweight.h.

Member Typedef Documentation

typedef std::map<double, double> nbw::skzpReweight::WeightMap_t
private

Definition at line 119 of file skzpReweight.h.

Constructor & Destructor Documentation

nbw::skzpReweight::skzpReweight ( std::string  fpath = "/nova/data/flux/SKZPdata/fluka05ptxf.root",
std::string  bpath = "/nova/data/flux/SKZPdata/IPNDhists.root",
int  flag = 2 
)

Definition at line 18 of file skzpReweight.cxx.

19  {
20  //default parameters specified by minos-doc-7146
21  fFPar.push_back(1.56);
22  fFPar.push_back(-6.42);
23  fFPar.push_back(1.11);
24  fFPar.push_back(0.13);
25  fFPar.push_back(1.00);
26  fFPar.push_back(1.25);
27  fFPar.push_back(3.50);
28  fFPar.push_back(4.83);
29  fFPar.push_back(1.51);
30  fFPar.push_back(0.29);
31  fFPar.push_back(0.97);
32  fFPar.push_back(2.16);
33  fFPar.push_back(1.04);
34  fFPar.push_back(-0.89);
35  fFPar.push_back(0.88);
36  fFPar.push_back(0.05);
37 
38  fBPar.push_back(-3.85);
39  fBPar.push_back(1.39);
40 
41  fFpath = fpath;
42  fBpath = bpath;
43  fBflag = flag;
44  FlukConfig();
45  if (fBflag>0)
46  BeamConfig();
47  }
std::vector< double > fFPar
Definition: skzpReweight.h:59
std::string fBpath
Definition: skzpReweight.h:116
std::string fFpath
Definition: skzpReweight.h:60
std::vector< double > fBPar
Definition: skzpReweight.h:115
nbw::skzpReweight::~skzpReweight ( )

Definition at line 50 of file skzpReweight.cxx.

51  {
52  //Deconstructs members for Fluk
53  for (std::vector<Conventions::ParticleType_t>::iterator itPlist=fPlist.begin();itPlist!=fPlist.end(); itPlist++)
54  {
55  delete fWeightHist[*itPlist];
56  delete fPTPZ[*itPlist];
57  delete fWeightedPTPZ[*itPlist];
58  }
59  fWeightHist.clear();
60  fPTPZ.clear();
61  fWeightedPTPZ.clear();
62  fMeanPT.clear();
63  fWeightedMeanPT.clear();
64  fN.clear();
65  fNWeighted.clear();
66  fPlist.clear();
67  fFPar.clear();
68  fNBinsX.clear();
69  fNBinsY.clear();
70  }
std::map< Conventions::ParticleType_t, int > fNBinsY
Definition: skzpReweight.h:69
intermediate_table::iterator iterator
std::vector< double > fFPar
Definition: skzpReweight.h:59
std::map< Conventions::ParticleType_t, double > fN
Definition: skzpReweight.h:66
std::map< Conventions::ParticleType_t, int > fNBinsX
Definition: skzpReweight.h:69
std::vector< Conventions::ParticleType_t > fPlist
Definition: skzpReweight.h:61
std::map< Conventions::ParticleType_t, double > fNWeighted
Definition: skzpReweight.h:67
std::map< Conventions::ParticleType_t, TH2F * > fWeightedPTPZ
Definition: skzpReweight.h:63
std::map< Conventions::ParticleType_t, TH2F * > fPTPZ
Definition: skzpReweight.h:62
std::map< Conventions::ParticleType_t, TH2F * > fWeightHist
Definition: skzpReweight.h:64
std::map< Conventions::ParticleType_t, double > fWeightedMeanPT
Definition: skzpReweight.h:68
std::map< Conventions::ParticleType_t, double > fMeanPT
Definition: skzpReweight.h:65

Member Function Documentation

void nbw::skzpReweight::BeamConfig ( )
private

Definition at line 254 of file skzpReweight.cxx.

255  {
256  std::cout <<"BeamConfig reading data from: "<<fBpath<<std::endl;
257 
258  bool FoundHist = false;
259  TDirectory *save = gDirectory;
260  fBeamSysFile = new TFile(fBpath.c_str());
261  if (fBeamSysFile->IsZombie()){
262  std::cout << "Don't recognize path: " << fBpath << std::endl;
263  return;
264  }
265 
266  int ntype[] ={ 56, 55 , 53 , 52};
267  for (int inu = 0; inu < 4; ++inu) {
268  //'End' enums are there so one can just change conventions when needed.
269  for (int eff=1;eff<Conventions::kBeamSysEnd;eff++) {
270  for (int beam=1;beam<Conventions::kBeamEnd;beam++) {
271  for (int det=1;det<Conventions::kDetEnd;det++) {
272  std::string hname = GetHname(inu,eff,beam,det);
273  //I know this is ugly, but it works and I don't know why the dynamic_cast does not work like it is supposed to. I'll fix it later.
274  TH1D* hist=dynamic_cast<TH1D*> (fBeamSysFile->Get(hname.c_str()));
275  if (hist)
276  {
277  FoundHist = true;
279  FillVector(hist,ntype[inu],eff,beam,det);
280  else //for far/near error bar internally use kUnknownDet
281  FillVector(hist,ntype[inu],eff,beam,Conventions::kUnknownDet);
282  }
283  TH1F* hist2=dynamic_cast<TH1F*> (fBeamSysFile->Get(hname.c_str()));
284  if (hist2)
285  {
286  FoundHist = true;
288  FillVector(hist2,ntype[inu],eff,beam,det);
289  else //for far/near error bar internally use kUnknownDet
290  FillVector(hist2,ntype[inu],eff,beam,Conventions::kUnknownDet);
291  }
292  }
293  }
294  }
295  }
296  if (!FoundHist)
297  std::cout<<"No histograms found in " << fBpath << " with name formatted as expected by flag: " << fBflag << std::endl;
298  fBeamSysFile->Close();
299  delete fBeamSysFile;
300  fBeamSysFile = 0;
301  save->cd();
302  return;
303  }
std::string string
Definition: nybbler.cc:12
std::string GetHname(int inu, int eff, int beam, int det)
std::string fBpath
Definition: skzpReweight.h:116
def save(obj, fname)
Definition: root.py:19
void FillVector(TH1D *hist, int ntype, int eff, int beam, int det)
QTextStream & endl(QTextStream &s)
std::string nbw::skzpReweight::BeamSysToString ( Conventions::BeamSys_t  bstype)
private

Definition at line 413 of file skzpReweight.cxx.

414  {
415  if(fBflag>=0 || fBflag<=2)
416  {
417  switch (bstype)
418  {
419  case Conventions::kHornIMiscal: return "HornIMiscal"; break;
420  case Conventions::kHornIDist : return "HornIDist"; break;
421  case Conventions::kUnknownEff : return "Unknown"; break;
422  default : return "Unknown"; break;
423  }
424  }
425  return "Unknown";
426  }
std::string nbw::skzpReweight::BeamTypeToString ( Conventions::BeamType_t  btype)
private

Definition at line 428 of file skzpReweight.cxx.

429  {
430  if(fBflag == 0)
431  {
432  switch (btype)
433  {
434  case Conventions::kLE : return "LE"; break;
435  case Conventions::kLE010z185i : return "LE010z185i"; break;
436  case Conventions::kLE100z200i : return "LE100z200i"; break;
437  case Conventions::kLE250z200i : return "LE250z200i"; break;
438  case Conventions::kLE010z185iL: return "LE010z185iL"; break;
439  case Conventions::kLE010z170i : return "LE010z170i"; break;
440  case Conventions::kLE010z200i : return "LE010z200i"; break;
441  case Conventions::kLE010z000i : return "LE010z000i"; break;
442  case Conventions::kLE150z200i : return "LE150z200i"; break;
443  case Conventions::kUnknownBeam: return "Unknown"; break;
444  default : return "Unknown"; break;
445  }
446  }
447  else if (fBflag == 1 || fBflag == 2)
448  {
449  switch (btype)
450  {
451  case Conventions::kLE : return "L"; break;
452  case Conventions::kLE010z185i : return "L010z185i"; break;
453  case Conventions::kLE100z200i : return "L100z200i"; break;
454  case Conventions::kLE250z200i : return "L250z200i"; break;
455  case Conventions::kLE010z185iL: return "L010z185i_lowint"; break;
456  case Conventions::kLE010z170i : return "L010z170i"; break;
457  case Conventions::kLE010z200i : return "L010z200i"; break;
458  case Conventions::kLE010z000i : return "L010z000i"; break;
459  case Conventions::kLE150z200i : return "L150z200i"; break;
460  case Conventions::kUnknownBeam: return "Unknown"; break;
461  default : return "Unknown"; break;
462  }
463  }
464  return "Unknown";
465  }
std::string nbw::skzpReweight::DetTypeToString ( Conventions::DetType_t  dtype)
private

Definition at line 467 of file skzpReweight.cxx.

468  {
469  if(fBflag == 1)
470  {
471  switch (dtype)
472  {
473  case Conventions::kMINOSnd : return "Near"; break;
474  case Conventions::kMINOSfd : return "Far"; break;
475  case Conventions::kMINOSrat : return "FN"; break;
476  case Conventions::kUnknownDet: return "Unknown"; break;
477  default : return "Unknown"; break;
478  }
479  }
480  else if (fBflag == 0 || fBflag == 2)
481  {
482  switch (dtype)
483  {
484  case Conventions::kNOvAnd : return "NOvAnd"; break;
485  case Conventions::kNOvAfd : return "NOvAfd"; break;
486  case Conventions::kIPND : return "IPND"; break;
487  case Conventions::kMINOSnd : return "MINOSnd"; break;
488  case Conventions::kMINOSfd : return "MINOSfd"; break;
489  case Conventions::kNOvArat : return "NOvArat"; break;
490  case Conventions::kMINOSrat : return "MINOSrat"; break;
491  case Conventions::kUnknownDet: return "Unknown"; break;
492  default : return "Unknown"; break;
493  }
494  }
495  return "Unknown";
496  }
void nbw::skzpReweight::FillVector ( TH1D *  hist,
int  ntype,
int  eff,
int  beam,
int  det 
)
private

Definition at line 323 of file skzpReweight.cxx.

324  {
325  mapkey dex;
326  dex.NuDex=ntype;
327  dex.EffDex=eff;
328  dex.BeamDex=beam;
329  dex.DetDex=det;
330 
332  if (dexit == fBeamSysMap.end())
333  {
334  WeightMap_t wmap;
335  //Endex is the upper-edge for the energy of each bin
336  double EnDex=0;
337  for (int i=1; i<=hist->GetNbinsX();i++)
338  {
339  EnDex+=hist->GetBinWidth(i);
340  wmap.insert(WeightMap_t::value_type(EnDex,hist->GetBinContent(i)));
341  }
342  fBeamSysMap.insert(std::map<mapkey, WeightMap_t, LessThan >::value_type(dex, wmap));
343  }
344  else
345  std::cout << "Already loaded this histogram!" << std::endl;
346 
347  return;
348  }
intermediate_table::iterator iterator
std::map< mapkey, WeightMap_t, LessThan > fBeamSysMap
Definition: skzpReweight.h:121
std::map< double, double > WeightMap_t
Definition: skzpReweight.h:119
QTextStream & endl(QTextStream &s)
void nbw::skzpReweight::FillVector ( TH1F *  hist,
int  ntype,
int  eff,
int  beam,
int  det 
)
private

Definition at line 351 of file skzpReweight.cxx.

352  {
353  mapkey dex;
354  dex.NuDex=ntype;
355  dex.EffDex=eff;
356  dex.BeamDex=beam;
357  dex.DetDex=det;
358 
360  if (dexit == fBeamSysMap.end())
361  {
362  WeightMap_t wmap;
363  //Endex is the upper-edge for the energy of each bin
364  double EnDex=0;
365  for (int i=1; i<hist->GetNbinsX()+1;i++)
366  {
367  EnDex+=hist->GetBinWidth(i);
368  wmap.insert(WeightMap_t::value_type(EnDex,hist->GetBinContent(i)));
369  }
370  fBeamSysMap.insert(std::map<mapkey, WeightMap_t, LessThan >::value_type(dex, wmap));
371  }
372  else
373  std::cout << "Already loaded this histogram!" << std::endl;
374 
375  return;
376  }
intermediate_table::iterator iterator
std::map< mapkey, WeightMap_t, LessThan > fBeamSysMap
Definition: skzpReweight.h:121
std::map< double, double > WeightMap_t
Definition: skzpReweight.h:119
QTextStream & endl(QTextStream &s)
void nbw::skzpReweight::FlukConfig ( )
private

Definition at line 146 of file skzpReweight.cxx.

147  {
148  std::cout <<"FlukConfig reading data from: "<<fFpath<<std::endl;
149 
150  //Builds maps and histograms for Fluk from file
151  //In fluka05ptxf.root, file and histogram names refer to xf, but they actually mean pz, which is proportional to xF:The conversion takes place in this code and not the histogram.
152 
153  TFile* hFile=new TFile(fFpath.c_str());
154 
155  fPlist.push_back(Conventions::kPiPlus);
156  fPlist.push_back(Conventions::kPiMinus);
157  fPlist.push_back(Conventions::kKPlus);
158  fPlist.push_back(Conventions::kKMinus);
159  fPlist.push_back(Conventions::kK0L);
160 
161  for (std::vector<Conventions::ParticleType_t>::iterator itPlist=fPlist.begin();itPlist!=fPlist.end(); itPlist++)
162  {
163  TH2F* hist=dynamic_cast <TH2F*> (hFile->Get(Form("hF05ptxf%s",PartEnumToString(*itPlist).c_str()))->Clone());
164  hist->SetDirectory(0);
165  TH2F* hist2=dynamic_cast <TH2F*> (hist->Clone(Form("hWeightedPTXF%s",PartEnumToString(*itPlist).c_str())));
166  hist2->SetDirectory(0);
167  hist2->SetTitle(Form("%s weighted pt-pz",PartEnumToString(*itPlist).c_str()));
168  fWeightHist[*itPlist]=dynamic_cast<TH2F*> (hist->Clone(Form("hWeight%s",PartEnumToString(*itPlist).c_str())));
169 
170  fWeightHist[*itPlist]->Divide(hist2);
171  fWeightHist[*itPlist]->SetDirectory(0);
172 
173  std::pair<Conventions::ParticleType_t, TH2F* > p(*itPlist, hist);
174  std::pair<Conventions::ParticleType_t, TH2F* > p2(*itPlist, hist2);
175  fPTPZ.insert(p);
176  fWeightedPTPZ.insert(p2);
177 
178  fPTPZ[*itPlist]->SetDirectory(0);
179  fWeightedPTPZ[*itPlist]->SetDirectory(0);
180 
181  fMeanPT[*itPlist]=fPTPZ[*itPlist]->ProjectionY()->GetMean()*1000.;
182  fWeightedMeanPT[*itPlist]=fWeightedPTPZ[*itPlist]->ProjectionY()->GetMean()*1000.;
183 
184  double N=fPTPZ[*itPlist]->ProjectionY()->GetSumOfWeights();
185  double wN=fWeightedPTPZ[*itPlist]->ProjectionY()->GetSumOfWeights();
186 
187  fN[*itPlist]=N;
188  fNWeighted[*itPlist]=wN;
189  }
190  hFile->Close();
191  delete hFile;
192  hFile=0;
193  //Updates maps and histograms for Fluk with GetFlukWeight (the weight for KOL depends on fN and fNWeighted for KPlus and KMinus).
194  for (std::vector<Conventions::ParticleType_t>::iterator itPlist=fPlist.begin(); itPlist!=fPlist.end(); itPlist++)
195  {
196  double pt,pz;
197  double meanpt(0.);
198  double sum(0.);
199  for (int i=0;i<fPTPZ[*itPlist]->GetNbinsX()+1;i++)
200  {
201  for (int j=0;j<fPTPZ[*itPlist]->GetNbinsY()+1;j++)
202  {
203  pz=fPTPZ[*itPlist]->GetXaxis()->GetBinCenter(i);
204  pt=fPTPZ[*itPlist]->GetYaxis()->GetBinCenter(j);
205  fWeightedPTPZ[*itPlist]->SetBinContent(i,j,(fPTPZ[*itPlist]->GetBinContent(i,j)*GetFlukWeight(*itPlist,pt,pz)));
206  fWeightHist[*itPlist]->SetBinContent(i,j,GetFlukWeight(*itPlist,pt,pz));
207  meanpt+=fWeightedPTPZ[*itPlist]->GetBinContent(i,j)*pt;
208  sum+=fWeightedPTPZ[*itPlist]->GetBinContent(i,j);
209  }
210  }
211  meanpt/=sum;
212  fWeightedMeanPT[*itPlist]=meanpt*1000.; //GeV to MeV
213  fNWeighted[*itPlist]=sum;
214  meanpt=0.;
215  sum=0.;
216  }
217  return;
218  }
intermediate_table::iterator iterator
std::map< Conventions::ParticleType_t, double > fN
Definition: skzpReweight.h:66
std::vector< Conventions::ParticleType_t > fPlist
Definition: skzpReweight.h:61
std::map< Conventions::ParticleType_t, double > fNWeighted
Definition: skzpReweight.h:67
std::string fFpath
Definition: skzpReweight.h:60
std::map< Conventions::ParticleType_t, TH2F * > fWeightedPTPZ
Definition: skzpReweight.h:63
std::map< Conventions::ParticleType_t, TH2F * > fPTPZ
Definition: skzpReweight.h:62
std::map< Conventions::ParticleType_t, TH2F * > fWeightHist
Definition: skzpReweight.h:64
std::map< Conventions::ParticleType_t, double > fWeightedMeanPT
Definition: skzpReweight.h:68
p
Definition: test.py:223
std::string PartEnumToString(Conventions::ParticleType_t ptype)
double GetFlukWeight(const simb::MCFlux *mcf)
Definition: skzpReweight.h:45
std::map< Conventions::ParticleType_t, double > fMeanPT
Definition: skzpReweight.h:65
QTextStream & endl(QTextStream &s)
Conventions::ParticleType_t nbw::skzpReweight::GeantToEnum ( int  ptype)
private

Definition at line 379 of file skzpReweight.cxx.

380  {
381  switch (ptype)
382  {
383  case 8 : return Conventions::kPiPlus; break;
384  case 211 : return Conventions::kPiPlus; break;
385  case 9 : return Conventions::kPiMinus; break;
386  case -211: return Conventions::kPiMinus; break;
387  case 11 : return Conventions::kKPlus; break;
388  case 321 : return Conventions::kKPlus; break;
389  case 12 : return Conventions::kKMinus; break;
390  case -321: return Conventions::kKMinus; break;
391  case 10 : return Conventions::kK0L; break;
392  case 130 : return Conventions::kK0L; break;
393  default : return Conventions::kUnknown; break;
394  }
395  return Conventions::kUnknown;
396  }
double nbw::skzpReweight::GetBeamWeight ( int  ntype,
double  Enu,
int  det = 1,
int  beam = 2 
)

Definition at line 221 of file skzpReweight.cxx.

222  {
223  double weight=1.;
224 
225  if (ntype == 14) ntype = 56;
226  else if (ntype == -14) ntype = 55;
227  else if (ntype == 12) ntype = 53;
228  else if (ntype == -12) ntype = 52;
229 
230  struct mapkey dex;
231  dex.NuDex=ntype;
232  dex.BeamDex=beam;
233  dex.DetDex=det;
234 
235  for(int sys = 1; sys <= 2; ++sys){//Loops through each beam-focusing correction applied
236  double w = 0.;
237  dex.EffDex=sys;
238 
240  if (dexit != fBeamSysMap.end()){
241  for (WeightMap_t::iterator EnDex = (dexit->second).begin(); EnDex!=(dexit->second).end(); EnDex++){
242  if(EnDex->first > Enu){//the first EnDex that goes over is mapped to the weight value
243  w = EnDex->second;
244  break;
245  }
246  }
247  }
248  weight *= std::abs(w)*fBPar[sys-1]+1.;
249  }
250  return weight;
251  }
intermediate_table::iterator iterator
decltype(auto) constexpr end(T &&obj)
ADL-aware version of std::end.
Definition: StdUtils.h:77
T abs(T value)
std::map< mapkey, WeightMap_t, LessThan > fBeamSysMap
Definition: skzpReweight.h:121
std::vector< double > fBPar
Definition: skzpReweight.h:115
weight
Definition: test.py:257
double nbw::skzpReweight::GetFlukWeight ( const simb::MCFlux mcf)
inline

Definition at line 45 of file skzpReweight.h.

45  {
46  double pt = sqrt(mcf->ftpx*mcf->ftpx + mcf->ftpy*mcf->ftpy);
47  return GetFlukWeight(mcf->ftptype,pt,mcf->ftpz); };
double ftpx
Definition: MCFlux.h:79
int ftptype
Definition: MCFlux.h:82
double ftpz
Definition: MCFlux.h:81
double GetFlukWeight(const simb::MCFlux *mcf)
Definition: skzpReweight.h:45
double ftpy
Definition: MCFlux.h:80
double nbw::skzpReweight::GetFlukWeight ( int  ptype,
double  pT,
double  pz 
)

Definition at line 73 of file skzpReweight.cxx.

74  {
75  double weight=1.;
76  double xF=pz/120.;
77  double A,B,C;
78  double Ap,Bp,Cp;
80 
81  // This is SKZP parameterization
82  if (xF>1.||xF<0.) return weight;
83  if (pT>1.||pT<0.) return weight;
84 
85  if (eptype==Conventions::kPiPlus || eptype==Conventions::kPiMinus)
86  {
87  //Calculate weight for pions
88  if (pT<0.03) // for low pT
89  pT=0.03;
90  // calculate the A, B and C in SKZP parameterization
91  // A, B and C are best fit to Fluka 05
92 
93  A=-0.00761*pow((1.-xF),4.045)*(1.+9620.*xF)*pow(xF,-2.975);
94  B=0.05465*pow((1.-xF),2.675)*(1.+69590.*xF)*pow(xF,-3.144);
95 
96  if (xF<0.22)
97  C=-7.058/pow(xF,-0.1419)+9.188;
98  else
99  C = (3.008/exp((xF-0.1984)*3.577)) + 2.616*xF+0.1225;
100 
101  // scale/skew A, B and C
102  Ap=(fFPar[0]+fFPar[1]*xF)*A;
103  Bp=(fFPar[2]+fFPar[3]*xF)*B;
104  Cp=(fFPar[4]+fFPar[5]*xF)*C;
105 
106  weight=(Ap+Bp*pT)/(A+B*pT)*exp(-(Cp-C)*pow(pT,3./2.));
107  if(eptype == Conventions::kPiMinus && pz > 4)
108  weight *= ( fFPar[12] + fFPar[13] * xF );
109  }
110  else if (eptype==Conventions::kKPlus || eptype==Conventions::kKMinus)
111  {
112  //Calculate weight for kaons
113  if (pT<0.05) // for low pT
114  pT=0.05;
115  // calculate the A, B and C in SKZP parameterization
116  // A, B and C are best fit to Fluka 05
117 
118  A=-0.005187*pow((1.-xF),4.119)*(1.+2170.*xF)*pow(xF,-2.767);
119  B=0.4918*pow((1.-xF),2.672)*(1.+1373.*xF)*pow(xF,-2.927);
120 
121  if (xF<0.22)
122  C=-16.10/pow(xF,-0.04582)+17.92;
123  else
124  C = (6.905/exp((xF+0.163)*6.718)) - 0.4257*xF + 2.486;
125 
126  // scale/skew A, B and C
127  Ap=(fFPar[6]+fFPar[7]*xF)*A;
128  Bp=(fFPar[8]+fFPar[9]*xF)*B;
129  Cp=(fFPar[10]+fFPar[11]*xF)*C;
130 
131  weight=(Ap+Bp*pT)/(A+B*pT)*exp(-(Cp-C)*pow(pT,3./2.));
132  if(eptype==Conventions::kKMinus)
133  weight*=(fFPar[14]+fFPar[15]*xF);
134  }
135  else if (eptype==Conventions::kK0L)
136  {
137  //N(K0L) approximately given with (N(K+)+3*N(K-))/4
139  }
140 
141  if (weight>10.) weight=10.;
142  return weight;
143  }
std::vector< double > fFPar
Definition: skzpReweight.h:59
std::map< Conventions::ParticleType_t, double > fN
Definition: skzpReweight.h:66
constexpr T pow(T x)
Definition: pow.h:72
std::map< Conventions::ParticleType_t, double > fNWeighted
Definition: skzpReweight.h:67
enum nbw::Conventions::EParticleType ParticleType_t
#define A
Definition: memgrp.cpp:38
weight
Definition: test.py:257
Conventions::ParticleType_t GeantToEnum(int ptype)
std::string nbw::skzpReweight::GetHname ( int  inu,
int  eff,
int  beam,
int  det 
)
private

Definition at line 306 of file skzpReweight.cxx.

307  {
308  std::string hname;
309  if (fBflag == 1 || fBflag == 2)
310  {
311  std::string nus[]={"NuMu","NuMuBar","NuE","NuEBar"};
312  hname=nus[inu]+"_";
313  }
314  hname=hname + BeamSysToString(Conventions::BeamSys_t(eff));
315  if (fBflag == 1 || fBflag == 2) hname = hname + "_";
317  if (fBflag == 1 || fBflag == 2) hname = hname + "_";
318  hname=hname + DetTypeToString(Conventions::DetType_t(det));
319  return hname;
320  }
enum nbw::Conventions::EDetType DetType_t
std::string string
Definition: nybbler.cc:12
std::string BeamTypeToString(Conventions::BeamType_t btype)
enum nbw::Conventions::EBeamType BeamType_t
enum nbw::Conventions::EBeamSys BeamSys_t
std::string BeamSysToString(Conventions::BeamSys_t bstype)
std::string DetTypeToString(Conventions::DetType_t dtype)
double nbw::skzpReweight::GetWeight ( const simb::MCFlux mcf,
double  Enu,
int  det,
int  beam 
)
inline

Definition at line 41 of file skzpReweight.h.

41  {
42  double pt = sqrt(mcf->ftpx*mcf->ftpx + mcf->ftpy*mcf->ftpy);
43  return GetFlukWeight(mcf->ftptype,pt,mcf->ftpz)*GetBeamWeight(mcf->fntype,Enu,beam,det); };
double ftpx
Definition: MCFlux.h:79
int ftptype
Definition: MCFlux.h:82
double GetBeamWeight(int ntype, double Enu, int det=1, int beam=2)
double ftpz
Definition: MCFlux.h:81
int fntype
Definition: MCFlux.h:51
double GetFlukWeight(const simb::MCFlux *mcf)
Definition: skzpReweight.h:45
double ftpy
Definition: MCFlux.h:80
std::string nbw::skzpReweight::PartEnumToString ( Conventions::ParticleType_t  ptype)
private

Definition at line 398 of file skzpReweight.cxx.

399  {
400  switch (ptype)
401  {
402  case Conventions::kPiPlus : return "PiPlus"; break;
403  case Conventions::kPiMinus: return "PiMinus"; break;
404  case Conventions::kKPlus : return "KPlus"; break;
405  case Conventions::kKMinus : return "KMinus"; break;
406  case Conventions::kK0L : return "K0L"; break;
407  case Conventions::kUnknown: return "Unknown"; break;
408  default : return "Unknown"; break;
409  }
410  return "Unknown";
411  }
void nbw::skzpReweight::SetBeamParams ( std::vector< double >  bpar)
inline

Definition at line 37 of file skzpReweight.h.

37  {
38  fBPar=bpar; BeamConfig();
39  return; };
std::vector< double > fBPar
Definition: skzpReweight.h:115
void nbw::skzpReweight::SetFlukParams ( std::vector< double >  fpar)
inline

Definition at line 34 of file skzpReweight.h.

34  {
35  fFPar=fpar; FlukConfig();
36  return; };
std::vector< double > fFPar
Definition: skzpReweight.h:59
void nbw::skzpReweight::SetParams ( std::vector< double >  fpar,
std::vector< double >  bpar 
)
inline

Definition at line 30 of file skzpReweight.h.

30  {
31  fFPar=fpar; FlukConfig();
32  fBPar=bpar; BeamConfig();
33  return; };
std::vector< double > fFPar
Definition: skzpReweight.h:59
std::vector< double > fBPar
Definition: skzpReweight.h:115

Member Data Documentation

TFile* nbw::skzpReweight::fBeamSysFile
private

Definition at line 117 of file skzpReweight.h.

std::map<mapkey, WeightMap_t, LessThan > nbw::skzpReweight::fBeamSysMap
private

Definition at line 121 of file skzpReweight.h.

int nbw::skzpReweight::fBflag
private

Definition at line 126 of file skzpReweight.h.

std::vector<double> nbw::skzpReweight::fBPar
private

Definition at line 115 of file skzpReweight.h.

std::string nbw::skzpReweight::fBpath
private

Definition at line 116 of file skzpReweight.h.

std::vector<double> nbw::skzpReweight::fFPar
private

Definition at line 59 of file skzpReweight.h.

std::string nbw::skzpReweight::fFpath
private

Definition at line 60 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, double > nbw::skzpReweight::fMeanPT
private

Definition at line 65 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, double > nbw::skzpReweight::fN
private

Definition at line 66 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, int > nbw::skzpReweight::fNBinsX
private

Definition at line 69 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, int > nbw::skzpReweight::fNBinsY
private

Definition at line 69 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, double > nbw::skzpReweight::fNWeighted
private

Definition at line 67 of file skzpReweight.h.

std::vector<Conventions::ParticleType_t> nbw::skzpReweight::fPlist
private

Definition at line 61 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, TH2F* > nbw::skzpReweight::fPTPZ
private

Definition at line 62 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, double > nbw::skzpReweight::fWeightedMeanPT
private

Definition at line 68 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, TH2F* > nbw::skzpReweight::fWeightedPTPZ
private

Definition at line 63 of file skzpReweight.h.

std::map<Conventions::ParticleType_t, TH2F* > nbw::skzpReweight::fWeightHist
private

Definition at line 64 of file skzpReweight.h.


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