| Geant4 Cross Reference |
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // Author: F. Poignant, floriane.poignant@gmail.com
27 //
28 // file STCyclotronSensitiveTarget.cc
29 //
30 #include "STCyclotronRun.hh"
31 #include "STCyclotronSensitiveTarget.hh"
32
33 #include "G4AnalysisManager.hh"
34 #include "G4RunManager.hh"
35 #include "G4HCofThisEvent.hh"
36 #include "G4UnitsTable.hh"
37 #include "G4Step.hh"
38 #include "G4SteppingManager.hh"
39 #include "G4ThreeVector.hh"
40 #include "G4SDManager.hh"
41 #include "G4ios.hh"
42 #include "G4SystemOfUnits.hh"
43 #include "G4ThreeVector.hh"
44 #include "G4Track.hh"
45 #include "G4ParticleDefinition.hh"
46 #include "G4DecayTable.hh"
47 #include "G4VDecayChannel.hh"
48 #include "G4TrackVector.hh"
49 #include "G4VProcess.hh"
50 #include "G4Tubs.hh"
51 #include <map>
52
53 STCyclotronSensitiveTarget::STCyclotronSensitiveTarget(G4String name,
54 STCyclotronDetectorConstruction* det)
55 : G4VSensitiveDetector(name),
56 fDet(det)
57 {
58 fTempTrack = 0;
59 fTempTrack1 = 0;
60 fTempEnergy = 0.;
61 fTempVector = G4ThreeVector(0.,0.,0.);
62 fTrack=0;
63 }
64
65 STCyclotronSensitiveTarget::~STCyclotronSensitiveTarget()
66 {
67 delete fTrack;
68 }
69
70 G4bool STCyclotronSensitiveTarget::ProcessHits(G4Step* aStep, G4TouchableHistory*)
71 {
72
73 STCyclotronRun* fRun = static_cast<STCyclotronRun*>(G4RunManager::GetRunManager()->GetNonConstCurrentRun());
74 fTrack = aStep->GetTrack();
75
76 auto analysisManager = G4AnalysisManager::Instance();
77
78 //----------------------------------------------
79 // Volume info
80 //----------------------------------------------
81 G4double targetHalfDiameter= (fDet->GetTargetDiameter())/2.;
82
83 //----------------------------------------------
84 // Step information
85 //----------------------------------------------
86
87 G4double edep = aStep->GetTotalEnergyDeposit();
88 G4double energy = aStep->GetPreStepPoint()->GetKineticEnergy();
89 G4ThreeVector momentumDirection = aStep->GetPreStepPoint()->GetMomentumDirection();
90 G4ThreeVector vectorPosition = aStep->GetPreStepPoint()->GetPosition();
91
92 //----------------------------------------------
93 // Track
94 //----------------------------------------------
95
96 G4ParticleDefinition* thePartDef = fTrack->GetDefinition();
97 G4String partType= fTrack->GetDefinition()->GetParticleType();
98 G4String name = fTrack->GetDefinition()->GetParticleName();
99 G4double timeLife = fTrack->GetDefinition()->GetPDGLifeTime(); //<---
100 const G4VProcess* process = fTrack->GetCreatorProcess();
101
102 //----------------------------------------------
103 // Collect general information concerning all of the particles
104 // Collect energy deposition ; separe decay case to beam case
105 //----------------------------------------------
106
107 fRun->EnergyDepositionTarget(edep);
108
109
110 //----------------------------------------------
111 //Collect information about protons and deuterons
112 //----------------------------------------------
113
114 if(name == "proton" || name == "deuteron")
115 {
116
117 if(fTrack->GetTrackID()!=fTempTrack && (momentumDirection.getZ()>0.) &&
118 vectorPosition.getX()<targetHalfDiameter &&
119 vectorPosition.getX()>-targetHalfDiameter &&
120 vectorPosition.getY()<targetHalfDiameter &&
121 vectorPosition.getY()>-targetHalfDiameter)
122 {
123 analysisManager->FillH2(0,vectorPosition.getX(),vectorPosition.getY());
124 analysisManager->FillH1(0,energy);
125 fRun->CountParticlesTarget();
126 fTempTrack = fTrack->GetTrackID();
127 }
128
129 if(fTempTrack1 == 0)
130 {
131 fTempTrack1 = fTrack->GetTrackID();
132 }
133
134 if(fTrack->GetTrackID()!=fTempTrack1 && (momentumDirection.getZ()>0.) &&
135 fTempVector.getX()<targetHalfDiameter &&
136 fTempVector.getX()>-targetHalfDiameter &&
137 fTempVector.getY()<targetHalfDiameter &&
138 fTempVector.getY()>-targetHalfDiameter )
139 {
140
141 analysisManager->FillH2(4,fTempVector.getX(),fTempVector.getY());
142 analysisManager->FillH1(2,fTempEnergy);
143 fTempTrack1 = fTrack->GetTrackID();
144 }
145
146 fTempVector = aStep->GetPostStepPoint()->GetPosition(); //vectorPosition;
147 fTempEnergy = aStep->GetPostStepPoint()->GetKineticEnergy(); //energy;
148
149 analysisManager->FillH2(3,vectorPosition.getZ(),energy);
150
151 }
152
153 //----------------------------------------------
154 // Store ID for particles that are
155 // not protons/electrons or deuterons
156 //----------------------------------------------
157
158 if((name != "proton") && (name != "e-") && (name != "deuteron"))
159 {
160 fRun->StoreIsotopeID(fTrack->GetTrackID(),name);
161 }
162
163 //----------------------------------------------
164 // Collect of information for unstable isotopes
165 // generated from an interaction with the target
166 //----------------------------------------------
167
168 if (name!="deuteron")
169 {
170 if (( partType == "nucleus") && !(thePartDef->GetPDGStable()) && (fTrack->GetCurrentStepNumber()==1) && timeLife!=0.)
171 {
172 //G4cout << "Saving unstable particles ..." << G4endl;
173 G4int Z=thePartDef->GetAtomicNumber();
174 G4int A=thePartDef->GetAtomicMass();
175 analysisManager->FillH2(2,Z,A);
176
177 //----------------------------------------------
178 // isotopes count
179 //----------------------------------------------
180
181 fRun->PrimaryIsotopeCountTarget(name,timeLife);
182 analysisManager->FillH1(4,fTrack->GetPosition().getZ());
183 //particle that created the nucleus
184 std::map<G4int,G4String> parentID = fRun->GetIsotopeID();
185 G4String nameParent = parentID[fTrack->GetParentID()];
186 fRun->ParticleParent(name, process->GetProcessName());
187
188 //G4cout << name << " : " << process->GetProcessName() << " with track ID " << fTrack->GetTrackID() << " and step ID " << fTrack->GetCurrentStepNumber() << G4endl;
189
190 }
191 }
192
193 //----------------------------------------------
194 // Collect of information for stable isotopes
195 // generated from an interaction with the target
196 //----------------------------------------------
197
198 if (name!="deuteron")
199 {
200 if (( partType == "nucleus") && (thePartDef->GetPDGStable()) && (process->GetProcessName() != "RadioactiveDecay") && (fTrack->GetCurrentStepNumber()==1) )
201 {
202 //----------------------------------------------
203 // isotopes count
204 //----------------------------------------------
205 fRun->CountStableIsotopes(name);
206 }
207 }
208
209
210 //----------------------------------------------
211 // Collect unstable isotopes from decay
212 //----------------------------------------------
213
214 if (( partType == "nucleus") && !(thePartDef->GetPDGStable()) && (process->GetProcessName() == "RadioactiveDecay") && (fTrack->GetCurrentStepNumber()==1) && timeLife!=0)
215 {
216 std::map<G4int,G4String>::iterator itbis;
217 std::map<G4int,G4String> parentID = fRun->GetIsotopeID();
218 G4String nameParent = parentID[fTrack->GetParentID()];
219 fRun->DecayIsotopeCountTarget(name,nameParent,timeLife);
220 }
221
222 //----------------------------------------------
223 // Collect any other particles emitted
224 //----------------------------------------------
225
226 if((partType!="nucleus")&&(name!="proton")&&(name!="deuteron"))
227 {
228 fRun->ParticleCountTarget(name);
229 //Condition so the particle will be counted for only one step
230 if((fTrack->GetCurrentStepNumber()==1))
231 {
232 if(process->GetProcessName() != "RadioactiveDecay")
233 {
234 if(name=="e+"){
235 analysisManager->FillH1(5,energy);
236 }
237 if(name=="e-"){
238 analysisManager->FillH1(6,energy);
239 }
240 if(name=="gamma"){
241 analysisManager->FillH1(7,energy);
242 }
243 if(name=="neutron"){
244 analysisManager->FillH1(8,energy);
245 }
246 }
247
248 if(process->GetProcessName() == "RadioactiveDecay")
249 {
250 if(name=="e+"){
251 analysisManager->FillH1(9,energy);
252 }
253 if(name=="e-"){
254 analysisManager->FillH1(10,energy);
255 }
256 if(name=="gamma"){
257 analysisManager->FillH1(11,energy);
258 }
259 if(name=="neutron"){
260 analysisManager->FillH1(12,energy);
261 }
262 if(name=="nu_e"){
263 analysisManager->FillH1(13,energy);
264 }
265 if(name=="anti_nu_e"){
266 analysisManager->FillH1(14,energy);
267 }
268 }
269 }
270 }
271
272
273 fRun->SetTargetVolume(fDet->GetTargetVolume());
274 fRun->SetTargetThickness(fDet->GetTargetThickness());
275 fRun->SetTargetDiameter(fDet->GetTargetDiameter());
276
277 return true;
278
279 }
280
281
282