G4ChannelingFastSimModel Class Reference

#include <G4ChannelingFastSimModel.hh>

Inheritance diagram for G4ChannelingFastSimModel:

Public Member Functions
	G4ChannelingFastSimModel (const G4String &, G4Region *)
	G4ChannelingFastSimModel (const G4String &)
	~G4ChannelingFastSimModel ()
G4bool	IsApplicable (const G4ParticleDefinition &) override
	– IsApplicable
G4bool	ModelTrigger (const G4FastTrack &) override
	– ModelTrigger
void	DoIt (const G4FastTrack &, G4FastStep &) override
	– User method DoIt
void	Input (const G4Material *crystal, const G4String &lattice)
	special functions
void	Input (const G4Material *crystal, const G4String &lattice, const G4String &filePath)
void	RadiationModelActivate ()
G4ChannelingFastSimCrystalData *	GetCrystalData ()
G4BaierKatkov *	GetRadiationModel ()
G4bool	GetIfRadiationModelActive ()
void	SetLowKineticEnergyLimit (G4double ekinetic, const G4String &particleName)
	set cuts
void	SetLindhardAngleNumberHighLimit (G4double angleNumber, const G4String &particleName)
void	SetHighAngleLimit (G4double anglemax, const G4String &particleName)
void	SetDefaultLowKineticEnergyLimit (G4double ekinetic)
void	SetDefaultLindhardAngleNumberHighLimit (G4double angleNumber)
void	SetDefaultHighAngleLimit (G4double anglemax)
void	SetMaxPhotonsProducedPerStep (G4double nPhotons)
G4double	GetLowKineticEnergyLimit (G4int particleDefinitionID)
	get cuts
G4double	GetLindhardAngleNumberHighLimit (G4int particleDefinitionID)
G4double	GetHighAngleLimit (G4int particleDefinitionID)
G4int	GetMaxPhotonsProducedPerStep ()
	get the maximal number of photons that can be produced per fastStep
Public Member Functions inherited from G4VFastSimulationModel
	G4VFastSimulationModel (const G4String &aName)
	G4VFastSimulationModel (const G4String &aName, G4Envelope *, G4bool IsUnique=FALSE)
virtual	~G4VFastSimulationModel ()=default
virtual G4bool	AtRestModelTrigger (const G4FastTrack &)
virtual void	AtRestDoIt (const G4FastTrack &, G4FastStep &)
virtual void	Flush ()
const G4String	GetName () const
G4bool	operator== (const G4VFastSimulationModel &) const

Detailed Description

Definition at line 54 of file G4ChannelingFastSimModel.hh.

Constructor & Destructor Documentation

G4ChannelingFastSimModel() [1/2]

G4ChannelingFastSimModel::G4ChannelingFastSimModel	(	const G4String &	modelName,
		G4Region *	envelope )

Definition at line 45 of file G4ChannelingFastSimModel.cc.

: G4VFastSimulationModel(modelName, envelope)
{
}

G4ChannelingFastSimModel() [2/2]

G4ChannelingFastSimModel::G4ChannelingFastSimModel

(

const G4String &

modelName

)

Definition at line 53 of file G4ChannelingFastSimModel.cc.

: G4VFastSimulationModel(modelName)
{
 
}

~G4ChannelingFastSimModel()

G4ChannelingFastSimModel::~G4ChannelingFastSimModel

(

)

Definition at line 61 of file G4ChannelingFastSimModel.cc.

{
}

Member Function Documentation

DoIt()

void G4ChannelingFastSimModel::DoIt ( const G4FastTrack & fastTrack, G4FastStep & fastStep )

overridevirtual

– User method DoIt

Implements G4VFastSimulationModel.

Definition at line 141 of file G4ChannelingFastSimModel.cc.

