PWG2/RESONANCES/AliRsnCutESD2010.cxx

   1 //
   2 // Class AliRsnCutESD2010
   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
  21 #include "AliESDpid.h"
  22 #include "AliTOFT0maker.h"
  23 #include "AliTOFcalib.h"
  24 #include "AliCDBManager.h"
  25 #include "AliITSPIDResponse.h"
  26
  27 #include "AliRsnEvent.h"
  28 #include "AliRsnDaughter.h"
  29 #include "AliRsnCutESD2010.h"
  30
  31 ClassImp(AliRsnCutESD2010)
  32
  33 //_________________________________________________________________________________________________
  34 AliRsnCutESD2010::AliRsnCutESD2010() :
  35   AliRsnCut(AliRsnCut::kDaughter),
  36   fIsMC(kFALSE),
  37   fCheckITS(kTRUE),
  38   fCheckTPC(kTRUE),
  39   fCheckTOF(kTRUE),
  40   fUseGlobal(kTRUE),
  41   fUseITSSA(kTRUE),
  42   fMaxITSband(1E6),
  43   fTPCpLimit(0.35),
  44   fMinTPCband(-1E6),
  45   fMaxTPCband( 1E6),
  46   fESDtrackCutsTPC(),
  47   fESDtrackCutsITS(),
  48   fESDpid(0x0),
  49   fTOFmaker(0x0),
  50   fTOFcalib(0x0),
  51   fTOFcalibrateESD(kFALSE),
  52   fTOFcorrectTExp(kFALSE),
  53   fTOFuseT0(kFALSE),
  54   fTOFtuneMC(kFALSE),
  55   fTOFresolution(0.0),
  56   fLastRun(-1)
  57 {
  58 //
  59 // Default constructor.
  60 //
  61 }
  62
  63 //_________________________________________________________________________________________________
  64 AliRsnCutESD2010::AliRsnCutESD2010
  65 (const char *name) :
  66   AliRsnCut(name, AliRsnCut::kDaughter, 0.0, 0.0),
  67   fIsMC(kFALSE),
  68   fCheckITS(kTRUE),
  69   fCheckTPC(kTRUE),
  70   fCheckTOF(kTRUE),
  71   fUseGlobal(kTRUE),
  72   fUseITSSA(kTRUE),
  73   fMaxITSband(1E6),
  74   fTPCpLimit(0.35),
  75   fMinTPCband(-1E6),
  76   fMaxTPCband( 1E6),
  77   fESDtrackCutsTPC(),
  78   fESDtrackCutsITS(),
  79   fESDpid(0x0),
  80   fTOFmaker(0x0),
  81   fTOFcalib(0x0),
  82   fTOFcalibrateESD(kFALSE),
  83   fTOFcorrectTExp(kFALSE),
  84   fTOFuseT0(kFALSE),
  85   fTOFtuneMC(kFALSE),
  86   fTOFresolution(0.0),
  87   fLastRun(-1)
  88 {
  89 //
  90 // Main constructor.
  91 //
  92 }
  93
  94 //_________________________________________________________________________________________________
  95 AliRsnCutESD2010::AliRsnCutESD2010
  96 (const AliRsnCutESD2010& copy) :
  97   AliRsnCut(copy),
  98   fIsMC(kFALSE),
  99   fCheckITS(copy.fCheckITS),
 100   fCheckTPC(copy.fCheckTPC),
 101   fCheckTOF(copy.fCheckTOF),
 102   fUseGlobal(copy.fUseGlobal),
 103   fUseITSSA(copy.fUseITSSA),
 104   fMaxITSband(copy.fMaxITSband),
 105   fTPCpLimit(copy.fTPCpLimit),
 106   fMinTPCband(copy.fMinTPCband),
 107   fMaxTPCband(copy.fMaxTPCband),
 108   fESDtrackCutsTPC(copy.fESDtrackCutsTPC),
 109   fESDtrackCutsITS(copy.fESDtrackCutsITS),
 110   fESDpid(0x0),
 111   fTOFmaker(0x0),
 112   fTOFcalib(0x0),
 113   fTOFcalibrateESD(copy.fTOFcalibrateESD),
 114   fTOFcorrectTExp(copy.fTOFcorrectTExp),
 115   fTOFuseT0(copy.fTOFuseT0),
 116   fTOFtuneMC(copy.fTOFtuneMC),
 117   fTOFresolution(copy.fTOFresolution),
 118   fLastRun(-1)
 119 {
 120 //
 121 // Main constructor.
 122 //
 123
 124   Initialize();
 125 }
 126
 127 //_________________________________________________________________________________________________
 128 void AliRsnCutESD2010::Initialize()
 129 {
 130 //
 131 // Main constructor.
 132 //
 133 }
 134
 135 //_________________________________________________________________________________________________
 136 void AliRsnCutESD2010::SetEvent(AliRsnEvent *event)
 137 {
 138   // don't do anything if the event is the same as before
 139   if (fEvent != 0x0 && fEvent == event) return;
 140
 141   // retrieve the ESD event
 142   AliESDEvent *esd = event->GetRefESD();
 143   if (!esd)
 144   {
 145     fEvent = 0x0;
 146     return;
 147   }
 148   else
 149   {
 150     fEvent = event;
 151   }
 152
 153   // if absent, initialize ESD pid responst
 154   if (!fESDpid)
 155   {
 156     fESDpid = new AliESDpid;
 157     fESDpid->GetTPCResponse().SetBetheBlochParameters(fTPCpar[0],fTPCpar[1],fTPCpar[2],fTPCpar[3],fTPCpar[4]);
 158   }
 159
 160   // initialize DB to current run
 161   Int_t run = esd->GetRunNumber();
 162   if (run != fLastRun)
 163   {
 164     cout << "Run = " << run << " -- LAST = " << fLastRun << endl;
 165     fLastRun = run;
 166
 167     // setup TOF maker & calibration
 168     if (!fTOFcalib) fTOFcalib = new AliTOFcalib;
 169     fTOFcalib->SetCorrectTExp(fTOFcorrectTExp);
 170     if (!fTOFmaker) fTOFmaker = new AliTOFT0maker(fESDpid, fTOFcalib);
 171     fTOFmaker->SetTimeResolution(fTOFresolution);
 172
 173     AliCDBManager *cdb = AliCDBManager::Instance();
 174     cdb->SetDefaultStorage("raw://");
 175     cdb->SetRun(run);
 176     fTOFcalib->SetCorrectTExp(fTOFcorrectTExp);
 177     fTOFcalib->Init();
 178   }
 179
 180   // if required, calibrate the TOF t0 maker with current event
 181   if (fTOFcalibrateESD) fTOFcalib->CalibrateESD(esd);
 182   if (fTOFtuneMC) fTOFmaker->TuneForMC(esd);
 183   if (fTOFuseT0)
 184   {
 185     fTOFmaker->ComputeT0TOF(esd);
 186     fTOFmaker->ApplyT0TOF(esd);
 187     fESDpid->MakePID(esd, kFALSE, 0.);
 188   }
 189 }
 190
 191 //_________________________________________________________________________________________________
 192 Bool_t AliRsnCutESD2010::IsSelected(TObject *obj1, TObject* /*obj2*/)
 193 {
 194 //
 195 // Cut checker.
 196 //
 197
 198   // coherence check: require an ESD track
 199   AliRsnDaughter *daughter = dynamic_cast<AliRsnDaughter*>(obj1);
 200   if (!daughter) return kFALSE;
 201   AliESDtrack *track = daughter->GetRefESDtrack();
 202   if (!track) return kFALSE;
 203
 204   // if no reference event, skip
 205   if (!fEvent) return kFALSE;
 206
 207   // ITS stuff #1 create the response function
 208   AliITSPIDResponse itsrsp(fIsMC);
 209
 210   // TOF: define fixed function for compatibility range
 211   Double_t a1 = 0.01, a2 = -0.03;
 212   Double_t b1 = 0.25, b2 =  0.25;
 213   Double_t c1 = 0.05, c2 = -0.03;
 214   Double_t ymax, ymin;
 215
 216   ULong_t  status;
 217   Int_t    k, nITS;
 218   Double_t times[10], tpcNSigma, tpcMaxNSigma, itsSignal, itsNSigma, mom, tofTime, tofSigma, tofRef, tofRel;
 219   Bool_t   okQuality, okTOF, okTPC, okITS, okTrack, isTPC, isITSSA;
 220   UChar_t  itsCluMap;
 221
 222   // get commonly used variables
 223   status  = (ULong_t)track->GetStatus();
 224   mom     = track->P();
 225   isTPC   = ((status & AliESDtrack::kTPCin)  != 0);
 226   isITSSA = ((status & AliESDtrack::kTPCin)  == 0 && (status & AliESDtrack::kITSrefit) != 0 && (status & AliESDtrack::kITSpureSA) == 0 && (status & AliESDtrack::kITSpid) != 0);
 227
 228   // accept only tracks which are TPC+ITS or ITS standalone
 229   if (!isTPC   && !isITSSA) return kFALSE;
 230   if ( isTPC   && !fUseGlobal) return kFALSE;
 231   if ( isITSSA && !fUseITSSA) return kFALSE;
 232
 233   if (isTPC)
 234   {
 235     // check standard ESD cuts
 236     okQuality = fESDtrackCutsTPC.IsSelected(track);
 237     //cout << "GLOBAL -- quality = " << (okQuality ? "GOOD" : "BAD") << endl;
 238     if (!okQuality) return kFALSE;
 239
 240     // check TPC dE/dx
 241     if (fCheckTPC)
 242     {
 243       tpcNSigma = TMath::Abs(fESDpid->NumberOfSigmasTPC(track, AliPID::kKaon));
 244       if (track->GetInnerParam()->P() > fTPCpLimit) tpcMaxNSigma = fMinTPCband; else tpcMaxNSigma = fMaxTPCband;
 245       okTPC = (tpcNSigma <= tpcMaxNSigma);
 246       //cout << "RSN -- TPC    -- nsigma = " << tpcNSigma << ", max = " << tpcMaxNSigma << " --> " << (okTPC ? "OK" : "FAILED") << endl;
 247       //cout << "RSNTPC -- " << fTPCpar[0] << ' ' << fTPCpar[1] << ' ' << fTPCpar[2] << ' ' << fTPCpar[3] << ' ' << fTPCpar[4] << endl;
 248     }
 249     else
 250     {
 251       okTPC = kTRUE;
 252     }
 253
 254     // check TOF (only if momentum is large than function asymptote and flags are OK)
 255     if (fCheckTOF)
 256     {
 257       if (((status & AliESDtrack::kTOFout) != 0) && ((status & AliESDtrack::kTIME) != 0) && mom > TMath::Max(b1, b2))
 258       {
 259         track->GetIntegratedTimes(times);
 260         tofTime  = (Double_t)track->GetTOFsignal();
 261         tofSigma = fTOFmaker->GetExpectedSigma(mom, times[AliPID::kKaon], AliPID::ParticleMass(AliPID::kKaon));
 262         tofRef   = times[AliPID::kKaon];
 263         if (tofRef > 0.0)
 264         {
 265           tofRel   = (tofTime - tofRef) / tofRef;
 266           ymax     = a1 / (mom - b1) + c1;
 267           ymin     = a2 / (mom - b2) + c2;
 268           okTOF    = (tofRel >= ymin && tofRel <= ymax);
 269           //cout << "TOF    -- diff = " << tofDiff << ", rel diff = " << tofRel << ", range = " << ymin << " to " << ymax << ", sigma = " << tofSigma << " --> " << (okTOF ? "OK" : "FAILED") << endl;
 270         }
 271         else
 272         {
 273           okTOF = kTRUE;
 274           //cout << "TOF    -- not checked due to ZERO reference time" << endl;
 275         }
 276       }
 277       else
 278       {
 279         okTOF = kTRUE;
 280         //cout << "TOF    -- not checked because TOF pid absent" << endl;
 281       }
 282     }
 283     else
 284     {
 285       okTOF = kTRUE;
 286     }
 287
 288     // combine checks
 289     okTrack = okQuality && okTPC && okTOF;
 290     //cout << "GLOBAL -- overall = " << (okTrack ? "ACCEPTED" : "REJECTED") << endl;
 291   }
 292   else
 293   {
 294     // check standard ESD cuts
 295     okQuality = fESDtrackCutsITS.IsSelected(track);
 296     //cout << "ITSSA  -- quality = " << (okQuality ? "GOOD" : "BAD") << endl;
 297     if (!okQuality) return kFALSE;
 298
 299     // check dE/dx
 300     if (fCheckITS)
 301     {
 302       itsSignal = track->GetITSsignal();
 303       itsCluMap = track->GetITSClusterMap();
 304       nITS      = 0;
 305       for(k = 2; k < 6; k++) if(itsCluMap & (1 << k)) ++nITS;
 306       if (nITS < 3)
 307       {
 308         okITS = kFALSE;
 309         //cout << "ITS    -- not checked due to too few PID clusters" << endl;
 310       }
 311       else
 312       {
 313         itsNSigma = itsrsp.GetNumberOfSigmas(mom, itsSignal, AliPID::kKaon, nITS, kTRUE);
 314         okITS = (TMath::Abs(itsNSigma) <= fMaxITSband);
 315         //cout << "ITS    -- nsigma = " << itsNSigma << ", max = " << fMaxITSband << " --> " << (okITS ? "OK" : "FAILED") << endl;
 316       }
 317     }
 318     else
 319     {
 320       okITS = kTRUE;
 321     }
 322
 323     okTrack = okQuality && okITS;
 324     //cout << "ITSSA  -- overall = " << (okTrack ? "ACCEPTED" : "REJECTED") << endl;
 325   }
 326
 327   return okTrack;
 328 }