--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2006, 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. *
+ **************************************************************************/
+
+/////////////////////////////////////////////////////////////
+//
+// Base class for AOD reconstructed decay
+//
+// Author: A.Dainese, andrea.dainese@lnl.infn.it
+/////////////////////////////////////////////////////////////
+
+#include <TDatabasePDG.h>
+#include <TVector3.h>
+#include "AliVirtualParticle.h"
+#include "AliAODRecoDecay.h"
+
+ClassImp(AliAODRecoDecay)
+
+//--------------------------------------------------------------------------
+AliAODRecoDecay::AliAODRecoDecay() :
+ AliVirtualParticle(),
+ fSecondaryVtx(0x0),
+ fCharge(0),
+ fNProngs(0), fNDCA(0), fNPID(0),
+ fPx(0x0), fPy(0x0), fPz(0x0),
+ fd0(0x0),
+ fDCA(0x0),
+ fPID(0x0),
+ fEventNumber(-1),fRunNumber(-1)
+{
+ //
+ // Default Constructor
+ //
+ fSecondaryVtx = new AliAODVertex();
+}
+//--------------------------------------------------------------------------
+AliAODRecoDecay::AliAODRecoDecay(AliAODVertex *vtx2,Int_t nprongs,
+ Short_t charge,
+ Double_t *px,Double_t *py,Double_t *pz,
+ Double_t *d0) :
+ AliVirtualParticle(),
+ fSecondaryVtx(vtx2),
+ fCharge(charge),
+ fNProngs(nprongs), fNDCA(0), fNPID(0),
+ fPx(0x0), fPy(0x0), fPz(0x0),
+ fd0(0x0),
+ fDCA(0x0),
+ fPID(0x0),
+ fEventNumber(-1),fRunNumber(-1)
+{
+ //
+ // Constructor with AliAODVertex for decay vertex
+ //
+
+ fPx = new Double_t[GetNProngs()];
+ fPy = new Double_t[GetNProngs()];
+ fPz = new Double_t[GetNProngs()];
+ fd0 = new Double_t[GetNProngs()];
+ for(Int_t i=0; i<GetNProngs(); i++) {
+ fPx[i] = px[i];
+ fPy[i] = py[i];
+ fPz[i] = pz[i];
+ fd0[i] = d0[i];
+ }
+}
+//--------------------------------------------------------------------------
+AliAODRecoDecay::AliAODRecoDecay(AliAODVertex *vtx2,Int_t nprongs,
+ Short_t charge,
+ Double_t *d0) :
+ AliVirtualParticle(),
+ fSecondaryVtx(vtx2),
+ fCharge(charge),
+ fNProngs(nprongs), fNDCA(0), fNPID(0),
+ fPx(0x0), fPy(0x0), fPz(0x0),
+ fd0(0x0),
+ fDCA(0x0),
+ fPID(0x0),
+ fEventNumber(-1),fRunNumber(-1)
+{
+ //
+ // Constructor with AliAODVertex for decay vertex and without prongs momenta
+ //
+
+ fd0 = new Double_t[GetNProngs()];
+ for(Int_t i=0; i<GetNProngs(); i++) fd0[i] = d0[i];
+}
+//--------------------------------------------------------------------------
+AliAODRecoDecay::AliAODRecoDecay(const AliAODRecoDecay &source) :
+ AliVirtualParticle(source),
+ fSecondaryVtx(source.fSecondaryVtx),
+ fCharge(source.fCharge),
+ fNProngs(source.fNProngs), fNDCA(source.fNDCA), fNPID(source.fNPID),
+ fPx(0x0), fPy(0x0), fPz(0x0),
+ fd0(0x0),
+ fDCA(0x0),
+ fPID(0x0),
+ fEventNumber(source.fEventNumber),fRunNumber(source.fRunNumber)
+{
+ //
+ // Copy constructor
+ //
+ if(source.GetNProngs()>0) {
+ fd0 = new Double_t[GetNProngs()];
+ memcpy(fd0,source.fd0,GetNProngs()*sizeof(Double_t));
+ if(source.fPx) {
+ fPx = new Double_t[GetNProngs()];
+ fPy = new Double_t[GetNProngs()];
+ fPz = new Double_t[GetNProngs()];
+ memcpy(fPx,source.fPx,GetNProngs()*sizeof(Double_t));
+ memcpy(fPy,source.fPy,GetNProngs()*sizeof(Double_t));
+ memcpy(fPz,source.fPz,GetNProngs()*sizeof(Double_t));
+ }
+ if(source.fPID) {
+ fPID = new Double_t[5*GetNProngs()];
+ memcpy(fPID,source.fPID,GetNProngs()*sizeof(Double_t));
+ }
+ if(source.fDCA) {
+ fDCA = new Float_t[GetNProngs()*(GetNProngs()-1)/2];
+ memcpy(fDCA,source.fDCA,(GetNProngs()*(GetNProngs()-1)/2)*sizeof(Float_t));
+ }
+ }
+}
+//--------------------------------------------------------------------------
+AliAODRecoDecay &AliAODRecoDecay::operator=(const AliAODRecoDecay &source)
+{
+ //
+ // assignment operator
+ //
+ if(&source == this) return *this;
+ fSecondaryVtx = source.fSecondaryVtx;
+ fCharge = source.fCharge;
+ fNProngs = source.fNProngs;
+ fNDCA = source.fNDCA;
+ fNPID = source.fNPID;
+ fEventNumber = source.fEventNumber;
+ fRunNumber = source.fRunNumber;
+ if(source.GetNProngs()>0) {
+ fd0 = new Double_t[GetNProngs()];
+ memcpy(fd0,source.fd0,GetNProngs()*sizeof(Double_t));
+ if(source.fPx) {
+ fPx = new Double_t[GetNProngs()];
+ fPy = new Double_t[GetNProngs()];
+ fPz = new Double_t[GetNProngs()];
+ memcpy(fPx,source.fPx,GetNProngs()*sizeof(Double_t));
+ memcpy(fPy,source.fPy,GetNProngs()*sizeof(Double_t));
+ memcpy(fPz,source.fPz,GetNProngs()*sizeof(Double_t));
+ }
+ if(source.fPID) {
+ fPID = new Double_t[5*GetNProngs()];
+ memcpy(fPID,source.fPID,GetNProngs()*sizeof(Double_t));
+ }
+ if(source.fDCA) {
+ fDCA = new Float_t[GetNProngs()*(GetNProngs()-1)/2];
+ memcpy(fDCA,source.fDCA,(GetNProngs()*(GetNProngs()-1)/2)*sizeof(Float_t));
+ }
+ }
+ return *this;
+}
+//--------------------------------------------------------------------------
+AliAODRecoDecay::~AliAODRecoDecay() {
+ //
+ // Default Destructor
+ //
+ if(fSecondaryVtx) delete fSecondaryVtx;
+ if(fPx) delete [] fPx;
+ if(fPy) delete [] fPy;
+ if(fPz) delete [] fPz;
+ if(fd0) delete [] fd0;
+ if(fPID) delete [] fPID;
+ if(fDCA) delete [] fDCA;
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::Alpha() const
+{
+ //
+ // Armenteros-Podolanski alpha for 2-prong decays
+ //
+ if(GetNProngs()!=2) {
+ printf("Can be called only for 2-prong decays");
+ return (Double_t)-99999.;
+ }
+ return 1.-2./(1.+QlProng(0)/QlProng(1));
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::DecayLength(Double_t point[3]) const
+{
+ //
+ // Decay length assuming it is produced at "point" [cm]
+ //
+ return TMath::Sqrt((point[0]-GetSecVtxX())
+ *(point[0]-GetSecVtxX())
+ +(point[1]-GetSecVtxY())
+ *(point[1]-GetSecVtxY())
+ +(point[2]-GetSecVtxZ())
+ *(point[2]-GetSecVtxZ()));
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::DecayLengthXY(Double_t point[3]) const
+{
+ //
+ // Decay length in XY assuming it is produced at "point" [cm]
+ //
+ return TMath::Sqrt((point[0]-GetSecVtxX())
+ *(point[0]-GetSecVtxX())
+ +(point[1]-GetSecVtxY())
+ *(point[1]-GetSecVtxY()));
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::CosPointingAngle(Double_t point[3]) const
+{
+ //
+ // Cosine of pointing angle in space assuming it is produced at "point"
+ //
+ TVector3 mom(Px(),Py(),Pz());
+ TVector3 fline(GetSecVtxX()-point[0],
+ GetSecVtxY()-point[1],
+ GetSecVtxZ()-point[2]);
+
+ Double_t pta = mom.Angle(fline);
+
+ return TMath::Cos(pta);
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::CosPointingAngleXY(Double_t point[3]) const
+{
+ //
+ // Cosine of pointing angle in transverse plane assuming it is produced
+ // at "point"
+ //
+ TVector3 momXY(Px(),Py(),0.);
+ TVector3 flineXY(GetSecVtxX()-point[0],
+ GetSecVtxY()-point[1],
+ 0.);
+
+ Double_t ptaXY = momXY.Angle(flineXY);
+
+ return TMath::Cos(ptaXY);
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::CosThetaStar(Int_t ip,UInt_t pdgvtx,UInt_t pdgprong0,UInt_t pdgprong1) const
+{
+ //
+ // Only for 2-prong decays:
+ // Cosine of decay angle (theta*) in the rest frame of the mother particle
+ // for prong ip (0 or 1) with mass hypotheses pdgvtx for mother particle,
+ // pdgprong0 for prong 0 and pdgprong1 for prong1
+ //
+ if(GetNProngs()!=2) {
+ printf("Can be called only for 2-prong decays");
+ return (Double_t)-99999.;
+ }
+ Double_t massvtx = TDatabasePDG::Instance()->GetParticle(pdgvtx)->Mass();
+ Double_t massp[2];
+ massp[0] = TDatabasePDG::Instance()->GetParticle(pdgprong0)->Mass();
+ massp[1] = TDatabasePDG::Instance()->GetParticle(pdgprong1)->Mass();
+
+ Double_t pStar = TMath::Sqrt(TMath::Power(massvtx*massvtx-massp[0]*massp[0]-massp[1]*massp[1],2.)-4.*massp[0]*massp[0]*massp[1]*massp[1])/(2.*massvtx);
+
+ Double_t beta = P()/E(pdgvtx);
+ Double_t gamma = E(pdgvtx)/massvtx;
+
+ Double_t cts = (QlProng(ip)/gamma-beta*TMath::Sqrt(pStar*pStar+massp[ip]*massp[ip]))/pStar;
+
+ return cts;
+}
+//---------------------------------------------------------------------------
+Double_t AliAODRecoDecay::Ct(UInt_t pdg,Double_t point[3]) const
+{
+ //
+ // Decay time * c assuming it is produced at "point" [cm]
+ //
+ Double_t mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
+ return DecayLength(point)*mass/P();
+}
+//---------------------------------------------------------------------------
+Double_t AliAODRecoDecay::E(UInt_t pdg) const
+{
+ //
+ // Energy
+ //
+ Double_t mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
+ return TMath::Sqrt(mass*mass+P()*P());
+}
+//---------------------------------------------------------------------------
+Double_t AliAODRecoDecay::EProng(Int_t ip,UInt_t pdg) const
+{
+ //
+ // Energy of ip-th prong
+ //
+ Double_t mass = TDatabasePDG::Instance()->GetParticle(pdg)->Mass();
+ return TMath::Sqrt(mass*mass+PProng(ip)*PProng(ip));
+}
+//---------------------------------------------------------------------------
+/*Int_t AliAODRecoDecay::GetIndexProng(Int_t ip) const
+{
+ //
+ // Index of prong ip
+ //
+ if(!GetNProngs()) return 999999;
+ UShort_t *indices = fSecondaryVtx->GetIndices();
+ return indices[ip];
+}*/
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::ImpParXY(Double_t point[3]) const
+{
+ //
+ // Impact parameter in the bending plane of the particle
+ // w.r.t. to "point"
+ //
+ Double_t k = -(GetSecVtxX()-point[0])*Px()-(GetSecVtxY()-point[1])*Py();
+ k /= Pt()*Pt();
+ Double_t dx = GetSecVtxX()-point[0]+k*Px();
+ Double_t dy = GetSecVtxY()-point[1]+k*Py();
+ Double_t absImpPar = TMath::Sqrt(dx*dx+dy*dy);
+ TVector3 mom(Px(),Py(),Pz());
+ TVector3 fline(GetSecVtxX()-point[0],
+ GetSecVtxY()-point[1],
+ GetSecVtxZ()-point[2]);
+ TVector3 cross = mom.Cross(fline);
+ return (cross.Z()>0. ? absImpPar : -absImpPar);
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::InvMass(Int_t npdg,UInt_t *pdg) const
+{
+ //
+ // Invariant mass for prongs mass hypotheses in pdg array
+ //
+ if(GetNProngs()!=npdg) {
+ printf("npdg != GetNProngs()");
+ return (Double_t)-99999.;
+ }
+ Double_t energysum = 0.;
+
+ for(Int_t i=0; i<GetNProngs(); i++) {
+ energysum += EProng(i,pdg[i]);
+ }
+
+ Double_t mass = TMath::Sqrt(energysum*energysum-P()*P());
+
+ return mass;
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::InvMass2Prongs(Int_t ip1,Int_t ip2,
+ UInt_t pdg1,UInt_t pdg2) const
+{
+ //
+ // 2-prong(ip1,ip2) invariant mass for prongs mass hypotheses in pdg1,2
+ //
+ Double_t energysum = EProng(ip1,pdg1) + EProng(ip2,pdg2);
+ Double_t psum2 = (PxProng(ip1)+PxProng(ip2))*(PxProng(ip1)+PxProng(ip2))
+ +(PyProng(ip1)+PyProng(ip2))*(PyProng(ip1)+PyProng(ip2))
+ +(PzProng(ip1)+PzProng(ip2))*(PzProng(ip1)+PzProng(ip2));
+ Double_t mass = TMath::Sqrt(energysum*energysum-psum2);
+
+ return mass;
+}
+//---------------------------------------------------------------------------
+void AliAODRecoDecay::Print(Option_t* /*option*/) const
+{
+ //
+ // Print some information
+ //
+ printf("AliAODRecoDecay with %d prongs\n",GetNProngs());
+ printf("Secondary Vertex: (%f, %f, %f)\n",GetSecVtxX(),GetSecVtxY(),GetSecVtxZ());
+
+ return;
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::ProngsRelAngle(Int_t ip1,Int_t ip2) const
+{
+ //
+ // Relative angle between two prongs
+ //
+ TVector3 momA(PxProng(ip1),PyProng(ip1),PzProng(ip1));
+ TVector3 momB(PxProng(ip2),PyProng(ip2),PzProng(ip2));
+
+ Double_t angle = momA.Angle(momB);
+
+ return angle;
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::QlProng(Int_t ip) const
+{
+ //
+ // Longitudinal momentum of prong w.r.t. to total momentum
+ //
+ TVector3 mom(PxProng(ip),PyProng(ip),PzProng(ip));
+ TVector3 momTot(Px(),Py(),Pz());
+
+ return mom.Dot(momTot)/momTot.Mag();
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::QtProng(Int_t ip) const
+{
+ //
+ // Transverse momentum of prong w.r.t. to total momentum
+ //
+ TVector3 mom(PxProng(ip),PyProng(ip),PzProng(ip));
+ TVector3 momTot(Px(),Py(),Pz());
+
+ return mom.Perp(momTot);
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::QlProngFlightLine(Int_t ip,Double_t point[3]) const
+{
+ //
+ // Longitudinal momentum of prong w.r.t. to flight line between "point"
+ // and fSecondaryVtx
+ //
+ TVector3 mom(PxProng(ip),PyProng(ip),PzProng(ip));
+ TVector3 fline(GetSecVtxX()-point[0],
+ GetSecVtxY()-point[1],
+ GetSecVtxZ()-point[2]);
+
+ return mom.Dot(fline)/fline.Mag();
+}
+//----------------------------------------------------------------------------
+Double_t AliAODRecoDecay::QtProngFlightLine(Int_t ip,Double_t point[3]) const
+{
+ //
+ // Transverse momentum of prong w.r.t. to flight line between "point" and
+ // fSecondaryVtx
+ //
+ TVector3 mom(PxProng(ip),PyProng(ip),PzProng(ip));
+ TVector3 fline(GetSecVtxX()-point[0],
+ GetSecVtxY()-point[1],
+ GetSecVtxZ()-point[2]);
+
+ return mom.Perp(fline);
+}
+//--------------------------------------------------------------------------
--- /dev/null
+#ifndef ALIAODRECODECAY_H
+#define ALIAODRECODECAY_H
+/* Copyright(c) 1998-2006, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//***********************************************************
+// Class AliAODRecoDecay
+// base class for AOD reconstructed decays
+// Author: A.Dainese, andrea.dainese@lnl.infn.it
+//***********************************************************
+
+#include <TMath.h>
+#include "AliAODVertex.h"
+#include "AliVirtualParticle.h"
+
+class AliAODRecoDecay : public AliVirtualParticle {
+
+ public:
+
+ AliAODRecoDecay();
+ AliAODRecoDecay(AliAODVertex *vtx2,Int_t nprongs,Short_t charge,
+ Double_t *px,Double_t *py,Double_t *pz,
+ Double_t *d0);
+ AliAODRecoDecay(AliAODVertex *vtx2,Int_t nprongs,Short_t charge,
+ Double_t *d0);
+ virtual ~AliAODRecoDecay();
+
+ AliAODRecoDecay(const AliAODRecoDecay& source);
+ AliAODRecoDecay& operator=(const AliAODRecoDecay& source);
+
+
+ // decay vertex
+ Double_t GetSecVtxX() const {return fSecondaryVtx->GetX();}
+ Double_t GetSecVtxY() const {return fSecondaryVtx->GetY();}
+ Double_t GetSecVtxZ() const {return fSecondaryVtx->GetZ();}
+ Double_t RadiusSecVtx() const;
+ void SetSecondaryVtx(AliAODVertex *vtx2) {fSecondaryVtx=vtx2;}
+ AliAODVertex* GetSecondaryVtx() const {return fSecondaryVtx;}
+ void GetSecondaryVtx(Double_t vtx[3]) const;
+ Double_t GetReducedChi2() const {return fSecondaryVtx->GetChi2perNDF();}
+ Short_t Charge() const {return fCharge;}
+ Short_t GetCharge() const {return fCharge;}
+ void SetCharge(Short_t charge=0) {fCharge=charge;}
+
+ // PID
+ void SetPID(Int_t nprongs,Double_t *pid);
+ Double_t *GetPID() const { return fPID; }
+ void GetPIDProng(Int_t ip,Double_t *pid) const;
+ virtual const Double_t *PID() const { return fPID; }
+
+ // prong-to-prong DCAs
+ void SetDCAs(Int_t nDCA,Float_t *dca);
+ void SetDCA(Float_t dca); // 2 prong
+ Float_t GetDCA(Int_t i=0) const {return fDCA[i];}
+
+ //event and run number
+ void SetEventRunNumbers(Int_t en,Int_t rn)
+ { fEventNumber=en; fRunNumber=rn; return; }
+ Int_t GetEventNumber() const { return fEventNumber; }
+ Int_t GetRunNumber() const { return fRunNumber; }
+
+ // kinematics & topology
+ Double_t Px() const;
+ Double_t Py() const;
+ Double_t Pz() const;
+ Double_t P() const {return TMath::Sqrt(Px()*Px()+Py()*Py()+Pz()*Pz());}
+ Double_t Pt() const {return TMath::Sqrt(Px()*Px()+Py()*Py());}
+ Double_t OneOverPt() const {return (Pt() ? 1./Pt() : 0.);}
+ Double_t Phi() const {return TMath::ATan2(Py(),Px());}
+ Double_t Theta() const {return 0.5*TMath::Pi()-TMath::ATan(Pz()/(Pt()+1.e-13));}
+ Double_t Eta() const {return 0.5*TMath::Log((P()+Pz())/(P()-Pz()+1.e-13));}
+ Double_t E(UInt_t pdg) const;
+ Double_t Y(UInt_t pdg) const {return 0.5*TMath::Log((E(pdg)+Pz())/(E(pdg)-Pz()+1.e-13));}
+ Double_t DecayLength(Double_t point[3]) const;
+ Double_t DecayLength(AliAODVertex *vtx1) const
+ {return fSecondaryVtx->DistanceToVertex(vtx1);}
+ Double_t DecayLengthError(AliAODVertex *vtx1) const
+ {return fSecondaryVtx->ErrorDistanceToVertex(vtx1);}
+ Double_t NormalizedDecayLength(AliAODVertex *vtx1) const
+ {return DecayLength(vtx1)/DecayLengthError(vtx1);}
+ Double_t DecayLengthXY(Double_t point[3]) const;
+ Double_t DecayLengthXY(AliAODVertex *vtx1) const
+ {return fSecondaryVtx->DistanceXYToVertex(vtx1);}
+ Double_t DecayLengthXYError(AliAODVertex *vtx1) const
+ {return fSecondaryVtx->ErrorDistanceXYToVertex(vtx1);}
+ Double_t NormalizedDecayLengthXY(AliAODVertex *vtx1) const
+ {return DecayLengthXY(vtx1)/DecayLengthXYError(vtx1);}
+ Double_t Ct(UInt_t pdg,Double_t point[3]) const;
+ Double_t Ct(UInt_t pdg,AliAODVertex *vtx1) const;
+ Double_t CosPointingAngle(Double_t point[3]) const;
+ Double_t CosPointingAngle(AliAODVertex *vtx1) const;
+ Double_t CosPointingAngleXY(Double_t point[3]) const;
+ Double_t CosPointingAngleXY(AliAODVertex *vtx1) const;
+ Double_t CosThetaStar(Int_t ip,UInt_t pdgvtx,UInt_t pdgprong0,UInt_t pdgprong1) const;
+ Double_t InvMass(Int_t npdg,UInt_t *pdg) const;
+ Double_t ImpParXY(Double_t point[3]) const;
+ Double_t ImpParXY(AliAODVertex *vtx1) const;
+
+ // prongs
+ //Int_t GetNProngs() const {return fSecondaryVtx->GetNDaughters();}
+ Int_t GetNProngs() const {return fNProngs;}
+
+ Short_t ChargeProng(Int_t ip) const;
+ Double_t Getd0Prong(Int_t ip) const {return fd0[ip];}
+ Double_t Prodd0d0(Int_t ip1=0,Int_t ip2=0) const {return fd0[ip1]*fd0[ip2];}
+ Double_t PxProng(Int_t ip) const {return fPx[ip];}
+ Double_t PyProng(Int_t ip) const {return fPy[ip];}
+ Double_t PzProng(Int_t ip) const {return fPz[ip];}
+ Double_t PtProng(Int_t ip) const;
+ Double_t PProng(Int_t ip) const;
+ Double_t PhiProng(Int_t ip) const
+ {return TMath::ATan2(PyProng(ip),PxProng(ip));}
+ Double_t ThetaProng(Int_t ip) const
+ {return 0.5*TMath::Pi()-TMath::ATan(PzProng(ip)/(PtProng(ip)+1.e-13));}
+ Double_t EtaProng(Int_t ip) const
+ {return -TMath::Log(TMath::Tan(0.5*ThetaProng(ip)));}
+ Double_t EProng(Int_t ip,UInt_t pdg) const;
+ Double_t YProng(Int_t ip,UInt_t pdg) const
+ {return 0.5*TMath::Log((EProng(ip,pdg)+PzProng(ip))/(EProng(ip,pdg)-PzProng(ip)+1.e-13));}
+ Double_t Alpha() const; // for Armenteros-Podolanski plot (V0's)
+ Double_t QlProng(Int_t ip) const;
+ Double_t QtProng(Int_t ip=0) const; // for Armenteros-Podolanski plot (V0's)
+ Double_t QlProngFlightLine(Int_t ip,Double_t point[3]) const;
+ Double_t QlProngFlightLine(Int_t ip,AliAODVertex *vtx1) const;
+ Double_t QtProngFlightLine(Int_t ip,Double_t point[3]) const;
+ Double_t QtProngFlightLine(Int_t ip,AliAODVertex *vtx1) const;
+ Double_t InvMass2Prongs(Int_t ip1,Int_t ip2,UInt_t pdg1,UInt_t pdg2) const;
+ Double_t ProngsRelAngle(Int_t ip1=0,Int_t ip2=1) const;
+
+ // relate to other objects
+ //Double_t DistanceToVertex(AliAODVertex *vtx) // distance to a AliAODVertex
+ //Double_t DistanceToTrack(AliAODTrack *trk) // distance to a AliAODTrack
+
+
+ // print
+ void Print(Option_t* option = "") const;
+ //void PrintIndices() const {fSecondaryVtx->PrintIndices();}
+
+ // dummy functions for inheritance from AliVirtualParticle
+ Double_t E() const
+ {printf("Dummy function; use AliAODRecoDecay::E(UInt_t pdg) instead"); return (Double_t)-999.;}
+ Double_t Y() const
+ {printf("Dummy function; use AliAODRecoDecay::Y(UInt_t pdg) instead"); return (Double_t)-999.;}
+ Double_t M() const
+ {printf("Dummy function"); return (Double_t)-999.;}
+
+ protected:
+
+ AliAODVertex *fSecondaryVtx; // decay vertex
+ Short_t fCharge; // charge, use this convention for prongs charges:
+ // if(charge== 0) even-index prongs are +
+ // odd-index prongs are -
+ // if(charge==+1) even-index prongs are +
+ // odd-index prongs are -
+ // if(charge==-1) even-index prongs are -
+ // odd-index prongs are +
+
+ // TEMPORARY, to be removed when we do analysis on AliAODEvent
+ Int_t fNProngs; // number of prongs
+ Int_t fNDCA; // number of dca's
+ Int_t fNPID; // number of PID probabilities
+ Double_t *fPx; //[fNProngs] px of tracks at the vertex [GeV/c]
+ Double_t *fPy; //[fNProngs] py of tracks at the vertex [GeV/c]
+ Double_t *fPz; //[fNProngs] pz of tracks at the vertex [GeV/c]
+ Double_t *fd0; //[fNProngs] rphi impact params w.r.t. Primary Vtx [cm]
+ Float_t *fDCA; //[fNDCA] prong-to-prong DCA [cm]
+ // convention:fDCA[0]=p0p1,fDCA[1]=p0p2,fDCA[2]=p1p2,...
+ Double_t *fPID; //[fNPID] combined pid
+ // (combined detector response probabilities)
+
+ // TEMPORARY, to be removed when we do analysis on AliAODEvent
+ Int_t fEventNumber;
+ Int_t fRunNumber;
+ // TO BE PUT IN SPECIAL MC CLASS
+ //Bool_t fSignal; // TRUE if signal, FALSE if background (for simulation)
+ //Int_t fEvent; // number of the event this candidate comes from
+ //Int_t fTrkNum[2]; // numbers of the two decay tracks
+ //Int_t fPdg[2]; // PDG codes of the two tracks (for sim.)
+ //Int_t fMum[2]; // PDG codes of the mothers (for sim.)
+
+ //
+
+ ClassDef(AliAODRecoDecay,1) // base class for AOD reconstructed decays
+};
+
+
+inline Short_t AliAODRecoDecay::ChargeProng(Int_t ip) const
+{
+ if(fCharge==0 || fCharge==+1) {
+ if(ip%2==0) {
+ return (Short_t)1;
+ } else {
+ return (Short_t)-1;
+ }
+ } else { // fCharge==-1
+ if(ip%2==0) {
+ return (Short_t)-1;
+ } else {
+ return (Short_t)1;
+ }
+ }
+}
+
+inline Double_t AliAODRecoDecay::RadiusSecVtx() const
+{
+ return TMath::Sqrt(GetSecVtxX()*GetSecVtxX()+GetSecVtxY()*GetSecVtxY());
+}
+
+inline void AliAODRecoDecay::GetSecondaryVtx(Double_t vtx[3]) const
+{
+ fSecondaryVtx->GetPosition(vtx);
+ return;
+}
+
+inline Double_t AliAODRecoDecay::Px() const
+{
+ Double_t px=0.;
+ for(Int_t i=0;i<GetNProngs();i++) px+=PxProng(i);
+ return px;
+}
+
+inline Double_t AliAODRecoDecay::Py() const
+{
+ Double_t py=0.;
+ for(Int_t i=0;i<GetNProngs();i++) py+=PyProng(i);
+ return py;
+}
+
+inline Double_t AliAODRecoDecay::Pz() const
+{
+ Double_t pz=0.;
+ for(Int_t i=0;i<GetNProngs();i++) pz+=PzProng(i);
+ return pz;
+}
+
+inline Double_t AliAODRecoDecay::Ct(UInt_t pdg,AliAODVertex *vtx1) const
+{
+ Double_t v[3];
+ vtx1->GetPosition(v);
+ return Ct(pdg,v);
+}
+
+inline Double_t AliAODRecoDecay::CosPointingAngle(AliAODVertex *vtx1) const
+{
+ Double_t v[3];
+ vtx1->GetPosition(v);
+ return CosPointingAngle(v);
+}
+
+inline Double_t AliAODRecoDecay::CosPointingAngleXY(AliAODVertex *vtx1) const
+{
+ Double_t v[3];
+ vtx1->GetPosition(v);
+ return CosPointingAngleXY(v);
+}
+
+inline Double_t AliAODRecoDecay::ImpParXY(AliAODVertex *vtx1) const
+{
+ Double_t v[3];
+ vtx1->GetPosition(v);
+ return ImpParXY(v);
+}
+
+inline Double_t AliAODRecoDecay::PtProng(Int_t ip) const
+{
+ return TMath::Sqrt(PxProng(ip)*PxProng(ip)+PyProng(ip)*PyProng(ip));
+}
+
+inline Double_t AliAODRecoDecay::PProng(Int_t ip) const
+{
+ return TMath::Sqrt(PtProng(ip)*PtProng(ip)+PzProng(ip)*PzProng(ip));
+}
+
+inline Double_t AliAODRecoDecay::QlProngFlightLine(Int_t ip,AliAODVertex *vtx1) const
+{
+ Double_t v[3];
+ vtx1->GetPosition(v);
+ return QlProngFlightLine(ip,v);
+}
+
+inline Double_t AliAODRecoDecay::QtProngFlightLine(Int_t ip,AliAODVertex *vtx1) const
+{
+ Double_t v[3];
+ vtx1->GetPosition(v);
+ return QtProngFlightLine(ip,v);
+}
+
+inline void AliAODRecoDecay::SetDCAs(Int_t nDCA,Float_t *dca)
+{
+ if(nDCA!=(GetNProngs()*(GetNProngs()-1)/2)) {
+ printf("Wrong number of DCAs, must be nProngs*(nProngs-1)/2");
+ return;
+ }
+ if(fDCA) delete [] fDCA;
+ fDCA = new Float_t[nDCA];
+ for(Int_t i=0;i<nDCA;i++)
+ fDCA[i] = dca[i];
+ return;
+}
+
+inline void AliAODRecoDecay::SetDCA(Float_t dca)
+{
+ Float_t ddca[1]; ddca[0]=dca;
+ SetDCAs(1,ddca);
+ return;
+}
+
+inline void AliAODRecoDecay::SetPID(Int_t nprongs,Double_t *pid)
+{
+ if(nprongs!=GetNProngs()) {
+ printf("Wrong number of prongs");
+ return;
+ }
+ if(fPID) delete [] fPID;
+ fPID = new Double_t[nprongs*5];
+ for(Int_t i=0;i<nprongs;i++)
+ for(Int_t j=0;j<5;j++)
+ fPID[i*5+j] = pid[i*5+j];
+ return;
+}
+
+inline void AliAODRecoDecay::GetPIDProng(Int_t ip,Double_t *pid) const
+{
+ for(Int_t j=0;j<5;j++)
+ pid[j] = fPID[ip*5+j];
+ return;
+}
+
+
+
+#endif
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff