G4BertiniEvaporation Class Reference

#include <G4BertiniEvaporation.hh>

Inheritance diagram for G4BertiniEvaporation:

Public Member Functions
	G4BertiniEvaporation ()
	~G4BertiniEvaporation ()
virtual G4FragmentVector *	BreakItUp (const G4Fragment &theNucleus)
G4FragmentVector *	BreakItUp (G4LayeredNucleus &nucleus)
void	setVerboseLevel (const G4int verbose)
Public Member Functions inherited from G4VEvaporation
	G4VEvaporation ()
virtual	~G4VEvaporation ()
virtual G4FragmentVector *	BreakItUp (const G4Fragment &theNucleus)=0
virtual void	Initialise ()
virtual void	SetPhotonEvaporation (G4VEvaporationChannel *ptr)
G4VEvaporationChannel *	GetPhotonEvaporation ()
void	SetOPTxs (G4int opt)
void	UseSICB (G4bool use)

Additional Inherited Members
Protected Attributes inherited from G4VEvaporation
G4VEvaporationChannel *	thePhotonEvaporation
G4int	OPTxs
G4bool	useSICB

Detailed Description

Definition at line 41 of file G4BertiniEvaporation.hh.

Constructor & Destructor Documentation

G4BertiniEvaporation()

G4BertiniEvaporation::G4BertiniEvaporation

(

)

Definition at line 48 of file G4BertiniEvaporation.cc.

{
    G4cout << "Info G4BertiniEvaporation: This is still very fresh code."<< G4endl;
    G4cout << "     G4BertiniEvaporation: feed-back for improvement is very wellcome."<< G4endl;
    verboseLevel = 0;
 
    channelVector.push_back( new G4BENeutronChannel );
    channelVector.push_back( new G4BEProtonChannel  );
    channelVector.push_back( new G4BEDeuteronChannel);
    channelVector.push_back( new G4BETritonChannel  );
    channelVector.push_back( new G4BEHe3Channel     );
    channelVector.push_back( new G4BEHe4Channel     );
}

~G4BertiniEvaporation()

G4BertiniEvaporation::~G4BertiniEvaporation

(

)

Definition at line 63 of file G4BertiniEvaporation.cc.

{
  while ( !channelVector.empty() )
    {
      delete channelVector.back();
      channelVector.pop_back();
    }
}

Member Function Documentation

BreakItUp() [1/2]

virtual G4FragmentVector * G4BertiniEvaporation::BreakItUp ( const G4Fragment &  theNucleus )

inlinevirtual

Implements G4VEvaporation.

Definition at line 48 of file G4BertiniEvaporation.hh.

  {
    G4LayeredNucleus aNuc( theNucleus.GetA(), theNucleus.GetZ() );
    aNuc.AddExcitationEnergy(theNucleus.GetExcitationEnergy());
    return BreakItUp(aNuc);
  }

Referenced by BreakItUp().

BreakItUp() [2/2]

G4FragmentVector * G4BertiniEvaporation::BreakItUp

(

G4LayeredNucleus &

nucleus

)

Definition at line 84 of file G4BertiniEvaporation.cc.

