#include <G4QuasiElasticChannel.hh>
Public Member Functions | |
G4QuasiElasticChannel () | |
~G4QuasiElasticChannel () | |
G4double | GetFraction (G4Nucleus &theNucleus, const G4DynamicParticle &thePrimary) |
G4KineticTrackVector * | Scatter (G4Nucleus &theNucleus, const G4DynamicParticle &thePrimary) |
Definition at line 52 of file G4QuasiElasticChannel.hh.
G4QuasiElasticChannel::G4QuasiElasticChannel | ( | ) |
Definition at line 61 of file G4QuasiElasticChannel.cc.
00062 : theQuasiElastic(G4QuasiElRatios::GetPointer()), 00063 the3DNucleus(new G4Fancy3DNucleus) {}
G4QuasiElasticChannel::~G4QuasiElasticChannel | ( | ) |
G4double G4QuasiElasticChannel::GetFraction | ( | G4Nucleus & | theNucleus, | |
const G4DynamicParticle & | thePrimary | |||
) |
Definition at line 70 of file G4QuasiElasticChannel.cc.
References G4cout, G4endl, G4Nucleus::GetA_asInt(), G4DynamicParticle::GetDefinition(), G4Nucleus::GetN_asInt(), G4ParticleDefinition::GetPDGEncoding(), G4QuasiElRatios::GetRatios(), G4DynamicParticle::GetTotalMomentum(), and G4Nucleus::GetZ_asInt().
Referenced by G4TheoFSGenerator::ApplyYourself().
00072 { 00073 #ifdef debug_scatter 00074 G4cout << "G4QuasiElasticChannel:: P=" << thePrimary.GetTotalMomentum() 00075 << ", pPDG=" << thePrimary.GetDefinition()->GetPDGEncoding() 00076 << ", Z = " << theNucleus.GetZ_asInt()) 00077 << ", N = " << theNucleus.GetN_asInt()) 00078 << ", A = " << theNucleus.GetA_asInt() << G4endl; 00079 #endif 00080 00081 std::pair<G4double,G4double> ratios; 00082 ratios=theQuasiElastic->GetRatios(thePrimary.GetTotalMomentum(), 00083 thePrimary.GetDefinition()->GetPDGEncoding(), 00084 theNucleus.GetZ_asInt(), 00085 theNucleus.GetN_asInt()); 00086 #ifdef debug_scatter 00087 G4cout << "G4QuasiElasticChannel::ratios " << ratios.first << " x " <<ratios.second 00088 << " = " << ratios.first*ratios.second << G4endl; 00089 #endif 00090 00091 return ratios.first*ratios.second; 00092 }
G4KineticTrackVector * G4QuasiElasticChannel::Scatter | ( | G4Nucleus & | theNucleus, | |
const G4DynamicParticle & | thePrimary | |||
) |
Definition at line 94 of file G4QuasiElasticChannel.cc.
References G4ParticleTable::FindIon(), G4cout, G4endl, G4lrint(), G4UniformRand, G4DynamicParticle::Get4Momentum(), G4Nucleus::GetA_asInt(), G4DynamicParticle::GetDefinition(), G4V3DNucleus::GetMass(), G4V3DNucleus::GetNucleons(), G4ParticleDefinition::GetParticleName(), G4ParticleTable::GetParticleTable(), G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGEncoding(), G4ParticleDefinition::GetPDGMass(), G4Nucleus::GetZ_asInt(), G4V3DNucleus::Init(), G4Neutron::Neutron(), G4QuasiElRatios::Scatter(), and sqr().
Referenced by G4TheoFSGenerator::ApplyYourself().
00096 { 00097 G4int A=theNucleus.GetA_asInt(); 00098 G4int Z=theNucleus.GetZ_asInt(); 00099 // build Nucleus and choose random nucleon to scatter with 00100 the3DNucleus->Init(theNucleus.GetA_asInt(),theNucleus.GetZ_asInt()); 00101 const std::vector<G4Nucleon>& nucleons=the3DNucleus->GetNucleons(); 00102 G4double targetNucleusMass=the3DNucleus->GetMass(); 00103 G4LorentzVector targetNucleus4Mom(0.,0.,0.,targetNucleusMass); 00104 G4int index; 00105 do { 00106 index=G4lrint((A-1)*G4UniformRand()); 00107 } while (index < 0 || index >= static_cast<G4int>(nucleons.size())); 00108 00109 G4ParticleDefinition * pDef= nucleons[index].GetDefinition(); 00110 00111 G4int resA=A - 1; 00112 G4int resZ=Z - static_cast<int>(pDef->GetPDGCharge()); 00113 G4ParticleDefinition* resDef; 00114 G4double residualNucleusMass; 00115 if(resZ) 00116 { 00117 resDef=G4ParticleTable::GetParticleTable()->FindIon(resZ,resA,0,resZ); 00118 residualNucleusMass=resDef->GetPDGMass(); 00119 } 00120 else { 00121 resDef=G4Neutron::Neutron(); 00122 residualNucleusMass=resA * G4Neutron::Neutron()->GetPDGMass(); 00123 } 00124 #ifdef debug_scatter 00125 G4cout<<"G4QElChan::Scatter: neutron - proton? A ="<<A<<", Z="<<Z<<", projName=" 00126 <<pDef->GetParticleName()<<G4endl; 00127 #endif 00128 00129 G4LorentzVector pNucleon=nucleons[index].Get4Momentum(); 00130 G4double residualNucleusEnergy=std::sqrt(sqr(residualNucleusMass) + 00131 pNucleon.vect().mag2()); 00132 pNucleon.setE(targetNucleusMass-residualNucleusEnergy); 00133 G4LorentzVector residualNucleus4Mom=targetNucleus4Mom-pNucleon; 00134 00135 std::pair<G4LorentzVector,G4LorentzVector> result; 00136 00137 result=theQuasiElastic->Scatter(pDef->GetPDGEncoding(),pNucleon, 00138 thePrimary.GetDefinition()->GetPDGEncoding(), 00139 thePrimary.Get4Momentum()); 00140 G4LorentzVector scatteredHadron4Mom; 00141 if (result.first.e() > 0.) 00142 scatteredHadron4Mom=result.second; 00143 else { //scatter failed 00144 //G4cout << "Warning - G4QuasiElasticChannel::Scatter no scattering" << G4endl; 00145 //return 0; //no scatter 00146 scatteredHadron4Mom=thePrimary.Get4Momentum(); 00147 residualNucleus4Mom=G4LorentzVector(0.,0.,0.,targetNucleusMass); 00148 resDef=G4ParticleTable::GetParticleTable()->FindIon(Z,A,0,Z); 00149 } 00150 00151 #ifdef debug_scatter 00152 G4LorentzVector EpConservation=pNucleon+thePrimary.Get4Momentum() 00153 - result.first - result.second; 00154 if ( (EpConservation.vect().mag2() > .01*MeV*MeV ) 00155 || (std::abs(EpConservation.e()) > 0.1 * MeV ) ) 00156 { 00157 G4cout << "Warning - G4QuasiElasticChannel::Scatter E-p non conservation : " 00158 << EpConservation << G4endl; 00159 } 00160 #endif 00161 00162 G4KineticTrackVector * ktv = new G4KineticTrackVector(); 00163 G4KineticTrack * sPrim=new G4KineticTrack(thePrimary.GetDefinition(), 00164 0.,G4ThreeVector(0), scatteredHadron4Mom); 00165 ktv->push_back(sPrim); 00166 if (result.first.e() > 0.) 00167 { 00168 G4KineticTrack * sNuc=new G4KineticTrack(pDef, 0.,G4ThreeVector(0), result.first); 00169 ktv->push_back(sNuc); 00170 } 00171 if(resZ || resA==1) // For the only neutron or for tnuclei with Z>0 00172 { 00173 G4KineticTrack * rNuc=new G4KineticTrack(resDef, 00174 0.,G4ThreeVector(0), residualNucleus4Mom); 00175 ktv->push_back(rNuc); 00176 } 00177 else // The residual nucleus consists of only neutrons 00178 { 00179 residualNucleus4Mom/=resA; // Split 4-mom of A*n system equally 00180 for(G4int in=0; in<resA; in++) // Loop over neutrons in A*n system. 00181 { 00182 G4KineticTrack* rNuc=new G4KineticTrack(resDef, 00183 0.,G4ThreeVector(0), residualNucleus4Mom); 00184 ktv->push_back(rNuc); 00185 } 00186 } 00187 #ifdef debug_scatter 00188 G4cout<<"G4QElC::Scat: Nucleon: "<<result.first <<" mass "<<result.first.mag() << G4endl; 00189 G4cout<<"G4QElC::Scat: Project: "<<result.second<<" mass "<<result.second.mag()<< G4endl; 00190 #endif 00191 return ktv; 00192 }