G4QGSDiffractiveExcitation Class Reference

#include <G4QGSDiffractiveExcitation.hh>

Inheritance diagram for G4QGSDiffractiveExcitation:

Public Member Functions
	G4QGSDiffractiveExcitation ()
virtual	~G4QGSDiffractiveExcitation ()
virtual G4bool	ExciteParticipants (G4VSplitableHadron *aPartner, G4VSplitableHadron *bPartner, G4bool ProjectileDiffraction=TRUE) const
virtual G4ExcitedString *	String (G4VSplitableHadron *aHadron, G4bool isProjectile) const

Detailed Description

Definition at line 50 of file G4QGSDiffractiveExcitation.hh.

Constructor & Destructor Documentation

G4QGSDiffractiveExcitation()

G4QGSDiffractiveExcitation::G4QGSDiffractiveExcitation

(

)

Definition at line 68 of file G4QGSDiffractiveExcitation.cc.

{
}

~G4QGSDiffractiveExcitation()

G4QGSDiffractiveExcitation::~G4QGSDiffractiveExcitation ( )

virtual

Definition at line 72 of file G4QGSDiffractiveExcitation.cc.

{
}

Member Function Documentation

ExciteParticipants()

G4bool G4QGSDiffractiveExcitation::ExciteParticipants ( G4VSplitableHadron * aPartner, G4VSplitableHadron * bPartner, G4bool ProjectileDiffraction = TRUE ) const

virtual

Reimplemented in G4SingleDiffractiveExcitation.

Definition at line 77 of file G4QGSDiffractiveExcitation.cc.

