File indexing completed on 2024-04-06 12:30:11
0001 #ifndef SimG4CMS_ShowerLibraryProducer_HFShowerG4Hit_h
0002 #define SimG4CMS_ShowerLibraryProducer_HFShowerG4Hit_h
0003
0004 #include "SimDataFormats/CaloHit/interface/HFShowerPhoton.h"
0005 #include "DataFormats/Math/interface/Point3D.h"
0006
0007 #include "G4VHit.hh"
0008 #include "G4THitsCollection.hh"
0009 #include "G4Allocator.hh"
0010 #include "G4ThreeVector.hh"
0011 #include "G4LogicalVolume.hh"
0012
0013 #include <vector>
0014
0015 class HFShowerG4Hit : public G4VHit {
0016 public:
0017 HFShowerG4Hit();
0018 HFShowerG4Hit(G4int hitId, G4int tkID, double edep, double time);
0019 ~HFShowerG4Hit() override;
0020 HFShowerG4Hit(const HFShowerG4Hit& right);
0021 const HFShowerG4Hit& operator=(const HFShowerG4Hit& right);
0022 G4int operator==(const HFShowerG4Hit& right) const;
0023
0024 inline void* operator new(size_t);
0025 inline void operator delete(void* aHit);
0026
0027 private:
0028 G4int theHitId;
0029 G4int theTrackId;
0030 G4double theEdep;
0031 G4double theTime;
0032 G4ThreeVector localPos;
0033 G4ThreeVector globalPos;
0034 G4ThreeVector momDir;
0035
0036 public:
0037 inline void setHitId(G4int hitId) { theHitId = hitId; }
0038 inline void setTrackId(G4int trackId) { theTrackId = trackId; }
0039 inline void setEnergy(G4double edep) { theEdep = edep; }
0040 inline void updateEnergy(G4double edep) { theEdep += edep; }
0041 inline void setTime(G4double t) { theTime = t; }
0042 inline void setLocalPos(const G4ThreeVector& xyz) { localPos = xyz; }
0043 inline void setGlobalPos(const G4ThreeVector& xyz) { globalPos = xyz; }
0044 inline void setPrimMomDir(const G4ThreeVector& xyz) { momDir = xyz; }
0045
0046 inline G4int hitId() const { return theHitId; }
0047 inline G4int trackId() const { return theTrackId; }
0048 inline G4double edep() const { return theEdep; };
0049 inline G4double time() const { return theTime; }
0050 inline G4ThreeVector localPosition() const { return localPos; }
0051 inline G4ThreeVector globalPosition() const { return globalPos; }
0052 inline G4ThreeVector primaryMomDir() const { return momDir; }
0053 };
0054
0055 typedef G4THitsCollection<HFShowerG4Hit> HFShowerG4HitsCollection;
0056
0057 extern G4ThreadLocal G4Allocator<HFShowerG4Hit>* fHFShowerG4HitAllocator;
0058
0059 inline void* HFShowerG4Hit::operator new(size_t) {
0060 if (!fHFShowerG4HitAllocator)
0061 fHFShowerG4HitAllocator = new G4Allocator<HFShowerG4Hit>;
0062 return (void*)fHFShowerG4HitAllocator->MallocSingle();
0063 }
0064
0065 inline void HFShowerG4Hit::operator delete(void* aHit) { fHFShowerG4HitAllocator->FreeSingle((HFShowerG4Hit*)aHit); }
0066 #endif