G4Quasmon Class Reference

#include <G4Quasmon.hh>


Public Member Functions

 G4Quasmon (G4QContent qQCont=G4QContent(0, 0, 0, 0, 0, 0), G4LorentzVector q4M=G4LorentzVector(0., 0., 0., 0.), G4LorentzVector ph4M=G4LorentzVector(0., 0., 0., 0.))
 G4Quasmon (const G4Quasmon &right)
 G4Quasmon (G4Quasmon *right)
 ~G4Quasmon ()
const G4Quasmonoperator= (const G4Quasmon &right)
G4bool operator== (const G4Quasmon &right) const
G4bool operator!= (const G4Quasmon &right) const
G4double GetTemper () const
G4double GetSOverU () const
G4double GetEtaSup () const
G4LorentzVector Get4Momentum () const
G4QContent GetQC () const
G4QPDGCode GetQPDG () const
G4int GetStatus () const
G4int GetCharge () const
G4int GetBaryonNumber () const
G4int GetStrangeness () const
void Set4Momentum (G4LorentzVector Q4M)
void SetQC (G4QContent QQC)
void Boost (const G4LorentzVector &theBoost)
void Boost (const G4ThreeVector &B)
G4QHadronVectorFragment (G4QNucleus &nucEnviron, G4int nQ=1)
G4QHadronVectorDecayQuasmon ()
G4QHadronVectorDecayQHadron (G4QHadron *hadron)
void ClearOutput ()
void InitQuasmon (const G4QContent &qQCont, const G4LorentzVector &q4M)
void IncreaseBy (const G4Quasmon *pQuasm)
void IncreaseBy (G4QContent &qQCont, const G4LorentzVector &q4M)
void ClearQuasmon ()
void KillQuasmon ()
G4int CalculateNumberOfQPartons (G4double qMass)

Static Public Member Functions

static void SetParameters (G4double temper=180., G4double ssin2g=.3, G4double etaetap=.3)
static void SetTemper (G4double temperature)
static void SetSOverU (G4double ssin2g)
static void SetEtaSup (G4double etaetap)
static void OpenElectromagneticDecays ()
static void CloseElectromagneticDecays ()


Detailed Description

Definition at line 59 of file G4Quasmon.hh.


Constructor & Destructor Documentation

G4Quasmon::G4Quasmon ( G4QContent  qQCont = G4QContent(0, 0, 0, 0, 0, 0),
G4LorentzVector  q4M = G4LorentzVector(0., 0., 0., 0.),
G4LorentzVector  ph4M = G4LorentzVector(0., 0., 0., 0.) 
)

Definition at line 74 of file G4Quasmon.cc.

References G4QContent::DecQAQ(), G4cout, and G4endl.

00075  : q4Mom(q4M), valQ(qQCont), theWorld(0), phot4M(ph4M), f2all(0), rEP(0.), rMo(0.)
00076 {
00077 #ifdef debug
00078   G4cout<<"G4Quasmon:Constructor:QC="<<qQCont<<",Q4M="<<q4M<<",photonE="<<ph4M.e()<<G4endl;
00079 #endif
00080 #ifdef pardeb
00081   G4cout<<"**>G4Q:Con:(1),T="<<Temperature<<",S="<<SSin2Gluons<<",E="<<EtaEtaprime<<G4endl;
00082 #endif
00083   if(phot4M.e()>0.) q4Mom+=phot4M; // InCaseOf CaptureByQuark it will be subtracted back
00084   valQ.DecQAQ(-1);
00085   status=4;                                   
00086 }

G4Quasmon::G4Quasmon ( const G4Quasmon right  ) 

Definition at line 88 of file G4Quasmon.cc.

References f2all, nBarClust, nOfQ, phot4M, q4Mom, rEP, rMo, status, theQCandidates, theQHadrons, theWorld, and valQ.

00088                                           :  totMass(0), bEn(0), mbEn(0)
00089 {
00090   q4Mom                 = right.q4Mom;
00091   valQ                  = right.valQ;
00092   //theEnvironment        = right.theEnvironment;
00093   status                = right.status;
00094   //theQHadrons (Vector)
00095   G4int nQH             = right.theQHadrons.size();
00096   if(nQH) for(G4int ih=0; ih<nQH; ih++)
00097   {
00098     G4QHadron* curQH    = new G4QHadron(right.theQHadrons[ih]);
00099     theQHadrons.push_back(curQH);
00100   }
00101   theWorld              = right.theWorld;
00102   phot4M                = right.phot4M;
00103   nBarClust             = right.nBarClust;
00104   nOfQ                  = right.nOfQ;
00105   //theQCandidates (Vector)
00106   G4int nQC             = right.theQCandidates.size();
00107   if(nQC) for(G4int iq=0; iq<nQC; iq++)
00108   {
00109     G4QCandidate* curQC = new G4QCandidate(right.theQCandidates[iq]);
00110     theQCandidates.push_back(curQC);
00111   }
00112   f2all                 = right.f2all;
00113   rEP                   = right.rEP;
00114   rMo                   = right.rMo;
00115 }

G4Quasmon::G4Quasmon ( G4Quasmon right  ) 

Definition at line 117 of file G4Quasmon.cc.

References f2all, G4cout, G4endl, nBarClust, nOfQ, phot4M, q4Mom, rEP, rMo, status, theEnvironment, theQCandidates, theQHadrons, theWorld, and valQ.

