| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 //
27 //
28
29 //....oooOO0OOooo........oooOO0OOooo........oo
30 //....oooOO0OOooo........oooOO0OOooo........oo
31
32 #include <iostream>
33
34 #include "FCALSteppingAction.hh"
35 #include "G4SteppingManager.hh"
36
37 #include "globals.hh"
38 #include "G4SystemOfUnits.hh"
39 #include "G4Track.hh"
40 #include "G4DynamicParticle.hh"
41 #include "G4Material.hh"
42
43 #include "G4LogicalVolume.hh"
44 #include "G4VPhysicalVolume.hh"
45 #include "G4VTouchable.hh"
46 #include "G4TouchableHistory.hh"
47
48 #include "G4Event.hh"
49
50 #include "G4ThreeVector.hh"
51
52 #include "G4ios.hh"
53
54 //....oooOO0OOooo........oooOO0OOooo........oo
55
56 FCALSteppingAction::FCALSteppingAction():IDold
57 {;}
58
59 //....oooOO0OOooo........oooOO0OOooo........oo
60
61 FCALSteppingAction::~FCALSteppingAction()
62 {;}
63
64 //....oooOO0OOooo........oooOO0OOooo........oo
65
66 void FCALSteppingAction::UserSteppingAction(co
67 {
68 // Get Edep
69 G4double Edep = astep->GetTotalEnergyDeposi
70
71 // Get Track
72 G4Track* aTrack = astep->GetTrack();
73
74 // Get Touchable History
75 G4TouchableHistory* theTouchable = (G4Touch
76
77 // Energy deposit in FCAL1 and FCAL2
78 if(Edep != 0.)
79 {
80 G4VPhysicalVolume* physVol = theTouchabl
81
82 if(strcmp(physVol->GetName(),"FCALEmModu
83 strcmp(physVol->GetName(),"F1LArGapPhysical
84 {
85 EdepFCALEm = EdepFCALEm + Edep;
86 };
87
88 if( (strcmp(physVol->GetName(), "FCALHad
89 (strcmp(physVol->GetName(), "CuPlateAPhysi
90 (strcmp(physVol->GetName(), "CuPlateBPhysi
91 (strcmp(physVol->GetName(), "WAbsorberPhys
92 (strcmp(physVol->GetName(), "F2RodPhysical
93 (strcmp(physVol->GetName(), "F2LArGapPhysi
94 {
95 EdepFCALHad = EdepFCALHad + Edep;
96 };
97 };
98
99 // Get Tracks properties
100 G4int TrackID = aTrack->GetTrackID();
101 G4int ParentID = aTrack->GetParentID();
102 // Get Associated particle
103 const G4DynamicParticle * aDynamicParticle =
104 G4ParticleDefinition * aParticle = aTrack->G
105 G4String ParticleName = aParticle->GetPartic
106
107 IDnow = EventNo + 10000*TrackID+ 100000000*P
108
109 if(IDnow != IDold)
110 {
111 IDold = IDnow;
112
113 // Get the primary particle
114 if(TrackID==1 && ParentID==0 && (aTrack-
115 {
116 PrimaryVertex = aTrack->GetVertexPositi
117 PrimaryDirection = aTrack->GetVertexMoment
118
119 NSecondaries = 1;
120 Secondaries[NSecondaries][1] = aParticle->
121 Secondaries[NSecondaries][2] = PrimaryVert
122 Secondaries[NSecondaries][3] = PrimaryVert
123 Secondaries[NSecondaries][4] = PrimaryVert
124 Secondaries[NSecondaries][5] = (aDynamicPa
125 Secondaries[NSecondaries][6] = (aDynamicPa
126 Secondaries[NSecondaries][7] = (aDynamicPa
127 Secondaries[NSecondaries][8] = aDynamicPar
128 Secondaries[NSecondaries][9] = aDynamicPar
129 Secondaries[NSecondaries][10] = aDynamicPa
130
131 G4cout << " **** Primary : " << EventNo <
132 G4cout << " Vertex : " << PrimaryVertex <<
133 }
134
135
136 // Get secondaries in air close to the p
137 G4double DCACut = 2.*mm;
138 G4String Material = aTrack->GetMaterial(
139 G4ThreeVector TrackPos = aTrack->GetVert
140
141 if(TrackID != 1 && ParentID == 1 && (str
142 {
143 SecondaryVertex = aTrack->GetVertexPositio
144 SecondaryDirection = aTrack->GetVertexMome
145
146 // calculate DCA of secondries to primary
147 Distance = PrimaryVertex - SecondaryVertex
148 VectorProduct = PrimaryDirection.cross(Sec
149 if(VectorProduct == G4ThreeVector() &&
150 PrimaryDirection != G4ThreeVector() &&
151 {
152 G4ThreeVector Temp = Distance.cross(Pr
153 VectorProduct = Temp.cross(PrimaryDire
154 };
155
156 VectorProductMagnitude = VectorProduct.mag();
157 if(VectorProductMagnitude == 0.)
158 {
159 VectorProductNorm = G4ThreeVector();
160 } else {
161 VectorProductNorm = (1./VectorProduct.
162 };
163 DistOfClosestApproach = Distance * VectorP
164
165 if(std::abs(DistOfClosestApproach) < DCACu
166 {
167 NSecondaries++;
168 Secondaries[0][0] = NSecondaries;
169 Secondaries[NSecondaries][1] = aPartic
170 Secondaries[NSecondaries][2] = (aTrack
171 Secondaries[NSecondaries][3] = (aTrack
172 Secondaries[NSecondaries][4] = (aTrack
173 Secondaries[NSecondaries][5] =(aDynami
174 Secondaries[NSecondaries][6] = (aDynam
175 Secondaries[NSecondaries][7] = (aDynam
176 Secondaries[NSecondaries][8] = aDynami
177 Secondaries[NSecondaries][9] = aDynami
178 Secondaries[NSecondaries][10] =aDynami
179 };
180 };
181 };
182
183
184 // Get the World leaving particle
185 if(aTrack->GetNextVolume() == 0) {
186 if(IDnow != IDout) {
187 IDout = IDnow;
188
189 NTracks++;
190
191 OutOfWorldTracksData[0][0] = NTracks;
192
193 OutOfWorldTracksData[NTracks][1] = aPart
194
195 OutOfWorldTracksData[NTracks][2] = (aTra
196 OutOfWorldTracksData[NTracks][3] = (aTra
197 OutOfWorldTracksData[NTracks][4] = (aTra
198
199 OutOfWorldTracksData[NTracks][5] = (aDyn
200 OutOfWorldTracksData[NTracks][6] = (aDyn
201 OutOfWorldTracksData[NTracks][7] = (aDyn
202
203 OutOfWorldTracksData[NTracks][8] = aDyna
204
205 OutOfWorldTracksData[NTracks][9] = aDyna
206
207 OutOfWorldTracksData[NTracks][10] = aDyn
208 };
209 };
210
211
212 }
213
214 void FCALSteppingAction::initialize(G4int Nev)
215 EventNo = Nev;
216 NTracks = 0;
217 NSecondaries = 0;
218 EdepFCALEm = EdepFCALHad = 0.;
219
220 for(G4int i=0; i<6000; i++)
221 {
222 for(G4int j=0; j<11; j++)
223 {
224 OutOfWorldTracksData[i][j] = 0.;
225 Secondaries[i][j] = 0.;
226 }
227 };
228 }
229
230 G4double FCALSteppingAction::GetOutOfWorldTrac
231 return OutOfWorldTracksData[i][j];
232 }
233
234 G4double FCALSteppingAction::GetSecondaries(G4
235 return Secondaries[i][j];
236 }
237
238 G4double FCALSteppingAction::GetEdepFCAL(G4Str
239 if(strcmp(FCAL,"FCALEm") == 0) {
240 return EdepFCALEm;
241 } else {
242 if(strcmp(FCAL,"FCALHad") == 0) {
243 return EdepFCALHad;}
244 }
245 return 0.0;
246 }
247
248
249 //....oooOO0OOooo........oooOO0OOooo........oo
250
251
252
253