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 **************************************************************************/
17 //_________________________________________________________________________
18 // Count events with different selections
20 // It produces a histogram with the number of events with 9 bins:
21 // 0: all events (that passed the physics selection if it was on)
22 // 1: same but cross check that event pointer did exist
23 // 2: passes vertex cut
24 // 3: passes track number cut, tracks for eta < 0.8
30 // 9: not pileup from SPD
35 // 14: 10 && 2 && 3 && 5
39 // Author: Gustavo Conesa Balbastre (LPSC)
41 //_________________________________________________________________________
44 #include "AliAODHeader.h"
45 #include "AliTriggerAnalysis.h"
46 #include "AliESDEvent.h"
47 #include "AliAODEvent.h"
48 #include "AliESDtrackCuts.h"
49 #include "AliAnalysisManager.h"
50 #include "AliInputEventHandler.h"
52 #include "AliAnalysisTaskCounter.h"
53 ClassImp(AliAnalysisTaskCounter)
55 //________________________________________________________________________
56 AliAnalysisTaskCounter::AliAnalysisTaskCounter(const char *name)
57 : AliAnalysisTaskSE(name),
58 fAcceptFastCluster(kTRUE),
60 fTrackMultEtaCut(0.8),
61 fCaloFilterPatch(kFALSE),
62 fOutputContainer(0x0),
63 fESDtrackCuts(AliESDtrackCuts::GetStandardITSTPCTrackCuts2010()),
64 fTriggerAnalysis (new AliTriggerAnalysis),
66 fhXVertex(0), fhYVertex(0), fhZVertex(0),
67 fhXGoodVertex(0),fhYGoodVertex(0),fhZGoodVertex(0)
70 DefineOutput(1, TList::Class());
73 //________________________________________________________________________
74 AliAnalysisTaskCounter::AliAnalysisTaskCounter()
75 : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisTaskCounter"),
76 fAcceptFastCluster(kTRUE),
78 fTrackMultEtaCut(0.8),
79 fCaloFilterPatch(kFALSE),
80 fOutputContainer(0x0),
81 fESDtrackCuts(AliESDtrackCuts::GetStandardITSTPCTrackCuts2010()),
82 fTriggerAnalysis (new AliTriggerAnalysis),
84 fhXVertex(0), fhYVertex(0), fhZVertex(0),
85 fhXGoodVertex(0),fhYGoodVertex(0),fhZGoodVertex(0)
88 DefineOutput(1, TList::Class());
91 //__________________________________________________
92 AliAnalysisTaskCounter::~AliAnalysisTaskCounter()
96 fOutputContainer->Delete() ;
97 delete fOutputContainer ;
100 if(fESDtrackCuts) delete fESDtrackCuts;
101 if(fTriggerAnalysis) delete fTriggerAnalysis;
106 //-------------------------------------------------------------------
107 void AliAnalysisTaskCounter::UserCreateOutputObjects()
111 fOutputContainer = new TList();
113 fhZVertex = new TH1F("hZVertex", " Z vertex distribution" , 200 , -50 , 50 ) ;
114 fhZVertex->SetXTitle("v_{z} (cm)");
115 fOutputContainer->Add(fhZVertex);
117 fhZGoodVertex = new TH1F("hZGoodVertex", " Good Z vertex distribution" , 200 , -50 , 50 ) ;
118 fhZGoodVertex->SetXTitle("v_{z} (cm)");
119 fOutputContainer->Add(fhZGoodVertex);
121 fhXVertex = new TH1F("hXVertex", " X vertex distribution" , 200 , -2 , 2 ) ;
122 fhXVertex->SetXTitle("v_{x} (cm)");
123 fOutputContainer->Add(fhXVertex);
125 fhXGoodVertex = new TH1F("hXGoodVertex", " Good X vertex distribution" , 200 , -2 , 2 ) ;
126 fhXGoodVertex->SetXTitle("v_{x} (cm)");
127 fOutputContainer->Add(fhXGoodVertex);
129 fhYVertex = new TH1F("hYVertex", " Y vertex distribution" , 200 , -2 , 2 ) ;
130 fhYVertex->SetXTitle("v_{y} (cm)");
131 fOutputContainer->Add(fhYVertex);
133 fhYGoodVertex = new TH1F("hYGoodVertex", " Good Y vertex distribution" , 200 , -2 , 2 ) ;
134 fhYGoodVertex->SetXTitle("v_{y} (cm)");
135 fOutputContainer->Add(fhYGoodVertex);
138 fhNEvents = new TH1I("hNEvents", "Number of analyzed events", 21, 0, 21) ;
139 fhNEvents->SetXTitle("Selection");
140 fhNEvents->SetYTitle("# events");
141 fhNEvents->GetXaxis()->SetBinLabel(1 ,"1 = PS");
142 fhNEvents->GetXaxis()->SetBinLabel(2 ,"2 = 1 & ESD");
143 fhNEvents->GetXaxis()->SetBinLabel(3 ,"3 = 2 & |Z|<10");
144 fhNEvents->GetXaxis()->SetBinLabel(4 ,"4 = 2 & |eta|<0.8");
145 fhNEvents->GetXaxis()->SetBinLabel(5 ,"5 = 3 & 4");
146 fhNEvents->GetXaxis()->SetBinLabel(6 ,"6 = 2 & V0AND");
147 fhNEvents->GetXaxis()->SetBinLabel(7 ,"7 = 3 & 6");
148 fhNEvents->GetXaxis()->SetBinLabel(8 ,"8 = 4 & 6");
149 fhNEvents->GetXaxis()->SetBinLabel(9 ,"9 = 5 & 6");
150 fhNEvents->GetXaxis()->SetBinLabel(10,"10 = 2 & not pileup");
151 fhNEvents->GetXaxis()->SetBinLabel(11,"11 = 2 & good vertex");
152 fhNEvents->GetXaxis()->SetBinLabel(12,"12 = 3 & 11");
153 fhNEvents->GetXaxis()->SetBinLabel(13,"13 = 4 & 11");
154 fhNEvents->GetXaxis()->SetBinLabel(14,"14 = 6 & 11");
155 fhNEvents->GetXaxis()->SetBinLabel(15,"15 = 9 & 11");
156 fhNEvents->GetXaxis()->SetBinLabel(16,"16 = 10 & 11");
157 fhNEvents->GetXaxis()->SetBinLabel(17,"17 = 6 & 10");
158 fhNEvents->GetXaxis()->SetBinLabel(17,"17 = 1 & |Z|<50");
159 fhNEvents->GetXaxis()->SetBinLabel(18,"18 = Reject EMCAL 1");
160 fhNEvents->GetXaxis()->SetBinLabel(19,"19 = 18 & 2");
161 fhNEvents->GetXaxis()->SetBinLabel(20,"20 = Reject EMCAL 2");
162 fhNEvents->GetXaxis()->SetBinLabel(21,"20 = 20 & 2");
164 fOutputContainer->Add(fhNEvents);
166 fOutputContainer->SetOwner(kTRUE);
168 PostData(1,fOutputContainer);
172 //________________________________________________________________________
173 void AliAnalysisTaskCounter::UserExec(Option_t *)
176 // Called for each event
178 //printf("___ Event __ %d __\n",(Int_t)Entry());
180 fhNEvents->Fill(0.5);
182 AliVEvent * event = InputEvent();
184 printf("AliAnalysisTaskCounter::UserExec() - ERROR: event not available \n");
188 AliESDEvent * esdevent = dynamic_cast<AliESDEvent*> (event);
189 AliAODEvent * aodevent = dynamic_cast<AliAODEvent*> (event);
191 TString triggerclasses = "";
192 if(esdevent) triggerclasses = esdevent->GetFiredTriggerClasses();
193 if(aodevent) triggerclasses = aodevent->GetFiredTriggerClasses();
195 if (triggerclasses.Contains("FAST") && !triggerclasses.Contains("ALL") && !fAcceptFastCluster) {
196 //printf("Do not count events from fast cluster, trigger name %s\n",triggerclasses.Data());
200 fhNEvents->Fill(1.5);
203 Bool_t bSelectVZ = kFALSE;
204 Bool_t bV0AND = kFALSE;
205 Bool_t bPileup = kFALSE;
206 Bool_t bGoodV = kFALSE;
207 Bool_t bSelectTrack = kFALSE;
210 //---------------------------------
211 //Get the primary vertex, cut on Z
212 //---------------------------------
214 event->GetPrimaryVertex()->GetXYZ(v) ;
215 fhXVertex->Fill(v[0]);
216 fhYVertex->Fill(v[1]);
217 fhZVertex->Fill(v[2]);
219 if(TMath::Abs(v[2]) < fZVertexCut) {
221 fhNEvents->Fill(2.5);
223 //else printf("Vertex out %f \n",v[2]);
226 //--------------------------------------------------
227 //Tweak for calorimeter only productions
228 //--------------------------------------------------
229 if(fCaloFilterPatch && !esdevent){
230 if(event->GetNumberOfCaloClusters() > 0) {
231 AliVCluster * calo = event->GetCaloCluster(0);
232 if(calo->GetNLabels() == 4){
233 Int_t * selection = calo->GetLabels();
234 bPileup = selection[0];
235 bGoodV = selection[1];
236 bV0AND = selection[2];
237 trackMult = selection[3];
238 //if(selection[0] || selection[1] || selection[2])
239 //printf(" pu %d, gv %d, v0 %d, track mult %d\n ", selection[0], selection[1], selection[2], selection[3]);
241 bSelectTrack = kFALSE;
243 //First filtered AODs, track multiplicity stored there.
244 trackMult = (Int_t) ((AliAODHeader*)fInputEvent->GetHeader())->GetCentrality();
246 }else{//at least one cluster
247 //printf("AliAnalysisTaskCounter::UserExec() - No clusters in event\n");
248 //Remove events with vertex (0,0,0), bad vertex reconstruction
249 if(TMath::Abs(v[0]) < 1.e-6 && TMath::Abs(v[1]) < 1.e-6 && TMath::Abs(v[2]) < 1.e-6) bGoodV = kFALSE;
251 //First filtered AODs, track multiplicity stored there.
252 trackMult = (Int_t) ((AliAODHeader*)fInputEvent->GetHeader())->GetCentrality();
256 //--------------------------------------------------
257 //Count tracks, cut on number of tracks in eta < 0.8
258 //--------------------------------------------------
259 Int_t nTracks = event->GetNumberOfTracks() ;
260 for (Int_t itrack = 0; itrack < nTracks; itrack++) {////////////// track loop
261 AliVTrack * track = (AliVTrack*)event->GetTrack(itrack) ; // retrieve track from esd
264 if(esdevent && !fESDtrackCuts->AcceptTrack((AliESDtrack*)track)) continue;
266 //Do not count tracks out of acceptance cut
267 if(TMath::Abs(track->Eta())< fTrackMultEtaCut) trackMult++;
271 //printf("AliAnalysisTaskCounter::UserExec() - Track Mult %d \n",trackMult);
273 //--------------------------------------------------
274 // At least one track
275 //--------------------------------------------------
277 bSelectTrack = kTRUE;
278 fhNEvents->Fill(3.5);
279 if(bSelectVZ) fhNEvents->Fill(4.5);
282 //---------------------------------
284 //---------------------------------
285 if(esdevent && !fCaloFilterPatch) bV0AND = fTriggerAnalysis->IsOfflineTriggerFired(esdevent, AliTriggerAnalysis::kV0AND);
286 //else if(aodevent && !fCaloFilterPatch) bV0AND = //FIXME FOR AODs
290 fhNEvents->Fill(5.5);
291 if (bSelectVZ) fhNEvents->Fill(6.5);
292 if (bSelectTrack) fhNEvents->Fill(7.5);
293 if (bSelectVZ && bSelectTrack) fhNEvents->Fill(8.5);
296 //---------------------------------
298 //---------------------------------
299 if(!fCaloFilterPatch)
300 bPileup = event->IsPileupFromSPD(3, 0.8, 3., 2., 5.); //Default values, if not it does not compile
301 //bPileup = event->IsPileupFromSPD();
303 fhNEvents->Fill(9.5);
304 if(bV0AND) fhNEvents->Fill(16.5);
307 //---------------------------------
309 //---------------------------------
310 if(esdevent && !fCaloFilterPatch) bGoodV = CheckForPrimaryVertex();
313 fhXGoodVertex->Fill(v[0]);
314 fhYGoodVertex->Fill(v[1]);
315 fhZGoodVertex->Fill(v[2]);
317 fhNEvents->Fill(10.5);
318 if(bSelectVZ) fhNEvents->Fill(11.5);
319 if(bSelectTrack) fhNEvents->Fill(12.5);
320 if(bV0AND) fhNEvents->Fill(13.5);
321 if(bSelectVZ && bSelectTrack && bV0AND)
322 fhNEvents->Fill(14.5);
323 if(!bPileup) fhNEvents->Fill(15.5);
326 //printf("AliAnalysisTaskCounter::UserExec() : z vertex %d, good vertex %d, v0and %d, pile up %d, track mult %d\n ", bSelectVZ, bGoodV, bV0AND, bPileup, trackMult);
328 // Events that could be rejected in EMCAL
329 // LHC11a, SM4 and some SM3 events cut with this
330 Bool_t bEMCALRejected = kFALSE;
331 for (Int_t i = 0; i < InputEvent()->GetNumberOfCaloClusters(); i++)
333 AliVCluster *clus = InputEvent()->GetCaloCluster(i);
335 if ((clus->E() > 500 && clus->GetNCells() > 200 ) || clus->GetNCells() > 200) {
337 //printf("Counter: Reject event with cluster: E %f, ncells %d\n",clus->E(),clus->GetNCells());
339 fhNEvents->Fill(17.5);
340 if(bSelectVZ) fhNEvents->Fill(18.5);
341 bEMCALRejected = kTRUE;
347 //LHC11a, 3 last runs, cut with this
349 // Count number of cells in SM3 with energy larger than 0.1, cut on this number
352 for(Int_t icell = 0; icell < event->GetEMCALCells()->GetNumberOfCells(); icell++){
353 if(event->GetEMCALCells()->GetAmplitude(icell) > 0.1 && event->GetEMCALCells()->GetCellNumber(icell)/(24*48)==3) ncellsSM3++;
354 if(event->GetEMCALCells()->GetAmplitude(icell) > 0.1 && event->GetEMCALCells()->GetCellNumber(icell)/(24*48)==4) ncellsSM4++;
358 if(triggerclasses.Contains("EMC")) ncellcut = 35;
360 if( ncellsSM3 >= ncellcut || ncellsSM4 >= 100 ){
362 //printf("Counter: reject event with ncells in SM3: ncells %d\n",ncells);
364 fhNEvents->Fill(19.5);
365 if(bSelectVZ) fhNEvents->Fill(20.5);
370 PostData(1,fOutputContainer);
374 //____________________________________________________________________________
375 Bool_t AliAnalysisTaskCounter::CheckForPrimaryVertex(){
376 //Check if the vertex was well reconstructed, copy from V0Reader of conversion group
377 //It only works for ESDs
379 AliESDEvent * event = dynamic_cast<AliESDEvent*> (InputEvent());
382 if(event->GetPrimaryVertexTracks()->GetNContributors() > 0) {
386 if(event->GetPrimaryVertexTracks()->GetNContributors() < 1) {
388 if(event->GetPrimaryVertexSPD()->GetNContributors() > 0) {
389 //cout<<"spd vertex type::"<< fESDEvent->GetPrimaryVertex()->GetName() << endl;
393 if(event->GetPrimaryVertexSPD()->GetNContributors() < 1) {
394 // cout<<"bad vertex type::"<< fESDEvent->GetPrimaryVertex()->GetName() << endl;
399 //return fInputEvent->GetPrimaryVertex()->GetNContributors()>0;
404 //_____________________________________________________
405 void AliAnalysisTaskCounter::FinishTaskOutput()
407 // Put in the output some event summary histograms
409 AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
410 AliInputEventHandler *inputH = dynamic_cast<AliInputEventHandler*>(am->GetInputEventHandler());
412 TH2F *histStat = dynamic_cast<TH2F*>(inputH->GetStatistics());
413 TH2F *histBin0 = dynamic_cast<TH2F*>(inputH->GetStatistics("BIN0"));
416 fOutputContainer->Add(histStat);
417 else if(DebugLevel() > 1)
418 printf("AliAnalysisTaskCounter::FinishTaskOutput() - Stat histogram not available check, \n if ESDs, that AliPhysicsSelection was on, \n if AODs, if EventStat_temp.root exists \n");
421 fOutputContainer->Add(histBin0);