Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
EventAction.cc
Go to the documentation of this file.
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 // Delage et al. PDB4DNA: implementation of DNA geometry from the Protein Data
31 // Bank (PDB) description for Geant4-DNA Monte-Carlo
32 // simulations (submitted to Comput. Phys. Commun.)
33 // The Geant4-DNA web site is available at http://geant4-dna.org
34 //
35 // $Id$
36 //
39 
40 #include "EventAction.hh"
41 
42 #include "Analysis.hh"
43 #include "EventActionMessenger.hh"
44 #include "G4Event.hh"
45 #include "G4UnitsTable.hh"
46 #include "G4SystemOfUnits.hh"
47 #include "Randomize.hh"
48 
49 #include <algorithm>
50 
51 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
52 
54 {
55  //default parameter values
56  //
57  fThresEdepForSSB=8.22*eV;
58  fThresDistForDSB=10;
59  fTotalEnergyDeposit=0;
60 
61  //create commands
62  //
63  fpEventMessenger = new EventActionMessenger(this);
64 }
65 
66 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
67 
69 {
70  delete fpEventMessenger;
71 }
72 
73 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
74 
76 {
77  // Initialization of parameters
78  //
79  fTotalEnergyDeposit=0.;
80  fEdepStrand1.clear();
81  fEdepStrand2.clear();
82 }
83 
84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
85 
87 {
88  // At the end of an event, compute the number of strand breaks
89  //
90  G4int sb[2] = {0,0};
91  ComputeStrandBreaks(sb);
92  // Fill histograms
93  //
94  G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
95 
96  if ( fTotalEnergyDeposit>0. )
97  {
98  analysisManager->FillH1(1,fTotalEnergyDeposit);
99  }
100  if ( sb[0]>0 )
101  {
102  analysisManager->FillH1(2,sb[0]);
103  }
104  if ( sb[1]>0 )
105  {
106  analysisManager->FillH1(3,sb[1]);
107  }
108 }
109 
110 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
111 
112 void EventAction::ComputeStrandBreaks(G4int* sb)
113 {
114  // sb quantities
115  //
116  G4int ssb1=0;
117  G4int ssb2=0;
118  G4int dsb=0;
119 
120  // nucleotide id and energy deposit for each strand
121  G4int nucl1;
122  G4int nucl2;
123  G4double edep1;
124  G4double edep2;
125 
126  //Read strand1
127  //
128  while ( !fEdepStrand1.empty() )
129  {
130  nucl1 = fEdepStrand1.begin()->first;
131  edep1 = fEdepStrand1.begin()->second;
132  fEdepStrand1.erase( fEdepStrand1.begin() );
133 
134  // SSB in strand1
135  //
136  if ( edep1 >= fThresEdepForSSB/eV )
137  {
138  ssb1++;
139  }
140 
141  // Look at strand2
142  //
143  if ( !fEdepStrand2.empty() )
144  {
145  do
146  {
147  nucl2 = fEdepStrand2.begin()->first;
148  edep2 = fEdepStrand2.begin()->second;
149  if ( edep2 >= fThresEdepForSSB/eV )
150  {
151  ssb2++;
152  }
153  fEdepStrand2.erase( fEdepStrand2.begin() );
154  } while ( ((nucl1-nucl2)>fThresDistForDSB) && (!fEdepStrand2.empty()) );
155 
156  // no dsb
157  //
158  if ( nucl2-nucl1 > fThresDistForDSB )
159  {
160  fEdepStrand2[nucl2]=edep2;
161  if ( edep2 >= fThresEdepForSSB/eV )
162  {
163  ssb2--;
164  }
165  }
166 
167  // one dsb
168  //
169  if ( std::abs(nucl2-nucl1) <= fThresDistForDSB )
170  {
171  if ( ( edep2 >= fThresEdepForSSB/eV ) &&
172  ( edep1 >= fThresEdepForSSB/eV ) )
173  {
174  ssb1--;
175  ssb2--;
176  dsb++;
177  }
178  }
179  }
180  }
181 
182  // End with not processed data
183  //
184  while ( !fEdepStrand1.empty() )
185  {
186  nucl1 = fEdepStrand1.begin()->first;
187  edep1 = fEdepStrand1.begin()->second;
188  if ( edep1 >= fThresEdepForSSB/eV )
189  {
190  ssb1++;
191  }
192  fEdepStrand1.erase( fEdepStrand1.begin() );
193  }
194 
195  while ( !fEdepStrand2.empty() )
196  {
197  nucl2 = fEdepStrand2.begin()->first;
198  edep2 = fEdepStrand2.begin()->second;
199  if ( edep2 >= fThresEdepForSSB/eV )
200  {
201  ssb2++;
202  }
203  fEdepStrand2.erase( fEdepStrand2.begin() );
204  }
205 
206  sb[0]=ssb1+ssb2;
207  sb[1]=dsb;
208 }
void BeginOfEventAction(const G4Event *)
Definition: EventAction.cc:61
int G4int
Definition: G4Types.hh:78
void EndOfEventAction(const G4Event *)
Definition: EventAction.cc:95
static constexpr double eV
Definition: G4SIunits.hh:215
double G4double
Definition: G4Types.hh:76
G4CsvAnalysisManager G4AnalysisManager
Definition: g4csv_defs.hh:77