#include <G4SingleDiffractiveExcitation.hh>
Inheritance diagram for G4SingleDiffractiveExcitation:
Public Member Functions | |
G4SingleDiffractiveExcitation (G4double sigmaPt=0.6 *CLHEP::GeV, G4double minExtraMass=250 *CLHEP::MeV, G4double x0mass=250 *CLHEP::MeV) | |
~G4SingleDiffractiveExcitation () | |
G4bool | ExciteParticipants (G4VSplitableHadron *aPartner, G4VSplitableHadron *bPartner) const |
Definition at line 50 of file G4SingleDiffractiveExcitation.hh.
G4SingleDiffractiveExcitation::G4SingleDiffractiveExcitation | ( | G4double | sigmaPt = 0.6 *CLHEP::GeV , |
|
G4double | minExtraMass = 250 *CLHEP::MeV , |
|||
G4double | x0mass = 250 *CLHEP::MeV | |||
) |
Definition at line 49 of file G4SingleDiffractiveExcitation.cc.
00050 : 00051 widthOfPtSquare(-2*sqr(sigmaPt)) , minExtraMass(minextraMass), 00052 minmass(x0mass) 00053 {}
G4SingleDiffractiveExcitation::~G4SingleDiffractiveExcitation | ( | ) |
G4bool G4SingleDiffractiveExcitation::ExciteParticipants | ( | G4VSplitableHadron * | aPartner, | |
G4VSplitableHadron * | bPartner | |||
) | const [virtual] |
Reimplemented from G4QGSDiffractiveExcitation.
Definition at line 59 of file G4SingleDiffractiveExcitation.cc.
References G4UniformRand, G4VSplitableHadron::Get4Momentum(), G4VSplitableHadron::GetDefinition(), G4ParticleDefinition::GetPDGMass(), G4VSplitableHadron::Set4Momentum(), and sqr().
Referenced by G4QGSParticipants::SelectInteractions().
00060 { 00061 00062 G4LorentzVector Pprojectile=projectile->Get4Momentum(); 00063 G4double Mprojectile2=sqr(projectile->GetDefinition()->GetPDGMass() + minExtraMass); 00064 00065 G4LorentzVector Ptarget=target->Get4Momentum(); 00066 G4double Mtarget2=sqr(target->GetDefinition()->GetPDGMass() + minExtraMass); 00067 // G4cout << "E proj, target :" << Pprojectile.e() << ", " << 00068 // Ptarget.e() << G4endl; 00069 00070 G4bool KeepProjectile= G4UniformRand() > 0.5; 00071 00072 // reset the min.mass of the non diffractive particle to its value, ( minus a bit for rounding...) 00073 if ( KeepProjectile ) 00074 { 00075 // cout << " Projectile fix" << G4endl; 00076 Mprojectile2 = sqr(projectile->GetDefinition()->GetPDGMass() * (1-perCent) ); 00077 } else { 00078 // cout << " Target fix" << G4endl; 00079 Mtarget2=sqr(target->GetDefinition()->GetPDGMass() * (1-perCent) ); 00080 } 00081 00082 // Transform momenta to cms and then rotate parallel to z axis; 00083 00084 G4LorentzVector Psum; 00085 Psum=Pprojectile+Ptarget; 00086 00087 G4LorentzRotation toCms(-1*Psum.boostVector()); 00088 00089 G4LorentzVector Ptmp=toCms*Pprojectile; 00090 00091 if ( Ptmp.pz() <= 0. ) 00092 { 00093 // "String" moving backwards in CMS, abort collision !! 00094 // G4cout << " abort Collision!! " << G4endl; 00095 return false; 00096 } 00097 00098 toCms.rotateZ(-1*Ptmp.phi()); 00099 toCms.rotateY(-1*Ptmp.theta()); 00100 00101 // G4cout << "Pprojectile be4 boost " << Pprojectile << G4endl; 00102 // G4cout << "Ptarget be4 boost : " << Ptarget << G4endl; 00103 00104 00105 00106 G4LorentzRotation toLab(toCms.inverse()); 00107 00108 Pprojectile.transform(toCms); 00109 Ptarget.transform(toCms); 00110 00111 G4LorentzVector Qmomentum; 00112 G4int whilecount=0; 00113 do { 00114 // Generate pt 00115 00116 G4double maxPtSquare=sqr(Ptarget.pz()); 00117 if (whilecount++ >= 500 && (whilecount%100)==0) 00118 // G4cout << "G4SingleDiffractiveExcitation::ExciteParticipants possibly looping" 00119 // << ", loop count/ maxPtSquare : " 00120 // << whilecount << " / " << maxPtSquare << G4endl; 00121 if (whilecount > 1000 ) 00122 { 00123 Qmomentum=G4LorentzVector(0.,0.,0.,0.); 00124 // G4cout << "G4SingleDiffractiveExcitation::ExciteParticipants: Aborting loop!" << G4endl; 00125 return false; // Ignore this interaction 00126 } 00127 Qmomentum=G4LorentzVector(GaussianPt(widthOfPtSquare,maxPtSquare),0); 00128 00129 00130 // Momentum transfer 00131 G4double Xmin = minmass / ( Pprojectile.e() + Ptarget.e() ); 00132 G4double Xmax=1.; 00133 G4double Xplus =ChooseX(Xmin,Xmax); 00134 G4double Xminus=ChooseX(Xmin,Xmax); 00135 00136 G4double pt2=G4ThreeVector(Qmomentum.vect()).mag2(); 00137 G4double Qplus =-1 * pt2 / Xminus/Ptarget.minus(); 00138 G4double Qminus= pt2 / Xplus /Pprojectile.plus(); 00139 00140 if ( KeepProjectile ) 00141 { 00142 Qminus = (sqr(projectile->GetDefinition()->GetPDGMass()) + pt2 ) 00143 / (Pprojectile.plus() + Qplus ) 00144 - Pprojectile.minus(); 00145 } else 00146 { 00147 Qplus = Ptarget.plus() 00148 - (sqr(target->GetDefinition()->GetPDGMass()) + pt2 ) 00149 / (Ptarget.minus() - Qminus ); 00150 } 00151 00152 Qmomentum.setPz( (Qplus-Qminus)/2 ); 00153 Qmomentum.setE( (Qplus+Qminus)/2 ); 00154 00155 // G4cout << "Qplus / Qminus " << Qplus << " / " << Qminus<<G4endl; 00156 // G4cout << "pt2 " << pt2 << G4endl; 00157 // G4cout << "Qmomentum " << Qmomentum << G4endl; 00158 // G4cout << " Masses (P/T) : " << (Pprojectile+Qmomentum).mag() << 00159 // " / " << (Ptarget-Qmomentum).mag() << G4endl; 00160 00161 } while ( (Ptarget-Qmomentum).mag2() <= Mtarget2 00162 || (Pprojectile+Qmomentum).mag2() <= Mprojectile2 00163 || (Ptarget-Qmomentum).e() < 0. 00164 || (Pprojectile+Qmomentum).e() < 0. ); 00165 00166 00167 // G4double Ecms=Pprojectile.e() + Ptarget.e(); 00168 00169 Pprojectile += Qmomentum; 00170 00171 Ptarget -= Qmomentum; 00172 00173 // G4cout << "Pprojectile.e() : " << Pprojectile.e() << G4endl; 00174 // G4cout << "Ptarget.e() : " << Ptarget.e() << G4endl; 00175 00176 // G4cout << "end event_______________________________________________"<<G4endl; 00177 // 00178 00179 00180 // G4cout << "Pprojectile with Q : " << Pprojectile << G4endl; 00181 // G4cout << "Ptarget with Q : " << Ptarget << G4endl; 00182 // G4cout << "Projectile back: " << toLab * Pprojectile << G4endl; 00183 // G4cout << "Target back: " << toLab * Ptarget << G4endl; 00184 00185 // Transform back and update SplitableHadron Participant. 00186 Pprojectile.transform(toLab); 00187 Ptarget.transform(toLab); 00188 00189 // G4cout << "G4SingleDiffractiveExcitation- Target mass " << Ptarget.mag() << G4endl; 00190 // G4cout << "G4SingleDiffractiveExcitation- Projectile mass " << Pprojectile.mag() << G4endl; 00191 00192 target->Set4Momentum(Ptarget); 00193 projectile->Set4Momentum(Pprojectile); 00194 00195 00196 return true; 00197 }