00117                                     :  totMass(0), bEn(0), mbEn(0)
00118 {
00119 #ifdef sdebug
00120   G4cout<<"G4Quasmon::Copy-Constructor: ***CALLED*** E="<<right->theEnvironment<<G4endl;
00121 #endif
00122   q4Mom                 = right->q4Mom;
00123   valQ                  = right->valQ;
00124   //theEnvironment        = right->theEnvironment;
00125   status                = right->status;
00126   //theQHadrons (Vector)
00127   G4int nQH             = right->theQHadrons.size();
00128 #ifdef sdebug
00129   G4cout<<"G4Quasmon::Copy-Constructor:nQH="<<nQH<<G4endl;
00130 #endif
00131   if(nQH) for(G4int ih=0; ih<nQH; ih++)
00132   {
00133 #ifdef debug
00134     G4cout<<"G4Quasmon:Copy-Constructor:H#"<<ih<<",QH="<<right->theQHadrons[ih]<<G4endl;
00135 #endif
00136     G4QHadron* curQH    = new G4QHadron(right->theQHadrons[ih]);
00137     theQHadrons.push_back(curQH);
00138   }
00139   theWorld              = right->theWorld;
00140   phot4M                = right->phot4M;
00141   nBarClust             = right->nBarClust;
00142   nOfQ                  = right->nOfQ;
00143   //theQCandidates (Vector)
00144   G4int nQC             = right->theQCandidates.size();
00145 #ifdef sdebug
00146   G4cout<<"G4Quasmon:Copy-Constructor: nCand="<<nQC<<G4endl;
00147 #endif
00148   if(nQC) for(G4int iq=0; iq<nQC; iq++)
00149   {
00150 #ifdef sdebug
00151     G4cout<<"G4Quasmon:Copy-Constructor:C#"<<iq<<",QC="<<right->theQCandidates[iq]<<G4endl;
00152 #endif
00153     G4QCandidate* curQC = new G4QCandidate(right->theQCandidates[iq]);
00154     theQCandidates.push_back(curQC);
00155   }
00156   f2all                 = right->f2all;
00157   rEP                   = right->rEP;
00158   rMo                   = right->rMo;
00159 #ifdef sdebug
00160   G4cout<<"G4Quasmon:Copy-Constructor: >->-> DONE <-<-<"<<G4endl;
00161 #endif
00162 }

G4Quasmon::~G4Quasmon (  ) 

Definition at line 164 of file G4Quasmon.cc.

References G4cout, and G4endl.

00165 {
00166 #ifdef sdebug
00167   G4cout<<"G4Quasmon::Destructor before theQCandidates delete"<<G4endl;
00168 #endif
00169   for_each(theQCandidates.begin(), theQCandidates.end(), DeleteQCandidate());
00170 #ifdef sdebug
00171   G4cout<<"G4Quasmon::Destructor before theQHadrons"<<G4endl;
00172 #endif
00173   for_each(theQHadrons.begin(), theQHadrons.end(), DeleteQHadron());
00174 #ifdef sdebug
00175   G4cout<<"G4Quasmon::Destructor === DONE ==="<<G4endl;
00176 #endif
00177 }


Member Function Documentation

void G4Quasmon::Boost ( const G4ThreeVector B  )  [inline]

Definition at line 99 of file G4Quasmon.hh.

00099 {q4Mom.boost(B);} // Boosts 4-Momentum using v/c

void G4Quasmon::Boost ( const G4LorentzVector theBoost  ) 

Definition at line 6202 of file G4Quasmon.cc.

06203 {  
06204   // see CERNLIB short writeup U101 for the algorithm
06205   G4double bm=boost4M.mag();
06206   G4double factor = (q4Mom.vect()*boost4M.vect()/(boost4M.e()+bm) - q4Mom.e())/bm;
06207   q4Mom.setE(q4Mom.dot(boost4M)/bm);
06208   q4Mom.setVect(factor*boost4M.vect() + q4Mom.vect());
06209 } // End of Boost

G4int G4Quasmon::CalculateNumberOfQPartons ( G4double  qMass  ) 

Definition at line 4423 of file G4Quasmon.cc.

References G4cout.

04424 {
04425   static const G4double mK0  = G4QPDGCode(311).GetMass();
04426   // @@ Temporary here. To have 3 quarks in Nucleon Temperature should be < M_N/4 (234 MeV)
04427   // M^2=4*n*(n-1)*T^2 => n = M/2T + 1/2 + T/4M + o(T^3/16M^3)
04428   // @@ Genius (better than 10**(-3) even for n=2!) but useless
04429   //G4double qMOver2T = qMass/(Temperature+Temperature);
04430   //G4double est = qMOver2T+1.+0.125/qMOver2T;
04431   // @@ Longer but exact
04432   G4double qMOverT = qMass/Temperature;
04433   G4int valc = valQ.GetTot();
04434   // .................................
04435   // --- Exponent, Double Split, Poisson 1 ----------------
04441   // --- Uncomment up to here =-----------=^^^^^^^^^
04442   // Exponent ------
04443   //else nOfQ=mq-2*mean*log(G4UniformRand());
04444   // Poisson 1 ------
04446   // Double Split ------
04447   //else
04448   //{
04449   //  G4int imean = static_cast<int>(mean);
04450   //  G4double dm = mean - imean;
04451   //  if(G4UniformRand()>dm) nOfQ=mq+imean+imean;
04452   //  else nOfQ=mq+imean+imean+2;
04453   //}
04454   // .........
04455   // Poisson 2 =-----------=
04456   //if(valc%2==0)nOfQ = 2*RandomPoisson((1.+sqrt(1.+qMOverT*qMOverT))/4.);// abs(b) is even
04457   //else   nOfQ = 1+2*RandomPoisson((1.+sqrt(1.+qMOverT*qMOverT))/4.-0.5);// abs(b) is odd
04458   // Poisson 3 =-----------=
04459   nOfQ = RandomPoisson((1.+sqrt(1.+qMOverT*qMOverT))/2.);
04460   G4int     ev = valc%2;
04461   if      (!ev && nOfQ<2) nOfQ=2; // #of valence quarks is even
04462   else if ( ev && nOfQ<3) nOfQ=3; // #of valence quarks is odd
04463   //
04464 #ifdef debug
04465   G4cout<<"G4Q::Calc#ofQP:QM="<<q4Mom<<qMass<<",T="<<Temperature<<",QC="<<valQ<<",n="<<nOfQ
04466          <<G4endl;
04467 #endif
04468   G4int absb = abs(valQ.GetBaryonNumber());
04469   G4int tabn = 0;
04470   if(absb)tabn=3*absb;      // Minimal QC for baryonic system fragmentation
04471   else    tabn=4;           // Minimal QC for mesonic system fragmentation (@@ ?)
04472   if (nOfQ<tabn) nOfQ=tabn;
04473   G4int nSeaPairs = (nOfQ-valc)/2;
04474   G4int stran = abs(valQ.GetS());
04475   G4int astra = abs(valQ.GetAS());
04476   if(astra>stran) stran=astra;
04477   G4int nMaxStrangeSea=static_cast<int>((qMass-stran*mK0)/(mK0+mK0));//KK is min for s-sea
04478   if (absb) nMaxStrangeSea=static_cast<int>((qMass-absb)/672.); //LambdaK is min for s-sea
04479 #ifdef debug
04480   G4cout<<"G4Q::Calc#ofQP:"<<valQ<<",INtot="<<valc<<",nOfQ="<<nOfQ<<",SeaPairs="<<nSeaPairs
04481         <<G4endl;
04482 #endif
04483   if (nSeaPairs)            // Add/subtract sea pairs to/from initial quark content
04484   {
04485 #ifdef debug
04486     G4int morDec=0;
04487 #endif
04488     if(nSeaPairs>0)valQ.IncQAQ(nSeaPairs,SSin2Gluons);
04489 #ifdef debug
04490     else    morDec=valQ.DecQAQ(-nSeaPairs);
04491     if(morDec) G4cout<<"G4Q::Calc#ofQP: "<<morDec<<" pairs can be reduced more"<<G4endl;
04492 #endif
04493     G4int sSea=valQ.GetS(); // Content of strange quarks
04494     G4int asSea=valQ.GetAS();
04495     if(asSea<sSea) sSea=asSea;
04496     if(sSea>nMaxStrangeSea) // @@@@@@@ Too many strange sea ??????
04497     {
04498 #ifdef debug
04499       G4cout<<"G4Q::Calc#ofQP:**Reduce** S="<<sSea<<",aS="<<asSea<<",maxS="<<nMaxStrangeSea
04500             <<G4endl;
04501 #endif
04502       sSea-=nMaxStrangeSea; // Strange sea excess
04503       valQ.DecS(sSea);      // Reduce strange sea to adoptable limit
04504       valQ.DecAS(sSea);
04505       valQ.IncQAQ(sSea,0.); // Add notstrange sea ????????
04506     }
04507   }
04508   // @@ Chocolate rule --- Temporary (?)
04509   //G4int nmin = valc+valc-2; // Chocolate
04510   //G4int nmin = valc+absb;   // String Junction
04511   //if(nOfQ<nmin) nOfQ=nmin;
04512   // --- End of Temporary
04513 #ifdef debug
04514   G4cout<<"G4Quasmon::Calc#ofQP: *** RESULT IN*** nQ="<<nOfQ<<", FinalQC="<<valQ<<G4endl;
04515 #endif
04516   return nOfQ;
04517 } // End of "CalculateNumberOfQPartons"

