Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/examples/extended/electromagnetic/TestEm5/src/TrackingAction.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 /// \file electromagnetic/TestEm5/src/TrackingAction.cc
 27 /// \brief Implementation of the TrackingAction class
 28 //
 29 //
 30 //
 31 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 32 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 33 
 34 #include "TrackingAction.hh"
 35 
 36 #include "DetectorConstruction.hh"
 37 #include "EventAction.hh"
 38 #include "HistoManager.hh"
 39 #include "Run.hh"
 40 
 41 #include "G4PhysicalConstants.hh"
 42 #include "G4RunManager.hh"
 43 #include "G4SystemOfUnits.hh"
 44 #include "G4Track.hh"
 45 
 46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 47 
 48 TrackingAction::TrackingAction(DetectorConstruction* det, EventAction* event)
 49   : fDetector(det), fEventAction(event)
 50 {}
 51 
 52 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 53 
 54 void TrackingAction::PreUserTrackingAction(const G4Track* aTrack)
 55 {
 56   // few initialisations
 57   //
 58   if (aTrack->GetTrackID() == 1) {
 59     fXstartAbs = fDetector->GetxstartAbs();
 60     fXendAbs = fDetector->GetxendAbs();
 61     fPrimaryCharge = aTrack->GetDefinition()->GetPDGCharge();
 62     fDirX = aTrack->GetMomentumDirection().x();
 63   }
 64 }
 65 
 66 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 67 
 68 void TrackingAction::PostUserTrackingAction(const G4Track* aTrack)
 69 {
 70   G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
 71 
 72   Run* run = static_cast<Run*>(G4RunManager::GetRunManager()->GetNonConstCurrentRun());
 73 
 74   G4ThreeVector position = aTrack->GetPosition();
 75   G4ThreeVector vertex = aTrack->GetVertexPosition();
 76   G4double charge = aTrack->GetDefinition()->GetPDGCharge();
 77   G4double energy = aTrack->GetKineticEnergy();
 78 
 79   G4bool transmit =
 80     (((fDirX >= 0.0 && position.x() >= fXendAbs) || (fDirX < 0.0 && position.x() <= fXstartAbs))
 81      && energy > 0.0);
 82   G4bool reflect =
 83     (((fDirX >= 0.0 && position.x() <= fXstartAbs) || (fDirX < 0.0 && position.x() <= fXendAbs))
 84      && energy > 0.0);
 85   G4bool notabsor = (transmit || reflect);
 86 
 87   // transmitted + reflected particles counter
 88   //
 89   G4int flag = 0;
 90   if (charge == fPrimaryCharge) flag = 1;
 91   if (aTrack->GetTrackID() == 1) flag = 2;
 92   if (transmit) fEventAction->SetTransmitFlag(flag);
 93   if (reflect) fEventAction->SetReflectFlag(flag);
 94 
 95   //
 96   // histograms
 97   //
 98   G4bool charged = (charge != 0.);
 99   G4bool neutral = !charged;
100 
101   // energy spectrum at exit
102   // zero energy charged particle are not taken into account
103   //
104   G4int id = 0;
105   if (transmit && charged)
106     id = 10;
107   else if (transmit && neutral)
108     id = 20;
109   else if (reflect && charged)
110     id = 30;
111   else if (reflect && neutral)
112     id = 40;
113 
114   if (id > 0) {
115     analysisManager->FillH1(id, energy);
116   }
117 
118   // energy leakage
119   //
120   if (notabsor) {
121     G4int trackID = aTrack->GetTrackID();
122     G4int index = 0;
123     if (trackID > 1) index = 1;  // primary=0, secondaries=1
124     G4double eleak = aTrack->GetKineticEnergy();
125     if ((aTrack->GetDefinition() == G4Positron::Positron()) && (trackID > 1))
126       eleak += 2 * electron_mass_c2;
127     run->AddEnergyLeak(eleak, index);
128   }
129 
130   // space angle distribution at exit : dN/dOmega
131   //
132   G4ThreeVector direction = aTrack->GetMomentumDirection();
133   id = 0;
134   if (transmit && charged)
135     id = 12;
136   else if (transmit && neutral)
137     id = 22;
138   else if (reflect && charged)
139     id = 32;
140   else if (reflect && neutral)
141     id = 42;
142 
143   if (id > 0 && std::abs(direction.x()) < 1.0) {
144     G4double theta = std::acos(direction.x());
145     if (theta > 0.0) {
146       G4double dteta = analysisManager->GetH1Width(id);
147       G4double unit = analysisManager->GetH1Unit(id);
148       G4double weight = unit / (twopi * std::sin(theta) * dteta);
149       /*
150       G4cout << "theta, dteta, unit, weight: "
151              << theta << "  "
152              << dteta << "  "
153              << unit << "  "
154              << weight << G4endl;
155       */
156       analysisManager->FillH1(id, theta, weight);
157     }
158   }
159 
160   // energy fluence at exit : dE(MeV)/dOmega
161   //
162   id = 0;
163   if (transmit && charged)
164     id = 11;
165   else if (reflect && charged)
166     id = 31;
167   else if (transmit && neutral)
168     id = 21;
169   else if (reflect && neutral)
170     id = 41;
171 
172   if (id > 0 && std::abs(direction.x()) < 1.0) {
173     G4double theta = std::acos(direction.x());
174     if (theta > 0.0) {
175       G4double dteta = analysisManager->GetH1Width(id);
176       G4double unit = analysisManager->GetH1Unit(id);
177       G4double weight = unit / (twopi * std::sin(theta) * dteta);
178       weight *= (aTrack->GetKineticEnergy() / MeV);
179       analysisManager->FillH1(id, theta, weight);
180     }
181   }
182 
183   // projected angles distribution at exit
184   //
185   id = 0;
186   if (transmit && charged)
187     id = 13;
188   else if (transmit && neutral)
189     id = 23;
190   else if (reflect && charged)
191     id = 33;
192   else if (reflect && neutral)
193     id = 43;
194 
195   if (id > 0) {
196     if (direction.x() != 0.0) {
197       G4double tet = std::atan(direction.y() / std::fabs(direction.x()));
198       analysisManager->FillH1(id, tet);
199       if (transmit && (flag == 2)) run->AddMscProjTheta(tet);
200 
201       tet = std::atan(direction.z() / std::fabs(direction.x()));
202       analysisManager->FillH1(id, tet);
203       if (transmit && (flag == 2)) run->AddMscProjTheta(tet);
204     }
205   }
206 
207   // projected position and radius at exit
208   //
209   if (transmit && energy > 0.0) {
210     G4double y = position.y(), z = position.z();
211     G4double r = std::sqrt(y * y + z * z);
212     analysisManager->FillH1(14, y);
213     analysisManager->FillH1(14, z);
214     analysisManager->FillH1(15, r);
215   }
216 
217   // x-vertex of charged secondaries
218   //
219   if ((aTrack->GetParentID() == 1) && charged) {
220     G4double xVertex = (aTrack->GetVertexPosition()).x();
221     analysisManager->FillH1(6, xVertex);
222     if (notabsor) analysisManager->FillH1(7, xVertex);
223   }
224 }
225 
226 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
227