Package upgrade.
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnDaughter.cxx
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 // Class AliRsnDaughter
18 //
19 //
20 // Light-weight 'track' object into an internal format used
21 // for further steps of resonance analysis.
22 // Provides converters from all kinds of input track type
23 // (ESD, AOD and MC).
24 // Contains also a facility to compute invariant mass of a pair.
25 //
26 // author: A. Pulvirenti --- email: alberto.pulvirenti@ct.infn.it
27 //=========================================================================
28
29 #include <Riostream.h>
30
31 #include <TParticle.h>
32 #include <TString.h>
33
34 #include "AliLog.h"
35 #include "AliESDtrack.h"
36 #include "AliAODTrack.h"
37 #include "AliMCParticle.h"
38
39 #include "AliRsnPID.h"
40 #include "AliRsnMCInfo.h"
41 #include "AliRsnDaughter.h"
42
43 ClassImp(AliRsnDaughter)
44
45 //_____________________________________________________________________________
46 AliRsnDaughter::AliRsnDaughter() :
47   AliVParticle(),
48   fIndex(-1),
49   fLabel(-1),
50   fCharge(0),
51   fFlags(0),
52   fPIDType(AliRsnPID::kUnknown),
53   fMass(0.0),
54   fMCInfo(0x0)
55 {
56 //
57 // Default constructor.
58 // Initializes all data-members with meaningless values.
59 //
60
61     Int_t i;
62     for (i = 0; i < AliRsnPID::kSpecies; i++) {
63         if (i < 3) {
64             fP[i] = 0.0;
65             fV[i] = 0.0;
66         }
67         fPIDWeight[i] = 0.0;
68         fPIDProb[i] = 0.0;
69     }
70 }
71
72 //_____________________________________________________________________________
73 AliRsnDaughter::AliRsnDaughter(const AliRsnDaughter &copy) :
74   AliVParticle(copy),
75   fIndex(copy.fIndex),
76   fLabel(copy.fLabel),
77   fCharge(copy.fCharge),
78   fFlags(copy.fFlags),
79   fPIDType(copy.fPIDType),
80   fMass(copy.fMass),
81   fMCInfo(0x0)
82 {
83 //
84 // Copy constructor.
85 //
86
87     Int_t i;
88     for (i = 0; i < AliRsnPID::kSpecies; i++) {
89         if (i < 3) {
90             fP[i] = copy.fP[i];
91             fV[i] = copy.fV[i];
92         }
93         fPIDWeight[i] = copy.fPIDWeight[i];
94         fPIDProb[i] = copy.fPIDProb[i];
95     }
96
97     // initialize particle object
98     // only if it is present in the template object
99     if (copy.fMCInfo) fMCInfo = new AliRsnMCInfo(*(copy.fMCInfo));
100 }
101
102 //_____________________________________________________________________________
103 AliRsnDaughter::AliRsnDaughter(AliESDtrack *track, Bool_t useTPC) :
104   AliVParticle(),
105   fIndex(-1),
106   fLabel(-1),
107   fCharge(0),
108   fFlags(0),
109   fPIDType(AliRsnPID::kUnknown),
110   fMass(0.0),
111   fMCInfo(0x0)
112 {
113 //
114 // Constructor to get data from an ESD track.
115 //
116
117     Int_t i;
118     for (i = 0; i < AliRsnPID::kSpecies; i++) fPIDProb[i] = 0.0;
119     Adopt(track, useTPC);
120 }
121
122 //_____________________________________________________________________________
123 AliRsnDaughter::AliRsnDaughter(AliAODTrack *track) :
124   AliVParticle(),
125   fIndex(-1),
126   fLabel(-1),
127   fCharge(0),
128   fFlags(0),
129   fPIDType(AliRsnPID::kUnknown),
130   fMass(0.0),
131   fMCInfo(0x0)
132 {
133 //
134 // Constructor to get data from an AOD track.
135 //
136
137     Int_t i;
138     for (i = 0; i < AliRsnPID::kSpecies; i++) fPIDProb[i] = 0.0;
139     Adopt(track);
140 }
141
142 //_____________________________________________________________________________
143 AliRsnDaughter::AliRsnDaughter(AliMCParticle *track) :
144   AliVParticle(),
145   fIndex(-1),
146   fLabel(-1),
147   fCharge(0),
148   fFlags(0),
149   fPIDType(AliRsnPID::kUnknown),
150   fMass(0.0),
151   fMCInfo(0x0)
152 {
153 //
154 // Constructor to get data from an MC track.
155 //
156
157     Int_t i;
158     for (i = 0; i < AliRsnPID::kSpecies; i++) fPIDProb[i] = 0.0;
159     Adopt(track);
160 }
161
162 //_____________________________________________________________________________
163 AliRsnDaughter& AliRsnDaughter::operator=(const AliRsnDaughter &copy)
164 {
165 //
166 // Assignment operator.
167 // Works like the copy constructor and returns a reference
168 // to the initialized object for which it is called.
169 //
170
171     fIndex  = copy.fIndex;
172     fLabel  = copy.fLabel;
173     fCharge = copy.fCharge;
174     fFlags  = copy.fFlags;
175
176     Int_t i;
177     for (i = 0; i < AliRsnPID::kSpecies; i++) {
178         if (i < 3) {
179             fP[i] = copy.fP[i];
180             fV[i] = copy.fV[i];
181         }
182         fPIDWeight[i] = copy.fPIDWeight[i];
183         fPIDProb[i] = copy.fPIDProb[i];
184     }
185
186     fPIDType = copy.fPIDType;
187     fMass    = copy.fMass;
188
189     // initialize particle object
190     // only if it is present in the template object;
191     // otherwise, it is just cleared and not replaced with anything
192     if (fMCInfo) {
193         delete fMCInfo;
194         fMCInfo = 0x0;
195     }
196     if (copy.fMCInfo) fMCInfo = new AliRsnMCInfo(*(copy.fMCInfo));
197
198     return (*this);
199 }
200
201 //_____________________________________________________________________________
202 AliRsnDaughter::~AliRsnDaughter()
203 {
204 //
205 // Destructor
206 //
207
208     if (fMCInfo) {
209         delete fMCInfo;
210         fMCInfo = 0;
211     }
212 }
213
214 //_____________________________________________________________________________
215 void AliRsnDaughter::SetPIDWeight(Int_t i, Double_t value)
216 {
217 //
218 // I the argument 'i' is in the correct range,
219 // sets the i-th PID weight to 'value'
220 //
221
222     if (i >= 0 && i < AliRsnPID::kSpecies) fPIDWeight[i] = value;
223     else {
224         AliError(Form("Cannot set a weight related to slot %d", i));
225     }
226 }
227
228 //_____________________________________________________________________________
229 void AliRsnDaughter::SetPIDProb(Int_t i, Double_t value)
230 {
231 //
232 // I the argument 'i' is in the correct range,
233 // sets the i-th PID probability to 'value'
234 //
235
236     if (i >= 0 && i < AliRsnPID::kSpecies) fPIDProb[i] = value;
237     else {
238         AliError(Form("Cannot set a weight related to slot %d", i));
239     }
240 }
241
242 //_____________________________________________________________________________
243 void AliRsnDaughter::SetPIDWeights(const Double_t *pid)
244 {
245 //
246 // Sets ALL PID weights at once.
247 // The passed array is supposed to have at least as many
248 // slots as the number of allowed particle species.
249 //
250
251    Int_t i;
252    for (i = 0; i < AliRsnPID::kSpecies; i++) fPIDWeight[i] = pid[i];
253 }
254
255
256 //_____________________________________________________________________________
257 Bool_t AliRsnDaughter::Adopt(AliESDtrack* esdTrack, Bool_t useTPCInnerParam)
258 {
259 //
260 // Copies data from an AliESDtrack into "this":
261 //  - charge sign
262 //  - momentum
263 //  - point of closest approach to primary vertex
264 //  - ESD pid weights
265 // In case of errors returns kFALSE, otherwise kTRUE.
266 //
267
268     if (!esdTrack) {
269         AliError("Passed NULL object: nothing done");
270         return kFALSE;
271     }
272
273     // copy momentum and vertex
274     if (!useTPCInnerParam) {
275         esdTrack->GetPxPyPz(fP);
276         esdTrack->GetXYZ(fV);
277     }
278     else {
279         if (!esdTrack->GetTPCInnerParam()) return kFALSE;
280         esdTrack->GetTPCInnerParam()->GetPxPyPz(fP);
281         esdTrack->GetTPCInnerParam()->GetXYZ(fV);
282     }
283
284     // copy PID weights
285     Int_t    i;
286     Double_t pid[5];
287     if (!useTPCInnerParam) {
288         esdTrack->GetESDpid(pid);
289     }
290     else {
291         esdTrack->GetTPCpid(pid);
292     }
293     for (i = 0; i < 5; i++) fPIDWeight[i] = pid[i];
294
295     // copy flags
296     fFlags = esdTrack->GetStatus();
297
298     // copy charge sign
299     fCharge = (Short_t)esdTrack->Charge();
300
301     return kTRUE;
302 }
303
304
305 //_____________________________________________________________________________
306 Bool_t AliRsnDaughter::Adopt(AliAODTrack* aodTrack)
307 {
308 //
309 // Copies data from an AliAODtrack into "this":
310 //  - charge sign
311 //  - momentum
312 //  - point of closest approach to primary vertex
313 //  - ESD pid weights
314 // In case of errors returns kFALSE, otherwise kTRUE.
315 //
316
317     if (!aodTrack) {
318         AliError("Passed NULL object: nothing done");
319         return kFALSE;
320     }
321
322     // copy momentum  and vertex
323     fP[0] = aodTrack->Px();
324     fP[1] = aodTrack->Py();
325     fP[2] = aodTrack->Pz();
326     fV[0] = aodTrack->Xv();
327     fV[1] = aodTrack->Yv();
328     fV[2] = aodTrack->Zv();
329
330     // copy PID weights
331     Int_t i;
332     for (i = 0; i < 5; i++) fPIDWeight[i] = aodTrack->PID()[i];
333
334     // copy sign
335     fCharge = aodTrack->Charge();
336
337     return kTRUE;
338 }
339
340
341 //_____________________________________________________________________________
342 Bool_t AliRsnDaughter::Adopt(AliMCParticle *mcParticle)
343 {
344 //
345 // Copies data from a MCParticle into "this":
346 //  - charge sign
347 //  - momentum
348 //  - point of closest approach to primary vertex
349 //  - ESD pid weights
350 //  - true PDG code
351 //  - mother
352 // In case of errors returns kFALSE, otherwise kTRUE.
353 //
354
355     if (!mcParticle) {
356         AliError("Passed NULL object: nothing done");
357         return kFALSE;
358     }
359
360         // retrieve the TParticle object from the argument
361         TParticle *particle = mcParticle->Particle();
362         if (!particle) {
363            AliError("AliMCParticle::Particle() returned NULL");
364            return kFALSE;
365     }
366
367     // copy momentum  and vertex
368     fP[0] = particle->Px();
369     fP[1] = particle->Py();
370     fP[2] = particle->Pz();
371     fV[0] = particle->Vx();
372     fV[1] = particle->Vy();
373     fV[2] = particle->Vz();
374
375     // recognize charge sign from PDG code sign
376     Int_t pdg = particle->GetPdgCode();
377     Int_t absPDG = TMath::Abs(pdg);
378     if (absPDG <= 15) {
379         if (pdg > 0) fCharge = -1; else fCharge = 1;
380     }
381     else if (absPDG < 3000) {
382         if (pdg > 0) fCharge = 1; else fCharge = -1;
383     }
384     else {
385         fCharge = 0;
386         return kFALSE;
387     }
388
389     // identify track perfectly using PDG code
390     fPIDType = AliRsnPID::InternalType(pdg);
391     fMass = AliRsnPID::ParticleMass(fPIDType);
392
393     // flags and PID weights make no sense with MC tracks
394     fFlags = 0;
395     for (pdg = 0; pdg < AliRsnPID::kSpecies; pdg++) fPIDWeight[pdg] = 0.0;
396     fPIDWeight[fPIDType] = 1.0;
397
398     // copy other MC info (mother PDG code cannot be retrieved here)
399     InitMCInfo(particle);
400
401     return kTRUE;
402 }
403
404 //_____________________________________________________________________________
405 void AliRsnDaughter::Print(Option_t *option) const
406 {
407 //
408 // Prints the values of data members, using the options:
409 // - P --> momentum
410 // - V --> DCA vertex
411 // - C --> electric charge
412 // - F --> flags
413 // - I --> identification (PID, probability and mass)
414 // - W --> PID weights
415 // - M --> Montecarlo (from AliRsnMCInfo)
416 // - ALL --> All oprions switched on
417 //
418 // Index and label are printed by default.
419 //
420
421     TString opt(option);
422     opt.ToUpper();
423
424     cout << ".......Index            : " << fIndex << endl;
425     cout << ".......Label            : " << fLabel << endl;
426
427     if (opt.Contains("P") || opt.Contains("ALL")) {
428         cout << ".......Px, Py, Pz, Pt   : " << Px() << ' ' << Py() << ' ' << Pz() << ' ' << Pt() << endl;
429     }
430     if (opt.Contains("V") || opt.Contains("ALL")) {
431         cout << ".......Vx, Vy, Vz       : " << Xv() << ' ' << Yv() << ' ' << Zv() << endl;
432     }
433     if (opt.Contains("C") || opt.Contains("ALL")) {
434         cout << ".......Charge           : " << fCharge << endl;
435     }
436     if (opt.Contains("F") || opt.Contains("ALL")) {
437         cout << ".......Flags            : " << fFlags << endl;
438     }
439     if (opt.Contains("I") || opt.Contains("ALL")) {
440         cout << ".......PID              : " << AliRsnPID::ParticleName(fPIDType) << endl;
441         if (fPIDType > 0 && fPIDType < AliRsnPID::kSpecies) {
442             cout << ".......PID probability  : " << fPIDProb[fPIDType] << endl;
443         }
444         cout << ".......Mass             : " << fMass << endl;
445     }
446     if (opt.Contains("W") || opt.Contains("ALL")) {
447         cout << ".......Weights          : ";
448         Int_t i;
449         for (i = 0; i < AliRsnPID::kSpecies; i++) cout << fPIDWeight[i] << ' ';
450         cout << endl;
451     }
452     if (opt.Contains("M") || opt.Contains("ALL")) {
453         if (fMCInfo) {
454             cout << ".......PDG code         : " << fMCInfo->PDG() << endl;
455             cout << ".......Mother (label)   : " << fMCInfo->Mother() << endl;
456             cout << ".......Mother (PDG code): " << fMCInfo->MotherPDG() << endl;
457         }
458         else {
459             cout << ".......MC info not present" << endl;
460         }
461     }
462 }
463
464 //_____________________________________________________________________________
465 void AliRsnDaughter::InitMCInfo()
466 {
467 //
468 // Initializes the particle object with default constructor.
469 //
470
471     fMCInfo = new AliRsnMCInfo;
472 }
473
474 //_____________________________________________________________________________
475 Bool_t AliRsnDaughter::InitMCInfo(TParticle *particle)
476 {
477 //
478 // Copies data from an MC particle into the object which
479 // contains all MC details taken from kinematics info.
480 // If requested by second argument, momentum and vertex
481 // of the Particle are copied into the 'fP' and 'fV'
482 // data members, to simulate a perfect reconstruction.
483 // If something goes wrong, returns kFALSE,
484 // otherwise returns kTRUE.
485 //
486
487     // retrieve the TParticle object pointed by this MC track
488     if (!particle) {
489         AliError("Passed NULL particle object");
490         return kFALSE;
491     }
492
493     // initialize object if not initialized yet
494     if (fMCInfo) delete fMCInfo;
495     fMCInfo = new AliRsnMCInfo;
496     fMCInfo->Adopt(particle);
497
498     return kTRUE;
499 }
500
501 //_____________________________________________________________________________
502 Bool_t AliRsnDaughter::InitMCInfo(AliMCParticle *mcParticle)
503 {
504 //
505 // Copies data from an MC particle into the object which
506 // contains all MC details taken from kinematics info.
507 // If requested by second argument, momentum and vertex
508 // of the Particle are copied into the 'fP' and 'fV'
509 // data members, to simulate a perfect reconstruction.
510 // If something goes wrong, returns kFALSE,
511 // otherwise returns kTRUE.
512 //
513
514     // retrieve the TParticle object pointed by this MC track
515     TParticle *particle = mcParticle->Particle();
516     return InitMCInfo(particle);
517 }