#include <G4FTFAnnihilation.hh>
Public Member Functions | |
| G4FTFAnnihilation () | |
| virtual | ~G4FTFAnnihilation () |
| virtual G4bool | Annihilate (G4VSplitableHadron *aPartner, G4VSplitableHadron *bPartner, G4VSplitableHadron *&AdditionalString, G4FTFParameters *theParameters) const |
Definition at line 47 of file G4FTFAnnihilation.hh.
| G4FTFAnnihilation::G4FTFAnnihilation | ( | ) |
| G4FTFAnnihilation::~G4FTFAnnihilation | ( | ) | [virtual] |
| G4bool G4FTFAnnihilation::Annihilate | ( | G4VSplitableHadron * | aPartner, | |
| G4VSplitableHadron * | bPartner, | |||
| G4VSplitableHadron *& | AdditionalString, | |||
| G4FTFParameters * | theParameters | |||
| ) | const [virtual] |
Definition at line 76 of file G4FTFAnnihilation.cc.
References G4cout, G4endl, G4UniformRand, G4VSplitableHadron::Get4Momentum(), G4FTFParameters::GetAveragePt2(), G4VSplitableHadron::GetDefinition(), G4ParticleDefinition::GetPDGEncoding(), G4ParticleDefinition::GetPDGMass(), G4VSplitableHadron::GetPosition(), G4VSplitableHadron::GetTimeOfCreation(), G4VSplitableHadron::IncrementCollisionCount(), G4INCL::Math::pi, G4VSplitableHadron::Set4Momentum(), G4VSplitableHadron::SetFirstParton(), G4VSplitableHadron::SetPosition(), G4VSplitableHadron::SetSecondParton(), G4VSplitableHadron::SetStatus(), G4VSplitableHadron::SetTimeOfCreation(), G4VSplitableHadron::SplitUp(), and sqr().
00080 { 00081 // -------------------- Projectile parameters ----------------------- 00082 G4LorentzVector Pprojectile=projectile->Get4Momentum(); 00083 //G4cout<<"---------------------------- Annihilation----------------"<<G4endl; 00084 //G4cout<<"Pprojectile "<<Pprojectile<<G4endl; 00085 //G4cout<<"Pprojectile.mag2 "<<Pprojectile.mag2()<<G4endl; 00086 G4int ProjectilePDGcode=projectile->GetDefinition()->GetPDGEncoding(); 00087 if(ProjectilePDGcode > 0) 00088 { 00089 target->SetStatus(2); 00090 return false; 00091 } 00092 00093 // G4double M0projectile = Pprojectile.mag(); 00094 // G4double M0projectile2= projectile->GetDefinition()->GetPDGMass()* 00095 // projectile->GetDefinition()->GetPDGMass(); 00096 G4double M0projectile2=Pprojectile.mag2(); 00097 // -------------------- Target parameters ------------------------- 00098 G4int TargetPDGcode=target->GetDefinition()->GetPDGEncoding(); 00099 00100 G4LorentzVector Ptarget=target->Get4Momentum(); 00101 //G4cout<<"Ptarget "<<Ptarget<<G4endl; 00102 //G4cout<<"Ptarget.mag2 "<<Ptarget.mag2()<<G4endl; 00103 00104 // G4double M0target = Ptarget.mag(); 00105 // G4double M0target2= target->GetDefinition()->GetPDGMass()* 00106 // target->GetDefinition()->GetPDGMass(); 00107 G4double M0target2=Ptarget.mag2(); 00108 00109 //G4cout<<"Annihilate "<<ProjectilePDGcode<<" "<<TargetPDGcode<<G4endl; 00110 //G4cout<<"Pprojec "<<Pprojectile<<" "<<Pprojectile.mag2()<<G4endl; 00111 //G4cout<<"Ptarget "<<Ptarget <<" "<<Ptarget.mag2() <<G4endl; 00112 //G4cout<<"M0 proj target "<<M0projectile<<" "<<M0target<<G4endl; 00113 00114 G4double AveragePt2=theParameters->GetAveragePt2(); 00115 00116 // Kinematical properties of the interactions -------------- 00117 G4LorentzVector Psum; // 4-momentum in CMS 00118 Psum=Pprojectile+Ptarget; 00119 G4double S=Psum.mag2(); 00120 //G4cout<<"Psum S"<<Psum<<" "<<S<<G4endl; 00121 // Transform momenta to cms and then rotate parallel to z axis; 00122 G4LorentzRotation toCms(-1*Psum.boostVector()); 00123 //G4cout<<"G4LorentzRotation toCms(-1*Psum.boostVector());"<<G4endl; 00124 G4LorentzVector Ptmp=toCms*Pprojectile; 00125 00126 /* // For anti-baryons it is not needed ! 00127 if ( Ptmp.pz() <= 0. ) 00128 { 00129 target->SetStatus(2); 00130 // "String" moving backwards in CMS, abort collision !! 00131 return false; 00132 } 00133 */ 00134 00135 toCms.rotateZ(-1*Ptmp.phi()); 00136 toCms.rotateY(-1*Ptmp.theta()); 00137 00138 G4LorentzRotation toLab(toCms.inverse()); 00139 00140 G4double SqrtS=std::sqrt(S); 00141 00142 G4double maxPtSquare; 00143 00144 //G4cout<<"M0projectile+M0target Sqrt(S) (GeV) "<<M0projectile2/GeV<<" "<<M0target2/GeV<<" "<<(M0projectile2+M0target2)/GeV<<" "<<SqrtS/GeV<<G4endl; 00145 00146 G4double X_a(0.), X_b(0.), X_c(0.), X_d(0.); 00147 G4double MesonProdThreshold=projectile->GetDefinition()->GetPDGMass()+ 00148 target->GetDefinition()->GetPDGMass()+ 00149 (2.*140.+16.)*MeV; // 2 Mpi +DeltaE 00150 00151 G4double Prel2= S*S + M0projectile2*M0projectile2 + M0target2*M0target2 - 00152 2.*S*M0projectile2 - 2.*S*M0target2 - 2.*M0projectile2*M0target2; 00153 Prel2/=S; 00154 00155 if(Prel2 < 0. ) // *MeV*MeV 1600. 00156 { // Annihilation at rest! Values are copied from Paratemets. 00157 X_a= 625.1; // mb // 3-shirt diagram 00158 X_b= 9.780; // mb // anti-quark-quark annihilation 00159 X_c= 49.989; // mb 00160 X_d= 6.614; // mb 00161 } 00162 else 00163 { // Annihilation in flight! 00164 G4double FlowF=1./std::sqrt(Prel2)*GeV; 00165 00166 //G4cout<<"Annig FlowF "<<FlowF<<" sqrt "<<SqrtS/GeV<<G4endl; 00167 00168 // Process cross sections --------------------------------------------------- 00169 X_a=25.*FlowF; // mb 3-shirt diagram 00170 00171 // mb anti-quark-quark annihilation 00172 if(SqrtS < MesonProdThreshold) 00173 { 00174 X_b=3.13+140.*std::pow((MesonProdThreshold - SqrtS)/GeV,2.5); 00175 } 00176 else 00177 { 00178 X_b=6.8*GeV/SqrtS; 00179 } 00180 if(projectile->GetDefinition()->GetPDGMass()+ 00181 target->GetDefinition()->GetPDGMass() > SqrtS) {X_b=0.;} 00182 // This can be in an interaction of low energy anti-baryon with off-shell nuclear nucleon 00183 00184 // ???????????????????????????????????????? 00185 X_c=2.*FlowF*sqr(projectile->GetDefinition()->GetPDGMass()+ 00186 target->GetDefinition()->GetPDGMass())/S; 00187 // mb re-arrangement of 2 quarks and 2 anti-quarks 00188 // ???????????????????????????????????????? 00189 X_d=23.3*GeV*GeV/S; // mb anti-quark-quark string creation 00190 } // end of if(Prel2 < 1600. ) // *MeV*MeV 00191 00192 //G4cout<<"Annih X a b c d "<<X_a<<" "<<X_b<<" "<<X_c<<" "<<X_d<<G4endl; 00193 00194 if((ProjectilePDGcode == -2212)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00195 {X_b*=5.; X_c*=5.; X_d*=6.;} // Pbar P 00196 else if((ProjectilePDGcode == -2212)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00197 {X_b*=4.; X_c*=4.; X_d*=4.;} // Pbar N 00198 else if((ProjectilePDGcode == -2112)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00199 {X_b*=4.; X_c*=4.; X_d*=4.;} // NeutrBar P 00200 else if((ProjectilePDGcode == -2112)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00201 {X_b*=5.; X_c*=5.; X_d*=6.;} // NeutrBar N 00202 else if((ProjectilePDGcode == -3122)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00203 {X_b*=3.; X_c*=3.; X_d*=2.;} // LambdaBar P 00204 else if((ProjectilePDGcode == -3122)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00205 {X_b*=3.; X_c*=3.; X_d*=2.;} // LambdaBar N 00206 else if((ProjectilePDGcode == -3112)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00207 {X_b*=2.; X_c*=2.; X_d*=0.;} // Sigma-Bar P 00208 else if((ProjectilePDGcode == -3112)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00209 {X_b*=4.; X_c*=4.; X_d*=2.;} // Sigma-Bar N 00210 else if((ProjectilePDGcode == -3212)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00211 {X_b*=3.; X_c*=3.; X_d*=2.;} // Sigma0Bar P 00212 else if((ProjectilePDGcode == -3212)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00213 {X_b*=3.; X_c*=3.; X_d*=2.;} // Sigma0Bar N 00214 else if((ProjectilePDGcode == -3222)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00215 {X_b*=4.; X_c*=4.; X_d*=2.;} // Sigma+Bar P 00216 else if((ProjectilePDGcode == -3222)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00217 {X_b*=2.; X_c*=2.; X_d*=0.;} // Sigma+Bar N 00218 else if((ProjectilePDGcode == -3312)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00219 {X_b*=1.; X_c*=1.; X_d*=0.;} // Xi-Bar P 00220 else if((ProjectilePDGcode == -3312)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00221 {X_b*=2.; X_c*=2.; X_d*=0.;} // Xi-Bar N 00222 else if((ProjectilePDGcode == -3322)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00223 {X_b*=2.; X_c*=2.; X_d*=0.;} // Xi0Bar P 00224 else if((ProjectilePDGcode == -3322)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00225 {X_b*=1.; X_c*=1.; X_d*=0.;} // Xi0Bar N 00226 else if((ProjectilePDGcode == -3334)&&((TargetPDGcode == 2212)||(TargetPDGcode == 2214))) 00227 {X_b*=0.; X_c*=0.; X_d*=0.;} // Omega-Bar P 00228 else if((ProjectilePDGcode == -3334)&&((TargetPDGcode == 2112)||(TargetPDGcode == 2114))) 00229 {X_b*=0.; X_c*=0.; X_d*=0.;} // Omega-Bar N 00230 else {G4cout<<"Unknown anti-baryon for FTF annihilation: PDGcodes - "<<ProjectilePDGcode<<" "<<TargetPDGcode<<G4endl;} 00231 00232 //G4cout<<"Annih X a b c d "<<X_a<<" "<<X_b<<" "<<X_c<<" "<<X_d<<G4endl; 00233 //========================================= 00234 //X_a=0.; 00235 //X_b=0.; 00236 //X_c=0.; 00237 //X_d=0.; 00238 //G4cout<<"Annih X a b c d "<<X_a<<" "<<X_b<<" "<<X_c<<" "<<X_d<<G4endl; 00239 //========================================= 00240 G4double Xannihilation=X_a+X_b+X_c+X_d; 00241 // ------------------------------------------------------ 00242 // ------ Projectile unpacking -------------------------- 00243 G4int AQ[3]; 00244 UnpackBaryon(ProjectilePDGcode, AQ[0], AQ[1], AQ[2]); 00245 00246 // ------ Target unpacking ------------------------------ 00247 G4int Q[3]; 00248 UnpackBaryon(TargetPDGcode, Q[0], Q[1], Q[2]); 00249 00250 // ------------------------------------------------------ 00251 G4double Ksi=G4UniformRand(); 00252 00253 if(Ksi < X_a/Xannihilation) 00254 { 00255 //============================================================ 00256 // Simulation of 3 anti-quark-quark strings creation 00257 // Sampling of anti-quark order in projectile 00258 //G4cout<<"Process a"<<G4endl; 00259 G4int SampledCase=CLHEP::RandFlat::shootInt(G4long(6)); 00260 00261 G4int Tmp1(0), Tmp2(0); 00262 if(SampledCase == 0) { } 00263 if(SampledCase == 1) {Tmp1=AQ[1]; AQ[1]=AQ[2]; AQ[2]=Tmp1;} 00264 if(SampledCase == 2) {Tmp1=AQ[0]; AQ[0]=AQ[1]; AQ[1]=Tmp1;} 00265 if(SampledCase == 3) {Tmp1=AQ[0]; Tmp2=AQ[1]; AQ[0]=AQ[2]; AQ[1]=Tmp1; AQ[2]=Tmp2;} 00266 if(SampledCase == 4) {Tmp1=AQ[0]; Tmp2=AQ[1]; AQ[0]=Tmp2; AQ[1]=AQ[2]; AQ[2]=Tmp1;} 00267 if(SampledCase == 5) {Tmp1=AQ[0]; Tmp2=AQ[1]; AQ[0]=AQ[2]; AQ[1]=Tmp2; AQ[2]=Tmp1;} 00268 00269 // --------------- Set the string properties --------------- 00270 //G4cout<<"String 1 "<<AQ[0]<<" "<<Q[0]<<G4endl; 00271 projectile->SplitUp(); 00272 00273 projectile->SetFirstParton(AQ[0]); 00274 projectile->SetSecondParton(Q[0]); 00275 projectile->SetStatus(1); 00276 00277 //G4cout<<"String 2 "<<Q[1]<<" "<<AQ[1]<<G4endl; 00278 target->SplitUp(); 00279 00280 target->SetFirstParton(Q[1]); 00281 target->SetSecondParton(AQ[1]); 00282 target->SetStatus(1); 00283 00284 //G4cout<<"String 3 "<<AQ[2]<<" "<<Q[2]<<G4endl; 00285 AdditionalString=new G4DiffractiveSplitableHadron(); 00286 AdditionalString->SplitUp(); 00287 AdditionalString->SetFirstParton(AQ[2]); 00288 AdditionalString->SetSecondParton(Q[2]); 00289 AdditionalString->SetStatus(1); 00290 //G4cout<<G4endl<<"*AdditionalString in Annih"<<AdditionalString<<G4endl; 00291 00292 // Sampling kinematical properties 00293 // 1 string AQ[0]-Q[0]// 2 string AQ[1]-Q[1]// 3 string AQ[2]-Q[2] 00294 00295 G4ThreeVector Quark_Mom[6]; 00296 G4double ModMom2[6]; //ModMom[6], 00297 00298 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00299 AveragePt2=200.*200.; maxPtSquare=S; 00300 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00301 00302 G4double SumMt(0.); 00303 00304 G4double MassQ2=0.; //100.*100.*MeV*MeV; 00305 00306 G4int NumberOfTries(0); 00307 G4double ScaleFactor(1.); 00308 do // while(SumMt >SqrtS) 00309 { 00310 NumberOfTries++; 00311 00312 if(NumberOfTries == 100*(NumberOfTries/100)) 00313 { // At large number of tries it would be better to reduce the values of <Pt^2> 00314 ScaleFactor/=2.; 00315 AveragePt2 *=ScaleFactor; 00316 } 00317 00318 G4ThreeVector PtSum(0.,0.,0.); 00319 for(G4int i=0; i<6; i++) 00320 { 00321 Quark_Mom[i]=GaussianPt(AveragePt2, maxPtSquare); 00322 PtSum+=Quark_Mom[i]; 00323 } 00324 00325 PtSum/=6.; 00326 00327 SumMt=0.; 00328 for(G4int i=0; i<6; i++) 00329 { 00330 Quark_Mom[i]-=PtSum; 00331 // ModMom[i] =Quark_Mom[i].mag(); 00332 ModMom2[i]=Quark_Mom[i].mag2(); 00333 SumMt+=std::sqrt(ModMom2[i]+MassQ2); 00334 } 00335 } while(SumMt > SqrtS); 00336 00337 G4double WminusTarget(0.), WplusProjectile(0.); 00338 /* //--------------------- Closed is variant with sampling of Xs at minimum 00339 G4double SumMod_anti=ModMom[0]+ModMom[1]+ModMom[2]; 00340 Quark_Mom[0].setZ(ModMom[0]/SumMod_anti); 00341 Quark_Mom[1].setZ(ModMom[1]/SumMod_anti); 00342 Quark_Mom[2].setZ(ModMom[2]/SumMod_anti); 00343 00344 G4double SumMod_bary=ModMom[3]+ModMom[4]+ModMom[5]; 00345 Quark_Mom[3].setZ(ModMom[3]/SumMod_bary); 00346 Quark_Mom[4].setZ(ModMom[4]/SumMod_bary); 00347 Quark_Mom[5].setZ(ModMom[5]/SumMod_bary); 00348 00349 G4double Alfa=SumMod_anti*SumMod_anti; 00350 G4double Beta=SumMod_bary*SumMod_bary; 00351 //------------------------------------ 00352 G4double DecayMomentum2=S*S + Alfa*Alfa + Beta*Beta 00353 - 2.*S*Alfa - 2.*S*Beta - 2.*Alfa*Beta; 00354 00355 WminusTarget=(S-Alfa+Beta+std::sqrt(DecayMomentum2))/2./SqrtS; 00356 WplusProjectile=SqrtS-Beta/WminusTarget; 00357 */ //--------------------- Closed is variant with sampling of Xs at minimum 00358 // 00359 // // ------------------------------------------------ 00360 // Sampling X's of anti-baryon ------- 00361 G4double Alfa_R=0.5; 00362 00363 NumberOfTries=0; 00364 ScaleFactor=1.; 00365 00366 G4bool Succes(true); 00367 do // while(!Succes) 00368 { 00369 Succes=true; 00370 NumberOfTries++; 00371 00372 if(NumberOfTries == 100*(NumberOfTries/100)) 00373 { // At large number of tries it would be better to reduce the values of Pt's 00374 ScaleFactor/=2.; 00375 } 00376 00377 if(Alfa_R == 1.) 00378 { 00379 G4double Xaq1=1.-std::sqrt(G4UniformRand()); 00380 G4double Xaq2=(1.-Xaq1)*G4UniformRand(); 00381 G4double Xaq3=1.-Xaq1-Xaq2; 00382 00383 Quark_Mom[0].setZ(Xaq1); Quark_Mom[1].setZ(Xaq2); Quark_Mom[2].setZ(Xaq3); 00384 } 00385 else 00386 { 00387 G4double Xaq1=sqr(G4UniformRand()); 00388 G4double Xaq2=(1.-Xaq1)*sqr(std::sin(pi/2.*G4UniformRand())); 00389 G4double Xaq3=1.-Xaq1-Xaq2; 00390 00391 Quark_Mom[0].setZ(Xaq1); Quark_Mom[1].setZ(Xaq2); Quark_Mom[2].setZ(Xaq3); 00392 } // end of if(Alfa_R == 0.) 00393 00394 // Sampling X's of baryon ------------ 00395 if(Alfa_R == 1.) 00396 { 00397 G4double Xq1=1.-std::sqrt(G4UniformRand()); 00398 G4double Xq2=(1.-Xq1)*G4UniformRand(); 00399 G4double Xq3=1.-Xq1-Xq2; 00400 00401 Quark_Mom[3].setZ(Xq1); Quark_Mom[4].setZ(Xq2); Quark_Mom[5].setZ(Xq3); 00402 } 00403 else 00404 { 00405 G4double Xq1=sqr(G4UniformRand()); 00406 G4double Xq2=(1.-Xq1)*sqr(std::sin(pi/2.*G4UniformRand())); 00407 G4double Xq3=1.-Xq1-Xq2; 00408 00409 Quark_Mom[3].setZ(Xq1); Quark_Mom[4].setZ(Xq2); Quark_Mom[5].setZ(Xq3); 00410 } // end of if(Alfa_R == 0.) 00411 // 00412 G4double Alfa(0.), Beta(0.); 00413 00414 for(G4int i=0; i<3; i++) // For Anti-baryon 00415 { 00416 if(Quark_Mom[i].getZ() != 0.) 00417 {Alfa+=(ScaleFactor*ModMom2[i]+MassQ2)/Quark_Mom[i].getZ();} 00418 else {Succes=false;} 00419 } 00420 00421 for(G4int i=3; i<6; i++) // For baryon 00422 { 00423 if(Quark_Mom[i].getZ() != 0.) 00424 {Beta+=(ScaleFactor*ModMom2[i]+MassQ2)/Quark_Mom[i].getZ();} 00425 else {Succes=false;} 00426 } 00427 00428 if(!Succes) continue; 00429 00430 if(std::sqrt(Alfa)+std::sqrt(Beta) > SqrtS) {Succes=false; continue;} 00431 00432 G4double DecayMomentum2=S*S + Alfa*Alfa + Beta*Beta 00433 - 2.*S*Alfa - 2.*S*Beta - 2.*Alfa*Beta; 00434 00435 WminusTarget=(S-Alfa+Beta+std::sqrt(DecayMomentum2))/2./SqrtS; 00436 WplusProjectile=SqrtS-Beta/WminusTarget; 00437 00438 } while(!Succes); 00439 // //-------------------------------------------------- 00440 00441 G4double SqrtScaleF=std::sqrt(ScaleFactor); 00442 00443 for(G4int i=0; i<3; i++) 00444 { 00445 G4double Pz=WplusProjectile*Quark_Mom[i].getZ()/2.- 00446 (ScaleFactor*ModMom2[i]+MassQ2)/(2.*WplusProjectile*Quark_Mom[i].getZ()); 00447 Quark_Mom[i].setZ(Pz); 00448 00449 if(ScaleFactor != 1.) 00450 { 00451 Quark_Mom[i].setX(SqrtScaleF*Quark_Mom[i].getX()); 00452 Quark_Mom[i].setY(SqrtScaleF*Quark_Mom[i].getY()); 00453 } 00454 } 00455 00456 for(G4int i=3; i<6; i++) 00457 { 00458 G4double Pz=-WminusTarget*Quark_Mom[i].getZ()/2.+ 00459 (ScaleFactor*ModMom2[i]+MassQ2)/(2.*WminusTarget*Quark_Mom[i].getZ()); 00460 Quark_Mom[i].setZ(Pz); 00461 00462 if(ScaleFactor != 1.) 00463 { 00464 Quark_Mom[i].setX(SqrtScaleF*Quark_Mom[i].getX()); 00465 Quark_Mom[i].setY(SqrtScaleF*Quark_Mom[i].getY()); 00466 } 00467 } 00468 //G4cout<<"Sum AQ "<<Quark_Mom[0]+Quark_Mom[1]+Quark_Mom[2]<<G4endl; 00469 //G4cout<<"Sum Q "<<Quark_Mom[3]+Quark_Mom[4]+Quark_Mom[5]<<G4endl; 00470 //------------------------------------- 00471 00472 G4ThreeVector tmp=Quark_Mom[0]+Quark_Mom[3]; 00473 G4LorentzVector Pstring1(tmp,std::sqrt(Quark_Mom[0].mag2()+MassQ2)+ 00474 std::sqrt(Quark_Mom[3].mag2()+MassQ2)); 00475 G4double Ystring1=Pstring1.rapidity(); 00476 /* 00477 G4cout<<"Mom 1 string "<<G4endl; 00478 G4cout<<Quark_Mom[0]<<G4endl; 00479 G4cout<<Quark_Mom[3]<<G4endl; 00480 G4cout<<tmp<<" "<<tmp.mag()<<G4endl; 00481 */ 00482 //G4cout<<"1 str "<<Pstring1<<" "<<Pstring1.mag()<<" "<<Ystring1<<G4endl; 00483 00484 tmp=Quark_Mom[1]+Quark_Mom[4]; 00485 G4LorentzVector Pstring2(tmp,std::sqrt(Quark_Mom[1].mag2()+MassQ2)+ 00486 std::sqrt(Quark_Mom[4].mag2()+MassQ2)); 00487 G4double Ystring2=Pstring2.rapidity(); 00488 /* 00489 G4cout<<"Mom 2 string "<<G4endl; 00490 G4cout<<Quark_Mom[1]<<G4endl; 00491 G4cout<<Quark_Mom[4]<<G4endl; 00492 G4cout<<tmp<<" "<<tmp.mag()<<G4endl; 00493 */ 00494 //G4cout<<"2 str "<<Pstring2<<" "<<Pstring2.mag()<<" "<<Ystring2<<G4endl; 00495 00496 tmp=Quark_Mom[2]+Quark_Mom[5]; 00497 G4LorentzVector Pstring3(tmp,std::sqrt(Quark_Mom[2].mag2()+MassQ2)+ 00498 std::sqrt(Quark_Mom[5].mag2()+MassQ2)); 00499 G4double Ystring3=Pstring3.rapidity(); 00500 /* 00501 G4cout<<"Mom 3 string "<<G4endl; 00502 G4cout<<Quark_Mom[2]<<G4endl; 00503 G4cout<<Quark_Mom[5]<<G4endl; 00504 G4cout<<tmp<<" "<<tmp.mag()<<G4endl; 00505 */ 00506 //G4cout<<"3 str "<<Pstring3<<" "<<Pstring3.mag()<<" "<<Ystring3<<G4endl; 00507 //G4cout<<"SumE "<<Pstring1.e()+Pstring2.e()+Pstring3.e()<<G4endl; 00508 //G4cout<<Pstring1.mag()<<" "<<Pstring2.mag()<<" "<<Pstring3.mag()<<G4endl; 00509 //G4int Uzhi; G4cin>>Uzhi; 00510 //-------------------------------- 00511 G4LorentzVector LeftString(0.,0.,0.,0.); 00512 //----- 00513 if((Ystring1 > Ystring2)&&(Ystring2 > Ystring3)) 00514 { 00515 Pprojectile=Pstring1; 00516 LeftString =Pstring2; 00517 Ptarget =Pstring3; 00518 } 00519 00520 if((Ystring1 > Ystring3)&&(Ystring3 > Ystring2)) 00521 { 00522 Pprojectile=Pstring1; 00523 LeftString =Pstring3; 00524 Ptarget =Pstring2; 00525 } 00526 //----- 00527 if((Ystring2 > Ystring1)&&(Ystring1 > Ystring3)) 00528 { 00529 Pprojectile=Pstring2; 00530 LeftString =Pstring1; 00531 Ptarget =Pstring3; 00532 } 00533 00534 if((Ystring2 > Ystring3)&&(Ystring3 > Ystring1)) 00535 { 00536 Pprojectile=Pstring2; 00537 LeftString =Pstring3; 00538 Ptarget =Pstring1; 00539 } 00540 //----- 00541 if((Ystring3 > Ystring1)&&(Ystring1 > Ystring2)) 00542 { 00543 Pprojectile=Pstring3; 00544 LeftString =Pstring1; 00545 Ptarget =Pstring2; 00546 } 00547 00548 if((Ystring3 > Ystring2)&&(Ystring2 > Ystring1)) 00549 { 00550 Pprojectile=Pstring3; 00551 LeftString =Pstring2; 00552 Ptarget =Pstring1; 00553 } 00554 00555 //------------------------------------------------------- 00556 //G4cout<<"SumP "<<Pprojectile+LeftString+Ptarget<<" "<<SqrtS<<G4endl; 00557 00558 Pprojectile.transform(toLab); 00559 LeftString.transform(toLab); 00560 Ptarget.transform(toLab); 00561 //G4cout<<"SumP "<<Pprojectile+LeftString+Ptarget<<" "<<SqrtS<<G4endl; 00562 00563 // Calculation of the creation time --------------------- 00564 projectile->SetTimeOfCreation(target->GetTimeOfCreation()); 00565 projectile->SetPosition(target->GetPosition()); 00566 00567 AdditionalString->SetTimeOfCreation(target->GetTimeOfCreation()); 00568 AdditionalString->SetPosition(target->GetPosition()); 00569 // Creation time and position of target nucleon were determined at 00570 // ReggeonCascade() of G4FTFModel 00571 // ------------------------------------------------------ 00572 00573 //G4cout<<"Mproj "<<Pprojectile.mag()<<G4endl; 00574 //G4cout<<"Mtarg "<<Ptarget.mag()<<G4endl; 00575 projectile->Set4Momentum(Pprojectile); 00576 AdditionalString->Set4Momentum(LeftString); 00577 target->Set4Momentum(Ptarget); 00578 00579 projectile->IncrementCollisionCount(1); 00580 AdditionalString->IncrementCollisionCount(1); 00581 target->IncrementCollisionCount(1); 00582 00583 return true; 00584 } 00585 00586 //============================================================ 00587 // Simulation of anti-diquark-diquark string creation 00588 // 00589 if(Ksi < (X_a+X_b)/Xannihilation) 00590 { 00591 //G4cout<<"Process b"<<G4endl; 00592 G4int CandidatsN(0), CandAQ[9][2], CandQ[9][2]; 00593 G4int LeftAQ1(0), LeftAQ2(0), LeftQ1(0), LeftQ2(0); 00594 //------------------------------------------------------------ 00595 for(G4int iAQ=0; iAQ<3; iAQ++) 00596 { 00597 for(G4int iQ=0; iQ<3; iQ++) 00598 { 00599 if(-AQ[iAQ] == Q[iQ]) 00600 { 00601 if(iAQ == 0) {CandAQ[CandidatsN][0]=1; CandAQ[CandidatsN][1]=2;} 00602 if(iAQ == 1) {CandAQ[CandidatsN][0]=0; CandAQ[CandidatsN][1]=2;} 00603 if(iAQ == 2) {CandAQ[CandidatsN][0]=0; CandAQ[CandidatsN][1]=1;} 00604 if(iQ == 0) {CandQ[CandidatsN][0] =1; CandQ[CandidatsN][1]=2;} 00605 if(iQ == 1) {CandQ[CandidatsN][0] =0; CandQ[CandidatsN][1]=2;} 00606 if(iQ == 2) {CandQ[CandidatsN][0] =0; CandQ[CandidatsN][1]=1;} 00607 CandidatsN++; 00608 } //end of if(-AQ[i] == Q[j]) 00609 } //end of cycle on targ. quarks 00610 } //end of cycle on proj. anti-quarks 00611 //------------------------------------------------------------ 00612 //G4cout<<"CandidatsN "<<CandidatsN<<G4endl; 00613 00614 if(CandidatsN != 0) 00615 { 00616 G4int SampledCase=CLHEP::RandFlat::shootInt(G4long(CandidatsN)); 00617 00618 LeftAQ1=AQ[CandAQ[SampledCase][0]]; 00619 LeftAQ2=AQ[CandAQ[SampledCase][1]]; 00620 00621 LeftQ1=Q[CandQ[SampledCase][0]]; 00622 LeftQ2=Q[CandQ[SampledCase][1]]; 00623 00624 // -------- Build anti-diquark and diquark 00625 G4int Anti_DQ(0), DQ(0); 00626 00627 if(std::abs(LeftAQ1) > std::abs(LeftAQ2)) 00628 { 00629 Anti_DQ=1000*LeftAQ1+100*LeftAQ2-3; // 1 00630 } else 00631 { 00632 Anti_DQ=1000*LeftAQ2+100*LeftAQ1-3; // 1 00633 } 00634 // if(G4UniformRand() > 0.5) Anti_DQ-=2; 00635 00636 if(std::abs(LeftQ1) > std::abs(LeftQ2)) 00637 { 00638 DQ=1000*LeftQ1+100*LeftQ2+3; // 1 00639 } else 00640 { 00641 DQ=1000*LeftQ2+100*LeftQ1+3; // 1 00642 } 00643 // if(G4UniformRand() > 0.5) DQ+=2; 00644 00645 // --------------- Set the string properties --------------- 00646 //G4cout<<"Left ADiQ DiQ "<<Anti_DQ<<" "<<DQ<<G4endl; 00647 00648 projectile->SplitUp(); 00649 00650 // projectile->SetFirstParton(Anti_DQ); 00651 // projectile->SetSecondParton(DQ); 00652 projectile->SetFirstParton(DQ); 00653 projectile->SetSecondParton(Anti_DQ); 00654 00655 projectile->SetStatus(1); 00656 target->SetStatus(3); // The target nucleon has annihilated 00657 00658 Pprojectile.setPx(0.); // VU Mar1 00659 Pprojectile.setPy(0.); // VU Mar1 00660 Pprojectile.setPz(0.); 00661 Pprojectile.setE(SqrtS); 00662 Pprojectile.transform(toLab); 00663 00664 // Calculation of the creation time --------------------- 00665 projectile->SetTimeOfCreation(target->GetTimeOfCreation()); 00666 projectile->SetPosition(target->GetPosition()); 00667 // Creation time and position of target nucleon were determined at 00668 // ReggeonCascade() of G4FTFModel 00669 // ------------------------------------------------------ 00670 00671 //G4cout<<"Mproj "<<Pprojectile.mag()<<G4endl; 00672 //G4cout<<"Mtarg "<<Ptarget.mag()<<G4endl; 00673 projectile->Set4Momentum(Pprojectile); 00674 00675 projectile->IncrementCollisionCount(1); 00676 00677 return true; 00678 } // end of if(CandidatsN != 0) 00679 } // if(Ksi < (X_a+X_b)/Xannihilation) 00680 00681 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00682 00683 if(Ksi < (X_a+X_b+X_c)/Xannihilation) 00684 { 00685 //============================================================ 00686 // Simulation of 2 anti-quark-quark strings creation 00687 //G4cout<<"Process c"<<G4endl; 00688 G4int CandidatsN(0), CandAQ[9][2], CandQ[9][2]; 00689 G4int LeftAQ1(0), LeftAQ2(0), LeftQ1(0), LeftQ2(0); 00690 //------------------------------------------------------------ 00691 for(G4int iAQ=0; iAQ<3; iAQ++) 00692 { 00693 for(G4int iQ=0; iQ<3; iQ++) 00694 { 00695 if(-AQ[iAQ] == Q[iQ]) 00696 { 00697 if(iAQ == 0) {CandAQ[CandidatsN][0]=1; CandAQ[CandidatsN][1]=2;} 00698 if(iAQ == 1) {CandAQ[CandidatsN][0]=0; CandAQ[CandidatsN][1]=2;} 00699 if(iAQ == 2) {CandAQ[CandidatsN][0]=0; CandAQ[CandidatsN][1]=1;} 00700 if(iQ == 0) {CandQ[CandidatsN][0] =1; CandQ[CandidatsN][1]=2;} 00701 if(iQ == 1) {CandQ[CandidatsN][0] =0; CandQ[CandidatsN][1]=2;} 00702 if(iQ == 2) {CandQ[CandidatsN][0] =0; CandQ[CandidatsN][1]=1;} 00703 CandidatsN++; 00704 } //end of if(-AQ[i] == Q[j]) 00705 } //end of cycle on targ. quarks 00706 } //end of cycle on proj. anti-quarks 00707 //------------------------------------------------------------ 00708 //G4cout<<"CandidatsN "<<CandidatsN<<G4endl; 00709 00710 if(CandidatsN != 0) 00711 { 00712 G4int SampledCase=CLHEP::RandFlat::shootInt(G4long(CandidatsN)); 00713 00714 LeftAQ1=AQ[CandAQ[SampledCase][0]]; 00715 LeftAQ2=AQ[CandAQ[SampledCase][1]]; 00716 00717 if(G4UniformRand() < 0.5) 00718 { 00719 LeftQ1=Q[CandQ[SampledCase][0]]; 00720 LeftQ2=Q[CandQ[SampledCase][1]]; 00721 } else 00722 { 00723 LeftQ2=Q[CandQ[SampledCase][0]]; 00724 LeftQ1=Q[CandQ[SampledCase][1]]; 00725 } 00726 00727 // --------------- Set the string properties --------------- 00728 //G4cout<<"String 1 "<<LeftAQ1<<" "<<LeftQ1<<G4endl; 00729 projectile->SplitUp(); 00730 00731 projectile->SetFirstParton(LeftAQ1); 00732 projectile->SetSecondParton(LeftQ1); 00733 projectile->SetStatus(1); 00734 00735 //G4cout<<"String 2 "<<LeftAQ2<<" "<<LeftQ2<<G4endl; 00736 target->SplitUp(); 00737 00738 target->SetFirstParton(LeftQ2); 00739 target->SetSecondParton(LeftAQ2); 00740 target->SetStatus(1); 00741 00742 // Sampling kinematical properties 00743 // 1 string LeftAQ1-LeftQ1// 2 string LeftAQ2-LeftQ2 00744 00745 G4ThreeVector Quark_Mom[4]; 00746 G4double ModMom2[4]; //ModMom[4], 00747 00748 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00749 AveragePt2=200.*200.; maxPtSquare=S; 00750 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00751 00752 G4double SumMt(0.); 00753 00754 G4double MassQ2=0.; //100.*100.*MeV*MeV; 00755 00756 G4int NumberOfTries(0); 00757 G4double ScaleFactor(1.); 00758 do // while(SumMt >SqrtS) 00759 { 00760 NumberOfTries++; 00761 00762 if(NumberOfTries == 100*(NumberOfTries/100)) 00763 { // At large number of tries it would be better to reduce the values of <Pt^2> 00764 ScaleFactor/=2.; 00765 AveragePt2 *=ScaleFactor; 00766 } 00767 00768 G4ThreeVector PtSum(0.,0.,0.); 00769 for(G4int i=0; i<4; i++) 00770 { 00771 Quark_Mom[i]=GaussianPt(AveragePt2, maxPtSquare); 00772 PtSum+=Quark_Mom[i]; 00773 } 00774 00775 PtSum/=4.; 00776 00777 SumMt=0.; 00778 for(G4int i=0; i<4; i++) 00779 { 00780 Quark_Mom[i]-=PtSum; 00781 // ModMom[i] =Quark_Mom[i].mag(); 00782 ModMom2[i]=Quark_Mom[i].mag2(); 00783 SumMt+=std::sqrt(ModMom2[i]+MassQ2); 00784 } 00785 } while(SumMt > SqrtS); 00786 00787 G4double WminusTarget(0.), WplusProjectile(0.); 00788 00789 // Sampling X's of anti-baryon ------- 00790 G4double Alfa_R=0.5; 00791 00792 NumberOfTries=0; 00793 ScaleFactor=1.; 00794 00795 G4bool Succes(true); 00796 do // while(!Succes) 00797 { 00798 Succes=true; 00799 NumberOfTries++; 00800 00801 if(NumberOfTries == 100*(NumberOfTries/100)) 00802 { // At large number of tries it would be better to reduce the values of Pt's 00803 ScaleFactor/=2.; 00804 } 00805 00806 if(Alfa_R == 1.) 00807 { 00808 G4double Xaq1=std::sqrt(G4UniformRand()); 00809 G4double Xaq2=1.-Xaq1; 00810 00811 Quark_Mom[0].setZ(Xaq1); Quark_Mom[1].setZ(Xaq2); 00812 } 00813 else 00814 { 00815 G4double Xaq1=sqr(std::sin(pi/2.*G4UniformRand())); 00816 G4double Xaq2=1.-Xaq1; 00817 00818 Quark_Mom[0].setZ(Xaq1); Quark_Mom[1].setZ(Xaq2); 00819 } // end of if(Alfa_R == 0.) 00820 00821 // Sampling X's of baryon ------------ 00822 if(Alfa_R == 1.) 00823 { 00824 G4double Xq1=1.-std::sqrt(G4UniformRand()); 00825 G4double Xq2=1.-Xq1; 00826 00827 Quark_Mom[2].setZ(Xq1); Quark_Mom[3].setZ(Xq2); 00828 } 00829 else 00830 { 00831 G4double Xq1=sqr(std::sin(pi/2.*G4UniformRand())); 00832 G4double Xq2=1.-Xq1; 00833 00834 Quark_Mom[2].setZ(Xq1); Quark_Mom[3].setZ(Xq2); 00835 } // end of if(Alfa_R == 0.) 00836 // 00837 G4double Alfa(0.), Beta(0.); 00838 00839 for(G4int i=0; i<2; i++) // For Anti-baryon 00840 { 00841 if(Quark_Mom[i].getZ() != 0.) 00842 {Alfa+=(ScaleFactor*ModMom2[i]+MassQ2)/Quark_Mom[i].getZ();} 00843 else {Succes=false;} 00844 } 00845 00846 for(G4int i=2; i<4; i++) // For baryon 00847 { 00848 if(Quark_Mom[i].getZ() != 0.) 00849 {Beta+=(ScaleFactor*ModMom2[i]+MassQ2)/Quark_Mom[i].getZ();} 00850 else {Succes=false;} 00851 } 00852 00853 if(!Succes) continue; 00854 00855 if(std::sqrt(Alfa)+std::sqrt(Beta) > SqrtS) {Succes=false; continue;} 00856 00857 G4double DecayMomentum2=S*S + Alfa*Alfa + Beta*Beta 00858 - 2.*S*Alfa - 2.*S*Beta - 2.*Alfa*Beta; 00859 00860 WminusTarget=(S-Alfa+Beta+std::sqrt(DecayMomentum2))/2./SqrtS; 00861 WplusProjectile=SqrtS-Beta/WminusTarget; 00862 00863 } while(!Succes); 00864 // //-------------------------------------------------- 00865 00866 G4double SqrtScaleF=std::sqrt(ScaleFactor); 00867 00868 for(G4int i=0; i<2; i++) 00869 { 00870 G4double Pz=WplusProjectile*Quark_Mom[i].getZ()/2.- 00871 (ScaleFactor*ModMom2[i]+MassQ2)/(2.*WplusProjectile*Quark_Mom[i].getZ()); 00872 Quark_Mom[i].setZ(Pz); 00873 00874 if(ScaleFactor != 1.) 00875 { 00876 Quark_Mom[i].setX(SqrtScaleF*Quark_Mom[i].getX()); 00877 Quark_Mom[i].setY(SqrtScaleF*Quark_Mom[i].getY()); 00878 } 00879 //G4cout<<"Anti Q "<<i<<" "<<Quark_Mom[i]<<G4endl; 00880 } 00881 00882 for(G4int i=2; i<4; i++) 00883 { 00884 G4double Pz=-WminusTarget*Quark_Mom[i].getZ()/2.+ 00885 (ScaleFactor*ModMom2[i]+MassQ2)/(2.*WminusTarget*Quark_Mom[i].getZ()); 00886 Quark_Mom[i].setZ(Pz); 00887 00888 if(ScaleFactor != 1.) 00889 { 00890 Quark_Mom[i].setX(SqrtScaleF*Quark_Mom[i].getX()); 00891 Quark_Mom[i].setY(SqrtScaleF*Quark_Mom[i].getY()); 00892 } 00893 //G4cout<<"Bary Q "<<i<<" "<<Quark_Mom[i]<<G4endl; 00894 } 00895 //G4cout<<"Sum AQ "<<Quark_Mom[0]+Quark_Mom[1]<<G4endl; 00896 //G4cout<<"Sum Q "<<Quark_Mom[2]+Quark_Mom[3]<<G4endl; 00897 //------------------------------------- 00898 00899 G4ThreeVector tmp=Quark_Mom[0]+Quark_Mom[2]; 00900 G4LorentzVector Pstring1(tmp,std::sqrt(Quark_Mom[0].mag2()+MassQ2)+ 00901 std::sqrt(Quark_Mom[2].mag2()+MassQ2)); 00902 G4double Ystring1=Pstring1.rapidity(); 00903 /* 00904 G4cout<<"Mom 1 string "<<G4endl; 00905 G4cout<<Quark_Mom[0]<<G4endl; 00906 G4cout<<Quark_Mom[2]<<G4endl; 00907 G4cout<<tmp<<" "<<tmp.mag()<<G4endl; 00908 //G4cout<<"1 str "<<Pstring1<<" "<<Pstring1.mag()<<" "<<Ystring1<<G4endl; 00909 */ 00910 00911 tmp=Quark_Mom[1]+Quark_Mom[3]; 00912 G4LorentzVector Pstring2(tmp,std::sqrt(Quark_Mom[1].mag2()+MassQ2)+ 00913 std::sqrt(Quark_Mom[3].mag2()+MassQ2)); 00914 G4double Ystring2=Pstring2.rapidity(); 00915 /* 00916 G4cout<<"Mom 2 string "<<G4endl; 00917 G4cout<<Quark_Mom[1]<<G4endl; 00918 G4cout<<Quark_Mom[3]<<G4endl; 00919 G4cout<<tmp<<" "<<tmp.mag()<<G4endl; 00920 G4cout<<"2 str "<<Pstring2<<" "<<Pstring2.mag()<<" "<<Ystring2<<G4endl; 00921 */ 00922 //-------------------------------- 00923 if(Ystring1 > Ystring2) 00924 { 00925 Pprojectile=Pstring1; 00926 Ptarget =Pstring2; 00927 } else 00928 { 00929 Pprojectile=Pstring2; 00930 Ptarget =Pstring1; 00931 } 00932 00933 //------------------------------------------------------- 00934 //G4cout<<"SumP CMS "<<Pprojectile+Ptarget<<" "<<SqrtS<<G4endl; 00935 00936 Pprojectile.transform(toLab); 00937 Ptarget.transform(toLab); 00938 //G4cout<<"SumP Lab "<<Pprojectile+Ptarget<<" "<<SqrtS<<G4endl; 00939 00940 // Calculation of the creation time --------------------- 00941 projectile->SetTimeOfCreation(target->GetTimeOfCreation()); 00942 projectile->SetPosition(target->GetPosition()); 00943 00944 // Creation time and position of target nucleon were determined at 00945 // ReggeonCascade() of G4FTFModel 00946 // ------------------------------------------------------ 00947 00948 //G4cout<<"Mproj "<<Pprojectile.mag()<<G4endl; 00949 //G4cout<<"Mtarg "<<Ptarget.mag()<<G4endl; 00950 projectile->Set4Momentum(Pprojectile); 00951 00952 target->Set4Momentum(Ptarget); 00953 00954 projectile->IncrementCollisionCount(1); 00955 target->IncrementCollisionCount(1); 00956 00957 return true; 00958 } // End of if(CandidatsN != 0) 00959 } 00960 00961 //============================================================ 00962 // Simulation of anti-quark-quark string creation 00963 // 00964 if(Ksi < (X_a+X_b+X_c+X_d)/Xannihilation) 00965 { 00966 //G4cout<<"Process d"<<G4endl; 00967 G4int CandidatsN(0), CandAQ[9], CandQ[9]; 00968 G4int LeftAQ(0), LeftQ(0); 00969 //------------------------------------------------------------ 00970 for(G4int iAQ1=0; iAQ1<3; iAQ1++) 00971 { 00972 for(G4int iAQ2=0; iAQ2<3; iAQ2++) 00973 { 00974 if(iAQ1 != iAQ2) 00975 { 00976 for(G4int iQ1=0; iQ1<3; iQ1++) 00977 { 00978 for(G4int iQ2=0; iQ2<3; iQ2++) 00979 { 00980 if(iQ1 != iQ2) 00981 { 00982 if((-AQ[iAQ1] == Q[iQ1]) && (-AQ[iAQ2] == Q[iQ2])) 00983 { 00984 if((iAQ1 == 0) && (iAQ2 == 1)){CandAQ[CandidatsN]=2;} 00985 if((iAQ1 == 1) && (iAQ2 == 0)){CandAQ[CandidatsN]=2;} 00986 00987 if((iAQ1 == 0) && (iAQ2 == 2)){CandAQ[CandidatsN]=1;} 00988 if((iAQ1 == 2) && (iAQ2 == 0)){CandAQ[CandidatsN]=1;} 00989 00990 if((iAQ1 == 1) && (iAQ2 == 2)){CandAQ[CandidatsN]=0;} 00991 if((iAQ1 == 2) && (iAQ2 == 1)){CandAQ[CandidatsN]=0;} 00992 //---------------------------------------------------------------- 00993 if((iQ1 == 0) && (iQ2 == 1)){CandQ[CandidatsN]=2;} 00994 if((iQ1 == 1) && (iQ2 == 0)){CandQ[CandidatsN]=2;} 00995 00996 if((iQ1 == 0) && (iQ2 == 2)){CandQ[CandidatsN]=1;} 00997 if((iQ1 == 2) && (iQ2 == 0)){CandQ[CandidatsN]=1;} 00998 00999 if((iQ1 == 1) && (iQ2 == 2)){CandQ[CandidatsN]=0;} 01000 if((iQ1 == 2) && (iQ2 == 1)){CandQ[CandidatsN]=0;} 01001 CandidatsN++; 01002 }//-------------------------- 01003 } //end of if(jQ1 != jQ2) 01004 } //end of for(G4int jQ2=0; j<3; j++) 01005 } //end of for(G4int jQ=0; j<3; j++) 01006 } //end of if(iAQ1 != iAQ2) 01007 } //end of for(G4int iAQ2=0; i<3; i++) 01008 } //end of for(G4int iAQ1=0; i<3; i++) 01009 //------------------------------------------------------------ 01010 01011 if(CandidatsN != 0) 01012 { 01013 G4int SampledCase=CLHEP::RandFlat::shootInt(G4long(CandidatsN)); 01014 01015 LeftAQ=AQ[CandAQ[SampledCase]]; 01016 01017 LeftQ =Q[CandQ[SampledCase]]; 01018 01019 // --------------- Set the string properties --------------- 01020 //G4cout<<"Left Aq Q "<<LeftAQ<<" "<<LeftQ<<G4endl; 01021 01022 projectile->SplitUp(); 01023 01024 // projectile->SetFirstParton(LeftAQ); 01025 // projectile->SetSecondParton(LeftQ); 01026 projectile->SetFirstParton(LeftQ); 01027 projectile->SetSecondParton(LeftAQ); 01028 01029 projectile->SetStatus(1); 01030 target->SetStatus(3); // The target nucleon has annihilated 01031 01032 Pprojectile.setPx(0.); // VU Mar1 01033 Pprojectile.setPy(0.); // Vu Mar1 01034 Pprojectile.setPz(0.); 01035 Pprojectile.setE(SqrtS); 01036 Pprojectile.transform(toLab); 01037 01038 // Calculation of the creation time --------------------- 01039 projectile->SetTimeOfCreation(target->GetTimeOfCreation()); 01040 projectile->SetPosition(target->GetPosition()); 01041 // Creation time and position of target nucleon were determined at 01042 // ReggeonCascade() of G4FTFModel 01043 // ------------------------------------------------------ 01044 01045 //G4cout<<"Mproj "<<Pprojectile.mag()<<G4endl; 01046 //G4cout<<"Mtarg "<<Ptarget.mag()<<G4endl; 01047 projectile->Set4Momentum(Pprojectile); 01048 01049 projectile->IncrementCollisionCount(1); 01050 return true; 01051 } // end of if(CandidatsN != 0) 01052 } // if(Ksi < (X_a+X_b+X_c+X_d/Xannihilation) 01053 01054 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 01055 //G4cout<<"Pr Y "<<Pprojectile.rapidity()<<" Tr Y "<<Ptarget.rapidity()<<G4endl; 01056 return true; 01057 }
1.4.7