#include <G4ElectroNuclearReaction.hh>
Inheritance diagram for G4ElectroNuclearReaction:
Public Member Functions | |
G4ElectroNuclearReaction () | |
~G4ElectroNuclearReaction () | |
G4HadFinalState * | ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &aTargetNucleus) |
void | Description () const |
Definition at line 61 of file G4ElectroNuclearReaction.hh.
G4ElectroNuclearReaction::G4ElectroNuclearReaction | ( | ) | [inline] |
Definition at line 64 of file G4ElectroNuclearReaction.hh.
References Description(), G4VPartonStringModel::SetFragmentationModel(), G4TheoFSGenerator::SetHighEnergyGenerator(), G4HadronicInteraction::SetMaxEnergy(), G4HadronicInteraction::SetMinEnergy(), and G4TheoFSGenerator::SetTransport().
00064 :G4HadronicInteraction("CHIPSElectroNuclear") 00065 { 00066 SetMinEnergy(0*CLHEP::GeV); 00067 SetMaxEnergy(100*CLHEP::TeV); 00068 00069 theHEModel = new G4TheoFSGenerator; 00070 theCascade = new G4GeneratorPrecompoundInterface; 00071 theHEModel->SetTransport(theCascade); 00072 theHEModel->SetHighEnergyGenerator(&theStringModel); 00073 theStringDecay = new G4ExcitedStringDecay(&theFragmentation); 00074 theStringModel.SetFragmentationModel(theStringDecay); 00075 theHEModel->SetMinEnergy(2.5*CLHEP::GeV); 00076 theHEModel->SetMaxEnergy(100*CLHEP::TeV); 00077 theElectronData = new G4ElectroNuclearCrossSection; 00078 thePhotonData = new G4PhotoNuclearCrossSection; 00079 Description(); 00080 }
G4ElectroNuclearReaction::~G4ElectroNuclearReaction | ( | ) | [inline] |
Definition at line 82 of file G4ElectroNuclearReaction.hh.
00083 { 00084 delete theHEModel; 00085 delete theCascade; 00086 delete theStringDecay; 00087 delete theElectronData; 00088 delete thePhotonData; 00089 }
G4HadFinalState * G4ElectroNuclearReaction::ApplyYourself | ( | const G4HadProjectile & | aTrack, | |
G4Nucleus & | aTargetNucleus | |||
) | [inline, virtual] |
Implements G4HadronicInteraction.
Definition at line 142 of file G4ElectroNuclearReaction.hh.
References G4HadFinalState::AddSecondaries(), G4TheoFSGenerator::ApplyYourself(), G4ChiralInvariantPhaseSpace::ApplyYourself(), G4HadFinalState::Clear(), G4Electron::Electron(), G4Electron::ElectronDefinition(), G4cerr, G4cout, G4endl, G4UniformRand, G4Gamma::GammaDefinition(), G4HadProjectile::Get4Momentum(), G4VCrossSectionDataSet::GetCrossSection(), G4HadProjectile::GetDefinition(), G4Element::GetElementTable(), G4ElectroNuclearCrossSection::GetEquivalentPhotonEnergy(), G4ElectroNuclearCrossSection::GetEquivalentPhotonQ2(), G4ParticleDefinition::GetPDGMass(), G4ElectroNuclearCrossSection::GetVirtualFactor(), G4Nucleus::GetZ_asInt(), isAlive, G4Neutron::Neutron(), position, G4Positron::PositronDefinition(), G4Proton::Proton(), G4HadFinalState::SetEnergyChange(), G4DynamicParticle::SetKineticEnergy(), G4HadFinalState::SetMomentumChange(), and G4HadFinalState::SetStatusChange().
00144 { 00145 theResult.Clear(); 00146 static const G4double dM=G4Proton::Proton()->GetPDGMass() + 00147 G4Neutron::Neutron()->GetPDGMass(); // MeanDoubleNucleon Mass = m_n+m_p (@@ no binding) 00148 static const G4double me=G4Electron::Electron()->GetPDGMass(); // electron mass 00149 static const G4double me2=me*me; // squared electron mass 00150 G4DynamicParticle theTempEl(const_cast<G4ParticleDefinition *>(aTrack.GetDefinition()), 00151 aTrack.Get4Momentum().vect()); 00152 const G4DynamicParticle* theElectron=&theTempEl; 00153 const G4ParticleDefinition* aD = theElectron->GetDefinition(); 00154 if((aD != G4Electron::ElectronDefinition()) && (aD != G4Positron::PositronDefinition())) 00155 throw G4HadronicException(__FILE__, __LINE__, 00156 "G4ElectroNuclearReaction::ApplyYourself called for neither electron or positron"); 00157 const G4ElementTable* aTab = G4Element::GetElementTable(); 00158 G4Element * anElement = 0; 00159 G4int aZ = static_cast<G4int>(aTargetNucleus.GetZ_asInt()+.1); 00160 for(size_t ii=0; ii<aTab->size(); ++ii) if ( std::abs((*aTab)[ii]->GetZ()-aZ) < .1) 00161 { 00162 anElement = (*aTab)[ii]; 00163 break; 00164 } 00165 if(0==anElement) 00166 { 00167 G4cerr<<"***G4ElectroNuclearReaction::ApplyYourself: element with Z=" 00168 <<aTargetNucleus.GetZ_asInt()<<" is not in the element table"<<G4endl; 00169 throw G4HadronicException(__FILE__, __LINE__, "Anomalous element error."); 00170 } 00171 00172 // Note: high energy gamma nuclear now implemented. 00173 G4double xSec = theElectronData->GetCrossSection(theElectron, anElement);// Check out XS 00174 if(xSec<=0.) 00175 { 00176 theResult.SetStatusChange(isAlive); 00177 theResult.SetEnergyChange(theElectron->GetKineticEnergy()); 00178 // new direction for the electron 00179 theResult.SetMomentumChange(theElectron->GetMomentumDirection()); 00180 return &theResult; // DO-NOTHING condition 00181 } 00182 G4double photonEnergy = theElectronData->GetEquivalentPhotonEnergy(); 00183 G4double theElectronKinEnergy=theElectron->GetKineticEnergy(); 00184 if( theElectronKinEnergy < photonEnergy ) 00185 { 00186 G4cout<<"G4ElectroNuclearReaction::ApplyYourself: photonEnergy is very high"<<G4endl; 00187 G4cout<<">>> If this condition persists, please contact Geant4 group"<<G4endl; 00188 theResult.SetStatusChange(isAlive); 00189 theResult.SetEnergyChange(theElectron->GetKineticEnergy()); 00190 // new direction for the electron 00191 theResult.SetMomentumChange(theElectron->GetMomentumDirection()); 00192 return &theResult; // DO-NOTHING condition 00193 } 00194 G4double photonQ2 = theElectronData->GetEquivalentPhotonQ2(photonEnergy); 00195 G4double W=photonEnergy-photonQ2/dM; // Hadronic energy flow from the virtual photon 00196 if(getenv("debug_G4ElectroNuclearReaction") ) 00197 { 00198 G4cout << "G4ElectroNuclearReaction: Equivalent Energy = "<<W<<G4endl; 00199 } 00200 if(W<0.) 00201 { 00202 theResult.SetStatusChange(isAlive); 00203 theResult.SetEnergyChange(theElectron->GetKineticEnergy()); 00204 // new direction for the electron 00205 theResult.SetMomentumChange(theElectron->GetMomentumDirection()); 00206 return &theResult; // DO-NOTHING condition 00207 // throw G4HadronicException(__FILE__, __LINE__, 00208 // "G4ElectroNuclearReaction::ApplyYourself: negative equivalent energy"); 00209 } 00210 G4DynamicParticle* theDynamicPhoton = new 00211 G4DynamicParticle(G4Gamma::GammaDefinition(), 00212 G4ParticleMomentum(1.,0.,0.), photonEnergy*CLHEP::MeV); //----->-* 00213 G4double sigNu=thePhotonData->GetCrossSection(theDynamicPhoton, anElement); // | 00214 theDynamicPhoton->SetKineticEnergy(W); // Redefine photon with equivalent energy | 00215 G4double sigK =thePhotonData->GetCrossSection(theDynamicPhoton, anElement); // | 00216 delete theDynamicPhoton; // <--------------------------------------------------------* 00217 G4double rndFraction = theElectronData->GetVirtualFactor(photonEnergy, photonQ2); 00218 if(sigNu*G4UniformRand()>sigK*rndFraction) 00219 { 00220 theResult.SetStatusChange(isAlive); 00221 theResult.SetEnergyChange(theElectron->GetKineticEnergy()); 00222 // new direction for the electron 00223 theResult.SetMomentumChange(theElectron->GetMomentumDirection()); 00224 return &theResult; // DO-NOTHING condition 00225 } 00226 theResult.SetStatusChange(isAlive); 00227 // Scatter an electron and make gamma+A reaction 00228 G4double iniE=theElectronKinEnergy+me; // Initial total energy of electron 00229 G4double finE=iniE-photonEnergy; // Final total energy of electron 00230 theResult.SetEnergyChange(std::max(0.,finE-me)); // Modifies the KINETIC ENERGY 00231 G4double EEm=iniE*finE-me2; // Just an intermediate value to avoid "2*" 00232 G4double iniP=std::sqrt(iniE*iniE-me2); // Initial momentum of the electron 00233 G4double finP=std::sqrt(finE*finE-me2); // Final momentum of the electron 00234 G4double cost=(EEm+EEm-photonQ2)/iniP/finP;// std::cos(theta) for the electron scattering 00235 if(cost> 1.) cost= 1.; 00236 if(cost<-1.) cost=-1.; 00237 G4ThreeVector dir=theElectron->GetMomentumDirection(); // Direction of primary electron 00238 G4ThreeVector ort=dir.orthogonal(); // Not normed orthogonal vector (!) (to dir) 00239 G4ThreeVector ortx = ort.unit(); // First unit vector orthogonal to the direction 00240 G4ThreeVector orty = dir.cross(ortx); // Second unit vector orthoganal to the direction 00241 G4double sint=std::sqrt(1.-cost*cost); // Perpendicular component 00242 G4double phi=CLHEP::twopi*G4UniformRand(); // phi of scattered electron 00243 G4double sinx=sint*std::sin(phi); // x-component 00244 G4double siny=sint*std::cos(phi); // y-component 00245 G4ThreeVector findir=cost*dir+sinx*ortx+siny*orty; 00246 theResult.SetMomentumChange(findir); // new direction for the electron 00247 G4ThreeVector photonMomentum=iniP*dir-finP*findir; 00248 G4DynamicParticle localGamma(G4Gamma::GammaDefinition(), photonEnergy, photonMomentum); 00249 //G4DynamicParticle localGamma(G4Gamma::GammaDefinition(),photonDirection, photonEnergy); 00250 //G4DynamicParticle localGamma(G4Gamma::GammaDefinition(), photonLorentzVector); 00251 G4ThreeVector position(0,0,0); 00252 G4HadProjectile localTrack(localGamma); 00253 G4HadFinalState * result; 00254 if (photonEnergy < 3*CLHEP::GeV) 00255 result = theLEModel.ApplyYourself(localTrack, aTargetNucleus, &theResult); 00256 else { 00257 // G4cout << "0) Getting a high energy electro-nuclear reaction"<<G4endl; 00258 G4HadFinalState * aResult = theHEModel->ApplyYourself(localTrack, aTargetNucleus); 00259 theResult.AddSecondaries(aResult); 00260 aResult->Clear(); 00261 result = &theResult; 00262 } 00263 return result; 00264 }
void G4ElectroNuclearReaction::Description | ( | ) | const [inline] |
Definition at line 94 of file G4ElectroNuclearReaction.hh.
References G4HadronicInteraction::GetModelName().
Referenced by G4ElectroNuclearReaction().
00095 { 00096 char* dirName = getenv("G4PhysListDocDir"); 00097 if (dirName) { 00098 std::ofstream outFile; 00099 G4String outFileName = GetModelName() + ".html"; 00100 G4String pathName = G4String(dirName) + "/" + outFileName; 00101 00102 outFile.open(pathName); 00103 outFile << "<html>\n"; 00104 outFile << "<head>\n"; 00105 00106 outFile << "<title>Description of CHIPS ElectroNuclear Model</title>\n"; 00107 outFile << "</head>\n"; 00108 outFile << "<body>\n"; 00109 00110 outFile << "G4ElectroNuclearReaction handles the inelastic scattering\n" 00111 << "of e- and e+ from nuclei using the Chiral Invariant Phase\n" 00112 << "Space (CHIPS) model of M. Kossov. This model uses the\n" 00113 << "Equivalent Photon Approximation in which the incoming\n" 00114 << "electron generates a virtual photon at the electromagnetic\n" 00115 << "vertex, and the virtual photon is converted to a real photon\n" 00116 << "before it interacts with the nucleus. The real photon\n" 00117 << "interacts with the hadrons in the target by producing\n" 00118 << "quasmons (or generalized excited hadrons) which then decay\n" 00119 << "into final state hadrons. This model is valid for e- and\n" 00120 << "e+ of all incident energies.\n"; 00121 00122 outFile << "</body>\n"; 00123 outFile << "</html>\n"; 00124 outFile.close(); 00125 } 00126 }