{  
  G4double etotal;//particle total energy
  G4double etotalPreStep;//etotal at the previous step
  G4double etotalToSetParticleProperties;//etotal value at which
                                         //SetParticleProperties is calculated
  G4double ekinetic = 0;//kinetic energy
  G4double eDeposited = 0.;//deposited energy along the trajectory
  G4double elossAccum = 0;// accumulate local energy loss (not radiation)
  G4double mass;  //particle mass
  G4double charge;//particle charge
  G4double tGlobal; //global time
  G4double tGlobalPreStep; //global time at the previous step
  G4ThreeVector xyz0;// the coordinates in the local reference system of the volume
  G4ThreeVector xyz0PreStep;// xyz at the previous step
  G4ThreeVector xyz;// the coordinates in the co-rotating reference system within
                    //a channel (elementary periodic cell)
  G4double x,y,z;   // the coordinates in the co-rotating reference system within
                    //a channel (elementary periodic cell)
  G4double tx0,ty0; // the angles in the local reference system of the volume
  G4double tx,ty;   // the angles in the co-rotating reference system within
                    //a channel (elementary periodic cell)
  G4double txPreStep,tyPreStep;// tx,ty at the previous step
  G4ThreeVector momentumDirection;
  G4ThreeVector scatteringAnglesAndEnergyLoss;//output of scattering functions
  G4double lindhardAngleNumberHighLimit0; //current high limit of the angle expressed in
                                          //[Lindhard angle] units
  G4double highAngleLimit0; //current absolute high limit of the angle expressed
 
  //coordinates in Runge-Kutta calculations
  G4double x1=0.,x2=0.,x3=0.,x4=0.,y1=0.,y2=0.,y3=0.,y4=0.;
  //angles in Runge-Kutta calculations
  G4double tx1=0.,tx2=0.,tx3=0.,tx4=0.,ty1=0.,ty2=0.,ty3=0.,ty4=0.;
  //variables in Runge-Kutta calculations
  G4double kvx1=0.,kvx2=0.,kvx3=0.,kvx4=0.,kvy1=0.,kvy2=0.,kvy3=0.,kvy4=0.;
  //simulation step along z (internal step of the model) and its parts
  G4double dz,dzd3,dzd8;//dzd3 = dz/3; dzd8 = dz/8;
  //simulation step along the momentum direction
  G4double momentumDirectionStep;
  //effective simulation step (taking into account nuclear density along the trajectory)
  G4double effectiveStep;
 
  // flag, if Inside(xyz0) switches to kInside
  G4bool inside = false;
 
  G4LogicalVolume* crystallogic = fastTrack.GetEnvelopeLogicalVolume();
  fCrystalData->SetGeometryParameters(crystallogic);
 
  //set the max number of secondaries (photons) that can be added at this fastStep
  if (fRad)
  {
      fastStep.SetNumberOfSecondaryTracks(fMaxPhotonsProducedPerStep);
      //reseting the BaierKatkov integral to start it with the new trajectory
      fBaierKatkov->ResetRadIntegral();//to avoid any memory from the previous trajectory
  }
 
  mass = fastTrack.GetPrimaryTrack()->GetParticleDefinition()->GetPDGMass();
  etotal = mass + fastTrack.GetPrimaryTrack()->GetKineticEnergy();
  charge = fastTrack.GetPrimaryTrack()->GetParticleDefinition()->GetPDGCharge();
 
  G4String particleName =
      fastTrack.GetPrimaryTrack()->GetParticleDefinition()->GetParticleName();
 
  lindhardAngleNumberHighLimit0 =
      GetLindhardAngleNumberHighLimit(fastTrack.GetPrimaryTrack()->
                                      GetParticleDefinition()->GetParticleDefinitionID());
  highAngleLimit0 = GetHighAngleLimit(fastTrack.GetPrimaryTrack()->
                                      GetParticleDefinition()->GetParticleDefinitionID());
 
  //set fCrystalData parameters depending on the particle parameters
  fCrystalData->SetParticleProperties(etotal, mass, charge, particleName);
 
  //global time
  tGlobal = fastTrack.GetPrimaryTrack()->GetGlobalTime();
 
  //coordinates in the co-rotating reference system within a channel
  xyz0= fastTrack.GetPrimaryTrackLocalPosition();
  xyz = fCrystalData->CoordinatesFromBoxToLattice(xyz0);
  x=xyz.x();
  y=xyz.y();
  z=xyz.z();
 
  momentumDirection=fastTrack.GetPrimaryTrackLocalDirection();
  //angle in the co-rotating reference system within a channel
  //(!!! ONLY FORWARD DIRECTION, momentumDirection.getZ()>0,
  //valid for high energies defined by the standard energy cuts)
  tx0 = std::atan(momentumDirection.x()/momentumDirection.z());
  ty0 = std::atan(momentumDirection.y()/momentumDirection.z());
 
  //angles in the co-rotating reference system within a channel
  tx = fCrystalData->AngleXFromBoxToLattice(tx0,z);
  ty = ty0;
 
  etotalToSetParticleProperties = etotal*0.999;
  G4bool inCrystal=true;//flag necessary to escape the cycle (at inCrystal=0;)
  //do calculations until the particle is inside the volume
  do
  {
      //remember the global time before the next step dz
      tGlobalPreStep=tGlobal;
      //remember the coordinates before the next step dz
      xyz0PreStep = xyz0;
      //remember the angles and the total energy before the step dz
      txPreStep = tx;
      tyPreStep = ty;
      etotalPreStep = etotal;
 
      dz = fCrystalData->GetSimulationStep(tx,ty);
      dzd3=dz/3;
      dzd8=dz/8;
 
      //trajectory calculation:
      //Runge-Cutt "3/8"
      //fCrystalData->GetCurv(z)*fCrystalData->GetCorrectionZ() is due to dependence of
      //the radius on x; GetCurv gets 1/R for the central ("central plane/axis")
 
      //first step
      kvx1=fCrystalData->Ex(x,y);
      x1=x+tx*dzd3;
      tx1=tx+(kvx1-fCrystalData->GetCurv(z)*fCrystalData->GetCorrectionZ())*dzd3;
      if (fCrystalData->GetModel()==2)
      {
         kvy1=fCrystalData->Ey(x,y);
         y1=y+ty*dzd3;
         ty1=ty+kvy1*dzd3;
      }
 
      //second step
      kvx2=fCrystalData->Ex(x1,y1);
      x2=x-tx*dzd3+tx1*dz;
      tx2=tx-(kvx1-fCrystalData->GetCurv(z)*fCrystalData->GetCorrectionZ())*dzd3+
              (kvx2-fCrystalData->GetCurv(z)*fCrystalData->GetCorrectionZ())*dz;
      if (fCrystalData->GetModel()==2)
      {
         kvy2=fCrystalData->Ey(x1,y1);
         y2=y-ty*dzd3+ty1*dz;
         ty2=ty-kvy1*dzd3+kvy2*dz;
      }
 
      //third step
      kvx3=fCrystalData->Ex(x2,y2);
      x3=x+(tx-tx1+tx2)*dz;
      tx3=tx+(kvx1-kvx2+kvx3-fCrystalData->GetCurv(z)*fCrystalData->GetCorrectionZ())*dz;
      if (fCrystalData->GetModel()==2)
      {
         kvy3=fCrystalData->Ey(x2,y2);
         y3=y+(ty-ty1+ty2)*dz;
         ty3=ty+(kvy1-kvy2+kvy3)*dz;
      }
 
      //fourth step
      kvx4=fCrystalData->Ex(x3,y3);
      x4=x+(tx+3.*tx1+3.*tx2+tx3)*dzd8;
      tx4=tx+(kvx1+3.*kvx2+3.*kvx3+kvx4)*dzd8-
              fCrystalData->GetCurv(z)*fCrystalData->GetCorrectionZ()*dz;
      if (fCrystalData->GetModel()==2)
      {
          kvy4=fCrystalData->Ey(x3,y3);
          y4=y+(ty+3.*ty1+3.*ty2+ty3)*dzd8;
          ty4=ty+(kvy1+3.*kvy2+3.*kvy3+kvy4)*dzd8;
      }
      else
      {
          y4 =y+ty*dz;
          ty4=ty;
      }
 
      x=x4;
      tx=tx4;
      y=y4;
      ty=ty4;
 
      z+=dz*fCrystalData->GetCorrectionZ();//motion along the z coordinate
                                          //("central plane/axis", no current plane/axis)
 
      xyz = fCrystalData->ChannelChange(x,y,z);
      x=xyz.x();
      y=xyz.y();
      z=xyz.z();
 
      //the coordinates in the local reference system of the volume
      //this vector will be used in the cycle escape condition and
      //in the radiation model (if activated)
      xyz0=fCrystalData->CoordinatesFromLatticeToBox(xyz);
 
      momentumDirectionStep=
              dz*std::sqrt(1+std::pow(std::tan(tx),2)+std::pow(std::tan(ty),2));
      tGlobal+=momentumDirectionStep/(fCrystalData->GetBeta())/CLHEP::c_light;
 
      //default scattering and energy loss 0
      scatteringAnglesAndEnergyLoss = G4ThreeVector(0.,0.,0.);
 
      //calculate separately for each element of the crystal
      for (G4int i = 0; i < fCrystalData->GetNelements(); i++)
      {
          //effective step taking into account nuclear density along the trajectory
          effectiveStep = momentumDirectionStep*fCrystalData->NuclearDensity(x,y,i);
          //Coulomb scattering on screened atomic potential (both multiple and single)
          scatteringAnglesAndEnergyLoss += fCrystalData->
                         CoulombAtomicScattering(effectiveStep,momentumDirectionStep,i);
 
          //Amorphous part of ionization energy losses
          elossAccum += fCrystalData->IonizationLosses(momentumDirectionStep, i);
      }
      //electron scattering and coherent part of ionization energy losses
      scatteringAnglesAndEnergyLoss += fCrystalData->CoulombElectronScattering(
                                                   fCrystalData->MinIonizationEnergy(x,y),
                                                   fCrystalData->ElectronDensity(x,y),
                                                   momentumDirectionStep);
      tx += scatteringAnglesAndEnergyLoss.x();
      ty += scatteringAnglesAndEnergyLoss.y();
      elossAccum += scatteringAnglesAndEnergyLoss.z();
 
      // recalculate the energy depended parameters
      //(only if the energy decreased enough, not at each step)
      if (etotalToSetParticleProperties>etotal)
      {
          fCrystalData->SetParticleProperties(etotal, mass, charge, particleName);
          etotalToSetParticleProperties = etotal*0.999;
      }
 
      //chain of conditions to escape the cycle
      // if Inside(xyz0)==kInside has been already true
      //(a particle has been inside the crystal)
      if (inside)
      {
          // if low energy
         if (etotal-mass<=GetLowKineticEnergyLimit(fastTrack.GetPrimaryTrack()->
                                                   GetParticleDefinition()->
                                                   GetParticleDefinitionID()))
             {inCrystal = false;}//escape the cycle
         //check if the angle w.r.t. the axes or planes is too high =>
         //return to standard Geant4:
         else if (fCrystalData->GetModel()==1) //1D model, field of planes
         {
             //if the angle w.r.t. the planes is too high
             if (std::abs(tx) >=
                 std::max(lindhardAngleNumberHighLimit0*
                              fCrystalData->GetLindhardAngle(),
                          highAngleLimit0))
                {inCrystal = false;}//escape the cycle
         }
         else if (fCrystalData->GetModel()==2) //2D model, field of axes
         {
             //if the angle w.r.t. the axes is too high
             if (std::sqrt(tx*tx+ty*ty) >=
                 std::max(lindhardAngleNumberHighLimit0*
                              fCrystalData->GetLindhardAngle(),
                          highAngleLimit0))
                {inCrystal = false;}//escape the cycle
         }
 
           //radiation production & radiation energy losses
           //works only if the radiation model is activated
           if (fRad)
           {
               //back to the local reference system of the volume
               tx0 = fCrystalData->AngleXFromLatticeToBox(tx,z);
               ty0 = ty;
               //xyz0 was calculated above
 
               //running the radiation model and checking if a photon has been emitted
               if(fBaierKatkov->DoRadiation(etotal,mass,
                                           tx0,ty0,
                                           scatteringAnglesAndEnergyLoss.x(),
                                           scatteringAnglesAndEnergyLoss.y(),
                                           momentumDirectionStep,tGlobal,xyz0,
                                           crystallogic->
                                           GetSolid()->
                                           Inside(xyz0)!=kInside&&inCrystal))
               // also it was checked if the particle is escaping the volume
               // calculate the radiation integral immidiately in this case
               {
                   //a photon has been emitted!
                   //shift the particle back into the radiation point
                   etotal = fBaierKatkov->GetParticleNewTotalEnergy();
                   tx0 = fBaierKatkov->GetParticleNewAngleX();
                   ty0 = fBaierKatkov->GetParticleNewAngleY();
                   tGlobal = fBaierKatkov->GetNewGlobalTime();
                   xyz0 = fBaierKatkov->GetParticleNewCoordinateXYZ();
 
                   //add secondary photon
                   fBaierKatkov->GeneratePhoton(fastStep);
 
                   //particle energy was changed
                   fCrystalData->SetParticleProperties(etotal, mass, charge, particleName);
 
                   //coordinates in the co-rotating reference system within a channel
                   xyz = fCrystalData->CoordinatesFromBoxToLattice(xyz0);
                   x=xyz.x();
                   y=xyz.y();
                   z=xyz.z();
 
                   //angles in the co-rotating reference system within a channel
                   tx = fCrystalData->AngleXFromBoxToLattice(tx0,z);
                   ty = ty0;
               }
           }
           else
           {
               //we calculate deposited energy and energy losses ONLY in absence
               //of radiation otherwise we do it only at the end of model
               etotal -= elossAccum;
               eDeposited += elossAccum;
               elossAccum=0;
               ekinetic = etotal-mass;
               if(ekinetic<1*keV)
               {
                   G4cout << "Warning in G4ChannelingFastSimModel: " <<
                   ekinetic << "<" << 1*keV << " !" << G4endl;
                   eDeposited-=(1*keV-ekinetic);
                   ekinetic = 1*keV;
                   G4cout << "Setting deposited energy=" <<
                   eDeposited << " & ekinetic=" << ekinetic << G4endl;
                   etotal = mass+ekinetic;
               }
           }
 
           //precise check if the particle is escaping the volume
           if (crystallogic->GetSolid()->
                   Inside(xyz0)!=kInside)
           {
               //one step back to remain inside the volume
               //after the escape of the volume
               tGlobal = tGlobalPreStep;
               xyz0 = xyz0PreStep;
               tx = txPreStep;
               ty = tyPreStep;
               etotal = etotalPreStep;
               z-=dz*fCrystalData->GetCorrectionZ();
               // change the flag => this particle will not enter
               // the model before escape this volume
 
               inCrystal = false; //escape the cycle
           }
      }
      else
      {
          // if Inside(xyz0)==kInside we can enable checking of particle escape
          if (crystallogic->GetSolid()->
                           Inside(xyz0)==kInside)
                 {inside = true;}
          // a very rare case, if a particle remains
          // on the boundary and escapes the crystal
          else if (crystallogic->GetSolid()->
                           Inside(xyz0)==kOutside)
                 {inCrystal = false;}//escape the cycle
      }
  }
  while (inCrystal);
 
  //the angles in the local reference system of the volume
  tx0 = fCrystalData->AngleXFromLatticeToBox(tx,z);
  ty0 = ty;
 
  //set global time
  fastStep.ProposePrimaryTrackFinalTime(tGlobal);
  //set final position
  fastStep.ProposePrimaryTrackFinalPosition(xyz0);
 
  //set deposited energy (due to ionization)
  etotal -= elossAccum;
  eDeposited += elossAccum;
  ekinetic = etotal-mass;
  if(ekinetic<1*keV)
  {
      G4cout << "Warning in G4ChannelingFastSimModel: " <<
      ekinetic << "<" << 1*keV << " !" << G4endl;
      eDeposited-=(1*keV-ekinetic);
      ekinetic = 1*keV;
      G4cout << "Setting deposited energy=" <<
      eDeposited << " & ekinetic=" << ekinetic << G4endl;
  }
  fastStep.ProposeTotalEnergyDeposited(eDeposited);
  //set final kinetic energy
  fastStep.ProposePrimaryTrackFinalKineticEnergy(ekinetic);
 
 
  //set final momentum direction
  G4double momentumDirectionZ =
          1./std::sqrt(1.+std::pow(std::tan(tx0),2)+std::pow(std::tan(ty0),2));
  fastStep.ProposePrimaryTrackFinalMomentumDirection(
              G4ThreeVector(momentumDirectionZ*std::tan(tx0),
                            momentumDirectionZ*std::tan(ty0),
                            momentumDirectionZ));
}

GetCrystalData()

G4ChannelingFastSimCrystalData * G4ChannelingFastSimModel::GetCrystalData ( )

inline

Definition at line 80 of file G4ChannelingFastSimModel.hh.

{return fCrystalData;}

Referenced by G4CoherentPairProductionPhysics::ConstructProcess().

GetHighAngleLimit()

G4double G4ChannelingFastSimModel::GetHighAngleLimit ( G4int particleDefinitionID )

inline

Definition at line 119 of file G4ChannelingFastSimModel.hh.

               {return (fHighAngleLimit.count(particleDefinitionID) == 1)
                              ? fHighAngleLimit[particleDefinitionID]
                              : fDefaultHighAngleLimit;}

Referenced by DoIt(), and ModelTrigger().

GetIfRadiationModelActive()

G4bool G4ChannelingFastSimModel::GetIfRadiationModelActive ( )

inline

Definition at line 84 of file G4ChannelingFastSimModel.hh.

{return fRad;}

GetLindhardAngleNumberHighLimit()

G4double G4ChannelingFastSimModel::GetLindhardAngleNumberHighLimit ( G4int particleDefinitionID )

inline

Definition at line 115 of file G4ChannelingFastSimModel.hh.

               {return (fLindhardAngleNumberHighLimit.count(particleDefinitionID) == 1)
                        ? fLindhardAngleNumberHighLimit[particleDefinitionID]
                        : fDefaultLindhardAngleNumberHighLimit;}

