QA ref defaut storage setter in sim and rec
[u/mrichter/AliRoot.git] / PWG1 / AliESDRecInfo.cxx
CommitLineData
6fc428f0 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///////////////////////////////////////////////////////////////////////////////
18// //
19// Time Projection Chamber //
5b662585 20// Comparison macro for reconstructed tracks //
6fc428f0 21// responsible:
22// marian.ivanov@cern.ch //
23//
24//
25
26
27
28
29
6fc428f0 30//ROOT includes
31#include "Rtypes.h"
32//
33//ALIROOT includes
34//
35#include "AliESDtrack.h"
5b662585 36#include "AliTracker.h"
6fc428f0 37#include "AliTPCParam.h"
38#include "AliTrackReference.h"
39#include "AliTPCParamSR.h"
6fc428f0 40#include "AliESDfriend.h"
41#include "AliESDtrack.h"
42#include "AliTPCseed.h"
43#include "AliITStrackMI.h"
4a50b4f7 44#include "AliTRDtrackV1.h"
76472f75 45#include "AliMCInfo.h"
6fc428f0 46#include "AliESDRecInfo.h"
47
48
49
50ClassImp(AliESDRecInfo)
51
52
53
54
55AliTPCParam * GetTPCParam(){
56 AliTPCParamSR * par = new AliTPCParamSR;
57 par->Update();
58 return par;
59}
60
61
62
63
64AliESDRecInfo::AliESDRecInfo():
65 fITSOn(0), // ITS refitted inward
66 fTRDOn(0), // ITS refitted inward
67 fDeltaP(0), //delta of momenta
68 fSign(0), // sign
69 fReconstructed(0), //flag if track was reconstructed
70 fFake(0), // fake track
71 fMultiple(0), // number of reconstructions
72 fTPCOn(0), // TPC refitted inward
73 fBestTOFmatch(0), //best matching between times
74 fESDtrack(0), // esd track
75 fTrackF(0), // friend track
76 fTPCtrack(0), // tpc track
77 fITStrack(0), // its track
18e1354a 78 fTRDtrack(0), // trd track
79 fTracks(0) // array of tracks with the same label
6fc428f0 80{
81 //
82 // default constructor
83 //
84}
85
86
87AliESDRecInfo::AliESDRecInfo(const AliESDRecInfo& recinfo):
cd875161 88 TObject(),
89 fITSOn(0), // ITS refitted inward
90 fTRDOn(0), // ITS refitted inward
91 fDeltaP(0), //delta of momenta
92 fSign(0), // sign
93 fReconstructed(0), //flag if track was reconstructed
94 fFake(0), // fake track
95 fMultiple(0), // number of reconstructions
96 fTPCOn(0), // TPC refitted inward
97 fBestTOFmatch(0), //best matching between times
98 fESDtrack(0), // esd track
99 fTrackF(0), // friend track
100 fTPCtrack(0), // tpc track
101 fITStrack(0), // its track
18e1354a 102 fTRDtrack(0), // trd track
103 fTracks(0) // array of tracks with the same label
6fc428f0 104{
105 //
106 //
107 //
108 memcpy(this,&recinfo, sizeof(recinfo));
109 fESDtrack=0; fTrackF=0; fTPCtrack=0;fITStrack=0;fTRDtrack=0;
18e1354a 110 fTracks=0;
6fc428f0 111 SetESDtrack(recinfo.GetESDtrack());
112}
113
114
cd875161 115AliESDRecInfo& AliESDRecInfo::operator=(const AliESDRecInfo& info) {
116 //
117 // Assignment operator
118 //
a1e6aa99 119 this->~AliESDRecInfo();
cd875161 120 new (this) AliESDRecInfo(info);
121 return *this;
122}
123
124
125
6fc428f0 126AliESDRecInfo::~AliESDRecInfo()
127
128{
129 //
130 // destructor
131 //
132 if (fESDtrack) { delete fESDtrack; fESDtrack=0;}
133 if (fTrackF) { delete fTrackF; fTrackF=0;}
134 if (fTPCtrack) { delete fTPCtrack; fTPCtrack=0;}
135 if (fITStrack) { delete fITStrack; fITStrack=0;}
136 if (fTRDtrack) { delete fTRDtrack; fTRDtrack=0;}
18e1354a 137 if (fTracks) { delete fTracks; fTracks=0;}
6fc428f0 138}
139
140
141
142void AliESDRecInfo::Reset()
143{
144 //
145 // reset info
146 //
147 fMultiple =0;
148 fFake =0;
149 fReconstructed=0;
150 if (fESDtrack) { delete fESDtrack; fESDtrack=0;}
151 if (fTrackF) { delete fTrackF; fTrackF=0;}
152 if (fTPCtrack) { delete fTPCtrack; fTPCtrack=0;}
153 if (fITStrack) { delete fITStrack; fITStrack=0;}
154 if (fTRDtrack) { delete fTRDtrack; fTRDtrack=0;}
18e1354a 155 if (fTracks) { delete fTracks; fTracks=0;}
6fc428f0 156}
157
158void AliESDRecInfo::SetESDtrack(const AliESDtrack *track){
159 //
160 //
161 //
162 if (fESDtrack) delete fESDtrack;
163 fESDtrack = (AliESDtrack*)track->Clone();
9ec43bae 164 if (track->GetFriendTrack()){
6fc428f0 165 if (fTrackF) delete fTrackF;
166 fTrackF = (AliESDfriendTrack*)track->GetFriendTrack()->Clone();
9ec43bae 167 Int_t icalib=0;
168 TObject *cobject=0;
169 //
170 while (fTrackF->GetCalibObject(icalib)){
171 cobject=fTrackF->GetCalibObject(icalib);
172 if (dynamic_cast<AliTPCseed*>(cobject)){
4a50b4f7 173 if (fTPCtrack) delete fTPCtrack;
174 fTPCtrack = (AliTPCseed*)(dynamic_cast<AliTPCseed*>(cobject))->Clone();
175 } else if (dynamic_cast<AliTRDtrackV1*>(cobject)){
176 if (fTRDtrack) delete fTRDtrack;
177 fTRDtrack = (AliTRDtrackV1*)(dynamic_cast<AliTRDtrackV1*>(cobject))->Clone();
9ec43bae 178 }
179 icalib++;
6fc428f0 180 }
181 }
182
183}
184
18e1354a 185
186
187void AliESDRecInfo::AddESDtrack(const AliESDtrack *track, AliMCInfo* info){
188 //
189 // Add ESDtrack
190 //
191 AliESDtrack *nctrack = (AliESDtrack*) track;
192 fMultiple++;
193 if (!fESDtrack) {
194 SetESDtrack(track);
195 Update(info,0,kTRUE);
196 return;
197 }
198 if (!fTracks) fTracks = new TClonesArray("AliESDtrack",10);
199 Int_t ntracks = fTracks->GetEntriesFast();
200 new ((*fTracks)[ntracks]) AliESDtrack(*track);
201 if (nctrack->GetKinkIndex(0)>0) return;
202 //
203 //
204 //
205 if (!nctrack->IsOn(AliESDtrack::kTPCrefit) && fStatus[1]==2) return;
206 if (!nctrack->IsOn(AliESDtrack::kITSin) && fStatus[0]>0) return;
207 if ( nctrack->GetParameter()[4]*info->GetCharge()<0) return; //returning track
208
209 Float_t dtheta = TMath::ATan(nctrack->GetTgl())-info->GetParticle().Theta()-TMath::Pi()-2;
210 if (TMath::Abs(dtheta)>0.1) return;
211
212 SetESDtrack(track);
213 Update(info,0,kTRUE);
214}
215
216
217
6fc428f0 218void AliESDRecInfo::UpdatePoints(AliESDtrack*track)
219{
220 //
221 //
222 Int_t iclusters[200];
223 Float_t density[160];
224 for (Int_t i=0;i<160;i++) density[i]=-1.;
225 fTPCPoints[0]= 160;
226 fTPCPoints[1] = -1;
227 //
228 if (fTPCPoints[0]<fTPCPoints[1]) return;
229 // Int_t nclusters=track->GetTPCclusters(iclusters);
230
231 Int_t ngood=0;
232 Int_t undeff=0;
233 Int_t nall =0;
234 Int_t range=20;
235 for (Int_t i=0;i<160;i++){
236 Int_t last = i-range;
237 if (nall<range) nall++;
238 if (last>=0){
239 if (iclusters[last]>0&& (iclusters[last]&0x8000)==0) ngood--;
240 if (iclusters[last]==-1) undeff--;
241 }
242 if (iclusters[i]>0&& (iclusters[i]&0x8000)==0) ngood++;
243 if (iclusters[i]==-1) undeff++;
244 if (nall==range &&undeff<range/2) density[i-range/2] = Float_t(ngood)/Float_t(nall-undeff);
245 }
246 Float_t maxdens=0;
247 Int_t indexmax =0;
248 for (Int_t i=0;i<160;i++){
249 if (density[i]<0) continue;
250 if (density[i]>maxdens){
251 maxdens=density[i];
252 indexmax=i;
253 }
254 }
255 //
256 //max dens point
257 fTPCPoints[3] = maxdens;
258 fTPCPoints[1] = indexmax;
259 //
260 // last point
261 for (Int_t i=indexmax;i<160;i++){
262 if (density[i]<0) continue;
263 if (density[i]<maxdens/2.) {
264 break;
265 }
266 fTPCPoints[2]=i;
267 }
268 //
269 // first point
270 for (Int_t i=indexmax;i>0;i--){
271 if (density[i]<0) continue;
272 if (density[i]<maxdens/2.) {
273 break;
274 }
275 fTPCPoints[0]=i;
276 }
277 //
278 // Density at the last 30 padrows
279 //
280 //
281 nall = 0;
282 ngood = 0;
283 for (Int_t i=159;i>0;i--){
284 if (iclusters[i]==-1) continue; //dead zone
285 nall++;
286 if (iclusters[i]>0) ngood++;
287 if (nall>20) break;
288 }
289 fTPCPoints[4] = Float_t(ngood)/Float_t(nall);
290 //
291 if ((track->GetStatus()&AliESDtrack::kITSrefit)>0) fTPCPoints[0]=-1;
292
293
294}
295
296//
297//
298void AliESDRecInfo::Update(AliMCInfo* info,AliTPCParam * /*par*/, Bool_t reconstructed)
299{
300 //
6fc428f0 301 //calculates derived variables
302 //
6fc428f0 303 UpdatePoints(fESDtrack);
5b662585 304 UpdateStatus(info,reconstructed);
305 UpdateITS(info);
306 UpdateTPC(info);
307 UpdateTOF(info);
308}
309
310
311void AliESDRecInfo::UpdateStatus(AliMCInfo* info, Bool_t reconstructed){
6fc428f0 312 //
5b662585 313 // Interpret bit mask flags
6fc428f0 314 //
315 for (Int_t i=0;i<4;i++) fStatus[i] =0;
316 fReconstructed = kFALSE;
317 fTPCOn = kFALSE;
318 fITSOn = kFALSE;
319 fTRDOn = kFALSE;
320 if (reconstructed==kFALSE) return;
321
322 fLabels[0] = info->fLabel;
323 fLabels[1] = info->fPrimPart;
324 fReconstructed = kTRUE;
325 fTPCOn = ((fESDtrack->GetStatus()&AliESDtrack::kTPCrefit)>0) ? kTRUE : kFALSE;
326 fITSOn = ((fESDtrack->GetStatus()&AliESDtrack::kITSrefit)>0) ? kTRUE : kFALSE;
327 fTRDOn = ((fESDtrack->GetStatus()&AliESDtrack::kTRDrefit)>0) ? kTRUE : kFALSE;
328 //
329 //
330 if ((fESDtrack->GetStatus()&AliESDtrack::kTPCrefit)>0){
331 fStatus[1] =3;
332 }
333 else{
334 if ((fESDtrack->GetStatus()&AliESDtrack::kTPCout)>0){
335 fStatus[1] =2;
336 }
337 else{
338 if ((fESDtrack->GetStatus()&AliESDtrack::kTPCin)>0)
339 fStatus[1]=1;
340 }
341 }
342 //
343 if ((fESDtrack->GetStatus()&AliESDtrack::kITSout)>0){
344 fStatus[0] =2;
345 }
346 else{
347 if ((fESDtrack->GetStatus()&AliESDtrack::kITSrefit)>0){
348 fStatus[0] =1;
349 }
350 else{
351 fStatus[0]=0;
352 }
353 }
6fc428f0 354 //
355 //
356 if ((fESDtrack->GetStatus()&AliESDtrack::kTRDrefit)>0){
357 fStatus[2] =2;
358 }
359 else{
360 if ((fESDtrack->GetStatus()&AliESDtrack::kTRDout)>0){
361 fStatus[2] =1;
362 }
363 }
364 if ((fESDtrack->GetStatus()&AliESDtrack::kTRDStop)>0){
365 fStatus[2] =10;
5b662585 366 }
367}
6fc428f0 368
5b662585 369void AliESDRecInfo::UpdateTOF(AliMCInfo* info){
6fc428f0 370 //
5b662585 371 // Update TOF related comparison information
372 //
373 fBestTOFmatch=1000;
6fc428f0 374 if (((fESDtrack->GetStatus()&AliESDtrack::kTOFout)>0)){
375 //
376 // best tof match
377 Double_t times[5];
378 fESDtrack->GetIntegratedTimes(times);
379 for (Int_t i=0;i<5;i++){
380 if ( TMath::Abs(fESDtrack->GetTOFsignal()-times[i]) <TMath::Abs(fBestTOFmatch) ){
381 fBestTOFmatch = fESDtrack->GetTOFsignal()-times[i];
382 }
383 }
384 Int_t toflabel[3];
385 fESDtrack->GetTOFLabel(toflabel);
386 Bool_t toffake=kTRUE;
387 Bool_t tofdaughter=kFALSE;
388 for (Int_t i=0;i<3;i++){
389 if (toflabel[i]<0) continue;
390 if (toflabel[i]== TMath::Abs(fESDtrack->GetLabel())) toffake=kFALSE;
391 if (toflabel[i]==info->fParticle.GetDaughter(0) || (toflabel[i]==info->fParticle.GetDaughter(1))) tofdaughter=kTRUE; // decay product of original particle
392 fStatus[3]=1;
393 }
394 if (toffake) fStatus[3] =3; //total fake
395 if (tofdaughter) fStatus[3]=2; //fake because of decay
396 }else{
397 fStatus[3]=0;
398 }
5b662585 399}
6fc428f0 400
401
5b662585 402
403void AliESDRecInfo::UpdateITS(AliMCInfo* info){
404 //
405 // Update ITS related comparison information
406 //
407 fITSinP0[0]=info->fParticle.Px();
408 fITSinP0[1]=info->fParticle.Py();
409 fITSinP0[2]=info->fParticle.Pz();
410 fITSinP0[3]=info->fParticle.Pt();
411 //
412 fITSinR0[0]=info->fParticle.Vx();
413 fITSinR0[1]=info->fParticle.Vy();
414 fITSinR0[2]=info->fParticle.Vz();
415 fITSinR0[3] = TMath::Sqrt(fITSinR0[0]*fITSinR0[0]+fITSinR0[1]*fITSinR0[1]);
416 fITSinR0[4] = TMath::ATan2(fITSinR0[1],fITSinR0[0]);
417 //
418 //
419 if (fITSinP0[3]>0.0000001){
420 fITSAngle0[0] = TMath::ATan2(fITSinP0[1],fITSinP0[0]);
421 fITSAngle0[1] = TMath::ATan(fITSinP0[2]/fITSinP0[3]);
6fc428f0 422 }
423 if (fITSOn){
424 // ITS
425 Double_t param[5],x;
426 fESDtrack->GetExternalParameters(x,param);
427 // fESDtrack->GetConstrainedExternalParameters(x,param);
428 Double_t cov[15];
429 fESDtrack->GetExternalCovariance(cov);
430 //fESDtrack->GetConstrainedExternalCovariance(cov);
431 if (TMath::Abs(param[4])<0.0000000001) return;
5b662585 432
6fc428f0 433 fESDtrack->GetXYZ(fITSinR1);
434 fESDtrack->GetPxPyPz(fITSinP1);
435 fITSinP1[3] = TMath::Sqrt(fITSinP1[0]*fITSinP1[0]+fITSinP1[1]*fITSinP1[1]);
436 //
437 fITSinR1[3] = TMath::Sqrt(fITSinR1[0]*fITSinR1[0]+fITSinR1[1]*fITSinR1[1]);
438 fITSinR1[4] = TMath::ATan2(fITSinR1[1],fITSinR1[0]);
439 //
440 //
441 if (fITSinP1[3]>0.0000001){
442 fITSAngle1[0] = TMath::ATan2(fITSinP1[1],fITSinP1[0]);
443 fITSAngle1[1] = TMath::ATan(fITSinP1[2]/fITSinP1[3]);
444 }
445 //
446 //
447 fITSDelta[0] = (fITSinR0[4]-fITSinR1[4])*fITSinR1[3]; //delta rfi
448 fITSPools[0] = fITSDelta[0]/TMath::Sqrt(cov[0]);
449 fITSDelta[1] = (fITSinR0[2]-fITSinR1[2]); //delta z
450 fITSPools[1] = fITSDelta[1]/TMath::Sqrt(cov[2]);
451 fITSDelta[2] = (fITSAngle0[0]-fITSAngle1[0]);
452 fITSPools[2] = fITSDelta[2]/TMath::Sqrt(cov[5]);
453 fITSDelta[3] = (TMath::Tan(fITSAngle0[1])-TMath::Tan(fITSAngle1[1]));
454 fITSPools[3] = fITSDelta[3]/TMath::Sqrt(cov[9]);
455 fITSDelta[4] = (fITSinP0[3]-fITSinP1[3]);
456 Double_t sign = (param[4]>0) ? 1:-1;
457 fSign = sign;
5b662585 458 fITSPools[4] = sign*(1./fITSinP0[3]-1./fITSinP1[3])/TMath::Sqrt(cov[14]);
6fc428f0 459 }
6fc428f0 460}
461
5b662585 462void AliESDRecInfo::UpdateTPC(AliMCInfo* info){
463 //
464 //
465 // Update TPC related information
466 //
467 AliTrackReference * ref = &(info->fTrackRef);
468 fTPCinR0[0] = info->fTrackRef.X();
469 fTPCinR0[1] = info->fTrackRef.Y();
470 fTPCinR0[2] = info->fTrackRef.Z();
471 fTPCinR0[3] = TMath::Sqrt(fTPCinR0[0]*fTPCinR0[0]+fTPCinR0[1]*fTPCinR0[1]);
472 fTPCinR0[4] = TMath::ATan2(fTPCinR0[1],fTPCinR0[0]);
473 //
474 fTPCinP0[0] = ref->Px();
475 fTPCinP0[1] = ref->Py();
476 fTPCinP0[2] = ref->Pz();
477 fTPCinP0[3] = ref->Pt();
478 fTPCinP0[4] = ref->P();
479 fDeltaP = (ref->P()-info->fParticle.P())/info->fParticle.P();
480 //
481 //
482 if (fTPCinP0[3]>0.0000001){
483 //
484 fTPCAngle0[0] = TMath::ATan2(fTPCinP0[1],fTPCinP0[0]);
485 fTPCAngle0[1] = TMath::ATan(fTPCinP0[2]/fTPCinP0[3]);
486 }
487 //
488 if (fStatus[1]>0 &&info->fNTPCRef>0&&TMath::Abs(fTPCinP0[3])>0.0001){
489 //TPC
490 fESDtrack->GetInnerXYZ(fTPCinR1);
491 fTPCinR1[3] = TMath::Sqrt(fTPCinR1[0]*fTPCinR1[0]+fTPCinR1[1]*fTPCinR1[1]);
492 fTPCinR1[4] = TMath::ATan2(fTPCinR1[1],fTPCinR1[0]);
493 fESDtrack->GetInnerPxPyPz(fTPCinP1);
494 fTPCinP1[3] = TMath::Sqrt(fTPCinP1[0]*fTPCinP1[0]+fTPCinP1[1]*fTPCinP1[1]);
495 fTPCinP1[4] = TMath::Sqrt(fTPCinP1[3]*fTPCinP1[3]+fTPCinP1[2]*fTPCinP1[2]);
496 //
497 //
498 if (fTPCinP1[3]>0.000000000000001){
499 fTPCAngle1[0] = TMath::ATan2(fTPCinP1[1],fTPCinP1[0]);
500 fTPCAngle1[1] = TMath::ATan(fTPCinP1[2]/fTPCinP1[3]);
501 }
502 Double_t cov[15], param[5],x, alpha;
503 fESDtrack->GetInnerExternalCovariance(cov);
504 fESDtrack->GetInnerExternalParameters(alpha, x,param);
505 if (x<50) return ;
506 //
507 fTPCDelta[0] = (fTPCinR0[4]-fTPCinR1[4])*fTPCinR1[3]; //delta rfi
508 fTPCPools[0] = fTPCDelta[0]/TMath::Sqrt(cov[0]);
509 fTPCDelta[1] = (fTPCinR0[2]-fTPCinR1[2]); //delta z
510 fTPCPools[1] = fTPCDelta[1]/TMath::Sqrt(cov[2]);
511 fTPCDelta[2] = (fTPCAngle0[0]-fTPCAngle1[0]);
512 fTPCPools[2] = fTPCDelta[2]/TMath::Sqrt(cov[5]);
513 fTPCDelta[3] = (TMath::Tan(fTPCAngle0[1])-TMath::Tan(fTPCAngle1[1]));
514 fTPCPools[3] = fTPCDelta[3]/TMath::Sqrt(cov[9]);
515 fTPCDelta[4] = (fTPCinP0[3]-fTPCinP1[3]);
516 Double_t sign = (param[4]>0)? 1.:-1;
517 fSign =sign;
518 fTPCPools[4] = sign*(1./fTPCinP0[3]-1./fTPCinP1[3])/TMath::Sqrt(TMath::Abs(cov[14]));
519 }
5b662585 520}
521
522
523
524
6fc428f0 525