Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4LossTableBuilder Class Reference

#include <G4LossTableBuilder.hh>

Public Member Functions

 G4LossTableBuilder (G4bool master=true)
 
 ~G4LossTableBuilder ()
 
void BuildDEDXTable (G4PhysicsTable *dedxTable, const std::vector< G4PhysicsTable * > &)
 
void BuildRangeTable (const G4PhysicsTable *dedxTable, G4PhysicsTable *rangeTable)
 
void BuildInverseRangeTable (const G4PhysicsTable *rangeTable, G4PhysicsTable *invRangeTable)
 
G4PhysicsTableBuildTableForModel (G4PhysicsTable *table, G4VEmModel *model, const G4ParticleDefinition *, G4double emin, G4double emax, G4bool spline)
 
void InitialiseBaseMaterials (const G4PhysicsTable *table=nullptr)
 
const std::vector< G4int > * GetCoupleIndexes () const
 
const std::vector< G4double > * GetDensityFactors () const
 
G4bool GetFlag (size_t idx)
 
G4bool GetBaseMaterialFlag ()
 
void SetSplineFlag (G4bool flag)
 
void SetInitialisationFlag (G4bool flag)
 
void SetBaseMaterialActive (G4bool flag)
 
G4LossTableBuilderoperator= (const G4LossTableBuilder &right)=delete
 
 G4LossTableBuilder (const G4LossTableBuilder &)=delete
 

Detailed Description

Definition at line 60 of file G4LossTableBuilder.hh.

Constructor & Destructor Documentation

◆ G4LossTableBuilder() [1/2]

G4LossTableBuilder::G4LossTableBuilder ( G4bool  master = true)

Definition at line 78 of file G4LossTableBuilder.cc.

78 : isMaster(master)
79{
80 theParameters = G4EmParameters::Instance();
81 if(nullptr == theFlag) {
82#ifdef G4MULTITHREADED
83 G4MUTEXLOCK(&ltbMutex);
84 if(nullptr == theFlag) {
85#endif
86 if(!isMaster) {
88 ed << "Initialisation called from a worker thread ";
89 G4Exception("G4LossTableBuilder: ", "em0001",
90 JustWarning, ed);
91 }
92 theDensityFactor = new std::vector<G4double>;
93 theDensityIdx = new std::vector<G4int>;
94 theFlag = new std::vector<G4bool>;
95#ifdef G4MULTITHREADED
96 }
97 G4MUTEXUNLOCK(&ltbMutex);
98#endif
99 }
100}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
#define G4MUTEXLOCK(mutex)
Definition: G4Threading.hh:251
#define G4MUTEXUNLOCK(mutex)
Definition: G4Threading.hh:254
static G4EmParameters * Instance()

◆ ~G4LossTableBuilder()

G4LossTableBuilder::~G4LossTableBuilder ( )

Definition at line 104 of file G4LossTableBuilder.cc.

105{
106 if(isMaster) {
107 delete theDensityFactor;
108 delete theDensityIdx;
109 delete theFlag;
110 theDensityFactor = nullptr;
111 theDensityIdx = nullptr;
112 theFlag = nullptr;
113 }
114}

◆ G4LossTableBuilder() [2/2]

G4LossTableBuilder::G4LossTableBuilder ( const G4LossTableBuilder )
delete

Member Function Documentation

◆ BuildDEDXTable()

void G4LossTableBuilder::BuildDEDXTable ( G4PhysicsTable dedxTable,
const std::vector< G4PhysicsTable * > &  list 
)

Definition at line 149 of file G4LossTableBuilder.cc.

