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
00038
00039 #include "G4DiffractiveSplitableHadron.hh"
00040
00041 #include "G4ParticleDefinition.hh"
00042 #include "Randomize.hh"
00043
00044 G4DiffractiveSplitableHadron::G4DiffractiveSplitableHadron()
00045
00046 {
00047 PartonIndex=-1;
00048 Parton[0] = new G4Parton(1);
00049 Parton[1] = new G4Parton(-1);
00050 }
00051
00052 G4DiffractiveSplitableHadron::G4DiffractiveSplitableHadron(const G4ReactionProduct & aPrimary)
00053 : G4VSplitableHadron(aPrimary)
00054 {
00055 PartonIndex=-2;
00056 Parton[0]=NULL;
00057 }
00058
00059 G4DiffractiveSplitableHadron::G4DiffractiveSplitableHadron(const G4Nucleon & aNucleon)
00060 : G4VSplitableHadron(aNucleon)
00061 {
00062 PartonIndex=-2;
00063 Parton[0]=NULL;
00064 }
00065
00066 G4DiffractiveSplitableHadron::G4DiffractiveSplitableHadron(const G4VKineticNucleon * aNucleon)
00067 : G4VSplitableHadron(aNucleon)
00068 {
00069 PartonIndex=-2;
00070 Parton[0]=NULL;
00071 }
00072
00073 G4DiffractiveSplitableHadron::~G4DiffractiveSplitableHadron()
00074 {
00075
00076
00077 }
00078
00079
00080 void G4DiffractiveSplitableHadron::SplitUp()
00081 {
00082
00083 if (IsSplit()) return;
00084 Splitting();
00085
00086 if (Parton[0] != NULL) return;
00087
00088
00089
00090 G4int PDGcode=GetDefinition()->GetPDGEncoding();
00091
00092 G4int stringStart, stringEnd;
00093 ChooseStringEnds(PDGcode, &stringStart,&stringEnd);
00094
00095 Parton[0] = new G4Parton(stringStart);
00096 Parton[1] = new G4Parton(stringEnd);
00097 PartonIndex=-1;
00098 }
00099
00100 G4Parton * G4DiffractiveSplitableHadron::GetNextParton()
00101 {
00102 ++PartonIndex;
00103 if ( PartonIndex > 1 || PartonIndex < 0 ) return NULL;
00104 G4int PartonInd(PartonIndex);
00105 if(PartonIndex == 1) PartonIndex=-1;
00106 return Parton[PartonInd];
00107
00108 }
00109
00110 G4Parton * G4DiffractiveSplitableHadron::GetNextAntiParton()
00111 {
00112 ++PartonIndex;
00113 if ( PartonIndex > 1 || PartonIndex < 0 ) return NULL;
00114 G4int PartonInd(PartonIndex);
00115 if(PartonIndex == 1) PartonIndex=-1;
00116 return Parton[PartonInd];
00117
00118 }
00119
00120 void G4DiffractiveSplitableHadron::SetFirstParton(G4int PDGcode)
00121 {
00122 delete Parton[0];
00123 Parton[0]=new G4Parton(PDGcode);
00124 }
00125
00126 void G4DiffractiveSplitableHadron::SetSecondParton(G4int PDGcode)
00127 {
00128 delete Parton[1];
00129 Parton[1]=new G4Parton(PDGcode);
00130 }
00131
00132
00133
00134 void G4DiffractiveSplitableHadron::ChooseStringEnds(G4int PDGcode,G4int * aEnd, G4int * bEnd) const
00135 {
00136 const G4double udspin1= 1./6.;
00137 const G4double uuspin1= 1./3.;
00138
00139
00140 G4int absPDGcode=std::abs(PDGcode);
00141
00142 if ( absPDGcode < 1000 )
00143 {
00144 G4int heavy= absPDGcode/ 100;
00145 G4int light= (absPDGcode %100)/10;
00146
00147
00148 G4int anti= 1 -2 * ( std::max( heavy, light ) % 2 );
00149 if (PDGcode < 0 ) anti *=-1;
00150
00151 heavy *= anti;
00152 light *= -1 * anti;
00153
00154 if ( G4UniformRand() < 0.5 )
00155 {
00156 *aEnd=heavy;
00157 *bEnd=light;
00158 }
00159 else
00160 {
00161 *aEnd=light;
00162 *bEnd=heavy;
00163 }
00164 }
00165 else
00166 {
00167 G4int j1000 = PDGcode/ 1000;
00168 G4int j100 = (PDGcode % 1000) / 100;
00169 G4int j10 = (PDGcode % 100) / 10;
00170
00171 G4double random= G4UniformRand();
00172
00173
00174 if ( std::abs(j100) >= std::abs(j10) )
00175 {
00176 if ( random < udspin1 )
00177 {
00178 *aEnd=j1000;
00179 *bEnd= Diquark( j100, j10, 1);
00180 } else if ( random < (udspin1 + uuspin1) )
00181 {
00182 *aEnd= j10;
00183 *bEnd= Diquark( j1000, j100, 1);
00184 } else
00185 {
00186 *aEnd=j100;
00187
00188
00189 *bEnd= Diquark( j1000, j10, j1000 != j100 ? 0 : 1);
00190 }
00191 } else
00192 {
00193
00194 if ( random < udspin1 )
00195 {
00196 *aEnd=j1000;
00197
00198 *bEnd= Diquark( j100, j10, j100 != j10 ? 0 : 10);
00199 } else if ( random < (udspin1 + uuspin1) )
00200 {
00201 *aEnd= j10;
00202 *bEnd= Diquark( j1000, j100, 1);
00203 } else
00204 {
00205 *aEnd=j100;
00206 *bEnd= Diquark( j1000, j10, 1);
00207 }
00208
00209 }
00210
00211 }
00212
00213
00214
00215 }
00216
00217 G4int G4DiffractiveSplitableHadron::Diquark(G4int aquark,G4int bquark,G4int Spin) const
00218 {
00219 G4int diquarkPDG;
00220
00221 diquarkPDG = std::max( std::abs(aquark), std::abs(bquark) )*1000 +
00222 std::min( std::abs(aquark), std::abs(bquark) )*100 +
00223 2*Spin + 1;
00224 return ( aquark > 0 && bquark > 0 ) ? diquarkPDG : -1*diquarkPDG;
00225 }