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
00040
00041
00042
00043 #include "G4QChipolino.hh"
00044 #include <cstdlib>
00045 using namespace std;
00046
00047 G4QChipolino::G4QChipolino(G4QContent& QCont)
00048 {
00049
00050 G4QPDGCode Pi0(111);
00051 G4double mPi0 = Pi0.GetMass();
00052 G4QContent Pi0QC = Pi0.GetQuarkContent();
00053 G4int ban =QCont.GetBaryonNumber();
00054 G4int tban=abs(3*ban);
00055 G4int tot=QCont.GetTot();
00056 G4int tod=tot%2;
00057 if ( (!tod && (tot < 4 || (ban && tot < tban) ) ) || (tod && tot < tban+2) )
00058 QCont.IncQAQ(1,0.);
00059
00060 G4QContent rQC=QCont;
00061 tot=rQC.GetTot();
00062 if (tot%2)rQC.DecQAQ(-tban-2);
00063 else if(ban)rQC.DecQAQ(-tban);
00064 else rQC.DecQAQ(-4);
00065 tot=rQC.GetTot();
00066 #ifdef debug
00067 cout<<"G4QChipolino is called with QC="<<QCont<<",rQC="<<rQC<<",tot="<<tot<<G4endl;
00068 #endif
00069 minM=1000000.;
00070 theQPDG1 = Pi0;
00071 theQPDG2 = Pi0;
00072 theQCont1 = Pi0QC;
00073 if (!tot)
00074 {
00075 G4cerr<<"***G4QChipolino: shouldn't be here 1 QC="<<rQC<<G4endl;
00076 }
00077 else if (tot==2 || tot==3)
00078 {
00079 G4cerr<<"***G4QChipolino: shouldn't be here 2 QC="<<rQC<<G4endl;
00080 theQCont1= rQC;
00081 theQPDG1.InitByQCont(rQC);
00082 theQCont = rQC+Pi0QC;
00083 }
00084 else if (tot==4)
00085 {
00086 G4QContent bQC=rQC.IndQ();
00087 #ifdef debug
00088 G4cout<<"G4QChipolino: tot=4,rQC="<<rQC<<",bQC="<<bQC<<G4endl;
00089 #endif
00090 for(int j=0; j<2; j++)
00091 {
00092 G4QContent aQC=rQC.IndAQ(j);
00093 G4QContent cQC=bQC+aQC;
00094 G4QPDGCode cQPDG(cQC);
00095 G4double M1=cQPDG.GetMass();
00096 if(cQPDG.GetPDGCode()==221) M1=mPi0;
00097 G4QContent oQC=rQC-cQC;
00098 #ifdef debug
00099 cout<<"G4QChipolino: aQC="<<aQC<<", cQC="<<cQC<<", oQC="<<oQC<<G4endl;
00100 #endif
00101 G4QPDGCode oQPDG(oQC);
00102 G4double M2=oQPDG.GetMass();
00103 if(oQPDG.GetPDGCode()==221) M2=mPi0;
00104 G4double m_value=M1+M2;
00105 #ifdef debug
00106 cout<<"G4QChipolino: c="<<cQPDG<<",cM="<<M1<<",o="<<oQPDG<<",oM="<<M2
00107 <<",cM+0M="<<m_value<<", curMinM="<<minM<<G4endl;
00108 #endif
00109 if(m_value<minM)
00110 {
00111 minM=m_value;
00112 theQPDG1 = cQPDG;
00113 theQCont1 = cQC;
00114 theQPDG2 = oQPDG;
00115 }
00116 }
00117 }
00118 else if (tot==5)
00119 {
00120 G4int nQ=rQC.GetQ();
00121 G4int nA=rQC.GetAQ();
00122 G4bool fl=nA>nQ;
00123 #ifdef pdebug
00124 cout<<"G4QChipolino: Baryon case nQ="<<nQ<<",nA="<<nA<<",QC="<<rQC
00125 <<",fl="<<fl<<G4endl;
00126 #endif
00127 G4QContent bQC;
00128 if (fl) bQC=rQC.IndQ();
00129 else bQC=rQC.IndAQ();
00130 for (int i=0; i<4; i++)
00131 {
00132 G4QContent cQC;
00133 if (fl) cQC=bQC+rQC.IndAQ(i);
00134 else cQC=bQC+rQC.IndQ(i);
00135 G4QPDGCode cQPDG(cQC);
00136 G4double M1=cQPDG.GetMass();
00137 if(cQPDG.GetPDGCode()==221) M1=mPi0;
00138 G4QContent oQC=rQC-cQC;
00139 G4QPDGCode oQPDG(oQC);
00140 G4double M2=oQPDG.GetMass();
00141 if(oQPDG.GetPDGCode()==221) M2=mPi0;
00142 G4double m_value=M1+M2;
00143 if(m_value<minM)
00144 {
00145 minM=m_value;
00146 theQPDG1 = cQPDG;
00147 theQCont1 = cQC;
00148 theQPDG2 = oQPDG;
00149 }
00150 }
00151 #ifdef pdebug
00152 cout<<"G4QChipolino: Baryon case minM="<<minM<<", M="<<theQCont1<<theQPDG1
00153 <<", B="<<theQPDG2<<G4endl;
00154 #endif
00155 }
00156 else if (tot==6)
00157 {
00158 if(ban)
00159 {
00160 G4int nQ=rQC.GetQ();
00161 G4int nA=rQC.GetAQ();
00162 G4bool fl=nA>nQ;
00163 #ifdef debug
00164 cout<<"G4QChipolino: Di-Bar. case nQ="<<nQ<<",nA="<<nA<<",QC="<<rQC<<",fl="<<fl<<G4endl;
00165 #endif
00166 for (int i=0; i<4; i++)
00167 {
00168 G4QContent aQC;
00169 if (fl) aQC=rQC.IndAQ(i);
00170 else aQC=rQC.IndQ(i);
00171 for (int j=i+1; j<5; j++)
00172 {
00173 G4QContent bQC;
00174 if (fl) bQC=aQC+rQC.IndAQ(j);
00175 else bQC=aQC+rQC.IndQ(j);
00176 for (int k=j+1; k<6; k++)
00177 {
00178 G4QContent cQC;
00179 if (fl) cQC=bQC+rQC.IndAQ(k);
00180 else cQC=bQC+rQC.IndQ(k);
00181 G4QPDGCode cQPDG(cQC);
00182 G4double M1=cQPDG.GetMass();
00183 if(cQPDG.GetPDGCode()==221) M1=mPi0;
00184 G4QContent oQC=rQC-cQC;
00185 G4QPDGCode oQPDG=(oQC);
00186 G4double M2=oQPDG.GetMass();
00187 if(oQPDG.GetPDGCode()==221) M2=mPi0;
00188 G4double m_value=M1+M2;
00189 if(m_value<minM)
00190 {
00191 minM=m_value;
00192 theQPDG1 = cQPDG;
00193 theQCont1 = cQC;
00194 theQPDG2 = oQPDG;
00195 }
00196 }
00197 }
00198 }
00199 }
00200 else
00201 {
00202 theQCont1 = rQC.IndQ(0)+rQC.IndQ(1)+rQC.IndQ(2);
00203 theQPDG1.InitByQCont(theQCont1);
00204 theQPDG2.InitByQCont(rQC.IndAQ(0)+rQC.IndAQ(1)+rQC.IndAQ(2));
00205 }
00206 }
00207 else if(((rQC.GetU() )>(rQC.GetS() -4) && (rQC.GetD() )>(rQC.GetS() -4)) ||
00208 ((rQC.GetAU())>(rQC.GetAS()-4) && (rQC.GetAD())>(rQC.GetAS()-4)) )
00209 {
00210 G4int kD=rQC.GetD();
00211 G4int kU=rQC.GetU();
00212 G4int kS=rQC.GetS();
00213 G4int mD=rQC.GetAD();
00214 G4int mU=rQC.GetAU();
00215 G4int mS=rQC.GetAS();
00216 G4int nQ=rQC.GetQ();
00217 G4int nA=rQC.GetAQ();
00218 G4bool fl=nA>nQ;
00219 #ifdef debug
00220 G4cout<<"G4QChipolino: NucFragment case nQ="<<nQ<<",nAQ="<<nA<<", QC="<<rQC<<",fl="<<fl
00221 <<G4endl;
00222 #endif
00223 if( (fl && kS>1) || (!fl && mS>1))
00224 {
00225 #ifdef debug
00226 G4cerr<<"***G4QChipolino: ***Overfowed by strange quarks*** rQC="<<rQC<<G4endl;
00227
00228 #endif
00229 }
00230 else if(fl)
00231 {
00232
00233
00234 if(!mS)
00235 {
00236 G4int nI=mU-mD;
00237 G4int nN=(mU+mD-nI*3)/6;
00238 if(!kS)
00239 {
00240 if((nI>=0&&nN>=0)||(nI<0&&nN>=-nI))
00241 {
00242 if(nI>0)
00243 {
00244 theQPDG1 = G4QPDGCode(-(90000000+1000*(nN+nI-1)+nN));
00245 theQPDG2 = G4QPDGCode(-2212);
00246 }
00247 else
00248 {
00249 theQPDG1 = G4QPDGCode(-(90000000+1000*(nN+nI)+nN-1));
00250 theQPDG2 = G4QPDGCode(-2112);
00251 }
00252 }
00253 else if((nI>=0&&nN>-2)||(nI<0&&nN>-nI-2))
00254 {
00255 if(nI>0)
00256 {
00257 theQPDG1=G4QPDGCode(-(90000000+1000*(nN+nI-2)+nN+1));
00258 theQPDG2=G4QPDGCode(-2224);
00259 }
00260 else
00261 {
00262 theQPDG1=G4QPDGCode(-(90000000+1000*(nN+nI+1)+nN-2));
00263 theQPDG2=G4QPDGCode(-1114);
00264 }
00265 }
00266 else
00267 {
00268 G4cerr<<"***G4QChipolino:**A**IsotopicAsymmetry (without S),rQC="<<rQC<<G4endl;
00269
00270 }
00271 }
00272 else if(kS<2)
00273 {
00274 nN =(mU+mD-4-nI*3)/6;
00275 if(nI>0)
00276 {
00277 nN+=1;
00278 theQPDG1 = G4QPDGCode(-(90000000+1000*(nN+nI-1)+nN));
00279 theQPDG2 = G4QPDGCode(-321);
00280 }
00281 else
00282 {
00283 theQPDG1 = G4QPDGCode(-(90000000+1000*(nN+nI+1)+nN));
00284 theQPDG2 = G4QPDGCode(-311);
00285 }
00286 }
00287 else
00288 {
00289 G4cerr<<"***G4QChipolino: ***Too many kaons are needed*** rQC="<<rQC<<G4endl;
00290
00291 }
00292 }
00293 else
00294 {
00295 if(mS<=mU&&mS<=mD)
00296 {
00297 G4int nI=mU-mD;
00298 G4int nN=(mU+mD-mS-mS-nI*3)/6;
00299 if((nI>=0&&nN>=0)||(nI<0&&nN>=-nI))
00300 {
00301 if(nI>0)
00302 {
00303 theQPDG1 = G4QPDGCode(-(90000000+1000*(kS*1000+nN+nI-1)+nN));
00304 theQPDG2 = G4QPDGCode(-2212);
00305 }
00306 else
00307 {
00308 theQPDG1 = G4QPDGCode(-(90000000+1000*(kS*1000+nN+nI)+nN-1));
00309 theQPDG2 = G4QPDGCode(-2112);
00310 }
00311 }
00312 else if((nI>=0&&nN>-2)||(nI<0&&nN>-nI-2))
00313 {
00314 if(nI>0)
00315 {
00316 theQPDG1=G4QPDGCode(-(90000000+1000*(kS*1000+nN+nI-2)+nN+1));
00317 theQPDG2=G4QPDGCode(-2224);
00318 }
00319 else
00320 {
00321 theQPDG1=G4QPDGCode(-(90000000+1000*(kS*1000+nN+nI+1)+nN-2));
00322 theQPDG2=G4QPDGCode(-1114);
00323 }
00324 }
00325 else
00326 {
00327 G4cerr<<"***G4QChipolino:**A**IsotopicAssimetry (with S)*** rQC="<<rQC<<G4endl;
00328
00329 }
00330 }
00331 else
00332 {
00333 G4int lam=mU;
00334 if (lam>mD) lam=mD;
00335 G4int lD=mD-lam;
00336 G4int lU=mU-lam;
00337 G4int lS=mS-lam;
00338 if(lD+lU+lS!=3||lD<0||lU<0||lS<0)
00339 {
00340 G4cerr<<"***G4QChipolino:*AntiFragment* rQC="<<rQC<<",s="<<lS<<",u="<<lU<<",d"
00341 <<lD<<G4endl;
00342
00343 }
00344 if ( !lD && lU==2) theQPDG2=G4QPDGCode(-3222);
00345 else if( !lU && lD==2) theQPDG2=G4QPDGCode(-3112);
00346 else if( !lD && lU==1) theQPDG2=G4QPDGCode(-3322);
00347 else if( !lU && lD==1) theQPDG2=G4QPDGCode(-3312);
00348 else theQPDG2=G4QPDGCode(-3334);
00349 theQPDG1=G4QPDGCode(-(90+lam)*1000000);
00350 }
00351 theQCont1 = rQC-theQPDG2.GetQuarkContent();
00352 theQCont = rQC;
00353 }
00354 }
00355 else
00356 {
00357 if(!kS)
00358 {
00359 G4int nI=kU-kD;
00360 G4int nN=(kU+kD-nI*3)/6;
00361 if(!mS)
00362 {
00363 if((nI>=0&&nN>=0)||(nI<0&&nN>=-nI))
00364 {
00365 if(nI>0)
00366 {
00367 theQPDG1 = G4QPDGCode(90000000+1000*(nN+nI-1)+nN);
00368 theQPDG2 = G4QPDGCode(2212);
00369 }
00370 else
00371 {
00372 theQPDG1 = G4QPDGCode(90000000+1000*(nN+nI)+nN-1);
00373 theQPDG2 = G4QPDGCode(2112);
00374 }
00375 }
00376 else if((nI>=0&&nN>-2)||(nI<0&&nN>-nI-2))
00377 {
00378 if(nI>0)
00379 {
00380 theQPDG1=G4QPDGCode(90000000+1000*(nN+nI-2)+nN+1);
00381 theQPDG2=G4QPDGCode(2224);
00382 }
00383 else
00384 {
00385 theQPDG1=G4QPDGCode(90000000+1000*(nN+nI+1)+nN-2);
00386 theQPDG2=G4QPDGCode(1114);
00387 }
00388 }
00389 else
00390 {
00391 G4cerr<<"***G4QChipolino:***Isotopic assimetry (without S), rQC="<<rQC<<G4endl;
00392
00393 }
00394 }
00395 else if(mS<2)
00396 {
00397 nN =(kU+kD-4-nI*3)/6;
00398 if(nI>0)
00399 {
00400 nN+=1;
00401 theQPDG1 = G4QPDGCode(90000000+1000*(nN+nI-1)+nN);
00402 theQPDG2 = G4QPDGCode(321);
00403 }
00404 else
00405 {
00406 theQPDG1 = G4QPDGCode(90000000+1000*(nN+nI+1)+nN);
00407 theQPDG2 = G4QPDGCode(311);
00408 }
00409 }
00410 else
00411 {
00412 G4cerr<<"***G4QChipolino: ***Too many kaons are needed*** rQC="<<rQC<<G4endl;
00413
00414 }
00415 }
00416 else
00417 {
00418 if(kS<=kU&&kS<=kD)
00419 {
00420 G4int nI=kU-kD;
00421 G4int nN=(kU+kD-kS-kS-nI*3)/6;
00422 if((nI>=0&&nN>=0)||(nI<0&&nN>=-nI))
00423 {
00424 if(nI>0)
00425 {
00426 theQPDG1 = G4QPDGCode(90000000+1000*(kS*1000+nN+nI-1)+nN);
00427 theQPDG2 = G4QPDGCode(2212);
00428 }
00429 else
00430 {
00431 theQPDG1 = G4QPDGCode(90000000+1000*(kS*1000+nN+nI)+nN-1);
00432 theQPDG2 = G4QPDGCode(2112);
00433 }
00434 }
00435 else if((nI>=0&&nN>-2)||(nI<0&&nN>-nI-2))
00436 {
00437 if(nI>0)
00438 {
00439 theQPDG1=G4QPDGCode(90000000+1000*(kS*1000+nN+nI-2)+nN+1);
00440 theQPDG2=G4QPDGCode(2224);
00441 }
00442 else
00443 {
00444 theQPDG1=G4QPDGCode(90000000+1000*(kS*1000+nN+nI+1)+nN-2);
00445 theQPDG2=G4QPDGCode(1114);
00446 }
00447 }
00448 else
00449 {
00450 G4cerr<<"***G4QChipolino: ***Isotopic assimetry (with S)*** rQC="<<rQC<<G4endl;
00451
00452 }
00453 }
00454 else
00455 {
00456 G4int lam=kU;
00457 if (lam>kD) lam=kD;
00458 G4int lD=kD-lam;
00459 G4int lU=kU-lam;
00460 G4int lS=kS-lam;
00461 if(lD+lU+lS!=3||lD<0||lU<0||lS<0)
00462 {
00463 G4cerr<<"***G4QChipolino:*Fragment*rQC="<<rQC<<",s="<<lS<<",u="<<lU<<",d"
00464 <<lD<<G4endl;
00465
00466 }
00467 if ( !lD && lU==2) theQPDG2=G4QPDGCode(3222);
00468 else if( !lU && lD==2) theQPDG2=G4QPDGCode(3112);
00469 else if( !lD && lU==1) theQPDG2=G4QPDGCode(3322);
00470 else if( !lU && lD==1) theQPDG2=G4QPDGCode(3312);
00471 else theQPDG2=G4QPDGCode(3334);
00472 theQPDG1=G4QPDGCode((90+lam)*1000000);
00473 }
00474 theQCont1 = rQC-theQPDG2.GetQuarkContent();
00475 theQCont = rQC;
00476 }
00477 }
00478 }
00479 else
00480 {
00481 G4cerr<<"***G4QChipolino: ***Exotics*** rQC="<<rQC<<G4endl;
00482
00483 }
00484 }
00485
00486 G4QChipolino::G4QChipolino(const G4QChipolino& right)
00487 {
00488 theQPDG1 = right.theQPDG1;
00489 theQPDG2 = right.theQPDG2;
00490 theQCont = right.theQCont;
00491 theQCont1 = right.theQCont1;
00492 minM = right.minM;
00493 }
00494
00495 G4QChipolino::G4QChipolino(G4QChipolino* right)
00496 {
00497 theQPDG1 = right->theQPDG1;
00498 theQPDG2 = right->theQPDG2;
00499 theQCont = right->theQCont;
00500 theQCont1 = right->theQCont1;
00501 minM = right->minM;
00502 }
00503
00504 const G4QChipolino& G4QChipolino::operator=(const G4QChipolino &right)
00505 {
00506 if(this != &right)
00507 {
00508 theQPDG1 = right.theQPDG1;
00509 theQPDG2 = right.theQPDG2;
00510 theQCont = right.theQCont;
00511 theQCont1 = right.theQCont1;
00512 minM = right.minM;
00513 }
00514 return *this;
00515 }
00516
00517 G4QChipolino::~G4QChipolino() {}
00518
00519
00520 ostream& operator<<(ostream& lhs, G4QChipolino& rhs)
00521 {
00522 lhs<<"{1="<<rhs.GetQPDG1()<<",2="<<rhs.GetQPDG2()<< "}";
00523 return lhs;
00524 }
00525
00526
00527
00528
00529
00530