G4S2Light Class Reference

#include <G4S2Light.hh>

Inheritance diagram for G4S2Light:

Public Member Functions
	G4S2Light (const G4String &processName="S2", G4ProcessType type=fElectromagnetic)
	~G4S2Light ()
G4bool	IsApplicable (const G4ParticleDefinition &aParticleType)
G4double	GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *)
G4double	GetMeanLifeTime (const G4Track &aTrack, G4ForceCondition *)
G4VParticleChange *	PostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
G4VParticleChange *	AtRestDoIt (const G4Track &aTrack, const G4Step &aStep)
void	SetTrackSecondariesFirst (const G4bool state)
G4bool	GetTrackSecondariesFirst () const
void	SetScintillationYieldFactor (const G4double yieldfactor)
G4double	GetScintillationYieldFactor () const
void	SetScintillationExcitationRatio (const G4double excitationratio)
G4double	GetScintillationExcitationRatio () const

Protected Attributes
G4bool	fTrackSecondariesFirst
G4double	YieldFactor
G4double	ExcitationRatio

Detailed Description

Definition at line 22 of file G4S2Light.hh.

Constructor & Destructor Documentation

G4S2Light::G4S2Light	(	const G4String &	processName = "S2",
		G4ProcessType	type = fElectromagnetic
	)

Definition at line 66 of file G4S2Light.cc.

      : G4VRestDiscreteProcess(processName, type)
{
  thr[18] = 0.190; E_eV[18] = 9.7; MolarMass[18] = 39.948;
  tau1[18] = 6*CLHEP::ns; tau3[18] = 1600*CLHEP::ns; ConvertEff[18]=0.7857;
  thr[54] = 0.474; E_eV[54] = 7.143; MolarMass[54] = 131.293;
  ConvertEff[54] = 1.01; //50% when LXe TPC is dirty
        //luxManager = LUXSimManager::GetManager();
        SetProcessSubType(fScintillation);
        fTrackSecondariesFirst = false;
        if (verboseLevel>0) {
          G4cout << GetProcessName() << " is created " << G4endl;
        }
}

G4S2Light::~G4S2Light

(

)

Definition at line 82 of file G4S2Light.cc.

{} //destructor

Member Function Documentation

G4VParticleChange * G4S2Light::AtRestDoIt	(	const G4Track &	aTrack,
		const G4Step &	aStep
	)

Definition at line 85 of file G4S2Light.cc.

{
  return G4S2Light::PostStepDoIt(aTrack, aStep);
}

G4double G4S2Light::GetMeanFreePath	(	const G4Track &	aTrack,
		G4double	,
		G4ForceCondition *	condition
	)

Definition at line 258 of file G4S2Light.cc.

{
  const G4Material* aMaterial = aTrack.GetMaterial();
  if ( !aMaterial ) return DBL_MIN;
  const G4ElementVector* theElementVector = aMaterial->GetElementVector();
  G4Element* Element = (*theElementVector)[0]; G4int z1;
  if ( Element ) z1 = (G4int)(Element->GetZ()); else z1 = -1;
  if ( (z1==2 || z1==10 || z1==18 || z1==36 || z1==54) &&
       aMaterial->GetState() == kStateLiquid ) {
    if ( aTrack.GetPosition()[2] < (GAT-5*CLHEP::mm) || aTrack.GetWeight() == 2 )
      return E_PURITY;
    else return 0.000*CLHEP::mm;
  }
  else if ((z1==2 || z1==10 || z1==18 || z1==36 || z1==54) &&
   aMaterial->GetState()== kStateGas && aMaterial->GetDensity()>0.0*(CLHEP::g/CLHEP::cm3)) {
    G4MaterialPropertiesTable* aMaterialPropertiesTable =
      aMaterial->GetMaterialPropertiesTable();
    if(!aMaterialPropertiesTable) return DBL_MIN;
    if ( aTrack.GetPosition()[2] < (GAT-5*CLHEP::mm) ) return DBL_MAX;
    G4double ElectricField = fabs(aMaterialPropertiesTable->
                                  GetConstProperty("ELECTRICFIELDANODE"));
    ElectricField = ElectricField/(CLHEP::volt/CLHEP::cm);
    G4double mfp, nDensity = 
      (aMaterial->GetDensity()/(CLHEP::g/CLHEP::cm3))/(MolarMass[z1])*AVO;
    if ( (1/E_eV[z1])*(ElectricField/nDensity)*1e17 == thr[z1] ) {
        mfp = 0*CLHEP::cm;
    }
    else { //denominator of formula OK
      mfp = 1 / 
        (nDensity*1e-17*((1/E_eV[z1])*(ElectricField/nDensity)*1e17-thr[z1]));
    }
    //equation 8 within Monteiro's paper JINST 2007 (where 140=1/7.14.. eV)
    mfp*=CLHEP::cm/ConvertEff[z1];
    if(mfp<0)mfp=0;
    if(mfp>GASGAP/2.7)mfp=GASGAP/exp(1);
    return mfp; //mean free path for 1 photon
  }
  else {
    *condition = NotForced;
    const G4Material* aMaterial = aTrack.GetMaterial();
    if ( aMaterial->GetDensity() == GRID_DENSITY )
      return DBL_MAX; //pass electrons through grid wires
    return 0*CLHEP::nm; //kill elsewhere
  }
}

