Skip to content
Snippets Groups Projects
Commit f74efb26 authored by Ewelina Maria Lobodzinska's avatar Ewelina Maria Lobodzinska Committed by Graeme Stewart
Browse files

fix inheritance of ReadHepEvtFromAscii (TruthIO-00-00-12)

        * inherit ReadHepEvtFromAscii from AthAlgorithm
        * tag as TruthIO-00-00-12

2015-04-17 E. M. Lobodzinska <ewelina@mail.desy.de>
        * add ReadHepEvtFromAscii (Dirk Sammel)
        * tag as TruthIO-00-00-11
parent ccb37a5b
No related branches found
No related tags found
No related merge requests found
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
#include "AthenaBaseComps/AthAlgorithm.h"
#include "HepMC/IO_GenEvent.h"
class StoreGateSvc;
class ReadHepEvtFromAscii:public AthAlgorithm {
public:
ReadHepEvtFromAscii(const std::string& name, ISvcLocator* pSvcLocator);
StatusCode initialize();
StatusCode execute();
StatusCode finalize();
private:
StoreGateSvc* m_sgSvc;
HepMC::IO_GenEvent* m_ascii_in;
// Setable Properties:-
std::string m_key;
std::string m_input_file;
};
......@@ -6,8 +6,10 @@ use AtlasPolicy AtlasPolicy-*
use GaudiInterface GaudiInterface-* External
use GeneratorModules GeneratorModules-* Generators
use AtlasHepMC AtlasHepMC-* External
use AthenaBaseComps AthenaBaseComps-* Control
private
use StoreGate StoreGate-* Control
use GeneratorObjects GeneratorObjects-* Generators
use EventInfo EventInfo-* Event
use HepPDT v* LCG_Interfaces
......@@ -15,10 +17,11 @@ end_private
private
#macro_append HepMC_linkopts " $(HepMC_IO_linkopts) "
macro_append TruthIO_shlibflags "$(HepMC_IO_linkopts)"
end_private
public
library TruthIO *.cxx -s=components *.cxx
library TruthIO *.cxx *.F -s=components *.cxx
apply_pattern component_library
apply_pattern declare_joboptions files="*.py"
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
#include <set>
#include "TruthIO/ReadHepEvtFromAscii.h"
#include "GeneratorObjects/McEventCollection.h"
#include "HepMC/GenEvent.h"
#include "HepMC/IO_HEPEVT.h"
#include "HepMC/HEPEVT_Wrapper.h"
#include "GaudiKernel/DataSvc.h"
#include "StoreGate/StoreGateSvc.h"
extern "C" {
void openfile_(int & ifile, const char* filename, int);
void readevt_(int& ifile);
void closefile_(int& ifile);
}
ReadHepEvtFromAscii::ReadHepEvtFromAscii(const std::string& name, ISvcLocator* pSvcLocator)
: AthAlgorithm(name, pSvcLocator),
m_sgSvc(0), m_ascii_in(0)
{
// Set users' request
declareProperty("McEventKey", m_key = "GEN_EVENT");
declareProperty("HepEvtFile", m_input_file = "events.ascii");
}
StatusCode ReadHepEvtFromAscii::initialize(){
msg(MSG::INFO) << ">>> ReadHepEvtFromAscii from Initialize" << endreq;
StatusCode sc = service("StoreGateSvc", m_sgSvc);
if (sc.isFailure()) {
msg(MSG::ERROR) << "Could not find StoreGateSvc" << endreq;
return sc;
}
HepMC::HEPEVT_Wrapper::set_sizeof_int(4);
HepMC::HEPEVT_Wrapper::set_sizeof_real(8);
HepMC::HEPEVT_Wrapper::set_max_number_entries(10000);
// Initialize input file
// m_ascii_in = new HepMC::IO_Ascii(m_input_file.c_str(), std::ios::in);
int ifile=5;
closefile_(ifile);
openfile_(ifile,m_input_file.c_str(),m_input_file.size());
// Initialization terminated
return StatusCode::SUCCESS;
}
StatusCode ReadHepEvtFromAscii::execute() {
msg(MSG::INFO) << ">>> ReadHepEvtFromAscii from execute" << endreq;
McEventCollection* mcEvtColl;
if ( m_sgSvc->contains<McEventCollection>(m_key) && m_sgSvc->retrieve(mcEvtColl, m_key).isSuccess() ) {
if (msgLvl(MSG::VERBOSE)) msg(MSG::VERBOSE) << "found an McEventCollecion in store" << endreq;
} else {
// McCollection doesn't exist. Create it (empty)
if (msgLvl(MSG::VERBOSE)) msg(MSG::VERBOSE) << "create new McEventCollecion in store" << endreq;
mcEvtColl = new McEventCollection;
StatusCode status = m_sgSvc->record( mcEvtColl, m_key );
if (status.isFailure()) {
msg(MSG::ERROR) << "Could not record McEventCollection" << endreq;
return status;
}
}
HepMC::GenEvent* evt = new HepMC::GenEvent();
int ifile=5;
readevt_(ifile);
HepMC::IO_HEPEVT hepio;
hepio.set_print_inconsistency_errors(0);
hepio.fill_next_event(evt);
// Convert GeV to MeV
for ( HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p )
{
HepMC::FourVector newMomentum(0.,0.,0.,0.);
newMomentum.setPx( (*p)->momentum().px() * 1000. );
newMomentum.setPy( (*p)->momentum().py() * 1000. );
newMomentum.setPz( (*p)->momentum().pz() * 1000. );
newMomentum.setE( (*p)->momentum().e() * 1000. );
(*p)->set_momentum(newMomentum);
}
mcEvtColl->push_back(evt);
// End of execution for each event
return StatusCode::SUCCESS;
}
StatusCode ReadHepEvtFromAscii::finalize() {
msg(MSG::INFO) << ">>> ReadHepEvtFromAscii from finalize" << endreq;
int ifile=5;
closefile_(ifile);
// End of finalization step
return StatusCode::SUCCESS;
}
......@@ -4,12 +4,14 @@
#include "TruthIO/WriteHepMC.h"
#include "TruthIO/PrintMC.h"
#include "TruthIO/PrintHijingPars.h"
#include "TruthIO/ReadHepEvtFromAscii.h"
DECLARE_ALGORITHM_FACTORY(DumpMC)
DECLARE_ALGORITHM_FACTORY(HepMCReadFromFile)
DECLARE_ALGORITHM_FACTORY(PrintMC)
DECLARE_ALGORITHM_FACTORY(WriteHepMC)
DECLARE_ALGORITHM_FACTORY(PrintHijingPars)
DECLARE_ALGORITHM_FACTORY(ReadHepEvtFromAscii)
DECLARE_FACTORY_ENTRIES(GeneratorModules) {
DECLARE_ALGORITHM(DumpMC)
......@@ -17,4 +19,5 @@ DECLARE_FACTORY_ENTRIES(GeneratorModules) {
DECLARE_ALGORITHM(PrintMC)
DECLARE_ALGORITHM(WriteHepMC)
DECLARE_ALGORITHM(PrintHijingPars)
DECLARE_ALGORITHM(ReadHepEvtFromAscii)
}
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
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment