]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/hfe/AliHFEextraCuts.cxx
committing on behalf of Silvia... new macros for HFE analysis on the train (Markus
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliHFEextraCuts.cxx
CommitLineData
809a4336 1/**************************************************************************
2* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
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**************************************************************************/
50685501 15//
16// Extra cuts implemented by the ALICE Heavy Flavour Electron Group
17// Cuts stored here:
18// - ITS pixels
19// - TPC cluster ratio
20// - TRD tracklets
21//
22// Authors:
23// Markus Fasel <M.Fasel@gsi.de>
24//
3a72645a 25#include <TBits.h>
809a4336 26#include <TClass.h>
27#include <TH1F.h>
28#include <TH2F.h>
29#include <TList.h>
30#include <TString.h>
75d81601 31#include <TMath.h>
809a4336 32
3a72645a 33#include "AliAODTrack.h"
34#include "AliAODPid.h"
35#include "AliESDEvent.h"
36#include "AliESDVertex.h"
809a4336 37#include "AliESDtrack.h"
70da6c5a 38#include "AliLog.h"
809a4336 39#include "AliMCParticle.h"
3a72645a 40#include "AliVEvent.h"
41#include "AliVTrack.h"
42#include "AliVParticle.h"
43#include "AliVertexerTracks.h"
44#include "AliVVertex.h"
809a4336 45
46#include "AliHFEextraCuts.h"
47
48ClassImp(AliHFEextraCuts)
49
50//______________________________________________________
51AliHFEextraCuts::AliHFEextraCuts(const Char_t *name, const Char_t *title):
52 AliCFCutBase(name, title),
3a72645a 53 fEvent(NULL),
809a4336 54 fCutCorrelation(0),
55 fRequirements(0),
3a72645a 56 fTPCiter1(kFALSE),
57 fMinNClustersTPC(0),
809a4336 58 fClusterRatioTPC(0.),
59 fMinTrackletsTRD(0),
60 fPixelITS(0),
70da6c5a 61 fCheck(kFALSE),
62 fQAlist(0x0) ,
dbe3abbe 63 fDebugLevel(0)
809a4336 64{
65 //
66 // Default Constructor
67 //
68 memset(fImpactParamCut, 0, sizeof(Float_t) * 4);
69}
70
71//______________________________________________________
72AliHFEextraCuts::AliHFEextraCuts(const AliHFEextraCuts &c):
73 AliCFCutBase(c),
3a72645a 74 fEvent(c.fEvent),
809a4336 75 fCutCorrelation(c.fCutCorrelation),
76 fRequirements(c.fRequirements),
3a72645a 77 fTPCiter1(c.fTPCiter1),
78 fMinNClustersTPC(c.fMinNClustersTPC),
809a4336 79 fClusterRatioTPC(c.fClusterRatioTPC),
80 fMinTrackletsTRD(c.fMinTrackletsTRD),
81 fPixelITS(c.fPixelITS),
75d81601 82 fCheck(c.fCheck),
dbe3abbe 83 fQAlist(0x0),
70da6c5a 84 fDebugLevel(0)
809a4336 85{
86 //
87 // Copy constructor
88 // Performs a deep copy
89 //
90 memcpy(fImpactParamCut, c.fImpactParamCut, sizeof(Float_t) * 4);
91 if(IsQAOn()){
92 fIsQAOn = kTRUE;
93 fQAlist = dynamic_cast<TList *>(c.fQAlist->Clone());
94 fQAlist->SetOwner();
95 }
96}
97
98//______________________________________________________
99AliHFEextraCuts &AliHFEextraCuts::operator=(const AliHFEextraCuts &c){
100 //
101 // Assignment operator
102 //
103 if(this != &c){
104 AliCFCutBase::operator=(c);
3a72645a 105 fEvent = c.fEvent;
809a4336 106 fCutCorrelation = c.fCutCorrelation;
107 fRequirements = c.fRequirements;
3a72645a 108 fTPCiter1 = c.fTPCiter1;
809a4336 109 fClusterRatioTPC = c.fClusterRatioTPC;
3a72645a 110 fMinNClustersTPC = c.fMinNClustersTPC;
809a4336 111 fMinTrackletsTRD = c.fMinTrackletsTRD;
112 fPixelITS = c.fPixelITS;
75d81601 113 fCheck = c.fCheck;
dbe3abbe 114 fDebugLevel = c.fDebugLevel;
809a4336 115 memcpy(fImpactParamCut, c.fImpactParamCut, sizeof(Float_t) * 4);
116 if(IsQAOn()){
117 fIsQAOn = kTRUE;
118 fQAlist = dynamic_cast<TList *>(c.fQAlist->Clone());
119 fQAlist->SetOwner();
120 }else fQAlist = 0x0;
121 }
122 return *this;
123}
124
125//______________________________________________________
126AliHFEextraCuts::~AliHFEextraCuts(){
127 //
128 // Destructor
129 //
3a72645a 130 if(fQAlist) delete fQAlist;
131}
132
133//______________________________________________________
134void AliHFEextraCuts::SetRecEventInfo(const TObject *event){
135 //
136 // Set Virtual event an make a copy
137 //
138 if (!event) {
139 AliError("Pointer to AliVEvent !");
140 return;
141 }
142 TString className(event->ClassName());
143 if (! (className.CompareTo("AliESDEvent")==0 || className.CompareTo("AliAODEvent")==0)) {
144 AliError("argument must point to an AliESDEvent or AliAODEvent !");
145 return ;
809a4336 146 }
3a72645a 147 fEvent = (AliVEvent*) event;
148
809a4336 149}
150
151//______________________________________________________
9eeae5d5 152Bool_t AliHFEextraCuts::IsSelected(TObject *o){
809a4336 153 //
154 // Steering function for the track selection
155 //
3a72645a 156 TString type = o->IsA()->GetName();
157 AliDebug(2, Form("Object type %s", type.Data()));
158 if(!type.CompareTo("AliESDtrack") || !type.CompareTo("AliAODTrack")){
159 return CheckRecCuts(dynamic_cast<AliVTrack *>(o));
809a4336 160 }
3a72645a 161 return CheckMCCuts(dynamic_cast<AliVParticle *>(o));
809a4336 162}
163
164//______________________________________________________
3a72645a 165Bool_t AliHFEextraCuts::CheckRecCuts(AliVTrack *track){
809a4336 166 //
167 // Checks cuts on reconstructed tracks
168 // returns true if track is selected
169 // QA histograms are filled before track selection and for
170 // selected tracks after track selection
171 //
70da6c5a 172 AliDebug(1, "Called");
809a4336 173 ULong64_t survivedCut = 0; // Bitmap for cuts which are passed by the track, later to be compared with fRequirements
3a72645a 174 if(IsQAOn()) FillQAhistosRec(track, kBeforeCuts);
809a4336 175 // Apply cuts
75d81601 176 Float_t impactR, impactZ, ratioTPC;
3a72645a 177 Double_t hfeimpactR, hfeimpactnsigmaR;
178 Double_t hfeimpactRcut, hfeimpactnsigmaRcut;
179 GetImpactParameters(track, impactR, impactZ);
180 if(TESTBIT(fRequirements, kMinHFEImpactParamR) || TESTBIT(fRequirements, kMinHFEImpactParamNsigmaR)){
181 // Protection for PbPb
182 GetHFEImpactParameterCuts(track, hfeimpactRcut, hfeimpactnsigmaRcut);
183 GetHFEImpactParameters(track, hfeimpactR, hfeimpactnsigmaR);
184 }
185 UInt_t nTPCf = GetTPCfindableClusters(track, fTPCiter1), nclsTPC = GetTPCncls(track, fTPCiter1);
809a4336 186 // printf("Check TPC findable clusters: %d, found Clusters: %d\n", track->GetTPCNclsF(), track->GetTPCNcls());
3a72645a 187 ratioTPC = nTPCf > 0. ? static_cast<Float_t>(nclsTPC)/static_cast<Float_t>(nTPCf) : 1.;
809a4336 188 UChar_t trdTracklets;
3a72645a 189 trdTracklets = GetTRDnTrackletsPID(track);
809a4336 190 UChar_t itsPixel = track->GetITSClusterMap();
3a72645a 191 Int_t status1 = GetITSstatus(track, 0);
192 Int_t status2 = GetITSstatus(track, 1);
70da6c5a 193 Bool_t statusL0 = CheckITSstatus(status1);
194 Bool_t statusL1 = CheckITSstatus(status2);
809a4336 195 if(TESTBIT(fRequirements, kMinImpactParamR)){
196 // cut on min. Impact Parameter in Radial direction
75d81601 197 if(TMath::Abs(impactR) >= fImpactParamCut[0]) SETBIT(survivedCut, kMinImpactParamR);
809a4336 198 }
199 if(TESTBIT(fRequirements, kMinImpactParamZ)){
200 // cut on min. Impact Parameter in Z direction
75d81601 201 if(TMath::Abs(impactZ) >= fImpactParamCut[1]) SETBIT(survivedCut, kMinImpactParamZ);
809a4336 202 }
203 if(TESTBIT(fRequirements, kMaxImpactParamR)){
204 // cut on max. Impact Parameter in Radial direction
75d81601 205 if(TMath::Abs(impactR) <= fImpactParamCut[2]) SETBIT(survivedCut, kMaxImpactParamR);
809a4336 206 }
207 if(TESTBIT(fRequirements, kMaxImpactParamZ)){
208 // cut on max. Impact Parameter in Z direction
75d81601 209 if(TMath::Abs(impactZ) <= fImpactParamCut[3]) SETBIT(survivedCut, kMaxImpactParamZ);
809a4336 210 }
3a72645a 211 if(TESTBIT(fRequirements, kMinHFEImpactParamR)){
212 // cut on min. HFE Impact Parameter in Radial direction
213 if(TMath::Abs(hfeimpactR) >= hfeimpactRcut) SETBIT(survivedCut, kMinHFEImpactParamR);
214 }
215 if(TESTBIT(fRequirements, kMinHFEImpactParamNsigmaR)){
216 // cut on max. HFE Impact Parameter n sigma in Radial direction
217 if(TMath::Abs(hfeimpactnsigmaR) >= hfeimpactnsigmaRcut) SETBIT(survivedCut, kMinHFEImpactParamNsigmaR);
218 }
809a4336 219 if(TESTBIT(fRequirements, kClusterRatioTPC)){
220 // cut on min ratio of found TPC clusters vs findable TPC clusters
221 if(ratioTPC >= fClusterRatioTPC) SETBIT(survivedCut, kClusterRatioTPC);
222 }
223 if(TESTBIT(fRequirements, kMinTrackletsTRD)){
224 // cut on minimum number of TRD tracklets
70da6c5a 225 AliDebug(1, Form("Min TRD cut: [%d|%d]\n", fMinTrackletsTRD, trdTracklets));
809a4336 226 if(trdTracklets >= fMinTrackletsTRD) SETBIT(survivedCut, kMinTrackletsTRD);
227 }
3a72645a 228 if(TESTBIT(fRequirements, kMinNClustersTPC)){
229 // cut on minimum number of TRD tracklets
230 AliDebug(1, Form("Min TPC cut: [%d|%d]\n", fMinNClustersTPC, nclsTPC));
231 if(nclsTPC >= fMinNClustersTPC) SETBIT(survivedCut, kMinNClustersTPC);
232 }
809a4336 233 if(TESTBIT(fRequirements, kPixelITS)){
234 // cut on ITS pixel layers
70da6c5a 235 AliDebug(1, "ITS cluster Map: ");
236 //PrintBitMap(itsPixel);
809a4336 237 switch(fPixelITS){
70da6c5a 238 case kFirst:
239 AliDebug(2, "First");
240 if(TESTBIT(itsPixel, 0))
241 SETBIT(survivedCut, kPixelITS);
242 else
243 if(fCheck && !statusL0)
244 SETBIT(survivedCut, kPixelITS);
dbe3abbe 245 break;
246 case kSecond:
70da6c5a 247 AliDebug(2, "Second");
248 if(TESTBIT(itsPixel, 1))
249 SETBIT(survivedCut, kPixelITS);
250 else
251 if(fCheck && !statusL1)
252 SETBIT(survivedCut, kPixelITS);
dbe3abbe 253 break;
254 case kBoth:
70da6c5a 255 AliDebug(2, "Both");
256 if(TESTBIT(itsPixel,0) && TESTBIT(itsPixel,1))
257 SETBIT(survivedCut, kPixelITS);
258 else
259 if(fCheck && !(statusL0 && statusL1))
260 SETBIT(survivedCut, kPixelITS);
261 break;
dbe3abbe 262 case kAny:
70da6c5a 263 AliDebug(2, "Any");
264 if(TESTBIT(itsPixel, 0) || TESTBIT(itsPixel, 1))
265 SETBIT(survivedCut, kPixelITS);
266 else
267 if(fCheck && !(statusL0 || statusL1))
268 SETBIT(survivedCut, kPixelITS);
dbe3abbe 269 break;
70da6c5a 270 default:
271 AliDebug(2, "None");
272 break;
809a4336 273 }
70da6c5a 274 AliDebug(1, Form("Survived Cut: %s\n", TESTBIT(survivedCut, kPixelITS) ? "YES" : "NO"));
809a4336 275 }
276 if(fRequirements == survivedCut){
277 //
278 // Track selected
279 //
3a72645a 280 AliDebug(2, "Track Survived cuts\n");
281 if(IsQAOn()) FillQAhistosRec(track, kAfterCuts);
809a4336 282 return kTRUE;
283 }
3a72645a 284 AliDebug(2, "Track cut");
809a4336 285 if(IsQAOn()) FillCutCorrelation(survivedCut);
286 return kFALSE;
287}
288
289//______________________________________________________
3a72645a 290Bool_t AliHFEextraCuts::CheckMCCuts(AliVParticle */*track*/) const {
809a4336 291 //
292 // Checks cuts on Monte Carlo tracks
293 // returns true if track is selected
294 // QA histograms are filled before track selection and for
295 // selected tracks after track selection
296 //
297 return kTRUE; // not yet implemented
298}
299
300//______________________________________________________
3a72645a 301void AliHFEextraCuts::FillQAhistosRec(AliVTrack *track, UInt_t when){
809a4336 302 //
303 // Fill the QA histograms for ESD tracks
304 // Function can be called before cuts or after cut application (second argument)
305 //
3a72645a 306 const Int_t kNhistos = 6;
75d81601 307 Float_t impactR, impactZ;
3a72645a 308 GetImpactParameters(track, impactR, impactZ);
309 Int_t nTPCf = GetTPCfindableClusters(track, fTPCiter1), nclsTPC = GetTPCncls(track, fTPCiter1);
310 (dynamic_cast<TH1F *>(fQAlist->At(0 + when * kNhistos)))->Fill(impactR);
311 (dynamic_cast<TH1F *>(fQAlist->At(1 + when * kNhistos)))->Fill(impactZ);
809a4336 312 // printf("TPC findable clusters: %d, found Clusters: %d\n", track->GetTPCNclsF(), track->GetTPCNcls());
3a72645a 313 (dynamic_cast<TH1F *>(fQAlist->At(2 + when * kNhistos)))->Fill(nTPCf > 0. ? static_cast<Float_t>(nclsTPC)/static_cast<Float_t>(nTPCf) : 1.);
314 (dynamic_cast<TH1F *>(fQAlist->At(3 + when * kNhistos)))->Fill(GetTRDnTrackletsPID(track));
315 (dynamic_cast<TH1F *>(fQAlist->At(4 + when * kNhistos)))->Fill(nclsTPC);
809a4336 316 UChar_t itsPixel = track->GetITSClusterMap();
3a72645a 317 TH1 *pixelHist = dynamic_cast<TH1F *>(fQAlist->At(5 + when * kNhistos));
e3fc062d 318 //Int_t firstEntry = pixelHist->GetXaxis()->GetFirst();
319 Double_t firstEntry = 0.5;
809a4336 320 if(!((itsPixel & BIT(0)) || (itsPixel & BIT(1))))
321 pixelHist->Fill(firstEntry + 3);
322 else{
323 if(itsPixel & BIT(0)){
324 pixelHist->Fill(firstEntry);
325 if(itsPixel & BIT(1)) pixelHist->Fill(firstEntry + 2);
326 else pixelHist->Fill(firstEntry + 4);
327 }
328 if(itsPixel & BIT(1)){
329 pixelHist->Fill(firstEntry + 1);
330 if(!(itsPixel & BIT(0))) pixelHist->Fill(firstEntry + 5);
331 }
332 }
333}
334
335// //______________________________________________________
336// void AliHFEextraCuts::FillQAhistosMC(AliMCParticle *track, UInt_t when){
337// //
338// // Fill the QA histograms for MC tracks
339// // Function can be called before cuts or after cut application (second argument)
340// // Not yet implemented
341// //
342// }
343
344//______________________________________________________
345void AliHFEextraCuts::FillCutCorrelation(ULong64_t survivedCut){
346 //
347 // Fill cut correlation histograms for tracks that didn't pass cuts
348 //
3a72645a 349 const Int_t kNhistos = 6;
350 TH2 *correlation = dynamic_cast<TH2F *>(fQAlist->At(2 * kNhistos));
809a4336 351 for(Int_t icut = 0; icut < kNcuts; icut++){
352 if(!TESTBIT(fRequirements, icut)) continue;
353 for(Int_t jcut = icut; jcut < kNcuts; jcut++){
354 if(!TESTBIT(fRequirements, jcut)) continue;
355 if(TESTBIT(survivedCut, icut) && TESTBIT(survivedCut, jcut))
3a72645a 356 correlation->Fill(icut, jcut);
809a4336 357 }
358 }
359}
360
361//______________________________________________________
362void AliHFEextraCuts::AddQAHistograms(TList *qaList){
363 //
364 // Add QA histograms
365 // For each cut a histogram before and after track cut is created
366 // Histos before respectively after cut are stored in different lists
367 // Additionally a histogram with the cut correlation is created and stored
368 // in the top directory
369 //
3a72645a 370
371 const Int_t kNhistos = 6;
809a4336 372 TH1 *histo1D = 0x0;
373 TH2 *histo2D = 0x0;
809a4336 374 TString cutstr[2] = {"before", "after"};
3a72645a 375
376 if(!fQAlist) fQAlist = new TList; // for internal representation, not owner
809a4336 377 for(Int_t icond = 0; icond < 2; icond++){
3a72645a 378 qaList->AddAt((histo1D = new TH1F(Form("%s_impactParamR%s",GetName(),cutstr[icond].Data()), "Radial Impact Parameter", 100, 0, 10)), 0 + icond * kNhistos);
379 fQAlist->AddAt(histo1D, 0 + icond * kNhistos);
809a4336 380 histo1D->GetXaxis()->SetTitle("Impact Parameter");
381 histo1D->GetYaxis()->SetTitle("Number of Tracks");
3a72645a 382 qaList->AddAt((histo1D = new TH1F(Form("%s_impactParamZ%s",GetName(),cutstr[icond].Data()), "Z Impact Parameter", 200, 0, 20)), 1 + icond * kNhistos);
383 fQAlist->AddAt(histo1D, 1 + icond * kNhistos);
809a4336 384 histo1D->GetXaxis()->SetTitle("Impact Parameter");
385 histo1D->GetYaxis()->SetTitle("Number of Tracks");
3a72645a 386 qaList->AddAt((histo1D = new TH1F(Form("%s_tpcClr%s",GetName(),cutstr[icond].Data()), "Cluster Ratio TPC", 10, 0, 1)), 2 + icond * kNhistos);
387 fQAlist->AddAt(histo1D, 2 + icond * kNhistos);
809a4336 388 histo1D->GetXaxis()->SetTitle("Cluster Ratio TPC");
389 histo1D->GetYaxis()->SetTitle("Number of Tracks");
3a72645a 390 qaList->AddAt((histo1D = new TH1F(Form("%s_trdTracklets%s",GetName(),cutstr[icond].Data()), "Number of TRD tracklets", 7, 0, 7)), 3 + icond * kNhistos);
391 fQAlist->AddAt(histo1D, 3 + icond * kNhistos);
809a4336 392 histo1D->GetXaxis()->SetTitle("Number of TRD Tracklets");
393 histo1D->GetYaxis()->SetTitle("Number of Tracks");
3a72645a 394 qaList->AddAt((histo1D = new TH1F(Form("%s_tpcClusters%s",GetName(),cutstr[icond].Data()), "Number of TPC clusters", 161, 0, 160)), 4 + icond * kNhistos);
395 fQAlist->AddAt(histo1D, 4 + icond * kNhistos);
396 histo1D->GetXaxis()->SetTitle("Number of TPC clusters");
397 histo1D->GetYaxis()->SetTitle("Number of Tracks");
398 qaList->AddAt((histo1D = new TH1F(Form("%s_itsPixel%s",GetName(),cutstr[icond].Data()), "ITS Pixel Hits", 6, 0, 6)), 5 + icond * kNhistos);
399 fQAlist->AddAt(histo1D, 5 + icond * kNhistos);
809a4336 400 histo1D->GetXaxis()->SetTitle("ITS Pixel");
401 histo1D->GetYaxis()->SetTitle("Number of Tracks");
402 Int_t first = histo1D->GetXaxis()->GetFirst();
403 TString binNames[6] = { "First", "Second", "Both", "None", "Exclusive First", "Exclusive Second"};
404 for(Int_t ilabel = 0; ilabel < 6; ilabel++)
405 histo1D->GetXaxis()->SetBinLabel(first + ilabel, binNames[ilabel].Data());
406 }
809a4336 407 // Add cut correlation
3a72645a 408 qaList->AddAt((histo2D = new TH2F(Form("%s_cutcorrelation",GetName()), "Cut Correlation", kNcuts, 0, kNcuts - 1, kNcuts, 0, kNcuts -1)), 2 * kNhistos);
409 fQAlist->AddAt(histo2D, 2 * kNhistos);
410 TString labels[kNcuts] = {"MinImpactParamR", "MaxImpactParamR", "MinImpactParamZ", "MaxImpactParamZ", "ClusterRatioTPC", "MinTrackletsTRD", "ITSpixel", "kMinHFEImpactParamR", "kMinHFEImpactParamNsigmaR", "TPC Number of clusters"};
809a4336 411 Int_t firstx = histo2D->GetXaxis()->GetFirst(), firsty = histo2D->GetYaxis()->GetFirst();
412 for(Int_t icut = 0; icut < kNcuts; icut++){
413 histo2D->GetXaxis()->SetBinLabel(firstx + icut, labels[icut].Data());
414 histo2D->GetYaxis()->SetBinLabel(firsty + icut, labels[icut].Data());
415 }
809a4336 416}
dbe3abbe 417
418//______________________________________________________
419void AliHFEextraCuts::PrintBitMap(Int_t bitmap){
420 for(Int_t ibit = 32; ibit--; )
421 printf("%d", bitmap & BIT(ibit) ? 1 : 0);
422 printf("\n");
423}
424
425//______________________________________________________
75d81601 426Bool_t AliHFEextraCuts::CheckITSstatus(Int_t itsStatus) const {
dbe3abbe 427 //
428 // Check whether ITS area is dead
429 //
430 Bool_t status;
431 switch(itsStatus){
432 case 2: status = kFALSE; break;
433 case 3: status = kFALSE; break;
434 case 7: status = kFALSE; break;
435 default: status = kTRUE;
436 }
437 return status;
438}
3a72645a 439
440//______________________________________________________
441Int_t AliHFEextraCuts::GetTRDnTrackletsPID(AliVTrack *track){
442 //
443 // Get Number of TRD tracklets
444 //
445 Int_t nTracklets = 0;
446 if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
447 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
448 nTracklets = esdtrack->GetTRDntrackletsPID();
449 } else if(!TString(track->IsA()->GetName()).CompareTo("AliAODTrack")){
450 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
451 AliAODPid *pidobject = aodtrack->GetDetPid();
452 // this is normally NOT the way to do this, but due to limitation in the
453 // AOD track it is not possible in a different way
454 if(pidobject){
455 Float_t *trdmom = pidobject->GetTRDmomentum();
456 for(Int_t ily = 0; ily < 6; ily++){
457 if(trdmom[ily] > -1) nTracklets++;
458 }
459 } else nTracklets = 6; // No Cut possible
460 }
461 return nTracklets;
462}
463
464//______________________________________________________
465Int_t AliHFEextraCuts::GetITSstatus(AliVTrack *track, Int_t layer){
466 //
467 // Check ITS layer status
468 //
469 Int_t status = 0;
470 if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
471 Int_t det;
472 Float_t xloc, zloc;
473 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
474 esdtrack->GetITSModuleIndexInfo(layer, det, status, xloc, zloc);
475 }
476 return status;
477}
478
479//______________________________________________________
480Int_t AliHFEextraCuts::GetTPCfindableClusters(AliVTrack *track, Bool_t iter1){
481 //
482 // Get Number of findable clusters in the TPC
483 //
484 AliDebug(1, Form("Using TPC clusters from iteration 1: %s", iter1 ? "Yes" : "No"));
485 Int_t nClusters = 159; // in case no Information available consider all clusters findable
486 if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
487 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
488 nClusters = esdtrack->GetTPCNclsF();
489 }
490 return nClusters;
491}
492
493//______________________________________________________
494Int_t AliHFEextraCuts::GetTPCncls(AliVTrack *track, Bool_t iter1){
495 //
496 // Get Number of findable clusters in the TPC
497 //
498 AliDebug(1, Form("Using TPC clusters from iteration 1: %s", iter1 ? "Yes" : "No"));
499 Int_t nClusters = 0; // in case no Information available consider all clusters findable
500 TString type = track->IsA()->GetName();
501 if(!type.CompareTo("AliESDtrack")){
502 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
503 if(iter1)
504 nClusters = esdtrack->GetTPCNclsIter1();
505 else
506 nClusters = esdtrack->GetTPCNcls();
507 }
508 else if(!type.CompareTo("AliAODTrack")){
509 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
510 const TBits &tpcmap = aodtrack->GetTPCClusterMap();
511 for(UInt_t ibit = 0; ibit < tpcmap.GetNbits(); ibit++)
512 if(tpcmap.TestBitNumber(ibit)) nClusters++;
513
514 }
515 return nClusters;
516}
517
518//______________________________________________________
519void AliHFEextraCuts::GetImpactParameters(AliVTrack *track, Float_t &radial, Float_t &z){
520 //
521 // Get impact parameter
522 //
523 TString type = track->IsA()->GetName();
524 if(!type.CompareTo("AliESDtrack")){
525 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
526 esdtrack->GetImpactParameters(radial, z);
527 }
528 else if(!type.CompareTo("AliAODTrack")){
529 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
530 Double_t xyz[3];
531 aodtrack->XYZAtDCA(xyz);
532 z = xyz[2];
533 radial = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]+xyz[1]);
534 }
535}
536
537//______________________________________________________
538void AliHFEextraCuts::GetHFEImpactParameters(AliVTrack *track, Double_t &dcaxy, Double_t &dcansigmaxy){
539 //
540 // Get HFE impact parameter (with recalculated primary vertex)
541 //
542 dcaxy=0;
543 dcansigmaxy=0;
544 if(!fEvent){
545 AliDebug(1, "No Input event available\n");
546 return;
547 }
548 const Double_t kBeampiperadius=3.;
549 TString type = track->IsA()->GetName();
550 Double_t dca[2]={-999.,-999.};
551 Double_t cov[3]={-999.,-999.,-999.};
552
553 // recalculate primary vertex
554 AliVertexerTracks vertexer(fEvent->GetMagneticField());
555 vertexer.SetITSMode();
556 vertexer.SetMinClusters(4);
557 Int_t skipped[2];
558 skipped[0] = track->GetID();
559 vertexer.SetSkipTracks(1,skipped);
560 AliVVertex *vtxESDSkip = vertexer.FindPrimaryVertex(fEvent);
561 vertexer.SetSkipTracks(1,skipped);
562 if(vtxESDSkip->GetNContributors()<2) return;
563
564 // Getting the DCA
565 // Propagation always done on a working copy to not disturb the track params of the original track
566 AliESDtrack *esdtrack = NULL;
567 if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
568 // Case ESD track: take copy constructor
569 esdtrack = new AliESDtrack(*dynamic_cast<AliESDtrack *>(track));
570 } else {
571 // Case AOD track: take different constructor
572 esdtrack = new AliESDtrack(track);
573 }
574 if(esdtrack->PropagateToDCA(vtxESDSkip, fEvent->GetMagneticField(), kBeampiperadius, dca, cov)){
575 // protection
576 dcaxy = dca[0];
577 if(cov[0]) dcansigmaxy = dcaxy/TMath::Sqrt(cov[0]);
578 if(!cov[0]) dcansigmaxy = -99.;
579 }
580 delete esdtrack;
581 delete vtxESDSkip;
582}
583
584
585//______________________________________________________
586void AliHFEextraCuts::GetHFEImpactParameterCuts(AliVTrack *track, Double_t &hfeimpactRcut, Double_t &hfeimpactnsigmaRcut){
587 //
588 // Get HFE impact parameter cut(pt dependent)
589 //
590
591 TString type = track->IsA()->GetName();
592 if(!type.CompareTo("AliESDtrack")){
593 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
594
595 Double_t pt = esdtrack->Pt();
596 //hfeimpactRcut=0.0064+0.078*exp(-0.56*pt); // used Carlo's old parameter
597 hfeimpactRcut=0.011+0.077*exp(-0.65*pt); // used Carlo's new parameter
598 hfeimpactnsigmaRcut=3; // 3 sigma trail cut
599 }
600}
601