Line data Source code
1 : //--------------------------------------------------------------------------
2 : //
3 : // GenEventStreamIO.cc
4 : // Author: Lynn Garren
5 : //
6 : // Implement operator >> and operator <<
7 : //
8 : // ----------------------------------------------------------------------
9 :
10 : #include <iostream>
11 : #include <ostream>
12 : #include <istream>
13 : #include <sstream>
14 :
15 : #include "HepMC/GenEvent.h"
16 : #include "HepMC/GenCrossSection.h"
17 : #include "HepMC/StreamInfo.h"
18 : #include "HepMC/StreamHelpers.h"
19 : #include "HepMC/Version.h"
20 : #include "HepMC/IO_Exception.h"
21 :
22 : namespace HepMC {
23 :
24 : // ------------------------- local methods ----------------
25 :
26 : /// This method is called by the stream destructor.
27 : /// It does cleanup on stored user data (StreamInfo)
28 : /// and is registered by the first call to get_stream_info().
29 : void HepMCStreamCallback(std::ios_base::event e, std::ios_base& b, int i)
30 : {
31 : // only clean up if the stream object is going away.
32 0 : if(i!=0 && e!= std::ios_base::erase_event) return;
33 :
34 : // retrieve the pointer to the object
35 0 : StreamInfo* hd = (StreamInfo*)b.pword(i);
36 0 : b.pword(i) = 0;
37 0 : b.iword(i) = 0;
38 : #ifdef HEPMC_DEBUG
39 : // the following line is just for sanity checking
40 : if(hd) std::cerr << "deleted StreamInfo " << hd->stream_id() << "\n";
41 : #endif
42 0 : delete hd;
43 0 : }
44 :
45 : // ------------------------- iomanip ----------------
46 :
47 : /// A custom iomanip that allows us to store and access user data (StreamInfo)
48 : /// associated with the stream.
49 : /// This method creates the StreamInfo object the first time it is called.
50 : template <class IO>
51 : StreamInfo& get_stream_info(IO& iost)
52 : {
53 0 : if(iost.iword(0) == 0)
54 : {
55 : // make sure we add the callback if this is the first time through
56 0 : iost.iword(0)=1;
57 0 : iost.register_callback(&HepMCStreamCallback, 0);
58 : // this is our special "context" record.
59 : // there is one of these at the head of each IO block.
60 : // allocate room for a StreamInfo in the userdata area
61 0 : iost.pword(0) = new StreamInfo;
62 : #ifdef HEPMC_DEBUG
63 : // the following line is just for sanity checking
64 : std::cerr << "created StreamInfo " << ((StreamInfo*)iost.pword(0))->stream_id() << "\n";
65 : #endif
66 0 : }
67 0 : return *(StreamInfo*)iost.pword(0);
68 0 : }
69 :
70 : // ------------------------- GenEvent member functions ----------------
71 :
72 : std::ostream& GenEvent::write( std::ostream& os )
73 : {
74 : /// Writes evt to an output stream.
75 :
76 : //
77 0 : StreamInfo & info = get_stream_info(os);
78 : //
79 : // if this is the first event, set precision
80 0 : if ( !info.finished_first_event() ) {
81 : // precision 16 (# digits following decimal point) is the minimum that
82 : // will capture the full information stored in a double
83 : // However, we let the user set precision, since that is the expected functionality
84 : // we use decimal to store integers, because it is smaller than hex!
85 0 : os.setf(std::ios::dec,std::ios::basefield);
86 0 : os.setf(std::ios::scientific,std::ios::floatfield);
87 : //
88 0 : info.set_finished_first_event(true);
89 0 : }
90 : //
91 : // output the event data including the number of primary vertices
92 : // and the total number of vertices
93 : //std::vector<long> random_states = random_states();
94 0 : os << 'E';
95 0 : detail::output( os, event_number() );
96 0 : detail::output( os, mpi() );
97 0 : detail::output( os, event_scale() );
98 0 : detail::output( os, alphaQCD() );
99 0 : detail::output( os, alphaQED() );
100 0 : detail::output( os, signal_process_id() );
101 0 : detail::output( os, ( signal_process_vertex() ?
102 0 : signal_process_vertex()->barcode() : 0 ) );
103 0 : detail::output( os, vertices_size() ); // total number of vertices.
104 0 : write_beam_particles( os, beam_particles() );
105 : // random state
106 0 : detail::output( os, (int)m_random_states.size() );
107 0 : for ( std::vector<long>::iterator rs = m_random_states.begin();
108 0 : rs != m_random_states.end(); ++rs ) {
109 0 : detail::output( os, *rs );
110 : }
111 : // weights
112 : // we need to iterate over the map so that the weights printed
113 : // here will be in the same order as the names printed next
114 0 : os << ' ' << (int)weights().size() ;
115 0 : for ( WeightContainer::const_map_iterator w = weights().map_begin();
116 0 : w != weights().map_end(); ++w ) {
117 0 : detail::output( os, m_weights[w->second] );
118 : }
119 0 : detail::output( os,'\n');
120 : // now add names for weights
121 : // note that this prints a new line if and only if the weight container
122 : // is not empty
123 0 : if ( ! weights().empty() ) {
124 0 : os << "N " << weights().size() << " " ;
125 0 : for ( WeightContainer::const_map_iterator w = weights().map_begin();
126 0 : w != weights().map_end(); ++w ) {
127 0 : detail::output( os,'"');
128 0 : os << w->first;
129 0 : detail::output( os,'"');
130 0 : detail::output( os,' ');
131 : }
132 0 : detail::output( os,'\n');
133 0 : }
134 : //
135 : // Units
136 0 : os << "U " << name(momentum_unit());
137 0 : os << " " << name(length_unit());
138 0 : detail::output( os,'\n');
139 : //
140 : // write GenCrossSection if it has been set
141 0 : if( m_cross_section ) m_cross_section->write(os);
142 : //
143 : // write HeavyIon and PdfInfo if they have been set
144 0 : if( m_heavy_ion ) os << heavy_ion() ;
145 0 : if( m_pdf_info ) os << pdf_info() ;
146 : //
147 : // Output all of the vertices - note there is no real order.
148 0 : for ( GenEvent::vertex_const_iterator v = vertices_begin();
149 0 : v != vertices_end(); ++v ) {
150 0 : write_vertex(os, *v);
151 : }
152 0 : return os;
153 0 : }
154 :
155 : std::istream& GenEvent::read( std::istream& is )
156 : {
157 : /// read a GenEvent from streaming input
158 : //
159 0 : StreamInfo & info = get_stream_info(is);
160 0 : clear();
161 : //
162 : // search for event listing key before first event only.
163 0 : if ( !info.finished_first_event() ) {
164 : //
165 0 : find_file_type(is);
166 0 : info.set_finished_first_event(true);
167 0 : }
168 : //
169 : // make sure the stream is good
170 0 : if ( !is ) {
171 0 : std::cerr << "streaming input: end of stream found "
172 0 : << "setting badbit." << std::endl;
173 0 : is.clear(std::ios::badbit);
174 0 : return is;
175 : }
176 :
177 : //
178 : // test to be sure the next entry is of type "E" then ignore it
179 0 : if ( is.peek()!='E' ) {
180 : // if the E is not the next entry, then check to see if it is
181 : // the end event listing key - if yes, search for another start key
182 0 : int ioendtype;
183 0 : find_end_key(is,ioendtype);
184 0 : if ( ioendtype == info.io_type() ) {
185 0 : find_file_type(is);
186 : // are we at the end of the file?
187 0 : if( !is ) return is;
188 0 : } else if ( ioendtype > 0 ) {
189 0 : std::cerr << "streaming input: end key does not match start key "
190 0 : << "setting badbit." << std::endl;
191 0 : is.clear(std::ios::badbit);
192 0 : return is;
193 0 : } else if ( !info.has_key() ) {
194 0 : find_file_type(is);
195 : // are we at the end of the file?
196 0 : if( !is ) return is;
197 : } else {
198 0 : std::cerr << "streaming input: end key not found "
199 0 : << "setting badbit." << std::endl;
200 0 : is.clear(std::ios::badbit);
201 0 : return is;
202 : }
203 0 : }
204 :
205 0 : int signal_process_vertex = 0;
206 0 : int num_vertices = 0, bp1 = 0, bp2 = 0;
207 : bool units_line = false;
208 : // OK - now ready to start reading the event, so set the header flag
209 0 : info.set_reading_event_header(true);
210 : // The flag will be set to false when we reach the end of the header
211 0 : while(info.reading_event_header()) {
212 0 : switch(is.peek()) {
213 : case 'E':
214 : { // deal with the event line
215 0 : process_event_line( is, num_vertices, bp1, bp2, signal_process_vertex );
216 0 : } break;
217 : case 'N':
218 : { // get weight names
219 0 : read_weight_names( is );
220 0 : } break;
221 : case 'U':
222 : { // get unit information if it exists
223 : units_line = true;
224 0 : if( info.io_type() == gen ) {
225 0 : read_units( is );
226 0 : }
227 : } break;
228 : case 'C':
229 : { // we have a GenCrossSection line
230 : // create cross section
231 0 : GenCrossSection xs;
232 : // check for invalid data
233 : try {
234 : // read the line
235 0 : xs.read(is);
236 0 : }
237 : catch (IO_Exception& e) {
238 0 : detail::find_event_end( is );
239 0 : }
240 0 : if(xs.is_set()) {
241 0 : set_cross_section( xs );
242 : }
243 0 : } break;
244 : case 'H':
245 : { // we have a HeavyIon line OR an unexpected HepMC... line
246 0 : if( info.io_type() == gen || info.io_type() == extascii ) {
247 : // get HeavyIon
248 0 : HeavyIon ion;
249 : // check for invalid data
250 : try {
251 0 : is >> &ion;
252 0 : }
253 : catch (IO_Exception& e) {
254 0 : detail::find_event_end( is );
255 0 : }
256 0 : if(ion.is_valid()) {
257 0 : set_heavy_ion( ion );
258 : }
259 0 : }
260 : } break;
261 : case 'F':
262 : { // we have a PdfInfo line
263 0 : if( info.io_type() == gen || info.io_type() == extascii ) {
264 : // get PdfInfo
265 0 : PdfInfo pdf;
266 : // check for invalid data
267 : try {
268 0 : is >> &pdf;
269 0 : }
270 : catch (IO_Exception& e) {
271 0 : detail::find_event_end( is );
272 0 : }
273 0 : if(pdf.is_valid()) {
274 0 : set_pdf_info( pdf );
275 : }
276 0 : }
277 : } break;
278 : case 'V':
279 : {
280 : // this should be the first vertex line - exit this loop
281 0 : info.set_reading_event_header(false);
282 0 : } break;
283 : case 'P':
284 : { // we should not find this line
285 0 : std::cerr << "streaming input: found unexpected line P" << std::endl;
286 0 : info.set_reading_event_header(false);
287 0 : } break;
288 : default:
289 : // ignore everything else
290 : break;
291 : } // switch on line type
292 : } // while reading_event_header
293 : // before proceeding - did we find a units line?
294 0 : if( !units_line ) {
295 0 : use_units( info.io_momentum_unit(),
296 0 : info.io_position_unit() );
297 0 : }
298 : //
299 : // the end vertices of the particles are not connected until
300 : // after the event is read --- we store the values in a map until then
301 0 : TempParticleMap particle_to_end_vertex;
302 : //
303 : // read in the vertices
304 0 : for ( int iii = 1; iii <= num_vertices; ++iii ) {
305 0 : GenVertex* v = new GenVertex();
306 : try {
307 0 : detail::read_vertex(is,particle_to_end_vertex,v);
308 0 : }
309 : catch (IO_Exception& e) {
310 0 : for( TempParticleMap::orderIterator it = particle_to_end_vertex.order_begin();
311 0 : it != particle_to_end_vertex.order_end(); ++it ) {
312 0 : GenParticle* p = it->second;
313 : // delete particles only if they are not already owned by a vertex
314 0 : if( p->production_vertex() ) {
315 0 : } else if( p->end_vertex() ) {
316 : } else {
317 0 : delete p;
318 : }
319 : }
320 0 : delete v;
321 0 : detail::find_event_end( is );
322 0 : }
323 0 : add_vertex( v );
324 : }
325 : // set the signal process vertex
326 0 : if ( signal_process_vertex ) {
327 0 : set_signal_process_vertex(
328 0 : barcode_to_vertex(signal_process_vertex) );
329 : }
330 : //
331 : // last connect particles to their end vertices
332 : GenParticle* beam1(0);
333 : GenParticle* beam2(0);
334 0 : for ( TempParticleMap::orderIterator pmap
335 0 : = particle_to_end_vertex.order_begin();
336 0 : pmap != particle_to_end_vertex.order_end(); ++pmap ) {
337 0 : GenParticle* p = pmap->second;
338 0 : int vtx = particle_to_end_vertex.end_vertex( p );
339 0 : GenVertex* itsDecayVtx = barcode_to_vertex(vtx);
340 0 : if ( itsDecayVtx ) itsDecayVtx->add_particle_in( p );
341 : else {
342 0 : std::cerr << "read_io_genevent: ERROR particle points"
343 0 : << " to null end vertex. " <<std::endl;
344 : }
345 : // also look for the beam particles
346 0 : if( p->barcode() == bp1 ) beam1 = p;
347 0 : if( p->barcode() == bp2 ) beam2 = p;
348 : }
349 0 : set_beam_particles(beam1,beam2);
350 : return is;
351 0 : }
352 :
353 : // ------------------------- operator << and operator >> ----------------
354 :
355 : std::ostream & operator << (std::ostream & os, GenEvent & evt)
356 : {
357 : /// Writes evt to an output stream.
358 0 : evt.write(os);
359 0 : return os;
360 : }
361 :
362 : std::istream & operator >> (std::istream & is, GenEvent & evt)
363 : {
364 0 : evt.read(is);
365 0 : return is;
366 : }
367 :
368 : // ------------------------- set units ----------------
369 :
370 : std::istream & set_input_units(std::istream & is,
371 : Units::MomentumUnit mom,
372 : Units::LengthUnit len )
373 : {
374 : //
375 0 : StreamInfo & info = get_stream_info(is);
376 0 : info.use_input_units( mom, len );
377 0 : return is;
378 : }
379 :
380 : // ------------------------- begin and end block lines ----------------
381 :
382 : std::ostream & write_HepMC_IO_block_begin(std::ostream & os )
383 : {
384 : //
385 0 : StreamInfo & info = get_stream_info(os);
386 :
387 0 : if( !info.finished_first_event() ) {
388 0 : os << "\n" << "HepMC::Version " << versionName();
389 0 : os << "\n";
390 0 : os << info.IO_GenEvent_Key() << "\n";
391 0 : }
392 0 : return os;
393 0 : }
394 :
395 : std::ostream & write_HepMC_IO_block_end(std::ostream & os )
396 : {
397 : //
398 0 : StreamInfo & info = get_stream_info(os);
399 :
400 0 : if( info.finished_first_event() ) {
401 0 : os << info.IO_GenEvent_End() << "\n";
402 0 : os << std::flush;
403 0 : }
404 0 : return os;
405 0 : }
406 :
407 : std::istream & GenEvent::process_event_line( std::istream & is,
408 : int & num_vertices,
409 : int & bp1, int & bp2,
410 : int & signal_process_vertex )
411 : {
412 : //
413 0 : if ( !is ) {
414 0 : std::cerr << "GenEvent::process_event_line setting badbit." << std::endl;
415 0 : is.clear(std::ios::badbit);
416 0 : return is;
417 : }
418 : //
419 0 : StreamInfo & info = get_stream_info(is);
420 0 : std::string line;
421 0 : std::getline(is,line);
422 0 : std::istringstream iline(line);
423 0 : std::string firstc;
424 0 : iline >> firstc;
425 : //
426 : // read values into temp variables, then fill GenEvent
427 0 : int event_number = 0, signal_process_id = 0,
428 0 : random_states_size = 0, nmpi = -1;
429 0 : double eventScale = 0, alpha_qcd = 0, alpha_qed = 0;
430 0 : iline >> event_number;
431 0 : if(!iline) detail::find_event_end( is );
432 0 : if( info.io_type() == gen || info.io_type() == extascii ) {
433 0 : iline >> nmpi;
434 0 : if(!iline) detail::find_event_end( is );
435 0 : set_mpi( nmpi );
436 0 : }
437 0 : iline >> eventScale ;
438 0 : if(!iline) detail::find_event_end( is );
439 0 : iline >> alpha_qcd ;
440 0 : if(!iline) detail::find_event_end( is );
441 0 : iline >> alpha_qed;
442 0 : if(!iline) detail::find_event_end( is );
443 0 : iline >> signal_process_id ;
444 0 : if(!iline) detail::find_event_end( is );
445 0 : iline >> signal_process_vertex;
446 0 : if(!iline) detail::find_event_end( is );
447 0 : iline >> num_vertices;
448 0 : if(!iline) detail::find_event_end( is );
449 0 : if( info.io_type() == gen || info.io_type() == extascii ) {
450 0 : iline >> bp1 ;
451 0 : if(!iline) detail::find_event_end( is );
452 0 : iline >> bp2;
453 0 : if(!iline) detail::find_event_end( is );
454 : }
455 0 : iline >> random_states_size;
456 0 : if(!iline) detail::find_event_end( is );
457 0 : std::vector<long> random_states(random_states_size);
458 0 : for ( int i = 0; i < random_states_size; ++i ) {
459 0 : iline >> random_states[i];
460 0 : if(!iline) detail::find_event_end( is );
461 : }
462 0 : WeightContainer::size_type weights_size = 0;
463 0 : iline >> weights_size;
464 0 : if(!iline) detail::find_event_end( is );
465 0 : std::vector<double> wgt(weights_size);
466 0 : for ( WeightContainer::size_type ii = 0; ii < weights_size; ++ii ) {
467 0 : iline >> wgt[ii];
468 0 : if(!iline) detail::find_event_end( is );
469 : }
470 : // weight names will be added later if they exist
471 0 : if( weights_size > 0 ) m_weights = wgt;
472 : //
473 : // fill signal_process_id, event_number, random_states, etc.
474 0 : set_signal_process_id( signal_process_id );
475 0 : set_event_number( event_number );
476 0 : set_random_states( random_states );
477 0 : set_event_scale( eventScale );
478 0 : set_alphaQCD( alpha_qcd );
479 0 : set_alphaQED( alpha_qed );
480 : //
481 : return is;
482 0 : }
483 :
484 : std::istream & GenEvent::read_weight_names( std::istream & is )
485 : {
486 : // now check for a named weight line
487 0 : if ( !is ) {
488 0 : std::cerr << "GenEvent::read_weight_names setting badbit." << std::endl;
489 0 : is.clear(std::ios::badbit);
490 0 : return is;
491 : }
492 : // Test to be sure the next entry is of type "N"
493 : // If we have no named weight line, this is not an error
494 : // releases prior to 2.06.00 do not have named weights
495 0 : if ( is.peek() !='N') {
496 0 : return is;
497 : }
498 : // now get this line and process it
499 0 : std::string line;
500 0 : std::getline(is,line);
501 0 : std::istringstream wline(line);
502 0 : std::string firstc;
503 0 : WeightContainer::size_type name_size = 0;
504 0 : wline >> firstc >> name_size;
505 0 : if(!wline) detail::find_event_end( is );
506 0 : if( firstc != "N") {
507 0 : std::cout << "debug: first character of named weights is " << firstc << std::endl;
508 0 : std::cout << "debug: We should never get here" << std::endl;
509 0 : is.clear(std::ios::badbit);
510 0 : return is;
511 : }
512 0 : if( m_weights.size() != name_size ) {
513 0 : std::cout << "debug: weight sizes do not match "<< std::endl;
514 0 : std::cout << "debug: weight vector size is " << m_weights.size() << std::endl;
515 0 : std::cout << "debug: weight name size is " << name_size << std::endl;
516 0 : is.clear(std::ios::badbit);
517 0 : return is;
518 : }
519 0 : std::string name;
520 0 : std::string::size_type i1 = line.find("\"");
521 : std::string::size_type i2;
522 0 : std::string::size_type len = line.size();
523 0 : WeightContainer namedWeight;
524 0 : for ( WeightContainer::size_type ii = 0; ii < name_size; ++ii ) {
525 : // weight names may contain blanks
526 0 : if(i1 >= len) {
527 0 : std::cout << "debug: attempting to read past the end of the named weight line " << std::endl;
528 0 : std::cout << "debug: We should never get here" << std::endl;
529 0 : std::cout << "debug: Looking for the end of this event" << std::endl;
530 0 : detail::find_event_end( is );
531 : }
532 0 : i2 = line.find("\"",i1+1);
533 0 : name = line.substr(i1+1,i2-i1-1);
534 0 : namedWeight[name] = m_weights[ii];
535 0 : i1 = line.find("\"",i2+1);
536 : }
537 0 : m_weights = namedWeight;
538 : return is;
539 0 : }
540 :
541 : std::istream & GenEvent::read_units( std::istream & is )
542 : {
543 : //
544 0 : if ( !is ) {
545 0 : std::cerr << "GenEvent::read_units setting badbit." << std::endl;
546 0 : is.clear(std::ios::badbit);
547 0 : return is;
548 : }
549 : //
550 0 : StreamInfo & info = get_stream_info(is);
551 : // test to be sure the next entry is of type "U" then ignore it
552 : // if we have no units, this is not an error
553 : // releases prior to 2.04.00 did not write unit information
554 0 : if ( is.peek() !='U') {
555 0 : use_units( info.io_momentum_unit(),
556 0 : info.io_position_unit() );
557 0 : return is;
558 : }
559 0 : is.ignore(); // ignore the first character in the line
560 0 : std::string mom, pos;
561 0 : is >> mom >> pos;
562 0 : is.ignore(1); // eat the extra whitespace
563 0 : use_units(mom,pos);
564 : //
565 : return is;
566 0 : }
567 :
568 : std::istream & GenEvent::find_file_type( std::istream & istr )
569 : {
570 : //
571 : // make sure the stream is good
572 0 : if ( !istr ) return istr;
573 :
574 : //
575 0 : StreamInfo & info = get_stream_info(istr);
576 :
577 : // if there is no input block line, then we assume this stream
578 : // is in the IO_GenEvent format
579 0 : if ( istr.peek()=='E' ) {
580 0 : info.set_io_type( gen );
581 0 : info.set_has_key(false);
582 0 : return istr;
583 : }
584 :
585 0 : std::string line;
586 0 : while ( std::getline(istr,line) ) {
587 : //
588 : // search for event listing key before first event only.
589 : //
590 0 : if( line == info.IO_GenEvent_Key() ) {
591 0 : info.set_io_type( gen );
592 0 : info.set_has_key(true);
593 0 : return istr;
594 0 : } else if( line == info.IO_Ascii_Key() ) {
595 0 : info.set_io_type( ascii );
596 0 : info.set_has_key(true);
597 0 : return istr;
598 0 : } else if( line == info.IO_ExtendedAscii_Key() ) {
599 0 : info.set_io_type( extascii );
600 0 : info.set_has_key(true);
601 0 : return istr;
602 0 : } else if( line == info.IO_Ascii_PDT_Key() ) {
603 0 : info.set_io_type( ascii_pdt );
604 0 : info.set_has_key(true);
605 0 : return istr;
606 0 : } else if( line == info.IO_ExtendedAscii_PDT_Key() ) {
607 0 : info.set_io_type( extascii_pdt );
608 0 : info.set_has_key(true);
609 0 : return istr;
610 : }
611 : }
612 0 : info.set_io_type( 0 );
613 0 : info.set_has_key(false);
614 0 : return istr;
615 0 : }
616 :
617 : std::istream & GenEvent::find_end_key( std::istream & istr, int & iotype )
618 : {
619 0 : iotype = 0;
620 : // peek at the first character before proceeding
621 0 : if( istr.peek()!='H' ) return istr;
622 : //
623 : // we only check the next line
624 0 : std::string line;
625 0 : std::getline(istr,line);
626 : //
627 0 : StreamInfo & info = get_stream_info(istr);
628 : //
629 : // check to see if this is an end key
630 0 : if( line == info.IO_GenEvent_End() ) {
631 0 : iotype = gen;
632 0 : } else if( line == info.IO_Ascii_End() ) {
633 0 : iotype = ascii;
634 0 : } else if( line == info.IO_ExtendedAscii_End() ) {
635 0 : iotype = extascii;
636 0 : } else if( line == info.IO_Ascii_PDT_End() ) {
637 0 : iotype = ascii_pdt;
638 0 : } else if( line == info.IO_ExtendedAscii_PDT_End() ) {
639 0 : iotype = extascii_pdt;
640 0 : }
641 0 : if( iotype != 0 && info.io_type() != iotype ) {
642 0 : std::cerr << "GenEvent::find_end_key: iotype keys have changed" << std::endl;
643 : } else {
644 0 : return istr;
645 : }
646 : //
647 : // if we get here, then something has gotten badly confused
648 0 : std::cerr << "GenEvent::find_end_key: MALFORMED INPUT" << std::endl;
649 0 : istr.clear(std::ios::badbit);
650 0 : return istr;
651 0 : }
652 :
653 : std::ostream & establish_output_stream_info( std::ostream & os )
654 : {
655 0 : StreamInfo & info = get_stream_info(os);
656 0 : if ( !info.finished_first_event() ) {
657 : // precision 16 (# digits following decimal point) is the minimum that
658 : // will capture the full information stored in a double
659 0 : os.precision(16);
660 : // we use decimal to store integers, because it is smaller than hex!
661 0 : os.setf(std::ios::dec,std::ios::basefield);
662 0 : os.setf(std::ios::scientific,std::ios::floatfield);
663 0 : }
664 0 : return os;
665 : }
666 :
667 : std::istream & establish_input_stream_info( std::istream & is )
668 : {
669 0 : StreamInfo & info = get_stream_info(is);
670 0 : if ( !info.finished_first_event() ) {
671 : // precision 16 (# digits following decimal point) is the minimum that
672 : // will capture the full information stored in a double
673 0 : is.precision(16);
674 : // we use decimal to store integers, because it is smaller than hex!
675 0 : is.setf(std::ios::dec,std::ios::basefield);
676 0 : is.setf(std::ios::scientific,std::ios::floatfield);
677 0 : }
678 0 : return is;
679 : }
680 :
681 :
682 : // ------------------------- helper functions ----------------
683 :
684 : namespace detail {
685 :
686 : // The functions defined here need to use get_stream_info
687 :
688 : std::istream & read_particle( std::istream & is,
689 : TempParticleMap & particle_to_end_vertex,
690 : GenParticle * p )
691 : {
692 : // get the next line
693 0 : std::string line;
694 0 : std::getline(is,line);
695 0 : std::istringstream iline(line);
696 0 : std::string firstc;
697 0 : iline >> firstc;
698 0 : if( firstc != "P" ) {
699 0 : std::cerr << "StreamHelpers::detail::read_particle invalid line type: "
700 0 : << firstc << std::endl;
701 0 : std::cerr << "StreamHelpers::detail::read_particle setting badbit."
702 0 : << std::endl;
703 0 : is.clear(std::ios::badbit);
704 0 : return is;
705 : }
706 : //
707 0 : StreamInfo & info = get_stream_info(is);
708 : //testHepMC.cc
709 : // declare variables to be read in to, and read everything except flow
710 0 : double px = 0., py = 0., pz = 0., e = 0., m = 0., theta = 0., phi = 0.;
711 0 : int bar_code = 0, id = 0, status = 0, end_vtx_code = 0, flow_size = 0;
712 : // check that the input stream is still OK after reading item
713 0 : iline >> bar_code ;
714 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
715 0 : iline >> id ;
716 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
717 0 : iline >> px ;
718 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
719 0 : iline >> py ;
720 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
721 0 : iline >> pz ;
722 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
723 0 : iline >> e ;
724 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
725 0 : if( info.io_type() != ascii ) {
726 0 : iline >> m ;
727 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
728 : }
729 0 : iline >> status ;
730 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
731 0 : iline >> theta ;
732 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
733 0 : iline >> phi ;
734 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
735 0 : iline >> end_vtx_code ;
736 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
737 0 : iline >> flow_size;
738 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
739 : //
740 : // read flow patterns if any exist
741 0 : Flow flow;
742 0 : int code_index, code;
743 0 : for ( int i = 1; i <= flow_size; ++i ) {
744 0 : iline >> code_index >> code;
745 0 : if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
746 0 : flow.set_icode( code_index,code);
747 : }
748 0 : p->set_momentum( FourVector(px,py,pz,e) );
749 0 : p->set_pdg_id( id );
750 0 : p->set_status( status );
751 0 : p->set_flow( flow );
752 0 : p->set_polarization( Polarization(theta,phi) );
753 0 : if( info.io_type() == ascii ) {
754 0 : p->set_generated_mass( p->momentum().m() );
755 : } else {
756 0 : p->set_generated_mass( m );
757 : }
758 0 : p->suggest_barcode( bar_code );
759 : //
760 : // all particles are connected to their end vertex separately
761 : // after all particles and vertices have been created - so we keep
762 : // a map of all particles that have end vertices
763 0 : if ( end_vtx_code != 0 ) {
764 0 : particle_to_end_vertex.addEndParticle(p,end_vtx_code);
765 : }
766 : return is;
767 0 : }
768 :
769 : std::ostream & establish_output_stream_info( std::ostream & os )
770 : {
771 0 : StreamInfo & info = get_stream_info(os);
772 0 : if ( !info.finished_first_event() ) {
773 : // precision 16 (# digits following decimal point) is the minimum that
774 : // will capture the full information stored in a double
775 0 : os.precision(16);
776 : // we use decimal to store integers, because it is smaller than hex!
777 0 : os.setf(std::ios::dec,std::ios::basefield);
778 0 : os.setf(std::ios::scientific,std::ios::floatfield);
779 0 : }
780 0 : return os;
781 : }
782 :
783 : std::istream & establish_input_stream_info( std::istream & is )
784 : {
785 0 : StreamInfo & info = get_stream_info(is);
786 0 : if ( !info.finished_first_event() ) {
787 : // precision 16 (# digits following decimal point) is the minimum that
788 : // will capture the full information stored in a double
789 0 : is.precision(16);
790 : // we use decimal to store integers, because it is smaller than hex!
791 0 : is.setf(std::ios::dec,std::ios::basefield);
792 0 : is.setf(std::ios::scientific,std::ios::floatfield);
793 0 : }
794 0 : return is;
795 : }
796 :
797 : } // detail
798 :
799 : } // HepMC
|
