/************************************************************************** * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * * * * Author: The ALICE Off-line Project. * * Contributors are mentioned in the code where appropriate. * * * * Permission to use, copy, modify and distribute this software and its * * documentation strictly for non-commercial purposes is hereby granted * * without fee, provided that the above copyright notice appears in all * * copies and that both the copyright notice and this permission notice * * appear in the supporting documentation. The authors make no claims * * about the suitability of this software for any purpose. It is * * provided "as is" without express or implied warranty. * **************************************************************************/ //////////////////////////////////////////////////////////////////////////////// // // This class implements a candidate resonance. It has two pointers to its // two candidate daughters, whose 4-momenta are combined to obtain the mother // invariant mass and other kinematical quantities. // This class contains also some methods used to compute kinematical relations // between the candidate resonance and other particles. // // authors: A. Pulvirenti (alberto.pulvirenti@ct.infn.it) // M. Vala (martin.vala@cern.ch) // //////////////////////////////////////////////////////////////////////////////// #include #include #include "AliAODMCParticle.h" #include "AliMCParticle.h" #include "AliRsnEvent.h" #include "AliRsnMother.h" ClassImp(AliRsnMother) //__________________________________________________________________________________________________ AliRsnMother::AliRsnMother(const AliRsnMother &obj) : TObject(obj), fRefEvent(obj.fRefEvent), fSum(obj.fSum), fSumMC(obj.fSumMC), fRef(obj.fRef), fRefMC(obj.fRefMC), fDCAproduct(obj.fDCAproduct) { // // Copy constructor. // Does not duplicate pointers. // Int_t i; for (i = 0; i < 2; i++) fDaughter[i] = obj.fDaughter[i]; } //__________________________________________________________________________________________________ AliRsnMother &AliRsnMother::operator=(const AliRsnMother &obj) { // // Assignment operator. // Does not duplicate pointers. // if (this == &obj) return *this; fSum = obj.fSum; fRef = obj.fRef; fSumMC = obj.fSumMC; fRefMC = obj.fRefMC; fRefEvent = obj.fRefEvent; fDaughter[0] = obj.fDaughter[0]; fDaughter[1] = obj.fDaughter[1]; fDCAproduct = obj.fDCAproduct; return (*this); } //__________________________________________________________________________________________________ AliRsnMother::~AliRsnMother() { // // Desctructor. // Does nothing, since pointers are not created in this class. // } //_______________________________________________________________________________________________________________________ void AliRsnMother::Reset() { // // Resets the mother, zeroing all data members. // Int_t i; for (i = 0; i < 2; i++) fDaughter[i] = 0x0; fRefEvent = 0x0; fSum.SetXYZM(0.0, 0.0, 0.0, 0.0); fRef.SetXYZM(0.0, 0.0, 0.0, 0.0); } //__________________________________________________________________________________________________ Int_t AliRsnMother::CommonMother() const { // // If MC info is available, checks if the two tracks in the pair have the same mother. // If the mother label is the same, the function returns the PDG code of mother, // otherwise it returns 0. // The two arguments passed by reference contain the GEANT labels of the mother // of the two particles to which the two daughters point. This is for being able // to check if they are really coming from a resonance (indexes >= 0) or not. // Int_t m1 = fDaughter[0]->GetMother(); Int_t m2 = fDaughter[1]->GetMother(); Int_t out = 0; // a true mother makes sense only if both mothers // are not-negative and equal if (m1 >= 0 && m2 >= 0 && m1 == m2) { out = TMath::Abs(fDaughter[0]->GetMotherPDG()); AliDebugClass(1, Form("Mothers are: %d %d --> EQUAL (PDG = %d)", m1, m2, out)); } return out; } //__________________________________________________________________________________________________ Double_t AliRsnMother::AngleToLeading(Bool_t &success) { // // Compute the angle betwee this and the leading particls // of the reference event (if this was set properly). // In case one of the two daughters is the leading, return // a meaningless value, in order to skip this pair. // if second argument is kTRUE, use MC values. // if (!fRefEvent) { success = kFALSE; return -99.0; } Int_t id1 = fDaughter[0]->GetID(); Int_t id2 = fDaughter[1]->GetID(); Int_t idL = fRefEvent->GetLeadingIndex(); if (id1 == idL || id2 == idL) { success = kFALSE; return -99.0; } AliRsnDaughter leading = fRefEvent->GetDaughter(idL, kFALSE); AliVParticle *ref = leading.GetRef(); Double_t angle = ref->Phi() - fSum.Phi(); //return angle w.r.t. leading particle in the range -pi/2, 3/2pi while (angle >= 1.5 * TMath::Pi()) angle -= 2 * TMath::Pi(); while (angle < -0.5 * TMath::Pi()) angle += 2 * TMath::Pi(); success = kTRUE; return angle; } //__________________________________________________________________________________________________ void AliRsnMother::ComputeSum(Double_t mass0, Double_t mass1, Double_t motherMass) { // // Sets the masses for the 4-momenta of the daughters and then // sums them, taking into account that the space part is set to // each of them when the reference object is set (see AliRsnDaughter::SetRef) // fDaughter[0]->FillP(mass0); fDaughter[1]->FillP(mass1); // sum fSum = fDaughter[0]->Prec() + fDaughter[1]->Prec(); fSumMC = fDaughter[0]->Psim() + fDaughter[1]->Psim(); // reference fRef.SetXYZM(fSum.X(), fSum.Y(), fSum.Z(), motherMass); fRefMC.SetXYZM(fSumMC.X(), fSumMC.Y(), fSumMC.Z(), motherMass); } //__________________________________________________________________________________________________ Double_t AliRsnMother::CosThetaStar(Bool_t first, Bool_t useMC) { // // Computes the cosine of theta*, which is the angle of one of the daughters // with respect to the total momentum of the resonance, in its rest frame. // The arguments are needed to choose which of the daughters one want to use // and if reconstructed or MC momentum must be used. // [Contribution from Z. Feckova] // TLorentzVector &mother = Sum(useMC); TLorentzVector &daughter0 = (first ? fDaughter[0]->P(useMC) : fDaughter[1]->P(useMC)); // TLorentzVector &daughter1 = (first ? fDaughter[1]->P(useMC) : fDaughter[0]->P(useMC)); TVector3 momentumM(mother.Vect()); TVector3 normal(mother.Y() / momentumM.Mag(), -mother.X() / momentumM.Mag(), 0.0); // // Computes first the invariant mass of the mother // Double_t mass0 = fDaughter[0]->P(useMC).M(); // Double_t mass1 = fDaughter[1]->P(useMC).M(); // Double_t p0 = daughter0.Vect().Mag(); // Double_t p1 = daughter1.Vect().Mag(); // Double_t E0 = TMath::Sqrt(mass0 * mass0 + p0 * p0); // Double_t E1 = TMath::Sqrt(mass1 * mass1 + p1 * p1); // Double_t MotherMass = TMath::Sqrt((E0 + E1) * (E0 + E1) - (p0 * p0 + 2.0 * daughter0.Vect().Dot(daughter1.Vect()) + p1 * p1)); // MotherMass = fSum.M(); // Computes components of beta Double_t betaX = -mother.X() / mother.E(); Double_t betaY = -mother.Y() / mother.E(); Double_t betaZ = -mother.Z() / mother.E(); // Computes Lorentz transformation of the momentum of the first daughter // into the rest frame of the mother and theta* daughter0.Boost(betaX, betaY, betaZ); TVector3 momentumD = daughter0.Vect(); Double_t cosThetaStar = normal.Dot(momentumD) / momentumD.Mag(); return cosThetaStar; } //__________________________________________________________________________________________________ Double_t AliRsnMother::DCAproduct() { // // returns product of DCA of the two daughters // AliRsnEvent *event1 = fDaughter[0]->GetOwnerEvent(); AliRsnEvent *event2 = fDaughter[1]->GetOwnerEvent(); if(event1 != event2){ AliError("Attempting to build pair with tracks coming from different events"); return 0.0; } if (event1->IsAOD()) { AliAODEvent *aodEvent = (AliAODEvent*)event1->GetRefAOD(); if (!aodEvent) return 0.0; AliAODTrack *track1 = (AliAODTrack*)fDaughter[0]->Ref2AODtrack(); AliAODTrack *track2 = (AliAODTrack*)fDaughter[1]->Ref2AODtrack(); AliVVertex *vertex = aodEvent->GetPrimaryVertex(); if (!vertex || !track1 || !track2) return 0.0; Double_t b1[2], cov1[3], b2[2], cov2[3]; if ( track1->PropagateToDCA(vertex, aodEvent->GetMagneticField(), kVeryBig, b1, cov1) && track2->PropagateToDCA(vertex, aodEvent->GetMagneticField(), kVeryBig, b2, cov2) ) fDCAproduct = b1[0]*b2[0]; else fDCAproduct = 0.0; } else { AliESDEvent *esdEvent = (AliESDEvent*)event1->GetRefESD(); if (!esdEvent) return 0.0; AliESDtrack *track1 = (AliESDtrack*)fDaughter[0]->Ref2ESDtrack(); AliESDtrack *track2 = (AliESDtrack*)fDaughter[1]->Ref2ESDtrack(); const AliVVertex *vertex = esdEvent->GetPrimaryVertex(); if (!vertex || !track1 || !track2) return 0.0; Double_t b1[2], cov1[3], b2[2], cov2[3]; if ( track1->PropagateToDCA(vertex, esdEvent->GetMagneticField(), kVeryBig, b1, cov1) && track2->PropagateToDCA(vertex, esdEvent->GetMagneticField(), kVeryBig, b2, cov2) ) fDCAproduct = b1[0]*b2[0]; else fDCAproduct = 0.0; } return fDCAproduct; } //__________________________________________________________________________________________________ void AliRsnMother::PrintInfo(const Option_t * /*option*/) const { // // Print some info of the pair. // The options are passed to the AliRsnDaughter::Print() method // AliInfo("======== BEGIN PAIR INFO ==========="); AliInfo("Track #1"); fDaughter[0]->Print(); AliInfo("Track #2"); fDaughter[1]->Print(); AliInfo("========= END PAIR INFO ==========="); } //__________________________________________________________________________________________________ Bool_t AliRsnMother::CheckPair(Bool_t checkMC) const { // // Checks that the pair is well initialized: // - both daughters are good pointers // - if MC is required, both daughters have a MC reference // if (!fDaughter[0] || !fDaughter[1]) { AliError("One of the two tracks is NULL in this pair!"); return kFALSE; } if (checkMC) { if (fDaughter[0]->GetRefMC() == 0x0 || fDaughter[1]->GetRefMC() == 0x0) { AliError("Required MC info but not all MC refs are available"); return kFALSE; } } return kTRUE; }