hepmc - Blame information for rev 523

Subversion Repositories:
Rev:
Rev Author Line No. Line
390 garren 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   if(i!=0 && e!= std::ios_base::erase_event) return;
33  
34   // retrieve the pointer to the object
35   StreamInfo* hd = (StreamInfo*)b.pword(i);
36   b.pword(i) = 0;
37   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   delete hd;
43 }
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   if(iost.iword(0) == 0)
54     {
55       // make sure we add the callback if this is the first time through
56       iost.iword(0)=1;
57       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       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     }
67   return *(StreamInfo*)iost.pword(0);
68 }
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     StreamInfo & info = get_stream_info(os);
78     //
79     // if this is the first event, set precision
80     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         os.setf(std::ios::dec,std::ios::basefield);
86         os.setf(std::ios::scientific,std::ios::floatfield);
87         //
88         info.set_finished_first_event(true);
89     }
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     os << 'E';
95     detail::output( os, event_number() );
96     detail::output( os, mpi() );
97     detail::output( os, event_scale() );
98     detail::output( os, alphaQCD() );
99     detail::output( os, alphaQED() );
100     detail::output( os, signal_process_id() );
101     detail::output( os,   ( signal_process_vertex() ?
102                 signal_process_vertex()->barcode() : 0 )   );
103     detail::output( os, vertices_size() ); // total number of vertices.
104     write_beam_particles( os, beam_particles() );
105     // random state
106     detail::output( os, (int)m_random_states.size() );
107     for ( std::vector<long>::iterator rs = m_random_states.begin();
108           rs != m_random_states.end(); ++rs ) {
109          detail::output( os, *rs );
110     }
111     // weights
432 garren 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
390 garren 114     os << ' ' << (int)weights().size() ;
432 garren 115     for ( WeightContainer::const_map_iterator w = weights().map_begin();
116           w != weights().map_end(); ++w ) {
117         detail::output( os, m_weights[w->second] );
390 garren 118     }
119     detail::output( os,'\n');
432 garren 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     if ( ! weights().empty() ) {
124         os << "N " << weights().size() << " " ;
125         for ( WeightContainer::const_map_iterator w = weights().map_begin();
126               w != weights().map_end(); ++w ) {
127             detail::output( os,'"');
128             os << w->first;
129             detail::output( os,'"');
130             detail::output( os,' ');
131         }
132         detail::output( os,'\n');
133     }
390 garren 134     //
135     // Units
136     os << "U " << name(momentum_unit());
137     os << " " << name(length_unit());
138     detail::output( os,'\n');
139     //
140     // write GenCrossSection if it has been set
141     if( m_cross_section ) m_cross_section->write(os);
142     //
432 garren 143     // write HeavyIon and PdfInfo if they have been set
144     if( m_heavy_ion ) os << heavy_ion() ;
145     if( m_pdf_info ) os << pdf_info() ;
390 garren 146     //
147     // Output all of the vertices - note there is no real order.
148     for ( GenEvent::vertex_const_iterator v = vertices_begin();
149           v != vertices_end(); ++v ) {
150         write_vertex(os, *v);
151     }
152     return os;
153 }
154  
155 std::istream& GenEvent::read( std::istream& is )
156 {
157     /// read a GenEvent from streaming input
158     //
159     StreamInfo & info = get_stream_info(is);
160     clear();
161     //
162     // search for event listing key before first event only.
163     if ( !info.finished_first_event() ) {
164         //
165         find_file_type(is);
166         info.set_finished_first_event(true);
167     }
168     //
169     // make sure the stream is good
170     if ( !is ) {
171         std::cerr << "streaming input: end of stream found "
172                   << "setting badbit." << std::endl;
173         is.clear(std::ios::badbit);
174         return is;
175     }
176  
177     //
178     // test to be sure the next entry is of type "E" then ignore it
179     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         int ioendtype;
183         find_end_key(is,ioendtype);
184         if ( ioendtype == info.io_type() ) {
185             find_file_type(is);
186             // are we at the end of the file?
187             if( !is ) return is;
188         } else if ( ioendtype > 0 ) {
189             std::cerr << "streaming input: end key does not match start key "
190                       << "setting badbit." << std::endl;
191             is.clear(std::ios::badbit);
192             return is;
193         } else if ( !info.has_key() ) {
194             find_file_type(is);
195             // are we at the end of the file?
196             if( !is ) return is;
197         } else {
198             std::cerr << "streaming input: end key not found "
199                       << "setting badbit." << std::endl;
200             is.clear(std::ios::badbit);
201             return is;
202         }
203     }
432 garren 204  
205     int signal_process_vertex = 0;
206     int num_vertices = 0, bp1 = 0, bp2 = 0;
523 garren 207     bool units_line = false;
432 garren 208     // OK - now ready to start reading the event, so set the header flag
209     info.set_reading_event_header(true);
210     // The flag will be set to false when we reach the end of the header
211     while(info.reading_event_header()) {
212         switch(is.peek()) {
213             case 'E':
214             {   // deal with the event line
215                 process_event_line( is, num_vertices, bp1, bp2, signal_process_vertex );
216             } break;
217             case 'N':
218             {   // get weight names
219                 read_weight_names( is );
220             } break;
221             case 'U':
222             {   // get unit information if it exists
523 garren 223                 units_line = true;
432 garren 224                 if( info.io_type() == gen ) {
225                     read_units( is );
226                 }
227             } break;
228             case 'C':
229             {   // we have a GenCrossSection line
230                 // create cross section
440 garren 231                 GenCrossSection xs;
432 garren 232                 // check for invalid data
233                 try {
234                     // read the line
440 garren 235                     xs.read(is);
432 garren 236                 }
237                 catch (IO_Exception& e) {
238                     detail::find_event_end( is );
239                 }
440 garren 240                 if(xs.is_set()) {
241                     set_cross_section( xs );
432 garren 242                 }
243             } break;
244             case 'H':
455 garren 245             {   // we have a HeavyIon line OR an unexpected HepMC... line
432 garren 246                 if( info.io_type() == gen || info.io_type() == extascii ) {
247                     // get HeavyIon
440 garren 248                     HeavyIon ion;
432 garren 249                     // check for invalid data
250                     try {
440 garren 251                         is >> &ion;
432 garren 252                     }
253                     catch (IO_Exception& e) {
254                         detail::find_event_end( is );
255                     }
440 garren 256                     if(ion.is_valid()) {
257                         set_heavy_ion( ion );
432 garren 258                     }
259                 }
260             } break;
261             case 'F':
262             {   // we have a PdfInfo line
263                 if( info.io_type() == gen || info.io_type() == extascii ) {
264                     // get PdfInfo
440 garren 265                     PdfInfo pdf;
432 garren 266                     // check for invalid data
267                     try {
440 garren 268                         is >> &pdf;
432 garren 269                     }
270                     catch (IO_Exception& e) {
271                         detail::find_event_end( is );
272                     }
440 garren 273                     if(pdf.is_valid()) {
274                         set_pdf_info( pdf );
432 garren 275                     }
276                 }
277             } break;
278             case 'V':
279             {  
280                 // this should be the first vertex line - exit this loop
281                 info.set_reading_event_header(false);
282             } break;
283             case 'P':
284             {   // we should not find this line
285                 std::cerr << "streaming input: found unexpected line P" << std::endl;
286                 info.set_reading_event_header(false);
287             } break;
288             default:
289                 // ignore everything else
290                 break;
291         } // switch on line type
292     } // while reading_event_header
523 garren 293     // before proceeding - did we find a units line?
294     if( !units_line ) {
295         use_units( info.io_momentum_unit(),
296                        info.io_position_unit() );
297     }
390 garren 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     TempParticleMap particle_to_end_vertex;
302     //
303     // read in the vertices
304     for ( int iii = 1; iii <= num_vertices; ++iii ) {
305         GenVertex* v = new GenVertex();
443 garren 306         try {
307             detail::read_vertex(is,particle_to_end_vertex,v);
308         }
309         catch (IO_Exception& e) {
444 garren 310             for( TempParticleMap::orderIterator it = particle_to_end_vertex.order_begin();
311                  it != particle_to_end_vertex.order_end(); ++it ) {
312                 GenParticle* p = it->second;
313                 // delete particles only if they are not already owned by a vertex
314                 if( p->production_vertex() ) {
315                 } else if( p->end_vertex() ) {
316                 } else {
317                      delete p;
318                 }
319             }
443 garren 320             delete v;
321             detail::find_event_end( is );
322         }
390 garren 323         add_vertex( v );
324     }
325     // set the signal process vertex
326     if ( signal_process_vertex ) {
327         set_signal_process_vertex(
328             barcode_to_vertex(signal_process_vertex) );
329     }
330     //
331     // last connect particles to their end vertices
332     GenParticle* beam1(0);
333     GenParticle* beam2(0);
444 garren 334     for ( TempParticleMap::orderIterator pmap
390 garren 335               = particle_to_end_vertex.order_begin();
336           pmap != particle_to_end_vertex.order_end(); ++pmap ) {
337         GenParticle* p =  pmap->second;
338         int vtx = particle_to_end_vertex.end_vertex( p );
339         GenVertex* itsDecayVtx = barcode_to_vertex(vtx);
340         if ( itsDecayVtx ) itsDecayVtx->add_particle_in( p );
341         else {
342             std::cerr << "read_io_genevent: ERROR particle points"
343                       << " to null end vertex. " <<std::endl;
344         }
345         // also look for the beam particles
346         if( p->barcode() == bp1 ) beam1 = p;
347         if( p->barcode() == bp2 ) beam2 = p;
348     }
349     set_beam_particles(beam1,beam2);
350     return is;
351 }
352  
353 // ------------------------- operator << and operator >> ----------------
354  
355 std::ostream & operator << (std::ostream & os, GenEvent & evt)
356 {
357     /// Writes evt to an output stream.
358     evt.write(os);
359     return os;
360 }
361  
362 std::istream & operator >> (std::istream & is, GenEvent & evt)
363 {
364     evt.read(is);
365     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     StreamInfo & info = get_stream_info(is);
376     info.use_input_units( mom, len );
377     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     StreamInfo & info = get_stream_info(os);
386  
387     if( !info.finished_first_event() ) {
388     os << "\n" << "HepMC::Version " << versionName();
389     os << "\n";
390     os << info.IO_GenEvent_Key() << "\n";
391     }
392     return os;
393 }
394  
395 std::ostream & write_HepMC_IO_block_end(std::ostream & os )
396 {
397     //
398     StreamInfo & info = get_stream_info(os);
399  
400     if( info.finished_first_event() ) {
401         os << info.IO_GenEvent_End() << "\n";
402         os << std::flush;
403     }
404     return os;
405 }
406  
432 garren 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     if ( !is ) {
414         std::cerr << "GenEvent::process_event_line setting badbit." << std::endl;
415         is.clear(std::ios::badbit);
416         return is;
417     }
418     //
419     StreamInfo & info = get_stream_info(is);
420     std::string line;
421     std::getline(is,line);
422     std::istringstream iline(line);
423     std::string firstc;
424     iline >> firstc;
425     //
426     // read values into temp variables, then fill GenEvent
427     int event_number = 0, signal_process_id = 0,
428         random_states_size = 0, nmpi = -1;
429     double eventScale = 0, alpha_qcd = 0, alpha_qed = 0;
430     iline >> event_number;
431     if(!iline) detail::find_event_end( is );
432     if( info.io_type() == gen || info.io_type() == extascii ) {
433         iline >> nmpi;
434         if(!iline) detail::find_event_end( is );
435         set_mpi( nmpi );
436     }
437     iline >> eventScale ;
438     if(!iline) detail::find_event_end( is );
439     iline >> alpha_qcd ;
440     if(!iline) detail::find_event_end( is );
441     iline >> alpha_qed;
442     if(!iline) detail::find_event_end( is );
443     iline >> signal_process_id ;
444     if(!iline) detail::find_event_end( is );
445     iline >> signal_process_vertex;
446     if(!iline) detail::find_event_end( is );
447     iline >> num_vertices;
448     if(!iline) detail::find_event_end( is );
449     if( info.io_type() == gen || info.io_type() == extascii ) {
450         iline >> bp1 ;
451         if(!iline) detail::find_event_end( is );
452         iline >> bp2;
453         if(!iline) detail::find_event_end( is );
454     }
455     iline >> random_states_size;
456     if(!iline) detail::find_event_end( is );
457     std::vector<long> random_states(random_states_size);
458     for ( int i = 0; i < random_states_size; ++i ) {
459         iline >> random_states[i];
460         if(!iline) detail::find_event_end( is );
461     }
462     WeightContainer::size_type weights_size = 0;
463     iline >> weights_size;
464     if(!iline) detail::find_event_end( is );
465     std::vector<double> wgt(weights_size);
466     for ( WeightContainer::size_type ii = 0; ii < weights_size; ++ii ) {
467         iline >> wgt[ii];
468         if(!iline) detail::find_event_end( is );
469     }
470     // weight names will be added later if they exist
471     if( weights_size > 0 ) m_weights = wgt;
472     //
473     // fill signal_process_id, event_number, random_states, etc.
474     set_signal_process_id( signal_process_id );
475     set_event_number( event_number );
476     set_random_states( random_states );
477     set_event_scale( eventScale );
478     set_alphaQCD( alpha_qcd );
479     set_alphaQED( alpha_qed );
480     //
481     return is;
482 }
483  
484 std::istream & GenEvent::read_weight_names( std::istream & is )
485 {
486     // now check for a named weight line
487     if ( !is ) {
488         std::cerr << "GenEvent::read_weight_names setting badbit." << std::endl;
489         is.clear(std::ios::badbit);
490         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     if ( is.peek() !='N') {
496         return is;
497     }
498     // now get this line and process it
499     std::string line;
500     std::getline(is,line);
501     std::istringstream wline(line);
502     std::string firstc;
503     WeightContainer::size_type name_size = 0;
504     wline >> firstc >> name_size;
505     if(!wline) detail::find_event_end( is );
506     if( firstc != "N") {
507         std::cout << "debug: first character of named weights is " << firstc << std::endl;
508         std::cout << "debug: We should never get here" << std::endl;
509         is.clear(std::ios::badbit);
510         return is;
511     }
512     if( m_weights.size() != name_size ) {
513         std::cout << "debug: weight sizes do not match "<< std::endl;
514         std::cout << "debug: weight vector size is " << m_weights.size() << std::endl;
515         std::cout << "debug: weight name size is " << name_size << std::endl;
516         is.clear(std::ios::badbit);
517         return is;
518     }
519     std::string name;
520     std::string::size_type i1 = line.find("\"");
521     std::string::size_type i2;
522     std::string::size_type len = line.size();
523     WeightContainer namedWeight;
524     for ( WeightContainer::size_type ii = 0; ii < name_size; ++ii ) {
525         // weight names may contain blanks
526         if(i1 >= len) {
527             std::cout << "debug: attempting to read past the end of the named weight line " << std::endl;
528             std::cout << "debug: We should never get here" << std::endl;
529             std::cout << "debug: Looking for the end of this event" << std::endl;
530             detail::find_event_end( is );
531         }
532         i2 = line.find("\"",i1+1);
533         name = line.substr(i1+1,i2-i1-1);
534         namedWeight[name] = m_weights[ii];
535         i1 = line.find("\"",i2+1);
536     }
537     m_weights = namedWeight;
538     return is;
539 }
540  
390 garren 541 std::istream & GenEvent::read_units( std::istream & is )
542 {
543     //
544     if ( !is ) {
432 garren 545         std::cerr << "GenEvent::read_units setting badbit." << std::endl;
390 garren 546         is.clear(std::ios::badbit);
547         return is;
548     }
549     //
550     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     if ( is.peek() !='U') {
555         use_units( info.io_momentum_unit(),
556                        info.io_position_unit() );
557         return is;
558     }
559     is.ignore();        // ignore the first character in the line
560     std::string mom, pos;
561     is >> mom >> pos;
562     is.ignore(1);      // eat the extra whitespace
563     use_units(mom,pos);
564     //
565     return is;
566 }
567  
568 std::istream & GenEvent::find_file_type( std::istream & istr )
569 {
570     //
571     // make sure the stream is good
572     if ( !istr ) return istr;
573  
574     //
575     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     if ( istr.peek()=='E' ) {
580         info.set_io_type( gen );
581         info.set_has_key(false);
582         return istr;
583     }
584  
585     std::string line;
586     while ( std::getline(istr,line) ) {
587         //
588         // search for event listing key before first event only.
589         //
590         if( line == info.IO_GenEvent_Key() ) {
591             info.set_io_type( gen );
592             info.set_has_key(true);
593             return istr;
594         } else if( line == info.IO_Ascii_Key() ) {
595             info.set_io_type( ascii );
596             info.set_has_key(true);
597             return istr;
598         } else if( line == info.IO_ExtendedAscii_Key() ) {
599             info.set_io_type( extascii );
600             info.set_has_key(true);
601             return istr;
602         } else if( line == info.IO_Ascii_PDT_Key() ) {
603             info.set_io_type( ascii_pdt );
604             info.set_has_key(true);
605             return istr;
606         } else if( line == info.IO_ExtendedAscii_PDT_Key() ) {
607             info.set_io_type( extascii_pdt );
608             info.set_has_key(true);
609             return istr;
610         }
611     }
612     info.set_io_type( 0 );
613     info.set_has_key(false);
614     return istr;
615 }
616  
617 std::istream & GenEvent::find_end_key( std::istream & istr, int & iotype )
618 {
619     iotype = 0;
620     // peek at the first character before proceeding
621     if( istr.peek()!='H' ) return istr;
622     //
623     // we only check the next line
624     std::string line;
625     std::getline(istr,line);
626     //
627     StreamInfo & info = get_stream_info(istr);
628     //
629     // check to see if this is an end key
630     if( line == info.IO_GenEvent_End() ) {
631         iotype = gen;
632     } else if( line == info.IO_Ascii_End() ) {
633         iotype = ascii;
634     } else if( line == info.IO_ExtendedAscii_End() ) {
635         iotype = extascii;
636     } else if( line == info.IO_Ascii_PDT_End() ) {
637         iotype = ascii_pdt;
638     } else if( line == info.IO_ExtendedAscii_PDT_End() ) {
639         iotype = extascii_pdt;
640     }
641     if( iotype != 0 && info.io_type() != iotype ) {
642         std::cerr << "GenEvent::find_end_key: iotype keys have changed" << std::endl;
643     } else {
644         return istr;
645     }
646     //
647     // if we get here, then something has gotten badly confused
648     std::cerr << "GenEvent::find_end_key: MALFORMED INPUT" << std::endl;
649     istr.clear(std::ios::badbit);
650     return istr;
651 }
652  
653 std::ostream & establish_output_stream_info( std::ostream & os )
654 {
655     StreamInfo & info = get_stream_info(os);
656     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         os.precision(16);
660         // we use decimal to store integers, because it is smaller than hex!
661         os.setf(std::ios::dec,std::ios::basefield);
662         os.setf(std::ios::scientific,std::ios::floatfield);
663     }
664     return os;
665 }
666  
667 std::istream & establish_input_stream_info( std::istream & is )
668 {
669     StreamInfo & info = get_stream_info(is);
670     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         is.precision(16);
674         // we use decimal to store integers, because it is smaller than hex!
675         is.setf(std::ios::dec,std::ios::basefield);
676         is.setf(std::ios::scientific,std::ios::floatfield);
677     }
678     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     std::string line;
694     std::getline(is,line);
695     std::istringstream iline(line);
696     std::string firstc;
697     iline >> firstc;
698     if( firstc != "P" ) {
699         std::cerr << "StreamHelpers::detail::read_particle invalid line type: "
700                   << firstc << std::endl;
701         std::cerr << "StreamHelpers::detail::read_particle setting badbit."
702                   << std::endl;
703         is.clear(std::ios::badbit);
704         return is;
705     }
706     //
707     StreamInfo & info = get_stream_info(is);
443 garren 708     //testHepMC.cc
390 garren 709     // declare variables to be read in to, and read everything except flow
710     double px = 0., py = 0., pz = 0., e = 0., m = 0., theta = 0., phi = 0.;
711     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     iline >> bar_code ;
443 garren 714     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 715     iline >> id ;
443 garren 716     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 717     iline >> px ;
443 garren 718     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 719     iline >> py ;
443 garren 720     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 721     iline >> pz ;
443 garren 722     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 723     iline >> e ;
443 garren 724     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 725     if( info.io_type() != ascii ) {
726         iline >> m ;
443 garren 727         if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 728     }
729     iline >> status ;
443 garren 730     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 731     iline >> theta ;
443 garren 732     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 733     iline >> phi ;
443 garren 734     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 735     iline >> end_vtx_code ;
443 garren 736     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 737     iline >> flow_size;
443 garren 738     if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 739     //
740     // read flow patterns if any exist
741     Flow flow;
742     int code_index, code;
743     for ( int i = 1; i <= flow_size; ++i ) {
744         iline >> code_index >> code;
443 garren 745         if(!iline) {  delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); }
390 garren 746         flow.set_icode( code_index,code);
747     }
748     p->set_momentum( FourVector(px,py,pz,e) );
749     p->set_pdg_id( id );
750     p->set_status( status );
751     p->set_flow( flow );
752     p->set_polarization( Polarization(theta,phi) );
753     if( info.io_type() == ascii ) {
754         p->set_generated_mass( p->momentum().m() );
755     } else {
756         p->set_generated_mass( m );
757     }
758     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     if ( end_vtx_code != 0 ) {
764         particle_to_end_vertex.addEndParticle(p,end_vtx_code);
765     }
766     return is;
767 }
768  
769 std::ostream & establish_output_stream_info( std::ostream & os )
770 {
771     StreamInfo & info = get_stream_info(os);
772     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         os.precision(16);
776         // we use decimal to store integers, because it is smaller than hex!
777         os.setf(std::ios::dec,std::ios::basefield);
778         os.setf(std::ios::scientific,std::ios::floatfield);
779     }
780     return os;
781 }
782  
783 std::istream & establish_input_stream_info( std::istream & is )
784 {
785     StreamInfo & info = get_stream_info(is);
786     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         is.precision(16);
790         // we use decimal to store integers, because it is smaller than hex!
791         is.setf(std::ios::dec,std::ios::basefield);
792         is.setf(std::ios::scientific,std::ios::floatfield);
793     }
794     return is;
795 }
796  
797 } // detail
798  
799 } // HepMC