Line data Source code
1 : // LHEF3.cc is a part of the PYTHIA event generator.
2 : // Copyright (C) 2015 Torbjorn Sjostrand.
3 : // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
4 : // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 :
6 : // This file is written by Stefan Prestel.
7 : // It contains the main class for LHEF 3.0 functionalities.
8 : // Function definitions.
9 :
10 : #include "Pythia8/LHEF3.h"
11 :
12 : namespace Pythia8 {
13 :
14 : //==========================================================================
15 :
16 : // The XMLTag struct is used to represent all information within an XML tag.
17 : // It contains the attributes as a map, any sub-tags as a vector of pointers
18 : // to other XMLTag objects, and any other information as a single string.
19 :
20 : //==========================================================================
21 :
22 : // The LHAweights struct.
23 :
24 : //--------------------------------------------------------------------------
25 :
26 : // Construct from XML tag.
27 :
28 0 : LHAweights::LHAweights(const XMLTag & tag) {
29 0 : for ( map<string,string>::const_iterator it = tag.attr.begin();
30 0 : it != tag.attr.end(); ++it ) {
31 0 : string v = it->second.c_str();
32 0 : attributes[it->first] = v;
33 0 : }
34 :
35 0 : contents = tag.contents;
36 :
37 0 : istringstream iss(tag.contents);
38 0 : double w;
39 0 : while ( iss >> w ) weights.push_back(w);
40 0 : }
41 :
42 : //--------------------------------------------------------------------------
43 :
44 : // Print out.
45 :
46 : void LHAweights::print(ostream & file) const {
47 0 : file << "<weights";
48 0 : for ( map<string,string>::const_iterator it = attributes.begin();
49 0 : it != attributes.end(); ++it )
50 0 : file << " " << it->first << "=\"" << it->second << "\"";
51 0 : file << " >";
52 0 : for ( int j = 0, M = weights.size(); j < M; ++j ) file << " " << weights[j];
53 0 : file << "</weights>" << endl;
54 0 : }
55 :
56 : //==========================================================================
57 :
58 : // The LHAscales struct: Collect different scales relevant for an event.
59 :
60 : //--------------------------------------------------------------------------
61 :
62 : // Construct from an XML-tag.
63 :
64 0 : LHAscales::LHAscales(const XMLTag & tag, double defscale)
65 0 : : muf(defscale), mur(defscale), mups(defscale), SCALUP(defscale) {
66 0 : for ( map<string,string>::const_iterator it = tag.attr.begin();
67 0 : it != tag.attr.end(); ++it ) {
68 0 : double v = atof(it->second.c_str());
69 0 : if ( it->first == "muf" ) muf = v;
70 0 : else if ( it->first == "mur" ) mur = v;
71 0 : else if ( it->first == "mups" ) mups = v;
72 0 : else attributes[it->first] = v;
73 : }
74 0 : contents = tag.contents;
75 0 : }
76 :
77 : //--------------------------------------------------------------------------
78 :
79 : // Print out the corresponding XML-tag.
80 :
81 : void LHAscales::print(ostream & file) const {
82 0 : file << "<scales";
83 0 : file << " muf=\"" << muf << "\"";
84 0 : file << " mur=\"" << mur << "\"";
85 0 : file << " mups=\"" << mups << "\"";
86 0 : for ( map<string,double>::const_iterator it = attributes.begin();
87 0 : it != attributes.end(); ++it )
88 0 : file << " " << it->first << "=\"" << it->second << "\"";
89 0 : file << contents;
90 0 : file << "</scales>" << endl;
91 0 : }
92 :
93 : //==========================================================================
94 :
95 : // The LHAgenerator struct: Collect generator information for an event file.
96 :
97 : //--------------------------------------------------------------------------
98 :
99 : // Construct from an XML-tag
100 :
101 0 : LHAgenerator::LHAgenerator(const XMLTag & tag, string defname)
102 0 : : name(defname), version(defname), contents(defname) {
103 0 : for ( map<string,string>::const_iterator it = tag.attr.begin();
104 0 : it != tag.attr.end(); ++it ) {
105 0 : string v = it->second.c_str();
106 0 : if ( it->first == "name" ) name = v;
107 0 : else if ( it->first == "version" ) version = v;
108 0 : else attributes[it->first] = v;
109 0 : }
110 0 : contents = tag.contents;
111 0 : }
112 :
113 : //--------------------------------------------------------------------------
114 :
115 : // Print out the corresponding XML-tag.
116 :
117 : void LHAgenerator::print(ostream & file) const {
118 0 : file << "<generator";
119 0 : if ( name != "" ) file << " name=\"" << name << "\"";
120 0 : if ( version != "" ) file << " version=\"" << version << "\"";
121 0 : for ( map<string,string>::const_iterator it = attributes.begin();
122 0 : it != attributes.end(); ++it )
123 0 : file << " " << it->first << "=\"" << it->second << "\"";
124 0 : file << " >";
125 0 : file << contents;
126 0 : file << "</generator>" << endl;
127 0 : }
128 :
129 : //==========================================================================
130 :
131 : // The LHAwgt struct: Collect the wgt information.
132 :
133 : //--------------------------------------------------------------------------
134 :
135 : // Construct from an XML-tag
136 :
137 0 : LHAwgt::LHAwgt(const XMLTag & tag, double defwgt)
138 0 : : id(""), contents(defwgt) {
139 0 : for ( map<string,string>::const_iterator it = tag.attr.begin();
140 0 : it != tag.attr.end(); ++it ) {
141 0 : string v = it->second;
142 0 : if ( it->first == "id" ) id = v;
143 0 : else attributes[it->first] = v;
144 0 : }
145 0 : contents = atof(tag.contents.c_str());
146 0 : }
147 :
148 : //--------------------------------------------------------------------------
149 :
150 : // Print out the corresponding XML-tag.
151 :
152 : void LHAwgt::print(ostream & file) const {
153 0 : file << "<wgt";
154 0 : if ( id != "" ) file << " id=\"" << id << "\"";
155 0 : for ( map<string,string>::const_iterator it = attributes.begin();
156 0 : it != attributes.end(); ++it )
157 0 : file << " " << it->first << "=\"" << it->second << "\"";
158 0 : file << " >";
159 0 : file << contents;
160 0 : file << "</wgt>" << endl;
161 0 : }
162 :
163 : //==========================================================================
164 :
165 : // The LHAweight struct: Collect the weight information.
166 :
167 : //--------------------------------------------------------------------------
168 :
169 : // Construct from an XML-tag.
170 :
171 0 : LHAweight::LHAweight(const XMLTag & tag, string defname)
172 0 : : id(defname), contents(defname) {
173 0 : for ( map<string,string>::const_iterator it = tag.attr.begin();
174 0 : it != tag.attr.end(); ++it ) {
175 0 : string v = it->second;
176 0 : if ( it->first == "id" ) id = v;
177 0 : else attributes[it->first] = v;
178 0 : }
179 0 : contents = tag.contents;
180 0 : }
181 :
182 : //--------------------------------------------------------------------------
183 :
184 : // Print out the corresponding XML-tag.
185 :
186 : void LHAweight::print(ostream & file) const {
187 0 : file << "<weight";
188 0 : if ( id != "" ) file << " id=\"" << id << "\"";
189 0 : for ( map<string,string>::const_iterator it = attributes.begin();
190 0 : it != attributes.end(); ++it )
191 0 : file << " " << it->first << "=\"" << it->second << "\"";
192 0 : file << " >";
193 0 : file << contents;
194 0 : file << "</weight>" << endl;
195 0 : }
196 :
197 : //==========================================================================
198 :
199 : // The LHAweightgroup struct: The LHAweightgroup assigns a group-name to a set
200 : // of LHAweight objects.
201 :
202 : //--------------------------------------------------------------------------
203 :
204 : // Construct a group of LHAweight objects from an XML tag and
205 : // insert them in the given vector.
206 :
207 0 : LHAweightgroup::LHAweightgroup(const XMLTag & tag) {
208 :
209 0 : for ( map<string,string>::const_iterator it = tag.attr.begin();
210 0 : it != tag.attr.end(); ++it ) {
211 0 : string v = it->second.c_str();
212 0 : if ( it->first == "name" ) name = v;
213 0 : else attributes[it->first] = v;
214 0 : }
215 :
216 0 : contents = tag.contents;
217 :
218 : // Now add the weight's step by step.
219 0 : string s;
220 0 : vector<XMLTag*> tags = XMLTag::findXMLTags(tag.contents, &s);
221 0 : for ( int i = 0, N = tags.size(); i < N; ++i ) {
222 0 : const XMLTag & tagnow = *tags[i];
223 0 : LHAweight wt(tagnow);
224 0 : weights.insert(make_pair(wt.id, wt));
225 0 : }
226 0 : for ( int i = 0, N = tag.tags.size(); i < N; ++i ) {
227 0 : const XMLTag & tagnow = *tag.tags[i];
228 0 : const LHAweight & wt(tagnow);
229 0 : weights.insert(make_pair(wt.id, wt));
230 0 : }
231 :
232 0 : for ( int i = 0, N = tags.size(); i < N; ++i ) if (tags[i]) delete tags[i];
233 :
234 0 : }
235 :
236 : //--------------------------------------------------------------------------
237 :
238 : // Print out the corresponding XML-tag.
239 :
240 : void LHAweightgroup::print(ostream & file) const {
241 0 : file << "<weightgroup";
242 0 : if ( name != "" ) file << " name=\"" << name << "\"";
243 0 : for ( map<string,string>::const_iterator it = attributes.begin();
244 0 : it != attributes.end(); ++it )
245 0 : file << " " << it->first << "=\"" << it->second << "\"";
246 0 : file << " >\n";
247 0 : for ( map<string,LHAweight>::const_iterator it = weights.begin();
248 0 : it != weights.end(); ++it ) it->second.print(file);
249 0 : file << "</weightgroup>" << endl;
250 0 : }
251 :
252 : //==========================================================================
253 :
254 : // The LHArwgt struct: Assigns a group-name to a set of LHAwgt objects.
255 :
256 : //--------------------------------------------------------------------------
257 :
258 : // Construct a group of LHAwgt objects from an XML tag and
259 : // insert them in the given vector.
260 :
261 0 : LHArwgt::LHArwgt(const XMLTag & tag) {
262 :
263 0 : for ( map<string,string>::const_iterator it = tag.attr.begin();
264 0 : it != tag.attr.end(); ++it ) {
265 0 : string v = it->second.c_str();
266 0 : attributes[it->first] = v;
267 0 : }
268 0 : contents = tag.contents;
269 :
270 : // Now add the wgt's step by step.
271 0 : string s;
272 0 : vector<XMLTag*> tags = XMLTag::findXMLTags(tag.contents, &s);
273 0 : for ( int i = 0, N = tags.size(); i < N; ++i ) {
274 0 : const XMLTag & tagnow = *tags[i];
275 0 : LHAwgt wt(tagnow);
276 0 : wgts.insert(make_pair(wt.id, wt));
277 0 : }
278 0 : for ( int i = 0, N = tag.tags.size(); i < N; ++i ) {
279 0 : const XMLTag & tagnow = *tag.tags[i];
280 0 : LHAwgt wt(tagnow);
281 0 : wgts.insert(make_pair(wt.id, wt));
282 0 : }
283 :
284 0 : for ( int i = 0, N = tags.size(); i < N; ++i ) if (tags[i]) delete tags[i];
285 :
286 0 : }
287 :
288 : //--------------------------------------------------------------------------
289 :
290 : // Print out the corresponding XML-tag.
291 :
292 : void LHArwgt::print(ostream & file) const {
293 0 : file << "<rwgt";
294 0 : for ( map<string,string>::const_iterator it = attributes.begin();
295 0 : it != attributes.end(); ++it )
296 0 : file << " " << it->first << "=\"" << it->second << "\"";
297 0 : file << " >\n";
298 0 : for ( map<string,LHAwgt>::const_iterator it = wgts.begin();
299 0 : it != wgts.end(); ++it ) it->second.print(file);
300 0 : file << "</rgwt>" << endl;
301 0 : }
302 :
303 : //==========================================================================
304 :
305 : // The LHAinitrwgt assigns a group-name to a set of LHAweightgroup objects.
306 :
307 : //--------------------------------------------------------------------------
308 :
309 : // Construct a group of LHAweightgroup objects from an XML tag and
310 : // insert them in the given vector.
311 :
312 0 : LHAinitrwgt::LHAinitrwgt(const XMLTag & tag) {
313 0 : for ( map<string,string>::const_iterator it = tag.attr.begin();
314 0 : it != tag.attr.end(); ++it ) {
315 0 : string v = it->second.c_str();
316 0 : attributes[it->first] = v;
317 0 : }
318 0 : contents = tag.contents;
319 :
320 : // Now add the wgt's step by step.
321 0 : string s;
322 0 : vector<XMLTag*> tags = XMLTag::findXMLTags(tag.contents, &s);
323 0 : for ( int i = 0, N = tags.size(); i < N; ++i ) {
324 0 : const XMLTag & tagnow = *tags[i];
325 0 : if ( tagnow.name == "weightgroup" ) {
326 0 : LHAweightgroup wgroup(tagnow);
327 0 : string wgname = wgroup.name;
328 0 : weightgroups.insert(make_pair(wgname, wgroup));
329 0 : string ss;
330 0 : vector<XMLTag*> tags2 = XMLTag::findXMLTags(tagnow.contents, &ss);
331 0 : for ( int k = 0, M = tags2.size(); k < M; ++k ) {
332 0 : const XMLTag & tagnow2 = *tags2[k];
333 0 : if ( tagnow2.name == "weight" ) {
334 0 : LHAweight wt(tagnow2);
335 0 : string wtname = wt.id;
336 0 : weights.insert(make_pair(wtname, wt));
337 0 : }
338 : }
339 0 : for ( int j = 0, M = tags2.size(); j < M; ++j )
340 0 : if (tags2[j]) delete tags2[j];
341 0 : } else if ( tagnow.name == "weight" ) {
342 0 : LHAweight wt(tagnow);
343 0 : string wtname = wt.id;
344 0 : weights.insert(make_pair(wtname, wt));
345 0 : }
346 : }
347 :
348 : // Now add the wgt's step by step.
349 0 : for ( int i = 0, N = tag.tags.size(); i < N; ++i ) {
350 0 : const XMLTag & tagnow = *tag.tags[i];
351 0 : if ( tagnow.name == "weightgroup" ) {
352 0 : LHAweightgroup wgroup(tagnow);
353 0 : string wgname = wgroup.name;
354 0 : weightgroups.insert(make_pair(wgname, wgroup));
355 0 : string ss;
356 0 : vector<XMLTag*> tags2 = XMLTag::findXMLTags(tagnow.contents, &ss);
357 0 : for ( int k = 0, M = tags2.size(); k < M; ++k ) {
358 0 : const XMLTag & tagnow2 = *tags2[k];
359 0 : if ( tagnow2.name == "weight" ) {
360 0 : LHAweight wt(tagnow2);
361 0 : string wtname = wt.id;
362 0 : weights.insert(make_pair(wtname, wt));
363 0 : }
364 : }
365 0 : for ( int k = 0, M = tagnow.tags.size(); k < M; ++k ) {
366 0 : const XMLTag & tagnow2 = *tagnow.tags[k];
367 0 : if ( tagnow2.name == "weight" ) {
368 0 : LHAweight wt(tagnow2);
369 0 : string wtname = wt.id;
370 0 : weights.insert(make_pair(wtname, wt));
371 0 : }
372 : }
373 0 : for ( int j = 0, M = tags2.size(); j < M; ++j )
374 0 : if (tags2[j]) delete tags2[j];
375 0 : } else if ( tagnow.name == "weight" ) {
376 0 : LHAweight wt(tagnow);
377 0 : string wtname = wt.id;
378 0 : weights.insert(make_pair(wtname, wt));
379 0 : }
380 : }
381 :
382 0 : for ( int i = 0, N = tags.size(); i < N; ++i ) if (tags[i]) delete tags[i];
383 :
384 0 : }
385 :
386 : //--------------------------------------------------------------------------
387 :
388 : // Print out the corresponding XML-tag.
389 :
390 : void LHAinitrwgt::print(ostream & file) const {
391 0 : file << "<initrwgt";
392 0 : for ( map<string,string>::const_iterator it = attributes.begin();
393 0 : it != attributes.end(); ++it )
394 0 : file << " " << it->first << "=\"" << it->second << "\"";
395 0 : file << " >\n";
396 0 : for ( map<string,LHAweightgroup>::const_iterator it = weightgroups.begin();
397 0 : it != weightgroups.end(); ++it ) it->second.print(file);
398 0 : for ( map<string,LHAweight>::const_iterator it = weights.begin();
399 0 : it != weights.end(); ++it ) it->second.print(file);
400 0 : file << "</initrgwt>" << endl;
401 0 : }
402 :
403 : //==========================================================================
404 :
405 : // The HEPEUP class is a simple container corresponding to the Les Houches
406 : // accord (<A HREF="http://arxiv.org/abs/hep-ph/0109068">hep-ph/0109068</A>)
407 : // common block with the same name. The members are named in the same
408 : // way as in the common block. However, fortran arrays are represented
409 : // by vectors, except for the arrays of length two which are
410 : // represented by pair objects.
411 :
412 : //--------------------------------------------------------------------------
413 :
414 : // Copy information from the given HEPEUP.
415 :
416 : HEPEUP & HEPEUP::setEvent(const HEPEUP & x) {
417 :
418 0 : NUP = x.NUP;
419 0 : IDPRUP = x.IDPRUP;
420 0 : XWGTUP = x.XWGTUP;
421 0 : XPDWUP = x.XPDWUP;
422 0 : SCALUP = x.SCALUP;
423 0 : AQEDUP = x.AQEDUP;
424 0 : AQCDUP = x.AQCDUP;
425 0 : IDUP = x.IDUP;
426 0 : ISTUP = x.ISTUP;
427 0 : MOTHUP = x.MOTHUP;
428 0 : ICOLUP = x.ICOLUP;
429 0 : PUP = x.PUP;
430 0 : VTIMUP = x.VTIMUP;
431 0 : SPINUP = x.SPINUP;
432 0 : heprup = x.heprup;
433 0 : scales = x.scales;
434 0 : weights = x.weights;
435 0 : weights_detailed = x.weights_detailed;
436 0 : weights_compressed = x.weights_compressed;
437 0 : rwgt = x.rwgt;
438 0 : attributes = x.attributes;
439 0 : return *this;
440 :
441 : }
442 :
443 : //--------------------------------------------------------------------------
444 :
445 : // Reset the HEPEUP object.
446 :
447 : void HEPEUP::reset() {
448 0 : NUP = 0;
449 0 : weights_detailed.clear();
450 0 : weights_compressed.clear();
451 0 : weights.clear();
452 0 : rwgt.clear();
453 0 : scales.clear();
454 0 : attributes.clear();
455 0 : }
456 :
457 : //--------------------------------------------------------------------------
458 :
459 : // Assuming the NUP variable, corresponding to the number of
460 : // particles in the current event, is correctly set, resize the
461 : // relevant vectors accordingly.
462 :
463 : void HEPEUP::resize() {
464 0 : IDUP.resize(NUP);
465 0 : ISTUP.resize(NUP);
466 0 : MOTHUP.resize(NUP);
467 0 : ICOLUP.resize(NUP);
468 0 : PUP.resize(NUP, vector<double>(5));
469 0 : VTIMUP.resize(NUP);
470 0 : SPINUP.resize(NUP);
471 0 : }
472 :
473 : //==========================================================================
474 :
475 : // The Reader class is initialized with a stream from which to read a
476 : // version 1/2 Les Houches Accord event file. In the constructor of
477 : // the Reader object the optional header information is read and then
478 : // the mandatory init is read. After this the whole header block
479 : // including the enclosing lines with tags are available in the public
480 : // headerBlock member variable. Also the information from the init
481 : // block is available in the heprup member variable and any additional
482 : // comment lines are available in initComments. After each successful
483 : // call to the readEvent() function the standard Les Houches Accord
484 : // information about the event is available in the hepeup member
485 : // variable and any additional comments in the eventComments
486 : // variable. A typical reading sequence would look as follows:
487 :
488 : //--------------------------------------------------------------------------
489 :
490 : // Used internally in the constructors to read header and init blocks.
491 : bool Reader::init() {
492 :
493 : bool readingHeader = false;
494 : bool readingInit = false;
495 :
496 : // Make sure we are reading a LHEF file:
497 0 : getLine();
498 :
499 0 : if ( currentLine.find("<LesHouchesEvents" ) == string::npos )
500 0 : return false;
501 0 : version = 0;
502 0 : if ( currentLine.find("version=\"1" ) != string::npos )
503 0 : version = 1;
504 0 : else if ( currentLine.find("version=\"2" ) != string::npos )
505 0 : version = 2;
506 0 : else if ( currentLine.find("version=\"3" ) != string::npos )
507 0 : version = 3;
508 : else
509 0 : return false;
510 :
511 : // Loop over all lines until we hit the </init> tag.
512 0 : while ( getLine() && currentLine.find("</init>") == string::npos ) {
513 0 : if ( currentLine.find("<header") != string::npos ) {
514 : // We have hit the header block, so we should dump this and
515 : // all following lines to headerBlock until we hit the end of
516 : // it.
517 : readingHeader = true;
518 0 : headerBlock = currentLine + "\n";
519 0 : }
520 0 : else if ( currentLine.find("<init>") != string::npos
521 0 : || currentLine.find("<init ") != string::npos ) {
522 : // We have hit the init block, so we should expect to find the
523 : // standard information in the following.
524 : readingInit = true;
525 :
526 : // The first line tells us how many lines to read next.
527 0 : getLine();
528 0 : istringstream iss(currentLine);
529 0 : if ( !( iss >> heprup.IDBMUP.first >> heprup.IDBMUP.second
530 0 : >> heprup.EBMUP.first >> heprup.EBMUP.second
531 0 : >> heprup.PDFGUP.first >> heprup.PDFGUP.second
532 0 : >> heprup.PDFSUP.first >> heprup.PDFSUP.second
533 0 : >> heprup.IDWTUP >> heprup.NPRUP ) ) {
534 0 : heprup.NPRUP = -42;
535 0 : return false;
536 : }
537 0 : heprup.resize();
538 :
539 0 : for ( int i = 0; i < heprup.NPRUP; ++i ) {
540 0 : getLine();
541 0 : istringstream isss(currentLine);
542 0 : if ( !( isss >> heprup.XSECUP[i] >> heprup.XERRUP[i]
543 0 : >> heprup.XMAXUP[i] >> heprup.LPRUP[i] ) ) {
544 0 : heprup.NPRUP = -42;
545 0 : return false;
546 : }
547 0 : }
548 0 : }
549 0 : else if ( currentLine.find("</header>") != string::npos ) {
550 : // The end of the header block. Dump this line as well to the
551 : // headerBlock and we're done.
552 : readingHeader = false;
553 0 : headerBlock += currentLine + "\n";
554 0 : }
555 0 : else if ( readingHeader ) {
556 : // We are in the process of reading the header block. Dump the
557 : // line to headerBlock.
558 0 : headerBlock += currentLine + "\n";
559 0 : headerComments += currentLine + "\n";
560 0 : }
561 0 : else if ( readingInit ) {
562 : // Here we found a comment line. Dump it to initComments.
563 0 : initComments += currentLine + "\n";
564 0 : }
565 : else {
566 : // We found some other stuff outside the standard tags.
567 0 : outsideBlock += currentLine + "\n";
568 : }
569 : }
570 :
571 0 : if ( file == NULL ) heprup.NPRUP = -42;
572 :
573 : // Scan the header block for XML tags
574 0 : string leftovers;
575 0 : vector<XMLTag*> tags1 = XMLTag::findXMLTags(headerComments, &leftovers);
576 0 : if ( leftovers.find_first_not_of(" \t\n") == string::npos )
577 0 : leftovers="";
578 :
579 0 : for ( int i = 0, N = tags1.size(); i < N; ++i ) {
580 0 : const XMLTag & tag = *tags1[i];
581 :
582 0 : if ( tag.name == "initrwgt" ) {
583 0 : LHAinitrwgt irwgt(tag);
584 0 : heprup.initrwgt = irwgt;
585 0 : for ( int j = 0, M = tag.tags.size(); j < M; ++j ) {
586 0 : XMLTag & ctag = *tag.tags[j];
587 0 : if ( ctag.name == "weightgroup" ) {
588 0 : LHAweightgroup wgroup(ctag);
589 0 : string wgname = wgroup.name;
590 0 : heprup.weightgroups.insert(make_pair(wgname, wgroup));
591 :
592 0 : string ss;
593 0 : vector<XMLTag*> tags2 = XMLTag::findXMLTags(ctag.contents, &ss);
594 0 : for ( int k = 0, O = tags2.size(); k < O; ++k ) {
595 0 : const XMLTag & tagnow2 = *tags2[k];
596 0 : if ( tagnow2.name == "weight" ) {
597 0 : LHAweight wt(tagnow2);
598 0 : string wtname = wt.id;
599 0 : heprup.weights.insert(make_pair(wtname, wt));
600 0 : }
601 : }
602 0 : for ( int k = 0, O = ctag.tags.size(); k < O; ++k ) {
603 0 : const XMLTag & tagnow2 = *ctag.tags[k];
604 0 : if ( tagnow2.name == "weight" ) {
605 0 : LHAweight wt(tagnow2);
606 0 : string wtname = wt.id;
607 0 : heprup.weights.insert(make_pair(wtname, wt));
608 0 : }
609 : }
610 0 : } else if ( ctag.name == "weight" ) {
611 0 : string tname = ctag.attr["id"];
612 0 : heprup.weights.insert(make_pair(tname, LHAweight(ctag)));
613 0 : }
614 : }
615 0 : }
616 : }
617 :
618 0 : heprup.generators.clear();
619 : // Scan the init block for XML tags
620 0 : leftovers="";
621 0 : vector<XMLTag*> tags2 = XMLTag::findXMLTags(initComments, &leftovers);
622 0 : if ( leftovers.find_first_not_of(" \t\n") == string::npos )
623 0 : leftovers="";
624 :
625 0 : for ( int i = 0, N = tags2.size(); i < N; ++i ) {
626 0 : const XMLTag & tag = *tags2[i];
627 0 : if ( tag.name == "generator" ) {
628 0 : heprup.generators.push_back(LHAgenerator(tag));
629 0 : }
630 : }
631 :
632 0 : for ( int i = 0, N = tags1.size(); i < N; ++i )
633 0 : if (tags1[i]) delete tags1[i];
634 0 : for ( int i = 0, N = tags2.size(); i < N; ++i )
635 0 : if (tags2[i]) delete tags2[i];
636 :
637 : // Done
638 : return true;
639 :
640 0 : }
641 :
642 : //--------------------------------------------------------------------------
643 :
644 : // Read an event from the file and store it in the hepeup
645 : // object. Optional comment lines are stored in the eventComments
646 : // member variable. return true if the read was successful.
647 :
648 : bool Reader::readEvent(HEPEUP * peup) {
649 :
650 0 : HEPEUP & eup = (peup? *peup: hepeup);
651 0 : eup.clear();
652 0 : eup.heprup = &heprup;
653 :
654 : // Check if the initialization was successful. Otherwise we will
655 : // not read any events.
656 0 : if ( heprup.NPRUP < 0 ) return false;
657 0 : eventComments = "";
658 0 : outsideBlock = "";
659 0 : eup.NUP = 0;
660 :
661 : // Keep reading lines until we hit the next event or the end of
662 : // the event block. Save any inbetween lines. Exit if we didn't
663 : // find an event.
664 0 : while ( getLine() && currentLine.find("<event") == string::npos )
665 0 : outsideBlock += currentLine + "\n";
666 :
667 : // Get event attributes.
668 0 : if (currentLine != "") {
669 0 : string eventLine(currentLine);
670 0 : eventLine += "</event>";
671 0 : vector<XMLTag*> evtags = XMLTag::findXMLTags(eventLine);
672 0 : XMLTag & evtag = *evtags[0];
673 0 : for ( map<string,string>::const_iterator it = evtag.attr.begin();
674 0 : it != evtag.attr.end(); ++it ) {
675 0 : string v = it->second.c_str();
676 0 : eup.attributes[it->first] = v;
677 0 : }
678 0 : for ( int i = 0, N = evtags.size(); i < N; ++i )
679 0 : if (evtags[i]) delete evtags[i];
680 0 : }
681 :
682 0 : if ( !getLine() ) return false;
683 :
684 : // We found an event. The first line determines how many
685 : // subsequent particle lines we have.
686 0 : istringstream iss(currentLine);
687 0 : if ( !( iss >> eup.NUP >> eup.IDPRUP >> eup.XWGTUP
688 0 : >> eup.SCALUP >> eup.AQEDUP >> eup.AQCDUP ) )
689 0 : return false;
690 0 : eup.resize();
691 :
692 : // Read all particle lines.
693 0 : for ( int i = 0; i < eup.NUP; ++i ) {
694 0 : if ( !getLine() ) return false;
695 0 : istringstream isss(currentLine);
696 0 : if ( !( isss >> eup.IDUP[i] >> eup.ISTUP[i]
697 0 : >> eup.MOTHUP[i].first >> eup.MOTHUP[i].second
698 0 : >> eup.ICOLUP[i].first >> eup.ICOLUP[i].second
699 0 : >> eup.PUP[i][0] >> eup.PUP[i][1] >> eup.PUP[i][2]
700 0 : >> eup.PUP[i][3] >> eup.PUP[i][4]
701 0 : >> eup.VTIMUP[i] >> eup.SPINUP[i] ) )
702 0 : return false;
703 0 : }
704 :
705 : // Now read any additional comments.
706 0 : while ( getLine() && currentLine.find("</event>") == string::npos )
707 0 : eventComments += currentLine + "\n";
708 :
709 0 : if ( file == NULL ) return false;
710 :
711 0 : eup.scales = LHAscales(eup.SCALUP);
712 :
713 : // Scan the init block for XML tags
714 0 : string leftovers;
715 0 : vector<XMLTag*> tags = XMLTag::findXMLTags(eventComments, &leftovers);
716 0 : if ( leftovers.find_first_not_of(" \t\n") == string::npos )
717 0 : leftovers="";
718 :
719 0 : eventComments = "";
720 0 : istringstream f(leftovers);
721 0 : string l;
722 0 : while (getline(f, l)) {
723 0 : size_t p = l.find_first_not_of(" \t");
724 0 : l.erase(0, p);
725 0 : p = l.find_last_not_of(" \t");
726 0 : if (string::npos != p) l.erase(p+1);
727 0 : if (l.find_last_not_of("\n") != string::npos)
728 0 : eventComments += l + "\n";
729 : }
730 :
731 0 : for ( int i = 0, N = tags.size(); i < N; ++i ) {
732 0 : XMLTag & tag = *tags[i];
733 :
734 0 : if ( tag.name == "weights" ) {
735 0 : LHAweights wts(tag);
736 0 : eup.weights = wts;
737 :
738 0 : for ( int k = 0, M = int(wts.weights.size()); k < M; ++k ) {
739 0 : eup.weights_compressed.push_back(wts.weights[k]);
740 : }
741 :
742 0 : }
743 0 : else if ( tag.name == "scales" ) {
744 0 : eup.scales = LHAscales(tag, eup.SCALUP);
745 0 : }
746 0 : else if ( tag.name == "rwgt" ) {
747 0 : LHArwgt rwgt(tag);
748 0 : eup.rwgt = rwgt;
749 0 : string s;
750 0 : vector<XMLTag*> tags2 = XMLTag::findXMLTags(rwgt.contents, &s);
751 0 : for ( int k = 0, M = tags2.size(); k < M; ++k ) {
752 0 : const XMLTag & tagnow = *tags2[k];
753 0 : if ( tagnow.name == "wgt" ) {
754 0 : LHAwgt wt(tagnow);
755 0 : eup.weights_detailed.insert(make_pair(wt.id, wt.contents));
756 0 : }
757 : }
758 0 : for ( int k = 0, M = tag.tags.size(); k < M; ++k ) {
759 0 : const XMLTag & tagnow = *tag.tags[k];
760 0 : if ( tagnow.name == "wgt" ) {
761 0 : LHAwgt wt(tagnow);
762 0 : eup.weights_detailed.insert(make_pair(wt.id, wt.contents));
763 0 : }
764 : }
765 0 : }
766 : }
767 :
768 0 : for ( int i = 0, N = tags.size(); i < N; ++i ) if (tags[i]) delete tags[i];
769 :
770 : return true;
771 :
772 0 : }
773 :
774 : //==========================================================================
775 :
776 : // The Writer class is initialized with a stream to which to write a
777 : // version 3.0 Les Houches Accord event file.
778 :
779 : //--------------------------------------------------------------------------
780 :
781 : // Write out an optional header block followed by the standard init
782 : // block information together with any comment lines.
783 :
784 : void Writer::init() {
785 :
786 : // Write out the standard XML tag for the event file.
787 0 : if ( version == 1 )
788 0 : file << "<LesHouchesEvents version=\"1.0\">" << endl;
789 : else
790 0 : file << "<LesHouchesEvents version=\"3.0\">" << endl;
791 :
792 0 : file << setprecision(8);
793 :
794 : // Print headercomments and header init information.
795 0 : file << "<header>" << endl;
796 0 : file << hashline(headerStream.str(),true) << std::flush;
797 0 : if ( version != 1 ) heprup.initrwgt.print(file);
798 0 : file << "</header>" << endl;
799 :
800 0 : file << "<init>"<< endl
801 0 : << " " << setw(8) << heprup.IDBMUP.first
802 0 : << " " << setw(8) << heprup.IDBMUP.second
803 0 : << " " << setw(14) << heprup.EBMUP.first
804 0 : << " " << setw(14) << heprup.EBMUP.second
805 0 : << " " << setw(4) << heprup.PDFGUP.first
806 0 : << " " << setw(4) << heprup.PDFGUP.second
807 0 : << " " << setw(4) << heprup.PDFSUP.first
808 0 : << " " << setw(4) << heprup.PDFSUP.second
809 0 : << " " << setw(4) << heprup.IDWTUP
810 0 : << " " << setw(4) << heprup.NPRUP << endl;
811 0 : heprup.resize();
812 0 : for ( int i = 0; i < heprup.NPRUP; ++i )
813 0 : file << " " << setw(14) << heprup.XSECUP[i]
814 0 : << " " << setw(14) << heprup.XERRUP[i]
815 0 : << " " << setw(14) << heprup.XMAXUP[i]
816 0 : << " " << setw(6) << heprup.LPRUP[i] << endl;
817 :
818 0 : if ( version == 1 ) {
819 0 : file << hashline(initStream.str(),true) << std::flush
820 0 : << "</init>" << endl;
821 0 : initStream.str("");
822 0 : return;
823 : }
824 :
825 0 : for ( int i = 0, N = heprup.generators.size(); i < N; ++i ) {
826 0 : heprup.generators[i].print(file);
827 : }
828 :
829 0 : file << hashline(initStream.str(),true) << std::flush
830 0 : << "</init>" << endl;
831 0 : initStream.str("");
832 0 : }
833 :
834 : //--------------------------------------------------------------------------
835 :
836 : // Write out the event stored in hepeup, followed by optional
837 : // comment lines.
838 :
839 : bool Writer::writeEvent(HEPEUP * peup) {
840 :
841 0 : HEPEUP & eup = (peup? *peup: hepeup);
842 :
843 0 : file << "<event";
844 0 : for ( map<string,string>::const_iterator it = eup.attributes.begin();
845 0 : it != eup.attributes.end(); ++it )
846 0 : file << " " << it->first << "=\"" << it->second << "\"";
847 0 : file << ">" << endl;
848 0 : file << " " << setw(4) << eup.NUP
849 0 : << " " << setw(6) << eup.IDPRUP
850 0 : << " " << setw(14) << eup.XWGTUP
851 0 : << " " << setw(14) << eup.SCALUP
852 0 : << " " << setw(14) << eup.AQEDUP
853 0 : << " " << setw(14) << eup.AQCDUP << endl;
854 0 : eup.resize();
855 :
856 : int pDigits = 18;
857 :
858 0 : for ( int i = 0; i < eup.NUP; ++i )
859 0 : file << " " << setw(8) << eup.IDUP[i]
860 0 : << " " << setw(2) << eup.ISTUP[i]
861 0 : << " " << setw(4) << eup.MOTHUP[i].first
862 0 : << " " << setw(4) << eup.MOTHUP[i].second
863 0 : << " " << setw(4) << eup.ICOLUP[i].first
864 0 : << " " << setw(4) << eup.ICOLUP[i].second
865 0 : << " " << setw(pDigits) << eup.PUP[i][0]
866 0 : << " " << setw(pDigits) << eup.PUP[i][1]
867 0 : << " " << setw(pDigits) << eup.PUP[i][2]
868 0 : << " " << setw(pDigits) << eup.PUP[i][3]
869 0 : << " " << setw(pDigits) << eup.PUP[i][4]
870 0 : << " " << setw(1) << eup.VTIMUP[i]
871 0 : << " " << setw(1) << eup.SPINUP[i] << endl;
872 :
873 : // Write event comments.
874 0 : file << hashline(eventStream.str()) << std::flush;
875 0 : eventStream.str("");
876 :
877 0 : if ( version != 1 ) {
878 0 : eup.rwgt.print(file);
879 0 : eup.weights.print(file);
880 0 : eup.scales.print(file);
881 0 : }
882 :
883 0 : file << "</event>" << endl;
884 :
885 0 : if ( !file ) return false;
886 :
887 0 : return true;
888 :
889 0 : }
890 :
891 : //--------------------------------------------------------------------------
892 :
893 : // Make sure that each line in the string s starts with a
894 : // #-character and that the string ends with a new-line.
895 :
896 : string Writer::hashline(string s, bool comment) {
897 0 : string ret;
898 0 : istringstream is(s);
899 0 : string ss;
900 0 : while ( getline(is, ss) ) {
901 0 : if ( comment )
902 0 : ss = "# " + ss;
903 0 : ret += ss + '\n';
904 : }
905 : return ret;
906 0 : }
907 :
908 : //==========================================================================
909 :
910 : }
|
