3 //------------------------------------------------------------------
4 // This class introduces the weight's calculation
5 // according to the Lednicky's algorithm.
6 // The detailed description of the algorithm can be found
7 // in comments to fortran code:
10 //------------------------------------------------------------------
16 #include "AliHBTLLWeights.h"
17 #include "AliHBTPair.h"
18 #include "AliHBTParticle.h"
19 #include "WLedCOMMONS.h"
21 /*******************************************************************/
22 /****** ROUTINES USED FOR COMMUNUCATION ********/
23 /******************** WITH FORTRAN ********************/
24 /*******************************************************************/
26 # define led_bldata led_bldata_
27 # define fsiini fsiini_
28 # define ltran12 ltran12_
32 # define led_bldata LED_BLDATA
33 # define fsiini FSIINI
34 # define ltran12 LTRAN12
36 # define type_of_call _stdcall
38 /****************************************************************/
39 extern "C" void type_of_call led_bldata();
40 extern "C" void type_of_call fsiini();
41 extern "C" void type_of_call ltran12();
42 extern "C" void type_of_call fsiw();
43 /**************************************************************/
45 ClassImp(AliHBTLLWeights)
47 AliHBTLLWeights* AliHBTLLWeights::fgLLWeights=NULL;
49 AliHBTLLWeights::AliHBTLLWeights()
51 // Default Constructor
55 SetColWithResidNuclSwitch();
56 SetStrongInterSwitch();
57 SetQuantumStatistics();
64 AliHBTLLWeights* AliHBTLLWeights::Instance()
66 // Instantiates new object or returns a pointer to already exitsing one
70 fgLLWeights = new AliHBTLLWeights();
76 Double_t AliHBTLLWeights::GetWeight(const AliHBTPair* partpair)
78 // Returns the weignt of the pair "partpair"
79 AliHBTParticle *part1 = partpair->Particle1();
80 AliHBTParticle *part2 = partpair->Particle2();
82 if ( (part1 == 0x0) || (part2 == 0x0))
84 Error("GetWeight","Null particle pointer");
89 Double_t part1Momentum[]={part1->Px(),part1->Py(),part1->Pz()};
90 Double_t part2Momentum[]={part2->Px(),part2->Py(),part2->Pz()};
92 if ( (part1->Px() == part2->Px()) &&
93 (part1->Py() == part2->Py()) &&
94 (part1->Pz() == part2->Pz()) )
100 if ((!fRandomPosition) &&
101 (part1->Vx() == part2->Vx()) && (part1->Vy() == part2->Vy())
102 && (part1->Vz() == part2->Vz()) )
109 FSI_MOM.P1X=part1Momentum[0];
110 FSI_MOM.P1Y=part1Momentum[1];
111 FSI_MOM.P1Z=part1Momentum[2];
113 FSI_MOM.P2X=part2Momentum[0];
114 FSI_MOM.P2Y=part2Momentum[1];
115 FSI_MOM.P2Z=part2Momentum[2];
117 if (fRandomPosition){
119 Double_t rxcm = fsigma*gRandom->Gaus();
120 Double_t rycm = fsigma*gRandom->Gaus();
121 Double_t rzcm = fsigma*gRandom->Gaus();
123 FSI_PRF.X=rxcm*fwcons;
124 FSI_PRF.Y=rycm*fwcons;
125 FSI_PRF.Z=rzcm*fwcons;
128 Double_t rps=rxcm*rxcm+rycm*rycm+rzcm*rzcm;
129 Double_t rp=TMath::Sqrt(rps);
138 return LEDWEIGHT.WEIN;
142 /************************************************************/
143 void AliHBTLLWeights::Init()
145 //---------------------------------------------------------------------
146 //initial parameters of model
148 FSI_NS.NS = fApproximationModel;
150 if(!ftest){LEDWEIGHT.ITEST=0;}
154 if(fColoumbSwitch){FSI_NS.ICH =1;}
156 if(fStrongInterSwitch){FSI_NS.ISI=1;}
158 if(fQuantStatSwitch){FSI_NS.IQS=1;}
160 if(fColWithResidNuclSwitch){FSI_NS.I3C=1;}
164 if(fRandomPosition){LEDWEIGHT.IRANPOS=1;}
165 else{LEDWEIGHT.IRANPOS=0;}
168 if ( (fPID1 == 0) || (fPID2 == 0) )
170 Fatal("Init","Particles types are not set");
174 FSI_NS.LL = GetPairCode(fPID1,fPID2);
178 Fatal("Init","Particles types are not supported");
183 TParticlePDG* tpart1 = TDatabasePDG::Instance()->GetParticle(fPID1);
186 Fatal("init","We can not find particle with ID=%d in our DataBase",fPID1);
190 FSI_POC.AM1=tpart1->Mass();
191 FSI_POC.C1=tpart1->Charge();
193 TParticlePDG* tpart2 = TDatabasePDG::Instance()->GetParticle(fPID2);
200 Fatal("init","We can not find particle with ID=%d in our DataBase",fPID2);
204 FSI_POC.AM2=tpart2->Mass();
205 FSI_POC.C1=tpart2->Charge();
211 //constants for radii simulation
214 fsigma =TMath::Sqrt(2.)*fRadius;
219 Int_t AliHBTLLWeights::GetPairCode(const AliHBTPair* partpair)
221 // Return the code of the pair "partpair"
222 return GetPairCode(partpair->Particle1()->GetPdgCode(),partpair->Particle2()->GetPdgCode());
225 Int_t AliHBTLLWeights::GetPairCode(Int_t pid1,Int_t pid2)
227 // pairCode 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
228 // hpid: n p n alfa pi+ pi0 pi+ n p pi+ pi+ pi+ pi- K+ K+ K+ K- d d t t K0 K0 d p p p n
229 // lpid: n p p alfa pi- pi0 pi+ d d K- K+ p p K- K+ p p d alfa t alfa K0 K0b t t alfa lambda lambda
230 // NS=1 y/n: + + + + + - - - - - - - - - - - - - - - - - - - - - - -
232 //alphas, deuterons and tyts are NOT supported here
234 Int_t chargefactor = 1;
235 Int_t hpid; //pid in higher row
236 Int_t lpid; //pid in lower row
237 Int_t code; //pairCode
241 //determine the order of selcetion in switch
242 if (TMath::Abs(pid1) < TMath::Abs(pid2) )
244 if (pid1<0) chargefactor=-1;
245 hpid=pid2*chargefactor;
246 lpid=pid1*chargefactor;
251 if (pid2<0) chargefactor=-1;
252 hpid=pid1*chargefactor;
253 lpid=pid2*chargefactor;
262 //Determine the pair code
263 switch (hpid) //switch on first particle id
269 code = 1; //neutron neutron
273 code = 3; //neutron proton
277 code = 28; //neutron lambda
281 return 0; //given pair not supported
290 code = 2; //proton proton
294 code = 27;//proton lambda
298 return 0; //given pair not supported
309 code = 7; //piplus piplus
313 code = 5; //piplus piminus
317 code = 10; //piplus Kminus
321 code = 11; //piplus Kplus
325 code = 12; //piplus proton
330 return 0; //given pair not supported
342 return 0; //given pair not supported
351 code = 14; //Kplus Kminus
355 code = 15; //Kplus Kplus
359 code = 16; //Kplus proton
363 return 0; //given pair not supported
372 code = 17; //Kminus proton
377 return 0; //given pair not supported
386 code = 2; //Kzero Kzero
390 code = 17; //Kzero KzeroBar
394 return 0; //given pair not supported