Fixed warnings
[u/mrichter/AliRoot.git] / STEER / AliAODVertex.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2007, 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 /* $Id$ */
17
18 //-------------------------------------------------------------------------
19 //     AOD track base class
20 //     Base class for Analysis Object Data
21 //     Generic version
22 //     Author: Markus Oldenburg, CERN
23 //     Inheritance from AliVVertex: A. Dainese
24 //-------------------------------------------------------------------------
25
26 #include "AliAODVertex.h"
27 #include "AliAODTrack.h"
28
29 ClassImp(AliAODVertex)
30
31 //______________________________________________________________________________
32 AliAODVertex::AliAODVertex() : 
33   AliVVertex(),
34   fChi2perNDF(-999.),
35   fID(-1),
36   fType(kUndef),
37   fNprong(0),
38   fIprong(0),
39   fNContributors(0),
40   fCovMatrix(NULL),
41   fParent(),
42   fDaughters(),
43   fProngs(NULL)
44   {
45   // default constructor
46
47   fPosition[0] = fPosition[1] = fPosition[2] = -999.;
48 }
49
50 //______________________________________________________________________________
51 AliAODVertex::AliAODVertex(const Double_t position[3], 
52                            const Double_t covMatrix[6],
53                            Double_t  chi2perNDF,
54                            TObject  *parent,
55                            Short_t id,
56                            Char_t vtype, 
57                            Int_t  nprong) :
58   AliVVertex(),
59   fChi2perNDF(chi2perNDF),
60   fID(id),
61   fType(vtype),
62   fNprong(nprong),
63   fIprong(0),
64   fNContributors(0),
65   fCovMatrix(NULL),
66   fParent(parent),
67   fDaughters(),
68   fProngs(0)
69 {
70   // constructor
71
72   SetPosition(position);
73   if (covMatrix) SetCovMatrix(covMatrix);
74   MakeProngs();
75 }
76
77 //______________________________________________________________________________
78 AliAODVertex::AliAODVertex(const Float_t position[3], 
79                            const Float_t  covMatrix[6],
80                            Double_t  chi2perNDF,
81                            TObject  *parent,
82                            Short_t id,
83                            Char_t vtype,
84                            Int_t nprong) :
85
86   AliVVertex(),
87   fChi2perNDF(chi2perNDF),
88   fID(id),
89   fType(vtype),
90   fNprong(nprong),
91   fIprong(0),
92   fNContributors(0),
93   fCovMatrix(NULL),
94   fParent(parent),
95   fDaughters(),
96   fProngs(0)
97 {
98   // constructor
99
100   SetPosition(position);
101   if (covMatrix) SetCovMatrix(covMatrix);
102   MakeProngs();
103 }
104
105 //______________________________________________________________________________
106 AliAODVertex::AliAODVertex(const Double_t position[3], 
107                            Double_t  chi2perNDF,
108                            Char_t vtype, 
109                            Int_t nprong) :
110   AliVVertex(),
111   fChi2perNDF(chi2perNDF),
112   fID(-1),
113   fType(vtype),
114   fNprong(nprong),
115   fIprong(0),
116   fNContributors(0),
117   fCovMatrix(NULL),
118   fParent(),
119   fDaughters(),
120   fProngs(0)
121 {
122   // constructor without covariance matrix
123
124   SetPosition(position);  
125   MakeProngs();
126 }
127
128 //______________________________________________________________________________
129 AliAODVertex::AliAODVertex(const Float_t position[3], 
130                            Double_t  chi2perNDF,
131                            Char_t vtype, Int_t nprong) :
132   AliVVertex(),
133   fChi2perNDF(chi2perNDF),
134   fID(-1),
135   fType(vtype),
136   fNprong(nprong),
137   fIprong(0),
138   fNContributors(0),
139   fCovMatrix(NULL),
140   fParent(),
141   fDaughters(),
142   fProngs(0)
143 {
144   // constructor without covariance matrix
145
146   SetPosition(position);  
147   MakeProngs();
148 }
149
150 //______________________________________________________________________________
151 AliAODVertex::~AliAODVertex() 
152 {
153   // Destructor
154
155   delete fCovMatrix;
156   if (fNprong > 0) delete[] fProngs;
157 }
158
159 //______________________________________________________________________________
160 AliAODVertex::AliAODVertex(const AliAODVertex& vtx) :
161   AliVVertex(vtx),
162   fChi2perNDF(vtx.fChi2perNDF),
163   fID(vtx.fID),
164   fType(vtx.fType),
165   fNprong(vtx.fNprong),
166   fIprong(vtx.fIprong),
167   fNContributors(vtx.fNContributors),
168   fCovMatrix(NULL),
169   fParent(vtx.fParent),
170   fDaughters(vtx.fDaughters),
171   fProngs(0)
172 {
173   // Copy constructor.
174   
175   for (int i = 0; i < 3; i++) 
176     fPosition[i] = vtx.fPosition[i];
177
178   if (vtx.fCovMatrix) fCovMatrix=new AliAODRedCov<3>(*vtx.fCovMatrix);
179   MakeProngs();
180   for (int i = 0; i < fNprong; i++) {
181       fProngs[i] = vtx.fProngs[i];
182   }
183 }
184
185 //______________________________________________________________________________
186 AliAODVertex& AliAODVertex::operator=(const AliAODVertex& vtx) 
187 {
188   // Assignment operator
189   if (this != &vtx) {
190
191     // name and type
192     AliVVertex::operator=(vtx);
193
194     //momentum
195     for (int i = 0; i < 3; i++) 
196       fPosition[i] = vtx.fPosition[i];
197     
198     fChi2perNDF = vtx.fChi2perNDF;
199     fID = vtx.fID;
200     fType = vtx.fType;
201
202     //covariance matrix
203     delete fCovMatrix;
204     fCovMatrix = NULL;   
205     if (vtx.fCovMatrix) fCovMatrix=new AliAODRedCov<3>(*vtx.fCovMatrix);
206     
207     //other stuff
208     fParent = vtx.fParent;
209     fDaughters = vtx.fDaughters;
210     fNprong    = vtx.fNprong;
211     fIprong    = vtx.fIprong;  
212
213     MakeProngs();
214     for (int i = 0; i < fNprong; i++) {
215         fProngs[i] = vtx.fProngs[i];
216     }
217   }
218   
219   return *this;
220 }
221
222 //______________________________________________________________________________
223 void AliAODVertex::AddDaughter(TObject *daughter)
224 {
225   // Add reference to daughter track
226     if (!fProngs) {
227         if (fDaughters.GetEntries()==0) {
228             TRefArray* arr = &fDaughters;
229             new(arr)TRefArray(TProcessID::GetProcessWithUID(daughter));         
230         }
231         fDaughters.Add(daughter);       
232     } else {
233         if (fIprong < fNprong) {
234             fProngs[fIprong++] = daughter;
235         } else {
236             AliWarning("Number of daughters out of range !\n");
237         }
238     }
239   return;
240 }
241
242
243 //______________________________________________________________________________
244 template <class T> void AliAODVertex::GetSigmaXYZ(T sigma[3]) const
245 {
246   // Return errors on vertex position in thrust frame
247   
248   if(fCovMatrix) {
249     sigma[0]=fCovMatrix[3]; //GetCovXZ
250     sigma[1]=fCovMatrix[4]; //GetCovYZ
251     sigma[2]=fCovMatrix[5]; //GetCovZZ
252   } else 
253     sigma[0]=sigma[1]=sigma[2]=-999.;
254
255   /*
256   for (int i = 0, j = 6; i < 3; i++) {
257     j -= i+1;
258     sigma[2-i] = fCovMatrix ? TMath::Sqrt(fCovMatrix[j]) : -999.;
259   }
260   */
261 }
262
263 //______________________________________________________________________________
264 Int_t AliAODVertex::GetNContributors() const 
265 {
266   // Returns the number of tracks used to fit this vertex.
267   Int_t cont = 0;
268
269   if (!strcmp(GetTitle(), "vertexer: 3D")) {
270     cont = fNContributors;
271   } else {
272     for (Int_t iDaug = 0; iDaug < GetNDaughters(); iDaug++) {
273       if (((AliAODTrack*)fDaughters.At(iDaug))->GetUsedForVtxFit()) cont++;
274     }
275   }
276   return cont;
277 }
278
279 //______________________________________________________________________________
280 Bool_t AliAODVertex::HasDaughter(TObject *daughter) const 
281 {
282   // Checks if the given daughter (particle) is part of this vertex.
283     if (!fProngs) {
284         TRefArrayIter iter(&fDaughters);
285         while (TObject *daugh = iter.Next()) {
286             if (daugh == daughter) return kTRUE;
287         }
288         return kFALSE;
289     } else {
290         Bool_t has = kFALSE;
291         for (int i = 0; i < fNprong; i++) {
292             if (fProngs[i].GetObject() == daughter) has = kTRUE;
293         }
294         return has;
295     }
296 }
297
298 //______________________________________________________________________________
299 Double_t AliAODVertex::RotatedCovMatrixXX(Double_t phi, Double_t theta) const
300 {
301   // XX term of covariance matrix after rotation by phi around z-axis
302   // and, then, by theta around new y-axis
303
304   if (!fCovMatrix) {
305     //AliFatal("Covariance matrix not set");
306     return -999.;
307   }
308
309   Double_t covMatrix[6];
310
311   GetCovMatrix(covMatrix);
312
313   Double_t cp = TMath::Cos(phi);
314   Double_t sp = TMath::Sin(phi);
315   Double_t ct = TMath::Cos(theta);
316   Double_t st = TMath::Sin(theta);
317   return
318      covMatrix[0]*cp*cp*ct*ct  // GetCovXX
319     +covMatrix[1]*2.*cp*sp*ct*ct  // GetCovXY
320     +covMatrix[3]*2.*cp*ct*st  // GetCovXZ
321     +covMatrix[2]*sp*sp*ct*ct  // GetCovYY
322     +covMatrix[4]*2.*sp*ct*st  // GetCovYZ
323     +covMatrix[5]*st*st;  // GetCovZZ
324 }
325
326 //______________________________________________________________________________
327 Double_t AliAODVertex::RotatedCovMatrixXY(Double_t phi, Double_t theta) const
328 {
329   // XY term of covariance matrix after rotation by phi around z-axis
330   // and, then, by theta around new y-axis
331
332   if (!fCovMatrix) {
333     //AliFatal("Covariance matrix not set");
334     return -999.;
335   }
336
337   Double_t covMatrix[6];
338
339   GetCovMatrix(covMatrix);
340
341   Double_t cp = TMath::Cos(phi);
342   Double_t sp = TMath::Sin(phi);
343   Double_t ct = TMath::Cos(theta);
344   Double_t st = TMath::Sin(theta);
345   return 
346     -covMatrix[0]*cp*sp*ct  // GetCovXX
347     +covMatrix[1]*ct*(cp*cp-sp*sp)  // GetCovXY
348     -covMatrix[3]*sp*st  // GetCovXZ
349     +covMatrix[2]*cp*sp*ct  // GetCovYY
350     +covMatrix[4]*cp*st;  // GetCovYZ
351 }
352
353 //______________________________________________________________________________
354 Double_t AliAODVertex::RotatedCovMatrixXZ(Double_t phi, Double_t theta) const
355 {
356   // XZ term of covariance matrix after rotation by phi around z-axis
357   // and, then, by theta around new y-axis
358
359   if (!fCovMatrix) {
360     //AliFatal("Covariance matrix not set");
361     return -999.;
362   }
363
364   Double_t covMatrix[6];
365
366   GetCovMatrix(covMatrix);
367
368   Double_t cp = TMath::Cos(phi);
369   Double_t sp = TMath::Sin(phi);
370   Double_t ct = TMath::Cos(theta);
371   Double_t st = TMath::Sin(theta);
372   return 
373     -covMatrix[0]*cp*cp*ct*st  // GetCovXX
374     -covMatrix[1]*2.*cp*sp*ct*st  // GetCovXY
375     +covMatrix[3]*cp*(ct*ct-st*st)  // GetCovXZ
376     -covMatrix[2]*sp*sp*ct*st  // GetCovYY
377     +covMatrix[4]*sp*(ct*ct-st*st)  // GetCovYZ
378     +covMatrix[5]*ct*st;  // GetCovZZ
379 }
380
381 //______________________________________________________________________________
382 Double_t AliAODVertex::RotatedCovMatrixYY(Double_t phi) const
383 {
384   // YY term of covariance matrix after rotation by phi around z-axis
385   // and, then, by theta around new y-axis
386
387   if (!fCovMatrix) {
388     //AliFatal("Covariance matrix not set");
389     return -999.;
390   }
391
392   Double_t covMatrix[6];
393
394   GetCovMatrix(covMatrix);
395
396   Double_t cp = TMath::Cos(phi);
397   Double_t sp = TMath::Sin(phi);
398   return
399      covMatrix[0]*sp*sp  // GetCovXX
400     -covMatrix[1]*2.*cp*sp  // GetCovXY
401     +covMatrix[2]*cp*cp;  // GetCovYY
402 }
403
404 //______________________________________________________________________________
405 Double_t AliAODVertex::RotatedCovMatrixYZ(Double_t phi, Double_t theta) const
406 {
407   // YZ term of covariance matrix after rotation by phi around z-axis
408   // and, then, by theta around new y-axis
409
410   if (!fCovMatrix) {
411     //AliFatal("Covariance matrix not set");
412     return -999.;
413   }
414
415   Double_t covMatrix[6];
416
417   GetCovMatrix(covMatrix);
418
419   Double_t cp = TMath::Cos(phi);
420   Double_t sp = TMath::Sin(phi);
421   Double_t ct = TMath::Cos(theta);
422   Double_t st = TMath::Sin(theta);
423   return 
424      covMatrix[0]*cp*sp*st  // GetCovXX
425     +covMatrix[1]*st*(sp*sp-cp*cp)  // GetCovXY
426     -covMatrix[3]*sp*ct  // GetCovXZ
427     -covMatrix[2]*cp*sp*st  // GetCovYY
428     +covMatrix[4]*cp*ct;  // GetCovYZ
429 }
430
431 //______________________________________________________________________________
432 Double_t AliAODVertex::RotatedCovMatrixZZ(Double_t phi, Double_t theta) const
433 {
434   // ZZ term of covariance matrix after rotation by phi around z-axis
435   // and, then, by theta around new y-axis
436
437   if (!fCovMatrix) {
438     //AliFatal("Covariance matrix not set");
439     return -999.;
440   }
441
442   Double_t covMatrix[6];
443
444   GetCovMatrix(covMatrix);
445
446   Double_t cp = TMath::Cos(phi);
447   Double_t sp = TMath::Sin(phi);
448   Double_t ct = TMath::Cos(theta);
449   Double_t st = TMath::Sin(theta);
450   return
451      covMatrix[0]*cp*cp*st*st  // GetCovXX
452     +covMatrix[1]*2.*cp*sp*st*st  // GetCovXY
453     -covMatrix[3]*2.*cp*ct*st  // GetCovXZ
454     +covMatrix[2]*sp*sp*st*st  // GetCovYY
455     -covMatrix[4]*2.*sp*sp*ct*st  // GetCovYZ
456     +covMatrix[5]*ct*ct;  // GetCovZZ
457 }
458
459 //______________________________________________________________________________
460 Double_t AliAODVertex::DistanceToVertex(AliAODVertex *vtx) const
461 {
462   // distance in 3D to another AliAODVertex
463
464   Double_t dx = GetX()-vtx->GetX();
465   Double_t dy = GetY()-vtx->GetY();
466   Double_t dz = GetZ()-vtx->GetZ();
467
468   return TMath::Sqrt(dx*dx+dy*dy+dz*dz);
469 }
470
471 //______________________________________________________________________________
472 Double_t AliAODVertex::DistanceXYToVertex(AliAODVertex *vtx) const
473 {
474   // distance in XY to another AliAODVertex
475
476   Double_t dx = GetX()-vtx->GetX();
477   Double_t dy = GetY()-vtx->GetY();
478
479   return TMath::Sqrt(dx*dx+dy*dy);
480 }
481
482 //______________________________________________________________________________
483 Double_t AliAODVertex::ErrorDistanceToVertex(AliAODVertex *vtx) const
484 {
485   // error on the distance in 3D to another AliAODVertex
486
487   Double_t phi,theta;
488   PhiAndThetaToVertex(vtx,phi,theta);
489   // error2 due to this vertex
490   Double_t error2 = RotatedCovMatrixXX(phi,theta);
491   // error2 due to vtx vertex
492   Double_t error2vtx = vtx->RotatedCovMatrixXX(phi,theta);
493
494   return TMath::Sqrt(error2+error2vtx);
495 }
496
497 //______________________________________________________________________________
498 Double_t AliAODVertex::ErrorDistanceXYToVertex(AliAODVertex *vtx) const
499 {
500   // error on the distance in XY to another AliAODVertex
501
502   Double_t phi,theta;
503   PhiAndThetaToVertex(vtx,phi,theta);
504   // error2 due to this vertex
505   Double_t error2 = RotatedCovMatrixXX(phi);
506   // error2 due to vtx vertex
507   Double_t error2vtx = vtx->RotatedCovMatrixXX(phi);
508
509   return TMath::Sqrt(error2+error2vtx);
510 }
511
512 //______________________________________________________________________________
513 template <class T, class P> 
514 void AliAODVertex::PhiAndThetaToVertex(AliAODVertex *vtx, P &phi, T &theta) const
515 {
516   // rotation angles around z-axis (phi) and around new y-axis (theta)
517   // with which vtx is seen (used by RotatedCovMatrix... methods)
518
519   phi = TMath::Pi()+TMath::ATan2(-vtx->GetY()+GetY(),-vtx->GetX()+GetX());
520   Double_t vtxxphi = vtx->GetX()*TMath::Cos(phi)+vtx->GetY()*TMath::Sin(phi);
521   Double_t xphi = GetX()*TMath::Cos(phi)+GetY()*TMath::Sin(phi);
522   theta = TMath::ATan2(vtx->GetZ()-GetZ(),vtxxphi-xphi);
523 }
524
525 //______________________________________________________________________________
526 void AliAODVertex::PrintIndices() const 
527 {
528   // Print indices of particles originating form this vertex
529
530   TRefArrayIter iter(&fDaughters);
531   while (TObject *daugh = iter.Next()) {
532     printf("Particle %p originates from this vertex.\n", static_cast<void*>(daugh));
533   }
534 }
535
536 //______________________________________________________________________________
537 void AliAODVertex::Print(Option_t* /*option*/) const 
538 {
539   // Print information of all data members
540
541   printf("Vertex position:\n");
542   printf("     x = %f\n", fPosition[0]);
543   printf("     y = %f\n", fPosition[1]);
544   printf("     z = %f\n", fPosition[2]);
545   printf(" parent particle: %p\n", static_cast<void*>(fParent.GetObject()));
546   printf(" origin of %d particles\n", fDaughters.GetEntriesFast());
547   printf(" vertex type %d\n", fType);
548   
549   /*
550   if (fCovMatrix) {
551     printf("Covariance matrix:\n");
552     printf(" %12.10f  %12.10f  %12.10f\n %12.10f  %12.10f  %12.10f\n %12.10f  %12.10f  %12.10f\n", 
553            fCovMatrix[0],
554            fCovMatrix[1],
555            fCovMatrix[3],
556            fCovMatrix[1],
557            fCovMatrix[2],
558            fCovMatrix[4],
559            fCovMatrix[3],
560            fCovMatrix[4],
561            fCovMatrix[5]); 
562            } */
563   printf(" Chi^2/NDF = %f\n", fChi2perNDF);
564 }
565