[u/mrichter/AliRoot.git] / PWGLF / RESONANCES / AliRsnCutTrackQuality.cxx

//
// Class AliRsnCutTrackQuality
//
// General implementation of a single cut strategy, which can be:
// - a value contained in a given interval  [--> IsBetween()   ]
// - a value equal to a given reference     [--> MatchesValue()]
//
// In all cases, the reference value(s) is (are) given as data members
// and each kind of cut requires a given value type (Int, UInt, Double),
// but the cut check procedure is then automatized and chosen thanks to
// an enumeration of the implemented cut types.
// At the end, the user (or any other point which uses this object) has
// to use the method IsSelected() to check if this cut has been passed.
//
// authors: Martin Vala (martin.vala@cern.ch)
//          Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
//

#include <Riostream.h>
#include <TFormula.h>
#include <TBits.h>

#include "AliLog.h"
#include "AliESDtrackCuts.h"

#include "AliRsnEvent.h"
#include "AliRsnDaughter.h"
#include "AliRsnCutTrackQuality.h"

ClassImp(AliRsnCutTrackQuality)

//_________________________________________________________________________________________________
AliRsnCutTrackQuality::AliRsnCutTrackQuality(const char *name) :
   AliRsnCut(name, AliRsnTarget::kDaughter, 0.0, 0.0),
   fFlagsOn(0x0),
   fFlagsOff(0x0),
   fRejectKinkDaughters(kTRUE),
   fDCARfixed(kTRUE),
   fDCARptFormula(""),
   fDCARmax(1E20),
   fDCAZfixed(kTRUE),
   fDCAZptFormula(""),
   fDCAZmax(1E20),
   fSPDminNClusters(0),
   fITSminNClusters(0),
   fITSmaxChi2(1E20),
   fTPCminNClusters(0),
   fTPCmaxChi2(1E20),
   fCutMaxChi2TPCConstrainedVsGlobal(1E20),
   fAODTestFilterBit(-1),
   fESDtrackCuts(0x0)
{
//
// Default constructor.
// Initializes all cuts in such a way that all of them are disabled.
//

   SetPtRange(0.0, 1E20);
   SetEtaRange(-1E20, 1E20);
}

//_________________________________________________________________________________________________
AliRsnCutTrackQuality::AliRsnCutTrackQuality(const AliRsnCutTrackQuality &copy) :
   AliRsnCut(copy),
   fFlagsOn(copy.fFlagsOn),
   fFlagsOff(copy.fFlagsOff),
   fRejectKinkDaughters(copy.fRejectKinkDaughters),
   fDCARfixed(copy.fDCARfixed),
   fDCARptFormula(copy.fDCARptFormula),
   fDCARmax(copy.fDCARmax),
   fDCAZfixed(copy.fDCAZfixed),
   fDCAZptFormula(copy.fDCAZptFormula),
   fDCAZmax(copy.fDCAZmax),
   fSPDminNClusters(copy.fSPDminNClusters),
   fITSminNClusters(copy.fITSminNClusters),
   fITSmaxChi2(copy.fITSmaxChi2),
   fTPCminNClusters(copy.fTPCminNClusters),
   fTPCmaxChi2(copy.fTPCmaxChi2),
   fCutMaxChi2TPCConstrainedVsGlobal(copy.fCutMaxChi2TPCConstrainedVsGlobal),
   fAODTestFilterBit(copy.fAODTestFilterBit),
   fESDtrackCuts(copy.fESDtrackCuts)
{
//
// Copy constructor.
// Just copy all data member values.
//

   SetPtRange(copy.fPt[0], copy.fPt[1]);
   SetEtaRange(copy.fEta[0], copy.fEta[1]);
}

//_________________________________________________________________________________________________
AliRsnCutTrackQuality &AliRsnCutTrackQuality::operator=(const AliRsnCutTrackQuality &copy)
{
//
// Assignment operator.
// Just copy all data member values.
//

   if (this == &copy)
      return *this;

   fFlagsOn = copy.fFlagsOn;
   fFlagsOff = copy.fFlagsOff;
   fRejectKinkDaughters = copy.fRejectKinkDaughters;
   fDCARfixed = copy.fDCARfixed;
   fDCARptFormula = copy.fDCARptFormula;
   fDCARmax = copy.fDCARmax;
   fDCAZfixed = copy.fDCAZfixed;
   fDCAZptFormula = copy.fDCAZptFormula;
   fDCAZmax = copy.fDCAZmax;
   fSPDminNClusters = copy.fSPDminNClusters;
   fITSminNClusters = copy.fITSminNClusters;
   fITSmaxChi2 = copy.fITSmaxChi2;
   fTPCminNClusters = copy.fTPCminNClusters;
   fTPCmaxChi2 = copy.fTPCmaxChi2;
   fAODTestFilterBit = copy.fAODTestFilterBit;
   fESDtrackCuts = copy.fESDtrackCuts;
   SetPtRange(copy.fPt[0], copy.fPt[1]);
   SetEtaRange(copy.fEta[0], copy.fEta[1]);

   return (*this);
}

//_________________________________________________________________________________________________
void AliRsnCutTrackQuality::DisableAll()
{
//
// Disable all cuts
//

   fFlagsOn = 0x0;
   fFlagsOff = 0x0;
   fRejectKinkDaughters = kFALSE;
   fDCARfixed = kTRUE;
   fDCARptFormula = "";
   fDCARmax = 1E20;
   fDCAZfixed = kTRUE;
   fDCAZptFormula = "";
   fDCAZmax = 1E20;
   fSPDminNClusters = 0;
   fITSminNClusters = 0;
   fITSmaxChi2 = 1E20;
   fTPCminNClusters = 0;
   fTPCmaxChi2 = 1E20;
   fAODTestFilterBit = -1;
   if (fESDtrackCuts) {
      const char *cutsName = fESDtrackCuts->GetName();
      const char *cutsTitle = fESDtrackCuts->GetTitle();
      delete fESDtrackCuts;
      fESDtrackCuts = new AliESDtrackCuts(cutsName,cutsTitle);
   }
   SetPtRange(0.0, 1E20);
   SetEtaRange(-1E20, 1E20);
}

//_________________________________________________________________________________________________
Bool_t AliRsnCutTrackQuality::IsSelected(TObject *object)
{
//
// Cut checker.
// Checks the type of object being evaluated
// and then calls the appropriate sub-function (for ESD or AOD)
//

   // coherence check
   if (!TargetOK(object)) return kFALSE;

   // status is checked in the same way for all tracks, using AliVTrack
   // as a convention, if a the collection of 'on' flags is '0x0', it
   // is assumed that no flags are required, and this check is skipped;
   // for the collection of 'off' flags this is not needed
   AliVTrack *vtrack = fDaughter->Ref2Vtrack();
   if (!vtrack) {
      AliDebug(AliLog::kDebug + 2, "This object is not either an ESD nor AOD track");
      return kFALSE;
   }
   ULong_t status   = (ULong_t)vtrack->GetStatus();
   ULong_t checkOn  = status & fFlagsOn;
   ULong_t checkOff = status & fFlagsOff;
   if (fFlagsOn != 0x0 && checkOn != fFlagsOn) {
      AliDebug(AliLog::kDebug + 2, Form("Failed flag check: required  %s", Binary(fFlagsOn)));
      AliDebug(AliLog::kDebug + 2, Form("                   track has %s", Binary(status  )));
      return kFALSE;
   }
   if (checkOff != 0) {
      AliDebug(AliLog::kDebug + 2, Form("Failed flag check: forbidden %s", Binary(fFlagsOff)));
      AliDebug(AliLog::kDebug + 2, Form("                   track has %s", Binary(status  )));
      return kFALSE;
   }
   AliDebug(AliLog::kDebug + 3, Form("Flag check OK: required  %s", Binary(fFlagsOn)));
   AliDebug(AliLog::kDebug + 3, Form("               forbidden %s", Binary(fFlagsOff)));
   AliDebug(AliLog::kDebug + 3, Form("               track has %s", Binary(status  )));

   // retrieve real object type
   AliESDtrack *esdTrack = fDaughter->Ref2ESDtrack();
   AliAODTrack *aodTrack = fDaughter->Ref2AODtrack();
   if (esdTrack) {
      AliDebug(AliLog::kDebug + 2, "Checking an ESD track");
      if (fESDtrackCuts)
         return fESDtrackCuts->IsSelected(esdTrack);
      else
         return CheckESD(esdTrack);
   } else if (aodTrack) {
      AliDebug(AliLog::kDebug + 2, "Checking an AOD track");
      return CheckAOD(aodTrack);
   } else {
      AliDebug(AliLog::kDebug + 2, Form("This object is not either an ESD nor AOD track, it is an %s", fDaughter->GetRef()->ClassName()));
      return kFALSE;
   }
}

//_________________________________________________________________________________________________
Bool_t AliRsnCutTrackQuality::CheckESD(AliESDtrack *track)
{
//
// Check an ESD track.
// This is done using the default track checker for ESD.
// It is declared static, not to recreate it every time.
//

   static AliESDtrackCuts cuts;

   // general acceptance/pt cuts
   cuts.SetPtRange(fPt[0], fPt[1]);
   cuts.SetEtaRange(fEta[0], fEta[1]);

   // transverse DCA cuts
   if (fDCARfixed)
      cuts.SetMaxDCAToVertexXY(fDCARmax);
   else
      cuts.SetMaxDCAToVertexXYPtDep(fDCARptFormula.Data());

   // longitudinal DCA cuts
   if (fDCAZfixed)
      cuts.SetMaxDCAToVertexZ(fDCAZmax);
   else
      cuts.SetMaxDCAToVertexZPtDep(fDCAZptFormula.Data());

   // these options are always disabled in current version
   cuts.SetDCAToVertex2D(kFALSE);
   cuts.SetRequireSigmaToVertex(kFALSE);

   // TPC related cuts for TPC+ITS tracks
   cuts.SetMinNClustersTPC(fTPCminNClusters);
   cuts.SetMaxChi2PerClusterTPC(fTPCmaxChi2);
   cuts.SetAcceptKinkDaughters(!fRejectKinkDaughters);
   cuts.SetMaxChi2TPCConstrainedGlobal(fCutMaxChi2TPCConstrainedVsGlobal);

   // ITS related cuts for TPC+ITS tracks
   if (fSPDminNClusters > 0)
      cuts.SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
   cuts.SetMaxChi2PerClusterITS(fITSmaxChi2);

   // now that all is initialized, do the check
   return cuts.IsSelected(track);
}

//_________________________________________________________________________________________________
Bool_t AliRsnCutTrackQuality::CheckAOD(AliAODTrack *track)
{
//
// Check an AOD track.
// This is done doing directly all checks, since there is not
// an equivalend checker for AOD tracks
//

   // if a test bit is used, check it and skip the following
   if (fAODTestFilterBit >= 0) {
      UInt_t bit = 1 << fAODTestFilterBit;
      AliDebugClass(2, Form("Required a test filter bit for AOD check: %u (result: %s)", bit, (track->TestFilterBit(bit) ? "accept" : "reject")));
      if (!track->TestFilterBit(bit))
         return kFALSE;
      else {
         if (track->Pt() < fPt[0] || track->Pt() > fPt[1]) return kFALSE;
         if (track->Eta() < fEta[0] || track->Eta() > fEta[1]) return kFALSE;
         return kTRUE;
      }
   }

   // try to retrieve the reference AOD event
   AliAODEvent *aodEvent = 0x0;
   if (fEvent) aodEvent = fEvent->GetRefAOD();
   if (!aodEvent) {
      AliError("AOD reference event is not initialized!");
      return kFALSE;
   }

   // step #0: check SPD and ITS clusters
   Int_t nSPD = 0;
   nSPD  = TESTBIT(track->GetITSClusterMap(), 0);
   nSPD += TESTBIT(track->GetITSClusterMap(), 1);
   if (nSPD < fSPDminNClusters) {
      AliDebug(AliLog::kDebug + 2, "Not enough SPD clusters in this track. Rejected");
      return kFALSE;
   }


   // step #1: check number of clusters in TPC
   if (track->GetTPCNcls() < fTPCminNClusters) {
      AliDebug(AliLog::kDebug + 2, "Too few TPC clusters. Rejected");
      return kFALSE;
   }
   if (track->GetITSNcls() < fITSminNClusters) {
      AliDebug(AliLog::kDebug + 2, "Too few ITS clusters. Rejected");
      return kFALSE;
   }

   // step #2: check chi square
   if (track->Chi2perNDF() > fTPCmaxChi2) {
      AliDebug(AliLog::kDebug + 2, "Bad chi2. Rejected");
      return kFALSE;
   }
   if (track->Chi2perNDF() > fITSmaxChi2) {
      AliDebug(AliLog::kDebug + 2, "Bad chi2. Rejected");
      return kFALSE;
   }

   // step #3: reject kink daughters
   AliAODVertex *vertex = track->GetProdVertex();
   if (vertex && fRejectKinkDaughters) {
      if (vertex->GetType() == AliAODVertex::kKink) {
         AliDebug(AliLog::kDebug + 2, "Kink daughter. Rejected");
         return kFALSE;
      }
   }

   // step #4: DCA cut (transverse)
   // --> reject all tracks not ITS refitted
   Double_t b[2], cov[3];
   vertex = aodEvent->GetPrimaryVertex();
   if (!vertex) {
      AliDebug(AliLog::kDebug + 2, "NULL vertex");
      return kFALSE;
   }
   if ((track->GetStatus() & AliESDtrack::kITSrefit) == 0) {
      AliDebug(AliLog::kDebug + 2, "Not ITS refitted");
      return kFALSE;
   }
   if (!track->PropagateToDCA(vertex, aodEvent->GetMagneticField(), kVeryBig, b, cov)) {
      AliDebug(AliLog::kDebug + 2, "Failed propagation to vertex");
      return kFALSE;
   }
   // if the DCA cut is not fixed, compute current value
   if (!fDCARfixed) {
      static TString str(fDCARptFormula);
      str.ReplaceAll("pt", "x");
      static const TFormula dcaXY(Form("%s_dcaXY", GetName()), str.Data());
      fDCARmax = dcaXY.Eval(track->Pt());
   }
   // check the cut
   if (TMath::Abs(b[0]) > fDCARmax) {
      AliDebug(AliLog::kDebug + 2, "Too large transverse DCA");
      return kFALSE;
   }

   // step #5: DCA cut (longitudinal)
   // the DCA has already been computed above
   // if the DCA cut is not fixed, compute current value
   if (!fDCAZfixed) {
      static TString str(fDCAZptFormula);
      str.ReplaceAll("pt", "x");
      static const TFormula dcaZ(Form("%s_dcaXY", GetName()), str.Data());
      fDCAZmax = dcaZ.Eval(track->Pt());
   }
   // check the cut
   if (TMath::Abs(b[1]) > fDCAZmax) {
      AliDebug(AliLog::kDebug + 2, "Too large longitudinal DCA");
      return kFALSE;
   }

   // step #6: check eta/pt range
   if (track->Eta() < fEta[0] || track->Eta() > fEta[1]) {
      AliDebug(AliLog::kDebug + 2, "Outside ETA acceptance");
      return kFALSE;
   }
   if (track->Pt() < fPt[0] || track->Pt() > fPt[1]) {
      AliDebug(AliLog::kDebug + 2, "Outside PT acceptance");
      return kFALSE;
   }

   // if we are here, all cuts were passed and no exit point was got
   AliDebug(AliLog::kDebug + 2, "============================= ACCEPTED TRACK =====================================================");
   return kTRUE;
}

