74 if ( Pprojectile.
z() < 0.0 )
return false;
75 G4bool PutOnMassShell(
false );
81 if ( M0projectile < MminProjectile ) {
82 PutOnMassShell =
true;
85 G4double M0projectile2 = M0projectile * M0projectile;
95 if ( M0target < MminTarget ) {
96 PutOnMassShell =
true;
99 G4double M0target2 = M0target * M0target;
103 Psum = Pprojectile + Ptarget;
106 if ( Ptmp.pz() <= 0.0 )
return false;
109 toCms.
rotateZ( -1*Ptmp.phi() );
110 toCms.
rotateY( -1*Ptmp.theta() );
119 if ( SqrtS < M0projectile + M0target )
return false;
121 PZcms2 = ( S*S +
sqr( M0projectile2 ) +
sqr( M0target2 )
122 - 2*S*M0projectile2 - 2*S*M0target2 - 2*M0projectile2*M0target2 ) / 4.0 / S;
124 if ( PZcms2 < 0.0 ) {
129 M0projectile2 = M0projectile * M0projectile;
130 PZcms2= ( S*S +
sqr( M0projectile2 ) +
sqr( M0target2 )
131 - 2*S*M0projectile2 - 2*S*M0target2 - 2*M0projectile2*M0target2 ) / 4.0 / S;
132 if ( PZcms2 < 0.0 ) {
return false; }
139 PZcms = std::sqrt( PZcms2 );
141 if ( PutOnMassShell ) {
142 if ( Pprojectile.z() > 0.0 ) {
143 Pprojectile.setPz( PZcms );
144 Ptarget.setPz( -PZcms );
146 Pprojectile.setPz( -PZcms );
147 Ptarget.setPz( PZcms );
149 Pprojectile.setE( std::sqrt( M0projectile2 + Pprojectile.x() * Pprojectile.x() +
150 Pprojectile.y() * Pprojectile.y() + PZcms2 ) );
151 Ptarget.setE( std::sqrt( M0target2 + Ptarget.x() * Ptarget.x() + Ptarget.y() * Ptarget.y() +
163 const G4int maxNumberOfLoops = 1000;
164 G4int loopCounter = 0;
168 ProjMassT2 = M0projectile2 + Pt2;
169 ProjMassT = std::sqrt( ProjMassT2 );
170 TargMassT2 = M0target2 + Pt2;
171 TargMassT = std::sqrt( TargMassT2 );
172 }
while ( ( SqrtS < ProjMassT + TargMassT ) &&
173 ++loopCounter < maxNumberOfLoops );
174 if ( loopCounter >= maxNumberOfLoops ) {
178 PZcms2 = ( S*S +
sqr( ProjMassT2 ) +
sqr( TargMassT2 )
179 - 2.0*S*ProjMassT2 - 2.0*S*TargMassT2 - 2.0*ProjMassT2*TargMassT2 ) / 4.0 / S;
181 if ( PZcms2 < 0.0 ) { PZcms2 = 0.0; };
182 PZcms = std::sqrt( PZcms2 );
183 Pprojectile.setPz( PZcms );
184 Ptarget.setPz( -PZcms );
185 Pprojectile += Qmomentum;
186 Ptarget -= Qmomentum;
215 if ( AveragePt2 <= 0.0 ) {
219 (
G4Exp( -maxPtSquare/AveragePt2 ) -1.0 ) );
223 return G4ThreeVector( Pt * std::cos( phi ), Pt * std::sin( phi ), 0.0 );
231 "G4ElasticHNScattering copy contructor not meant to be called" );
244 "G4ElasticHNScattering = operator not meant to be called" );
252 "G4ElasticHNScattering == operator not meant to be called" );
260 "G4ElasticHNScattering != operator not meant to be called" );
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
const G4ParticleDefinition * GetDefinition() const
HepLorentzVector & transform(const HepRotation &)
G4double GetAvaragePt2ofElasticScattering()
CLHEP::Hep3Vector G4ThreeVector
virtual ~G4ElasticHNScattering()
G4ThreeVector GaussianPt(G4double AveragePt2, G4double maxPtSquare) const
virtual G4bool ElasticScattering(G4VSplitableHadron *aPartner, G4VSplitableHadron *bPartner, G4FTFParameters *theParameters) const
int operator==(const G4ElasticHNScattering &right) const
G4double GetPDGMass() const
G4double G4Log(G4double x)
int operator!=(const G4ElasticHNScattering &right) const
void SetTimeOfCreation(G4double aTime)
const G4ThreeVector & GetPosition() const
static constexpr double twopi
HepLorentzVector & rotateY(double)
HepLorentzVector & rotateZ(double)
void SetPosition(const G4ThreeVector &aPosition)
const G4ElasticHNScattering & operator=(const G4ElasticHNScattering &right)
Hep3Vector boostVector() const
const G4LorentzVector & Get4Momentum() const
CLHEP::HepLorentzVector G4LorentzVector
HepLorentzRotation & transform(const HepBoost &b)
G4double GetTimeOfCreation()
void Set4Momentum(const G4LorentzVector &a4Momentum)
void IncrementCollisionCount(G4int aCount)
G4double GetMinimumMass(const G4ParticleDefinition *p) const