gEvGenHadronNucleus.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*!
3 
4 \program gevgen_hadron
5 
6 \brief Generates hadron + nucleus interactions using GENIE's INTRANUKE
7  Similar to NEUGEN's pitest (S.Dytman & H.Gallagher)
8 
9  Syntax :
10  gevgen_hadron [-n nev] -p probe -t tgt [-r run#] -k KE
11  [-f flux] [-o prefix] [-m mode]
12  [--seed random_number_seed]
13  [--message-thresholds xml_file]
14  [--event-record-print-level level]
15  [--mc-job-status-refresh-rate rate]
16 
17  Options :
18  [] Denotes an optional argument
19  -n
20  Specifies the number of events to generate (default: 10000)
21  -p
22  Specifies the incoming hadron PDG code
23  -t
24  Specifies the nuclear target PDG code (10LZZZAAAI)
25  -r
26  Specifies the MC run number (default: 0)
27  -k
28  Specifies the incoming hadron's kinetic energy (in GeV).
29  The same option can be use to specify a kinetic energy range as
30  a comma-separated set of numbers (eg 0.1,1.2).
31  If no flux is specified then the hadrons will be fired with a
32  uniform kinetic energy distribution within the specified range.
33  If a kinetic energy spectrum is supplied then the hadrons will be
34  fired with a kinetic energy following the input spectrum within
35  the specified range.
36  -f
37  Specifies the incoming hadron's kinetic energy spectrum -
38  it can be either a function, eg 'x*x+4*exp(-x)' or a text file
39  containing 2 columns corresponding to (kinetic energy {GeV}, 'flux').
40  -o
41  Output filename prefix
42  -m
43  INTRANUKE mode <hA, hN> (default: hA)
44  --seed
45  Random number seed.
46  --message-thresholds
47  Allows users to customize the message stream thresholds.
48  The thresholds are specified using an XML file.
49  See $GENIE/config/Messenger.xml for the XML schema.
50  --event-record-print-level
51  Allows users to set the level of information shown when the event
52  record is printed in the screen. See GHepRecord::Print().
53  --mc-job-status-refresh-rate
54  Allows users to customize the refresh rate of the status file.
55 
56  Examples:
57 
58  (1) Generate 100k pi^{+}+Fe56 events with a pi^{+} kinetic energy
59  of 165 MeV:
60  % gevgen_hadron -n 100000 -p 211 -t 1000260560 -k 0.165
61 
62  (2) Generate 100k pi^{+}+Fe56 events with the pi^{+} kinetic energy
63  distributed uniformly in the [165 MeV, 1200 MeV] range:
64  % gevgen_hadron -n 100000 -p 211 -t 1000260560 -k 0.165,1.200
65 
66  (3) Generate 100k pi^{+}+Fe56 events with the pi^{+} kinetic energy
67  distributed as f(KE) = 1/KE in the [165 MeV, 1200 MeV] range:
68  % ghAevgen -n gevgen_hadron -p 211 -t 1000260560 -k 0.165,1.200 -f '1/x'
69 
70 \authors Steve Dytman, Minsuk Kim and Aaron Meyer
71  University of Pittsburgh
72 
73  Costas Andreopoulos,
74  University of Liverpool & STFC Rutherford Appleton Lab
75 
76 \version 1.3
77 
78 \created May 1, 2007
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 
89 #include <TSystem.h>
90 #include <TFile.h>
91 #include <TTree.h>
92 #include <TH1D.h>
93 #include <TF1.h>
94 
111 #include "Framework/Utils/AppInit.h"
115 #include "Framework/Utils/RunOpt.h"
117 
120 
121 using namespace genie;
122 using namespace genie::controls;
123 
124 using namespace genie::utils::intranuke;
125 
126 // Function prototypes
127 void GetCommandLineArgs (int argc, char ** argv);
128 const EventRecordVisitorI * GetIntranuke (void);
129 double GenProbeKineticEnergy (void);
131 void BuildSpectrum (void);
132 void PrintSyntax (void);
133 
134 // Default options
135 int kDefOptNevents = 10000; // n-events to generate
136 Long_t kDefOptRunNu = 0; // default run number
137 string kDefOptEvFilePrefix = "gntp.inuke"; // default output file prefix
138 string kDefOptMode = "hA"; // default mode
139 
140 // User-specified options:
141 string gOptMode; // mode variable
142 Long_t gOptRunNu; // run number
143 int gOptNevents; // n-events to generate
144 int gOptProbePdgCode; // probe PDG code
145 int gOptTgtPdgCode; // target PDG code
146 double gOptProbeKE; // incoming hadron kinetic enegy (GeV) - for monoenergetic probes
147 double gOptProbeKEmin; // incoming hadron kinetic enegy (GeV) - if using flux
148 double gOptProbeKEmax; // incoming hadron kinetic enegy (GeV) - if using flux
149 string gOptFlux; // input flux (function or flux file)
150 string gOptEvFilePrefix; // event file prefix
151 bool gOptUsingFlux=false; // using kinetic energy distribution?
152 long int gOptRanSeed ; // random number seed
153 
154 TH1D * gSpectrum = 0;
155 
156 //____________________________________________________________________________
157 int main(int argc, char ** argv)
158 {
159  // Parse command line arguments
160  GetCommandLineArgs(argc,argv);
161 
162  if ( ! RunOpt::Instance()->Tune() ) {
163  LOG("gmkspl", pFATAL) << " No TuneId in RunOption";
164  exit(-1);
165  }
167 
168  // Init random number generator generator with user-specified seed number,
169  // set user-specified mesg thresholds, set user-specified GHEP print-level
170  utils::app_init::MesgThresholds(RunOpt::Instance()->MesgThresholdFiles());
172  GHepRecord::SetPrintLevel(RunOpt::Instance()->EventRecordPrintLevel());
173 
174  // Build the incident hadron kinetic energy spectrum, if required
175  BuildSpectrum();
176 
177  LOG("gevgen_hadron",pNOTICE) << "finish setup";
178 
179  // Get the specified INTRANUKE model
180  const EventRecordVisitorI * intranuke = GetIntranuke();
181 
182  // Initialize an Ntuple Writer to save GHEP records into a ROOT tree
183  NtpWriter ntpw(kNFGHEP, gOptRunNu);
185  ntpw.Initialize();
186 
187  // Create an MC job monitor
188  GMCJMonitor mcjmonitor(gOptRunNu);
189 
190  LOG("gevgen_hadron",pNOTICE) << "ready to generate events";
191 
192  //
193  // Generate events
194  //
195 
196  int ievent = 0;
197  while (ievent < gOptNevents) {
198  LOG("gevgen_hadron", pNOTICE)
199  << " *** Generating event............ " << ievent;
200 
201  // initialize
202  EventRecord * evrec = InitializeEvent();
203 
204  // generate full h+A event
205  intranuke->ProcessEventRecord(evrec);
206 
207  // print generated events
208  LOG("gevgen_hadron", pNOTICE ) << *evrec;
209 
210  // add event at the output ntuple
211  ntpw.AddEventRecord(ievent, evrec);
212 
213  // refresh the mc job monitor
214  mcjmonitor.Update(ievent,evrec);
215 
216  ievent++;
217  delete evrec;
218 
219  } // end loop events
220 
221  // Save the generated MC events
222  ntpw.Save();
223 
224  // Clean-up
225  if(gSpectrum) {
226  delete gSpectrum;
227  gSpectrum = 0;
228  }
229 
230  return 0;
231 }
232 //____________________________________________________________________________
234 {
235 // get the requested INTRANUKE module
236 
237  string sname = "";
238  string sconf = "";
239 
240  if(gOptMode.compare("hA")==0) {
241  sname = "genie::HAIntranuke";
242  sconf = "Default";
243  }
244  else
245  if(gOptMode.compare("hN")==0) {
246  sname = "genie::HNIntranuke";
247  sconf = "Default";
248  }
249  else
250  if(gOptMode.compare("hA2019")==0) {
251  sname = "genie::HAIntranuke2019";
252  sconf = "Default";
253  }
254  else
255  if(gOptMode.compare("hN2019")==0) {
256  sname = "genie::HNIntranuke2019";
257  sconf = "Default";
258  }
259  else
260  if(gOptMode.compare("hA2018")==0) {
261  sname = "genie::HAIntranuke2018";
262  sconf = "Default";
263  }
264  else
265  if(gOptMode.compare("hN2018")==0) {
266  sname = "genie::HNIntranuke2018";
267  sconf = "Default";
268  }
269  else {
270  LOG("gevgen_hadron", pFATAL) << "Invalid Intranuke mode - Exiting";
271  gAbortingInErr = true;
272  exit(1);
273  }
274 
275  AlgFactory * algf = AlgFactory::Instance();
276  const EventRecordVisitorI * intranuke =
277  dynamic_cast<const EventRecordVisitorI *> (algf->GetAlgorithm(sname,sconf));
278  assert(intranuke);
279 
280  return intranuke;
281 }
282 //____________________________________________________________________________
284 {
285 // Initialize event record. Inserting the probe and target particles.
286 
287  EventRecord * evrec = new EventRecord();
288  Interaction * interaction = new Interaction;
289  evrec->AttachSummary(interaction);
290 
291  // dummy vertex position
292  TLorentzVector x4null(0.,0.,0.,0.);
293 
294  // incident hadron & target nucleon masses
295  PDGLibrary * pdglib = PDGLibrary::Instance();
296  double mh = pdglib -> Find (gOptProbePdgCode) -> Mass();
297  double M = pdglib -> Find (gOptTgtPdgCode ) -> Mass();
298 
299  // incident hadron kinetic energy
300  double ke = GenProbeKineticEnergy();
301 
302  // form incident hadron and target 4-momenta
303  double Eh = mh + ke;
304  double pzh = TMath::Sqrt(TMath::Max(0.,Eh*Eh-mh*mh));
305  TLorentzVector p4h (0.,0.,pzh,Eh);
306  TLorentzVector p4tgt (0.,0.,0., M);
307 
308  // insert probe and target entries
310  evrec->AddParticle(gOptProbePdgCode, ist, -1,-1,-1,-1, p4h, x4null);
311  evrec->AddParticle(gOptTgtPdgCode, ist, -1,-1,-1,-1, p4tgt, x4null);
312 
313  return evrec;
314 }
315 //____________________________________________________________________________
317 {
318  if(gOptUsingFlux) return gSpectrum->GetRandom(); // spectrum
319  else return gOptProbeKE; // mono-energetic
320 }
321 //____________________________________________________________________________
322 void BuildSpectrum(void)
323 {
324 // Create kinetic energy spectrum from input function
325 //
326 
327  if(!gOptUsingFlux) return;
328 
329  if(gSpectrum) {
330  delete gSpectrum;
331  gSpectrum = 0;
332  }
333 
334  LOG("gevgen_hadron", pNOTICE)
335  << "Generating a flux histogram ... ";
336 
337  int flux_bins = 300;
338  int flux_entries = 1000000;
339  double ke_min = gOptProbeKEmin;
340  double ke_max = gOptProbeKEmax;
341  double dke = ke_max - ke_min;
342  assert(dke>0 && ke_min>=0.);
343 
344  // kinetic energy spectrum
345  gSpectrum = new TH1D(
346  "spectrum","hadron kinetic energy spectrum", flux_bins, ke_min, ke_max);
347 
348  // check whether the input flux is a file or a functional form
349  bool input_is_file = ! gSystem->AccessPathName(gOptFlux.c_str());
350  if(input_is_file) {
351  Spline * input_flux = new Spline(gOptFlux.c_str());
352 
353  // generate the flux hisogram from the input flux file
354  int n = 100;
355  double ke_step = (ke_max-ke_min)/(n-1);
356  double ymax = -1, ry = -1, gy = -1, ke = -1;
357  for(int i=0; i<n; i++) {
358  ke = ke_min + i*ke_step;
359  ymax = TMath::Max(ymax, input_flux->Evaluate(ke));
360  }
361  ymax *= 1.3;
362 
364 
365  for(int ientry=0; ientry<flux_entries; ientry++) {
366  bool accept = false;
367  unsigned int iter=0;
368  while(!accept) {
369  iter++;
370  if(iter > kRjMaxIterations) {
371  LOG("gevgen_hadron", pFATAL)
372  << "Couldn't generate a flux histogram";
373  gAbortingInErr = true;
374  exit(1);
375  }
376  ke = ke_min + dke * r->RndGen().Rndm();
377  gy = ymax * r->RndGen().Rndm();
378  ry = input_flux->Evaluate(ke);
379  accept = gy < ry;
380  if(accept) gSpectrum->Fill(ke);
381  }
382  }
383  delete input_flux;
384 
385  } else {
386  // generate the flux histogram from the input functional form
387  TF1 * input_func = new TF1("input_func", gOptFlux.c_str(), ke_min, ke_max);
388  gSpectrum->FillRandom("input_func", flux_entries);
389  delete input_func;
390  }
391  TFile f("./input-hadron-flux.root","recreate");
392  gSpectrum->Write();
393  f.Close();
394 }
395 //____________________________________________________________________________
396 void GetCommandLineArgs(int argc, char ** argv)
397 {
398  LOG("gevgen_hadron", pINFO) << "Parsing command line arguments";
399 
400  // Common run options.
402 
403  // Parse run options for this app
404 
405  CmdLnArgParser parser(argc,argv);
406 
407  // number of events
408  if( parser.OptionExists('n') ) {
409  LOG("gevgen_hadron", pINFO) << "Reading number of events to generate";
410  gOptNevents = parser.ArgAsInt('n');
411  } else {
412  LOG("gevgen_hadron", pINFO)
413  << "Unspecified number of events to generate - Using default";
415  }
416 
417  // run number
418  if( parser.OptionExists('r') ) {
419  LOG("gevgen_hadron", pINFO) << "Reading MC run number";
420  gOptRunNu = parser.ArgAsLong('r');
421  } else {
422  LOG("gevgen_hadron", pINFO) << "Unspecified run number - Using default";
424  }
425 
426  // incoming hadron PDG code
427  if( parser.OptionExists('p') ) {
428  LOG("gevgen_hadron", pINFO) << "Reading rescattering particle PDG code";
429  gOptProbePdgCode = parser.ArgAsInt('p');
430  } else {
431  LOG("gevgen_hadron", pFATAL) << "Unspecified PDG code - Exiting";
432  PrintSyntax();
433  gAbortingInErr = true;
434  exit(1);
435  }
436 
437  // target PDG code
438  if( parser.OptionExists('t') ) {
439  LOG("gevgen_hadron", pINFO) << "Reading target PDG code";
440  gOptTgtPdgCode = parser.ArgAsInt('t');
441  } else {
442  LOG("gevgen_hadron", pFATAL) << "Unspecified target PDG code - Exiting";
443  PrintSyntax();
444  gAbortingInErr = true;
445  exit(1);
446  }
447 
448  // target PDG code
449  if( parser.OptionExists('m') ) {
450  LOG("gevgen_hadron", pINFO) << "Reading mode";
451  gOptMode = parser.ArgAsString('m');
452  } else {
453  LOG("gevgen_hadron", pFATAL) << "Unspecified mode - Using default";
455  }
456 
457  // flux functional form or flux file
458  if( parser.OptionExists('f') ) {
459  LOG("gevgen_hadron", pINFO) << "Reading hadron's kinetic energy spectrum";
460  gOptFlux = parser.ArgAsString('f');
461  gOptUsingFlux = true;
462  }
463 
464  // incoming hadron kinetic energy (or kinetic energy range, if using flux)
465  if( parser.OptionExists('k') ) {
466  LOG("gevgen_hadron", pINFO) << "Reading probe kinetic energy";
467  string ke = parser.ArgAsString('k');
468  // is it just a value or a range (comma separated set of values)
469  if(ke.find(",") != string::npos) {
470  // split the comma separated list
471  vector<string> kerange = utils::str::Split(ke, ",");
472  assert(kerange.size() == 2);
473  double kemin = atof(kerange[0].c_str());
474  double kemax = atof(kerange[1].c_str());
475  assert(kemax>kemin && kemin>0);
476  gOptProbeKE = -1;
479  // if no flux was specified, generate uniformly within given range
480  if(!gOptUsingFlux) {
481  gOptFlux = "1";
482  gOptUsingFlux = true;
483  }
484  } else {
485  gOptProbeKE = atof(ke.c_str());
486  gOptProbeKEmin = -1;
487  gOptProbeKEmax = -1;
488  if(gOptUsingFlux) {
489  LOG("gevgen_hadron", pFATAL)
490  << "You specified an input flux without a kinetic energy range";
491  PrintSyntax();
492  gAbortingInErr = true;
493  exit(1);
494  }
495  }
496  } else {
497  LOG("gevgen_hadron", pFATAL) << "Unspecified kinetic energy - Exiting";
498  PrintSyntax();
499  gAbortingInErr = true;
500  exit(1);
501  }
502 
503  // event file prefix
504  if( parser.OptionExists('o') ) {
505  LOG("gevgen_hadron", pINFO) << "Reading the event filename prefix";
506  gOptEvFilePrefix = parser.ArgAsString('o');
507  } else {
508  LOG("gevgen_hadron", pDEBUG)
509  << "Will set the default event filename prefix";
511  } //-o
512 
513  // INTRANUKE mode
514  if( parser.OptionExists('m') ) {
515  LOG("gevgen_hadron", pINFO) << "Reading mode";
516  gOptMode = parser.ArgAsString('m');
517  } else {
518  LOG("gevgen_hadron", pDEBUG)
519  << "Unspecified mode - Using default";
521  }
522 
523  // random number seed
524  if( parser.OptionExists("seed") ) {
525  LOG("gevgen_hadron", pINFO) << "Reading random number seed";
526  gOptRanSeed = parser.ArgAsLong("seed");
527  } else {
528  LOG("gevgen_hadron", pINFO) << "Unspecified random number seed - Using default";
529  gOptRanSeed = -1;
530  }
531 
532 
533  LOG("gevgen_hadron", pNOTICE)
534  << "\n"
535  << utils::print::PrintFramedMesg("gevgen_hadron job configuration");
536 
537  LOG("gevgen_hadron", pNOTICE) << "MC Run Number = " << gOptRunNu;
538  LOG("gevgen_hadron", pNOTICE) << "Random number seed = " << gOptRanSeed;
539  LOG("gevgen_hadron", pNOTICE) << "Mode = " << gOptMode;
540  LOG("gevgen_hadron", pNOTICE) << "Number of events = " << gOptNevents;
541  LOG("gevgen_hadron", pNOTICE) << "Probe PDG code = " << gOptProbePdgCode;
542  LOG("gevgen_hadron", pNOTICE) << "Target PDG code = " << gOptTgtPdgCode;
543  if(gOptProbeKEmin<0 && gOptProbeKEmax<0) {
544  LOG("gevgen_hadron", pNOTICE)
545  << "Hadron input KE = " << gOptProbeKE;
546  } else {
547  LOG("gevgen_hadron", pNOTICE)
548  << "Hadron input KE range = ["
549  << gOptProbeKEmin << ", " << gOptProbeKEmax << "]";
550  }
551  if(gOptUsingFlux) {
552  LOG("gevgen_hadron", pNOTICE)
553  << "Input flux = "
554  << gOptFlux;
555  }
556 
557  LOG("gevgen_hadron", pNOTICE) << "\n";
558  LOG("gevgen_hadron", pNOTICE) << *RunOpt::Instance();
559 }
560 //____________________________________________________________________________
561 void PrintSyntax(void)
562 {
563  LOG("gevgen_hadron", pNOTICE)
564  << "\n\n"
565  << "Syntax:" << "\n"
566  << " gevgen_hadron [-r run] [-n nev] -p hadron_pdg -t tgt_pdg -k KE [-m mode] "
567  << " [-f flux] "
568  << " [--seed random_number_seed]"
569  << " [--message-thresholds xml_file]"
570  << " [--event-record-print-level level]"
571  << " [--mc-job-status-refresh-rate rate]"
572  << "\n";
573 }
574 //____________________________________________________________________________
void RandGen(long int seed)
Definition: AppInit.cxx:31
static void SetPrintLevel(int print_level)
Definition: GHepRecord.cxx:992
int gOptNevents
long ArgAsLong(char opt)
void GetCommandLineArgs(int argc, char **argv)
Long_t kDefOptRunNu
string ArgAsString(char opt)
THE MAIN GENIE PROJECT NAMESPACE
Definition: GeneratorBase.h:8
virtual void ProcessEventRecord(GHepRecord *event_rec) const =0
double gOptProbeKEmax
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
Defines the EventRecordVisitorI interface. Concrete implementations of this interface use the &#39;Visito...
long int gOptRanSeed
Long_t gOptRunNu
void ReadFromCommandLine(int argc, char **argv)
Definition: RunOpt.cxx:107
enum genie::EGHepStatus GHepStatus_t
A numeric analysis tool class for interpolating 1-D functions.
Definition: Spline.h:47
#define pFATAL
Definition: Messenger.h:57
string kDefOptMode
const double kemax
void Update(int iev, const EventRecord *event)
Definition: GMCJMonitor.cxx:58
double GenProbeKineticEnergy(void)
double Mass(Resonance_t res)
resonance mass (GeV)
TH1D * gSpectrum
double Evaluate(double x) const
Definition: Spline.cxx:362
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
virtual void AttachSummary(Interaction *interaction)
Definition: GHepRecord.cxx:143
int kDefOptNevents
EventRecord * InitializeEvent(void)
Double_t ymax
Definition: plot.C:25
void BuildSpectrum(void)
Summary information for an interaction.
Definition: Interaction.h:56
Simple class to create & update MC job status files and env. vars. This is used to be able to keep tr...
Definition: GMCJMonitor.h:30
bool gOptUsingFlux
const char * Find(const char *fname)
Definition: Icons.cxx:12
int main(int argc, char **argv)
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
void PrintSyntax(void)
const Algorithm * GetAlgorithm(const AlgId &algid)
Definition: AlgFactory.cxx:86
void Save(void)
get the even tree
Definition: NtpWriter.cxx:214
#define pINFO
Definition: Messenger.h:63
int gOptProbePdgCode
string gOptFlux
void AddEventRecord(int ievent, const EventRecord *ev_rec)
save the event tree
Definition: NtpWriter.cxx:69
Misc GENIE control constants.
string gOptMode
void BuildTune()
build tune and inform XSecSplineList
Definition: RunOpt.cxx:100
Generated Event Record. It is a GHepRecord object that can accept / be visited by EventRecordVisitorI...
Definition: EventRecord.h:38
string gOptEvFilePrefix
void CustomizeFilenamePrefix(string prefix)
Definition: NtpWriter.cxx:123
void Initialize(void)
add event
Definition: NtpWriter.cxx:95
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:43
static RunOpt * Instance(void)
Definition: RunOpt.cxx:62
vector< string > Split(string input, string delim)
Definition: StringUtils.cxx:42
static AlgFactory * Instance()
Definition: AlgFactory.cxx:75
A utility class to facilitate creating the GENIE MC Ntuple from the output GENIE GHEP event records...
Definition: NtpWriter.h:40
Singleton class to load & serve a TDatabasePDG.
Definition: PDGLibrary.h:30
Var Sqrt(const Var &v)
Use to take sqrt of a var.
Definition: Var.cxx:326
TRandom3 & RndGen(void) const
rnd number generator for generic usage
Definition: RandomGen.h:81
string kDefOptEvFilePrefix
exit(0)
assert(nhit_max >=nhit_nbins)
TRandom3 r(0)
string PrintFramedMesg(string mesg, unsigned int nl=1, const char f='*')
Definition: PrintUtils.cxx:171
virtual void AddParticle(const GHepParticle &p)
Definition: GHepRecord.cxx:535
void MesgThresholds(string inpfile)
Definition: AppInit.cxx:100
double gOptProbeKEmin
Command line argument parser.
#define pNOTICE
Definition: Messenger.h:62
The GENIE Algorithm Factory.
Definition: AlgFactory.h:40
static const unsigned int kRjMaxIterations
Definition: Controls.h:27
double gOptProbeKE
const double kemin
bool gAbortingInErr
Definition: Messenger.cxx:56
Most commonly used PDG codes. A set of utility functions to handle PDG codes is provided in PDGUtils...
const EventRecordVisitorI * GetIntranuke(void)
bool OptionExists(char opt)
was option set?
#define pDEBUG
Definition: Messenger.h:64
int gOptTgtPdgCode