//_________________________________________________________________________________________________
void AliRsnCutTrackQuality::Print(const Option_t *) const
{
//
// Print information on this cut
//

   AliInfo(Form("Cut name                : %s", GetName()));
   AliInfo(Form("Required flags (off, on): %lx %lx", fFlagsOn, fFlagsOff));
   AliInfo(Form("Ranges in eta, pt       : %.2f - %.2f, %.2f - %.2f", fEta[0], fEta[1], fPt[0], fPt[1]));
   AliInfo(Form("Kink daughters are      : %s", (fRejectKinkDaughters ? "rejected" : "accepted")));
   AliInfo(Form("TPC requirements        : min. cluster = %d, max chi2 = %f", fTPCminNClusters, fTPCmaxChi2));
   AliInfo(Form("ITS requirements        : min. cluster = %d (all), %d (SPD), max chi2 = %f", fITSminNClusters, fSPDminNClusters, fITSmaxChi2));

   if (fDCARfixed) {
      AliInfo(Form("DCA r cut               : fixed to %f cm", fDCARmax));
   } else {
      AliInfo(Form("DCA r cut formula       : %s", fDCARptFormula.Data()));
   }

   if (fDCAZfixed) {
      AliInfo(Form("DCA z cut               : fixed to %f cm", fDCAZmax));
   } else {
      AliInfo(Form("DCA z cut formula       : %s", fDCAZptFormula.Data()));
   }
}
//__________________________________________________________________________________________________
void AliRsnCutTrackQuality::SetDefaults2010()
{
//
// Default settings for cuts used in 2010
//
   AddStatusFlag(AliESDtrack::kTPCin   , kTRUE);
   AddStatusFlag(AliESDtrack::kTPCrefit, kTRUE);
   AddStatusFlag(AliESDtrack::kITSrefit, kTRUE);
   SetPtRange(0.15, 1E+20);
   SetEtaRange(-0.8, 0.8);
   SetDCARPtFormula("0.0182+0.0350/pt^1.01");
   SetDCAZmax(2.0);
   SetSPDminNClusters(1);
   SetITSminNClusters(0);
   SetITSmaxChi2(36);
   SetMaxChi2TPCConstrainedGlobal(36);
   SetTPCminNClusters(70);
   SetTPCmaxChi2(4.0);
   SetRejectKinkDaughters();
   SetAODTestFilterBit(5);
}

//__________________________________________________________________________________________________
const char *AliRsnCutTrackQuality::Binary(UInt_t number)
{
//
// Convert an integer in binary
//

   static char b[50];
   b[0] = '\0';

   UInt_t z;
   for (z = 512; z > 0; z >>= 1)
      strncat(b, ((number & z) == z) ? "1" : "0", 1);

   return b;
}
Commit	Line	Data
35e49ca5	1	//
	2	// Class AliRsnCutTrackQuality
	3	//
	4	// General implementation of a single cut strategy, which can be:
	5	// - a value contained in a given interval [--> IsBetween() ]
	6	// - a value equal to a given reference [--> MatchesValue()]
	7	//
	8	// In all cases, the reference value(s) is (are) given as data members
	9	// and each kind of cut requires a given value type (Int, UInt, Double),
	10	// but the cut check procedure is then automatized and chosen thanks to
	11	// an enumeration of the implemented cut types.
	12	// At the end, the user (or any other point which uses this object) has
	13	// to use the method IsSelected() to check if this cut has been passed.
	14	//
	15	// authors: Martin Vala (martin.vala@cern.ch)
	16	// Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
	17	//
	18
	19	#include <Riostream.h>
	20	#include <TFormula.h>
	21	#include <TBits.h>
	22
	23	#include "AliLog.h"
35e49ca5	24	#include "AliESDtrackCuts.h"
	25
	26	#include "AliRsnEvent.h"
	27	#include "AliRsnDaughter.h"
	28	#include "AliRsnCutTrackQuality.h"
	29
	30	ClassImp(AliRsnCutTrackQuality)
	31
	32	//_________________________________________________________________________________________________
	33	AliRsnCutTrackQuality::AliRsnCutTrackQuality(const char *name) :
f34f960b	34	AliRsnCut(name, AliRsnTarget::kDaughter, 0.0, 0.0),
2a1c7696	35	fFlagsOn(0x0),
	36	fFlagsOff(0x0),
	37	fRejectKinkDaughters(kTRUE),
	38	fDCARfixed(kTRUE),
	39	fDCARptFormula(""),
