G4GeneratorPrecompoundInterface Class Reference

#include <G4GeneratorPrecompoundInterface.hh>

Inheritance diagram for G4GeneratorPrecompoundInterface:

G4VIntraNuclearTransportModel G4HadronicInteraction

Public Member Functions

 G4GeneratorPrecompoundInterface (G4VPreCompoundModel *p=0)
virtual ~G4GeneratorPrecompoundInterface ()
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
virtual G4ReactionProductVectorPropagate (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus)
void SetCaptureThreshold (G4double)
virtual void PropagateModelDescription (std::ostream &) const

Detailed Description

Definition at line 58 of file G4GeneratorPrecompoundInterface.hh.


Constructor & Destructor Documentation

G4GeneratorPrecompoundInterface::G4GeneratorPrecompoundInterface ( G4VPreCompoundModel p = 0  ) 

Definition at line 59 of file G4GeneratorPrecompoundInterface.cc.

References G4HadronicInteractionRegistry::FindModel(), G4HadronicInteractionRegistry::Instance(), G4Neutron::Neutron(), G4Proton::Proton(), and G4VIntraNuclearTransportModel::SetDeExcitation().

00060 : CaptureThreshold(10*MeV)
00061 {
00062     proton = G4Proton::Proton();
00063     neutron = G4Neutron::Neutron();
00064     if(preModel) { SetDeExcitation(preModel); }
00065     else  { 
00066       G4HadronicInteraction* hadi =  
00067         G4HadronicInteractionRegistry::Instance()->FindModel("PRECO");
00068       G4VPreCompoundModel* pre = static_cast<G4VPreCompoundModel*>(hadi);
00069       if(!pre) { pre = new G4PreCompoundModel(); }
00070       SetDeExcitation(pre); 
00071     }
00072 }

G4GeneratorPrecompoundInterface::~G4GeneratorPrecompoundInterface (  )  [virtual]

Definition at line 74 of file G4GeneratorPrecompoundInterface.cc.

00075 {
00076 }


Member Function Documentation

G4HadFinalState * G4GeneratorPrecompoundInterface::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
) [virtual]

Implements G4HadronicInteraction.

Definition at line 235 of file G4GeneratorPrecompoundInterface.cc.

References G4cout, and G4endl.

00236 {
00237     G4cout << "G4GeneratorPrecompoundInterface: ApplyYourself interface called stand-allone."
00238             << G4endl;
00239     G4cout << "This class is only a mediator between generator and precompound"<<G4endl;
00240     G4cout << "Please remove from your physics list."<<G4endl;
00241     throw G4HadronicException(__FILE__, __LINE__, "SEVERE: G4GeneratorPrecompoundInterface model interface called stand-allone.");
00242     return new G4HadFinalState;
00243 }

G4ReactionProductVector * G4GeneratorPrecompoundInterface::Propagate ( G4KineticTrackVector theSecondaries,
G4V3DNucleus theNucleus 
) [virtual]

Implements G4VIntraNuclearTransportModel.

Definition at line 85 of file G4GeneratorPrecompoundInterface.cc.

References G4Alpha::AlphaDefinition(), G4Nucleon::AreYouHit(), G4VPreCompoundModel::DeExcite(), G4Deuteron::DeuteronDefinition(), G4cout, G4endl, G4HadProjectile::Get4Momentum(), G4Nucleon::Get4Momentum(), G4Nucleon::GetBindingEnergy(), G4V3DNucleus::GetCharge(), G4Nucleon::GetDefinition(), G4IonTable::GetIon(), G4ParticleTable::GetIonTable(), G4V3DNucleus::GetMassNumber(), G4V3DNucleus::GetNextNucleon(), G4NucleiProperties::GetNuclearMass(), G4V3DNucleus::GetNuclearRadius(), G4ParticleTable::GetParticleTable(), G4ParticleDefinition::GetPDGCharge(), G4VIntraNuclearTransportModel::GetPrimaryProjectile(), G4He3::He3Definition(), G4Neutron::NeutronDefinition(), G4Proton::ProtonDefinition(), G4ReactionProduct::SetMomentum(), G4Fragment::SetNumberOfCharged(), G4Fragment::SetNumberOfHoles(), G4Fragment::SetNumberOfParticles(), G4ReactionProduct::SetTotalEnergy(), G4V3DNucleus::StartLoop(), G4VIntraNuclearTransportModel::theDeExcitation, and G4Triton::TritonDefinition().

00086 {
00087     G4ReactionProductVector * theTotalResult = new G4ReactionProductVector;
00088 
00089     // decay the strong resonances
00090     G4DecayKineticTracks decay(theSecondaries);
00091 
00092     // prepare the fragment
00093     G4int anA=theNucleus->GetMassNumber();
00094     G4int aZ=theNucleus->GetCharge();
00095     G4int numberOfEx = 0;
00096     G4int numberOfCh = 0;
00097     G4int numberOfHoles = 0;
00098     G4double exEnergy = 0.0;
00099     G4double R = theNucleus->GetNuclearRadius();
00100     G4ThreeVector exciton3Momentum(0.,0.,0.);
00101     G4ThreeVector captured3Momentum(0.,0.,0.);
00102     G4ThreeVector wounded3Momentum(0.,0.,0.);
00103 
00104     // loop over secondaries
00105 #ifdef exactExcitationEnergy
00106     G4LorentzVector secondary4Momemtum(0,0,0,0);
00107 #endif
00108     G4KineticTrackVector::iterator iter;
00109     for(iter=theSecondaries->begin(); iter !=theSecondaries->end(); ++iter)
00110     {
00111         G4ParticleDefinition* part = (*iter)->GetDefinition();
00112         G4double e = (*iter)->Get4Momentum().e();
00113         G4double mass = (*iter)->Get4Momentum().mag();
00114         G4ThreeVector mom = (*iter)->Get4Momentum().vect();
00115         if((part != proton && part != neutron) ||
00116                 (e > mass + CaptureThreshold) ||
00117                 ((*iter)->GetPosition().mag() > R)) {
00118             G4ReactionProduct * theNew = new G4ReactionProduct(part);
00119             theNew->SetMomentum(mom);
00120             theNew->SetTotalEnergy(e);
00121             theTotalResult->push_back(theNew);
00122 #ifdef exactExcitationEnergy
00123             secondary4Momemtum += (*iter)->Get4Momentum();
00124 #endif
00125         } else {
00126             // within the nucleus, neutron or proton
00127             // now calculate  A, Z of the fragment, momentum, number of exciton states
00128             ++anA;
00129             ++numberOfEx;
00130             G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
00131             aZ += Z;
00132             numberOfCh += Z;
00133             captured3Momentum += mom;
00134             exEnergy += (e - mass);
00135         }
00136         delete (*iter);
00137     }
00138     delete theSecondaries;
00139 
00140     // loop over wounded nucleus
00141     G4Nucleon * theCurrentNucleon =
00142             theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
00143     while(theCurrentNucleon) {
00144         if(theCurrentNucleon->AreYouHit()) {
00145             ++numberOfHoles;
00146             ++numberOfEx;
00147             --anA;
00148             aZ -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/eplus + 0.1);
00149             wounded3Momentum += theCurrentNucleon->Get4Momentum().vect();
00150             //G4cout << "hit nucleon " << theCurrentNucleon->Get4Momentum() << G4endl;
00151             exEnergy += theCurrentNucleon->GetBindingEnergy();
00152         }
00153         theCurrentNucleon = theNucleus->GetNextNucleon();
00154     }
00155     exciton3Momentum = captured3Momentum - wounded3Momentum;
00156 
00157     if(anA>0 && aZ>0) {
00158         G4double fMass =  G4NucleiProperties::GetNuclearMass(anA, aZ);
00159 #ifdef exactExcitationEnergy
00160         // recalculate exEnergy from Energy balance....
00161         const G4HadProjectile * primary = GetPrimaryProjectile();
00162         G4double Einitial= primary->Get4Momentum().e()
00163                                 + G4NucleiProperties::GetNuclearMass(theNucleus->GetMassNumber(),theNucleus->GetCharge());
00164         G4double Efinal = fMass + secondary4Momemtum.e();
00165         if ( (Einitial - Efinal) > 0 ) {
00166             // G4cout << "G4GPI::Propagate() : positive exact excitation Energy "
00167             //  << (Einitial - Efinal)/MeV << " MeV, exciton estimate " << exEnergy/MeV << " MeV" << G4endl;
00168             exEnergy=Einitial - Efinal;
00169         }
00170         else {
00171             //  G4cout << "G4GeneratorPrecompoundInterface::Propagate() : negative exact excitation Energy "
00172             //  << (Einitial - Efinal)/MeV << " MeV, setting  excitation to 0 MeV" << G4endl;
00173             exEnergy=0.;
00174         }
00175 #endif
00176         fMass += exEnergy;
00177         G4ThreeVector balance=primary->Get4Momentum().vect() - secondary4Momemtum.vect() - exciton3Momentum;
00178         #ifdef G4GPI_debug_excitation
00179         G4cout << "momentum balance init/final  " << balance << " value " << balance.mag() << G4endl
00180                 << "primary / secondaries "<< primary->Get4Momentum() << " / "
00181                 << secondary4Momemtum << " captured/wounded: " << captured3Momentum << " / " <<  wounded3Momentum
00182                 << "  exciton " << exciton3Momentum << G4endl
00183                 << secondary4Momemtum.vect() + exciton3Momentum << G4endl;
00184         #endif
00185 #ifdef exactExcitationEnergy
00186         G4LorentzVector exciton4Momentum(exciton3Momentum, fMass);
00187 #else
00188         G4LorentzVector exciton4Momentum(exciton3Momentum,
00189                 std::sqrt(exciton3Momentum.mag2() + fMass*fMass));
00190 #endif
00191         if ( exEnergy > 0.0 ) {  // Need to de-excite the remnant nucleus only if excitation energy > 0.
00192             G4Fragment anInitialState(anA, aZ, exciton4Momentum);
00193             anInitialState.SetNumberOfParticles(numberOfEx-numberOfHoles);
00194             anInitialState.SetNumberOfCharged(numberOfCh);
00195             anInitialState.SetNumberOfHoles(numberOfHoles);
00196 
00197             G4ReactionProductVector * aPrecoResult = theDeExcitation->DeExcite(anInitialState);
00198             // fill pre-compound part into the result, and return
00199             theTotalResult->insert(theTotalResult->end(),aPrecoResult->begin(),aPrecoResult->end() );
00200             delete aPrecoResult;
00201 
00202         } else {  // No/negative excitation energy, we only need to create the remnant nucleus
00203             //  energy is not conserved, ignore exciton momentum, i.e. remnant nucleus will be at rest
00204             G4ParticleDefinition* theKindOfFragment = 0;
00205             if (anA == 1 && aZ == 0) {
00206                 theKindOfFragment = G4Neutron::NeutronDefinition();
00207             } else if (anA == 1 && aZ == 1) {
00208                 theKindOfFragment = G4Proton::ProtonDefinition();
00209             } else if (anA == 2 && aZ == 1) {
00210                 theKindOfFragment = G4Deuteron::DeuteronDefinition();
00211             } else if (anA == 3 && aZ == 1) {
00212                 theKindOfFragment = G4Triton::TritonDefinition();
00213             } else if (anA == 3 && aZ == 2) {
00214                 theKindOfFragment = G4He3::He3Definition();
00215             } else if (anA == 4 && aZ == 2) {
00216                 theKindOfFragment = G4Alpha::AlphaDefinition();;
00217             } else {
00218                 theKindOfFragment =
00219                         G4ParticleTable::GetParticleTable()->GetIonTable()->GetIon(aZ,anA,0.0);
00220             }
00221             if (theKindOfFragment != 0) {
00222                 G4ReactionProduct * theNew = new G4ReactionProduct(theKindOfFragment);
00223                 theNew->SetMomentum(G4ThreeVector(0.,0.,0.));
00224                 theNew->SetTotalEnergy(fMass);
00225                 //theNew->SetFormationTime(??0.??);
00226                 theTotalResult->push_back(theNew);
00227             }
00228         }
00229     }
00230 
00231     return theTotalResult;
00232 }

void G4GeneratorPrecompoundInterface::PropagateModelDescription ( std::ostream &   )  const [virtual]

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 244 of file G4GeneratorPrecompoundInterface.cc.

00245 {
00246     outFile << "G4GeneratorPrecompoundInterface interfaces a high\n"
00247             << "energy model through the wounded nucleus to precompound de-excition.\n"
00248             << "Low energy protons and neutron present among secondaries produced by \n"
00249             << "the high energy generator and within the nucleus are captured. The wounded\n"
00250             << "nucleus and the captured particles form an excited nuclear fragment. This\n"
00251             << "fragment is passed to the Geant4 pre-compound model for de-excitation.\n"
00252             << "Nuclear de-excitation:\n";
00253     // preco
00254 
00255 }

void G4GeneratorPrecompoundInterface::SetCaptureThreshold ( G4double   )  [inline]

Definition at line 88 of file G4GeneratorPrecompoundInterface.hh.

00089 {
00090   CaptureThreshold=value;
00091 }


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:04 2013 for Geant4 by  doxygen 1.4.7