void G4Quasmon::ClearOutput (  )  [inline]

Definition at line 168 of file G4Quasmon.hh.

Referenced by KillQuasmon().

00169       {std::for_each(theQHadrons.begin(), theQHadrons.end(), DeleteQHadron());
00170        theQHadrons.clear();
00171       }

void G4Quasmon::ClearQuasmon (  )  [inline]

Definition at line 200 of file G4Quasmon.hh.

Referenced by KillQuasmon().

00201 {
00202   static const G4QContent zeroQC(0,0,0,0,0,0);
00203   static const G4LorentzVector nothing(0.,0.,0.,0.);
00204   phot4M= nothing;
00205   valQ  = zeroQC;
00206   q4Mom = nothing;
00207   status= 0;
00208   std::for_each(theQCandidates.begin(), theQCandidates.end(), DeleteQCandidate());
00209   theQCandidates.clear();
00210 
00211 }

void G4Quasmon::CloseElectromagneticDecays (  )  [static]

Definition at line 189 of file G4Quasmon.cc.

00189 {ElMaDecays=false;}

G4QHadronVector * G4Quasmon::DecayQHadron ( G4QHadron hadron  ) 

Definition at line 5779 of file G4Quasmon.cc.

References G4QHadron::DecayIn2(), G4QHadron::DecayIn3(), G4cerr, G4cout, G4UniformRand, G4QHadron::Get4Momentum(), G4QParticle::GetDecayVector(), G4QHadron::GetMass(), G4QPDGCode::GetPDGCode(), G4QCHIPSWorld::GetQParticle(), G4QHadron::GetQPDG(), G4QHadron::MakeAntiHadron(), G4QParticle::MinMassOfFragm(), and G4QHadron::Set4Momentum().

Referenced by G4QIonIonCollision::Breeder(), and G4QFragmentation::Breeder().

05780 {
05781   G4QHadronVector* theFragments = new G4QHadronVector;  // user is responsible to delete!
05782   G4QPDGCode theQPDG  = qH->GetQPDG();
05783   G4int      thePDG   = theQPDG.GetPDGCode();    // Get the PDG code of decaying hadron
05784   G4int        pap = 0;                          // --- particle
05785   if(thePDG<0) pap = 1;                          // --- anti-particle
05786   G4LorentzVector t = qH->Get4Momentum();        // Get 4-momentum of decaying hadron
05787   G4double m_value = t.m();                      // Get the mass value of decaying Hadron
05788   // --- Randomize a channel of decay
05789   G4QDecayChanVector decV = theWorld->GetQParticle(theQPDG)->GetDecayVector();
05790   G4int nChan = decV.size();
05791 #ifdef debug
05792   G4cout<<"G4Quasm::DecQHadron: PDG="<<thePDG<<",m="<<m_value<<",("<<nChan<<" channels)"<<G4endl;
05793 #endif
05794   if(nChan)
05795   {
05796     G4int i=0;
05797     if(nChan>1)
05798     {
05799       G4double rnd = G4UniformRand();            // Random value to select a Decay Channel
05800       for(i=0; i<nChan; i++)
05801       {
05802         G4QDecayChan* dC = decV[i];              // The i-th Decay Channel
05803 #ifdef debug
05804         G4cout<<"G4Quasmon::DecaQHadr:i="<<i<<",r="<<rnd<<"<dl="<<dC->GetDecayChanLimit()
05805               <<", mm="<<dC->GetMinMass()<<G4endl;
05806 #endif
05807         if(rnd<dC->GetDecayChanLimit() && m_value>dC->GetMinMass()) break;
05808       }
05809       if(i>nChan-1) i=nChan-1;
05810     }
05811     G4QPDGCodeVector cV=decV[i]->GetVecOfSecHadrons();// PDGVector of theSelectedDecChannel
05812     G4int nPart=cV.size();                         // A#of particles to decay in
05813 #ifdef debug
05814     G4cout<<"G4Quasmon::DecayQHadron: resi="<<i<<",nP="<<nPart<<":"<<cV[0]->GetPDGCode()
05815           <<","<<cV[1]->GetPDGCode();
05816     if(nPart>2) G4cout<<","<<cV[2]->GetPDGCode();
05817     G4cout<<G4endl;
05818 #endif
05819     if(nPart<2||nPart>3)
05820     {
05821       G4cerr<<"---Warning---G4Q::DecayQHadr:n="<<nPart<<",ch#"<<i<<",PDG="<<thePDG<<G4endl;
05822       theFragments->push_back(qH);                    // Fill as it is (del.equiv.)
05823       return theFragments;
05824     }
05825 #ifdef debug
05826     G4cout<<"G4Q::DecQH:Decay("<<ElMaDecays<<") PDG="<<thePDG<<t<<m_value<<",nP="<<nPart<<G4endl;
05827 #endif
05828     if(nPart==2)
05829     {
05830       G4QHadron* fHadr;
05831       G4QHadron* sHadr;
05832       G4int fPDG=cV[0]->GetPDGCode();
05833       G4int sPDG=cV[1]->GetPDGCode();
05834       // Radiative decays In2 (eta, eta', Sigma0) are closed if the ElMaDecays=false
05835       if ( (fPDG != 22 && sPDG != 22) || ElMaDecays) {
05836 #ifdef debug
05837         G4cout<<"G4Q::DecQH:Yes2,fPDG="<<fPDG<<",sPDG="<<sPDG<<",EMF="<<ElMaDecays<<G4endl;
05838 #endif
05839         if(cV[0]->GetWidth()==0.)
05840         { // Randomize only the second Hardon or none
05841           fHadr = new G4QHadron(cV[0]->GetPDGCode());   // the First Hadron is created *1*
05842           if(cV[1]->GetWidth()==0.)sHadr = new G4QHadron(sPDG);//theSecondHadron is created
05843           else
05844           {
05845             G4QParticle* sPart=theWorld->GetQParticle(cV[1]);// Pt for theSecondHadron
05846             G4double sdm = m_value - fHadr->GetMass();  // MaxMassLimit for the 2-nd Hadron
05847             sHadr = new G4QHadron(sPart,sdm);           // the Second Hadron is created *2*
05848             if(sPDG<0) sHadr->MakeAntiHadron();
05849           }
05850         }
05851         else                                              // Randomize masses ofBothHadrons
05852         {
05853           G4QParticle* sPart=theWorld->GetQParticle(cV[1]);// Pt for theSecondHadron
05854           G4double mim = sPart->MinMassOfFragm();         // MinMassLimit for theSecondHadr
05855           G4double fdm = m_value - mim;                   // MaxMassLimit for theFirstHadr
05856           G4QParticle* fPart=theWorld->GetQParticle(cV[0]);// Pt for the First Hadron
05857           fHadr = new G4QHadron(fPart,fdm);               // the 1-st Hadron is initialized
05858           if(fPDG<0) fHadr->MakeAntiHadron();
05859           G4double fm=fHadr->GetMass();                   // Mass of the first hadron
05860           G4double sdm = m_value - fm;                    // MaxMassLimit for theSecondHadr
05861           sHadr = new G4QHadron(sPart,sdm);               // the 2-nd Hadron is initialized
05862           if(sPDG<0) sHadr->MakeAntiHadron();
05863 #ifdef debug
05864           G4cout<<"G4Q::DQH:M="<<m_value<<",mi="<<mim<<",fd="<<fdm<<",fm="<<fm<<",sd="<<sdm
05865                 <<",sm="<<sHadr->GetMass()<<G4endl;
05866 #endif
05867         }
05868 #ifdef debug
05869         G4cout<<"G4Q::DQH:(DecayIn2)1="<<fHadr->GetMass()<<",2="<<sHadr->GetMass()<<G4endl;
05870 #endif
05871         if(pap)
05872         {
05873           fHadr->MakeAntiHadron();
05874           sHadr->MakeAntiHadron();
05875         }
05876         G4LorentzVector f4Mom = fHadr->Get4Momentum();    // Get First Hadron 4Mom (mass) 
05877         G4LorentzVector s4Mom = sHadr->Get4Momentum();    // Get Second Hadron 4Mom (mass) 
05878         if(!qH->DecayIn2(f4Mom,s4Mom))                    // Error in DecayIn2
05879         {
05880           delete fHadr;                                   // Delete "new fHadr"
05881           delete sHadr;                                   // Delete "new sHadr"
05882 #ifdef debug
05883           G4cerr<<"---Warning---G4Q::DecayQHadron:in2,PDGC="<<thePDG<<", ch#"<<i<<": 4M="
05884                 <<qH->Get4Momentum()<<"("<<qH->GetMass()<<")->"<<f4Mom<<"+"<<s4Mom<<G4endl;
05885           //throw G4QException("***Exception***G4Q::DecayQHadron: Failed to decay in 2");
05886 #endif
05887           theFragments->push_back(qH);                    // Fill as it is (del.equiv.)
05888           return theFragments;
05889         }
05890         else
05891         {
05892           //qH->SetNFragments(2);
05893           //theFragments.push_back(qH);                   // Fill with NFr=2 (del.equiv.)
05894           // Instead
05895           delete qH;                                      // Delete it (without History)
05896           //
05897           fHadr->Set4Momentum(f4Mom);             // Put the randomized 4Mom to 1-st Hadron
05898           G4QHadronVector* theTmpQHV=DecayQHadron(fHadr); // Try to decay
05899           G4int nProd=theTmpQHV->size();
05900 #ifdef debug
05901           G4cout<<"G4Q::DecayQHadr:(DecayIn2) nOfProdForQH1="<<nProd<<G4endl;
05902 #endif
05903           if(nProd==1) theFragments->push_back((*theTmpQHV)[0]);// Final = no Further Decay
05904           else for(G4int ip1=0; ip1<nProd; ip1++)
05905           {
05906             G4QHadronVector* intTmpQHV = DecayQHadron((*theTmpQHV)[ip1]);
05907             G4int tmpS=intTmpQHV->size();
05908             if(tmpS==1)theFragments->push_back((*intTmpQHV)[0]);// Final = no Further Decay
05909             else
05910             {
05911               theFragments->resize(tmpS+theFragments->size());// Resize theFragments length
05912               copy(intTmpQHV->begin(), intTmpQHV->end(), theFragments->end()-tmpS);
05913             }
05914 #ifdef debug
05915             G4cout<<"G4Q::DecayQHadr:(DecayIn2) Copy Sec11 nProd="<<tmpS<<G4endl;
05916 #endif
05917             intTmpQHV->clear();
05918             delete intTmpQHV;
05919           }
05920           theTmpQHV->clear();
05921           delete theTmpQHV;
05922           sHadr->Set4Momentum(s4Mom);             // Put the randomized 4Mom to 2-nd Hadron
05923           theTmpQHV=DecayQHadron(sHadr);          // Try to decay
05924           nProd=theTmpQHV->size();
05925 #ifdef debug
05926           G4cout<<"G4Q::DecayQHadr:(DecayIn2) nOfProdForQH2="<<nProd<<G4endl;
05927 #endif
05928           if(nProd==1) theFragments->push_back((*theTmpQHV)[0]);// Final = no Further Decay
05929           else for(G4int ip1=0; ip1<nProd; ip1++)
05930           {
05931             G4QHadronVector* intTmpQHV = DecayQHadron((*theTmpQHV)[ip1]);
05932             G4int tmpS=intTmpQHV->size();
05933             if(tmpS==1)theFragments->push_back((*intTmpQHV)[0]);// Final = no Further Decay
05934             else
05935             {
05936               theFragments->resize(tmpS+theFragments->size());// Resize theFragments length
05937               copy(intTmpQHV->begin(), intTmpQHV->end(), theFragments->end()-tmpS);
05938             }
05939 #ifdef debug
05940             G4cout<<"G4Q::DecayQHadr:(DecayIn2) Copy Sec12 nProd="<<tmpS<<G4endl;
05941 #endif
05942             intTmpQHV->clear();
05943             delete intTmpQHV;
05944           }
05945           theTmpQHV->clear();
05946           delete theTmpQHV;
05947         }
05948 #ifdef debug
05949         G4cout<<"G4Q::DecQHadr: DecayIn2 is made with nH="<<theFragments->size()<<G4endl;
05950 #endif
05951       }
05952       else
05953       {
05954 #ifdef debug
05955         if(thePDG==89999003||thePDG==90002999)G4cerr<<"*G4Q::DQH:8999003/90002999"<<G4endl;
05956 #endif
05957         theFragments->push_back(qH);               // Fill hadron as it is (del.equivalent)
05958       }
05959     }
05960     else if(nPart==3)
05961     {
05962       G4QHadron* fHadr;
05963       G4QHadron* sHadr;
05964       G4QHadron* tHadr;
05965       G4int fPDG=cV[0]->GetPDGCode();
05966       G4int sPDG=cV[1]->GetPDGCode();
05967       G4int tPDG=cV[2]->GetPDGCode();
05968       //The radiative decays of the GS hadrons In3 are closed if ElMaDecays=false
05969       if ( (fPDG != 22 && sPDG != 22 && tPDG != 22) || ElMaDecays)
05970       {
05971 #ifdef debug
05972         G4cout<<"G4Q::DQH:Y,f="<<fPDG<<",s="<<sPDG<<",t="<<tPDG<<",F="<<ElMaDecays<<G4endl;
05973 #endif
05974         if(cV[0]->GetWidth()==0.)                        // Don't randomize theFirstHardon
05975         {
05976           fHadr = new G4QHadron(fPDG);                   // theFirst Hadron is created  *1*
05977           if(cV[1]->GetWidth()==0.)
05978           {
05979             sHadr = new G4QHadron(sPDG);                 // theSecond Hadron is created *2*
05980             if(cV[2]->GetWidth()==0.)tHadr = new G4QHadron(tPDG);//theThirdHadron isCreated
05981             else
05982             {
05983               G4QParticle* tPart=theWorld->GetQParticle(cV[2]);// Pt for the3-rdH
05984               G4double tdm = m_value-fHadr->GetMass()-sHadr->GetMass();// MaxMass for the 2d Hadr
05985               tHadr = new G4QHadron(tPart,tdm);                  //the3rdHadron is created
05986               if(tPDG<0) tHadr->MakeAntiHadron();
05987             }
05988           }
05989           else                                              // Randomize 2nd & 3rd Hadrons
05990           {
05991             m_value-=fHadr->GetMass();                       // Reduce the residual MaxMass
05992             G4QParticle* tPart=theWorld->GetQParticle(cV[2]);// Pt for the 3-rd Hadron
05993             G4double mim = tPart->MinMassOfFragm();         // MinMassLimit for the 3rd Hd
05994             G4double sdm = m_value - mim;                   // MaxMassLimit for the 2nd Hd
05995             G4QParticle* sPart=theWorld->GetQParticle(cV[1]);// Pt for the 2-nd Hadron
05996             sHadr = new G4QHadron(sPart,sdm);               // theSecondHadron is created
05997             if(sPDG<0) sHadr->MakeAntiHadron();
05998             G4double tdm = m_value - sHadr->GetMass();      // MaxMassLimit for the 3-rd H
05999             tHadr = new G4QHadron(tPart,tdm);               // the Third Hadron is created
06000             if(tPDG<0) tHadr->MakeAntiHadron();
06001           }
06002         }
06003         else  // Randomize masses of all three Hadrons
06004         {
06005           G4QParticle* sPart=theWorld->GetQParticle(cV[1]); // Pt for theSecondHadr
06006           G4double smim = sPart->MinMassOfFragm();          // MinMassLim for SecondHadron
06007           G4QParticle* tPart=theWorld->GetQParticle(cV[2]); // Pt for the Third Hadron
06008           G4double tmim = tPart->MinMassOfFragm();          // MinMassLimit for theThirdHd
06009           G4double fdm = m_value - smim - tmim;             // MaxMassLimit for theFirstHd
06010           G4QParticle* fPart=theWorld->GetQParticle(cV[0]); // Pt for the First Hadron
06011           fHadr = new G4QHadron(fPart,fdm);                 // the First Hadron is created
06012           if(fPDG<0) fHadr->MakeAntiHadron();
06013           m_value-=fHadr->GetMass();                        // Reduce the residual MaxMass
06014           G4double  sdm = m_value - tmim;                   // MaxMassLimit for theSecondH
06015           sHadr = new G4QHadron(sPart,sdm);                 // theSecondHadron is created
06016           if(sPDG<0) sHadr->MakeAntiHadron();
06017           G4double  tdm = m_value - sHadr->GetMass();       // MaxMassLimit for theThird H
06018           tHadr = new G4QHadron(tPart,tdm);                 // the Third Hadron is created
06019           if(tPDG<0) tHadr->MakeAntiHadron();
06020         }     
06021 #ifdef debug
06022         G4cout<<"G4Quasmon::DecayQHadron:3Dec. m1="<<fHadr->GetMass()
06023               <<",m2="<<sHadr->GetMass()<<",m3="<<tHadr->GetMass()<<G4endl;
06024 #endif
06025         if(pap)
06026         {
06027           fHadr->MakeAntiHadron();
06028           sHadr->MakeAntiHadron();
06029           tHadr->MakeAntiHadron();
06030         }
06031         G4LorentzVector f4Mom = fHadr->Get4Momentum(); // Get 4M of the First Hadron (mass)
06032         G4LorentzVector s4Mom = sHadr->Get4Momentum(); // Get 4M of the SecondHadron (mass)
06033         G4LorentzVector t4Mom = tHadr->Get4Momentum(); // Get 4M of the Third Hadron (mass)
06034         if(!qH->DecayIn3(f4Mom,s4Mom,t4Mom))
06035         {
06036           delete fHadr;                                // Delete "new fHadr"
06037           delete sHadr;                                // Delete "new sHadr"
06038           delete tHadr;                                // Delete "new tHadr"
06039           G4cerr<<"---Warning---G4Q::DecayQHadron:in3,PDGC="<<thePDG<<", ch#"<<i<<G4endl;
06040           theFragments->push_back(qH);             // Fill as it is (delete equivalent)
06041           return theFragments;
06042         }
06043         else
06044         {
06045           //qH->SetNFragments(3);
06046           //theFragments.push_back(q);              // Fill with NFr=3 (del.equiv.)
06047           // Instead
06048           delete qH;
06049           //
06050           fHadr->Set4Momentum(f4Mom);             // Put the randomized 4Mom to 1-st Hadron
06051           G4QHadronVector* theTmpQHV=DecayQHadron(fHadr); // Try to decay
06052           G4int nProd=theTmpQHV->size();
06053 #ifdef debug
06054           G4cout<<"G4Q::DecayQHadr:(DecayIn3) nOfProdForQH1="<<nProd<<G4endl;
06055 #endif
06056           if(nProd==1) theFragments->push_back((*theTmpQHV)[0]);// Final = no Further Decay
06057           else for(G4int ip1=0; ip1<nProd; ip1++)
06058           {
06059             G4QHadronVector* intTmpQHV = DecayQHadron((*theTmpQHV)[ip1]);
06060             G4int tmpS=intTmpQHV->size();
06061             if(tmpS==1)theFragments->push_back((*intTmpQHV)[0]);// Final = no Further Decay
06062             else
06063             {
06064               theFragments->resize(tmpS+theFragments->size());// Resize theFragments length
06065               copy(intTmpQHV->begin(), intTmpQHV->end(), theFragments->end()-tmpS);
06066             }
06067 #ifdef debug
06068             G4cout<<"G4Q::DecayQHadr:(DecayIn3) Copy Sec11 nProd="<<tmpS<<G4endl;
06069 #endif
06070             intTmpQHV->clear();
06071             delete intTmpQHV;
06072           }
06073           theTmpQHV->clear();
06074           delete theTmpQHV;
06075 
06076           sHadr->Set4Momentum(s4Mom);             // Put the randomized 4Mom to 2-nd Hadron
06077           theTmpQHV=DecayQHadron(sHadr);          // Try to decay
06078           nProd=theTmpQHV->size();
06079 #ifdef debug
06080           G4cout<<"G4Q::DecayQHadr:(DecayIn3) nOfProdForQH2="<<nProd<<G4endl;
06081 #endif
06082           if(nProd==1) theFragments->push_back((*theTmpQHV)[0]);// Final = no Further Decay
06083           else for(G4int ip1=0; ip1<nProd; ip1++)
06084           {
06085             G4QHadronVector* intTmpQHV = DecayQHadron((*theTmpQHV)[ip1]);
06086             G4int tmpS=intTmpQHV->size();
06087             if(tmpS==1)theFragments->push_back((*intTmpQHV)[0]);// Final = no Further Decay
06088             else
06089             {
06090               theFragments->resize(tmpS+theFragments->size());// Resize theFragments length
06091               copy(intTmpQHV->begin(), intTmpQHV->end(), theFragments->end()-tmpS);
06092             }
06093 #ifdef debug
06094             G4cout<<"G4Q::DecayQHadr:(DecayIn3) Copy Sec12 nProd="<<tmpS<<G4endl;
06095 #endif
06096             intTmpQHV->clear();
06097             delete intTmpQHV;
06098           }
06099           theTmpQHV->clear();
06100           delete theTmpQHV;
06101 
06102           tHadr->Set4Momentum(t4Mom);             // Put the randomized 4Mom to 3-rd Hadron
06103           theTmpQHV=DecayQHadron(tHadr);          // Try to decay
06104           nProd=theTmpQHV->size();
06105 #ifdef debug
06106           G4cout<<"G4Q::DecayQHadr:(DecayIn3) nOfProdForQH3="<<nProd<<G4endl;
06107 #endif
06108           if(nProd==1) theFragments->push_back((*theTmpQHV)[0]);// Final = no Further Decay
06109           else for(G4int ip1=0; ip1<nProd; ip1++)
06110           {
06111             G4QHadronVector* intTmpQHV = DecayQHadron((*theTmpQHV)[ip1]);
06112             G4int tmpS=intTmpQHV->size();
06113             if(tmpS==1)theFragments->push_back((*intTmpQHV)[0]);// Final = no Further Decay
06114             else
06115             {
06116               theFragments->resize(tmpS+theFragments->size());// Resize theFragments length
06117               copy(intTmpQHV->begin(), intTmpQHV->end(), theFragments->end()-tmpS);
06118             }
06119 #ifdef debug
06120             G4cout<<"G4Q::DecayQHadr:(DecayIn3) Copy Sec13 nProd="<<tmpS<<G4endl;
06121 #endif
06122             intTmpQHV->clear();
06123             delete intTmpQHV;
06124           }
06125           theTmpQHV->clear();
06126           delete theTmpQHV;
06127 
06128         }
06129 #ifdef debug
06130         G4cout<<"G4Q::DecQHadr: DecayIn3 is made with nH="<<theFragments->size()<<G4endl;
06131 #endif
06132       }
06133       else theFragments->push_back(qH);            // Fill hadron as it is (del.equivalent)
06134     }
06135   }
06136   else
06137   {
06138 #ifdef debug
06139     G4cout<<"G4Quas::DecQHadr:Fill PDG= "<<thePDG<<t<<m_value<<" as it is ***0***>>"<<G4endl;
06140 #endif
06141     if(thePDG==89999003||thePDG==90002999)G4cerr<<"-War-G4Q::DQH:8999003/90002999"<<G4endl;
06142     theFragments->push_back(qH);                   // Fill as it is (delete equivalent)
06143   }
06144 #ifdef debug
06145   G4cout<<"G4Q::DecQHadr:=-= HADRON IS DECAYED =-= with nH="<<theFragments->size()<<G4endl;
06146 #endif
06147   return theFragments;
06148 } // End of "DecayOutHadron"

G4QHadronVector * G4Quasmon::DecayQuasmon (  ) 

Definition at line 3769 of file G4Quasmon.cc.

References G4cerr, and G4cout.

03770 {
03771   G4QHadron* thisQuasmon = new G4QHadron(valQ,q4Mom);  // create a Hadron for this Quasmon
03772   FillHadronVector(thisQuasmon);                       // Fill it as a hadron
03773   G4QHadronVector* theFragments = new G4QHadronVector; // user is responsible to delete!
03774   G4int nHadrs=theQHadrons.size();
03775 #ifdef debug
03776   G4cout<<"G4Q::DecayQuasmon:After decay (FillHadronVector byItself) nH="<<nHadrs<<G4endl;
03777 #endif
03778   if(nHadrs) for (int hadron=0; hadron<nHadrs; hadron++)
03779   {
03780     G4QHadron* curHadr = new G4QHadron(theQHadrons[hadron]);
03781     theFragments->push_back(curHadr);                  // (user must delete)
03782   }
03783 #ifdef pdebug
03784   else G4cerr<<"*******G4Quasmon::DecayQuasmon: *** Nothing is in the output ***"<<G4endl;
03785 #endif
03786   valQ=G4QContent(0,0,0,0,0,0);                        // Wipe the Quasmon out
03787   q4Mom=G4LorentzVector(0.,0.,0.,0.);                  // ... with its 4-momentum
03788   return theFragments;
03789 } // End of "DecayQuasmon"

G4QHadronVector * G4Quasmon::Fragment ( G4QNucleus nucEnviron,
G4int  nQ = 1 
)

Definition at line 6178 of file G4Quasmon.cc.

References G4cerr, G4cout, and G4QNucleus::GetProbability().

Referenced by G4QInelastic::PostStepDoIt(), and G4QDiffractionRatio::ProjFragment().

06179 {
06180 #ifdef debug
06181   G4cout<<"G4Quasmon::Fragment called E="<<nucEnviron<<nucEnviron.GetProbability()<<G4endl;
06182 #endif
06183   G4int nQs=nQ;
06184   HadronizeQuasmon(nucEnviron,nQs);
06185   G4int nHadrs=theQHadrons.size();
06186 #ifdef debug
06187   G4cout<<"G4Quasm::Fragm:after HadronizeQuasmon nH="<<nHadrs<<",Env="<<nucEnviron<<G4endl;
06188 #endif
06189   G4QHadronVector* theFragments = new G4QHadronVector;// user is responsible for delition !
06190   if(nHadrs) for (int hadron=0; hadron<nHadrs; hadron++)
06191   {
06192     G4QHadron* curHadr = new G4QHadron(theQHadrons[hadron]);
06193     theFragments->push_back(curHadr);         // (delete equivalent - user)
06194   }
06195 #ifdef pdebug
06196   else G4cerr<<"*******G4Quasmon::Fragment *** Nothing is in the output ***"<<G4endl;
06197 #endif
06198   return theFragments;
06199 } // End of "Fragment"