d7712d44	40	fDCARmax(1E20),
2a1c7696	41	fDCAZfixed(kTRUE),
2a1c7696	42	fDCAZptFormula(""),
d7712d44	43	fDCAZmax(1E20),
2a1c7696	44	fSPDminNClusters(0),
2a1c7696	45	fITSminNClusters(0),
d7712d44	46	fITSmaxChi2(1E20),
2a1c7696	47	fTPCminNClusters(0),
f34f960b	48	fTPCmaxChi2(1E20),
74d60285	49	fCutMaxChi2TPCConstrainedVsGlobal(1E20),
	50	fAODTestFilterBit(-1),
	51	fESDtrackCuts(0x0)
35e49ca5	52	{
	53	//
	54	// Default constructor.
	55	// Initializes all cuts in such a way that all of them are disabled.
	56	//
	57
d7712d44	58	SetPtRange(0.0, 1E20);
d7712d44	59	SetEtaRange(-1E20, 1E20);
35e49ca5	60	}
	61
	62	//_________________________________________________________________________________________________
	63	AliRsnCutTrackQuality::AliRsnCutTrackQuality(const AliRsnCutTrackQuality &copy) :
2a1c7696	64	AliRsnCut(copy),
	65	fFlagsOn(copy.fFlagsOn),
	66	fFlagsOff(copy.fFlagsOff),
	67	fRejectKinkDaughters(copy.fRejectKinkDaughters),
	68	fDCARfixed(copy.fDCARfixed),
	69	fDCARptFormula(copy.fDCARptFormula),
	70	fDCARmax(copy.fDCARmax),
	71	fDCAZfixed(copy.fDCAZfixed),
	72	fDCAZptFormula(copy.fDCAZptFormula),
	73	fDCAZmax(copy.fDCAZmax),
	74	fSPDminNClusters(copy.fSPDminNClusters),
	75	fITSminNClusters(copy.fITSminNClusters),
	76	fITSmaxChi2(copy.fITSmaxChi2),
	77	fTPCminNClusters(copy.fTPCminNClusters),
f34f960b	78	fTPCmaxChi2(copy.fTPCmaxChi2),
74d60285	79	fCutMaxChi2TPCConstrainedVsGlobal(copy.fCutMaxChi2TPCConstrainedVsGlobal),
5f17a44d	80	fAODTestFilterBit(copy.fAODTestFilterBit),
5f17a44d	81	fESDtrackCuts(copy.fESDtrackCuts)
35e49ca5	82	{
	83	//
	84	// Copy constructor.
	85	// Just copy all data member values.
	86	//
	87
2a1c7696	88	SetPtRange(copy.fPt[0], copy.fPt[1]);
2a1c7696	89	SetEtaRange(copy.fEta[0], copy.fEta[1]);
35e49ca5	90	}
	91
	92	//_________________________________________________________________________________________________
61f275d1	93	AliRsnCutTrackQuality &AliRsnCutTrackQuality::operator=(const AliRsnCutTrackQuality &copy)
35e49ca5	94	{
	95	//
	96	// Assignment operator.
	97	// Just copy all data member values.
	98	//
	99
61f275d1	100	if (this == &copy)
	101	return *this;
	102
2a1c7696	103	fFlagsOn = copy.fFlagsOn;
	104	fFlagsOff = copy.fFlagsOff;
	105	fRejectKinkDaughters = copy.fRejectKinkDaughters;
	106	fDCARfixed = copy.fDCARfixed;
	107	fDCARptFormula = copy.fDCARptFormula;
	108	fDCARmax = copy.fDCARmax;
	109	fDCAZfixed = copy.fDCAZfixed;
	110	fDCAZptFormula = copy.fDCAZptFormula;
	111	fDCAZmax = copy.fDCAZmax;
	112	fSPDminNClusters = copy.fSPDminNClusters;
	113	fITSminNClusters = copy.fITSminNClusters;
	114	fITSmaxChi2 = copy.fITSmaxChi2;
	115	fTPCminNClusters = copy.fTPCminNClusters;
	116	fTPCmaxChi2 = copy.fTPCmaxChi2;
f34f960b	117	fAODTestFilterBit = copy.fAODTestFilterBit;
5f17a44d	118	fESDtrackCuts = copy.fESDtrackCuts;
2a1c7696	119	SetPtRange(copy.fPt[0], copy.fPt[1]);
	120	SetEtaRange(copy.fEta[0], copy.fEta[1]);
	121
	122	return (*this);
35e49ca5	123	}
	124
	125	//_________________________________________________________________________________________________
	126	void AliRsnCutTrackQuality::DisableAll()
	127	{
	128	//
	129	// Disable all cuts
	130	//
2a1c7696	131
	132	fFlagsOn = 0x0;
	133	fFlagsOff = 0x0;
	134	fRejectKinkDaughters = kFALSE;
	135	fDCARfixed = kTRUE;
	136	fDCARptFormula = "";
d7712d44	137	fDCARmax = 1E20;
2a1c7696	138	fDCAZfixed = kTRUE;
2a1c7696	139	fDCAZptFormula = "";
d7712d44	140	fDCAZmax = 1E20;
2a1c7696	141	fSPDminNClusters = 0;
2a1c7696	142	fITSminNClusters = 0;
d7712d44	143	fITSmaxChi2 = 1E20;
2a1c7696	144	fTPCminNClusters = 0;
d7712d44	145	fTPCmaxChi2 = 1E20;
f34f960b	146	fAODTestFilterBit = -1;
5f17a44d	147	if (fESDtrackCuts) {
74d60285	148	const char *cutsName = fESDtrackCuts->GetName();
	149	const char *cutsTitle = fESDtrackCuts->GetTitle();
	150	delete fESDtrackCuts;
	151	fESDtrackCuts = new AliESDtrackCuts(cutsName,cutsTitle);
5f17a44d	152	}
d7712d44	153	SetPtRange(0.0, 1E20);
d7712d44	154	SetEtaRange(-1E20, 1E20);
35e49ca5	155	}
	156
	157	//_________________________________________________________________________________________________
	158	Bool_t AliRsnCutTrackQuality::IsSelected(TObject *object)
	159	{
	160	//
	161	// Cut checker.
	162	// Checks the type of object being evaluated
	163	// and then calls the appropriate sub-function (for ESD or AOD)
	164	//
	165
2a1c7696	166	// coherence check
	167	if (!TargetOK(object)) return kFALSE;
	168
	169	// status is checked in the same way for all tracks, using AliVTrack
	170	// as a convention, if a the collection of 'on' flags is '0x0', it
	171	// is assumed that no flags are required, and this check is skipped;
	172	// for the collection of 'off' flags this is not needed
f34f960b	173	AliVTrack *vtrack = fDaughter->Ref2Vtrack();
2a1c7696	174	if (!vtrack) {
f34f960b	175	AliDebug(AliLog::kDebug + 2, "This object is not either an ESD nor AOD track");
2a1c7696	176	return kFALSE;
	177	}
	178	ULong_t status = (ULong_t)vtrack->GetStatus();
	179	ULong_t checkOn = status & fFlagsOn;
	180	ULong_t checkOff = status & fFlagsOff;
	181	if (fFlagsOn != 0x0 && checkOn != fFlagsOn) {
f34f960b	182	AliDebug(AliLog::kDebug + 2, Form("Failed flag check: required %s", Binary(fFlagsOn)));
f34f960b	183	AliDebug(AliLog::kDebug + 2, Form(" track has %s", Binary(status )));
2a1c7696	184	return kFALSE;
	185	}
	186	if (checkOff != 0) {
f34f960b	187	AliDebug(AliLog::kDebug + 2, Form("Failed flag check: forbidden %s", Binary(fFlagsOff)));
f34f960b	188	AliDebug(AliLog::kDebug + 2, Form(" track has %s", Binary(status )));
2a1c7696	189	return kFALSE;
2a1c7696	190	}
f34f960b	191	AliDebug(AliLog::kDebug + 3, Form("Flag check OK: required %s", Binary(fFlagsOn)));
	192	AliDebug(AliLog::kDebug + 3, Form(" forbidden %s", Binary(fFlagsOff)));
	193	AliDebug(AliLog::kDebug + 3, Form(" track has %s", Binary(status )));
2a1c7696	194
2a1c7696	195	// retrieve real object type
f34f960b	196	AliESDtrack *esdTrack = fDaughter->Ref2ESDtrack();
f34f960b	197	AliAODTrack *aodTrack = fDaughter->Ref2AODtrack();
2a1c7696	198	if (esdTrack) {
2a1c7696	199	AliDebug(AliLog::kDebug + 2, "Checking an ESD track");
74d60285	200	if (fESDtrackCuts)
	201	return fESDtrackCuts->IsSelected(esdTrack);
	202	else
	203	return CheckESD(esdTrack);
2a1c7696	204	} else if (aodTrack) {
	205	AliDebug(AliLog::kDebug + 2, "Checking an AOD track");
	206	return CheckAOD(aodTrack);
	207	} else {
	208	AliDebug(AliLog::kDebug + 2, Form("This object is not either an ESD nor AOD track, it is an %s", fDaughter->GetRef()->ClassName()));
	209	return kFALSE;
	210	}
35e49ca5	211	}
	212
	213	//_________________________________________________________________________________________________
	214	Bool_t AliRsnCutTrackQuality::CheckESD(AliESDtrack *track)
	215	{
	216	//
	217	// Check an ESD track.
	218	// This is done using the default track checker for ESD.
	219	// It is declared static, not to recreate it every time.
	220	//
	221
2a1c7696	222	static AliESDtrackCuts cuts;
	223
	224	// general acceptance/pt cuts
	225	cuts.SetPtRange(fPt[0], fPt[1]);
	226	cuts.SetEtaRange(fEta[0], fEta[1]);
	227
	228	// transverse DCA cuts
	229	if (fDCARfixed)
	230	cuts.SetMaxDCAToVertexXY(fDCARmax);
	231	else
	232	cuts.SetMaxDCAToVertexXYPtDep(fDCARptFormula.Data());
	233
	234	// longitudinal DCA cuts
	235	if (fDCAZfixed)
	236	cuts.SetMaxDCAToVertexZ(fDCAZmax);
	237	else
	238	cuts.SetMaxDCAToVertexZPtDep(fDCAZptFormula.Data());
	239
	240	// these options are always disabled in current version
	241	cuts.SetDCAToVertex2D(kFALSE);
	242	cuts.SetRequireSigmaToVertex(kFALSE);
	243
	244	// TPC related cuts for TPC+ITS tracks
	245	cuts.SetMinNClustersTPC(fTPCminNClusters);
	246	cuts.SetMaxChi2PerClusterTPC(fTPCmaxChi2);
	247	cuts.SetAcceptKinkDaughters(!fRejectKinkDaughters);
5f17a44d	248	cuts.SetMaxChi2TPCConstrainedGlobal(fCutMaxChi2TPCConstrainedVsGlobal);
74d60285	249
2a1c7696	250	// ITS related cuts for TPC+ITS tracks
2a1c7696	251	if (fSPDminNClusters > 0)
74d60285	252	cuts.SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
5f17a44d	253	cuts.SetMaxChi2PerClusterITS(fITSmaxChi2);
2a1c7696	254
	255	// now that all is initialized, do the check
	256	return cuts.IsSelected(track);
35e49ca5	257	}
	258
	259	//_________________________________________________________________________________________________
	260	Bool_t AliRsnCutTrackQuality::CheckAOD(AliAODTrack *track)
	261	{
	262	//
	263	// Check an AOD track.
	264	// This is done doing directly all checks, since there is not
	265	// an equivalend checker for AOD tracks
	266	//
2a1c7696	267
f34f960b	268	// if a test bit is used, check it and skip the following
f34f960b	269	if (fAODTestFilterBit >= 0) {
e5ebff7e	270	UInt_t bit = 1 << fAODTestFilterBit;
f34f960b	271	AliDebugClass(2, Form("Required a test filter bit for AOD check: %u (result: %s)", bit, (track->TestFilterBit(bit) ? "accept" : "reject")));
61f275d1	272	if (!track->TestFilterBit(bit))
f34f960b	273	return kFALSE;
	274	else {
	275	if (track->Pt() < fPt[0] \|\| track->Pt() > fPt[1]) return kFALSE;
	276	if (track->Eta() < fEta[0] \|\| track->Eta() > fEta[1]) return kFALSE;
	277	return kTRUE;
	278	}
	279	}
61f275d1	280
d7712d44	281	// try to retrieve the reference AOD event
	282	AliAODEvent *aodEvent = 0x0;
	283	if (fEvent) aodEvent = fEvent->GetRefAOD();
	284	if (!aodEvent) {
	285	AliError("AOD reference event is not initialized!");
	286	return kFALSE;
	287	}
	288
2a1c7696	289	// step #0: check SPD and ITS clusters
	290	Int_t nSPD = 0;
	291	nSPD = TESTBIT(track->GetITSClusterMap(), 0);
	292	nSPD += TESTBIT(track->GetITSClusterMap(), 1);
	293	if (nSPD < fSPDminNClusters) {
	294	AliDebug(AliLog::kDebug + 2, "Not enough SPD clusters in this track. Rejected");
35e49ca5	295	return kFALSE;
2a1c7696	296	}
2a1c7696	297
5f17a44d	298
2a1c7696	299	// step #1: check number of clusters in TPC
	300	if (track->GetTPCNcls() < fTPCminNClusters) {
	301	AliDebug(AliLog::kDebug + 2, "Too few TPC clusters. Rejected");
	302	return kFALSE;
	303	}
	304	if (track->GetITSNcls() < fITSminNClusters) {
	305	AliDebug(AliLog::kDebug + 2, "Too few ITS clusters. Rejected");
	306	return kFALSE;
	307	}
	308
	309	// step #2: check chi square
	310	if (track->Chi2perNDF() > fTPCmaxChi2) {
	311	AliDebug(AliLog::kDebug + 2, "Bad chi2. Rejected");
	312	return kFALSE;
	313	}
	314	if (track->Chi2perNDF() > fITSmaxChi2) {
	315	AliDebug(AliLog::kDebug + 2, "Bad chi2. Rejected");
	316	return kFALSE;
	317	}
	318
	319	// step #3: reject kink daughters
	320	AliAODVertex *vertex = track->GetProdVertex();
	321	if (vertex && fRejectKinkDaughters) {
	322	if (vertex->GetType() == AliAODVertex::kKink) {
	323	AliDebug(AliLog::kDebug + 2, "Kink daughter. Rejected");
	324	return kFALSE;
	325	}
	326	}
	327
	328	// step #4: DCA cut (transverse)
c865cb1d	329	// --> reject all tracks not ITS refitted
2a1c7696	330	Double_t b[2], cov[3];
d7712d44	331	vertex = aodEvent->GetPrimaryVertex();
2a1c7696	332	if (!vertex) {
	333	AliDebug(AliLog::kDebug + 2, "NULL vertex");
	334	return kFALSE;
	335	}
c865cb1d	336	if ((track->GetStatus() & AliESDtrack::kITSrefit) == 0) {
	337	AliDebug(AliLog::kDebug + 2, "Not ITS refitted");
	338	return kFALSE;
	339	}
d7712d44	340	if (!track->PropagateToDCA(vertex, aodEvent->GetMagneticField(), kVeryBig, b, cov)) {
2a1c7696	341	AliDebug(AliLog::kDebug + 2, "Failed propagation to vertex");
	342	return kFALSE;
	343	}
	344	// if the DCA cut is not fixed, compute current value
	345	if (!fDCARfixed) {
	346	static TString str(fDCARptFormula);
	347	str.ReplaceAll("pt", "x");
	348	static const TFormula dcaXY(Form("%s_dcaXY", GetName()), str.Data());
	349	fDCARmax = dcaXY.Eval(track->Pt());
	350	}
	351	// check the cut
	352	if (TMath::Abs(b[0]) > fDCARmax) {
	353	AliDebug(AliLog::kDebug + 2, "Too large transverse DCA");
	354	return kFALSE;
	355	}
	356
	357	// step #5: DCA cut (longitudinal)
	358	// the DCA has already been computed above
	359	// if the DCA cut is not fixed, compute current value
	360	if (!fDCAZfixed) {
	361	static TString str(fDCAZptFormula);
	362	str.ReplaceAll("pt", "x");
	363	static const TFormula dcaZ(Form("%s_dcaXY", GetName()), str.Data());
	364	fDCAZmax = dcaZ.Eval(track->Pt());
	365	}
	366	// check the cut
	367	if (TMath::Abs(b[1]) > fDCAZmax) {
	368	AliDebug(AliLog::kDebug + 2, "Too large longitudinal DCA");
	369	return kFALSE;
	370	}
	371
	372	// step #6: check eta/pt range
	373	if (track->Eta() < fEta[0] \|\| track->Eta() > fEta[1]) {
	374	AliDebug(AliLog::kDebug + 2, "Outside ETA acceptance");
	375	return kFALSE;
	376	}
	377	if (track->Pt() < fPt[0] \|\| track->Pt() > fPt[1]) {
	378	AliDebug(AliLog::kDebug + 2, "Outside PT acceptance");
	379	return kFALSE;
	380	}
	381
	382	// if we are here, all cuts were passed and no exit point was got
c865cb1d	383	AliDebug(AliLog::kDebug + 2, "============================= ACCEPTED TRACK =====================================================");
2a1c7696	384	return kTRUE;
35e49ca5	385	}
a909ffad	386
	387	//_________________________________________________________________________________________________
	388	void AliRsnCutTrackQuality::Print(const Option_t *) const
	389	{
	390	//
	391	// Print information on this cut
	392	//
	393
2a1c7696	394	AliInfo(Form("Cut name : %s", GetName()));
	395	AliInfo(Form("Required flags (off, on): %lx %lx", fFlagsOn, fFlagsOff));
	396	AliInfo(Form("Ranges in eta, pt : %.2f - %.2f, %.2f - %.2f", fEta[0], fEta[1], fPt[0], fPt[1]));
	397	AliInfo(Form("Kink daughters are : %s", (fRejectKinkDaughters ? "rejected" : "accepted")));
	398	AliInfo(Form("TPC requirements : min. cluster = %d, max chi2 = %f", fTPCminNClusters, fTPCmaxChi2));
	399	AliInfo(Form("ITS requirements : min. cluster = %d (all), %d (SPD), max chi2 = %f", fITSminNClusters, fSPDminNClusters, fITSmaxChi2));
	400
	401	if (fDCARfixed) {
	402	AliInfo(Form("DCA r cut : fixed to %f cm", fDCARmax));
	403	} else {
	404	AliInfo(Form("DCA r cut formula : %s", fDCARptFormula.Data()));
	405	}
	406
	407	if (fDCAZfixed) {
	408	AliInfo(Form("DCA z cut : fixed to %f cm", fDCAZmax));
	409	} else {
	410	AliInfo(Form("DCA z cut formula : %s", fDCAZptFormula.Data()));
	411	}
a909ffad	412	}
5f17a44d	413	//__________________________________________________________________________________________________
	414	void AliRsnCutTrackQuality::SetDefaults2010()
	415	{
	416	//
	417	// Default settings for cuts used in 2010
	418	//
	419	AddStatusFlag(AliESDtrack::kTPCin , kTRUE);
	420	AddStatusFlag(AliESDtrack::kTPCrefit, kTRUE);
	421	AddStatusFlag(AliESDtrack::kITSrefit, kTRUE);
	422	SetPtRange(0.15, 1E+20);
	423	SetEtaRange(-0.8, 0.8);
	424	SetDCARPtFormula("0.0182+0.0350/pt^1.01");
	425	SetDCAZmax(2.0);
	426	SetSPDminNClusters(1);
	427	SetITSminNClusters(0);
	428	SetITSmaxChi2(36);
	429	SetMaxChi2TPCConstrainedGlobal(36);
	430	SetTPCminNClusters(70);
	431	SetTPCmaxChi2(4.0);
	432	SetRejectKinkDaughters();
	433	SetAODTestFilterBit(5);
	434	}
	435
	436	//__________________________________________________________________________________________________
	437	const char *AliRsnCutTrackQuality::Binary(UInt_t number)
	438	{
	439	//
	440	// Convert an integer in binary
	441	//
	442
	443	static char b[50];
	444	b[0] = '\0';
	445
	446	UInt_t z;
	447	for (z = 512; z > 0; z >>= 1)
	448	strncat(b, ((number & z) == z) ? "1" : "0", 1);
	449
	450	return b;
	451	}