1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 ////////////////////////////////////////////////////////////////////////////////
18 // This class implements a candidate resonance. It has two pointers to its
19 // two candidate daughters, whose 4-momenta are combined to obtain the mother
20 // invariant mass and other kinematical quantities.
21 // This class contains also some methods used to compute kinematical relations
22 // between the candidate resonance and other particles.
24 // authors: A. Pulvirenti (alberto.pulvirenti@ct.infn.it)
25 // M. Vala (martin.vala@cern.ch)
27 ////////////////////////////////////////////////////////////////////////////////
29 #include <Riostream.h>
32 #include "AliAODMCParticle.h"
33 #include "AliMCParticle.h"
34 #include "AliRsnEvent.h"
36 #include "AliRsnMother.h"
38 ClassImp(AliRsnMother)
40 //__________________________________________________________________________________________________
41 AliRsnMother::AliRsnMother(const AliRsnMother &obj) :
43 fRefEvent(obj.fRefEvent),
48 fDCAproduct(obj.fDCAproduct)
52 // Does not duplicate pointers.
56 for (i = 0; i < 2; i++) fDaughter[i] = obj.fDaughter[i];
59 //__________________________________________________________________________________________________
60 AliRsnMother &AliRsnMother::operator=(const AliRsnMother &obj)
63 // Assignment operator.
64 // Does not duplicate pointers.
73 fRefEvent = obj.fRefEvent;
74 fDaughter[0] = obj.fDaughter[0];
75 fDaughter[1] = obj.fDaughter[1];
76 fDCAproduct = obj.fDCAproduct;
81 //__________________________________________________________________________________________________
82 AliRsnMother::~AliRsnMother()
86 // Does nothing, since pointers are not created in this class.
90 //_______________________________________________________________________________________________________________________
91 void AliRsnMother::Reset()
94 // Resets the mother, zeroing all data members.
98 for (i = 0; i < 2; i++) fDaughter[i] = 0x0;
100 fSum.SetXYZM(0.0, 0.0, 0.0, 0.0);
101 fRef.SetXYZM(0.0, 0.0, 0.0, 0.0);
104 //__________________________________________________________________________________________________
105 Int_t AliRsnMother::CommonMother() const
108 // If MC info is available, checks if the two tracks in the pair have the same mother.
109 // If the mother label is the same, the function returns the PDG code of mother,
110 // otherwise it returns 0.
111 // The two arguments passed by reference contain the GEANT labels of the mother
112 // of the two particles to which the two daughters point. This is for being able
113 // to check if they are really coming from a resonance (indexes >= 0) or not.
116 Int_t m1 = fDaughter[0]->GetMother();
117 Int_t m2 = fDaughter[1]->GetMother();
120 // a true mother makes sense only if both mothers
121 // are not-negative and equal
122 if (m1 >= 0 && m2 >= 0 && m1 == m2) {
123 out = TMath::Abs(fDaughter[0]->GetMotherPDG());
124 AliDebugClass(1, Form("Mothers are: %d %d --> EQUAL (PDG = %d)", m1, m2, out));
130 //__________________________________________________________________________________________________
131 Double_t AliRsnMother::AngleToLeading(Bool_t &success)
134 // Compute the angle betwee this and the leading particls
135 // of the reference event (if this was set properly).
136 // In case one of the two daughters is the leading, return
137 // a meaningless value, in order to skip this pair.
138 // if second argument is kTRUE, use MC values.
146 Int_t id1 = fDaughter[0]->GetID();
147 Int_t id2 = fDaughter[1]->GetID();
148 Int_t idL = fRefEvent->GetLeadingIndex();
150 if (id1 == idL || id2 == idL) {
155 AliRsnDaughter leading = fRefEvent->GetDaughter(idL, kFALSE);
156 AliVParticle *ref = leading.GetRef();
157 Double_t angle = ref->Phi() - fSum.Phi();
159 //return angle w.r.t. leading particle in the range -pi/2, 3/2pi
160 while (angle >= 1.5 * TMath::Pi()) angle -= 2 * TMath::Pi();
161 while (angle < -0.5 * TMath::Pi()) angle += 2 * TMath::Pi();
167 //__________________________________________________________________________________________________
168 void AliRsnMother::ComputeSum(Double_t mass0, Double_t mass1, Double_t motherMass)
171 // Sets the masses for the 4-momenta of the daughters and then
172 // sums them, taking into account that the space part is set to
173 // each of them when the reference object is set (see AliRsnDaughter::SetRef)
176 fDaughter[0]->FillP(mass0);
177 fDaughter[1]->FillP(mass1);
180 fSum = fDaughter[0]->Prec() + fDaughter[1]->Prec();
181 fSumMC = fDaughter[0]->Psim() + fDaughter[1]->Psim();
184 fRef.SetXYZM(fSum.X(), fSum.Y(), fSum.Z(), motherMass);
185 fRefMC.SetXYZM(fSumMC.X(), fSumMC.Y(), fSumMC.Z(), motherMass);
188 //__________________________________________________________________________________________________
189 Double_t AliRsnMother::CosThetaStar(Bool_t first, Bool_t useMC)
192 // Computes the cosine of theta*, which is the angle of one of the daughters
193 // with respect to the total momentum of the resonance, in its rest frame.
194 // The arguments are needed to choose which of the daughters one want to use
195 // and if reconstructed or MC momentum must be used.
196 // [Contribution from Z. Feckova]
199 TLorentzVector &mother = Sum(useMC);
200 TLorentzVector &daughter0 = (first ? fDaughter[0]->P(useMC) : fDaughter[1]->P(useMC));
201 // TLorentzVector &daughter1 = (first ? fDaughter[1]->P(useMC) : fDaughter[0]->P(useMC));
202 TVector3 momentumM(mother.Vect());
203 TVector3 normal(mother.Y() / momentumM.Mag(), -mother.X() / momentumM.Mag(), 0.0);
205 // // Computes first the invariant mass of the mother
206 // Double_t mass0 = fDaughter[0]->P(useMC).M();
207 // Double_t mass1 = fDaughter[1]->P(useMC).M();
208 // Double_t p0 = daughter0.Vect().Mag();
209 // Double_t p1 = daughter1.Vect().Mag();
210 // Double_t E0 = TMath::Sqrt(mass0 * mass0 + p0 * p0);
211 // Double_t E1 = TMath::Sqrt(mass1 * mass1 + p1 * p1);
212 // Double_t MotherMass = TMath::Sqrt((E0 + E1) * (E0 + E1) - (p0 * p0 + 2.0 * daughter0.Vect().Dot(daughter1.Vect()) + p1 * p1));
213 // MotherMass = fSum.M();
215 // Computes components of beta
216 Double_t betaX = -mother.X() / mother.E();
217 Double_t betaY = -mother.Y() / mother.E();
218 Double_t betaZ = -mother.Z() / mother.E();
220 // Computes Lorentz transformation of the momentum of the first daughter
221 // into the rest frame of the mother and theta*
222 daughter0.Boost(betaX, betaY, betaZ);
223 TVector3 momentumD = daughter0.Vect();
225 Double_t cosThetaStar = normal.Dot(momentumD) / momentumD.Mag();
230 //__________________________________________________________________________________________________
231 Double_t AliRsnMother::DCAproduct()
234 // returns product of DCA of the two daughters
236 AliRsnEvent *event1 = fDaughter[0]->GetOwnerEvent();
237 AliRsnEvent *event2 = fDaughter[1]->GetOwnerEvent();
239 if(event1 != event2){
240 AliError("Attempting to build pair with tracks coming from different events");
244 if (event1->IsAOD()) {
245 AliAODEvent *aodEvent = (AliAODEvent*)event1->GetRefAOD();
246 if (!aodEvent) return 0.0;
247 AliAODTrack *track1 = (AliAODTrack*)fDaughter[0]->Ref2AODtrack();
248 AliAODTrack *track2 = (AliAODTrack*)fDaughter[1]->Ref2AODtrack();
249 AliVVertex *vertex = aodEvent->GetPrimaryVertex();
250 if (!vertex || !track1 || !track2) return 0.0;
252 Double_t b1[2], cov1[3], b2[2], cov2[3];
253 if ( track1->PropagateToDCA(vertex, aodEvent->GetMagneticField(), kVeryBig, b1, cov1) &&
254 track2->PropagateToDCA(vertex, aodEvent->GetMagneticField(), kVeryBig, b2, cov2) )
255 fDCAproduct = b1[0]*b2[0];
256 else fDCAproduct = 0.0;
258 AliESDEvent *esdEvent = (AliESDEvent*)event1->GetRefESD();
259 if (!esdEvent) return 0.0;
260 AliESDtrack *track1 = (AliESDtrack*)fDaughter[0]->Ref2ESDtrack();
261 AliESDtrack *track2 = (AliESDtrack*)fDaughter[1]->Ref2ESDtrack();
262 const AliVVertex *vertex = esdEvent->GetPrimaryVertex();
263 if (!vertex || !track1 || !track2) return 0.0;
265 Double_t b1[2], cov1[3], b2[2], cov2[3];
266 if ( track1->PropagateToDCA(vertex, esdEvent->GetMagneticField(), kVeryBig, b1, cov1) &&
267 track2->PropagateToDCA(vertex, esdEvent->GetMagneticField(), kVeryBig, b2, cov2) )
268 fDCAproduct = b1[0]*b2[0];
269 else fDCAproduct = 0.0;
274 //__________________________________________________________________________________________________
275 void AliRsnMother::PrintInfo(const Option_t * /*option*/) const
278 // Print some info of the pair.
279 // The options are passed to the AliRsnDaughter::Print() method
282 AliInfo("======== BEGIN PAIR INFO ===========");
284 fDaughter[0]->Print();
286 fDaughter[1]->Print();
287 AliInfo("========= END PAIR INFO ===========");
290 //__________________________________________________________________________________________________
291 Bool_t AliRsnMother::CheckPair(Bool_t checkMC) const
294 // Checks that the pair is well initialized:
295 // - both daughters are good pointers
296 // - if MC is required, both daughters have a MC reference
299 if (!fDaughter[0] || !fDaughter[1]) {
300 AliError("One of the two tracks is NULL in this pair!");
305 if (fDaughter[0]->GetRefMC() == 0x0 || fDaughter[1]->GetRefMC() == 0x0) {
306 AliError("Required MC info but not all MC refs are available");