Referenced by DoIt(), and ModelTrigger().

GetLowKineticEnergyLimit()

G4double G4ChannelingFastSimModel::GetLowKineticEnergyLimit ( G4int particleDefinitionID )

inline

get cuts

Definition at line 111 of file G4ChannelingFastSimModel.hh.

               {return (fLowEnergyLimit.count(particleDefinitionID) == 1)
                        ? fLowEnergyLimit[particleDefinitionID]
                        : fDefaultLowEnergyLimit;}

Referenced by DoIt(), and ModelTrigger().

GetMaxPhotonsProducedPerStep()

G4int G4ChannelingFastSimModel::GetMaxPhotonsProducedPerStep ( )

inline

get the maximal number of photons that can be produced per fastStep

Definition at line 125 of file G4ChannelingFastSimModel.hh.

{return fMaxPhotonsProducedPerStep;}

GetRadiationModel()

G4BaierKatkov * G4ChannelingFastSimModel::GetRadiationModel ( )

inline

Definition at line 82 of file G4ChannelingFastSimModel.hh.

{return fBaierKatkov;}

Input() [1/2]

void G4ChannelingFastSimModel::Input ( const G4Material * crystal, const G4String & lattice )

inline

special functions

Definition at line 70 of file G4ChannelingFastSimModel.hh.

            {Input(crystal,lattice,"");}

