G4LossTableBuilder Class Reference

#include <G4LossTableBuilder.hh>

Public Member Functions
	G4LossTableBuilder ()
virtual	~G4LossTableBuilder ()
void	BuildDEDXTable (G4PhysicsTable *dedxTable, const std::vector< G4PhysicsTable * > &)
void	BuildRangeTable (const G4PhysicsTable *dedxTable, G4PhysicsTable *rangeTable, G4bool isIonisation=false)
void	BuildInverseRangeTable (const G4PhysicsTable *rangeTable, G4PhysicsTable *invRangeTable, G4bool isIonisation=false)
G4PhysicsTable *	BuildTableForModel (G4PhysicsTable *table, G4VEmModel *model, const G4ParticleDefinition *, G4double emin, G4double emax, G4bool spline)
void	InitialiseBaseMaterials (G4PhysicsTable *table)
const std::vector< G4int > *	GetCoupleIndexes ()
const std::vector< G4double > *	GetDensityFactors ()
G4bool	GetFlag (size_t idx) const
void	SetSplineFlag (G4bool flag)
void	SetInitialisationFlag (G4bool flag)

---

Detailed Description

Definition at line 61 of file G4LossTableBuilder.hh.

---

Constructor & Destructor Documentation

G4LossTableBuilder::G4LossTableBuilder

(

)

Definition at line 72 of file G4LossTableBuilder.cc.

{
  splineFlag = true;
  isInitialized = false;

  theDensityFactor = new std::vector<G4double>;
  theDensityIdx = new std::vector<G4int>;
  theFlag = new std::vector<G4bool>;
}

G4LossTableBuilder::~G4LossTableBuilder

(

)

[virtual]

Definition at line 84 of file G4LossTableBuilder.cc.

{
  delete theDensityFactor;
  delete theDensityIdx;
  delete theFlag;
}

---

Member Function Documentation

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

Definition at line 94 of file G4LossTableBuilder.cc.

References G4PhysicsVector::FillSecondDerivatives(), G4PhysicsVector::PutValue(), G4PhysicsTableHelper::SetPhysicsVector(), and G4PhysicsVector::SetSpline().

{
  size_t n_processes = list.size();
  //G4cout << "Nproc= " << n_processes << " Ncoup= " << dedxTable->size() << G4endl;
  if(1 >= n_processes) { return; }

  size_t nCouples = dedxTable->size();
  if(0 >= nCouples) { return; }

  for (size_t i=0; i<nCouples; ++i) {
    //    if ((*theFlag)[i]) {
    G4PhysicsLogVector* pv0 = static_cast<G4PhysicsLogVector*>((*(list[0]))[i]);
    if(pv0) {
      size_t npoints = pv0->GetVectorLength();
      G4PhysicsLogVector* pv = new G4PhysicsLogVector(*pv0);
      //    pv = new G4PhysicsLogVector(elow, ehigh, npoints-1);
      pv->SetSpline(splineFlag);
      for (size_t j=0; j<npoints; ++j) {
        G4double dedx = 0.0;
        for (size_t k=0; k<n_processes; ++k) {
          G4PhysicsVector* pv1   = (*(list[k]))[i];
          dedx += (*pv1)[j];
        }
        pv->PutValue(j, dedx);
      }
      if(splineFlag) { pv->FillSecondDerivatives(); }
      G4PhysicsTableHelper::SetPhysicsVector(dedxTable, i, pv);
    }
  }
}

void G4LossTableBuilder::BuildInverseRangeTable	(	const G4PhysicsTable *	rangeTable,
		G4PhysicsTable *	invRangeTable,
		G4bool	isIonisation = false
	)

Definition at line 213 of file G4LossTableBuilder.cc.

References G4PhysicsVector::Energy(), G4PhysicsVector::FillSecondDerivatives(), G4PhysicsVector::GetVectorLength(), G4LPhysicsFreeVector::PutValues(), G4PhysicsTableHelper::SetPhysicsVector(), and G4PhysicsVector::SetSpline().

