Geant4: G4MesonAbsorption Class Reference

G4MesonAbsorption::G4MesonAbsorption ( ) [inline]

Definition at line 37 of file G4MesonAbsorption.hh.

00037 {}

virtual G4MesonAbsorption::~G4MesonAbsorption ( ) [inline, virtual]

Definition at line 38 of file G4MesonAbsorption.hh.

00038 {}

G4CollisionInitialState* G4MesonAbsorption::GetCollision	(	G4KineticTrack *	projectile,
		std::vector< G4KineticTrack * >	targets
	)

const std::vector< G4CollisionInitialState * > & G4MesonAbsorption::GetCollisions	(	G4KineticTrack *	aProjectile,
		std::vector< G4KineticTrack * > &	someCandidates,
		G4double	aCurrentTime
	)			`[virtual]`

Implements G4BCAction.

Definition at line 44 of file G4MesonAbsorption.cc.

References DBL_MAX.

00047 {
00048   theCollisions.clear();
00049   if(someCandidates.size() >1)
00050   {
00051     std::vector<G4KineticTrack *>::iterator j=someCandidates.begin();
00052     for(; j != someCandidates.end(); ++j)
00053     {
00054       G4double collisionTime = GetTimeToAbsorption(*aProjectile, **j);
00055       if(collisionTime == DBL_MAX)
00056       {
00057         continue;
00058       }
00059       G4KineticTrackVector aTarget;
00060       aTarget.push_back(*j);
00061       FindAndFillCluster(aTarget, aProjectile, someCandidates);
00062       if(aTarget.size()>=2)
00063       {
00064         theCollisions.push_back(
00065         new G4CollisionInitialState(collisionTime+aCurrentTime, aProjectile, aTarget, this) );
00066       }
00067     }
00068   }
00069   return theCollisions;
00070 }

G4KineticTrackVector * G4MesonAbsorption::GetFinalState	(	G4KineticTrack *	aProjectile,
		std::vector< G4KineticTrack * > &	theTargets
	)			`[virtual]`

Implements G4BCAction.

Definition at line 104 of file G4MesonAbsorption.cc.

References G4endl, G4lrint(), G4UniformRand, G4KineticTrack::Get4Momentum(), G4ParticleDefinition::GetBaryonNumber(), G4KineticTrack::GetDefinition(), G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGMass(), G4Neutron::NeutronDefinition(), G4INCL::Math::pi, G4Proton::ProtonDefinition(), and sqr().

00106 {
00107   // G4cout << "We have an absorption !!!!!!!!!!!!!!!!!!!!!!"<<G4endl;
00108   // Only 2-body absorption for the time being.
00109   // If insufficient, use 2-body absorption and clusterization to emulate 3-,4-,etc body absorption.
00110   G4LorentzVector thePro = projectile->Get4Momentum();
00111   G4LorentzVector theT1 = targets[0]->Get4Momentum();
00112   G4LorentzVector theT2 = targets[1]->Get4Momentum();
00113   G4LorentzVector incoming = thePro + theT1 + theT2;
00114   G4double energyBalance = incoming.t();
00115   G4int chargeBalance = G4lrint(projectile->GetDefinition()->GetPDGCharge()
00116                        + targets[0]->GetDefinition()->GetPDGCharge()
00117                        + targets[1]->GetDefinition()->GetPDGCharge());
00118 
00119   G4int baryonBalance = projectile->GetDefinition()->GetBaryonNumber()
00120                        + targets[0]->GetDefinition()->GetBaryonNumber()
00121                        + targets[1]->GetDefinition()->GetBaryonNumber();
00122 
00123 
00124   // boost all to MMS
00125   G4LorentzRotation toSPS((-1)*(thePro + theT1 + theT2).boostVector());
00126   theT1 = toSPS * theT1;
00127   theT2 = toSPS * theT2;
00128   thePro = toSPS * thePro;
00129   G4LorentzRotation fromSPS(toSPS.inverse());
00130 
00131   // rotate to z
00132   G4LorentzRotation toZ;
00133   G4LorentzVector Ptmp=projectile->Get4Momentum();
00134   toZ.rotateZ(-1*Ptmp.phi());
00135   toZ.rotateY(-1*Ptmp.theta());
00136   theT1 = toZ * theT1;
00137   theT2 = toZ * theT2;
00138   thePro = toZ * thePro;
00139   G4LorentzRotation toLab(toZ.inverse());
00140 
00141   // Get definitions
00142   G4ParticleDefinition * d1 = targets[0]->GetDefinition();
00143   G4ParticleDefinition * d2 = targets[1]->GetDefinition();
00144   if(0.5>G4UniformRand())
00145   {
00146     G4ParticleDefinition * temp;
00147     temp=d1;d1=d2;d2=temp;
00148   }
00149   G4ParticleDefinition * dp = projectile->GetDefinition();
00150   if(dp->GetPDGCharge()<-.5)
00151   {
00152     if(d1->GetPDGCharge()>.5)
00153     {
00154       if(d2->GetPDGCharge()>.5 && 0.5>G4UniformRand())
00155       {
00156         d2 = G4Neutron::NeutronDefinition();
00157       }
00158       else
00159       {
00160         d1 = G4Neutron::NeutronDefinition();
00161       }
00162     }
00163     else if(d2->GetPDGCharge()>.5)
00164     {
00165       d2 = G4Neutron::NeutronDefinition();
00166     }
00167   }
00168   else if(dp->GetPDGCharge()>.5)
00169   {
00170     if(d1->GetPDGCharge()<.5)
00171     {
00172       if(d2->GetPDGCharge()<.5 && 0.5>G4UniformRand())
00173       {
00174         d2 = G4Proton::ProtonDefinition();
00175       }
00176       else
00177       {
00178         d1 = G4Proton::ProtonDefinition();
00179       }
00180     }
00181     else if(d2->GetPDGCharge()<.5)
00182     {
00183       d2 = G4Proton::ProtonDefinition();
00184     }
00185   }
00186 
00187   // calculate the momenta.
00188   G4double M_sq  = (thePro+theT1+theT2).mag2();
00189   G4double m1_sq = sqr(d1->GetPDGMass());
00190   G4double m2_sq = sqr(d2->GetPDGMass());
00191   G4double m_sq  = M_sq-m1_sq-m2_sq;
00192   G4double p = std::sqrt((m_sq*m_sq - 4.*m1_sq * m2_sq)/(4.*M_sq));
00193   G4double costh = 2.*G4UniformRand()-1.;
00194   G4double phi = 2.*pi*G4UniformRand();
00195   G4ThreeVector pFinal(p*std::sin(std::acos(costh))*std::cos(phi), p*std::sin(std::acos(costh))*std::sin(phi), p*costh);
00196 
00197   // G4cout << "testing p "<<p-pFinal.mag()<<G4endl;
00198   // construct the final state particles lorentz momentum.
00199   G4double eFinal1 = std::sqrt(m1_sq+pFinal.mag2());
00200   G4LorentzVector final1(pFinal, eFinal1);
00201   G4double eFinal2 = std::sqrt(m2_sq+pFinal.mag2());
00202   G4LorentzVector final2(-1.*pFinal, eFinal2);
00203 
00204   // rotate back.
00205   final1 = toLab * final1;
00206   final2 = toLab * final2;
00207 
00208   // boost back to LAB
00209   final1 = fromSPS * final1;
00210   final2 = fromSPS * final2;
00211 
00212   // make the final state
00213   G4KineticTrack * f1 =
00214       new G4KineticTrack(d1, 0., targets[0]->GetPosition(), final1);
00215   G4KineticTrack * f2 =
00216       new G4KineticTrack(d2, 0., targets[1]->GetPosition(), final2);
00217   G4KineticTrackVector * result = new G4KineticTrackVector;
00218   result->push_back(f1);
00219   result->push_back(f2);
00220 
00221   for(size_t hpw=0; hpw<result->size(); hpw++)
00222   {
00223     energyBalance-=result->operator[](hpw)->Get4Momentum().t();
00224     chargeBalance-=G4lrint(result->operator[](hpw)->GetDefinition()->GetPDGCharge());
00225     baryonBalance-=result->operator[](hpw)->GetDefinition()->GetBaryonNumber();
00226   }
00227   if(getenv("AbsorptionEnergyBalanceCheck"))
00228     std::cout << "DEBUGGING energy balance B: "
00229               <<energyBalance<<" "
00230               <<chargeBalance<<" "
00231               <<baryonBalance<<" "
00232               <<G4endl;
00233   return result;
00234 }


Public Member Functions
	G4MesonAbsorption ()
virtual	~G4MesonAbsorption ()
virtual const std::vector< G4CollisionInitialState * > &	GetCollisions (G4KineticTrack aProjectile, std::vector< G4KineticTrack > &someCandidates, G4double aCurrentTime)
virtual G4KineticTrackVector *	GetFinalState (G4KineticTrack aProjectile, std::vector< G4KineticTrack > &theTargets)
G4CollisionInitialState *	GetCollision (G4KineticTrack projectile, std::vector< G4KineticTrack > targets)

G4MesonAbsorption Class Reference

Public Member Functions

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation