]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITStrackMI.cxx
Adding MUON HLT code to the repository.
[u/mrichter/AliRoot.git] / ITS / AliITStrackMI.cxx
CommitLineData
15dd636f 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//-------------------------------------------------------------------------
17// Implementation of the ITS track class
18//
19// Origin: Marian Ivanov, CERN, Marian.Ivanov@cern.ch
20// dEdx analysis by: Boris Batyunya, JINR, Boris.Batiounia@cern.ch
21//-------------------------------------------------------------------------
22
23#include <TMatrixD.h>
24
25#include <TMath.h>
26
27#include "AliCluster.h"
28#include "AliESDtrack.h"
29#include "AliITStrackMI.h"
30
31ClassImp(AliITStrackMI)
32
33const Int_t kWARN=5;
34
35//____________________________________________________________________________
36AliITStrackMI::AliITStrackMI():AliITStrackV2(),
37 fNUsed(0),
38 fNSkipped(0),
39 fNDeadZone(0),
40 fDeadZoneProbability(0),
41 fReconstructed(kFALSE),
42 fConstrain(kFALSE)
43{
44 for(Int_t i=0; i<kMaxLayer; i++) fClIndex[i]=-1;
45 for(Int_t i=0; i<6; i++) { fNy[i]=0; fNz[i]=0; fNormQ[i]=0; fNormChi2[i]=1000;}
46 for(Int_t i=0; i<12; i++) {fDy[i]=0; fDz[i]=0; fSigmaY[i]=0; fSigmaZ[i]=0; fChi2MIP[i]=0;}
47 fD[0]=0; fD[1]=0;
48 fExpQ=40;
49 fdEdxMismatch=0;
50 fChi22=0;
628e7bb0 51 fGoldV0 = kFALSE;
15dd636f 52}
53
54//____________________________________________________________________________
55AliITStrackMI::AliITStrackMI(AliESDtrack& t,Bool_t c) throw (const Char_t *) :
56AliITStrackV2(t,c) {
57 //------------------------------------------------------------------
58 // Conversion ESD track -> ITS track.
59 // If c==kTRUE, create the ITS track out of the constrained params.
60 //------------------------------------------------------------------
61 fNUsed = 0;
62 fReconstructed = kFALSE;
63 fNSkipped =0;
64 fNDeadZone = 0;
65 fDeadZoneProbability = 0;
66 for(Int_t i=0; i<6; i++) {fClIndex[i]=-1; fNy[i]=0; fNz[i]=0; fNormQ[i]=0; fNormChi2[i]=1000;}
67 for(Int_t i=0; i<12; i++) {fDy[i]=0; fDz[i]=0; fSigmaY[i]=0; fSigmaZ[i]=0;fChi2MIP[i]=0;}
68 fD[0]=0; fD[1]=0;
69 fExpQ=40;
70 fConstrain = kFALSE;
71 fdEdxMismatch=0;
72 fChi22 =0;
628e7bb0 73 fGoldV0 = kFALSE;
15dd636f 74 //if (!Invariant()) throw "AliITStrackV2: conversion failed !\n";
75
76}
77
78void AliITStrackMI::UpdateESDtrack(ULong_t flags) {
79 fESDtrack->UpdateTrackParams(this,flags);
80 if (flags == AliESDtrack::kITSin) fESDtrack->SetITSChi2MIP(fChi2MIP);
81}
82
83//____________________________________________________________________________
84AliITStrackMI::AliITStrackMI(const AliITStrackMI& t) : AliITStrackV2(t) {
85 //------------------------------------------------------------------
86 //Copy constructor
87 //------------------------------------------------------------------
88 fNUsed = t.fNUsed;
89 fReconstructed = t.fReconstructed;
90 fNSkipped = t.fNSkipped;
91 fNDeadZone = t.fNDeadZone;
92 fDeadZoneProbability = t.fDeadZoneProbability;
93 fLab = t.fLab;
94 fFakeRatio = t.fFakeRatio;
95 fdEdxMismatch = t.fdEdxMismatch;
96 fChi22 = t.fChi22;
628e7bb0 97 fGoldV0 = t.fGoldV0;;
15dd636f 98
99 fD[0]=t.fD[0]; fD[1]=t.fD[1];
628e7bb0 100 fDnorm[0] = t.fDnorm[0]; fDnorm[1]=t.fDnorm[1];
15dd636f 101 fExpQ= t.fExpQ;
102 for(Int_t i=0; i<6; i++) {
103 fClIndex[i]= t.fClIndex[i]; fNy[i]=t.fNy[i]; fNz[i]=t.fNz[i]; fNormQ[i]=t.fNormQ[i]; fNormChi2[i] = t.fNormChi2[i];
104 }
105 for(Int_t i=0; i<12; i++) {fDy[i]=t.fDy[i]; fDz[i]=t.fDz[i];
106 fSigmaY[i]=t.fSigmaY[i]; fSigmaZ[i]=t.fSigmaZ[i];fChi2MIP[i]=t.fChi2MIP[i];}
107 fConstrain = t.fConstrain;
108 //memcpy(fDy,t.fDy,6*sizeof(Float_t));
109 //memcpy(fDz,t.fDz,6*sizeof(Float_t));
110 //memcpy(fSigmaY,t.fSigmaY,6*sizeof(Float_t));
111 //memcpy(fSigmaZ,t.fSigmaZ,6*sizeof(Float_t));
112 //memcpy(fChi2MIP,t.fChi2MIP,12*sizeof(Float_t));
113}
114
115//_____________________________________________________________________________
116Int_t AliITStrackMI::Compare(const TObject *o) const {
117 //-----------------------------------------------------------------
118 // This function compares tracks according to the their curvature
119 //-----------------------------------------------------------------
120 AliITStrackMI *t=(AliITStrackMI*)o;
121 //Double_t co=TMath::Abs(t->Get1Pt());
122 //Double_t c =TMath::Abs(Get1Pt());
123 Double_t co=t->GetSigmaY2()*t->GetSigmaZ2()*(0.5+TMath::Sqrt(0.5*t->fD[0]*t->fD[0]+t->fD[1]*t->fD[1]));
124 Double_t c =GetSigmaY2()*GetSigmaZ2()*(0.5+TMath::Sqrt(0.5*fD[0]*fD[0]+fD[1]*fD[1]));
125 if (c>co) return 1;
126 else if (c<co) return -1;
127 return 0;
128}
129
130
131Double_t AliITStrackMI::GetPredictedChi2MI(Double_t cy, Double_t cz, Double_t cerry, Double_t cerrz) const
132{
133 //-----------------------------------------------------------------
134 // This function calculates a predicted chi2 increment.
135 //-----------------------------------------------------------------
136 Double_t r00=cerry*cerry, r01=0., r11=cerrz*cerrz;
137 r00+=fC00; r01+=fC10; r11+=fC11;
138 //
139 Double_t det=r00*r11 - r01*r01;
140 if (TMath::Abs(det) < 1.e-30) {
141 Int_t n=GetNumberOfClusters();
142 if (n>kWARN)
143 Warning("GetPredictedChi2","Singular matrix (%d) !\n",n);
144 return 1e10;
145 }
146 Double_t tmp=r00; r00=r11; r11=tmp; r01=-r01;
147
148 Double_t dy=cy - fP0, dz=cz - fP1;
149
150 return (dy*r00*dy + 2*r01*dy*dz + dz*r11*dz)/det;
151}
152
153//____________________________________________________________________________
154Int_t AliITStrackMI::CorrectForMaterial(Double_t d, Double_t x0) {
155 //------------------------------------------------------------------
156 //This function corrects the track parameters for crossed material
157 //------------------------------------------------------------------
158 // Double_t p2=(1.+ GetTgl()*GetTgl())/(Get1Pt()*Get1Pt());
159 Double_t p2=(1.+ fP3*fP3)/(Get1Pt()*Get1Pt());
160 Double_t et = p2 + GetMass()*GetMass();
161 Double_t beta2=p2/et;
162 et = sqrt(et);
163 d*=TMath::Sqrt((1.+ fP3*fP3)/(1.- fP2*fP2));
164 //d*=TMath::Sqrt(1.+ fP3*fP3 +fP2*fP2/(1.- fP2*fP2));
165
166 //Multiple scattering******************
167 if (d!=0) {
168 Double_t theta2=14.1*14.1/(beta2*p2*1e6)*TMath::Abs(d);
169 //Double_t theta2=1.0259e-6*14*14/28/(beta2*p2)*TMath::Abs(d)*9.36*2.33;
170 fC22 += theta2*(1.- fP2*fP2)*(1. + fP3*fP3);
171 fC33 += theta2*(1. + fP3*fP3)*(1. + fP3*fP3);
172 fC43 += theta2*fP3*fP4*(1. + fP3*fP3);
173 fC44 += theta2*fP3*fP4*fP3*fP4;
174 }
175
176 //Energy losses************************
177 if (x0!=0.) {
178 d*=x0;
179 // Double_t dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2)) - beta2)*d;
180 //Double_t dE=0;
181 Double_t dE = 0.265*0.153e-3*(39.2-55.6*beta2+28.7*beta2*beta2+27.41/beta2)*d;
182 /*
183 if (beta2/(1-beta2)>3.5*3.5){
184 dE=0.153e-3/beta2*(log(3.5*5940)+0.5*log(beta2/(1-beta2)) - beta2)*d;
185 }
186 else{
187 dE=0.153e-3/beta2*(log(5940*beta2/(1-beta2)) - beta2)*d;
188 dE+=0.06e-3/(beta2*beta2)*d;
189 }
190 */
191 fP4*=(1.- et/p2*dE);
192 Double_t delta44 = (dE*fP4*et/p2);
193 delta44*=delta44;
194 fC44+= delta44/400.;
195 }
196
197 if (!Invariant()) return 0;
198
199 return 1;
200}
201
202//____________________________________________________________________________
203Int_t AliITStrackMI::UpdateMI(Double_t cy, Double_t cz, Double_t cerry, Double_t cerrz, Double_t chi2,UInt_t index) {
204 //------------------------------------------------------------------
205 //This function updates track parameters
206 //------------------------------------------------------------------
207 Double_t p0=fP0,p1=fP1,p2=fP2,p3=fP3,p4=fP4;
208 Double_t c00=fC00;
209 Double_t c10=fC10, c11=fC11;
210 Double_t c20=fC20, c21=fC21, c22=fC22;
211 Double_t c30=fC30, c31=fC31, c32=fC32, c33=fC33;
212 Double_t c40=fC40, c41=fC41, c42=fC42, c43=fC43, c44=fC44;
213
214
215 Double_t r00=cerry*cerry, r01=0., r11=cerrz*cerrz;
216 r00+=fC00; r01+=fC10; r11+=fC11;
217 Double_t det=r00*r11 - r01*r01;
218 Double_t tmp=r00; r00=r11/det; r11=tmp/det; r01=-r01/det;
219
220
221 Double_t k00=fC00*r00+fC10*r01, k01=fC00*r01+fC10*r11;
222 Double_t k10=fC10*r00+fC11*r01, k11=fC10*r01+fC11*r11;
223 Double_t k20=fC20*r00+fC21*r01, k21=fC20*r01+fC21*r11;
224 Double_t k30=fC30*r00+fC31*r01, k31=fC30*r01+fC31*r11;
225 Double_t k40=fC40*r00+fC41*r01, k41=fC40*r01+fC41*r11;
226
227 Double_t dy=cy - fP0, dz=cz - fP1;
228 Int_t layer = (index & 0xf0000000) >> 28;
229 fDy[layer] = dy;
230 fDz[layer] = dz;
231 fSigmaY[layer] = TMath::Sqrt(cerry*cerry+fC00);
232 fSigmaZ[layer] = TMath::Sqrt(cerrz*cerrz+fC11);
233
234 Double_t sf=fP2 + k20*dy + k21*dz;
235
236 fP0 += k00*dy + k01*dz;
237 fP1 += k10*dy + k11*dz;
238 fP2 = sf;
239 fP3 += k30*dy + k31*dz;
240 fP4 += k40*dy + k41*dz;
241
242 Double_t c01=fC10, c02=fC20, c03=fC30, c04=fC40;
243 Double_t c12=fC21, c13=fC31, c14=fC41;
244
245 fC00-=k00*fC00+k01*fC10; fC10-=k00*c01+k01*fC11;
246 fC20-=k00*c02+k01*c12; fC30-=k00*c03+k01*c13;
247 fC40-=k00*c04+k01*c14;
248
249 fC11-=k10*c01+k11*fC11;
250 fC21-=k10*c02+k11*c12; fC31-=k10*c03+k11*c13;
251 fC41-=k10*c04+k11*c14;
252
253 fC22-=k20*c02+k21*c12; fC32-=k20*c03+k21*c13;
254 fC42-=k20*c04+k21*c14;
255
256 fC33-=k30*c03+k31*c13;
257 fC43-=k30*c04+k31*c14;
258
259 fC44-=k40*c04+k41*c14;
260
261 if (!Invariant()) {
262 fP0=p0; fP1=p1; fP2=p2; fP3=p3; fP4=p4;
263 fC00=c00;
264 fC10=c10; fC11=c11;
265 fC20=c20; fC21=c21; fC22=c22;
266 fC30=c30; fC31=c31; fC32=c32; fC33=c33;
267 fC40=c40; fC41=c41; fC42=c42; fC43=c43; fC44=c44;
268 return 0;
269 }
270
271 if (chi2<0) return 1;
272 Int_t n=GetNumberOfClusters();
273 fIndex[n]=index;
274 SetNumberOfClusters(n+1);
275 SetChi2(GetChi2()+chi2);
276
277 return 1;
278}
279
280Int_t AliITStrackMI::GetProlongationFast(Double_t alp, Double_t xk,Double_t &y, Double_t &z)
281{
282 //-----------------------------------------------------------------------------
283 //get fast prolongation
284 //-----------------------------------------------------------------------------
285 Double_t ca=TMath::Cos(alp-fAlpha), sa=TMath::Sin(alp-fAlpha);
286 Double_t cf=TMath::Sqrt(1.- fP2*fP2);
287 // **** rotation **********************
288 y= -fX*sa + fP0*ca;
289 // **** translation ******************
290 Double_t dx = xk- fX*ca - fP0*sa;
291 Double_t f1=fP2*ca - cf*sa, f2=f1 + fP4*dx;
292 if (TMath::Abs(f2) >= 0.9999) {
293 return 0;
294 }
295 Double_t r1=TMath::Sqrt(1.- f1*f1), r2=TMath::Sqrt(1.- f2*f2);
296 y += dx*(f1+f2)/(r1+r2);
297 z = fP1+dx*(f1+f2)/(f1*r2 + f2*r1)*fP3;
298 return 1;
299}
300
628e7bb0 301
302Bool_t AliITStrackMI::IsGoldPrimary()
303{
304 //
305 // Indicates gold pimary track
306 //
307 Bool_t isGold=kTRUE;
308 if (!fConstrain) return kFALSE; //
309 if (fNDeadZone+fNDeadZone<5.5) isGold = kFALSE; // short track
310 //
311 if (fChi2/Float_t(fN)>2.){
312 if (fChi2MIP[0]+fNUsed>3.5) isGold = kFALSE;
313 }
314 if (fChi2MIP[2]>4.5) isGold = kFALSE; //back propagation chi2
315 //
316 if (fDnorm[0]>0&&fDnorm[1]>0){
317 const Float_t distcut2 =2.5*2.5; //normalize distance cut
318 Float_t dist2 = fD[0]*fD[0]/(fDnorm[0]*fDnorm[0])+fD[1]*fD[1]/(fDnorm[1]*fDnorm[1]); //normalize distance to the vertex (pools)
319 if (dist2>distcut2) isGold = kFALSE;
320 }
321 return isGold;
322}