-/**************************************************************************\r
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
- * *\r
- * Author: The ALICE Off-line Project. *\r
- * Contributors are mentioned in the code where appropriate. *\r
- * *\r
- * Permission to use, copy, modify and distribute this software and its *\r
- * documentation strictly for non-commercial purposes is hereby granted *\r
- * without fee, provided that the above copyright notice appears in all *\r
- * copies and that both the copyright notice and this permission notice *\r
- * appear in the supporting documentation. The authors make no claims *\r
- * about the suitability of this software for any purpose. It is *\r
- * provided "as is" without express or implied warranty. *\r
- **************************************************************************/\r
-\r
-/*\r
- * author: M.Kalisky@gsi.de\r
- * 08/Dec/2010\r
- *\r
- * Description: This class allows with purely kinematical cuts\r
- * to select clean samples of electrons, pions and protons from the\r
- * V0 online finder ESD V0 candidates for PID and dectector resonse\r
- * studies.\r
- */\r
-\r
-#include <TVector3.h>\r
-#include <TDatabasePDG.h>\r
-\r
-#include "AliESDv0.h"\r
-#include "AliESDtrack.h"\r
-#include "AliESDEvent.h"\r
-#include "AliVEvent.h"\r
-#include "AliLog.h"\r
-#include "AliKFParticle.h"\r
-#include "AliVTrack.h"\r
-#include "AliKFVertex.h"\r
-\r
-#include "AliESDv0KineCuts.h"\r
-\r
-ClassImp(AliESDv0KineCuts)\r
-\r
-//____________________________________________________________________\r
-AliESDv0KineCuts::AliESDv0KineCuts() :\r
- fEvent(0x0)\r
- , fPrimaryVertex(0x0)\r
- , fType(0)\r
- , fMode(0)\r
- , fTPCNcls(1)\r
- , fTPCrefit(kTRUE)\r
- , fTPCchi2perCls(4.0)\r
- , fTPCclsRatio(0.6)\r
- , fNoKinks(kTRUE)\r
- , fGcutChi2NDF(10)\r
- , fGcutInvMass(0.05)\r
- , fK0cutChi2NDF(10)\r
- , fLcutChi2NDF(10)\r
- , fUseExternalVertex(kFALSE)\r
- , fDeleteVertex(kFALSE)\r
-{\r
- //\r
- // Default constructor\r
- //\r
-\r
- // default single track cuts\r
- fTPCNcls = 1; // minimal number of the TPC clusters\r
- fTPCrefit = kTRUE; // TPC refit\r
- fTPCchi2perCls = 4.0; // chi2 per TPC cluster\r
- fTPCclsRatio = 0.6; // minimal foun/findable TPC cluster ratio\r
- fNoKinks = kTRUE; // kinks - no [kTRUE] or do not care [kFalse]\r
-\r
-\r
- // default gamma cuts values\r
- fGcutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle gamma\r
- fGcutCosPoint[0] = 0; // cos of the pointing angle [min, max]\r
- fGcutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]\r
- fGcutDCA[0] = 0.; // DCA between the daughter tracks [min, max]\r
- fGcutDCA[1] = 0.25; // DCA between the daughter tracks [min, max]\r
- fGcutVertexR[0] = 3.; // radius of the conversion point [min, max]\r
- fGcutVertexR[1] = 90.; // radius of the conversion point [min, max]\r
- fGcutPsiPair[0] = 0.; // value of the psi pair cut [min, max]\r
- fGcutPsiPair[1] = 0.05; // value of the psi pair cut [min, max]\r
- fGcutInvMass = 0.05; // upper value on the gamma invariant mass\r
- // default K0 cuts\r
- fK0cutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle K0\r
- fK0cutCosPoint[0] = 0.; // cos of the pointing angle [min, max]\r
- fK0cutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]\r
- fK0cutDCA[0] = 0.; // DCA between the daughter tracks [min, max]\r
- fK0cutDCA[1] = 0.2; // DCA between the daughter tracks [min, max]\r
- fK0cutVertexR[0] = 2.0; // radius of the decay point [min, max]\r
- fK0cutVertexR[1] = 30.0; // radius of the decay point [min, max]\r
- fK0cutInvMass[0] = 0.486; // invariant mass window\r
- fK0cutInvMass[1] = 0.508; // invariant mass window\r
- // Lambda & anti-Lambda cut values\r
- fLcutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle K0\r
- fLcutCosPoint[0] = 0.; // cos of the pointing angle [min, max]\r
- fLcutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]\r
- fLcutDCA[0] = 0.; // DCA between the daughter tracks [min, max]\r
- fLcutDCA[1] = 0.2; // DCA between the daughter tracks [min, max]\r
- fLcutVertexR[0] = 2.0; // radius of the decay point [min, max]\r
- fLcutVertexR[1] = 40.0; // radius of the decay point [min, max]\r
- fLcutInvMass[0] = 1.11; // invariant mass window\r
- fLcutInvMass[1] = 1.12; // invariant mass window\r
- \r
-}\r
-//____________________________________________________________________\r
-AliESDv0KineCuts::~AliESDv0KineCuts(){\r
- //\r
- // Destructor\r
- //\r
-\r
-\r
-}\r
-//____________________________________________________________________\r
-AliESDv0KineCuts::AliESDv0KineCuts(const AliESDv0KineCuts &ref):\r
- TObject(ref)\r
- , fEvent(0x0)\r
- , fPrimaryVertex(0x0)\r
- , fType(0)\r
- , fMode(0)\r
- , fTPCNcls(1)\r
- , fTPCrefit(kTRUE)\r
- , fTPCchi2perCls(4.0)\r
- , fTPCclsRatio(0.6)\r
- , fNoKinks(kTRUE)\r
- , fGcutChi2NDF(10)\r
- , fGcutInvMass(0.05)\r
- , fK0cutChi2NDF(10)\r
- , fLcutChi2NDF(10)\r
- , fUseExternalVertex(kFALSE)\r
- , fDeleteVertex(kFALSE)\r
-{\r
- //\r
- // Copy operator\r
- //\r
-\r
- ref.Copy(*this);\r
-}\r
-//____________________________________________________________________\r
-AliESDv0KineCuts &AliESDv0KineCuts::operator=(const AliESDv0KineCuts &ref){\r
- //\r
- // assignment operator\r
- //\r
- if(this != &ref)\r
- ref.Copy(*this);\r
- return *this; \r
-}\r
-//____________________________________________________________________\r
-void AliESDv0KineCuts::Copy(TObject &ref) const {\r
- //\r
- // Performs the copying of the object\r
- //\r
-\r
- TObject::Copy(ref);\r
-\r
- AliESDv0KineCuts &target = dynamic_cast<AliESDv0KineCuts &>(ref);\r
-\r
- // default single track cuts\r
- target.fTPCNcls = fTPCNcls;\r
- target.fTPCrefit = fTPCrefit;\r
- target.fTPCchi2perCls = fTPCchi2perCls;\r
- target.fTPCclsRatio = fTPCclsRatio;\r
- target.fNoKinks = fNoKinks;\r
- target.fUseExternalVertex = fUseExternalVertex; //added december 2nd 2011\r
- target.fDeleteVertex = fDeleteVertex; //added december 2nd 2011\r
-\r
- // default gamma cuts values\r
- target.fGcutChi2NDF = fGcutChi2NDF;\r
- memcpy(target.fGcutCosPoint, fGcutCosPoint, sizeof(Float_t) * 2);\r
- memcpy(target.fGcutDCA, fGcutDCA, sizeof(Float_t) * 2); \r
- memcpy(target.fGcutVertexR, fGcutVertexR, sizeof(Float_t) * 2);\r
- memcpy(target.fGcutPsiPair, fGcutPsiPair, sizeof(Float_t) * 2);\r
- target.fGcutInvMass = fGcutInvMass;\r
- // default K0 cuts\r
- target.fK0cutChi2NDF = fK0cutChi2NDF;\r
- memcpy(target.fK0cutCosPoint, fK0cutCosPoint, sizeof(Float_t) * 2);\r
- memcpy(target.fK0cutDCA, fK0cutDCA, sizeof(Float_t) * 2);\r
- memcpy(target.fK0cutVertexR, fK0cutVertexR, sizeof(Float_t) * 2);\r
- memcpy(target.fK0cutInvMass, fK0cutInvMass, sizeof(Float_t) * 2);\r
- // Lambda & anti-Lambda cut values\r
- target.fLcutChi2NDF = fLcutChi2NDF;\r
- memcpy(target.fLcutCosPoint, fLcutCosPoint, sizeof(Float_t) * 2);\r
- memcpy(target.fLcutDCA, fLcutDCA, sizeof(Float_t) * 2);\r
- memcpy(target.fLcutVertexR, fLcutVertexR, sizeof(Float_t) * 2);\r
- memcpy(target.fLcutInvMass, fLcutInvMass, sizeof(Float_t) * 2);\r
- \r
-}\r
-//____________________________________________________________________\r
-Bool_t AliESDv0KineCuts::ProcessV0(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const \r
-{\r
- //\r
- // main user function\r
- //\r
-\r
- if(!v0) return kFALSE;\r
- if(!fEvent){\r
- AliErrorClass("No valid Event pointer available, provide it first");\r
- return kFALSE;\r
- }\r
-\r
- if(!V0CutsCommon(v0)) return kFALSE;\r
-\r
- const Int_t id = PreselectV0(v0);\r
-\r
- if(!SingleTrackCuts(v0)) return kFALSE;\r
-\r
- switch(id){\r
- case kUndef:\r
- return kFALSE;\r
- case kGamma:\r
- return CaseGamma(v0, pdgV0, pdgP, pdgN);\r
- case kK0:\r
- return CaseK0(v0, pdgV0, pdgP, pdgN);\r
- case kLambda:\r
- return CaseLambda(v0, pdgV0, pdgP, pdgN, 0);\r
- case kALambda:\r
- return CaseLambda(v0, pdgV0, pdgP, pdgN, 1);\r
- default:\r
- return kFALSE; \r
- }\r
-\r
- return kFALSE;\r
-}\r
-//____________________________________________________________________\r
-Bool_t AliESDv0KineCuts::ProcessV0(AliESDv0* const v0, Int_t &pdgP, Int_t &pdgN) const \r
-{\r
- //\r
- // main user function, simplified if the V0 identity is not necessary\r
- //\r
-\r
- if(!v0) return kFALSE;\r
- if(!fEvent){\r
- AliErrorClass("No valid Event pointer available, provide it first");\r
- return kFALSE;\r
- }\r
-\r
- Int_t idV0 = -1;\r
- return ProcessV0(v0, idV0, pdgP, pdgN);\r
-\r
-}\r
-//____________________________________________________________________\r
-Int_t AliESDv0KineCuts::PreselectV0(AliESDv0* const v0) const \r
-{\r
- //\r
- // Make a preselection (exclusive) of the V0 cadidates based on\r
- // Armenteros plot\r
- // the armenteros cut values are currently fixed and user is not able to set them via\r
- // set funcions. The reason is that these cuts are optimized and furneter changes should \r
- // not be necessary. To prove otherwise please study in detail before changing the values\r
- //\r
- \r
- Float_t ap[2] = {-1., -1.};\r
- Armenteros(v0, ap);\r
- // for clarity\r
- const Float_t alpha = ap[0];\r
- const Float_t qt = ap[1];\r
-\r
- // selection cuts \r
- // - the reagions for different candidates must not overlap \r
-\r
- // Gamma cuts\r
- const Double_t cutAlphaG = 0.35; \r
- const Double_t cutQTG = 0.05;\r
- const Double_t cutAlphaG2[2] = {0.6, 0.8};\r
- const Double_t cutQTG2 = 0.04;\r
-\r
- // K0 cuts\r
- const Float_t cutQTK0[2] = {0.1075, 0.215};\r
- const Float_t cutAPK0[2] = {0.199, 0.8}; // parameters for curved QT cut\r
- \r
- // Lambda & A-Lambda cuts\r
- const Float_t cutQTL = 0.03;\r
- const Float_t cutAlphaL[2] = {0.35, 0.7};\r
- const Float_t cutAlphaAL[2] = {-0.7, -0.35};\r
- const Float_t cutAPL[3] = {0.107, -0.69, 0.5}; // parameters fir curved QT cut\r
-\r
-\r
- if(kPurity == fMode){\r
- // Check for Gamma candidates\r
- if(qt < cutQTG){\r
- if( (TMath::Abs(alpha) < cutAlphaG) ) return kGamma;\r
- }\r
- // additional region - should help high pT gammas\r
- if(qt < cutQTG2){\r
- if( (TMath::Abs(alpha) > cutAlphaG2[0]) && (TMath::Abs(alpha) < cutAlphaG2[1]) ) return kGamma;\r
- }\r
- }\r
- if(kEffGamma == fMode){\r
- if(qt < cutQTG) return kGamma;\r
- }\r
-\r
- \r
- // Check for K0 candidates\r
- Float_t q = cutAPK0[0] * TMath::Sqrt(TMath::Abs(1 - alpha*alpha/(cutAPK0[1]*cutAPK0[1])));\r
- if( (qt > cutQTK0[0]) && (qt < cutQTK0[1]) && (qt > q) ){\r
- return kK0;\r
- }\r
-\r
- // Check for Lambda candidates\r
- q = cutAPL[0] * TMath::Sqrt(TMath::Abs(1 - ( (alpha + cutAPL[1]) * (alpha + cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) ));\r
- if( (alpha > cutAlphaL[0]) && (alpha < cutAlphaL[1]) && (qt > cutQTL) && (qt < q) ){\r
- return kLambda;\r
- }\r
-\r
- // Check for A-Lambda candidates\r
- q = cutAPL[0] * TMath::Sqrt(TMath::Abs(1 - ( (alpha - cutAPL[1]) * (alpha - cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) ));\r
- if( (alpha > cutAlphaAL[0]) && (alpha < cutAlphaAL[1]) && (qt > cutQTL) && (qt < q) ){\r
- return kALambda;\r
- }\r
- \r
- return kUndef;\r
-}\r
-//____________________________________________________________________\r
-Bool_t AliESDv0KineCuts::SingleTrackCuts(AliESDv0 * const v0) const \r
-{\r
- //\r
- // apply single track cuts\r
- // correct sign not relevat here\r
- //\r
-\r
- if(!v0) return kFALSE;\r
- \r
- Int_t pIndex = 0, nIndex = 0;\r
- pIndex = v0->GetPindex();\r
- nIndex = v0->GetNindex();\r
- AliESDtrack* d[2];\r
- d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));\r
- d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));\r
- \r
- for(Int_t i=0; i<2; ++i){\r
- if(!d[i]) return kFALSE;\r
- \r
- // status word\r
- ULong_t status = d[i]->GetStatus();\r
-\r
- // No. of TPC clusters leave to the users\r
- if(d[i]->GetTPCNcls() < 1) return kFALSE;\r
-\r
- // TPC refit\r
- if(!(status & AliESDtrack::kTPCrefit)) return kFALSE;\r
- \r
- // Chi2 per TPC cluster\r
- Int_t nTPCclusters = d[i]->GetTPCNcls();\r
- Float_t chi2perTPCcluster = d[i]->GetTPCchi2()/Float_t(nTPCclusters);\r
- if(chi2perTPCcluster > 4) return kFALSE;\r
-\r
- // TPC cluster ratio\r
- Float_t cRatioTPC = d[i]->GetTPCNclsF() > 0. ? static_cast<Float_t>(d[i]->GetTPCNcls())/static_cast<Float_t> (d[i]->GetTPCNclsF()) : 1.;\r
- if(cRatioTPC < 0.6) return kFALSE;\r
- \r
- // kinks\r
- if(d[i]->GetKinkIndex(0) != 0) return kFALSE;\r
- \r
- }\r
-\r
- return kTRUE;\r
-}\r
-//____________________________________________________________________\r
-Bool_t AliESDv0KineCuts::CaseGamma(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const \r
-{\r
- //\r
- // process the gamma conversion candidate\r
- //\r
-\r
- if(!v0) return kFALSE;\r
-\r
- AliVTrack* daughter[2];\r
- Int_t pIndex = 0, nIndex = 0;\r
-\r
- Bool_t sign = CheckSigns(v0);\r
- if(sign){\r
- pIndex = v0->GetPindex();\r
- nIndex = v0->GetNindex();\r
- }\r
- else{\r
- pIndex = v0->GetNindex();\r
- nIndex = v0->GetPindex(); \r
- }\r
- daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));\r
- daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));\r
- if(!daughter[0] || !daughter[1]) return kFALSE;\r
-\r
- AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kElectron), TMath::Abs(kElectron));\r
- if(!kfMother) return kFALSE;\r
-\r
- AliESDtrack* d[2];\r
- d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));\r
- d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));\r
-\r
- Float_t iMass = v0->GetEffMass(0, 0);\r
-\r
- // cos pointing angle\r
- Double_t cosPoint = v0->GetV0CosineOfPointingAngle();\r
- cosPoint = TMath::ACos(cosPoint);\r
-\r
- // DCA between daughters\r
- Double_t dca = v0->GetDcaV0Daughters();\r
-\r
- // Production vertex\r
- Double_t x, y, z; \r
- v0->GetXYZ(x,y,z);\r
- Double_t r = TMath::Sqrt(x*x + y*y);\r
-\r
- Double_t xy[2];\r
- Double_t r2 = -1.;\r
- if ( GetConvPosXY(d[0], d[1], xy) ){\r
- r2 = TMath::Sqrt(xy[0]*xy[0] + xy[1]*xy[1]);\r
- }\r
-\r
- // psi pair \r
- Double_t psiPair = PsiPair(v0);\r
- \r
- // V0 chi2/ndf\r
- Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();\r
-\r
- if(kfMother) delete kfMother; \r
- \r
- // apply the cuts\r
-\r
- if(iMass > fGcutInvMass) return kFALSE;\r
-\r
- if(chi2ndf > fGcutChi2NDF) return kFALSE;\r
-\r
- if(cosPoint < fGcutCosPoint[0] || cosPoint > fGcutCosPoint[1]) return kFALSE;\r
-\r
- if(dca < fGcutDCA[0] || dca > fGcutDCA[1]) return kFALSE;\r
-\r
- if(r < fGcutVertexR[0] || r > fGcutVertexR[1]) return kFALSE;\r
-\r
- if(psiPair < fGcutPsiPair[0] || psiPair > fGcutPsiPair[1]) return kFALSE;\r
- \r
- // all cuts passed\r
-\r
- pdgV0 = 22;\r
- if(sign){\r
- pdgP = -11;\r
- pdgN = 11;\r
- }\r
- else{\r
- pdgP = 11;\r
- pdgN = -11;\r
- }\r
-\r
- return kTRUE;\r
-}\r
-//____________________________________________________________________\r
-Bool_t AliESDv0KineCuts::CaseK0(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const {\r
- //\r
- // process the K0 candidate\r
- //\r
-\r
- if(!v0) return kFALSE;\r
- \r
- AliVTrack* daughter[2];\r
- Int_t pIndex = 0, nIndex = 0;\r
- Bool_t sign = CheckSigns(v0);\r
- if(sign){\r
- pIndex = v0->GetPindex();\r
- nIndex = v0->GetNindex();\r
- }\r
- else{\r
- pIndex = v0->GetNindex();\r
- nIndex = v0->GetPindex(); \r
- }\r
- \r
- daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));\r
- daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));\r
- if(!daughter[0] || !daughter[1]) return kFALSE;\r
-\r
- AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kPiPlus));\r
- if(!kfMother) return kFALSE;\r
-\r
- AliESDtrack* d[2];\r
- d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));\r
- d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));\r
-\r
- Float_t iMass = v0->GetEffMass(2, 2);\r
-\r
- // cos pointing angle\r
- Double_t cosPoint = v0->GetV0CosineOfPointingAngle();\r
- cosPoint = TMath::ACos(cosPoint);\r
-\r
- // DCA between daughters\r
- Double_t dca = v0->GetDcaV0Daughters();\r
-\r
- // Production vertex\r
- Double_t x, y, z; \r
- v0->GetXYZ(x,y,z);\r
-\r
- Double_t r = TMath::Sqrt(x*x + y*y); \r
-\r
- // V0 chi2/ndf\r
- Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();\r
- \r
- if(kfMother) delete kfMother; \r
-\r
- //\r
- // apply the cuts\r
- //\r
- if(iMass < fK0cutInvMass[0] || iMass > fK0cutInvMass[1]) return kFALSE;\r
-\r
- if(chi2ndf > fK0cutChi2NDF) return kFALSE;\r
-\r
- if(cosPoint < fK0cutCosPoint[0] || cosPoint > fK0cutCosPoint[1]) return kFALSE;\r
-\r
- if(dca < fK0cutDCA[0] || dca > fK0cutDCA[1]) return kFALSE;\r
-\r
- if(r < fK0cutVertexR[0] || r > fK0cutVertexR[1]) return kFALSE;\r
-\r
- // all cuts passed\r
- pdgV0 = 310;\r
- if(sign){\r
- pdgP = 211;\r
- pdgN = -211;\r
- }\r
- else{\r
- pdgP = -211;\r
- pdgN = 211;\r
- }\r
-\r
- return kTRUE;\r
-}\r
-//____________________________________________________________________\r
-Bool_t AliESDv0KineCuts::CaseLambda(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN, Int_t id) const {\r
- //\r
- // process teh Lambda and Anti-Lambda candidate\r
- //\r
- \r
- if(!v0) return kFALSE;\r
-\r
- const Double_t cL0mass=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(); // PDG lambda mass\r
-\r
- AliVTrack* daughter[2];\r
- Int_t pIndex = 0, nIndex = 0;\r
- Float_t mMass[2] = {-1., -1.};\r
- Bool_t sign = CheckSigns(v0);\r
- if(sign){\r
- pIndex = v0->GetPindex();\r
- nIndex = v0->GetNindex();\r
- mMass[0] = v0->GetEffMass(4, 2);\r
- mMass[1] = v0->GetEffMass(2, 4);\r
- }\r
- else{\r
- pIndex = v0->GetNindex();\r
- nIndex = v0->GetPindex(); \r
- mMass[0] = v0->GetEffMass(2, 4);\r
- mMass[1] = v0->GetEffMass(4, 2);\r
- }\r
- \r
- daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));\r
- daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));\r
- if(!daughter[0] || !daughter[1]) return kFALSE;\r
-\r
- AliKFParticle *kfMother[2] = {0x0, 0x0};\r
- // Lambda\r
- kfMother[0] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kProton), TMath::Abs(kPiPlus));\r
- if(!kfMother[0]) return kFALSE;\r
- \r
- // Anti-Lambda\r
- kfMother[1] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kProton));\r
- if(!kfMother[1]) return kFALSE;\r
-\r
- Float_t dMass[2] = {TMath::Abs(mMass[0] - cL0mass), TMath::Abs(mMass[1] - cL0mass)};\r
- \r
- AliESDtrack* d[2];\r
- d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));\r
- d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));\r
- if(!d[0] || !d[1]) return kFALSE;\r
- \r
- Float_t p[2] = {d[0]->GetP(), d[1]->GetP()}; \r
-\r
- // check the 3 lambda - antilambda variables\r
- Int_t check[2] = {-1, -1}; // 0 : lambda, 1 : antilambda\r
- // 1) momentum of the daughter particles - proton is expected to have higher momentum than pion\r
- check[0] = (p[0] > p[1]) ? 0 : 1;\r
- // 2) mass of the mother particle\r
- check[1] = (dMass[0] < dMass[1]) ? 0 : 1;\r
- \r
- // require positive correlation of (1) and (2)\r
- if(check[0] != check[1]){\r
- if(kfMother[0]) delete kfMother[0]; \r
- if(kfMother[1]) delete kfMother[1]; \r
- return kFALSE;\r
- }\r
-\r
- // now that the check[0] == check[1]\r
- const Int_t type = check[0];\r
-\r
- // require that the input armenteros preselection agree:\r
- if(type != id) return kFALSE;\r
-\r
- Float_t iMass =0.;\r
- if(sign){\r
- iMass = (type == 0) ? v0->GetEffMass(4, 2) : v0->GetEffMass(2, 4);\r
- }\r
- else{\r
- iMass = (type == 0) ? v0->GetEffMass(2, 4) : v0->GetEffMass(4, 2);\r
- }\r
-\r
- // cos pointing angle\r
- Double_t cosPoint = v0->GetV0CosineOfPointingAngle();\r
- cosPoint = TMath::ACos(cosPoint);\r
-\r
- // DCA between daughters\r
- Double_t dca = v0->GetDcaV0Daughters();\r
- \r
- // Production vertex\r
- Double_t x, y, z; \r
- v0->GetXYZ(x,y,z);\r
- Double_t r = TMath::Sqrt(x*x + y*y);\r
-\r
- // proton - pion indices\r
- Int_t ix[2] = {0, 1};\r
- if(1 == type){\r
- ix[0] = 1;\r
- ix[1] = 0;\r
- }\r
-\r
- // V0 chi2/ndf\r
- Double_t chi2ndf = kfMother[type]->GetChi2()/kfMother[type]->GetNDF();\r
-\r
- if(kfMother[0]) delete kfMother[0]; \r
- if(kfMother[1]) delete kfMother[1]; \r
-\r
- //\r
- // apply the cuts\r
- //\r
-\r
- if(iMass < fLcutInvMass[0] || iMass > fLcutInvMass[1]) return kFALSE;\r
-\r
- if(chi2ndf > fLcutChi2NDF) return kFALSE;\r
-\r
- if(cosPoint < fLcutCosPoint[0] || cosPoint > fLcutCosPoint[1]) return kFALSE;\r
-\r
- if(dca < fLcutDCA[0] || dca > fLcutDCA[1]) return kFALSE;\r
-\r
- if(r < fLcutVertexR[0] || r > fLcutVertexR[1]) return kFALSE;\r
-\r
- // all cuts passed\r
-\r
- if(0 == type){\r
- pdgV0 = 3122;\r
- if(sign){\r
- pdgP = 2212;\r
- pdgN = -211;\r
- }\r
- else{\r
- pdgP = -211;\r
- pdgN = 2212;\r
- }\r
- }\r
- else{\r
- pdgV0 = -3122;\r
- if(sign){\r
- pdgP = 211;\r
- pdgN = -2212;\r
- }\r
- else{\r
- pdgP = -2212;\r
- pdgN = 211;\r
- }\r
- }\r
-\r
- return kTRUE;\r
-}\r
-//____________________________________________________________________\r
-Bool_t AliESDv0KineCuts::V0CutsCommon(const AliESDv0 * const v0) const \r
-{\r
- //\r
- // V0 cuts common to all V0s\r
- //\r
-\r
- AliESDtrack* dN, *dP; \r
- \r
- dP = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(v0->GetPindex()));\r
- dN = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(v0->GetNindex())); \r
- \r
- if(!dN || !dP) return kFALSE;\r
-\r
- Int_t qP = dP->Charge();\r
- Int_t qN = dN->Charge();\r
-\r
- if((qP*qN) != -1) return kFALSE;\r
-\r
- return kTRUE;\r
-}\r
-//____________________________________________________________________\r
-void AliESDv0KineCuts::Armenteros(AliESDv0* const v0, Float_t val[2]) const \r
-{\r
- //\r
- // computes the Armenteros variables for given V0\r
- // fills the histogram\r
- // returns the values via "val"\r
- //\r
- \r
- Double_t mn[3] = {0,0,0};\r
- Double_t mp[3] = {0,0,0}; \r
- Double_t mm[3] = {0,0,0}; \r
-\r
- if(CheckSigns(v0)){\r
- v0->GetNPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter\r
- v0->GetPPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter\r
- }\r
- else{\r
- v0->GetPPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter\r
- v0->GetNPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter\r
- }\r
- v0->GetPxPyPz(mm[0],mm[1],mm[2]); //reconstructed cartesian momentum components of mother\r
-\r
- TVector3 vecN(mn[0],mn[1],mn[2]);\r
- TVector3 vecP(mp[0],mp[1],mp[2]);\r
- TVector3 vecM(mm[0],mm[1],mm[2]);\r
- \r
- Double_t thetaP = acos((vecP * vecM)/(vecP.Mag() * vecM.Mag()));\r
- Double_t thetaN = acos((vecN * vecM)/(vecN.Mag() * vecM.Mag()));\r
- \r
- Double_t alfa = ((vecP.Mag())*cos(thetaP)-(vecN.Mag())*cos(thetaN))/\r
- ((vecP.Mag())*cos(thetaP)+(vecN.Mag())*cos(thetaN)) ;\r
- Double_t qt = vecP.Mag()*sin(thetaP);\r
-\r
- val[0] = alfa;\r
- val[1] = qt;\r
-}\r
-//____________________________________________________________________\r
-Bool_t AliESDv0KineCuts::CheckSigns(AliESDv0* const v0) const \r
-{\r
- //\r
- // check wheter the sign was correctly applied to \r
- // V0 daughter tracks\r
- //\r
- \r
- Bool_t correct = kFALSE;\r
-\r
- Int_t pIndex = 0, nIndex = 0;\r
- pIndex = v0->GetPindex();\r
- nIndex = v0->GetNindex();\r
- \r
- AliESDtrack* d[2];\r
- d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));\r
- d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));\r
-\r
- Int_t sign[2];\r
- sign[0] = (int)d[0]->GetSign();\r
- sign[1] = (int)d[1]->GetSign();\r
- \r
- if(-1 == sign[0] && 1 == sign[1]){\r
- correct = kFALSE;\r
- }\r
- else{\r
- correct = kTRUE;\r
- }\r
- \r
- return correct;\r
-}\r
-//________________________________________________________________\r
-Double_t AliESDv0KineCuts::PsiPair(AliESDv0* const v0) const \r
-{\r
- //\r
- // Angle between daughter momentum plane and plane \r
- // \r
-\r
- if(!fEvent) return -1.;\r
-\r
- Float_t magField = fEvent->GetMagneticField();\r
-\r
- Int_t pIndex = -1;\r
- Int_t nIndex = -1;\r
- if(CheckSigns(v0)){\r
- pIndex = v0->GetPindex();\r
- nIndex = v0->GetNindex();\r
- }\r
- else{\r
- pIndex = v0->GetNindex();\r
- nIndex = v0->GetPindex(); \r
- }\r
- \r
-\r
- AliESDtrack* daughter[2];\r
-\r
- daughter[0] = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(pIndex));\r
- daughter[1] = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(nIndex));\r
-\r
- Double_t x, y, z;\r
- v0->GetXYZ(x,y,z);//Reconstructed coordinates of V0; to be replaced by Markus Rammler's method in case of conversions!\r
- \r
- Double_t mn[3] = {0,0,0};\r
- Double_t mp[3] = {0,0,0};\r
- \r
-\r
- v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;\r
- v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter; \r
-\r
-\r
- Double_t deltat = 1.;\r
- deltat = TMath::ATan(mp[2]/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1])+1.e-13)) - TMath::ATan(mn[2]/(TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1])+1.e-13));//difference of angles of the two daughter tracks with z-axis\r
-\r
- Double_t radiussum = TMath::Sqrt(x*x + y*y) + 50;//radius to which tracks shall be propagated\r
-\r
- Double_t momPosProp[3];\r
- Double_t momNegProp[3];\r
- \r
- AliExternalTrackParam pt(*daughter[0]), nt(*daughter[1]);\r
- \r
- Double_t psiPair = 4.;\r
-\r
- if(nt.PropagateTo(radiussum,magField) == 0)//propagate tracks to the outside\r
- psiPair = -5.;\r
- if(pt.PropagateTo(radiussum,magField) == 0)\r
- psiPair = -5.;\r
- pt.GetPxPyPz(momPosProp);//Get momentum vectors of tracks after propagation\r
- nt.GetPxPyPz(momNegProp);\r
- \r
- Double_t pEle =\r
- TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter\r
- Double_t pPos =\r
- TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter\r
- \r
- Double_t scalarproduct =\r
- momPosProp[0]*momNegProp[0]+momPosProp[1]*momNegProp[1]+momPosProp[2]*momNegProp[2];//scalar product of propagated positive and negative daughters' momenta\r
- \r
- Double_t chipair = TMath::ACos(scalarproduct/(pEle*pPos));//Angle between propagated daughter tracks\r
-\r
- psiPair = TMath::Abs(TMath::ASin(deltat/chipair)); \r
-\r
- return psiPair; \r
-}\r
-//___________________________________________________________________\r
-Bool_t AliESDv0KineCuts::GetConvPosXY(AliESDtrack * const ptrack, AliESDtrack * const ntrack, Double_t convpos[2]) const\r
-{\r
- //\r
- // recalculate the gamma conversion XY postition\r
- //\r
-\r
- const Double_t b = fEvent->GetMagneticField();\r
-\r
- Double_t helixcenterpos[2];\r
- GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos);\r
-\r
- Double_t helixcenterneg[2];\r
- GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg);\r
-\r
- Double_t poshelix[6];\r
- ptrack->GetHelixParameters(poshelix,b);\r
- Double_t posradius = TMath::Abs(1./poshelix[4]);\r
-\r
- Double_t neghelix[6];\r
- ntrack->GetHelixParameters(neghelix,b);\r
- Double_t negradius = TMath::Abs(1./neghelix[4]);\r
-\r
- Double_t xpos = helixcenterpos[0];\r
- Double_t ypos = helixcenterpos[1];\r
- Double_t xneg = helixcenterneg[0];\r
- Double_t yneg = helixcenterneg[1];\r
-\r
- convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);\r
- convpos[1] = (ypos*negradius+ yneg*posradius)/(negradius+posradius);\r
-\r
- return 1;\r
-}\r
-//___________________________________________________________________\r
-Bool_t AliESDv0KineCuts::GetHelixCenter(AliESDtrack * const track, Double_t b,Int_t charge, Double_t center[2]) const\r
-{\r
- //\r
- // computes the center of the track helix\r
- //\r
- \r
- Double_t pi = TMath::Pi();\r
- \r
- Double_t helix[6];\r
- track->GetHelixParameters(helix,b);\r
- \r
- Double_t xpos = helix[5];\r
- Double_t ypos = helix[0];\r
- Double_t radius = TMath::Abs(1./helix[4]);\r
- Double_t phi = helix[2];\r
-\r
- if(phi < 0){\r
- phi = phi + 2*pi;\r
- }\r
-\r
- phi -= pi/2.;\r
- Double_t xpoint = radius * TMath::Cos(phi);\r
- Double_t ypoint = radius * TMath::Sin(phi);\r
-\r
- if(b<0){\r
- if(charge > 0){\r
- xpoint = - xpoint;\r
- ypoint = - ypoint;\r
- }\r
- /* avoid self assignment\r
- if(charge < 0){\r
- xpoint = xpoint;\r
- ypoint = ypoint;\r
- }\r
- */\r
- }\r
- if(b>0){\r
- /* avoid self assignment\r
- if(charge > 0){\r
- xpoint = xpoint;\r
- ypoint = ypoint;\r
- }\r
- */\r
- if(charge < 0){\r
- xpoint = - xpoint;\r
- ypoint = - ypoint;\r
- }\r
- }\r
- center[0] = xpos + xpoint;\r
- center[1] = ypos + ypoint;\r
-\r
- return 1;\r
-}\r
-//___________________________________________________________________\r
-AliKFParticle *AliESDv0KineCuts::CreateMotherParticle(const AliVTrack* const pdaughter, const AliVTrack* const ndaughter, Int_t pspec, Int_t nspec) const\r
-{\r
- //\r
- // Creates a mother particle\r
- //\r
- AliKFParticle pkfdaughter(*pdaughter, pspec);\r
- AliKFParticle nkfdaughter(*ndaughter, nspec);\r
- \r
- \r
- // Create the mother particle \r
- AliKFParticle *m = new AliKFParticle(pkfdaughter, nkfdaughter);\r
- m->SetField(fEvent->GetMagneticField());\r
- if(TMath::Abs(kElectron) == pspec && TMath::Abs(kElectron) == nspec) m->SetMassConstraint(0, 0.001);\r
- else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(), 0.);\r
- else if(TMath::Abs(kProton) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);\r
- else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kProton) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);\r
- else{\r
- AliErrorClass("Wrong daughter ID - mass constraint can not be set");\r
- }\r
-\r
- AliKFVertex improvedVertex = *fPrimaryVertex;\r
- improvedVertex += *m;\r
- m->SetProductionVertex(improvedVertex);\r
- \r
- // update 15/06/2010\r
- // mother particle will not be added to primary vertex but only to its copy \r
- // as this confilcts with calling\r
- // m->SetPrimaryVertex() function and\r
- // subsequently removing the mother particle afterwards\r
- // Source: Sergey Gorbunov\r
-\r
- return m;\r
-}\r
-//____________________________________________________________________\r
-void AliESDv0KineCuts::SetEvent(AliESDEvent* const event){\r
- //\r
- // direct setter of ESD event\r
- //\r
- fEvent = event;\r
- if(!fEvent){\r
- AliErrorClass("Invalid input event pointer");\r
- return;\r
- }\r
-if (fUseExternalVertex) return;\r
-else{\r
- if(fPrimaryVertex && fDeleteVertex){\r
- delete fPrimaryVertex;\r
- fPrimaryVertex=0x0;\r
- }\r
- fPrimaryVertex = new AliKFVertex(*(fEvent->GetPrimaryVertex()));\r
- fDeleteVertex=kTRUE;\r
- }\r
-\r
-\r
-\r
-}\r
-//____________________________________________________________________\r
-void AliESDv0KineCuts::SetEvent(AliVEvent* const event){\r
- //\r
- // direct setter of ESD event\r
- //\r
-\r
- fEvent = dynamic_cast<AliESDEvent*>(event);\r
- if(!fEvent){\r
- AliErrorClass("Invalid input event pointer");\r
- return;\r
- }\r
- \r
- if (fUseExternalVertex) return;\r
- else{\r
- if(fPrimaryVertex && fDeleteVertex){\r
- delete fPrimaryVertex;\r
- fPrimaryVertex=0x0;\r
- }\r
- fPrimaryVertex = new AliKFVertex(*(fEvent->GetPrimaryVertex()));\r
- fDeleteVertex=kTRUE;\r
- }\r
-}\r
-\r
-\r
-//________________________________________________________________\r
-void AliESDv0KineCuts::UseExternalVertex(Bool_t use_external){\r
- //\r
- // Reenable primary Vertex from ESD event\r
- //\r
- if (use_external) fUseExternalVertex =kTRUE;\r
- else fUseExternalVertex =kFALSE;\r
-}\r
-\r
-\r
-\r
-\r
-//________________________________________________________________\r
-void AliESDv0KineCuts::SetPrimaryVertex(AliKFVertex* const v){\r
- //\r
- // set the primary vertex of the event\r
- //\r
- if(fPrimaryVertex && fDeleteVertex){ \r
- delete fPrimaryVertex;\r
- fPrimaryVertex =0x0;\r
- fDeleteVertex = kFALSE;\r
- } \r
- fUseExternalVertex=kTRUE; \r
- fPrimaryVertex = v; // set primary Vertex\r
- if(!fPrimaryVertex){\r
- AliErrorClass("Failed to initialize the primary vertex");\r
- return;\r
- }\r
-}\r
-//___________________________________________________________________\r
-void AliESDv0KineCuts::SetMode(Int_t mode, Int_t type){\r
- //\r
- // this function allows the user to select (prior running the 'ProcessV0' function)\r
- // to select different approaches to V0 selection - the 'mode'\r
- // - and -\r
- // different systems (pp, PbPb) - 'type' \r
- //\r
- // To see the cut values for different modes please refer to the\r
- // function SetCuts()\r
- //\r
- // Important notice: based on the parameters particular sets of cuts will\r
- // be activated for teh V0 selection. If some additional changes to single\r
- // cuts are needed please us the SetXXXcut function (see the header file)\r
- // \r
-\r
- switch(mode){\r
- case kPurity:\r
- fMode = kPurity; // used to obtain highest purity possible - the efficiency may be low\r
- break;\r
- case kEffGamma:\r
- fMode = kEffGamma; // used to obtain highes efficiency possible - the purity may be worse\r
- break;\r
- default:\r
- AliError("V0 selection mode not recognozed, setting 'kPurity'");\r
- fMode = kPurity;\r
- }\r
-\r
- switch(type){\r
- case kPP:\r
- fType = kPP; // cuts optimized for low multiplicity \r
- break;\r
- case kPbPb:\r
- fType = kPbPb; // cuts optimized for high multiplicity\r
- break;\r
- }\r
- \r
- // setup the cut values for selected mode & type\r
- SetCuts();\r
-\r
-}\r
-//___________________________________________________________________\r
-void AliESDv0KineCuts::SetMode(Int_t mode, const char* type){\r
- //\r
- // overloaded function - please see above\r
- // \r
- \r
- Int_t t = -1;\r
-\r
- if(!strcmp("pp", type)) t = kPP;\r
- else if(!(strcmp("PbPb", type))) t = kPbPb;\r
- else{\r
- AliError("data type not recognized, setting 'pp'");\r
- t = kPP; \r
- }\r
-\r
- SetMode(mode, t);\r
-\r
-}\r
-//___________________________________________________________________\r
-void AliESDv0KineCuts::SetCuts(){\r
- //\r
- // this funciton sets the default cut values based on the selected\r
- // fMode and fType.\r
- // please note that only the cuts that have different values than the default\r
- // cuts are updated here\r
- //\r
- \r
- // last update: 14/02/2011\r
- // as a very preliminary - the only change to default cuts is to apply\r
- // less restricting gamma conversion selection in PreselectV0() function\r
- \r
-\r
- \r
-}\r
+/**************************************************************************
+ * 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. *
+ **************************************************************************/
+
+/*
+ * author: M.Kalisky@gsi.de
+ * 08/Dec/2010
+ *
+ * Description: This class allows with purely kinematical cuts
+ * to select clean samples of electrons, pions and protons from the
+ * V0 online finder ESD V0 candidates for PID and dectector resonse
+ * studies.
+ */
+
+#include <TVector3.h>
+#include <TDatabasePDG.h>
+
+#include "AliESDv0.h"
+#include "AliESDtrack.h"
+#include "AliESDEvent.h"
+#include "AliVEvent.h"
+#include "AliLog.h"
+#include "AliKFParticle.h"
+#include "AliVTrack.h"
+#include "AliKFVertex.h"
+
+#include "AliESDv0KineCuts.h"
+
+ClassImp(AliESDv0KineCuts)
+
+//____________________________________________________________________
+AliESDv0KineCuts::AliESDv0KineCuts() :
+ fEvent(0x0)
+ , fPrimaryVertex(0x0)
+ , fType(0)
+ , fMode(0)
+ , fTPCNcls(1)
+ , fTPCrefit(kTRUE)
+ , fTPCchi2perCls(4.0)
+ , fTPCclsRatio(0.6)
+ , fNoKinks(kTRUE)
+ , fGcutChi2NDF(10)
+ , fGcutInvMass(0.05)
+ , fK0cutChi2NDF(10)
+ , fLcutChi2NDF(10)
+ , fUseExternalVertex(kFALSE)
+ , fDeleteVertex(kFALSE)
+{
+ //
+ // Default constructor
+ //
+
+ // default single track cuts
+ fTPCNcls = 1; // minimal number of the TPC clusters
+ fTPCrefit = kTRUE; // TPC refit
+ fTPCchi2perCls = 4.0; // chi2 per TPC cluster
+ fTPCclsRatio = 0.6; // minimal foun/findable TPC cluster ratio
+ fNoKinks = kTRUE; // kinks - no [kTRUE] or do not care [kFalse]
+
+
+ // default gamma cuts values
+ fGcutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle gamma
+ fGcutCosPoint[0] = 0; // cos of the pointing angle [min, max]
+ fGcutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
+ fGcutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
+ fGcutDCA[1] = 0.25; // DCA between the daughter tracks [min, max]
+ fGcutVertexR[0] = 3.; // radius of the conversion point [min, max]
+ fGcutVertexR[1] = 90.; // radius of the conversion point [min, max]
+ fGcutPsiPair[0] = 0.; // value of the psi pair cut [min, max]
+ fGcutPsiPair[1] = 0.05; // value of the psi pair cut [min, max]
+ fGcutInvMass = 0.05; // upper value on the gamma invariant mass
+ // default K0 cuts
+ fK0cutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle K0
+ fK0cutCosPoint[0] = 0.; // cos of the pointing angle [min, max]
+ fK0cutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
+ fK0cutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
+ fK0cutDCA[1] = 0.2; // DCA between the daughter tracks [min, max]
+ fK0cutVertexR[0] = 2.0; // radius of the decay point [min, max]
+ fK0cutVertexR[1] = 30.0; // radius of the decay point [min, max]
+ fK0cutInvMass[0] = 0.486; // invariant mass window
+ fK0cutInvMass[1] = 0.508; // invariant mass window
+ // Lambda & anti-Lambda cut values
+ fLcutChi2NDF = 10; // Chi2NF cut value for the AliKFparticle K0
+ fLcutCosPoint[0] = 0.; // cos of the pointing angle [min, max]
+ fLcutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
+ fLcutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
+ fLcutDCA[1] = 0.2; // DCA between the daughter tracks [min, max]
+ fLcutVertexR[0] = 2.0; // radius of the decay point [min, max]
+ fLcutVertexR[1] = 40.0; // radius of the decay point [min, max]
+ fLcutInvMass[0] = 1.11; // invariant mass window
+ fLcutInvMass[1] = 1.12; // invariant mass window
+
+}
+//____________________________________________________________________
+AliESDv0KineCuts::~AliESDv0KineCuts(){
+ //
+ // Destructor
+ //
+
+
+}
+//____________________________________________________________________
+AliESDv0KineCuts::AliESDv0KineCuts(const AliESDv0KineCuts &ref):
+ TObject(ref)
+ , fEvent(0x0)
+ , fPrimaryVertex(0x0)
+ , fType(0)
+ , fMode(0)
+ , fTPCNcls(1)
+ , fTPCrefit(kTRUE)
+ , fTPCchi2perCls(4.0)
+ , fTPCclsRatio(0.6)
+ , fNoKinks(kTRUE)
+ , fGcutChi2NDF(10)
+ , fGcutInvMass(0.05)
+ , fK0cutChi2NDF(10)
+ , fLcutChi2NDF(10)
+ , fUseExternalVertex(kFALSE)
+ , fDeleteVertex(kFALSE)
+{
+ //
+ // Copy operator
+ //
+
+ ref.Copy(*this);
+}
+//____________________________________________________________________
+AliESDv0KineCuts &AliESDv0KineCuts::operator=(const AliESDv0KineCuts &ref){
+ //
+ // assignment operator
+ //
+ if(this != &ref)
+ ref.Copy(*this);
+ return *this;
+}
+//____________________________________________________________________
+void AliESDv0KineCuts::Copy(TObject &ref) const {
+ //
+ // Performs the copying of the object
+ //
+
+ TObject::Copy(ref);
+
+ AliESDv0KineCuts &target = dynamic_cast<AliESDv0KineCuts &>(ref);
+
+ // default single track cuts
+ target.fTPCNcls = fTPCNcls;
+ target.fTPCrefit = fTPCrefit;
+ target.fTPCchi2perCls = fTPCchi2perCls;
+ target.fTPCclsRatio = fTPCclsRatio;
+ target.fNoKinks = fNoKinks;
+ target.fUseExternalVertex = fUseExternalVertex; //added december 2nd 2011
+ target.fDeleteVertex = fDeleteVertex; //added december 2nd 2011
+
+ // default gamma cuts values
+ target.fGcutChi2NDF = fGcutChi2NDF;
+ memcpy(target.fGcutCosPoint, fGcutCosPoint, sizeof(Float_t) * 2);
+ memcpy(target.fGcutDCA, fGcutDCA, sizeof(Float_t) * 2);
+ memcpy(target.fGcutVertexR, fGcutVertexR, sizeof(Float_t) * 2);
+ memcpy(target.fGcutPsiPair, fGcutPsiPair, sizeof(Float_t) * 2);
+ target.fGcutInvMass = fGcutInvMass;
+ // default K0 cuts
+ target.fK0cutChi2NDF = fK0cutChi2NDF;
+ memcpy(target.fK0cutCosPoint, fK0cutCosPoint, sizeof(Float_t) * 2);
+ memcpy(target.fK0cutDCA, fK0cutDCA, sizeof(Float_t) * 2);
+ memcpy(target.fK0cutVertexR, fK0cutVertexR, sizeof(Float_t) * 2);
+ memcpy(target.fK0cutInvMass, fK0cutInvMass, sizeof(Float_t) * 2);
+ // Lambda & anti-Lambda cut values
+ target.fLcutChi2NDF = fLcutChi2NDF;
+ memcpy(target.fLcutCosPoint, fLcutCosPoint, sizeof(Float_t) * 2);
+ memcpy(target.fLcutDCA, fLcutDCA, sizeof(Float_t) * 2);
+ memcpy(target.fLcutVertexR, fLcutVertexR, sizeof(Float_t) * 2);
+ memcpy(target.fLcutInvMass, fLcutInvMass, sizeof(Float_t) * 2);
+
+}
+//____________________________________________________________________
+Bool_t AliESDv0KineCuts::ProcessV0(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const
+{
+ //
+ // main user function
+ //
+
+ if(!v0) return kFALSE;
+ if(!fEvent){
+ AliErrorClass("No valid Event pointer available, provide it first");
+ return kFALSE;
+ }
+
+ if(!V0CutsCommon(v0)) return kFALSE;
+
+ const Int_t id = PreselectV0(v0);
+
+ if(!SingleTrackCuts(v0)) return kFALSE;
+
+ switch(id){
+ case kUndef:
+ return kFALSE;
+ case kGamma:
+ return CaseGamma(v0, pdgV0, pdgP, pdgN);
+ case kK0:
+ return CaseK0(v0, pdgV0, pdgP, pdgN);
+ case kLambda:
+ return CaseLambda(v0, pdgV0, pdgP, pdgN, 0);
+ case kALambda:
+ return CaseLambda(v0, pdgV0, pdgP, pdgN, 1);
+ default:
+ return kFALSE;
+ }
+
+ return kFALSE;
+}
+//____________________________________________________________________
+Bool_t AliESDv0KineCuts::ProcessV0(AliESDv0* const v0, Int_t &pdgP, Int_t &pdgN) const
+{
+ //
+ // main user function, simplified if the V0 identity is not necessary
+ //
+
+ if(!v0) return kFALSE;
+ if(!fEvent){
+ AliErrorClass("No valid Event pointer available, provide it first");
+ return kFALSE;
+ }
+
+ Int_t idV0 = -1;
+ return ProcessV0(v0, idV0, pdgP, pdgN);
+
+}
+//____________________________________________________________________
+Int_t AliESDv0KineCuts::PreselectV0(AliESDv0* const v0) const
+{
+ //
+ // Make a preselection (exclusive) of the V0 cadidates based on
+ // Armenteros plot
+ // the armenteros cut values are currently fixed and user is not able to set them via
+ // set funcions. The reason is that these cuts are optimized and furneter changes should
+ // not be necessary. To prove otherwise please study in detail before changing the values
+ //
+
+ Float_t ap[2] = {-1., -1.};
+ Armenteros(v0, ap);
+ // for clarity
+ const Float_t alpha = ap[0];
+ const Float_t qt = ap[1];
+
+ // selection cuts
+ // - the reagions for different candidates must not overlap
+
+ // Gamma cuts
+ const Double_t cutAlphaG = 0.35;
+ const Double_t cutQTG = 0.05;
+ const Double_t cutAlphaG2[2] = {0.6, 0.8};
+ const Double_t cutQTG2 = 0.04;
+
+ // K0 cuts
+ const Float_t cutQTK0[2] = {0.1075, 0.215};
+ const Float_t cutAPK0[2] = {0.199, 0.8}; // parameters for curved QT cut
+
+ // Lambda & A-Lambda cuts
+ const Float_t cutQTL = 0.03;
+ const Float_t cutAlphaL[2] = {0.35, 0.7};
+ const Float_t cutAlphaAL[2] = {-0.7, -0.35};
+ const Float_t cutAPL[3] = {0.107, -0.69, 0.5}; // parameters fir curved QT cut
+
+
+ if(kPurity == fMode){
+ // Check for Gamma candidates
+ if(qt < cutQTG){
+ if( (TMath::Abs(alpha) < cutAlphaG) ) return kGamma;
+ }
+ // additional region - should help high pT gammas
+ if(qt < cutQTG2){
+ if( (TMath::Abs(alpha) > cutAlphaG2[0]) && (TMath::Abs(alpha) < cutAlphaG2[1]) ) return kGamma;
+ }
+ }
+ if(kEffGamma == fMode){
+ if(qt < cutQTG) return kGamma;
+ }
+
+
+ // Check for K0 candidates
+ Float_t q = cutAPK0[0] * TMath::Sqrt(TMath::Abs(1 - alpha*alpha/(cutAPK0[1]*cutAPK0[1])));
+ if( (qt > cutQTK0[0]) && (qt < cutQTK0[1]) && (qt > q) ){
+ return kK0;
+ }
+
+ // Check for Lambda candidates
+ q = cutAPL[0] * TMath::Sqrt(TMath::Abs(1 - ( (alpha + cutAPL[1]) * (alpha + cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) ));
+ if( (alpha > cutAlphaL[0]) && (alpha < cutAlphaL[1]) && (qt > cutQTL) && (qt < q) ){
+ return kLambda;
+ }
+
+ // Check for A-Lambda candidates
+ q = cutAPL[0] * TMath::Sqrt(TMath::Abs(1 - ( (alpha - cutAPL[1]) * (alpha - cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) ));
+ if( (alpha > cutAlphaAL[0]) && (alpha < cutAlphaAL[1]) && (qt > cutQTL) && (qt < q) ){
+ return kALambda;
+ }
+
+ return kUndef;
+}
+//____________________________________________________________________
+Bool_t AliESDv0KineCuts::SingleTrackCuts(AliESDv0 * const v0) const
+{
+ //
+ // apply single track cuts
+ // correct sign not relevat here
+ //
+
+ if(!v0) return kFALSE;
+
+ Int_t pIndex = 0, nIndex = 0;
+ pIndex = v0->GetPindex();
+ nIndex = v0->GetNindex();
+ AliESDtrack* d[2];
+ d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
+ d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
+
+ for(Int_t i=0; i<2; ++i){
+ if(!d[i]) return kFALSE;
+
+ // status word
+ ULong_t status = d[i]->GetStatus();
+
+ // No. of TPC clusters leave to the users
+ if(d[i]->GetTPCNcls() < 1) return kFALSE;
+
+ // TPC refit
+ if(!(status & AliESDtrack::kTPCrefit)) return kFALSE;
+
+ // Chi2 per TPC cluster
+ Int_t nTPCclusters = d[i]->GetTPCNcls();
+ Float_t chi2perTPCcluster = d[i]->GetTPCchi2()/Float_t(nTPCclusters);
+ if(chi2perTPCcluster > 4) return kFALSE;
+
+ // TPC cluster ratio
+ Float_t cRatioTPC = d[i]->GetTPCNclsF() > 0. ? static_cast<Float_t>(d[i]->GetTPCNcls())/static_cast<Float_t> (d[i]->GetTPCNclsF()) : 1.;
+ if(cRatioTPC < 0.6) return kFALSE;
+
+ // kinks
+ if(d[i]->GetKinkIndex(0) != 0) return kFALSE;
+
+ }
+
+ return kTRUE;
+}
+//____________________________________________________________________
+Bool_t AliESDv0KineCuts::CaseGamma(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const
+{
+ //
+ // process the gamma conversion candidate
+ //
+
+ if(!v0) return kFALSE;
+
+ AliVTrack* daughter[2];
+ Int_t pIndex = 0, nIndex = 0;
+
+ Bool_t sign = CheckSigns(v0);
+ if(sign){
+ pIndex = v0->GetPindex();
+ nIndex = v0->GetNindex();
+ }
+ else{
+ pIndex = v0->GetNindex();
+ nIndex = v0->GetPindex();
+ }
+ daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
+ daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
+ if(!daughter[0] || !daughter[1]) return kFALSE;
+
+ AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kElectron), TMath::Abs(kElectron));
+ if(!kfMother) return kFALSE;
+
+ AliESDtrack* d[2];
+ d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
+ d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
+
+ Float_t iMass = v0->GetEffMass(0, 0);
+
+ // cos pointing angle
+ Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
+ cosPoint = TMath::ACos(cosPoint);
+
+ // DCA between daughters
+ Double_t dca = v0->GetDcaV0Daughters();
+
+ // Production vertex
+ Double_t x, y, z;
+ v0->GetXYZ(x,y,z);
+ Double_t r = TMath::Sqrt(x*x + y*y);
+
+ Double_t xy[2];
+ Double_t r2 = -1.;
+ if ( GetConvPosXY(d[0], d[1], xy) ){
+ r2 = TMath::Sqrt(xy[0]*xy[0] + xy[1]*xy[1]);
+ }
+
+ // psi pair
+ Double_t psiPair = PsiPair(v0);
+
+ // V0 chi2/ndf
+ Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();
+
+ if(kfMother) delete kfMother;
+
+ // apply the cuts
+
+ if(iMass > fGcutInvMass) return kFALSE;
+
+ if(chi2ndf > fGcutChi2NDF) return kFALSE;
+
+ if(cosPoint < fGcutCosPoint[0] || cosPoint > fGcutCosPoint[1]) return kFALSE;
+
+ if(dca < fGcutDCA[0] || dca > fGcutDCA[1]) return kFALSE;
+
+ if(r < fGcutVertexR[0] || r > fGcutVertexR[1]) return kFALSE;
+
+ if(psiPair < fGcutPsiPair[0] || psiPair > fGcutPsiPair[1]) return kFALSE;
+
+ // all cuts passed
+
+ pdgV0 = 22;
+ if(sign){
+ pdgP = -11;
+ pdgN = 11;
+ }
+ else{
+ pdgP = 11;
+ pdgN = -11;
+ }
+
+ return kTRUE;
+}
+//____________________________________________________________________
+Bool_t AliESDv0KineCuts::CaseK0(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN) const {
+ //
+ // process the K0 candidate
+ //
+
+ if(!v0) return kFALSE;
+
+ AliVTrack* daughter[2];
+ Int_t pIndex = 0, nIndex = 0;
+ Bool_t sign = CheckSigns(v0);
+ if(sign){
+ pIndex = v0->GetPindex();
+ nIndex = v0->GetNindex();
+ }
+ else{
+ pIndex = v0->GetNindex();
+ nIndex = v0->GetPindex();
+ }
+
+ daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
+ daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
+ if(!daughter[0] || !daughter[1]) return kFALSE;
+
+ AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kPiPlus));
+ if(!kfMother) return kFALSE;
+
+ AliESDtrack* d[2];
+ d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
+ d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
+
+ Float_t iMass = v0->GetEffMass(2, 2);
+
+ // cos pointing angle
+ Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
+ cosPoint = TMath::ACos(cosPoint);
+
+ // DCA between daughters
+ Double_t dca = v0->GetDcaV0Daughters();
+
+ // Production vertex
+ Double_t x, y, z;
+ v0->GetXYZ(x,y,z);
+
+ Double_t r = TMath::Sqrt(x*x + y*y);
+
+ // V0 chi2/ndf
+ Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();
+
+ if(kfMother) delete kfMother;
+
+ //
+ // apply the cuts
+ //
+ if(iMass < fK0cutInvMass[0] || iMass > fK0cutInvMass[1]) return kFALSE;
+
+ if(chi2ndf > fK0cutChi2NDF) return kFALSE;
+
+ if(cosPoint < fK0cutCosPoint[0] || cosPoint > fK0cutCosPoint[1]) return kFALSE;
+
+ if(dca < fK0cutDCA[0] || dca > fK0cutDCA[1]) return kFALSE;
+
+ if(r < fK0cutVertexR[0] || r > fK0cutVertexR[1]) return kFALSE;
+
+ // all cuts passed
+ pdgV0 = 310;
+ if(sign){
+ pdgP = 211;
+ pdgN = -211;
+ }
+ else{
+ pdgP = -211;
+ pdgN = 211;
+ }
+
+ return kTRUE;
+}
+//____________________________________________________________________
+Bool_t AliESDv0KineCuts::CaseLambda(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN, Int_t id) const {
+ //
+ // process teh Lambda and Anti-Lambda candidate
+ //
+
+ if(!v0) return kFALSE;
+
+ const Double_t cL0mass=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(); // PDG lambda mass
+
+ AliVTrack* daughter[2];
+ Int_t pIndex = 0, nIndex = 0;
+ Float_t mMass[2] = {-1., -1.};
+ Bool_t sign = CheckSigns(v0);
+ if(sign){
+ pIndex = v0->GetPindex();
+ nIndex = v0->GetNindex();
+ mMass[0] = v0->GetEffMass(4, 2);
+ mMass[1] = v0->GetEffMass(2, 4);
+ }
+ else{
+ pIndex = v0->GetNindex();
+ nIndex = v0->GetPindex();
+ mMass[0] = v0->GetEffMass(2, 4);
+ mMass[1] = v0->GetEffMass(4, 2);
+ }
+
+ daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
+ daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
+ if(!daughter[0] || !daughter[1]) return kFALSE;
+
+ AliKFParticle *kfMother[2] = {0x0, 0x0};
+ // Lambda
+ kfMother[0] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kProton), TMath::Abs(kPiPlus));
+ if(!kfMother[0]) return kFALSE;
+
+ // Anti-Lambda
+ kfMother[1] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kProton));
+ if(!kfMother[1]) return kFALSE;
+
+ Float_t dMass[2] = {TMath::Abs(mMass[0] - cL0mass), TMath::Abs(mMass[1] - cL0mass)};
+
+ AliESDtrack* d[2];
+ d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
+ d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
+ if(!d[0] || !d[1]) return kFALSE;
+
+ Float_t p[2] = {d[0]->GetP(), d[1]->GetP()};
+
+ // check the 3 lambda - antilambda variables
+ Int_t check[2] = {-1, -1}; // 0 : lambda, 1 : antilambda
+ // 1) momentum of the daughter particles - proton is expected to have higher momentum than pion
+ check[0] = (p[0] > p[1]) ? 0 : 1;
+ // 2) mass of the mother particle
+ check[1] = (dMass[0] < dMass[1]) ? 0 : 1;
+
+ // require positive correlation of (1) and (2)
+ if(check[0] != check[1]){
+ if(kfMother[0]) delete kfMother[0];
+ if(kfMother[1]) delete kfMother[1];
+ return kFALSE;
+ }
+
+ // now that the check[0] == check[1]
+ const Int_t type = check[0];
+
+ // require that the input armenteros preselection agree:
+ if(type != id) return kFALSE;
+
+ Float_t iMass =0.;
+ if(sign){
+ iMass = (type == 0) ? v0->GetEffMass(4, 2) : v0->GetEffMass(2, 4);
+ }
+ else{
+ iMass = (type == 0) ? v0->GetEffMass(2, 4) : v0->GetEffMass(4, 2);
+ }
+
+ // cos pointing angle
+ Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
+ cosPoint = TMath::ACos(cosPoint);
+
+ // DCA between daughters
+ Double_t dca = v0->GetDcaV0Daughters();
+
+ // Production vertex
+ Double_t x, y, z;
+ v0->GetXYZ(x,y,z);
+ Double_t r = TMath::Sqrt(x*x + y*y);
+
+ // proton - pion indices
+ Int_t ix[2] = {0, 1};
+ if(1 == type){
+ ix[0] = 1;
+ ix[1] = 0;
+ }
+
+ // V0 chi2/ndf
+ Double_t chi2ndf = kfMother[type]->GetChi2()/kfMother[type]->GetNDF();
+
+ if(kfMother[0]) delete kfMother[0];
+ if(kfMother[1]) delete kfMother[1];
+
+ //
+ // apply the cuts
+ //
+
+ if(iMass < fLcutInvMass[0] || iMass > fLcutInvMass[1]) return kFALSE;
+
+ if(chi2ndf > fLcutChi2NDF) return kFALSE;
+
+ if(cosPoint < fLcutCosPoint[0] || cosPoint > fLcutCosPoint[1]) return kFALSE;
+
+ if(dca < fLcutDCA[0] || dca > fLcutDCA[1]) return kFALSE;
+
+ if(r < fLcutVertexR[0] || r > fLcutVertexR[1]) return kFALSE;
+
+ // all cuts passed
+
+ if(0 == type){
+ pdgV0 = 3122;
+ if(sign){
+ pdgP = 2212;
+ pdgN = -211;
+ }
+ else{
+ pdgP = -211;
+ pdgN = 2212;
+ }
+ }
+ else{
+ pdgV0 = -3122;
+ if(sign){
+ pdgP = 211;
+ pdgN = -2212;
+ }
+ else{
+ pdgP = -2212;
+ pdgN = 211;
+ }
+ }
+
+ return kTRUE;
+}
+//____________________________________________________________________
+Bool_t AliESDv0KineCuts::V0CutsCommon(const AliESDv0 * const v0) const
+{
+ //
+ // V0 cuts common to all V0s
+ //
+
+ AliESDtrack* dN, *dP;
+
+ dP = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(v0->GetPindex()));
+ dN = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(v0->GetNindex()));
+
+ if(!dN || !dP) return kFALSE;
+
+ Int_t qP = dP->Charge();
+ Int_t qN = dN->Charge();
+
+ if((qP*qN) != -1) return kFALSE;
+
+ return kTRUE;
+}
+//____________________________________________________________________
+void AliESDv0KineCuts::Armenteros(AliESDv0* const v0, Float_t val[2]) const
+{
+ //
+ // computes the Armenteros variables for given V0
+ // fills the histogram
+ // returns the values via "val"
+ //
+
+ Double_t mn[3] = {0,0,0};
+ Double_t mp[3] = {0,0,0};
+ Double_t mm[3] = {0,0,0};
+
+ if(CheckSigns(v0)){
+ v0->GetNPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
+ v0->GetPPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
+ }
+ else{
+ v0->GetPPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
+ v0->GetNPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
+ }
+ v0->GetPxPyPz(mm[0],mm[1],mm[2]); //reconstructed cartesian momentum components of mother
+
+ TVector3 vecN(mn[0],mn[1],mn[2]);
+ TVector3 vecP(mp[0],mp[1],mp[2]);
+ TVector3 vecM(mm[0],mm[1],mm[2]);
+
+ Double_t thetaP = acos((vecP * vecM)/(vecP.Mag() * vecM.Mag()));
+ Double_t thetaN = acos((vecN * vecM)/(vecN.Mag() * vecM.Mag()));
+
+ Double_t alfa = ((vecP.Mag())*cos(thetaP)-(vecN.Mag())*cos(thetaN))/
+ ((vecP.Mag())*cos(thetaP)+(vecN.Mag())*cos(thetaN)) ;
+ Double_t qt = vecP.Mag()*sin(thetaP);
+
+ val[0] = alfa;
+ val[1] = qt;
+}
+//____________________________________________________________________
+Bool_t AliESDv0KineCuts::CheckSigns(AliESDv0* const v0) const
+{
+ //
+ // check wheter the sign was correctly applied to
+ // V0 daughter tracks
+ //
+
+ Bool_t correct = kFALSE;
+
+ Int_t pIndex = 0, nIndex = 0;
+ pIndex = v0->GetPindex();
+ nIndex = v0->GetNindex();
+
+ AliESDtrack* d[2];
+ d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
+ d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
+
+ Int_t sign[2];
+ sign[0] = (int)d[0]->GetSign();
+ sign[1] = (int)d[1]->GetSign();
+
+ if(-1 == sign[0] && 1 == sign[1]){
+ correct = kFALSE;
+ }
+ else{
+ correct = kTRUE;
+ }
+
+ return correct;
+}
+//________________________________________________________________
+Double_t AliESDv0KineCuts::PsiPair(AliESDv0* const v0) const
+{
+ //
+ // Angle between daughter momentum plane and plane
+ //
+
+ if(!fEvent) return -1.;
+
+ Float_t magField = fEvent->GetMagneticField();
+
+ Int_t pIndex = -1;
+ Int_t nIndex = -1;
+ if(CheckSigns(v0)){
+ pIndex = v0->GetPindex();
+ nIndex = v0->GetNindex();
+ }
+ else{
+ pIndex = v0->GetNindex();
+ nIndex = v0->GetPindex();
+ }
+
+
+ AliESDtrack* daughter[2];
+
+ daughter[0] = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(pIndex));
+ daughter[1] = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(nIndex));
+
+ Double_t x, y, z;
+ v0->GetXYZ(x,y,z);//Reconstructed coordinates of V0; to be replaced by Markus Rammler's method in case of conversions!
+
+ Double_t mn[3] = {0,0,0};
+ Double_t mp[3] = {0,0,0};
+
+
+ v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
+ v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
+
+
+ Double_t deltat = 1.;
+ deltat = TMath::ATan(mp[2]/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1])+1.e-13)) - TMath::ATan(mn[2]/(TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1])+1.e-13));//difference of angles of the two daughter tracks with z-axis
+
+ Double_t radiussum = TMath::Sqrt(x*x + y*y) + 50;//radius to which tracks shall be propagated
+
+ Double_t momPosProp[3];
+ Double_t momNegProp[3];
+
+ AliExternalTrackParam pt(*daughter[0]), nt(*daughter[1]);
+
+ Double_t psiPair = 4.;
+
+ if(nt.PropagateTo(radiussum,magField) == 0)//propagate tracks to the outside
+ psiPair = -5.;
+ if(pt.PropagateTo(radiussum,magField) == 0)
+ psiPair = -5.;
+ pt.GetPxPyPz(momPosProp);//Get momentum vectors of tracks after propagation
+ nt.GetPxPyPz(momNegProp);
+
+ Double_t pEle =
+ TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
+ Double_t pPos =
+ TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
+
+ Double_t scalarproduct =
+ momPosProp[0]*momNegProp[0]+momPosProp[1]*momNegProp[1]+momPosProp[2]*momNegProp[2];//scalar product of propagated positive and negative daughters' momenta
+
+ Double_t chipair = TMath::ACos(scalarproduct/(pEle*pPos));//Angle between propagated daughter tracks
+
+ psiPair = TMath::Abs(TMath::ASin(deltat/chipair));
+
+ return psiPair;
+}
+//___________________________________________________________________
+Bool_t AliESDv0KineCuts::GetConvPosXY(AliESDtrack * const ptrack, AliESDtrack * const ntrack, Double_t convpos[2]) const
+{
+ //
+ // recalculate the gamma conversion XY postition
+ //
+
+ const Double_t b = fEvent->GetMagneticField();
+
+ Double_t helixcenterpos[2];
+ GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos);
+
+ Double_t helixcenterneg[2];
+ GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg);
+
+ Double_t poshelix[6];
+ ptrack->GetHelixParameters(poshelix,b);
+ Double_t posradius = TMath::Abs(1./poshelix[4]);
+
+ Double_t neghelix[6];
+ ntrack->GetHelixParameters(neghelix,b);
+ Double_t negradius = TMath::Abs(1./neghelix[4]);
+
+ Double_t xpos = helixcenterpos[0];
+ Double_t ypos = helixcenterpos[1];
+ Double_t xneg = helixcenterneg[0];
+ Double_t yneg = helixcenterneg[1];
+
+ convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);
+ convpos[1] = (ypos*negradius+ yneg*posradius)/(negradius+posradius);
+
+ return 1;
+}
+//___________________________________________________________________
+Bool_t AliESDv0KineCuts::GetHelixCenter(AliESDtrack * const track, Double_t b,Int_t charge, Double_t center[2]) const
+{
+ //
+ // computes the center of the track helix
+ //
+
+ Double_t pi = TMath::Pi();
+
+ Double_t helix[6];
+ track->GetHelixParameters(helix,b);
+
+ Double_t xpos = helix[5];
+ Double_t ypos = helix[0];
+ Double_t radius = TMath::Abs(1./helix[4]);
+ Double_t phi = helix[2];
+
+ if(phi < 0){
+ phi = phi + 2*pi;
+ }
+
+ phi -= pi/2.;
+ Double_t xpoint = radius * TMath::Cos(phi);
+ Double_t ypoint = radius * TMath::Sin(phi);
+
+ if(b<0){
+ if(charge > 0){
+ xpoint = - xpoint;
+ ypoint = - ypoint;
+ }
+ /* avoid self assignment
+ if(charge < 0){
+ xpoint = xpoint;
+ ypoint = ypoint;
+ }
+ */
+ }
+ if(b>0){
+ /* avoid self assignment
+ if(charge > 0){
+ xpoint = xpoint;
+ ypoint = ypoint;
+ }
+ */
+ if(charge < 0){
+ xpoint = - xpoint;
+ ypoint = - ypoint;
+ }
+ }
+ center[0] = xpos + xpoint;
+ center[1] = ypos + ypoint;
+
+ return 1;
+}
+//___________________________________________________________________
+AliKFParticle *AliESDv0KineCuts::CreateMotherParticle(const AliVTrack* const pdaughter, const AliVTrack* const ndaughter, Int_t pspec, Int_t nspec) const
+{
+ //
+ // Creates a mother particle
+ //
+ AliKFParticle pkfdaughter(*pdaughter, pspec);
+ AliKFParticle nkfdaughter(*ndaughter, nspec);
+
+
+ // Create the mother particle
+ AliKFParticle *m = new AliKFParticle(pkfdaughter, nkfdaughter);
+ m->SetField(fEvent->GetMagneticField());
+ if(TMath::Abs(kElectron) == pspec && TMath::Abs(kElectron) == nspec) m->SetMassConstraint(0, 0.001);
+ else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(), 0.);
+ else if(TMath::Abs(kProton) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);
+ else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kProton) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);
+ else{
+ AliErrorClass("Wrong daughter ID - mass constraint can not be set");
+ }
+
+ AliKFVertex improvedVertex = *fPrimaryVertex;
+ improvedVertex += *m;
+ m->SetProductionVertex(improvedVertex);
+
+ // update 15/06/2010
+ // mother particle will not be added to primary vertex but only to its copy
+ // as this confilcts with calling
+ // m->SetPrimaryVertex() function and
+ // subsequently removing the mother particle afterwards
+ // Source: Sergey Gorbunov
+
+ return m;
+}
+//____________________________________________________________________
+void AliESDv0KineCuts::SetEvent(AliESDEvent* const event){
+ //
+ // direct setter of ESD event
+ //
+ fEvent = event;
+ if(!fEvent){
+ AliErrorClass("Invalid input event pointer");
+ return;
+ }
+if (fUseExternalVertex) return;
+else{
+ if(fPrimaryVertex && fDeleteVertex){
+ delete fPrimaryVertex;
+ fPrimaryVertex=0x0;
+ }
+ fPrimaryVertex = new AliKFVertex(*(fEvent->GetPrimaryVertex()));
+ fDeleteVertex=kTRUE;
+ }
+
+
+
+}
+//____________________________________________________________________
+void AliESDv0KineCuts::SetEvent(AliVEvent* const event){
+ //
+ // direct setter of ESD event
+ //
+
+ fEvent = dynamic_cast<AliESDEvent*>(event);
+ if(!fEvent){
+ AliErrorClass("Invalid input event pointer");
+ return;
+ }
+
+ if (fUseExternalVertex) return;
+ else{
+ if(fPrimaryVertex && fDeleteVertex){
+ delete fPrimaryVertex;
+ fPrimaryVertex=0x0;
+ }
+ fPrimaryVertex = new AliKFVertex(*(fEvent->GetPrimaryVertex()));
+ fDeleteVertex=kTRUE;
+ }
+}
+
+
+//________________________________________________________________
+void AliESDv0KineCuts::UseExternalVertex(Bool_t use_external){
+ //
+ // Reenable primary Vertex from ESD event
+ //
+ if (use_external) fUseExternalVertex =kTRUE;
+ else fUseExternalVertex =kFALSE;
+}
+
+
+
+
+//________________________________________________________________
+void AliESDv0KineCuts::SetPrimaryVertex(AliKFVertex* const v){
+ //
+ // set the primary vertex of the event
+ //
+ if(fPrimaryVertex && fDeleteVertex){
+ delete fPrimaryVertex;
+ fPrimaryVertex =0x0;
+ fDeleteVertex = kFALSE;
+ }
+ fUseExternalVertex=kTRUE;
+ fPrimaryVertex = v; // set primary Vertex
+ if(!fPrimaryVertex){
+ AliErrorClass("Failed to initialize the primary vertex");
+ return;
+ }
+}
+//___________________________________________________________________
+void AliESDv0KineCuts::SetMode(Int_t mode, Int_t type){
+ //
+ // this function allows the user to select (prior running the 'ProcessV0' function)
+ // to select different approaches to V0 selection - the 'mode'
+ // - and -
+ // different systems (pp, PbPb) - 'type'
+ //
+ // To see the cut values for different modes please refer to the
+ // function SetCuts()
+ //
+ // Important notice: based on the parameters particular sets of cuts will
+ // be activated for teh V0 selection. If some additional changes to single
+ // cuts are needed please us the SetXXXcut function (see the header file)
+ //
+
+ switch(mode){
+ case kPurity:
+ fMode = kPurity; // used to obtain highest purity possible - the efficiency may be low
+ break;
+ case kEffGamma:
+ fMode = kEffGamma; // used to obtain highes efficiency possible - the purity may be worse
+ break;
+ default:
+ AliError("V0 selection mode not recognozed, setting 'kPurity'");
+ fMode = kPurity;
+ }
+
+ switch(type){
+ case kPP:
+ fType = kPP; // cuts optimized for low multiplicity
+ break;
+ case kPbPb:
+ fType = kPbPb; // cuts optimized for high multiplicity
+ break;
+ }
+
+ // setup the cut values for selected mode & type
+ SetCuts();
+
+}
+//___________________________________________________________________
+void AliESDv0KineCuts::SetMode(Int_t mode, const char* type){
+ //
+ // overloaded function - please see above
+ //
+
+ Int_t t = -1;
+
+ if(!strcmp("pp", type)) t = kPP;
+ else if(!(strcmp("PbPb", type))) t = kPbPb;
+ else{
+ AliError("data type not recognized, setting 'pp'");
+ t = kPP;
+ }
+
+ SetMode(mode, t);
+
+}
+//___________________________________________________________________
+void AliESDv0KineCuts::SetCuts(){
+ //
+ // this funciton sets the default cut values based on the selected
+ // fMode and fType.
+ // please note that only the cuts that have different values than the default
+ // cuts are updated here
+ //
+
+ // last update: 14/02/2011
+ // as a very preliminary - the only change to default cuts is to apply
+ // less restricting gamma conversion selection in PreselectV0() function
+
+
+
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff