1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 // The class AliCFTrackKineCuts is designed to select both generated
17 // and reconstructed tracks of a given range in momentum space,
18 // electric charge and azimuthal emission angle phi
19 // and to provide corresponding QA histograms.
20 // This class inherits from the Analysis' Framework abstract base class
21 // AliAnalysisCuts and is a part of the Correction Framework.
22 // This class acts on single, generated and reconstructed tracks, it is
23 // applicable on ESD and AOD data.
24 // It mainly consists of a IsSelected function that returns a boolean.
25 // This function checks whether the considered track passes a set of cuts:
37 // The cut values for these cuts are set with the corresponding set functions.
38 // All cut classes provided by the correction framework are supposed to be
39 // added in the Analysis Framwork's class AliAnalysisFilter and applied by
40 // the filter via a loop.
42 // author: I. Kraus (Ingrid.Kraus@cern.ch)
44 // AliESDtrackCuts writte by Jan Fiete Grosse-Oetringhaus and
45 // AliRsnDaughterCut class written by A. Pulvirenti.
48 #include <TDirectory.h>
52 #include <AliVParticle.h>
54 #include "AliCFTrackKineCuts.h"
56 ClassImp(AliCFTrackKineCuts)
58 //__________________________________________________________________________________
59 AliCFTrackKineCuts::AliCFTrackKineCuts() :
91 fhBinLimMomentum(0x0),
97 fhBinLimRapidity(0x0),
102 // Default constructor
106 //__________________________________________________________________________________
107 AliCFTrackKineCuts::AliCFTrackKineCuts(Char_t* name, Char_t* title) :
108 AliCFCutBase(name,title),
126 fRequireIsCharged(0),
139 fhBinLimMomentum(0x0),
145 fhBinLimRapidity(0x0),
154 //__________________________________________________________________________________
155 AliCFTrackKineCuts::AliCFTrackKineCuts(const AliCFTrackKineCuts& c) :
157 fMomentumMin(c.fMomentumMin),
158 fMomentumMax(c.fMomentumMax),
169 fRapidityMin(c.fRapidityMin),
170 fRapidityMax(c.fRapidityMax),
174 fRequireIsCharged(c.fRequireIsCharged),
175 fhCutStatistics(c.fhCutStatistics),
176 fhCutCorrelation(c.fhCutCorrelation),
178 fhNBinsMomentum(c.fhNBinsMomentum),
179 fhNBinsPt(c.fhNBinsPt),
180 fhNBinsPx(c.fhNBinsPx),
181 fhNBinsPy(c.fhNBinsPy),
182 fhNBinsPz(c.fhNBinsPz),
183 fhNBinsEta(c.fhNBinsEta),
184 fhNBinsRapidity(c.fhNBinsRapidity),
185 fhNBinsPhi(c.fhNBinsPhi),
186 fhNBinsCharge(c.fhNBinsCharge),
187 fhBinLimMomentum(c.fhBinLimMomentum),
188 fhBinLimPt(c.fhBinLimPt),
189 fhBinLimPx(c.fhBinLimPx),
190 fhBinLimPy(c.fhBinLimPy),
191 fhBinLimPz(c.fhBinLimPz),
192 fhBinLimEta(c.fhBinLimEta),
193 fhBinLimRapidity(c.fhBinLimRapidity),
194 fhBinLimPhi(c.fhBinLimPhi),
195 fhBinLimCharge(c.fhBinLimCharge)
200 ((AliCFTrackKineCuts &) c).Copy(*this);
202 //__________________________________________________________________________________
203 AliCFTrackKineCuts& AliCFTrackKineCuts::operator=(const AliCFTrackKineCuts& c)
206 // Assignment operator
209 AliCFCutBase::operator=(c) ;
210 fMomentumMin = c.fMomentumMin ;
211 fMomentumMax = c.fMomentumMax ;
220 fEtaMin = c.fEtaMin ;
221 fEtaMax = c.fEtaMax ;
222 fRapidityMin = c.fRapidityMin ;
223 fRapidityMax = c.fRapidityMax ;
224 fPhiMin = c.fPhiMin ;
225 fPhiMax = c.fPhiMax ;
226 fCharge = c.fCharge ;
227 fRequireIsCharged = c.fRequireIsCharged ;
228 fhCutStatistics = c.fhCutStatistics ;
229 fhCutCorrelation = c.fhCutCorrelation ;
231 fhNBinsMomentum = c.fhNBinsMomentum;
232 fhNBinsPt = c.fhNBinsPt;
233 fhNBinsPx = c.fhNBinsPx;
234 fhNBinsPy = c.fhNBinsPy;
235 fhNBinsPz = c.fhNBinsPz;
236 fhNBinsEta = c.fhNBinsEta;
237 fhNBinsRapidity = c.fhNBinsRapidity;
238 fhNBinsPhi = c.fhNBinsPhi;
239 fhNBinsCharge = c.fhNBinsCharge;
240 fhBinLimMomentum = c.fhBinLimMomentum;
241 fhBinLimPt = c.fhBinLimPt;
242 fhBinLimPx = c.fhBinLimPx;
243 fhBinLimPy = c.fhBinLimPy;
244 fhBinLimPz = c.fhBinLimPz;
245 fhBinLimEta = c.fhBinLimEta;
246 fhBinLimRapidity = c.fhBinLimRapidity;
247 fhBinLimPhi = c.fhBinLimPhi;
248 fhBinLimCharge = c.fhBinLimCharge;
250 for (Int_t i=0; i<c.kNHist; i++){
251 for (Int_t j=0; j<c.kNStepQA; j++){
252 if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
255 ((AliCFTrackKineCuts &) c).Copy(*this);
259 //__________________________________________________________________________________
260 AliCFTrackKineCuts::~AliCFTrackKineCuts()
265 if (fhCutStatistics) delete fhCutStatistics;
266 if (fhCutCorrelation) delete fhCutCorrelation;
268 for (Int_t i=0; i<kNHist; i++){
269 for (Int_t j=0; j<kNStepQA; j++){
270 if(fhQA[i][j]) delete fhQA[i][j];
273 if(fBitmap) delete fBitmap;
274 if(fhBinLimMomentum) delete fhBinLimMomentum;
275 if(fhBinLimPt) delete fhBinLimPt;
276 if(fhBinLimPx) delete fhBinLimPx;
277 if(fhBinLimPy) delete fhBinLimPy;
278 if(fhBinLimPz) delete fhBinLimPz;
279 if(fhBinLimEta) delete fhBinLimEta;
280 if(fhBinLimRapidity) delete fhBinLimRapidity;
281 if(fhBinLimPhi) delete fhBinLimPhi;
282 if(fhBinLimCharge) delete fhBinLimCharge;
284 //__________________________________________________________________________________
285 void AliCFTrackKineCuts::Initialise()
288 // sets everything to zero
307 fRequireIsCharged = 0;
319 SetRequireIsCharged();
321 for (Int_t i=0; i<kNHist; i++){
322 for (Int_t j=0; j<kNStepQA; j++) {
328 fhCutCorrelation = 0;
329 fBitmap=new TBits(0);
331 //set default bining for QA histograms
332 SetHistogramBins(kCutP,200,0.,20.);
333 SetHistogramBins(kCutPt,200,0.,20.);
334 SetHistogramBins(kCutPx,400,-20.,20.);
335 SetHistogramBins(kCutPy,400,-20.,20.);
336 SetHistogramBins(kCutPz,400,-20.,20.);
337 SetHistogramBins(kCutRapidity,200,-10.,10.);
338 SetHistogramBins(kCutEta,200,-10.,10.);
339 SetHistogramBins(kCutPhi,38,-0.6,7.);
340 SetHistogramBins(kCutCharge,21,-10.5,10.5);
342 //__________________________________________________________________________________
343 void AliCFTrackKineCuts::Copy(TObject &c) const
348 AliCFTrackKineCuts& target = (AliCFTrackKineCuts &) c;
352 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
353 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
355 for (Int_t i=0; i<kNHist; i++){
356 for (Int_t j=0; j<kNStepQA; j++){
357 if(fhQA[i][j]) target.fhQA[i][j] = (TH1F*)fhQA[i][j]->Clone();
362 //__________________________________________________________________________________
363 void AliCFTrackKineCuts::GetBitMap(TObject* obj, TBits *bitmap) {
365 // retrieve the pointer to the bitmap
367 TBits *bm = SelectionBitMap(obj);
370 //__________________________________________________________________________________
371 TBits* AliCFTrackKineCuts::SelectionBitMap(TObject* obj) {
373 // test if the track passes the single cuts
374 // and store the information in a bitmap
377 // bitmap stores the decision of each single cut
378 for(Int_t i=0; i<kNCuts; i++)fBitmap->SetBitNumber(i,kFALSE);
380 // cast TObject into VParticle
381 AliVParticle* particle = dynamic_cast<AliVParticle *>(obj);
382 if ( !particle ) return fBitmap ;
384 for(Int_t i=0; i<kNCuts; i++)fBitmap->SetBitNumber(i,kTRUE);
387 if((particle->P() < fMomentumMin) || (particle->P() > fMomentumMax))
388 fBitmap->SetBitNumber(iCutBit,kFALSE);
390 if ((particle->Pt() < fPtMin) || (particle->Pt() > fPtMax))
391 fBitmap->SetBitNumber(iCutBit,kFALSE);
393 if ((particle->Px() < fPxMin) || (particle->Px() > fPxMax))
394 fBitmap->SetBitNumber(iCutBit,kFALSE);
396 if ((particle->Py() < fPyMin) || (particle->Py() > fPyMax))
397 fBitmap->SetBitNumber(iCutBit,kFALSE);
399 if ((particle->Pz() < fPzMin) || (particle->Pz() > fPzMax))
400 fBitmap->SetBitNumber(iCutBit,kFALSE);
402 if ((particle->Eta() < fEtaMin) || (particle->Eta() > fEtaMax))
403 fBitmap->SetBitNumber(iCutBit,kFALSE);
405 if ((particle->Y() < fRapidityMin) || (particle->Y() > fRapidityMax))
406 fBitmap->SetBitNumber(iCutBit,kFALSE);
408 if ((particle->Phi() < fPhiMin) || (particle->Phi() > fPhiMax))
409 fBitmap->SetBitNumber(iCutBit,kFALSE);
411 if (fCharge < 10 && particle->Charge() != fCharge)
412 fBitmap->SetBitNumber(iCutBit,kFALSE);
414 if (fRequireIsCharged && particle->Charge()==0)
415 fBitmap->SetBitNumber(iCutBit,kFALSE);
419 //__________________________________________________________________________________
420 Bool_t AliCFTrackKineCuts::IsSelected(TObject* obj) {
422 // loops over decisions of single cuts and returns if the track is accepted
424 TBits* bitmap = SelectionBitMap(obj);
426 Bool_t isSelected = kTRUE;
428 for (UInt_t icut=0; icut<bitmap->GetNbits();icut++)
429 if(!bitmap->TestBitNumber(icut)) isSelected = kFALSE;
433 //__________________________________________________________________________________
434 void AliCFTrackKineCuts::Init() {
436 // initialises all histograms and the TList which holds the histograms
441 //__________________________________________________________________________________
442 void AliCFTrackKineCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t *bins)
451 fhNBinsMomentum=nbins+1;
452 fhBinLimMomentum=new Double_t[nbins+1];
453 for(Int_t i=0;i<nbins+1;i++)fhBinLimMomentum[i]=bins[i];
458 fhBinLimPt=new Double_t[nbins+1];
459 for(Int_t i=0;i<nbins+1;i++)fhBinLimPt[i]=bins[i];
464 fhBinLimPx=new Double_t[nbins+1];
465 for(Int_t i=0;i<nbins+1;i++)fhBinLimPx[i]=bins[i];
470 fhBinLimPy=new Double_t[nbins+1];
471 for(Int_t i=0;i<nbins+1;i++)fhBinLimPy[i]=bins[i];
476 fhBinLimPz=new Double_t[nbins+1];
477 for(Int_t i=0;i<nbins+1;i++)fhBinLimPz[i]=bins[i];
481 fhNBinsRapidity=nbins+1;
482 fhBinLimRapidity=new Double_t[nbins+1];
483 for(Int_t i=0;i<nbins+1;i++)fhBinLimRapidity[i]=bins[i];
488 fhBinLimEta=new Double_t[nbins+1];
489 for(Int_t i=0;i<nbins+1;i++)fhBinLimEta[i]=bins[i];
494 fhBinLimPhi=new Double_t[nbins+1];
495 for(Int_t i=0;i<nbins+1;i++)fhBinLimPhi[i]=bins[i];
499 fhNBinsCharge=nbins+1;
500 fhBinLimCharge=new Double_t[nbins+1];
501 for(Int_t i=0;i<nbins+1;i++)fhBinLimCharge[i]=bins[i];
505 //__________________________________________________________________________________
506 void AliCFTrackKineCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t xmin, Double_t xmax)
513 fhNBinsMomentum=nbins+1;
514 fhBinLimMomentum=new Double_t[nbins+1];
515 for(Int_t i=0;i<nbins+1;i++)fhBinLimMomentum[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
520 fhBinLimPt=new Double_t[nbins+1];
521 for(Int_t i=0;i<nbins+1;i++)fhBinLimPt[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
526 fhBinLimPx=new Double_t[nbins+1];
527 for(Int_t i=0;i<nbins+1;i++)fhBinLimPx[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
532 fhBinLimPy=new Double_t[nbins+1];
533 for(Int_t i=0;i<nbins+1;i++)fhBinLimPy[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
538 fhBinLimPz=new Double_t[nbins+1];
539 for(Int_t i=0;i<nbins+1;i++)fhBinLimPz[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
543 fhNBinsRapidity=nbins+1;
544 fhBinLimRapidity=new Double_t[nbins+1];
545 for(Int_t i=0;i<nbins+1;i++)fhBinLimRapidity[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
550 fhBinLimEta=new Double_t[nbins+1];
551 for(Int_t i=0;i<nbins+1;i++)fhBinLimEta[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
556 fhBinLimPhi=new Double_t[nbins+1];
557 for(Int_t i=0;i<nbins+1;i++)fhBinLimPhi[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
561 fhNBinsCharge=nbins+1;
562 fhBinLimCharge=new Double_t[nbins+1];
563 for(Int_t i=0;i<nbins+1;i++)fhBinLimCharge[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
567 //__________________________________________________________________________________
568 void AliCFTrackKineCuts::DefineHistograms() {
570 // histograms for cut variables, cut statistics and cut correlations
574 // book cut statistics and cut correlation histograms
575 fhCutStatistics = new TH1F(Form("%s_cut_statistics",GetName()), Form("%s cut statistics",GetName()), kNCuts,0.5,kNCuts+0.5);
576 fhCutStatistics->SetLineWidth(2);
577 fhCutStatistics->GetXaxis()->SetBinLabel(1,"p");
578 fhCutStatistics->GetXaxis()->SetBinLabel(2,"pt");
579 fhCutStatistics->GetXaxis()->SetBinLabel(3,"px");
580 fhCutStatistics->GetXaxis()->SetBinLabel(4,"py");
581 fhCutStatistics->GetXaxis()->SetBinLabel(5,"pz");
582 fhCutStatistics->GetXaxis()->SetBinLabel(6,"eta");
583 fhCutStatistics->GetXaxis()->SetBinLabel(7,"y");
584 fhCutStatistics->GetXaxis()->SetBinLabel(8,"phi");
585 fhCutStatistics->GetXaxis()->SetBinLabel(9,"charge");
586 fhCutStatistics->GetXaxis()->SetBinLabel(10,"is charged");
588 fhCutCorrelation = new TH2F(Form("%s_cut_correlation",GetName()), Form("%s cut correlation",GetName()), kNCuts,0.5,kNCuts+0.5,kNCuts,0.5,kNCuts+0.5);
589 fhCutCorrelation->SetLineWidth(2);
590 fhCutCorrelation->GetXaxis()->SetBinLabel(1,"p");
591 fhCutCorrelation->GetXaxis()->SetBinLabel(2,"pt");
592 fhCutCorrelation->GetXaxis()->SetBinLabel(3,"px");
593 fhCutCorrelation->GetXaxis()->SetBinLabel(4,"py");
594 fhCutCorrelation->GetXaxis()->SetBinLabel(5,"pz");
595 fhCutCorrelation->GetXaxis()->SetBinLabel(6,"eta");
596 fhCutCorrelation->GetXaxis()->SetBinLabel(7,"y");
597 fhCutCorrelation->GetXaxis()->SetBinLabel(8,"phi");
598 fhCutCorrelation->GetXaxis()->SetBinLabel(9,"charge");
599 fhCutCorrelation->GetXaxis()->SetBinLabel(10,"is charged");
601 fhCutCorrelation->GetYaxis()->SetBinLabel(1,"p");
602 fhCutCorrelation->GetYaxis()->SetBinLabel(2,"pt");
603 fhCutCorrelation->GetYaxis()->SetBinLabel(3,"px");
604 fhCutCorrelation->GetYaxis()->SetBinLabel(4,"py");
605 fhCutCorrelation->GetYaxis()->SetBinLabel(5,"pz");
606 fhCutCorrelation->GetYaxis()->SetBinLabel(6,"eta");
607 fhCutCorrelation->GetYaxis()->SetBinLabel(7,"y");
608 fhCutCorrelation->GetYaxis()->SetBinLabel(8,"phi");
609 fhCutCorrelation->GetYaxis()->SetBinLabel(9,"charge");
610 fhCutCorrelation->GetYaxis()->SetBinLabel(10,"is charged");
613 // book QA histograms
615 for (Int_t i=0; i<kNStepQA; i++) {
616 if (i==0) sprintf(str," ");
617 else sprintf(str,"_cut");
619 fhQA[kCutP][i] = new TH1F(Form("%s_momentum%s",GetName(),str), "",fhNBinsMomentum-1,fhBinLimMomentum);
620 fhQA[kCutPt][i] = new TH1F(Form("%s_transverse_momentum%s",GetName(),str),"",fhNBinsPt-1,fhBinLimPt);
621 fhQA[kCutPx][i] = new TH1F(Form("%s_px%s",GetName(),str), "",fhNBinsPx-1,fhBinLimPx);
622 fhQA[kCutPy][i] = new TH1F(Form("%s_py%s",GetName(),str), "",fhNBinsPy-1,fhBinLimPy);
623 fhQA[kCutPz][i] = new TH1F(Form("%s_pz%s",GetName(),str), "",fhNBinsPz-1,fhBinLimPz);
624 fhQA[kCutRapidity][i]=new TH1F(Form("%s_rapidity%s",GetName(),str), "",fhNBinsRapidity-1,fhBinLimRapidity);
625 fhQA[kCutEta][i] = new TH1F(Form("%s_eta%s",GetName(),str), "",fhNBinsEta-1,fhBinLimEta);
626 fhQA[kCutPhi][i] = new TH1F(Form("%s_phi%s",GetName(),str), "",fhNBinsPhi-1,fhBinLimPhi);
627 fhQA[kCutCharge][i] = new TH1F(Form("%s_charge%s",GetName(),str), "",fhNBinsCharge-1,fhBinLimCharge);
629 fhQA[kCutP][i] ->SetXTitle("momentum p (GeV/c)");
630 fhQA[kCutPt][i] ->SetXTitle("p_{T} (GeV/c)");
631 fhQA[kCutPx][i] ->SetXTitle("p_{x} (GeV/c)");
632 fhQA[kCutPy][i] ->SetXTitle("p_{y} (GeV/c)");
633 fhQA[kCutPz][i] ->SetXTitle("p_{z} (GeV/c)");
634 fhQA[kCutRapidity][i]->SetXTitle("rapidity y");
635 fhQA[kCutEta][i] ->SetXTitle("pseudo rapidity #eta");
636 fhQA[kCutPhi][i] ->SetXTitle("azimuth #phi (rad)");
637 fhQA[kCutCharge][i] ->SetXTitle("charge");
639 for(Int_t i=0; i<kNHist; i++) fhQA[i][1]->SetLineColor(color);
641 //__________________________________________________________________________________
642 void AliCFTrackKineCuts::FillHistograms(TObject* obj, Bool_t b)
645 // fill the QA histograms
649 // cast TObject into VParticle
650 AliVParticle* particle = dynamic_cast<AliVParticle *>(obj);
651 if ( !particle ) return;
653 // index = 0: fill histograms before cuts
654 // index = 1: fill histograms after cuts
656 index = ((b) ? 1 : 0);
658 fhQA[kCutP][index]->Fill(particle->P());
659 fhQA[kCutPt][index]->Fill(particle->Pt());
660 fhQA[kCutPx][index]->Fill(particle->Px());
661 fhQA[kCutPy][index]->Fill(particle->Py());
662 fhQA[kCutPz][index]->Fill(particle->Pz());
663 fhQA[kCutRapidity][index]->Fill(particle->Y());
664 fhQA[kCutEta][index]->Fill(particle->Eta());
665 fhQA[kCutPhi][index]->Fill(particle->Phi());
666 fhQA[kCutCharge][index]->Fill((float)particle->Charge());
668 // fill cut statistics and cut correlation histograms with information from the bitmap
672 TBits* bitmap = SelectionBitMap(obj);
674 // Number of single cuts in this class
675 UInt_t ncuts = bitmap->GetNbits();
676 for(UInt_t bit=0; bit<ncuts;bit++) {
677 if (!bitmap->TestBitNumber(bit)) {
678 fhCutStatistics->Fill(bit+1);
679 for (UInt_t bit2=bit; bit2<ncuts;bit2++) {
680 if (!bitmap->TestBitNumber(bit2))
681 fhCutCorrelation->Fill(bit+1,bit2+1);
686 //__________________________________________________________________________________
687 void AliCFTrackKineCuts::SaveHistograms(const Char_t* dir) {
689 // saves the histograms in a directory (dir)
696 gDirectory->mkdir(dir);
699 gDirectory->mkdir("before_cuts");
700 gDirectory->mkdir("after_cuts");
702 fhCutStatistics->Write();
703 fhCutCorrelation->Write();
705 for (Int_t j=0; j<kNStepQA; j++) {
707 gDirectory->cd("before_cuts");
709 gDirectory->cd("after_cuts");
711 for(Int_t i=0; i<kNHist; i++) fhQA[i][j]->Write();
713 gDirectory->cd("../");
715 gDirectory->cd("../");
717 //__________________________________________________________________________________
718 void AliCFTrackKineCuts::DrawHistograms(Bool_t drawLogScale)
721 // draws some histograms
726 Float_t right = 0.03;
727 Float_t left = 0.175;
729 Float_t bottom = 0.175;
731 TCanvas* canvas1 = new TCanvas("Track_QA_Kinematics_1", "Track QA Kinematics 1", 800, 500);
732 canvas1->Divide(2, 1);
735 fhCutStatistics->SetStats(kFALSE);
736 fhCutStatistics->LabelsOption("v");
737 gPad->SetLeftMargin(left);
738 gPad->SetBottomMargin(0.25);
739 gPad->SetRightMargin(right);
740 gPad->SetTopMargin(0.1);
741 fhCutStatistics->Draw();
744 fhCutCorrelation->SetStats(kFALSE);
745 fhCutCorrelation->LabelsOption("v");
746 gPad->SetLeftMargin(0.30);
747 gPad->SetRightMargin(bottom);
748 gPad->SetTopMargin(0.1);
749 gPad->SetBottomMargin(0.25);
750 fhCutCorrelation->Draw("COLZ");
752 canvas1->SaveAs(Form("%s.eps", canvas1->GetName()));
753 canvas1->SaveAs(Form("%s.ps", canvas1->GetName()));
757 TCanvas* canvas2 = new TCanvas("Track_QA_Kinematics_2", "Track QA Kinematics 2", 1600, 800);
758 canvas2->Divide(4, 2);
761 fhQA[kCutP][0]->SetStats(kFALSE);
762 if(drawLogScale) gPad->SetLogy();
763 gPad->SetRightMargin(right);
764 gPad->SetLeftMargin(left);
765 gPad->SetTopMargin(top);
766 gPad->SetBottomMargin(bottom);
767 fhQA[kCutP][0]->Draw();
768 fhQA[kCutP][1]->Draw("same");
771 fhQA[kCutPt][0]->SetStats(kFALSE);
772 if(drawLogScale) gPad->SetLogy();
773 gPad->SetRightMargin(right);
774 gPad->SetLeftMargin(left);
775 gPad->SetTopMargin(top);
776 gPad->SetBottomMargin(bottom);
777 fhQA[kCutPt][0]->Draw();
778 fhQA[kCutPt][1]->Draw("same");
781 fhQA[kCutRapidity][0]->SetStats(kFALSE);
782 if(drawLogScale) gPad->SetLogy();
783 gPad->SetRightMargin(right);
784 gPad->SetLeftMargin(left);
785 gPad->SetTopMargin(top);
786 gPad->SetBottomMargin(bottom);
787 fhQA[kCutRapidity][0]->Draw();
788 fhQA[kCutRapidity][1]->Draw("same");
791 fhQA[kCutEta][0]->SetStats(kFALSE);
792 if(drawLogScale) gPad->SetLogy();
793 gPad->SetRightMargin(right);
794 gPad->SetLeftMargin(left);
795 gPad->SetTopMargin(top);
796 gPad->SetBottomMargin(bottom);
797 fhQA[kCutEta][0]->Draw();
798 fhQA[kCutEta][1]->Draw("same");
801 fhQA[kCutPx][0]->SetStats(kFALSE);
802 if(drawLogScale) gPad->SetLogy();
803 gPad->SetRightMargin(right);
804 gPad->SetLeftMargin(left);
805 gPad->SetTopMargin(top);
806 gPad->SetBottomMargin(bottom);
807 fhQA[kCutPx][0]->Draw();
808 fhQA[kCutPx][1]->Draw("same");
811 fhQA[kCutPy][0]->SetStats(kFALSE);
812 if(drawLogScale) gPad->SetLogy();
813 gPad->SetRightMargin(right);
814 gPad->SetLeftMargin(left);
815 gPad->SetTopMargin(top);
816 gPad->SetBottomMargin(bottom);
817 fhQA[kCutPy][0]->Draw();
818 fhQA[kCutPy][1]->Draw("same");
821 fhQA[kCutPz][0]->SetStats(kFALSE);
822 if(drawLogScale) gPad->SetLogy();
823 gPad->SetRightMargin(right);
824 gPad->SetLeftMargin(left);
825 gPad->SetTopMargin(top);
826 gPad->SetBottomMargin(bottom);
827 fhQA[kCutPz][0]->Draw();
828 fhQA[kCutPz][1]->Draw("same");
830 canvas2->SaveAs(Form("%s.eps", canvas2->GetName()));
831 canvas2->SaveAs(Form("%s.ps", canvas2->GetName()));
835 TCanvas* canvas3 = new TCanvas("Track_QA_Kinematics_3", "Track QA Kinematics 3", 800, 400);
836 canvas3->Divide(2, 1);
839 gPad->SetRightMargin(right);
840 gPad->SetLeftMargin(left);
841 gPad->SetTopMargin(top);
842 gPad->SetBottomMargin(bottom);
843 fhQA[kCutPhi][0]->SetStats(kFALSE);
844 fhQA[kCutPhi][0]->Draw();
845 fhQA[kCutPhi][1]->Draw("same");
848 gPad->SetRightMargin(right);
849 gPad->SetLeftMargin(left);
850 gPad->SetTopMargin(top);
851 gPad->SetBottomMargin(bottom);
852 fhQA[kCutCharge][0]->SetStats(kFALSE);
853 fhQA[kCutCharge][0]->Draw();
854 fhQA[kCutCharge][1]->Draw("same");
856 canvas3->SaveAs(Form("%s.eps", canvas3->GetName()));
857 canvas3->SaveAs(Form("%s.ps", canvas3->GetName()));
859 //__________________________________________________________________________________
860 void AliCFTrackKineCuts::AddQAHistograms(TList *qaList) const {
862 // saves the histograms in a TList
866 qaList->Add(fhCutStatistics);
867 qaList->Add(fhCutCorrelation);
869 for (Int_t j=0; j<kNStepQA; j++) {
870 for(Int_t i=0; i<kNHist; i++)
871 qaList->Add(fhQA[i][j]);