Referenced by Input().

Input() [2/2]

void G4ChannelingFastSimModel::Input	(	const G4Material *	crystal,
		const G4String &	lattice,
		const G4String &	filePath )

Definition at line 531 of file G4ChannelingFastSimModel.cc.

{
   //initializing the class with containing all
   //the crystal material and crystal lattice data and
   //Channeling scattering and ionization processes
   fCrystalData = new G4ChannelingFastSimCrystalData();
   //setting all the crystal material and lattice data
   fCrystalData->SetMaterialProperties(crystal,lattice,filePath);
 
   //setting default low energy cuts for kinetic energy
   SetLowKineticEnergyLimit(1*GeV,"proton");
   SetLowKineticEnergyLimit(1*GeV,"anti_proton");
   SetLowKineticEnergyLimit(200*MeV,"e-");
   SetLowKineticEnergyLimit(200*MeV,"e+");
 
   //set the model high limit of the angle expressed in [Lindhard angle] units
   SetLindhardAngleNumberHighLimit(100.,"proton");
   SetLindhardAngleNumberHighLimit(100.,"anti_proton");
   SetLindhardAngleNumberHighLimit(100.,"e-");
   SetLindhardAngleNumberHighLimit(100.,"e+");
}

IsApplicable()

G4bool G4ChannelingFastSimModel::IsApplicable ( const G4ParticleDefinition & particleType )

overridevirtual

– IsApplicable

Implements G4VFastSimulationModel.

Definition at line 67 of file G4ChannelingFastSimModel.cc.

{
  return std::abs(particleType.GetPDGCharge())>DBL_EPSILON;
}

ModelTrigger()

G4bool G4ChannelingFastSimModel::ModelTrigger ( const G4FastTrack & fastTrack )

overridevirtual

– ModelTrigger

Implements G4VFastSimulationModel.

Definition at line 74 of file G4ChannelingFastSimModel.cc.

{
  //default output
  G4bool modelTrigger = false;
 
  G4int particleDefinitionID =
          fastTrack.GetPrimaryTrack()->GetParticleDefinition()->GetParticleDefinitionID();
  //kinetic energy
  G4double ekinetic = fastTrack.GetPrimaryTrack()->GetKineticEnergy();
 
  //energy cut, at the beginning, to not check everything else
  if(ekinetic > GetLowKineticEnergyLimit(particleDefinitionID))
  {
      //current logical volume
      G4LogicalVolume* crystallogic = fastTrack.GetEnvelopeLogicalVolume();
      fCrystalData->SetGeometryParameters(crystallogic);
 
      G4ThreeVector momentumDirection = fastTrack.GetPrimaryTrackLocalDirection();
      // the particle angle vs crystal plane or axis
      G4double angle = std::atan(momentumDirection.x()/momentumDirection.z());
      //recalculate angle into the lattice reference system
      angle = fCrystalData->
              AngleXFromBoxToLattice(angle,
                                     (fCrystalData->CoordinatesFromBoxToLattice(
                                          fastTrack.GetPrimaryTrackLocalPosition())).z());
      if (fCrystalData->GetModel()==2)
      {
          angle = std::sqrt(angle*angle+
                            std::pow(std::atan(momentumDirection.y()/
                                               momentumDirection.z()),2));
      }
 
      //particle mass
      G4double mass = fastTrack.GetPrimaryTrack()->GetParticleDefinition()->GetPDGMass();
      //particle total energy
      G4double etotal = mass + ekinetic;
      //particle charge
      G4double charge = fastTrack.GetPrimaryTrack()->
                        GetParticleDefinition()->GetPDGCharge();
 
      //Particle position
      G4ThreeVector xyz0 = fastTrack.GetPrimaryTrackLocalPosition();
      //Step estimate
      G4double dz0 = fCrystalData->GetMaxSimulationStep(etotal,mass,charge);
      xyz0 += 2*dz0*momentumDirection;//overestimated particle shift on the next step
                                      //in channeling
 
      //Applies the parameterisation not at the last step, only forward local direction
      //above low energy limit and below angular limit
 
      modelTrigger = (crystallogic->GetSolid()->
                      Inside(xyz0)==kInside) &&
                      momentumDirection.z()>0. &&
                      std::abs(angle) <
                      std::max(
                          GetLindhardAngleNumberHighLimit(particleDefinitionID) *
                          fCrystalData->GetLindhardAngle(etotal,
                                                         mass,
                                                         charge),
                          GetHighAngleLimit(particleDefinitionID));
  }
 
  return modelTrigger;
}

