Updated definition of theta* and corrected a bug in selection of true pairs from MC
[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 "AliRsnDaughter.h"
15 #include "AliRsnPairDef.h"
16 #include "AliRsnMother.h"
17
18 ClassImp(AliRsnMother)
19
20 //_____________________________________________________________________________
21 AliRsnMother::AliRsnMother() : 
22   fUseMC(kFALSE),
23   fDefaultMass(0.0),
24   fSum(),
25   fSumMC(),
26   fRef(),
27   fRefMC()
28 {
29 //
30 // Constructor.
31 // Initializes all variables to meaningless values.
32 //
33
34   Int_t i;
35   for (i = 0; i < 2; i++) fDaughter[i] = 0x0;
36 }
37
38 //_____________________________________________________________________________
39 AliRsnMother::AliRsnMother(const AliRsnMother &obj) : 
40   TObject(obj), 
41   fUseMC(obj.fUseMC),
42   fDefaultMass(obj.fDefaultMass),
43   fSum(obj.fSum),
44   fSumMC(obj.fSumMC),
45   fRef(obj.fRef),
46   fRefMC(obj.fRefMC)
47 {
48 //
49 // Copy constructor.
50 // Initializes all variables to copy values.
51 // Does not duplicate pointers.
52 //
53
54   Int_t i;
55   for (i = 0; i < 2; i++) fDaughter[i] = obj.fDaughter[i];
56 }
57
58 //_____________________________________________________________________________
59 AliRsnMother& AliRsnMother::operator=(const AliRsnMother &obj)
60 {
61 //
62 // Assignment operator.
63 // Initializes all variables to copy values.
64 // Does not duplicate pointers.
65 //
66
67   Int_t i;
68   
69   fDefaultMass = obj.fDefaultMass;
70   fSum = obj.fSum;
71   fRef = obj.fRef;
72   fSumMC = obj.fSumMC;
73   fRefMC = obj.fRefMC;
74   
75   for (i = 0; i < 2; i++) fDaughter[i] = obj.fDaughter[i];
76
77   return (*this);
78 }
79
80 //_____________________________________________________________________________
81 AliRsnMother::~AliRsnMother()
82 {
83 //
84 // Desctructor.
85 // Does nothing, since pointers are not created in this class.
86 //
87 }
88
89 //_____________________________________________________________________________
90 Int_t AliRsnMother::CommonMother() const
91 {
92 //
93 // Checks if the two tracks in the pair have the same mother.
94 // This can be known if MC info is present.
95 // If the mother label is the same, rhe PDG code of the mother is returned,
96 // otherwise the method returns 0.
97 //
98
99   // if MC info is not available, the pairs is not true by default
100   if (!fDaughter[0]->GetRefMC() || !fDaughter[1]->GetRefMC()) 
101   {
102     AliInfo("Cannot know if the pairs is true or not because MC Info is not present");
103     return 0;
104   }
105
106   // check that labels are the same
107   if (fDaughter[0]->GetParticle()->GetFirstMother() != fDaughter[1]->GetParticle()->GetFirstMother())
108     return 0;
109
110   // if we reach this point, the two tracks have the same mother
111   // let's check now the PDG code of this common mother
112   return TMath::Abs(fDaughter[0]->GetMotherPDG());
113 }
114
115 //_____________________________________________________________________________
116 void AliRsnMother::SetDaughters
117 (AliRsnDaughter *d0, Double_t mass0, AliRsnDaughter *d1, Double_t mass1)
118 {
119 //
120 // Sets the pair defined in this usind tso passed daughters and two masses
121 // which will be assigned to them, in order to recompute their 4-momenta
122 // and sum them into the datamembers of this object.
123 //
124
125   if (d0) fDaughter[0] = d0;
126   if (d1) fDaughter[1] = d1;
127   
128   if (!d0 || !d1) return;
129   
130   fDaughter[0]->SetMass(mass0);
131   fDaughter[1]->SetMass(mass1);
132   
133   fSum   = fDaughter[0]->P(kFALSE) + fDaughter[1]->P(kFALSE);
134   fSumMC = fDaughter[0]->P(kTRUE)  + fDaughter[1]->P(kTRUE);
135   
136   fRef  .SetXYZM(fSum  .X(), fSum  .Y(), fSum  .Z(), fDefaultMass);
137   fRefMC.SetXYZM(fSumMC.X(), fSumMC.Y(), fSumMC.Z(), fDefaultMass);
138 }
139
140 //_____________________________________________________________________________
141 void AliRsnMother::ResetPair()
142 {
143 //
144 // Resets the mother, nullifying all data members
145 //
146
147   Int_t i;
148   for (i = 0; i < 2; i++) fDaughter[i] = 0x0;
149   
150   fSum  .SetXYZM(0.0, 0.0, 0.0, 0.0);
151   fRef  .SetXYZM(0.0, 0.0, 0.0, 0.0);
152   fSumMC.SetXYZM(0.0, 0.0, 0.0, 0.0);
153   fRefMC.SetXYZM(0.0, 0.0, 0.0, 0.0);
154 }
155
156 //_____________________________________________________________________________
157 Double_t AliRsnMother::ThetaStar(Bool_t first, Bool_t useMC)
158 {
159 //
160 // Returns the theta* as the angle of the first daughter
161 // w.r. to the mother momentum, in its rest frame
162 //
163
164   TLorentzVector &mother   = (useMC ? fSumMC : fSum);
165   TLorentzVector &daughter = (first ? fDaughter[0]->P() : fDaughter[1]->P());
166   
167   Double_t theta = mother.Vect().Theta();
168   Double_t phi   = mother.Vect().Phi();
169   Double_t beta  = mother.Beta();
170   
171   Double_t betax = beta * TMath::Sin(theta) * TMath::Cos(phi);
172   Double_t betay = beta * TMath::Sin(theta) * TMath::Sin(phi);
173   Double_t betaz = beta * TMath::Cos(theta);
174   Double_t gamx  = 1.0 / TMath::Sqrt(1.0 - betax*betax);
175   Double_t gamy  = 1.0 / TMath::Sqrt(1.0 - betay*betay);
176   Double_t gamz  = 1.0 / TMath::Sqrt(1.0 - betaz*betaz);
177   
178   Double_t pstarx = gamx * (daughter.Vect().X() - betax * daughter.E());
179   Double_t pstary = gamy * (daughter.Vect().Y() - betay * daughter.E());
180   Double_t pstarz = gamz * (daughter.Vect().Z() - betaz * daughter.E());
181   
182   TVector3 dstar(pstarx, pstary, pstarz);
183   TVector3 vnorm(mother.Vect().Y() / mother.Perp(), mother.Vect().X() / mother.Perp(), 0.0);
184   
185   Double_t cosThetaStar = dstar.Dot(vnorm) / dstar.Mag();
186   
187   return cosThetaStar;
188 }
189
190 //_____________________________________________________________________________
191 void AliRsnMother::PrintInfo(const Option_t * /*option*/) const
192 {
193 //
194 // Print some info of the pair.
195 // The options are passed to the AliRsnDaughter::Print() method
196 //
197
198   AliInfo("======== BEGIN PAIR INFO ===========");
199   AliInfo("Track #1");
200   fDaughter[0]->Print();
201   AliInfo("Track #2");
202   fDaughter[1]->Print();
203   AliInfo("========= END PAIR INFO ===========");
204 }
205
206 //_____________________________________________________________________________
207 Bool_t AliRsnMother::CheckPair() const
208 {
209 //
210 // Checks that the pair is well initialized:
211 // - both daughters are good pointers
212 // - if MC is required, both daughters have a MC reference
213 //
214
215   if (!fDaughter[0] || !fDaughter[1]) 
216   {
217     AliError("One of the two tracks is NULL in this pair!");
218     return kFALSE;
219   }
220   
221   if (fUseMC)
222   {
223     if (fDaughter[0]->GetRefMC() == 0x0 || fDaughter[1]->GetRefMC() == 0x0)
224     {
225       AliError("Required MC info but not all MC refs are available");
226       return kFALSE;
227     }
228   }
229   
230   return kTRUE;
231 }
232
233 //_____________________________________________________________________________
234 Bool_t AliRsnMother::MatchesDef(AliRsnPairDef *def)
235 {
236 //
237 // Checks if the daughters, in any order, do match a given decay channel,
238 // using the specified identification method, which can be the 'true' one
239 // or the 'realistic' one only.
240 //
241
242   if (!def) return kFALSE;
243   if (!fDaughter[0]->GetRefMC()) return kFALSE;
244   if (!fDaughter[1]->GetRefMC()) return kFALSE;
245
246   Bool_t decayMatch = kFALSE;
247   Int_t  pdg[2], ref[2];
248   pdg[0] = fDaughter[0]->GetRefMC()->Particle()->GetPdgCode();
249   pdg[1] = fDaughter[1]->GetRefMC()->Particle()->GetPdgCode();
250   ref[0] = TMath::Abs(AliPID::ParticleCode(def->GetPID(0)));
251   ref[1] = TMath::Abs(AliPID::ParticleCode(def->GetPID(1)));
252
253   // check #1:
254   // if first member of pairDef has same sign as first member of this,
255   // daughter[0] perfect PID must match first member of pairDef
256   // daughter[1] perfect PID must march second member of pairDef
257   if (fDaughter[0]->IsSign(def->GetCharge(0)) && fDaughter[1]->IsSign(def->GetCharge(1))) 
258   {
259     decayMatch = (pdg[0] == ref[0] && pdg[1] == ref[1]);
260   }
261
262   // check #2:
263   // if first member of pairDef has same sign as second member of this,
264   // daughter[0] perfect PID must match second member of pairDef
265   // daughter[1] perfect PID must march first member of pairDef
266   if (fDaughter[1]->IsSign(def->GetCharge(0)) && fDaughter[0]->IsSign(def->GetCharge(1))) 
267   {
268     decayMatch = (pdg[0] == ref[1] && pdg[1] == ref[0]);
269   }
270
271   return (decayMatch && (CommonMother() == def->GetMotherPDG()));
272 }