#include <TString.h>
#include <TMath.h>
+#include "AliAODEvent.h"
#include "AliAODTrack.h"
#include "AliAODPid.h"
+#include "AliAODVertex.h"
+#include "AliAODTrack.h"
#include "AliESDEvent.h"
#include "AliESDVertex.h"
#include "AliESDtrack.h"
#include "AliVParticle.h"
#include "AliVertexerTracks.h"
#include "AliVVertex.h"
+#include "AliExternalTrackParam.h"
#include "AliHFEextraCuts.h"
ClassImp(AliHFEextraCuts)
-const Int_t AliHFEextraCuts::fgkNQAhistos = 8;
+const Int_t AliHFEextraCuts::fgkNQAhistos = 10;
//______________________________________________________
AliHFEextraCuts::AliHFEextraCuts(const Char_t *name, const Char_t *title):
fMinNClustersTPCPID(0),
fClusterRatioTPC(0.),
fMinTrackletsTRD(0),
+ fMaxChi2TRD(5.0),
+ fMinNbITScls(0),
fTRDtrackletsExact(0),
fPixelITS(0),
+ fDriftITS(0),
fTPCclusterDef(0),
fTPCclusterRatioDef(0),
fFractionTPCShared(-1.0),
+ fAbsHFEImpactParamNsigmaR(kTRUE),
+ fTOFsignalDx(1.0),
+ fTOFsignalDz(1.0),
fCheck(kFALSE),
fQAlist(0x0) ,
fDebugLevel(0)
//
// Default Constructor
//
+ //printf("Set the number of min ITS clusters %d\n",(Int_t)fMinNbITScls);
memset(fImpactParamCut, 0, sizeof(Float_t) * 4);
memset(fIPcutParam, 0, sizeof(Float_t) * 4);
}
fMinNClustersTPCPID(c.fMinNClustersTPCPID),
fClusterRatioTPC(c.fClusterRatioTPC),
fMinTrackletsTRD(c.fMinTrackletsTRD),
+ fMaxChi2TRD(c.fMaxChi2TRD),
+ fMinNbITScls(c.fMinNbITScls),
fTRDtrackletsExact(c.fTRDtrackletsExact),
fPixelITS(c.fPixelITS),
+ fDriftITS(c.fDriftITS),
fTPCclusterDef(c.fTPCclusterDef),
fTPCclusterRatioDef(c.fTPCclusterRatioDef),
fFractionTPCShared(c.fFractionTPCShared),
+ fAbsHFEImpactParamNsigmaR(c.fAbsHFEImpactParamNsigmaR),
+ fTOFsignalDx(c.fTOFsignalDx),
+ fTOFsignalDz(c.fTOFsignalDz),
fCheck(c.fCheck),
fQAlist(0x0),
fDebugLevel(0)
fMinNClustersTPC = c.fMinNClustersTPC;
fMinNClustersTPCPID = c.fMinNClustersTPCPID;
fMinTrackletsTRD = c.fMinTrackletsTRD;
+ fMaxChi2TRD = c.fMaxChi2TRD;
+ fMinNbITScls = c.fMinNbITScls;
fTRDtrackletsExact = c.fTRDtrackletsExact;
fPixelITS = c.fPixelITS;
+ fDriftITS = c.fDriftITS;
fTPCclusterDef = c.fTPCclusterDef;
fTPCclusterRatioDef = c.fTPCclusterRatioDef;
fFractionTPCShared = c.fFractionTPCShared;
+ fAbsHFEImpactParamNsigmaR = c.fAbsHFEImpactParamNsigmaR;
+ fTOFsignalDx = c.fTOFsignalDx;
+ fTOFsignalDz = c.fTOFsignalDz;
fCheck = c.fCheck;
fDebugLevel = c.fDebugLevel;
memcpy(fImpactParamCut, c.fImpactParamCut, sizeof(Float_t) * 4);
//
// Steering function for the track selection
//
- TString type = o->IsA()->GetName();
- AliDebug(2, Form("Object type %s", type.Data()));
- if(!type.CompareTo("AliESDtrack") || !type.CompareTo("AliAODTrack")){
+ TClass *type = o->IsA();
+ AliDebug(2, Form("Object type %s", type->GetName()));
+ if((type == AliESDtrack::Class()) || (type == AliAODTrack::Class())){
return CheckRecCuts(dynamic_cast<AliVTrack *>(o));
}
return CheckMCCuts(dynamic_cast<AliVParticle *>(o));
// QA histograms are filled before track selection and for
// selected tracks after track selection
//
- AliDebug(1, "Called");
+ AliDebug(1, Form("%s: Called", GetName()));
ULong64_t survivedCut = 0; // Bitmap for cuts which are passed by the track, later to be compared with fRequirements
if(IsQAOn()) FillQAhistosRec(track, kBeforeCuts);
// Apply cuts
- Float_t impactR, impactZ;
+ Float_t impactR = -999.;
+ Float_t impactZ = -999.;
Double_t hfeimpactR, hfeimpactnsigmaR;
Double_t hfeimpactRcut, hfeimpactnsigmaRcut;
Double_t maximpactRcut;
GetHFEImpactParameterCuts(track, hfeimpactRcut, hfeimpactnsigmaRcut);
GetHFEImpactParameters(track, hfeimpactR, hfeimpactnsigmaR);
}
+ Int_t nclsITS = GetITSNbOfcls(track);
UInt_t nclsTPC = GetTPCncls(track);
- UInt_t nclsTPCPID = GetTPCnclusdEdx(track);
+ UInt_t nclsTPCPID = track->GetTPCsignalN();
// printf("Check TPC findable clusters: %d, found Clusters: %d\n", track->GetTPCNclsF(), track->GetTPCNcls());
Float_t fractionSharedClustersTPC = GetTPCsharedClustersRatio(track);
Double_t ratioTPC = GetTPCclusterRatio(track);
UChar_t trdTracklets;
trdTracklets = track->GetTRDntrackletsPID();
+ Float_t trdchi2=-999.;
+ trdchi2=GetTRDchi(track);
UChar_t itsPixel = track->GetITSClusterMap();
Int_t status1 = GetITSstatus(track, 0);
Int_t status2 = GetITSstatus(track, 1);
Bool_t statusL0 = CheckITSstatus(status1);
Bool_t statusL1 = CheckITSstatus(status2);
+ Double_t tofsignalDx = 0.0;
+ Double_t tofsignalDz = 0.0;
+ GetTOFsignalDxDz(track,tofsignalDx,tofsignalDz);
+
if(TESTBIT(fRequirements, kTPCfractionShared)) {
// cut on max fraction of shared TPC clusters
if(TMath::Abs(fractionSharedClustersTPC) < fFractionTPCShared) SETBIT(survivedCut, kTPCfractionShared);
}
if(TESTBIT(fRequirements, kMinHFEImpactParamNsigmaR)){
// cut on max. HFE Impact Parameter n sigma in Radial direction
- if(TMath::Abs(hfeimpactnsigmaR) >= hfeimpactnsigmaRcut) SETBIT(survivedCut, kMinHFEImpactParamNsigmaR);
+ if(fAbsHFEImpactParamNsigmaR) {
+// if((TMath::Abs(hfeimpactnsigmaR) >= hfeimpactnsigmaRcut) && (TMath::Abs(hfeimpactnsigmaR) < 8)) { //mj debug
+ if(TMath::Abs(hfeimpactnsigmaR) >= hfeimpactnsigmaRcut) {
+ SETBIT(survivedCut, kMinHFEImpactParamNsigmaR);
+// printf("0: hfeimpactnsigmaR= %lf hfeimpactnsigmaRcut= %lf\n",hfeimpactnsigmaR,hfeimpactnsigmaRcut);
+ }
+ }
+ else {
+ if(hfeimpactnsigmaRcut>0){
+ if(hfeimpactnsigmaR >= hfeimpactnsigmaRcut) {
+ SETBIT(survivedCut, kMinHFEImpactParamNsigmaR);
+ //printf("1: hfeimpactnsigmaR= %lf hfeimpactnsigmaRcut= %lf\n",hfeimpactnsigmaR,hfeimpactnsigmaRcut);
+ }
+ }
+ else{
+ if(hfeimpactnsigmaR <= hfeimpactnsigmaRcut) {
+ SETBIT(survivedCut, kMinHFEImpactParamNsigmaR);
+ //printf("2: hfeimpactnsigmaR= %lf hfeimpactnsigmaRcut= %lf\n",hfeimpactnsigmaR,hfeimpactnsigmaRcut);
+ }
+ }
+ }
}
if(TESTBIT(fRequirements, kClusterRatioTPC)){
// cut on min ratio of found TPC clusters vs findable TPC clusters
}
if(TESTBIT(fRequirements, kMinTrackletsTRD)){
// cut on minimum number of TRD tracklets
- AliDebug(1, Form("Min TRD cut: [%d|%d], Require exact number [%s]\n", fMinTrackletsTRD, trdTracklets, fTRDtrackletsExact ? "Yes" : "No"));
+ AliDebug(1, Form("%s: Min TRD cut: [%d|%d], Require exact number [%s]\n", GetName(), fMinTrackletsTRD, trdTracklets, fTRDtrackletsExact ? "Yes" : "No"));
if(fTRDtrackletsExact){
if(trdTracklets == fMinTrackletsTRD) {
SETBIT(survivedCut, kMinTrackletsTRD);
- AliDebug(1, "Track Selected");
+ AliDebug(1, Form("%s: Track Selected", GetName()));
}
}else{
if(trdTracklets >= fMinTrackletsTRD){
SETBIT(survivedCut, kMinTrackletsTRD);
- AliDebug(1, "Track Selected");
+ AliDebug(1, Form("%s: Track Selected", GetName()));
}
//printf("Min number of tracklets %d\n",fMinTrackletsTRD);
}
}
+
+ if(TESTBIT(fRequirements, kMaxTRDChi2)){
+ // cut on TRD chi2
+ AliDebug(1, Form("%s: Cut on TRD chi2: [%f|%f]\n", GetName(),fMaxChi2TRD, trdchi2));
+ if(trdchi2 < fMaxChi2TRD) {
+ SETBIT(survivedCut, kMaxTRDChi2);
+ AliDebug(1,Form("%s: survived %f\n",GetName(),trdchi2));
+ }
+ }
+
+ if(TESTBIT(fRequirements, kMinNbITScls)){
+ // cut on minimum number of ITS clusters
+ //printf(Form("Min ITS clusters: [%d|%d]\n", (Int_t)fMinNbITScls, nclsITS));
+ AliDebug(1, Form("%s: Min ITS clusters: [%d|%d]\n", GetName(), fMinNbITScls, nclsITS));
+ if(nclsITS >= ((Int_t)fMinNbITScls)) SETBIT(survivedCut, kMinNbITScls);
+ }
+
if(TESTBIT(fRequirements, kMinNClustersTPC)){
- // cut on minimum number of TRD tracklets
- AliDebug(1, Form("Min TPC cut: [%d|%d]\n", fMinNClustersTPC, nclsTPC));
+ // cut on minimum number of TPC tracklets
+ //printf(Form("Min TPC cut: [%d|%d]\n", fMinNClustersTPC, nclsTPC));
+ AliDebug(1, Form("%s: Min TPC cut: [%d|%d]\n", GetName(), fMinNClustersTPC, nclsTPC));
if(nclsTPC >= fMinNClustersTPC) SETBIT(survivedCut, kMinNClustersTPC);
}
if(TESTBIT(fRequirements, kMinNClustersTPCPID)){
- AliDebug(1, Form("Min TPC PID cut: [%d|%d]\n", fMinNClustersTPCPID, nclsTPCPID));
+ AliDebug(1, Form("%s: Min TPC PID cut: [%d|%d]\n", GetName(), fMinNClustersTPCPID, nclsTPCPID));
if(nclsTPCPID >= fMinNClustersTPCPID) SETBIT(survivedCut, kMinNClustersTPCPID);
}
+ if(TESTBIT(fRequirements, kDriftITS)){
+ //printf("Require drift\n");
+ switch(fDriftITS){
+ case kFirstD:
+ if(TESTBIT(itsPixel, 2)) SETBIT(survivedCut, kDriftITS);
+ break;
+ default:
+ SETBIT(survivedCut, kDriftITS);
+ break;
+ }
+ }
if(TESTBIT(fRequirements, kPixelITS)){
// cut on ITS pixel layers
- AliDebug(1, "ITS cluster Map: ");
+ AliDebug(1, Form("%s: ITS cluster Map: ", GetName()));
//PrintBitMap(itsPixel);
switch(fPixelITS){
case kFirst:
AliDebug(2, "First");
+ //printf("First\n");
if(TESTBIT(itsPixel, 0))
SETBIT(survivedCut, kPixelITS);
else
SETBIT(survivedCut, kPixelITS);
break;
case kSecond:
+ //printf("Second\n");
AliDebug(2, "Second");
if(TESTBIT(itsPixel, 1))
SETBIT(survivedCut, kPixelITS);
SETBIT(survivedCut, kPixelITS);
break;
case kBoth:
+ //printf("Both\n");
AliDebug(2, "Both");
if(TESTBIT(itsPixel,0) && TESTBIT(itsPixel,1))
SETBIT(survivedCut, kPixelITS);
SETBIT(survivedCut, kPixelITS);
break;
case kAny:
+ //printf("Any\n");
AliDebug(2, "Any");
if(TESTBIT(itsPixel, 0) || TESTBIT(itsPixel, 1))
SETBIT(survivedCut, kPixelITS);
SETBIT(survivedCut, kPixelITS);
}
break;
+ case kNone:
+ // No cut applied, set as survived
+ SETBIT(survivedCut, kPixelITS);
+ break;
default:
+ // default, defined as no cut applied
AliDebug(2, "None");
+ SETBIT(survivedCut, kPixelITS);
break;
}
- AliDebug(1, Form("Survived Cut: %s\n", TESTBIT(survivedCut, kPixelITS) ? "YES" : "NO"));
+ AliDebug(1, Form("%s: Survived Cut: %s\n", GetName(), TESTBIT(survivedCut, kPixelITS) ? "YES" : "NO"));
}
if(TESTBIT(fRequirements, kTOFPID)){
if(TESTBIT(fRequirements, kTPCPIDCleanUp)){
// cut on TPC PID cleanup
- Int_t fClusterdEdx=GetTPCnclusdEdx(track);
Bool_t fBitsAboveThreshold=GetTPCCountSharedMapBitsAboveThreshold(track);
- if((fBitsAboveThreshold==0)&&(fClusterdEdx>80)) SETBIT(survivedCut, kTPCPIDCleanUp);
+ if((fBitsAboveThreshold==0)&&(nclsTPCPID>80)) SETBIT(survivedCut, kTPCPIDCleanUp);
}
if(TESTBIT(fRequirements, kEMCALmatch)){
if(track->GetStatus() && AliESDtrack::kEMCALmatch) SETBIT(survivedCut, kEMCALmatch);
}
+ if(TESTBIT(fRequirements, kRejectKinkDaughter)){
+ //printf("test daughter\n");
+ if(!IsKinkDaughter(track)) SETBIT(survivedCut, kRejectKinkDaughter);
+ //else printf("Found kink daughter\n");
+ }
+
+ if(TESTBIT(fRequirements, kRejectKinkMother)){
+ //printf("test mother\n");
+ if(!IsKinkMother(track)) SETBIT(survivedCut, kRejectKinkMother);
+ //else printf("Found kink mother\n");
+ }
+ if(TESTBIT(fRequirements, kTOFsignalDxy)){
+ // cut on TOF matching cluster
+ if((TMath::Abs(tofsignalDx) <= fTOFsignalDx) && (TMath::Abs(tofsignalDz) <= fTOFsignalDz)) SETBIT(survivedCut, kTOFsignalDxy);
+ }
+ if(TESTBIT(fRequirements, kITSpattern)){
+ // cut on ITS pattern (every layer with a working ITS module must have an ITS cluster)
+ if(CheckITSpattern(track)) SETBIT(survivedCut, kITSpattern);
+ }
+
if(fRequirements == survivedCut){
- //
+ //
// Track selected
//
- AliDebug(2, "Track Survived cuts\n");
+ AliDebug(2, Form("%s: Track Survived cuts\n", GetName()));
if(IsQAOn()) FillQAhistosRec(track, kAfterCuts);
return kTRUE;
}
- AliDebug(2, "Track cut");
+ AliDebug(2, Form("%s: Track cut", GetName()));
if(IsQAOn()) FillCutCorrelation(survivedCut);
return kFALSE;
}
if(track->GetStatus() && AliESDtrack::kEMCALmatch) hStatusBits->Fill(3);
if(GetTPCCountSharedMapBitsAboveThreshold(track)==0) hStatusBits->Fill(4);
}
- if((htmp = dynamic_cast<TH1F *>(fQAlist->At(7 + when * fgkNQAhistos)))) htmp->Fill(GetTPCnclusdEdx(track));
+ if((htmp = dynamic_cast<TH1F *>(fQAlist->At(7 + when * fgkNQAhistos)))) htmp->Fill(track->GetTPCsignalN());
+
+ if((htmp = dynamic_cast<TH1F *>(fQAlist->At(8 + when * fgkNQAhistos)))) htmp->Fill(GetTRDchi(track));
+
+ if((htmp = dynamic_cast<TH1F *>(fQAlist->At(9 + when * fgkNQAhistos)))) htmp->Fill(GetITSNbOfcls(track));
}
// //______________________________________________________
fQAlist->AddAt(histo1D, 1 + icond * fgkNQAhistos);
histo1D->GetXaxis()->SetTitle("Impact Parameter");
histo1D->GetYaxis()->SetTitle("Number of Tracks");
- qaList->AddAt((histo1D = new TH1F(Form("%s_tpcClr%s",GetName(),cutstr[icond].Data()), "Cluster Ratio TPC", 100, 0, 1)), 2 + icond * fgkNQAhistos);
+ qaList->AddAt((histo1D = new TH1F(Form("%s_tpcClr%s",GetName(),cutstr[icond].Data()), "Cluster Ratio TPC", 100, 0, 1.3)), 2 + icond * fgkNQAhistos);
fQAlist->AddAt(histo1D, 2 + icond * fgkNQAhistos);
histo1D->GetXaxis()->SetTitle("Cluster Ratio TPC");
histo1D->GetYaxis()->SetTitle("Number of Tracks");
fQAlist->AddAt(histo1D, 7 + icond * fgkNQAhistos);
histo1D->GetXaxis()->SetTitle("Number of TPC clusters for dEdx calculation");
histo1D->GetYaxis()->SetTitle("counts");
+ qaList->AddAt((histo1D = new TH1F(Form("%s_trdchi2perTracklet%s",GetName(),cutstr[icond].Data()), "chi2 per TRD tracklet", 100, 0, 10)), 8 + icond * fgkNQAhistos);
+ fQAlist->AddAt(histo1D, 8 + icond * fgkNQAhistos);
+ histo1D->GetXaxis()->SetTitle("Chi2 per TRD Tracklet");
+ histo1D->GetYaxis()->SetTitle("Number of Tracks");
+ qaList->AddAt((histo1D = new TH1F(Form("%s_ITScls%s",GetName(),cutstr[icond].Data()), "ITScls", 6, 0, 6)), 9 + icond * fgkNQAhistos);
+ fQAlist->AddAt(histo1D, 9 + icond * fgkNQAhistos);
+ histo1D->GetXaxis()->SetTitle("ITScls");
+ histo1D->GetYaxis()->SetTitle("Number of ITS cls");
+
}
// Add cut correlation
qaList->AddAt((histo2D = new TH2F(Form("%s_cutcorrelation",GetName()), "Cut Correlation", kNcuts, 0, kNcuts - 1, kNcuts, 0, kNcuts -1)), 2 * fgkNQAhistos);
}
//______________________________________________________
-Int_t AliHFEextraCuts::GetITSstatus(AliVTrack *track, Int_t layer){
+Int_t AliHFEextraCuts::GetITSstatus(const AliVTrack * const track, Int_t layer) const {
//
// Check ITS layer status
//
if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
Int_t det;
Float_t xloc, zloc;
- AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+ const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
if(esdtrack) esdtrack->GetITSModuleIndexInfo(layer, det, status, xloc, zloc);
}
return status;
// Checks if number of shared bits is above 1
//
Int_t nsharebit = 1;
- TString type = track->IsA()->GetName();
- if(!type.CompareTo("AliESDtrack")){
- AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- if(esdtrack){ // coverity
- TBits shared = esdtrack->GetTPCSharedMap();
- if(shared.CountBits() >= 2) nsharebit=1;
- else nsharebit=0;
- }
- }
+ const TBits *shared;
+ if(track->IsA() == AliESDtrack::Class())
+ shared = &((static_cast<AliESDtrack *>(track))->GetTPCSharedMap());
+ else if(track->IsA() == AliAODTrack::Class())
+ shared = &((static_cast<AliAODTrack *>(track))->GetTPCSharedMap());
+ else
+ return 0;
+
+ if(shared->CountBits() >= 2) nsharebit=1;
+ else nsharebit=0;
+
return nsharebit;
}
-
-
//______________________________________________________
UInt_t AliHFEextraCuts::GetTPCncls(AliVTrack *track){
- //
- // Get Number of findable clusters in the TPC
- //
- Int_t nClusters = 0; // in case no Information available consider all clusters findable
- TString type = track->IsA()->GetName();
- if(!type.CompareTo("AliESDtrack")){
- AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- if(esdtrack){ // coverity
- if(TESTBIT(fTPCclusterDef, kFoundIter1)){
- nClusters = esdtrack->GetTPCNclsIter1();
- AliDebug(2, ("Using def kFoundIter1"));
- } else if(TESTBIT(fTPCclusterDef, kCrossedRows)){
- AliDebug(2, ("Using def kCrossedRows"));
- nClusters = static_cast<UInt_t>(esdtrack->GetTPCClusterInfo(2,1));
- } else{
- AliDebug(2, ("Using def kFound"));
- nClusters = esdtrack->GetTPCNcls();
- }
- }
- }
- else if(!type.CompareTo("AliAODTrack")){
- AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
- if(aodtrack){
- const TBits &tpcmap = aodtrack->GetTPCClusterMap();
- for(UInt_t ibit = 0; ibit < tpcmap.GetNbits(); ibit++)
- if(tpcmap.TestBitNumber(ibit)) nClusters++;
- }
- }
- return nClusters;
-}
-
-//______________________________________________________
-UInt_t AliHFEextraCuts::GetTPCnclusdEdx(AliVTrack *track){
//
- // Get number of TPC cluster used to calculate dEdx
+ // Get Number of findable clusters in the TPC
//
- Int_t nClustersdEdx = 0;
- TString type = track->IsA()->GetName();
- if(!type.CompareTo("AliESDtrack")){
- AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- if(esdtrack){ // coverity
- nClustersdEdx = esdtrack->GetTPCsignalN();
+ Int_t nClusters = 0; // in case no Information available consider all clusters findable
+ TClass *type = track->IsA();
+ if(TESTBIT(fTPCclusterDef, kFoundIter1)){
+ if(type == AliESDtrack::Class()){
+ AliDebug(2, ("Using def kFoundIter1"));
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ nClusters = esdtrack->GetTPCNclsIter1();
+ } else {
+ AliDebug(2, ("Number of clusters in iteration 1 not available for AOD tracks"));
}
+ } else if(TESTBIT(fTPCclusterDef, kCrossedRows)){
+ AliDebug(2, ("Using def kCrossedRows"));
+ nClusters = static_cast<UInt_t>(track->GetTPCClusterInfo(2,1));
+ } else if(TESTBIT(fTPCclusterDef, kFound)){
+ AliDebug(2, ("Using def kFound"));
+ nClusters = track->GetTPCNcls();
}
- return nClustersdEdx;
+ else {
+ AliDebug(2, ("Using def kFoundAll"));
+ if(type == AliESDtrack::Class()){
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ const TBits &clusterTPC = esdtrack->GetTPCClusterMap();
+ nClusters = clusterTPC.CountBits();
+ }
+ else {
+ AliAODTrack *aodtrack = static_cast<AliAODTrack *>(track);
+ const TBits &clusterTPC = aodtrack->GetTPCClusterMap();
+ nClusters = clusterTPC.CountBits();
+ }
+ }
+ return nClusters;
}
-
-
//______________________________________________________
Float_t AliHFEextraCuts::GetTPCsharedClustersRatio(AliVTrack *track){
//
// Get fraction of shared TPC clusters
//
Float_t fracClustersTPCShared = 0.0;
- TString type = track->IsA()->GetName();
- if(!type.CompareTo("AliESDtrack")){
- AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- if(esdtrack){ // coverity
- Float_t nClustersTPC = esdtrack->GetTPCclusters(0);
- Int_t nClustersTPCShared = esdtrack->GetTPCnclsS();
- if (nClustersTPC!=0) {
- fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClustersTPC);
- }
- }
+ Int_t nClustersTPC = track->GetTPCNcls();
+ Int_t nClustersTPCShared = 0;
+ TClass *type = track->IsA();
+ if(type == AliESDtrack::Class()){
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ nClustersTPCShared = esdtrack->GetTPCnclsS();
+ } else if(type == AliAODTrack::Class()){
+ AliAODTrack *aodtrack = static_cast<AliAODTrack *>(track);
+ const TBits &shared = aodtrack->GetTPCSharedMap();
+ nClustersTPCShared = shared.CountBits(0) - shared.CountBits(159);
+ }
+ if (nClustersTPC!=0) {
+ fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClustersTPC);
}
return fracClustersTPCShared;
}
// Only implemented for ESD tracks
//
Double_t clusterRatio = 1.; // in case no Information available consider all clusters findable
- TString type = track->IsA()->GetName();
- if(!type.CompareTo("AliESDtrack")){
- AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- if(esdtrack){ // coverity
- if(TESTBIT(fTPCclusterRatioDef, kCROverFindable)){
- AliDebug(2, "Using ratio def kCROverFindable");
- clusterRatio = esdtrack->GetTPCNclsF() ? esdtrack->GetTPCClusterInfo(2,1)/static_cast<Double_t>(esdtrack->GetTPCNclsF()) : 1.; // crossed rows/findable
- } else if(TESTBIT(fTPCclusterRatioDef, kFoundOverCR)){
- AliDebug(2, "Using ratio def kFoundOverCR");
- clusterRatio = esdtrack->GetTPCClusterInfo(2,0); // found/crossed rows
- } else if(TESTBIT(fTPCclusterRatioDef, kFoundOverFindableIter1)){
- AliDebug(2, "Using ratio def kFoundOverFindableIter1");
- clusterRatio = esdtrack->GetTPCNclsFIter1() ? static_cast<Double_t>(esdtrack->GetTPCNclsIter1())/static_cast<Double_t>(esdtrack->GetTPCNclsFIter1()) : 1.; // crossed
- } else {
- AliDebug(2, "Using ratio def kFoundOverFindable");
- clusterRatio = esdtrack->GetTPCNclsF() ? static_cast<Double_t>(esdtrack->GetTPCNcls())/static_cast<Double_t>(esdtrack->GetTPCNclsF()) : 1.; // found/findable
- }
+ TClass *type = track->IsA();
+ if(TESTBIT(fTPCclusterRatioDef, kCROverFindable)){
+ AliDebug(2, "Using ratio def kCROverFindable");
+ clusterRatio = track->GetTPCNclsF() ? track->GetTPCClusterInfo(2,1)/static_cast<Double_t>(track->GetTPCNclsF()) : 1.; // crossed rows/findable
+ } else if(TESTBIT(fTPCclusterRatioDef, kFoundOverCR)){
+ AliDebug(2, "Using ratio def kFoundOverCR");
+ clusterRatio = track->GetTPCClusterInfo(2,0); // found/crossed rows
+ } else if(TESTBIT(fTPCclusterRatioDef, kFoundOverFindableIter1)){
+ if(type == AliESDtrack::Class()){
+ AliDebug(2, "Using ratio def kFoundOverFindableIter1");
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ clusterRatio = esdtrack->GetTPCNclsFIter1() ? static_cast<Double_t>(esdtrack->GetTPCNclsIter1())/static_cast<Double_t>(esdtrack->GetTPCNclsFIter1()) : 1.; // crossed
+ } else {
+ AliDebug(2, "Cluster ratio after iteration 1 not possible for AOD tracks");
+ clusterRatio = 1.;
}
+ } else if(TESTBIT(fTPCclusterRatioDef, kFoundOverFindable)) {
+ AliDebug(2, "Using ratio def kFoundOverFindable");
+ clusterRatio = track->GetTPCNclsF() ? static_cast<Double_t>(track->GetTPCNcls())/static_cast<Double_t>(track->GetTPCNclsF()) : 1.; // found/findable
+ }
+ else {
+ AliDebug(2, "Using ratio def kFoundAllOverFindable");
+ Int_t allclusters = 0;
+ if(type == AliESDtrack::Class()){
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ const TBits &clusterTPC = esdtrack->GetTPCClusterMap();
+ allclusters = clusterTPC.CountBits();
+ }
+ else {
+ AliAODTrack *aodtrack = static_cast<AliAODTrack *>(track);
+ const TBits &clusterTPC = aodtrack->GetTPCClusterMap();
+ allclusters = clusterTPC.CountBits();
+ }
+ clusterRatio = track->GetTPCNclsF() ? static_cast<Double_t>(allclusters)/static_cast<Double_t>(track->GetTPCNclsF()) : 1.; // foundall/findable
}
return clusterRatio;
+
}
+//___________________________________________________
+void AliHFEextraCuts::GetImpactParameters(AliVTrack *track, Float_t &radial, Float_t &z){
+ //
+ // Get impact parameter
+ //
+ const Double_t kBeampiperadius=3.;
+ Double_t dcaD[2]={-999.,-999.},
+ covD[3]={-999.,-999.,-999.};
+ TClass *type = track->IsA();
+ if(type == AliESDtrack::Class()){
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ esdtrack->GetImpactParameters(radial, z);
+ }
+ else if(type == AliAODTrack::Class()){
+
+ //case of AOD tracks
+ AliAODEvent *aodevent = dynamic_cast<AliAODEvent *>(fEvent);
+ if(!aodevent) {
+ AliDebug(1, "No aod event available\n");
+ return;
+ }
+
+ //Case ESD track: take copy constructor
+ AliAODTrack *aodtrack = NULL;
+ AliAODTrack *tmptrack = dynamic_cast<AliAODTrack *>(track);
+ if(tmptrack) aodtrack = new AliAODTrack(*tmptrack);
+
+ AliAODVertex *vtxAODSkip = aodevent->GetPrimaryVertex();
+ if(!vtxAODSkip) return;
+ AliExternalTrackParam etp; etp.CopyFromVTrack(aodtrack);
+ if(etp.PropagateToDCA(vtxAODSkip,aodevent->GetMagneticField(), kBeampiperadius, dcaD, covD)) {
+ radial = dcaD[0];
+ z = dcaD[1];
+ }
+ //if(vtxAODSkip) delete vtxAODSkip;
+ if(aodtrack) delete aodtrack;
+ }
+}
//______________________________________________________
-void AliHFEextraCuts::GetImpactParameters(AliVTrack *track, Float_t &radial, Float_t &z){
- //
- // Get impact parameter
- //
- TString type = track->IsA()->GetName();
- if(!type.CompareTo("AliESDtrack")){
- AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- if(esdtrack) esdtrack->GetImpactParameters(radial, z);
- }
- else if(!type.CompareTo("AliAODTrack")){
- AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
+Bool_t AliHFEextraCuts::IsKinkDaughter(AliVTrack *track){
+ //
+ // Is kink Daughter
+ //
+ TClass *type = track->IsA();
+ if(type == AliESDtrack::Class()){
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ if(!esdtrack) return kFALSE;
+ if(esdtrack->GetKinkIndex(0)>0) return kTRUE;
+ else return kFALSE;
+
+ }
+ else if(type == AliAODTrack::Class()){
+ AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
if(aodtrack){
- Double_t xyz[3];
- aodtrack->XYZAtDCA(xyz);
- z = xyz[2];
- radial = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]+xyz[1]);
+ AliAODVertex *aodvertex = aodtrack->GetProdVertex();
+ if(!aodvertex) return kFALSE;
+ if(aodvertex->GetType()==AliAODVertex::kKink) return kTRUE;
+ else return kFALSE;
}
- }
+ else return kFALSE;
+ }
+
+ return kFALSE;
+}
+//______________________________________________________
+Bool_t AliHFEextraCuts::IsKinkMother(AliVTrack *track){
+ //
+ // Is kink Mother: only for ESD since need to loop over vertices for AOD
+ //
+ //
+
+ TClass *type = track->IsA();
+ if(type == AliESDtrack::Class()){
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ if(!esdtrack) return kFALSE;
+ if(esdtrack->GetKinkIndex(0)!=0) return kTRUE;
+ else return kFALSE;
+ }
+
+ return kFALSE;
+
+}
+
+//______________________________________________________
+Float_t AliHFEextraCuts::GetTRDchi(AliVTrack *track){
+ //
+ // Get TRDchi2
+ //
+ Int_t ntls(0);
+ TClass *type = track->IsA();
+ if(type == AliESDtrack::Class()){
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ ntls = esdtrack->GetTRDntracklets();
+ return ntls ? esdtrack->GetTRDchi2()/ntls : -999;
+ }
+ else if(type == AliAODTrack::Class()){
+ AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
+ if(aodtrack){
+ return 999.;
+ }
+ }
+
+ return 999.;
+
}
+//______________________________________________________
+Int_t AliHFEextraCuts::GetITSNbOfcls(AliVTrack *track){
+ //
+ // Get ITS nb of clusters
+ //
+ TClass *type = track->IsA();
+ if(type == AliESDtrack::Class()){
+ AliESDtrack *esdtrack = static_cast<AliESDtrack *>(track);
+ return esdtrack->GetITSclusters(0);
+
+ }
+ else if(type == AliAODTrack::Class()){
+ AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
+ if(aodtrack){
+ return aodtrack->GetITSNcls();
+ }
+ }
+
+ return 0;
+}
//______________________________________________________
void AliHFEextraCuts::GetHFEImpactParameters(AliVTrack *track, Double_t &dcaxy, Double_t &dcansigmaxy){
+ //
+ // Get HFE impact parameter (with recalculated primary vertex)
+ //
+ dcaxy=0;
+ dcansigmaxy=0;
+ if(!fEvent){
+ AliDebug(1, "No Input event available\n");
+ return;
+ }
+ TString type = track->IsA()->GetName();
+ const Double_t kBeampiperadius=3.;
+ Double_t dcaD[2]={-999.,-999.},
+ covD[3]={-999.,-999.,-999.};
+
+ if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
+ //case of ESD tracks
+ AliESDEvent *esdevent = dynamic_cast<AliESDEvent *>(fEvent);
+ if(!esdevent) {
+ AliDebug(1, "No esd event available\n");
+ return;
+ }
+
+ //case ESD track: take copy constructor
+ AliESDtrack *esdtrack = NULL;
+ AliESDtrack *tmptrack = dynamic_cast<AliESDtrack *>(track);
+ if(tmptrack) esdtrack = new AliESDtrack(*tmptrack);
+
+ const AliVVertex *vtxESDSkip = esdevent->GetPrimaryVertex();
+ if(!vtxESDSkip) return;
+ if( vtxESDSkip->GetNContributors() < 30){ // if vertex contributor is smaller than 30, recalculate the primary vertex
+ vtxESDSkip = RemoveDaughtersFromPrimaryVtx(esdevent, track);
+ }
+
+ if(!vtxESDSkip) return;
+ if(esdtrack && esdtrack->PropagateToDCA(vtxESDSkip, fEvent->GetMagneticField(), kBeampiperadius, dcaD, covD)){
+ dcaxy = dcaD[0];
+ if(covD[0]) dcansigmaxy = dcaxy/TMath::Sqrt(covD[0]);
+ }
+ //delete vtxESDSkip;
+ delete esdtrack;
+ }
+ else {
+ //case of AOD tracks
+ AliAODEvent *aodevent = dynamic_cast<AliAODEvent *>(fEvent);
+ if(!aodevent) {
+ AliDebug(1, "No aod event available\n");
+ return;
+ }
+
+ //Case ESD track: take copy constructor
+ AliAODTrack *aodtrack = NULL;
+ AliAODTrack *tmptrack = dynamic_cast<AliAODTrack *>(track);
+ if(tmptrack) aodtrack = new AliAODTrack(*tmptrack);
+
+ AliAODVertex *vtxAODSkip = aodevent->GetPrimaryVertex();
+ if(!vtxAODSkip) return;
+ if(vtxAODSkip->GetNContributors() < 30){ // if vertex contributor is smaller than 30, recalculate the primary vertex
+ vtxAODSkip = RemoveDaughtersFromPrimaryVtx(aodevent, track);
+ }
+ if(!vtxAODSkip) return;
+ AliExternalTrackParam etp; etp.CopyFromVTrack(aodtrack);
+ if(etp.PropagateToDCA(vtxAODSkip,aodevent->GetMagneticField(), kBeampiperadius, dcaD, covD)) {
+ dcaxy = dcaD[0];
+ if(covD[0]) dcansigmaxy = dcaxy/TMath::Sqrt(covD[0]);
+ }
+ //if(vtxAODSkip) delete vtxAODSkip;
+ if(aodtrack) delete aodtrack;
+ }
+}
+
+//______________________________________________________
+void AliHFEextraCuts::GetHFEImpactParameters(AliVTrack *track, Double_t dcaD[2], Double_t covD[3]){
//
// Get HFE impact parameter (with recalculated primary vertex)
//
- dcaxy=0;
- dcansigmaxy=0;
if(!fEvent){
AliDebug(1, "No Input event available\n");
return;
}
const Double_t kBeampiperadius=3.;
TString type = track->IsA()->GetName();
- Double_t dcaD[2]={-999.,-999.};
- Double_t covD[3]={-999.,-999.,-999.};
- Float_t dcaF[2]={-999.,-999.};
- Float_t covF[3]={-999.,-999.,-999.};
- AliESDEvent *esdevent = dynamic_cast<AliESDEvent *>(fEvent);
- if(!esdevent) {
- AliDebug(1, "No esd event available\n");
- return;
- }
- const AliVVertex *vtxESDSkip = esdevent->GetPrimaryVertex();
- if( vtxESDSkip->GetNContributors() < 30){ // if vertex contributor is smaller than 30, recalculate the primary vertex
- // recalculate primary vertex for peri. and pp
- AliVertexerTracks vertexer(fEvent->GetMagneticField());
- vertexer.SetITSMode();
- vertexer.SetMinClusters(4);
- Int_t skipped[2];
- skipped[0] = track->GetID();
- vertexer.SetSkipTracks(1,skipped);
-
- //----diamond constraint-----------------------------
- vertexer.SetConstraintOn();
- Float_t diamondcovxy[3];
- esdevent->GetDiamondCovXY(diamondcovxy);
- Double_t pos[3]={esdevent->GetDiamondX(),esdevent->GetDiamondY(),0.};
- Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
- AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1);
- vertexer.SetVtxStart(diamond);
- delete diamond; diamond=NULL;
- //----------------------------------------------------
-
- vtxESDSkip = vertexer.FindPrimaryVertex(fEvent);
-
- // Getting the DCA
- // Propagation always done on a working copy to not disturb the track params of the original track
- AliESDtrack *esdtrack = NULL;
- if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
- // Case ESD track: take copy constructor
- AliESDtrack *tmptrack = dynamic_cast<AliESDtrack *>(track);
- if(tmptrack) esdtrack = new AliESDtrack(*tmptrack);
- } else {
- // Case AOD track: take different constructor
- esdtrack = new AliESDtrack(track);
- }
- if(esdtrack && esdtrack->PropagateToDCA(vtxESDSkip, fEvent->GetMagneticField(), kBeampiperadius, dcaD, covD)){
- // protection
- dcaxy = dcaD[0];
- if(covD[0]) dcansigmaxy = dcaxy/TMath::Sqrt(covD[0]);
- if(!covD[0]) dcansigmaxy = -49.;
- }
- delete esdtrack;
- delete vtxESDSkip;
- }
- else{
- AliESDtrack *esdtrack = NULL;
- if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
- // Case ESD track: take copy constructor
+ if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
+ //case of ESD tracks
+ AliESDEvent *esdevent = dynamic_cast<AliESDEvent *>(fEvent);
+ if(!esdevent) {
+ AliDebug(1, "No esd event available\n");
+ return;
+ }
+
+ //case ESD track: take copy constructor
+ AliESDtrack *esdtrack = NULL;
AliESDtrack *tmptrack = dynamic_cast<AliESDtrack *>(track);
if(tmptrack) esdtrack = new AliESDtrack(*tmptrack);
- } else {
- // Case AOD track: take different constructor
- esdtrack = new AliESDtrack(track);
- }
- if(esdtrack) esdtrack->GetImpactParameters(dcaF, covF);
- dcaxy = dcaF[0];
- if(covF[0]) dcansigmaxy = dcaxy/TMath::Sqrt(covF[0]);
- if(!covF[0]) dcansigmaxy = -49.;
- delete esdtrack;
+
+ const AliVVertex *vtxESDSkip = esdevent->GetPrimaryVertex();
+ if(!vtxESDSkip) return;
+ if( vtxESDSkip->GetNContributors() < 30){ // if vertex contributor is smaller than 30, recalculate the primary vertex
+ vtxESDSkip = RemoveDaughtersFromPrimaryVtx(esdevent, track);
+ }
+ if(!vtxESDSkip) return;
+
+ if(esdtrack) esdtrack->PropagateToDCA(vtxESDSkip, fEvent->GetMagneticField(), kBeampiperadius, dcaD, covD);
+
+ //delete vtxESDSkip;
+ delete esdtrack;
}
-}
+ else {
+ //case of AOD tracks
+ AliAODEvent *aodevent = dynamic_cast<AliAODEvent *>(fEvent);
+ if(!aodevent) {
+ AliDebug(1, "No aod event available\n");
+ return;
+ }
+
+ //Case ESD track: take copy constructor
+ AliAODTrack *aodtrack = NULL;
+ AliAODTrack *tmptrack = dynamic_cast<AliAODTrack *>(track);
+ if(tmptrack) aodtrack = new AliAODTrack(*tmptrack);
+ AliAODVertex *vtxAODSkip = aodevent->GetPrimaryVertex();
+ if(!vtxAODSkip) return;
+ if(vtxAODSkip->GetNContributors() < 30){ // if vertex contributor is smaller than 30, recalculate the primary vertex
+ vtxAODSkip = RemoveDaughtersFromPrimaryVtx(aodevent, track);
+ }
+ if(!vtxAODSkip) return;
+ AliExternalTrackParam etp; etp.CopyFromVTrack(aodtrack);
+ etp.PropagateToDCA(vtxAODSkip,aodevent->GetMagneticField(), kBeampiperadius, dcaD, covD);
+
+ //delete vtxAODSkip;
+ delete aodtrack;
+ }
+}
//______________________________________________________
void AliHFEextraCuts::GetHFEImpactParameterCuts(AliVTrack *track, Double_t &hfeimpactRcut, Double_t &hfeimpactnsigmaRcut){
else maximpactRcut = 9999999999.0;
}
}
+//______________________________________________________
+void AliHFEextraCuts::GetTOFsignalDxDz(AliVTrack *track, Double_t &tofsignalDx, Double_t &tofsignalDz){
+ //
+ // TOF matching
+ //
+
+ TString type = track->IsA()->GetName();
+ if(!type.CompareTo("AliESDtrack")){
+ AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+ if(!esdtrack) return;
+ tofsignalDx = esdtrack->GetTOFsignalDx();
+ tofsignalDz = esdtrack->GetTOFsignalDz();
+ }
+
+}
+
+//______________________________________________________
+Bool_t AliHFEextraCuts::CheckITSpattern(const AliVTrack *const track) const {
+ //
+ // Check if every ITS layer, which has a module which is alive, also
+ // has an ITS cluster
+ //
+ Bool_t patternOK(kTRUE);
+ Int_t status(0);
+ for(Int_t ily = 0; ily < 6; ily++){
+ status = GetITSstatus(track, ily);
+ if(CheckITSstatus(status)){
+ // pixel alive, check whether layer has a cluster
+ if(!TESTBIT(track->GetITSClusterMap(),ily)){
+ // No cluster even though pixel is alive - reject track
+ patternOK = kFALSE;
+ break;
+ }
+ }
+ }
+ return patternOK;
+}
+
+//---------------------------------------------------------------------------
+const AliVVertex* AliHFEextraCuts::RemoveDaughtersFromPrimaryVtx(AliESDEvent *esdevent, AliVTrack *track) {
+ //
+ // This method returns a primary vertex without the daughter tracks of the
+ // candidate and it recalculates the impact parameters and errors for ESD tracks.
+ //
+ // The output vertex is created with "new".
+ //
+
+ const AliVVertex *vtxESDSkip = esdevent->GetPrimaryVertex();
+
+ AliVertexerTracks vertexer(fEvent->GetMagneticField());
+ vertexer.SetITSMode();
+ vertexer.SetMinClusters(4);
+ Int_t skipped[2];
+ skipped[0] = track->GetID();
+ vertexer.SetSkipTracks(1,skipped);
+
+ //diamond constraint
+ vertexer.SetConstraintOn();
+ Float_t diamondcovxy[3];
+ esdevent->GetDiamondCovXY(diamondcovxy);
+ Double_t pos[3]={esdevent->GetDiamondX(),esdevent->GetDiamondY(),0.};
+ Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
+ AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1);
+ vertexer.SetVtxStart(diamond);
+ delete diamond; diamond=NULL;
+
+ vtxESDSkip = vertexer.FindPrimaryVertex(fEvent);
+
+ return vtxESDSkip;
+}
+
+//---------------------------------------------------------------------------
+AliAODVertex* AliHFEextraCuts::RemoveDaughtersFromPrimaryVtx(AliAODEvent *aod, AliVTrack *track) {
+ //
+ // This method returns a primary vertex without the daughter tracks of the
+ // candidate and it recalculates the impact parameters and errors for AOD tracks.
+ // The output vertex is created with "new".
+ //
+
+ AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
+ if(!vtxAOD) return 0;
+ TString title=vtxAOD->GetTitle();
+ if(!title.Contains("VertexerTracks")) return 0;
+
+ AliVertexerTracks *vertexer = new AliVertexerTracks(aod->GetMagneticField());
+
+ vertexer->SetITSMode();
+ vertexer->SetMinClusters(3);
+ vertexer->SetConstraintOff();
+
+ if(title.Contains("WithConstraint")) {
+ Float_t diamondcovxy[3];
+ aod->GetDiamondCovXY(diamondcovxy);
+ Double_t pos[3]={aod->GetDiamondX(),aod->GetDiamondY(),0.};
+ Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
+ AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1);
+ vertexer->SetVtxStart(diamond);
+ delete diamond; diamond=NULL;
+ }
+ Int_t skipped[2]; for(Int_t i=0;i<2;i++) skipped[i]=-1;
+ Int_t id = (Int_t)track->GetID();
+ if(!(id<0)) skipped[0] = id;
+
+ /*// exclude tracks with global constrained parameters
+ Int_t nTracks=aod->GetNumberOfTracks();
+ for(Int_t i=0; i<nTracks; i++){
+ t = aod->GetTrack(i);
+ if(t->TestFilterMask(512)){
+ id = (Int_t)t->GetID();
+ skipped[nTrksToSkip++] = id;
+ }
+ }*/
+
+ vertexer->SetSkipTracks(1,skipped);
+ AliESDVertex *vtxESDNew = vertexer->FindPrimaryVertex(aod);
+
+ delete vertexer; vertexer=NULL;
+
+ if(!vtxESDNew) return 0;
+ if(vtxESDNew->GetNContributors()<=0) {
+ delete vtxESDNew; vtxESDNew=NULL;
+ return 0;
+ }
+
+ // convert to AliAODVertex
+ Double_t pos[3],cov[6],chi2perNDF;
+ vtxESDNew->GetXYZ(pos); // position
+ vtxESDNew->GetCovMatrix(cov); //covariance matrix
+ chi2perNDF = vtxESDNew->GetChi2toNDF();
+ delete vtxESDNew; vtxESDNew=NULL;
+
+ AliAODVertex *vtxAODNew = new AliAODVertex(pos,cov,chi2perNDF);
+
+ return vtxAODNew;
+}