G4LorentzVector G4Quasmon::Get4Momentum (  )  const [inline]

Definition at line 161 of file G4Quasmon.hh.

Referenced by G4QEnvironment::AddQuasmon(), G4QFragmentation::Fragment(), and IncreaseBy().

00161 {return q4Mom;}

G4int G4Quasmon::GetBaryonNumber (  )  const [inline]

Definition at line 164 of file G4Quasmon.hh.

References G4QContent::GetBaryonNumber().

Referenced by G4QEnvironment::AddQuasmon().

00164 {return valQ.GetBaryonNumber();}

G4int G4Quasmon::GetCharge (  )  const [inline]

Definition at line 163 of file G4Quasmon.hh.

References G4QContent::GetCharge().

Referenced by G4QEnvironment::AddQuasmon().

00163 {return valQ.GetCharge();}

G4double G4Quasmon::GetEtaSup (  )  const [inline]

Definition at line 160 of file G4Quasmon.hh.

00160 {return EtaEtaprime;}

G4QContent G4Quasmon::GetQC (  )  const [inline]

Definition at line 162 of file G4Quasmon.hh.

Referenced by IncreaseBy().

00162 {return valQ;}

G4QPDGCode G4Quasmon::GetQPDG (  )  const [inline]

Definition at line 166 of file G4Quasmon.hh.

00166 {return G4QPDGCode(valQ);}

G4double G4Quasmon::GetSOverU (  )  const [inline]

Definition at line 159 of file G4Quasmon.hh.

00159 {return SSin2Gluons;}

G4int G4Quasmon::GetStatus (  )  const [inline]

Definition at line 167 of file G4Quasmon.hh.

00167 {return status;}

