git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Update by Ludmila Malinina
[u/mrichter/AliRoot.git] / HBTAN / AliHBTLLWeights.cxx
CommitLineData
7f92929e 1#include "AliHBTLLWeights.h"
2#include "AliPDG.h"
3#include "AliHBTPair.h"
4#include "AliHBTParticle.h"
5#include <TList.h>
6#include <TRandom.h>
7#include <TMath.h>
8
9/*******************************************************************/
10/****** ROUTINES USED FOR COMMUNUCATION ********/
11/******************** WITH FORTRAN ********************/
12/*******************************************************************/
13#ifndef WIN32
14# define led_bldata led_bldata_
15# define fsiini fsiini_
16# define ltran12 ltran12_
17# define fsiw fsiw_
18# define type_of_call
19#else
20# define led_bldata LED_BLDATA
21# define fsiini FSIINI
22# define ltran12 LTRAN12
23# define fsiw FSIW
24# define type_of_call _stdcall
25#endif
26/****************************************************************/
27extern "C" void type_of_call led_bldata();
28extern "C" void type_of_call fsiini();
29extern "C" void type_of_call ltran12();
30extern "C" void type_of_call fsiw();
31/**************************************************************/
32
33ClassImp(AliHBTLLWeights)
7f92929e 34
35AliHBTLLWeights* AliHBTLLWeights::fgLLWeights=NULL;
36
37AliHBTLLWeights::AliHBTLLWeights()
38{
39// Default Constructor
40 fPID1 = 0;
41 fPID2 = 0;
42 SetRandomPosition();
43 SetColWithResidNuclSwitch();
44 SetStrongInterSwitch();
45 SetQuantumStatistics();
46 SetColoumb();
47 SetTest();
2f8eea63 48
7f92929e 49}
50
51
52 AliHBTLLWeights* AliHBTLLWeights::Instance()
53{
2f8eea63 54 if (fgLLWeights) {
55 return fgLLWeights;
56 } else {
57 fgLLWeights = new AliHBTLLWeights();
58 return fgLLWeights;
7f92929e 59 }
60}
61
62
63Double_t AliHBTLLWeights::GetWeight(const AliHBTPair* partpair)
64{
65 AliHBTParticle *part1 = partpair->Particle1();
66 AliHBTParticle *part2 = partpair->Particle2();
2f8eea63 67
7f92929e 68 if ( (part1 == 0x0) || (part2 == 0x0))
2f8eea63 69 {
7f92929e 70 Error("GetWeight","Null particle pointer");
71 return 0.0;
2f8eea63 72 }
73
7f92929e 74
75 Double_t part1Momentum[]={part1->Px(),part1->Py(),part1->Pz()};
76 Double_t part2Momentum[]={part2->Px(),part2->Py(),part2->Pz()};
77
2f8eea63 78 if ( (part1->Px() == part2->Px()) &&
79 (part1->Py() == part2->Py()) &&
80 (part1->Pz() == part2->Pz()) )
81 {
82 return 0.0;
83 }
7f92929e 84
85
2f8eea63 86 if ((!fRandomPosition) &&
87 (part1->Vx() == part2->Vx()) && (part1->Vy() == part2->Vy())
88 && (part1->Vz() == part2->Vz()) )
89 {
90 return 0.0;
91 }
7f92929e 92
7f92929e 93
2f8eea63 94
95 FSI_MOM.P1X=part1Momentum[0];
96 FSI_MOM.P1Y=part1Momentum[1];
97 FSI_MOM.P1Z=part1Momentum[2];
7f92929e 98
2f8eea63 99 FSI_MOM.P2X=part2Momentum[0];
100 FSI_MOM.P2Y=part2Momentum[1];
101 FSI_MOM.P2Z=part2Momentum[2];
102
103 if (fRandomPosition){
104
105 Double_t rxcm = fsigma*gRandom->Gaus();
106 Double_t rycm = fsigma*gRandom->Gaus();
107 Double_t rzcm = fsigma*gRandom->Gaus();
7f92929e 108
2f8eea63 109 FSI_PRF.X=rxcm*fwcons;
110 FSI_PRF.Y=rycm*fwcons;
111 FSI_PRF.Z=rzcm*fwcons;
112 FSI_PRF.T=0.;
113
114 Double_t rps=rxcm*rxcm+rycm*rycm+rzcm*rzcm;
115 Double_t rp=TMath::Sqrt(rps);
116 FSI_PRF.RP=rp;
117 FSI_PRF.RPS=rps;
118
119 }
120
121 ltran12();
122 fsiw();
123
124 return LEDWEIGHT.WEIN;
125
7f92929e 126 }
127
128/************************************************************/
129void AliHBTLLWeights::Init()
130 {
131//---------------------------------------------------------------------
132
133//initial parameters of model
134
2f8eea63 135 FSI_NS.NS = approximationModel;
136
137 if(!ftest){LEDWEIGHT.ITEST=0;}
7f92929e 138
2f8eea63 139 if(ftest){
140 LEDWEIGHT.ITEST=1;
141 if(fColoumbSwitch){FSI_NS.ICH =1;}
142 else{FSI_NS.ICH=0;}
143 if(fStrongInterSwitch){FSI_NS.ISI=1;}
144 else{FSI_NS.ISI=0;}
145 if(fQuantStatSwitch){FSI_NS.IQS=1;}
146 else{FSI_NS.IQS=0;}
147 if(fColWithResidNuclSwitch){FSI_NS.I3C=1;}
148 else{FSI_NS.I3C=0;}
149 }
150
151 if(fRandomPosition){LEDWEIGHT.IRANPOS=1;}
152 else{LEDWEIGHT.IRANPOS=0;}
7f92929e 153
7f92929e 154
2f8eea63 155 if ( (fPID1 == 0) || (fPID2 == 0) )
156 {
157 Fatal("Init","Particles types are not set");
158 return;//pro forma
159 }
7f92929e 160
2f8eea63 161 FSI_NS.LL = GetPairCode(fPID1,fPID2);
162
163 if (FSI_NS.LL == 0)
164 {
165 Fatal("Init","Particles types are not supported");
166 return;//pro forma
167 }
7f92929e 168
7f92929e 169
2f8eea63 170 TParticlePDG* tpart1 = TDatabasePDG::Instance()->GetParticle(fPID1);
171 if (tpart1 == 0x0)
172 {
173 Fatal("init","We can not find particle with ID=%d in our DataBase",fPID1);
174 return;
175 }
176
177 FSI_POC.AM1=tpart1->Mass();
178 FSI_POC.C1=tpart1->Charge();
179
180 TParticlePDG* tpart2 = TDatabasePDG::Instance()->GetParticle(fPID2);
181//mlv
182
183
184
185 if (tpart2 == 0x0)
186 {
187 Fatal("init","We can not find particle with ID=%d in our DataBase",fPID2);
188 return;
189 }
190
191 FSI_POC.AM2=tpart2->Mass();
192 FSI_POC.C1=tpart2->Charge();
7f92929e 193
2f8eea63 194 led_bldata();
195 fsiini();
7f92929e 196
197
198//constants for radii simulation
199
2f8eea63 200 if(fRandomPosition){
201 fsigma =TMath::Sqrt(2.)*fRadius;
202 fwcons =FSI_CONS.W;
203 }
7f92929e 204}
205
206Int_t AliHBTLLWeights::GetPairCode(const AliHBTPair* partpair)
207{
208 return GetPairCode(partpair->Particle1()->GetPdgCode(),partpair->Particle2()->GetPdgCode());
209}
210
211Int_t AliHBTLLWeights::GetPairCode(Int_t pid1,Int_t pid2)
212{
213// 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
214// 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
215// 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
216// NS=1 y/n: + + + + + - - - - - - - - - - - - - - - - - - - - - - -
217
218//alphas, deuterons and tyts are NOT supported here
219
220 Int_t chargefactor = 1;
221 Int_t hpid; //pid in higher row
222 Int_t lpid; //pid in lower row
223 Int_t code; //pairCode
224
225 Bool_t swap;
226
227//determine the order of selcetion in switch
228 if (TMath::Abs(pid1) < TMath::Abs(pid2) )
229 {
230 if (pid1<0) chargefactor=-1;
231 hpid=pid2*chargefactor;
232 lpid=pid1*chargefactor;
233 swap = kFALSE;
234 }
235 else
236 {
237 if (pid2<0) chargefactor=-1;
238 hpid=pid1*chargefactor;
239 lpid=pid2*chargefactor;
240 swap = kTRUE;
241 }
242
2f8eea63 243//mlv
244 hpid=pid1;
245 lpid=pid2;
246
247
7f92929e 248//Determine the pair code
249 switch (hpid) //switch on first particle id
250 {
251 case kNeutron:
252 switch (lpid)
253 {
254 case kNeutron:
255 code = 1; //neutron neutron
256 break;
257
258 case kProton:
259 code = 3; //neutron proton
260 break;
261
262 case kLambda0:
263 code = 28; //neutron lambda
264 break;
265
266 default:
267 return 0; //given pair not supported
268 break;
269 }
270 break;
271
272 case kProton:
273 switch (lpid)
274 {
275 case kProton:
276 code = 2; //proton proton
277 break;
278
279 case kLambda0:
280 code = 27;//proton lambda
281 break;
282
283 default:
284 return 0; //given pair not supported
285 break;
286
287 }
288 break;
289
290 case kPiPlus:
2f8eea63 291
7f92929e 292 switch (lpid)
293 {
294 case kPiPlus:
295 code = 7; //piplus piplus
296 break;
297
298 case kPiMinus:
299 code = 5; //piplus piminus
300 break;
301
302 case kKMinus:
303 code = 10; //piplus Kminus
304 break;
305
306 case kKPlus:
307 code = 11; //piplus Kplus
308 break;
309
310 case kProton:
311 code = 12; //piplus proton
312 chargefactor*=-1;
313 break;
314
315 default:
316 return 0; //given pair not supported
317 break;
318 }
319 break;
320 case kPi0:
321 switch (lpid)
322 {
323 case kPi0:
324 code = 6;
325 break;
326
327 default:
328 return 0; //given pair not supported
329 break;
330 }
331 break;
332
333 case kKPlus:
334 switch (lpid)
335 {
336 case kKMinus:
337 code = 14; //Kplus Kminus
338 break;
339
340 case kKPlus:
341 code = 15; //Kplus Kplus
342 break;
343
344 case kProton:
345 code = 16; //Kplus proton
346 break;
347
348 default:
349 return 0; //given pair not supported
350 break;
351 }
352 break;
353
354 case kKMinus:
355 switch (lpid)
356 {
357 case kProton:
358 code = 17; //Kminus proton
359 chargefactor*=1;
360 break;
361
362 default:
363 return 0; //given pair not supported
364 break;
365 }
366 break;
367
368 case kK0:
369 switch (lpid)
370 {
371 case kK0:
372 code = 2; //Kzero Kzero
373 break;
374
375 case kK0Bar:
376 code = 17; //Kzero KzeroBar
377 break;
378
379 default:
380 return 0; //given pair not supported
381 break;
382 }
383 break;
384
385 default: return 0;
386 }
387 return code;
388}
389