Changeset 571 in svn
- Timestamp:
- Nov 2, 2011, 5:06:22 PM (13 years ago)
- Location:
- trunk/Utilities/HepMC/src
- Files:
-
- 9 added
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/Utilities/HepMC/src/CompareGenEvent.cc
r349 r571 91 91 92 92 bool compareWeights( GenEvent* e1, GenEvent* e2 ) { 93 WeightContainer w1 = e1->weights(); 94 WeightContainer w2 = e2->weights(); 95 if( w1.size() != w2.size() ) { 96 std::cerr << "compareWeights: size of weight container differs " << std::endl; 97 return false; 98 } 99 for( int i=0; i<w1.size(); ++i ) { 100 if( w1[i] != w2[i] ) { 101 std::cerr << "compareWeights: weight container entry " 102 << i << " differs" << std::endl; 103 return false; 104 } 105 } 106 return true; 93 if( e1->weights() == e2->weights() ) return true; 94 std::cerr << "compareWeights: weight containers differ " << std::endl; 95 return false; 107 96 } 108 97 -
trunk/Utilities/HepMC/src/GenEvent.cc
r349 r571 14 14 15 15 #include "GenEvent.h" 16 #include "GenCrossSection.h" 16 17 #include "Version.h" 18 #include "StreamHelpers.h" 17 19 18 20 namespace HepMC { … … 38 40 m_vertex_barcodes(), 39 41 m_particle_barcodes(), 42 m_cross_section(0), 40 43 m_heavy_ion(0), 41 44 m_pdf_info(0), … … 71 74 m_vertex_barcodes(), 72 75 m_particle_barcodes(), 76 m_cross_section(0), 73 77 m_heavy_ion( new HeavyIon(ion) ), 74 78 m_pdf_info( new PdfInfo(pdf) ), … … 102 106 m_vertex_barcodes(), 103 107 m_particle_barcodes(), 108 m_cross_section(0), 104 109 m_heavy_ion(0), 105 110 m_pdf_info(0), … … 136 141 m_vertex_barcodes(), 137 142 m_particle_barcodes(), 143 m_cross_section(0), 138 144 m_heavy_ion( new HeavyIon(ion) ), 139 145 m_pdf_info( new PdfInfo(pdf) ), … … 162 168 m_vertex_barcodes ( /* inevent.m_vertex_barcodes */ ), 163 169 m_particle_barcodes ( /* inevent.m_particle_barcodes */ ), 170 m_cross_section ( inevent.cross_section() ? new GenCrossSection(*inevent.cross_section()) : 0 ), 164 171 m_heavy_ion ( inevent.heavy_ion() ? new HeavyIon(*inevent.heavy_ion()) : 0 ), 165 172 m_pdf_info ( inevent.pdf_info() ? new PdfInfo(*inevent.pdf_info()) : 0 ), … … 168 175 { 169 176 /// deep copy - makes a copy of all vertices! 170 //++s_counter;171 177 // 172 178 … … 234 240 m_vertex_barcodes.swap( other.m_vertex_barcodes ); 235 241 m_particle_barcodes.swap( other.m_particle_barcodes ); 242 std::swap(m_cross_section , other.m_cross_section ); 236 243 std::swap(m_heavy_ion , other.m_heavy_ion ); 237 244 std::swap(m_pdf_info , other.m_pdf_info ); … … 252 259 { 253 260 /// Deep destructor. 254 /// deletes all vertices/particles in this evt255 /// 261 /// deletes all vertices/particles in this GenEvent 262 /// deletes the associated HeavyIon and PdfInfo 256 263 delete_all_vertices(); 264 delete m_cross_section; 257 265 delete m_heavy_ion; 258 266 delete m_pdf_info; 259 //--s_counter;260 267 } 261 268 … … 281 288 : 0 ) 282 289 << "\n"; 283 //ostr << " Current Memory Usage: "284 // << GenEvent::counter() << " events, "285 // << GenVertex::counter() << " vertices, "286 // << GenParticle::counter() << " particles.\n";287 290 write_units( ostr ); 291 write_cross_section(ostr); 288 292 ostr << " Entries this event: " << vertices_size() << " vertices, " 289 293 << particles_size() << " particles.\n"; … … 302 306 // Weights 303 307 ostr << " Wgts(" << weights().size() << ")="; 304 for ( WeightContainer::const_iterator wgt = weights().begin(); 305 wgt != weights().end(); ++wgt ) { ostr << *wgt << " "; } 306 ostr << "\n"; 308 weights().print(ostr); 307 309 ostr << " EventScale " << event_scale() 308 310 << " [energy] \t alphaQCD=" << alphaQCD() … … 367 369 /// 368 370 delete_all_vertices(); 371 // remove existing objects and set pointers to null 372 delete m_cross_section; 373 m_cross_section = 0; 369 374 delete m_heavy_ion; 375 m_heavy_ion = 0; 370 376 delete m_pdf_info; 377 m_pdf_info = 0; 371 378 m_signal_process_id = 0; 372 379 m_beam_particle_1 = 0; 373 380 m_beam_particle_2 = 0; 374 381 m_event_number = 0; 382 m_mpi = -1; 375 383 m_event_scale = -1; 376 384 m_alphaQCD = -1; … … 378 386 m_weights = std::vector<double>(); 379 387 m_random_states = std::vector<long>(); 380 //m_momentum_unit = ; 381 //m_position_unit = ; 388 // resetting unit information 389 m_momentum_unit = Units::default_momentum_unit(); 390 m_position_unit = Units::default_length_unit(); 382 391 // error check just to be safe 383 392 if ( m_vertex_barcodes.size() != 0 … … 388 397 << m_vertex_barcodes.size() << " " 389 398 << m_particle_barcodes.size() << std::endl; 390 //std::cerr << "Total Number vtx,particle in memory "391 // << "after method called = "392 // << GenVertex::counter() << "\t"393 // << GenParticle::counter() << std::endl;394 399 } 395 400 return; … … 420 425 << m_vertex_barcodes.size() << " " 421 426 << m_particle_barcodes.size() << std::endl; 422 //std::cerr << "Total Number vtx,particle in memory "423 // << "after method called = "424 // << GenVertex::counter() << "\t"425 // << GenParticle::counter() << std::endl;426 427 } 427 428 } … … 602 603 } 603 604 605 void GenEvent::write_cross_section( std::ostream& os ) const 606 { 607 // write the GenCrossSection information if the cross section was set 608 if( !cross_section() ) return; 609 if( cross_section()->is_set() ) { 610 os << " Cross Section: " << cross_section()->cross_section() ; 611 os << " +/- " << cross_section()->cross_section_error() ; 612 os << std::endl; 613 } 614 } 615 604 616 bool GenEvent::use_momentum_unit( Units::MomentumUnit newunit ) { 605 617 // currently not exception-safe. … … 651 663 } 652 664 665 bool GenEvent::is_valid() const { 666 /// A GenEvent is presumed valid if it has both associated 667 /// particles and vertices. No other information is checked. 668 if ( vertices_empty() ) return false; 669 if ( particles_empty() ) return false; 670 return true; 671 } 672 673 std::ostream & GenEvent::write_beam_particles(std::ostream & os, 674 std::pair<HepMC::GenParticle *,HepMC::GenParticle *> pr ) 675 { 676 GenParticle* p = pr.first; 677 if(!p) { 678 detail::output( os, 0 ); 679 } else { 680 detail::output( os, p->barcode() ); 681 } 682 p = pr.second; 683 if(!p) { 684 detail::output( os, 0 ); 685 } else { 686 detail::output( os, p->barcode() ); 687 } 688 689 return os; 690 } 691 692 std::ostream & GenEvent::write_vertex(std::ostream & os, GenVertex const * v) 693 { 694 if ( !v || !os ) { 695 std::cerr << "GenEvent::write_vertex !v||!os, " 696 << "v="<< v << " setting badbit" << std::endl; 697 os.clear(std::ios::badbit); 698 return os; 699 } 700 // First collect info we need 701 // count the number of orphan particles going into v 702 int num_orphans_in = 0; 703 for ( GenVertex::particles_in_const_iterator p1 704 = v->particles_in_const_begin(); 705 p1 != v->particles_in_const_end(); ++p1 ) { 706 if ( !(*p1)->production_vertex() ) ++num_orphans_in; 707 } 708 // 709 os << 'V'; 710 detail::output( os, v->barcode() ); // v's unique identifier 711 detail::output( os, v->id() ); 712 detail::output( os, v->position().x() ); 713 detail::output( os, v->position().y() ); 714 detail::output( os, v->position().z() ); 715 detail::output( os, v->position().t() ); 716 detail::output( os, num_orphans_in ); 717 detail::output( os, (int)v->particles_out_size() ); 718 detail::output( os, (int)v->weights().size() ); 719 for ( WeightContainer::const_iterator w = v->weights().begin(); 720 w != v->weights().end(); ++w ) { 721 detail::output( os, *w ); 722 } 723 detail::output( os,'\n'); 724 // incoming particles 725 for ( GenVertex::particles_in_const_iterator p2 726 = v->particles_in_const_begin(); 727 p2 != v->particles_in_const_end(); ++p2 ) { 728 if ( !(*p2)->production_vertex() ) { 729 write_particle( os, *p2 ); 730 } 731 } 732 // outgoing particles 733 for ( GenVertex::particles_out_const_iterator p3 734 = v->particles_out_const_begin(); 735 p3 != v->particles_out_const_end(); ++p3 ) { 736 write_particle( os, *p3 ); 737 } 738 return os; 739 } 740 741 std::ostream & GenEvent::write_particle( std::ostream & os, GenParticle const * p ) 742 { 743 if ( !p || !os ) { 744 std::cerr << "GenEvent::write_particle !p||!os, " 745 << "p="<< p << " setting badbit" << std::endl; 746 os.clear(std::ios::badbit); 747 return os; 748 } 749 os << 'P'; 750 detail::output( os, p->barcode() ); 751 detail::output( os, p->pdg_id() ); 752 detail::output( os, p->momentum().px() ); 753 detail::output( os, p->momentum().py() ); 754 detail::output( os, p->momentum().pz() ); 755 detail::output( os, p->momentum().e() ); 756 detail::output( os, p->generated_mass() ); 757 detail::output( os, p->status() ); 758 detail::output( os, p->polarization().theta() ); 759 detail::output( os, p->polarization().phi() ); 760 // since end_vertex is oftentimes null, this CREATES a null vertex 761 // in the map 762 detail::output( os, ( p->end_vertex() ? p->end_vertex()->barcode() : 0 ) ); 763 os << ' ' << p->flow() << "\n"; 764 765 return os; 766 } 767 653 768 } // HepMC -
trunk/Utilities/HepMC/src/GenParticle.cc
r349 r571 39 39 m_pdg_id( inparticle.pdg_id() ), 40 40 m_status( inparticle.status() ), 41 m_flow( this),41 m_flow(inparticle.flow()), 42 42 m_polarization( inparticle.polarization() ), 43 43 m_production_vertex(0), … … 188 188 /// Dump this particle's full info to ostr 189 189 std::ostream& operator<<( std::ostream& ostr, const GenParticle& part ) { 190 // find the current stream state 191 std::ios_base::fmtflags orig = ostr.flags(); 192 std::streamsize prec = ostr.precision(); 190 193 ostr << " "; 191 194 ostr.width(9); … … 199 202 ostr << part.momentum().px() << ","; 200 203 ostr.width(9); 201 ostr.precision(2);202 204 ostr << part.momentum().py() << ","; 203 205 ostr.width(9); 204 ostr.precision(2);205 206 ostr << part.momentum().pz() << ","; 206 207 ostr.width(9); 207 ostr.precision(2);208 208 ostr << part.momentum().e() << " "; 209 209 ostr.setf(std::ios::fmtflags(0), std::ios::floatfield); … … 227 227 ostr << (void*)part.end_vertex(); 228 228 } 229 // restore the stream state 230 ostr.flags(orig); 231 ostr.precision(prec); 229 232 return ostr; 230 233 } -
trunk/Utilities/HepMC/src/GenVertex.cc
r349 r571 144 144 145 145 void GenVertex::print( std::ostream& ostr ) const { 146 // find the current stream state 147 std::ios_base::fmtflags orig = ostr.flags(); 148 std::streamsize prec = ostr.precision(); 146 149 if ( barcode()!=0 ) { 147 150 if ( position() != FourVector(0,0,0,0) ) { … … 159 162 ostr << position().x() << ","; 160 163 ostr.width(9); 161 ostr.precision(2);162 164 ostr << position().y() << ","; 163 165 ostr.width(9); 164 ostr.precision(2);165 166 ostr << position().z() << ","; 166 167 ostr.width(9); 167 ostr.precision(2);168 168 ostr << position().t(); 169 169 ostr.setf(std::ios::fmtflags(0), std::ios::floatfield); … … 198 198 ostr << position().x(); 199 199 ostr.width(9); 200 ostr.precision(2);201 200 ostr << position().y(); 202 201 ostr.width(9); 203 ostr.precision(2);204 202 ostr << position().z(); 205 203 ostr.width(9); 206 ostr.precision(2);207 204 ostr << position().t(); 208 205 ostr.setf(std::ios::fmtflags(0), std::ios::floatfield); … … 249 246 ostr << **part2 << std::endl; 250 247 } 248 // restore the stream state 249 ostr.flags(orig); 250 ostr.precision(prec); 251 251 } 252 252 -
trunk/Utilities/HepMC/src/IO_AsciiParticles.cc
r349 r571 11 11 #include "IO_AsciiParticles.h" 12 12 #include "GenEvent.h" 13 #include "ParticleDataTable.h"14 13 #include "Version.h" 15 14 -
trunk/Utilities/HepMC/src/IO_GenEvent.cc
r349 r571 8 8 9 9 #include "IO_GenEvent.h" 10 #include "IO_Exception.h" 10 11 #include "GenEvent.h" 11 #include "ParticleDataTable.h" 12 #include "HeavyIon.h" 13 #include "PdfInfo.h" 14 #include "CommonIO.h" 15 #include "Version.h" 12 #include "StreamHelpers.h" 16 13 17 14 namespace HepMC { … … 23 20 m_istr(0), 24 21 m_iostr(0), 25 m_finished_first_event_io(false),26 22 m_have_file(false), 27 m_common_io() 23 m_error_type(IO_Exception::OK), 24 m_error_message() 28 25 { 29 26 if ( (m_mode&std::ios::out && m_mode&std::ios::in) || 30 27 (m_mode&std::ios::app && m_mode&std::ios::in) ) { 31 std::cerr << "IO_GenEvent::IO_GenEvent Error, open of file requested " 32 << "of input AND output type. Not allowed. Closing file."33 28 m_error_type = IO_Exception::InputAndOutput; 29 m_error_message ="IO_GenEvent::IO_GenEvent Error, open of file requested of input AND output type. Not allowed. Closing file."; 30 std::cerr << m_error_message << std::endl; 34 31 m_file.close(); 35 32 return; 36 33 } 37 // precision 16 (# digits following decimal point) is the minimum that38 // will capture the full information stored in a double39 m_file.precision(16);40 // we use decimal to store integers, because it is smaller than hex!41 m_file.setf(std::ios::dec,std::ios::basefield);42 m_file.setf(std::ios::scientific,std::ios::floatfield);43 34 // now we set the streams 44 35 m_iostr = &m_file; … … 46 37 m_istr = &m_file; 47 38 m_ostr = NULL; 39 detail::establish_input_stream_info(m_file); 48 40 } 49 41 if ( m_mode&std::ios::out ) { 50 42 m_ostr = &m_file; 51 43 m_istr = NULL; 44 detail::establish_output_stream_info(m_file); 52 45 } 53 46 m_have_file = true; … … 59 52 m_istr(&istr), 60 53 m_iostr(&istr), 61 m_finished_first_event_io(false),62 54 m_have_file(false), 63 m_common_io() 64 { } 55 m_error_type(IO_Exception::OK), 56 m_error_message() 57 { 58 detail::establish_input_stream_info( istr ); 59 } 65 60 66 61 IO_GenEvent::IO_GenEvent( std::ostream & ostr ) … … 68 63 m_istr(0), 69 64 m_iostr(&ostr), 70 m_finished_first_event_io(false),71 65 m_have_file(false), 72 m_common_io() 73 { 74 // precision 16 (# digits following decimal point) is the minimum that 75 // will capture the full information stored in a double 76 m_ostr->precision(16); 77 // we use decimal to store integers, because it is smaller than hex! 78 m_ostr->setf(std::ios::dec,std::ios::basefield); 79 m_ostr->setf(std::ios::scientific,std::ios::floatfield); 80 } 66 m_error_type(IO_Exception::OK), 67 m_error_message() 68 { 69 detail::establish_output_stream_info( ostr ); 70 } 81 71 82 72 IO_GenEvent::~IO_GenEvent() { 83 write_end_listing(); 73 if ( m_ostr != NULL ) { 74 write_HepMC_IO_block_end(*m_ostr); 75 } 84 76 if(m_have_file) m_file.close(); 85 77 } 86 78 87 void IO_GenEvent::use_input_units( Units::MomentumUnit mom, Units::LengthUnit len ) { 88 m_common_io.use_input_units(mom,len); 79 void IO_GenEvent::use_input_units( Units::MomentumUnit mom, 80 Units::LengthUnit len ) { 81 if( m_istr != NULL ) { 82 set_input_units( *m_istr, mom, len ); 83 } 89 84 } 90 85 … … 99 94 } 100 95 96 void IO_GenEvent::precision( int size ) { 97 if( size > 16 ) { 98 std::cerr << "IO_GenEvent::precision Error, " 99 << "precision is greater than 16. " 100 << "Not allowed. Using default precision of 16." 101 << std::endl; 102 size = 16; 103 } 104 if(m_ostr) { 105 m_ostr->precision(size); 106 } 107 } 108 109 bool IO_GenEvent::fill_next_event( GenEvent* evt ){ 110 // 111 // reset error type 112 m_error_type = IO_Exception::OK; 113 // 114 // test that evt pointer is not null 115 if ( !evt ) { 116 m_error_type = IO_Exception::NullEvent; 117 m_error_message = "IO_GenEvent::fill_next_event error - passed null event."; 118 std::cerr << m_error_message << std::endl; 119 return false; 120 } 121 // make sure the stream is good, and that it is in input mode 122 if ( !(*m_istr) ) return false; 123 if ( !m_istr ) { 124 m_error_type = IO_Exception::WrongFileType; 125 m_error_message = "HepMC::IO_GenEvent::fill_next_event attempt to read from output file."; 126 std::cerr << m_error_message << std::endl; 127 return false; 128 } 129 // use streaming input 130 try { 131 *m_istr >> *evt; 132 } 133 catch (IO_Exception& e) { 134 m_error_type = IO_Exception::InvalidData; 135 m_error_message = e.what(); 136 evt->clear(); 137 return false; 138 } 139 if( evt->is_valid() ) return true; 140 return false; 141 } 142 101 143 void IO_GenEvent::write_event( const GenEvent* evt ) { 102 144 /// Writes evt to output stream. It does NOT delete the event after writing. … … 105 147 if ( !evt ) return; 106 148 if ( m_ostr == NULL ) { 107 std::cerr << "HepMC::IO_GenEvent::write_event " 108 << " attempt to write to input file." << std::endl; 149 m_error_type = IO_Exception::WrongFileType; 150 m_error_message = "HepMC::IO_GenEvent::write_event attempt to write to input file."; 151 std::cerr << m_error_message << std::endl; 109 152 return; 110 153 } 111 154 // 112 155 // write event listing key before first event only. 113 if ( !m_finished_first_event_io ) { 114 m_finished_first_event_io = 1; 115 *m_ostr << "\n" << "HepMC::Version " << versionName(); 116 *m_ostr << "\n"; 117 m_common_io.write_IO_GenEvent_Key(*m_ostr); 118 } 119 // 120 // output the event data including the number of primary vertices 121 // and the total number of vertices 122 std::vector<long> random_states = evt->random_states(); 123 *m_ostr << 'E'; 124 output( evt->event_number() ); 125 output( evt->mpi() ); 126 output( evt->event_scale() ); 127 output( evt->alphaQCD() ); 128 output( evt->alphaQED() ); 129 output( evt->signal_process_id() ); 130 output( ( evt->signal_process_vertex() ? 131 evt->signal_process_vertex()->barcode() : 0 ) ); 132 output( evt->vertices_size() ); // total number of vertices. 133 write_beam_particles( evt->beam_particles() ); 134 output( (int)random_states.size() ); 135 for ( std::vector<long>::iterator rs = random_states.begin(); 136 rs != random_states.end(); ++rs ) { 137 output( *rs ); 138 } 139 output( (int)evt->weights().size() ); 140 for ( WeightContainer::const_iterator w = evt->weights().begin(); 141 w != evt->weights().end(); ++w ) { 142 output( *w ); 143 } 144 output('\n'); 145 write_unit_info( evt ); 146 write_heavy_ion( evt->heavy_ion() ); 147 write_pdf_info( evt->pdf_info() ); 148 // 149 // Output all of the vertices - note there is no real order. 150 for ( GenEvent::vertex_const_iterator v = evt->vertices_begin(); 151 v != evt->vertices_end(); ++v ) { 152 write_vertex( *v ); 153 } 154 } 155 156 bool IO_GenEvent::fill_next_event( GenEvent* evt ){ 157 // 158 // 159 // test that evt pointer is not null 160 if ( !evt ) { 161 std::cerr 162 << "IO_GenEvent::fill_next_event error - passed null event." 163 << std::endl; 164 return false; 165 } 166 // make sure the stream is good, and that it is in input mode 167 if ( !(*m_istr) ) return false; 168 if ( !m_istr ) { 169 std::cerr << "HepMC::IO_GenEvent::fill_next_event " 170 << " attempt to read from output file." << std::endl; 171 return false; 172 } 173 // 174 // search for event listing key before first event only. 175 // 176 // skip through the file just after first occurence of the start_key 177 int iotype = 0; 178 if ( !m_finished_first_event_io ) { 179 iotype = m_common_io.find_file_type(*m_istr); 180 if( iotype <= 0 ) { 181 std::cerr << "IO_GenEvent::fill_next_event start key not found " 182 << "setting badbit." << std::endl; 183 m_istr->clear(std::ios::badbit); 184 return false; 185 } 186 m_finished_first_event_io = 1; 187 } 188 // 189 // test to be sure the next entry is of type "E" then ignore it 190 if ( !(*m_istr) ) { 191 std::cerr << "IO_GenEvent::fill_next_event end of stream found " 192 << "setting badbit." << std::endl; 193 m_istr->clear(std::ios::badbit); 194 return false; 195 } 196 if ( !(*m_istr) || m_istr->peek()!='E' ) { 197 // if the E is not the next entry, then check to see if it is 198 // the end event listing key - if yes, search for another start key 199 int ioendtype = m_common_io.find_end_key(*m_istr); 200 if ( ioendtype == iotype ) { 201 iotype = m_common_io.find_file_type(*m_istr); 202 if( iotype <= 0 ) { 203 // this is the only case where we set an EOF state 204 m_istr->clear(std::ios::eofbit); 205 return false; 206 } 207 } else if ( ioendtype > 0 ) { 208 std::cerr << "IO_GenEvent::fill_next_event end key does not match start key " 209 << "setting badbit." << std::endl; 210 m_istr->clear(std::ios::badbit); 211 return false; 212 } else { 213 std::cerr << "IO_GenEvent::fill_next_event end key not found " 214 << "setting badbit." << std::endl; 215 m_istr->clear(std::ios::badbit); 216 return false; 217 } 218 } 219 m_istr->ignore(); 220 // call the appropriate read method 221 if( m_common_io.io_type() == gen ) { 222 return m_common_io.read_io_genevent(m_istr, evt); 223 } else if( m_common_io.io_type() == ascii ) { 224 return m_common_io.read_io_ascii(m_istr, evt); 225 } else if( m_common_io.io_type() == extascii ) { 226 return m_common_io.read_io_extendedascii(m_istr, evt); 227 } else if( m_common_io.io_type() == ascii_pdt ) { 228 } else if( m_common_io.io_type() == extascii_pdt ) { 229 } 230 // should not get to this statement 231 return false; 156 write_HepMC_IO_block_begin(*m_ostr); 157 // explicit cast is necessary 158 GenEvent e = *evt; 159 *m_ostr << e ; 232 160 } 233 161 … … 236 164 if ( !(*m_ostr) ) return; 237 165 if ( m_ostr == NULL ) { 238 std::cerr << "HepMC::IO_GenEvent::write_comment " 239 << " attempt to write to input file." << std::endl; 166 m_error_type = IO_Exception::WrongFileType; 167 m_error_message = "HepMC::IO_GenEvent::write_event attempt to write to input file."; 168 std::cerr << m_error_message << std::endl; 240 169 return; 241 170 } 242 171 // write end of event listing key if events have already been written 243 write_ end_listing();172 write_HepMC_IO_block_end(*m_ostr); 244 173 // insert the comment key before the comment 245 174 *m_ostr << "\n" << "HepMC::IO_GenEvent-COMMENT\n"; 246 175 *m_ostr << comment << std::endl; 247 176 } 248 249 void IO_GenEvent::write_vertex( GenVertex* v ) {250 // assumes mode has already been checked251 if ( !v || !(*m_ostr) ) {252 std::cerr << "IO_GenEvent::write_vertex !v||!(*m_ostr), "253 << "v="<< v << " setting badbit" << std::endl;254 m_ostr->clear(std::ios::badbit);255 return;256 }257 // First collect info we need258 // count the number of orphan particles going into v259 int num_orphans_in = 0;260 for ( GenVertex::particles_in_const_iterator p1261 = v->particles_in_const_begin();262 p1 != v->particles_in_const_end(); ++p1 ) {263 if ( !(*p1)->production_vertex() ) ++num_orphans_in;264 }265 //266 *m_ostr << 'V';267 output( v->barcode() ); // v's unique identifier268 output( v->id() );269 output( v->position().x() );270 output( v->position().y() );271 output( v->position().z() );272 output( v->position().t() );273 output( num_orphans_in );274 output( (int)v->particles_out_size() );275 output( (int)v->weights().size() );276 for ( WeightContainer::iterator w = v->weights().begin();277 w != v->weights().end(); ++w ) {278 output( *w );279 }280 output('\n');281 for ( GenVertex::particles_in_const_iterator p2282 = v->particles_in_const_begin();283 p2 != v->particles_in_const_end(); ++p2 ) {284 if ( !(*p2)->production_vertex() ) {285 write_particle( *p2 );286 }287 }288 for ( GenVertex::particles_out_const_iterator p3289 = v->particles_out_const_begin();290 p3 != v->particles_out_const_end(); ++p3 ) {291 write_particle( *p3 );292 }293 }294 295 void IO_GenEvent::write_beam_particles(296 std::pair<HepMC::GenParticle *,HepMC::GenParticle *> pr ) {297 GenParticle* p = pr.first;298 //m_file << 'B';299 if(!p) {300 output( 0 );301 } else {302 output( p->barcode() );303 }304 p = pr.second;305 if(!p) {306 output( 0 );307 } else {308 output( p->barcode() );309 }310 }311 312 void IO_GenEvent::write_heavy_ion( HeavyIon const * ion ) {313 // assumes mode has already been checked314 if ( !(*m_ostr) ) {315 std::cerr << "IO_GenEvent::write_heavy_ion !(*m_ostr), "316 << " setting badbit" << std::endl;317 m_ostr->clear(std::ios::badbit);318 return;319 }320 *m_ostr << 'H';321 // HeavyIon* is set to 0 by default322 if ( !ion ) {323 output( 0 );324 output( 0 );325 output( 0 );326 output( 0 );327 output( 0 );328 output( 0 );329 output( 0 );330 output( 0 );331 output( 0 );332 output( 0. );333 output( 0. );334 output( 0. );335 output( 0. );336 output('\n');337 return;338 }339 //340 output( ion->Ncoll_hard() );341 output( ion->Npart_proj() );342 output( ion->Npart_targ() );343 output( ion->Ncoll() );344 output( ion->spectator_neutrons() );345 output( ion->spectator_protons() );346 output( ion->N_Nwounded_collisions() );347 output( ion->Nwounded_N_collisions() );348 output( ion->Nwounded_Nwounded_collisions() );349 output( ion->impact_parameter() );350 output( ion->event_plane_angle() );351 output( ion->eccentricity() );352 output( ion->sigma_inel_NN() );353 output('\n');354 }355 356 void IO_GenEvent::write_pdf_info( PdfInfo const * pdf ) {357 // assumes mode has already been checked358 if ( !(*m_ostr) ) {359 std::cerr << "IO_GenEvent::write_pdf_info !(*m_ostr), "360 << " setting badbit" << std::endl;361 m_ostr->clear(std::ios::badbit);362 return;363 }364 *m_ostr << 'F';365 // PdfInfo* is set to 0 by default366 if ( !pdf ) {367 output( 0 );368 output( 0 );369 output( 0. );370 output( 0. );371 output( 0. );372 output( 0. );373 output( 0. );374 output( 0 );375 output( 0 );376 output('\n');377 return;378 }379 //380 output( pdf->id1() );381 output( pdf->id2() );382 output( pdf->x1() );383 output( pdf->x2() );384 output( pdf->scalePDF() );385 output( pdf->pdf1() );386 output( pdf->pdf2() );387 output( pdf->pdf_id1() );388 output( pdf->pdf_id2() );389 output('\n');390 }391 392 void IO_GenEvent::write_unit_info( const GenEvent* evt ) {393 if ( !(*m_ostr) ) {394 std::cerr << "IO_GenEvent::write_unit_info !(*m_ostr), "395 << " setting badbit" << std::endl;396 m_ostr->clear(std::ios::badbit);397 return;398 }399 // could write enums here, but strings are more readable400 *m_ostr << "U " << name(evt->momentum_unit());401 *m_ostr << " " << name(evt->length_unit());402 output('\n');403 }404 405 void IO_GenEvent::write_particle( GenParticle* p ) {406 // assumes mode has already been checked407 if ( !p || !(*m_ostr) ) {408 std::cerr << "IO_GenEvent::write_particle !p||!(*m_ostr), "409 << "p="<< p << " setting badbit" << std::endl;410 m_ostr->clear(std::ios::badbit);411 return;412 }413 *m_ostr << 'P';414 output( p->barcode() );415 output( p->pdg_id() );416 output( p->momentum().px() );417 output( p->momentum().py() );418 output( p->momentum().pz() );419 output( p->momentum().e() );420 output( p->generated_mass() );421 output( p->status() );422 output( p->polarization().theta() );423 output( p->polarization().phi() );424 // since end_vertex is oftentimes null, this CREATES a null vertex425 // in the map426 output( ( p->end_vertex() ? p->end_vertex()->barcode() : 0 ) );427 *m_ostr << ' ' << p->flow() << "\n";428 }429 430 void IO_GenEvent::write_particle_data( const ParticleData* pdata ) {431 // assumes mode has already been checked432 if ( !pdata || !(*m_ostr) ) {433 std::cerr << "IO_GenEvent::write_particle_data !pdata||!(*m_ostr), "434 << "pdata="<< pdata << " setting badbit" << std::endl;435 m_ostr->clear(std::ios::badbit);436 return;437 }438 *m_ostr << 'D';439 output( pdata->pdg_id() );440 output( pdata->charge() );441 output( pdata->mass() );442 output( pdata->clifetime() );443 output( (int)(pdata->spin()*2.+.1) );444 // writes the first 21 characters starting with 0445 *m_ostr << " " << pdata->name().substr(0,21) << "\n";446 }447 448 HeavyIon* IO_GenEvent::read_heavy_ion()449 {450 // assumes mode has already been checked451 //452 // test to be sure the next entry is of type "H" then ignore it453 if ( !(*m_istr) || m_istr->peek()!='H' ) {454 std::cerr << "IO_GenEvent::read_heavy_ion setting badbit." << std::endl;455 m_istr->clear(std::ios::badbit);456 return false;457 }458 m_istr->ignore();459 // read values into temp variables, then create a new HeavyIon object460 int nh =0, np =0, nt =0, nc =0,461 neut = 0, prot = 0, nw =0, nwn =0, nwnw =0;462 float impact = 0., plane = 0., xcen = 0., inel = 0.;463 *m_istr >> nh >> np >> nt >> nc >> neut >> prot464 >> nw >> nwn >> nwnw >> impact >> plane >> xcen >> inel;465 m_istr->ignore(2,'\n');466 if( nh == 0 ) return false;467 HeavyIon* ion = new HeavyIon(nh, np, nt, nc, neut, prot,468 nw, nwn, nwnw,469 impact, plane, xcen, inel );470 //471 return ion;472 }473 474 GenParticle* IO_GenEvent::read_particle(475 TempParticleMap& particle_to_end_vertex ){476 // assumes mode has already been checked477 //478 // test to be sure the next entry is of type "P" then ignore it479 if ( !(*m_istr) || m_istr->peek()!='P' ) {480 std::cerr << "IO_GenEvent::read_particle setting badbit."481 << std::endl;482 m_istr->clear(std::ios::badbit);483 return false;484 }485 m_istr->ignore();486 //487 // declare variables to be read in to, and read everything except flow488 double px = 0., py = 0., pz = 0., e = 0., m = 0., theta = 0., phi = 0.;489 int bar_code = 0, id = 0, status = 0, end_vtx_code = 0, flow_size = 0;490 *m_istr >> bar_code >> id >> px >> py >> pz >> e >> m >> status491 >> theta >> phi >> end_vtx_code >> flow_size;492 //493 // read flow patterns if any exist494 Flow flow;495 int code_index, code;496 for ( int i = 1; i <= flow_size; ++i ) {497 *m_istr >> code_index >> code;498 flow.set_icode( code_index,code);499 }500 m_istr->ignore(2,'\n'); // '\n' at end of entry501 GenParticle* p = new GenParticle( FourVector(px,py,pz,e),502 id, status, flow,503 Polarization(theta,phi) );504 p->set_generated_mass( m );505 p->suggest_barcode( bar_code );506 //507 // all particles are connected to their end vertex separately508 // after all particles and vertices have been created - so we keep509 // a map of all particles that have end vertices510 if ( end_vtx_code != 0 ) {511 particle_to_end_vertex.addEndParticle(p,end_vtx_code);512 }513 return p;514 }515 516 ParticleData* IO_GenEvent::read_particle_data( ParticleDataTable* pdt ) {517 // assumes mode has already been checked518 //519 // test to be sure the next entry is of type "D" then ignore it520 if ( !(*m_istr) || m_istr->peek()!='D' ) return false;521 m_istr->ignore();522 //523 // read values into temp variables then create new ParticleData object524 char its_name[22];525 int its_id = 0, its_spin = 0;526 double its_charge = 0, its_mass = 0, its_clifetime = 0;527 *m_istr >> its_id >> its_charge >> its_mass528 >> its_clifetime >> its_spin;529 m_istr->ignore(1); // eat the " "530 m_istr->getline( its_name, 22, '\n' );531 ParticleData* pdata = new ParticleData( its_name, its_id, its_charge,532 its_mass, its_clifetime,533 double(its_spin)/2.);534 pdt->insert(pdata);535 return pdata;536 }537 538 bool IO_GenEvent::write_end_listing() {539 if ( m_finished_first_event_io && (m_ostr != NULL) ) {540 m_common_io.write_IO_GenEvent_End(*m_ostr);541 *m_ostr << std::flush;542 m_finished_first_event_io = 0;543 return true;544 }545 return false;546 }547 177 548 178 } // HepMC -
trunk/Utilities/HepMC/src/Polarization.cc
r349 r571 9 9 namespace HepMC { 10 10 11 Polarization::Polarization( ) 12 : m_theta( 0. ), 13 m_phi( 0. ), 14 m_defined( false ) 15 { } 16 11 17 Polarization::Polarization( double theta, double phi ) 12 18 : m_theta( valid_theta(theta) ), 13 m_phi ( valid_phi(phi) ) 19 m_phi ( valid_phi(phi) ), 20 m_defined( true ) 14 21 { } 15 22 16 23 Polarization::Polarization( const Polarization& inpolar ) 17 24 : m_theta( valid_theta( inpolar.theta() ) ), 18 m_phi ( valid_phi( inpolar.phi() ) ) 25 m_phi ( valid_phi( inpolar.phi() ) ), 26 m_defined( inpolar.is_defined() ) 19 27 { } 20 28 21 29 Polarization::Polarization( const ThreeVector& vec3in ) 22 30 : m_theta( valid_theta( vec3in.theta() ) ), 23 m_phi ( valid_phi( vec3in.phi() ) ) 31 m_phi ( valid_phi( vec3in.phi() ) ), 32 m_defined( true ) 24 33 { } 25 34 … … 28 37 std::swap( m_theta, other.m_theta ); 29 38 std::swap( m_phi, other.m_phi ); 39 std::swap( m_defined, other.m_defined ); 30 40 } 31 41 … … 64 74 return m_phi = valid_phi( phi ); 65 75 } 76 77 bool Polarization::is_defined( ) const { 78 return m_defined; 79 } 80 81 void Polarization::set_undefined() { 82 m_defined = false; 83 m_theta = 0.; 84 m_phi = 0.; 85 } 66 86 67 87 void Polarization::set_theta_phi( double theta, double phi ) { 68 88 set_theta( theta ); 69 89 set_phi( phi ) ; 90 m_defined = true; 70 91 } 71 92 … … 73 94 set_theta( vec3in.theta() ); 74 95 set_phi( vec3in.phi() ); 96 m_defined = true; 75 97 return vec3in; 76 98 } -
trunk/Utilities/HepMC/src/Units.cc
r349 r571 52 52 // if this function fails to compile, rerun configure using --with-length_units 53 53 LengthUnit default_length_unit() { 54 //return @HEPMC_DEFAULT_LEN_UNIT@ ;55 54 return MM ; 56 55 } … … 58 57 // if this function fails to compile, rerun configure using --with-momentum_units 59 58 MomentumUnit default_momentum_unit() { 60 //return @HEPMC_DEFAULT_MOM_UNIT@ ; 61 return GEV; 59 return GEV ; 62 60 } 63 61
Note:
See TracChangeset
for help on using the changeset viewer.