+//--------------------------------------------------------------------------
+void AliESDtrackCuts::SetPtDepDCACuts(Double_t pt) {
+ //
+ // set the pt-dependent DCA cuts
+ //
+
+ if(f1CutMaxDCAToVertexXYPtDep) {
+ fCutMaxDCAToVertexXY=f1CutMaxDCAToVertexXYPtDep->Eval(pt);
+ }
+
+ if(f1CutMaxDCAToVertexZPtDep) {
+ fCutMaxDCAToVertexZ=f1CutMaxDCAToVertexZPtDep->Eval(pt);
+ }
+
+ if(f1CutMinDCAToVertexXYPtDep) {
+ fCutMinDCAToVertexXY=f1CutMinDCAToVertexXYPtDep->Eval(pt);
+ }
+
+ if(f1CutMinDCAToVertexZPtDep) {
+ fCutMinDCAToVertexZ=f1CutMinDCAToVertexZPtDep->Eval(pt);
+ }
+
+
+ return;
+}
+
+
+
+//--------------------------------------------------------------------------
+Bool_t AliESDtrackCuts::CheckPtDepDCA(TString dist,Bool_t print) const {
+ //
+ // Check the correctness of the string syntax
+ //
+ Bool_t retval=kTRUE;
+
+ if(!dist.Contains("pt")) {
+ if(print) AliError("string must contain \"pt\"");
+ retval= kFALSE;
+ }
+ return retval;
+}
+
+ void AliESDtrackCuts::SetMaxDCAToVertexXYPtDep(const char *dist){
+
+ if(f1CutMaxDCAToVertexXYPtDep){
+ delete f1CutMaxDCAToVertexXYPtDep;
+ // resetiing both
+ f1CutMaxDCAToVertexXYPtDep = 0;
+ fCutMaxDCAToVertexXYPtDep = "";
+ }
+ if(!CheckPtDepDCA(dist,kTRUE)){
+ return;
+ }
+ fCutMaxDCAToVertexXYPtDep = dist;
+ TString tmp(dist);
+ tmp.ReplaceAll("pt","x");
+ f1CutMaxDCAToVertexXYPtDep = new TFormula("f1CutMaxDCAToVertexXYPtDep",tmp.Data());
+
+}
+
+ void AliESDtrackCuts::SetMaxDCAToVertexZPtDep(const char *dist){
+
+
+ if(f1CutMaxDCAToVertexZPtDep){
+ delete f1CutMaxDCAToVertexZPtDep;
+ // resetiing both
+ f1CutMaxDCAToVertexZPtDep = 0;
+ fCutMaxDCAToVertexZPtDep = "";
+ }
+ if(!CheckPtDepDCA(dist,kTRUE))return;
+
+ fCutMaxDCAToVertexZPtDep = dist;
+ TString tmp(dist);
+ tmp.ReplaceAll("pt","x");
+ f1CutMaxDCAToVertexZPtDep = new TFormula("f1CutMaxDCAToVertexZPtDep",tmp.Data());
+
+
+}
+
+
+ void AliESDtrackCuts::SetMinDCAToVertexXYPtDep(const char *dist){
+
+
+ if(f1CutMinDCAToVertexXYPtDep){
+ delete f1CutMinDCAToVertexXYPtDep;
+ // resetiing both
+ f1CutMinDCAToVertexXYPtDep = 0;
+ fCutMinDCAToVertexXYPtDep = "";
+ }
+ if(!CheckPtDepDCA(dist,kTRUE))return;
+
+ fCutMinDCAToVertexXYPtDep = dist;
+ TString tmp(dist);
+ tmp.ReplaceAll("pt","x");
+ f1CutMinDCAToVertexXYPtDep = new TFormula("f1CutMinDCAToVertexXYPtDep",tmp.Data());
+
+}
+
+
+ void AliESDtrackCuts::SetMinDCAToVertexZPtDep(const char *dist){
+
+
+
+ if(f1CutMinDCAToVertexZPtDep){
+ delete f1CutMinDCAToVertexZPtDep;
+ // resetiing both
+ f1CutMinDCAToVertexZPtDep = 0;
+ fCutMinDCAToVertexZPtDep = "";
+ }
+ if(!CheckPtDepDCA(dist,kTRUE))return;
+ fCutMinDCAToVertexZPtDep = dist;
+ TString tmp(dist);
+ tmp.ReplaceAll("pt","x");
+ f1CutMinDCAToVertexZPtDep = new TFormula("f1CutMinDCAToVertexZPtDep",tmp.Data());
+}
+
+AliESDtrackCuts* AliESDtrackCuts::GetMultEstTrackCuts(MultEstTrackCuts cut)
+{
+ // returns the multiplicity estimator track cuts objects to allow for user configuration
+ // upon first call the objects are created
+ //
+ // the cut defined here correspond to GetStandardITSTPCTrackCuts2010 (apart from the one for without SPD)
+
+ if (!fgMultEstTrackCuts[kMultEstTrackCutGlobal])
+ {
+ // quality cut on ITS+TPC tracks
+ fgMultEstTrackCuts[kMultEstTrackCutGlobal] = new AliESDtrackCuts();
+ fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetMinNClustersTPC(70);
+ fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetMaxChi2PerClusterTPC(4);
+ fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetAcceptKinkDaughters(kFALSE);
+ fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetRequireTPCRefit(kTRUE);
+ fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetRequireITSRefit(kTRUE);
+ //multiplicity underestimate if we use global tracks with |eta| > 0.9
+ fgMultEstTrackCuts[kMultEstTrackCutGlobal]->SetEtaRange(-0.9, 0.9);
+
+ // quality cut on ITS_SA tracks (complementary to ITS+TPC)
+ fgMultEstTrackCuts[kMultEstTrackCutITSSA] = new AliESDtrackCuts();
+ fgMultEstTrackCuts[kMultEstTrackCutITSSA]->SetRequireITSRefit(kTRUE);
+
+ // primary selection for tracks with SPD hits
+ fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD] = new AliESDtrackCuts();
+ fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
+ fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
+ fgMultEstTrackCuts[kMultEstTrackCutDCAwSPD]->SetMaxDCAToVertexZ(2);
+
+ // primary selection for tracks w/o SPD hits
+ fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD] = new AliESDtrackCuts();
+ fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kNone);
+ fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetMaxDCAToVertexXYPtDep("1.5*(0.0182+0.0350/pt^1.01)");
+ fgMultEstTrackCuts[kMultEstTrackCutDCAwoSPD]->SetMaxDCAToVertexZ(2);
+ }
+
+ return fgMultEstTrackCuts[cut];
+}
+
+Int_t AliESDtrackCuts::GetReferenceMultiplicity(const AliESDEvent* esd, MultEstTrackType trackType, Float_t etaRange)
+{
+ // Get multiplicity estimate based on TPC/ITS tracks and tracklets
+ // Adapted for AliESDtrackCuts from a version developed by: Ruben Shahoyan, Anton Alkin, Arvinder Palaha
+ //
+ // Returns a negative value if no reliable estimate can be provided:
+ // -1 SPD vertex missing
+ // -2 SPD VertexerZ dispersion too large
+ // -3 Track vertex missing (not checked for kTracklets)
+ // -4 Distance between SPD and track vertices too large (not checked for kTracklets)
+ //
+ // WARNING This functions does not cut on the z vtx. Depending on the eta range requested, you need to restrict your z vertex range!
+ //
+ // Strategy for combined estimators
+ // 1. Count global tracks and flag them
+ // 2. Count ITSSA as complementaries for ITSTPC+ or as main tracks
+ // 3. Count the complementary SPD tracklets
+
+ const AliESDVertex* vertices[2];
+ vertices[0] = esd->GetPrimaryVertexSPD();
+ vertices[1] = esd->GetPrimaryVertexTracks();
+
+ if (!vertices[0]->GetStatus())
+ {
+ AliDebugClass(AliLog::kDebug, "No SPD vertex. Not able to make a reliable multiplicity estimate.");
+ return -1;
+ }
+
+ if (vertices[0]->IsFromVertexerZ() && vertices[0]->GetDispersion() > 0.02)
+ {
+ AliDebugClass(AliLog::kDebug, "Vertexer z dispersion > 0.02. Not able to make a reliable multiplicity estimate.");
+ return -2;
+ }
+
+ Int_t multiplicityEstimate = 0;
+
+ // SPD tracklet-only estimate
+ if (trackType == kTracklets)
+ {
+ const AliMultiplicity* spdmult = esd->GetMultiplicity(); // spd multiplicity object
+ for (Int_t i=0; i<spdmult->GetNumberOfTracklets(); ++i)
+ {
+ if (TMath::Abs(spdmult->GetEta(i)) > etaRange)
+ continue; // eta selection for tracklets
+ multiplicityEstimate++;
+ }
+ return multiplicityEstimate;
+ }
+
+ if (!vertices[1]->GetStatus())
+ {
+ AliDebugClass(AliLog::kDebug, "No track vertex. Not able to make a reliable multiplicity estimate.");
+ return -3;
+ }
+
+ // TODO value of displacement to be studied
+ const Float_t maxDisplacement = 0.5;
+ //check for displaced vertices
+ Double_t displacement = TMath::Abs(vertices[0]->GetZ() - vertices[1]->GetZ());
+ if (displacement > maxDisplacement)
+ {
+ AliDebugClass(AliLog::kDebug, Form("Displaced vertices %f > %f",displacement,maxDisplacement));
+ return -4;
+ }
+
+ // update eta range in track cuts
+ GetMultEstTrackCuts(kMultEstTrackCutITSSA)->SetEtaRange(-etaRange, etaRange);
+ GetMultEstTrackCuts(kMultEstTrackCutDCAwSPD)->SetEtaRange(-etaRange, etaRange);
+ GetMultEstTrackCuts(kMultEstTrackCutDCAwoSPD)->SetEtaRange(-etaRange, etaRange);
+
+ //*******************************************************************************************************
+ //set counters to initial zeros
+ Int_t tracksITSTPC = 0; //number of global tracks for a given event
+ Int_t tracksITSSA = 0; //number of ITS standalone tracks for a given event
+ Int_t tracksITSTPCSA_complementary = 0; //number of ITS standalone tracks complementary to TPC for a given event
+ Int_t trackletsITSTPC_complementary = 0;//number of SPD tracklets complementary to global/ITSSA tracks for a given events
+ Int_t trackletsITSSA_complementary = 0; //number of SPD tracklets complementary to ITSSA tracks for a given event
+
+ const Int_t nESDTracks = esd->GetNumberOfTracks();
+ Int_t highestID = 0;
+
+ // flags for secondary and rejected tracks
+ const Int_t kRejBit = BIT(15); // set this bit in global tracks if it is rejected by a cut
+ const Int_t kSecBit = BIT(16); // set this bit in global tracks if it is secondary according to a cut
+
+ for(Int_t itracks=0; itracks < nESDTracks; itracks++) {
+ if(esd->GetTrack(itracks)->GetLabel() > highestID) highestID = esd->GetTrack(itracks)->GetLabel();
+ esd->GetTrack(itracks)->ResetBit(kSecBit|kRejBit); //reset bits used for flagging secondaries and rejected tracks in case they were changed before this analysis
+ }
+ const Int_t maxid = highestID+1; // used to define bool array for check multiple associations of tracklets to one track. array starts at 0.
+
+ // bit mask for esd tracks, to check if multiple tracklets are associated to it
+ Bool_t globalBits[maxid], pureITSBits[maxid];
+ for(Int_t i=0; i<maxid; i++){ // set all bools to false
+ globalBits[i]=kFALSE;
+ pureITSBits[i]=kFALSE;
+ }