QwGeant4
QweakSimTarget_DetectorSD.cc
Go to the documentation of this file.
1 
2 // QweakSimTarget_DetectorSD.cc
3 // Martin McHugh
4 // 2013-07-20
5 
6 /////// --------------------------------------------------------------------
7 
9 
10 //--- user classes
14 
15 /////// --------------------------------------------------------------------
16 
18 : G4VSensitiveDetector(name)
19 {
20  collectionName.insert("TargetCollection");
23 }
24 
25 /////// --------------------------------------------------------------------
26 
27 void QweakSimTarget_DetectorSD::Initialize(G4HCofThisEvent* HCE)
28 {
29  Target_DetectorHitsCollection = new QweakSimTarget_DetectorHitsCollection(SensitiveDetectorName,collectionName[0]);
30  Target_CollectionID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
31  HCE->AddHitsCollection(Target_CollectionID, Target_DetectorHitsCollection);
32 }
33 
34 /////// --------------------------------------------------------------------
35 
36 G4bool QweakSimTarget_DetectorSD::ProcessHits(G4Step* aStep, G4TouchableHistory* )
37 {
38  //--- Not Entering Geometry
39  if (aStep->GetPreStepPoint()->GetStepStatus() != fGeomBoundary ) return false;
40 
41  G4TouchableHandle theTouchable = aStep->GetPreStepPoint()->GetTouchableHandle();
42 
43  //G4double parentID = aStep->GetTrack()->GetParentID();
44  G4double trackID = aStep->GetTrack()->GetTrackID();
45 
46  G4ParticleDefinition* fpParticleDefinition = aStep->GetTrack()->GetDefinition();
47  G4String particleName = fpParticleDefinition->GetParticleName();
48  //G4double pdgCharge = fpParticleDefinition->GetPDGCharge();
49  G4int pdgEncoding = fpParticleDefinition->GetPDGEncoding();
50 
51  //GlobalTimeOfHit
52  G4double globalTime = aStep->GetPreStepPoint()->GetGlobalTime();
53 
54  G4ThreeVector worldPos = aStep->GetPreStepPoint()->GetPosition();
55  G4ThreeVector localPos = theTouchable->GetHistory()->GetTopTransform().TransformPoint(worldPos);
56 
57  // Local Origin (where did the hit come from)
58  G4ThreeVector originVertexPosition = aStep->GetTrack()->GetVertexPosition();
59 
60  G4ThreeVector worldMomentum = aStep->GetPreStepPoint()->GetMomentum();
61  G4ThreeVector localMomentum = theTouchable->GetHistory()->GetTopTransform().TransformPoint(worldMomentum);
62  G4ThreeVector currentMomentumDirection = aStep->GetTrack()->GetMomentumDirection();
63  G4ThreeVector originVertexMomentumDirection = aStep->GetTrack()->GetVertexMomentumDirection();
64 
65  G4double originVertexKineticEnergy = aStep->GetTrack()->GetVertexKineticEnergy();
66  G4double originVertexTotalEnergy = aStep->GetTrack()->GetVertexKineticEnergy();//GetVertexTotalEnergy();
67 
68  G4double currentKineticEnergy = aStep->GetTrack()->GetKineticEnergy();
69  G4double currentTotalEnergy = aStep->GetTrack()->GetTotalEnergy();
70 
71  //--- Target deposited energy
72  G4double depositedEnergy = aStep->GetTotalEnergyDeposit();
73 
74 
75  //----------------- Store Hit information
76 
78 
79  aHit -> StoreTrackID(trackID);
80 
81  aHit -> StoreParticleName(particleName);
82  aHit -> StoreParticleType(pdgEncoding);
83 
84  aHit -> StoreGlobalTime(globalTime);
85 
86  aHit -> StoreWorldPosition(worldPos);
87  aHit -> StoreLocalPosition(localPos);
88  aHit -> StoreOriginVertexPosition(originVertexPosition);
89 
90  aHit -> StoreWorldMomentum(worldMomentum);
91  aHit -> StoreLocalMomentum(localMomentum);
92  aHit -> StoreMomentumDirection(currentMomentumDirection);
93  aHit -> StoreOriginVertexMomentumDirection(originVertexMomentumDirection);
94 
95  aHit -> StoreOriginVertexKineticEnergy(originVertexKineticEnergy);
96  aHit -> StoreOriginVertexTotalEnergy(originVertexTotalEnergy); // kinetic energy
97 
98  aHit -> StoreKineticEnergy(currentKineticEnergy);
99  aHit -> StoreTotalEnergy(currentTotalEnergy);
100 
101  //--- Target deposited energy
102  aHit -> StoreDepositedEnergy(depositedEnergy);
103 
104 
105  //----------------- check if it is first touch
106 
107  if ( !(aHit -> GetLogVolume()) )
108  {
109  //--- store translation and rotation matrix of the drift cell
110  //--- for the sake of drawing the hit
111  aHit -> StoreLogVolume(theTouchable->GetVolume()->GetLogicalVolume());
112  G4AffineTransform aTrans = theTouchable->GetHistory()->GetTopTransform();
113  aTrans.Invert();
114  aHit -> StoreCellRotation(aTrans.NetRotation());
115  aHit -> StoreCellPosition(aTrans.NetTranslation());
116  }
117 
118  Target_DetectorHitsCollection->insert(aHit);
119 
120  return true;
121 }
122 
123 /////// --------------------------------------------------------------------
124 
126 {
127  //G4int NbDCHits = DC_hitsCollection->entries();
128 
129  //G4cout << "\n-------->Hits Collection: in this event they are " << NbDCHits
130  // << " hits in the Drift Cells : " << G4endl;
131  //for (G4int i=0;i<NbDCHits;i++) (
132  // *DC_hitsCollection)[i]->Print();
133 }
134 
135 /////// --------------------------------------------------------------------
QweakSimTarget_DetectorHitsCollection * Target_DetectorHitsCollection
G4bool ProcessHits(G4Step *aStep, G4TouchableHistory *ROhist)
G4THitsCollection< QweakSimTarget_DetectorHit > QweakSimTarget_DetectorHitsCollection
void Initialize(G4HCofThisEvent *HCE)
void EndOfEvent(G4HCofThisEvent *HCE)