]>
Commit | Line | Data |
---|---|---|
6eb67451 | 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. * | |
70521312 | 14 | ************************************************************************** |
15 | ************************************************************************** | |
6eb67451 | 16 | * * |
17 | * This class builds AliTPCtrack objects from generated tracks to feed * | |
18 | * ITS tracking (V2). The AliTPCtrack is built from its first hit in * | |
19 | * the TPC. The track is assigned a Kalman-like covariance matrix * | |
20 | * depending on its pT and pseudorapidity and track parameters are * | |
70521312 | 21 | * smeared according to this covariance matrix. * |
6eb67451 | 22 | * Output file contains sorted tracks, ready for matching with ITS * |
23 | * * | |
24 | * Test macro is: AliBarrelRec_TPCparam.C * | |
25 | * * | |
26 | * Origin: Andrea Dainese, Padova - e-mail: andrea.dainese@pd.infn.it * | |
27 | * * | |
28 | **************************************************************************/ | |
29 | #include "AliTPCtrackerParam.h" | |
30 | #include "alles.h" | |
31 | #include "AliMagF.h" | |
32 | #include "AliTPCtrack.h" | |
33 | #include "TMatrixD.h" | |
34 | #include "AliKalmanTrack.h" | |
70521312 | 35 | #include "AliMagFCM.h" |
36 | #include "AliGausCorr.h" | |
37 | ||
38 | ||
39 | ||
6eb67451 | 40 | |
41 | ClassImp(AliTPCtrackerParam) | |
42 | ||
43 | //----------------------------------------------------------------- | |
44 | AliTPCtrackerParam::AliTPCtrackerParam(const Int_t coll,const Double_t Bz) | |
45 | { | |
46 | fColl = coll; // collision code (0: PbPb6000) | |
47 | fBz = Bz; // value of the z component of L3 field (Tesla) | |
48 | ||
49 | } | |
50 | //----------------------------------------------------------------- | |
51 | AliTPCtrackerParam::~AliTPCtrackerParam() | |
52 | {} | |
53 | //----------------------------------------------------------------- | |
54 | ||
55 | Int_t AliTPCtrackerParam::BuildTPCtracks(const TFile *inp, TFile *out, Int_t n) | |
56 | { | |
57 | ||
58 | if(fColl!=0) { | |
f3d60f45 | 59 | cerr<<"AliTPCtrackerParam::BuildTPCtracks: Invalid collision!\n"; |
60 | cerr<<" Available: 0 -> PbPb6000"<<endl; return 0; | |
6eb67451 | 61 | } |
62 | if(fBz!=0.4) { | |
f3d60f45 | 63 | cerr<<"AliTPCtrackerParam::BuildTPCtracks: Invalid field!\n"; |
64 | cerr<<" Available: 0.4"<<endl; return 0; | |
6eb67451 | 65 | } |
66 | ||
67 | TFile *infile=(TFile*)inp; | |
68 | ||
6eb67451 | 69 | // Get gAlice object from file |
70 | if(!(gAlice=(AliRun*)infile->Get("gAlice"))) { | |
71 | cerr<<"gAlice has not been found on galice.root !\n"; | |
72 | return 1; | |
73 | } | |
74 | ||
70521312 | 75 | AliMagFCM *fiel = (AliMagFCM*)gAlice->Field(); |
76 | Double_t fieval=(Double_t)fiel->SolenoidField()/10.; | |
77 | printf("Magnetic field is %6.2f Tesla\n",fieval); | |
78 | if(fBz!=fieval) { | |
79 | cerr<<"AliTPCtrackerParam::BuildTPCtracks: Invalid field!"<<endl; | |
80 | cerr<<"Field selected is: "<<fBz<<" T\n"; | |
81 | cerr<<"Field found on file is: "<<fieval<<" T\n"; | |
82 | return 0; | |
83 | } | |
84 | ||
85 | AliKalmanTrack::SetConvConst(100/0.299792458/fBz); | |
6eb67451 | 86 | |
87 | ||
88 | // loop over first n events in file | |
89 | for(Int_t evt=0; evt<n; evt++){ | |
90 | cerr<<"+++\n+++ Processing event "<<evt<<"\n+++\n"; | |
91 | ||
92 | AliTPCtrack *tpctrack=0; | |
93 | ||
94 | // tree for TPC tracks | |
95 | Char_t tname[100]; | |
96 | sprintf(tname,"TreeT_TPC_%d",evt); | |
97 | TTree *tracktree = new TTree(tname,"Tree with TPC tracks"); | |
98 | tracktree->Branch("tracks","AliTPCtrack",&tpctrack,20000,0); | |
99 | ||
100 | // array for TPC tracks | |
101 | TObjArray tarray(20000); | |
102 | ||
103 | ||
104 | Int_t nparticles=gAlice->GetEvent(evt); | |
105 | ||
106 | Bool_t done[500000]; | |
107 | for(Int_t l=0; l<500000; l++) { done[l]=kFALSE; } | |
108 | ||
109 | ||
110 | // Get TPC detector | |
111 | AliTPC *TPC=(AliTPC*)gAlice->GetDetector("TPC"); | |
112 | Int_t ver = TPC->IsVersion(); | |
113 | cerr<<"+++ TPC version "<<ver<<" has been found !\n"; | |
114 | AliTPCParam *digp=(AliTPCParam*)infile->Get("75x40_100x60"); | |
115 | if(!digp) { cerr<<"TPC parameters have not been found !\n"; return 1; } | |
116 | TPC->SetParam(digp); | |
117 | ||
118 | // Get TreeH with hits | |
119 | TTree *TH=gAlice->TreeH(); | |
120 | Int_t ntracks=(Int_t)TH->GetEntries(); | |
121 | cerr<<"+++\n+++ Number of particles in event "<<evt<<": "<<nparticles<<"\n+++\n+++ Number of \"primary tracks\" (entries in TreeH): "<<ntracks<<"\n+++\n\n"; | |
122 | ||
123 | TParticle* Part; | |
124 | AliTPChit* tpcHit; | |
125 | Double_t hPx,hPy,hPz,hPt,xg,yg,zg,xl,yl,zl; | |
126 | Double_t Alpha; | |
127 | Float_t cosAlpha,sinAlpha; | |
128 | Int_t label,Pdg,charge,bin; | |
129 | Int_t tracks=0; | |
130 | //Int_t N1=0,N2=0; | |
131 | ||
132 | // loop over entries in TreeH | |
133 | for(Int_t i=0; i<ntracks; i++) { | |
134 | if(i%1000==0) cerr<<" --- Processing primary track "<<i<<" of "<<ntracks<<" ---\r"; | |
135 | TPC->ResetHits(); | |
136 | TH->GetEvent(i); | |
137 | // Get FirstHit | |
138 | tpcHit=(AliTPChit*)TPC->FirstHit(-1); | |
139 | for( ; tpcHit; tpcHit=(AliTPChit*)TPC->NextHit() ) { | |
140 | if(tpcHit->fQ !=0.) continue; | |
141 | // Get particle momentum at hit | |
142 | hPx=tpcHit->X(); hPy=tpcHit->Y(); hPz=tpcHit->Z(); | |
143 | hPt=TMath::Sqrt(hPx*hPx+hPy*hPy); | |
144 | // reject hits with Pt<mag*0.45 GeV/c | |
145 | if(hPt<(fBz*0.45)) continue; | |
146 | ||
147 | // Get track label | |
148 | label=tpcHit->Track(); | |
149 | // check if this track has already been processed | |
150 | if(done[label]) continue; | |
151 | // electric charge | |
152 | Part = gAlice->Particle(label); | |
153 | Pdg = Part->GetPdgCode(); | |
154 | if(Pdg>200 || Pdg==-11 || Pdg==-13) { charge=1; } | |
155 | else if(Pdg<-200 || Pdg==11 || Pdg==13) { charge=-1; } | |
156 | else continue; | |
157 | ||
158 | ||
159 | if((tpcHit=(AliTPChit*)TPC->NextHit())==0) break; | |
160 | if(tpcHit->fQ != 0.) continue; | |
161 | // Get global coordinates of hit | |
162 | xg=tpcHit->X(); yg=tpcHit->Y(); zg=tpcHit->Z(); | |
163 | if(TMath::Sqrt(xg*xg+yg*yg)>90.) continue; | |
164 | ||
165 | // Get TPC sector, Alpha angle and local coordinates | |
166 | // printf("Sector %d\n",tpcHit->fSector); | |
167 | digp->AdjustCosSin(tpcHit->fSector,cosAlpha,sinAlpha); | |
168 | Alpha = TMath::ATan2(sinAlpha,cosAlpha); | |
169 | xl = xg*cosAlpha + yg*sinAlpha; | |
170 | yl =-xg*sinAlpha + yg*cosAlpha; | |
171 | zl = zg; | |
172 | //printf("Alpha %f xl %f yl %f zl %f\n",Alpha,xl,yl,zl); | |
173 | ||
174 | // reject tracks which are not in the TPC acceptance | |
175 | if(TMath::Abs(zl+(244.-xl)*hPz/hPt)>252.) continue; | |
176 | ||
177 | // Get bin in pT,eta | |
178 | bin = GetBin(hPt,Part->Eta()); | |
179 | ||
180 | // Apply selection according to TPC efficiency | |
181 | //if(abs(Pdg)==211) N1++; | |
182 | if(!SelectedTrack(Pdg,hPt,Part->Eta())) continue; | |
183 | //if(abs(Pdg)==211) N2++; | |
184 | ||
185 | // Mark track as "done" | |
186 | done[label]=kTRUE; tracks++; | |
187 | ||
188 | // create AliTPCtrack object | |
189 | tpctrack = BuildTrack(Alpha,xl,yl,zl,hPx,hPy,hPz,hPt,charge,label); | |
190 | ||
191 | // put track in array | |
192 | tarray.AddLast(tpctrack); | |
193 | ||
194 | } | |
195 | ||
196 | } // loop over entries in TreeH | |
197 | ||
198 | TObjArray newtarray(20000); | |
199 | ||
200 | // assing covariance matrixes and smear track parameters | |
201 | CookTracks(tarray,newtarray); | |
202 | ||
203 | // sort array with TPC tracks (decreasing pT) | |
204 | newtarray.Sort(); | |
205 | ||
206 | ||
207 | Int_t arrentr = newtarray.GetEntriesFast(); | |
208 | //printf("\n %d \n\n",arrentr); | |
209 | for(Int_t l=0; l<arrentr; l++) { | |
210 | tpctrack=(AliTPCtrack*)newtarray.UncheckedAt(l); | |
211 | tracktree->Fill(); | |
212 | } | |
213 | ||
214 | // write the tree with tracks in the output file | |
215 | out->cd(); | |
216 | tracktree->Write(); | |
217 | ||
218 | delete tracktree; | |
219 | ||
220 | printf("\n\n+++\n+++ Number of TPC tracks: %d\n+++\n",tracks); | |
221 | //printf("Average Eff: %f\n",(Float_t)N2/N1); | |
222 | ||
223 | } // loop on events | |
224 | ||
225 | return 0; | |
226 | } | |
227 | ||
228 | //----------------------------------------------------------------- | |
229 | AliTPCtrack* AliTPCtrackerParam::BuildTrack(Double_t alpha,Double_t x, | |
230 | Double_t y,Double_t z,Double_t Px, | |
231 | Double_t Py,Double_t Pz,Double_t Pt, | |
232 | Int_t ch,Int_t lab) | |
233 | { | |
234 | Double_t xref = x; | |
235 | Double_t xx[5],cc[15]; | |
236 | cc[0]=cc[2]=cc[5]=cc[9]=cc[14]=10.; | |
237 | cc[1]=cc[3]=cc[4]=cc[6]=cc[7]=cc[8]=cc[10]=cc[11]=cc[12]=cc[13]=0.; | |
238 | ||
239 | // Magnetic field | |
240 | TVector3 Bfield(0.,0.,fBz); | |
241 | ||
242 | ||
243 | // radius [cm] of track projection in (x,y) | |
244 | Double_t rho = Pt*100/0.299792458/Bfield.Z(); | |
245 | // center of track projection in local reference frame | |
246 | TVector3 MomH,PosH; | |
247 | ||
248 | ||
249 | // position (local) and momentum (local) at the hit | |
250 | // in the bending plane (z=0) | |
251 | PosH.SetXYZ(x,y,0.); | |
252 | MomH.SetXYZ(Px*TMath::Cos(alpha)+Py*TMath::Sin(alpha),-Px*TMath::Sin(alpha)+Py*TMath::Cos(alpha),0.); | |
253 | TVector3 Vrho = MomH.Cross(Bfield); | |
254 | Vrho *= ch; | |
255 | Vrho.SetMag(rho); | |
256 | ||
257 | TVector3 Vcenter = PosH+Vrho; | |
258 | ||
259 | Double_t x0 = Vcenter.X(); | |
260 | ||
261 | // fX = xref X-coordinate of this track (reference plane) | |
262 | // fAlpha = Alpha Rotation angle the local (TPC sector) | |
263 | // fP0 = YL Y-coordinate of a track | |
264 | // fP1 = ZG Z-coordinate of a track | |
265 | // fP2 = C*x0 x0 is center x in rotated frame | |
266 | // fP3 = Tgl tangent of the track momentum dip angle | |
267 | // fP4 = C track curvature | |
268 | xx[0] = y; | |
269 | xx[1] = z; | |
270 | xx[3] = Pz/Pt; | |
271 | xx[4] = -ch/rho; | |
272 | xx[2] = xx[4]*x0; | |
273 | ||
274 | // create the object AliTPCtrack | |
275 | AliTPCtrack* track = new AliTPCtrack(0,xx,cc,xref,alpha); | |
276 | // set the label | |
277 | track->SetLabel(lab); | |
278 | ||
279 | return track; | |
280 | } | |
281 | ||
282 | //----------------------------------------------------------------- | |
283 | Bool_t AliTPCtrackerParam::SelectedTrack(Int_t pdg,Double_t pt,Double_t eta) | |
284 | { | |
285 | Double_t eff=0.; | |
286 | ||
287 | //eff computed with | zl+(244-xl)*pz/pt | < 252 | |
288 | Double_t effPi[27] = {0.724587,0.743389,0.619273,0.798477,0.812036,0.823195,0.771437,0.775826,0.784136,0.809071,0.762001,0.774576,0.848834,0.787201,0.792548,0.942089,0.951631,0.951085,0.960885,0.971451,0.969103,0.983245,0.978939,0.988706,0.990852,0.985679,0.993606}; | |
289 | Double_t effK[18] = {0.377934,0.363962,0.321721,0.518784,0.547459,0.517878,0.612704,0.619101,0.620894,0.733411,0.732128,0.750373,0.790630,0.806565,0.791353,0.967486,0.970483,0.974527}; | |
290 | Double_t effP[15] = {0.131173,0.165114,0.229658,0.365357,0.412989,0.483297,0.454614,0.505173,0.658615,0.694753,0.730661,0.815680,0.873461,0.887227,0.899324}; | |
291 | Double_t effEl[15] = {0.835549,0.853746,0.718207,0.835230,0.831489,0.862222,0.757783,0.747301,0.824096,0.867949,0.871891,0.808480,0.890625,0.911765,0.973684}; | |
292 | Double_t effMu[15] = {0.553486,0.641392,0.609932,0.591126,0.706729,0.750755,0.747952,0.729051,0.760849,0.898810,0.737500,0.830357,0.735294,0.800000,0.882353}; | |
293 | ||
294 | ||
295 | if(TMath::Abs(pdg)==211) eff = LinearInterpolation(9,effPi,pt,eta); | |
296 | if(TMath::Abs(pdg)==321) eff = LinearInterpolation(6,effK,pt,eta); | |
297 | if(TMath::Abs(pdg)==2212) eff = LinearInterpolation(5,effP,pt,eta); | |
298 | if(TMath::Abs(pdg)==11) eff = LinearInterpolation(5,effEl,pt,eta); | |
299 | if(TMath::Abs(pdg)==13) eff = LinearInterpolation(5,effMu,pt,eta); | |
300 | ||
301 | if(gRandom->Rndm() < eff) return kTRUE; | |
302 | ||
303 | return kFALSE; | |
304 | } | |
305 | ||
306 | //----------------------------------------------------------------- | |
307 | Double_t AliTPCtrackerParam::LinearInterpolation(Int_t ptBins,Double_t *value, | |
308 | Double_t pt,Double_t eta) | |
309 | { | |
310 | Double_t kValue=0,kValue1=0,kValue2=0; | |
311 | Int_t kEtaSide = (TMath::Abs(eta)<.45 ? 0 : 1); | |
312 | Double_t kEta[3]={0.15,0.45,0.75}; | |
313 | Double_t kPt[9]={0.244,0.390,0.676,1.190,2.36,4.,6.,10.,20.}; | |
314 | if(ptBins==6) kPt[5]=10.; | |
315 | ||
316 | for(Int_t i=0; i<ptBins; i++) { | |
317 | if(pt<kPt[i]) { | |
318 | if(i==0) i=1; | |
319 | kValue1 = value[3*i-3+kEtaSide]+(value[3*i-3+kEtaSide]-value[3*i+kEtaSide])/(kPt[i-1]-kPt[i])*(pt-kPt[i-1]); | |
320 | kValue2 = value[3*i-3+kEtaSide+1]+(value[3*i-3+kEtaSide+1]-value[3*i+kEtaSide+1])/(kPt[i-1]-kPt[i])*(pt-kPt[i-1]); | |
321 | ||
322 | kValue = kValue1+(kValue1-kValue2)/(kEta[kEtaSide]-kEta[kEtaSide+1])*(eta-kEta[kEtaSide]); | |
323 | break; | |
324 | } | |
325 | if(i==ptBins-1) { | |
326 | i=ptBins-2; | |
327 | kValue1 = value[3*i-3+kEtaSide]+(value[3*i-3+kEtaSide]-value[3*i+kEtaSide])/(kPt[i-1]-kPt[i])*(pt-kPt[i-1]); | |
328 | kValue2 = value[3*i-3+kEtaSide+1]+(value[3*i-3+kEtaSide+1]-value[3*i+kEtaSide+1])/(kPt[i-1]-kPt[i])*(pt-kPt[i-1]); | |
329 | ||
330 | kValue = kValue1+(kValue1-kValue2)/(kEta[kEtaSide]-kEta[kEtaSide+1])*(eta-kEta[kEtaSide]); | |
331 | break; | |
332 | } | |
333 | } | |
334 | ||
335 | return kValue; | |
336 | } | |
337 | ||
338 | //----------------------------------------------------------------- | |
339 | Int_t AliTPCtrackerParam::GetBin(Double_t pt,Double_t eta) { | |
340 | ||
341 | if(TMath::Abs(eta)<0.3) { | |
342 | if(pt<0.3) return 0; | |
343 | if(pt>=0.3 && pt<0.5) return 3; | |
344 | if(pt>=0.5 && pt<1.) return 6; | |
345 | if(pt>=1. && pt<1.5) return 9; | |
346 | if(pt>=1.5 && pt<3.) return 12; | |
347 | if(pt>=3. && pt<5.) return 15; | |
348 | if(pt>=5. && pt<7.) return 18; | |
349 | if(pt>=7. && pt<15.) return 21; | |
350 | if(pt>=15.) return 24; | |
351 | } | |
352 | if(TMath::Abs(eta)>=0.3 && TMath::Abs(eta)<0.6) { | |
353 | if(pt<0.3) return 1; | |
354 | if(pt>=0.3 && pt<0.5) return 4; | |
355 | if(pt>=0.5 && pt<1.) return 7; | |
356 | if(pt>=1. && pt<1.5) return 10; | |
357 | if(pt>=1.5 && pt<3.) return 13; | |
358 | if(pt>=3. && pt<5.) return 16; | |
359 | if(pt>=5. && pt<7.) return 19; | |
360 | if(pt>=7. && pt<15.) return 22; | |
361 | if(pt>=15.) return 25; | |
362 | } | |
363 | if(TMath::Abs(eta)>=0.6) { | |
364 | if(pt<0.3) return 2; | |
365 | if(pt>=0.3 && pt<0.5) return 5; | |
366 | if(pt>=0.5 && pt<1.) return 8; | |
367 | if(pt>=1. && pt<1.5) return 11; | |
368 | if(pt>=1.5 && pt<3.) return 14; | |
369 | if(pt>=3. && pt<5.) return 17; | |
370 | if(pt>=5. && pt<7.) return 20; | |
371 | if(pt>=7. && pt<15.) return 23; | |
372 | if(pt>=15.) return 26; | |
373 | } | |
374 | ||
375 | return -1; | |
376 | ||
377 | } | |
378 | ||
379 | //----------------------------------------------------------------- | |
380 | ||
381 | TMatrixD AliTPCtrackerParam::GetSmearingMatrix(Double_t* cc,Double_t pt,Double_t eta) | |
382 | { | |
383 | TMatrixD CovMat(5,5); | |
384 | ||
385 | CovMat(0,0)=cc[0]; | |
386 | CovMat(1,0)=cc[1]; CovMat(0,1)=CovMat(1,0); | |
387 | CovMat(1,1)=cc[2]; | |
388 | CovMat(2,0)=cc[3]; CovMat(0,2)=CovMat(2,0); | |
389 | CovMat(2,1)=cc[4]; CovMat(1,2)=CovMat(2,1); | |
390 | CovMat(2,2)=cc[5]; | |
391 | CovMat(3,0)=cc[6]; CovMat(0,3)=CovMat(3,0); | |
392 | CovMat(3,1)=cc[7]; CovMat(1,3)=CovMat(3,1); | |
393 | CovMat(3,2)=cc[8]; CovMat(2,3)=CovMat(3,2); | |
394 | CovMat(3,3)=cc[9]; | |
395 | CovMat(4,0)=cc[10]; CovMat(0,4)=CovMat(4,0); | |
396 | CovMat(4,1)=cc[11]; CovMat(1,4)=CovMat(4,1); | |
397 | CovMat(4,2)=cc[12]; CovMat(2,4)=CovMat(4,2); | |
398 | CovMat(4,3)=cc[13]; CovMat(3,4)=CovMat(4,3); | |
399 | CovMat(4,4)=cc[14]; | |
400 | ||
401 | TMatrixD ScaleMat(5,5); | |
402 | for(Int_t k=0;k<5;k++) { | |
403 | for(Int_t l=0;l<5;l++) { | |
404 | ScaleMat(k,l)=0.; | |
405 | } | |
406 | } | |
407 | ||
408 | ||
409 | Double_t s00[27]={1.3553,1.2973,1.2446,1.3428,1.3031,1.2345,1.3244,1.2681,1.2046,1.3046,1.2430,1.1606,1.2462,1.2104,1.1082,1.2207,1.1189,1.0789,1.2079,1.1300,1.0426,1.1502,1.1122,1.0559,1.1419,1.1072,1.0208}; | |
410 | Double_t s11[27]={1.0890,1.1463,1.2313,1.1149,1.1597,1.2280,1.1225,1.1584,1.2329,1.1149,1.1550,1.2369,1.1095,1.1353,1.2050,1.0649,1.0858,1.1546,1.0663,1.0672,1.1340,1.0416,1.0738,1.0945,1.0663,1.0654,1.0909}; | |
411 | Double_t s22[27]={1.1709,1.1367,1.0932,1.2247,1.1832,1.1247,1.2470,1.2017,1.1441,1.2202,1.1653,1.1050,1.1819,1.1269,1.0583,1.1546,1.0621,0.9928,1.1305,1.0512,0.9576,1.0995,1.0445,0.9884,1.0968,1.0368,0.9574}; | |
412 | Double_t s33[27]={0.9720,0.9869,1.0271,1.0030,1.0223,1.0479,1.0164,1.0305,1.0559,1.0339,1.0450,1.0686,1.0450,1.0507,1.0784,1.0334,1.0208,1.0863,1.0252,1.0114,1.0835,0.9854,1.0144,1.0507,1.0124,1.0159,1.0464}; | |
413 | Double_t s44[27]={1.1104,1.0789,1.0479,1.1597,1.1234,1.0728,1.2087,1.1602,1.1041,1.1942,1.1428,1.0831,1.1572,1.1036,1.0431,1.1296,1.0498,0.9844,1.1145,1.0266,0.9489,1.0959,1.0450,0.9875,1.0775,1.0266,0.9406}; | |
414 | ||
415 | ||
416 | ScaleMat(0,0) = LinearInterpolation(9,s00,pt,eta); | |
417 | ScaleMat(1,1) = LinearInterpolation(9,s11,pt,eta); | |
418 | ScaleMat(2,2) = LinearInterpolation(9,s22,pt,eta); | |
419 | ScaleMat(3,3) = LinearInterpolation(9,s33,pt,eta); | |
420 | ScaleMat(4,4) = LinearInterpolation(9,s44,pt,eta); | |
421 | ||
422 | TMatrixD Mat(ScaleMat,TMatrixD::kMult,CovMat); | |
423 | TMatrixD CovMatSmear(Mat,TMatrixD::kMult,ScaleMat); | |
424 | ||
425 | ||
426 | return CovMatSmear; | |
427 | } | |
428 | ||
429 | //----------------------------------------------------------------- | |
430 | void AliTPCtrackerParam::SmearTrack(Double_t* xx,Double_t* xxsm,TMatrixD cov) { | |
431 | ||
70521312 | 432 | AliGausCorr *corgen = new AliGausCorr(cov,5); |
433 | TArrayD corr(5); | |
434 | corgen->GetGaussN(corr); | |
435 | delete corgen; | |
436 | corgen = 0; | |
6eb67451 | 437 | |
438 | for(Int_t l=0;l<5;l++) { | |
439 | xxsm[l] = xx[l]+corr[l]; | |
440 | } | |
441 | ||
442 | return; | |
443 | } | |
444 | ||
445 | //----------------------------------------------------------------- | |
446 | void AliTPCtrackerParam::CookTracks(TObjArray& tarray,TObjArray& newtarray) | |
447 | { | |
448 | TString* s = new TString("$ALICE_ROOT/TPC/CovMatrixDB_"); | |
449 | if(fColl==0) s->Append("PbPb6000"); | |
450 | if(fBz==0.4) s->Append("_B0.4T.root"); | |
451 | // open file with matrixes DB | |
452 | TFile* DBfile = new TFile(s->Data()); | |
453 | ||
70521312 | 454 | AliTPCtrack* track = 0; |
6eb67451 | 455 | Int_t entr = (Int_t)tarray.GetEntriesFast(); |
456 | ||
457 | for(Int_t k=0; k<entr; k++) { | |
458 | track=(AliTPCtrack*)tarray.UncheckedAt(k); | |
459 | ||
460 | Double_t pt = 1/TMath::Abs(track->Get1Pt()); | |
461 | Double_t eta = -TMath::Log(TMath::Tan(0.25*TMath::Pi()-0.5*TMath::ATan(track->GetTgl()))); | |
462 | ||
463 | Int_t bin = GetBin(pt,eta); | |
464 | ||
465 | // string with tree name | |
466 | TString str("Tree_bin"); | |
467 | str+=bin; | |
468 | ||
469 | // get the right tree | |
470 | TTree* CovTree = (TTree*)DBfile->Get(str.Data()); | |
471 | TArrayF* matrix = new TArrayF; | |
472 | CovTree->SetBranchAddress("matrixes",&matrix); | |
473 | ||
474 | // get random entry from the tree | |
475 | Int_t Entries = (Int_t)CovTree->GetEntries(); | |
476 | CovTree->GetEvent(gRandom->Integer(Entries)); | |
477 | ||
478 | Double_t xref,alpha,xx[5],xxsm[5],cc[15]; | |
479 | Int_t lab; | |
480 | // get P and Cosl from track | |
481 | Double_t Cosl = TMath::Cos(TMath::ATan(track->GetTgl())); | |
482 | Double_t P = 1/TMath::Abs(track->Get1Pt())/Cosl; | |
483 | ||
484 | // get covariance matrix from regularized matrix | |
485 | cc[0]=matrix->At(0)*(1.588e-3+1.911e-4/TMath::Power(P,1.5)); | |
486 | cc[1]=matrix->At(1); | |
487 | cc[2]=matrix->At(2)*(1.195e-3+0.8102e-3/P); | |
488 | cc[3]=matrix->At(3)*(7.275e-5+1.181e-5/TMath::Power(P,1.5)); | |
489 | cc[4]=matrix->At(4); | |
490 | cc[5]=matrix->At(5)*(5.163e-6+1.138e-6/TMath::Power(P,2.)/Cosl); | |
491 | cc[6]=matrix->At(6); | |
492 | cc[7]=matrix->At(7)*(1.176e-5+1.175e-5/TMath::Power(P,1.5)/Cosl/Cosl/Cosl); | |
493 | cc[8]=matrix->At(8); | |
494 | cc[9]=matrix->At(9)*(1.289e-7+4.636e-7/TMath::Power(P,1.7)/Cosl/Cosl/Cosl/Cosl); | |
495 | cc[10]=matrix->At(10)*(4.056e-7+4.379e-8/TMath::Power(P,1.5)); | |
496 | cc[11]=matrix->At(11); | |
497 | cc[12]=matrix->At(12)*(3.049e-8+8.244e-9/TMath::Power(P,2.)/Cosl); | |
498 | cc[13]=matrix->At(13); | |
499 | cc[14]=matrix->At(14)*(1.847e-10+5.822e-11/TMath::Power(P,2.)/Cosl/Cosl/Cosl); | |
500 | ||
501 | TMatrixD CovMatSmear(5,5); | |
502 | ||
503 | CovMatSmear = GetSmearingMatrix(cc,pt,eta); | |
504 | ||
505 | // get track original parameters | |
506 | xref=track->GetX(); | |
507 | alpha=track->GetAlpha(); | |
508 | xx[0]=track->GetY(); | |
509 | xx[1]=track->GetZ(); | |
510 | xx[2]=track->GetX()*track->GetC()-track->GetSnp(); | |
511 | xx[3]=track->GetTgl(); | |
512 | xx[4]=track->GetC(); | |
513 | lab=track->GetLabel(); | |
514 | ||
515 | // use smearing matrix to smear the original parameters | |
516 | SmearTrack(xx,xxsm,CovMatSmear); | |
517 | ||
518 | AliTPCtrack* tpctrack = new AliTPCtrack(0,xxsm,cc,xref,alpha); | |
519 | tpctrack->SetLabel(lab); | |
520 | ||
521 | // fill the array | |
522 | newtarray.AddLast(tpctrack); | |
70521312 | 523 | |
524 | delete matrix; | |
525 | ||
6eb67451 | 526 | } |
527 | ||
528 | DBfile->Close(); | |
529 | ||
70521312 | 530 | delete s; |
531 | delete DBfile; | |
532 | ||
533 | ||
6eb67451 | 534 | return; |
535 | } | |
536 | //----------------------------------------------------------------- | |
537 | ||
538 | ||
539 | ||
540 | ||
541 | ||
542 | ||
543 | ||
544 | ||
545 | ||
546 | ||
547 | ||
548 |