151{
152 InitialiseBaseMaterials(dedxTable);
153 std::size_t n_processes = list.size();
154 if(1 >= n_processes) { return; }
155
156 std::size_t nCouples = dedxTable->size();
157 //G4cout << "Nproc= " << n_processes << " nCouples=" << nCouples << " Nv= "
158 // << dedxTable->size() << G4endl;
159 if(0 >= nCouples) { return; }
160
161 for (std::size_t i=0; i<nCouples; ++i) {
162 auto pv0 = static_cast<G4PhysicsLogVector*>((*(list[0]))[i]);
163 if(pv0 == nullptr) { continue; }
164 std::size_t npoints = pv0->GetVectorLength();
165 auto pv = new G4PhysicsLogVector(*pv0);
166 for (std::size_t j=0; j<npoints; ++j) {
167 G4double dedx = 0.0;
168 for (std::size_t k=0; k<n_processes; ++k) {
169 const G4PhysicsVector* pv1 = (*(list[k]))[i];
170 dedx += (*pv1)[j];
171 }
172 pv->PutValue(j, dedx);
173 }
174 if(splineFlag) { pv->FillSecondDerivatives(); }
176 }
177 //G4cout << "### G4LossTableBuilder::BuildDEDXTable " << G4endl;
178 //G4cout << *dedxTable << G4endl;
179}
double G4double
Definition: G4Types.hh:83
void InitialiseBaseMaterials(const G4PhysicsTable *table=nullptr)
static void SetPhysicsVector(G4PhysicsTable *physTable, std::size_t idx, G4PhysicsVector *vec)
void PutValue(const std::size_t index, const G4double value)
std::size_t GetVectorLength() const

◆ BuildInverseRangeTable()

void G4LossTableBuilder::BuildInverseRangeTable ( const G4PhysicsTable rangeTable,
G4PhysicsTable invRangeTable 
)

Definition at line 263 of file G4LossTableBuilder.cc.

266{
267 std::size_t nCouples = rangeTable->size();
268 if(0 >= nCouples) { return; }
269
270 for (std::size_t i=0; i<nCouples; ++i) {
271 G4PhysicsVector* pv = (*rangeTable)[i];
272 if((pv == nullptr) || (isBaseMatActive && !(*theFlag)[i])) { continue; }
273 std::size_t npoints = pv->GetVectorLength();
274
275 delete (*invRangeTable)[i];
276 auto v = new G4PhysicsFreeVector(npoints,splineFlag);
277
278 for (std::size_t j=0; j<npoints; ++j) {
279 G4double e = pv->Energy(j);
280 G4double r = (*pv)[j];
281 v->PutValues(j,r,e);
282 }
283 if(splineFlag) { v->FillSecondDerivatives(); }
284
285 G4PhysicsTableHelper::SetPhysicsVector(invRangeTable, i, v);
286 }
287 //G4cout << "### Inverse range table" << G4endl;
288 //G4cout << *invRangeTable << G4endl;
289}
G4double Energy(const std::size_t index) const

Referenced by G4TablesForExtrapolator::Initialisation().

◆ BuildRangeTable()

void G4LossTableBuilder::BuildRangeTable ( const G4PhysicsTable dedxTable,
G4PhysicsTable rangeTable 
)

Definition at line 183 of file G4LossTableBuilder.cc.