G4int G4Quasmon::GetStrangeness (  )  const [inline]

Definition at line 165 of file G4Quasmon.hh.

References G4QContent::GetStrangeness().

00165 {return valQ.GetStrangeness();}

G4double G4Quasmon::GetTemper (  )  const [inline]

Definition at line 158 of file G4Quasmon.hh.

00158 {return Temperature;}

void G4Quasmon::IncreaseBy ( G4QContent qQCont,
const G4LorentzVector q4M 
) [inline]

Definition at line 186 of file G4Quasmon.hh.

00187 {
00188   valQ  += qQCont;
00189   q4Mom += q4M;
00190   status= 3;
00191 }

void G4Quasmon::IncreaseBy ( const G4Quasmon pQuasm  )  [inline]

Definition at line 180 of file G4Quasmon.hh.

References Get4Momentum(), and GetQC().

00181 {
00182   valQ  += pQuasm->GetQC();
00183   q4Mom += pQuasm->Get4Momentum();
00184   status= 3;
00185 }

void G4Quasmon::InitQuasmon ( const G4QContent qQCont,
const G4LorentzVector q4M 
) [inline]

Definition at line 193 of file G4Quasmon.hh.

00194 {
00195   valQ  = qQCont;
00196   q4Mom = q4M;
00197   status= 3;
00198 }

void G4Quasmon::KillQuasmon (  )  [inline]

Definition at line 213 of file G4Quasmon.hh.

References ClearOutput(), and ClearQuasmon().

00214 {
00215   ClearQuasmon();
00216   ClearOutput();
00217 }

void G4Quasmon::OpenElectromagneticDecays (  )  [static]

Definition at line 186 of file G4Quasmon.cc.

00186 {ElMaDecays=true;}

G4bool G4Quasmon::operator!= ( const G4Quasmon right  )  const [inline]

Definition at line 157 of file G4Quasmon.hh.

00157 {return this != &rhs;}

const G4Quasmon & G4Quasmon::operator= ( const G4Quasmon right  ) 

Definition at line 202 of file G4Quasmon.cc.

References f2all, nBarClust, nOfQ, phot4M, q4Mom, rEP, rMo, status, theQCandidates, theQHadrons, theWorld, and valQ.

00203 {
00204   if(this != &right)                          // Beware of self assignment
00205   {
00206     q4Mom                 = right.q4Mom;
00207     valQ                  = right.valQ;
00208     //theEnvironment        = right.theEnvironment;
00209     status                = right.status;
00210     //theQHadrons (Vector)
00211     G4int iQH             = theQHadrons.size();
00212     if(iQH) for(G4int jh=0; jh<iQH; jh++) delete theQHadrons[jh];
00213     theQHadrons.clear();
00214     G4int nQH             = right.theQHadrons.size();
00215     if(nQH) for(G4int ih=0; ih<nQH; ih++)
00216     {
00217       G4QHadron* curQH    = new G4QHadron(right.theQHadrons[ih]);
00218       theQHadrons.push_back(curQH);
00219     }
00220     theWorld              = right.theWorld;
00221     phot4M                = right.phot4M;
00222     nBarClust             = right.nBarClust;
00223     nOfQ                  = right.nOfQ;
00224     //theQCandidates (Vector)
00225     G4int iQC             = theQCandidates.size();
00226     if(iQC) for(G4int jq=0; jq<iQC; jq++) delete theQCandidates[jq];
00227     theQCandidates.clear();
00228     G4int nQC             = right.theQCandidates.size();
00229     if(nQC) for(G4int iq=0; iq<nQC; iq++)
00230     {
00231       G4QCandidate* curQC = new G4QCandidate(right.theQCandidates[iq]);
00232       theQCandidates.push_back(curQC);
00233     }
00234     f2all                 = right.f2all;
00235     rEP                   = right.rEP;
00236     rMo                   = right.rMo;
00237   }
00238   return *this;
00239 } // End of "="

G4bool G4Quasmon::operator== ( const G4Quasmon right  )  const [inline]

Definition at line 156 of file G4Quasmon.hh.

00156 {return this == &rhs;}

void G4Quasmon::Set4Momentum ( G4LorentzVector  Q4M  )  [inline]

Definition at line 96 of file G4Quasmon.hh.

Referenced by G4QIonIonCollision::Fragment(), and G4QFragmentation::Fragment().

00096 {q4Mom=Q4M;} // Set new value for the Quasmon 4mom

void G4Quasmon::SetEtaSup ( G4double  etaetap  )  [static]

Definition at line 200 of file G4Quasmon.cc.

00200 {EtaEtaprime=etaetap;}

void G4Quasmon::SetParameters ( G4double  temper = 180.,
G4double  ssin2g = .3,
G4double  etaetap = .3 
) [static]

Definition at line 192 of file G4Quasmon.cc.

Referenced by G4ChiralInvariantPhaseSpace::ApplyYourself(), G4QAtomicElectronScattering::G4QAtomicElectronScattering(), G4QCaptureAtRest::G4QCaptureAtRest(), G4QInelastic::G4QInelastic(), G4StringChipsParticleLevelInterface::Propagate(), G4StringChipsInterface::Propagate(), G4QStringChipsParticleLevelInterface::Propagate(), G4QInelastic::SetParameters(), G4QCaptureAtRest::SetParameters(), and G4QAtomicElectronScattering::SetParameters().

00193 {
00194   Temperature=temperature;
00195   SSin2Gluons=ssin2g;
00196   EtaEtaprime=etaetap;
00197 }

void G4Quasmon::SetQC ( G4QContent  QQC  )  [inline]

Definition at line 97 of file G4Quasmon.hh.

Referenced by G4QIonIonCollision::Fragment(), and G4QFragmentation::Fragment().

00097 {valQ=QQC;}  // Set new Quark Cont for the Quasmon

void G4Quasmon::SetSOverU ( G4double  ssin2g  )  [static]

Definition at line 199 of file G4Quasmon.cc.

00199 {SSin2Gluons=ssin2g;}

void G4Quasmon::SetTemper ( G4double  temperature  )  [static]

Definition at line 198 of file G4Quasmon.cc.

00198 {Temperature=temperature;}


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:53:15 2013 for Geant4 by  doxygen 1.4.7