00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037 #ifndef G4QNucleus_h
00038 #define G4QNucleus_h 1
00039
00040 #include <utility>
00041 #include <vector>
00042 #include <CLHEP/Units/SystemOfUnits.h>
00043
00044 #include "globals.hh"
00045 #include "G4RandomDirection.hh"
00046 #include "G4QCandidateVector.hh"
00047 #include "G4QHadronVector.hh"
00048 #include "G4LorentzRotation.hh"
00049 #include "G4QChipolino.hh"
00050
00051 class G4QNucleus : public G4QHadron
00052 {
00053 public:
00054 G4QNucleus();
00055 G4QNucleus(G4int nucPDG);
00056 G4QNucleus(G4LorentzVector p, G4int nucPDG);
00057 G4QNucleus(G4QContent nucQC);
00058 G4QNucleus(G4QContent nucQC, G4LorentzVector p);
00059 G4QNucleus(G4int z, G4int n, G4int s=0);
00060 G4QNucleus(G4int z, G4int n, G4int s, G4LorentzVector p);
00061 G4QNucleus(G4QNucleus* right, G4bool cop3D = false);
00062 G4QNucleus(const G4QNucleus &right, G4bool cop3D=false);
00063 ~G4QNucleus();
00064
00065 const G4QNucleus& operator=(const G4QNucleus& right);
00066 G4bool operator==(const G4QNucleus &right) const {return this==&right;}
00067 G4bool operator!=(const G4QNucleus &right) const {return this!=&right;}
00068
00069 G4int GetPDG() const {return 90000000+1000*(1000*S+Z)+N;}
00070 G4int GetZ() const {return Z;}
00071 G4int GetN() const {return N;}
00072 G4int GetS() const {return S;}
00073 G4int GetA() const {return Z+N+S;}
00074 G4int GetDZ() const {return dZ;}
00075 G4int GetDN() const {return dN;}
00076 G4int GetDS() const {return dS;}
00077 G4int GetDA() const {return dZ+dN+dS;}
00078 G4int GetMaxClust() const {return maxClust;}
00079 G4double GetProbability(G4int bn=0) const {return probVect[bn];}
00080 G4double GetMZNS() const {return GetQPDG().GetNuclMass(Z,N,S);}
00081 G4double GetTbIntegral();
00082 G4double GetGSMass() const {return GetQPDG().GetMass();}
00083 G4QContent GetQCZNS() const
00084 {
00085 if(S>=0) return G4QContent(Z+N+N+S,Z+Z+N+S,S,0,0,0);
00086 else return G4QContent(Z+N+N+S,Z+Z+N+S,0,0,0,-S);
00087 }
00088 G4int GetNDefMesonC() const{return nDefMesonC;};
00089 G4int GetNDefBaryonC()const{return nDefBaryonC;};
00090 G4double GetDensity(const G4ThreeVector&aPos) {return rho0*GetRelativeDensity(aPos);}
00091 G4double GetRho0() {return rho0;}
00092 G4double GetRelativeDensity(const G4ThreeVector& aPosition);
00093 G4double GetRelWSDensity(const G4double& r)
00094 {return 1./(1.+std::exp((r-radius)/WoodSaxonSurf));}
00095 G4double GetRelOMDensity(const G4double& r2){return std::exp(-r2/radius);}
00096 G4double GetRadius(const G4double maxRelativeDenisty=0.5);
00097 G4double GetOuterRadius();
00098 G4double GetDeriv(const G4ThreeVector& point);
00099 G4double GetFermiMomentum(G4double density);
00100 G4QHadron* GetNextNucleon()
00101 {
00102
00103 return (currentNucleon>=0&¤tNucleon<GetA()) ? theNucleons[currentNucleon++] : 0;
00104 }
00105 void SubtractNucleon(G4QHadron* pNucleon);
00106 void DeleteNucleons();
00107 G4LorentzVector GetNucleons4Momentum()
00108 {
00109 G4LorentzVector sum(0.,0.,0.,0.);
00110 for(unsigned i=0; i<theNucleons.size(); i++) sum += theNucleons[i]->Get4Momentum();
00111 sum.setE(std::sqrt(sqr(GetGSMass())+sum.v().mag2()));
00112 return sum;
00113 }
00114 std::vector<G4double> const* GetBThickness() {return &Tb;}
00115
00116
00117 G4bool EvaporateBaryon(G4QHadron* h1,G4QHadron* h2);
00118 void EvaporateNucleus(G4QHadron* hA, G4QHadronVector* oHV);
00119
00120 void DecayDibaryon(G4QHadron* dB, G4QHadronVector* oHV);
00121 void DecayAntiDibaryon(G4QHadron* dB, G4QHadronVector* oHV);
00122 void DecayIsonucleus(G4QHadron* dB, G4QHadronVector* oHV);
00123 void DecayMultyBaryon(G4QHadron* dB, G4QHadronVector* oHV);
00124 void DecayAntiStrange(G4QHadron* dB, G4QHadronVector* oHV);
00125 void DecayAlphaBar(G4QHadron* dB, G4QHadronVector* oHV);
00126 void DecayAlphaDiN(G4QHadron* dB, G4QHadronVector* oHV);
00127 void DecayAlphaAlpha(G4QHadron* dB, G4QHadronVector* oHV);
00128 G4int SplitBaryon();
00129 G4int HadrToNucPDG(G4int hPDG);
00130 G4int NucToHadrPDG(G4int nPDG);
00131 G4bool Split2Baryons();
00132 void ActivateBThickness();
00133 G4double GetBThickness(G4double b);
00134 G4double GetThickness(G4double b);
00135 void InitByPDG(G4int newPDG);
00136 void InitByQC(G4QContent newQC)
00137 {G4int PDG=G4QPDGCode(newQC).GetPDGCode(); InitByPDG(PDG);}
00138 void IncProbability(G4int bn);
00139 void Increase(G4int PDG, G4LorentzVector LV = G4LorentzVector(0.,0.,0.,0.));
00140 void Increase(G4QContent QC, G4LorentzVector LV = G4LorentzVector(0.,0.,0.,0.));
00141 void Reduce(G4int PDG);
00142 void CalculateMass() {Set4Momentum(G4LorentzVector(0.,0.,0.,GetGSMass()));}
00143 void SetMaxClust(G4int maxC){maxClust=maxC;}
00144 void InitCandidateVector(G4QCandidateVector& theQCandidates,
00145 G4int nM=45, G4int nB=72, G4int nC=117);
00146 void PrepareCandidates(G4QCandidateVector& theQCandidates, G4bool piF=false, G4bool
00147 gaF=false, G4LorentzVector LV=G4LorentzVector(0.,0.,0.,0.));
00148 G4int UpdateClusters(G4bool din);
00149 G4QNucleus operator+=(const G4QNucleus& rhs);
00150 G4QNucleus operator-=(const G4QNucleus& rhs);
00151 G4QNucleus operator*=(const G4int& rhs);
00152 G4bool StartLoop();
00153 G4bool ReduceSum(G4ThreeVector* vectors, G4ThreeVector sum);
00154 void SimpleSumReduction(G4ThreeVector* vectors, G4ThreeVector sum);
00155 void DoLorentzBoost(const G4LorentzVector& theBoost)
00156 {
00157 theMomentum.boost(theBoost);
00158 for(unsigned i=0; i<theNucleons.size(); i++) theNucleons[i]->Boost(theBoost);
00159 }
00160 void DoLorentzRotation(const G4LorentzRotation& theLoRot)
00161 {
00162 theMomentum=theLoRot*theMomentum;
00163 for(unsigned i=0; i<theNucleons.size(); i++) theNucleons[i]->LorentzRotate(theLoRot);
00164 }
00165 void DoLorentzBoost(const G4ThreeVector& theBeta)
00166 {
00167 theMomentum.boost(theBeta);
00168 for(unsigned i=0; i<theNucleons.size(); i++) theNucleons[i]->Boost(theBeta);
00169 }
00170 void DoLorentzContraction(const G4LorentzVector&B){DoLorentzContraction(B.vect()/B.e());}
00171 void DoLorentzContraction(const G4ThreeVector& theBeta);
00172 void DoTranslation(const G4ThreeVector& theShift);
00173
00174
00175 static void SetParameters(G4double fN=.1,G4double fD=.05, G4double cP=4., G4double mR=1.,
00176 G4double nD=.8*CLHEP::fermi);
00177
00178
00179 G4int RandomizeBinom(G4double p,G4int N);
00180 G4double CoulombBarGen(const G4double& rZ, const G4double& rA, const G4double& cZ,
00181 const G4double& cA);
00182 G4double CoulombBarrier(const G4double& cZ=1, const G4double& cA=1, G4double dZ=0.,
00183 G4double dA=0.);
00184 G4double FissionCoulombBarrier(const G4double& cZ, const G4double& cA, G4double dZ=0.,
00185 G4double dA=0.);
00186 G4double BindingEnergy(const G4double& cZ=0, const G4double& cA=0, G4double dZ=0.,
00187 G4double dA=0.);
00188 G4double CoulBarPenProb(const G4double& CB, const G4double& E, const G4int& C,
00189 const G4int& B);
00190 std::pair<G4double, G4double> ChooseImpactXandY(G4double maxImpact);
00191 void ChooseNucleons();
00192 void ChoosePositions();
00193 void ChooseFermiMomenta();
00194 void InitDensity();
00195 void Init3D();
00196 private:
00197
00198 void SetZNSQC(G4int z, G4int n, G4int s);
00199 G4QNucleus GetThis() const {return G4QNucleus(Z,N,S);}
00200
00201
00202 private:
00203
00204 static const G4int nDefMesonC =45;
00205 static const G4int nDefBaryonC=72;
00206
00207 static G4double freeNuc;
00208 static G4double freeDib;
00209 static G4double clustProb;
00210 static G4double mediRatio;
00211 static G4double nucleonDistance;
00212 static G4double WoodSaxonSurf;
00213
00214 G4int Z;
00215 G4int N;
00216 G4int S;
00217
00218 G4int dZ;
00219 G4int dN;
00220 G4int dS;
00221 G4int maxClust;
00222 G4double probVect[256];
00223
00224 G4QHadronVector theNucleons;
00225 G4int currentNucleon;
00226 G4double rho0;
00227 G4double radius;
00228 std::vector<G4double> Tb;
00229 G4bool TbActive;
00230 G4bool RhoActive;
00231 };
00232
00233 std::ostream& operator<<(std::ostream& lhs, G4QNucleus& rhs);
00234 std::ostream& operator<<(std::ostream& lhs, const G4QNucleus& rhs);
00235
00236 #endif