{
  G4int nucleusA;
  G4int nucleusZ;
  G4int i;
  G4double totalProbability;
  G4double excE;
  G4double newExcitation;
  G4double nucleusTotalMomentum;
  G4double nucleusKineticEnergy;
  G4double mRes; // Mass of residual nucleus.
  G4ThreeVector nucleusMomentumVector;
  G4DynamicParticle *pEmittedParticle;
  std::vector< G4DynamicParticle * > secondaryParticleVector;
  G4FragmentVector * result = new G4FragmentVector;
  
  // Read properties of the nucleus.
  nucleusA = nucleus.GetA_asInt();
  nucleusZ = nucleus.GetZ_asInt();
  excE                  = nucleus.GetEnergyDeposit();
  nucleusMomentumVector = nucleus.GetMomentum();
 
  // Move to CMS frame, save initial velocity of the nucleus to boostToLab vector.
  G4ThreeVector boostToLab( ( 1/G4NucleiProperties::GetNuclearMass( nucleusA, nucleusZ ) ) 
                * nucleusMomentumVector ); // xx mass ok?
 
  if ( verboseLevel >= 10 )
    G4cout << " G4BertiniEvaporation : initial kinematics : boostToLab vector = " << boostToLab << G4endl
       << "                     excitation energy  : " << excE<< G4endl;
 
  if ( nucleusA == 8 && nucleusZ == 4 )
    {
      splitBe8( excE, boostToLab, secondaryParticleVector );
      fillResult( secondaryParticleVector, result);
      return result;
    }
  
  // Initialize evaporation channels and calculate sum of emission
  // probabilities.
  std::vector< G4BertiniEvaporationChannel * >::iterator iChannel = channelVector.begin();  
  totalProbability = 0;
  for ( ; iChannel != channelVector.end() ; *iChannel++ )                  
     {
       ( *iChannel )->setNucleusA( nucleusA );
       ( *iChannel )->setNucleusZ( nucleusZ );
       ( *iChannel )->setExcitationEnergy( excE );
       ( *iChannel )->setPairingCorrection( 1 );
       ( *iChannel )->calculateProbability();
       totalProbability += ( *iChannel )->getProbability();
     }
 
  // If no evaporation process is possible, try without pairing energy.
  if ( totalProbability == 0 )
    {
      if ( verboseLevel >= 4 )
    G4cout << "G4BertiniEvaporation : no emission possible with pairing correction, trying without it" << G4endl;
      for ( iChannel = channelVector.begin() ; iChannel != channelVector.end() ; *iChannel++ )
    {
      ( *iChannel )->setPairingCorrection(0);
      ( *iChannel )->calculateProbability();
      totalProbability += ( *iChannel )->getProbability();
    }
      if ( verboseLevel >= 4 )
    G4cout << "                       probability without correction " << totalProbability << G4endl;
 
    }
 
  // Normal evaporation cycle follows.
  // Particles are evaporated till all probabilities are zero.
  while ( totalProbability > 0 )
    {
      G4BertiniEvaporationChannel *pSelectedChannel;
 
      // Sample active emission channel.
      G4double sampledProbability = G4UniformRand() * totalProbability;
      if ( verboseLevel >= 10 ) G4cout << "          RandomProb : " << sampledProbability << G4endl;
      for ( iChannel = channelVector.begin() ; iChannel != channelVector.end() ; *iChannel++ )
    {
      sampledProbability -= ( *iChannel )->getProbability();
      if ( sampledProbability < 0 ) break;
    }
      pSelectedChannel = ( *iChannel );
      if ( iChannel == channelVector.end() )
    throw G4HadronicException(__FILE__, __LINE__,  "G4BertiniEvaporation : Error while sampling evaporation particle" );
 
      if ( verboseLevel >= 4 ) 
    G4cout << "G4BertiniEvaporation : particle " << pSelectedChannel->getName() << " selected " << G4endl;
 
      // Max 10 tries to get a physically acceptable particle energy
      // in this emission channel.
      i = 0;
 
      do
    {
      pEmittedParticle = pSelectedChannel->emit();
      // This loop checks that particle with too large energy is not emitted.
      // CMS frame is considered in this loop. Nonrelativistic treatment. xxx
 
 
      const G4int zRes = nucleusZ - pSelectedChannel->getParticleZ(); 
      const G4int aRes  = nucleusA - pSelectedChannel->getParticleA(); 
      // const G4double eBind = G4NucleiProperties::GetBindingEnergy( aRes, zRes );  // Binding energy of the nucleus.
      mRes  = G4NucleiProperties::GetNuclearMass( aRes, zRes ); // Mass of the target nucleus
      //      In HETC88:
      //      eBind = Z * (-0.78244) + A * 8.36755 - cameron ( A , Z );
      //      mRes =  zRes * 938.79304 + ( aRes - zRes ) * 939.57548 - eBind; 
      //      where cameron ( A, Z ) is the mass excess by Cameron, see Canadian Journal of Physics,
      //      vol. 35, 1957, p.1022
 
      nucleusTotalMomentum = pEmittedParticle->GetTotalMomentum(); // CMS frame
      nucleusKineticEnergy = std::pow( nucleusTotalMomentum, 2 ) / ( 2 * mRes );
      newExcitation = excE - pEmittedParticle->GetKineticEnergy() - nucleusKineticEnergy - pSelectedChannel->getQ();
 
      if ( verboseLevel >= 10)
        G4cout << "G4BertiniEvaporation : Kinematics " << G4endl
           << "                    part kinetic E in CMS = " 
           << pEmittedParticle->GetKineticEnergy() << G4endl
           << "                    new excitation E      =  " 
           << newExcitation << G4endl;
 
      i++;
      if ( !( newExcitation > 0 ) && i < 10) delete pEmittedParticle;
    } while ( !( newExcitation > 0 )  &&  i < 10 );
 
      if ( i >= 10 ) 
    {
      // No appropriate particle energy found.
      // Set probability of this channel to zero 
      // and try to sample another particle on the 
      // next round.
      delete pEmittedParticle;
 
      if ( verboseLevel >= 4 )
        G4cout << "G4BertiniEvaporation : No appropriate energy for particle " 
           << pSelectedChannel->getName() << " found." << G4endl;
 
      pSelectedChannel->setProbability( 0 );
 
      totalProbability = 0;
      for ( ; iChannel != channelVector.end() ; *iChannel++ )                  
        totalProbability += ( *iChannel )->getProbability();
    } // Treatment of physically unacceptable particle ends.
      else
    {
      // Good particle found. 
 
      // Transform particle properties to lab frame and save it.
      G4LorentzVector particleFourVector = pEmittedParticle->Get4Momentum();
      G4LorentzVector nucleusFourVector(  - pEmittedParticle->GetMomentum(), // CMS Frame.
                          nucleusKineticEnergy + mRes ); // Total energy.
      particleFourVector.boost( boostToLab );
      nucleusFourVector.boost(  boostToLab );
      G4DynamicParticle *pPartLab  = new G4DynamicParticle( pEmittedParticle->GetDefinition(),
                                particleFourVector.e(), // Total energy
                                particleFourVector.vect() ); // momentum vector
      secondaryParticleVector.push_back( pPartLab );
 
      // Update residual nucleus and boostToLab vector.
      nucleusA = nucleusA - pSelectedChannel->getParticleA();
      nucleusZ = nucleusZ - pSelectedChannel->getParticleZ();
      excE = newExcitation;
      boostToLab = G4ThreeVector ( ( 1/mRes ) * nucleusFourVector.vect() );
 
      if ( verboseLevel >= 10 )
        G4cout << "  particle mom in cms " << pEmittedParticle->GetMomentum() << G4endl 
           << "  particle mom in cms " << pEmittedParticle->GetTotalMomentum() << G4endl 
           << "          Etot in cms " << pEmittedParticle->GetTotalEnergy() << G4endl
           << "          Ekin in cms " << pEmittedParticle->GetKineticEnergy() << G4endl
           << "        particle mass " << pEmittedParticle->GetMass() << G4endl
           << "          boost  vect " << boostToLab << G4endl
           << "          boosted 4v  " << particleFourVector << G4endl
           << "          nucleus 4v  " << nucleusFourVector << G4endl
           << "          nucl cm mom " << nucleusTotalMomentum << G4endl
           << " particle k.e. in lab " << pPartLab->GetKineticEnergy() << G4endl
           << "     new boost vector " << boostToLab << G4endl;
 
      delete pEmittedParticle;
 
      // If the residual nucleus is Be8, split it.
      if ( nucleusA == 8 && nucleusZ == 4 )
        {
          splitBe8( excE, boostToLab, secondaryParticleVector );
              fillResult( secondaryParticleVector, result);
          return result;
        }
 
      // Update evaporation channels and 
      // total emission probability.
      totalProbability = 0;
      for ( iChannel = channelVector.begin() ; iChannel != channelVector.end() ; *iChannel++ )  
        {
          ( *iChannel )->setNucleusA( nucleusA ); 
          ( *iChannel )->setNucleusZ( nucleusZ );
          ( *iChannel )->setExcitationEnergy( excE ); 
          ( *iChannel )->calculateProbability();
          totalProbability += ( *iChannel )->getProbability();
        }
    } // Treatment of physically acceptable particle ends.
    } // End of evaporation cycle.
 
  // Gamma de-excitation.
  G4BEGammaDeexcitation *pGammaChannel = new G4BEGammaDeexcitation;
  pGammaChannel->setNucleusA( nucleusA );
  pGammaChannel->setNucleusZ( nucleusZ );
  pGammaChannel->setVerboseLevel( verboseLevel );
  pGammaChannel->setExcitationEnergy( excE );
 
  // Emit equally sampled photons until all the excitation energy is
  // used; the last photon gets the remaining de-excitation energy.
  // Change of momentum of the nucleus is neglected.
  G4double gammaEnergy;
  G4double totalDeExcEnergy = 0;
 
  while ( excE > 0 )
    {
      pEmittedParticle = pGammaChannel->emit();
      gammaEnergy = pEmittedParticle->GetKineticEnergy();
      
      if ( ( totalDeExcEnergy + gammaEnergy ) > excE ) 
    {
      // All the remaining energy is used here.
      gammaEnergy = excE - totalDeExcEnergy;
      excE = 0;
    }
 
      totalDeExcEnergy += gammaEnergy; 
 
      if ( verboseLevel >= 10 )
    G4cout << " G4BertiniEvaporation : gamma de-excitation, g of " << gammaEnergy << " MeV " << G4endl;
      
      G4LorentzVector gammaFourVector( pEmittedParticle->GetMomentum(),
                       pEmittedParticle->GetKineticEnergy() );
      gammaFourVector.boost( boostToLab );
      pEmittedParticle->SetKineticEnergy( gammaFourVector.e() );
      pEmittedParticle->SetMomentumDirection( gammaFourVector.vect().unit() );
 
      secondaryParticleVector.push_back( pEmittedParticle );
 
      } 
 
  delete pGammaChannel;
 
  // Residual nucleus is not returned.
 
  fillResult( secondaryParticleVector, result);
  
  return result;
}

setVerboseLevel()

void G4BertiniEvaporation::setVerboseLevel

(

const G4int

verbose

)

Definition at line 73 of file G4BertiniEvaporation.cc.

{
  verboseLevel = verbose;
 
  // Update verbose level to all evaporation channels.
  std::vector< G4BertiniEvaporationChannel * >::iterator iChannel = channelVector.begin();  
  for ( ; iChannel != channelVector.end() ; *iChannel++ )                  
    ( *iChannel )->setVerboseLevel( verboseLevel );
}

---

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

• G4BertiniEvaporation.hh
• G4BertiniEvaporation.cc