File indexing completed on 2023-03-17 11:24:27
0001
0002 #ifndef FP420_FP420SD_h
0003 #define FP420_FP420SD_h
0004
0005
0006 #include "SimG4Core/Notification/interface/Observer.h"
0007 #include "SimG4Core/SensitiveDetector/interface/SensitiveTkDetector.h"
0008
0009 #include "SimG4Core/Notification/interface/BeginOfRun.h"
0010 #include "SimG4Core/Notification/interface/BeginOfEvent.h"
0011 #include "SimG4Core/Notification/interface/EndOfEvent.h"
0012
0013 #include "SimG4Core/Notification/interface/BeginOfTrack.h"
0014 #include "SimG4Core/Notification/interface/BeginOfJob.h"
0015
0016 #include "DataFormats/GeometryVector/interface/LocalPoint.h"
0017
0018 #include "SimG4CMS/FP420/interface/FP420G4Hit.h"
0019 #include "SimG4CMS/FP420/interface/FP420G4HitCollection.h"
0020 #include "SimG4CMS/FP420/interface/FP420NumberingScheme.h"
0021
0022 #include "FWCore/Framework/interface/Frameworkfwd.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSetfwd.h"
0024
0025 #include "G4Step.hh"
0026 #include "G4StepPoint.hh"
0027 #include "G4Track.hh"
0028 #include "G4VPhysicalVolume.hh"
0029
0030 #include <string>
0031
0032 class TrackingSlaveSD;
0033
0034 class FrameRotation;
0035 class TrackInformation;
0036 class SimTrackManager;
0037 class TrackingSlaveSD;
0038 class UpdatablePSimHit;
0039 class G4ProcessTypeEnumerator;
0040 class G4TrackToParticleID;
0041
0042
0043
0044 class FP420SD : public SensitiveTkDetector,
0045 public Observer<const BeginOfRun*>,
0046 public Observer<const BeginOfEvent*>,
0047 public Observer<const EndOfEvent*> {
0048 public:
0049 FP420SD(const std::string&, const SensitiveDetectorCatalog&, edm::ParameterSet const&, const SimTrackManager*);
0050
0051 ~FP420SD() override;
0052
0053 bool ProcessHits(G4Step*, G4TouchableHistory*) override;
0054 uint32_t setDetUnitId(const G4Step*) override;
0055
0056 void Initialize(G4HCofThisEvent* HCE) override;
0057 void EndOfEvent(G4HCofThisEvent* eventHC) override;
0058 void clear() override;
0059 void DrawAll() override;
0060 void PrintAll() override;
0061
0062 virtual double getEnergyDeposit(G4Step* step);
0063 void fillHits(edm::PSimHitContainer&, const std::string&) override;
0064 void clearHits() override;
0065
0066 private:
0067 void update(const BeginOfRun*) override;
0068 void update(const BeginOfEvent*) override;
0069 void update(const ::EndOfEvent*) override;
0070
0071 G4ThreeVector SetToLocal(const G4ThreeVector& global);
0072 G4ThreeVector SetToLocalExit(const G4ThreeVector& globalPoint);
0073 void GetStepInfo(G4Step* aStep);
0074 G4bool HitExists();
0075 void CreateNewHit();
0076 void UpdateHit();
0077 void StoreHit(FP420G4Hit*);
0078 void ResetForNewPrimary();
0079 void Summarize();
0080
0081 private:
0082
0083 TrackingSlaveSD* slave;
0084 FP420NumberingScheme* numberingScheme;
0085
0086 G4ThreeVector entrancePoint, exitPoint;
0087 G4ThreeVector theEntryPoint;
0088 G4ThreeVector theExitPoint;
0089
0090 float incidentEnergy;
0091
0092 G4int hcID;
0093 FP420G4HitCollection* theHC;
0094 const SimTrackManager* theManager;
0095
0096 G4int tsID;
0097 FP420G4Hit* currentHit;
0098 G4Track* theTrack;
0099 G4VPhysicalVolume* currentPV;
0100
0101 uint32_t unitID, previousUnitID;
0102 G4int tSliceID;
0103 unsigned int primaryID, primID;
0104
0105 G4double tSlice;
0106
0107 G4StepPoint* preStepPoint;
0108 G4StepPoint* postStepPoint;
0109 float edeposit;
0110
0111 G4ThreeVector hitPoint;
0112 G4ThreeVector hitPointExit;
0113 G4ThreeVector hitPointLocal;
0114 G4ThreeVector hitPointLocalExit;
0115 float Pabs;
0116 float Tof;
0117 float Eloss;
0118 short ParticleType;
0119
0120 float ThetaAtEntry;
0121 float PhiAtEntry;
0122
0123 int ParentId;
0124 float Vx, Vy, Vz;
0125 float X, Y, Z;
0126
0127
0128
0129
0130 int eventno;
0131
0132 protected:
0133 float edepositEM, edepositHAD;
0134 G4int emPDG;
0135 G4int epPDG;
0136 G4int gammaPDG;
0137 };
0138
0139 #endif