80 #ifdef debugDoubleDiffraction
81 G4cout<<
G4endl<<
"G4QGSDiffractiveExcitation::ExciteParticipants - Double diffraction."<<
G4endl;
94 if(M0projectile < projectile->GetDefinition()->GetPDGMass())
105 if(M0target < target->GetDefinition()->GetPDGMass())
115 if (SqrtS < M0projectile + M0target) {
return false;}
118 G4double Mprojectile2 = M0projectile * M0projectile;
119 G4double Mtarget2 = M0target * M0target;
127 if ( Ptmp.pz() <= 0. )
135 toCms.
rotateY(-1*Ptmp.theta());
142 G4double PZcms2=(S*S+Mprojectile2*Mprojectile2+Mtarget2*Mtarget2-
143 2*S*Mprojectile2-2*S*Mtarget2-2*Mprojectile2*Mtarget2)/4./S;
145 if (PZcms2 < 0) {
return false;}
151 if (Pprojectile.z() > 0.)
153 Pprojectile.setPz( PZcms);
154 Ptarget.setPz( -PZcms);
158 Pprojectile.setPz(-PZcms);
159 Ptarget.setPz( PZcms);
162 Pprojectile.setE(std::sqrt(Mprojectile2+
sqr(Pprojectile.x())+
sqr(Pprojectile.y())+PZcms2));
163 Ptarget.setE( std::sqrt( Mtarget2+
sqr( Ptarget.x())+
sqr( Ptarget.y())+PZcms2));
168 #ifdef debugDoubleDiffraction
169 G4cout <<
"Pprojectile after boost to CMS: " << Pprojectile <<
" "<<Pprojectile.mag()<<
G4endl;
170 G4cout <<
"Ptarget after boost to CMS: " << Ptarget <<
" "<<Ptarget.mag() <<
G4endl;
178 if (M0projectile <= projectile->GetDefinition()->GetPDGMass())
180 if( absPrPDGcode > 1000 )
182 MinPrDiffMass = 1.16;
185 else if( absPrPDGcode == 211 || PrPDGcode == 111)
190 else if( absPrPDGcode == 321 || absPrPDGcode == 130 || absPrPDGcode == 310)
197 MinPrDiffMass = 1.16;
203 MinPrDiffMass = M0projectile + 220.0*
MeV;
207 MinPrDiffMass = MinPrDiffMass *
GeV;
208 AveragePt2 = AveragePt2 * GeV*
GeV;
212 if (SqrtS < MinPrDiffMass + MinTrDiffMass) {
return false;}
214 G4double MinPrDiffMass2 = MinPrDiffMass * MinPrDiffMass;
215 G4double MinTrDiffMass2 = MinTrDiffMass * MinTrDiffMass;
227 if (whilecount++ >= 500 && (whilecount%100)==0)
228 if (whilecount > 1000 ) {
237 ProjMassT2=MinPrDiffMass2+Pt2;
238 ProjMassT =std::sqrt(ProjMassT2);
240 TargMassT2=MinTrDiffMass2+Pt2;
241 TargMassT =std::sqrt(TargMassT2);
243 if (SqrtS < ProjMassT + TargMassT)
continue;
245 PZcms2=(S*S+ProjMassT2*ProjMassT2+TargMassT2*TargMassT2-
246 2.*S*ProjMassT2-2.*S*TargMassT2-2.*ProjMassT2*TargMassT2)/4./S;
247 if (PZcms2 < 0 ) {PZcms2=0;};
248 PZcms =std::sqrt(PZcms2);
250 G4double PMinusMin=std::sqrt(ProjMassT2+PZcms2)-PZcms;
253 PMinusNew=
ChooseP(PMinusMin,PMinusMax);
254 Qminus=PMinusNew-Pprojectile.minus();
256 G4double TPlusMin=std::sqrt(TargMassT2+PZcms2)-PZcms;
259 TPlusNew=
ChooseP(TPlusMin, TPlusMax);
260 Qplus=-(TPlusNew-Ptarget.plus());
262 Qmomentum.
setPz( (Qplus-Qminus)/2 );
263 Qmomentum.
setE( (Qplus+Qminus)/2 );
265 }
while ( (Pprojectile+Qmomentum).mag2() < MinPrDiffMass2 ||
266 (Ptarget -Qmomentum).mag2() < MinTrDiffMass2 );
268 Pprojectile += Qmomentum;
269 Ptarget -= Qmomentum;
275 #ifdef debugDoubleDiffraction
276 G4cout <<
"Pprojectile after boost to Lab: " << Pprojectile <<
" "<<Pprojectile.mag()<<
G4endl;
277 G4cout <<
"Ptarget after boost to Lab: " << Ptarget <<
" "<<Ptarget.mag() <<
G4endl;
293 {
G4cout <<
" G4QGSDiffractiveExcitation::String() Error:No start parton found"<<
G4endl;
298 {
G4cout <<
" G4QGSDiffractiveExcitation::String() Error:No end parton found"<<
G4endl;
331 G4int Sign= isProjectile ? -1 : 1;
333 G4double endMinus = 0.5 * (tm1 + Sign*tm2);
339 Pstart.
setPz(0.5*(startPlus - startMinus));
340 Pstart.
setE(0.5*(startPlus + startMinus));
342 Pend.
setPz(0.5*(endPlus - endMinus));
343 Pend.
setE(0.5*(endPlus + endMinus));
348 #ifdef debugQGSdiffExictation
369 if ( Pmin <= 0. || range <=0. )
371 G4cout <<
" Pmin, range : " << Pmin <<
" , " << range <<
G4endl;
372 throw G4HadronicException(__FILE__, __LINE__,
"G4QGSDiffractiveExcitation::ChooseP : Invalid arguments ");
392 return G4ThreeVector (Pt*std::cos(phi), Pt*std::sin(phi), 0.);