gMakeSplines.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*!
3 
4 \program gmkspl
5 
6 \brief GENIE utility program building XML cross section splines that can
7  be loaded into GENIE to speed-up event generation.
8  The list of neutrino PDG codes is passed from the command line.
9  The list of nuclear target PDG codes is either passed from the
10  command line or extracted from the input ROOT/GEANT geometry.
11 
12  Syntax :
13  gmkspl -p nupdg
14  <-t target_pdg_codes,
15  -f geometry_file>
16  <-o | --output-cross-sections> output_xml_xsec_file
17  [-n nknots]
18  [-e max_energy]
19  [--no-copy]
20  [--seed random_number_seed]
21  [--input-cross-sections xml_file]
22  [--event-generator-list list_name]
23  [--tune genie_tune]
24  [--message-thresholds xml_file]
25  [--xml-path config_xml_dir]
26 
27  Note :
28  [] marks optional arguments.
29  <> marks a list of arguments out of which only one can be
30  selected at any given time.
31 
32  Options :
33  -p
34  A comma separated list of nu PDG codes.
35  -t
36  A comma separated list of tgt PDG codes.
37  PDG code format: 10LZZZAAAI
38  -f
39  A ROOT file containing a ROOT/GEANT geometry description.
40  -o, --output-cross-sections
41  Name of output XML file containing computed cross-section data.
42  Default: `xsec_splines.xml'.
43  -n
44  Number of knots per spline.
45  Default: 15 knots per decade of energy range with a minimum
46  of 30 knots totally.
47  -e
48  Maximum energy in spline.
49  Default: The max energy in the validity range of the spline
50  generating thread.
51  --no-copy
52  Does not write out the input cross-sections in the output file
53  --seed
54  Random number seed.
55  --input-cross-sections
56  Name (incl. full path) of an XML file with pre-computed
57  free-nucleon cross-section values. If loaded, it can speed-up
58  cross-section calculation for nuclear targets.
59  --event-generator-list
60  List of event generators to load in event generation drivers.
61  [default: "Default"].
62  --tune
63  Specifies a GENIE comprehensive neutrino interaction model tune.
64  [default: "Default"].
65  --message-thresholds
66  Allows users to customize the message stream thresholds.
67  The thresholds are specified using an XML file.
68  See $GENIE/config/Messenger.xml for the XML schema.
69  --xml-path
70  A directory to load XML files from - overrides $GXMLPATH, and $GENIE/config
71 
72  *** See the User Manual for more details and examples. ***
73 
74 
75 \author Costas Andreopoulos <costas.andreopoulos \at stfc.ac.uk>
76  University of Liverpool & STFC Rutherford Appleton Lab
77 
78 \created September 27, 2005
79 
80 \cpright Copyright (c) 2003-2019, The GENIE Collaboration
81  For the full text of the license visit http://copyright.genie-mc.org
82  or see $GENIE/LICENSE
83 */
84 //____________________________________________________________________________
85 
86 #include <cassert>
87 #include <cstdlib>
88 #include <string>
89 #include <vector>
90 
91 #if defined(HAVE_FENV_H) && defined(HAVE_FEENABLEEXCEPT)
92 #include <fenv.h> // for `feenableexcept`
93 #endif
94 
95 #include <TSystem.h>
96 
103 #include "Framework/Utils/RunOpt.h"
104 #include "Framework/Utils/AppInit.h"
106 //#include "Framework/Utils/SystemUtils.h"
110 
111 #ifdef __GENIE_GEOM_DRIVERS_ENABLED__
113 #endif
114 
115 using std::string;
116 using std::vector;
117 
118 using namespace genie;
119 
120 #ifdef __GENIE_GEOM_DRIVERS_ENABLED__
121 using namespace genie::geometry;
122 #endif
123 
124 // Prototypes:
125 void GetCommandLineArgs (int argc, char ** argv);
126 void PrintSyntax (void);
128 PDGCodeList * GetTargetCodes (void);
129 
130 // User-specified options:
131 string gOptNuPdgCodeList = "";
132 string gOptTgtPdgCodeList = "";
133 string gOptGeomFilename = "";
134 int gOptNKnots = -1;
135 double gOptMaxE = -1.;
136 bool gOptNoCopy = false;
137 long int gOptRanSeed = -1; // random number seed
138 string gOptInpXSecFile = ""; // input cross-section file
139 string gOptOutXSecFile = ""; // output cross-section file
140 
141 //____________________________________________________________________________
142 int main(int argc, char ** argv)
143 {
144  // Parse command line arguments
145  GetCommandLineArgs(argc,argv);
146 
147  if ( ! RunOpt::Instance()->Tune() ) {
148  LOG("gmkspl", pFATAL) << " No TuneId in RunOption";
149  exit(-1);
150  }
152 
153  // throw on NaNs and Infs...
154 #if defined(HAVE_FENV_H) && defined(HAVE_FEENABLEEXCEPT)
155  feenableexcept(FE_DIVBYZERO | FE_INVALID | FE_OVERFLOW);
156 #endif
157 
158  // Init
159  utils::app_init::MesgThresholds(RunOpt::Instance()->MesgThresholdFiles());
162 
163  // Get list of neutrinos and nuclear targets
164 
165  PDGCodeList * neutrinos = GetNeutrinoCodes();
166  PDGCodeList * targets = GetTargetCodes();
167 
168  if(!neutrinos || neutrinos->size() == 0 ) {
169  LOG("gmkspl", pFATAL) << "Empty neutrino PDG code list";
170  PrintSyntax();
171  exit(2);
172  }
173  if(!targets || targets->size() == 0 ) {
174  LOG("gmkspl", pFATAL) << "Empty target PDG code list";
175  PrintSyntax();
176  exit(3);
177  }
178 
179  LOG("gmkspl", pINFO) << "Neutrinos: " << *neutrinos;
180  LOG("gmkspl", pINFO) << "Targets: " << *targets;
181 
182  // Loop over all possible input init states and ask the GEVGDriver
183  // to build splines for all the interactions that its loaded list
184  // of event generators can generate.
185 
186  PDGCodeList::const_iterator nuiter;
187  PDGCodeList::const_iterator tgtiter;
188  for(nuiter = neutrinos->begin(); nuiter != neutrinos->end(); ++nuiter) {
189  for(tgtiter = targets->begin(); tgtiter != targets->end(); ++tgtiter) {
190  int nupdgc = *nuiter;
191  int tgtpdgc = *tgtiter;
192  InitialState init_state(tgtpdgc, nupdgc);
193  GEVGDriver driver;
195  driver.Configure(init_state);
197  }
198  }
199 
200  // Save the splines at the requested XML file
202  bool save_init = !gOptNoCopy;
203  xspl->SaveAsXml(gOptOutXSecFile, save_init);
204 
205  delete neutrinos;
206  delete targets;
207 
208  return 0;
209 }
210 //____________________________________________________________________________
211 void GetCommandLineArgs(int argc, char ** argv)
212 {
213  LOG("gmkspl", pINFO) << "Parsing command line arguments";
214 
215  // Common run options. Set defaults and read.
218 
219  // Parse run options for this app
220 
221  CmdLnArgParser parser(argc,argv);
222 
223  // output XML file name
224  if( parser.OptionExists('o') ||
225  parser.OptionExists("output-cross-sections") )
226  {
227  LOG("gmkspl", pINFO) << "Reading output filename";
228  if( parser.OptionExists('o') ) {
229  gOptOutXSecFile = parser.ArgAsString('o');
230  }
231  else {
232  gOptOutXSecFile = parser.ArgAsString("output-cross-sections");
233  }
234  } else {
235  LOG("gmkspl", pINFO) << "Unspecified filename - Using default";
236  gOptOutXSecFile = "xsec_splines.xml";
237  }
238 
239  // number of knots
240  if( parser.OptionExists('n') ) {
241  LOG("gmkspl", pINFO) << "Reading number of knots/spline";
242  gOptNKnots = parser.ArgAsInt('n');
243  } else {
244  LOG("gmkspl", pINFO)
245  << "Unspecified number of knots - Using default";
246  gOptNKnots = -1;
247  }
248 
249  // max spline energy (if < max of validity range)
250  if( parser.OptionExists('e') ) {
251  LOG("gmkspl", pINFO) << "Reading maximum spline energy";
252  gOptMaxE = parser.ArgAsDouble('e');
253  } else {
254  LOG("gmkspl", pINFO)
255  << "Unspecified maximum spline energy - Using default";
256  gOptMaxE = -1;
257  }
258 
259  // write out input splines?
260  if( parser.OptionExists("no-copy") ) {
261  LOG("gmkspl", pINFO) << "Not copying input splines to output";
262  gOptNoCopy = true;
263  }
264 
265  // comma-separated neutrino PDG code list
266  if( parser.OptionExists('p') ) {
267  LOG("gmkspl", pINFO) << "Reading neutrino PDG codes";
268  gOptNuPdgCodeList = parser.ArgAsString('p');
269  } else {
270  LOG("gmkspl", pFATAL)
271  << "Unspecified neutrino PDG code list - Exiting";
272  PrintSyntax();
273  exit(1);
274  }
275 
276  // comma-separated target PDG code list or input geometry file
277  bool tgt_cmd = true;
278  if( parser.OptionExists('t') ) {
279  LOG("gmkspl", pINFO) << "Reading target nuclei PDG codes";
280  gOptTgtPdgCodeList = parser.ArgAsString('t');
281  } else {
282  LOG("gmkspl", pINFO) << "No code list specified from the command line";
283  tgt_cmd = false;
284  }
285 
286  bool tgt_geom = true;
287  if( parser.OptionExists('f') ) {
288  LOG("gmkspl", pINFO) << "Reading ROOT geometry filename";
289  gOptGeomFilename = parser.ArgAsString('f');
290  } else {
291  LOG("gmkspl", pINFO) << "No geometry file was specified";
292  tgt_cmd = false;
293  }
294 
295  bool both = tgt_geom && tgt_cmd;
296  bool none = !tgt_geom && !tgt_cmd;
297  if(none) {
298  LOG("gmkspl", pFATAL)
299  << "No geom file or cmd line target list was specified - Exiting";
300  PrintSyntax();
301  exit(1);
302  }
303  if(both) {
304  LOG("gmkspl", pFATAL)
305  << "You specified both a geom file and a cmd line target list "
306  << "- Exiting confused";
307  PrintSyntax();
308  exit(1);
309  }
310 
311  // random number seed
312  if( parser.OptionExists("seed") ) {
313  LOG("gmkspl", pINFO) << "Reading random number seed";
314  gOptRanSeed = parser.ArgAsLong("seed");
315  } else {
316  LOG("gmkspl", pINFO) << "Unspecified random number seed - Using default";
317  gOptRanSeed = -1;
318  }
319 
320  // input cross-section file
321  if( parser.OptionExists("input-cross-sections") ) {
322  LOG("gmkspl", pINFO) << "Reading cross-section file";
323  gOptInpXSecFile = parser.ArgAsString("input-cross-sections");
324  } else {
325  LOG("gmkspl", pINFO) << "Unspecified input cross-section file";
326  gOptInpXSecFile = "";
327  }
328 
329  //
330  // print the command-line options
331  //
332  LOG("gmkspl", pNOTICE)
333  << "\n"
334  << utils::print::PrintFramedMesg("gmkspl job configuration")
335  << "\n Neutrino PDG codes : " << gOptNuPdgCodeList
336  << "\n Target PDG codes : " << gOptTgtPdgCodeList
337  << "\n Input ROOT geometry : " << gOptGeomFilename
338  << "\n Output cross-section file : " << gOptOutXSecFile
339  << "\n Input cross-section file : " << gOptInpXSecFile
340  << "\n Random number seed : " << gOptRanSeed
341  << "\n";
342 
343  LOG("gmkspl", pNOTICE) << *RunOpt::Instance();
344 }
345 //____________________________________________________________________________
346 void PrintSyntax(void)
347 {
348  LOG("gmkspl", pNOTICE)
349  << "\n\n" << "Syntax:" << "\n"
350  << " gmkspl -p nupdg <-t tgtpdg, -f geomfile> "
351  << " <-o | --output-cross-section> xsec_xml_file_name"
352  << " [-n nknots] [-e max_energy] "
353  << " [--seed seed_number]"
354  << " [--input-cross-section xml_file]"
355  << " [--event-generator-list list_name]"
356  << " [--xml-path config_xml_dir]"
357  << " [--message-thresholds xml_file]\n\n";
358 }
359 //____________________________________________________________________________
361 {
362  // split the comma separated list
363  vector<string> nuvec = utils::str::Split(gOptNuPdgCodeList, ",");
364 
365  // fill in the PDG code list
366  PDGCodeList * list = new PDGCodeList;
367  vector<string>::const_iterator iter;
368  for(iter = nuvec.begin(); iter != nuvec.end(); ++iter) {
369  list->push_back( atoi(iter->c_str()) );
370  }
371  return list;
372 }
373 //____________________________________________________________________________
375 {
376  bool from_geom_file = ( gOptGeomFilename.size() > 0 );
377  bool from_cmd_line = ( gOptTgtPdgCodeList.size() > 0 );
378 
379  if (from_cmd_line) {
380  // split the comma separated list
381  vector<string> tgtvec = utils::str::Split(gOptTgtPdgCodeList, ",");
382 
383  // fill in the PDG code list
384  PDGCodeList * list = new PDGCodeList;
385  vector<string>::const_iterator iter;
386  for(iter = tgtvec.begin(); iter != tgtvec.end(); ++iter) {
387  list->push_back( atoi(iter->c_str()) );
388  }
389  return list;
390  }
391 
392  if (from_geom_file) {
393 #ifdef __GENIE_GEOM_DRIVERS_ENABLED__
394  // create/configure a geometry driver
395  LOG("gmkspl", pINFO) << "Creating/configuring a ROOT geom. driver";
397 
399 
400  delete geom;
401  return list;
402 #else
403  LOG("gmkspl", pFATAL)
404  << "To read-in a ROOT geometry you need to enable the geometry drivers!";
405  gAbortingInErr = true;
406  exit(1);
407  return 0;
408 #endif
409 
410  }
411  return 0;
412 }
413 //____________________________________________________________________________
void RandGen(long int seed)
Definition: AppInit.cxx:31
PDGCodeList * GetTargetCodes(void)
void XSecTable(string inpfile, bool require_table)
Definition: AppInit.cxx:39
int gOptNKnots
long ArgAsLong(char opt)
PDGCodeList * GetNeutrinoCodes(void)
double ArgAsDouble(char opt)
string ArgAsString(char opt)
THE MAIN GENIE PROJECT NAMESPACE
Definition: GeneratorBase.h:8
string gOptNuPdgCodeList
long int gOptRanSeed
void ReadFromCommandLine(int argc, char **argv)
Definition: RunOpt.cxx:107
string gOptInpXSecFile
#define pFATAL
Definition: Messenger.h:57
GENIE geometry drivers.
static XSecSplineList * Instance()
virtual const PDGCodeList & ListOfTargetNuclei(void)
implement the GeomAnalyzerI interface
int main(int argc, char **argv)
A list of PDG codes.
Definition: PDGCodeList.h:33
string gOptOutXSecFile
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
void SetEventGeneratorList(string listname)
Definition: GEVGDriver.cxx:349
string gOptTgtPdgCodeList
A ROOT/GEANT4 geometry driver.
GENIE Event Generation Driver. A minimalist user interface object for generating neutrino interaction...
Definition: GEVGDriver.h:55
#define pINFO
Definition: Messenger.h:63
void BuildTune()
build tune and inform XSecSplineList
Definition: RunOpt.cxx:100
void SaveAsXml(const string &filename, bool save_init=true) const
::xsd::cxx::tree::string< char, simple_type > string
Definition: Database.h:154
static RunOpt * Instance(void)
Definition: RunOpt.cxx:62
vector< string > Split(string input, string delim)
Definition: StringUtils.cxx:42
void PrintSyntax(void)
void geom(int which=0)
Definition: geom.C:163
exit(0)
void Configure(int nu_pdgc, int Z, int A)
Definition: GEVGDriver.cxx:138
double gOptMaxE
string PrintFramedMesg(string mesg, unsigned int nl=1, const char f='*')
Definition: PrintUtils.cxx:171
A vector of EventGeneratorI objects.
void CreateSplines(int nknots=-1, double emax=-1, bool inLogE=true)
Definition: GEVGDriver.cxx:578
void GetCommandLineArgs(int argc, char **argv)
void MesgThresholds(string inpfile)
Definition: AppInit.cxx:100
Command line argument parser.
#define pNOTICE
Definition: Messenger.h:62
bool gOptNoCopy
bool gAbortingInErr
Definition: Messenger.cxx:56
List of cross section vs energy splines.
bool OptionExists(char opt)
was option set?
void EnableBareXSecPreCalc(bool flag)
Definition: RunOpt.h:61
void push_back(int pdg_code)
Definition: PDGCodeList.cxx:67
string gOptGeomFilename
Initial State information.
Definition: InitialState.h:49