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