Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // 26 // 27 // Geant4 class G4EmTableUtil 27 // Geant4 class G4EmTableUtil 28 // 28 // 29 // Author V.Ivanchenko 14.03.2022 29 // Author V.Ivanchenko 14.03.2022 30 // 30 // 31 31 32 #include "G4EmTableUtil.hh" 32 #include "G4EmTableUtil.hh" 33 #include "G4RegionStore.hh" 33 #include "G4RegionStore.hh" 34 #include "G4ProductionCutsTable.hh" 34 #include "G4ProductionCutsTable.hh" 35 #include "G4EmParameters.hh" 35 #include "G4EmParameters.hh" 36 #include "G4EmUtility.hh" 36 #include "G4EmUtility.hh" 37 #include "G4LossTableManager.hh" 37 #include "G4LossTableManager.hh" 38 #include "G4EmTableType.hh" 38 #include "G4EmTableType.hh" 39 #include "G4PhysicsModelCatalog.hh" 39 #include "G4PhysicsModelCatalog.hh" 40 #include "G4LowEnergyEmProcessSubType.hh" 40 #include "G4LowEnergyEmProcessSubType.hh" 41 #include "G4PhysicsTableHelper.hh" 41 #include "G4PhysicsTableHelper.hh" 42 #include "G4PhysicsLogVector.hh" 42 #include "G4PhysicsLogVector.hh" 43 #include "G4ProcessManager.hh" 43 #include "G4ProcessManager.hh" 44 #include "G4UIcommand.hh" 44 #include "G4UIcommand.hh" 45 #include "G4GenericIon.hh" 45 #include "G4GenericIon.hh" 46 #include <iostream> 46 #include <iostream> 47 47 48 //....oooOO0OOooo........oooOO0OOooo........oo 48 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 49 49 50 const G4DataVector* 50 const G4DataVector* 51 G4EmTableUtil::PrepareEmProcess(G4VEmProcess* 51 G4EmTableUtil::PrepareEmProcess(G4VEmProcess* proc, 52 const G4Partic 52 const G4ParticleDefinition* part, 53 const G4Partic 53 const G4ParticleDefinition* secPart, 54 G4EmModelManager* modelManager, 54 G4EmModelManager* modelManager, 55 const G4double 55 const G4double& maxKinEnergy, 56 G4int& secID, G4int& tripletID, 56 G4int& secID, G4int& tripletID, 57 G4int& mainSec 57 G4int& mainSec, const G4int& verb, 58 const G4bool& 58 const G4bool& master) 59 { 59 { 60 G4EmParameters* param = G4EmParameters::Inst 60 G4EmParameters* param = G4EmParameters::Instance(); 61 61 62 // initialisation of models 62 // initialisation of models 63 G4double plimit = param->MscThetaLimit(); 63 G4double plimit = param->MscThetaLimit(); 64 G4int nModels = modelManager->NumberOfModels 64 G4int nModels = modelManager->NumberOfModels(); 65 for(G4int i=0; i<nModels; ++i) { 65 for(G4int i=0; i<nModels; ++i) { 66 G4VEmModel* mod = modelManager->GetModel(i 66 G4VEmModel* mod = modelManager->GetModel(i); 67 if(nullptr == mod) { continue; } 67 if(nullptr == mod) { continue; } 68 mod->SetPolarAngleLimit(plimit); 68 mod->SetPolarAngleLimit(plimit); >> 69 mod->SetMasterThread(master); 69 if(mod->HighEnergyLimit() > maxKinEnergy) 70 if(mod->HighEnergyLimit() > maxKinEnergy) { 70 mod->SetHighEnergyLimit(maxKinEnergy); 71 mod->SetHighEnergyLimit(maxKinEnergy); 71 } 72 } 72 proc->SetEmModel(mod); 73 proc->SetEmModel(mod); 73 } 74 } 74 75 75 // defined ID of secondary particles and ver 76 // defined ID of secondary particles and verbosity 76 G4int stype = proc->GetProcessSubType(); 77 G4int stype = proc->GetProcessSubType(); 77 if(stype == fAnnihilation) { 78 if(stype == fAnnihilation) { 78 secID = _Annihilation; 79 secID = _Annihilation; 79 tripletID = _TripletGamma; 80 tripletID = _TripletGamma; 80 } else if(stype == fGammaConversion) { 81 } else if(stype == fGammaConversion) { 81 secID = _PairProduction; 82 secID = _PairProduction; 82 mainSec = 2; 83 mainSec = 2; 83 } else if(stype == fPhotoElectricEffect) { 84 } else if(stype == fPhotoElectricEffect) { 84 secID = _PhotoElectron; 85 secID = _PhotoElectron; 85 } else if(stype == fComptonScattering) { 86 } else if(stype == fComptonScattering) { 86 secID = _ComptonElectron; 87 secID = _ComptonElectron; 87 } else if(stype >= fLowEnergyElastic) { 88 } else if(stype >= fLowEnergyElastic) { 88 secID = fDNAUnknownModel; 89 secID = fDNAUnknownModel; 89 } 90 } 90 if(master) { 91 if(master) { 91 proc->SetVerboseLevel(param->Verbose()); 92 proc->SetVerboseLevel(param->Verbose()); 92 } else { 93 } else { 93 proc->SetVerboseLevel(param->WorkerVerbose 94 proc->SetVerboseLevel(param->WorkerVerbose()); 94 } 95 } 95 96 96 // model initialisation 97 // model initialisation 97 const G4DataVector* cuts = modelManager->Ini 98 const G4DataVector* cuts = modelManager->Initialise(part, secPart, verb); 98 99 99 if(1 < verb) { 100 if(1 < verb) { 100 G4cout << "### G4EmTableUtil::PreparePhysi 101 G4cout << "### G4EmTableUtil::PreparePhysicsTable() done for " 101 << proc->GetProcessName() 102 << proc->GetProcessName() 102 << " and particle " << part->GetPar 103 << " and particle " << part->GetParticleName() 103 << G4endl; 104 << G4endl; 104 } 105 } 105 return cuts; 106 return cuts; 106 } 107 } 107 108 108 //....oooOO0OOooo........oooOO0OOooo........oo 109 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 109 110 110 void G4EmTableUtil::BuildEmProcess(G4VEmProces 111 void G4EmTableUtil::BuildEmProcess(G4VEmProcess* proc, 111 const G4VEm 112 const G4VEmProcess* masterProc, 112 const G4ParticleDefinition* f 113 const G4ParticleDefinition* firstPart, 113 const G4ParticleDefinition* p 114 const G4ParticleDefinition* part, 114 const G4int nModels, const G4 115 const G4int nModels, const G4int verb, 115 const G4boo 116 const G4bool master, const G4bool isLocked, 116 const G4boo 117 const G4bool toBuild, G4bool& baseMat) 117 { 118 { 118 G4String num = part->GetParticleName(); 119 G4String num = part->GetParticleName(); 119 if(1 < verb) { 120 if(1 < verb) { 120 G4cout << "### G4EmTableUtil::BuildPhysics 121 G4cout << "### G4EmTableUtil::BuildPhysicsTable() for " 121 << proc->GetProcessName() << " and 122 << proc->GetProcessName() << " and particle " << num 122 << " buildLambdaTable=" << toBuild 123 << " buildLambdaTable=" << toBuild << " master= " << master 123 << G4endl; 124 << G4endl; 124 } 125 } 125 126 126 if(firstPart == part) { 127 if(firstPart == part) { 127 128 128 // worker initialisation 129 // worker initialisation 129 if(!master) { 130 if(!master) { 130 proc->SetLambdaTable(masterProc->LambdaT 131 proc->SetLambdaTable(masterProc->LambdaTable()); 131 proc->SetLambdaTablePrim(masterProc->Lam 132 proc->SetLambdaTablePrim(masterProc->LambdaTablePrim()); 132 proc->SetCrossSectionType(masterProc->Cr 133 proc->SetCrossSectionType(masterProc->CrossSectionType()); 133 proc->SetEnergyOfCrossSectionMax(masterP 134 proc->SetEnergyOfCrossSectionMax(masterProc->EnergyOfCrossSectionMax()); 134 135 135 // local initialisation of models 136 // local initialisation of models 136 baseMat = masterProc->UseBaseMaterial(); 137 baseMat = masterProc->UseBaseMaterial(); 137 G4bool printing = true; 138 G4bool printing = true; 138 for(G4int i=0; i<nModels; ++i) { 139 for(G4int i=0; i<nModels; ++i) { 139 G4VEmModel* mod = proc->GetModelByInde 140 G4VEmModel* mod = proc->GetModelByIndex(i, printing); 140 G4VEmModel* mod0= masterProc->GetModel 141 G4VEmModel* mod0= masterProc->GetModelByIndex(i, printing); 141 mod->SetUseBaseMaterials(baseMat); 142 mod->SetUseBaseMaterials(baseMat); 142 mod->InitialiseLocal(part, mod0); 143 mod->InitialiseLocal(part, mod0); 143 } 144 } 144 // master thread 145 // master thread 145 } else { 146 } else { 146 if(toBuild) { proc->BuildLambdaTable(); 147 if(toBuild) { proc->BuildLambdaTable(); } 147 auto fXSType = proc->CrossSectionType(); 148 auto fXSType = proc->CrossSectionType(); 148 auto v = proc->EnergyOfCrossSectionMax() 149 auto v = proc->EnergyOfCrossSectionMax(); 149 delete v; 150 delete v; 150 v = nullptr; 151 v = nullptr; 151 if(fXSType == fEmOnePeak) { 152 if(fXSType == fEmOnePeak) { 152 auto table = proc->LambdaTable(); 153 auto table = proc->LambdaTable(); 153 if(nullptr == table) { 154 if(nullptr == table) { 154 v = G4EmUtility::FindCrossSectionMax(proc, 155 v = G4EmUtility::FindCrossSectionMax(proc, part); 155 } else { 156 } else { 156 v = G4EmUtility::FindCrossSectionMax(table 157 v = G4EmUtility::FindCrossSectionMax(table); 157 } 158 } 158 if(nullptr == v) { proc->SetCrossSecti 159 if(nullptr == v) { proc->SetCrossSectionType(fEmIncreasing); } 159 } 160 } 160 proc->SetEnergyOfCrossSectionMax(v); 161 proc->SetEnergyOfCrossSectionMax(v); 161 } 162 } 162 } 163 } 163 // protection against double printout 164 // protection against double printout 164 if(isLocked) { return; } 165 if(isLocked) { return; } 165 166 166 // explicitly defined printout by particle n 167 // explicitly defined printout by particle name 167 if(1 < verb || (0 < verb && (num == "gamma" 168 if(1 < verb || (0 < verb && (num == "gamma" || num == "e-" || 168 num == "e+" || num == "mu+" || 169 num == "e+" || num == "mu+" || 169 num == "mu-" || num == "proton" 170 num == "mu-" || num == "proton"|| 170 num == "pi+" || num == "pi-" || 171 num == "pi+" || num == "pi-" || 171 num == "kaon+" || num == "kaon-" 172 num == "kaon+" || num == "kaon-" || 172 num == "alpha" || num == "anti_pr 173 num == "alpha" || num == "anti_proton" || 173 num == "GenericIon" || 174 num == "GenericIon" || 174 num == "alpha+" || num == "helium 175 num == "alpha+" || num == "helium" || 175 num == "hydrogen"))) { 176 num == "hydrogen"))) { 176 proc->StreamInfo(G4cout, *part); 177 proc->StreamInfo(G4cout, *part); 177 } 178 } 178 179 179 if(1 < verb) { 180 if(1 < verb) { 180 G4cout << "### G4EmTableUtil::BuildPhysics 181 G4cout << "### G4EmTableUtil::BuildPhysicsTable() done for " 181 << proc->GetProcessName() << " and 182 << proc->GetProcessName() << " and particle " << num 182 << " baseMat=" << baseMat << G4endl 183 << " baseMat=" << baseMat << G4endl; 183 } 184 } 184 } 185 } 185 186 186 //....oooOO0OOooo........oooOO0OOooo........oo 187 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 187 188 188 void G4EmTableUtil::BuildLambdaTable(G4VEmProc 189 void G4EmTableUtil::BuildLambdaTable(G4VEmProcess* proc, 189 const G4P 190 const G4ParticleDefinition* part, 190 G4EmModel 191 G4EmModelManager* modelManager, 191 G4LossTableBuilder* bld, 192 G4LossTableBuilder* bld, 192 G4Physics 193 G4PhysicsTable* theLambdaTable, 193 G4Physics 194 G4PhysicsTable* theLambdaTablePrim, 194 const G4d 195 const G4double minKinEnergy, 195 const G4d 196 const G4double minKinEnergyPrim, 196 const G4d 197 const G4double maxKinEnergy, 197 const G4d 198 const G4double scale, 198 const G4int verboseLevel, 199 const G4int verboseLevel, 199 const G4b 200 const G4bool startFromNull, 200 const G4b 201 const G4bool splineFlag) 201 { 202 { 202 if(1 < verboseLevel) { 203 if(1 < verboseLevel) { 203 G4cout << "G4EmTableUtil::BuildLambdaTable 204 G4cout << "G4EmTableUtil::BuildLambdaTable() for process " 204 << proc->GetProcessName() << " and 205 << proc->GetProcessName() << " and particle " 205 << part->GetParticleName() << G4end 206 << part->GetParticleName() << G4endl; 206 } 207 } 207 208 208 // Access to materials 209 // Access to materials 209 const G4ProductionCutsTable* theCoupleTable= 210 const G4ProductionCutsTable* theCoupleTable= 210 G4ProductionCutsTable::GetProductionCu 211 G4ProductionCutsTable::GetProductionCutsTable(); 211 std::size_t numOfCouples = theCoupleTable->G 212 std::size_t numOfCouples = theCoupleTable->GetTableSize(); 212 213 213 G4PhysicsLogVector* aVector = nullptr; 214 G4PhysicsLogVector* aVector = nullptr; 214 G4PhysicsLogVector* aVectorPrim = nullptr; 215 G4PhysicsLogVector* aVectorPrim = nullptr; 215 G4PhysicsLogVector* bVectorPrim = nullptr; 216 G4PhysicsLogVector* bVectorPrim = nullptr; 216 217 217 G4double emax1 = std::min(maxKinEnergy, minK 218 G4double emax1 = std::min(maxKinEnergy, minKinEnergyPrim); 218 219 219 for(std::size_t i=0; i<numOfCouples; ++i) { 220 for(std::size_t i=0; i<numOfCouples; ++i) { 220 221 221 if (bld->GetFlag(i)) { 222 if (bld->GetFlag(i)) { 222 // create physics vector and fill it 223 // create physics vector and fill it 223 const G4MaterialCutsCouple* couple = 224 const G4MaterialCutsCouple* couple = 224 theCoupleTable->GetMaterialCutsCouple( 225 theCoupleTable->GetMaterialCutsCouple((G4int)i); 225 226 226 // build main table 227 // build main table 227 if(nullptr != theLambdaTable) { 228 if(nullptr != theLambdaTable) { 228 delete (*theLambdaTable)[i]; 229 delete (*theLambdaTable)[i]; 229 230 230 // if start from zero then change the 231 // if start from zero then change the scale 231 G4double emin = minKinEnergy; 232 G4double emin = minKinEnergy; 232 G4bool startNull = false; 233 G4bool startNull = false; 233 if(startFromNull) { 234 if(startFromNull) { 234 G4double e = proc->MinPrimaryEnergy( 235 G4double e = proc->MinPrimaryEnergy(part, couple->GetMaterial()); 235 if(e >= emin) { 236 if(e >= emin) { 236 emin = e; 237 emin = e; 237 startNull = true; 238 startNull = true; 238 } 239 } 239 } 240 } 240 G4double emax = emax1; 241 G4double emax = emax1; 241 if(emax <= emin) { emax = 2*emin; } 242 if(emax <= emin) { emax = 2*emin; } 242 G4int bin = G4lrint(scale*G4Log(emax/e 243 G4int bin = G4lrint(scale*G4Log(emax/emin)); 243 bin = std::max(bin, 5); 244 bin = std::max(bin, 5); 244 aVector = new G4PhysicsLogVector(emin, 245 aVector = new G4PhysicsLogVector(emin, emax, bin, splineFlag); 245 modelManager->FillLambdaVector(aVector 246 modelManager->FillLambdaVector(aVector, couple, startNull); 246 if(splineFlag) { aVector->FillSecondDe 247 if(splineFlag) { aVector->FillSecondDerivatives(); } 247 G4PhysicsTableHelper::SetPhysicsVector 248 G4PhysicsTableHelper::SetPhysicsVector(theLambdaTable, i, aVector); 248 } 249 } 249 // build high energy table 250 // build high energy table 250 if(nullptr != theLambdaTablePrim && minK << 251 if(nullptr != theLambdaTablePrim) { 251 delete (*theLambdaTablePrim)[i]; 252 delete (*theLambdaTablePrim)[i]; 252 253 253 // start not from zero and always use 254 // start not from zero and always use spline 254 if(nullptr == bVectorPrim) { 255 if(nullptr == bVectorPrim) { 255 G4int bin = G4lrint(scale*G4Log(maxK 256 G4int bin = G4lrint(scale*G4Log(maxKinEnergy/minKinEnergyPrim)); 256 bin = std::max(bin, 5); 257 bin = std::max(bin, 5); 257 aVectorPrim = 258 aVectorPrim = 258 new G4PhysicsLogVector(minKinEnerg 259 new G4PhysicsLogVector(minKinEnergyPrim, maxKinEnergy, bin, true); 259 bVectorPrim = aVectorPrim; 260 bVectorPrim = aVectorPrim; 260 } else { 261 } else { 261 aVectorPrim = new G4PhysicsLogVector 262 aVectorPrim = new G4PhysicsLogVector(*bVectorPrim); 262 } 263 } 263 modelManager->FillLambdaVector(aVector 264 modelManager->FillLambdaVector(aVectorPrim, couple, false, 264 fIsCros 265 fIsCrossSectionPrim); 265 aVectorPrim->FillSecondDerivatives(); 266 aVectorPrim->FillSecondDerivatives(); 266 G4PhysicsTableHelper::SetPhysicsVector 267 G4PhysicsTableHelper::SetPhysicsVector(theLambdaTablePrim, i, 267 268 aVectorPrim); 268 } 269 } 269 } 270 } 270 } 271 } 271 272 272 if(1 < verboseLevel) { 273 if(1 < verboseLevel) { 273 G4cout << "Lambda table is built for " << 274 G4cout << "Lambda table is built for " << part->GetParticleName() << G4endl; 274 } 275 } 275 } 276 } 276 277 277 //....oooOO0OOooo........oooOO0OOooo........oo 278 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 278 279 279 void G4EmTableUtil::BuildLambdaTable(G4VEnerg 280 void G4EmTableUtil::BuildLambdaTable(G4VEnergyLossProcess* proc, 280 const G4P 281 const G4ParticleDefinition* part, 281 G4EmModel 282 G4EmModelManager* modelManager, 282 G4LossTableBuilder* bld, 283 G4LossTableBuilder* bld, 283 G4Physics 284 G4PhysicsTable* theLambdaTable, 284 const G4D 285 const G4DataVector* theCuts, 285 const G4d 286 const G4double minKinEnergy, 286 const G4d 287 const G4double maxKinEnergy, 287 const G4d 288 const G4double scale, 288 const G4int verboseLevel, 289 const G4int verboseLevel, 289 const G4b 290 const G4bool splineFlag) 290 { 291 { 291 if(1 < verboseLevel) { 292 if(1 < verboseLevel) { 292 G4cout << "G4EmTableUtil::BuildLambdaTable 293 G4cout << "G4EmTableUtil::BuildLambdaTable() for process " 293 << proc->GetProcessName() << " and 294 << proc->GetProcessName() << " and particle " 294 << part->GetParticleName() << G4end 295 << part->GetParticleName() << G4endl; 295 } 296 } 296 297 297 const G4ProductionCutsTable* theCoupleTable= 298 const G4ProductionCutsTable* theCoupleTable= 298 G4ProductionCutsTable::GetProductionCu 299 G4ProductionCutsTable::GetProductionCutsTable(); 299 std::size_t numOfCouples = theCoupleTable->G 300 std::size_t numOfCouples = theCoupleTable->GetTableSize(); 300 301 301 G4PhysicsLogVector* aVector = nullptr; 302 G4PhysicsLogVector* aVector = nullptr; 302 for(std::size_t i=0; i<numOfCouples; ++i) { 303 for(std::size_t i=0; i<numOfCouples; ++i) { 303 if (bld->GetFlag(i)) { 304 if (bld->GetFlag(i)) { 304 // create physics vector and fill it 305 // create physics vector and fill it 305 const G4MaterialCutsCouple* couple = 306 const G4MaterialCutsCouple* couple = 306 theCoupleTable->GetMaterialCutsCouple( 307 theCoupleTable->GetMaterialCutsCouple((G4int)i); 307 308 308 delete (*theLambdaTable)[i]; 309 delete (*theLambdaTable)[i]; 309 G4bool startNull = true; 310 G4bool startNull = true; 310 G4double emin = 311 G4double emin = 311 proc->MinPrimaryEnergy(part, couple->G 312 proc->MinPrimaryEnergy(part, couple->GetMaterial(), (*theCuts)[i]); 312 if(minKinEnergy > emin) { 313 if(minKinEnergy > emin) { 313 emin = minKinEnergy; 314 emin = minKinEnergy; 314 startNull = false; 315 startNull = false; 315 } 316 } 316 317 317 G4double emax = maxKinEnergy; 318 G4double emax = maxKinEnergy; 318 if(emax <= emin) { emax = 2*emin; } 319 if(emax <= emin) { emax = 2*emin; } 319 G4int bin = G4lrint(scale*G4Log(emax/emi 320 G4int bin = G4lrint(scale*G4Log(emax/emin)); 320 bin = std::max(bin, 5); 321 bin = std::max(bin, 5); 321 aVector = new G4PhysicsLogVector(emin, e 322 aVector = new G4PhysicsLogVector(emin, emax, bin, splineFlag); 322 modelManager->FillLambdaVector(aVector, 323 modelManager->FillLambdaVector(aVector, couple, startNull, fRestricted); 323 if(splineFlag) { aVector->FillSecondDeri 324 if(splineFlag) { aVector->FillSecondDerivatives(); } 324 G4PhysicsTableHelper::SetPhysicsVector(t 325 G4PhysicsTableHelper::SetPhysicsVector(theLambdaTable, i, aVector); 325 } 326 } 326 } 327 } 327 328 328 if(1 < verboseLevel) { 329 if(1 < verboseLevel) { 329 G4cout << "Lambda table is built for " << 330 G4cout << "Lambda table is built for " << part->GetParticleName() << G4endl; 330 } 331 } 331 } 332 } 332 333 333 //....oooOO0OOooo........oooOO0OOooo........oo 334 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 334 335 335 const G4ParticleDefinition* 336 const G4ParticleDefinition* 336 G4EmTableUtil::CheckIon(G4VEnergyLossProcess* 337 G4EmTableUtil::CheckIon(G4VEnergyLossProcess* proc, 337 const G4ParticleDefini 338 const G4ParticleDefinition* part, 338 const G4ParticleDefini 339 const G4ParticleDefinition* partLocal, 339 const G4int verb, G4bo 340 const G4int verb, G4bool& isIon) 340 { 341 { 341 if(1 < verb) { 342 if(1 < verb) { 342 G4cout << "G4EmTableUtil::CheckIon for " 343 G4cout << "G4EmTableUtil::CheckIon for " 343 << proc->GetProcessName() << " for 344 << proc->GetProcessName() << " for " << part->GetParticleName() 344 << " should be called from G4VEnerg 345 << " should be called from G4VEnergyLossProcess::PreparePhysicsTable" 345 << G4endl; 346 << G4endl; 346 } 347 } 347 const G4ParticleDefinition* particle = partL 348 const G4ParticleDefinition* particle = partLocal; 348 349 349 // Are particle defined? 350 // Are particle defined? 350 if(nullptr == particle) { particle = part; } 351 if(nullptr == particle) { particle = part; } 351 if(part->GetParticleType() == "nucleus") { 352 if(part->GetParticleType() == "nucleus") { 352 G4String pname = part->GetParticleName(); 353 G4String pname = part->GetParticleName(); 353 if(pname != "deuteron" && pname != "triton 354 if(pname != "deuteron" && pname != "triton" && 354 pname != "He3" && pname != "alpha+" && << 355 pname != "alpha+" && pname != "alpha") { 355 356 356 const G4ParticleDefinition* theGIon = G4 357 const G4ParticleDefinition* theGIon = G4GenericIon::GenericIon(); 357 isIon = true; 358 isIon = true; 358 359 359 // this is a loop to compare pointers of 360 // this is a loop to compare pointers of G4GenericIon processes in order 360 // to confirm that for given particle th 361 // to confirm that for given particle the G4GenericIon physics is used 361 if(particle != theGIon) { 362 if(particle != theGIon) { 362 G4ProcessManager* pm = theGIon->GetPro 363 G4ProcessManager* pm = theGIon->GetProcessManager(); 363 G4ProcessVector* v = pm->GetAlongStepP 364 G4ProcessVector* v = pm->GetAlongStepProcessVector(); 364 G4int n = (G4int)v->size(); 365 G4int n = (G4int)v->size(); 365 for(G4int j=0; j<n; ++j) { 366 for(G4int j=0; j<n; ++j) { 366 if((*v)[j] == proc) { 367 if((*v)[j] == proc) { 367 particle = theGIon; 368 particle = theGIon; 368 break; 369 break; 369 } 370 } 370 } 371 } 371 } 372 } 372 } 373 } 373 } 374 } 374 return particle; 375 return particle; 375 } 376 } 376 377 377 //....oooOO0OOooo........oooOO0OOooo........oo 378 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 378 379 379 void G4EmTableUtil::UpdateModels(G4VEnergyLoss 380 void G4EmTableUtil::UpdateModels(G4VEnergyLossProcess* proc, 380 G4EmModelManager* modelManager, 381 G4EmModelManager* modelManager, 381 const G4doubl 382 const G4double maxKinEnergy, 382 const G4int n 383 const G4int nModels, 383 G4int& secID, 384 G4int& secID, G4int& biasID, 384 G4int& mainSe 385 G4int& mainSec, const G4bool baseMat, 385 const G4bool, << 386 const G4bool isMaster, const G4bool useAGen) 386 { 387 { 387 // defined ID of secondary particles 388 // defined ID of secondary particles 388 G4int stype = proc->GetProcessSubType(); 389 G4int stype = proc->GetProcessSubType(); 389 if(stype == fBremsstrahlung) { 390 if(stype == fBremsstrahlung) { 390 secID = _Bremsstrahlung; 391 secID = _Bremsstrahlung; 391 biasID = _SplitBremsstrahlung; 392 biasID = _SplitBremsstrahlung; 392 } else if(stype == fPairProdByCharged) { 393 } else if(stype == fPairProdByCharged) { 393 secID = _PairProduction; 394 secID = _PairProduction; 394 mainSec = 2; 395 mainSec = 2; 395 } 396 } 396 397 397 // initialisation of models 398 // initialisation of models 398 for(G4int i=0; i<nModels; ++i) { 399 for(G4int i=0; i<nModels; ++i) { 399 G4VEmModel* mod = modelManager->GetModel(i 400 G4VEmModel* mod = modelManager->GetModel(i); >> 401 mod->SetMasterThread(isMaster); 400 mod->SetAngularGeneratorFlag(useAGen); 402 mod->SetAngularGeneratorFlag(useAGen); 401 if(mod->HighEnergyLimit() > maxKinEnergy) 403 if(mod->HighEnergyLimit() > maxKinEnergy) { 402 mod->SetHighEnergyLimit(maxKinEnergy); 404 mod->SetHighEnergyLimit(maxKinEnergy); 403 } 405 } 404 mod->SetUseBaseMaterials(baseMat); 406 mod->SetUseBaseMaterials(baseMat); 405 } 407 } 406 } 408 } 407 409 408 //....oooOO0OOooo........oooOO0OOooo........oo 410 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 409 411 410 void 412 void 411 G4EmTableUtil::BuildLocalElossProcess(G4VEnerg 413 G4EmTableUtil::BuildLocalElossProcess(G4VEnergyLossProcess* proc, 412 const G4VEnergyLossProcess* mast 414 const G4VEnergyLossProcess* masterProc, 413 const G4ParticleDefinition* part 415 const G4ParticleDefinition* part, 414 const G4int nModels) 416 const G4int nModels) 415 { 417 { 416 // copy table pointers from master thread 418 // copy table pointers from master thread 417 proc->SetDEDXTable(masterProc->DEDXTable(),f 419 proc->SetDEDXTable(masterProc->DEDXTable(),fRestricted); 418 proc->SetDEDXTable(masterProc->DEDXunRestric 420 proc->SetDEDXTable(masterProc->DEDXunRestrictedTable(),fTotal); 419 proc->SetDEDXTable(masterProc->IonisationTab 421 proc->SetDEDXTable(masterProc->IonisationTable(),fIsIonisation); 420 proc->SetRangeTableForLoss(masterProc->Range 422 proc->SetRangeTableForLoss(masterProc->RangeTableForLoss()); 421 proc->SetCSDARangeTable(masterProc->CSDARang 423 proc->SetCSDARangeTable(masterProc->CSDARangeTable()); 422 proc->SetInverseRangeTable(masterProc->Inver 424 proc->SetInverseRangeTable(masterProc->InverseRangeTable()); 423 proc->SetLambdaTable(masterProc->LambdaTable 425 proc->SetLambdaTable(masterProc->LambdaTable()); 424 proc->SetCrossSectionType(masterProc->CrossS 426 proc->SetCrossSectionType(masterProc->CrossSectionType()); 425 proc->SetEnergyOfCrossSectionMax(masterProc- 427 proc->SetEnergyOfCrossSectionMax(masterProc->EnergyOfCrossSectionMax()); 426 proc->SetTwoPeaksXS(masterProc->TwoPeaksXS() 428 proc->SetTwoPeaksXS(masterProc->TwoPeaksXS()); 427 proc->SetIonisation(masterProc->IsIonisation 429 proc->SetIonisation(masterProc->IsIonisationProcess()); 428 G4bool baseMat = masterProc->UseBaseMaterial 430 G4bool baseMat = masterProc->UseBaseMaterial(); 429 431 430 // local initialisation of models 432 // local initialisation of models 431 G4bool printing = true; 433 G4bool printing = true; 432 for(G4int i=0; i<nModels; ++i) { 434 for(G4int i=0; i<nModels; ++i) { 433 G4VEmModel* mod = proc->GetModelByIndex(i, 435 G4VEmModel* mod = proc->GetModelByIndex(i, printing); 434 G4VEmModel* mod0= masterProc->GetModelByIn 436 G4VEmModel* mod0= masterProc->GetModelByIndex(i, printing); 435 mod->SetUseBaseMaterials(baseMat); 437 mod->SetUseBaseMaterials(baseMat); 436 mod->InitialiseLocal(part, mod0); 438 mod->InitialiseLocal(part, mod0); 437 } 439 } 438 } 440 } 439 441 440 //....oooOO0OOooo........oooOO0OOooo........oo 442 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 441 443 442 void G4EmTableUtil::BuildDEDXTable(G4VEnergyLo 444 void G4EmTableUtil::BuildDEDXTable(G4VEnergyLossProcess* proc, 443 const G4ParticleDefinition* part, 445 const G4ParticleDefinition* part, 444 G4EmModelManager* modelManager, 446 G4EmModelManager* modelManager, 445 G4LossTableBuilder* bld, 447 G4LossTableBuilder* bld, 446 G4PhysicsTable* table, 448 G4PhysicsTable* table, 447 const G4double emin, 449 const G4double emin, 448 const G4double emax, 450 const G4double emax, 449 const G4int nbins, 451 const G4int nbins, 450 const G4int verbose, 452 const G4int verbose, 451 const G4EmTableType tType, 453 const G4EmTableType tType, 452 const G4bool spline) 454 const G4bool spline) 453 { 455 { 454 // Access to materials 456 // Access to materials 455 const G4ProductionCutsTable* theCoupleTable= 457 const G4ProductionCutsTable* theCoupleTable= 456 G4ProductionCutsTable::GetProductionCu 458 G4ProductionCutsTable::GetProductionCutsTable(); 457 std::size_t numOfCouples = theCoupleTable->G 459 std::size_t numOfCouples = theCoupleTable->GetTableSize(); 458 460 459 if(1 < verbose) { 461 if(1 < verbose) { 460 G4cout << numOfCouples << " couples" << " 462 G4cout << numOfCouples << " couples" << " minKinEnergy(MeV)= " << emin 461 << " maxKinEnergy(MeV)= " << emax < 463 << " maxKinEnergy(MeV)= " << emax << " nbins= " << nbins << G4endl; 462 } 464 } 463 G4PhysicsLogVector* aVector = nullptr; 465 G4PhysicsLogVector* aVector = nullptr; 464 G4PhysicsLogVector* bVector = nullptr; 466 G4PhysicsLogVector* bVector = nullptr; 465 467 466 for(std::size_t i=0; i<numOfCouples; ++i) { 468 for(std::size_t i=0; i<numOfCouples; ++i) { 467 469 468 if(1 < verbose) { 470 if(1 < verbose) { 469 G4cout << "G4EmTableUtil::BuildDEDXVecto 471 G4cout << "G4EmTableUtil::BuildDEDXVector idx= " << i 470 << " flagTable=" << table->GetFl 472 << " flagTable=" << table->GetFlag(i) 471 << " flagBuilder=" << bld->GetFla 473 << " flagBuilder=" << bld->GetFlag(i) << G4endl; 472 } 474 } 473 if(bld->GetFlag(i)) { 475 if(bld->GetFlag(i)) { 474 476 475 // create physics vector and fill it 477 // create physics vector and fill it 476 const G4MaterialCutsCouple* couple = 478 const G4MaterialCutsCouple* couple = 477 theCoupleTable->GetMaterialCutsCouple( 479 theCoupleTable->GetMaterialCutsCouple((G4int)i); 478 delete (*table)[i]; 480 delete (*table)[i]; 479 if(nullptr != bVector) { 481 if(nullptr != bVector) { 480 aVector = new G4PhysicsLogVector(*bVec 482 aVector = new G4PhysicsLogVector(*bVector); 481 } else { 483 } else { 482 bVector = new G4PhysicsLogVector(emin, 484 bVector = new G4PhysicsLogVector(emin, emax, nbins, spline); 483 aVector = bVector; 485 aVector = bVector; 484 } 486 } 485 487 486 modelManager->FillDEDXVector(aVector, co 488 modelManager->FillDEDXVector(aVector, couple, tType); 487 if(spline) { aVector->FillSecondDerivati 489 if(spline) { aVector->FillSecondDerivatives(); } 488 490 489 // Insert vector for this material into 491 // Insert vector for this material into the table 490 G4PhysicsTableHelper::SetPhysicsVector(t 492 G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector); 491 } 493 } 492 } 494 } 493 495 494 if(1 < verbose) { 496 if(1 < verbose) { 495 G4cout << "G4EmTableUtil::BuildDEDXTable() 497 G4cout << "G4EmTableUtil::BuildDEDXTable(): table is built for " 496 << part->GetParticleName() 498 << part->GetParticleName() 497 << " and process " << proc->GetProc 499 << " and process " << proc->GetProcessName() 498 << G4endl; 500 << G4endl; 499 if(2 < verbose) G4cout << (*table) << G4en 501 if(2 < verbose) G4cout << (*table) << G4endl; 500 } 502 } 501 } 503 } 502 504 503 //....oooOO0OOooo........oooOO0OOooo........oo 505 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 504 506 505 void G4EmTableUtil::PrepareMscProcess(G4VMulti 507 void G4EmTableUtil::PrepareMscProcess(G4VMultipleScattering* proc, 506 const G4ParticleDefinition& part 508 const G4ParticleDefinition& part, 507 G4EmModelManager* modelManager, 509 G4EmModelManager* modelManager, 508 G4MscStepLimitType& stepLimit, 510 G4MscStepLimitType& stepLimit, 509 G4double 511 G4double& facrange, 510 G4bool& latDisplacement, G4bool& 512 G4bool& latDisplacement, G4bool& master, 511 G4bool& isIon, G4bool& baseMat) 513 G4bool& isIon, G4bool& baseMat) 512 { 514 { 513 auto param = G4EmParameters::Instance(); 515 auto param = G4EmParameters::Instance(); 514 G4int verb = (master) ? param->Verbose() : p 516 G4int verb = (master) ? param->Verbose() : param->WorkerVerbose(); 515 proc->SetVerboseLevel(verb); 517 proc->SetVerboseLevel(verb); 516 518 517 if(part.GetPDGMass() > CLHEP::GeV || 519 if(part.GetPDGMass() > CLHEP::GeV || 518 part.GetParticleName() == "GenericIon") { 520 part.GetParticleName() == "GenericIon") { isIon = true; } 519 521 520 if(1 < verb) { 522 if(1 < verb) { 521 G4cout << "### G4EmTableUtil::PrepearPhysi 523 G4cout << "### G4EmTableUtil::PrepearPhysicsTable() for " 522 << proc->GetProcessName() 524 << proc->GetProcessName() 523 << " and particle " << part.GetPart 525 << " and particle " << part.GetParticleName() 524 << " isIon: " << isIon << " isMaste 526 << " isIon: " << isIon << " isMaster: " << master 525 << G4endl; 527 << G4endl; 526 } 528 } 527 529 528 // initialise process 530 // initialise process 529 proc->InitialiseProcess(&part); 531 proc->InitialiseProcess(&part); 530 532 531 // heavy particles 533 // heavy particles 532 if(part.GetPDGMass() > CLHEP::MeV) { 534 if(part.GetPDGMass() > CLHEP::MeV) { 533 stepLimit = param->MscMuHadStepLimitType() 535 stepLimit = param->MscMuHadStepLimitType(); 534 facrange = param->MscMuHadRangeFactor(); 536 facrange = param->MscMuHadRangeFactor(); 535 latDisplacement = param->MuHadLateralDispl 537 latDisplacement = param->MuHadLateralDisplacement(); 536 } else { 538 } else { 537 stepLimit = param->MscStepLimitType(); 539 stepLimit = param->MscStepLimitType(); 538 facrange = param->MscRangeFactor(); 540 facrange = param->MscRangeFactor(); 539 latDisplacement = param->LateralDisplaceme 541 latDisplacement = param->LateralDisplacement(); 540 } 542 } 541 543 542 // initialisation of models 544 // initialisation of models 543 auto numberOfModels = modelManager->NumberOf 545 auto numberOfModels = modelManager->NumberOfModels(); 544 for(G4int i=0; i<numberOfModels; ++i) { 546 for(G4int i=0; i<numberOfModels; ++i) { 545 G4VMscModel* msc = proc->GetModelByIndex(i 547 G4VMscModel* msc = proc->GetModelByIndex(i); 546 msc->SetIonisation(nullptr, &part); 548 msc->SetIonisation(nullptr, &part); >> 549 msc->SetMasterThread(master); 547 msc->SetPolarAngleLimit(param->MscThetaLim 550 msc->SetPolarAngleLimit(param->MscThetaLimit()); 548 G4double emax = std::min(msc->HighEnergyLi 551 G4double emax = std::min(msc->HighEnergyLimit(),param->MaxKinEnergy()); 549 msc->SetHighEnergyLimit(emax); 552 msc->SetHighEnergyLimit(emax); 550 msc->SetUseBaseMaterials(baseMat); 553 msc->SetUseBaseMaterials(baseMat); 551 } 554 } 552 modelManager->Initialise(&part, nullptr, ver 555 modelManager->Initialise(&part, nullptr, verb); 553 } 556 } 554 557 555 //....oooOO0OOooo........oooOO0OOooo........oo 558 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo..... 556 559 557 void G4EmTableUtil::BuildMscProcess(G4VMultipl 560 void G4EmTableUtil::BuildMscProcess(G4VMultipleScattering* proc, 558 const G4VM 561 const G4VMultipleScattering* masterProc, 559 const G4ParticleDefini 562 const G4ParticleDefinition& part, 560 const G4ParticleDefini 563 const G4ParticleDefinition* firstPart, 561 G4int nModels, G4bool master) 564 G4int nModels, G4bool master) 562 { 565 { 563 auto param = G4EmParameters::Instance(); 566 auto param = G4EmParameters::Instance(); 564 G4int verb = param->Verbose(); 567 G4int verb = param->Verbose(); 565 568 566 if (firstPart == &part) { << 569 if(!master && firstPart == &part) { 567 G4LossTableBuilder* bld = G4LossTableManag << 570 // initialisation of models 568 G4bool baseMat = bld->GetBaseMaterialFlag( << 571 G4bool baseMat = masterProc->UseBaseMaterial(); 569 if (master) { << 572 for(G4int i=0; i<nModels; ++i) { 570 for (G4int i=0; i<nModels; ++i) { << 573 G4VMscModel* msc = proc->GetModelByIndex(i); 571 G4VMscModel* msc = proc->GetModelByIndex(i); << 574 G4VMscModel* msc0 = masterProc->GetModelByIndex(i); 572 msc->SetUseBaseMaterials(baseMat); << 575 msc->SetUseBaseMaterials(baseMat); 573 // table is always built for low mass partic << 576 msc->SetCrossSectionTable(msc0->GetCrossSectionTable(), false); 574 if (part.GetParticleName() != "Generic << 577 msc->InitialiseLocal(&part, msc0); 575 (part.GetPDGMass() < CLHEP::GeV || msc-> << 576 G4double emin = << 577 std::max(msc->LowEnergyLimit(), msc->Low << 578 G4double emax = << 579 std::min(msc->HighEnergyLimit(), msc->Hi << 580 emin = std::max(emin, param->MinKinEnergy( << 581 emax = std::min(emax, param->MaxKinEnergy( << 582 if (emin < emax) { << 583 auto table = bld->BuildTableForModel(msc << 584 msc, &part, emin, emax, true); << 585 msc->SetCrossSectionTable(table, true); << 586 } << 587 } << 588 } << 589 } else { << 590 for (G4int i=0; i<nModels; ++i) { << 591 G4VMscModel* msc = proc->GetModelByIndex(i); << 592 G4VMscModel* msc0 = masterProc->GetModelByIn << 593 msc->SetUseBaseMaterials(baseMat); << 594 msc->SetCrossSectionTable(msc0->GetCrossSect << 595 msc->InitialiseLocal(&part, msc0); << 596 } << 597 } 578 } 598 } 579 } 599 if(!param->IsPrintLocked()) { 580 if(!param->IsPrintLocked()) { 600 const G4String& num = part.GetParticleName 581 const G4String& num = part.GetParticleName(); 601 582 602 // explicitly defined printout by particle 583 // explicitly defined printout by particle name 603 if(1 < verb || (0 < verb && (num == "e-" | 584 if(1 < verb || (0 < verb && (num == "e-" || 604 num == "e+" || num == "mu+" || 585 num == "e+" || num == "mu+" || 605 num == "mu-" || num == "proton"|| 586 num == "mu-" || num == "proton"|| 606 num == "pi+" || num == "pi-" || 587 num == "pi+" || num == "pi-" || 607 num == "kaon+" || num == "kaon-" || 588 num == "kaon+" || num == "kaon-" || 608 num == "alpha" || num == "anti_proton 589 num == "alpha" || num == "anti_proton" || 609 num == "GenericIon" || num == "alpha+ 590 num == "GenericIon" || num == "alpha+" || 610 num == "alpha" ))) { 591 num == "alpha" ))) { 611 proc->StreamInfo(G4cout, part); 592 proc->StreamInfo(G4cout, part); 612 } 593 } 613 } 594 } 614 if(1 < verb) { 595 if(1 < verb) { 615 G4cout << "### G4EmTableUtil::BuildPhysics 596 G4cout << "### G4EmTableUtil::BuildPhysicsTable() done for " 616 << proc->GetProcessName() 597 << proc->GetProcessName() 617 << " and particle " << part.GetParticleNa 598 << " and particle " << part.GetParticleName() << G4endl; 618 } 599 } 619 } 600 } 620 601 621 //....oooOO0OOooo........oooOO0OOooo........oo 602 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 622 603 623 G4bool G4EmTableUtil::StoreMscTable(G4VMultipl 604 G4bool G4EmTableUtil::StoreMscTable(G4VMultipleScattering* proc, 624 const G4Pa 605 const G4ParticleDefinition* part, 625 const G4St 606 const G4String& dir, 626 const G4in 607 const G4int nModels, const G4int verb, 627 const G4bool ascii) 608 const G4bool ascii) 628 { 609 { 629 G4bool ok = true; 610 G4bool ok = true; 630 for(G4int i=0; i<nModels; ++i) { 611 for(G4int i=0; i<nModels; ++i) { 631 G4VMscModel* msc = proc->GetModelByIndex(i 612 G4VMscModel* msc = proc->GetModelByIndex(i); 632 G4PhysicsTable* table = msc->GetCrossSecti 613 G4PhysicsTable* table = msc->GetCrossSectionTable(); 633 if (nullptr != table) { 614 if (nullptr != table) { 634 G4String ss = G4UIcommand::ConvertToStri 615 G4String ss = G4UIcommand::ConvertToString(i); 635 G4String name = 616 G4String name = 636 proc->GetPhysicsTableFileName(part, di 617 proc->GetPhysicsTableFileName(part, dir, "LambdaMod"+ss, ascii); 637 G4bool yes = table->StorePhysicsTable(na 618 G4bool yes = table->StorePhysicsTable(name,ascii); 638 619 639 if ( yes ) { 620 if ( yes ) { 640 if ( verb > 0 ) { 621 if ( verb > 0 ) { 641 G4cout << "Physics table are stored 622 G4cout << "Physics table are stored for " 642 << part->GetParticleName() 623 << part->GetParticleName() 643 << " and process " << proc->G 624 << " and process " << proc->GetProcessName() 644 << " with a name <" << name < 625 << " with a name <" << name << "> " << G4endl; 645 } 626 } 646 } else { 627 } else { 647 G4cout << "Fail to store Physics Table 628 G4cout << "Fail to store Physics Table for " 648 << part->GetParticleName() 629 << part->GetParticleName() 649 << " and process " << proc->Get 630 << " and process " << proc->GetProcessName() 650 << " in the directory <" << dir 631 << " in the directory <" << dir 651 << "> " << G4endl; 632 << "> " << G4endl; 652 ok = false; 633 ok = false; 653 } 634 } 654 } 635 } 655 } 636 } 656 return ok; 637 return ok; 657 } 638 } 658 639 659 //....oooOO0OOooo........oooOO0OOooo........oo 640 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 660 641 661 G4bool G4EmTableUtil::StoreTable(G4VProcess* p 642 G4bool G4EmTableUtil::StoreTable(G4VProcess* ptr, 662 const G4Parti 643 const G4ParticleDefinition* part, 663 G4PhysicsTabl 644 G4PhysicsTable* aTable, 664 const G4Strin 645 const G4String& dir, 665 const G4Strin 646 const G4String& tname, 666 const G4int v 647 const G4int verb, const G4bool ascii) 667 { 648 { 668 G4bool res = true; 649 G4bool res = true; 669 if (nullptr != aTable) { 650 if (nullptr != aTable) { 670 const G4String& name = 651 const G4String& name = 671 ptr->GetPhysicsTableFileName(part, dir, 652 ptr->GetPhysicsTableFileName(part, dir, tname, ascii); 672 if ( aTable->StorePhysicsTable(name, ascii 653 if ( aTable->StorePhysicsTable(name, ascii) ) { 673 if (1 < verb) G4cout << "Stored: " << na 654 if (1 < verb) G4cout << "Stored: " << name << G4endl; 674 } else { 655 } else { 675 res = false; 656 res = false; 676 G4cout << "G4EmTableUtil::StoreTable fai 657 G4cout << "G4EmTableUtil::StoreTable fail to store: " << name << G4endl; 677 } 658 } 678 } 659 } 679 return res; 660 return res; 680 } 661 } 681 662 682 //....oooOO0OOooo........oooOO0OOooo........oo 663 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 683 664 684 G4bool G4EmTableUtil::RetrieveTable(G4VProcess 665 G4bool G4EmTableUtil::RetrieveTable(G4VProcess* ptr, 685 const G4Pa 666 const G4ParticleDefinition* part, 686 G4PhysicsT 667 G4PhysicsTable* aTable, 687 const G4St 668 const G4String& dir, const G4String& tname, 688 const G4in 669 const G4int verb, const G4bool ascii, 689 const G4bo 670 const G4bool spline) 690 { 671 { 691 G4bool res = true; 672 G4bool res = true; 692 if (nullptr == aTable) { return res; } 673 if (nullptr == aTable) { return res; } 693 if (1 < verb) { << 674 if (0 < verb) { 694 G4cout << tname << " table for " << part-> 675 G4cout << tname << " table for " << part->GetParticleName() 695 << " will be retrieved " << G4endl; 676 << " will be retrieved " << G4endl; 696 } 677 } 697 const G4String& name = 678 const G4String& name = 698 ptr->GetPhysicsTableFileName(part, dir, tn 679 ptr->GetPhysicsTableFileName(part, dir, tname, ascii); 699 if(G4PhysicsTableHelper::RetrievePhysicsTabl 680 if(G4PhysicsTableHelper::RetrievePhysicsTable(aTable, name, ascii, spline)) { 700 if(spline) { 681 if(spline) { 701 for(auto & v : *aTable) { 682 for(auto & v : *aTable) { 702 if(nullptr != v) { v->FillSecondDerivatives( 683 if(nullptr != v) { v->FillSecondDerivatives(); } 703 } 684 } 704 } 685 } 705 if (0 < verb) { 686 if (0 < verb) { 706 G4cout << tname << " table for " << part 687 G4cout << tname << " table for " << part->GetParticleName() 707 << " is retrieved from <" << name << "> << 688 << " is Retrieved from <" << name << ">" 708 << G4endl; 689 << G4endl; 709 } 690 } 710 } else { 691 } else { 711 res = false; 692 res = false; 712 G4cout << "G4EmTableUtil::RetrieveTable fa 693 G4cout << "G4EmTableUtil::RetrieveTable fail to retrieve: " << tname 713 << " from " << name << " for " << p 694 << " from " << name << " for " << part->GetParticleName() << G4endl; 714 } 695 } 715 return res; 696 return res; 716 } 697 } 717 698 718 //....oooOO0OOooo........oooOO0OOooo........oo 699 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 719 700 720 701