{
  size_t nCouples = rangeTable->size();
  if(0 >= nCouples) { return; }

  for (size_t i=0; i<nCouples; ++i) {

    if(isIonisation) {
      if( !(*theFlag)[i] ) { continue; }
    }
    G4PhysicsVector* pv = (*rangeTable)[i];
    size_t npoints = pv->GetVectorLength();
    G4double rlow  = (*pv)[0];
    G4double rhigh = (*pv)[npoints-1];
      
    delete (*invRangeTable)[i];
    G4LPhysicsFreeVector* v = new G4LPhysicsFreeVector(npoints,rlow,rhigh);
    v->SetSpline(splineFlag);

    for (size_t j=0; j<npoints; ++j) {
      G4double e  = pv->Energy(j);
      G4double r  = (*pv)[j];
      v->PutValues(j,r,e);
    }
    if(splineFlag) { v->FillSecondDerivatives(); }

    G4PhysicsTableHelper::SetPhysicsVector(invRangeTable, i, v);
  }
}

void G4LossTableBuilder::BuildRangeTable	(	const G4PhysicsTable *	dedxTable,
		G4PhysicsTable *	rangeTable,
		G4bool	isIonisation = false
	)

Definition at line 128 of file G4LossTableBuilder.cc.

References G4PhysicsVector::Energy(), G4PhysicsVector::FillSecondDerivatives(), G4PhysicsVector::GetVectorLength(), CLHEP::detail::n, G4PhysicsVector::PutValue(), G4PhysicsTableHelper::SetPhysicsVector(), G4PhysicsVector::SetSpline(), and G4PhysicsVector::Value().

{
  size_t nCouples = dedxTable->size();
  if(0 >= nCouples) { return; }

  size_t n = 100;
  G4double del = 1.0/(G4double)n;

  for (size_t i=0; i<nCouples; ++i) {
    if(isIonisation) {
      if( !(*theFlag)[i] ) { continue; }
    }
    G4PhysicsLogVector* pv = static_cast<G4PhysicsLogVector*>((*dedxTable)[i]);
    size_t npoints = pv->GetVectorLength();
    size_t bin0    = 0;
    G4double elow  = pv->Energy(0);
    G4double ehigh = pv->Energy(npoints-1);
    G4double dedx1 = (*pv)[0];

    //G4cout << "i= " << i << "npoints= " << npoints << " dedx1= " << dedx1 << G4endl;

    // protection for specific cases dedx=0
    if(dedx1 == 0.0) {
      for (size_t k=1; k<npoints; ++k) {
        bin0++;
        elow  = pv->Energy(k);
        dedx1 = (*pv)[k];
        if(dedx1 > 0.0) { break; }
      }
      npoints -= bin0;
    }
    //G4cout<<"New Range vector" << G4endl;
    //G4cout<<"nbins= "<<npoints-1<<" elow= "<<elow<<" ehigh= "<<ehigh
    //      <<" bin0= " << bin0 <<G4endl;

    // initialisation of a new vector
    if(npoints < 2) { npoints = 2; }

    delete (*rangeTable)[i];
    G4PhysicsLogVector* v;
    if(0 == bin0) { v = new G4PhysicsLogVector(*pv); }
    else { v = new G4PhysicsLogVector(elow, ehigh, npoints-1); }

    // dedx is exact zero cannot build range table
    if(2 == npoints) {
      v->PutValue(0,1000.);
      v->PutValue(1,2000.);
      G4PhysicsTableHelper::SetPhysicsVector(rangeTable, i, v);
      return;
    }
    v->SetSpline(splineFlag);

    // assumed dedx proportional to beta
    G4double energy1 = v->Energy(0);
    G4double range   = 2.*energy1/dedx1;
    //G4cout << "range0= " << range << G4endl;
    v->PutValue(0,range);

    for (size_t j=1; j<npoints; ++j) {

      G4double energy2 = v->Energy(j);
      G4double de      = (energy2 - energy1) * del;
      G4double energy  = energy2 + de*0.5;
      G4double sum = 0.0;
      //G4cout << "j= " << j << " e1= " << energy1 << " e2= " << energy2 
      //       << " n= " << n << G4endl;
      for (size_t k=0; k<n; ++k) {
        energy -= de;
        dedx1 = pv->Value(energy);
        if(dedx1 > 0.0) { sum += de/dedx1; }
      }
      range += sum;
      v->PutValue(j,range);
      energy1 = energy2;
    }
    if(splineFlag) { v->FillSecondDerivatives(); }
    G4PhysicsTableHelper::SetPhysicsVector(rangeTable, i, v);
  }
}

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

