d1/d95/_g4_hadron_elastic_process_8cc_source.html

 //

 // ********************************************************************

 // * License and Disclaimer                                           *

 // *                                                                  *

 // * The  Geant4 software  is  copyright of the Copyright Holders  of *

 // * the Geant4 Collaboration.  It is provided  under  the terms  and *

 // * conditions of the Geant4 Software License,  included in the file *

 // * LICENSE and available at  http://cern.ch/geant4/license .  These *

 // * include a list of copyright holders.                             *

 // *                                                                  *

 // * Neither the authors of this software system, nor their employing *

 // * institutes,nor the agencies providing financial support for this *

 // * work  make  any representation or  warranty, express or implied, *

 // * regarding  this  software system or assume any liability for its *

 // * use.  Please see the license in the file  LICENSE  and URL above *

 // * for the full disclaimer and the limitation of liability.         *

 // *                                                                  *

 // * This  code  implementation is the result of  the  scientific and *

 // * technical work of the GEANT4 collaboration.                      *

 // * By using,  copying,  modifying or  distributing the software (or *

 // * any work based  on the software)  you  agree  to acknowledge its *

 // * use  in  resulting  scientific  publications,  and indicate your *

 // * acceptance of all terms of the Geant4 Software license.          *

 // ********************************************************************

 //

 // $Id: G4HadronElasticProcess.cc 66499 2012-12-19 09:16:35Z gcosmo $

 //

 // Geant4 Hadron Elastic Scattering Process

 //

 // Created 26 July 2012 V.Ivanchenko from G4WHadronElasticProcess

 //

 // Modified:

 // 14-Sep-12 M.Kelsey -- Pass subType code to base ctor


 #include <iostream>

 #include <typeinfo>


 #include "G4HadronElasticProcess.hh"

 #include "G4SystemOfUnits.hh"

 #include "G4Nucleus.hh"

 #include "G4ProcessManager.hh"

 #include "G4CrossSectionDataStore.hh"

 #include "G4HadronElasticDataSet.hh"

 #include "G4ProductionCutsTable.hh"

 #include "G4HadronicException.hh"

 #include "G4HadronicDeprecate.hh"


 G4HadronElasticProcess::G4HadronElasticProcess(const G4String& pName)

   : G4HadronicProcess(pName, fHadronElastic), isInitialised(false)

 {

   AddDataSet(new G4HadronElasticDataSet);

   lowestEnergy = 1.*keV;

 }


 G4HadronElasticProcess::~G4HadronElasticProcess()

 {}


 void G4HadronElasticProcess::Description() const

 {

   char* dirName = getenv("G4PhysListDocDir");

   if (dirName) {

     std::ofstream outFile;

     G4String outFileName = GetProcessName() + ".html";

     G4String pathName = G4String(dirName) + "/" + outFileName;

     outFile.open(pathName);

     outFile << "<html>\n";

     outFile << "<head>\n";


     outFile << "<title>Description of G4HadronElasticProcess</title>\n";

     outFile << "</head>\n";

     outFile << "<body>\n";


     outFile << "G4HadronElasticProcess handles the elastic scattering of\n"

             << "hadrons by invoking one or more hadronic models and one or\n"

             << "more hadronic cross sections.\n";


     outFile << "</body>\n";

     outFile << "</html>\n";

     outFile.close();

   }

 }


 G4VParticleChange*

 G4HadronElasticProcess::PostStepDoIt(const G4Track& track,

                                      const G4Step& /*step*/)

 {

   theTotalResult->Clear();

   theTotalResult->Initialize(track);

   G4double weight = track.GetWeight();

   theTotalResult->ProposeWeight(weight);


   // For elastic scattering, _any_ result is considered an interaction

   ClearNumberOfInteractionLengthLeft();


   G4double kineticEnergy = track.GetKineticEnergy();

   const G4DynamicParticle* dynParticle = track.GetDynamicParticle();

   const G4ParticleDefinition* part = dynParticle->GetDefinition();


   // NOTE:  Very low energy scatters were causing numerical (FPE) errors

   //        in earlier releases; these limits have not been changed since.

   if (kineticEnergy <= lowestEnergy)   return theTotalResult;


   G4Material* material = track.GetMaterial();

   G4Nucleus* targNucleus = GetTargetNucleusPointer();


   // Select element

   G4Element* elm = 0;

   try

     {

       elm = GetCrossSectionDataStore()->SampleZandA(dynParticle, material,

                                                     *targNucleus);

     }

   catch(G4HadronicException & aR)

     {

       G4ExceptionDescription ed;

       aR.Report(ed);

       DumpState(track,"SampleZandA",ed);

       ed << " PostStepDoIt failed on element selection" << G4endl;

       G4Exception("G4HadronElasticProcess::PostStepDoIt", "had003",

                   FatalException, ed);

     }

   G4HadronicInteraction* hadi = 0;

   try

     {

       hadi = ChooseHadronicInteraction( kineticEnergy, material, elm );

     }

   catch(G4HadronicException & aE)

     {

       G4ExceptionDescription ed;

       aE.Report(ed);

       ed << "Target element "<< elm->GetName()<<"  Z= "

          << targNucleus->GetZ_asInt() << "  A= "

          << targNucleus->GetA_asInt() << G4endl;

       DumpState(track,"ChooseHadronicInteraction",ed);

       ed << " No HadronicInteraction found out" << G4endl;

       G4Exception("G4HadronElasticProcess::PostStepDoIt", "had005",

                   FatalException, ed);

     }


   size_t idx = track.GetMaterialCutsCouple()->GetIndex();

   G4double tcut = (*(G4ProductionCutsTable::GetProductionCutsTable()

                      ->GetEnergyCutsVector(3)))[idx];

   hadi->SetRecoilEnergyThreshold(tcut);


   // Initialize the hadronic projectile from the track

   //  G4cout << "track " << track.GetDynamicParticle()->Get4Momentum()<<G4endl;

   G4HadProjectile theProj(track);

   if(verboseLevel>1) {

     G4cout << "G4HadronElasticProcess::PostStepDoIt for "

            << part->GetParticleName()

            << " in " << material->GetName()

            << " Target Z= " << targNucleus->GetZ_asInt()

            << " A= " << targNucleus->GetA_asInt() << G4endl;

   }


   G4HadFinalState* result = 0;

   try

     {

       result = hadi->ApplyYourself( theProj, *targNucleus);

     }

   catch(G4HadronicException aR)

     {

       G4ExceptionDescription ed;

       aR.Report(ed);

       ed << "Call for " << hadi->GetModelName() << G4endl;

       ed << "Target element "<< elm->GetName()<<"  Z= "

          << targNucleus->GetZ_asInt()

          << "  A= " << targNucleus->GetA_asInt() << G4endl;

       DumpState(track,"ApplyYourself",ed);

       ed << " ApplyYourself failed" << G4endl;

       G4Exception("G4HadronElasticProcess::PostStepDoIt", "had006",

                   FatalException, ed);

     }


   // Check the result for catastrophic energy non-conservation

   // cannot be applied because is not guranteed that recoil

   // nucleus is created

   // result = CheckResult(theProj, targNucleus, result);


   // directions

   G4ThreeVector indir = track.GetMomentumDirection();

   G4double phi = CLHEP::twopi*G4UniformRand();

   G4ThreeVector it(0., 0., 1.);

   G4ThreeVector outdir = result->GetMomentumChange();


   if(verboseLevel>1) {

     G4cout << "Efin= " << result->GetEnergyChange()

            << " de= " << result->GetLocalEnergyDeposit()

            << " nsec= " << result->GetNumberOfSecondaries()

            << " dir= " << outdir

            << G4endl;

   }


   // energies

   G4double edep = result->GetLocalEnergyDeposit();

   G4double efinal = result->GetEnergyChange();

   if(efinal < 0.0) { efinal = 0.0; }

   if(edep < 0.0)   { edep = 0.0; }


   // NOTE:  Very low energy scatters were causing numerical (FPE) errors

   //        in earlier releases; these limits have not been changed since.

   if(efinal <= lowestEnergy) {

     edep += efinal;

     efinal = 0.0;

   }


   // primary change

   theTotalResult->ProposeEnergy(efinal);


   G4TrackStatus status = track.GetTrackStatus();

   if(efinal > 0.0) {

     outdir.rotate(phi, it);

     outdir.rotateUz(indir);

     theTotalResult->ProposeMomentumDirection(outdir);

   } else {

     if(part->GetProcessManager()->GetAtRestProcessVector()->size() > 0)

          { status = fStopButAlive; }

     else { status = fStopAndKill; }

     theTotalResult->ProposeTrackStatus(status);

   }


   //G4cout << "Efinal= " << efinal << "  TrackStatus= " << status << G4endl;


   theTotalResult->SetNumberOfSecondaries(0);


   // recoil

   if(result->GetNumberOfSecondaries() > 0) {

     G4DynamicParticle* p = result->GetSecondary(0)->GetParticle();


     if(p->GetKineticEnergy() > tcut) {

       theTotalResult->SetNumberOfSecondaries(1);

       G4ThreeVector pdir = p->GetMomentumDirection();

       // G4cout << "recoil " << pdir << G4endl;

       pdir.rotate(phi, it);

       pdir.rotateUz(indir);

       // G4cout << "recoil rotated " << pdir << G4endl;

       p->SetMomentumDirection(pdir);


       // in elastic scattering time and weight are not changed

       G4Track* t = new G4Track(p, track.GetGlobalTime(),

                                track.GetPosition());

       t->SetWeight(weight);

       t->SetTouchableHandle(track.GetTouchableHandle());

       theTotalResult->AddSecondary(t);


     } else {

       edep += p->GetKineticEnergy();

       delete p;

     }

   }

   theTotalResult->ProposeLocalEnergyDeposit(edep);

   theTotalResult->ProposeNonIonizingEnergyDeposit(edep);

   result->Clear();


   return theTotalResult;

 }


 void

 G4HadronElasticProcess::PreparePhysicsTable(const G4ParticleDefinition& part)

 {

   if(!isInitialised) {

     isInitialised = true;

     if(G4Neutron::Neutron() == &part) { lowestEnergy = 1.e-6*eV; }

   }

   G4HadronicProcess::PreparePhysicsTable(part);

 }


 void G4HadronElasticProcess::SetLowestEnergy(G4double val)

 {

   lowestEnergy = val;

 }


 void

 G4HadronElasticProcess::SetLowestEnergyNeutron(G4double val)

 {

   lowestEnergy = val;

   G4HadronicDeprecate("G4HadronElasticProcess::SetLowestEnergyNeutron()");

 }


