PWGGA/GammaConv/ConvCorrelations/AliConversionTrackCuts.cxx

   1 /**************************************************************************
   2  * This file is property of and copyright by the ALICE HLT Project        *
   3  * ALICE Experiment at CERN, All rights reserved.                         *
   4  *                                                                        *
   5  * Primary Author: Svein Lindal <slindal@fys.uio.no>                      *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /// @file   AliConversionTrackCuts.cxx
  17 /// @author Svein Lindal
  18 /// @brief  Base class for analysation of conversion particle - track correlations
  19
  20
  21 #include "AliConversionTrackCuts.h"
  22 //#include "AliAODTrack.h"
  23 #include "AliAODEvent.h"
  24 #include <TFormula.h>
  25 #include <iostream>
  26 #include "TH2F.h"
  27
  28
  29 using namespace std;
  30 ClassImp(AliConversionTrackCuts)
  31
  32 //________________________________________________________________________
  33 AliConversionTrackCuts::AliConversionTrackCuts() :
  34 AliAnalysisCuts(),
  35   fFlagsOn(0x0),
  36   fFlagsOff(0x0),
  37   fRejectKinkDaughters(kTRUE),
  38   fDCARfixed(kTRUE),
  39   fDCARptFormula(""),
  40   fDCARmax(1E20),
  41   fDCAZfixed(kTRUE),
  42   fDCAZptFormula(""),
  43   fDCAZmax(1E20),
  44   fSPDminNClusters(0),
  45   fITSminNClusters(0),
  46   fITSmaxChi2(1E20),
  47   fTPCminNClusters(0),
  48   fTPCmaxChi2(1E20),
  49   fAODTestFilterBit(-1),
  50   fhPhi(NULL),
  51   fHistograms(NULL)
  52 {
  53   //Constructor
  54 }
  55 //________________________________________________________________________
  56 AliConversionTrackCuts::AliConversionTrackCuts(TString name, TString title = "title") :
  57   AliAnalysisCuts(name, title),
  58   fFlagsOn(0x0),
  59   fFlagsOff(0x0),
  60   fRejectKinkDaughters(kTRUE),
  61   fDCARfixed(kTRUE),
  62   fDCARptFormula(""),
  63   fDCARmax(1E20),
  64   fDCAZfixed(kTRUE),
  65   fDCAZptFormula(""),
  66   fDCAZmax(1E20),
  67   fSPDminNClusters(0),
  68   fITSminNClusters(0),
  69   fITSmaxChi2(1E20),
  70   fTPCminNClusters(0),
  71   fTPCmaxChi2(1E20),
  72   fAODTestFilterBit(-1),
  73   fhPhi(NULL),
  74   fHistograms(NULL)
  75 {
  76   //Constructor
  77 }
  78
  79
  80 //________________________________________________________________________________
  81  AliConversionTrackCuts::~AliConversionTrackCuts() {
  82    ///destructor
  83    // if(fHistograms)
  84    //    delete fHistograms;
  85    // fHistograms = NULL;
  86
  87 }
  88
  89 TList * AliConversionTrackCuts::CreateHistograms() {
  90
  91   if(!fHistograms) fHistograms = new TList();
  92
  93   fHistograms->SetOwner(kTRUE);
  94   fHistograms->SetName("trackCuts");
  95
  96   fhPhi = new TH2F("phi", "phi", 20, -0.5, 19.5, 128, 0, TMath::TwoPi());
  97   fHistograms->Add(fhPhi);
  98
  99   fhPhiPt = new TH2F("phipt", "phipt", 80, 0, 100, 128, 0, TMath::TwoPi());
 100   fHistograms->Add(fhPhiPt);
 101
 102
 103   return fHistograms;
 104 }
 105
 106
 107 void AliConversionTrackCuts::FillHistograms(Int_t cutIndex, AliVTrack * track) {
 108
 109   fhPhi->Fill(cutIndex, track->Phi());
 110   fhPhiPt->Fill(track->Pt(), track->Phi());
 111 }
 112
 113 Bool_t AliConversionTrackCuts::AcceptTrack(AliAODTrack * track, AliAODEvent * aodEvent) {
 114   // Check an AOD track.
 115 // This is done doing directly all checks, since there is not
 116 // an equivalend checker for AOD tracks
 117 //
 118
 119    // try to retrieve the reference AOD event
 120    // AliAODEvent *aodEvent = 0x0;
 121    // if (fEvent) aodEvent = fEvent->GetRefAOD();
 122    // if (!aodEvent) {
 123    //    AliError("AOD reference event is not initialized!");
 124    //    return kFALSE;
 125    // }
 126
 127   Int_t cutIndex = 0;
 128
 129   FillHistograms(cutIndex, track);
 130   cutIndex++;
 131   // step #0: check SPD and ITS clusters
 132   Int_t nSPD = 0;
 133   nSPD  = TESTBIT(track->GetITSClusterMap(), 0);
 134   nSPD += TESTBIT(track->GetITSClusterMap(), 1);
 135   if (nSPD < fSPDminNClusters) {
 136   FillHistograms(cutIndex, track);
 137         AliDebug(AliLog::kDebug + 2, "Not enough SPD clusters in this track. Rejected");
 138         return kFALSE;
 139   }
 140   cutIndex++;
 141
 142   // step #1: check number of clusters in TPC
 143   if (track->GetTPCNcls() < fTPCminNClusters) {
 144   FillHistograms(cutIndex, track);
 145         AliDebug(AliLog::kDebug + 2, "Too few TPC clusters. Rejected");
 146         return kFALSE;
 147   }
 148   cutIndex++;
 149
 150   if (track->GetITSNcls() < fITSminNClusters) {
 151   FillHistograms(cutIndex, track);
 152         AliDebug(AliLog::kDebug + 2, "Too few ITS clusters. Rejected");
 153         return kFALSE;
 154   }
 155   cutIndex++;
 156
 157   // step #2: check chi square
 158   if (track->Chi2perNDF() > fTPCmaxChi2) {
 159   FillHistograms(cutIndex, track);
 160         AliDebug(AliLog::kDebug + 2, "Bad chi2. Rejected");
 161         return kFALSE;
 162   }
 163   cutIndex++;
 164
 165   if (track->Chi2perNDF() > fITSmaxChi2) {
 166   FillHistograms(cutIndex, track);
 167         AliDebug(AliLog::kDebug + 2, "Bad chi2. Rejected");
 168         return kFALSE;
 169   }
 170   cutIndex++;
 171
 172
 173    // step #3: reject kink daughters
 174    AliAODVertex *vertex = track->GetProdVertex();
 175    if (vertex && fRejectKinkDaughters) {
 176       if (vertex->GetType() == AliAODVertex::kKink) {
 177   FillHistograms(cutIndex, track);
 178          AliDebug(AliLog::kDebug + 2, "Kink daughter. Rejected");
 179          return kFALSE;
 180       }
 181    }
 182   cutIndex++;
 183
 184
 185    // step #4: DCA cut (transverse)
 186    Double_t b[2], cov[3];
 187    vertex = aodEvent->GetPrimaryVertex();
 188    if (!vertex) {
 189   FillHistograms(cutIndex, track);
 190       AliDebug(AliLog::kDebug + 2, "NULL vertex");
 191       return kFALSE;
 192    }
 193   cutIndex++;
 194
 195    if (!track->PropagateToDCA(vertex, aodEvent->GetMagneticField(), kVeryBig, b, cov)) {
 196       AliDebug(AliLog::kDebug + 2, "Failed propagation to vertex");
 197   FillHistograms(cutIndex, track);
 198       return kFALSE;
 199    }
 200   cutIndex++;
 201
 202    // if the DCA cut is not fixed, compute current value
 203    if (!fDCARfixed) {
 204          FillHistograms(cutIndex, track);
 205       static TString str(fDCARptFormula);
 206       str.ReplaceAll("pt", "x");
 207       static const TFormula dcaXY(Form("%s_dcaXY", GetName()), str.Data());
 208       fDCARmax = dcaXY.Eval(track->Pt());
 209    }
 210   cutIndex++;
 211
 212    // check the cut
 213    if (TMath::Abs(b[0]) > fDCARmax) {
 214        FillHistograms(cutIndex, track);
 215            AliDebug(AliLog::kDebug + 2, "Too large transverse DCA");
 216       return kFALSE;
 217    }
 218   cutIndex++;
 219
 220
 221    // step #5: DCA cut (longitudinal)
 222    // the DCA has already been computed above
 223    // if the DCA cut is not fixed, compute current value
 224    if (!fDCAZfixed) {
 225         FillHistograms(cutIndex, track);
 226                 static TString str(fDCAZptFormula);
 227       str.ReplaceAll("pt", "x");
 228       static const TFormula dcaZ(Form("%s_dcaXY", GetName()), str.Data());
 229       fDCAZmax = dcaZ.Eval(track->Pt());
 230    }
 231   cutIndex++;
 232
 233    // check the cut
 234   if (TMath::Abs(b[1]) > fDCAZmax) {
 235        FillHistograms(cutIndex, track);
 236            AliDebug(AliLog::kDebug + 2, "Too large longitudinal DCA");
 237       return kFALSE;
 238    }
 239
 240   cutIndex++;
 241
 242    // step #6: check eta/pt range
 243    if (track->Eta() < fEta[0] || track->Eta() > fEta[1]) {
 244       FillHistograms(cutIndex, track);
 245           AliDebug(AliLog::kDebug + 2, "Outside ETA acceptance");
 246       return kFALSE;
 247    }
 248   cutIndex++;
 249
 250    // if (track->Pt() < fPt[0] || track->Pt() > fPt[1]) {
 251    //    AliDebug(AliLog::kDebug + 2, "Outside PT acceptance");
 252    //    return kFALSE;
 253    // }
 254
 255    // if we are here, all cuts were passed and no exit point was got
 256    return kTRUE;
 257 }
 258
 259 //_________________________________________________________________________________________________
 260 void AliConversionTrackCuts::Print(const Option_t *) const
 261 {
 262 //
 263 // Print information on this cut
 264 //
 265
 266    AliInfo(Form("Cut name                : %s", GetName()));
 267    AliInfo(Form("Required flags (off, on): %lx %lx", fFlagsOn, fFlagsOff));
 268    AliInfo(Form("Ranges in eta, pt       : %.2f - %.2f, %.2f - %.2f", fEta[0], fEta[1], fPt[0], fPt[1]));
 269    AliInfo(Form("Kink daughters are      : %s", (fRejectKinkDaughters ? "rejected" : "accepted")));
 270    AliInfo(Form("TPC requirements        : min. cluster = %d, max chi2 = %f", fTPCminNClusters, fTPCmaxChi2));
 271    AliInfo(Form("ITS requirements        : min. cluster = %d (all), %d (SPD), max chi2 = %f", fITSminNClusters, fSPDminNClusters, fITSmaxChi2));
 272
 273    if (fDCARfixed) {
 274          AliInfo(Form("DCA r cut               : fixed to %f cm", fDCARmax));
 275    } else {
 276          AliInfo(Form("DCA r cut formula       : %s", fDCARptFormula.Data()));
 277    }
 278
 279    if (fDCAZfixed) {
 280          AliInfo(Form("DCA z cut               : fixed to %f cm", fDCAZmax));
 281    } else {
 282          AliInfo(Form("DCA z cut formula       : %s", fDCAZptFormula.Data()));
 283    }
 284
 285
 286 }
 287
 288
 289
 290
 291
 292