| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 // 27 //
28 28
29 // ------------------------------------------- 29 // ------------------------------------------------------------
30 // GEANT 4 class implementation file 30 // GEANT 4 class implementation file
31 // 31 //
32 // ---------------- G4FTFParticipants---- 32 // ---------------- G4FTFParticipants----------------
33 // by Gunter Folger, June 1998. 33 // by Gunter Folger, June 1998.
34 // class finding colliding particles in 34 // class finding colliding particles in FTFPartonStringModel
35 // Changed in a part by V. Uzhinsky in oder t 35 // Changed in a part by V. Uzhinsky in oder to put in correcpondence
36 // with original FRITIOF mode. November 36 // with original FRITIOF mode. November - December 2006.
37 // Ajusted for (anti) nucleus - nucleus inter 37 // Ajusted for (anti) nucleus - nucleus interactions by V. Uzhinsky.
38 // (February 2011) 38 // (February 2011)
39 // ------------------------------------------- 39 // ------------------------------------------------------------
40 40
41 #include <utility> 41 #include <utility>
42 #include <vector> 42 #include <vector>
43 #include <algorithm> 43 #include <algorithm>
44 44
45 #include "G4FTFParticipants.hh" 45 #include "G4FTFParticipants.hh"
46 #include "G4ios.hh" 46 #include "G4ios.hh"
47 #include "Randomize.hh" 47 #include "Randomize.hh"
48 #include "G4SystemOfUnits.hh" 48 #include "G4SystemOfUnits.hh"
49 #include "G4FTFParameters.hh" 49 #include "G4FTFParameters.hh"
50 #include "G4DiffractiveSplitableHadron.hh" 50 #include "G4DiffractiveSplitableHadron.hh"
51 #include "G4VSplitableHadron.hh" 51 #include "G4VSplitableHadron.hh"
52 #include "G4PhysicalConstants.hh" 52 #include "G4PhysicalConstants.hh"
53 53
54 54
55 //============================================ 55 //============================================================================
56 56
57 //#define debugFTFparticipant 57 //#define debugFTFparticipant
58 58
59 59
60 //============================================ 60 //============================================================================
61 61
62 G4FTFParticipants::G4FTFParticipants() : Bimpa 62 G4FTFParticipants::G4FTFParticipants() : Bimpact( 0.0 ), BinInterval( false ),
63 Bmin2 63 Bmin2( -1.0 ), Bmax2( -1.0 ),
64 curre 64 currentInteraction( -1 )
65 {} 65 {}
66 66
67 67
68 //============================================ 68 //============================================================================
69 69
70 G4FTFParticipants::~G4FTFParticipants() {} 70 G4FTFParticipants::~G4FTFParticipants() {}
71 71
72 //============================================ 72 //============================================================================
73 73
74 void G4FTFParticipants::GetList( const G4React 74 void G4FTFParticipants::GetList( const G4ReactionProduct& thePrimary,
75 G4FTFParamete 75 G4FTFParameters* theParameters ) {
76 76
77 #ifdef debugFTFparticipant 77 #ifdef debugFTFparticipant
78 G4cout << "Participants::GetList" << G4endl 78 G4cout << "Participants::GetList" << G4endl
79 << "thePrimary " << thePrimary.GetMom 79 << "thePrimary " << thePrimary.GetMomentum() << G4endl << G4endl;
80 #endif 80 #endif
81 81
82 G4double betta_z = thePrimary.GetMomentum(). 82 G4double betta_z = thePrimary.GetMomentum().z() / thePrimary.GetTotalEnergy();
83 if ( betta_z < 1.0e-10 ) betta_z = 1.0e-10; 83 if ( betta_z < 1.0e-10 ) betta_z = 1.0e-10;
84 84
85 StartLoop(); // reset Loop over Interactions 85 StartLoop(); // reset Loop over Interactions
86 86
87 for ( unsigned int i = 0; i < theInteraction 87 for ( unsigned int i = 0; i < theInteractions.size(); i++ ) delete theInteractions[i];
88 theInteractions.clear(); 88 theInteractions.clear();
89 89
90 G4double deltaxy = 2.0 * fermi; // Extra nu 90 G4double deltaxy = 2.0 * fermi; // Extra nuclear radius
91 91
92 if ( theProjectileNucleus == nullptr ) { // 92 if ( theProjectileNucleus == nullptr ) { // Hadron-nucleus or anti-baryon-nucleus interactions
93 93
94 G4double impactX( 0.0 ), impactY( 0.0 ); 94 G4double impactX( 0.0 ), impactY( 0.0 );
95 G4double B( 0.0 ), B2( 0.0 ); 95 G4double B( 0.0 ), B2( 0.0 );
96 96
97 G4VSplitableHadron* primarySplitable = new 97 G4VSplitableHadron* primarySplitable = new G4DiffractiveSplitableHadron( thePrimary );
98 98
99 #ifdef debugFTFparticipant 99 #ifdef debugFTFparticipant
100 G4cout << "Hadron-nucleus or anti-baryon-n 100 G4cout << "Hadron-nucleus or anti-baryon-nucleus interactions" << G4endl;
101 #endif 101 #endif
102 102
103 G4double xyradius = theNucleus->GetOuterRa 103 G4double xyradius = theNucleus->GetOuterRadius() + deltaxy; // Range of impact parameter sampling
104 104
105 const G4int maxNumberOfLoops = 1000; 105 const G4int maxNumberOfLoops = 1000;
106 G4int loopCounter = 0; 106 G4int loopCounter = 0;
107 do { 107 do {
108 108
109 std::pair< G4double, G4double > theImpac 109 std::pair< G4double, G4double > theImpactParameter;
110 if ( SampleBinInterval() ) { 110 if ( SampleBinInterval() ) {
111 B2 = GetBmin2() + G4UniformRand() * ( 111 B2 = GetBmin2() + G4UniformRand() * ( GetBmax2() - GetBmin2() );
112 B = B2 > 0.0 ? std::sqrt( B2 ) : 0.0; 112 B = B2 > 0.0 ? std::sqrt( B2 ) : 0.0;
113 G4double Phi = twopi * G4UniformRand() 113 G4double Phi = twopi * G4UniformRand();
114 impactX = B * std::cos( Phi ); 114 impactX = B * std::cos( Phi );
115 impactY = B * std::sin( Phi ); 115 impactY = B * std::sin( Phi );
116 SetImpactParameter( B ); 116 SetImpactParameter( B );
117 } else { 117 } else {
118 theImpactParameter = theNucleus->Choos 118 theImpactParameter = theNucleus->ChooseImpactXandY( xyradius );
119 impactX = theImpactParameter.first; 119 impactX = theImpactParameter.first;
120 impactY = theImpactParameter.second; 120 impactY = theImpactParameter.second;
121 SetImpactParameter( std::sqrt( sqr(imp 121 SetImpactParameter( std::sqrt( sqr(impactX) + sqr(impactY) ) );
122 } 122 }
123 123
124 #ifdef debugFTFparticipant 124 #ifdef debugFTFparticipant
125 G4cout << "New interaction list," << " b 125 G4cout << "New interaction list," << " b[fm]= "
126 << std::sqrt( sqr(impactX ) + sqr 126 << std::sqrt( sqr(impactX ) + sqr( impactY ) )/fermi << G4endl;
127 #endif 127 #endif
128 128
129 G4ThreeVector thePosition( impactX, impa 129 G4ThreeVector thePosition( impactX, impactY, 0.0 );
130 primarySplitable->SetPosition( thePositi 130 primarySplitable->SetPosition( thePosition );
131 131
132 theNucleus->StartLoop(); 132 theNucleus->StartLoop();
133 G4Nucleon* nucleon; 133 G4Nucleon* nucleon;
134 134
135 #ifdef debugFTFparticipant 135 #ifdef debugFTFparticipant
136 G4int TrN( 0 ); 136 G4int TrN( 0 );
137 #endif 137 #endif
138 138
139 while ( ( nucleon = theNucleus->GetNextN 139 while ( ( nucleon = theNucleus->GetNextNucleon() ) ) { /* Loop checking, 10.08.2015, A.Ribon */
140 140
141 G4double impact2 = sqr( impactX - nucl 141 G4double impact2 = sqr( impactX - nucleon->GetPosition().x() ) +
142 sqr( impactY - nucl 142 sqr( impactY - nucleon->GetPosition().y() );
143 143
144 if ( theParameters->GetProbabilityOfIn 144 if ( theParameters->GetProbabilityOfInteraction( impact2/fermi/fermi ) >
145 G4UniformRand() ) { 145 G4UniformRand() ) {
146 primarySplitable->SetStatus( 1 ); / 146 primarySplitable->SetStatus( 1 ); // It takes part in the interaction
147 G4VSplitableHadron* targetSplitable 147 G4VSplitableHadron* targetSplitable = 0;
148 if ( ! nucleon->AreYouHit() ) { 148 if ( ! nucleon->AreYouHit() ) {
149 targetSplitable = new G4Diffractiv 149 targetSplitable = new G4DiffractiveSplitableHadron( *nucleon );
150 nucleon->Hit( targetSplitable ); 150 nucleon->Hit( targetSplitable );
151 targetSplitable->SetStatus( 1 ); / 151 targetSplitable->SetStatus( 1 ); // It takes part in the interaction
152 152
153 #ifdef debugFTFparticipant 153 #ifdef debugFTFparticipant
154 G4cout << "Participated nucleons # 154 G4cout << "Participated nucleons #, " << TrN << " " << "Splitable Pr* Tr* "
155 << primarySplitable << " " 155 << primarySplitable << " " << targetSplitable << G4endl;
156 #endif 156 #endif
157 157
158 } 158 }
159 G4InteractionContent* aInteraction = 159 G4InteractionContent* aInteraction = new G4InteractionContent( primarySplitable );
160 G4Nucleon* PrNucleon = 0; 160 G4Nucleon* PrNucleon = 0;
161 aInteraction->SetProjectileNucleon( 161 aInteraction->SetProjectileNucleon( PrNucleon );
162 aInteraction->SetTarget( targetSplit 162 aInteraction->SetTarget( targetSplitable );
163 aInteraction->SetTargetNucleon( nucl 163 aInteraction->SetTargetNucleon( nucleon );
164 aInteraction->SetStatus( 1 ); 164 aInteraction->SetStatus( 1 );
165 aInteraction->SetInteractionTime( ( 165 aInteraction->SetInteractionTime( ( primarySplitable->GetPosition().z() +
166 166 nucleon->GetPosition().z() ) / betta_z );
167 theInteractions.push_back( aInteract 167 theInteractions.push_back( aInteraction );
168 } 168 }
169 169
170 #ifdef debugFTFparticipant 170 #ifdef debugFTFparticipant
171 TrN++; 171 TrN++;
172 #endif 172 #endif
173 173
174 } 174 }
175 175
176 } while ( ( theInteractions.size() == 0 ) 176 } while ( ( theInteractions.size() == 0 ) &&
177 ++loopCounter < maxNumberOfLoops 177 ++loopCounter < maxNumberOfLoops ); /* Loop checking, 10.08.2015, A.Ribon */
178 if ( loopCounter >= maxNumberOfLoops ) { 178 if ( loopCounter >= maxNumberOfLoops ) {
179 #ifdef debugFTFparticipant 179 #ifdef debugFTFparticipant
180 G4cout << "BAD situation: forced exit fr 180 G4cout << "BAD situation: forced exit from the while loop!" << G4endl;
181 #endif 181 #endif
182 return; 182 return;
183 } 183 }
184 184
185 #ifdef debugFTFparticipant 185 #ifdef debugFTFparticipant
186 G4cout << "Number of Hit nucleons " << the 186 G4cout << "Number of Hit nucleons " << theInteractions.size() << "\t Bx[fm] " << impactX/fermi
187 << "\t By[fm] " << impactY/fermi << 187 << "\t By[fm] " << impactY/fermi << "\t B[fm] "
188 << std::sqrt( sqr( impactX ) + sqr( 188 << std::sqrt( sqr( impactX ) + sqr( impactY ) )/fermi << G4endl << G4endl;
189 #endif 189 #endif
190 190
191 //SortInteractionsIncT(); // Not need beca 191 //SortInteractionsIncT(); // Not need because nucleons are sorted in increasing z-coordinates.
192 ShiftInteractionTime(); // To put correct 192 ShiftInteractionTime(); // To put correct times and z-coordinates
193 return; 193 return;
194 194
195 } // end of if ( theProjectileNucleus == 0 195 } // end of if ( theProjectileNucleus == 0 )
196 196
197 // Projectile and target are nuclei 197 // Projectile and target are nuclei
198 198
199 #ifdef debugFTFparticipant 199 #ifdef debugFTFparticipant
200 G4cout << "Projectile and target are nuclei" 200 G4cout << "Projectile and target are nuclei" << G4endl;
201 #endif 201 #endif
202 202
203 //G4cout<<theProjectileNucleus->GetOuterRadi 203 //G4cout<<theProjectileNucleus->GetOuterRadius()/fermi<<" "<<theNucleus->GetOuterRadius()/fermi<<" "<<deltaxy/fermi<<G4endl;
204 204
205 // Range of impact parameter sampling 205 // Range of impact parameter sampling
206 G4double xyradius = theProjectileNucleus->Ge 206 G4double xyradius = theProjectileNucleus->GetOuterRadius() + theNucleus->GetOuterRadius() + deltaxy;
207 207
208 G4double impactX( 0.0 ), impactY( 0.0 ); 208 G4double impactX( 0.0 ), impactY( 0.0 );
209 G4double B( 0.0 ), B2( 0.0 ); 209 G4double B( 0.0 ), B2( 0.0 );
210 210
211 const G4int maxNumberOfLoops = 1000; 211 const G4int maxNumberOfLoops = 1000;
212 G4int loopCounter = 0; 212 G4int loopCounter = 0;
213 do { 213 do {
214 214
215 std::pair< G4double, G4double > theImpactP 215 std::pair< G4double, G4double > theImpactParameter;
216 if ( SampleBinInterval() ) { 216 if ( SampleBinInterval() ) {
217 B2 = GetBmin2() + G4UniformRand() * ( Ge 217 B2 = GetBmin2() + G4UniformRand() * ( GetBmax2() - GetBmin2() );
218 B = B2 > 0.0 ? std::sqrt( B2 ) : 0.0; 218 B = B2 > 0.0 ? std::sqrt( B2 ) : 0.0; // In G4 internal units (mm)
219 G4double Phi = twopi * G4UniformRand(); 219 G4double Phi = twopi * G4UniformRand();
220 impactX = B * std::cos( Phi ); 220 impactX = B * std::cos( Phi );
221 impactY = B * std::sin( Phi ); 221 impactY = B * std::sin( Phi );
222 SetImpactParameter( B ); 222 SetImpactParameter( B );
223 } else { 223 } else {
224 theImpactParameter = theNucleus->ChooseI 224 theImpactParameter = theNucleus->ChooseImpactXandY( xyradius );
225 impactX = theImpactParameter.first; 225 impactX = theImpactParameter.first;
226 impactY = theImpactParameter.second; 226 impactY = theImpactParameter.second;
227 SetImpactParameter( std::sqrt( sqr(impac 227 SetImpactParameter( std::sqrt( sqr(impactX) + sqr(impactY) ) );
228 } 228 }
229 229
230 #ifdef debugFTFparticipant 230 #ifdef debugFTFparticipant
231 G4cout << "New interaction list, " << "b[f 231 G4cout << "New interaction list, " << "b[fm] "
232 << std::sqrt( sqr( impactX ) + sqr( 232 << std::sqrt( sqr( impactX ) + sqr( impactY ) )/fermi << G4endl;
233 #endif 233 #endif
234 234
235 G4ThreeVector theBeamPosition( impactX, im 235 G4ThreeVector theBeamPosition( impactX, impactY, 0.0 );
236 236
237 theProjectileNucleus->StartLoop(); 237 theProjectileNucleus->StartLoop();
238 G4Nucleon* ProjectileNucleon; 238 G4Nucleon* ProjectileNucleon;
239 239
240 #ifdef debugFTFparticipant 240 #ifdef debugFTFparticipant
241 G4int PrNuclN( 0 ); 241 G4int PrNuclN( 0 );
242 #endif 242 #endif
243 243
244 while ( ( ProjectileNucleon = theProjectil 244 while ( ( ProjectileNucleon = theProjectileNucleus->GetNextNucleon() ) ) { /* Loop checking, 10.08.2015, A.Ribon */
245 245
246 G4VSplitableHadron* ProjectileSplitable 246 G4VSplitableHadron* ProjectileSplitable = 0;
247 theNucleus->StartLoop(); 247 theNucleus->StartLoop();
248 G4Nucleon* TargetNucleon = 0; 248 G4Nucleon* TargetNucleon = 0;
249 249
250 #ifdef debugFTFparticipant 250 #ifdef debugFTFparticipant
251 G4int TrNuclN( 0 ); 251 G4int TrNuclN( 0 );
252 #endif 252 #endif
253 253
254 while ( ( TargetNucleon = theNucleus->Ge 254 while ( ( TargetNucleon = theNucleus->GetNextNucleon() ) ) { /* Loop checking, 10.08.2015, A.Ribon */
255 255
256 G4double impact2 = sqr( impactX + Proj 256 G4double impact2 = sqr( impactX + ProjectileNucleon->GetPosition().x() -
257 TargetNucleon- 257 TargetNucleon->GetPosition().x() ) +
258 sqr( impactY + Proj 258 sqr( impactY + ProjectileNucleon->GetPosition().y() -
259 TargetNucleon- 259 TargetNucleon->GetPosition().y() );
260 G4VSplitableHadron* TargetSplitable = 260 G4VSplitableHadron* TargetSplitable = 0;
261 if ( theParameters->GetProbabilityOfIn 261 if ( theParameters->GetProbabilityOfInteraction( impact2/fermi/fermi ) >
262 G4UniformRand() ) { // An intera 262 G4UniformRand() ) { // An interaction has happend!
263 263
264 #ifdef debugFTFparticipant 264 #ifdef debugFTFparticipant
265 G4cout << G4endl << "An Interaction 265 G4cout << G4endl << "An Interaction has happend" << G4endl << "Proj N mom " << PrNuclN
266 << " " << ProjectileNucleon-> 266 << " " << ProjectileNucleon->Get4Momentum() << "-------------" << G4endl
267 << "Targ N mom " << TrNuclN < 267 << "Targ N mom " << TrNuclN << " " << TargetNucleon->Get4Momentum() << G4endl
268 << "PrN TrN Z coords [fm]" << 268 << "PrN TrN Z coords [fm]" << ProjectileNucleon->GetPosition().z()/fermi
269 << " " << TargetNucleon->GetP 269 << " " << TargetNucleon->GetPosition().z()/fermi
270 << " " << ProjectileNucleon-> 270 << " " << ProjectileNucleon->GetPosition().z()/fermi +
271 TargetNucleon->GetP 271 TargetNucleon->GetPosition().z()/fermi << G4endl;
272 #endif 272 #endif
273 273
274 if ( ! ProjectileNucleon->AreYouHit( 274 if ( ! ProjectileNucleon->AreYouHit() ) {
275 // Projectile nucleon was not invo 275 // Projectile nucleon was not involved until now.
276 ProjectileSplitable = new G4Diffra 276 ProjectileSplitable = new G4DiffractiveSplitableHadron( *ProjectileNucleon );
277 ProjectileNucleon->Hit( Projectile 277 ProjectileNucleon->Hit( ProjectileSplitable );
278 ProjectileSplitable->SetStatus( 1 278 ProjectileSplitable->SetStatus( 1 ); // It takes part in the interaction
279 } else { // Projectile nucleon was 279 } else { // Projectile nucleon was involved before.
280 ProjectileSplitable = ProjectileNu 280 ProjectileSplitable = ProjectileNucleon->GetSplitableHadron();
281 } 281 }
282 282
283 if ( ! TargetNucleon->AreYouHit() ) 283 if ( ! TargetNucleon->AreYouHit() ) { // Target nucleon was not involved until now
284 TargetSplitable = new G4Diffractiv 284 TargetSplitable = new G4DiffractiveSplitableHadron( *TargetNucleon );
285 TargetNucleon->Hit( TargetSplitabl 285 TargetNucleon->Hit( TargetSplitable );
286 TargetSplitable->SetStatus( 1 ); 286 TargetSplitable->SetStatus( 1 ); // It takes part in the interaction
287 } else { // Target nucleon was invo 287 } else { // Target nucleon was involved before.
288 TargetSplitable = TargetNucleon->G 288 TargetSplitable = TargetNucleon->GetSplitableHadron();
289 } 289 }
290 290
291 G4InteractionContent* anInteraction 291 G4InteractionContent* anInteraction = new G4InteractionContent( ProjectileSplitable );
292 anInteraction->SetTarget( TargetSpli 292 anInteraction->SetTarget( TargetSplitable );
293 anInteraction->SetProjectileNucleon( 293 anInteraction->SetProjectileNucleon( ProjectileNucleon );
294 anInteraction->SetTargetNucleon( Tar 294 anInteraction->SetTargetNucleon( TargetNucleon );
295 anInteraction->SetInteractionTime( ( 295 anInteraction->SetInteractionTime( ( ProjectileNucleon->GetPosition().z() +
296 296 TargetNucleon->GetPosition().z() ) / betta_z );
297 anInteraction->SetStatus( 1 ); 297 anInteraction->SetStatus( 1 );
298 298
299 #ifdef debugFTFparticipant 299 #ifdef debugFTFparticipant
300 G4cout << "Part anInteraction->GetIn 300 G4cout << "Part anInteraction->GetInteractionTime() [fm] "
301 << anInteraction->GetInteract 301 << anInteraction->GetInteractionTime()/fermi << G4endl
302 << "Splitable Pr* Tr* " << Pr 302 << "Splitable Pr* Tr* " << ProjectileSplitable << " "
303 << TargetSplitable << G4endl; 303 << TargetSplitable << G4endl;
304 #endif 304 #endif
305 305
306 theInteractions.push_back( anInterac 306 theInteractions.push_back( anInteraction );
307 307
308 } // End of an Interaction has happend 308 } // End of an Interaction has happend!
309 309
310 #ifdef debugFTFparticipant 310 #ifdef debugFTFparticipant
311 TrNuclN++; 311 TrNuclN++;
312 #endif 312 #endif
313 313
314 } // End of while ( ( TargetNucleon = t 314 } // End of while ( ( TargetNucleon = theNucleus->GetNextNucleon() ) )
315 315
316 #ifdef debugFTFparticipant 316 #ifdef debugFTFparticipant
317 PrNuclN++; 317 PrNuclN++;
318 #endif 318 #endif
319 319
320 } // End of while ( ( ProjectileNucleon = 320 } // End of while ( ( ProjectileNucleon = theProjectileNucleus->GetNextNucleon() ) )
321 321
322 if ( theInteractions.size() != 0 ) theProj 322 if ( theInteractions.size() != 0 ) theProjectileNucleus->DoTranslation( theBeamPosition );
323 323
324 } while ( ( theInteractions.size() == 0 ) && 324 } while ( ( theInteractions.size() == 0 ) &&
325 ++loopCounter < maxNumberOfLoops ) 325 ++loopCounter < maxNumberOfLoops ); /* Loop checking, 10.08.2015, A.Ribon */
326 if ( loopCounter >= maxNumberOfLoops ) { 326 if ( loopCounter >= maxNumberOfLoops ) {
327 #ifdef debugFTFparticipant 327 #ifdef debugFTFparticipant
328 G4cout << "BAD situation: forced exit from 328 G4cout << "BAD situation: forced exit from the while loop!" << G4endl;
329 #endif 329 #endif
330 return; 330 return;
331 } 331 }
332 332
333 SortInteractionsIncT(); 333 SortInteractionsIncT();
334 ShiftInteractionTime(); 334 ShiftInteractionTime();
335 335
336 #ifdef debugFTFparticipant 336 #ifdef debugFTFparticipant
337 G4cout << G4endl << "Number of primary colli 337 G4cout << G4endl << "Number of primary collisions " << theInteractions.size()
338 << "\t Bx[fm] " << impactX/fermi << " 338 << "\t Bx[fm] " << impactX/fermi << "\t By[fm] " << impactY/fermi
339 << "\t B[fm] " << std::sqrt( sqr( imp 339 << "\t B[fm] " << std::sqrt( sqr( impactX ) + sqr( impactY ) )/fermi << G4endl
340 << "FTF participant End. ############ 340 << "FTF participant End. #######################" << G4endl << G4endl;
341 #endif 341 #endif
342 return; 342 return;
343 } 343 }
344 344
345 345
346 //============================================ 346 //============================================================================
347 347
348 bool G4FTFPartHelperForSortInT( const G4Intera 348 bool G4FTFPartHelperForSortInT( const G4InteractionContent* Int1,
349 const G4Intera 349 const G4InteractionContent* Int2 ) {
350 return Int1->GetInteractionTime() < Int2->Ge 350 return Int1->GetInteractionTime() < Int2->GetInteractionTime();
351 } 351 }
352 352
353 353
354 //============================================ 354 //============================================================================
355 355
356 void G4FTFParticipants::SortInteractionsIncT() 356 void G4FTFParticipants::SortInteractionsIncT() { // on increased T
357 if ( theInteractions.size() < 2 ) return; / 357 if ( theInteractions.size() < 2 ) return; // Avoid unnecesary work
358 std::sort( theInteractions.begin(), theInter 358 std::sort( theInteractions.begin(), theInteractions.end(), G4FTFPartHelperForSortInT );
359 } 359 }
360 360
361 361
362 //============================================ 362 //============================================================================
363 363
364 void G4FTFParticipants::ShiftInteractionTime() 364 void G4FTFParticipants::ShiftInteractionTime() {
365 G4double InitialTime = theInteractions[0]->G 365 G4double InitialTime = theInteractions[0]->GetInteractionTime();
366 for ( unsigned int i = 1; i < theInteraction 366 for ( unsigned int i = 1; i < theInteractions.size(); i++ ) {
367 G4double InterTime = theInteractions[i]->G 367 G4double InterTime = theInteractions[i]->GetInteractionTime() - InitialTime;
368 theInteractions[i]->SetInteractionTime( In 368 theInteractions[i]->SetInteractionTime( InterTime );
369 G4InteractionContent* aCollision = theInte 369 G4InteractionContent* aCollision = theInteractions[i];
370 G4VSplitableHadron* projectile = aCollisio 370 G4VSplitableHadron* projectile = aCollision->GetProjectile();
371 G4VSplitableHadron* target = aCollision->G 371 G4VSplitableHadron* target = aCollision->GetTarget();
372 G4ThreeVector prPosition = projectile->Get 372 G4ThreeVector prPosition = projectile->GetPosition();
373 prPosition.setZ( target->GetPosition().z() 373 prPosition.setZ( target->GetPosition().z() );
374 projectile->SetPosition( prPosition ); 374 projectile->SetPosition( prPosition );
375 projectile->SetTimeOfCreation( InterTime ) 375 projectile->SetTimeOfCreation( InterTime );
376 target->SetTimeOfCreation( InterTime ); 376 target->SetTimeOfCreation( InterTime );
377 } 377 }
378 return; 378 return;
379 } 379 }
380 380
381 381
382 //============================================ 382 //============================================================================
383 383
384 void G4FTFParticipants::Clean() { 384 void G4FTFParticipants::Clean() {
385 for ( size_t i = 0; i < theInteractions.size 385 for ( size_t i = 0; i < theInteractions.size(); i++ ) {
386 if ( theInteractions[ i ] ) { 386 if ( theInteractions[ i ] ) {
387 delete theInteractions[ i ]; 387 delete theInteractions[ i ];
388 theInteractions[ i ] = 0; 388 theInteractions[ i ] = 0;
389 } 389 }
390 } 390 }
391 theInteractions.clear(); 391 theInteractions.clear();
392 currentInteraction = -1; 392 currentInteraction = -1;
393 } 393 }
394 394
395 395