G4Nucleus::GetA_asInt
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109

G4HadFinalState
Definition: G4HadFinalState.hh:45

G4ProductionCutsTable::GetEnergyCutsVector
const std::vector< G4double > * GetEnergyCutsVector(size_t pcIdx) const
Definition: G4ProductionCutsTable.hh:250

G4HadronElasticProcess::G4HadronElasticProcess
G4HadronElasticProcess(const G4String &procName="hadElastic")
Definition: G4HadronElasticProcess.cc:48

G4Nucleus
Definition: G4Nucleus.hh:50

G4VProcess::verboseLevel
G4int verboseLevel
Definition: G4VProcess.hh:368

G4ExceptionDescription
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76

G4HadFinalState::GetSecondary
G4HadSecondary * GetSecondary(size_t i)
Definition: G4HadFinalState.cc:109

G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const

G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42

outFile
std::ofstream outFile
Definition: GammaRayTel.cc:68

fStopButAlive
Definition: G4TrackStatus.hh:54

G4Track::GetDynamicParticle
const G4DynamicParticle * GetDynamicParticle() const

fHadronElastic
Definition: G4HadronicProcessType.hh:44

G4HadronElasticProcess::lowestEnergy
G4double lowestEnergy
Definition: G4HadronElasticProcess.hh:77