{
  #ifdef debugDoubleDiffraction
    G4cout<<G4endl<<"G4QGSDiffractiveExcitation::ExciteParticipants - Double diffraction."<<G4endl;
    G4cout<<"Proj Targ "<<projectile->GetDefinition()->GetParticleName()<<" "<<target->GetDefinition()->GetParticleName()<<G4endl;
    G4cout<<"Proj 4 Mom "<<projectile->Get4Momentum()<<" "<<projectile->Get4Momentum().mag()<<G4endl;
    G4cout<<"Targ 4 Mom "<<target->Get4Momentum()    <<" "<<target->Get4Momentum().mag()    <<G4endl;
  #endif
 
  G4LorentzVector Pprojectile=projectile->Get4Momentum();
 
  // -------------------- Projectile parameters -----------------------------------
  G4bool PutOnMassShell=0;
 
  G4double M0projectile = Pprojectile.mag();                        // Without de-excitation
 
  if(M0projectile < projectile->GetDefinition()->GetPDGMass())
  {
    PutOnMassShell=1;
    M0projectile=projectile->GetDefinition()->GetPDGMass();
  }
 
  // -------------------- Target parameters ----------------------------------------------
  G4LorentzVector Ptarget=target->Get4Momentum();
 
  G4double M0target = Ptarget.mag();
 
  if(M0target < target->GetDefinition()->GetPDGMass())
  {
    PutOnMassShell=1;
    M0target=target->GetDefinition()->GetPDGMass();
  }
 
  G4LorentzVector Psum=Pprojectile+Ptarget;
  G4double S=Psum.mag2();
  G4double SqrtS=std::sqrt(S);
 
  if (SqrtS < M0projectile + M0target) {return false;}  // The model cannot work for pp-interactions
                                                    // at Plab < 1.3 GeV/c.
 
  G4double Mprojectile2 = M0projectile * M0projectile;
  G4double Mtarget2     = M0target     * M0target;    //Ptarget.mag2(); // for AA-inter.
 
  // Transform momenta to cms and then rotate parallel to z axis;
 
  G4LorentzRotation toCms(-1*Psum.boostVector());
 
  G4LorentzVector Ptmp=toCms*Pprojectile;
 
  if ( Ptmp.pz() <= 0. )
  {
    // "String" moving backwards in  CMS, abort collision !!
    //G4cout << " abort Collision!! " << G4endl;
    return false;
  }
 
  toCms.rotateZ(-1*Ptmp.phi());
  toCms.rotateY(-1*Ptmp.theta());
 
  G4LorentzRotation toLab(toCms.inverse());
 
  Pprojectile.transform(toCms);
  Ptarget.transform(toCms);
 
  G4double PZcms2=(S*S+Mprojectile2*Mprojectile2+Mtarget2*Mtarget2-
           2*S*Mprojectile2-2*S*Mtarget2-2*Mprojectile2*Mtarget2)/4./S;
 
  if (PZcms2 < 0) {return false;}   // It can be in an interaction with off-shell nuclear nucleon
 
  G4double PZcms = std::sqrt(PZcms2);
 
  if (PutOnMassShell)
  {
    if (Pprojectile.z() > 0.)
    {
      Pprojectile.setPz( PZcms);
      Ptarget.setPz(    -PZcms);
    }
    else
    {
      Pprojectile.setPz(-PZcms);
      Ptarget.setPz(     PZcms);
    };
 
    Pprojectile.setE(std::sqrt(Mprojectile2+sqr(Pprojectile.x())+sqr(Pprojectile.y())+PZcms2));
    Ptarget.setE(    std::sqrt(    Mtarget2+sqr(    Ptarget.x())+sqr(    Ptarget.y())+PZcms2));
  }
 
  G4double maxPtSquare = PZcms2;
 
  #ifdef debugDoubleDiffraction
    G4cout << "Pprojectile after boost to CMS: " << Pprojectile <<" "<<Pprojectile.mag()<<G4endl;
    G4cout << "Ptarget     after boost to CMS: " << Ptarget     <<" "<<Ptarget.mag()    <<G4endl;
  #endif
 
  G4int    PrPDGcode=projectile->GetDefinition()->GetPDGEncoding();
  G4int    absPrPDGcode=std::abs(PrPDGcode);
  G4double MinPrDiffMass(0.);
  G4double AveragePt2(0.);
 
  if (M0projectile <= projectile->GetDefinition()->GetPDGMass())
  { // Normal projectile 
    if( absPrPDGcode > 1000 )                         //------Projectile is baryon --------
    {
      if ( absPrPDGcode > 4000 && absPrPDGcode < 6000 )  // Projectile is a charm or bottom baryon
      {
        MinPrDiffMass = projectile->GetDefinition()->GetPDGMass()/CLHEP::GeV + 0.25;  // GeV
        AveragePt2 = 0.3;                                                             // GeV^2
      }
      else
      {
        MinPrDiffMass = 1.16;              // GeV
        AveragePt2 = 0.3;                           // GeV^2
      }
    }
    else if( absPrPDGcode == 211 || PrPDGcode ==  111) //------Projectile is Pion -----------
    {
      MinPrDiffMass = 1.0;                       // GeV
      AveragePt2 = 0.3;                          // GeV^2
    }
    else if( absPrPDGcode == 321 || absPrPDGcode == 130 || absPrPDGcode == 310) //-Projectile is Kaon-
    {
      MinPrDiffMass = 1.1;                        // GeV
      AveragePt2 = 0.3;                           // GeV^2
    }
    else if( absPrPDGcode > 400 && absPrPDGcode < 600)  // Projectile is a charm or bottom meson
    {
      MinPrDiffMass = projectile->GetDefinition()->GetPDGMass()/CLHEP::GeV + 0.25;  // GeV
      AveragePt2 = 0.3;                                                             // GeV^2
    }
    else                                              //------Projectile is undefined, Nucleon assumed
    {
      MinPrDiffMass = 1.16;                       // GeV
      AveragePt2 = 0.3;                           // GeV^2
    }
  }
  else
  { // Excited projectile
    MinPrDiffMass = M0projectile + 220.0*MeV;
    AveragePt2 = 0.3; 
  }
 
  MinPrDiffMass = MinPrDiffMass * GeV;
  AveragePt2 = AveragePt2 * GeV*GeV;
  //---------------------------------------------
  G4double MinTrDiffMass = 1.16*GeV;
 
  if (SqrtS < MinPrDiffMass + MinTrDiffMass) {return false;}   // The model cannot work at low energy
 
  G4double MinPrDiffMass2 = MinPrDiffMass * MinPrDiffMass;
  G4double MinTrDiffMass2 = MinTrDiffMass * MinTrDiffMass;
 
  G4double Pt2;
  G4double ProjMassT2, ProjMassT;
  G4double TargMassT2, TargMassT;
  G4double PMinusNew, TPlusNew;
 
  G4LorentzVector Qmomentum;
  G4double Qminus, Qplus;
 
  G4int whilecount=0;
  do {
    if (whilecount++ >= 500 && (whilecount%100)==0)
      if (whilecount > 1000 ) {
    Qmomentum=G4LorentzVector(0.,0.,0.,0.);
    return false;     //  Ignore this interaction
      }
 
    //  Generate pt
    Qmomentum=G4LorentzVector(GaussianPt(AveragePt2,maxPtSquare),0);
 
    Pt2=G4ThreeVector(Qmomentum.vect()).mag2();
    ProjMassT2=MinPrDiffMass2+Pt2;
    ProjMassT =std::sqrt(ProjMassT2);
 
    TargMassT2=MinTrDiffMass2+Pt2;
    TargMassT =std::sqrt(TargMassT2);
 
    if (SqrtS < ProjMassT + TargMassT) continue;
 
    PZcms2=(S*S+ProjMassT2*ProjMassT2+TargMassT2*TargMassT2-
        2.*S*ProjMassT2-2.*S*TargMassT2-2.*ProjMassT2*TargMassT2)/4./S;
    if (PZcms2 < 0 ) {PZcms2=0;};
    PZcms =std::sqrt(PZcms2);
 
    G4double PMinusMin=std::sqrt(ProjMassT2+PZcms2)-PZcms;
    G4double PMinusMax=SqrtS-TargMassT;
 
    PMinusNew=ChooseP(PMinusMin,PMinusMax);
    Qminus=PMinusNew-Pprojectile.minus();
 
    G4double TPlusMin=std::sqrt(TargMassT2+PZcms2)-PZcms;
    G4double TPlusMax=SqrtS-ProjMassT;
 
    TPlusNew=ChooseP(TPlusMin, TPlusMax);
    Qplus=-(TPlusNew-Ptarget.plus());
 
    Qmomentum.setPz( (Qplus-Qminus)/2 );
    Qmomentum.setE(  (Qplus+Qminus)/2 );
 
  } while ( (Pprojectile+Qmomentum).mag2() <  MinPrDiffMass2 ||
        (Ptarget    -Qmomentum).mag2() <  MinTrDiffMass2   ); 
 
  Pprojectile += Qmomentum;
  Ptarget     -= Qmomentum;
 
  // Transform back and update SplitableHadron Participant.
  Pprojectile.transform(toLab);
  Ptarget.transform(toLab);
 
  #ifdef debugDoubleDiffraction
    G4cout << "Pprojectile after boost to Lab: " << Pprojectile <<" "<<Pprojectile.mag()<<G4endl;
    G4cout << "Ptarget     after boost to Lab: " << Ptarget     <<" "<<Ptarget.mag()    <<G4endl;
  #endif
 
  target->Set4Momentum(Ptarget);
  projectile->Set4Momentum(Pprojectile);
 
  return true;
}

