hepmc - Blame information for rev 442
Subversion Repositories:
| 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; | ||
| 207 | // OK - now ready to start reading the event, so set the header flag | ||
| 208 | info.set_reading_event_header(true); | ||
| 209 | // The flag will be set to false when we reach the end of the header | ||
| 210 | while(info.reading_event_header()) { | ||
| 211 | switch(is.peek()) { | ||
| 212 | case 'E': | ||
| 213 | { // deal with the event line | ||
| 214 | process_event_line( is, num_vertices, bp1, bp2, signal_process_vertex ); | ||
| 215 | } break; | ||
| 216 | case 'N': | ||
| 217 | { // get weight names | ||
| 218 | read_weight_names( is ); | ||
| 219 | } break; | ||
| 220 | case 'U': | ||
| 221 | { // get unit information if it exists | ||
| 222 | if( info.io_type() == gen ) { | ||
| 223 | read_units( is ); | ||
| 224 | } | ||
| 225 | } break; | ||
| 226 | case 'C': | ||
| 227 | { // we have a GenCrossSection line | ||
| 228 | // create cross section | ||
| 440 | garren | 229 | GenCrossSection xs; |
| 432 | garren | 230 | // check for invalid data |
| 231 | try { | ||
| 232 | // read the line | ||
| 440 | garren | 233 | xs.read(is); |
| 432 | garren | 234 | } |
| 235 | catch (IO_Exception& e) { | ||
| 236 | detail::find_event_end( is ); | ||
| 237 | } | ||
| 440 | garren | 238 | if(xs.is_set()) { |
| 239 | set_cross_section( xs ); | ||
| 432 | garren | 240 | } |
| 241 | } break; | ||
| 242 | case 'H': | ||
| 243 | { // we have a HeavyIon line | ||
| 244 | if( info.io_type() == gen || info.io_type() == extascii ) { | ||
| 245 | // get HeavyIon | ||
| 440 | garren | 246 | HeavyIon ion; |
| 432 | garren | 247 | // check for invalid data |
| 248 | try { | ||
| 440 | garren | 249 | is >> &ion; |
| 432 | garren | 250 | } |
| 251 | catch (IO_Exception& e) { | ||
| 252 | detail::find_event_end( is ); | ||
| 253 | } | ||
| 440 | garren | 254 | if(ion.is_valid()) { |
| 255 | set_heavy_ion( ion ); | ||
| 432 | garren | 256 | } |
| 257 | } | ||
| 258 | } break; | ||
| 259 | case 'F': | ||
| 260 | { // we have a PdfInfo line | ||
| 261 | if( info.io_type() == gen || info.io_type() == extascii ) { | ||
| 262 | // get PdfInfo | ||
| 440 | garren | 263 | PdfInfo pdf; |
| 432 | garren | 264 | // check for invalid data |
| 265 | try { | ||
| 440 | garren | 266 | is >> &pdf; |
| 432 | garren | 267 | } |
| 268 | catch (IO_Exception& e) { | ||
| 269 | detail::find_event_end( is ); | ||
| 270 | } | ||
| 440 | garren | 271 | if(pdf.is_valid()) { |
| 272 | set_pdf_info( pdf ); | ||
| 432 | garren | 273 | } |
| 274 | } | ||
| 275 | } break; | ||
| 276 | case 'V': | ||
| 277 | { | ||
| 278 | // this should be the first vertex line - exit this loop | ||
| 279 | info.set_reading_event_header(false); | ||
| 280 | } break; | ||
| 281 | case 'P': | ||
| 282 | { // we should not find this line | ||
| 283 | std::cerr << "streaming input: found unexpected line P" << std::endl; | ||
| 284 | info.set_reading_event_header(false); | ||
| 285 | } break; | ||
| 286 | default: | ||
| 287 | // ignore everything else | ||
| 288 | break; | ||
| 289 | } // switch on line type | ||
| 290 | } // while reading_event_header | ||
| 390 | garren | 291 | // |
| 292 | // the end vertices of the particles are not connected until | ||
| 293 | // after the event is read --- we store the values in a map until then | ||
| 294 | TempParticleMap particle_to_end_vertex; | ||
| 295 | // | ||
| 296 | // read in the vertices | ||
| 297 | for ( int iii = 1; iii <= num_vertices; ++iii ) { | ||
| 298 | GenVertex* v = new GenVertex(); | ||
| 299 | detail::read_vertex(is,particle_to_end_vertex,v); | ||
| 300 | add_vertex( v ); | ||
| 301 | } | ||
| 302 | // set the signal process vertex | ||
| 303 | if ( signal_process_vertex ) { | ||
| 304 | set_signal_process_vertex( | ||
| 305 | barcode_to_vertex(signal_process_vertex) ); | ||
| 306 | } | ||
| 307 | // | ||
| 308 | // last connect particles to their end vertices | ||
| 309 | GenParticle* beam1(0); | ||
| 310 | GenParticle* beam2(0); | ||
| 311 | for ( std::map<int,GenParticle*>::iterator pmap | ||
| 312 | = particle_to_end_vertex.order_begin(); | ||
| 313 | pmap != particle_to_end_vertex.order_end(); ++pmap ) { | ||
| 314 | GenParticle* p = pmap->second; | ||
| 315 | int vtx = particle_to_end_vertex.end_vertex( p ); | ||
| 316 | GenVertex* itsDecayVtx = barcode_to_vertex(vtx); | ||
| 317 | if ( itsDecayVtx ) itsDecayVtx->add_particle_in( p ); | ||
| 318 | else { | ||
| 319 | std::cerr << "read_io_genevent: ERROR particle points" | ||
| 320 | << " to null end vertex. " <<std::endl; | ||
| 321 | } | ||
| 322 | // also look for the beam particles | ||
| 323 | if( p->barcode() == bp1 ) beam1 = p; | ||
| 324 | if( p->barcode() == bp2 ) beam2 = p; | ||
| 325 | } | ||
| 326 | set_beam_particles(beam1,beam2); | ||
| 327 | return is; | ||
| 328 | } | ||
| 329 | |||
| 330 | // ------------------------- operator << and operator >> ---------------- | ||
| 331 | |||
| 332 | std::ostream & operator << (std::ostream & os, GenEvent & evt) | ||
| 333 | { | ||
| 334 | /// Writes evt to an output stream. | ||
| 335 | evt.write(os); | ||
| 336 | return os; | ||
| 337 | } | ||
| 338 | |||
| 339 | std::istream & operator >> (std::istream & is, GenEvent & evt) | ||
| 340 | { | ||
| 341 | evt.read(is); | ||
| 342 | return is; | ||
| 343 | } | ||
| 344 | |||
| 345 | // ------------------------- set units ---------------- | ||
| 346 | |||
| 347 | std::istream & set_input_units(std::istream & is, | ||
| 348 | Units::MomentumUnit mom, | ||
| 349 | Units::LengthUnit len ) | ||
| 350 | { | ||
| 351 | // | ||
| 352 | StreamInfo & info = get_stream_info(is); | ||
| 353 | info.use_input_units( mom, len ); | ||
| 354 | return is; | ||
| 355 | } | ||
| 356 | |||
| 357 | // ------------------------- begin and end block lines ---------------- | ||
| 358 | |||
| 359 | std::ostream & write_HepMC_IO_block_begin(std::ostream & os ) | ||
| 360 | { | ||
| 361 | // | ||
| 362 | StreamInfo & info = get_stream_info(os); | ||
| 363 | |||
| 364 | if( !info.finished_first_event() ) { | ||
| 365 | os << "\n" << "HepMC::Version " << versionName(); | ||
| 366 | os << "\n"; | ||
| 367 | os << info.IO_GenEvent_Key() << "\n"; | ||
| 368 | } | ||
| 369 | return os; | ||
| 370 | } | ||
| 371 | |||
| 372 | std::ostream & write_HepMC_IO_block_end(std::ostream & os ) | ||
| 373 | { | ||
| 374 | // | ||
| 375 | StreamInfo & info = get_stream_info(os); | ||
| 376 | |||
| 377 | if( info.finished_first_event() ) { | ||
| 378 | os << info.IO_GenEvent_End() << "\n"; | ||
| 379 | os << std::flush; | ||
| 380 | } | ||
| 381 | return os; | ||
| 382 | } | ||
| 383 | |||
| 432 | garren | 384 | std::istream & GenEvent::process_event_line( std::istream & is, |
| 385 | int & num_vertices, | ||
| 386 | int & bp1, int & bp2, | ||
| 387 | int & signal_process_vertex ) | ||
| 388 | { | ||
| 389 | // | ||
| 390 | if ( !is ) { | ||
| 391 | std::cerr << "GenEvent::process_event_line setting badbit." << std::endl; | ||
| 392 | is.clear(std::ios::badbit); | ||
| 393 | return is; | ||
| 394 | } | ||
| 395 | // | ||
| 396 | StreamInfo & info = get_stream_info(is); | ||
| 397 | std::string line; | ||
| 398 | std::getline(is,line); | ||
| 399 | std::istringstream iline(line); | ||
| 400 | std::string firstc; | ||
| 401 | iline >> firstc; | ||
| 402 | // | ||
| 403 | // read values into temp variables, then fill GenEvent | ||
| 404 | int event_number = 0, signal_process_id = 0, | ||
| 405 | random_states_size = 0, nmpi = -1; | ||
| 406 | double eventScale = 0, alpha_qcd = 0, alpha_qed = 0; | ||
| 407 | iline >> event_number; | ||
| 408 | if(!iline) detail::find_event_end( is ); | ||
| 409 | if( info.io_type() == gen || info.io_type() == extascii ) { | ||
| 410 | iline >> nmpi; | ||
| 411 | if(!iline) detail::find_event_end( is ); | ||
| 412 | set_mpi( nmpi ); | ||
| 413 | } | ||
| 414 | iline >> eventScale ; | ||
| 415 | if(!iline) detail::find_event_end( is ); | ||
| 416 | iline >> alpha_qcd ; | ||
| 417 | if(!iline) detail::find_event_end( is ); | ||
| 418 | iline >> alpha_qed; | ||
| 419 | if(!iline) detail::find_event_end( is ); | ||
| 420 | iline >> signal_process_id ; | ||
| 421 | if(!iline) detail::find_event_end( is ); | ||
| 422 | iline >> signal_process_vertex; | ||
| 423 | if(!iline) detail::find_event_end( is ); | ||
| 424 | iline >> num_vertices; | ||
| 425 | if(!iline) detail::find_event_end( is ); | ||
| 426 | if( info.io_type() == gen || info.io_type() == extascii ) { | ||
| 427 | iline >> bp1 ; | ||
| 428 | if(!iline) detail::find_event_end( is ); | ||
| 429 | iline >> bp2; | ||
| 430 | if(!iline) detail::find_event_end( is ); | ||
| 431 | } | ||
| 432 | iline >> random_states_size; | ||
| 433 | if(!iline) detail::find_event_end( is ); | ||
| 434 | std::vector<long> random_states(random_states_size); | ||
| 435 | for ( int i = 0; i < random_states_size; ++i ) { | ||
| 436 | iline >> random_states[i]; | ||
| 437 | if(!iline) detail::find_event_end( is ); | ||
| 438 | } | ||
| 439 | WeightContainer::size_type weights_size = 0; | ||
| 440 | iline >> weights_size; | ||
| 441 | if(!iline) detail::find_event_end( is ); | ||
| 442 | std::vector<double> wgt(weights_size); | ||
| 443 | for ( WeightContainer::size_type ii = 0; ii < weights_size; ++ii ) { | ||
| 444 | iline >> wgt[ii]; | ||
| 445 | if(!iline) detail::find_event_end( is ); | ||
| 446 | } | ||
| 447 | // weight names will be added later if they exist | ||
| 448 | if( weights_size > 0 ) m_weights = wgt; | ||
| 449 | // | ||
| 450 | // fill signal_process_id, event_number, random_states, etc. | ||
| 451 | set_signal_process_id( signal_process_id ); | ||
| 452 | set_event_number( event_number ); | ||
| 453 | set_random_states( random_states ); | ||
| 454 | set_event_scale( eventScale ); | ||
| 455 | set_alphaQCD( alpha_qcd ); | ||
| 456 | set_alphaQED( alpha_qed ); | ||
| 457 | // | ||
| 458 | return is; | ||
| 459 | } | ||
| 460 | |||
| 461 | std::istream & GenEvent::read_weight_names( std::istream & is ) | ||
| 462 | { | ||
| 463 | // now check for a named weight line | ||
| 464 | if ( !is ) { | ||
| 465 | std::cerr << "GenEvent::read_weight_names setting badbit." << std::endl; | ||
| 466 | is.clear(std::ios::badbit); | ||
| 467 | return is; | ||
| 468 | } | ||
| 469 | // Test to be sure the next entry is of type "N" | ||
| 470 | // If we have no named weight line, this is not an error | ||
| 471 | // releases prior to 2.06.00 do not have named weights | ||
| 472 | if ( is.peek() !='N') { | ||
| 473 | return is; | ||
| 474 | } | ||
| 475 | // now get this line and process it | ||
| 476 | std::string line; | ||
| 477 | std::getline(is,line); | ||
| 478 | std::istringstream wline(line); | ||
| 479 | std::string firstc; | ||
| 480 | WeightContainer::size_type name_size = 0; | ||
| 481 | wline >> firstc >> name_size; | ||
| 482 | if(!wline) detail::find_event_end( is ); | ||
| 483 | if( firstc != "N") { | ||
| 484 | std::cout << "debug: first character of named weights is " << firstc << std::endl; | ||
| 485 | std::cout << "debug: We should never get here" << std::endl; | ||
| 486 | is.clear(std::ios::badbit); | ||
| 487 | return is; | ||
| 488 | } | ||
| 489 | if( m_weights.size() != name_size ) { | ||
| 490 | std::cout << "debug: weight sizes do not match "<< std::endl; | ||
| 491 | std::cout << "debug: weight vector size is " << m_weights.size() << std::endl; | ||
| 492 | std::cout << "debug: weight name size is " << name_size << std::endl; | ||
| 493 | is.clear(std::ios::badbit); | ||
| 494 | return is; | ||
| 495 | } | ||
| 496 | std::string name; | ||
| 497 | std::string::size_type i1 = line.find("\""); | ||
| 498 | std::string::size_type i2; | ||
| 499 | std::string::size_type len = line.size(); | ||
| 500 | WeightContainer namedWeight; | ||
| 501 | for ( WeightContainer::size_type ii = 0; ii < name_size; ++ii ) { | ||
| 502 | // weight names may contain blanks | ||
| 503 | if(i1 >= len) { | ||
| 504 | std::cout << "debug: attempting to read past the end of the named weight line " << std::endl; | ||
| 505 | std::cout << "debug: We should never get here" << std::endl; | ||
| 506 | std::cout << "debug: Looking for the end of this event" << std::endl; | ||
| 507 | detail::find_event_end( is ); | ||
| 508 | } | ||
| 509 | i2 = line.find("\"",i1+1); | ||
| 510 | name = line.substr(i1+1,i2-i1-1); | ||
| 511 | namedWeight[name] = m_weights[ii]; | ||
| 512 | i1 = line.find("\"",i2+1); | ||
| 513 | } | ||
| 514 | m_weights = namedWeight; | ||
| 515 | return is; | ||
| 516 | } | ||
| 517 | |||
| 390 | garren | 518 | std::istream & GenEvent::read_units( std::istream & is ) |
| 519 | { | ||
| 520 | // | ||
| 521 | if ( !is ) { | ||
| 432 | garren | 522 | std::cerr << "GenEvent::read_units setting badbit." << std::endl; |
| 390 | garren | 523 | is.clear(std::ios::badbit); |
| 524 | return is; | ||
| 525 | } | ||
| 526 | // | ||
| 527 | StreamInfo & info = get_stream_info(is); | ||
| 528 | // test to be sure the next entry is of type "U" then ignore it | ||
| 529 | // if we have no units, this is not an error | ||
| 530 | // releases prior to 2.04.00 did not write unit information | ||
| 531 | if ( is.peek() !='U') { | ||
| 532 | use_units( info.io_momentum_unit(), | ||
| 533 | info.io_position_unit() ); | ||
| 534 | return is; | ||
| 535 | } | ||
| 536 | is.ignore(); // ignore the first character in the line | ||
| 537 | std::string mom, pos; | ||
| 538 | is >> mom >> pos; | ||
| 539 | is.ignore(1); // eat the extra whitespace | ||
| 540 | use_units(mom,pos); | ||
| 541 | // | ||
| 542 | return is; | ||
| 543 | } | ||
| 544 | |||
| 545 | std::istream & GenEvent::find_file_type( std::istream & istr ) | ||
| 546 | { | ||
| 547 | // | ||
| 548 | // make sure the stream is good | ||
| 549 | if ( !istr ) return istr; | ||
| 550 | |||
| 551 | // | ||
| 552 | StreamInfo & info = get_stream_info(istr); | ||
| 553 | |||
| 554 | // if there is no input block line, then we assume this stream | ||
| 555 | // is in the IO_GenEvent format | ||
| 556 | if ( istr.peek()=='E' ) { | ||
| 557 | info.set_io_type( gen ); | ||
| 558 | info.set_has_key(false); | ||
| 559 | return istr; | ||
| 560 | } | ||
| 561 | |||
| 562 | std::string line; | ||
| 563 | while ( std::getline(istr,line) ) { | ||
| 564 | // | ||
| 565 | // search for event listing key before first event only. | ||
| 566 | // | ||
| 567 | if( line == info.IO_GenEvent_Key() ) { | ||
| 568 | info.set_io_type( gen ); | ||
| 569 | info.set_has_key(true); | ||
| 570 | return istr; | ||
| 571 | } else if( line == info.IO_Ascii_Key() ) { | ||
| 572 | info.set_io_type( ascii ); | ||
| 573 | info.set_has_key(true); | ||
| 574 | return istr; | ||
| 575 | } else if( line == info.IO_ExtendedAscii_Key() ) { | ||
| 576 | info.set_io_type( extascii ); | ||
| 577 | info.set_has_key(true); | ||
| 578 | return istr; | ||
| 579 | } else if( line == info.IO_Ascii_PDT_Key() ) { | ||
| 580 | info.set_io_type( ascii_pdt ); | ||
| 581 | info.set_has_key(true); | ||
| 582 | return istr; | ||
| 583 | } else if( line == info.IO_ExtendedAscii_PDT_Key() ) { | ||
| 584 | info.set_io_type( extascii_pdt ); | ||
| 585 | info.set_has_key(true); | ||
| 586 | return istr; | ||
| 587 | } | ||
| 588 | } | ||
| 589 | info.set_io_type( 0 ); | ||
| 590 | info.set_has_key(false); | ||
| 591 | return istr; | ||
| 592 | } | ||
| 593 | |||
| 594 | std::istream & GenEvent::find_end_key( std::istream & istr, int & iotype ) | ||
| 595 | { | ||
| 596 | iotype = 0; | ||
| 597 | // peek at the first character before proceeding | ||
| 598 | if( istr.peek()!='H' ) return istr; | ||
| 599 | // | ||
| 600 | // we only check the next line | ||
| 601 | std::string line; | ||
| 602 | std::getline(istr,line); | ||
| 603 | // | ||
| 604 | StreamInfo & info = get_stream_info(istr); | ||
| 605 | // | ||
| 606 | // check to see if this is an end key | ||
| 607 | if( line == info.IO_GenEvent_End() ) { | ||
| 608 | iotype = gen; | ||
| 609 | } else if( line == info.IO_Ascii_End() ) { | ||
| 610 | iotype = ascii; | ||
| 611 | } else if( line == info.IO_ExtendedAscii_End() ) { | ||
| 612 | iotype = extascii; | ||
| 613 | } else if( line == info.IO_Ascii_PDT_End() ) { | ||
| 614 | iotype = ascii_pdt; | ||
| 615 | } else if( line == info.IO_ExtendedAscii_PDT_End() ) { | ||
| 616 | iotype = extascii_pdt; | ||
| 617 | } | ||
| 618 | if( iotype != 0 && info.io_type() != iotype ) { | ||
| 619 | std::cerr << "GenEvent::find_end_key: iotype keys have changed" << std::endl; | ||
| 620 | } else { | ||
| 621 | return istr; | ||
| 622 | } | ||
| 623 | // | ||
| 624 | // if we get here, then something has gotten badly confused | ||
| 625 | std::cerr << "GenEvent::find_end_key: MALFORMED INPUT" << std::endl; | ||
| 626 | istr.clear(std::ios::badbit); | ||
| 627 | return istr; | ||
| 628 | } | ||
| 629 | |||
| 630 | std::ostream & establish_output_stream_info( std::ostream & os ) | ||
| 631 | { | ||
| 632 | StreamInfo & info = get_stream_info(os); | ||
| 633 | if ( !info.finished_first_event() ) { | ||
| 634 | // precision 16 (# digits following decimal point) is the minimum that | ||
| 635 | // will capture the full information stored in a double | ||
| 636 | os.precision(16); | ||
| 637 | // we use decimal to store integers, because it is smaller than hex! | ||
| 638 | os.setf(std::ios::dec,std::ios::basefield); | ||
| 639 | os.setf(std::ios::scientific,std::ios::floatfield); | ||
| 640 | } | ||
| 641 | return os; | ||
| 642 | } | ||
| 643 | |||
| 644 | std::istream & establish_input_stream_info( std::istream & is ) | ||
| 645 | { | ||
| 646 | StreamInfo & info = get_stream_info(is); | ||
| 647 | if ( !info.finished_first_event() ) { | ||
| 648 | // precision 16 (# digits following decimal point) is the minimum that | ||
| 649 | // will capture the full information stored in a double | ||
| 650 | is.precision(16); | ||
| 651 | // we use decimal to store integers, because it is smaller than hex! | ||
| 652 | is.setf(std::ios::dec,std::ios::basefield); | ||
| 653 | is.setf(std::ios::scientific,std::ios::floatfield); | ||
| 654 | } | ||
| 655 | return is; | ||
| 656 | } | ||
| 657 | |||
| 658 | |||
| 659 | // ------------------------- helper functions ---------------- | ||
| 660 | |||
| 661 | namespace detail { | ||
| 662 | |||
| 663 | // The functions defined here need to use get_stream_info | ||
| 664 | |||
| 665 | std::istream & read_particle( std::istream & is, | ||
| 666 | TempParticleMap & particle_to_end_vertex, | ||
| 667 | GenParticle * p ) | ||
| 668 | { | ||
| 669 | // get the next line | ||
| 670 | std::string line; | ||
| 671 | std::getline(is,line); | ||
| 672 | std::istringstream iline(line); | ||
| 673 | std::string firstc; | ||
| 674 | iline >> firstc; | ||
| 675 | if( firstc != "P" ) { | ||
| 676 | std::cerr << "StreamHelpers::detail::read_particle invalid line type: " | ||
| 677 | << firstc << std::endl; | ||
| 678 | std::cerr << "StreamHelpers::detail::read_particle setting badbit." | ||
| 679 | << std::endl; | ||
| 680 | is.clear(std::ios::badbit); | ||
| 681 | return is; | ||
| 682 | } | ||
| 683 | // | ||
| 684 | StreamInfo & info = get_stream_info(is); | ||
| 685 | // | ||
| 686 | // declare variables to be read in to, and read everything except flow | ||
| 687 | double px = 0., py = 0., pz = 0., e = 0., m = 0., theta = 0., phi = 0.; | ||
| 688 | int bar_code = 0, id = 0, status = 0, end_vtx_code = 0, flow_size = 0; | ||
| 689 | // check that the input stream is still OK after reading item | ||
| 690 | iline >> bar_code ; | ||
| 442 | garren | 691 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 692 | iline >> id ; |
| 442 | garren | 693 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 694 | iline >> px ; |
| 442 | garren | 695 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 696 | iline >> py ; |
| 442 | garren | 697 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 698 | iline >> pz ; |
| 442 | garren | 699 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 700 | iline >> e ; |
| 442 | garren | 701 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 702 | if( info.io_type() != ascii ) { |
| 703 | iline >> m ; | ||
| 442 | garren | 704 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 705 | } |
| 706 | iline >> status ; | ||
| 442 | garren | 707 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 708 | iline >> theta ; |
| 442 | garren | 709 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 710 | iline >> phi ; |
| 442 | garren | 711 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 712 | iline >> end_vtx_code ; |
| 442 | garren | 713 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 714 | iline >> flow_size; |
| 442 | garren | 715 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 716 | // |
| 717 | // read flow patterns if any exist | ||
| 718 | Flow flow; | ||
| 719 | int code_index, code; | ||
| 720 | for ( int i = 1; i <= flow_size; ++i ) { | ||
| 721 | iline >> code_index >> code; | ||
| 442 | garren | 722 | if(!iline) { delete p; particle_to_end_vertex.clear_temp_map(); detail::find_event_end( is ); } |
| 390 | garren | 723 | flow.set_icode( code_index,code); |
| 724 | } | ||
| 725 | p->set_momentum( FourVector(px,py,pz,e) ); | ||
| 726 | p->set_pdg_id( id ); | ||
| 727 | p->set_status( status ); | ||
| 728 | p->set_flow( flow ); | ||
| 729 | p->set_polarization( Polarization(theta,phi) ); | ||
| 730 | if( info.io_type() == ascii ) { | ||
| 731 | p->set_generated_mass( p->momentum().m() ); | ||
| 732 | } else { | ||
| 733 | p->set_generated_mass( m ); | ||
| 734 | } | ||
| 735 | p->suggest_barcode( bar_code ); | ||
| 736 | // | ||
| 737 | // all particles are connected to their end vertex separately | ||
| 738 | // after all particles and vertices have been created - so we keep | ||
| 739 | // a map of all particles that have end vertices | ||
| 740 | if ( end_vtx_code != 0 ) { | ||
| 741 | particle_to_end_vertex.addEndParticle(p,end_vtx_code); | ||
| 742 | } | ||
| 743 | return is; | ||
| 744 | } | ||
| 745 | |||
| 746 | std::istream & read_units( std::istream & is, GenEvent & evt ) | ||
| 747 | { | ||
| 748 | // | ||
| 749 | if ( !is ) { | ||
| 750 | std::cerr << "StreamHelpers read_units setting badbit." << std::endl; | ||
| 751 | is.clear(std::ios::badbit); | ||
| 752 | return is; | ||
| 753 | } | ||
| 754 | // | ||
| 755 | StreamInfo & info = get_stream_info(is); | ||
| 756 | // test to be sure the next entry is of type "U" then ignore it | ||
| 757 | // if we have no units, this is not an error | ||
| 758 | // releases prior to 2.04.00 did not write unit information | ||
| 759 | if ( is.peek() !='U') { | ||
| 760 | evt.use_units( info.io_momentum_unit(), | ||
| 761 | info.io_position_unit() ); | ||
| 762 | return is; | ||
| 763 | } | ||
| 764 | is.ignore(); // ignore the first character in the line | ||
| 765 | std::string mom, pos; | ||
| 766 | is >> mom >> pos; | ||
| 767 | is.ignore(1); // eat the extra whitespace | ||
| 768 | evt.use_units(mom,pos); | ||
| 769 | // | ||
| 770 | return is; | ||
| 771 | } | ||
| 772 | |||
| 773 | std::ostream & establish_output_stream_info( std::ostream & os ) | ||
| 774 | { | ||
| 775 | StreamInfo & info = get_stream_info(os); | ||
| 776 | if ( !info.finished_first_event() ) { | ||
| 777 | // precision 16 (# digits following decimal point) is the minimum that | ||
| 778 | // will capture the full information stored in a double | ||
| 779 | os.precision(16); | ||
| 780 | // we use decimal to store integers, because it is smaller than hex! | ||
| 781 | os.setf(std::ios::dec,std::ios::basefield); | ||
| 782 | os.setf(std::ios::scientific,std::ios::floatfield); | ||
| 783 | } | ||
| 784 | return os; | ||
| 785 | } | ||
| 786 | |||
| 787 | std::istream & establish_input_stream_info( std::istream & is ) | ||
| 788 | { | ||
| 789 | StreamInfo & info = get_stream_info(is); | ||
| 790 | if ( !info.finished_first_event() ) { | ||
| 791 | // precision 16 (# digits following decimal point) is the minimum that | ||
| 792 | // will capture the full information stored in a double | ||
| 793 | is.precision(16); | ||
| 794 | // we use decimal to store integers, because it is smaller than hex! | ||
| 795 | is.setf(std::ios::dec,std::ios::basefield); | ||
| 796 | is.setf(std::ios::scientific,std::ios::floatfield); | ||
| 797 | } | ||
| 798 | return is; | ||
| 799 | } | ||
| 800 | |||
| 801 | } // detail | ||
| 802 | |||
| 803 | } // HepMC |
