Geant4 Cross Reference |
1 // 1 2 // ******************************************* 3 // * License and Disclaimer 4 // * 5 // * The Geant4 software is copyright of th 6 // * the Geant4 Collaboration. It is provided 7 // * conditions of the Geant4 Software License 8 // * LICENSE and available at http://cern.ch/ 9 // * include a list of copyright holders. 10 // * 11 // * Neither the authors of this software syst 12 // * institutes,nor the agencies providing fin 13 // * work make any representation or warran 14 // * regarding this software system or assum 15 // * use. Please see the license in the file 16 // * for the full disclaimer and the limitatio 17 // * 18 // * This code implementation is the result 19 // * technical work of the GEANT4 collaboratio 20 // * By using, copying, modifying or distri 21 // * any work based on the software) you ag 22 // * use in resulting scientific publicati 23 // * acceptance of all terms of the Geant4 Sof 24 // ******************************************* 25 // 26 // ------------------------------------------- 27 // 28 // GEANT4 Class file 29 // 30 // 31 // File name: G4BraggIonModel 32 // 33 // Author: Vladimir Ivanchenko 34 // 35 // Creation date: 13.10.2004 36 // 37 // Modifications: 38 // 11-05-05 Major optimisation of internal int 39 // 29-11-05 Do not use G4Alpha class (V.Ivantc 40 // 15-02-06 ComputeCrossSectionPerElectron, Co 41 // 25-04-06 Add stopping data from ASTAR (V.Iv 42 // 23-10-06 Reduce lowestKinEnergy to 0.25 keV 43 // 12-08-08 Added methods GetParticleCharge, G 44 // CorrectionsAlongStep needed for io 45 // 46 47 // Class Description: 48 // 49 // Implementation of energy loss and delta-ele 50 // slow charged heavy particles 51 52 // ------------------------------------------- 53 // 54 55 //....oooOO0OOooo........oooOO0OOooo........oo 56 //....oooOO0OOooo........oooOO0OOooo........oo 57 58 #include "G4BraggIonModel.hh" 59 #include "G4PhysicalConstants.hh" 60 #include "G4SystemOfUnits.hh" 61 #include "Randomize.hh" 62 #include "G4Electron.hh" 63 #include "G4ParticleChangeForLoss.hh" 64 #include "G4EmCorrections.hh" 65 #include "G4DeltaAngle.hh" 66 #include "G4ICRU90StoppingData.hh" 67 #include "G4ASTARStopping.hh" 68 #include "G4PSTARStopping.hh" 69 #include "G4NistManager.hh" 70 #include "G4Log.hh" 71 #include "G4Exp.hh" 72 #include "G4AutoLock.hh" 73 74 //....oooOO0OOooo........oooOO0OOooo........oo 75 76 G4ASTARStopping* G4BraggIonModel::fASTAR = nul 77 78 namespace 79 { 80 G4Mutex alphaMutex = G4MUTEX_INITIALIZER; 81 } 82 83 G4BraggIonModel::G4BraggIonModel(const G4Parti 84 const G4Strin 85 : G4BraggModel(p, nam) 86 { 87 HeMass = 3.727417*CLHEP::GeV; 88 massFactor = 1000.*CLHEP::amu_c2/HeMass; 89 } 90 91 //....oooOO0OOooo........oooOO0OOooo........oo 92 93 G4BraggIonModel::~G4BraggIonModel() 94 { 95 if(isFirstAlpha) { 96 delete fASTAR; 97 fASTAR = nullptr; 98 } 99 } 100 101 //....oooOO0OOooo........oooOO0OOooo........oo 102 103 void G4BraggIonModel::Initialise(const G4Parti 104 const G4DataV 105 { 106 G4BraggModel::Initialise(p, ref); 107 const G4String& pname = particle->GetParticl 108 if(pname == "alpha") { isAlpha = true; } 109 if(isAlpha && fASTAR == nullptr) { 110 G4AutoLock l(&alphaMutex); 111 if(fASTAR == nullptr) { 112 isFirstAlpha = true; 113 fASTAR = new G4ASTARStopping(); 114 } 115 l.unlock(); 116 } 117 if(isFirstAlpha) { 118 fASTAR->Initialise(); 119 } 120 } 121 122 123 //....oooOO0OOooo........oooOO0OOooo........oo 124 125 G4double G4BraggIonModel::GetChargeSquareRatio 126 127 128 { 129 // this method is called only for ions, so n 130 if(isAlpha) { return 1.0; } 131 return G4BraggModel::GetChargeSquareRatio(p, 132 } 133 134 //....oooOO0OOooo........oooOO0OOooo........oo 135 136 G4double G4BraggIonModel::ComputeCrossSectionP 137 con 138 139 140 141 142 { 143 G4double sigma = 144 Z*ComputeCrossSectionPerElectron(p,kinEner 145 if(isAlpha) { 146 sigma *= (HeEffChargeSquare(Z, kinEnergy/C 147 } 148 return sigma; 149 } 150 151 //....oooOO0OOooo........oooOO0OOooo........oo 152 153 G4double G4BraggIonModel::CrossSectionPerVolum 154 con 155 con 156 157 158 159 { 160 G4double sigma = material->GetElectronDensit 161 ComputeCrossSectionPerElectron(p,kinEnergy 162 if(isAlpha) { 163 const G4double zeff = material->GetTotNbOf 164 material->GetTotNbOfAtomsPerVolume(); 165 sigma *= (HeEffChargeSquare(zeff, kinEnerg 166 } 167 return sigma; 168 } 169 170 //....oooOO0OOooo........oooOO0OOooo........oo 171 172 G4double G4BraggIonModel::ComputeDEDXPerVolume 173 174 175 176 { 177 const G4double tmax = MaxSecondaryEnergy(p, 178 const G4double tlim = lowestKinEnergy*massRa 179 const G4double tmin = std::max(std::min(cut, 180 G4double dedx = 0.0; 181 182 if(kineticEnergy < tlim) { 183 dedx = HeDEDX(material, tlim)*std::sqrt(ki 184 } else { 185 dedx = HeDEDX(material, kineticEnergy); 186 187 if (tmin < tmax) { 188 const G4double tau = kineticEnergy/mass; 189 const G4double x = tmin/tmax; 190 191 G4double del = 192 (G4Log(x)*(tau + 1.)*(tau + 1.)/(tau * 193 CLHEP::twopi_mc2_rcl2*material->GetElectronD 194 if(isAlpha) { 195 const G4double zeff = material->GetTotNbOfEl 196 material->GetTotNbOfAtomsPerVolume(); 197 heChargeSquare = HeEffChargeSquare(zeff, kin 198 del *= heChargeSquare; 199 } 200 dedx += del; 201 } 202 } 203 dedx = std::max(dedx, 0.0); 204 /* 205 G4cout << "BraggIon: " << material->GetNam 206 << " E(MeV)=" << kineticEnergy/MeV 207 << " Tmin(MeV)=" << tmin << " dedx( 208 << dedx*gram/(MeV*cm2*material->Get 209 << " q2=" << chargeSquare << G4endl 210 */ 211 return dedx; 212 } 213 214 //....oooOO0OOooo........oooOO0OOooo........oo 215 216 void G4BraggIonModel::CorrectionsAlongStep(con 217 con 218 con 219 G4d 220 { 221 // no correction for alpha 222 if(isAlpha) { return; } 223 224 // no correction at a small step at the last 225 const G4double preKinEnergy = dp->GetKinetic 226 if(eloss >= preKinEnergy || eloss < preKinEn 227 228 // corrections only for ions 229 const G4ParticleDefinition* p = dp->GetDefin 230 if(p != particle) { SetParticle(p); } 231 232 // effective energy and charge at a step 233 const G4Material* mat = couple->GetMaterial( 234 const G4double e = std::max(preKinEnergy - e 235 const G4double q20 = corr->EffectiveChargeSq 236 const G4double q2 = corr->EffectiveChargeSqu 237 const G4double qfactor = q2/q20; 238 /* 239 G4cout << "G4BraggIonModel::CorrectionsAlo 240 << preKinEnergy << " Eeff(MeV)=" << e 241 << " eloss=" << eloss << " elossnew=" << e 242 << " qfactor=" << qfactor << " Qpre=" << q 243 << p->GetParticleName() <<G4endl; 244 */ 245 eloss *= qfactor; 246 } 247 248 //....oooOO0OOooo........oooOO0OOooo........oo 249 250 G4int G4BraggIonModel::HasMaterialForHe(const 251 { 252 const G4String& chFormula = mat->GetChemical 253 if(chFormula.empty()) { return -1; } 254 255 // ICRU Report N49, 1993. Ziegler model for 256 257 static const G4int numberOfMolecula = 11; 258 static const G4String molName[numberOfMolecu 259 "CaF_2", "Cellulose_Nitrate", "LiF", "Po 260 "(C_2H_4)_N-Polyethylene", "(C_2H_4)_N-Po 261 "Polysterene", "SiO_2", "NaI", "H_2O", 262 "Graphite" }; 263 264 // Search for the material in the table 265 for (G4int i=0; i<numberOfMolecula; ++i) { 266 if (chFormula == molName[i]) { 267 return i; 268 } 269 } 270 return -1; 271 } 272 273 //....oooOO0OOooo........oooOO0OOooo........oo 274 275 G4double G4BraggIonModel::HeStoppingPower(cons 276 { 277 G4double ionloss = 0.0; 278 if (iMolecula >= 0) { 279 280 // The data and the fit from: 281 // ICRU Report N49, 1993. Ziegler's model 282 // He energy in internal units of parametr 283 // Input scaled energy of a proton or Gene 284 G4double T = kineticEnergy/(massRate*CLHEP 285 286 static const G4float a[11][5] = { 287 {9.43672f, 0.54398f, 84.341f, 1.3705f, 288 {67.1503f, 0.41409f, 404.512f, 148.97f, 289 {5.11203f, 0.453f, 36.718f, 50.6f, 290 {61.793f, 0.48445f, 361.537f, 57.889f, 291 {7.83464f, 0.49804f, 160.452f, 3.192f, 292 {19.729f, 0.52153f, 162.341f, 58.35f, 293 {26.4648f, 0.50112f, 188.913f, 30.079f, 294 {7.8655f, 0.5205f, 63.96f, 51.32f, 295 {8.8965f, 0.5148f, 339.36f, 1.7205f, 296 {2.959f, 0.53255f, 34.247f, 60.655f, 297 {3.80133f, 0.41590f, 12.9966f, 117.83f, 298 299 static const G4double atomicWeight[11] = { 300 101.96128f, 44.0098f, 16.0426f, 28.0536 301 104.1512f, 44.665f, 60.0843f, 18.0152 302 303 const G4int i = iMolecula; 304 305 G4double slow = (G4double)(a[i][0]); 306 307 G4double x1 = (G4double)(a[i][1]); 308 G4double x2 = (G4double)(a[i][2]); 309 G4double x3 = (G4double)(a[i][3]); 310 G4double x4 = (G4double)(a[i][4]); 311 312 // Free electron gas model 313 if ( T < 0.001 ) { 314 G4double shigh = G4Log( 1.0 + x3*1000.0 315 ionloss = slow*shigh / (slow + shigh) ; 316 ionloss *= std::sqrt(T*1000.0) ; 317 318 // Main parametrisation 319 } else { 320 slow *= G4Exp(G4Log(T*1000.0)*x1) ; 321 G4double shigh = G4Log( 1.0 + x3/T + x4* 322 ionloss = slow*shigh / (slow + shigh) ; 323 /* 324 G4cout << "## " << i << ". T= " << T 325 << " a0= " << a[i][0] << " a1= " << a 326 << " shigh= " << shigh 327 << " dedx= " << ionloss << " q^2= " < 328 << G4endl; 329 */ 330 } 331 ionloss = std::max(ionloss, 0.0) * atomicW 332 } 333 return ionloss; 334 } 335 336 //....oooOO0OOooo........oooOO0OOooo........oo 337 338 G4double G4BraggIonModel::HeElectronicStopping 339 const G4double kinet 340 { 341 G4double ionloss ; 342 G4int i = std::min(z-1, 91); // index of at 343 //G4cout << "ElectronicStoppingPower z=" << 344 // << " E=" << kineticEnergy << G4endl; 345 // The data and the fit from: 346 // ICRU Report 49, 1993. Ziegler's type of p 347 // Proton kinetic energy for parametrisation 348 // He energy in internal units of parametris 349 //G4double T = kineticEnergy*rateMassHe2p/CL 350 G4double T = kineticEnergy/CLHEP::MeV; 351 352 static const G4float a[92][5] = { 353 { 0.35485f, 0.6456f, 6.01525f, 20.8933f, 354 },{ 0.58f, 0.59f, 6.3f, 130.0f, 355 },{ 1.42f, 0.49f, 12.25f, 32.0f, 356 },{ 2.206f, 0.51f, 15.32f, 0.25f, 357 // },{ 2.1895f, 0.47183,7.2362f, 134 358 },{ 3.691f, 0.4128f, 18.48f, 50.72f, 359 },{ 3.83523f, 0.42993f,12.6125f, 227.41f, 360 // },{ 1.9259f, 0.5550f, 27.15125f, 26 361 },{ 1.9259f, 0.5550f, 27.1513f, 26.0665f, 362 },{ 2.81015f, 0.4759f, 50.0253f, 10.556f, 363 },{ 1.533f, 0.531f, 40.44f, 18.41f, 364 },{ 2.303f, 0.4861f, 37.01f, 37.96f, 365 // Z= 11-20 366 },{ 9.894f, 0.3081f, 23.65f, 0.384f, 367 },{ 4.3f, 0.47f, 34.3f, 3.3f, 368 },{ 2.5f, 0.625f, 45.7f, 0.1f, 369 },{ 2.1f, 0.65f, 49.34f, 1.788f, 370 },{ 1.729f, 0.6562f, 53.41f, 2.405f, 371 },{ 1.402f, 0.6791f, 58.98f, 3.528f, 372 },{ 1.117f, 0.7044f, 69.69f, 3.705f, 373 },{ 2.291f, 0.6284f, 73.88f, 4.478f, 374 },{ 8.554f, 0.3817f, 83.61f, 11.84f, 375 },{ 6.297f, 0.4622f, 65.39f, 10.14f, 376 // Z= 21-30 377 },{ 5.307f, 0.4918f, 61.74f, 12.4f, 378 },{ 4.71f, 0.5087f, 65.28f, 8.806f, 379 },{ 6.151f, 0.4524f, 83.0f, 18.31f, 380 },{ 6.57f, 0.4322f, 84.76f, 15.53f, 381 },{ 5.738f, 0.4492f, 84.6f, 14.18f, 382 },{ 5.013f, 0.4707f, 85.8f, 16.55f, 383 },{ 4.32f, 0.4947f, 76.14f, 10.85f, 384 },{ 4.652f, 0.4571f, 80.73f, 22.0f, 385 },{ 3.114f, 0.5236f, 76.67f, 7.62f, 386 },{ 3.114f, 0.5236f, 76.67f, 7.62f, 387 // Z= 31-40 388 },{ 3.114f, 0.5236f, 76.67f, 7.62f, 389 },{ 5.746f, 0.4662f, 79.24f, 1.185f, 390 },{ 2.792f, 0.6346f, 106.1f, 0.2986f, 391 },{ 4.667f, 0.5095f, 124.3f, 2.102f, 392 },{ 2.44f, 0.6346f, 105.0f, 0.83f, 393 },{ 1.413f, 0.7377f, 147.9f, 1.466f, 394 },{ 11.72f, 0.3826f, 102.8f, 9.231f, 395 },{ 7.126f, 0.4804f, 119.3f, 5.784f, 396 },{ 11.61f, 0.3955f, 146.7f, 7.031f, 397 },{ 10.99f, 0.41f, 163.9f, 7.1f, 398 // Z= 41-50 399 },{ 9.241f, 0.4275f, 163.1f, 7.954f, 400 },{ 9.276f, 0.418f, 157.1f, 8.038f, 401 },{ 3.999f, 0.6152f, 97.6f, 1.297f, 402 },{ 4.306f, 0.5658f, 97.99f, 5.514f, 403 },{ 3.615f, 0.6197f, 86.26f, 0.333f, 404 },{ 5.8f, 0.49f, 147.2f, 6.903f, 405 },{ 5.6f, 0.49f, 130.0f, 10.0f, 406 },{ 3.55f, 0.6068f, 124.7f, 1.112f, 407 },{ 3.6f, 0.62f, 105.8f, 0.1692f, 408 },{ 5.4f, 0.53f, 103.1f, 3.931f, 409 // Z= 51-60 410 },{ 3.97f, 0.6459f, 131.8f, 0.2233f, 411 },{ 3.65f, 0.64f, 126.8f, 0.6834f, 412 },{ 3.118f, 0.6519f, 164.9f, 1.208f, 413 },{ 3.949f, 0.6209f, 200.5f, 1.878f, 414 },{ 14.4f, 0.3923f, 152.5f, 8.354f, 415 },{ 10.99f, 0.4599f, 138.4f, 4.811f, 416 },{ 16.6f, 0.3773f, 224.1f, 6.28f, 417 },{ 10.54f, 0.4533f, 159.3f, 4.832f, 418 },{ 10.33f, 0.4502f, 162.0f, 5.132f, 419 },{ 10.15f, 0.4471f, 165.6f, 5.378f, 420 // Z= 61-70 421 },{ 9.976f, 0.4439f, 168.0f, 5.721f, 422 },{ 9.804f, 0.4408f, 176.2f, 5.675f, 423 },{ 14.22f, 0.363f, 228.4f, 7.024f, 424 },{ 9.952f, 0.4318f, 233.5f, 5.065f, 425 },{ 9.272f, 0.4345f, 210.0f, 4.911f, 426 },{ 10.13f, 0.4146f, 225.7f, 5.525f, 427 },{ 8.949f, 0.4304f, 213.3f, 5.071f, 428 },{ 11.94f, 0.3783f, 247.2f, 6.655f, 429 },{ 8.472f, 0.4405f, 195.5f, 4.051f, 430 },{ 8.301f, 0.4399f, 203.7f, 3.667f, 431 // Z= 71-80 432 },{ 6.567f, 0.4858f, 193.0f, 2.65f, 433 },{ 5.951f, 0.5016f, 196.1f, 2.662f, 434 },{ 7.495f, 0.4523f, 251.4f, 3.433f, 435 },{ 6.335f, 0.4825f, 255.1f, 2.834f, 436 },{ 4.314f, 0.5558f, 214.8f, 2.354f, 437 },{ 4.02f, 0.5681f, 219.9f, 2.402f, 438 },{ 3.836f, 0.5765f, 210.2f, 2.742f, 439 },{ 4.68f, 0.5247f, 244.7f, 2.749f, 440 },{ 2.892f, 0.6204f, 208.6f, 2.415f, 441 // },{ 3.223f, 0.5883f, 232.7f, 2.9 442 },{ 2.892f, 0.6204f, 208.6f, 2.415f, 443 // Z= 81-90 444 },{ 4.728f, 0.5522f, 217.0f, 3.091f, 445 },{ 6.18f, 0.52f, 170.0f, 4.0f, 446 },{ 9.0f, 0.47f, 198.0f, 3.8f, 447 },{ 2.324f, 0.6997f, 216.0f, 1.599f, 448 },{ 1.961f, 0.7286f, 223.0f, 1.621f, 449 },{ 1.75f, 0.7427f, 350.1f, 0.9789f, 450 },{ 10.31f, 0.4613f, 261.2f, 4.738f, 451 },{ 7.962f, 0.519f, 235.7f, 4.347f, 452 },{ 6.227f, 0.5645f, 231.9f, 3.961f, 453 },{ 5.246f, 0.5947f, 228.6f, 4.027f, 454 // Z= 91-92 455 },{ 5.408f, 0.5811f, 235.7f, 3.961f, 456 },{ 5.218f, 0.5828f, 245.0f, 3.838f, 457 }; 458 459 G4double slow = (G4double)(a[i][0]); 460 461 G4double x1 = (G4double)(a[i][1]); 462 G4double x2 = (G4double)(a[i][2]); 463 G4double x3 = (G4double)(a[i][3]); 464 G4double x4 = (G4double)(a[i][4]); 465 466 // Free electron gas model 467 if ( T < 0.001 ) { 468 G4double shigh = G4Log( 1.0 + x3*1000.0 + 469 ionloss = slow*shigh*std::sqrt(T*1000.0) 470 471 // Main parametrisation 472 } else { 473 slow *= G4Exp(G4Log(T*1000.0)*x1); 474 G4double shigh = G4Log( 1.0 + x3/T + x4*T 475 ionloss = slow*shigh / (slow + shigh) ; 476 /* 477 G4cout << "## " << i << ". T= " << T << " 478 << " a0= " << a[i][0] << " a1= " << 479 << " shigh= " << shigh 480 << " dedx= " << ionloss << " q^2= " 481 << G4endl; 482 */ 483 } 484 ionloss = std::max(ionloss, 0.0); 485 return ionloss; 486 } 487 488 //....oooOO0OOooo........oooOO0OOooo........oo 489 490 G4double G4BraggIonModel::HeDEDX(const G4Mater 491 const G4double 492 { 493 // aEnergy is energy of alpha 494 G4double eloss = 0.0; 495 // check DB 496 if(material != currentMaterial) { 497 currentMaterial = material; 498 baseMaterial = material->GetBaseMaterial() 499 ? material->GetBaseMaterial() : material 500 iPSTAR = -1; 501 iASTAR = -1; 502 iMolecula = -1; 503 iICRU90 = (nullptr != fICRU90) ? fICRU90-> 504 505 if(iICRU90 < 0) { 506 if(isAlpha) { 507 iASTAR = fASTAR->GetIndex(baseMaterial); 508 if(iASTAR < 0) { iMolecula = HasMaterialForH 509 } else { 510 iPSTAR = fPSTAR->GetIndex(baseMaterial); 511 } 512 } 513 /* 514 G4cout << "%%% " <<material->GetName() << 515 << iMolecula << " iASTAR= " << iAS 516 << " iICRU90= " << iICRU90<< G4end 517 */ 518 } 519 // ICRU90 520 if(iICRU90 >= 0) { 521 eloss = (isAlpha) 522 ? fICRU90->GetElectronicDEDXforAlpha(iIC 523 : fICRU90->GetElectronicDEDXforProton(iI 524 if(eloss > 0.0) { return eloss*material->G 525 } 526 // PSTAR parameterisation 527 if( iPSTAR >= 0 ) { 528 return fPSTAR->GetElectronicDEDX(iPSTAR, a 529 *material->GetDensity(); 530 } 531 // ASTAR 532 if( iASTAR >= 0 ) { 533 eloss = fASTAR->GetElectronicDEDX(iASTAR, 534 /* 535 G4cout << "ASTAR: E=" << aEnergy 536 << " dedx=" << eloss*material->GetDensity 537 << " " << particle->GetParticleName() << 538 */ 539 if(eloss > 0.0) { return eloss*material->G 540 } 541 542 const std::size_t numberOfElements = materia 543 const G4ElementVector* theElmVector = materi 544 const G4double* theAtomicNumDensityVector = 545 material->GetAtomicNumDensityVector(); 546 547 // molecular data use proton stopping power 548 // element data from ICRU49 include data for 549 if(iMolecula >= 0) { 550 const G4double zeff = material->GetTotNbOf 551 material->GetTotNbOfAtomsPerVolume(); 552 heChargeSquare = HeEffChargeSquare(zeff, a 553 eloss = HeStoppingPower(aEnergy)*heChargeS 554 555 // pure material 556 } else if(1 == numberOfElements) { 557 558 const G4Element* element = (*theElmVector) 559 eloss = HeElectronicStoppingPower(element- 560 * (material->GetTotNbOfAtomsPerVolume()) 561 562 // Brugg's rule calculation 563 } else { 564 // loop for the elements in the material 565 for (std::size_t i=0; i<numberOfElements; 566 const G4Element* element = (*theElmVecto 567 eloss += HeElectronicStoppingPower(eleme 568 * theAtomicNumDensityVector[i]; 569 } 570 } 571 return eloss*theZieglerFactor; 572 } 573 574 //....oooOO0OOooo........oooOO0OOooo........oo 575 576 G4double 577 G4BraggIonModel::HeEffChargeSquare(const G4dou 578 const G4dou 579 { 580 // The aproximation of He effective charge f 581 // J.F.Ziegler, J.P. Biersack, U. Littmark 582 // The Stopping and Range of Ions in Matter, 583 // Vol.1, Pergamon Press, 1985 584 585 static const G4double c[6] = {0.2865, 0.126 586 0.02402,-0.011 587 588 G4double e = std::max(0.0, G4Log(kinEnergyHe 589 G4double x = c[0] ; 590 G4double y = 1.0 ; 591 for (G4int i=1; i<6; ++i) { 592 y *= e; 593 x += y * c[i]; 594 } 595 596 G4double w = 7.6 - e ; 597 w = 1.0 + (0.007 + 0.00005*z) * G4Exp( -w*w 598 w = 4.0 * (1.0 - G4Exp(-x)) * w * w ; 599 600 return w; 601 } 602 603 //....oooOO0OOooo........oooOO0OOooo........oo 604 605