]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITStrack.cxx
New version from Boris and Enrico with improved ghost removal
[u/mrichter/AliRoot.git] / ITS / AliITStrack.cxx
CommitLineData
2a309a77 1#include <iostream.h>
2#include <TMath.h>
3#include <TVector.h>
4#include <TMatrix.h>
5#include <TObjArray.h>
6
7#include "AliITS.h"
8#include "AliRun.h"
9#include "AliITStrack.h"
10#include "AliGenerator.h"
8af13b4b 11#include "AliITSRad.h"
2a309a77 12
13
14ClassImp(AliITStrack)
15
16AliITStrack::AliITStrack() {
17//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
18// default constructor
19
53e7090a 20 fC00=fC10=fC11=fC20=fC21=fC22=fC30=fC31=fC32=fC33=fC40=fC41=fC42=fC43=fC44=0.;
2a309a77 21 flistCluster = new TObjArray;
22 fNumClustInTrack =0;
23 fChi2=-1;
24 flabel =0;
25 fVertex.ResizeTo(3);
26 fErrorVertex.ResizeTo(3);
27 fLayer = -1;
28 ClusterInTrack = new TMatrix(6,9);
29 Int_t i;
30 for(i=0; i<6; i++) (*ClusterInTrack)(i,6)=(*ClusterInTrack)(i,7)=
31 (*ClusterInTrack)(i,8)=-1.;
32 rtrack=0.;
2a309a77 33 d2.ResizeTo(6);
34 tgl2.ResizeTo(6);
35 dtgl.ResizeTo(6);
36}
37
38
39
40AliITStrack::AliITStrack(const AliITStrack &cobj) {
41//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
42
2a309a77 43 ClusterInTrack = new TMatrix(6,9);
44 Int_t i;
45 for(i=0; i<6; i++) (*ClusterInTrack)(i,6)=(*ClusterInTrack)(i,7)=
46 (*ClusterInTrack)(i,8)=-1.;
47 flistCluster = new TObjArray;
48 fVertex.ResizeTo(3);
49 fErrorVertex.ResizeTo(3);
50 fVertex = cobj.fVertex;
51 fErrorVertex = cobj.fErrorVertex;
52 flabel = cobj.flabel;
53 fLayer=cobj.fLayer;
54 fTPCtrack = cobj.fTPCtrack;
55 fNumClustInTrack = cobj.fNumClustInTrack;
56 fChi2= cobj.fChi2;
53e7090a 57 fX0=cobj.fX0; fX1=cobj.fX1; fX2=cobj.fX2; fX3=cobj.fX3; fX4=cobj.fX4;
2a309a77 58 rtrack=cobj.rtrack;
59 Dv=cobj.Dv;
60 Zv=cobj.Zv;
61 sigmaDv=cobj.sigmaDv;
62 sigmaZv=cobj.sigmaZv;
63 d2.ResizeTo(6);
64 tgl2.ResizeTo(6);
65 dtgl.ResizeTo(6);
66 d2=cobj.d2;
67 tgl2=cobj.tgl2;
68 dtgl=cobj.dtgl;
8af13b4b 69
53e7090a 70 fC00=cobj.fC00; fC10=cobj.fC10; fC11=cobj.fC11; fC20=cobj.fC20; fC21=cobj.fC21;
71 fC22=cobj.fC22; fC30=cobj.fC30; fC31=cobj.fC31; fC32=cobj.fC32; fC33=cobj.fC33;
72 fC40=cobj.fC40; fC41=cobj.fC41; fC42=cobj.fC42; fC43=cobj.fC43; fC44=cobj.fC44;
2a309a77 73
74 *ClusterInTrack = *cobj.ClusterInTrack;
75
76 for(i=0; i<cobj.flistCluster->GetSize(); i++)
77 flistCluster->AddLast(cobj.flistCluster->At(i));
78
79}
80
81AliITStrack::AliITStrack(AliTPCtrack &obj)
82{
83//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
84
85 fTPCtrack = &obj;
2a309a77 86 fVertex.ResizeTo(3);
87 fErrorVertex.ResizeTo(3);
88
89 d2.ResizeTo(6);
90 tgl2.ResizeTo(6);
91 dtgl.ResizeTo(6);
92 AliGenerator *gener = gAlice->Generator();
93 Float_t Vxg,Vyg,Vzg;
94 gener->GetOrigin(Vxg,Vyg,Vzg);
95
96
97 fVertex(0)=(Double_t)Vxg;
98 fVertex(1)=(Double_t)Vyg;
99 fVertex(2)=(Double_t)Vzg;
100
101 fLayer = 7;
53e7090a 102 //fmCovariance = new TMatrix(5,5);
2a309a77 103 ClusterInTrack = new TMatrix(6,9);
104
105 Int_t i;
106 for(i=0; i<6; i++) (*ClusterInTrack)(i,6)=(*ClusterInTrack)(i,7)=
107 (*ClusterInTrack)(i,8)=-1.;
108 flistCluster = new TObjArray;
109 fNumClustInTrack = 0;
110 LmTPC();
111
112}
113
53e7090a 114
2a309a77 115AliITStrack::~AliITStrack() {
116
117 //destructor
53e7090a 118
2a309a77 119 if(flistCluster) delete flistCluster;
120 if(ClusterInTrack) delete ClusterInTrack;
121
122}
123
53e7090a 124void AliITStrack::PutCElements(Double_t C00, Double_t C10, Double_t C11, Double_t C20, Double_t C21,
125Double_t C22, Double_t C30, Double_t C31, Double_t C32, Double_t C33, Double_t C40,
126Double_t C41, Double_t C42, Double_t C43, Double_t C44){
127
128 fC00=C00; fC10=C10; fC11=C11; fC20=C20; fC21=C21; fC22=C22; fC30=C30; fC31=C31; fC32=C32; fC33=C33;
129 fC40=C40; fC41=C41; fC42=C42; fC43=C43; fC44=C44;
130}
131
132void AliITStrack::GetCElements(Double_t &C00, Double_t &C10, Double_t &C11, Double_t &C20, Double_t &C21,
133Double_t &C22, Double_t &C30, Double_t &C31, Double_t &C32, Double_t &C33, Double_t &C40,
134Double_t &C41, Double_t &C42, Double_t &C43, Double_t &C44){
2a309a77 135
53e7090a 136 C00=fC00; C10=fC10; C11=fC11; C20=fC20; C21=fC21; C22=fC22; C30=fC30; C31=fC31; C32=fC32; C33=fC33;
137 C40=fC40; C41=fC41; C42=fC42; C43=fC43; C44=fC44;
138
139}
140
141void AliITStrack::GetXElements(Double_t &X0, Double_t &X1, Double_t &X2, Double_t &X3, Double_t &X4) {
142 X0=fX0; X1=fX1; X2=fX2; X3=fX3; X4=fX4;
143}
144
145void AliITStrack::PutXElements(Double_t X0, Double_t X1, Double_t X2, Double_t X3, Double_t X4){
146 fX0=X0; fX1=X1; fX2=X2; fX3=X3; fX4=X4;
147}
2a309a77 148
149void AliITStrack::LmTPC() {
150
151// Transform the TPC state vector from TPC-local to master and build a new state vector ITS-type
152// The covariance matrix is also modified accordingly
153//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
154
53e7090a 155
2a309a77 156 Double_t Alpha = fTPCtrack->GetAlpha();
157
158 //printf("LmTPC: Alpha %f\n",Alpha);
53e7090a 159
160 Double_t YTPC = fTPCtrack->GetY();
161 Double_t ZTPC = fTPCtrack->GetZ();
162 Double_t CTPC = fTPCtrack->GetC();
163 Double_t EtaTPC = fTPCtrack->GetEta();
164 Double_t TglTPC = fTPCtrack->GetTgl();
165
2a309a77 166
167 Double_t xm, ym, zm;
168 Double_t sina = TMath::Sin(Alpha);
169 Double_t cosa = TMath::Cos(Alpha);
170 Double_t xl= fTPCtrack->GetX();
53e7090a 171 xm = xl * cosa - YTPC*sina;
172 ym = xl * sina + YTPC*cosa;
173 zm = ZTPC;
2a309a77 174 //cout<<" xl e alpha = "<<xl<<" "<<Alpha<<"\n"; getchar();
175
176 Double_t x0m,y0m;
177
178 ///////////////////////////////////// determine yo //////////////////////////////////////////////////
53e7090a 179
2a309a77 180 Double_t Vxl=fVertex(0)*cosa+fVertex(1)*sina;
181 Double_t Vyl= -fVertex(0)*sina+fVertex(1)*cosa;
182 Double_t Xo,Yo, signy;
53e7090a 183 Double_t R = 1./CTPC;
184 Xo = EtaTPC / CTPC;
2a309a77 185 xoTPC=Xo;
186 Double_t Yo1, Yo2, diffsq1, diffsq2;
53e7090a 187 Yo1 = YTPC + TMath::Sqrt(R*R - (xl-Xo)*(xl-Xo));
188 Yo2 = YTPC - TMath::Sqrt(R*R - (xl-Xo)*(xl-Xo));
2a309a77 189 diffsq1=TMath::Abs((Yo1-Vyl)*(Yo1-Vyl)+(Xo-Vxl)*(Xo-Vxl)-R*R);
190 diffsq2=TMath::Abs((Yo2-Vyl)*(Yo2- Vyl)+(Xo-Vxl)*(Xo-Vxl)-R*R);
191 if(diffsq1<diffsq2) {Yo=Yo1; signy=1.;} else {Yo=Yo2; signy=-1.;};
53e7090a 192
2a309a77 193 ////////////////////////////////////////////////////////////////////////////////////////////////////
194
195 x0m = Xo * cosa - Yo * sina;
196 y0m = Xo * sina + Yo * cosa;
197
198 rtrack=TMath::Sqrt(xm*xm+ym*ym);
199 Double_t pigre=TMath::Pi();
200 Double_t phi=0.0;
201 if(ym == 0. || xm == 0.) {
202 if(ym == 0. && xm == 0.) {cout << " Error in AliITStrack::LmTPC x=0 and y=0 \n"; getchar();}
203 if(ym ==0. && xm>0.) phi=0.;
204 if(ym==0. && xm<0.) phi=pigre;
205 if(xm==0 && ym>0.) phi=pigre/2.;
206 if(xm==0 && ym<0.) phi=1.5*pigre;
207 }
208 else {
209 if (xm>0. && ym>0.) phi=TMath::ATan(ym/xm);
210 if (xm<0. && ym>0.) phi=pigre+TMath::ATan(ym/xm);
211 if (xm<0. && ym<0.) phi=pigre+TMath::ATan(ym/xm);
212 if (xm>0. && ym<0.) phi=(2*pigre)+TMath::ATan(ym/xm);
213 };
214 if(phi<0. || phi>(2*pigre)) {cout<<"attention error on phi in AliITStrack:LmTPC \n"; getchar();}
215
53e7090a 216 fX0=phi;
217 fX1=zm;
218 fX3=TglTPC;
219 fX4=CTPC;
8af13b4b 220
2a309a77 221
222 Double_t dd=TMath::Sqrt((x0m-fVertex(0))*(x0m-fVertex(0))+(y0m-fVertex(1))*(y0m-fVertex(1)));
223 Double_t signdd;
224 if (R>0) signdd=1.; else signdd=-1.;
53e7090a 225 fX2=signdd*dd-R;
2a309a77 226 //cout<<" fvertex = "<<fVertex(0)<<" "<<fVertex(1)<<" "<<fVertex(2)<<"\n";
53e7090a 227
2a309a77 228 Double_t cov[15];
229 fTPCtrack->GetCovariance(cov);
8af13b4b 230
2a309a77 231 Double_t dfidy, dDdy, dDdC, dDdeta;
232
233 dfidy=(xm*cosa+ym*sina)/(rtrack*rtrack);
234 dDdy=signdd*((y0m-fVertex(1))*cosa-(x0m-fVertex(0))*sina)/dd;
53e7090a 235 Double_t dyodr=signy*(R+(xl-Xo)*EtaTPC)/TMath::Sqrt(R*R-(xl-Xo)*(xl-Xo));
236 Double_t dyomdr=sina*EtaTPC+cosa*dyodr;
237 Double_t dxomdr=cosa*EtaTPC-sina*dyodr;
2a309a77 238 Double_t ddddR=((x0m-fVertex(0))*dxomdr+(y0m-fVertex(1))*dyomdr)/dd;
239 dDdC=-R*R*(signdd*ddddR-1.);
240 Double_t dyoldxol=signy*(xl-Xo)/TMath::Sqrt(R*R-(xl-Xo)*(xl-Xo));
241 Double_t dxomdeta=R*(cosa-sina*dyoldxol);
242 Double_t dyomdeta=R*(sina+cosa*dyoldxol);
243 dDdeta=signdd*((x0m-fVertex(0))*dxomdeta+(y0m-fVertex(1))*dyomdeta)/dd;
53e7090a 244
245 Double_t F00=dfidy;
246 Double_t F20=dDdy;
247 Double_t F22=dDdC;
248 Double_t F23=dDdeta;
249
250 Double_t T00=cov[0]*F00;
251 Double_t T02=cov[0]*F20+cov[6]*F22+cov[3]*F23;
252 Double_t T20=cov[6]*F00;
253 Double_t T22=cov[6]*F20+cov[9]*F22+cov[8]*F23;
254
255 fC00=F00*T00;
256 fC10=cov[1]*F00;
257 fC11=cov[2];
258 fC20=F20*T00+F22*T20+F23*cov[3]*F00;
259 fC21=F20*cov[1]+F22*cov[7]+F23*cov[4];
260 fC22=F20*T02+F22*T22+F23*(cov[3]*F20+cov[8]*F22+cov[5]*F23);
261 fC30=cov[10]*F00;
262 fC31=cov[11];
263 fC32=cov[10]*F20+cov[13]*F22+cov[12]*F23;
264 fC33=cov[14];
265 fC40=T20;
266 fC41=cov[7];
267 fC42=T22;
268 fC43=cov[13];
269 fC44=cov[9];
270
271 //cout<<" C32 e C44 = "<<fC32<<" "<<fC44<<"\n"; getchar();
272
2a309a77 273}
274
275
276AliITStrack &AliITStrack::operator=(AliITStrack obj) {
277//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
278
2a309a77 279 delete flistCluster;
280 delete ClusterInTrack;
2a309a77 281 ClusterInTrack = new TMatrix(6,9);
282 flistCluster = new TObjArray;
283 flabel = obj.flabel;
53e7090a 284 fTPCtrack = obj.fTPCtrack;
2a309a77 285 fNumClustInTrack = obj.fNumClustInTrack;
286 fChi2= obj.fChi2;
287 fVertex=obj.fVertex;
288 fErrorVertex=obj.fErrorVertex;
53e7090a 289 fX0=obj.fX0; fX1=obj.fX1; fX2=obj.fX2; fX3=obj.fX3; fX4=obj.fX4;
2a309a77 290 fLayer=obj.fLayer;
291 rtrack=obj.rtrack;
292 Dv=obj.Dv;
293 Zv=obj.Zv;
294 sigmaDv=obj.sigmaDv;
295 sigmaZv=obj.sigmaZv;
296 d2=obj.d2;
297 tgl2=obj.tgl2;
298 dtgl=obj.dtgl;
53e7090a 299
300 fC00=obj.fC00; fC10=obj.fC10; fC11=obj.fC11; fC20=obj.fC20; fC21=obj.fC21;
301 fC22=obj.fC22; fC30=obj.fC30; fC31=obj.fC31; fC32=obj.fC32; fC33=obj.fC33;
302 fC40=obj.fC40; fC41=obj.fC41; fC42=obj.fC42; fC43=obj.fC43; fC44=obj.fC44;
303
304
2a309a77 305 *ClusterInTrack = *obj.ClusterInTrack;
306 Int_t i;
307 for(i=0; i<obj.flistCluster->GetSize(); i++) flistCluster->AddLast(obj.flistCluster->At(i));
308
309 return *this;
310
311}
312
313void AliITStrack::PutCluster(Int_t layerc, TVector vecclust) {
314//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
315
316 (*ClusterInTrack)(layerc,0) = vecclust(0);
317 (*ClusterInTrack)(layerc,1) = vecclust(1);
318 (*ClusterInTrack)(layerc,2) = vecclust(2);
319 (*ClusterInTrack)(layerc,3) = vecclust(3);
320 (*ClusterInTrack)(layerc,4) = vecclust(4);
321 (*ClusterInTrack)(layerc,5) = vecclust(5);
322 (*ClusterInTrack)(layerc,6) = vecclust(6);
323 (*ClusterInTrack)(layerc,7) = vecclust(7);
324 (*ClusterInTrack)(layerc,8) = vecclust(8);
325
326}
327
328
329void AliITStrack::GetClusters() {
330
331 TMatrix A(*ClusterInTrack);
332 TMatrix B(6,3);
333 Int_t i;
334 for(i=0;i<6; i++){
335 B(i,0)=A(i,6); B(i,1)=A(i,7); B(i,2)=A(i,8);
336 }
337 A.Print();
338 // B.Print();
339
340}
341
342
343TVector AliITStrack::GetLabTrack(Int_t lay) {
344 TVector VecLabel(3);
345 VecLabel(0)=( (Float_t) (*ClusterInTrack)(lay,6) );
346 VecLabel(1)=( (Float_t) (*ClusterInTrack)(lay,7) );
347 VecLabel(2)=( (Float_t) (*ClusterInTrack)(lay,8) );
348 return VecLabel;
349}
350
351void AliITStrack::Search(TVector VecTotLabref, Long_t &labref, Int_t &freq){
352//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
353// define label
354
355 Int_t vecfreq[18];
356
357 Int_t i,j;
358 for(i=0; i<18; i++) vecfreq[i]=0;
359
360 for(i=0; i<18; i++) {
361 for(j=0; j<18; j++) {
53e7090a 362 // if(VecTotLabref(i) == 0.) VecTotLabref(i)=-3.; //commentato il 5-3-2001
2a309a77 363 if( (VecTotLabref(i)>=0.) && (VecTotLabref(i)==VecTotLabref(j)) ) vecfreq[i]++;
364 }
365 }
366 Int_t imax=-1000;
367 Long_t labdefault= (Long_t)1000000.;
368 freq=0;
369 for(i=0; i<18; i++) {
370 if(vecfreq[i]>freq) {freq=vecfreq[i]; imax=i;}
371 }
372 if(imax<0) labref=labdefault; else labref=(Long_t) VecTotLabref(imax);
373}
374
375
376void AliITStrack::Propagation(Double_t rk) {
377//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
378//Propagation of track
379 Double_t duepi=2.*TMath::Pi();
380 Double_t rkm1=rtrack;
381//cout<<" rk e rkm1 dentro Propagation "<<rk<<" "<<rkm1<<"\n";
382
383 //
384 Double_t Ak=argA(rk), Akm1=argA(rkm1);
385 Double_t ak=arga(rk), akm1=arga(rkm1);
53e7090a 386 fX0+=TMath::ASin(Ak)-TMath::ASin(Akm1);
2a309a77 387
53e7090a 388 if(fX0>duepi) fX0-=duepi;
389 if(fX0<0.) fX0+=duepi;
2a309a77 390
53e7090a 391 Double_t tgl=fX3;
392 Double_t C=fX4;
393 Double_t D=fX2;
2a309a77 394 Double_t Cy=C/2;
53e7090a 395 fX1+=tgl/Cy*(TMath::ASin(ak)-TMath::ASin(akm1));
2a309a77 396 rtrack=rk;
53e7090a 397
2a309a77 398
399 Double_t Bk=argB(rk), Bkm1=argB(rkm1);
400 Double_t Ck=argC(rk), Ckm1=argC(rkm1);
2a309a77 401
53e7090a 402 Double_t F02=Ck/TMath::Sqrt(1.-Ak*Ak) - Ckm1/TMath::Sqrt(1.-Akm1*Akm1);
403 Double_t F04=Bk/TMath::Sqrt(1.-Ak*Ak) - Bkm1/TMath::Sqrt(1.-Akm1*Akm1);
404 Double_t F12=tgl*D*(1./rk - 1./rkm1);
405 Double_t F13=rk - rkm1;
406
407
408 Double_t C00=fC00;
409 Double_t C10=fC10;
410 Double_t C11=fC11;
411 Double_t C20=fC20;
412 Double_t C21=fC21;
413 Double_t C22=fC22;
414 Double_t C30=fC30;
415 Double_t C31=fC31; // provare se si puo' fare a meno
416 Double_t C32=fC32;
417 Double_t C33=fC33;
418 Double_t C40=fC40;
419 Double_t C41=fC41;
420 Double_t C42=fC42;
421 Double_t C43=fC43;
422 Double_t C44=fC44;
423
424 Double_t R10=C10+C21*F02+C41*F04;
425 Double_t R20=C20+C22*F02+C42*F04;
426 Double_t R30=C30+C32*F02+C43*F04;
427 Double_t R40=C40+C42*F02+C44*F04;
428 Double_t R21=C21+C22*F12+C32*F13;
429 Double_t R31=C31+C32*F12+C33*F13;
430 Double_t R41=C41+C42*F12+fC43*F13;
431
432 fC00=C00+C20*F02+C40*F04+F02*R20+F04*R40;
433 fC10=R10+F12*R20+F13*R30;
434 fC11=C11+C21*F12+C31*F13+F12*R21+F13*R31;
435 fC20=R20;
436 fC21=R21;
437 fC30=R30;
438 fC31=R31;
439 fC40=R40;
440 fC41=R41;
2a309a77 441
2a309a77 442}
53e7090a 443
03454f1c 444void AliITStrack::AddEL(AliITSRad *rl, Double_t signdE, Bool_t flagtot, Double_t mass) {
8af13b4b 445//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
446// add energy loss
447
448 TVector s(6);
b3e91712 449 s(0)=0.0026+0.00283; s(1)=0.018; s(2)=0.0094; s(3)=0.0095; s(4)=0.0091; s(5)=0.0087;
8af13b4b 450 //0.00277 is added in the first layer to take into account the energy loss in the beam pipe
2a309a77 451
8af13b4b 452 //for(int k=0; k<6; k++) cout<<s(k)<<" "; cout<<"\n";
b3e91712 453 for(int k=0; k<6; k++) s(k)=s(k)*1.6;
454
8af13b4b 455
53e7090a 456 Double_t phi=fX0;
b3e91712 457
458 if(phi<0.174 ) s(5)=s(5)+0.012;
459 if(phi>6.1 ) s(5)=s(5)+0.012; // to take into account rail
460 if(phi>2.96 && phi<3.31 ) s(5)=s(5)+0.012;
461
462
53e7090a 463 Double_t tgl=fX3;
8af13b4b 464 Double_t theta=((TMath::Pi())/2.)-TMath::ATan(tgl);
465 //phi*=180./TMath::Pi();
466 //theta*=180./TMath::Pi();
467 //Double_t rad90=(TMath::Pi())/2.;
468 Double_t rad40=(TMath::Pi())*40./180.;
469 Double_t rad100=(TMath::Pi())*100/180;
470 Double_t rad360=(TMath::Pi())*2.;
471 Int_t imax=rl->Getimax();
472 Int_t jmax=rl->Getjmax();
473 Int_t i=(Int_t) ( (theta-rad40)/rad100*imax);
474 Int_t j=(Int_t) ( phi/rad360*jmax );
475 //Int_t i=(Int_t)( ((theta-((TMath::Pi())/4.))/((TMath::Pi())/2.))*imax );
476 //Int_t j=(Int_t)( (phi/((TMath::Pi())*2.))*jmax );
477 if(i<0) i=0;
478 if(i>=imax) i=imax-1;
479 if(j<0) j=0;
480 if(j>=jmax) j=jmax-1;
b3e91712 481 /*
8af13b4b 482 s(0) = 0.0028/TMath::Sin(theta)+( rl->GetRadMatrix1() )(i,j); // 0.0028 takes into account the beam pipe
483 s(1) = ( rl->GetRadMatrix2() )(i,j);
484 s(2) = ( rl->GetRadMatrix3() )(i,j);
485 s(3) = ( rl->GetRadMatrix4() )(i,j);
486 s(4) = ( rl->GetRadMatrix5() )(i,j);
487 s(5) = ( rl->GetRadMatrix6() )(i,j);
488
b3e91712 489
490 */
491
8af13b4b 492 //for(int k=0; k<6; k++) cout<<s(k)<<" "; getchar();
493
494 //if(phi>60) {cout<<" phi = "<<phi<<"\n"; getchar();}
495 //if(theta<45 || theta>135) {cout<<" theta = "<<theta<<"\n"; getchar();}
496 //cout<<" dentro AddEl: phi, theta = "<<phi<<" "<<theta<<"\n"; getchar();
497
53e7090a 498 Double_t cl=1.+fX3*fX3; // cl=1/(cosl)**2 = 1 + (tgl)**2
8af13b4b 499 Double_t sqcl=TMath::Sqrt(cl);
500 Double_t pt=GetPt();
501
502 Double_t p2=pt*pt*cl;
503 Double_t E=TMath::Sqrt(p2+mass*mass);
504 Double_t beta2=p2/(p2+mass*mass);
505
506 Double_t dE;
507 if(flagtot) {
508 Double_t stot=s(0)+s(1)+s(2)+s(3)+s(4)+s(5);
509 dE=0.153/beta2*(log(5940*beta2/(1-beta2)) - beta2)*stot*21.82*sqcl;
510 } else {
511 dE=0.153/beta2*(log(5940*beta2/(1-beta2)) - beta2)*s(fLayer-1)*21.82*sqcl;
512 }
513 dE=signdE*dE/1000.;
514
515 E+=dE;
516 Double_t p=TMath::Sqrt(E*E-mass*mass);
517 Double_t sign=1.;
53e7090a 518 if(fX4 < 0.) sign=-1.;
8af13b4b 519 pt=sign*p/sqcl;
8ab89103 520 Double_t CC=(0.299792458*0.2)/(pt*100.);
b3e91712 521 //Double_t CC=(0.299792458*0.5)/(pt*100.);
53e7090a 522 fX4=CC;
8af13b4b 523
8af13b4b 524}
53e7090a 525
2a309a77 526void AliITStrack::Correct(Double_t rk) {
527//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
528// correct track to take into account real geometry detector
529
530 Double_t duepi=2.*TMath::Pi();
531 Double_t rkm1=rtrack;
532 Double_t Ak=argA(rk), Akm1=argA(rkm1);
533 Double_t ak=arga(rk), akm1=arga(rkm1);
534
53e7090a 535 fX0+=TMath::ASin(Ak)-TMath::ASin(Akm1);
536 if(fX0>duepi) fX0-=duepi;
537 if(fX0<0.) fX0+=duepi;
2a309a77 538
53e7090a 539 Double_t tgl=fX3;
540 Double_t C=fX4;
2a309a77 541 Double_t Cy=C/2;
53e7090a 542 fX1+=tgl/Cy*(TMath::ASin(ak)-TMath::ASin(akm1));
2a309a77 543 rtrack=rk;
544
545}
2a309a77 546
8af13b4b 547void AliITStrack::AddMS(AliITSRad *rl) {
548
549////////// Modification of the covariance matrix to take into account multiple scattering ///////////
550//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
551
552 TVector s(6);
553
b3e91712 554 s(0)=0.0026+0.00283; s(1)=0.018; s(2)=0.0094; s(3)=0.0095; s(4)=0.0091; s(5)=0.0087;
8af13b4b 555//0.00277 is added in the first layer to take into account the energy loss in the beam pipe
556
b3e91712 557 for(int k=0; k<6; k++) s(k)=s(k)*1.6;
558
53e7090a 559 Double_t phi=fX0;
b3e91712 560 if(phi<0.174 ) s(5)=s(5)+0.012; //aggiunta provvisoria
561 if(phi>6.1 ) s(5)=s(5)+0.012; //aggiunta provvisoria
562 if(phi>2.96 && phi< 3.31) s(5)=s(5)+0.012; //aggiunta provvisoria
563
53e7090a 564 Double_t tgl=fX3;
8af13b4b 565 Double_t theta=((TMath::Pi())/2.)-TMath::ATan(tgl);
53e7090a 566 Double_t rad40=(TMath::Pi())*40./180.; // rivedere
8af13b4b 567 Double_t rad100=(TMath::Pi())*100/180;
568 Double_t rad360=(TMath::Pi())*2.;
569 Int_t imax=rl->Getimax();
570 Int_t jmax=rl->Getjmax();
571 Int_t i=(Int_t) ( (theta-rad40)/rad100*imax);
572 Int_t j=(Int_t) ( phi/rad360*jmax);
8af13b4b 573
574 if(i<0) i=0;
575 if(i>=imax) i=imax-1;
576 if(j<0) j=0;
577 if(j>=jmax) j=jmax-1;
b3e91712 578 /*
8af13b4b 579 s(0) = 0.0028/TMath::Sin(theta)+( rl->GetRadMatrix1() )(i,j); // 0.0028 takes into account the beam pipe
580 s(1) = ( rl->GetRadMatrix2() )(i,j);
581 s(2) = ( rl->GetRadMatrix3() )(i,j);
582 s(3) = ( rl->GetRadMatrix4() )(i,j);
583 s(4) = ( rl->GetRadMatrix5() )(i,j);
584 s(5) = ( rl->GetRadMatrix6() )(i,j);
b3e91712 585 */
8af13b4b 586 Double_t mass=0.1396;
587 Int_t layer=(Int_t)GetLayer();
588
8af13b4b 589 Double_t cosl=TMath::Cos(TMath::ATan(tgl));
53e7090a 590 Double_t D=fX2;
591 Double_t C=fX4;
8af13b4b 592 Double_t Cy=C/2.;
593 Double_t Q20=1./(cosl*cosl);
594 Double_t Q30=C*tgl;
595
596 Double_t Q40=Cy*(rtrack*rtrack-D*D)/(1.+ 2.*Cy*D);
597 Double_t dd=D+Cy*D*D-Cy*rtrack*rtrack;
53e7090a 598 Double_t dprova=rtrack*rtrack - dd*dd;
599 Double_t Q41=0.;
600 if(dprova>0.) Q41=-1./cosl*TMath::Sqrt(dprova)/(1.+ 2.*Cy*D);
601
2a309a77 602 Double_t p2=(GetPt()*GetPt())/(cosl*cosl);
603 Double_t beta2=p2/(p2+mass*mass);
8af13b4b 604 Double_t theta2=14.1*14.1/(beta2*p2*1.e6)*(s(layer-1)/cosl);
53e7090a 605
606 fC22+=theta2*(Q40*Q40+Q41*Q41);
607 fC32+=theta2*Q20*Q40;
608 fC33+=theta2*Q20*Q20;
609 fC42+=theta2*Q30*Q40;
610 fC43+=theta2*Q30*Q20;
611 fC44+=theta2*Q30*Q30;
612
2a309a77 613}
8af13b4b 614void AliITStrack::PrimaryTrack(AliITSRad *rl) {
2a309a77 615//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
616// calculation of part of covariance matrix for vertex constraint
617
618 Double_t Rlayer[6];
619
620 Rlayer[0]=4.; Rlayer[1]=7.; Rlayer[2]=14.9; Rlayer[3]=23.8;
621 Rlayer[4]=39.1; Rlayer[5]=43.6;
622
53e7090a 623 Double_t Cy=fX4/2.;
624 Double_t tgl=(fX1-Zv)*Cy/TMath::ASin(Cy*rtrack);
2a309a77 625 Double_t rtrack=1.;
53e7090a 626 fX0=0.;
627 fX1=rtrack*tgl;
628 fX2=Dv;
629 fX3=tgl;
630
631 fC00=fC10=fC11=fC20=fC21=fC22=fC30=fC31=fC32=fC33=fC40=fC41=fC42=fC43=0.;
632
8af13b4b 633 AddEL(rl,1.,1);
2a309a77 634 fLayer=0;
635 Int_t i;
636 for (i=0; i<6; i++) {
637 Propagation(Rlayer[i]);
638 fLayer++;
53e7090a 639 d2(i)=fC22;
640 tgl2(i)=fC33;
641 dtgl(i)=fC32;
8af13b4b 642 AddMS(rl);
643 AddEL(rl,-1,0);
2a309a77 644 }
645}
646
647
648Int_t AliITStrack::DoNotCross(Double_t rk) const{
53e7090a 649 Double_t C=fX4;
650 Double_t D=fX2;
2a309a77 651 Double_t Cy=C/2.;
652 return (TMath::Abs((Cy*rk+(1.+Cy*D)*D/rk)/(1.+2.*Cy*D))>=1.)?1:0;
653}
654
655
656Double_t AliITStrack::argA(Double_t rk) const {
53e7090a 657 Double_t C=fX4;
658 Double_t D=fX2;
2a309a77 659 Double_t Cy=C/2.;
660 Double_t arg=(Cy*rk + (1 + Cy*D)*D/rk)/(1.+ 2.*Cy*D);
661 if (TMath::Abs(arg) < 1.) return arg;
662 //cout<<"class AliITSTrack: argA out of range !\n";/* getchar();*/
663 return (arg>0) ? 0.99999999999 : -0.9999999999;
664}
665
666Double_t AliITStrack::arga(Double_t rk) const {
53e7090a 667 Double_t C=fX4;
668 Double_t D=fX2;
2a309a77 669 Double_t Cy=C/2.;
670 Double_t arg=(rk*rk - D*D)/(1.+ 2.*Cy*D);
671 if (arg<0.) {/*cout<<"class AliITSTrack: arga out of range !\n";*/ arg=0.;}
672 return Cy*TMath::Sqrt(arg);
673}
674
675Double_t AliITStrack::argB(Double_t rk) const {
53e7090a 676 Double_t C=fX4;
677 Double_t D=fX2;
2a309a77 678 Double_t Cy=C/2.;
679 return (rk*rk - D*D)/(rk*(1.+ 2.*Cy*D)*(1.+ 2.*Cy*D));
680}
681
682Double_t AliITStrack::argC(Double_t rk) const {
53e7090a 683 Double_t C=fX4;
684 Double_t D=fX2;
2a309a77 685 Double_t Cy=C/2.;
686 return (1./rk - 2.*Cy*argA(rk)/(1.+ 2.*Cy*D));
687}
53e7090a 688
2a309a77 689Double_t AliITStrack::PhiDef(Double_t x, Double_t y){
690 Double_t pigre= TMath::Pi();
53e7090a 691 Double_t phi=10000.;
2a309a77 692 if(y == 0. || x == 0.) {
693 if(y == 0. && x == 0.) {
694 cout << " Error in AliITStracking::PhiDef x=0 and y=0 \n"; getchar();
695 }
696 if(y==0. && x>0.) phi=0.;
697 if(y==0. && x<0.) phi=pigre;
698 if(x==0 && y>0.) phi=pigre/2.;
699 if(x==0 && y<0.) phi=1.5*pigre;
700 }
701 else {
702 if (x>0. && y>0.) phi=TMath::ATan(y/x);
703 if (x<0. && y>0.) phi=pigre+TMath::ATan(y/x);
704 if (x<0. && y<0.) phi=pigre+TMath::ATan(y/x);
705 if (x>0. && y<0.) phi=(2.*pigre)+TMath::ATan(y/x);
706 }
707 if(phi<0. || phi>(2*pigre)) {
708 cout<<" Error on phi in AliITStracking::PhiDef \n"; getchar();
709 }
710 return phi;
711}
53e7090a 712