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 **************************************************************************/
15 // To select on the Event Class: look at the Trigger mask and the ZDC info.
16 // Only pp-running trigger types implemented so far
17 // handles all masks for the trigger description
18 // and some general combinations like MB1,MB2,MB3,MB4 and MB5.
19 // The argument of IsSelected member function (passed object) is cast into
20 // an AliVEvent, but cuts have a true meaning only for AliESD(AOD)Event
22 // The class derives from AliCFCutBase
23 // Author:S.Arcelli Silvia.Arcelli@cern.ch
30 #include "AliVEvent.h"
31 #include "AliCFEventClassCuts.h"
32 ClassImp(AliCFEventClassCuts)
33 //____________________________________________________________________
34 AliCFEventClassCuts::AliCFEventClassCuts() :
38 fZDCN1EnergyMin(-1.e99),
39 fZDCP1EnergyMin(-1.e99),
40 fZDCN2EnergyMin(-1.e99),
41 fZDCP2EnergyMin(-1.e99),
42 fZDCEM1EnergyMin(-1.e99),
43 fZDCEM2EnergyMin(-1.e99),
44 fZDCN1EnergyMax(1.e99),
45 fZDCP1EnergyMax(1.e99),
46 fZDCN2EnergyMax(1.e99),
47 fZDCP2EnergyMax(1.e99),
48 fZDCEM1EnergyMax(1.e99),
49 fZDCEM2EnergyMax(1.e99),
74 //____________________________________________________________________
75 AliCFEventClassCuts::AliCFEventClassCuts(Char_t* name, Char_t* title) :
76 AliCFCutBase(name,title),
79 fZDCN1EnergyMin(-1.e99),
80 fZDCP1EnergyMin(-1.e99),
81 fZDCN2EnergyMin(-1.e99),
82 fZDCP2EnergyMin(-1.e99),
83 fZDCEM1EnergyMin(-1.e99),
84 fZDCEM2EnergyMin(-1.e99),
85 fZDCN1EnergyMax(1.e99),
86 fZDCP1EnergyMax(1.e99),
87 fZDCN2EnergyMax(1.e99),
88 fZDCP2EnergyMax(1.e99),
89 fZDCEM1EnergyMax(1.e99),
90 fZDCEM2EnergyMax(1.e99),
110 fBitMap=new TBits(0);
114 //_____________________________________________________________________________
115 AliCFEventClassCuts::AliCFEventClassCuts(const AliCFEventClassCuts& c) :
117 fTriggerType(c.fTriggerType),
118 fTriggerAND(c.fTriggerAND),
119 fZDCN1EnergyMin(c.fZDCN1EnergyMin),
120 fZDCP1EnergyMin(c.fZDCP1EnergyMin),
121 fZDCN2EnergyMin(c.fZDCN2EnergyMin),
122 fZDCP2EnergyMin(c.fZDCP2EnergyMin),
123 fZDCEM1EnergyMin(c.fZDCEM1EnergyMin),
124 fZDCEM2EnergyMin(c.fZDCEM2EnergyMin),
125 fZDCN1EnergyMax(c.fZDCN1EnergyMax),
126 fZDCP1EnergyMax(c.fZDCP1EnergyMax),
127 fZDCN2EnergyMax(c.fZDCN2EnergyMax),
128 fZDCP2EnergyMax(c.fZDCP2EnergyMax),
129 fZDCEM1EnergyMax(c.fZDCEM1EnergyMax),
130 fZDCEM2EnergyMax(c.fZDCEM2EnergyMax),
132 fhNBinsTrigger(c.fhNBinsTrigger),
133 fhBinLimTrigger(c.fhBinLimTrigger ),
134 fhNBinsZDCEnN1(c.fhNBinsZDCEnN1),
135 fhBinLimZDCEnN1(c.fhBinLimZDCEnN1),
136 fhNBinsZDCEnP1(c.fhNBinsZDCEnP1),
137 fhBinLimZDCEnP1(c.fhBinLimZDCEnP1),
138 fhNBinsZDCEnN2(c.fhNBinsZDCEnN2),
139 fhBinLimZDCEnN2(c.fhBinLimZDCEnN2),
140 fhNBinsZDCEnP2(c.fhNBinsZDCEnP2),
141 fhBinLimZDCEnP2(c.fhBinLimZDCEnP2),
142 fhNBinsZDCEnEM1(c.fhNBinsZDCEnEM1),
143 fhBinLimZDCEnEM1(c.fhBinLimZDCEnEM1),
144 fhNBinsZDCEnEM2(c.fhNBinsZDCEnEM2),
145 fhBinLimZDCEnEM2(c.fhBinLimZDCEnEM2)
153 //_____________________________________________________________________________
154 AliCFEventClassCuts& AliCFEventClassCuts::operator=(const AliCFEventClassCuts& c){
156 // Assignment operator
159 AliCFCutBase::operator=(c) ;
160 fTriggerType = c.fTriggerType ;
161 fTriggerAND = c.fTriggerAND ;
162 fZDCN1EnergyMin = c.fZDCN1EnergyMin;
163 fZDCP1EnergyMin = c.fZDCP1EnergyMin;
164 fZDCN2EnergyMin = c.fZDCN2EnergyMin;
165 fZDCP2EnergyMin = c.fZDCP2EnergyMin;
166 fZDCEM1EnergyMin = c.fZDCEM1EnergyMin;
167 fZDCEM2EnergyMin = c.fZDCEM2EnergyMin;
168 fZDCN1EnergyMax = c.fZDCN1EnergyMax;
169 fZDCP1EnergyMax = c.fZDCP1EnergyMax;
170 fZDCN2EnergyMax = c.fZDCN2EnergyMax;
171 fZDCP2EnergyMax = c.fZDCP2EnergyMax;
172 fZDCEM1EnergyMax = c.fZDCEM1EnergyMax;
173 fZDCEM2EnergyMax = c.fZDCEM2EnergyMax;
175 fhNBinsTrigger = c.fhNBinsTrigger;
176 fhBinLimTrigger = c.fhBinLimTrigger ;
177 fhNBinsZDCEnN1 = c.fhNBinsZDCEnN1;
178 fhBinLimZDCEnN1 = c.fhBinLimZDCEnN1;
179 fhNBinsZDCEnP1 = c.fhNBinsZDCEnP1;
180 fhBinLimZDCEnP1 = c.fhBinLimZDCEnP1;
181 fhNBinsZDCEnN2 = c.fhNBinsZDCEnN2;
182 fhBinLimZDCEnN2 = c.fhBinLimZDCEnN2;
183 fhNBinsZDCEnP2 = c.fhNBinsZDCEnP2;
184 fhBinLimZDCEnP2 = c.fhBinLimZDCEnP2;
185 fhNBinsZDCEnEM1 = c.fhNBinsZDCEnEM1;
186 fhBinLimZDCEnEM1 = c.fhBinLimZDCEnEM1;
187 fhNBinsZDCEnEM2 = c.fhNBinsZDCEnEM2;
188 fhBinLimZDCEnEM2 = c.fhBinLimZDCEnEM2;
192 for (Int_t i=0; i<c.kNCuts; i++){
193 for (Int_t j=0; j<c.kNStepQA; j++){
194 if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
201 //_____________________________________________________________________________
202 AliCFEventClassCuts::~AliCFEventClassCuts()
207 for (Int_t i=0; i<kNCuts; i++){
208 for (Int_t j=0; j<kNStepQA; j++){
209 if(fhQA[i][j]) delete fhQA[i][j];
213 if(fBitMap)delete fBitMap;
214 if(fhBinLimTrigger)delete fhBinLimTrigger;
215 if(fhBinLimZDCEnN1)delete fhBinLimZDCEnN1;
216 if(fhBinLimZDCEnP1)delete fhBinLimZDCEnP1;
217 if(fhBinLimZDCEnN2)delete fhBinLimZDCEnN2;
218 if(fhBinLimZDCEnP2)delete fhBinLimZDCEnP2;
219 if(fhBinLimZDCEnEM1)delete fhBinLimZDCEnEM1;
220 if(fhBinLimZDCEnEM2)delete fhBinLimZDCEnEM2;
224 //_____________________________________________________________________________
225 void AliCFEventClassCuts::Init() {
227 // initialises all QA histograms
233 //_____________________________________________________________________________
234 void AliCFEventClassCuts::Initialise()
242 // sets pointers to histos to zero
245 for(Int_t i=0; i<kNCuts; i++){
246 for(Int_t j =0; j<kNStepQA; j++){
251 //set default bin number/ranges for QA histograms
253 SetHistogramBins(kTrigger,23,-0.5,22.5);
254 SetHistogramBins(kZDCEnergyN1,800,-500,7500);
255 SetHistogramBins(kZDCEnergyP1,800,-500,7500);
256 SetHistogramBins(kZDCEnergyN2,800,-500,7500);
257 SetHistogramBins(kZDCEnergyP2,800,-500,7500);
258 SetHistogramBins(kZDCEnergyEM1,800,-500,7500);
259 SetHistogramBins(kZDCEnergyEM2,800,-500,7500);
263 //____________________________________________________________________
264 Bool_t AliCFEventClassCuts::IsSelected(TObject* obj) {
266 //Check if the requested cuts are passed
269 TBits* bitmap = SelectionBitMap(obj);
271 Bool_t isSelected = kTRUE;
273 for (UInt_t icut=0; icut<bitmap->GetNbits();icut++)
274 if(!bitmap->TestBitNumber(icut)) isSelected = kFALSE;
279 //____________________________________________________________________
280 TBits *AliCFEventClassCuts::SelectionBitMap(TObject* obj) {
282 //cut on trigger type (just pp running trigger types implemented so far)
283 //and on the energy observed in the ZDC. The argument is cast into
284 //an AliVEvent, but has true meaning only for AliESDEvent type objects.
285 //Check if the requested cuts are passed and return a bitmap
288 for(Int_t j=0;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE);
289 AliVEvent* esd = dynamic_cast<AliVEvent *>(obj);
290 if (!esd ) return fBitMap ;
293 //now start checking the cuts
294 //first assume the event will be accepted:
295 for(Int_t j=0;j<kNCuts;j++)fBitMap->SetBitNumber(j,kTRUE);
300 //look at the Trigger mask in current event
301 TBits *triggerBitMap=new TBits(0);
302 TriggerBitMap(esd,triggerBitMap);
303 //now compare to what was requested as a Trigger:
304 if(fTriggerType.GetNbits()>0)fBitMap->SetBitNumber(0,kFALSE); //trigger required, initialize to false
305 for(Int_t j=0;j<kNTriggers+kNTriggersMB;j++){
306 if(fTriggerType.TestBitNumber(j)){
308 if(triggerBitMap->TestBitNumber(j) == fTriggerType.TestBitNumber(j)){
309 fBitMap->SetBitNumber(0,kTRUE);
311 break;// @least one requested bit fired, ok
314 if(!triggerBitMap->TestBitNumber(j)){
321 delete triggerBitMap;
322 //Then, cut on the energy observed in the ZDC
324 if( esd->GetZDCN1Energy()<fZDCN1EnergyMin || esd->GetZDCN1Energy()>fZDCN1EnergyMax)fBitMap->SetBitNumber(1,kFALSE);
325 if( esd->GetZDCP1Energy()<fZDCP1EnergyMin || esd->GetZDCP1Energy()>fZDCP1EnergyMax)fBitMap->SetBitNumber(2,kFALSE);
326 if( esd->GetZDCN2Energy()<fZDCN2EnergyMin || esd->GetZDCN2Energy()>fZDCN2EnergyMax)fBitMap->SetBitNumber(3,kFALSE);
327 if( esd->GetZDCP2Energy()<fZDCP2EnergyMin || esd->GetZDCP2Energy()>fZDCP2EnergyMax)fBitMap->SetBitNumber(4,kFALSE);
328 if( esd->GetZDCEMEnergy(0)<fZDCEM1EnergyMin || esd->GetZDCEMEnergy(0)>fZDCEM1EnergyMax)fBitMap->SetBitNumber(5,kFALSE);
329 if( esd->GetZDCEMEnergy(1)<fZDCEM2EnergyMin || esd->GetZDCEMEnergy(1)>fZDCEM2EnergyMax)fBitMap->SetBitNumber(6,kFALSE);
334 //_____________________________________________________________________________
335 Bool_t AliCFEventClassCuts::IsTriggered(AliVEvent* ev, TriggerType trigger) {
337 //look at the Trigger mask in current event
338 TBits *triggerBitMap=new TBits(0);
339 TriggerBitMap(ev,triggerBitMap);
340 Bool_t isTriggered=kFALSE;
341 if(triggerBitMap->TestBitNumber(trigger))isTriggered=kTRUE;
342 delete triggerBitMap;
347 //_____________________________________________________________________________
348 void AliCFEventClassCuts::TriggerBitMap(AliVEvent* ev, TBits *bitmapT ) {
351 for(Int_t itrig=0;itrig<kNTriggers+kNTriggersMB;itrig++)bitmapT->SetBitNumber(itrig,kFALSE);
354 ULong64_t triggerMask = ev->GetTriggerMask();
355 //run over the different triggers in the mask, and check which bits have fired
356 for(Int_t itrig=0;itrig<kNTriggers;itrig++){
357 bitmapT->SetBitNumber(itrig,kFALSE);
358 if (triggerMask&(0x1 <<itrig)){
359 bitmapT->SetBitNumber(itrig,kTRUE);
363 //Trigger combinations, Minimum bias triggers
365 //MB1 case: (GFO || V0OR) && !BG
366 if((bitmapT->TestBitNumber(5) || (bitmapT->TestBitNumber(0) || bitmapT->TestBitNumber(1))) && !bitmapT->TestBitNumber(2)) bitmapT->SetBitNumber(17,kTRUE);
368 //MB2 case: (GFO && V0OR) && !BG
369 if((bitmapT->TestBitNumber(5) && (bitmapT->TestBitNumber(0) || bitmapT->TestBitNumber(1))) && !bitmapT->TestBitNumber(2)) bitmapT->SetBitNumber(18,kTRUE);
371 //MB3 case : (GFO && V0AND) && !BG
372 if((bitmapT->TestBitNumber(5) && (bitmapT->TestBitNumber(0) && bitmapT->TestBitNumber(1))) && !bitmapT->TestBitNumber(2)) bitmapT->SetBitNumber(19,kTRUE);
374 //MB4 case: (GFO || V0AND) && !BG
375 if((bitmapT->TestBitNumber(5) || (bitmapT->TestBitNumber(0) && bitmapT->TestBitNumber(1))) && !bitmapT->TestBitNumber(2)) bitmapT->SetBitNumber(20,kTRUE);
377 //MB5 case:: GFO && !BG
378 if(bitmapT->TestBitNumber(5) && !bitmapT->TestBitNumber(2)) bitmapT->SetBitNumber(21,kTRUE);
383 //_____________________________________________________________________________
384 void AliCFEventClassCuts::GetBitMap(TObject* obj, TBits* bitmap){
386 // retrieve the pointer to the bitmap
389 bitmap = SelectionBitMap(obj);
393 //_____________________________________________________________________________
394 void AliCFEventClassCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t *bins)
397 // QA histogram axis parameters
398 // variable bin size:user inputs nbins and the vector of bin limits
403 fhNBinsTrigger=nbins;
404 fhBinLimTrigger=new Double_t[nbins+1];
405 for(Int_t i=0;i<nbins+1;i++)fhBinLimTrigger[i]=bins[i];
409 fhNBinsZDCEnN1=nbins;
410 fhBinLimZDCEnN1=new Double_t[nbins+1];
411 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnN1[i]=bins[i];
415 fhNBinsZDCEnP1=nbins;
416 fhBinLimZDCEnP1=new Double_t[nbins+1];
417 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnP1[i]=bins[i];
421 fhNBinsZDCEnN2=nbins;
422 fhBinLimZDCEnN2=new Double_t[nbins+1];
423 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnN2[i]=bins[i];
427 fhNBinsZDCEnP2=nbins;
428 fhBinLimZDCEnP2=new Double_t[nbins+1];
429 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnP2[i]=bins[i];
433 fhNBinsZDCEnEM1=nbins;
434 fhBinLimZDCEnEM1=new Double_t[nbins+1];
435 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnEM1[i]=bins[i];
439 fhNBinsZDCEnEM2=nbins;
440 fhBinLimZDCEnEM2=new Double_t[nbins+1];
441 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnEM2[i]=bins[i];
448 //_____________________________________________________________________________
449 void AliCFEventClassCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t xmin, Double_t xmax)
452 // QA histogram axis parameters
453 // fixed bin size: user inputs nbins, xmin and xmax
457 fhNBinsTrigger=nbins;
458 fhBinLimTrigger=new Double_t[nbins+1];
459 for(Int_t i=0;i<nbins+1;i++)fhBinLimTrigger[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
463 fhNBinsZDCEnN1=nbins;
464 fhBinLimZDCEnN1=new Double_t[nbins+1];
465 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnN1[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
469 fhNBinsZDCEnP1=nbins;
470 fhBinLimZDCEnP1=new Double_t[nbins+1];
471 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnP1[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
475 fhNBinsZDCEnN2=nbins;
476 fhBinLimZDCEnN2=new Double_t[nbins+1];
477 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnN2[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
481 fhNBinsZDCEnP2=nbins;
482 fhBinLimZDCEnP2=new Double_t[nbins+1];
483 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnP2[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
487 fhNBinsZDCEnEM1=nbins;
488 fhBinLimZDCEnEM1=new Double_t[nbins+1];
489 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnEM1[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
493 fhNBinsZDCEnEM2=nbins;
494 fhBinLimZDCEnEM2=new Double_t[nbins+1];
495 for(Int_t i=0;i<nbins+1;i++)fhBinLimZDCEnEM2[i]=xmin+i*(xmax-xmin)/Double_t(nbins);
500 //_____________________________________________________________________________
501 void AliCFEventClassCuts::DefineHistograms() {
503 // histograms for cut variables
508 AliInfo(Form("Nn QA histos requested, Please first set the QA flag on!"));
512 // book QA histograms
515 for (Int_t i=0; i<kNStepQA; i++) {
516 if (i==0) sprintf(str," ");
517 else sprintf(str,"_cut");
519 fhQA[kTrigger][i] = new TH1F(Form("%s_TriggerBits%s",GetName(),str), "",fhNBinsTrigger,fhBinLimTrigger);
520 fhQA[kZDCEnergyN1][i] = new TH1F(Form("%s_ZDC_Energy_N1%s",GetName(),str), "",fhNBinsZDCEnN1,fhBinLimZDCEnN1);
521 fhQA[kZDCEnergyP1][i] = new TH1F(Form("%s_ZDC_Energy_P1%s",GetName(),str), "",fhNBinsZDCEnP1,fhBinLimZDCEnP1);
522 fhQA[kZDCEnergyN2][i] = new TH1F(Form("%s_ZDC_Energy_N2%s",GetName(),str), "",fhNBinsZDCEnN2,fhBinLimZDCEnN2);
523 fhQA[kZDCEnergyP2][i] = new TH1F(Form("%s_ZDC_Energy_P2%s",GetName(),str), "",fhNBinsZDCEnP2,fhBinLimZDCEnP2);
524 fhQA[kZDCEnergyEM1][i] = new TH1F(Form("%s_ZDC_Energy_EM1%s",GetName(),str), "",fhNBinsZDCEnEM1,fhBinLimZDCEnEM1);
525 fhQA[kZDCEnergyEM2][i] = new TH1F(Form("%s_ZDC_Energy_EM2%s",GetName(),str), "",fhNBinsZDCEnEM2,fhBinLimZDCEnEM2);
528 fhQA[kTrigger][i] ->SetXTitle("Trigger Bits");
529 fhQA[kZDCEnergyN1][i] ->SetXTitle("ZDC Energy N1 (GeV)");
530 fhQA[kZDCEnergyP1][i] ->SetXTitle("ZDC Energy P1 (GeV)");
531 fhQA[kZDCEnergyN2][i] ->SetXTitle("ZDC Energy N2 (GeV)");
532 fhQA[kZDCEnergyP2][i] ->SetXTitle("ZDC Energy P2 (GeV)");
533 fhQA[kZDCEnergyEM1][i] ->SetXTitle("ZDC Energy EM1 (GeV)");
534 fhQA[kZDCEnergyEM2][i] ->SetXTitle("ZDC Energy EM2 (GeV)");
538 for(Int_t i=0; i<kNCuts; i++) fhQA[i][1]->SetLineColor(color);
542 //_____________________________________________________________________________
543 void AliCFEventClassCuts::FillHistograms(TObject* obj, Bool_t b)
546 // fill the QA histograms
550 // cast TObject into VParticle
551 AliVEvent* esd = dynamic_cast<AliVEvent *>(obj);
554 // index = 0: fill histograms before cuts
555 // index = 1: fill histograms after cuts
557 index = ((b) ? 1 : 0);
560 //look at the Trigger mask in current event
561 TBits *triggerBitMap=new TBits(0);
562 TriggerBitMap(esd, triggerBitMap);
565 for(Int_t itrig=0;itrig<kNTriggers+kNTriggersMB;itrig++){
566 if(triggerBitMap->TestBitNumber(itrig)){
567 fhQA[kTrigger][index]->Fill(itrig);
571 delete triggerBitMap;
574 fhQA[kZDCEnergyN1][index]->Fill(esd->GetZDCN1Energy());
575 fhQA[kZDCEnergyP1][index]->Fill(esd->GetZDCP1Energy());
576 fhQA[kZDCEnergyN2][index]->Fill(esd->GetZDCN2Energy());
577 fhQA[kZDCEnergyP2][index]->Fill(esd->GetZDCP2Energy());
578 fhQA[kZDCEnergyEM1][index]->Fill(esd->GetZDCEMEnergy(0));
579 fhQA[kZDCEnergyEM2][index]->Fill(esd->GetZDCEMEnergy(1));
583 //_____________________________________________________________________________
584 void AliCFEventClassCuts::AddQAHistograms(TList *list) const {
586 // saves the histograms in a TList
590 for (Int_t j=0; j<kNStepQA; j++) {
591 for(Int_t i=0; i<kNCuts; i++)
592 list->Add(fhQA[i][j]);