ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4INCLStore.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4INCLStore.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 // INCL++ intra-nuclear cascade model
27 // Alain Boudard, CEA-Saclay, France
28 // Joseph Cugnon, University of Liege, Belgium
29 // Jean-Christophe David, CEA-Saclay, France
30 // Pekka Kaitaniemi, CEA-Saclay, France, and Helsinki Institute of Physics, Finland
31 // Sylvie Leray, CEA-Saclay, France
32 // Davide Mancusi, CEA-Saclay, France
33 //
34 #define INCLXX_IN_GEANT4_MODE 1
35 
36 #include "globals.hh"
37 
38 #include "G4INCLStore.hh"
39 #include <fstream>
40 #include "G4INCLLogger.hh"
41 #include "G4INCLCluster.hh"
42 
43 namespace G4INCL {
44 
45  Store::Store(Config const * const config) :
46  loadedA(0),
47  loadedZ(0),
48  loadedStoppingTime(0.),
49  theConfig(config)
50  {
51  }
52 
54  theBook.reset();
55  clear();
56  }
57 
59  inside.push_back(p);
60  }
61 
62  void Store::add(ParticleList const &pL) {
63  inside.insert(inside.end(), pL.begin(), pL.end());
64  }
65 
67  // Add the avatar to the avatar map
68  avatarList.push_back(a);
69 
70  ParticleList pList = a->getParticles();
71  for(ParticleIter i=pList.begin(), e=pList.end(); i!=e; ++i) {
72  addIncomingParticle((*i));
73  // Connect each particle to the avatar
75  }
76  }
77 
79  for(IAvatarIter a=al.begin(), e=al.end(); a!=e; ++a)
81  }
82 
83  void Store::add(IAvatar *a) {
84  // Add the avatar to the avatar map
85  avatarList.push_back(a);
86 
87  ParticleList pList = a->getParticles();
88  for(ParticleIter i=pList.begin(), e=pList.end(); i!=e; ++i) {
89  // Connect each particle to the avatar
91  }
92 
93  }
94 
96  incoming.push_back(p);
97  }
98 
100  particleAvatarConnections.insert(PAPair(p,a));
101  }
102 
104  PAIterPair iterPair = particleAvatarConnections.equal_range(p);
105  for(PAIter i=iterPair.first, last=iterPair.second; i!=last; ++i) {
106  if(i->second==a) {
107  particleAvatarConnections.erase(i);
108  return;
109  }
110  }
111  INCL_WARN("Loop in Store::disconnectAvatarFromParticle fell through." << std::endl
112  << "This indicates an inconsistent state of the particleAvatarConnections map." << std::endl);
113  }
114 
115  void Store::removeAvatar(IAvatar * const avatar) {
116  // Disconnect the avatar from particles
117  ParticleList particlesRelatedToAvatar = avatar->getParticles();
118  for(ParticleIter particleIter = particlesRelatedToAvatar.begin(), e = particlesRelatedToAvatar.end();
119  particleIter != e; ++particleIter) {
120  disconnectAvatarFromParticle(avatar, *particleIter);
121  }
122 
123  // Remove the avatar itself
124  avatarList.remove(avatar);
125  }
126 
128  PAIterPair iterPair = particleAvatarConnections.equal_range(particle);
129  for(PAIter i=iterPair.first, last=iterPair.second; i!=last; ++i) {
130  avatarsToBeRemoved.insert(i->second);
131  }
132  }
133 
135  for(ASIter a=avatarsToBeRemoved.begin(), e=avatarsToBeRemoved.end(); a!=e; ++a) {
136  removeAvatar(*a);
137  delete *a;
138  }
139  avatarsToBeRemoved.clear();
140  }
141 
143  if(avatarList.empty()) return NULL;
144 
145 #ifdef INCL_AVATAR_SEARCH_FullSort
146 
147  /* Full sort algorithm.
148  *
149  * Simple, but guaranteed to work.
150  */
152  IAvatar *avatar = avatarList.front();
153 
154 #elif defined(INCL_AVATAR_SEARCH_MinElement)
155 
156  /* Algorithm provided by the C++ stdlib. */
157  IAvatar *avatar = *(std::min_element(avatarList.begin(), avatarList.end(),
159 
160 #else
161 #error Unrecognized INCL_AVATAR_SEARCH. Allowed values are: FullSort, MinElement.
162 #endif
163 
164  removeAvatar(avatar);
165  return avatar;
166  }
167 
169  for(ParticleIter particleIter = inside.begin(), particleEnd=inside.end();
170  particleIter != particleEnd; ++particleIter) {
171  (*particleIter)->propagate(step);
172  }
173  }
174 
177  // The particle will be destroyed when destroying the Store
178  inside.remove(p);
179  }
180 
183  // Have to destroy the particle here, the Store will forget about it
184  inside.remove(p);
185  delete p;
186  }
187 
189  removeFromIncoming(particle);
190  add(particle);
191  }
192 
194  for(IAvatarIter iter = avatarList.begin(), e = avatarList.end(); iter != e; ++iter) {
195  delete *iter;
196  }
197 
199  avatarList.clear();
200  avatarsToBeRemoved.clear();
201  }
202 
203  void Store::clear() {
204  clearAvatars();
205 
206  clearInside();
207  clearOutgoing();
208 
209  if( incoming.size() != 0 ) {
210  INCL_WARN("Incoming list is not empty when Store::clear() is called" << '\n');
211  }
212  incoming.clear();
213 
214  }
215 
217  for(ParticleIter iter=inside.begin(), e=inside.end(); iter!=e; ++iter) {
218  delete *iter;
219  }
220  inside.clear();
221  }
222 
224  for(ParticleIter iter=outgoing.begin(), e=outgoing.end(); iter!=e; ++iter) {
225  if((*iter)->isCluster()) {
226  Cluster *c = dynamic_cast<Cluster *>(*iter);
227 // assert(c);
228 #ifdef INCLXX_IN_GEANT4_MODE
229  if(!c)
230  continue;
231 #endif
232  c->deleteParticles();
233  }
234  delete (*iter);
235  }
236  outgoing.clear();
237  }
238 
239  void Store::loadParticles(std::string const &filename) {
240  clear();
241  G4int projectileA, projectileZ, A, Z;
242  G4double stoppingTime, cutNN;
243  G4int ID, type, isParticipant;
244  G4double x, y, z;
245  G4double px, py, pz, E, v;
246 
247  std::ifstream in(filename.c_str());
248  in >> projectileA >> projectileZ >> A >> Z >> stoppingTime >> cutNN;
249  loadedA = A;
250  loadedZ = Z;
251  loadedStoppingTime = stoppingTime;
252 
253  G4int readA = 0;
254  G4int readZ = 0;
255  while(1) { /* Loop checking, 10.07.2015, D.Mancusi */
256  in >> ID >> type >> isParticipant >> x >> y >> z >> px >> py >> pz >> E >> v;
257  if(!in.good()) break;
258  ParticleType t;
259  if(type == 1) {
260  t = Proton;
261  readZ++;
262  readA++;
263  }
264  else if(type == -1) {
265  t = Neutron;
266  readA++;
267  }
268  else {
269  INCL_FATAL("Unrecognized particle type while loading particles; type=" << type << '\n');
270  t = UnknownParticle;
271  }
272 
273  Particle *p = new Particle(t, E, ThreeVector(px, py, pz),
274  ThreeVector(x, y, z));
275  p->setPotentialEnergy(v);
276  if(isParticipant == 1) {
277  p->makeParticipant();
279  }
280  add(p);
281  }
282 
283  in.close();
284  }
285 
287  std::stringstream ss;
288  G4int A = 0, Z = 0;
289  for(ParticleIter i=inside.begin(), e=inside.end(); i!=e; ++i) {
290  if((*i)->getType() == Proton) {
291  A++;
292  Z++;
293  }
294  if((*i)->getType() == Neutron) {
295  A++;
296  }
297  }
298  // Note: Projectile A and Z are set to 0 (we don't really know
299  // anything about them at this point).
300  ss << "0 0 " << A << " " << Z << " "
301  << "100.0" << " "
302  << "0.0" << '\n';
303 
304  for(ParticleIter i=inside.begin(), e=inside.end(); i!=e; ++i) {
305  G4int ID = (*i)->getID();
306  G4int type = 0;
307  if((*i)->getType() == Proton) {
308  type = 1;
309  }
310  if((*i)->getType() == Neutron) {
311  type = -1;
312  }
313 
314  G4int isParticipant = 0;
315  if((*i)->isParticipant()) {
316  isParticipant = 1;
317  }
318 
319  G4double x = (*i)->getPosition().getX();
320  G4double y = (*i)->getPosition().getY();
321  G4double z = (*i)->getPosition().getZ();
322  G4double E = (*i)->getEnergy();
323  G4double px = (*i)->getMomentum().getX();
324  G4double py = (*i)->getMomentum().getY();
325  G4double pz = (*i)->getMomentum().getZ();
326  G4double V = (*i)->getPotentialEnergy();
327 
328  ss << ID << " " << type << " " << isParticipant << " "
329  << x << " " << y << " " << z << " "
330  << px << " " << py << " " << pz << " "
331  << E << " " << V << '\n';
332  }
333 
334  return ss.str();
335  }
336 
337  void Store::writeParticles(std::string const &filename) {
338  std::ofstream out(filename.c_str());
340  out.close();
341  }
342 
343  std::string Store::printAvatars() {
344  std::stringstream ss;
345  for(IAvatarIter i = avatarList.begin(), e = avatarList.end(); i != e; ++i) {
346  ss << (*i)->toString() << '\n';
347  }
348  return ss.str();
349  }
350 
352  for(IAvatarIter i = avatarList.begin(), e = avatarList.end(); i != e; ++i)
353  if((*i)->getType()==CollisionAvatarType) return true;
354  return false;
355  }
356 
357 }