diff --git a/Generators/TruthIO/CMakeLists.txt b/Generators/TruthIO/CMakeLists.txt index 401a857bb29a813e356bdc46cf30fe9eabe9bffd..bd91ec10b475e6d27a587cd53c87d17ba905e0cc 100644 --- a/Generators/TruthIO/CMakeLists.txt +++ b/Generators/TruthIO/CMakeLists.txt @@ -22,7 +22,6 @@ find_package( HepPDT ) # Component(s) in the package: atlas_add_component( TruthIO src/*.cxx - src/*.F src/components/*.cxx INCLUDE_DIRS ${HEPPDT_INCLUDE_DIRS} LINK_LIBRARIES ${HEPPDT_LIBRARIES} AtlasHepMCLib AtlasHepMCfioLib AthenaBaseComps GaudiKernel GeneratorModulesLib StoreGateLib SGtests EventInfo GeneratorObjects ) diff --git a/Generators/TruthIO/TruthIO/ReadHepEvtFromAscii.h b/Generators/TruthIO/TruthIO/ReadHepEvtFromAscii.h index 87992b5d260b8173d2b6c29c80314c8d11d01e08..689a2a56df3912f43d1187be390dc2c7d9238fb5 100644 --- a/Generators/TruthIO/TruthIO/ReadHepEvtFromAscii.h +++ b/Generators/TruthIO/TruthIO/ReadHepEvtFromAscii.h @@ -21,5 +21,8 @@ private: // Setable Properties:- std::string m_key; std::string m_input_file; + std::ifstream m_file; + bool read_hepevt_particle( int i); + bool read_hepevt_event_header(); }; diff --git a/Generators/TruthIO/src/HepEvt.cxx b/Generators/TruthIO/src/HepEvt.cxx new file mode 100644 index 0000000000000000000000000000000000000000..726b2fd747cbcad1546c812569c83d3dc7fec003 --- /dev/null +++ b/Generators/TruthIO/src/HepEvt.cxx @@ -0,0 +1,13 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + + +extern "C" { +const unsigned int hepevt_bytes_allocation_ATLAS = + sizeof(long int) * ( 2 + 6 * 10000 ) + + sizeof(double) * ( 9 * 10000 ); +struct { + char data[hepevt_bytes_allocation_ATLAS]; + } hepevt_; +} \ No newline at end of file diff --git a/Generators/TruthIO/src/ReadHepEvtFromAscii.cxx b/Generators/TruthIO/src/ReadHepEvtFromAscii.cxx index 06831e5d4a22c5c1d296edea493defe0feb9fa65..795d640bce94d4e01aa5449d6758e8e56ef9fa3e 100644 --- a/Generators/TruthIO/src/ReadHepEvtFromAscii.cxx +++ b/Generators/TruthIO/src/ReadHepEvtFromAscii.cxx @@ -14,13 +14,63 @@ #include "StoreGate/StoreGateSvc.h" -extern "C" { - void openfile_(int & ifile, const char* filename, int); - void readevt_(int& ifile); - void closefile_(int& ifile); -} +bool ReadHepEvtFromAscii::read_hepevt_event_header() +{ + const size_t max_e_buffer_size=512; + char buf_e[max_e_buffer_size]; + bool eventline=false; + int im=0, pm=0; + while(!eventline) + { + m_file.getline(buf_e,max_e_buffer_size); + if( strlen(buf_e) == 0 ) return false; + std::stringstream st_e(buf_e); + char attr=' '; + eventline=false; + while (!eventline) + { + if (!(st_e>>attr)) break; + if (attr==' ') continue; + else eventline=false; + if (attr=='E') + { + eventline=static_cast<bool>(st_e>>im>>pm); + } + } + } + HepMC::HEPEVT_Wrapper::set_event_number(im); + HepMC::HEPEVT_Wrapper::set_number_entries(pm); + return eventline; +} +bool ReadHepEvtFromAscii::read_hepevt_particle( int i) +{ + const size_t max_p_buffer_size=512; + const size_t max_v_buffer_size=512; + char buf_p[max_p_buffer_size]; + char buf_v[max_v_buffer_size]; + int intcodes[6]; + double fltcodes1[5]; + double fltcodes2[4]; + m_file.getline(buf_p,max_p_buffer_size); + if( strlen(buf_p) == 0 ) return false; + m_file.getline(buf_v,max_v_buffer_size); + if( strlen(buf_v) == 0 ) return false; + std::stringstream st_p(buf_p); + std::stringstream st_v(buf_v); + if (!static_cast<bool>(st_p>>intcodes[0]>>intcodes[1]>>intcodes[2]>>intcodes[3]>>intcodes[4]>>intcodes[5]>>fltcodes1[0]>>fltcodes1[1]>>fltcodes1[2]>>fltcodes1[3]>>fltcodes1[4])) { return false;} + if (!static_cast<bool>(st_v>>fltcodes2[0]>>fltcodes2[1]>>fltcodes2[2]>>fltcodes2[3])) { return false;} + HepMC::HEPEVT_Wrapper::set_status(i,intcodes[0]); + HepMC::HEPEVT_Wrapper::set_id(i,intcodes[1]); + HepMC::HEPEVT_Wrapper::set_parents(i,intcodes[2],std::max(intcodes[2],intcodes[3]));/* Pythia writes second mother 0*/ + HepMC::HEPEVT_Wrapper::set_children(i,intcodes[4],intcodes[5]); + HepMC::HEPEVT_Wrapper::set_momentum(i,fltcodes1[0],fltcodes1[1],fltcodes1[2],fltcodes1[3]); + HepMC::HEPEVT_Wrapper::set_mass(i,fltcodes1[4]); + HepMC::HEPEVT_Wrapper::set_position(i,fltcodes2[0],fltcodes2[1],fltcodes2[2],fltcodes2[3]); + return true; +} + ReadHepEvtFromAscii::ReadHepEvtFromAscii(const std::string& name, ISvcLocator* pSvcLocator) : AthAlgorithm(name, pSvcLocator), m_sgSvc(0) @@ -46,10 +96,8 @@ StatusCode ReadHepEvtFromAscii::initialize(){ HepMC::HEPEVT_Wrapper::set_max_number_entries(10000); // Initialize input file - int ifile=5; - closefile_(ifile); - openfile_(ifile,m_input_file.c_str(),m_input_file.size()); - + m_file=std::ifstream(m_input_file.c_str()); + if( !m_file.is_open() ) return StatusCode::FAILURE; // Initialization terminated return StatusCode::SUCCESS; } @@ -73,8 +121,10 @@ StatusCode ReadHepEvtFromAscii::execute() { } HepMC::GenEvent* evt = new HepMC::GenEvent(); - int ifile=5; - readevt_(ifile); + HepMC::HEPEVT_Wrapper::zero_everything(); + bool fileok=read_hepevt_event_header(); + for (int i=1; (i<=HepMC::HEPEVT_Wrapper::number_entries())&&fileok; i++) fileok=read_hepevt_particle(i); + if (!fileok) return StatusCode::FAILURE; HepMC::IO_HEPEVT hepio; hepio.set_print_inconsistency_errors(0); @@ -100,8 +150,7 @@ StatusCode ReadHepEvtFromAscii::finalize() { msg(MSG::INFO) << ">>> ReadHepEvtFromAscii from finalize" << endmsg; - int ifile=5; - closefile_(ifile); + if(m_file.is_open()) m_file.close(); // End of finalization step return StatusCode::SUCCESS; diff --git a/Generators/TruthIO/src/readEvt.F b/Generators/TruthIO/src/readEvt.F deleted file mode 100644 index 1e337806f881d90a252e03f2741633c4112ceffc..0000000000000000000000000000000000000000 --- a/Generators/TruthIO/src/readEvt.F +++ /dev/null @@ -1,37 +0,0 @@ - subroutine openFile(ifile,filename) - integer ifile - character*(*) filename - logical ex - - inquire(FILE=filename,EXIST=ex) - if(.not. ex) write(*,*) 'File does not exist:',filename - open(ifile,FILE=filename,status='old') - end - - subroutine readEvt(ifile) - integer ifile -c...Hepevt common block - integer nmxhep,nevhep,nhep,isthep,idhep,jmohep,jdahep - double precision phep,vhep - PARAMETER (NMXHEP=10000) - COMMON/HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP), - & JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP) - - integer i - - read(ifile,*) nevhep,nhep - do i=1,nhep - read(ifile,*) dummy,isthep(i),idhep(i),jmohep(1,i), - & jmohep(2,i),jdahep(1,i),jdahep(2,i) - read(ifile,*) phep(1,i),phep(2,i),phep(3,i),phep(4,i),phep(5,i) - read(ifile,*) vhep(1,i),vhep(2,i),vhep(3,i),vhep(4,i) - enddo - - end - - subroutine closeFile(ifile) - integer ifile - - close(ifile) - - end