186{
187 //G4cout << "### G4LossTableBuilder::BuildRangeTable: DEDX table" << G4endl;
188 //G4cout << *const_cast<G4PhysicsTable*>(dedxTable) << G4endl;
189 const std::size_t nCouples = dedxTable->size();
190 if(0 >= nCouples) { return; }
191
192 const std::size_t n = 100;
193 const G4double del = 1.0/(G4double)n;
194
195 for (std::size_t i=0; i<nCouples; ++i) {
196 auto pv = static_cast<G4PhysicsLogVector*>((*dedxTable)[i]);
197 if((pv == nullptr) || (isBaseMatActive && !(*theFlag)[i])) { continue; }
198 std::size_t npoints = pv->GetVectorLength();
199 std::size_t bin0 = 0;
200 G4double elow = pv->Energy(0);
201 G4double ehigh = pv->Energy(npoints-1);
202 G4double dedx1 = (*pv)[0];
203
204 //G4cout << "i= " << i << "npoints= " << npoints << " dedx1= "
205 //<< dedx1 << G4endl;
206
207 // protection for specific cases dedx=0
208 if(dedx1 == 0.0) {
209 for (std::size_t k=1; k<npoints; ++k) {
210 ++bin0;
211 elow = pv->Energy(k);
212 dedx1 = (*pv)[k];
213 if(dedx1 > 0.0) { break; }
214 }
215 npoints -= bin0;
216 }
217 //G4cout<<"New Range vector" << G4endl;
218 //G4cout<<"nbins= "<<npoints-1<<" elow= "<<elow<<" ehigh= "<<ehigh
219 // <<" bin0= " << bin0 <<G4endl;
220
221 // initialisation of a new vector
222 if(npoints < 3) { npoints = 3; }
223
224 delete (*rangeTable)[i];
226 if(0 == bin0) { v = new G4PhysicsLogVector(*pv); }
227 else { v = new G4PhysicsLogVector(elow, ehigh, npoints-1, splineFlag); }
228
229 // assumed dedx proportional to beta
230 G4double energy1 = v->Energy(0);
231 G4double range = 2.*energy1/dedx1;
232 //G4cout << "range0= " << range << G4endl;
233 v->PutValue(0,range);
234
235 for (std::size_t j=1; j<npoints; ++j) {
236
237 G4double energy2 = v->Energy(j);
238 G4double de = (energy2 - energy1) * del;
239 G4double energy = energy2 + de*0.5;
240 G4double sum = 0.0;
241 std::size_t idx = j - 1;
242 //G4cout << "j= " << j << " e1= " << energy1 << " e2= " << energy2
243 // << " n= " << n << G4endl;
244 for (std::size_t k=0; k<n; ++k) {
245 energy -= de;
246 dedx1 = pv->Value(energy, idx);
247 if(dedx1 > 0.0) { sum += de/dedx1; }
248 }
249 range += sum;
250 v->PutValue(j,range);
251 energy1 = energy2;
252 }
253 if(splineFlag) { v->FillSecondDerivatives(); }
255 }
256 //G4cout << "### Range table" << G4endl;
257 //G4cout << *rangeTable << G4endl;
258}
void FillSecondDerivatives(const G4SplineType=G4SplineType::Base, const G4double dir1=0.0, const G4double dir2=0.0)
G4double energy(const ThreeVector &p, const G4double m)

Referenced by G4TablesForExtrapolator::Initialisation().

◆ BuildTableForModel()

G4PhysicsTable * G4LossTableBuilder::BuildTableForModel ( G4PhysicsTable table,
G4VEmModel model,
const G4ParticleDefinition part,
G4double  emin,
G4double  emax,
G4bool  spline 
)

Definition at line 384 of file G4LossTableBuilder.cc.

