Changeset 1055 for trunk/source/processes/electromagnetic/utils/src
- Timestamp:
- May 28, 2009, 4:26:57 PM (15 years ago)
- Location:
- trunk/source/processes/electromagnetic/utils/src
- Files:
-
- 15 edited
-
G4DummyModel.cc (modified) (2 diffs)
-
G4EmCalculator.cc (modified) (2 diffs)
-
G4EmConfigurator.cc (modified) (4 diffs)
-
G4EmElementSelector.cc (modified) (7 diffs)
-
G4EmModelManager.cc (modified) (16 diffs)
-
G4EmProcessOptions.cc (modified) (3 diffs)
-
G4EnergyLossMessenger.cc (modified) (4 diffs)
-
G4LossTableBuilder.cc (modified) (3 diffs)
-
G4LossTableManager.cc (modified) (4 diffs)
-
G4VEmFluctuationModel.cc (modified) (2 diffs)
-
G4VEmModel.cc (modified) (6 diffs)
-
G4VEmProcess.cc (modified) (12 diffs)
-
G4VEnergyLossProcess.cc (modified) (23 diffs)
-
G4VMscModel.cc (modified) (5 diffs)
-
G4VMultipleScattering.cc (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/source/processes/electromagnetic/utils/src/G4DummyModel.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4DummyModel.cc,v 1. 3 2007/05/22 17:31:58vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4DummyModel.cc,v 1.4 2009/04/07 18:39:47 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 52 52 53 53 G4DummyModel::G4DummyModel(const G4String& nam) 54 : G4V EmModel(nam)54 : G4VMscModel(nam) 55 55 {} 56 56 -
trunk/source/processes/electromagnetic/utils/src/G4EmCalculator.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4EmCalculator.cc,v 1.4 4 2008/08/03 18:47:15vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4EmCalculator.cc,v 1.46 2009/02/24 09:56:03 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 755 755 if(currentProcess) currentProcessName = currentProcess->GetProcessName(); 756 756 757 if(p->GetParticleType() == "nucleus" && 758 currentParticleName != "deuteron" && currentParticleName != "triton") { 757 if(p->GetParticleType() == "nucleus" 758 && currentParticleName != "deuteron" 759 && currentParticleName != "triton" 760 && currentParticleName != "alpha+" 761 && currentParticleName != "helium" 762 && currentParticleName != "hydrogen" 763 ) { 759 764 baseParticle = theGenericIon; 760 765 massRatio = baseParticle->GetPDGMass()/p->GetPDGMass(); -
trunk/source/processes/electromagnetic/utils/src/G4EmConfigurator.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4EmConfigurator.cc,v 1. 3 2008/11/21 12:30:29vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4EmConfigurator.cc,v 1.4 2009/02/26 11:33:33 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 117 117 G4String fname = ""; 118 118 if(fm) fname = fm->GetName(); 119 AddModelForRegion(particleName, processName, mod->GetName(), regionName, 119 G4String mname = ""; 120 if(mod) mname = mod->GetName(); 121 AddModelForRegion(particleName, processName, mname, regionName, 120 122 emin, emax, fname); 121 123 } … … 203 205 204 206 for(G4int i=0; i<nm; i++) { 205 if(modelName == modelList[i]->GetName() && 207 G4String mname = ""; 208 if(modelList[i]) mname = modelList[i]->GetName(); 209 G4String fname = ""; 210 if(flucModelList[i]) fname = flucModelList[i]->GetName(); 211 if(modelName == mname && flucModelName == fname && 206 212 (particleList[i] == "" || particleList[i] == particleName) ) { 207 213 mod = modelList[i]; … … 214 220 215 221 if(!mod) { 216 G4cout << "### G4EmConfigurator WARNING: fails to find a model <" 217 << modelName << "> for process <" 218 << processName << "> and " << particleName 219 << G4endl; 220 if(flucModelName != "") 221 G4cout << " fluctuation model <" 222 << flucModelName << G4endl; 222 223 // set fluctuation model for ionisation processes 224 if(fluc && ptype == eloss) { 225 G4VEnergyLossProcess* p = reinterpret_cast<G4VEnergyLossProcess*>(proc); 226 p->SetFluctModel(fluc); 227 228 } else { 229 G4cout << "### G4EmConfigurator WARNING: fails to find a model <" 230 << modelName << "> for process <" 231 << processName << "> and " << particleName 232 << G4endl; 233 if(flucModelName != "") { 234 G4cout << " fluctuation model <" 235 << flucModelName << G4endl; 236 } 237 } 223 238 } else { 224 239 -
trunk/source/processes/electromagnetic/utils/src/G4EmElementSelector.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4EmElementSelector.cc,v 1. 4 2008/08/21 18:53:32vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4EmElementSelector.cc,v 1.10 2009/05/26 16:59:35 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 58 58 G4double emin, 59 59 G4double emax, 60 G4bool spline):60 G4bool /*spline*/): 61 61 model(mod), material(mat), nbins(bins), cutEnergy(-1.0), 62 62 lowEnergy(emin), highEnergy(emax) … … 65 65 nElmMinusOne = n - 1; 66 66 theElementVector = material->GetElementVector(); 67 element = (*theElementVector)[0]; 67 68 if(nElmMinusOne > 0) { 68 for(G4int i=0; i<nElmMinusOne; i++) { 69 xSections.reserve(n); 70 for(G4int i=0; i<n; i++) { 69 71 G4PhysicsLogVector* v = new G4PhysicsLogVector(lowEnergy,highEnergy,nbins); 70 v->SetSpline(spline);72 //v->SetSpline(spline); 71 73 xSections.push_back(v); 72 74 } … … 80 82 { 81 83 if(nElmMinusOne > 0) { 82 for(G4int i=0; i< nElmMinusOne; i++) {84 for(G4int i=0; i<=nElmMinusOne; i++) { 83 85 delete xSections[i]; 84 86 } … … 101 103 102 104 G4int i; 103 G4int n = nElmMinusOne + 1;104 G4double* xsec = new G4double[n];105 105 106 106 // loop over bins … … 110 110 cross = 0.0; 111 111 //G4cout << "j= " << j << " e(MeV)= " << e/MeV << G4endl; 112 for (i=0; i< n; i++) {112 for (i=0; i<=nElmMinusOne; i++) { 113 113 cross += theAtomNumDensityVector[i]* 114 114 model->ComputeCrossSectionPerAtom(part, (*theElementVector)[i], e, 115 115 cutEnergy, e); 116 xsec[i] = cross; 117 } 118 if(DBL_MIN >= cross) cross = 1.0; 119 // normalise cross section sum 120 for (i=0; i<nElmMinusOne; i++) { 121 xSections[i]->PutValue(j, xsec[i]/cross); 122 //G4cout << "i= " << i << " xs= " << xsec[i]/cross << G4endl; 116 xSections[i]->PutValue(j, cross); 123 117 } 124 118 } 125 delete [] xsec; 119 120 // xSections start from null, so use probabilities from the next bin 121 if(DBL_MIN >= (*xSections[nElmMinusOne])[0]) { 122 for (i=0; i<=nElmMinusOne; i++) { 123 xSections[i]->PutValue(0, (*xSections[i])[1]); 124 } 125 } 126 // xSections ends with null, so use probabilities from the previous bin 127 if(DBL_MIN >= (*xSections[nElmMinusOne])[nbins-1]) { 128 for (i=0; i<=nElmMinusOne; i++) { 129 xSections[i]->PutValue(nbins-1, (*xSections[i])[nbins-2]); 130 } 131 } 132 // perform normalization 133 for(G4int j=0; j<nbins; j++) { 134 cross = (*xSections[nElmMinusOne])[j]; 135 // only for positive X-section 136 if(cross > DBL_MIN) { 137 for (i=0; i<nElmMinusOne; i++) { 138 xSections[i]->PutValue(j, (*xSections[i])[j]/cross); 139 } 140 } 141 } 126 142 } 127 143 … … 139 155 } 140 156 } 141 G4cout << "Last Element in element vector "157 G4cout << "Last Element in element vector " 142 158 << (*theElementVector)[nElmMinusOne]->GetName() 143 159 << G4endl; -
trunk/source/processes/electromagnetic/utils/src/G4EmModelManager.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4EmModelManager.cc,v 1.4 6 2008/10/13 14:56:56vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4EmModelManager.cc,v 1.49 2009/04/17 10:35:32 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 61 61 // 15-03-07 Add maxCutInRange (V.Ivanchenko) 62 62 // 12-04-07 Add verbosity at destruction (V.Ivanchenko) 63 // 08-04-08 Fixed and simplified initialisation of G4RegionModel (VI) 63 64 // 64 65 // Class Description: … … 134 135 theGamma = G4Gamma::Gamma(); 135 136 thePositron = G4Positron::Positron(); 137 models.reserve(4); 138 flucModels.reserve(4); 139 regions.reserve(4); 140 orderOfModels.reserve(4); 141 isUsed.reserve(4); 136 142 } 137 143 … … 178 184 regions.push_back(r); 179 185 orderOfModels.push_back(num); 186 isUsed.push_back(0); 180 187 p->DefineForRegion(r); 181 188 nEmModels++; … … 187 194 G4double emin, G4double emax) 188 195 { 189 if (nEmModels ) {196 if (nEmModels > 0) { 190 197 for(G4int i=0; i<nEmModels; i++) { 191 198 if(nam == models[i]->GetName()) { … … 225 232 { 226 233 verboseLevel = val; 234 G4String partname = p->GetParticleName(); 227 235 if(1 < verboseLevel) { 228 236 G4cout << "G4EmModelManager::Initialise() for " 229 << p->GetParticleName() 230 << G4endl; 237 << partname << G4endl; 231 238 } 232 239 // Are models defined? 233 240 if(!nEmModels) { 234 G4Exception("G4EmModelManager::Initialise without any model defined ");241 G4Exception("G4EmModelManager::Initialise without any model defined for "+partname); 235 242 } 236 243 particle = p; … … 271 278 G4ProductionCutsTable::GetProductionCutsTable(); 272 279 G4int numOfCouples = theCoupleTable->GetTableSize(); 273 idxOfRegionModels = new G4int[numOfCouples+1]; 274 idxOfRegionModels[numOfCouples] = 0; 280 if(nRegions > 1) idxOfRegionModels = new G4int[numOfCouples]; 275 281 setOfRegionModels = new G4RegionModels*[nRegions]; 276 282 277 283 std::vector<G4int> modelAtRegion(nEmModels); 278 284 std::vector<G4int> modelOrd(nEmModels); 279 G4DataVector eLow(nEmModels );285 G4DataVector eLow(nEmModels+1); 280 286 G4DataVector eHigh(nEmModels); 281 G4int nmax = nEmModels;282 287 283 288 // Order models for regions … … 305 310 if (region) G4cout << region->GetName(); 306 311 G4cout << "> " 307 308 309 << "; order= " << ord310 312 << " tmin(MeV)= " << tmin/MeV 313 << "; tmax(MeV)= " << tmax/MeV 314 << "; order= " << ord 315 << G4endl; 311 316 } 312 317 313 if (n == 0) n++; 314 else { 318 if(n > 0) { 319 320 // extend energy range to previous models 315 321 tmin = std::min(tmin, eHigh[n-1]); 316 if(tmin >= tmax) push = false; 317 else { 318 319 // high energy model 320 if(tmin == eHigh[n-1] && tmax > eHigh[n-1]) n++; 321 else if (tmax > eHigh[n-1]) { 322 323 // compare order of models 324 for(G4int k = n-1; k>=0; k--) { 325 326 if (ord >= modelOrd[k]) { 327 tmin = std::max(tmin, eHigh[k]); 328 if(k < n-1) n = k + 2; 329 break; 330 } else if (tmin > eLow[k]) { 331 eHigh[k] = tmin; 332 n = k + 2; 333 break; 334 } else if (tmin == eLow[k]) { 335 n = k + 1; 336 break; 322 tmax = std::max(tmax, eLow[0]); 323 //G4cout << "tmin= " << tmin << " tmax= " 324 // << tmax << " ord= " << ord <<G4endl; 325 // empty energy range 326 if( tmax - tmin <= eV) push = false; 327 // low-energy model 328 else if (tmax == eLow[0]) { 329 push = false; 330 insert = true; 331 idx = 0; 332 // resolve intersections 333 } else if(tmin < eHigh[n-1]) { 334 // compare order 335 for(G4int k=0; k<n; k++) { 336 // new model has lower application 337 if(ord >= modelOrd[k]) { 338 if(tmin < eHigh[k] && tmin >= eLow[k]) tmin = eHigh[k]; 339 if(tmax <= eHigh[k] && tmax > eLow[k]) tmax = eLow[k]; 340 if(tmax > eHigh[k] && tmin < eLow[k]) { 341 if(tmax - eHigh[k] > eLow[k] - tmin) tmin = eHigh[k]; 342 else tmax = eLow[k]; 343 } 344 if( tmax - tmin <= eV) { 345 push = false; 346 break; 337 347 } 338 348 } 339 if(tmin < eLow[0]) n = 1; 340 idx = n - 1; 341 342 // low energy model 343 } else { 344 345 tmax = std::max(tmax, eLow[0]); 346 insert = true; 347 push = false; 348 idx = 0; 349 if(tmax <= eLow[0]) tmax = eLow[0]; 350 else { 351 352 for(G4int k=0; k<n; k++) { 353 354 if (ord >= modelOrd[k]) { 355 if(k == 0) { 356 if(tmin < eLow[0]) tmax = eLow[0]; 357 else insert = false; 358 break; 359 } else { 360 insert = false; 361 break; 362 } 363 } else if(tmax < eHigh[k]) { 364 idx = k; 365 if(k > 0) tmin = eLow[k]; 366 eLow[k] = tmax; 367 break; 368 } else if(tmax == eHigh[k]) { 369 insert = false; 370 push = true; 371 idx = k; 372 if(k > 0) tmin = eLow[k]; 373 else tmin = std::min(tmin,eLow[0]); 374 break; 375 } else { 376 modelAtRegion[k] = ii; 377 modelOrd[k] = ord; 378 if(k == 0) eLow[idx] = std::min(tmin,eLow[0]); 349 } 350 //G4cout << "tmin= " << tmin << " tmax= " 351 // << tmax << " push= " << push << " idx= " << idx <<G4endl; 352 if(push) { 353 if (tmax == eLow[0]) { 354 push = false; 355 insert = true; 356 idx = 0; 357 // continue resolve intersections 358 } else if(tmin < eHigh[n-1]) { 359 // last energy interval 360 if(tmin > eLow[n-1] && tmax >= eHigh[n-1]) { 361 eHigh[n-1] = tmin; 362 // first energy interval 363 } else if(tmin <= eLow[0] && tmax < eHigh[0]) { 364 eLow[0] = tmax; 365 push = false; 366 insert = true; 367 idx = 0; 368 } else { 369 // find energy interval to replace 370 for(G4int k=0; k<n; k++) { 371 if(tmin <= eLow[k] && tmax >= eHigh[k]) { 372 push = false; 373 modelAtRegion[k] = ii; 374 modelOrd[k] = ord; 375 isUsed[ii] = 1; 376 } 379 377 } 380 378 } 381 379 } 382 if(insert && idx < n) n++;383 else insert = false;384 380 } 385 381 } 386 382 } 387 if(n > nmax) {388 nmax = n;389 modelAtRegion.resize(nmax);390 modelOrd.resize(nmax);391 eLow.resize(nmax);392 eHigh.resize(nmax);393 }394 383 if(insert) { 395 for(G4int k=n- 2; k>=idx; k--) {384 for(G4int k=n-1; k>=idx; k--) { 396 385 modelAtRegion[k+1] = modelAtRegion[k]; 397 386 modelOrd[k+1] = modelOrd[k]; … … 400 389 } 401 390 } 391 //G4cout << "push= " << push << " insert= " << insert 392 //<< " idx= " << idx <<G4endl; 402 393 if (push || insert) { 394 n++; 403 395 modelAtRegion[idx] = ii; 404 396 modelOrd[idx] = ord; 405 397 eLow[idx] = tmin; 406 398 eHigh[idx] = tmax; 399 isUsed[ii] = 1; 407 400 } 408 401 } … … 416 409 } 417 410 eLow[0] = 0.0; 418 if(n >= nmax) eLow.resize(nmax+1);419 411 eLow[n] = eHigh[n-1]; 420 412 … … 430 422 } 431 423 424 currRegionModel = setOfRegionModels[0]; 425 432 426 // Access to materials and build cuts 433 434 427 for(G4int i=0; i<numOfCouples; i++) { 435 428 … … 439 432 const G4ProductionCuts* pcuts = couple->GetProductionCuts(); 440 433 441 G4int reg = nRegions; 442 do {reg--;} while (reg>0 && pcuts != (setr[reg]->GetProductionCuts())); 443 idxOfRegionModels[i] = reg; 444 434 G4int reg = 0; 435 if(nRegions > 1) { 436 reg = nRegions; 437 do {reg--;} while (reg>0 && pcuts != (setr[reg]->GetProductionCuts())); 438 idxOfRegionModels[i] = reg; 439 } 445 440 if(1 < verboseLevel) { 446 441 G4cout << "G4EmModelManager::Initialise() for " 447 448 449 450 442 << material->GetName() 443 << " indexOfCouple= " << i 444 << " indexOfRegion= " << reg 445 << G4endl; 451 446 } 452 447 … … 494 489 495 490 for(G4int jj=0; jj<nEmModels; jj++) { 496 models[jj]->Initialise(particle, theCuts); 497 if(flucModels[jj]) flucModels[jj]->InitialiseMe(particle); 498 } 499 491 if(1 == isUsed[jj]) { 492 models[jj]->Initialise(particle, theCuts); 493 if(flucModels[jj]) flucModels[jj]->InitialiseMe(particle); 494 } 495 } 500 496 501 497 if(1 < verboseLevel) { … … 534 530 } 535 531 536 G4int reg = idxOfRegionModels[i]; 532 G4int reg = 0; 533 if(nRegions > 1) reg = idxOfRegionModels[i]; 537 534 const G4RegionModels* regModels = setOfRegionModels[reg]; 538 535 G4int nmod = regModels->NumberOfModels(); … … 683 680 } 684 681 685 G4int reg = idxOfRegionModels[i]; 682 G4int reg = 0; 683 if(nRegions > 1) reg = idxOfRegionModels[i]; 686 684 const G4RegionModels* regModels = setOfRegionModels[reg]; 687 685 G4int nmod = regModels->NumberOfModels(); … … 799 797 G4RegionModels* r = setOfRegionModels[i]; 800 798 const G4Region* reg = r->Region(); 801 if(verb >1 || nRegions > 1) {802 }799 // if(verb > -1 || nRegions > 1) { 800 // } 803 801 G4int n = r->NumberOfModels(); 804 if(verb > 1 || n > 0) {802 if(verb > -1 || n > 0) { 805 803 G4cout << " ===== EM models for the G4Region " << reg->GetName() 806 804 << " ======" << G4endl;; -
trunk/source/processes/electromagnetic/utils/src/G4EmProcessOptions.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4EmProcessOptions.cc,v 1.2 4 2008/04/17 10:33:27 vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4EmProcessOptions.cc,v 1.26 2009/02/18 14:43:27 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 59 59 #include "G4VMultipleScattering.hh" 60 60 #include "G4Region.hh" 61 #include "G4RegionStore.hh" 61 62 62 63 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... … … 392 393 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 393 394 394 void G4EmProcessOptions::ActivateDeexcitation(G4bool val, const G4Region* r) 395 { 396 const std::vector<G4VEnergyLossProcess*>& v = 397 theManager->GetEnergyLossProcessVector(); 398 std::vector<G4VEnergyLossProcess*>::const_iterator itr; 399 for(itr = v.begin(); itr != v.end(); itr++) { 400 G4VEnergyLossProcess* p = *itr; 401 if(p) p->ActivateDeexcitation(val,r); 402 } 403 const std::vector<G4VEmProcess*>& w = 404 theManager->GetEmProcessVector(); 405 std::vector<G4VEmProcess*>::const_iterator itp; 406 for(itp = w.begin(); itp != w.end(); itp++) { 407 G4VEmProcess* q = *itp; 408 if(q) q->ActivateDeexcitation(val,r); 395 void G4EmProcessOptions::ActivateDeexcitation(const G4String& pname, 396 G4bool val, 397 const G4String& reg) 398 { 399 G4RegionStore* regionStore = G4RegionStore::GetInstance(); 400 const G4Region* r = 0; 401 if(reg == "" || reg == "World") { 402 r = regionStore->GetRegion("DefaultRegionForTheWorld", false); 403 } else { 404 r = regionStore->GetRegion(reg, false); 405 } 406 if(!r) { 407 G4cout << "G4EmProcessOptions::ActivateDeexcitation ERROR: G4Region <" 408 << reg << "> not found, the command ignored" << G4endl; 409 return; 410 } 411 412 const std::vector<G4VEnergyLossProcess*>& v = 413 theManager->GetEnergyLossProcessVector(); 414 std::vector<G4VEnergyLossProcess*>::const_iterator itr; 415 for(itr = v.begin(); itr != v.end(); itr++) { 416 G4VEnergyLossProcess* p = *itr; 417 if(p) { 418 if(pname == p->GetProcessName()) p->ActivateDeexcitation(val,r); 419 } 420 } 421 const std::vector<G4VEmProcess*>& w = 422 theManager->GetEmProcessVector(); 423 std::vector<G4VEmProcess*>::const_iterator itp; 424 for(itp = w.begin(); itp != w.end(); itp++) { 425 G4VEmProcess* q = *itp; 426 if(q) { 427 if(pname == q->GetProcessName()) q->ActivateDeexcitation(val,r); 428 } 409 429 } 410 430 } -
trunk/source/processes/electromagnetic/utils/src/G4EnergyLossMessenger.cc
r1007 r1055 25 25 // 26 26 // 27 // $Id: G4EnergyLossMessenger.cc,v 1.3 5 2008/10/20 13:27:45vnivanch Exp $28 // GEANT4 tag $Name: geant4-09-0 2$27 // $Id: G4EnergyLossMessenger.cc,v 1.37 2009/02/18 14:43:27 vnivanch Exp $ 28 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 29 29 // 30 30 // ------------------------------------------------------------------- … … 178 178 aplCmd->AvailableForStates(G4State_PreInit,G4State_Idle); 179 179 180 deexCmd = new G4UIcommand("/process/em/deexcitation",this); 181 deexCmd->SetGuidance("Set deexcitation flag per process and G4Region."); 182 deexCmd->SetGuidance(" procName : process name"); 183 deexCmd->SetGuidance(" flag : flag"); 184 deexCmd->SetGuidance(" regName : G4Region name"); 185 186 G4UIparameter* pName = new G4UIparameter("pName",'s',false); 187 deexCmd->SetParameter(pName); 188 189 G4UIparameter* flag = new G4UIparameter("flag",'s',false); 190 deexCmd->SetParameter(flag); 191 192 G4UIparameter* regName = new G4UIparameter("regName",'s',false); 193 deexCmd->SetParameter(regName); 194 180 195 dedxCmd = new G4UIcmdWithAnInteger("/process/eLoss/binsDEDX",this); 181 196 dedxCmd->SetGuidance("Set number of bins for DEDX tables"); … … 259 274 delete MinSubSecCmd; 260 275 delete StepFuncCmd; 276 delete deexCmd; 261 277 delete eLossDirectory; 262 278 delete mscDirectory; … … 320 336 } 321 337 338 if (command == deexCmd) { 339 G4String s1 (""), s2(""), s3(""); 340 G4bool b = false; 341 std::istringstream is(newValue); 342 is >> s1 >> s2 >> s3; 343 if(s2 == "true") b = true; 344 opt->ActivateDeexcitation(s1,b,s3); 345 } 346 322 347 if (command == mscCmd) { 323 348 if(newValue == "Minimal") -
trunk/source/processes/electromagnetic/utils/src/G4LossTableBuilder.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4LossTableBuilder.cc,v 1.2 7 2008/07/22 15:55:15vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4LossTableBuilder.cc,v 1.28 2009/02/18 16:24:47 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 132 132 G4double dedx1 = pv->GetValue(elow, b); 133 133 134 //G4cout << "nbins= " << nbins << " dedx1= " << dedx1 << G4endl; 135 134 136 // protection for specific cases dedx=0 135 137 if(dedx1 == 0.0) { 136 for (size_t k=1; k<nbins ; k++) {138 for (size_t k=1; k<nbins-1; k++) { 137 139 bin0++; 138 140 elow = pv->GetLowEdgeEnergy(k); … … 143 145 } 144 146 147 //G4cout << "nbins= " << nbins << " elow= " << elow << " ehigh= " << ehigh << G4endl; 145 148 // initialisation of a new vector 146 149 G4PhysicsLogVector* v = new G4PhysicsLogVector(elow, ehigh, nbins); 150 // dedx is exect zero 151 if(nbins <= 2) { 152 v->PutValue(0,1000.); 153 v->PutValue(1,2000.); 154 G4PhysicsTableHelper::SetPhysicsVector(rangeTable, i, v); 155 return; 156 } 147 157 v->SetSpline(splineFlag); 148 158 -
trunk/source/processes/electromagnetic/utils/src/G4LossTableManager.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4LossTableManager.cc,v 1.9 5 2008/11/13 18:23:39 schaelicExp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4LossTableManager.cc,v 1.96 2009/04/09 16:10:57 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 93 93 #include "G4EmCorrections.hh" 94 94 #include "G4EmSaturation.hh" 95 #include "G4EmConfigurator.hh" 95 96 #include "G4EmTableType.hh" 96 97 #include "G4LossTableBuilder.hh" … … 164 165 emCorrections= new G4EmCorrections(); 165 166 emSaturation = new G4EmSaturation(); 167 emConfigurator = new G4EmConfigurator(); 166 168 integral = true; 167 169 integralActive = false; … … 932 934 } 933 935 936 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 937 938 G4EmConfigurator* G4LossTableManager::EmConfigurator() 939 { 940 return emConfigurator; 941 } 942 934 943 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... -
trunk/source/processes/electromagnetic/utils/src/G4VEmFluctuationModel.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4VEmFluctuationModel.cc,v 1. 3 2008/07/15 16:56:39vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4VEmFluctuationModel.cc,v 1.4 2009/02/19 11:25:50 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 65 65 } 66 66 67 void G4VEmFluctuationModel::InitialiseMe(const G4ParticleDefinition*) 68 {} 67 69 70 void G4VEmFluctuationModel::SetParticleAndCharge(const G4ParticleDefinition*, 71 G4double) 72 {} 73 74 -
trunk/source/processes/electromagnetic/utils/src/G4VEmModel.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4VEmModel.cc,v 1.2 0 2008/11/13 23:13:18 schaelicExp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4VEmModel.cc,v 1.27 2009/05/26 15:00:49 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 41 41 // 25.10.2005 Set default highLimit=100.TeV (V.Ivanchenko) 42 42 // 06.02.2006 add method ComputeMeanFreePath() (mma) 43 // 16.02.2009 Move implementations of virtual methods to source (VI) 43 44 // 44 45 // … … 53 54 #include "G4LossTableManager.hh" 54 55 #include "G4ProductionCutsTable.hh" 56 #include "G4ParticleChangeForLoss.hh" 57 #include "G4ParticleChangeForGamma.hh" 55 58 56 59 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... … … 60 63 fluc(0), name(nam), lowLimit(0.1*keV), highLimit(100.0*TeV), 61 64 polarAngleLimit(0.0),secondaryThreshold(DBL_MAX),theLPMflag(false), 62 pParticleChange(0),nuclearStopping(false),nsec(5) 65 pParticleChange(0),nuclearStopping(false), 66 currentCouple(0),currentElement(0), 67 nsec(5),flagDeexcitation(false) 63 68 { 64 69 xsec.resize(nsec); … … 78 83 } 79 84 } 85 } 86 87 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 88 89 G4ParticleChangeForLoss* G4VEmModel::GetParticleChangeForLoss() 90 { 91 G4ParticleChangeForLoss* p = 0; 92 if (pParticleChange) { 93 p = static_cast<G4ParticleChangeForLoss*>(pParticleChange); 94 } else { 95 p = new G4ParticleChangeForLoss(); 96 } 97 return p; 98 } 99 100 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 101 102 G4ParticleChangeForGamma* G4VEmModel::GetParticleChangeForGamma() 103 { 104 G4ParticleChangeForGamma* p = 0; 105 if (pParticleChange) { 106 p = static_cast<G4ParticleChangeForGamma*>(pParticleChange); 107 } else { 108 p = new G4ParticleChangeForGamma(); 109 } 110 return p; 111 } 112 113 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 114 115 void G4VEmModel::InitialiseElementSelectors(const G4ParticleDefinition* p, 116 const G4DataVector& cuts) 117 { 118 // initialise before run 119 flagDeexcitation = false; 120 121 G4int nbins = G4int(std::log10(highLimit/lowLimit) + 0.5); 122 if(nbins < 3) nbins = 3; 123 G4bool spline = G4LossTableManager::Instance()->SplineFlag(); 124 125 G4ProductionCutsTable* theCoupleTable= 126 G4ProductionCutsTable::GetProductionCutsTable(); 127 G4int numOfCouples = theCoupleTable->GetTableSize(); 128 129 // prepare vector 130 if(numOfCouples > nSelectors) elmSelectors.reserve(numOfCouples); 131 132 // initialise vector 133 for(G4int i=0; i<numOfCouples; i++) { 134 const G4MaterialCutsCouple* couple = 135 theCoupleTable->GetMaterialCutsCouple(i); 136 const G4Material* material = couple->GetMaterial(); 137 G4int idx = couple->GetIndex(); 138 139 // selector already exist check if should be deleted 140 G4bool create = true; 141 if(i < nSelectors) { 142 if(elmSelectors[i]) { 143 if(material == elmSelectors[i]->GetMaterial()) create = false; 144 else delete elmSelectors[i]; 145 } 146 } else { 147 nSelectors++; 148 elmSelectors.push_back(0); 149 } 150 if(create) { 151 elmSelectors[i] = new G4EmElementSelector(this,material,nbins, 152 lowLimit,highLimit,spline); 153 } 154 elmSelectors[i]->Initialise(p, cuts[idx]); 155 //elmSelectors[i]->Dump(p); 156 } 157 } 158 159 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 160 161 G4double G4VEmModel::ComputeDEDXPerVolume(const G4Material*, 162 const G4ParticleDefinition*, 163 G4double,G4double) 164 { 165 return 0.0; 80 166 } 81 167 … … 107 193 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 108 194 109 G4double G4VEmModel::ComputeMeanFreePath(const G4ParticleDefinition* p, 110 G4double ekin, 111 const G4Material* material, 112 G4double emin, 113 G4double emax) 114 { 115 G4double mfp = DBL_MAX; 116 G4double cross = CrossSectionPerVolume(material,p,ekin,emin,emax); 117 if (cross > DBL_MIN) mfp = 1./cross; 118 return mfp; 119 } 120 121 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 122 123 void G4VEmModel::InitialiseElementSelectors(const G4ParticleDefinition* p, 124 const G4DataVector& cuts) 125 { 126 G4int nbins = G4int(std::log10(highLimit/lowLimit) + 0.5); 127 if(nbins < 3) nbins = 3; 128 G4bool spline = G4LossTableManager::Instance()->SplineFlag(); 129 130 G4ProductionCutsTable* theCoupleTable= 131 G4ProductionCutsTable::GetProductionCutsTable(); 132 G4int numOfCouples = theCoupleTable->GetTableSize(); 133 134 // prepare vector 135 if(numOfCouples > nSelectors) { 136 elmSelectors.resize(numOfCouples); 137 nSelectors = numOfCouples; 138 } 139 140 // initialise vector 141 for(G4int i=0; i<numOfCouples; i++) { 142 const G4MaterialCutsCouple* couple = 143 theCoupleTable->GetMaterialCutsCouple(i); 144 const G4Material* material = couple->GetMaterial(); 145 G4int idx = couple->GetIndex(); 146 147 // selector already exist check if should be deleted 148 G4bool create = true; 149 if(elmSelectors[i]) { 150 if(material == elmSelectors[i]->GetMaterial()) create = false; 151 else delete elmSelectors[i]; 152 } 153 if(create) { 154 elmSelectors[i] = new G4EmElementSelector(this,material,nbins, 155 lowLimit,highLimit,spline); 156 } 157 elmSelectors[i]->Initialise(p, cuts[idx]); 158 //elmSelectors[i]->Dump(p); 159 } 160 } 161 162 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 163 164 195 G4double G4VEmModel::ComputeCrossSectionPerAtom(const G4ParticleDefinition*, 196 G4double, G4double, G4double, 197 G4double, G4double) 198 { 199 return 0.0; 200 } 201 202 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 203 204 G4double G4VEmModel::MinEnergyCut(const G4ParticleDefinition*, 205 const G4MaterialCutsCouple*) 206 { 207 return 0.0; 208 } 209 210 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 211 212 G4double G4VEmModel::GetChargeSquareRatio(const G4ParticleDefinition* p, 213 const G4Material*, G4double) 214 { 215 G4double q = p->GetPDGCharge()/CLHEP::eplus; 216 return q*q; 217 } 218 219 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 220 221 G4double G4VEmModel::GetParticleCharge(const G4ParticleDefinition* p, 222 const G4Material*, G4double) 223 { 224 return p->GetPDGCharge(); 225 } 226 227 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 228 229 void G4VEmModel::CorrectionsAlongStep(const G4MaterialCutsCouple*, 230 const G4DynamicParticle*, 231 G4double&,G4double&,G4double) 232 {} 233 234 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 235 236 void G4VEmModel::SampleDeexcitationAlongStep(const G4Material*, 237 const G4Track&, 238 G4double& ) 239 {} 240 241 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 242 243 G4double G4VEmModel::MaxSecondaryEnergy(const G4ParticleDefinition*, 244 G4double kineticEnergy) 245 { 246 return kineticEnergy; 247 } 248 249 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 250 251 void G4VEmModel::DefineForRegion(const G4Region*) 252 {} 253 254 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 255 256 void G4VEmModel::SetupForMaterial(const G4ParticleDefinition*, 257 const G4Material*, G4double) 258 {} 259 260 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... -
trunk/source/processes/electromagnetic/utils/src/G4VEmProcess.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4VEmProcess.cc,v 1.6 0 2008/10/17 14:46:16vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4VEmProcess.cc,v 1.66 2009/04/17 10:35:32 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 78 78 #include "G4Positron.hh" 79 79 #include "G4PhysicsTableHelper.hh" 80 #include "G4EmConfigurator.hh" 80 81 81 82 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... … … 91 92 applyCuts(false), 92 93 startFromNull(true), 93 nRegions(0), 94 selectedModel(0), 94 useDeexcitation(false), 95 nDERegions(0), 96 idxDERegions(0), 97 currentModel(0), 95 98 particle(0), 96 99 currentCouple(0) … … 100 103 // Size of tables assuming spline 101 104 minKinEnergy = 0.1*keV; 102 maxKinEnergy = 10 0.0*TeV;103 nLambdaBins = 84;105 maxKinEnergy = 10.0*TeV; 106 nLambdaBins = 77; 104 107 105 108 // default lambda factor … … 139 142 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 140 143 144 void G4VEmProcess::Clear() 145 { 146 delete [] theEnergyOfCrossSectionMax; 147 delete [] theCrossSectionMax; 148 delete [] idxDERegions; 149 theEnergyOfCrossSectionMax = 0; 150 theCrossSectionMax = 0; 151 idxDERegions = 0; 152 currentCouple = 0; 153 preStepLambda = 0.0; 154 mfpKinEnergy = DBL_MAX; 155 deRegions.clear(); 156 nDERegions = 0; 157 } 158 159 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 160 161 void G4VEmProcess::AddEmModel(G4int order, G4VEmModel* p, 162 const G4Region* region) 163 { 164 G4VEmFluctuationModel* fm = 0; 165 modelManager->AddEmModel(order, p, fm, region); 166 if(p) p->SetParticleChange(pParticleChange); 167 } 168 169 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 170 171 void G4VEmProcess::SetModel(G4VEmModel* p, G4int index) 172 { 173 G4int n = emModels.size(); 174 if(index >= n) for(G4int i=n; i<index+1; i++) {emModels.push_back(0);} 175 emModels[index] = p; 176 } 177 178 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 179 180 G4VEmModel* G4VEmProcess::Model(G4int index) 181 { 182 G4VEmModel* p = 0; 183 if(index >= 0 && index < G4int(emModels.size())) p = emModels[index]; 184 return p; 185 } 186 187 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 188 189 void G4VEmProcess::UpdateEmModel(const G4String& nam, 190 G4double emin, G4double emax) 191 { 192 modelManager->UpdateEmModel(nam, emin, emax); 193 } 194 195 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 196 197 G4VEmModel* G4VEmProcess::GetModelByIndex(G4int idx, G4bool ver) 198 { 199 return modelManager->GetModel(idx, ver); 200 } 201 202 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 203 141 204 void G4VEmProcess::PreparePhysicsTable(const G4ParticleDefinition& part) 142 205 { … … 149 212 << G4endl; 150 213 } 214 215 (G4LossTableManager::Instance())->EmConfigurator()->AddModels(); 151 216 152 217 if(particle == &part) { … … 159 224 theCutsElectron = theCoupleTable->GetEnergyCutsVector(idxG4ElectronCut); 160 225 theCutsPositron = theCoupleTable->GetEnergyCutsVector(idxG4PositronCut); 161 if(buildLambdaTable) 226 if(buildLambdaTable){ 162 227 theLambdaTable = G4PhysicsTableHelper::PreparePhysicsTable(theLambdaTable); 163 } 164 } 165 166 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 167 168 void G4VEmProcess::Clear() 169 { 170 if(theEnergyOfCrossSectionMax) delete [] theEnergyOfCrossSectionMax; 171 if(theCrossSectionMax) delete [] theCrossSectionMax; 172 theEnergyOfCrossSectionMax = 0; 173 theCrossSectionMax = 0; 174 currentCouple = 0; 175 preStepLambda = 0.0; 176 mfpKinEnergy = DBL_MAX; 228 } 229 } 230 // Sub Cutoff and Deexcitation 231 if (nDERegions>0) { 232 233 const G4ProductionCutsTable* theCoupleTable= 234 G4ProductionCutsTable::GetProductionCutsTable(); 235 size_t numOfCouples = theCoupleTable->GetTableSize(); 236 237 idxDERegions = new G4bool[numOfCouples]; 238 239 for (size_t j=0; j<numOfCouples; j++) { 240 241 const G4MaterialCutsCouple* couple = 242 theCoupleTable->GetMaterialCutsCouple(j); 243 const G4ProductionCuts* pcuts = couple->GetProductionCuts(); 244 G4bool reg = false; 245 for(G4int i=0; i<nDERegions; i++) { 246 if(deRegions[i]) { 247 if(pcuts == deRegions[i]->GetProductionCuts()) reg = true; 248 } 249 } 250 idxDERegions[j] = reg; 251 } 252 } 253 if (1 < verboseLevel && nDERegions>0) { 254 G4cout << " Deexcitation is activated for regions: " << G4endl; 255 for (G4int i=0; i<nDERegions; i++) { 256 const G4Region* r = deRegions[i]; 257 G4cout << " " << r->GetName() << G4endl; 258 } 259 } 177 260 } 178 261 … … 237 320 G4cout << *theLambdaTable << G4endl; 238 321 } 322 } 323 } 324 325 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 326 327 void G4VEmProcess::PrintInfoDefinition() 328 { 329 if(verboseLevel > 0) { 330 G4cout << G4endl << GetProcessName() << ": for " 331 << particle->GetParticleName(); 332 if(integral) G4cout << ", integral: 1 "; 333 if(applyCuts) G4cout << ", applyCuts: 1 "; 334 G4cout << " SubType= " << GetProcessSubType() << G4endl; 335 if(buildLambdaTable) { 336 G4cout << " Lambda tables from " 337 << G4BestUnit(minKinEnergy,"Energy") 338 << " to " 339 << G4BestUnit(maxKinEnergy,"Energy") 340 << " in " << nLambdaBins << " bins, spline: " 341 << (G4LossTableManager::Instance())->SplineFlag() 342 << G4endl; 343 } 344 PrintInfo(); 345 modelManager->DumpModelList(verboseLevel); 346 } 347 348 if(verboseLevel > 2 && buildLambdaTable) { 349 G4cout << " LambdaTable address= " << theLambdaTable << G4endl; 350 if(theLambdaTable) G4cout << (*theLambdaTable) << G4endl; 239 351 } 240 352 } … … 304 416 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 305 417 306 G4double G4VEmProcess::GetMeanFreePath(const G4Track& track,307 G4double,308 G4ForceCondition* condition)309 {310 *condition = NotForced;311 return G4VEmProcess::MeanFreePath(track);312 }313 314 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....315 316 418 G4VParticleChange* G4VEmProcess::PostStepDoIt(const G4Track& track, 317 419 const G4Step&) … … 342 444 } 343 445 344 G4VEmModel* currentModel = SelectModel(finalT); 345 446 SelectModel(finalT); 447 if(useDeexcitation) { 448 currentModel->SetDeexcitationFlag(idxDERegions[currentMaterialIndex]); 449 } 346 450 /* 347 451 if(0 < verboseLevel) { … … 404 508 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 405 509 406 void G4VEmProcess::PrintInfoDefinition()407 {408 if(verboseLevel > 0) {409 G4cout << G4endl << GetProcessName() << ": for "410 << particle->GetParticleName();411 if(integral) G4cout << ", integral: 1 ";412 if(applyCuts) G4cout << ", applyCuts: 1 ";413 G4cout << " SubType= " << GetProcessSubType() << G4endl;414 if(buildLambdaTable) {415 G4cout << " Lambda tables from "416 << G4BestUnit(minKinEnergy,"Energy")417 << " to "418 << G4BestUnit(maxKinEnergy,"Energy")419 << " in " << nLambdaBins << " bins, spline: "420 << (G4LossTableManager::Instance())->SplineFlag()421 << G4endl;422 }423 PrintInfo();424 modelManager->DumpModelList(verboseLevel);425 }426 427 if(verboseLevel > 2 && buildLambdaTable) {428 G4cout << " LambdaTable address= " << theLambdaTable << G4endl;429 if(theLambdaTable) G4cout << (*theLambdaTable) << G4endl;430 }431 }432 433 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....434 435 G4double G4VEmProcess::CrossSectionPerVolume(G4double kineticEnergy,436 const G4MaterialCutsCouple* couple)437 {438 // Cross section per atom is calculated439 DefineMaterial(couple);440 G4double cross = 0.0;441 G4bool b;442 if(theLambdaTable) {443 cross = (((*theLambdaTable)[currentMaterialIndex])->444 GetValue(kineticEnergy, b));445 } else {446 G4VEmModel* model = SelectModel(kineticEnergy);447 cross =448 model->CrossSectionPerVolume(currentMaterial,particle,kineticEnergy);449 }450 451 return cross;452 }453 454 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....455 456 510 G4bool G4VEmProcess::StorePhysicsTable(const G4ParticleDefinition* part, 457 511 const G4String& directory, … … 521 575 522 576 return yes; 577 } 578 579 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 580 581 void G4VEmProcess::ActivateDeexcitation(G4bool val, const G4Region* r) 582 { 583 G4RegionStore* regionStore = G4RegionStore::GetInstance(); 584 const G4Region* reg = r; 585 if (!reg) {reg = regionStore->GetRegion("DefaultRegionForTheWorld", false);} 586 587 // the region is in the list 588 if (nDERegions) { 589 for (G4int i=0; i<nDERegions; i++) { 590 if (reg == deRegions[i]) { 591 if(!val) deRegions[i] = 0; 592 return; 593 } 594 } 595 } 596 597 // new region 598 if(val) { 599 useDeexcitation = true; 600 deRegions.push_back(reg); 601 nDERegions++; 602 } else { 603 useDeexcitation = false; 604 } 605 } 606 607 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 608 609 G4double G4VEmProcess::CrossSectionPerVolume(G4double kineticEnergy, 610 const G4MaterialCutsCouple* couple) 611 { 612 // Cross section per atom is calculated 613 DefineMaterial(couple); 614 G4double cross = 0.0; 615 G4bool b; 616 if(theLambdaTable) { 617 cross = (((*theLambdaTable)[currentMaterialIndex])-> 618 GetValue(kineticEnergy, b)); 619 } else { 620 SelectModel(kineticEnergy); 621 cross = currentModel->CrossSectionPerVolume(currentMaterial, 622 particle,kineticEnergy); 623 } 624 625 return cross; 626 } 627 628 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 629 630 G4double G4VEmProcess::GetMeanFreePath(const G4Track& track, 631 G4double, 632 G4ForceCondition* condition) 633 { 634 *condition = NotForced; 635 return G4VEmProcess::MeanFreePath(track); 523 636 } 524 637 -
trunk/source/processes/electromagnetic/utils/src/G4VEnergyLossProcess.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4VEnergyLossProcess.cc,v 1.14 3 2008/10/17 14:46:16vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4VEnergyLossProcess.cc,v 1.149 2009/04/17 10:35:32 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 142 142 #include "G4SafetyHelper.hh" 143 143 #include "G4TransportationManager.hh" 144 #include "G4EmConfigurator.hh" 144 145 145 146 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... … … 150 151 secondaryParticle(0), 151 152 nSCoffRegions(0), 153 nDERegions(0), 152 154 idxSCoffRegions(0), 155 idxDERegions(0), 153 156 nProcesses(0), 154 157 theDEDXTable(0), … … 176 179 isIonisation(true), 177 180 useSubCutoff(false), 181 useDeexcitation(false), 178 182 particle(0), 179 183 currentCouple(0), … … 188 192 // Size of tables assuming spline 189 193 minKinEnergy = 0.1*keV; 190 maxKinEnergy = 10 0.0*TeV;191 nBins = 84;194 maxKinEnergy = 10.0*TeV; 195 nBins = 77; 192 196 maxKinEnergyCSDA = 1.0*GeV; 193 197 nBinsCSDA = 35; … … 237 241 << G4endl; 238 242 delete vstrag; 239 Clea r();243 Clean(); 240 244 241 245 if ( !baseParticle ) { … … 291 295 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 292 296 293 void G4VEnergyLossProcess::Clea r()297 void G4VEnergyLossProcess::Clean() 294 298 { 295 299 if(1 < verboseLevel) { … … 302 306 delete [] theCrossSectionMax; 303 307 delete [] idxSCoffRegions; 308 delete [] idxDERegions; 304 309 305 310 theDEDXAtMaxEnergy = 0; … … 312 317 scProcesses.clear(); 313 318 nProcesses = 0; 319 } 320 321 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 322 323 G4double G4VEnergyLossProcess::MinPrimaryEnergy(const G4ParticleDefinition*, 324 const G4Material*, 325 G4double cut) 326 { 327 return cut; 328 } 329 330 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 331 332 void G4VEnergyLossProcess::AddEmModel(G4int order, G4VEmModel* p, 333 G4VEmFluctuationModel* fluc, 334 const G4Region* region) 335 { 336 modelManager->AddEmModel(order, p, fluc, region); 337 if(p) p->SetParticleChange(pParticleChange, fluc); 338 } 339 340 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 341 342 void G4VEnergyLossProcess::UpdateEmModel(const G4String& nam, 343 G4double emin, G4double emax) 344 { 345 modelManager->UpdateEmModel(nam, emin, emax); 346 } 347 348 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 349 350 void G4VEnergyLossProcess::SetEmModel(G4VEmModel* p, G4int index) 351 { 352 G4int n = emModels.size(); 353 if(index >= n) for(G4int i=n; i<index+1; i++) {emModels.push_back(0);} 354 emModels[index] = p; 355 } 356 357 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 358 359 G4VEmModel* G4VEnergyLossProcess::EmModel(G4int index) 360 { 361 G4VEmModel* p = 0; 362 if(index >= 0 && index < G4int(emModels.size())) p = emModels[index]; 363 return p; 364 } 365 366 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 367 368 G4VEmModel* G4VEnergyLossProcess::GetModelByIndex(G4int idx, G4bool ver) 369 { 370 return modelManager->GetModel(idx, ver); 371 } 372 373 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 374 375 G4int G4VEnergyLossProcess::NumberOfModels() 376 { 377 return modelManager->NumberOfModels(); 314 378 } 315 379 … … 364 428 } 365 429 366 Clear(); 430 Clean(); 431 lManager->EmConfigurator()->AddModels(); 367 432 368 433 // Base particle and set of models can be defined here … … 407 472 minSubRange, verboseLevel); 408 473 409 // Sub Cutoff Regime410 if (nSCoffRegions>0 ) {474 // Sub Cutoff and Deexcitation 475 if (nSCoffRegions>0 || nDERegions>0) { 411 476 theSubCuts = modelManager->SubCutoff(); 412 477 … … 414 479 G4ProductionCutsTable::GetProductionCutsTable(); 415 480 size_t numOfCouples = theCoupleTable->GetTableSize(); 416 idxSCoffRegions = new G4int[numOfCouples]; 481 482 if(nSCoffRegions>0) idxSCoffRegions = new G4bool[numOfCouples]; 483 if(nDERegions>0) idxDERegions = new G4bool[numOfCouples]; 417 484 418 485 for (size_t j=0; j<numOfCouples; j++) { … … 421 488 theCoupleTable->GetMaterialCutsCouple(j); 422 489 const G4ProductionCuts* pcuts = couple->GetProductionCuts(); 423 G4int reg = 0; 424 for(G4int i=0; i<nSCoffRegions; i++) { 425 if( pcuts == scoffRegions[i]->GetProductionCuts()) reg = 1; 426 } 427 idxSCoffRegions[j] = reg; 490 491 if(nSCoffRegions>0) { 492 G4bool reg = false; 493 for(G4int i=0; i<nSCoffRegions; i++) { 494 if( pcuts == scoffRegions[i]->GetProductionCuts()) reg = true; 495 } 496 idxSCoffRegions[j] = reg; 497 } 498 if(nDERegions>0) { 499 G4bool reg = false; 500 for(G4int i=0; i<nDERegions; i++) { 501 if( pcuts == deRegions[i]->GetProductionCuts()) reg = true; 502 } 503 idxDERegions[j] = reg; 504 } 428 505 } 429 506 } … … 445 522 } 446 523 } 524 if (nDERegions) { 525 G4cout << " Deexcitation is ON for regions: " << G4endl; 526 for (G4int i=0; i<nDERegions; i++) { 527 const G4Region* r = deRegions[i]; 528 G4cout << " " << r->GetName() << G4endl; 529 } 530 } 447 531 } 448 532 } … … 479 563 if(isIonisation) G4cout << " isIonisation flag = 1"; 480 564 G4cout << G4endl; 481 }482 }483 484 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....485 486 void G4VEnergyLossProcess::ActivateSubCutoff(G4bool val, const G4Region* r)487 {488 G4RegionStore* regionStore = G4RegionStore::GetInstance();489 if(val) {490 useSubCutoff = true;491 if (!r) {r = regionStore->GetRegion("DefaultRegionForTheWorld", false);}492 if (nSCoffRegions) {493 for (G4int i=0; i<nSCoffRegions; i++) {494 if (r == scoffRegions[i]) return;495 }496 }497 scoffRegions.push_back(r);498 nSCoffRegions++;499 } else {500 useSubCutoff = false;501 565 } 502 566 } … … 640 704 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 641 705 642 G4double G4VEnergyLossProcess::GetContinuousStepLimit( 643 const G4Track&, 644 G4double, G4double, G4double&) 645 { 646 return DBL_MAX; 706 void G4VEnergyLossProcess::PrintInfoDefinition() 707 { 708 if(0 < verboseLevel) { 709 G4cout << G4endl << GetProcessName() << ": for " 710 << particle->GetParticleName() 711 << " SubType= " << GetProcessSubType() 712 << G4endl 713 << " dE/dx and range tables from " 714 << G4BestUnit(minKinEnergy,"Energy") 715 << " to " << G4BestUnit(maxKinEnergy,"Energy") 716 << " in " << nBins << " bins" << G4endl 717 << " Lambda tables from threshold to " 718 << G4BestUnit(maxKinEnergy,"Energy") 719 << " in " << nBins << " bins, spline: " 720 << (G4LossTableManager::Instance())->SplineFlag() 721 << G4endl; 722 if(theRangeTableForLoss && isIonisation) { 723 G4cout << " finalRange(mm)= " << finalRange/mm 724 << ", dRoverRange= " << dRoverRange 725 << ", integral: " << integral 726 << ", fluct: " << lossFluctuationFlag 727 << ", linLossLimit= " << linLossLimit 728 << G4endl; 729 } 730 PrintInfo(); 731 modelManager->DumpModelList(verboseLevel); 732 if(theCSDARangeTable && isIonisation) { 733 G4cout << " CSDA range table up" 734 << " to " << G4BestUnit(maxKinEnergyCSDA,"Energy") 735 << " in " << nBinsCSDA << " bins" << G4endl; 736 } 737 if(nSCoffRegions>0 && isIonisation) { 738 G4cout << " Subcutoff sampling in " << nSCoffRegions 739 << " regions" << G4endl; 740 } 741 if(2 < verboseLevel) { 742 G4cout << " DEDXTable address= " << theDEDXTable << G4endl; 743 if(theDEDXTable && isIonisation) G4cout << (*theDEDXTable) << G4endl; 744 G4cout << "non restricted DEDXTable address= " 745 << theDEDXunRestrictedTable << G4endl; 746 if(theDEDXunRestrictedTable && isIonisation) { 747 G4cout << (*theDEDXunRestrictedTable) << G4endl; 748 } 749 if(theDEDXSubTable && isIonisation) { 750 G4cout << (*theDEDXSubTable) << G4endl; 751 } 752 G4cout << " CSDARangeTable address= " << theCSDARangeTable 753 << G4endl; 754 if(theCSDARangeTable && isIonisation) { 755 G4cout << (*theCSDARangeTable) << G4endl; 756 } 757 G4cout << " RangeTableForLoss address= " << theRangeTableForLoss 758 << G4endl; 759 if(theRangeTableForLoss && isIonisation) { 760 G4cout << (*theRangeTableForLoss) << G4endl; 761 } 762 G4cout << " InverseRangeTable address= " << theInverseRangeTable 763 << G4endl; 764 if(theInverseRangeTable && isIonisation) { 765 G4cout << (*theInverseRangeTable) << G4endl; 766 } 767 G4cout << " LambdaTable address= " << theLambdaTable << G4endl; 768 if(theLambdaTable && isIonisation) { 769 G4cout << (*theLambdaTable) << G4endl; 770 } 771 G4cout << " SubLambdaTable address= " << theSubLambdaTable << G4endl; 772 if(theSubLambdaTable && isIonisation) { 773 G4cout << (*theSubLambdaTable) << G4endl; 774 } 775 } 776 } 777 } 778 779 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 780 781 void G4VEnergyLossProcess::ActivateSubCutoff(G4bool val, const G4Region* r) 782 { 783 G4RegionStore* regionStore = G4RegionStore::GetInstance(); 784 const G4Region* reg = r; 785 if (!reg) {reg = regionStore->GetRegion("DefaultRegionForTheWorld", false);} 786 787 // the region is in the list 788 if (nSCoffRegions) { 789 for (G4int i=0; i<nSCoffRegions; i++) { 790 if (reg == scoffRegions[i]) { 791 if(!val) deRegions[i] = 0; 792 return; 793 } 794 } 795 } 796 797 // new region 798 if(val) { 799 useSubCutoff = true; 800 scoffRegions.push_back(reg); 801 nSCoffRegions++; 802 } else { 803 useSubCutoff = false; 804 } 805 } 806 807 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 808 809 void G4VEnergyLossProcess::ActivateDeexcitation(G4bool val, const G4Region* r) 810 { 811 G4RegionStore* regionStore = G4RegionStore::GetInstance(); 812 const G4Region* reg = r; 813 if (!reg) {reg = regionStore->GetRegion("DefaultRegionForTheWorld", false);} 814 815 // the region is in the list 816 if (nDERegions) { 817 for (G4int i=0; i<nDERegions; i++) { 818 if (reg == deRegions[i]) { 819 if(!val) deRegions[i] = 0; 820 return; 821 } 822 } 823 } 824 825 // new region 826 if(val) { 827 useDeexcitation = true; 828 deRegions.push_back(reg); 829 nDERegions++; 830 } else { 831 useDeexcitation = false; 832 } 647 833 } 648 834 … … 674 860 // <<" stepLimit= "<<x<<G4endl; 675 861 return x; 676 }677 678 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....679 680 G4double G4VEnergyLossProcess::GetMeanFreePath(681 const G4Track& track,682 G4double,683 G4ForceCondition* condition)684 685 {686 *condition = NotForced;687 return MeanFreePath(track);688 862 } 689 863 … … 806 980 if (length >= fRange) { 807 981 eloss = preStepKinEnergy; 808 currentModel->CorrectionsAlongStep(currentCouple, dynParticle, 809 eloss, esecdep, length); 982 if (useDeexcitation) { 983 if(idxDERegions[currentMaterialIndex]) { 984 currentModel->SampleDeexcitationAlongStep(currentMaterial, track, eloss); 985 if(eloss < 0.0) eloss = 0.0; 986 } 987 } 810 988 fParticleChange.SetProposedKineticEnergy(0.0); 811 fParticleChange.ProposeLocalEnergyDeposit( preStepKinEnergy);989 fParticleChange.ProposeLocalEnergyDeposit(eloss); 812 990 return &fParticleChange; 813 991 } … … 934 1112 G4double tmax = 935 1113 std::min(currentModel->MaxSecondaryKinEnergy(dynParticle),cut); 1114 G4double emean = eloss; 936 1115 eloss = fluc->SampleFluctuations(currentMaterial,dynParticle, 937 tmax,length,e loss);1116 tmax,length,emean); 938 1117 /* 939 1118 if(-1 < verboseLevel) … … 950 1129 eloss += esecdep; 951 1130 if(eloss < 0.0) eloss = 0.0; 1131 1132 // deexcitation 1133 else if (useDeexcitation) { 1134 if(idxDERegions[currentMaterialIndex]) { 1135 currentModel->SampleDeexcitationAlongStep(currentMaterial, track, eloss); 1136 if(eloss < 0.0) eloss = 0.0; 1137 } 1138 } 952 1139 953 1140 // Energy balanse … … 1106 1293 1107 1294 SelectModel(postStepScaledEnergy); 1295 if(useDeexcitation) { 1296 currentModel->SetDeexcitationFlag(idxDERegions[currentMaterialIndex]); 1297 } 1298 1108 1299 const G4DynamicParticle* dynParticle = track.GetDynamicParticle(); 1109 1300 G4double tcut = (*theCuts)[currentMaterialIndex]; … … 1140 1331 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1141 1332 1142 void G4VEnergyLossProcess::PrintInfoDefinition() 1143 { 1144 if(0 < verboseLevel) { 1145 G4cout << G4endl << GetProcessName() << ": for " 1146 << particle->GetParticleName() 1147 << " SubType= " << GetProcessSubType() 1148 << G4endl 1149 << " dE/dx and range tables from " 1150 << G4BestUnit(minKinEnergy,"Energy") 1151 << " to " << G4BestUnit(maxKinEnergy,"Energy") 1152 << " in " << nBins << " bins" << G4endl 1153 << " Lambda tables from threshold to " 1154 << G4BestUnit(maxKinEnergy,"Energy") 1155 << " in " << nBins << " bins, spline: " 1156 << (G4LossTableManager::Instance())->SplineFlag() 1333 G4bool G4VEnergyLossProcess::StorePhysicsTable( 1334 const G4ParticleDefinition* part, const G4String& directory, 1335 G4bool ascii) 1336 { 1337 G4bool res = true; 1338 if ( baseParticle || part != particle ) return res; 1339 1340 if(!StoreTable(part,theDEDXTable,ascii,directory,"DEDX")) 1341 {res = false;} 1342 1343 if(!StoreTable(part,theDEDXunRestrictedTable,ascii,directory,"DEDXnr")) 1344 {res = false;} 1345 1346 if(!StoreTable(part,theDEDXSubTable,ascii,directory,"SubDEDX")) 1347 {res = false;} 1348 1349 if(!StoreTable(part,theIonisationTable,ascii,directory,"Ionisation")) 1350 {res = false;} 1351 1352 if(!StoreTable(part,theIonisationSubTable,ascii,directory,"SubIonisation")) 1353 {res = false;} 1354 1355 if(isIonisation && 1356 !StoreTable(part,theCSDARangeTable,ascii,directory,"CSDARange")) 1357 {res = false;} 1358 1359 if(isIonisation && 1360 !StoreTable(part,theRangeTableForLoss,ascii,directory,"Range")) 1361 {res = false;} 1362 1363 if(isIonisation && 1364 !StoreTable(part,theInverseRangeTable,ascii,directory,"InverseRange")) 1365 {res = false;} 1366 1367 if(!StoreTable(part,theLambdaTable,ascii,directory,"Lambda")) 1368 {res = false;} 1369 1370 if(!StoreTable(part,theSubLambdaTable,ascii,directory,"SubLambda")) 1371 {res = false;} 1372 1373 if ( res ) { 1374 if(0 < verboseLevel) { 1375 G4cout << "Physics tables are stored for " << particle->GetParticleName() 1376 << " and process " << GetProcessName() 1377 << " in the directory <" << directory 1378 << "> " << G4endl; 1379 } 1380 } else { 1381 G4cout << "Fail to store Physics Tables for " 1382 << particle->GetParticleName() 1383 << " and process " << GetProcessName() 1384 << " in the directory <" << directory 1385 << "> " << G4endl; 1386 } 1387 return res; 1388 } 1389 1390 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 1391 1392 G4bool G4VEnergyLossProcess::RetrievePhysicsTable( 1393 const G4ParticleDefinition* part, const G4String& directory, 1394 G4bool ascii) 1395 { 1396 G4bool res = true; 1397 const G4String particleName = part->GetParticleName(); 1398 1399 if(1 < verboseLevel) { 1400 G4cout << "G4VEnergyLossProcess::RetrievePhysicsTable() for " 1401 << particleName << " and process " << GetProcessName() 1402 << "; tables_are_built= " << tablesAreBuilt 1157 1403 << G4endl; 1158 if(theRangeTableForLoss && isIonisation) { 1159 G4cout << " finalRange(mm)= " << finalRange/mm 1160 << ", dRoverRange= " << dRoverRange 1161 << ", integral: " << integral 1162 << ", fluct: " << lossFluctuationFlag 1163 << ", linLossLimit= " << linLossLimit 1164 << G4endl; 1165 } 1166 PrintInfo(); 1167 modelManager->DumpModelList(verboseLevel); 1168 if(theCSDARangeTable && isIonisation) { 1169 G4cout << " CSDA range table up" 1170 << " to " << G4BestUnit(maxKinEnergyCSDA,"Energy") 1171 << " in " << nBinsCSDA << " bins" << G4endl; 1172 } 1173 if(nSCoffRegions>0 && isIonisation) { 1174 G4cout << " Subcutoff sampling in " << nSCoffRegions 1175 << " regions" << G4endl; 1176 } 1177 if(2 < verboseLevel) { 1178 G4cout << " DEDXTable address= " << theDEDXTable << G4endl; 1179 if(theDEDXTable && isIonisation) G4cout << (*theDEDXTable) << G4endl; 1180 G4cout << "non restricted DEDXTable address= " 1181 << theDEDXunRestrictedTable << G4endl; 1182 if(theDEDXunRestrictedTable && isIonisation) { 1183 G4cout << (*theDEDXunRestrictedTable) << G4endl; 1184 } 1185 if(theDEDXSubTable && isIonisation) { 1186 G4cout << (*theDEDXSubTable) << G4endl; 1187 } 1188 G4cout << " CSDARangeTable address= " << theCSDARangeTable 1404 } 1405 if(particle == part) { 1406 1407 // G4bool yes = true; 1408 if ( !baseParticle ) { 1409 1410 G4bool fpi = true; 1411 if(!RetrieveTable(part,theDEDXTable,ascii,directory,"DEDX",fpi)) 1412 {fpi = false;} 1413 1414 if(!RetrieveTable(part,theIonisationTable,ascii,directory,"Ionisation",false)) 1415 {fpi = false;} 1416 1417 if(!RetrieveTable(part,theRangeTableForLoss,ascii,directory,"Range",fpi)) 1418 {res = false;} 1419 1420 if(!RetrieveTable(part,theDEDXunRestrictedTable,ascii,directory,"DEDXnr",false)) 1421 {res = false;} 1422 1423 if(!RetrieveTable(part,theCSDARangeTable,ascii,directory,"CSDARange",false)) 1424 {res = false;} 1425 1426 if(!RetrieveTable(part,theInverseRangeTable,ascii,directory,"InverseRange",fpi)) 1427 {res = false;} 1428 1429 if(!RetrieveTable(part,theLambdaTable,ascii,directory,"Lambda",true)) 1430 {res = false;} 1431 1432 G4bool yes = false; 1433 if(nSCoffRegions > 0) {yes = true;} 1434 1435 if(!RetrieveTable(part,theDEDXSubTable,ascii,directory,"SubDEDX",yes)) 1436 {res = false;} 1437 1438 if(!RetrieveTable(part,theSubLambdaTable,ascii,directory,"SubLambda",yes)) 1439 {res = false;} 1440 1441 if(!fpi) yes = false; 1442 if(!RetrieveTable(part,theIonisationSubTable,ascii,directory,"SubIonisation",yes)) 1443 {res = false;} 1444 } 1445 } 1446 1447 return res; 1448 } 1449 1450 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 1451 1452 G4bool G4VEnergyLossProcess::StoreTable(const G4ParticleDefinition* part, 1453 G4PhysicsTable* aTable, G4bool ascii, 1454 const G4String& directory, 1455 const G4String& tname) 1456 { 1457 G4bool res = true; 1458 if ( aTable ) { 1459 const G4String name = GetPhysicsTableFileName(part,directory,tname,ascii); 1460 if( !aTable->StorePhysicsTable(name,ascii)) res = false; 1461 } 1462 return res; 1463 } 1464 1465 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 1466 1467 G4bool G4VEnergyLossProcess::RetrieveTable(const G4ParticleDefinition* part, 1468 G4PhysicsTable* aTable, G4bool ascii, 1469 const G4String& directory, 1470 const G4String& tname, 1471 G4bool mandatory) 1472 { 1473 G4bool res = true; 1474 G4String filename = GetPhysicsTableFileName(part,directory,tname,ascii); 1475 G4bool yes = aTable->ExistPhysicsTable(filename); 1476 if(yes) { 1477 yes = G4PhysicsTableHelper::RetrievePhysicsTable(aTable,filename,ascii); 1478 if((G4LossTableManager::Instance())->SplineFlag()) { 1479 size_t n = aTable->length(); 1480 for(size_t i=0; i<n; i++) {(*aTable)[i]->SetSpline(true);} 1481 } 1482 } 1483 if(yes) { 1484 if (0 < verboseLevel) { 1485 G4cout << tname << " table for " << part->GetParticleName() 1486 << " is Retrieved from <" << filename << ">" 1189 1487 << G4endl; 1190 if(theCSDARangeTable && isIonisation) { 1191 G4cout << (*theCSDARangeTable) << G4endl; 1192 } 1193 G4cout << " RangeTableForLoss address= " << theRangeTableForLoss 1488 } 1489 } else { 1490 if(mandatory) res = false; 1491 if(mandatory || 1 < verboseLevel) { 1492 G4cout << tname << " table for " << part->GetParticleName() 1493 << " from file <" 1494 << filename << "> is not Retrieved" 1194 1495 << G4endl; 1195 if(theRangeTableForLoss && isIonisation) { 1196 G4cout << (*theRangeTableForLoss) << G4endl; 1197 } 1198 G4cout << " InverseRangeTable address= " << theInverseRangeTable 1199 << G4endl; 1200 if(theInverseRangeTable && isIonisation) { 1201 G4cout << (*theInverseRangeTable) << G4endl; 1202 } 1203 G4cout << " LambdaTable address= " << theLambdaTable << G4endl; 1204 if(theLambdaTable && isIonisation) { 1205 G4cout << (*theLambdaTable) << G4endl; 1206 } 1207 G4cout << " SubLambdaTable address= " << theSubLambdaTable << G4endl; 1208 if(theSubLambdaTable && isIonisation) { 1209 G4cout << (*theSubLambdaTable) << G4endl; 1210 } 1496 } 1497 } 1498 return res; 1499 } 1500 1501 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1502 1503 G4double G4VEnergyLossProcess::GetDEDXDispersion( 1504 const G4MaterialCutsCouple *couple, 1505 const G4DynamicParticle* dp, 1506 G4double length) 1507 { 1508 DefineMaterial(couple); 1509 G4double ekin = dp->GetKineticEnergy(); 1510 SelectModel(ekin*massRatio); 1511 G4double tmax = currentModel->MaxSecondaryKinEnergy(dp); 1512 tmax = std::min(tmax,(*theCuts)[currentMaterialIndex]); 1513 G4double d = 0.0; 1514 G4VEmFluctuationModel* fm = currentModel->GetModelOfFluctuations(); 1515 if(fm) d = fm->Dispersion(currentMaterial,dp,tmax,length); 1516 return d; 1517 } 1518 1519 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1520 1521 G4double G4VEnergyLossProcess::CrossSectionPerVolume( 1522 G4double kineticEnergy, const G4MaterialCutsCouple* couple) 1523 { 1524 // Cross section per volume is calculated 1525 DefineMaterial(couple); 1526 G4double cross = 0.0; 1527 G4bool b; 1528 if(theLambdaTable) { 1529 cross = 1530 ((*theLambdaTable)[currentMaterialIndex])->GetValue(kineticEnergy, b); 1531 } else { 1532 SelectModel(kineticEnergy); 1533 cross = 1534 currentModel->CrossSectionPerVolume(currentMaterial, 1535 particle, kineticEnergy, 1536 (*theCuts)[currentMaterialIndex]); 1537 } 1538 1539 return cross; 1540 } 1541 1542 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1543 1544 G4double G4VEnergyLossProcess::MeanFreePath(const G4Track& track) 1545 { 1546 DefineMaterial(track.GetMaterialCutsCouple()); 1547 preStepLambda = GetLambdaForScaledEnergy(track.GetKineticEnergy()*massRatio); 1548 G4double x = DBL_MAX; 1549 if(DBL_MIN < preStepLambda) x = 1.0/preStepLambda; 1550 return x; 1551 } 1552 1553 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1554 1555 G4double G4VEnergyLossProcess::ContinuousStepLimit(const G4Track& track, 1556 G4double x, G4double y, 1557 G4double& z) 1558 { 1559 G4GPILSelection sel; 1560 return AlongStepGetPhysicalInteractionLength(track, x, y, z, &sel); 1561 } 1562 1563 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1564 1565 G4double G4VEnergyLossProcess::GetMeanFreePath( 1566 const G4Track& track, 1567 G4double, 1568 G4ForceCondition* condition) 1569 1570 { 1571 *condition = NotForced; 1572 return MeanFreePath(track); 1573 } 1574 1575 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1576 1577 G4double G4VEnergyLossProcess::GetContinuousStepLimit( 1578 const G4Track&, 1579 G4double, G4double, G4double&) 1580 { 1581 return DBL_MAX; 1582 } 1583 1584 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1585 1586 G4PhysicsVector* G4VEnergyLossProcess::LambdaPhysicsVector( 1587 const G4MaterialCutsCouple* couple, G4double cut) 1588 { 1589 // G4double cut = (*theCuts)[couple->GetIndex()]; 1590 // G4int nbins = nLambdaBins; 1591 G4double tmin = 1592 std::max(MinPrimaryEnergy(particle, couple->GetMaterial(), cut), 1593 minKinEnergy); 1594 if(tmin >= maxKinEnergy) tmin = 0.5*maxKinEnergy; 1595 G4PhysicsVector* v = new G4PhysicsLogVector(tmin, maxKinEnergy, nBins); 1596 v->SetSpline((G4LossTableManager::Instance())->SplineFlag()); 1597 1598 return v; 1599 } 1600 1601 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1602 1603 void G4VEnergyLossProcess::AddCollaborativeProcess( 1604 G4VEnergyLossProcess* p) 1605 { 1606 G4bool add = true; 1607 if(p->GetProcessName() != "eBrem") add = false; 1608 if(add && nProcesses > 0) { 1609 for(G4int i=0; i<nProcesses; i++) { 1610 if(p == scProcesses[i]) { 1611 add = false; 1612 break; 1613 } 1614 } 1615 } 1616 if(add) { 1617 scProcesses.push_back(p); 1618 nProcesses++; 1619 if (1 < verboseLevel) { 1620 G4cout << "### The process " << p->GetProcessName() 1621 << " is added to the list of collaborative processes of " 1622 << GetProcessName() << G4endl; 1211 1623 } 1212 1624 } … … 1377 1789 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 1378 1790 1379 G4PhysicsVector* G4VEnergyLossProcess::LambdaPhysicsVector(1380 const G4MaterialCutsCouple* couple, G4double cut)1381 {1382 // G4double cut = (*theCuts)[couple->GetIndex()];1383 // G4int nbins = nLambdaBins;1384 G4double tmin =1385 std::max(MinPrimaryEnergy(particle, couple->GetMaterial(), cut),1386 minKinEnergy);1387 if(tmin >= maxKinEnergy) tmin = 0.5*maxKinEnergy;1388 G4PhysicsVector* v = new G4PhysicsLogVector(tmin, maxKinEnergy, nBins);1389 v->SetSpline((G4LossTableManager::Instance())->SplineFlag());1390 1391 return v;1392 }1393 1394 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....1395 1396 G4double G4VEnergyLossProcess::CrossSectionPerVolume(1397 G4double kineticEnergy, const G4MaterialCutsCouple* couple)1398 {1399 // Cross section per volume is calculated1400 DefineMaterial(couple);1401 G4double cross = 0.0;1402 G4bool b;1403 if(theLambdaTable) {1404 cross =1405 ((*theLambdaTable)[currentMaterialIndex])->GetValue(kineticEnergy, b);1406 } else {1407 SelectModel(kineticEnergy);1408 cross =1409 currentModel->CrossSectionPerVolume(currentMaterial,1410 particle, kineticEnergy,1411 (*theCuts)[currentMaterialIndex]);1412 }1413 1414 return cross;1415 }1416 1417 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....1418 1419 G4bool G4VEnergyLossProcess::StorePhysicsTable(1420 const G4ParticleDefinition* part, const G4String& directory,1421 G4bool ascii)1422 {1423 G4bool res = true;1424 if ( baseParticle || part != particle ) return res;1425 1426 if(!StoreTable(part,theDEDXTable,ascii,directory,"DEDX"))1427 {res = false;}1428 1429 if(!StoreTable(part,theDEDXunRestrictedTable,ascii,directory,"DEDXnr"))1430 {res = false;}1431 1432 if(!StoreTable(part,theDEDXSubTable,ascii,directory,"SubDEDX"))1433 {res = false;}1434 1435 if(!StoreTable(part,theIonisationTable,ascii,directory,"Ionisation"))1436 {res = false;}1437 1438 if(!StoreTable(part,theIonisationSubTable,ascii,directory,"SubIonisation"))1439 {res = false;}1440 1441 if(isIonisation &&1442 !StoreTable(part,theCSDARangeTable,ascii,directory,"CSDARange"))1443 {res = false;}1444 1445 if(isIonisation &&1446 !StoreTable(part,theRangeTableForLoss,ascii,directory,"Range"))1447 {res = false;}1448 1449 if(isIonisation &&1450 !StoreTable(part,theInverseRangeTable,ascii,directory,"InverseRange"))1451 {res = false;}1452 1453 if(!StoreTable(part,theLambdaTable,ascii,directory,"Lambda"))1454 {res = false;}1455 1456 if(!StoreTable(part,theSubLambdaTable,ascii,directory,"SubLambda"))1457 {res = false;}1458 1459 if ( res ) {1460 if(0 < verboseLevel) {1461 G4cout << "Physics tables are stored for " << particle->GetParticleName()1462 << " and process " << GetProcessName()1463 << " in the directory <" << directory1464 << "> " << G4endl;1465 }1466 } else {1467 G4cout << "Fail to store Physics Tables for "1468 << particle->GetParticleName()1469 << " and process " << GetProcessName()1470 << " in the directory <" << directory1471 << "> " << G4endl;1472 }1473 return res;1474 }1475 1476 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....1477 1478 G4bool G4VEnergyLossProcess::StoreTable(const G4ParticleDefinition* part,1479 G4PhysicsTable* aTable, G4bool ascii,1480 const G4String& directory,1481 const G4String& tname)1482 {1483 G4bool res = true;1484 if ( aTable ) {1485 const G4String name = GetPhysicsTableFileName(part,directory,tname,ascii);1486 if( !aTable->StorePhysicsTable(name,ascii)) res = false;1487 }1488 return res;1489 }1490 1491 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....1492 1493 G4bool G4VEnergyLossProcess::RetrievePhysicsTable(1494 const G4ParticleDefinition* part, const G4String& directory,1495 G4bool ascii)1496 {1497 G4bool res = true;1498 const G4String particleName = part->GetParticleName();1499 1500 if(1 < verboseLevel) {1501 G4cout << "G4VEnergyLossProcess::RetrievePhysicsTable() for "1502 << particleName << " and process " << GetProcessName()1503 << "; tables_are_built= " << tablesAreBuilt1504 << G4endl;1505 }1506 if(particle == part) {1507 1508 // G4bool yes = true;1509 if ( !baseParticle ) {1510 1511 G4bool fpi = true;1512 if(!RetrieveTable(part,theDEDXTable,ascii,directory,"DEDX",fpi))1513 {fpi = false;}1514 1515 if(!RetrieveTable(part,theIonisationTable,ascii,directory,"Ionisation",false))1516 {fpi = false;}1517 1518 if(!RetrieveTable(part,theRangeTableForLoss,ascii,directory,"Range",fpi))1519 {res = false;}1520 1521 if(!RetrieveTable(part,theDEDXunRestrictedTable,ascii,directory,"DEDXnr",false))1522 {res = false;}1523 1524 if(!RetrieveTable(part,theCSDARangeTable,ascii,directory,"CSDARange",false))1525 {res = false;}1526 1527 if(!RetrieveTable(part,theInverseRangeTable,ascii,directory,"InverseRange",fpi))1528 {res = false;}1529 1530 if(!RetrieveTable(part,theLambdaTable,ascii,directory,"Lambda",true))1531 {res = false;}1532 1533 G4bool yes = false;1534 if(nSCoffRegions > 0) {yes = true;}1535 1536 if(!RetrieveTable(part,theDEDXSubTable,ascii,directory,"SubDEDX",yes))1537 {res = false;}1538 1539 if(!RetrieveTable(part,theSubLambdaTable,ascii,directory,"SubLambda",yes))1540 {res = false;}1541 1542 if(!fpi) yes = false;1543 if(!RetrieveTable(part,theIonisationSubTable,ascii,directory,"SubIonisation",yes))1544 {res = false;}1545 }1546 }1547 1548 return res;1549 }1550 1551 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....1552 1553 G4bool G4VEnergyLossProcess::RetrieveTable(const G4ParticleDefinition* part,1554 G4PhysicsTable* aTable, G4bool ascii,1555 const G4String& directory,1556 const G4String& tname,1557 G4bool mandatory)1558 {1559 G4bool res = true;1560 G4String filename = GetPhysicsTableFileName(part,directory,tname,ascii);1561 G4bool yes = aTable->ExistPhysicsTable(filename);1562 if(yes) {1563 yes = G4PhysicsTableHelper::RetrievePhysicsTable(aTable,filename,ascii);1564 if((G4LossTableManager::Instance())->SplineFlag()) {1565 size_t n = aTable->length();1566 for(size_t i=0; i<n; i++) {(*aTable)[i]->SetSpline(true);}1567 }1568 }1569 if(yes) {1570 if (0 < verboseLevel) {1571 G4cout << tname << " table for " << part->GetParticleName()1572 << " is Retrieved from <" << filename << ">"1573 << G4endl;1574 }1575 } else {1576 if(mandatory) res = false;1577 if(mandatory || 1 < verboseLevel) {1578 G4cout << tname << " table for " << part->GetParticleName()1579 << " from file <"1580 << filename << "> is not Retrieved"1581 << G4endl;1582 }1583 }1584 return res;1585 }1586 1587 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....1588 1589 void G4VEnergyLossProcess::AddCollaborativeProcess(1590 G4VEnergyLossProcess* p)1591 {1592 G4bool add = true;1593 if(p->GetProcessName() != "eBrem") add = false;1594 if(add && nProcesses > 0) {1595 for(G4int i=0; i<nProcesses; i++) {1596 if(p == scProcesses[i]) {1597 add = false;1598 break;1599 }1600 }1601 }1602 if(add) {1603 scProcesses.push_back(p);1604 nProcesses++;1605 if (1 < verboseLevel) {1606 G4cout << "### The process " << p->GetProcessName()1607 << " is added to the list of collaborative processes of "1608 << GetProcessName() << G4endl;1609 }1610 }1611 }1612 1613 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....1614 1615 G4double G4VEnergyLossProcess::GetDEDXDispersion(1616 const G4MaterialCutsCouple *couple,1617 const G4DynamicParticle* dp,1618 G4double length)1619 {1620 DefineMaterial(couple);1621 G4double ekin = dp->GetKineticEnergy();1622 SelectModel(ekin*massRatio);1623 G4double tmax = currentModel->MaxSecondaryKinEnergy(dp);1624 tmax = std::min(tmax,(*theCuts)[currentMaterialIndex]);1625 G4double d = 0.0;1626 G4VEmFluctuationModel* fm = currentModel->GetModelOfFluctuations();1627 if(fm) d = fm->Dispersion(currentMaterial,dp,tmax,length);1628 return d;1629 }1630 1631 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....1632 1633 void G4VEnergyLossProcess::ActivateDeexcitation(G4bool, const G4Region*)1634 {}1635 1636 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....1637 -
trunk/source/processes/electromagnetic/utils/src/G4VMscModel.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4VMscModel.cc,v 1. 4 2008/03/10 18:39:45vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4VMscModel.cc,v 1.12 2009/05/10 19:36:19 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 49 49 50 50 #include "G4VMscModel.hh" 51 #include "G4ParticleChangeForMSC.hh" 52 #include "G4TransportationManager.hh" 51 53 52 54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... … … 55 57 G4VMscModel::G4VMscModel(const G4String& nam): 56 58 G4VEmModel(nam), 57 facrange(0.02), 59 safetyHelper(0), 60 facrange(0.04), 58 61 facgeom(2.5), 59 62 facsafety(0.25), … … 61 64 dtrl(0.05), 62 65 lambdalimit(mm), 66 geommax(1.e50*mm), 63 67 steppingAlgorithm(fUseSafety), 64 68 samplez(false), … … 72 76 73 77 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 78 79 G4ParticleChangeForMSC* G4VMscModel::GetParticleChangeForMSC() 80 { 81 G4ParticleChangeForMSC* p = 0; 82 if (pParticleChange) { 83 p = static_cast<G4ParticleChangeForMSC*>(pParticleChange); 84 } else { 85 p = new G4ParticleChangeForMSC(); 86 } 87 return p; 88 } 89 90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 91 92 void G4VMscModel::InitialiseSafetyHelper() 93 { 94 if(!safetyHelper) { 95 safetyHelper = G4TransportationManager::GetTransportationManager() 96 ->GetSafetyHelper(); 97 safetyHelper->InitialiseHelper(); 98 } 99 } 100 101 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 102 103 void G4VMscModel::ComputeDisplacement(G4ParticleChangeForMSC* fParticleChange, 104 const G4ThreeVector& dir, 105 G4double displacement, 106 G4double postsafety) 107 { 108 const G4ThreeVector* pos = fParticleChange->GetProposedPosition(); 109 G4double r = displacement; 110 if(r > postsafety) { 111 G4double newsafety = safetyHelper->ComputeSafety(*pos); 112 if(r > newsafety) r = newsafety; 113 } 114 if(r > 0.) { 115 116 // compute new endpoint of the Step 117 G4ThreeVector newPosition = *pos + r*dir; 118 119 // definitely not on boundary 120 if(displacement == r) { 121 safetyHelper->ReLocateWithinVolume(newPosition); 122 123 } else { 124 // check safety after displacement 125 G4double postsafety = safetyHelper->ComputeSafety(newPosition); 126 127 // displacement to boundary 128 if(postsafety <= 0.0) { 129 safetyHelper->Locate(newPosition, 130 *fParticleChange->GetProposedMomentumDirection()); 131 132 // not on the boundary 133 } else { 134 safetyHelper->ReLocateWithinVolume(newPosition); 135 } 136 } 137 fParticleChange->ProposePosition(newPosition); 138 } 139 } 140 141 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 142 143 void G4VMscModel::SampleScattering(const G4DynamicParticle*, G4double) 144 {} 145 146 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 147 148 G4double G4VMscModel::ComputeTruePathLengthLimit(const G4Track&, 149 G4PhysicsTable*, 150 G4double) 151 { 152 return DBL_MAX; 153 } 154 155 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 156 157 G4double G4VMscModel::ComputeGeomPathLength(G4double truePathLength) 158 { 159 return truePathLength; 160 } 161 162 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 163 164 G4double G4VMscModel::ComputeTrueStepLength(G4double geomPathLength) 165 { 166 return geomPathLength; 167 } 168 169 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 170 171 void G4VMscModel::SampleSecondaries(std::vector<G4DynamicParticle*>*, 172 const G4MaterialCutsCouple*, 173 const G4DynamicParticle*, 174 G4double, G4double) 175 {} 176 177 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... -
trunk/source/processes/electromagnetic/utils/src/G4VMultipleScattering.cc
r1007 r1055 24 24 // ******************************************************************** 25 25 // 26 // $Id: G4VMultipleScattering.cc,v 1.6 0 2008/11/20 20:32:40vnivanch Exp $27 // GEANT4 tag $Name: geant4-09-0 2$26 // $Id: G4VMultipleScattering.cc,v 1.66 2009/05/27 11:55:02 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-03-beta-cand-01 $ 28 28 // 29 29 // ------------------------------------------------------------------- … … 83 83 #include "G4GenericIon.hh" 84 84 #include "G4Electron.hh" 85 #include "G4EmConfigurator.hh" 85 86 86 87 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... … … 94 95 stepLimit(fUseSafety), 95 96 skin(3.0), 96 facrange(0.0 2),97 facrange(0.04), 97 98 facgeom(2.5), 98 99 latDisplasment(true), … … 105 106 // Size of tables assuming spline 106 107 minKinEnergy = 0.1*keV; 107 maxKinEnergy = 10 0.0*TeV;108 nBins = 84;108 maxKinEnergy = 10.0*TeV; 109 nBins = 77; 109 110 110 111 // default limit on polar angle … … 132 133 } 133 134 (G4LossTableManager::Instance())->DeRegister(this); 135 } 136 137 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 138 139 void G4VMultipleScattering::AddEmModel(G4int order, G4VEmModel* p, 140 const G4Region* region) 141 { 142 G4VEmFluctuationModel* fm = 0; 143 modelManager->AddEmModel(order, p, fm, region); 144 if(p) p->SetParticleChange(pParticleChange); 145 } 146 147 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 148 149 G4VEmModel* G4VMultipleScattering::GetModelByIndex(G4int idx, G4bool ver) 150 { 151 return modelManager->GetModel(idx, ver); 134 152 } 135 153 … … 208 226 << G4endl; 209 227 } 228 229 (G4LossTableManager::Instance())->EmConfigurator()->AddModels(); 210 230 211 231 if(firstParticle == &part) {
Note: See TracChangeset
for help on using the changeset viewer.
