]>
Commit | Line | Data |
---|---|---|
0ca57c2f | 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 | |
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 | os << ' ' << (int)weights().size() ; | |
115 | for ( WeightContainer::const_map_iterator w = weights().map_begin(); | |
116 | w != weights().map_end(); ++w ) { | |
117 | detail::output( os, m_weights[w->second] ); | |
118 | } | |
119 | 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 | 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 | } | |
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 | // | |
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() ; | |
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 | } | |
204 | ||
205 | int signal_process_vertex = 0; | |
206 | 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 | 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 | |
223 | units_line = true; | |
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 | |
231 | GenCrossSection xs; | |
232 | // check for invalid data | |
233 | try { | |
234 | // read the line | |
235 | xs.read(is); | |
236 | } | |
237 | catch (IO_Exception& e) { | |
238 | detail::find_event_end( is ); | |
239 | } | |
240 | if(xs.is_set()) { | |
241 | set_cross_section( xs ); | |
242 | } | |
243 | } break; | |
244 | case 'H': | |
245 | { // we have a HeavyIon line OR an unexpected HepMC... line | |
246 | if( info.io_type() == gen || info.io_type() == extascii ) { | |
247 | // get HeavyIon | |
248 | HeavyIon ion; | |
249 | // check for invalid data | |
250 | try { | |
251 | is >> &ion; | |
252 | } | |
253 | catch (IO_Exception& e) { | |
254 | detail::find_event_end( is ); | |
255 | } | |
256 | if(ion.is_valid()) { | |
257 | set_heavy_ion( ion ); | |
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 | |
265 | PdfInfo pdf; | |
266 | // check for invalid data | |
267 | try { | |
268 | is >> &pdf; | |
269 | } | |
270 | catch (IO_Exception& e) { | |
271 | detail::find_event_end( is ); | |
272 | } | |
273 | if(pdf.is_valid()) { | |
274 | set_pdf_info( pdf ); | |
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 | |
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 | } | |
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(); | |
306 | try { | |
307 | detail::read_vertex(is,particle_to_end_vertex,v); | |
308 | } | |
309 | catch (IO_Exception& e) { | |
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 | } | |
320 | delete v; | |
321 | detail::find_event_end( is ); | |
322 | } | |
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); | |
334 | for ( TempParticleMap::orderIterator pmap | |
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 | ||
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 | ||
541 | std::istream & GenEvent::read_units( std::istream & is ) | |
542 | { | |
543 | // | |
544 | if ( !is ) { | |
545 | std::cerr << "GenEvent::read_units setting badbit." << std::endl; | |
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); | |
708 | //testHepMC.cc | |
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 ; | |
714 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
715 | iline >> id ; | |
716 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
717 | iline >> px ; | |
718 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
719 | iline >> py ; | |
720 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
721 | iline >> pz ; | |
722 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
723 | iline >> e ; | |
724 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
725 | if( info.io_type() != ascii ) { | |
726 | iline >> m ; | |
727 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
728 | } | |
729 | iline >> status ; | |
730 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
731 | iline >> theta ; | |
732 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
733 | iline >> phi ; | |
734 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
735 | iline >> end_vtx_code ; | |
736 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
737 | iline >> flow_size; | |
738 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
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; | |
745 | if(!iline) { delete p; throw IO_Exception("read_particle input stream encounterd invalid data"); } | |
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 |