G4Material::GetName
const G4String & GetName() const
Definition: G4Material.hh:176

G4HadFinalState::Clear
void Clear()
Definition: G4HadFinalState.cc:89

G4Track::GetPosition
const G4ThreeVector & GetPosition() const

G4Track::GetTrackStatus
G4TrackStatus GetTrackStatus() const

G4Material
Definition: G4Material.hh:118

G4Track::GetMaterialCutsCouple
const G4MaterialCutsCouple * GetMaterialCutsCouple() const

G4DynamicParticle
Definition: G4DynamicParticle.hh:73

G4Element
Definition: G4Element.hh:97

G4HadProjectile
Definition: G4HadProjectile.hh:39

G4HadronicProcess::PreparePhysicsTable
virtual void PreparePhysicsTable(const G4ParticleDefinition &)
Definition: G4HadronicProcess.cc:158

G4Track::SetTouchableHandle
void SetTouchableHandle(const G4TouchableHandle &apValue)

G4HadronicException.hh

G4HadronElasticProcess::SetLowestEnergy
virtual void SetLowestEnergy(G4double)
Definition: G4HadronElasticProcess.cc:269

G4HadronElasticProcess::~G4HadronElasticProcess
virtual ~G4HadronElasticProcess()
Definition: G4HadronElasticProcess.cc:55

G4DynamicParticle::SetMomentumDirection
void SetMomentumDirection(const G4ThreeVector &aDirection)

G4HadFinalState::GetEnergyChange
G4double GetEnergyChange() const
Definition: G4HadFinalState.cc:52

G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const

G4MaterialCutsCouple::GetIndex
G4int GetIndex() const
Definition: G4MaterialCutsCouple.hh:111

G4HadronicInteraction::GetModelName
const G4String & GetModelName() const
Definition: G4HadronicInteraction.hh:119

G4ParticleDefinition
Definition: G4ParticleDefinition.hh:111

G4HadronicDeprecate
#define G4HadronicDeprecate(name)
Definition: G4HadronicDeprecate.hh:33

G4VProcess::ClearNumberOfInteractionLengthLeft
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:447

G4ParticleDefinition::GetProcessManager
G4ProcessManager * GetProcessManager() const

G4HadronicInteraction
Definition: G4HadronicInteraction.hh:63

G4ParticleDefinition::GetParticleName
const G4String & GetParticleName() const
Definition: G4ParticleDefinition.hh:159

G4Track::SetWeight
void SetWeight(G4double aValue)

G4VParticleChange::ProposeLocalEnergyDeposit
void ProposeLocalEnergyDeposit(G4double anEnergyPart)

G4CrossSectionDataStore.hh

G4VParticleChange::ProposeWeight
void ProposeWeight(G4double finalWeight)

G4HadronicProcess::DumpState
void DumpState(const G4Track &, const G4String &, G4ExceptionDescription &)
Definition: G4HadronicProcess.cc:823

G4ProcessManager::GetAtRestProcessVector
G4ProcessVector * GetAtRestProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const

G4Track::GetKineticEnergy
G4double GetKineticEnergy() const

G4HadronicProcess::AddDataSet
void AddDataSet(G4VCrossSectionDataSet *aDataSet)
Definition: G4HadronicProcess.hh:117

G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:87

G4cout
G4GLOB_DLL std::ostream G4cout

G4HadronicProcess
Definition: G4HadronicProcess.hh:69

G4HadronicProcess::theTotalResult
G4ParticleChange * theTotalResult
Definition: G4HadronicProcess.hh:222

G4HadronElasticDataSet.hh

G4VParticleChange::ProposeNonIonizingEnergyDeposit
void ProposeNonIonizingEnergyDeposit(G4double anEnergyPart)

G4HadronElasticDataSet
Definition: G4HadronElasticDataSet.hh:44

G4DynamicParticle::GetMomentumDirection
const G4ThreeVector & GetMomentumDirection() const

G4Nucleus.hh

G4HadronicProcess::ChooseHadronicInteraction
G4HadronicInteraction * ChooseHadronicInteraction(G4double kineticEnergy, G4Material *aMaterial, G4Element *anElement)
Definition: G4HadronicProcess.hh:142

G4HadronicProcess::GetCrossSectionDataStore
G4CrossSectionDataStore * GetCrossSectionDataStore()
Definition: G4HadronicProcess.hh:170

G4HadronicDeprecate.hh

G4ProcessManager.hh

G4Step
Definition: G4Step.hh:76

G4Neutron::Neutron
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104

G4Track::GetGlobalTime
G4double GetGlobalTime() const

G4VProcess::GetProcessName
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408

G4Track::GetTouchableHandle
const G4TouchableHandle & GetTouchableHandle() const

G4Track::GetMaterial
G4Material * GetMaterial() const

G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

G4ProcessVector::size
G4int size() const

G4ParticleChange::Initialize
virtual void Initialize(const G4Track &)
Definition: G4ParticleChange.cc:228

G4ProductionCutsTable::GetProductionCutsTable
static G4ProductionCutsTable * GetProductionCutsTable()
Definition: G4ProductionCutsTable.cc:63

eV
static const double eV
Definition: G4SIunits.hh:194

G4Track::GetMomentumDirection
const G4ThreeVector & GetMomentumDirection() const

G4ProductionCutsTable.hh

G4HadronicProcess::GetTargetNucleusPointer
G4Nucleus * GetTargetNucleusPointer()
Definition: G4HadronicProcess.hh:149

G4HadFinalState::GetMomentumChange
const G4ThreeVector & GetMomentumChange() const
Definition: G4HadFinalState.cc:68

G4VParticleChange::SetNumberOfSecondaries
void SetNumberOfSecondaries(G4int totSecondaries)

G4HadSecondary::GetParticle
G4DynamicParticle * GetParticle()
Definition: G4HadSecondary.hh:39

G4Nucleus::GetZ_asInt
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115

FatalException
Definition: G4ExceptionSeverity.hh:60

G4ParticleChange::ProposeEnergy
void ProposeEnergy(G4double finalEnergy)

fStopAndKill
Definition: G4TrackStatus.hh:56

G4Track
Definition: G4Track.hh:73

G4ParticleChange::AddSecondary
void AddSecondary(G4Track *aSecondary)
Definition: G4ParticleChange.cc:218

G4Track::GetWeight
G4double GetWeight() const

G4HadronElasticProcess::Description
virtual void Description() const
Definition: G4HadronElasticProcess.cc:58

G4endl
#define G4endl
Definition: G4ios.hh:61

keV
static const double keV
Definition: G4SIunits.hh:195

G4HadronicInteraction::SetRecoilEnergyThreshold
void SetRecoilEnergyThreshold(G4double val)
Definition: G4HadronicInteraction.hh:142

G4HadronicException
Definition: G4HadronicException.hh:33

G4ParticleChange::ProposeMomentumDirection
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)

G4CrossSectionDataStore::SampleZandA
G4Element * SampleZandA(const G4DynamicParticle *, const G4Material *, G4Nucleus &target)
Definition: G4CrossSectionDataStore.cc:234

G4double
double G4double
Definition: G4Types.hh:76

G4VParticleChange::ProposeTrackStatus
void ProposeTrackStatus(G4TrackStatus status)

G4Element::GetName
const G4String & GetName() const
Definition: G4Element.hh:127

G4SystemOfUnits.hh

G4HadronElasticProcess::PostStepDoIt
virtual G4VParticleChange * PostStepDoIt(const G4Track &aTrack, const G4Step &aStep)
Definition: G4HadronElasticProcess.cc:84

G4VParticleChange
Definition: G4VParticleChange.hh:94

G4HadronElasticProcess.hh

G4TrackStatus
G4TrackStatus
Definition: G4TrackStatus.hh:49

G4HadronicException::Report
void Report(std::ostream &aS)
Definition: G4HadronicException.hh:53

G4HadFinalState::GetNumberOfSecondaries
G4int GetNumberOfSecondaries() const
Definition: G4HadFinalState.cc:40

G4HadronicInteraction::ApplyYourself
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0

G4HadFinalState::GetLocalEnergyDeposit
G4double GetLocalEnergyDeposit() const
Definition: G4HadFinalState.cc:131

G4VParticleChange::Clear
void Clear()

G4HadronElasticProcess::PreparePhysicsTable
virtual void PreparePhysicsTable(const G4ParticleDefinition &)
Definition: G4HadronElasticProcess.cc:260

G4HadronElasticProcess::SetLowestEnergyNeutron
virtual void SetLowestEnergyNeutron(G4double)
Definition: G4HadronElasticProcess.cc:275

G4String
Definition: G4String.hh:45

G4HadronElasticProcess::isInitialised
G4bool isInitialised
Definition: G4HadronElasticProcess.hh:78