Definition at line 400 of file G4LossTableBuilder.cc.

References G4PhysicsTable::clearAndDestroy(), G4PhysicsVector::Energy(), G4PhysicsVector::FillSecondDerivatives(), GetFlag(), G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4LossTableManager::GetNumberOfBinsPerDecade(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), InitialiseBaseMaterials(), G4LossTableManager::Instance(), G4VEmModel::MinPrimaryEnergy(), CLHEP::detail::n, G4PhysicsTableHelper::PreparePhysicsTable(), G4PhysicsVector::PutValue(), G4PhysicsTableHelper::SetPhysicsVector(), G4PhysicsVector::SetSpline(), and G4VEmModel::Value().

{
  // check input
  G4PhysicsTable* table = G4PhysicsTableHelper::PreparePhysicsTable(aTable);
  if(!table) { return table; }
  if(emin >= emax) { 
    table->clearAndDestroy();
    delete table;
    table = 0;
    return table; 
  }
  InitialiseBaseMaterials(table);

  G4int nbins = G4int(std::log10(emax/emin) + 0.5)
    *G4LossTableManager::Instance()->GetNumberOfBinsPerDecade();
  if(nbins < 3) { nbins = 3; }

  // Access to materials
  const G4ProductionCutsTable* theCoupleTable=
        G4ProductionCutsTable::GetProductionCutsTable();
  size_t numOfCouples = theCoupleTable->GetTableSize();

  G4PhysicsLogVector* aVector = 0;
  G4PhysicsLogVector* bVector = 0;

  for(size_t i=0; i<numOfCouples; ++i) {

    //G4cout<< "i= " << i << " Flag=  " << GetFlag(i) << G4endl;

    if (GetFlag(i)) {

      // create physics vector and fill it
      const G4MaterialCutsCouple* couple = 
        theCoupleTable->GetMaterialCutsCouple(i);
      delete (*table)[i];

      // if start from zero then change the scale

      const G4Material* mat = couple->GetMaterial();

      G4double tmin = std::max(emin,model->MinPrimaryEnergy(mat,part));
      if(0.0 >= tmin) { tmin = eV; }
      G4int n = nbins + 1;

      if(tmin >= emax) {
        aVector = 0;
      } else if(tmin > emin) {
        G4int bin = nbins*G4int(std::log10(emax/tmin) + 0.5);
        if(bin < 3) { bin = 3; }
        n = bin + 1;
        aVector = new G4PhysicsLogVector(tmin, emax, bin);

      } else if(!bVector) {
        aVector = new G4PhysicsLogVector(emin, emax, nbins);
        bVector = aVector;

      } else {
        aVector = new G4PhysicsLogVector(*bVector);
      }

      if(aVector) {
        aVector->SetSpline(spline);
        for(G4int j=0; j<n; ++j) {
          aVector->PutValue(j, model->Value(couple, part, aVector->Energy(j)));
        }
        if(spline) { aVector->FillSecondDerivatives(); }
      }
      G4PhysicsTableHelper::SetPhysicsVector(table, i, aVector);
    }
  }
  /*
  G4cout << "G4LossTableBuilder::BuildTableForModel done for "
         << part->GetParticleName() << " and "<< model->GetName()
         << "  " << table << G4endl;
  */
  return table; 
}

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

(

)

[inline]

Definition at line 123 of file G4LossTableBuilder.hh.

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

{
  if(theDensityIdx->size() == 0) { InitialiseCouples(); }
  return theDensityIdx;
}

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

(

)

[inline]

Definition at line 130 of file G4LossTableBuilder.hh.

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

{
  if(theDensityIdx->size() == 0) { InitialiseCouples(); }
  return theDensityFactor;
}

G4bool G4LossTableBuilder::GetFlag

(

size_t

idx

)

const [inline]

Definition at line 136 of file G4LossTableBuilder.hh.

Referenced by G4VEnergyLossProcess::BuildDEDXTable(), G4VEnergyLossProcess::BuildLambdaTable(), and BuildTableForModel().