389{
390 // check input
392 if(nullptr == table) { return table; }
393 if(emin >= emax) {
394 table->clearAndDestroy();
395 delete table;
396 table = nullptr;
397 return table;
398 }
400 G4int nbins = theParameters->NumberOfBinsPerDecade();
401
402 // Access to materials
403 const G4ProductionCutsTable* theCoupleTable=
405 std::size_t numOfCouples = theCoupleTable->GetTableSize();
406
407 G4PhysicsLogVector* aVector = nullptr;
408
409 for(G4int i=0; i<(G4int)numOfCouples; ++i) {
410 if ((*theFlag)[i]) {
411
412 // create physics vector and fill it
413 auto couple = theCoupleTable->GetMaterialCutsCouple(i);
414 delete (*table)[i];
415
416 // if start from zero then change the scale
417
418 const G4Material* mat = couple->GetMaterial();
419
420 G4double tmin = std::max(emin,model->MinPrimaryEnergy(mat,part));
421 if(0.0 >= tmin) { tmin = CLHEP::eV; }
422 G4int n = nbins;
423
424 if(tmin >= emax) {
425 aVector = nullptr;
426 } else {
427 n *= G4lrint(std::log10(emax/tmin));
428 n = std::max(n, 3);
429 aVector = new G4PhysicsLogVector(tmin, emax, n, spline);
430 }
431
432 if(nullptr != aVector) {
433 //G4cout << part->GetParticleName() << " in " << mat->GetName()
434 // << " tmin= " << tmin << G4endl;
435 for(G4int j=0; j<=n; ++j) {
436 aVector->PutValue(j, model->Value(couple, part,
437 aVector->Energy(j)));
438 }
439 if(spline) { aVector->FillSecondDerivatives(); }
440 }
442 }
443 }
444 /*
445 G4cout << "G4LossTableBuilder::BuildTableForModel done for "
446 << part->GetParticleName() << " and "<< model->GetName()
447 << " " << table << G4endl;
448 */
449 //G4cout << *table << G4endl;
450 return table;
451}
int G4int
Definition: G4Types.hh:85
G4int NumberOfBinsPerDecade() const
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:691
static G4PhysicsTable * PreparePhysicsTable(G4PhysicsTable *physTable)
void clearAndDestroy()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
std::size_t GetTableSize() const
static G4ProductionCutsTable * GetProductionCutsTable()
virtual G4double MinPrimaryEnergy(const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
Definition: G4VEmModel.cc:358
virtual G4double Value(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:349
int G4lrint(double ad)
Definition: templates.hh:134

Referenced by G4VMscModel::GetParticleChangeForMSC().

◆ GetBaseMaterialFlag()

◆ GetCoupleIndexes()

const std::vector< G4int > * G4LossTableBuilder::GetCoupleIndexes ( ) const

◆ GetDensityFactors()

const std::vector< G4double > * G4LossTableBuilder::GetDensityFactors ( ) const

Definition at line 125 of file G4LossTableBuilder.cc.

126{
127 return theDensityFactor;
128}

Referenced by G4VEmModel::G4VEmModel(), G4VEmProcess::G4VEmProcess(), and G4VEnergyLossProcess::G4VEnergyLossProcess().

◆ GetFlag()

G4bool G4LossTableBuilder::GetFlag ( size_t  idx)

Definition at line 132 of file G4LossTableBuilder.cc.

133{
134 if(theFlag->empty()) { InitialiseBaseMaterials(); }
135 return (idx < theFlag->size()) ? (*theFlag)[idx] : false;
136}

Referenced by G4EmTableUtil::BuildDEDXTable(), G4EmTableUtil::BuildLambdaTable(), and G4GammaGeneralProcess::BuildPhysicsTable().

◆ InitialiseBaseMaterials()

void G4LossTableBuilder::InitialiseBaseMaterials ( const G4PhysicsTable table = nullptr)

Definition at line 293 of file G4LossTableBuilder.cc.

294{
295 if(!isMaster) { return; }
296 const G4ProductionCutsTable* theCoupleTable=
298 std::size_t nCouples = theCoupleTable->GetTableSize();
299 std::size_t nFlags = theFlag->size();
300 /*
301 G4cout << "### InitialiseBaseMaterials: nCouples=" << nCouples
302 << " nFlags=" << nFlags << " isInit:" << isInitialized
303 << " baseMat:" << baseMatFlag << G4endl;
304 */
305 // define base material flag
306 if(isBaseMatActive && !baseMatFlag) {
307 for(G4int i=0; i<(G4int)nCouples; ++i) {
308 if(nullptr != theCoupleTable->GetMaterialCutsCouple(i)->GetMaterial()->GetBaseMaterial()) {
309 baseMatFlag = true;
310 isInitialized = false;
311 break;
312 }
313 }
314 }
315
316 if(nFlags != nCouples) { isInitialized = false; }
317 if(isInitialized) { return; }
318
319 // reserve memory
320 theFlag->resize(nCouples, true);
321 if(nullptr == table) { return; }
322
323 if(baseMatFlag) {
324 theDensityFactor->resize(nCouples,1.0);
325 theDensityIdx->resize(nCouples);
326 }
327
328 // define default flag and index of used material cut couple
329 for(G4int i=0; i<(G4int)nCouples; ++i) {
330 (*theFlag)[i] = table->GetFlag(i);
331 if(baseMatFlag) { (*theDensityIdx)[i] = i; }
332 }
333 isInitialized = true;
334 if(baseMatFlag) {
335 // use base materials
336 for(G4int i=0; i<(G4int)nCouples; ++i) {
337 // base material is needed only for a couple which is not
338 // initialised and for which tables will be computed
339 auto couple = theCoupleTable->GetMaterialCutsCouple(i);
340 auto pcuts = couple->GetProductionCuts();
341 auto mat = couple->GetMaterial();
342 auto bmat = mat->GetBaseMaterial();
343
344 // base material exists - find it and check if it can be reused
345 if(nullptr != bmat) {
346 for(G4int j=0; j<(G4int)nCouples; ++j) {
347 if(j == i) { continue; }
348 auto bcouple = theCoupleTable->GetMaterialCutsCouple(j);
349
350 if(bcouple->GetMaterial() == bmat &&
351 bcouple->GetProductionCuts() == pcuts) {
352
353 // based couple exist in the same region
354 (*theDensityFactor)[i] = mat->GetDensity()/bmat->GetDensity();
355 (*theDensityIdx)[i] = j;
356 (*theFlag)[i] = false;
357
358 // ensure that there will no double initialisation
359 (*theDensityFactor)[j] = 1.0;
360 (*theDensityIdx)[j] = j;
361 (*theFlag)[j] = true;
362 break;
363 }
364 }
365 }
366 }
367 }
368 /*
369 G4cout << "### G4LossTableBuilder::InitialiseBaseMaterials: flag="
370 << baseMatFlag << G4endl;
371 for(std::size_t i=0; i<nCouples; ++i) {
372 G4cout << "CoupleIdx=" << i << " Flag= " << (*theFlag)[i] << " "
373 << theCoupleTable->GetMaterialCutsCouple(i)->GetMaterial()->GetName()
374 << " TableFlag= " << table->GetFlag(i)
375 << " " << (*table)[i]
376 << G4endl;
377 }
378 */
379}
const G4Material * GetMaterial() const
G4ProductionCuts * GetProductionCuts() const
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:228
G4bool GetFlag(std::size_t i) const

Referenced by BuildDEDXTable(), BuildTableForModel(), GetBaseMaterialFlag(), GetFlag(), G4VEmProcess::PreparePhysicsTable(), and G4VEnergyLossProcess::PreparePhysicsTable().

◆ operator=()

G4LossTableBuilder & G4LossTableBuilder::operator= ( const G4LossTableBuilder right)
delete

◆ SetBaseMaterialActive()

void G4LossTableBuilder::SetBaseMaterialActive ( G4bool  flag)
inline

Definition at line 137 of file G4LossTableBuilder.hh.

138{
139 isBaseMatActive = flag;
140 if(!flag) {
141 baseMatFlag = false;
142 isInitialized = false;
143 }
144}

Referenced by G4TablesForExtrapolator::Initialisation().

◆ SetInitialisationFlag()

void G4LossTableBuilder::SetInitialisationFlag ( G4bool  flag)
inline

Definition at line 132 of file G4LossTableBuilder.hh.

133{
134 isInitialized = flag;
135}

Referenced by G4LossTableManager::ResetParameters().

◆ SetSplineFlag()

void G4LossTableBuilder::SetSplineFlag ( G4bool  flag)
inline

Definition at line 127 of file G4LossTableBuilder.hh.

128{
129 splineFlag = flag;
130}

The documentation for this class was generated from the following files: