Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4MesonAbsorption Class Reference

#include <G4MesonAbsorption.hh>

Inheritance diagram for G4MesonAbsorption:
Collaboration diagram for G4MesonAbsorption:

Public Member Functions

 G4MesonAbsorption ()
 
virtual ~G4MesonAbsorption ()
 
virtual const std::vector
< G4CollisionInitialState * > & 
GetCollisions (G4KineticTrack *aProjectile, std::vector< G4KineticTrack * > &someCandidates, G4double aCurrentTime)
 
virtual G4KineticTrackVectorGetFinalState (G4KineticTrack *aProjectile, std::vector< G4KineticTrack * > &theTargets)
 
G4CollisionInitialStateGetCollision (G4KineticTrack *projectile, std::vector< G4KineticTrack * > targets)
 
- Public Member Functions inherited from G4BCAction
 G4BCAction ()
 
virtual ~G4BCAction ()
 

Detailed Description

Definition at line 34 of file G4MesonAbsorption.hh.

Constructor & Destructor Documentation

G4MesonAbsorption::G4MesonAbsorption ( )
inline

Definition at line 37 of file G4MesonAbsorption.hh.

37 {}
virtual G4MesonAbsorption::~G4MesonAbsorption ( )
inlinevirtual

Definition at line 38 of file G4MesonAbsorption.hh.

38 {}

Member Function Documentation

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.

47 {
48  theCollisions.clear();
49  if(someCandidates.size() >1)
50  {
51  std::vector<G4KineticTrack *>::iterator j=someCandidates.begin();
52  for(; j != someCandidates.end(); ++j)
53  {
54  G4double collisionTime = GetTimeToAbsorption(*aProjectile, **j);
55  if(collisionTime == DBL_MAX)
56  {
57  continue;
58  }
59  G4KineticTrackVector aTarget;
60  aTarget.push_back(*j);
61  FindAndFillCluster(aTarget, aProjectile, someCandidates);
62  if(aTarget.size()>=2)
63  {
64  theCollisions.push_back(
65  new G4CollisionInitialState(collisionTime+aCurrentTime, aProjectile, aTarget, this) );
66  }
67  }
68  }
69  return theCollisions;
70 }
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
G4KineticTrackVector * G4MesonAbsorption::GetFinalState ( G4KineticTrack aProjectile,
std::vector< G4KineticTrack * > &  theTargets 
)
virtual

Implements G4BCAction.

Definition at line 104 of file G4MesonAbsorption.cc.

106 {
107  // G4cout << "We have an absorption !!!!!!!!!!!!!!!!!!!!!!"<<G4endl;
108  // Only 2-body absorption for the time being.
109  // If insufficient, use 2-body absorption and clusterization to emulate 3-,4-,etc body absorption.
110  G4LorentzVector thePro = projectile->Get4Momentum();
111  G4LorentzVector theT1 = targets[0]->Get4Momentum();
112  G4LorentzVector theT2 = targets[1]->Get4Momentum();
113  G4LorentzVector incoming = thePro + theT1 + theT2;
114  G4double energyBalance = incoming.t();
115  G4int chargeBalance = G4lrint(projectile->GetDefinition()->GetPDGCharge()
116  + targets[0]->GetDefinition()->GetPDGCharge()
117  + targets[1]->GetDefinition()->GetPDGCharge());
118 
119  G4int baryonBalance = projectile->GetDefinition()->GetBaryonNumber()
120  + targets[0]->GetDefinition()->GetBaryonNumber()
121  + targets[1]->GetDefinition()->GetBaryonNumber();
122 
123 
124  // boost all to MMS
125  G4LorentzRotation toSPS((-1)*(thePro + theT1 + theT2).boostVector());
126  theT1 = toSPS * theT1;
127  theT2 = toSPS * theT2;
128  thePro = toSPS * thePro;
129  G4LorentzRotation fromSPS(toSPS.inverse());
130 
131  // rotate to z
132  G4LorentzRotation toZ;
133  G4LorentzVector Ptmp=projectile->Get4Momentum();
134  toZ.rotateZ(-1*Ptmp.phi());
135  toZ.rotateY(-1*Ptmp.theta());
136  theT1 = toZ * theT1;
137  theT2 = toZ * theT2;
138  thePro = toZ * thePro;
139  G4LorentzRotation toLab(toZ.inverse());
140 
141  // Get definitions
142  const G4ParticleDefinition * d1 = targets[0]->GetDefinition();
143  const G4ParticleDefinition * d2 = targets[1]->GetDefinition();
144  if(0.5>G4UniformRand())
145  {
146  const G4ParticleDefinition * temp;
147  temp=d1;d1=d2;d2=temp;
148  }
149  const G4ParticleDefinition * dp = projectile->GetDefinition();
150  if(dp->GetPDGCharge()<-.5)
151  {
152  if(d1->GetPDGCharge()>.5)
153  {
154  if(d2->GetPDGCharge()>.5 && 0.5>G4UniformRand())
155  {
157  }
158  else
159  {
161  }
162  }
163  else if(d2->GetPDGCharge()>.5)
164  {
166  }
167  }
168  else if(dp->GetPDGCharge()>.5)
169  {
170  if(d1->GetPDGCharge()<.5)
171  {
172  if(d2->GetPDGCharge()<.5 && 0.5>G4UniformRand())
173  {
175  }
176  else
177  {
179  }
180  }
181  else if(d2->GetPDGCharge()<.5)
182  {
184  }
185  }
186 
187  // calculate the momenta.
188  G4double M_sq = (thePro+theT1+theT2).mag2();
189  G4double m1_sq = sqr(d1->GetPDGMass());
190  G4double m2_sq = sqr(d2->GetPDGMass());
191  G4double m_sq = M_sq-m1_sq-m2_sq;
192  G4double p = std::sqrt((m_sq*m_sq - 4.*m1_sq * m2_sq)/(4.*M_sq));
193  G4double costh = 2.*G4UniformRand()-1.;
194  G4double phi = 2.*pi*G4UniformRand();
195  G4ThreeVector pFinal(p*std::sin(std::acos(costh))*std::cos(phi), p*std::sin(std::acos(costh))*std::sin(phi), p*costh);
196 
197  // G4cout << "testing p "<<p-pFinal.mag()<<G4endl;
198  // construct the final state particles lorentz momentum.
199  G4double eFinal1 = std::sqrt(m1_sq+pFinal.mag2());
200  G4LorentzVector final1(pFinal, eFinal1);
201  G4double eFinal2 = std::sqrt(m2_sq+pFinal.mag2());
202  G4LorentzVector final2(-1.*pFinal, eFinal2);
203 
204  // rotate back.
205  final1 = toLab * final1;
206  final2 = toLab * final2;
207 
208  // boost back to LAB
209  final1 = fromSPS * final1;
210  final2 = fromSPS * final2;
211 
212  // make the final state
213  G4KineticTrack * f1 =
214  new G4KineticTrack(d1, 0., targets[0]->GetPosition(), final1);
215  G4KineticTrack * f2 =
216  new G4KineticTrack(d2, 0., targets[1]->GetPosition(), final2);
218  result->push_back(f1);
219  result->push_back(f2);
220 
221  for(size_t hpw=0; hpw<result->size(); hpw++)
222  {
223  energyBalance-=result->operator[](hpw)->Get4Momentum().t();
224  chargeBalance-=G4lrint(result->operator[](hpw)->GetDefinition()->GetPDGCharge());
225  baryonBalance-=result->operator[](hpw)->GetDefinition()->GetBaryonNumber();
226  }
227  if(getenv("AbsorptionEnergyBalanceCheck"))
228  std::cout << "DEBUGGING energy balance B: "
229  <<energyBalance<<" "
230  <<chargeBalance<<" "
231  <<baryonBalance<<" "
232  <<G4endl;
233  return result;
234 }
G4double G4ParticleHPJENDLHEData::G4double result
const char * p
Definition: xmltok.h:285
static const G4double d2
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
HepLorentzVector & rotateZ(double)
int G4int
Definition: G4Types.hh:78
double phi() const
#define G4UniformRand()
Definition: Randomize.hh:97
double theta() const
static const G4double d1
G4double GetPDGMass() const
int G4lrint(double ad)
Definition: templates.hh:163
HepLorentzVector & rotateY(double)
#define G4endl
Definition: G4ios.hh:61
static constexpr double pi
Definition: G4SIunits.hh:75
T sqr(const T &x)
Definition: templates.hh:145
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const
static G4Neutron * NeutronDefinition()
Definition: G4Neutron.cc:99

Here is the call graph for this function:


The documentation for this class was generated from the following files: