Added reader class for HepMC and updated library dependencies of EVGEN and those...
authorbthorsbr <brian.peter.thorsbro@cern.ch>
Mon, 14 Jul 2014 13:44:50 +0000 (15:44 +0200)
committerhristov <Peter.Hristov@cern.ch>
Fri, 8 Aug 2014 14:37:53 +0000 (16:37 +0200)
ALIROOT/CMakebinaliroot.pkg
EVE/CMakebinalieve.pkg
EVGEN/AliGenReaderHepMC.cxx [new file with mode: 0644]
EVGEN/AliGenReaderHepMC.h [new file with mode: 0644]
EVGEN/CMakelibEVGEN.pkg
EVGEN/EVGENLinkDef.h

index dfcd4db..15eaee9 100644 (file)
@@ -34,7 +34,7 @@ set ( ELIBSDIR )
 set ( ELIBS  MUONevaluation    MUONmapping MUONshuttle MUONgraphics MUONsim MUONrec MUONgeometry MUONcalib MUONbase MUONraw MUONtrigger MUONcore TPCbase TPCsim TPCrec TPCutil  ITSbase ITSsim ITSrec PMDbase PMDsim PMDrec TRDbase TRDsim TRDrec FMDbase FMDsim FMDrec TOFbase TOFrec TOFsim PHOSUtils PHOSbase PHOSsim PHOSrec ADbase ADsim ADrec ACORDEbase ACORDEsim ACORDErec HMPIDbase HMPIDrec HMPIDsim ZDCbase ZDCsim ZDCrec VZERObase VZEROsim VZEROrec MFTbase MFTsim MFTrec EMCALUtils EMCALbase EMCALsim EMCALrec EMCALraw BCM STRUCT T0base T0sim T0rec FASTSIM microcern HLTbase HLTshuttle TRIGGERbase STEER STAT CDB AOD  STEERBase ESD ANALYSIS RAWDatasim RAWDatarec RAWDatabase)
 
 if(PYTHIA6)
-  list(APPEND ELIBS EVGEN)
+  list(APPEND ELIBS EVGEN HepMC HepMCParser)
 endif(PYTHIA6)
 
 if( ALICE_TARGET STREQUAL "macosx")
index 38c81d5..dcfc9ac 100644 (file)
@@ -35,7 +35,7 @@ set ( ELIBSDIR )
 set ( ELIBS  EveBase EveDet MUONmapping MUONevaluation MUONsim MUONrec MUONgeometry MUONcalib MUONbase MUONraw MUONtrigger MUONcore TPCsim TPCrec TPCbase TPCutil ITSbase ITSsim ITSrec ITSUpgradeBase ITSUpgradeSim PMDbase PMDsim PMDrec TRDbase TRDsim TRDrec FMDbase FMDsim FMDrec TOFbase TOFrec TOFsim PHOSUtils PHOSbase PHOSsim PHOSrec ACORDEbase ACORDEsim ACORDErec HMPIDbase HMPIDrec HMPIDsim ZDCbase ZDCsim ZDCrec VZERObase VZEROsim VZEROrec EMCALUtils EMCALbase EMCALsim EMCALrec EMCALraw BCM STRUCT T0base T0sim T0rec STEER CDB AOD ESD STEERBase ANALYSISalice ANALYSIS THijing hijing THbtp TEPEMGEN FASTSIM microcern RAWDatasim RAWDatarec RAWDatabase HLTbase XMLParser STAT OADB)
 
 if(PYTHIA6)
-  list (APPEND ELIBS EVGEN)
+  list (APPEND ELIBS EVGEN HepMC HepMCParser)
 endif(PYTHIA6)
 
 set ( PACKBLIBS  ${ROOTCLIBS} ${ROOTPLIBS} -lTreePlayer -lGeomPainter -lGed -lRGL -lEve ${SYSLIBS} )
diff --git a/EVGEN/AliGenReaderHepMC.cxx b/EVGEN/AliGenReaderHepMC.cxx
new file mode 100644 (file)
index 0000000..0bd0367
--- /dev/null
@@ -0,0 +1,69 @@
+#include <TVirtualMC.h>
+#include <TDatabasePDG.h>
+#include <TParticle.h>
+
+#include "AliGenReaderHepMC.h"
+#include "AliRun.h"
+#include "AliStack.h"
+
+
+
+
+ClassImp(AliGenReaderHepMC)
+
+
+AliGenReaderHepMC& AliGenReaderHepMC::operator=(const  AliGenReaderHepMC& rhs)
+{
+   // Assignment operator
+   rhs.Copy(*this);
+   return *this;
+}
+
+void AliGenReaderHepMC::Copy(TObject&) const
+{
+   //
+   // Copy
+   //
+   Fatal("Copy","Not implemented!\n");
+}
+
+void AliGenReaderHepMC::Init()
+{
+   // Initialisation
+   fEventsHandle = new HepMC::IO_GenEvent(fFileName, std::ios::in);
+   fParticleArray = new TClonesArray("TParticle");
+   fParticleIterator = new TIter(fParticleArray);
+
+
+
+}
+
+Int_t AliGenReaderHepMC::NextEvent()
+{
+   // Read the next event
+   if ((fGenEvent = fEventsHandle->read_next_event())) {
+      THepMCParser::ParseGenEvent2TCloneArray(fGenEvent,fParticleArray,false);
+      fParticleIterator->Reset();
+
+
+      // implement header... somewhere
+
+      return fGenEvent->particles_size();
+   }
+   return 0;
+}
+
+TParticle* AliGenReaderHepMC::NextParticle()
+{
+   // Read next particle
+   TParticle * particle = (TParticle*)fParticleIterator->Next();
+   if (particle && particle->GetStatusCode()==1) {
+      particle->SetBit(kTransportBit);
+   }
+   return particle;
+}
+
+void AliGenReaderHepMC::RewindEvent()
+{
+   fParticleIterator->Reset();
+}
diff --git a/EVGEN/AliGenReaderHepMC.h b/EVGEN/AliGenReaderHepMC.h
new file mode 100644 (file)
index 0000000..88b22c4
--- /dev/null
@@ -0,0 +1,49 @@
+#ifndef ALIGENREADERHEPMC_H
+#define ALIGENREADERHEPMC_H
+// Realisations of the AliGenReader interface to be used with AliGenExFile.
+// NextEvent() loops over events
+// and NextParticle() loops over particles.
+// This implementation reads HepMC output formats
+// Author: brian.peter.thorsbro@cern.ch, brian@thorsbro.dk
+// Based on AliGenReaderSL by andreas.morsch@cern.ch
+
+#include <TClonesArray.h>
+
+#include "AliGenReader.h"
+#include "AliGenEventHeader.h"
+#include "THepMCParser.h"
+#include "HepMC/IO_BaseClass.h"
+#include "HepMC/GenEvent.h"
+
+class TParticle;
+
+class AliGenReaderHepMC : public AliGenReader
+{
+public:
+   inline AliGenReaderHepMC():fEventsHandle(0), fGenEvent(0), fParticleArray(0), fParticleIterator(0), fGenEventHeader(0) {;}
+   AliGenReaderHepMC(const AliGenReaderHepMC &reader)
+   :AliGenReader(reader), fEventsHandle(0), fGenEvent(0), fParticleArray(0), fParticleIterator(0), fGenEventHeader(0) {reader.Copy(*this);}
+   inline virtual ~AliGenReaderHepMC(){ delete fEventsHandle; delete fGenEvent; delete fParticleArray; delete fParticleIterator;} // not deleting fGenEventHeader as it is returned out
+   AliGenEventHeader * GetGenEventHeader() const {return fGenEventHeader;};
+   virtual void Init();
+   virtual Int_t NextEvent();
+   virtual TParticle* NextParticle();
+   virtual void RewindEvent();
+   AliGenReaderHepMC & operator=(const AliGenReaderHepMC & rhs);
+
+protected:
+   HepMC::IO_BaseClass * fEventsHandle;   // pointer to the HepMC file handler
+   HepMC::GenEvent * fGenEvent;           // pointer to a generated event
+   TClonesArray * fParticleArray;         // pointer to array containing particles of current event
+   TIter * fParticleIterator;             // iterator coupled to the array
+   AliGenEventHeader * fGenEventHeader;   // AliGenEventHeader
+
+private:
+   void Copy(TObject&) const;
+
+   ClassDef(AliGenReaderHepMC, 1) //Generate particles from external file
+};
+#endif
+
+
+
index 987c42e..52cf3e6 100644 (file)
@@ -86,7 +86,8 @@ set ( SRCS
   AliGenMUONLMR.cxx
   AliGenLcLib.cxx
   AliGenITSULib.cxx
-  AliGenTunedOnPbPb.cxx  
+  AliGenTunedOnPbPb.cxx
+  AliGenReaderHepMC.cxx  
   )
 
 string ( REPLACE ".cxx" ".h" HDRS "${SRCS}" )
@@ -101,4 +102,6 @@ set ( EXPORT
   AliOmegaDalitz.h
   )
 
-set ( EINCLUDE FASTSIM THijing PYTHIA6 LHAPDF STEER/STEER STEER/ESD STEER/STEERBase)
+set ( EINCLUDE FASTSIM THijing PYTHIA6 LHAPDF STEER/STEER STEER/ESD STEER/STEERBase TEvtGen)
+
+# set ( ELIBS ${ROOTCLIBS} -lHepMC -lHepMCParser)
index 8a140da..2969cd7 100644 (file)
@@ -71,4 +71,6 @@
 #pragma link C++ class  AliGenLcLib+;
 #pragma link C++ class  AliGenITSULib;
 #pragma link C++ class  AliGenTunedOnPbPb;
+#pragma link C++ class  AliGenReaderHepMC++;
+
 #endif