Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4MoleculeCounter.hh
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 // Author: Mathieu Karamitros
27 
28 // The code is developed in the framework of the ESA AO7146
29 //
30 // We would be very happy hearing from you, send us your feedback! :)
31 //
32 // In order for Geant4-DNA to be maintained and still open-source,
33 // article citations are crucial.
34 // If you use Geant4-DNA chemistry and you publish papers about your software,
35 // in addition to the general paper on Geant4-DNA:
36 //
37 // Int. J. Model. Simul. Sci. Comput. 1 (2010) 157–178
38 //
39 // we would be very happy if you could please also cite the following
40 // reference papers on chemistry:
41 //
42 // J. Comput. Phys. 274 (2014) 841-882
43 // Prog. Nucl. Sci. Tec. 2 (2011) 503-508
44 
45 #ifndef G4MoleculeCounter_h
46 #define G4MoleculeCounter_h
47 
48 #include "G4VMoleculeCounter.hh"
49 #include <map>
50 #include <memory>
51 #include <set>
52 #include <vector>
53 
54 //------------------------------------------------------------------------------
55 
57  bool operator()(const double& a, const double& b) const{
58  if (std::fabs(a - b) < fPrecision){
59  return false;
60  }
61  else{
62  return a < b;
63  }
64  }
65 
66  static G4ThreadLocal double fPrecision;
67 };
68 
69 //------------------------------------------------------------------------------
70 
71 typedef std::map<G4double, G4int, compDoubleWithPrecision>
73 
74 typedef std::unique_ptr<std::set<G4double> > RecordedTimes;
75 typedef std::set<G4double>::iterator RecordedTimesIterator;
76 
77 //------------------------------------------------------------------------------
78 
80 {
81  //----------------------------------------------------------------------------
82 public:
83  typedef std::map<G4MolecularConfiguration*,
85  typedef std::unique_ptr<std::vector<G4MolecularConfiguration*> >
87 
88  //----------------------------------------------------------------------------
89 protected:
91  virtual ~G4MoleculeCounter();
92 
94  std::map<const G4MoleculeDefinition*, G4bool> fDontRegister;
95 
98 
99  struct Search
100  {
102  fLowerBoundSet = false;
103  }
104  CounterMapType::iterator fLastMoleculeSearched;
105  NbMoleculeAgainstTime::iterator fLowerBoundTime;
107  };
108 
109  std::unique_ptr<Search> fpLastSearch;
110 
111 #ifdef MOLECULE_COUNTER_TESTING
112 public:
113 #else
114 protected:
115 #endif
116 
117  friend class G4Molecule;
118  friend class G4VMoleculeCounter;
119 
120  //----------------------------------------------------------------------------
121 public:
122  static G4MoleculeCounter* Instance();
123 
124  //----------------------------------------------------------------------------
125 public:
126  void Initialize() override;
127 
128  void ResetCounter() override;
129 
131  G4double time,
132  const G4ThreeVector* position = nullptr,
133  int number = 1) override;
135  G4double time,
136  const G4ThreeVector* position = nullptr,
137  int number = 1) override;
138 
139  /* The dynamics of the given molecule won't be saved into memory.*/
140  inline void DontRegister(const G4MoleculeDefinition*) override;
141  inline bool IsRegistered(const G4MoleculeDefinition*) override;
142  inline void RegisterAll() override;
143 
144  //----------------------------------------------------------------------------
145 public:
147  double time);
148  inline const NbMoleculeAgainstTime&
150 
153 
154  inline void SetVerbose(G4int);
155  inline G4int GetVerbose();
156 
157  /*
158  * It sets the min time difference in between two time slices.
159  */
160  void SetTimeSlice(double);
161 
162 
163  void Dump();
164 
166  {
168  }
169 
170  inline void CheckTimeForConsistency(G4bool flag)
171  {
173  }
174 
175  //----------------------------------------------------------------------------
176 protected:
178  int SearchUpperBoundTime(double time, bool sameTypeOfMolecule);
179 };
180 
181 //------------------------------------------------------------------------------
182 
184 {
185  if(fVerbose)
186  {
187  G4cout << " ---> G4MoleculeCounter::ResetCounter" << G4endl;
188  }
189  fCounterMap.clear();
190  fpLastSearch.reset(0);
191 }
192 
193 //------------------------------------------------------------------------------
194 
195 inline const NbMoleculeAgainstTime&
197 {
198  return fCounterMap[molecule];
199 }
200 
201 //------------------------------------------------------------------------------
202 
204 {
205  fVerbose = level;
206 }
207 
208 //------------------------------------------------------------------------------
209 
211 {
212  return fVerbose;
213 }
214 
215 //------------------------------------------------------------------------------
216 
217 inline void
219 {
220  fDontRegister[molDef] = true;
221 }
222 
223 //------------------------------------------------------------------------------
224 
225 inline bool
227 {
228  if(fDontRegister.find(molDef) == fDontRegister.end()) return true;
229  return false;
230 }
231 
232 //------------------------------------------------------------------------------
233 
235 {
236  fDontRegister.clear();
237 }
238 
239 #endif
void AddAMoleculeAtTime(G4MolecularConfiguration *, G4double time, const G4ThreeVector *position=nullptr, int number=1) override
void CheckTimeForConsistency(G4bool flag)
int SearchUpperBoundTime(double time, bool sameTypeOfMolecule)
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
void ResetCounter() override
CounterMapType fCounterMap
std::set< G4double >::iterator RecordedTimesIterator
bool IsRegistered(const G4MoleculeDefinition *) override
#define G4ThreadLocal
Definition: tls.hh:89
std::unique_ptr< std::vector< G4MolecularConfiguration * > > RecordedMolecules
void RemoveAMoleculeAtTime(G4MolecularConfiguration *, G4double time, const G4ThreeVector *position=nullptr, int number=1) override
int G4int
Definition: G4Types.hh:78
G4bool SearchTimeMap(G4MolecularConfiguration *molecule)
std::map< const G4MoleculeDefinition *, G4bool > fDontRegister
int GetNMoleculesAtTime(G4MolecularConfiguration *molecule, double time)
std::unique_ptr< Search > fpLastSearch
tuple b
Definition: test.py:12
void RegisterAll() override
RecordedMolecules GetRecordedMolecules()
const NbMoleculeAgainstTime & GetNbMoleculeAgainstTime(G4MolecularConfiguration *molecule)
G4GLOB_DLL std::ostream G4cout
static G4MoleculeCounter * Instance()
std::map< G4MolecularConfiguration *, NbMoleculeAgainstTime > CounterMapType
bool G4bool
Definition: G4Types.hh:79
void Initialize() override
G4bool IsTimeCheckedForConsistency() const
bool operator()(const double &a, const double &b) const
CounterMapType::iterator fLastMoleculeSearched
std::unique_ptr< std::set< G4double > > RecordedTimes
static G4ThreadLocal double fPrecision
NbMoleculeAgainstTime::iterator fLowerBoundTime
void DontRegister(const G4MoleculeDefinition *) override
std::map< G4double, G4int, compDoubleWithPrecision > NbMoleculeAgainstTime
#define G4endl
Definition: G4ios.hh:61
G4bool fCheckTimeIsConsistentWithScheduler
RecordedTimes GetRecordedTimes()
double G4double
Definition: G4Types.hh:76
virtual ~G4MoleculeCounter()
void SetTimeSlice(double)