StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
LHEF3.h
1 // LHEF3.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2018 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 
6 // This file is written by Stefan Prestel. The code evolved from
7 // a LHEF 2.0 reader supplied by Leif Lonnblad.
8 // LHEF3.h contains the main class for LHEF 3.0 functionalities.
9 // Header file.
10 
11 #ifndef Pythia8_LHEF3_H
12 #define Pythia8_LHEF3_H
13 
14 #include "Pythia8/PythiaStdlib.h"
15 #include "Pythia8/Streams.h"
16 #include <stdexcept>
17 
18 namespace Pythia8 {
19 
20 //==========================================================================
21 
22 // The XMLTag struct is used to represent all information within an XML tag.
23 // It contains the attributes as a map, any sub-tags as a vector of pointers
24 // to other XMLTag objects, and any other information as a single string.
25 // The XMLTag struct written by Leif Lonnblad.
26 
27 struct XMLTag {
28 
29  // Convenient typdef.
30  typedef string::size_type pos_t;
31 
32  // Convenient alias for npos.
33  static const pos_t end;
34 
35  // The destructor also destroys any sub-tags.
36  ~XMLTag() {
37  for ( int i = 0, N = tags.size(); i < N; ++i )
38  if (tags[i]) delete tags[i];
39  }
40 
41  // The name of this tag.
42  string name;
43 
44  // The attributes of this tag.
45  map<string,string> attr;
46 
47  // A vector of sub-tags.
48  vector<XMLTag*> tags;
49 
50  // The contents of this tag.
51  string contents;
52 
53  // Find an attribute named n and set the double variable v to
54  // the corresponding value. Return false if no attribute was found.
55  bool getattr(string n, double & v) const {
56  map<string,string>::const_iterator it = attr.find(n);
57  if ( it == attr.end() ) return false;
58  v = atof(it->second.c_str());
59  return true;
60  }
61 
62  // Find an attribute named n and set the bool variable v to true if the
63  // corresponding value is "yes". Return false if no attribute was found.
64  bool getattr(string n, bool & v) const {
65  map<string,string>::const_iterator it = attr.find(n);
66  if ( it == attr.end() ) return false;
67  if ( it->second == "yes" ) v = true;
68  return true;
69  }
70 
71  // Find an attribute named n and set the long variable v to the
72  // corresponding value. Return false if no attribute was found.
73  bool getattr(string n, long & v) const {
74  map<string,string>::const_iterator it = attr.find(n);
75  if ( it == attr.end() ) return false;
76  v = atoi(it->second.c_str());
77  return true;
78  }
79 
80  // Find an attribute named n and set the long variable v to the
81  // corresponding value. Return false if no attribute was found.
82  bool getattr(string n, int & v) const {
83  map<string,string>::const_iterator it = attr.find(n);
84  if ( it == attr.end() ) return false;
85  v = int(atoi(it->second.c_str()));
86  return true;
87  }
88 
89  // Find an attribute named n and set the string variable v to the
90  // corresponding value. Return false if no attribute was found.
91  bool getattr(string n, string & v) const {
92  map<string,string>::const_iterator it = attr.find(n);
93  if ( it == attr.end() ) return false;
94  v = it->second;
95  return true;
96  }
97 
98  // Scan the given string and return all XML tags found as a vector
99  // of pointers to XMLTag objects.
100  static vector<XMLTag*> findXMLTags(string str,
101  string * leftover = 0) {
102  vector<XMLTag*> tags;
103  pos_t curr = 0;
104 
105  while ( curr != end ) {
106 
107  // Find the first tag.
108  pos_t begin = str.find("<", curr);
109 
110  // Skip comments.
111  if ( begin != end && str.find("<!--", curr) == begin ) {
112  pos_t endcom = str.find("-->", begin);
113  if ( endcom == end ) {
114  if ( leftover ) *leftover += str.substr(curr);
115  return tags;
116  }
117  if ( leftover ) *leftover += str.substr(curr, endcom - curr);
118  curr = endcom;
119  continue;
120  }
121 
122  // Also skip CDATA statements.
123  // Used for text data that should not be parsed by the XML parser.
124  // (e.g., JavaScript code contains a lot of "<" or "&" characters
125  // which XML would erroneously interpret as the start of a new
126  // element or the start of a character entity, respectively.)
127  // See eg http://www.w3schools.com/xml/xml_cdata.asp
128  if ( str.find("<![CDATA[", curr) == begin ) {
129  pos_t endcom = str.find("]]>", begin);
130  if ( endcom == end ) {
131  if ( leftover ) *leftover += str.substr(curr);
132  return tags;
133  }
134  if ( leftover ) *leftover += str.substr(curr, endcom - curr);
135  curr = endcom;
136  continue;
137  }
138 
139  if ( leftover ) *leftover += str.substr(curr, begin - curr);
140  if ( begin == end || begin > str.length() - 3 || str[begin + 1] == '/' )
141  return tags;
142 
143  pos_t close = str.find(">", curr);
144  if ( close == end ) return tags;
145 
146  // Find the tag name.
147  curr = str.find_first_of(" \t\n/>", begin);
148  tags.push_back(new XMLTag());
149  tags.back()->name = str.substr(begin + 1, curr - begin - 1);
150 
151  while ( true ) {
152 
153  // Now skip some white space to see if we can find an attribute.
154  curr = str.find_first_not_of(" \t\n", curr);
155  if ( curr == end || curr >= close ) break;
156 
157  pos_t tend = str.find_first_of("= \t\n", curr);
158  if ( tend == end || tend >= close ) break;
159 
160  string name = str.substr(curr, tend - curr);
161  curr = str.find("=", curr) + 1;
162 
163  // OK now find the beginning and end of the atribute.
164  curr = str.find("\"", curr);
165  if ( curr == end || curr >= close ) break;
166  pos_t bega = ++curr;
167  curr = str.find("\"", curr);
168  while ( curr != end && str[curr - 1] == '\\' )
169  curr = str.find("\"", curr + 1);
170 
171  string value = str.substr(bega, curr == end? end: curr - bega);
172 
173  tags.back()->attr[name] = value;
174 
175  ++curr;
176 
177  }
178 
179  curr = close + 1;
180  if ( str[close - 1] == '/' ) continue;
181 
182  pos_t endtag = str.find("</" + tags.back()->name + ">", curr);
183  if ( endtag == end ) {
184  tags.back()->contents = str.substr(curr);
185  curr = endtag;
186  } else {
187  tags.back()->contents = str.substr(curr, endtag - curr);
188  curr = endtag + tags.back()->name.length() + 3;
189  }
190 
191  string leftovers;
192  tags.back()->tags = findXMLTags(tags.back()->contents, &leftovers);
193  if ( leftovers.find_first_not_of(" \t\n") == end ) leftovers="";
194  tags.back()->contents = leftovers;
195 
196  }
197 
198  return tags;
199 
200  }
201 
202  // Print out this tag to a stream.
203  void list(ostream & os) const {
204  os << "<" << name;
205  for ( map<string,string>::const_iterator it = attr.begin();
206  it != attr.end(); ++it )
207  os << " " << it->first << "=\"" << it->second << "\"";
208  if ( contents.empty() && tags.empty() ) {
209  os << "/>" << endl;
210  return;
211  }
212  os << ">" << endl;
213  for ( int i = 0, N = tags.size(); i < N; ++i )
214  tags[i]->list(os);
215 
216  os << "````" << contents << "''''</" << name << ">" << endl;
217  }
218 
219 };
220 
221 //==========================================================================
222 
223 // The LHAweights struct represents the information in a weights tag.
224 
225 struct LHAweights {
226 
227  // Initialize default values.
228  LHAweights() {}
229 
230  // Construct from XML tag
231  LHAweights(const XMLTag & tag);
232 
233  // Print out an XML tag.
234  void list(ostream & file) const;
235 
236  // Function to reset this object.
237  void clear() {
238  contents="";
239  weights.clear();
240  attributes.clear();
241  }
242 
243  // The weights of this event.
244  vector<double> weights;
245 
246  // Any other attributes.
247  map<string,string> attributes;
248 
249  // The contents of the tag.
250  string contents;
251 
252  // Return number of weights.
253  int size() { return int(weights.size()); }
254 
255 };
256 
257 //==========================================================================
258 
259 // Collect different scales relevant for an event.
260 
261 struct LHAscales {
262 
263  // Empty constructor.
264  LHAscales(double defscale = -1.0)
265  : muf(defscale), mur(defscale), mups(defscale), SCALUP(defscale) {}
266 
267  // Construct from an XML-tag
268  LHAscales(const XMLTag & tag, double defscale = -1.0);
269 
270  // Print out the corresponding XML-tag.
271  void list(ostream & file) const;
272 
273  // Function to reset this object.
274  void clear() {
275  contents="";
276  muf=mur=mups=SCALUP;
277  attributes.clear();
278  }
279 
280  // The factorization scale used for this event.
281  double muf;
282 
283  // The renormalization scale used for this event.
284  double mur;
285 
286  // The starting scale for the parton shower as suggested by the
287  // matrix element generator.
288  double mups;
289 
290  // Any other scales reported by the matrix element generator.
291  map<string,double> attributes;
292 
293  // The default scale in this event.
294  double SCALUP;
295 
296  // The contents of the tag.
297  string contents;
298 
299 };
300 
301 //==========================================================================
302 
303 // Collect generator information for an event file.
304 
305 struct LHAgenerator {
306 
307  // Empty constructor.
308  LHAgenerator()
309  : name(""), version(""), contents("") {}
310 
311  // Construct from an XML-tag
312  LHAgenerator(const XMLTag & tag, string defname = "");
313 
314  // Print out the corresponding XML-tag.
315  void list(ostream & file) const;
316 
317  // Function to reset this object.
318  void clear() {
319  contents="";
320  name="";
321  version="";
322  attributes.clear();
323  }
324 
325  // The generator name used for this file.
326  string name;
327 
328  // The generator version used for this file.
329  string version;
330 
331  // Any other attributes.
332  map<string,string> attributes;
333 
334  // The contents of the tag.
335  string contents;
336 
337 };
338 
339 //==========================================================================
340 
341 // Collect the wgt information.
342 
343 struct LHAwgt {
344 
345  // Empty constructor.
346  LHAwgt(double defwgt = 1.0)
347  : id(""), contents(defwgt) {}
348 
349  // Construct from an XML-tag
350  LHAwgt(const XMLTag & tag, double defwgt = 1.0);
351 
352  // Print out the corresponding XML-tag.
353  void list(ostream & file) const;
354 
355  // Function to reset this object.
356  void clear() {
357  contents=0.0;
358  id="";
359  attributes.clear();
360  }
361 
362  // The identification number of this wgt tag.
363  string id;
364 
365  // Any other attributes.
366  map<string,string> attributes;
367 
368  // The weight associated to this tag.
369  double contents;
370 
371 };
372 
373 //==========================================================================
374 
375 // Collect the wgt information.
376 
377 struct LHAweight {
378 
379  // Empty constructor.
380  LHAweight(string defname = "")
381  : id(defname), contents(defname) {}
382 
383  // Construct from an XML-tag
384  LHAweight(const XMLTag & tag, string defname = "");
385 
386  // Print out the corresponding XML-tag.
387  void list(ostream & file) const;
388 
389  // Function to reset this object.
390  void clear() {
391  contents="";
392  id="";
393  attributes.clear();
394  }
395 
396  // The identification number of this weight tag.
397  string id;
398 
399  // Any other attributes.
400  map<string,string> attributes;
401 
402  // The weight description associated to this tag.
403  string contents;
404 
405 };
406 
407 //==========================================================================
408 
409 // The LHAweightgroup assigns a group-name to a set of LHAweight objects.
410 
412 
413  // Default constructor;
414  LHAweightgroup() {}
415 
416  // Construct a group of LHAweight objects from an XML tag and
417  // insert them in the given vector.
418  LHAweightgroup(const XMLTag & tag);
419 
420  // Print out the corresponding XML-tag.
421  void list(ostream & file) const;
422 
423  // Function to reset this object.
424  void clear() {
425  contents="";
426  name="";
427  weights.clear();
428  attributes.clear();
429  }
430 
431  // The contents of the tag.
432  string contents;
433 
434  // The name.
435  string name;
436 
437  // The vector of weights.
438  map<string, LHAweight> weights;
439  vector<string> weightsKeys;
440 
441  // Any other attributes.
442  map<string,string> attributes;
443 
444  // Return number of weights.
445  int size() { return int(weights.size()); }
446 
447 };
448 
449 //==========================================================================
450 
451 // The LHArwgt assigns a group-name to a set of LHAwgt objects.
452 
453 struct LHArwgt {
454 
455  // Default constructor;
456  LHArwgt() {}
457 
458  // Construct a group of LHAwgt objects from an XML tag and
459  // insert them in the given vector.
460  LHArwgt(const XMLTag & tag);
461 
462  // Print out the corresponding XML-tag.
463  void list(ostream & file) const;
464 
465  // Function to reset this object.
466  void clear() {
467  contents="";
468  wgts.clear();
469  attributes.clear();
470  }
471 
472  // The contents of the tag.
473  string contents;
474 
475  // The map of weights.
476  map<string, LHAwgt> wgts;
477  vector<string> wgtsKeys;
478 
479  // Any other attributes.
480  map<string,string> attributes;
481 
482  // Return number of weights.
483  int size() { return int(wgts.size()); }
484 
485 };
486 
487 //==========================================================================
488 
489 // The LHAinitrwgt assigns a group-name to a set of LHAweightgroup objects.
490 
491 struct LHAinitrwgt {
492 
493  // Default constructor;
494  LHAinitrwgt() {}
495 
496  // Construct a group of LHAweightgroup objects from an XML tag and
497  // insert them in the given vector.
498  LHAinitrwgt(const XMLTag & tag);
499 
500  // Print out the corresponding XML-tag.
501  void list(ostream & file) const;
502 
503  // Function to reset this object.
504  void clear() {
505  contents="";
506  weights.clear();
507  weightgroups.clear();
508  attributes.clear();
509  }
510 
511  // The contents of the tag.
512  string contents;
513 
514  // The vector of weight's.
515  map<string, LHAweight> weights;
516  vector<string> weightsKeys;
517 
518  // The vector of weightgroup's.
519  map<string, LHAweightgroup> weightgroups;
520  vector<string> weightgroupsKeys;
521 
522  // Any other attributes.
523  map<string,string> attributes;
524 
525  // Return number of weights.
526  int size() { return int(weights.size());}
527 
528  // Return number of weights.
529  int sizeWeightGroups() { return int(weightgroups.size()); }
530 
531 };
532 
533 //==========================================================================
534 
535 // The HEPRUP class is a simple container corresponding to the Les Houches
536 // accord (<A HREF="http://arxiv.org/abs/hep-ph/0109068">hep-ph/0109068</A>)
537 // common block with the same name. The members are named in the same
538 // way as in the common block. However, fortran arrays are represented
539 // by vectors, except for the arrays of length two which are
540 // represented by pair objects.
541 
542 class HEPRUP {
543 
544 public:
545 
546  // Default constructor.
547  HEPRUP() : IDWTUP(0), NPRUP(0) {}
548 
549  // Assignment operator.
550  HEPRUP & operator=(const HEPRUP & x) {
551  IDBMUP = x.IDBMUP;
552  EBMUP = x.EBMUP;
553  PDFGUP = x.PDFGUP;
554  PDFSUP = x.PDFSUP;
555  IDWTUP = x.IDWTUP;
556  NPRUP = x.NPRUP;
557  XSECUP = x.XSECUP;
558  XERRUP = x.XERRUP;
559  XMAXUP = x.XMAXUP;
560  LPRUP = x.LPRUP;
561  initrwgt = x.initrwgt;
562  generators = x.generators;
563  weightgroups = x.weightgroups;
564  weights = x.weights;
565  return *this;
566  }
567 
568  // Destructor.
569  ~HEPRUP() {}
570 
571  // Set the NPRUP variable, corresponding to the number of
572  // sub-processes, to \a nrup, and resize all relevant vectors
573  // accordingly.
574  void resize(int nrup) {
575  NPRUP = nrup;
576  resize();
577  }
578 
579  // Assuming the NPRUP variable, corresponding to the number of
580  // sub-processes, is correctly set, resize the relevant vectors
581  // accordingly.
582  void resize() {
583  XSECUP.resize(NPRUP);
584  XERRUP.resize(NPRUP);
585  XMAXUP.resize(NPRUP);
586  LPRUP.resize(NPRUP);
587  }
588 
589  // Clear all members.
590  void clear();
591 
592  // PDG id's of beam particles. (first/second is in +/-z direction).
593  pair<long,long> IDBMUP;
594 
595  // Energy of beam particles given in GeV.
596  pair<double,double> EBMUP;
597 
598  // The author group for the PDF used for the beams according to the
599  // PDFLib specification.
600  pair<int,int> PDFGUP;
601 
602  // The id number the PDF used for the beams according to the
603  // PDFLib specification.
604  pair<int,int> PDFSUP;
605 
606  // Master switch indicating how the ME generator envisages the
607  // events weights should be interpreted according to the Les Houches
608  // accord.
609  int IDWTUP;
610 
611  // The number of different subprocesses in this file.
612  int NPRUP;
613 
614  // The cross sections for the different subprocesses in pb.
615  vector<double> XSECUP;
616 
617  // The statistical error in the cross sections for the different
618  // subprocesses in pb.
619  vector<double> XERRUP;
620 
621  // The maximum event weights (in HEPEUP::XWGTUP) for different
622  // subprocesses.
623  vector<double> XMAXUP;
624 
625  // The subprocess code for the different subprocesses.
626  vector<int> LPRUP;
627 
628  // Contents of the LHAinitrwgt tag
629  LHAinitrwgt initrwgt;
630 
631  // Contents of the LHAgenerator tags.
632  vector<LHAgenerator> generators;
633 
634  // A map of the LHAweightgroup tags, indexed by name.
635  map<string,LHAweightgroup> weightgroups;
636 
637  // A map of the LHAweight tags, indexed by name.
638  map<string,LHAweight> weights;
639 
640 };
641 
642 //==========================================================================
643 
644 // The HEPEUP class is a simple container corresponding to the Les Houches
645 // accord (<A HREF="http://arxiv.org/abs/hep-ph/0109068">hep-ph/0109068</A>)
646 // common block with the same name. The members are named in the same
647 // way as in the common block. However, fortran arrays are represented
648 // by vectors, except for the arrays of length two which are
649 // represented by pair objects.
650 
651 class HEPEUP {
652 
653 public:
654 
655  // Default constructor.
656  HEPEUP()
657  : NUP(0), IDPRUP(0), XWGTUP(0.0), XPDWUP(0.0, 0.0),
658  SCALUP(0.0), AQEDUP(0.0), AQCDUP(0.0), heprup(0) {}
659 
660  // Copy constructor
661  HEPEUP(const HEPEUP & x) {
662  operator=(x);
663  }
664 
665  // Copy information from the given HEPEUP.
666  HEPEUP & setEvent(const HEPEUP & x);
667 
668  // Assignment operator.
669  HEPEUP & operator=(const HEPEUP & x) {
670  clear();
671  setEvent(x);
672  return *this;
673  }
674 
675  // Destructor.
676  ~HEPEUP() {
677  clear();
678  };
679 
680  // Reset the HEPEUP object.
681  void reset();
682 
683  // Clear the HEPEUP object.
684  void clear() {
685  reset();
686  }
687 
688  // Set the NUP variable, corresponding to the number of particles in
689  // the current event, to \a nup, and resize all relevant vectors
690  // accordingly.
691  void resize(int nup) {
692  NUP = nup;
693  resize();
694  }
695 
696  // Return the main weight for this event.
697  double weight() const {
698  return XWGTUP;
699  }
700 
701  // Assuming the NUP variable, corresponding to the number of
702  // particles in the current event, is correctly set, resize the
703  // relevant vectors accordingly.
704  void resize();
705 
706  // The number of particle entries in the current event.
707  int NUP;
708 
709  // The subprocess code for this event (as given in LPRUP).
710  int IDPRUP;
711 
712  // The weight for this event.
713  double XWGTUP;
714 
715  // The PDF weights for the two incoming partons. Note that this
716  // variable is not present in the current LesHouches accord
717  // (<A HREF="http://arxiv.org/abs/hep-ph/0109068">hep-ph/0109068</A>),
718  // hopefully it will be present in a future accord.
719  pair<double,double> XPDWUP;
720 
721  // The scale in GeV used in the calculation of the PDF's in this
722  // event.
723  double SCALUP;
724 
725  // The value of the QED coupling used in this event.
726  double AQEDUP;
727 
728  // The value of the QCD coupling used in this event.
729  double AQCDUP;
730 
731  // The PDG id's for the particle entries in this event.
732  vector<long> IDUP;
733 
734  // The status codes for the particle entries in this event.
735  vector<int> ISTUP;
736 
737  // Indices for the first and last mother for the particle entries in
738  // this event.
739  vector< pair<int,int> > MOTHUP;
740 
741  // The colour-line indices (first(second) is (anti)colour) for the
742  // particle entries in this event.
743  vector< pair<int,int> > ICOLUP;
744 
745  // Lab frame momentum (Px, Py, Pz, E and M in GeV) for the particle
746  // entries in this event.
747  vector< vector<double> > PUP;
748 
749  // Invariant lifetime (c*tau, distance from production to decay in
750  // mm) for the particle entries in this event.
751  vector<double> VTIMUP;
752 
753  // Spin info for the particle entries in this event given as the
754  // cosine of the angle between the spin vector of a particle and the
755  // 3-momentum of the decaying particle, specified in the lab frame.
756  vector<double> SPINUP;
757 
758  // A pointer to the current HEPRUP object.
759  HEPRUP * heprup;
760 
761  // The weights associated with this event, as given by the LHAwgt tags.
762  map<string,double> weights_detailed;
763 
764  // The weights associated with this event, as given by the LHAweights tags.
765  vector<double> weights_compressed;
766 
767  // Contents of the LHAscales tag
768  LHAscales scalesSave;
769 
770  // Contents of the LHAweights tag (compressed format)
771  LHAweights weightsSave;
772 
773  // Contents of the LHArwgt tag (detailed format)
774  LHArwgt rwgtSave;
775 
776  // Any other attributes.
777  map<string,string> attributes;
778 
779 };
780 
781 //==========================================================================
782 
783 // The Reader class is initialized with a stream from which to read a
784 // version 1/3 Les Houches Accord event file. In the constructor of
785 // the Reader object the optional header information is read and then
786 // the mandatory init is read. After this the whole header block
787 // including the enclosing lines with tags are available in the public
788 // headerBlock member variable. Also the information from the init
789 // block is available in the heprup member variable and any additional
790 // comment lines are available in initComments. After each successful
791 // call to the readEvent() function the standard Les Houches Accord
792 // information about the event is available in the hepeup member
793 // variable and any additional comments in the eventComments
794 // variable. A typical reading sequence would look as follows:
795 
796 class Reader {
797 
798 public:
799 
800  // Initialize the Reader with a filename from which to read an event
801  // file. After the constructor is called the whole header block
802  // including the enclosing lines with tags are available in the
803  // public headerBlock member variable. Also the information from the
804  // init block is available in the heprup member variable and any
805  // additional comment lines are available in initComments.
806  //
807  // filename: the name of the file to read from.
808  //
809  Reader(string filenameIn)
810  : filename(filenameIn), intstream(NULL), file(NULL) {
811  intstream = new igzstream(filename.c_str());
812  file = intstream;
813  isGood = init();
814  }
815 
816  Reader(istream* is)
817  : filename(""), intstream(NULL), file(is) {
818  isGood = init();
819  }
820 
821  // Clean up
822  ~Reader() {
823  if (intstream) delete intstream;
824  }
825 
826 
827 
828  // (Re)initialize the Reader with a filename from which to read an event
829  // file. After this, all information from the header and init block is
830  // available.
831  //
832  // filename: File name (not used as the input file stream is given)
833  // isIn : Name of the input file stream.
834  bool setup(string filenameIn) {
835  filename = filenameIn;
836  if (intstream) delete intstream;
837  intstream = new igzstream(filename.c_str());
838  file = intstream;
839  isGood = init();
840  return isGood;
841  }
842 
843 private:
844 
845  // Used internally in the constructors to read header and init blocks.
846  bool init();
847 
848 public:
849 
850  // Read an event from the file and store it in the hepeup
851  // object. Optional comment lines are stored in the eventComments
852  // member variable.
853  bool readEvent(HEPEUP * peup = 0);
854 
855  // Reset values of all event-related members to their defaults.
856  void clearEvent() {
857  currentLine = "";
858  hepeup.clear();
859  eventComments = "";
860  weights_detailed_vec.resize(0);
861  }
862 
863 protected:
864 
865  // Used internally to read a single line from the stream.
866  bool getLine() {
867  currentLine = "";
868  if(!getline(*file, currentLine)) return false;
869  // Replace single by double quotes
870  replace(currentLine.begin(),currentLine.end(),'\'','\"');
871  return true;
872  }
873 
874 protected:
875 
876  // Name of file-to-be-read.
877  string filename;
878 
879  // A local stream which is unused if a stream is supplied from the
880  // outside.
881  igzstream* intstream;
882 
883  // The stream we are reading from. This may be a pointer to an
884  // external stream or the internal intstream.
885  istream * file;
886 
887  // The last line read in from the stream in getline().
888  string currentLine;
889 
890 public:
891 
892  // Save if the initialisation worked.
893  bool isGood;
894 
895  // XML file version
896  int version;
897 
898  // All lines (since the last readEvent()) outside the header, init
899  // and event tags.
900  string outsideBlock;
901 
902  // All lines from the header block.
903  string headerBlock;
904  string headerComments;
905 
906  // The standard init information.
907  HEPRUP heprup;
908 
909  // Additional comments found in the init block.
910  string initComments;
911 
912  // The standard information about the last read event.
913  HEPEUP hepeup;
914 
915  // Additional comments found with the last read event.
916  string eventComments;
917 
918  // The detailed weights associated with this event, linearized to a vector.
919  vector<double> weights_detailed_vec;
920  vector<double> weights_detailed_vector() { return weights_detailed_vec; }
921 
922 private:
923 
924  // The default constructor should never be used.
925  Reader();
926 
927  // The copy constructor should never be used.
928  Reader(const Reader &);
929 
930  // The Reader cannot be assigned to.
931  Reader & operator=(const Reader &);
932 
933 };
934 
935 //==========================================================================
936 
937 // The Writer class is initialized with a stream to which to write a
938 // version 1.0 or 3.0 Les Houches Accord event file. In the init() function of
939 // the Writer object the main XML tag, header and init blocks are written,
940 // with the corresponding end tag is written by list_end_tag().
941 // After a Writer object (in the following called "writer") has been created,
942 // it is possible to assign version (3 by default) information by
943 //
944 // writer.version = <value>;
945 //
946 // The header block (called "someHeaderString" below) is assigned by
947 //
948 // writer.headerBlock() << someHeaderString;
949 //
950 // and the init block comments (called "someInitString" below) are assigned via
951 //
952 // writer.initComments() << someInitString;
953 //
954 // The standard init information (including amendments for LHEF 3.0) can
955 // be assigned by the heprup member variable:
956 //
957 // writer.heprup = heprup;
958 //
959 // where heprup is an object of type HEPRUP. All of the above information
960 // will be writen by calling the init() function.
961 //
962 // Before each event is written out with the writeEvent() function,
963 // the standard event information can be assigned to the hepeup
964 // variable by
965 //
966 // writer.hepeup = hepeup;
967 //
968 // where hepeup is of type HEPEUP. Event comments (called
969 // "someCommentString" below) can be assigned through
970 //
971 // writer.eventComments() << someCommentString;
972 //
973 // All of this event information is written by the writeEvent() function.
974 
975 class Writer {
976 
977 public:
978 
979  // Create a Writer object giving a stream to write to.
980  // @param os the stream where the event file is written.
981  Writer(ostream & os)
982  : file(os), version(3) {}
983 
984  // Create a Writer object giving a filename to write to.
985  // @param filename the name of the event file to be written.
986  Writer(string filename)
987  : intstream(filename.c_str()), file(intstream), version(3) {}
988 
989  // The destructor.
990  ~Writer() {}
991 
992  // Add header lines consisting of XML code with this stream.
993  ostream & headerBlock() {
994  return headerStream;
995  }
996 
997  // Add comment lines to the init block with this stream.
998  ostream & initComments() {
999  return initStream;
1000  }
1001 
1002  // Add comment lines to the next event to be written out with this stream.
1003  ostream & eventComments() {
1004  return eventStream;
1005  }
1006 
1007  // Write out the final XML end-tag.
1008  void list_end_tag() {
1009  file << "</LesHouchesEvents>" << endl;
1010  }
1011 
1012  // Write out an optional header block followed by the standard init
1013  // block information together with any comment lines.
1014  void init();
1015 
1016  // Write out the event stored in hepeup, followed by optional
1017  // comment lines.
1018  bool writeEvent(HEPEUP * peup = 0, int pDigits = 15);
1019  // Write out an event as a string.
1020  string getEventString(HEPEUP * peup = 0);
1021 
1022 protected:
1023 
1024  // Make sure that each line in the string \a s starts with a
1025  // #-character and that the string ends with a new-line.
1026  string hashline(string s, bool comment = false);
1027 
1028 protected:
1029 
1030  // A local stream which is unused if a stream is supplied from the
1031  // outside.
1032  ofstream intstream;
1033 
1034  // The stream we are writing to. This may be a reference to an
1035  // external stream or the internal intstream.
1036  ostream & file;
1037 
1038 public:
1039 
1040  // Stream to add all lines in the header block.
1041  ostringstream headerStream;
1042 
1043  // The standard init information.
1044  HEPRUP heprup;
1045 
1046  // Stream to add additional comments to be put in the init block.
1047  ostringstream initStream;
1048 
1049  // The standard information about the event we will write next.
1050  HEPEUP hepeup;
1051 
1052  // Stream to add additional comments to be written together the next event.
1053  ostringstream eventStream;
1054 
1055  // XML file version
1056  int version;
1057 
1058 private:
1059 
1060  // The default constructor should never be used.
1061  Writer();
1062 
1063  // The copy constructor should never be used.
1064  Writer(const Writer &);
1065 
1066  // The Writer cannot be assigned to.
1067  Writer & operator=(const Writer &);
1068 
1069 };
1070 
1071 //==========================================================================
1072 
1073 } // end namespace Pythia8
1074 
1075 #endif // end Pythia8_LHEF3_H