GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
HEDISStrucFunc.h
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*!
3 
4 \class genie::HEDISStrcuFunc
5 
6 \brief Singleton class to load Structure Functions used in HEDIS.
7 
8 \author Alfonso Garcia <alfonsog \at nikhef.nl>
9  NIKHEF
10 
11 \created August 28, 2019
12 
13 \cpright Copyright (c) 2003-2024, The GENIE Collaboration
14  For the full text of the license visit http://copyright.genie-mc.org
15  or see $GENIE/LICENSE
16 */
17 //____________________________________________________________________________
18 
19 #ifndef _HEDIS_STRUC_FUNC_H_
20 #define _HEDIS_STRUC_FUNC_H_
21 
24 
25 #include <map>
26 #include <vector>
27 #include <string>
28 #include <iostream>
29 #include <fstream>
30 
31 
32 using std::map;
33 using std::vector;
34 using std::string;
35 using std::to_string;
36 
37 namespace genie {
38 
39  struct SF_info {
40 
41  string LHAPDFset;
43  bool IsNLO;
44  string Scheme;
45  int QrkThrs;
46  int NGridX;
47  int NGridQ2;
48  double XGridMin;
50  double MassW, MassZ;
51  double Rho, Sin2ThW;
52  double Vud,Vus,Vub;
53  double Vcd,Vcs,Vcb;
54  double Vtd,Vts,Vtb;
55 
56  };
57 
58  inline bool operator== (const SF_info& a, const SF_info& b) {
59 
60  if ( a.LHAPDFmember != b.LHAPDFmember ) return false;
61  if ( a.LHAPDFset != b.LHAPDFset ) return false;
62  if ( a.IsNLO != b.IsNLO ) return false;
63  if ( a.Scheme != b.Scheme ) return false;
64  if ( a.QrkThrs != b.QrkThrs ) return false;
65  if ( a.NGridX != b.NGridX ) return false;
66  if ( a.NGridQ2 != b.NGridQ2 ) return false;
67  if ( abs(a.XGridMin-b.XGridMin)>1e-10 ) return false;
68  if ( abs(a.Q2GridMin-b.Q2GridMin)>1e-10 ) return false;
69  if ( abs(a.Q2GridMax-b.Q2GridMax)>1e-10 ) return false;
70  if ( abs(a.MassW-b.MassW)>1e-10 ) return false;
71  if ( abs(a.MassZ-b.MassZ)>1e-10 ) return false;
72  if ( abs(a.Rho-b.Rho)>1e-10 ) return false;
73  if ( abs(a.Sin2ThW-b.Sin2ThW)>1e-10 ) return false;
74  if ( abs(a.Vud-b.Vud)>1e-10 ) return false;
75  if ( abs(a.Vus-b.Vus)>1e-10 ) return false;
76  if ( abs(a.Vub-b.Vub)>1e-10 ) return false;
77  if ( abs(a.Vcd-b.Vcd)>1e-10 ) return false;
78  if ( abs(a.Vcs-b.Vcs)>1e-10 ) return false;
79  if ( abs(a.Vcb-b.Vcb)>1e-10 ) return false;
80  if ( abs(a.Vtd-b.Vtd)>1e-10 ) return false;
81  if ( abs(a.Vts-b.Vts)>1e-10 ) return false;
82  if ( abs(a.Vtb-b.Vtb)>1e-10 ) return false;
83  return true;
84 
85  }
86 
87  inline std::istream & operator>> (std::istream& is, SF_info& a) {
88 
89  string saux;
90  std::getline (is,saux); //# Header
91  std::getline (is,saux); //# Header
92  std::getline (is,saux); //# Header
93  std::getline (is,saux); //# Header
94  std::getline (is,saux); //# LHAPDF set
95  std::getline (is,saux); a.LHAPDFset=saux.c_str();
96  std::getline (is,saux); //# LHAPDF member
97  std::getline (is,saux); a.LHAPDFmember=atoi(saux.c_str());
98  std::getline (is,saux); //#NLO
99  std::getline (is,saux); a.IsNLO=atoi(saux.c_str());
100  std::getline (is,saux); //# Mass scheme
101  std::getline (is,saux); a.Scheme=saux.c_str();
102  std::getline (is,saux); //# Quark threshold
103  std::getline (is,saux); a.QrkThrs=atof(saux.c_str());
104  std::getline (is,saux); //# NGridX
105  std::getline (is,saux); a.NGridX=atoi(saux.c_str());
106  std::getline (is,saux); //# NGridQ2
107  std::getline (is,saux); a.NGridQ2=atoi(saux.c_str());
108  std::getline (is,saux); //# XGridMin
109  std::getline (is,saux); a.XGridMin=atof(saux.c_str());
110  std::getline (is,saux); //# Q2min
111  std::getline (is,saux); a.Q2GridMin=atof(saux.c_str());
112  std::getline (is,saux); //# Q2max
113  std::getline (is,saux); a.Q2GridMax=atof(saux.c_str());
114  std::getline (is,saux); //# Mass W
115  std::getline (is,saux); a.MassW=atof(saux.c_str());
116  std::getline (is,saux); //# Mass Z
117  std::getline (is,saux); a.MassZ=atof(saux.c_str());
118  std::getline (is,saux); //# Rho
119  std::getline (is,saux); a.Rho=atof(saux.c_str());
120  std::getline (is,saux); //# Sin2ThW
121  std::getline (is,saux); a.Sin2ThW=atof(saux.c_str());
122  std::getline (is,saux); //# CKM
123  std::getline (is,saux); a.Vud=atof(saux.c_str());
124  std::getline (is,saux); a.Vus=atof(saux.c_str());
125  std::getline (is,saux); a.Vub=atof(saux.c_str());
126  std::getline (is,saux); a.Vcd=atof(saux.c_str());
127  std::getline (is,saux); a.Vcs=atof(saux.c_str());
128  std::getline (is,saux); a.Vcb=atof(saux.c_str());
129  std::getline (is,saux); a.Vtd=atof(saux.c_str());
130  std::getline (is,saux); a.Vts=atof(saux.c_str());
131  std::getline (is,saux); a.Vtb=atof(saux.c_str());
132  return is;
133 
134  }
135 
136  inline std::ostream & operator<< (std::ostream& os, const SF_info& a) {
137 
138  return os << '\n'
139  << "#--------------------------------------------------------------------------------" << '\n'
140  << "# Metafile that stores information used to generate Structure Functions for HEDIS" << '\n'
141  << "#--------------------------------------------------------------------------------" << '\n'
142  << "# LHAPDF set" << '\n'
143  << a.LHAPDFset << '\n'
144  << "# LHAPDF member" << '\n'
145  << a.LHAPDFmember << '\n'
146  << "# NLO" << '\n'
147  << a.IsNLO << '\n'
148  << "# Mass Scheme" << '\n'
149  << a.Scheme << '\n'
150  << "# Quark threshold" << '\n'
151  << a.QrkThrs << '\n'
152  << "# NX" << '\n'
153  << a.NGridX << '\n'
154  << "# NQ2" << '\n'
155  << a.NGridQ2 << '\n'
156  << "# Xmin" << '\n'
157  << a.XGridMin << '\n'
158  << "# Q2min" << '\n'
159  << a.Q2GridMin << '\n'
160  << "# Q2max" << '\n'
161  << a.Q2GridMax << '\n'
162  << "# Mass W" << '\n'
163  << a.MassW << '\n'
164  << "# Mass Z" << '\n'
165  << a.MassZ << '\n'
166  << "# Rho" << '\n'
167  << a.Rho << '\n'
168  << "# Sin2ThW" << '\n'
169  << a.Sin2ThW << '\n'
170  << "# CKM" << '\n'
171  << a.Vud << '\n'
172  << a.Vus << '\n'
173  << a.Vub << '\n'
174  << a.Vcd << '\n'
175  << a.Vcs << '\n'
176  << a.Vcb << '\n'
177  << a.Vtd << '\n'
178  << a.Vts << '\n'
179  << a.Vtb << '\n';
180 
181  }
182 
183  struct SF_xQ2 {
184 
185  double F1;
186  double F2;
187  double F3;
188 
189  };
190 
192  {
193  public:
194 
195  // ................................................................
196  // HEDIS structure functions type
197  //
198 
199  typedef enum StrucFuncType {
206 
207  // ................................................................
208  // HEDIS form factor type
209  //
210 
212  {
213  public:
215  ~HEDISStrucFuncTable() { /* note: should delete the grids! */ }
216  map< HEDISStrucFunc::HEDISStrucFuncType_t, genie::BLI2DNonUnifGrid * > Table;
217  };
218 
219  // ................................................................
220 
221  static HEDISStrucFunc * Instance(SF_info sfinfo);
222 
223  // method to return values of the SF for a particular channel in x and Q2
224  SF_xQ2 EvalQrkSFLO ( const Interaction * in, double x, double Q2 );
225  SF_xQ2 EvalNucSFLO ( const Interaction * in, double x, double Q2 );
226  SF_xQ2 EvalNucSFNLO ( const Interaction * in, double x, double Q2 );
227 
228  private:
229 
230  // Ctors & dtor
231  HEDISStrucFunc(SF_info sfinfo);
233  ~HEDISStrucFunc();
234 
235  void CreateQrkSF ( const Interaction * in, string sfFile );
236  void CreateNucSF ( const Interaction * in, string sfFile );
237 
238  string QrkSFName ( const Interaction * in );
239  string NucSFName ( const Interaction * in ) ;
240  int QrkSFCode ( const Interaction * in );
241  int NucSFCode ( const Interaction * in ) ;
242 
243  // Self
245 
246  // These map holds all SF tables (interaction channel is the key)
247  map<int, HEDISStrucFuncTable> fQrkSFLOTables;
248  map<int, HEDISStrucFuncTable> fNucSFLOTables;
249  map<int, HEDISStrucFuncTable> fNucSFNLOTables;
250 
252  vector<double> sf_x_array;
253  vector<double> sf_q2_array;
254 
255  // singleton cleaner
256  struct Cleaner {
262  }
263  }
264  };
265  friend struct Cleaner;
266  };
267 
268 } // genie namespace
269 
270 #endif // _HEDIS_STRUC_FUNC_H_
static HEDISStrucFunc * fgInstance
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1077
SF_xQ2 EvalNucSFLO(const Interaction *in, double x, double Q2)
HEDISStrucFunc(SF_info sfinfo)
vector< double > sf_q2_array
map< int, HEDISStrucFuncTable > fQrkSFLOTables
int NucSFCode(const Interaction *in)
int QrkSFCode(const Interaction *in)
map< HEDISStrucFunc::HEDISStrucFuncType_t, genie::BLI2DNonUnifGrid * > Table
string QrkSFName(const Interaction *in)
enum genie::HEDISStrucFunc::StrucFuncType HEDISStrucFuncType_t
static constexpr double b
Definition: Units.h:78
static HEDISStrucFunc * Instance(SF_info sfinfo)
Summary information for an interaction.
Definition: Interaction.h:56
const double e
void CreateNucSF(const Interaction *in, string sfFile)
const double a
map< int, HEDISStrucFuncTable > fNucSFLOTables
std::istream & operator>>(std::istream &is, SF_info &a)
bool operator==(const TuneId &id1, const TuneId &id2)
Definition: TuneId.cxx:41
map< int, HEDISStrucFuncTable > fNucSFNLOTables
ostream & operator<<(ostream &stream, const AlgConfigPool &config_pool)
void CreateQrkSF(const Interaction *in, string sfFile)
string NucSFName(const Interaction *in)
SF_xQ2 EvalQrkSFLO(const Interaction *in, double x, double Q2)
vector< double > sf_x_array
SF_xQ2 EvalNucSFNLO(const Interaction *in, double x, double Q2)