| 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 << "### G4VEmProcess::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 << "### G4VEmProcess::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 << "### G4VEmProcess::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 << "G4EmProcess::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 << "G4EnergyLossProcess::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 << "G4VEnergyLossProcess::PreparePhysicsTable for "
343 << proc->GetProcessName() << " for 344 << proc->GetProcessName() << " for " << part->GetParticleName()
344 << " should be called from G4VEnerg <<
345 << G4endl; 345 << G4endl;
346 } 346 }
347 const G4ParticleDefinition* particle = partL 347 const G4ParticleDefinition* particle = partLocal;
348 348
349 // Are particle defined? 349 // Are particle defined?
350 if(nullptr == particle) { particle = part; } 350 if(nullptr == particle) { particle = part; }
351 if(part->GetParticleType() == "nucleus") { 351 if(part->GetParticleType() == "nucleus") {
352 G4String pname = part->GetParticleName(); 352 G4String pname = part->GetParticleName();
353 if(pname != "deuteron" && pname != "triton 353 if(pname != "deuteron" && pname != "triton" &&
354 pname != "He3" && pname != "alpha+" && << 354 pname != "alpha+" && pname != "alpha") {
355 355
356 const G4ParticleDefinition* theGIon = G4 356 const G4ParticleDefinition* theGIon = G4GenericIon::GenericIon();
357 isIon = true; 357 isIon = true;
358 358
359 // this is a loop to compare pointers of <<
360 // to confirm that for given particle th <<
361 if(particle != theGIon) { 359 if(particle != theGIon) {
362 G4ProcessManager* pm = theGIon->GetPro 360 G4ProcessManager* pm = theGIon->GetProcessManager();
363 G4ProcessVector* v = pm->GetAlongStepP 361 G4ProcessVector* v = pm->GetAlongStepProcessVector();
364 G4int n = (G4int)v->size(); 362 G4int n = (G4int)v->size();
365 for(G4int j=0; j<n; ++j) { 363 for(G4int j=0; j<n; ++j) {
366 if((*v)[j] == proc) { 364 if((*v)[j] == proc) {
367 particle = theGIon; 365 particle = theGIon;
368 break; 366 break;
369 } 367 }
370 } 368 }
371 } 369 }
372 } 370 }
373 } 371 }
374 return particle; 372 return particle;
375 } 373 }
376 374
377 //....oooOO0OOooo........oooOO0OOooo........oo 375 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
378 376
379 void G4EmTableUtil::UpdateModels(G4VEnergyLoss 377 void G4EmTableUtil::UpdateModels(G4VEnergyLossProcess* proc,
380 G4EmModelManager* modelManager, 378 G4EmModelManager* modelManager,
381 const G4doubl 379 const G4double maxKinEnergy,
382 const G4int n 380 const G4int nModels,
383 G4int& secID, 381 G4int& secID, G4int& biasID,
384 G4int& mainSe 382 G4int& mainSec, const G4bool baseMat,
385 const G4bool, << 383 const G4bool isMaster, const G4bool useAGen)
386 { 384 {
387 // defined ID of secondary particles 385 // defined ID of secondary particles
388 G4int stype = proc->GetProcessSubType(); 386 G4int stype = proc->GetProcessSubType();
389 if(stype == fBremsstrahlung) { 387 if(stype == fBremsstrahlung) {
390 secID = _Bremsstrahlung; 388 secID = _Bremsstrahlung;
391 biasID = _SplitBremsstrahlung; 389 biasID = _SplitBremsstrahlung;
392 } else if(stype == fPairProdByCharged) { 390 } else if(stype == fPairProdByCharged) {
393 secID = _PairProduction; 391 secID = _PairProduction;
394 mainSec = 2; 392 mainSec = 2;
395 } 393 }
396 394
397 // initialisation of models 395 // initialisation of models
398 for(G4int i=0; i<nModels; ++i) { 396 for(G4int i=0; i<nModels; ++i) {
399 G4VEmModel* mod = modelManager->GetModel(i 397 G4VEmModel* mod = modelManager->GetModel(i);
>> 398 mod->SetMasterThread(isMaster);
400 mod->SetAngularGeneratorFlag(useAGen); 399 mod->SetAngularGeneratorFlag(useAGen);
401 if(mod->HighEnergyLimit() > maxKinEnergy) 400 if(mod->HighEnergyLimit() > maxKinEnergy) {
402 mod->SetHighEnergyLimit(maxKinEnergy); 401 mod->SetHighEnergyLimit(maxKinEnergy);
403 } 402 }
404 mod->SetUseBaseMaterials(baseMat); 403 mod->SetUseBaseMaterials(baseMat);
405 } 404 }
406 } 405 }
407 406
408 //....oooOO0OOooo........oooOO0OOooo........oo 407 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
409 408
410 void 409 void
411 G4EmTableUtil::BuildLocalElossProcess(G4VEnerg 410 G4EmTableUtil::BuildLocalElossProcess(G4VEnergyLossProcess* proc,
412 const G4VEnergyLossProcess* mast 411 const G4VEnergyLossProcess* masterProc,
413 const G4ParticleDefinition* part 412 const G4ParticleDefinition* part,
414 const G4int nModels) 413 const G4int nModels)
415 { 414 {
416 // copy table pointers from master thread 415 // copy table pointers from master thread
417 proc->SetDEDXTable(masterProc->DEDXTable(),f 416 proc->SetDEDXTable(masterProc->DEDXTable(),fRestricted);
418 proc->SetDEDXTable(masterProc->DEDXunRestric 417 proc->SetDEDXTable(masterProc->DEDXunRestrictedTable(),fTotal);
419 proc->SetDEDXTable(masterProc->IonisationTab 418 proc->SetDEDXTable(masterProc->IonisationTable(),fIsIonisation);
420 proc->SetRangeTableForLoss(masterProc->Range 419 proc->SetRangeTableForLoss(masterProc->RangeTableForLoss());
421 proc->SetCSDARangeTable(masterProc->CSDARang 420 proc->SetCSDARangeTable(masterProc->CSDARangeTable());
422 proc->SetInverseRangeTable(masterProc->Inver 421 proc->SetInverseRangeTable(masterProc->InverseRangeTable());
423 proc->SetLambdaTable(masterProc->LambdaTable 422 proc->SetLambdaTable(masterProc->LambdaTable());
424 proc->SetCrossSectionType(masterProc->CrossS 423 proc->SetCrossSectionType(masterProc->CrossSectionType());
425 proc->SetEnergyOfCrossSectionMax(masterProc- 424 proc->SetEnergyOfCrossSectionMax(masterProc->EnergyOfCrossSectionMax());
426 proc->SetTwoPeaksXS(masterProc->TwoPeaksXS() 425 proc->SetTwoPeaksXS(masterProc->TwoPeaksXS());
427 proc->SetIonisation(masterProc->IsIonisation 426 proc->SetIonisation(masterProc->IsIonisationProcess());
428 G4bool baseMat = masterProc->UseBaseMaterial 427 G4bool baseMat = masterProc->UseBaseMaterial();
429 428
430 // local initialisation of models 429 // local initialisation of models
431 G4bool printing = true; 430 G4bool printing = true;
432 for(G4int i=0; i<nModels; ++i) { 431 for(G4int i=0; i<nModels; ++i) {
433 G4VEmModel* mod = proc->GetModelByIndex(i, 432 G4VEmModel* mod = proc->GetModelByIndex(i, printing);
434 G4VEmModel* mod0= masterProc->GetModelByIn 433 G4VEmModel* mod0= masterProc->GetModelByIndex(i, printing);
435 mod->SetUseBaseMaterials(baseMat); 434 mod->SetUseBaseMaterials(baseMat);
436 mod->InitialiseLocal(part, mod0); 435 mod->InitialiseLocal(part, mod0);
437 } 436 }
438 } 437 }
439 438
440 //....oooOO0OOooo........oooOO0OOooo........oo 439 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
441 440
442 void G4EmTableUtil::BuildDEDXTable(G4VEnergyLo 441 void G4EmTableUtil::BuildDEDXTable(G4VEnergyLossProcess* proc,
443 const G4ParticleDefinition* part, 442 const G4ParticleDefinition* part,
444 G4EmModelManager* modelManager, 443 G4EmModelManager* modelManager,
445 G4LossTableBuilder* bld, 444 G4LossTableBuilder* bld,
446 G4PhysicsTable* table, 445 G4PhysicsTable* table,
447 const G4double emin, 446 const G4double emin,
448 const G4double emax, 447 const G4double emax,
449 const G4int nbins, 448 const G4int nbins,
450 const G4int verbose, 449 const G4int verbose,
451 const G4EmTableType tType, 450 const G4EmTableType tType,
452 const G4bool spline) 451 const G4bool spline)
453 { 452 {
454 // Access to materials 453 // Access to materials
455 const G4ProductionCutsTable* theCoupleTable= 454 const G4ProductionCutsTable* theCoupleTable=
456 G4ProductionCutsTable::GetProductionCu 455 G4ProductionCutsTable::GetProductionCutsTable();
457 std::size_t numOfCouples = theCoupleTable->G 456 std::size_t numOfCouples = theCoupleTable->GetTableSize();
458 457
459 if(1 < verbose) { 458 if(1 < verbose) {
460 G4cout << numOfCouples << " couples" << " 459 G4cout << numOfCouples << " couples" << " minKinEnergy(MeV)= " << emin
461 << " maxKinEnergy(MeV)= " << emax < 460 << " maxKinEnergy(MeV)= " << emax << " nbins= " << nbins << G4endl;
462 } 461 }
463 G4PhysicsLogVector* aVector = nullptr; 462 G4PhysicsLogVector* aVector = nullptr;
464 G4PhysicsLogVector* bVector = nullptr; 463 G4PhysicsLogVector* bVector = nullptr;
465 464
466 for(std::size_t i=0; i<numOfCouples; ++i) { 465 for(std::size_t i=0; i<numOfCouples; ++i) {
467 466
468 if(1 < verbose) { 467 if(1 < verbose) {
469 G4cout << "G4EmTableUtil::BuildDEDXVecto << 468 G4cout << "G4VEnergyLossProcess::BuildDEDXVector idx= " << i
470 << " flagTable=" << table->GetFl 469 << " flagTable=" << table->GetFlag(i)
471 << " flagBuilder=" << bld->GetFla 470 << " flagBuilder=" << bld->GetFlag(i) << G4endl;
472 } 471 }
473 if(bld->GetFlag(i)) { 472 if(bld->GetFlag(i)) {
474 473
475 // create physics vector and fill it 474 // create physics vector and fill it
476 const G4MaterialCutsCouple* couple = 475 const G4MaterialCutsCouple* couple =
477 theCoupleTable->GetMaterialCutsCouple( 476 theCoupleTable->GetMaterialCutsCouple((G4int)i);
478 delete (*table)[i]; 477 delete (*table)[i];
479 if(nullptr != bVector) { 478 if(nullptr != bVector) {
480 aVector = new G4PhysicsLogVector(*bVec 479 aVector = new G4PhysicsLogVector(*bVector);
481 } else { 480 } else {
482 bVector = new G4PhysicsLogVector(emin, 481 bVector = new G4PhysicsLogVector(emin, emax, nbins, spline);
483 aVector = bVector; 482 aVector = bVector;
484 } 483 }
485 484
486 modelManager->FillDEDXVector(aVector, co 485 modelManager->FillDEDXVector(aVector, couple, tType);
487 if(spline) { aVector->FillSecondDerivati 486 if(spline) { aVector->FillSecondDerivatives(); }
488 487
489 // Insert vector for this material into 488 // Insert vector for this material into the table
490 G4PhysicsTableHelper::SetPhysicsVector(t 489 G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector);
491 } 490 }
492 } 491 }
493 492
494 if(1 < verbose) { 493 if(1 < verbose) {
495 G4cout << "G4EmTableUtil::BuildDEDXTable() << 494 G4cout << "G4VEnergyLossProcess::BuildDEDXTable(): table is built for "
496 << part->GetParticleName() 495 << part->GetParticleName()
497 << " and process " << proc->GetProc 496 << " and process " << proc->GetProcessName()
498 << G4endl; 497 << G4endl;
499 if(2 < verbose) G4cout << (*table) << G4en 498 if(2 < verbose) G4cout << (*table) << G4endl;
500 } 499 }
501 } 500 }
502 501
503 //....oooOO0OOooo........oooOO0OOooo........oo 502 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
504 503
505 void G4EmTableUtil::PrepareMscProcess(G4VMulti 504 void G4EmTableUtil::PrepareMscProcess(G4VMultipleScattering* proc,
506 const G4ParticleDefinition& part 505 const G4ParticleDefinition& part,
507 G4EmModelManager* modelManager, 506 G4EmModelManager* modelManager,
508 G4MscStepLimitType& stepLimit, 507 G4MscStepLimitType& stepLimit,
509 G4double 508 G4double& facrange,
510 G4bool& latDisplacement, G4bool& 509 G4bool& latDisplacement, G4bool& master,
511 G4bool& isIon, G4bool& baseMat) 510 G4bool& isIon, G4bool& baseMat)
512 { 511 {
513 auto param = G4EmParameters::Instance(); 512 auto param = G4EmParameters::Instance();
514 G4int verb = (master) ? param->Verbose() : p 513 G4int verb = (master) ? param->Verbose() : param->WorkerVerbose();
515 proc->SetVerboseLevel(verb); 514 proc->SetVerboseLevel(verb);
516 515
517 if(part.GetPDGMass() > CLHEP::GeV || 516 if(part.GetPDGMass() > CLHEP::GeV ||
518 part.GetParticleName() == "GenericIon") { 517 part.GetParticleName() == "GenericIon") { isIon = true; }
519 518
520 if(1 < verb) { 519 if(1 < verb) {
521 G4cout << "### G4EmTableUtil::PrepearPhysi << 520 G4cout << "### G4VMultipleScattering::PrepearPhysicsTable() for "
522 << proc->GetProcessName() 521 << proc->GetProcessName()
523 << " and particle " << part.GetPart 522 << " and particle " << part.GetParticleName()
524 << " isIon: " << isIon << " isMaste 523 << " isIon: " << isIon << " isMaster: " << master
525 << G4endl; 524 << G4endl;
526 } 525 }
527 526
528 // initialise process 527 // initialise process
529 proc->InitialiseProcess(&part); 528 proc->InitialiseProcess(&part);
530 529
531 // heavy particles 530 // heavy particles
532 if(part.GetPDGMass() > CLHEP::MeV) { 531 if(part.GetPDGMass() > CLHEP::MeV) {
533 stepLimit = param->MscMuHadStepLimitType() 532 stepLimit = param->MscMuHadStepLimitType();
534 facrange = param->MscMuHadRangeFactor(); 533 facrange = param->MscMuHadRangeFactor();
535 latDisplacement = param->MuHadLateralDispl 534 latDisplacement = param->MuHadLateralDisplacement();
536 } else { 535 } else {
537 stepLimit = param->MscStepLimitType(); 536 stepLimit = param->MscStepLimitType();
538 facrange = param->MscRangeFactor(); 537 facrange = param->MscRangeFactor();
539 latDisplacement = param->LateralDisplaceme 538 latDisplacement = param->LateralDisplacement();
540 } 539 }
541 540
542 // initialisation of models 541 // initialisation of models
543 auto numberOfModels = modelManager->NumberOf 542 auto numberOfModels = modelManager->NumberOfModels();
544 for(G4int i=0; i<numberOfModels; ++i) { 543 for(G4int i=0; i<numberOfModels; ++i) {
545 G4VMscModel* msc = proc->GetModelByIndex(i 544 G4VMscModel* msc = proc->GetModelByIndex(i);
546 msc->SetIonisation(nullptr, &part); 545 msc->SetIonisation(nullptr, &part);
>> 546 msc->SetMasterThread(master);
547 msc->SetPolarAngleLimit(param->MscThetaLim 547 msc->SetPolarAngleLimit(param->MscThetaLimit());
548 G4double emax = std::min(msc->HighEnergyLi 548 G4double emax = std::min(msc->HighEnergyLimit(),param->MaxKinEnergy());
549 msc->SetHighEnergyLimit(emax); 549 msc->SetHighEnergyLimit(emax);
550 msc->SetUseBaseMaterials(baseMat); 550 msc->SetUseBaseMaterials(baseMat);
551 } 551 }
552 modelManager->Initialise(&part, nullptr, ver 552 modelManager->Initialise(&part, nullptr, verb);
553 } 553 }
554 554
555 //....oooOO0OOooo........oooOO0OOooo........oo 555 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
556 556
557 void G4EmTableUtil::BuildMscProcess(G4VMultipl 557 void G4EmTableUtil::BuildMscProcess(G4VMultipleScattering* proc,
558 const G4VM 558 const G4VMultipleScattering* masterProc,
559 const G4ParticleDefini 559 const G4ParticleDefinition& part,
560 const G4ParticleDefini 560 const G4ParticleDefinition* firstPart,
561 G4int nModels, G4bool master) 561 G4int nModels, G4bool master)
562 { 562 {
563 auto param = G4EmParameters::Instance(); 563 auto param = G4EmParameters::Instance();
564 G4int verb = param->Verbose(); 564 G4int verb = param->Verbose();
565 565
566 if (firstPart == &part) { << 566 if(!master && firstPart == &part) {
567 G4LossTableBuilder* bld = G4LossTableManag << 567 // initialisation of models
568 G4bool baseMat = bld->GetBaseMaterialFlag( << 568 G4bool baseMat = masterProc->UseBaseMaterial();
569 if (master) { << 569 for(G4int i=0; i<nModels; ++i) {
570 for (G4int i=0; i<nModels; ++i) { << 570 G4VMscModel* msc = proc->GetModelByIndex(i);
571 G4VMscModel* msc = proc->GetModelByIndex(i); << 571 G4VMscModel* msc0 = masterProc->GetModelByIndex(i);
572 msc->SetUseBaseMaterials(baseMat); << 572 msc->SetUseBaseMaterials(baseMat);
573 // table is always built for low mass partic << 573 msc->SetCrossSectionTable(msc0->GetCrossSectionTable(), false);
574 if (part.GetParticleName() != "Generic << 574 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 } 575 }
598 } 576 }
599 if(!param->IsPrintLocked()) { 577 if(!param->IsPrintLocked()) {
600 const G4String& num = part.GetParticleName 578 const G4String& num = part.GetParticleName();
601 579
602 // explicitly defined printout by particle 580 // explicitly defined printout by particle name
603 if(1 < verb || (0 < verb && (num == "e-" | 581 if(1 < verb || (0 < verb && (num == "e-" ||
604 num == "e+" || num == "mu+" || 582 num == "e+" || num == "mu+" ||
605 num == "mu-" || num == "proton"|| 583 num == "mu-" || num == "proton"||
606 num == "pi+" || num == "pi-" || 584 num == "pi+" || num == "pi-" ||
607 num == "kaon+" || num == "kaon-" || 585 num == "kaon+" || num == "kaon-" ||
608 num == "alpha" || num == "anti_proton 586 num == "alpha" || num == "anti_proton" ||
609 num == "GenericIon" || num == "alpha+ 587 num == "GenericIon" || num == "alpha+" ||
610 num == "alpha" ))) { 588 num == "alpha" ))) {
611 proc->StreamInfo(G4cout, part); 589 proc->StreamInfo(G4cout, part);
612 } 590 }
613 } 591 }
614 if(1 < verb) { 592 if(1 < verb) {
615 G4cout << "### G4EmTableUtil::BuildPhysics << 593 G4cout << "### G4VMultipleScattering::BuildPhysicsTable() done for "
616 << proc->GetProcessName() 594 << proc->GetProcessName()
617 << " and particle " << part.GetParticleNa 595 << " and particle " << part.GetParticleName() << G4endl;
618 } 596 }
619 } 597 }
620 598
621 //....oooOO0OOooo........oooOO0OOooo........oo 599 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
622 600
623 G4bool G4EmTableUtil::StoreMscTable(G4VMultipl 601 G4bool G4EmTableUtil::StoreMscTable(G4VMultipleScattering* proc,
624 const G4Pa 602 const G4ParticleDefinition* part,
625 const G4St 603 const G4String& dir,
626 const G4in 604 const G4int nModels, const G4int verb,
627 const G4bool ascii) 605 const G4bool ascii)
628 { 606 {
629 G4bool ok = true; 607 G4bool ok = true;
630 for(G4int i=0; i<nModels; ++i) { 608 for(G4int i=0; i<nModels; ++i) {
631 G4VMscModel* msc = proc->GetModelByIndex(i 609 G4VMscModel* msc = proc->GetModelByIndex(i);
632 G4PhysicsTable* table = msc->GetCrossSecti 610 G4PhysicsTable* table = msc->GetCrossSectionTable();
633 if (nullptr != table) { 611 if (nullptr != table) {
634 G4String ss = G4UIcommand::ConvertToStri 612 G4String ss = G4UIcommand::ConvertToString(i);
635 G4String name = 613 G4String name =
636 proc->GetPhysicsTableFileName(part, di 614 proc->GetPhysicsTableFileName(part, dir, "LambdaMod"+ss, ascii);
637 G4bool yes = table->StorePhysicsTable(na 615 G4bool yes = table->StorePhysicsTable(name,ascii);
638 616
639 if ( yes ) { 617 if ( yes ) {
640 if ( verb > 0 ) { 618 if ( verb > 0 ) {
641 G4cout << "Physics table are stored 619 G4cout << "Physics table are stored for "
642 << part->GetParticleName() 620 << part->GetParticleName()
643 << " and process " << proc->G 621 << " and process " << proc->GetProcessName()
644 << " with a name <" << name < 622 << " with a name <" << name << "> " << G4endl;
645 } 623 }
646 } else { 624 } else {
647 G4cout << "Fail to store Physics Table 625 G4cout << "Fail to store Physics Table for "
648 << part->GetParticleName() 626 << part->GetParticleName()
649 << " and process " << proc->Get 627 << " and process " << proc->GetProcessName()
650 << " in the directory <" << dir 628 << " in the directory <" << dir
651 << "> " << G4endl; 629 << "> " << G4endl;
652 ok = false; 630 ok = false;
653 } 631 }
654 } 632 }
655 } 633 }
656 return ok; 634 return ok;
657 } 635 }
658 636
659 //....oooOO0OOooo........oooOO0OOooo........oo 637 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
660 638
661 G4bool G4EmTableUtil::StoreTable(G4VProcess* p 639 G4bool G4EmTableUtil::StoreTable(G4VProcess* ptr,
662 const G4Parti 640 const G4ParticleDefinition* part,
663 G4PhysicsTabl 641 G4PhysicsTable* aTable,
664 const G4Strin 642 const G4String& dir,
665 const G4Strin 643 const G4String& tname,
666 const G4int v 644 const G4int verb, const G4bool ascii)
667 { 645 {
668 G4bool res = true; 646 G4bool res = true;
669 if (nullptr != aTable) { 647 if (nullptr != aTable) {
670 const G4String& name = 648 const G4String& name =
671 ptr->GetPhysicsTableFileName(part, dir, 649 ptr->GetPhysicsTableFileName(part, dir, tname, ascii);
672 if ( aTable->StorePhysicsTable(name, ascii 650 if ( aTable->StorePhysicsTable(name, ascii) ) {
673 if (1 < verb) G4cout << "Stored: " << na 651 if (1 < verb) G4cout << "Stored: " << name << G4endl;
674 } else { 652 } else {
675 res = false; 653 res = false;
676 G4cout << "G4EmTableUtil::StoreTable fai << 654 G4cout << "Fail to store: " << name << G4endl;
677 } 655 }
678 } 656 }
679 return res; 657 return res;
680 } 658 }
681 659
682 //....oooOO0OOooo........oooOO0OOooo........oo 660 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
683 661
684 G4bool G4EmTableUtil::RetrieveTable(G4VProcess 662 G4bool G4EmTableUtil::RetrieveTable(G4VProcess* ptr,
685 const G4Pa 663 const G4ParticleDefinition* part,
686 G4PhysicsT 664 G4PhysicsTable* aTable,
687 const G4St 665 const G4String& dir, const G4String& tname,
688 const G4in 666 const G4int verb, const G4bool ascii,
689 const G4bo 667 const G4bool spline)
690 { 668 {
691 G4bool res = true; 669 G4bool res = true;
692 if (nullptr == aTable) { return res; } 670 if (nullptr == aTable) { return res; }
693 if (1 < verb) { << 671 G4cout << tname << " table for " << part->GetParticleName()
694 G4cout << tname << " table for " << part-> << 672 << " will be retrieved " << G4endl;
695 << " will be retrieved " << G4endl; <<
696 } <<
697 const G4String& name = 673 const G4String& name =
698 ptr->GetPhysicsTableFileName(part, dir, tn 674 ptr->GetPhysicsTableFileName(part, dir, tname, ascii);
699 if(G4PhysicsTableHelper::RetrievePhysicsTabl 675 if(G4PhysicsTableHelper::RetrievePhysicsTable(aTable, name, ascii, spline)) {
700 if(spline) { 676 if(spline) {
701 for(auto & v : *aTable) { 677 for(auto & v : *aTable) {
702 if(nullptr != v) { v->FillSecondDerivatives( 678 if(nullptr != v) { v->FillSecondDerivatives(); }
703 } 679 }
704 } 680 }
705 if (0 < verb) { 681 if (0 < verb) {
706 G4cout << tname << " table for " << part 682 G4cout << tname << " table for " << part->GetParticleName()
707 << " is retrieved from <" << name << "> << 683 << " is Retrieved from <" << name << ">"
708 << G4endl; 684 << G4endl;
709 } 685 }
710 } else { 686 } else {
711 res = false; 687 res = false;
712 G4cout << "G4EmTableUtil::RetrieveTable fa << 688 G4cout << "Fail to retrieve: " << tname << " from " << name << " for "
713 << " from " << name << " for " << p << 689 << part->GetParticleName() << G4endl;
714 } 690 }
715 return res; 691 return res;
716 } 692 }
717 693
718 //....oooOO0OOooo........oooOO0OOooo........oo 694 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
719 695
720 696