G4double G4S2Light::GetMeanLifeTime	(	const G4Track &	aTrack,
		G4ForceCondition *	condition
	)

Definition at line 308 of file G4S2Light.cc.

{
  *condition = InActivated;
  return DBL_MAX;
}

G4double G4S2Light::GetScintillationExcitationRatio ( ) const

inline

Definition at line 128 of file G4S2Light.hh.

{
        return ExcitationRatio;
}

G4double G4S2Light::GetScintillationYieldFactor ( ) const

inline

Definition at line 116 of file G4S2Light.hh.

{
        return YieldFactor;
}

G4bool G4S2Light::GetTrackSecondariesFirst ( ) const

inline

Definition at line 104 of file G4S2Light.hh.

{
        return fTrackSecondariesFirst;
}

G4bool G4S2Light::IsApplicable ( const G4ParticleDefinition &  aParticleType )

inline

Definition at line 90 of file G4S2Light.hh.

{
       if (aParticleType.GetParticleName() != "thermalelectron") return false;
       
       return true;
}

G4VParticleChange * G4S2Light::PostStepDoIt	(	const G4Track &	aTrack,
		const G4Step &	aStep
	)

Definition at line 91 of file G4S2Light.cc.

{
// The function where all the action happens! Not yet well commented, as this
// is but an alpha version of the electroluminescence code
        YieldFactor = 1.000;
        aParticleChange.Initialize(aTrack);
        const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
        const G4Material* aMaterial = aTrack.GetMaterial();
        if ( !aMaterial ) {
          aParticleChange.ProposeTrackStatus(fStopAndKill);
          return G4VRestDiscreteProcess::PostStepDoIt(aTrack, aStep);
        }
        G4MaterialPropertiesTable* aMaterialPropertiesTable = 
          aMaterial->GetMaterialPropertiesTable();
        if ( !aMaterialPropertiesTable ) {
          aParticleChange.ProposeTrackStatus(fStopAndKill);
          return G4VRestDiscreteProcess::PostStepDoIt(aTrack, aStep);
        }
        const G4ElementVector* theElementVector =
          aMaterial->GetElementVector();
        G4Element* Element = (*theElementVector)[0]; G4int z1;
        if ( Element ) z1 = (G4int)(Element->GetZ()); else z1 = -1;
        G4double z_anode = BORDER + GASGAP;
        G4double z_start = BORDER;
        //if(luxManager->GetGasRun()) z_start=GAT;
        tau1[54] = G4RandGauss::shoot(5.18*CLHEP::ns,1.55*CLHEP::ns);
        tau3[54] = G4RandGauss::shoot(100.1*CLHEP::ns,7.9*CLHEP::ns);
        
        G4StepPoint* pPreStepPoint = aStep.GetPreStepPoint();
        G4StepPoint* pPostStepPoint = aStep.GetPostStepPoint();
        G4ThreeVector x0 = pPreStepPoint->GetPosition();
        G4ThreeVector x1 = pPostStepPoint->GetPosition();
        G4double      t1 = pPostStepPoint->GetGlobalTime();
        G4double Density = aMaterial->GetDensity()/(CLHEP::g/CLHEP::cm3);
        G4double nDensity = (Density/MolarMass[z1])*AVO;
        G4int Phase = aMaterial->GetState();
        
        if ( Phase == kStateLiquid && x0[2] > (GAT-5*CLHEP::mm) && GAT != 0 ) {
          aParticleChange.ProposeEnergy(GetLiquidElectronDriftSpeed(
          aMaterial->GetTemperature(),fabs(aMaterialPropertiesTable->
                                           GetConstProperty("ELECTRICFIELDSURFACE")/(CLHEP::volt/CLHEP::cm)),MillerDriftSpeed,z1));
          aParticleChange.ProposeWeight(2.0);
          return G4VRestDiscreteProcess::PostStepDoIt(aTrack, aStep);
        }
        
        // absorb the electron if you are in liquid and there is non-perfect
        // purity, you turned S2 off with YieldFactor, the e- is above the
        // anode in the gas already, or it's not present in a noble element
        if ( Phase == kStateLiquid || !YieldFactor || x0[2] > z_anode ||
             x1[2] > z_anode || fabs(x1[2]-x0[2]) < 1e-7*CLHEP::nm ||
             (z1!=2 && z1!=10 && z1!=18 && z1!=36 && z1!=54) ) {
          G4double KE = aParticle->GetKineticEnergy();
          aParticleChange.ProposeTrackStatus(fStopAndKill);
          aParticleChange.ProposeEnergy(0.);
          aParticleChange.ProposeLocalEnergyDeposit(KE);
          return G4VRestDiscreteProcess::PostStepDoIt(aTrack, aStep);
        }
        if ( x0[2] <= z_start && x1[2] <= z_start )
          return G4VRestDiscreteProcess::PostStepDoIt(aTrack, aStep);
        
        G4double ElectricField;
        if(!WIN && !TOP && !ANE && !SRF && !GAT && !CTH && !BOT && !PMT)
          ElectricField = aMaterialPropertiesTable->
            GetConstProperty("ELECTRICFIELD");
        else
          ElectricField = aMaterialPropertiesTable->
            GetConstProperty("ELECTRICFIELDANODE");
        ElectricField = 
            fabs(ElectricField/(CLHEP::kilovolt/CLHEP::cm));
        
        if ( Density > 0. && Phase == kStateGas && (x0[2] > z_start ||
             x1[2] > z_start) && aParticleChange.GetWeight() >= 1.0 ) {
          aParticleChange.ProposeWeight(0);
          G4double ExtractEff = -0.0012052 + //Aprile 2004 IEEE No.5
            0.1638*ElectricField-0.0063782*pow(ElectricField,2.);
          if ( ElectricField > 10.0 ) ExtractEff = 1.00;
          if ( ElectricField <= 0.0 ) ExtractEff = 0.00;
          if ( ExtractEff < 0 ) ExtractEff = 0; 
          if ( ExtractEff > 1 ) ExtractEff = 1;
          //if ( luxManager->GetGasRun() ) ExtractEff = 1;
          if (G4UniformRand() > ExtractEff || !ElectricField ||
              (1/E_eV[z1])*1e17*((ElectricField*1e3)/nDensity)<=thr[z1]) {
            aParticleChange.ProposeTrackStatus(fStopAndKill);
            return G4VRestDiscreteProcess::PostStepDoIt(aTrack, aStep);
          }
          else {
            G4double eDrift =
              GetGasElectronDriftSpeed(1000.*ElectricField,nDensity);
            aParticleChange.ProposeEnergy(eDrift);
          }
        }
        
        if ( G4UniformRand () >= QE_EFF )
          return G4VRestDiscreteProcess::PostStepDoIt(aTrack, aStep);
        aParticleChange.SetNumberOfSecondaries(G4int(floor(YieldFactor)));
        if ( verboseLevel > 0 ) {
          G4cout << "\n Exiting from G4S2Light::DoIt -- "
                 << "NumberOfSecondaries = "
                 << aParticleChange.GetNumberOfSecondaries() << G4endl;
        }
        
        // start particle creation
        G4double sampledEnergy;
        G4DynamicParticle* aQuantum;
        
        // Generate random direction
        G4double cost = 1. - 2.*G4UniformRand();
        G4double sint = std::sqrt((1.-cost)*(1.+cost));
        G4double phi = CLHEP::twopi*G4UniformRand();
        G4double sinp = std::sin(phi);
        G4double cosp = std::cos(phi);
        G4double px = sint*cosp; G4double py = sint*sinp;
        G4double pz = cost;
        
        // Create momentum direction vector
        G4ParticleMomentum photonMomentum(px, py, pz);
        
        // Determine polarization of new photon
        G4double sx = cost*cosp; G4double sy = cost*sinp;
        G4double sz = -sint;
        G4ThreeVector photonPolarization(sx, sy, sz);
        G4ThreeVector perp = photonMomentum.cross(photonPolarization);
        phi = CLHEP::twopi*G4UniformRand();
        sinp = std::sin(phi); cosp = std::cos(phi);
        photonPolarization = cosp * photonPolarization + sinp * perp;
        photonPolarization = photonPolarization.unit();
        
        // Generate a new photon:
        sampledEnergy = G4RandGauss::shoot(E_eV[z1]*CLHEP::eV,0.2*CLHEP::eV);
        if ( z1==54 && sampledEnergy>8.5*CLHEP::eV ) sampledEnergy = 8.5*CLHEP::eV;
        aQuantum = 
          new G4DynamicParticle(G4OpticalPhoton::OpticalPhoton(),
                                photonMomentum);
        aQuantum->SetPolarization(photonPolarization.x(),
                                  photonPolarization.y(),
                                  photonPolarization.z());
        
        //assign energy to make particle real
        aQuantum->SetKineticEnergy(sampledEnergy);
        if (verboseLevel>1)
          G4cout << "sampledEnergy = " << sampledEnergy << G4endl;
        
        // electroluminesence emission time distribution
        G4double aSecondaryTime = t1, 
          SingTripRatio=.1; //guess: revisit
        if(G4UniformRand()<SingTripRatio/(1+SingTripRatio))
          aSecondaryTime -= tau1[z1]*log(G4UniformRand());
        else aSecondaryTime -= 
               tau3[z1]*log(G4UniformRand());
        
        //the position of the new secondary particle
        if ( x1[2] < BORDER + GASGAP / 2. )
          x1[2] += 0.001*CLHEP::mm;
        if ( x1[2] > BORDER + GASGAP / 2. )
          x1[2] -= 0.001*CLHEP::mm;
        G4ThreeVector aSecondaryPosition = x1;
        
        // GEANT4 business: stuff you need to make a new track
        G4Track* aSecondaryTrack = 
          new G4Track(aQuantum,aSecondaryTime,aSecondaryPosition);
        aParticleChange.AddSecondary(aSecondaryTrack);
        
        return G4VRestDiscreteProcess::PostStepDoIt(aTrack, aStep);
}

void G4S2Light::SetScintillationExcitationRatio ( const G4double  excitationratio )

inline

Definition at line 122 of file G4S2Light.hh.

{
        ExcitationRatio = excitationratio;
}

void G4S2Light::SetScintillationYieldFactor ( const G4double  yieldfactor )

inline

Definition at line 110 of file G4S2Light.hh.

{
  YieldFactor = yieldfactor;
}

void G4S2Light::SetTrackSecondariesFirst ( const G4bool  state )

inline

Definition at line 98 of file G4S2Light.hh.

{ 
        fTrackSecondariesFirst = state;
}

Member Data Documentation

G4double G4S2Light::ExcitationRatio

protected

Definition at line 78 of file G4S2Light.hh.

G4bool G4S2Light::fTrackSecondariesFirst

protected

Definition at line 76 of file G4S2Light.hh.

G4double G4S2Light::YieldFactor

protected

Definition at line 77 of file G4S2Light.hh.

---

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

• edep-sim/src/NESTVersion098/G4S2Light.hh
• edep-sim/src/NESTVersion098/G4S2Light.cc