{
  return (*theFlag)[idx];
}

void G4LossTableBuilder::InitialiseBaseMaterials

(

G4PhysicsTable *

table

)

Definition at line 249 of file G4LossTableBuilder.cc.

References G4Material::GetBaseMaterial(), G4Material::GetDensity(), G4PhysicsTable::GetFlag(), G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4MaterialCutsCouple::GetProductionCuts(), and G4ProductionCutsTable::GetProductionCutsTable().

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

{
  size_t nCouples = table->size();
  size_t nFlags = theFlag->size();

  if(nCouples == nFlags && isInitialized) { return; }

  isInitialized = true;

  //G4cout << "%%%%%% G4LossTableBuilder::InitialiseBaseMaterials Ncouples= " 
  //     << nCouples << "  FlagSize= " << nFlags << G4endl; 

  // variable density check
  const G4ProductionCutsTable* theCoupleTable=
    G4ProductionCutsTable::GetProductionCutsTable();

  /*
  for(size_t i=0; i<nFlags; ++i) {
    G4cout << "CoupleIdx= " << i << "  Flag= " <<  (*theFlag)[i] 
           << " tableFlag= " << table->GetFlag(i) << "  "
           << theCoupleTable->GetMaterialCutsCouple(i)->GetMaterial()->GetName()
           << G4endl;
  }
  */

  // expand vectors
  if(nFlags < nCouples) {
    for(size_t i=nFlags; i<nCouples; ++i) { theDensityFactor->push_back(1.0); }
    for(size_t i=nFlags; i<nCouples; ++i) { theDensityIdx->push_back(-1); }
    for(size_t i=nFlags; i<nCouples; ++i) { theFlag->push_back(true); }
  }
  for(size_t i=0; i<nCouples; ++i) {

    // base material is needed only for a couple which is not
    // initialised and for which tables will be computed
    (*theFlag)[i] = table->GetFlag(i);
    if ((*theDensityIdx)[i] < 0) {
      (*theDensityIdx)[i] = i;
      const G4MaterialCutsCouple* couple =
        theCoupleTable->GetMaterialCutsCouple(i);
      const G4ProductionCuts* pcuts = couple->GetProductionCuts();
      const G4Material* mat  = couple->GetMaterial();
      const G4Material* bmat = mat->GetBaseMaterial();

      // base material exists - find it and check if it can be reused
      if(bmat) {
        for(size_t j=0; j<nCouples; ++j) {

          if(j == i) { continue; }
          const G4MaterialCutsCouple* bcouple =
              theCoupleTable->GetMaterialCutsCouple(j);

          if(bcouple->GetMaterial() == bmat && 
             bcouple->GetProductionCuts() == pcuts) {

            // based couple exist in the same region
            (*theDensityIdx)[i] = j;
            (*theDensityFactor)[i] = mat->GetDensity()/bmat->GetDensity();
            (*theFlag)[i] = false;

            // ensure that there will no double initialisation
            (*theDensityIdx)[j] = j;
            (*theDensityFactor)[j] = 1.0;
            (*theFlag)[j] = true;
            break;
          }
        }
      }
    }
  }
  /*
  for(size_t i=0; i<nCouples; ++i) {
    G4cout << "CoupleIdx= " << i << "  Flag= " <<  (*theFlag)[i] 
           << "  TableFlag= " << table->GetFlag(i) << "  "
           << theCoupleTable->GetMaterialCutsCouple(i)->GetMaterial()->GetName()
           << G4endl;
  }
  G4cout << "%%%%%% G4LossTableBuilder::InitialiseBaseMaterials end" 
         << G4endl; 
  */
}

void G4LossTableBuilder::SetInitialisationFlag

(

G4bool

flag

)

[inline]

Definition at line 146 of file G4LossTableBuilder.hh.

Referenced by G4LossTableManager::PreparePhysicsTable().

{
  isInitialized = flag;
}

void G4LossTableBuilder::SetSplineFlag

(

G4bool

flag

)

[inline]

Definition at line 141 of file G4LossTableBuilder.hh.

Referenced by G4LossTableManager::SetSplineFlag().

{
  splineFlag = flag;
}

---

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

• G4LossTableBuilder.hh
• G4LossTableBuilder.cc

---