| 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 // * 20 // * *
21 // * Parts of this code which have been devel 21 // * Parts of this code which have been developed by Abdel-Waged *
22 // * et al under contract (31-465) to the King 22 // * et al under contract (31-465) to the King Abdul-Aziz City for *
23 // * Science and Technology (KACST), the Natio 23 // * Science and Technology (KACST), the National Centre of *
24 // * Mathematics and Physics (NCMP), Saudi Ara 24 // * Mathematics and Physics (NCMP), Saudi Arabia. *
25 // * 25 // * *
26 // * By using, copying, modifying or distri 26 // * By using, copying, modifying or distributing the software (or *
27 // * any work based on the software) you ag 27 // * any work based on the software) you agree to acknowledge its *
28 // * use in resulting scientific publicati 28 // * use in resulting scientific publications, and indicate your *
29 // * acceptance of all terms of the Geant4 Sof 29 // * acceptance of all terms of the Geant4 Software license. *
30 // ******************************************* 30 // ********************************************************************
31 // 31 //
32 /// \file hadronic/Hadr02/src/G4UrQMD1_3Model. 32 /// \file hadronic/Hadr02/src/G4UrQMD1_3Model.cc
33 /// \brief Implementation of the G4UrQMD1_3Mod 33 /// \brief Implementation of the G4UrQMD1_3Model class
34 // 34 //
35 // 35 //
36 ////////////////////////////////////////////// 36 ///////////////////////////////////////////////////////////////////////////////
37 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 37 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
38 // 38 //
39 // MODULE: G4UrQMD1_3Model.hh 39 // MODULE: G4UrQMD1_3Model.hh
40 // 40 //
41 // Version: 0.B 41 // Version: 0.B
42 // Date: 25/01/12 42 // Date: 25/01/12
43 // Authors: Kh. Abdel-Waged and Nuha Fe 43 // Authors: Kh. Abdel-Waged and Nuha Felemban
44 // Revised by: V.V. Uzhinskii << 44 // Revised by: V.V. Uzhinskii
45 // SPONSERED BY 45 // SPONSERED BY
46 // Customer: KAUST/NCMP 46 // Customer: KAUST/NCMP
47 // Contract: 31-465 47 // Contract: 31-465
48 // 48 //
49 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 49 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
50 // << 50 //
51 #ifdef G4_USE_URQMD 51 #ifdef G4_USE_URQMD
52 52
53 # include "G4UrQMD1_3Model.hh" << 53 #include "G4UrQMD1_3Model.hh"
54 << 54 #include "G4UrQMD1_3Interface.hh"
55 # include "G4UrQMD1_3Interface.hh" <<
56 //------------------------------- 55 //-------------------------------
57 # include "G4CollisionOutput.hh" << 56 #include "globals.hh"
58 # include "G4DynamicParticle.hh" << 57 #include "G4DynamicParticle.hh"
59 # include "G4IonTable.hh" << 58 #include "G4IonTable.hh"
60 # include "G4LorentzRotation.hh" << 59 #include "G4CollisionOutput.hh"
61 # include "G4Nucleus.hh" << 60 #include "G4V3DNucleus.hh"
62 # include "G4ParticleDefinition.hh" << 61 #include "G4Track.hh"
63 # include "G4ParticleTable.hh" << 62 #include "G4Nucleus.hh"
64 # include "G4PhysicalConstants.hh" << 63 #include "G4LorentzRotation.hh"
65 # include "G4SystemOfUnits.hh" << 64
66 # include "G4Track.hh" << 65 #include "G4ParticleDefinition.hh"
67 # include "G4V3DNucleus.hh" << 66 #include "G4ParticleTable.hh"
68 # include "globals.hh" << 67 #include "G4PhysicalConstants.hh"
69 << 68 #include "G4SystemOfUnits.hh"
70 // AND-> << 69
71 # include "G4Version.hh" << 70 //AND->
72 // AND<- << 71 #include "G4Version.hh"
>> 72 //AND<-
73 //----------------new_anti 73 //----------------new_anti
74 # include "G4AntiAlpha.hh" << 74 #include "G4AntiHe3.hh"
75 # include "G4AntiDeuteron.hh" << 75 #include "G4AntiDeuteron.hh"
76 # include "G4AntiHe3.hh" << 76 #include "G4AntiTriton.hh"
77 # include "G4AntiTriton.hh" << 77 #include "G4AntiAlpha.hh"
78 //--------------------------- 78 //---------------------------
79 # include <fstream> << 79 #include <fstream>
80 # include <string> << 80 #include <string>
>> 81
>> 82 /////////////////////////////////////////////////////////////////////////////////
81 83
82 //....oooOO0OOooo........oooOO0OOooo........oo <<
83 // 84 //
84 G4UrQMD1_3Model::G4UrQMD1_3Model(const G4Strin 85 G4UrQMD1_3Model::G4UrQMD1_3Model(const G4String& nam)
85 : G4VIntraNuclearTransportModel(nam), verbos << 86 :G4VIntraNuclearTransportModel(nam), verbose(0)
86 { 87 {
>> 88
>> 89
87 if (verbose > 3) { 90 if (verbose > 3) {
88 G4cout << " >>> G4UrQMD1_3Model default co 91 G4cout << " >>> G4UrQMD1_3Model default constructor" << G4endl;
89 } 92 }
90 93
91 // <<
92 // Set the minimum and maximum range for the <<
93 94
94 // SetMinEnergy(100.0*MeV); <<
95 // SetMaxEnergy(200.0*GeV); <<
96 95
97 // << 96 //
>> 97 // Set the minimum and maximum range for the UrQMD model
>> 98
>> 99 // SetMinEnergy(100.0*MeV);
>> 100 // SetMaxEnergy(200.0*GeV);
>> 101
>> 102 //
98 103
99 // << 104 //
100 WelcomeMessage(); 105 WelcomeMessage();
101 // << 106 //
102 CurrentEvent = 0; << 107 CurrentEvent=0;
103 // << 108 //
104 109
105 InitialiseDataTables(); << 110 InitialiseDataTables();
106 111
107 // <<
108 } <<
109 112
110 //....oooOO0OOooo........oooOO0OOooo........oo << 113 //
>> 114 }
>> 115 ////////////////////////////////////////////////////////////////////////////////
>> 116 //
111 // Destructor 117 // Destructor
112 // 118 //
113 G4UrQMD1_3Model::~G4UrQMD1_3Model() {} << 119 G4UrQMD1_3Model::~G4UrQMD1_3Model (){}
114 << 120 ////////////////////////////////////////////////////////////////////////////////
115 //....oooOO0OOooo........oooOO0OOooo........oo <<
116 121
117 G4ReactionProductVector* G4UrQMD1_3Model::Prop << 122 G4ReactionProductVector* G4UrQMD1_3Model::Propagate(G4KineticTrackVector* ,
118 { << 123 G4V3DNucleus* ) {
119 return 0; 124 return 0;
120 } 125 }
121 126
122 //....oooOO0OOooo........oooOO0OOooo........oo << 127 ////////////////////////////////////////////////////////////////////////////////
123 // 128 //
124 // ApplyYourself 129 // ApplyYourself
125 // 130 //
126 // Member function to process an event, and ge 131 // Member function to process an event, and get information about the products.
127 132
128 G4HadFinalState* G4UrQMD1_3Model::ApplyYoursel << 133
129 << 134 G4HadFinalState *G4UrQMD1_3Model::ApplyYourself (
>> 135 const G4HadProjectile &theTrack, G4Nucleus &theTarget)
130 { 136 {
131 // anti_new << 137 //anti_new
132 // ------------------define anti_light_nuc << 138 // ------------------define anti_light_nucleus
133 const G4ParticleDefinition* anti_deu = G4Ant << 139 const G4ParticleDefinition* anti_deu =
134 << 140 G4AntiDeuteron::AntiDeuteron();
135 const G4ParticleDefinition* anti_he3 = G4Ant << 141
136 << 142 const G4ParticleDefinition* anti_he3=
137 const G4ParticleDefinition* anti_tri = G4Ant << 143 G4AntiHe3::AntiHe3();
138 << 144
139 const G4ParticleDefinition* anti_alp = G4Ant << 145 const G4ParticleDefinition* anti_tri=
140 //------------------------------------------ << 146 G4AntiTriton::AntiTriton();
141 // << 147
142 // The secondaries will be returned in G4Had << 148 const G4ParticleDefinition* anti_alp=
143 // initialise this. The original track will << 149 G4AntiAlpha::AntiAlpha();
144 // secondaries followed. << 150 //---------------------------------------------------
145 // << 151 //
>> 152 // The secondaries will be returned in G4HadFinalState &theResult -
>> 153 // initialise this. The original track will always be discontinued and
>> 154 // secondaries followed.
>> 155 //
146 theResult.Clear(); 156 theResult.Clear();
147 theResult.SetStatusChange(stopAndKill); 157 theResult.SetStatusChange(stopAndKill);
148 158
149 G4DynamicParticle* cascadeParticle = 0; << 159 G4DynamicParticle* cascadeParticle=0;
150 // << 160 //
151 // << 161 //
152 // Get relevant information about the projec << 162 // Get relevant information about the projectile and target (A, Z, energy/nuc,
153 // momentum, etc). << 163 // momentum, etc).
154 // << 164 //
155 << 165
156 const G4ParticleDefinition* definitionP = th << 166
157 const G4double AP = definitionP->GetBaryonNu << 167 const G4ParticleDefinition *definitionP = theTrack.GetDefinition();
158 const G4double ZP = definitionP->GetPDGCharg << 168 const G4double AP = definitionP->GetBaryonNumber();
159 G4double AT = theTarget.GetN(); << 169 const G4double ZP = definitionP->GetPDGCharge();
160 G4double ZT = theTarget.GetZ(); << 170 G4double AT = theTarget.GetN();
161 // ---------------------------------------- << 171 G4double ZT = theTarget.GetZ();
162 G4int id = definitionP->GetPDGEncoding(); / << 172 // -----------------------------------------------
163 // ----------------------------------------- << 173 G4int id=definitionP->GetPDGEncoding(); //get particle encoding
>> 174 // ------------------------------------------------
164 G4int AP1 = G4lrint(AP); 175 G4int AP1 = G4lrint(AP);
165 G4int ZP1 = G4lrint(ZP); 176 G4int ZP1 = G4lrint(ZP);
166 G4int AT1 = G4lrint(AT); 177 G4int AT1 = G4lrint(AT);
167 G4int ZT1 = G4lrint(ZT); 178 G4int ZT1 = G4lrint(ZT);
168 // G4cout<<"------ap1--=="<<AP1<<"---zp1--- << 179 // G4cout<<"------ap1--=="<<AP1<<"---zp1---=="<<ZP1<<"---id-=="<<id<< G4endl;
169 // << 180 //
170 // ***************************************** << 181 // ****************************************************************************
171 // The following is the parameters necessary << 182 // The following is the parameters necessary to initiate Uinit() and UrQMD()
172 // ----------------------------------------- << 183 // ----------------------------------------------------------------------------
173 urqmdparams_.u_sptar = 0; //! 0= normal pro << 184 urqmdparams_.u_sptar=0; //!0= normal proj/target, 1=special proj/tar
174 urqmdparams_.u_spproj = 1; // projectile is << 185 urqmdparams_.u_spproj=1; // projectile is a special particle
175 << 186
176 // new_anti << 187 //new_anti
>> 188
>> 189 if (AP1>1 ||definitionP==anti_deu ||definitionP==anti_he3
>> 190 ||definitionP==anti_tri ||definitionP==anti_alp)
>> 191 {
>> 192
>> 193 urqmdparams_.u_ap=AP1;
>> 194 urqmdparams_.u_zp=ZP1;
>> 195
>> 196 urqmdparams_.u_spproj=0;
>> 197 } else if (id==2212) { //!proton
>> 198 urqmdparams_.u_ap=1;
>> 199 urqmdparams_.u_zp=1;
>> 200
>> 201
>> 202 } else if(id==-2212){ //! anti-proton
>> 203 urqmdparams_.u_ap=-1;
>> 204 urqmdparams_.u_zp=-1;
>> 205 } else if(id==2112){ //! neutron
>> 206 urqmdparams_.u_ap=1;
>> 207 urqmdparams_.u_zp=-1;
>> 208
>> 209 } else if(id==-2112){ //! anti-neutron
>> 210 urqmdparams_.u_ap=-1;
>> 211 urqmdparams_.u_zp=1;
>> 212
>> 213 } else if(id==211) { //! pi+
>> 214 urqmdparams_.u_ap=101;
>> 215 urqmdparams_.u_zp=2;
>> 216 } else if(id==-211) { //! pi-
>> 217 urqmdparams_.u_ap=101;
>> 218 urqmdparams_.u_zp=-2;
>> 219 } else if(id==321) { //! K+
>> 220 urqmdparams_.u_ap=106;
>> 221 urqmdparams_.u_zp=1;
>> 222 } else if(id==-321) { //! K-
>> 223 urqmdparams_.u_ap=-106;
>> 224 urqmdparams_.u_zp=-1;
>> 225 } else if(id==130 || id==310) { // ! K0
>> 226 urqmdparams_.u_ap=106;
>> 227 urqmdparams_.u_zp=-1;
>> 228 } else if(id==-130 || id==-310){ // ! K0bar
>> 229 urqmdparams_.u_ap=-106;
>> 230 urqmdparams_.u_zp=1;
>> 231 } else {
>> 232
>> 233 G4cout << " Sorry, No definition for particle for UrQMD::"
>> 234 <<id<< "found" << G4endl;
>> 235
>> 236 //AND->
>> 237 #if G4VERSION_NUMBER>=950
>> 238 //New signature (9.5) for G4Exception
>> 239 //Using G4HadronicException
>> 240 throw G4HadronicException(__FILE__,__LINE__,
>> 241 "Sorry, no definition for particle for UrQMD");
>> 242 #else
>> 243 G4Exception(" ");
>> 244 #endif
>> 245 //AND<-
>> 246 } //end if id
>> 247 //-------------------------------------------------------
>> 248
>> 249 urqmdparams_.u_at=AT1; // Target identified
>> 250 urqmdparams_.u_zt=ZT1;
>> 251 //----------------------------------------------------
>> 252 // identify Energy
>> 253 //
>> 254 G4ThreeVector Pbefore = theTrack.Get4Momentum().vect();
>> 255 G4double T = theTrack.GetKineticEnergy();
>> 256 G4double E = theTrack.GetTotalEnergy();
>> 257 G4double TotalEbefore = E*AP1 +
>> 258 theTarget.AtomicMass(AT1, ZT1) + theTarget.GetEnergyDeposit();
>> 259 // -----------------------------------------------------------------
>> 260
>> 261 if (AP1>1) {
>> 262 urqmdparams_.u_elab=T/(AP1*GeV); // Units are GeV/nuc for UrQMD
>> 263
>> 264 E = E/AP1; // Units are GeV/nuc
>> 265
>> 266 } else {
>> 267
>> 268 urqmdparams_.u_elab=T/GeV; //units are GeV
>> 269
>> 270 TotalEbefore = E +
>> 271 theTarget.AtomicMass(AT1, ZT1) + theTarget.GetEnergyDeposit();
>> 272 }
>> 273
>> 274
>> 275 //------------------------------------------------------------
>> 276 // identify impact parameter
>> 277 urqmdparams_.u_imp=-(1.1 * std::pow(G4double(AT1),(1./3.)));
>> 278 //units are in fm for UrQMD;
>> 279 //------------------------------------------------------------
>> 280 ///////////////////////// initialise/////////////////////
177 281
178 if (AP1 > 1 || definitionP == anti_deu || de << 282 if (CurrentEvent==0)
179 || definitionP == anti_alp) << 283 {
180 { << 284 G4cout << "\n creation of table, wait-------"<<G4endl;
181 urqmdparams_.u_ap = AP1; <<
182 urqmdparams_.u_zp = ZP1; <<
183 285
184 urqmdparams_.u_spproj = 0; << 286 G4cout << "\n"<<G4endl;
185 } <<
186 else if (id == 2212) { //! proton <<
187 urqmdparams_.u_ap = 1; <<
188 urqmdparams_.u_zp = 1; <<
189 } <<
190 else if (id == -2212) { //! anti-proton <<
191 urqmdparams_.u_ap = -1; <<
192 urqmdparams_.u_zp = -1; <<
193 } <<
194 else if (id == 2112) { //! neutron <<
195 urqmdparams_.u_ap = 1; <<
196 urqmdparams_.u_zp = -1; <<
197 } <<
198 else if (id == -2112) { //! anti-neutron <<
199 urqmdparams_.u_ap = -1; <<
200 urqmdparams_.u_zp = 1; <<
201 } <<
202 else if (id == 211) { //! pi+ <<
203 urqmdparams_.u_ap = 101; <<
204 urqmdparams_.u_zp = 2; <<
205 } <<
206 else if (id == -211) { //! pi- <<
207 urqmdparams_.u_ap = 101; <<
208 urqmdparams_.u_zp = -2; <<
209 } <<
210 else if (id == 321) { //! K+ <<
211 urqmdparams_.u_ap = 106; <<
212 urqmdparams_.u_zp = 1; <<
213 } <<
214 else if (id == -321) { //! K- <<
215 urqmdparams_.u_ap = -106; <<
216 urqmdparams_.u_zp = -1; <<
217 } <<
218 else if (id == 130 || id == 310) { // ! K0 <<
219 urqmdparams_.u_ap = 106; <<
220 urqmdparams_.u_zp = -1; <<
221 } <<
222 else if (id == -130 || id == -310) { // ! K <<
223 urqmdparams_.u_ap = -106; <<
224 urqmdparams_.u_zp = 1; <<
225 } <<
226 else { <<
227 G4cout << " Sorry, No definition for parti <<
228 287
229 // AND-> << 288 G4int io=0;
230 # if G4VERSION_NUMBER >= 950 <<
231 // New signature (9.5) for G4Exception <<
232 // Using G4HadronicException <<
233 throw G4HadronicException(__FILE__, __LINE <<
234 # else <<
235 G4Exception(" "); <<
236 # endif <<
237 // AND<- <<
238 } // end if id <<
239 //------------------------------------------ <<
240 <<
241 urqmdparams_.u_at = AT1; // Target identifi <<
242 urqmdparams_.u_zt = ZT1; <<
243 //------------------------------------------ <<
244 // identify Energy <<
245 // <<
246 G4ThreeVector Pbefore = theTrack.Get4Momentu <<
247 G4double T = theTrack.GetKineticEnergy(); <<
248 G4double E = theTrack.GetTotalEnergy(); <<
249 G4double TotalEbefore = E * AP1 + theTarget. <<
250 // -------------------------------------- <<
251 289
252 if (AP1 > 1) { << 290 uinit_ (&io);
253 urqmdparams_.u_elab = T / (AP1 * GeV); // <<
254 291
255 E = E / AP1; // Units are GeV/nuc <<
256 } <<
257 else { <<
258 urqmdparams_.u_elab = T / GeV; // units a <<
259 292
260 TotalEbefore = E + theTarget.AtomicMass(AT << 293 G4cout << "\n end to create table "<<G4endl;
261 } <<
262 294
263 //------------------------------------------ << 295 CurrentEvent=1;
264 // identify impact parameter << 296 }
265 urqmdparams_.u_imp = -(1.1 * std::pow(G4doub << 297 ////////////////////////////////////////////////////////
266 // units are in fm for UrQMD; <<
267 //------------------------------------------ <<
268 ///////////////////////// initialise//////// <<
269 298
270 if (CurrentEvent == 0) { << 299 //#ifdef debug_G4UrQMD1_3Model
271 G4cout << "\n creation of table, wait----- <<
272 300
273 G4cout << "\n" << G4endl; << 301 G4cout <<"UrQMDModel running-------------" <<G4endl;
274 302
275 G4int io = 0; << 303 urqmd_ ();
276 304
277 uinit_(&io); << 305 //#endif
278 306
279 G4cout << "\n end to create table " << G4 << 307 //G4cout <<"Number of produced particles: " <<sys_.npart<<G4endl;
280 308
281 CurrentEvent = 1; << 309 G4int n = sys_.npart; //no of produced particles
282 } << 310 if (n<2)
283 //////////////////////////////////////////// << 311 {
>> 312 G4cout <<"===============Warning================"<<G4endl;
>> 313 G4cout <<"======================================"<<G4endl;
284 314
285 // #ifdef debug_G4UrQMD1_3Model << 315 G4cout <<"Number of produced particles is very low: " <<sys_.npart<<G4endl;
>> 316 G4cout <<"============================================"<<G4endl;
286 317
287 G4cout << "UrQMDModel running-------------" << 318 //AND->
>> 319 #if G4VERSION_NUMBER>=950
>> 320 //New signature (9.5) for G4Exception
>> 321 //Using G4HadronicException instead of base class
>> 322 throw G4HadronicException(__FILE__,__LINE__,
>> 323 "Number of produced particle is very low");
>> 324 #else
>> 325 G4Exception(" "); //stop
>> 326 #endif
>> 327 //AND<-
>> 328 } else {
>> 329 for (G4int i=0; i<n; i++)
>> 330 {
>> 331
288 332
289 urqmd_(); << 333 G4int pid=pdgid_ (&isys_.ityp[i], &isys_.iso3[i]);
290 334
291 // #endif << 335 // Particle is a final state secondary and not a nucleus.
>> 336 // Determine what this secondary particle is, and if valid, load dynamic
>> 337 // parameters.
>> 338 //
292 339
293 // G4cout <<"Number of produced particles: <<
294 340
295 G4int n = sys_.npart; // no of produced par << 341 G4ParticleDefinition* pd=
296 if (n < 2) { << 342 G4ParticleTable::GetParticleTable()->FindParticle(pid);
297 G4cout << "===============Warning========= <<
298 G4cout << "=============================== <<
299 343
300 G4cout << "Number of produced particles is << 344 if (pd)
301 G4cout << "=============================== << 345 {
>> 346 G4double px = (coor_.px[i]+ffermi_.ffermpx[i])* GeV;
>> 347 //units are in MeV/c for G4
>> 348 G4double py = (coor_.py[i]+ffermi_.ffermpy[i])* GeV;
>> 349 G4double pz = (coor_.pz[i]+ffermi_.ffermpz[i])* GeV;
302 350
303 // AND-> << 351 G4double et = (coor_.p0[i]) *GeV;
304 # if G4VERSION_NUMBER >= 950 <<
305 // New signature (9.5) for G4Exception <<
306 // Using G4HadronicException instead of ba <<
307 throw G4HadronicException(__FILE__, __LINE <<
308 # else <<
309 G4Exception(" "); // stop <<
310 # endif <<
311 // AND<- <<
312 } <<
313 else { <<
314 for (G4int i = 0; i < n; i++) { <<
315 G4int pid = pdgid_(&isys_.ityp[i], &isys <<
316 352
317 // Particle is a final state secondary a <<
318 // Determine what this secondary particl <<
319 // parameters. <<
320 // <<
321 353
322 G4ParticleDefinition* pd = G4ParticleTab << 354 // ------------------------------Use only "Lorentz vector"----------
323 355
324 if (pd) { << 356 G4LorentzVector lorenzvec = G4LorentzVector(px,py,pz,et);
325 G4double px = (coor_.px[i] + ffermi_.f <<
326 // units are in MeV/c for G4 <<
327 G4double py = (coor_.py[i] + ffermi_.f <<
328 G4double pz = (coor_.pz[i] + ffermi_.f <<
329 357
330 G4double et = (coor_.p0[i]) * GeV; << 358 cascadeParticle = new G4DynamicParticle(pd, lorenzvec); //
331 359
332 // ------------------------------Us << 360 theResult.AddSecondary(cascadeParticle);
333 361
334 G4LorentzVector lorenzvec = G4LorentzV << 362 //======================================================================
335 363
336 cascadeParticle = new G4DynamicParticl <<
337 364
338 theResult.AddSecondary(cascadeParticle <<
339 365
340 //==================================== << 366 } //if
>> 367 } //for
341 368
342 } // if << 369 } //if warning
343 } // for <<
344 370
345 } // if warning << 371 //=======================================================================
>> 372 if (verbose >= 3) {
346 373
347 //========================================== << 374 //
348 if (verbose >= 3) { <<
349 // <<
350 G4double TotalEafter = 0.0; 375 G4double TotalEafter = 0.0;
351 G4ThreeVector TotalPafter; 376 G4ThreeVector TotalPafter;
352 G4double charge = 0.0; << 377 G4double charge = 0.0;
353 G4int baryon = 0; << 378 G4int baryon = 0;
354 G4int nSecondaries = theResult.GetNumberOf << 379 G4int nSecondaries = theResult.GetNumberOfSecondaries();
355 << 380
356 for (G4int j = 0; j < nSecondaries; j++) { << 381 for (G4int j=0; j<nSecondaries; j++) {
357 TotalEafter += theResult.GetSecondary(j) << 382 TotalEafter += theResult.GetSecondary(j)->
>> 383 GetParticle()->GetTotalEnergy();
358 384
359 TotalPafter += theResult.GetSecondary(j) << 385 TotalPafter += theResult.GetSecondary(j)->
>> 386 GetParticle()->GetMomentum();
360 387
361 G4ParticleDefinition* pd = theResult.Get << 388 G4ParticleDefinition *pd = theResult.GetSecondary(j)->
>> 389 GetParticle()->GetDefinition();
362 390
363 charge += pd->GetPDGCharge(); 391 charge += pd->GetPDGCharge();
364 baryon += pd->GetBaryonNumber(); 392 baryon += pd->GetBaryonNumber();
365 393
366 } // for secondaries << 394 } //for secondaries
367 395
368 G4cout << "------------------------------- << 396 G4cout <<"----------------------------------------"
369 << "------------------------------- << 397 <<"----------------------------------------"
370 G4cout << "Total energy before collision << 398 <<G4endl;
371 << " MeV" << G4endl; << 399 G4cout <<"Total energy before collision = " <<TotalEbefore ///MeV
372 G4cout << "Total energy after collision << 400 <<" MeV" <<G4endl;
373 << " MeV" << G4endl; << 401 G4cout <<"Total energy after collision = " <<TotalEafter //MeV
374 << 402 <<" MeV" <<G4endl;
375 G4cout << "------------------------------- << 403
376 << 404 G4cout <<"----------------------------------------"<<G4endl;
377 G4cout << "Total momentum before collision << 405
378 << " MeV/c" << G4endl; << 406 G4cout <<"Total momentum before collision = " <<Pbefore //MeV
379 G4cout << "Total momentum after collision << 407 <<" MeV/c" <<G4endl;
380 << " MeV/c" << G4endl; << 408 G4cout <<"Total momentum after collision = " <<TotalPafter //MeV
381 G4cout << "------------------------------- << 409 <<" MeV/c" <<G4endl;
>> 410 G4cout <<"----------------------------------------"<<G4endl;
382 411
383 if (verbose >= 4) { 412 if (verbose >= 4) {
384 G4cout << "Total charge before collision << 413 G4cout <<"Total charge before collision = " <<(ZP+ZT)*eplus
385 G4cout << "Total charge after collision << 414 <<G4endl;
>> 415 G4cout <<"Total charge after collision = " <<charge
>> 416 <<G4endl;
386 417
387 G4cout << "----------------------------- << 418 G4cout <<"----------------------------------------"<<G4endl;
388 419
389 G4cout << "Total baryon number before co << 420 G4cout <<"Total baryon number before collision = "<<AP+AT
390 G4cout << "Total baryon number after col << 421 <<G4endl;
391 G4cout << "----------------------------- << 422 G4cout <<"Total baryon number after collision = "<<baryon
>> 423 <<G4endl;
>> 424 G4cout <<"----------------------------------------"<<G4endl;
392 425
393 } // if verbose4 << 426 } //if verbose4
394 427
395 G4cout << "------------------------------- << 428 G4cout <<"----------------------------------------"
396 << "------------------------------- << 429 <<"----------------------------------------"
>> 430 <<G4endl;
397 431
398 } // if verbose3 << 432 } //if verbose3
399 433
400 return &theResult; <<
401 } // G4hadfinal <<
402 434
403 //....oooOO0OOooo........oooOO0OOooo........oo << 435 return &theResult;
>> 436 } //G4hadfinal
>> 437
>> 438
>> 439 //---------------------------------------------------------------------
>> 440
>> 441 //---------------------------------------------------------------------
404 // 442 //
405 // WelcomeMessage 443 // WelcomeMessage
406 // 444 //
407 void G4UrQMD1_3Model::WelcomeMessage() const << 445 void G4UrQMD1_3Model::WelcomeMessage () const
408 { 446 {
409 G4cout << G4endl; << 447 G4cout <<G4endl;
410 G4cout << " ******************************** << 448 G4cout <<" *****************************************************************"
411 G4cout << " Interface to G4UrQMD_1.3 << 449 <<G4endl;
412 G4cout << " Version number : 00.00.0B << 450 G4cout <<" Interface to G4UrQMD_1.3 activated"
413 G4cout << " (Interface written by Kh. Abdel- << 451 <<G4endl;
414 G4cout << G4endl; << 452 G4cout <<" Version number : 00.00.0B File date : 25/01/12" <<G4endl;
415 G4cout << " ******************************** << 453 G4cout <<" (Interface written by Kh. Abdel-Waged et al. for the KACST/NCMP)"
>> 454 <<G4endl;
>> 455 G4cout <<G4endl;
>> 456 G4cout <<" *****************************************************************"
>> 457 <<G4endl;
416 G4cout << G4endl; 458 G4cout << G4endl;
417 459
418 return; 460 return;
419 } 461 }
420 462
421 //....oooOO0OOooo........oooOO0OOooo........oo << 463 void G4UrQMD1_3Model::InitialiseDataTables ()
422 <<
423 void G4UrQMD1_3Model::InitialiseDataTables() <<
424 { 464 {
425 // << 465 //
426 // << 466 //
427 // The next line is to make sure the block d << 467 // The next line is to make sure the block data statements are
428 // executed. << 468 // executed.
429 // << 469 //
>> 470
>> 471 g4urqmdblockdata_ ();
>> 472
430 473
431 g4urqmdblockdata_(); << 474 ///////////////////////////////////////////////////
>> 475 /////// Dynamic seed //////////////////////////////
>> 476 //G4int ranseed=-time_ ();
>> 477 // Fixed seed ///////////////////////////
432 478
433 //////////////////////////////////////////// << 479 G4int ranseed=1097569630;
434 /////// Dynamic seed /////////////////////// <<
435 // G4int ranseed=-time_ (); <<
436 // Fixed seed ///////////////////////// <<
437 480
438 G4int ranseed = 1097569630; << 481 G4cout <<"\n seed: "<<ranseed<<G4endl;
439 482
440 G4cout << "\n seed: " << ranseed << G4endl; << 483 sseed_ (&ranseed);
441 484
442 sseed_(&ranseed); << 485 loginit_();
443 486
444 loginit_(); <<
445 } 487 }
446 488
447 #endif // G4_USE_URQMD << 489 #endif //G4_USE_URQMD
448 490