[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),
   fCheckOnlyFilterBit(kTRUE),
   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),
   fCheckOnlyFilterBit(copy.fCheckOnlyFilterBit),
   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;
   fCheckOnlyFilterBit = copy.fCheckOnlyFilterBit;
   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;
         if (fCheckOnlyFilterBit) 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) {
      TString str(fDCARptFormula);
      str.ReplaceAll("pt", "x");
      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) {
      TString str(fDCAZptFormula);
      str.ReplaceAll("pt", "x");
      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()));
   }

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