RadiationModelActivate()

void G4ChannelingFastSimModel::RadiationModelActivate

(

)

Definition at line 557 of file G4ChannelingFastSimModel.cc.

{
    fRad = true;
    //activate the Baier-Katkov radiation model
    fBaierKatkov = new G4BaierKatkov();
}

SetDefaultHighAngleLimit()

void G4ChannelingFastSimModel::SetDefaultHighAngleLimit ( G4double anglemax )

inline

Definition at line 101 of file G4ChannelingFastSimModel.hh.

           {fDefaultHighAngleLimit=anglemax;}

SetDefaultLindhardAngleNumberHighLimit()

void G4ChannelingFastSimModel::SetDefaultLindhardAngleNumberHighLimit ( G4double angleNumber )

inline

Definition at line 99 of file G4ChannelingFastSimModel.hh.

           {fDefaultLindhardAngleNumberHighLimit=angleNumber;}

SetDefaultLowKineticEnergyLimit()

void G4ChannelingFastSimModel::SetDefaultLowKineticEnergyLimit ( G4double ekinetic )

inline

Definition at line 97 of file G4ChannelingFastSimModel.hh.

           {fDefaultLowEnergyLimit=ekinetic;}

SetHighAngleLimit()

void G4ChannelingFastSimModel::SetHighAngleLimit ( G4double anglemax, const G4String & particleName )

inline

Definition at line 93 of file G4ChannelingFastSimModel.hh.

   {fHighAngleLimit[particleTable->FindParticle(particleName)->
                       GetParticleDefinitionID()] = anglemax;}

SetLindhardAngleNumberHighLimit()

void G4ChannelingFastSimModel::SetLindhardAngleNumberHighLimit ( G4double angleNumber, const G4String & particleName )

inline

Definition at line 90 of file G4ChannelingFastSimModel.hh.

   {fLindhardAngleNumberHighLimit[particleTable->FindParticle(particleName)->
              GetParticleDefinitionID()]=angleNumber;}

Referenced by Input().

SetLowKineticEnergyLimit()

void G4ChannelingFastSimModel::SetLowKineticEnergyLimit ( G4double ekinetic, const G4String & particleName )

inline

set cuts

Definition at line 87 of file G4ChannelingFastSimModel.hh.

   {fLowEnergyLimit[particleTable->FindParticle(particleName)->
              GetParticleDefinitionID()] = ekinetic;}

Referenced by Input().

SetMaxPhotonsProducedPerStep()

void G4ChannelingFastSimModel::SetMaxPhotonsProducedPerStep ( G4double nPhotons )

inline

get the maximal number of photons that can be produced per fastStep Caution: is redundant, if the radiation model is not activated

Definition at line 107 of file G4ChannelingFastSimModel.hh.

           {fMaxPhotonsProducedPerStep=nPhotons;}

---

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

• G4ChannelingFastSimModel.hh
• G4ChannelingFastSimModel.cc