Referenced by G4QGSParticipants::PerformDiffractiveCollisions().

String()

G4ExcitedString * G4QGSDiffractiveExcitation::String ( G4VSplitableHadron * aHadron, G4bool isProjectile ) const

virtual

Definition at line 300 of file G4QGSDiffractiveExcitation.cc.

{
  hadron->SplitUp();
  G4Parton *start= hadron->GetNextParton();
  if ( start==NULL) 
  { G4cout << " G4QGSDiffractiveExcitation::String() Error:No start parton found"<< G4endl;
    return NULL;
  }
  G4Parton *end  = hadron->GetNextParton();
  if ( end==NULL) 
  { G4cout << " G4QGSDiffractiveExcitation::String() Error:No end parton found"<< G4endl;
    return NULL;
  }
 
  G4ExcitedString * string;
  if ( isProjectile ) 
  {
    string= new G4ExcitedString(end,start, +1);
  } else {
    string= new G4ExcitedString(start,end, -1);
  }
 
  string->SetPosition(hadron->GetPosition());
 
  // momenta of string ends
 
  G4double maxAvailMomentumSquared=sqr(hadron->Get4Momentum().mag()/2.);
 
  G4double widthOfPtSquare = 0.5*sqr(GeV);
  G4ThreeVector pt=GaussianPt(widthOfPtSquare,maxAvailMomentumSquared);
 
  G4LorentzVector Pstart(G4LorentzVector(pt,0.));
  G4LorentzVector Pend;
  Pend.setPx(hadron->Get4Momentum().px() - pt.x());
  Pend.setPy(hadron->Get4Momentum().py() - pt.y());
 
  G4double tm1=hadron->Get4Momentum().minus() +
           ( Pend.perp2()-Pstart.perp2() ) / hadron->Get4Momentum().plus();
 
  G4double tm2= std::sqrt( std::max(0., sqr(tm1) -
               4. * Pend.perp2() * hadron->Get4Momentum().minus()
               /  hadron->Get4Momentum().plus() ));
 
  G4int Sign= isProjectile ? -1 : 1;
 
  G4double endMinus  = 0.5 * (tm1 + Sign*tm2);
  G4double startMinus= hadron->Get4Momentum().minus() - endMinus;
 
  G4double startPlus= Pstart.perp2() /  startMinus;
  G4double endPlus  = hadron->Get4Momentum().plus() - startPlus;
 
  Pstart.setPz(0.5*(startPlus - startMinus));
  Pstart.setE(0.5*(startPlus + startMinus));
 
  Pend.setPz(0.5*(endPlus - endMinus));
  Pend.setE(0.5*(endPlus + endMinus));
 
  start->Set4Momentum(Pstart);
  end->Set4Momentum(Pend);
 
  #ifdef debugQGSdiffExictation
    G4cout << " generated string flavors          " << start->GetPDGcode()   << " / " << end->GetPDGcode() << G4endl;
    G4cout << " generated string momenta:   quark " << start->Get4Momentum() << "mass : " <<start->Get4Momentum().mag()<< G4endl;
    G4cout << " generated string momenta: Diquark " << end ->Get4Momentum()  << "mass : " <<end->Get4Momentum().mag()<< G4endl;
    G4cout << " sum of ends                       " << Pstart+Pend << G4endl;
    G4cout << " Original                          " << hadron->Get4Momentum() << G4endl;
  #endif
 
  return string;    
}

---

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

• G4QGSDiffractiveExcitation.hh
• G4QGSDiffractiveExcitation.cc