ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
SteppingAction.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file SteppingAction.cc
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 // This example is provided by the Geant4-DNA collaboration
27 // Any report or published results obtained using the Geant4-DNA software
28 // shall cite the following Geant4-DNA collaboration publication:
29 // Med. Phys. 37 (2010) 4692-4708
30 // The Geant4-DNA web site is available at http://geant4-dna.org
31 //
32 //
35 
36 #include "Analysis.hh"
37 
38 #include "SteppingAction.hh"
39 #include "RunAction.hh"
40 #include "DetectorConstruction.hh"
41 #include "PrimaryGeneratorAction.hh"
42 
43 #include "G4SystemOfUnits.hh"
44 #include "G4SteppingManager.hh"
45 
46 #include "G4Electron.hh"
47 #include "G4Proton.hh"
48 #include "G4Gamma.hh"
49 #include "G4Alpha.hh"
50 #include "G4DNAGenericIonsManager.hh"
51 #include "G4EventManager.hh"
52 #include "G4Event.hh"
53 
54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
55 
56 SteppingAction::SteppingAction(): G4UserSteppingAction()
57 {}
58 
59 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
60 
61 SteppingAction::~SteppingAction()
62 {}
63 
64 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
65 
66 void SteppingAction::UserSteppingAction(const G4Step* step)
67 {
68  // Protection
69 
70  if (!step->GetPostStepPoint()) return;
71  if (!step->GetPostStepPoint()->GetProcessDefinedStep()) return;
72 
73  //
74 
75  G4double flagParticle=-1.;
76  G4double flagProcess=-1.;
77  G4double x,y,z,xp,yp,zp;
78 
79  // Particle identification
80 
81  // The following method avoids the usage of string comparison
82 
83  G4ParticleDefinition * partDef =
84  step->GetTrack()->GetDynamicParticle()->GetDefinition();
85 
86  if (partDef == G4Gamma::GammaDefinition())
87  flagParticle = 0;
88 
89  if (partDef == G4Electron::ElectronDefinition())
90  flagParticle = 1;
91 
92  if (partDef == G4Proton::ProtonDefinition())
93  flagParticle = 2;
94 
95  if (partDef == G4Alpha::AlphaDefinition())
96  flagParticle = 4;
97 
98  G4DNAGenericIonsManager *instance;
99  instance = G4DNAGenericIonsManager::Instance();
100 
101  // Usage example
102  /*
103  G4ParticleDefinition* protonDef = G4Proton::ProtonDefinition();
104  G4ParticleDefinition* hydrogenDef = instance->GetIon("hydrogen");
105  G4ParticleDefinition* alphaPlusPlusDef = instance->GetIon("alpha++");
106  G4ParticleDefinition* alphaPlusDef = instance->GetIon("alpha+");
107  G4ParticleDefinition* heliumDef = instance->GetIon("helium");
108  */
109 
110  if (partDef == instance->GetIon("hydrogen"))
111  flagParticle = 3;
112 
113  if (partDef == instance->GetIon("alpha+"))
114  flagParticle = 5;
115 
116  if (partDef == instance->GetIon("helium"))
117  flagParticle = 6;
118 
119  // Alternative method (based on string comparison)
120 
121  /*
122  const G4String& particleName = step->GetTrack()->GetDynamicParticle()->
123  GetDefinition()->GetParticleName();
124 
125  if (particleName == "gamma") flagParticle = 0;
126  else if (particleName == "e-") flagParticle = 1;
127  else if (particleName == "proton") flagParticle = 2;
128  else if (particleName == "hydrogen") flagParticle = 3;
129  else if (particleName == "alpha") flagParticle = 4;
130  else if (particleName == "alpha+") flagParticle = 5;
131  else if (particleName == "helium") flagParticle = 6;
132  */
133 
134  // Process identification
135 
136  // Process sub-types are listed in G4PhysicsListHelper.cc
137  // or in Geant4-DNA process class implementation files (*.cc)
138 
139  G4StepPoint * preStep = step->GetPreStepPoint();
140  G4StepPoint * postStep = step->GetPostStepPoint();
141  G4int procID = postStep->GetProcessDefinedStep()->GetProcessSubType();
142 
143  const G4String& processName = postStep->
144  GetProcessDefinedStep()->GetProcessName();
145 
146  if (processName=="Capture") flagProcess =1;
147  // (no subType and procID exists at the moment for this process)
148 
149  else if (flagParticle == 1)
150 
151  {
152  if (procID==58) flagProcess =10;
153  else if (procID==51) flagProcess =11;
154  else if (procID==52) flagProcess =12;
155  else if (procID==53) flagProcess =13;
156  else if (procID==55) flagProcess =14;
157  else if (procID==54) flagProcess =15;
158  else if (procID==10) flagProcess =110;
159  else if (procID==1) flagProcess =120;
160  else if (procID==2) flagProcess =130;
161  }
162 
163  else if (flagParticle == 2)
164 
165  {
166  if (procID==51) flagProcess =21;
167  else if (procID==52) flagProcess =22;
168  else if (procID==53) flagProcess =23;
169  else if (procID==56) flagProcess =24;
170  else if (procID==10) flagProcess =210;
171  else if (procID==1) flagProcess =220;
172  else if (procID==2) flagProcess =230;
173  else if (procID==8) flagProcess =240;
174  }
175 
176  else if (flagParticle == 3)
177 
178  {
179  if (procID==51) flagProcess =31;
180  else if (procID==52) flagProcess =32;
181  else if (procID==53) flagProcess =33;
182  else if (procID==57) flagProcess =35;
183  }
184 
185  else if (flagParticle == 4)
186 
187  {
188  if (procID==51) flagProcess =41;
189  else if (procID==52) flagProcess =42;
190  else if (procID==53) flagProcess =43;
191  else if (procID==56) flagProcess =44;
192  else if (procID==10) flagProcess =410;
193  else if (procID==1) flagProcess =420;
194  else if (procID==2) flagProcess =430;
195  else if (procID==8) flagProcess =440;
196  }
197 
198  else if (flagParticle == 5)
199 
200  {
201  if (procID==51) flagProcess =51;
202  else if (procID==52) flagProcess =52;
203  else if (procID==53) flagProcess =53;
204  else if (procID==56) flagProcess =54;
205  else if (procID==57) flagProcess =55;
206  else if (procID==10) flagProcess =510;
207  else if (procID==1) flagProcess =520;
208  else if (procID==2) flagProcess =530;
209  else if (procID==8) flagProcess =540;
210  }
211 
212  else if (flagParticle == 6)
213 
214  {
215  if (procID==51) flagProcess =61;
216  else if (procID==52) flagProcess =62;
217  else if (procID==53) flagProcess =63;
218  else if (procID==57) flagProcess =65;
219 
220  }
221 
222  else if (processName=="GenericIon_G4DNAIonisation") flagProcess =73;
223  else if (processName=="msc") flagProcess =710;
224  else if (processName=="CoulombScat") flagProcess =720;
225  else if (processName=="ionIoni") flagProcess =730;
226  else if (processName=="nuclearStopping") flagProcess =740;
227  // (for all GenericIons)
228 
229  // Alternatively, using process names
230 
231  /*
232  else if (processName=="e-_G4DNAElectronSolvation") flagProcess =10;
233  else if (processName=="e-_G4DNAElastic") flagProcess =11;
234  else if (processName=="e-_G4DNAExcitation") flagProcess =12;
235  else if (processName=="e-_G4DNAIonisation") flagProcess =13;
236  else if (processName=="e-_G4DNAAttachment") flagProcess =14;
237  else if (processName=="e-_G4DNAVibExcitation") flagProcess =15;
238 
239  else if (processName=="proton_G4DNAElastic") flagProcess =21;
240  else if (processName=="proton_G4DNAExcitation") flagProcess =22;
241  else if (processName=="proton_G4DNAIonisation") flagProcess =23;
242  else if (processName=="proton_G4DNAChargeDecrease") flagProcess =24;
243 
244  else if (processName=="hydrogen_G4DNAElastic") flagProcess =31;
245  else if (processName=="hydrogen_G4DNAExcitation") flagProcess =32;
246  else if (processName=="hydrogen_G4DNAIonisation") flagProcess =33;
247  else if (processName=="hydrogen_G4DNAChargeIncrease") flagProcess =35;
248 
249  else if (processName=="alpha_G4DNAElastic") flagProcess =41;
250  else if (processName=="alpha_G4DNAExcitation") flagProcess =42;
251  else if (processName=="alpha_G4DNAIonisation") flagProcess =43;
252  else if (processName=="alpha_G4DNAChargeDecrease") flagProcess =44;
253 
254  else if (processName=="alpha+_G4DNAElastic") flagProcess =51;
255  else if (processName=="alpha+_G4DNAExcitation") flagProcess =52;
256  else if (processName=="alpha+_G4DNAIonisation") flagProcess =53;
257  else if (processName=="alpha+_G4DNAChargeDecrease") flagProcess =54;
258  else if (processName=="alpha+_G4DNAChargeIncrease") flagProcess =55;
259 
260  else if (processName=="helium_G4DNAElastic") flagProcess =61;
261  else if (processName=="helium_G4DNAExcitation") flagProcess =62;
262  else if (processName=="helium_G4DNAIonisation") flagProcess =63;
263  else if (processName=="helium_G4DNAChargeIncrease") flagProcess =65;
264 
265  else if (processName=="GenericIon_G4DNAIonisation") flagProcess =73;
266 
267  */
268 
269  if (processName!="Transportation")
270  {
271  x=preStep->GetPosition().x()/nanometer;
272  y=preStep->GetPosition().y()/nanometer;
273  z=preStep->GetPosition().z()/nanometer;
274 
275  xp=postStep->GetPosition().x()/nanometer;
276  yp=postStep->GetPosition().y()/nanometer;
277  zp=postStep->GetPosition().z()/nanometer;
278 
279  // get analysis manager
280 
281  G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
282 
283  // fill ntuple
284  analysisManager->FillNtupleDColumn(0, flagParticle);
285  analysisManager->FillNtupleDColumn(1, flagProcess);
286  analysisManager->FillNtupleDColumn(2, x);
287  analysisManager->FillNtupleDColumn(3, y);
288  analysisManager->FillNtupleDColumn(4, z);
289  analysisManager->FillNtupleDColumn(5,
290  step->GetTotalEnergyDeposit()/eV);
291 
292  analysisManager->FillNtupleDColumn(6,
293  std::sqrt((x-xp)*(x-xp)+
294  (y-yp)*(y-yp)+(z-zp)*(z-zp)));
295 
296  analysisManager->FillNtupleDColumn(7,
297  (preStep->
298  GetKineticEnergy() -
299  postStep->
300  GetKineticEnergy())/eV );
301 
302  analysisManager->FillNtupleDColumn(8, preStep->
303  GetKineticEnergy()/eV);
304 
305  analysisManager->FillNtupleDColumn(9,
306  preStep->GetMomentumDirection()
307  *postStep->GetMomentumDirection() );
308 
309  analysisManager->FillNtupleIColumn(10,
310  G4EventManager::GetEventManager()->
311  GetConstCurrentEvent()->GetEventID());
312 
313  analysisManager->FillNtupleIColumn(11, step->GetTrack()->GetTrackID());
314 
315  analysisManager->FillNtupleIColumn(12, step->GetTrack()->GetParentID());
316 
317  analysisManager->FillNtupleIColumn(13,
318  step->GetTrack()->GetCurrentStepNumber());
319 
320  analysisManager->AddNtupleRow();
321  }
322 }