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 ///////////////////////////////////////////////////////////////////////////
17 // ---- CORRECTION FRAMEWORK ----
18 // AliCFAcceptanceCuts implementation
19 // Class to cut on the number of AliTrackReference's
21 ///////////////////////////////////////////////////////////////////////////
22 // author : R. Vernet (renaud.vernet@cern.ch)
23 ///////////////////////////////////////////////////////////////////////////
26 #include "AliMCParticle.h"
27 #include "AliCFAcceptanceCuts.h"
33 ClassImp(AliCFAcceptanceCuts)
35 //______________________________
36 AliCFAcceptanceCuts::AliCFAcceptanceCuts() :
45 fhCutCorrelation(0x0),
51 for (int i=0; i<kNCuts; i++) for (int j=0; j<kNStepQA; j++) fhQA[i][j]=0x0;
54 //______________________________
55 AliCFAcceptanceCuts::AliCFAcceptanceCuts(const Char_t* name, const Char_t* title) :
56 AliCFCutBase(name,title),
64 fhCutCorrelation(0x0),
70 for (int i=0; i<kNCuts; i++) for (int j=0; j<kNStepQA; j++) fhQA[i][j]=0x0;
73 //______________________________
74 AliCFAcceptanceCuts::AliCFAcceptanceCuts(const AliCFAcceptanceCuts& c) :
77 fMinNHitITS(c.fMinNHitITS),
78 fMinNHitTPC(c.fMinNHitTPC),
79 fMinNHitTRD(c.fMinNHitTRD),
80 fMinNHitTOF(c.fMinNHitTOF),
81 fMinNHitMUON(c.fMinNHitMUON),
82 fhCutStatistics(c.fhCutStatistics),
83 fhCutCorrelation(c.fhCutCorrelation),
89 for (int i=0; i<kNCuts; i++) for (int j=0; j<kNStepQA; j++) fhQA[i][j]=c.fhQA[i][j];
92 //______________________________
93 AliCFAcceptanceCuts& AliCFAcceptanceCuts::operator=(const AliCFAcceptanceCuts& c)
96 // Assignment operator
99 AliCFCutBase::operator=(c) ;
101 fMinNHitITS=c.fMinNHitITS;
102 fMinNHitTPC=c.fMinNHitTPC;
103 fMinNHitTRD=c.fMinNHitTRD;
104 fMinNHitTOF=c.fMinNHitTOF;
105 fMinNHitMUON=c.fMinNHitMUON;
106 fhCutStatistics = c.fhCutStatistics ;
107 fhCutCorrelation = c.fhCutCorrelation ;
108 fBitmap = c.fBitmap ;
109 for (int i=0; i<kNCuts; i++) for (int j=0; j<kNStepQA; j++) fhQA[i][j]=c.fhQA[i][j];
114 //______________________________________________________________
115 Bool_t AliCFAcceptanceCuts::IsSelected(TObject *obj) {
117 // check if selections on 'obj' are passed
118 // 'obj' must be an AliMCParticle
121 SelectionBitMap(obj);
123 if (fIsQAOn) FillHistograms(obj,kFALSE);
125 for (UInt_t icut=0; icut<fBitmap->GetNbits(); icut++)
126 if (!fBitmap->TestBitNumber(icut)) return kFALSE ;
128 if (fIsQAOn) FillHistograms(obj,kTRUE);
132 //______________________________
133 void AliCFAcceptanceCuts::SelectionBitMap(TObject* obj) {
135 // checks the number of track references associated to 'obj'
136 // 'obj' must be an AliMCParticle
139 for (UInt_t i=0; i<kNCuts; i++) fBitmap->SetBitNumber(i,kFALSE);
142 TString className(obj->ClassName());
143 if (className.CompareTo("AliMCParticle") != 0) {
144 AliError("obj must point to an AliMCParticle !");
148 AliMCParticle * part = dynamic_cast<AliMCParticle*>(obj);
151 Int_t nHitsITS=0, nHitsTPC=0, nHitsTRD=0, nHitsTOF=0, nHitsMUON=0 ;
152 for (Int_t iTrackRef=0; iTrackRef<part->GetNumberOfTrackReferences(); iTrackRef++) {
153 AliTrackReference * trackRef = part->GetTrackReference(iTrackRef);
155 Int_t detectorId = trackRef->DetectorId();
157 case AliTrackReference::kITS : nHitsITS++ ; break ;
158 case AliTrackReference::kTPC : nHitsTPC++ ; break ;
159 case AliTrackReference::kTRD : nHitsTRD++ ; break ;
160 case AliTrackReference::kTOF : nHitsTOF++ ; break ;
161 case AliTrackReference::kMUON : nHitsMUON++ ; break ;
169 if (nHitsITS >= fMinNHitITS ) fBitmap->SetBitNumber(iCutBit,kTRUE);
171 if (nHitsTPC >= fMinNHitTPC ) fBitmap->SetBitNumber(iCutBit,kTRUE);
173 if (nHitsTRD >= fMinNHitTRD ) fBitmap->SetBitNumber(iCutBit,kTRUE);
175 if (nHitsTOF >= fMinNHitTOF ) fBitmap->SetBitNumber(iCutBit,kTRUE);
177 if (nHitsMUON >= fMinNHitMUON) fBitmap->SetBitNumber(iCutBit,kTRUE);
181 void AliCFAcceptanceCuts::SetMCEventInfo(const TObject* mcInfo) {
183 // Sets pointer to MC event information (AliMCEvent)
187 AliError("Pointer to AliMCEvent !");
191 TString className(mcInfo->ClassName());
192 if (className.CompareTo("AliMCEvent") != 0) {
193 AliError("argument must point to an AliMCEvent !");
197 fMCInfo = (AliMCEvent*) mcInfo ;
201 //__________________________________________________________________________________
202 void AliCFAcceptanceCuts::AddQAHistograms(TList *qaList) {
204 // saves the histograms in a TList
209 qaList->Add(fhCutStatistics);
210 qaList->Add(fhCutCorrelation);
212 for (Int_t j=0; j<kNStepQA; j++) {
213 for(Int_t i=0; i<kNCuts; i++)
214 qaList->Add(fhQA[i][j]);
218 //__________________________________________________________________________________
219 void AliCFAcceptanceCuts::DefineHistograms() {
221 // histograms for cut variables, cut statistics and cut correlations
225 // book cut statistics and cut correlation histograms
226 fhCutStatistics = new TH1F(Form("%s_cut_statistics",GetName()),"",kNCuts,0.5,kNCuts+0.5);
227 fhCutStatistics->SetLineWidth(2);
229 fhCutStatistics->GetXaxis()->SetBinLabel(k,"hits ITS") ; k++;
230 fhCutStatistics->GetXaxis()->SetBinLabel(k,"hits TPC") ; k++;
231 fhCutStatistics->GetXaxis()->SetBinLabel(k,"hits TRD") ; k++;
232 fhCutStatistics->GetXaxis()->SetBinLabel(k,"hits TOF") ; k++;
233 fhCutStatistics->GetXaxis()->SetBinLabel(k,"hits MUON"); k++;
236 fhCutCorrelation = new TH2F(Form("%s_cut_correlation",GetName()),"",kNCuts,0.5,kNCuts+0.5,kNCuts,0.5,kNCuts+0.5);
237 fhCutCorrelation->SetLineWidth(2);
238 for (k=1; k<=kNCuts; k++) {
239 fhCutCorrelation->GetXaxis()->SetBinLabel(k,fhCutStatistics->GetXaxis()->GetBinLabel(k));
240 fhCutCorrelation->GetYaxis()->SetBinLabel(k,fhCutStatistics->GetXaxis()->GetBinLabel(k));
244 for (int i=0; i<kNStepQA; i++) {
245 if (i==0) sprintf(str," ");
246 else sprintf(str,"_cut");
247 fhQA[kCutHitsITS] [i] = new TH1F(Form("%s_HitsITS%s" ,GetName(),str),"",10,0,10);
248 fhQA[kCutHitsTPC] [i] = new TH1F(Form("%s_HitsTPC%s" ,GetName(),str),"",5,0,5);
249 fhQA[kCutHitsTRD] [i] = new TH1F(Form("%s_HitsTRD%s" ,GetName(),str),"",20,0,20);
250 fhQA[kCutHitsTOF] [i] = new TH1F(Form("%s_HitsTOF%s" ,GetName(),str),"",5,0,5);
251 fhQA[kCutHitsMUON][i] = new TH1F(Form("%s_HitsMUON%s" ,GetName(),str),"",5,0,5);
253 for(Int_t i=0; i<kNCuts; i++) fhQA[i][1]->SetLineColor(color);
256 //__________________________________________________________________________________
257 void AliCFAcceptanceCuts::FillHistograms(TObject* obj, Bool_t afterCuts)
260 // fill the QA histograms
262 AliMCParticle* part = dynamic_cast<AliMCParticle *>(obj);
264 Int_t nHitsITS=0, nHitsTPC=0, nHitsTRD=0, nHitsTOF=0, nHitsMUON=0 ;
265 for (Int_t iTrackRef=0; iTrackRef<part->GetNumberOfTrackReferences(); iTrackRef++) {
266 AliTrackReference * trackRef = part->GetTrackReference(iTrackRef);
268 Int_t detectorId = trackRef->DetectorId();
270 case AliTrackReference::kITS : nHitsITS++ ; break ;
271 case AliTrackReference::kTPC : nHitsTPC++ ; break ;
272 case AliTrackReference::kTRD : nHitsTRD++ ; break ;
273 case AliTrackReference::kTOF : nHitsTOF++ ; break ;
274 case AliTrackReference::kMUON : nHitsMUON++ ; break ;
280 fhQA[kCutHitsITS ][afterCuts]->Fill(nHitsITS );
281 fhQA[kCutHitsTPC ][afterCuts]->Fill(nHitsTPC );
282 fhQA[kCutHitsTRD ][afterCuts]->Fill(nHitsTRD );
283 fhQA[kCutHitsTOF ][afterCuts]->Fill(nHitsTOF );
284 fhQA[kCutHitsMUON][afterCuts]->Fill(nHitsMUON);
286 // fill cut statistics and cut correlation histograms with information from the bitmap
287 if (afterCuts) return;
289 // Number of single cuts in this class
290 UInt_t ncuts = fBitmap->GetNbits();
291 for(UInt_t bit=0; bit<ncuts;bit++) {
292 if (!fBitmap->TestBitNumber(bit)) {
293 fhCutStatistics->Fill(bit+1);
294 for (UInt_t bit2=bit; bit2<ncuts;bit2++) {
295 if (!fBitmap->TestBitNumber(bit2))
296 fhCutCorrelation->Fill(bit+1,bit2+1);