Generalized the MC reference object in order to cope with AOD format
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnMother.cxx
1 //
2 // Class AliRsnMother
3 //
4 // Implementation of a pair of tracks, for several purposes
5 // - computing the total 4-momentum & inv. mass for output histos filling
6 // - evaluating cut checks on the pair of particles
7 // - evaluating any kind of kinematic value over their sum
8 //
9 // authors: Martin Vala (martin.vala@cern.ch)
10 //          Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
11 //
12 #include <Riostream.h>
13 #include <TVector3.h>
14 #include "AliAODMCParticle.h"
15 #include "AliMCParticle.h"
16 #include "AliRsnDaughter.h"
17 #include "AliRsnPairDef.h"
18 #include "AliRsnMother.h"
19
20 ClassImp(AliRsnMother)
21
22 //_____________________________________________________________________________
23 AliRsnMother::AliRsnMother() : 
24   fUseMC(kFALSE),
25   fDefaultMass(0.0),
26   fSum(),
27   fSumMC(),
28   fRef(),
29   fRefMC()
30 {
31 //
32 // Constructor.
33 // Initializes all variables to meaningless values.
34 //
35
36   Int_t i;
37   for (i = 0; i < 2; i++) fDaughter[i] = 0x0;
38 }
39
40 //_____________________________________________________________________________
41 AliRsnMother::AliRsnMother(const AliRsnMother &obj) : 
42   TObject(obj), 
43   fUseMC(obj.fUseMC),
44   fDefaultMass(obj.fDefaultMass),
45   fSum(obj.fSum),
46   fSumMC(obj.fSumMC),
47   fRef(obj.fRef),
48   fRefMC(obj.fRefMC)
49 {
50 //
51 // Copy constructor.
52 // Initializes all variables to copy values.
53 // Does not duplicate pointers.
54 //
55
56   Int_t i;
57   for (i = 0; i < 2; i++) fDaughter[i] = obj.fDaughter[i];
58 }
59
60 //_____________________________________________________________________________
61 AliRsnMother& AliRsnMother::operator=(const AliRsnMother &obj)
62 {
63 //
64 // Assignment operator.
65 // Initializes all variables to copy values.
66 // Does not duplicate pointers.
67 //
68
69   Int_t i;
70   
71   fDefaultMass = obj.fDefaultMass;
72   fSum = obj.fSum;
73   fRef = obj.fRef;
74   fSumMC = obj.fSumMC;
75   fRefMC = obj.fRefMC;
76   
77   for (i = 0; i < 2; i++) fDaughter[i] = obj.fDaughter[i];
78
79   return (*this);
80 }
81
82 //_____________________________________________________________________________
83 AliRsnMother::~AliRsnMother()
84 {
85 //
86 // Desctructor.
87 // Does nothing, since pointers are not created in this class.
88 //
89 }
90
91 //_____________________________________________________________________________
92 Int_t AliRsnMother::CommonMother() const
93 {
94 //
95 // Checks if the two tracks in the pair have the same mother.
96 // This can be known if MC info is present.
97 // If the mother label is the same, rhe PDG code of the mother is returned,
98 // otherwise the method returns 0.
99 //
100
101   // if MC info is not available, the pairs is not true by default
102   if (!fDaughter[0]->GetRefMC() || !fDaughter[1]->GetRefMC()) 
103   {
104     AliInfo("Cannot know if the pairs is true or not because MC Info is not present");
105     return 0;
106   }
107
108   // check that labels are the same
109   if (fDaughter[0]->GetParticle()->GetFirstMother() != fDaughter[1]->GetParticle()->GetFirstMother())
110     return 0;
111
112   // if we reach this point, the two tracks have the same mother
113   // let's check now the PDG code of this common mother
114   return TMath::Abs(fDaughter[0]->GetMotherPDG());
115 }
116
117 //_____________________________________________________________________________
118 void AliRsnMother::SetDaughters
119 (AliRsnDaughter *d0, Double_t mass0, AliRsnDaughter *d1, Double_t mass1)
120 {
121 //
122 // Sets the pair defined in this usind tso passed daughters and two masses
123 // which will be assigned to them, in order to recompute their 4-momenta
124 // and sum them into the datamembers of this object.
125 //
126
127   if (d0) fDaughter[0] = d0;
128   if (d1) fDaughter[1] = d1;
129   
130   if (!d0 || !d1) return;
131   
132   fDaughter[0]->SetMass(mass0);
133   fDaughter[1]->SetMass(mass1);
134   
135   fSum   = fDaughter[0]->P(kFALSE) + fDaughter[1]->P(kFALSE);
136   fSumMC = fDaughter[0]->P(kTRUE)  + fDaughter[1]->P(kTRUE);
137   
138   fRef  .SetXYZM(fSum  .X(), fSum  .Y(), fSum  .Z(), fDefaultMass);
139   fRefMC.SetXYZM(fSumMC.X(), fSumMC.Y(), fSumMC.Z(), fDefaultMass);
140 }
141
142 //_____________________________________________________________________________
143 void AliRsnMother::ResetPair()
144 {
145 //
146 // Resets the mother, nullifying all data members
147 //
148
149   Int_t i;
150   for (i = 0; i < 2; i++) fDaughter[i] = 0x0;
151   
152   fSum  .SetXYZM(0.0, 0.0, 0.0, 0.0);
153   fRef  .SetXYZM(0.0, 0.0, 0.0, 0.0);
154   fSumMC.SetXYZM(0.0, 0.0, 0.0, 0.0);
155   fRefMC.SetXYZM(0.0, 0.0, 0.0, 0.0);
156 }
157
158 //_____________________________________________________________________________
159 Double_t AliRsnMother::CosThetaStar(Bool_t first, Bool_t useMC)
160 {
161   TLorentzVector mother    = (useMC ? fSumMC : fSum);
162   TLorentzVector daughter0 = (first ? fDaughter[0]->P() : fDaughter[1]->P());
163   TLorentzVector daughter1 = (first ? fDaughter[1]->P() : fDaughter[0]->P());
164   TVector3 momentumM          (mother.Vect());
165   //cout << mother.Vect().X() << ' ' << mother.Vect().Y() << ' ' << mother.Vect().Z() << endl;
166   TVector3 normal             (mother.Y()/momentumM.Mag(), -mother.X()/momentumM.Mag(), 0.0);
167   //cout << momentumM.Mag() << ' ' << normal.X() << ' ' << normal.Y() << ' ' << normal.Z() << endl;
168
169 // Computes first the invariant mass of the mother
170
171   Double_t mass0            = fDaughter[0]->P().M();
172   Double_t mass1            = fDaughter[1]->P().M();
173   Double_t p0               = daughter0.Vect().Mag();
174   Double_t p1               = daughter1.Vect().Mag();
175   Double_t E0               = TMath::Sqrt(mass0*mass0+p0*p0);
176   Double_t E1               = TMath::Sqrt(mass1*mass1+p1*p1);
177   Double_t MotherMass       = TMath::Sqrt((E0+E1)*(E0+E1)-(p0*p0+2.0*daughter0.Vect().Dot(daughter1.Vect())+p1*p1));
178   MotherMass = fSum.M();
179
180 // Computes components of beta
181
182   Double_t betaX = -mother.X()/mother.E(); //TMath::Sqrt(momentumM.Mag()*momentumM.Mag()+MotherMass*MotherMass);
183   Double_t betaY = -mother.Y()/mother.E(); //TMath::Sqrt(momentumM.Mag()*momentumM.Mag()+MotherMass*MotherMass);
184   Double_t betaZ = -mother.Z()/mother.E(); //TMath::Sqrt(momentumM.Mag()*momentumM.Mag()+MotherMass*MotherMass);
185
186 // Computes Lorentz transformation of the momentum of the first daughter
187 // into the rest frame of the mother and theta*
188   //cout << "Beta = " << betaX << ' ' << betaY << ' ' << betaZ << endl;
189
190   daughter0.Boost(betaX,betaY,betaZ);
191   TVector3 momentumD = daughter0.Vect();
192
193   Double_t cosThetaStar = normal.Dot(momentumD)/momentumD.Mag();
194
195   return cosThetaStar;
196
197 }
198
199 //_____________________________________________________________________________
200 void AliRsnMother::PrintInfo(const Option_t * /*option*/) const
201 {
202 //
203 // Print some info of the pair.
204 // The options are passed to the AliRsnDaughter::Print() method
205 //
206
207   AliInfo("======== BEGIN PAIR INFO ===========");
208   AliInfo("Track #1");
209   fDaughter[0]->Print();
210   AliInfo("Track #2");
211   fDaughter[1]->Print();
212   AliInfo("========= END PAIR INFO ===========");
213 }
214
215 //_____________________________________________________________________________
216 Bool_t AliRsnMother::CheckPair() const
217 {
218 //
219 // Checks that the pair is well initialized:
220 // - both daughters are good pointers
221 // - if MC is required, both daughters have a MC reference
222 //
223
224   if (!fDaughter[0] || !fDaughter[1]) 
225   {
226     AliError("One of the two tracks is NULL in this pair!");
227     return kFALSE;
228   }
229   
230   if (fUseMC)
231   {
232     if (fDaughter[0]->GetRefMC() == 0x0 || fDaughter[1]->GetRefMC() == 0x0)
233     {
234       AliError("Required MC info but not all MC refs are available");
235       return kFALSE;
236     }
237   }
238   
239   return kTRUE;
240 }
241
242 //_____________________________________________________________________________
243 Bool_t AliRsnMother::MatchesDef(AliRsnPairDef *def)
244 {
245 //
246 // Checks if the daughters, in any order, do match a given decay channel,
247 // using the specified identification method, which can be the 'true' one
248 // or the 'realistic' one only.
249 //
250
251   if (!def) return kFALSE;
252   if (!fDaughter[0]->GetRefMC()) return kFALSE;
253   if (!fDaughter[1]->GetRefMC()) return kFALSE;
254
255   Bool_t decayMatch = kFALSE;
256   Int_t  pdg[2], ref[2];
257   pdg[0] = fDaughter[0]->GetPDG();
258   pdg[1] = fDaughter[1]->GetPDG();
259   ref[0] = TMath::Abs(AliPID::ParticleCode(def->GetPID(0)));
260   ref[1] = TMath::Abs(AliPID::ParticleCode(def->GetPID(1)));
261
262   // check #1:
263   // if first member of pairDef has same sign as first member of this,
264   // daughter[0] perfect PID must match first member of pairDef
265   // daughter[1] perfect PID must march second member of pairDef
266   if (fDaughter[0]->IsSign(def->GetCharge(0)) && fDaughter[1]->IsSign(def->GetCharge(1))) 
267   {
268     decayMatch = (pdg[0] == ref[0] && pdg[1] == ref[1]);
269   }
270
271   // check #2:
272   // if first member of pairDef has same sign as second member of this,
273   // daughter[0] perfect PID must match second member of pairDef
274   // daughter[1] perfect PID must march first member of pairDef
275   if (fDaughter[1]->IsSign(def->GetCharge(0)) && fDaughter[0]->IsSign(def->GetCharge(1))) 
276   {
277     decayMatch = (pdg[0] == ref[1] && pdg[1] == ref[0]);
278   }
279
280   return (decayMatch && (CommonMother() == def->GetMotherPDG()));
281 }