AliPhysics  2c6b7ad (2c6b7ad)
AliAnalysisTaskEmcalRun2QA.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2016, 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  **************************************************************************/
15 
16 #include <cstring>
17 #include <Riostream.h>
18 
19 #include <TClonesArray.h>
20 #include <TH1F.h>
21 #include <TH2F.h>
22 #include <TH3F.h>
23 #include <THnSparse.h>
24 #include <THashList.h>
25 
26 #include "AliParticleContainer.h"
27 #include "AliClusterContainer.h"
28 #include "AliTrackContainer.h"
29 #include "AliCentrality.h"
30 #include "AliVCluster.h"
31 #include "AliVParticle.h"
32 #include "AliVTrack.h"
33 #include "AliLog.h"
34 #include "AliEMCALGeometry.h"
35 #include "AliEMCALGeoParams.h"
36 #include "AliPicoTrack.h"
37 #include "AliVVZERO.h"
38 #include "AliESDUtils.h"
39 
41 
42 using std::cout;
43 using std::endl;
44 
46 
47 //________________________________________________________________________
50 fCellEnergyCut(0.05),
51 fParticleLevel(kFALSE),
52 fIsMC(kFALSE),
53 fCentMethod2(""),
54 fCentMethod3(""),
55 fDoV0QA(0),
56 fDoEPQA(0),
57 fDoLeadingObjectPosition(0),
58 fMaxCellsInCluster(50),
59 fPtBinWidth(0.5),
60 fMaxPt(150),
61 fSeparateEMCalDCal(kTRUE),
62 fIsEmbedded(kFALSE),
63 fCent2(0),
64 fCent3(0),
65 fVZERO(0),
66 fV0ATotMult(0),
67 fV0CTotMult(0),
68 fLeadingTrack(),
69 //fGeoName("EMCAL_FIRSTYEARV1"),//fGeomEMCal(0),
70 fCaloUtils(0),
71 
72 fHistCellIDvsE(0),fHistCellIDvsELow(0),fHistCellEtaPhi(0),
73 fHistClusterIDvsE(0),fHistClusterIDvsELow(0),fHistClusterEtaPhi(0),fHistNrCellsInCluster(0),
74 
75 fOutputList_Event(), fOutputList_Cell(), fOutputList_Cluster(),
76 fHistManager("AliAnalysisTaskEmcalRun2QA")
77 {
78  InitConstants();
79 }
80 //________________________________________________________________________
82 AliAnalysisTaskEmcalLight(name, kTRUE),
83 fCellEnergyCut(0.05),
84 fParticleLevel(kFALSE),
85 fIsMC(kFALSE),
86 fCentMethod2(""),
87 fCentMethod3(""),
88 fDoV0QA(0),
89 fDoEPQA(0),
90 fDoLeadingObjectPosition(0),
91 fMaxCellsInCluster(50),
92 fPtBinWidth(0.5),
93 fMaxPt(150),
94 fSeparateEMCalDCal(kTRUE),
95 fIsEmbedded(kFALSE),
96 fCent2(0),
97 fCent3(0),
98 fVZERO(0),
99 fV0ATotMult(0),
100 fV0CTotMult(0),
101 //fGeoName("EMCAL_FIRSTYEARV1"),//fGeomEMCal(0),//
102 fCaloUtils(0),
103 
104 fHistCellIDvsE(0),fHistCellIDvsELow(0),fHistCellEtaPhi(0),
105 fHistClusterIDvsE(0),fHistClusterIDvsELow(0),fHistClusterEtaPhi(0),fHistNrCellsInCluster(0),
106 
107 fOutputList_Event(), fOutputList_Cell(), fOutputList_Cluster(),
108 fHistManager(name)
109 {
110  InitConstants();
111 }
112 //________________________________________________________________________
114 {
115  //..Standard
116  fRtoD=180.0/TMath::Pi();
117  memset(fNTotClusters, 0, sizeof(Int_t)*2);
118  fDebug=0; //printout some debug lines
119 
121 
122  //..super ugly but I can not retrieve the info directly from the
123  //..AliEMCALGeometry because it is initialized later than the OutputObjects (where it is needed)
124  //number of modules, 12 for EMCal + 8 for DCAL
125  //2013 fNumOfSuperMod=12; //HARD CODED eeeeewwwhh!
126  fNumOfSuperMod=20; //2015 HARD CODED eeeeewwwhh!
128 }
129 //________________________________________________________________________
131 {
132  //..Destructor
133  if (fCaloUtils) delete fCaloUtils ;
134 }
135 
136 //________________________________________________________________________
138 {
139  if(fDebug==1)cout<<"Inside of: AliAnalysisTaskEmcalRun2QA::UserCreateOutputObjects()"<<endl;
140 
141  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
142  // General Event Histograms
143  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
145 
146  AliEmcalContainer* cont = 0;
147 
148  TString histname;
149  TString title;
150 
151  Int_t nPtBins = TMath::CeilNint(fMaxPt / fPtBinWidth);
152 
153  fOutputList_Event = new TList();
154  fOutputList_Event ->SetOwner();
155  fOutputList_Event ->SetName("EventQA");
156 
157  fOutputList_Cell = new TList();
158  fOutputList_Cell ->SetOwner();
159  fOutputList_Cell ->SetName("CellQA");
160 
161  fOutputList_Cluster = new TList();
162  fOutputList_Cluster ->SetOwner();
163  fOutputList_Cluster ->SetName("ClusterQA");
164 
168 
169 
170  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
171  // Settings for Histograms
172  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
173 
174  //..per module!
175  Int_t fNMaxCols = AliEMCALGeoParams::fgkEMCALCols; //Eta direction = 48
176  Int_t fNMaxRows = AliEMCALGeoParams::fgkEMCALRows; //Phi direction = 24
177 
178  TString Calo="EMCal";
179 
180  //cout<<"Number of EMCal modules : "<<AliEMCALGeoParams::fgkEMCALModules<<", fNModules: "<<fNModules<<", fGeom: "<< fNumOfSuperMod<<endl;
181  //Int_t NumSupMod = fGeom->GetNumberOfSuperModules(); //GetNumberOfSuperModules()
182  //cout<<"Number of Supermodules: "<<NumSupMod<<endl;
183 
184  if(Calo=="EMCal")
185  {
186  fNMaxColsAbs = 2*fNMaxCols; //two supermodules next to each other
187  //..as the modules come in rows of 2 the maximum length is half the number of modules x row in a module
188  //..what about smaller modules?????
189  fNMaxRowsAbs = Int_t(fNumOfSuperMod/2)*fNMaxRows;
190  }
191  else
192  {
193  fNMaxColsAbs=fNumOfSuperMod*fNMaxRows;
194  }
195 
196  cout<<"Info1: "<<endl;
197  //number of rows
198  Int_t NRows = AliEMCALGeoParams::fgkEMCALRows;
199  cout<<"A"<<endl;
200  Int_t NCol = AliEMCALGeoParams::fgkEMCALCols;
201  cout<<"B"<<endl;
202  Int_t NSMod = AliEMCALGeoParams::fgkEMCALModules;
203 // cout<<" Number of super modules :"<<fGeom->GetNumberOfSuperModules()<<", "<<fCaloUtils->GetNumberOfSuperModulesUsed()<<"(fGeom, fCaloUtils)"<<endl;
204  cout<<"number of mo from geoparams: "<<NSMod<<endl;
205  cout<<"number of col from geoparams: "<<NCol<<endl;
206  cout<<"number of row from geoparams: "<<NRows<<endl;
207  cout<<"Number of rows: "<<NRows<<", rowsMax: "<<fNMaxRowsAbs<<", Number of collumns:"<<NCol<<", colMax: "<<fNMaxColsAbs<<endl;
208 
209  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
210  // Cluster Histograms
211  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
212 
213  for (auto cont_it : fClusterCollArray) {
214  cont = cont_it.second;
215 
216  fHistManager.CreateHistoGroup(cont->GetArrayName());
217  for (Int_t i = 0; i < GetNCentBins(); i++)
218  {
219  histname = TString::Format("%s/fHistRejectionReason_%d", cont->GetArrayName().Data(), i);
220  title = histname + ";Rejection reason;#it{E}_{cluster} (GeV);counts";
221  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, 40, 0, 100);
222  SetRejectionReasonLabels(hist->GetXaxis());
223 
224  histname = TString::Format("%s/fHistClusPosition_%d", cont->GetArrayName().Data(), i);
225  title = histname + ";#it{x} (cm);#it{y} (cm);#it{z} (cm)";
226  fHistManager.CreateTH3(histname.Data(), title.Data(), 50, -500, 500, 50, -500, 500, 50, -500, 500);
227 
228  histname = TString::Format("%s/fHistClusPhiEtaEnergy_%d", cont->GetArrayName().Data(), i);
229  title = histname + ";#eta;#phi;#it{E}_{cluster} (GeV)";
230  fHistManager.CreateTH3(histname.Data(), title.Data(), 50, -1, 1, 50, 0, TMath::TwoPi(), nPtBins, 0, fMaxPt);
231 
232  //Three different energies
233  histname = TString::Format("%s/fHistClusEnergy_%d", cont->GetArrayName().Data(), i);
234  title = histname + ";#it{E}_{cluster} (GeV);counts";
235  fHistManager.CreateTH1(histname.Data(), title.Data(), nPtBins, 0, fMaxPt);
236 
237  histname = TString::Format("%s/fHistClusNonLinCorrEnergy_%d", cont->GetArrayName().Data(), i);
238  title = histname + ";#it{E}_{cluster} (GeV);counts";
239  fHistManager.CreateTH1(histname.Data(), title.Data(), nPtBins, 0, fMaxPt);
240 
241  histname = TString::Format("%s/fHistClusHadCorrEnergy_%d", cont->GetArrayName().Data(), i);
242  title = histname + ";#it{E}_{cluster} (GeV);counts";
243  fHistManager.CreateTH1(histname.Data(), title.Data(), nPtBins, 0, fMaxPt);
244  }
245 
246  cout<<"create histograms for cluster name: "<<cont->GetArrayName()<<endl;
247 
252 
253  for (Int_t i = 0; i < fNumOfSuperMod+1; i++)
254  {
255  fHistClusterIDvsE[i] = new TH2F(Form("fHistClusterIDvsE_%d",i),Form("fHistClusterIDvsE_%d",i),18000,0,18000,100,0,50);
256  if(i==fNumOfSuperMod)fHistClusterIDvsE[i]->GetXaxis()->SetTitle("Leading cell ID (all SuperModules");
257  else fHistClusterIDvsE[i]->GetXaxis()->SetTitle(Form("Leading cell ID (SuperMod No. %d)",i));
258  fHistClusterIDvsE[i]->GetYaxis()->SetTitle("Energy probably in GeV");
260 
261  fHistClusterIDvsELow[i] = new TH2F(Form("fHistClusterIDvsELow_%d",i),Form("fHistClusterIDvsELow_%d",i),18000,0,18000,200,0,2);
262  if(i==fNumOfSuperMod)fHistClusterIDvsELow[i]->GetXaxis()->SetTitle("Leading cell ID (all SuperModules");
263  else fHistClusterIDvsELow[i]->GetXaxis()->SetTitle(Form("Leading cell ID (SuperMod No. %d)",i));
264  fHistClusterIDvsELow[i]->GetYaxis()->SetTitle("Energy probably in GeV");
266 
267 /* // fHistClusterEtaPhi[i] = new TH2F(Form("fHistClusterEtaPhi_%d",i),Form("fHistClusterEtaPhi_%d",i),1120,60,130,400,-1,1); //0.012625, 0.0125, 0.0136
268  fHistClusterEtaPhi[i] = new TH2F(Form("fHistClusterEtaPhi_%d",i),Form("fHistClusterEtaPhi_%d",i),5600,60,130,400,-1,1);
269  if(i==fNumOfSuperMod)fHistClusterEtaPhi[i]->GetXaxis()->SetTitle("Leading cell #phi (all SuperModules");
270  else fHistClusterEtaPhi[i]->GetXaxis()->SetTitle(Form("Leading cell #phi (SuperMod No. %d)",i));
271  fHistClusterEtaPhi[i]->GetYaxis()->SetTitle("Leading cell #eta");
272  fOutputList_Cluster->Add(fHistClusterEtaPhi[i]);
273 */
274  //For row and collumn
275  fHistClusterEtaPhi[i] = new TH2F(Form("fHistClusterEtaPhi_%d",i),Form("fHistClusterEtaPhi_%d",i),fNMaxColsAbs+2,-1.5,fNMaxColsAbs+0.5, fNMaxRowsAbs+2,-1.5,fNMaxRowsAbs+0.5);
276  if(i==fNumOfSuperMod)fHistClusterEtaPhi[i]->GetXaxis()->SetTitle("Leading cell column (#eta direction)(all SuperModules)");
277  else fHistClusterEtaPhi[i]->GetXaxis()->SetTitle(Form("Leading cell column (#eta direction)(SuperMod No. %d)",i));
278  fHistClusterEtaPhi[i]->GetYaxis()->SetTitle("Leading cell row (#phi direction)");
280 
281 
282  fHistNrCellsInCluster[i] = new TH1F(Form("fHistNrCellsInCluster_%d",i),Form("fHistNrCellsInCluster_%d",i),30,0,30);
283  if(i==fNumOfSuperMod)fHistNrCellsInCluster[i]->GetXaxis()->SetTitle("Number of cells (all SuperModules)");
284  else fHistNrCellsInCluster[i]->GetXaxis()->SetTitle(Form("Number of cells (SuperMod No. %d)",i));
285  fHistNrCellsInCluster[i]->GetYaxis()->SetTitle("N_{Cluster}");
287  }
288  }
289 
290  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
291  // Cells Histograms
292  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
293  if (!fCaloCellsName.IsNull())
294  {
295 
299 
300  for (Int_t i = 0; i < fNumOfSuperMod+1; i++)
301  {
302  fHistCellIDvsE[i] = new TH2F(Form("fHistCellIDvsE_%d",i),Form("fHistCellIDvsE_%d",i),18000,0,18000,100,0,50);
303  if(i==fNumOfSuperMod) fHistCellIDvsE[i]->GetXaxis()->SetTitle("Cell ID (all SuperModules");
304  else fHistCellIDvsE[i]->GetXaxis()->SetTitle(Form("Cell ID (SuperMod No. %d)",i));
305  fHistCellIDvsE[i]->GetYaxis()->SetTitle("Energy probably in GeV");
307 
308  fHistCellIDvsELow[i] = new TH2F(Form("fHistCellIDvsELow_%d",i),Form("fHistCellIDvsELow_%d",i),18000,0,18000,200,0,2);
309  if(i==fNumOfSuperMod) fHistCellIDvsELow[i]->GetXaxis()->SetTitle("Cell ID (all SuperModules");
310  else fHistCellIDvsELow[i]->GetXaxis()->SetTitle(Form("Cell ID (SuperMod No. %d)",i));
311  fHistCellIDvsELow[i]->GetYaxis()->SetTitle("Energy probably in GeV");
313 
314  /*
315  //for eta phi
316  if(i==fNumOfSuperMod)fHistCellEtaPhi[i] = new TH2F(Form("fHistCellEtaPhi_%d",i),Form("fHistCellEtaPhi_%d",i),700,60,130,200,-1,1); //superfine (not necessary)
317  else fHistCellEtaPhi[i] = new TH2F(Form("fHistCellEtaPhi_%d",i),Form("fHistCellEtaPhi_%d",i),5600,60,130,400,-1,1); //superfine (not necessary)
318  if(i==fNumOfSuperMod)fHistCellEtaPhi[i]->GetXaxis()->SetTitle("Cell #phi (all SuperModules)");
319  else fHistCellEtaPhi[i]->GetXaxis()->SetTitle(Form("Cell #phi (SuperMod No. %d)",i));
320  fHistCellEtaPhi[i]->GetYaxis()->SetTitle("Cell #eta");
321  fOutputList_Cell->Add(fHistCellEtaPhi[i]);
322  */
323 
324  //For row and collumn
325  fHistCellEtaPhi[i] = new TH2F(Form("fHistCellEtaPhi_%d",i),Form("fHistCellEtaPhi_%d",i),fNMaxColsAbs+2,-1.5,fNMaxColsAbs+0.5, fNMaxRowsAbs+2,-1.5,fNMaxRowsAbs+0.5);
326  if(i==fNumOfSuperMod)fHistCellEtaPhi[i]->GetXaxis()->SetTitle("column (#eta direction)(all SuperModules)");
327  else fHistCellEtaPhi[i]->GetXaxis()->SetTitle(Form("column (#eta direction)(SuperMod No. %d)",i));
328  fHistCellEtaPhi[i]->GetYaxis()->SetTitle("row (#phi direction)");
330  }
331  }
332 
333  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
334  // Other Stuff
335  // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
336 
337  Int_t dim = 0;
338  TString axistitle[40];
339  Int_t nbins[40] = {0};
340  Double_t min[40] = {0};
341  Double_t max[40] = {0};
342 
344  {
345  axistitle[dim] = "Centrality %";
346  nbins[dim] = 101;
347  min[dim] = 0;
348  max[dim] = 101;
349  dim++;
350 
351  if (!fCentMethod2.IsNull())
352  {
353  axistitle[dim] = Form("Centrality %s %%", fCentMethod2.Data());
354  nbins[dim] = 101;
355  min[dim] = 0;
356  max[dim] = 101;
357  dim++;
358  }
359 
360  if (!fCentMethod3.IsNull())
361  {
362  axistitle[dim] = Form("Centrality %s %%", fCentMethod3.Data());
363  nbins[dim] = 101;
364  min[dim] = 0;
365  max[dim] = 101;
366  dim++;
367  }
368 
369  if (fDoV0QA==1)
370  {
371  axistitle[dim] = "V0A total multiplicity";
372  nbins[dim] = 200;
373  min[dim] = 0;
374  max[dim] = 20000;
375  dim++;
376 
377  axistitle[dim] = "V0C total multiplicity";
378  nbins[dim] = 200;
379  min[dim] = 0;
380  max[dim] = 20000;
381  dim++;
382  }
383  else if (fDoV0QA==2)
384  {
385  axistitle[dim] = "V0A+V0C total multiplicity";
386  nbins[dim] = 300;
387  min[dim] = 0;
388  max[dim] = 30000;
389  dim++;
390  }
391 
392  if (fDoEPQA)
393  {
394  axistitle[dim] = "#psi_{EP}";
395  nbins[dim] = 200;
396  min[dim] = -TMath::Pi();
397  max[dim] = TMath::Pi();
398  dim++;
399  }
400  }
401 
402  if (fParticleCollArray.size() > 0)
403  {
404  axistitle[dim] = "No. of tracks";
406  {
407  nbins[dim] = 3000;
408  min[dim] = -0.5;
409  max[dim] = 6000-0.5;
410  }
411  else {
412  nbins[dim] = 200;
413  min[dim] = 0;
414  max[dim] = 200;
415  }
416  dim++;
417 
418  axistitle[dim] = "#it{p}_{T,track}^{leading} (GeV/c)";
419  nbins[dim] = nPtBins;
420  min[dim] = 0;
421  max[dim] = fMaxPt;
422  dim++;
423 
425  {
426  axistitle[dim] = "#eta_{track}^{leading}";
427  nbins[dim] = 100;
428  min[dim] = -1;
429  max[dim] = 1;
430  dim++;
431 
432  axistitle[dim] = "#phi_{track}^{leading}";
433  nbins[dim] = 100;
434  min[dim] = 0;
435  max[dim] = TMath::TwoPi();
436  dim++;
437  }
438  }
439 
440  if (fClusterCollArray.size() > 0)
441  {
442  axistitle[dim] = "No. of clusters";
443 
445  {
446  nbins[dim] = 2000;
447  min[dim] = -0.5;
448  max[dim] = 4000-0.5;
449  }
450  else
451  {
452  nbins[dim] = 200;
453  min[dim] = 0;
454  max[dim] = 200;
455  }
456  dim++;
457 
458  if (fSeparateEMCalDCal)
459  {
460  axistitle[dim] = "#it{E}_{EMCal cluster}^{leading} (GeV)";
461  nbins[dim] = nPtBins;
462  min[dim] = 0;
463  max[dim] = fMaxPt;
464  dim++;
465 
466  axistitle[dim] = "#it{E}_{DCal cluster}^{leading} (GeV)";
467  nbins[dim] = nPtBins;
468  min[dim] = 0;
469  max[dim] = fMaxPt;
470  dim++;
471 
473  {
474  axistitle[dim] = "#eta_{EMCal cluster}^{leading}";
475  nbins[dim] = 100;
476  min[dim] = -1;
477  max[dim] = 1;
478  dim++;
479 
480  axistitle[dim] = "#phi_{EMCal cluster}^{leading}";
481  nbins[dim] = 100;
482  min[dim] = 0;
483  max[dim] = TMath::TwoPi();
484  dim++;
485 
486  axistitle[dim] = "#eta_{DCal cluster}^{leading}";
487  nbins[dim] = 100;
488  min[dim] = -1;
489  max[dim] = 1;
490  dim++;
491 
492  axistitle[dim] = "#phi_{DCal cluster}^{leading}";
493  nbins[dim] = 100;
494  min[dim] = 0;
495  max[dim] = TMath::TwoPi();
496  dim++;
497  }
498  }
499  else
500  {
501  axistitle[dim] = "#it{E}_{cluster}^{leading} (GeV)";
502  nbins[dim] = nPtBins;
503  min[dim] = 0;
504  max[dim] = fMaxPt;
505  dim++;
506 
508  {
509  axistitle[dim] = "#eta_{cluster}^{leading}";
510  nbins[dim] = 100;
511  min[dim] = -1;
512  max[dim] = 1;
513  dim++;
514 
515  axistitle[dim] = "#phi_{cluster}^{leading}";
516  nbins[dim] = 100;
517  min[dim] = 0;
518  max[dim] = TMath::TwoPi();
519  dim++;
520  }
521  }
522  }
523 
524  if (!fCaloCellsName.IsNull())
525  {
526  axistitle[dim] = "No. of cells";
527 
529  {
530  nbins[dim] = 5000;
531  min[dim] = -0.5;
532  max[dim] = 10000-0.5;
533  }
534  else
535  {
536  nbins[dim] = 500;
537  min[dim] = -0.5;
538  max[dim] = 500-0.5;
539  }
540 
541  dim++;
542  }
543 
544  THnSparse* hn = fHistManager.CreateTHnSparse("fHistEventQA","fHistEventQA",dim,nbins,min,max);
545  for (Int_t i = 0; i < dim; i++)
546  hn->GetAxis(i)->SetTitle(axistitle[i]);
547 
549 
550  PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
551 }
552 
553 //________________________________________________________________________
555 {
556  if(fDebug==1)cout<<"Inside of: AliAnalysisTaskEmcalRun2QA::ExecOnce()"<<endl;
557 
558  if (fClusterCollArray.size() == 0 && fCaloCellsName.IsNull())
559  {
560  fNeedEmcalGeom = kFALSE;
561  }
562  else
563  {
564  fNeedEmcalGeom = kTRUE;
565  }
566 
567  //produces the fGeom Object, ...
569 
570  if (fDoV0QA)
571  {
572  fVZERO = InputEvent()->GetVZEROData();
573  if (!fVZERO)
574  {
575  AliError("AliVVZERO not available");
576  }
577  }
578  //Initialize the calo utils object with the current settings of the event
579  if(fCaloUtils)
580  {
581  cout<<"- - - - - fCaloUtils initialized - - - "<<endl;
582  // Set geometry matrices before filling arrays, in case recalibration/position calculation etc is needed
583  fCaloUtils->AccessGeometry(InputEvent()); // InputEvent()->GetRunNumber()
584  // Set the AODB calibration, bad channels etc. parameters at least once
585  fCaloUtils->AccessOADB(InputEvent());
586  //apparently not initialized correctly like eg in AliEMCALGeometry!
587  fCaloUtils->SetNumberOfSuperModulesUsed(fGeom->GetNumberOfSuperModules());
588  cout<<"- - - - - fCaloUtils initialized end- - - "<<endl;
589  }
590  else
591  {
592  cout<<"- - - - - no fCaloUtils initialized - - - "<<endl;
593  }
594 
595  if(fNumOfSuperMod != fGeom->GetNumberOfSuperModules())
596  {
597  cout<<"Will cause major error in histogram arrays - change hard coded value of fNumOfSuperMod!"<<endl;
598  cout<<"Hard coded value is: "<<fNumOfSuperMod<<", value from fGeom is: "<<fGeom->GetNumberOfSuperModules()<<endl;
599  }
600 
601  //Modules - super modules max rows max cells
602  cout<<"Info: "<<endl;
603  //number of rows
604  Int_t NRows = AliEMCALGeoParams::fgkEMCALRows;
605  Int_t NCol = AliEMCALGeoParams::fgkEMCALCols;
606  Int_t NSMod = AliEMCALGeoParams::fgkEMCALModules;
607  cout<<" Number of super modules :"<<fGeom->GetNumberOfSuperModules()<<", "<<fCaloUtils->GetNumberOfSuperModulesUsed()<<"(fGeom, fCaloUtils)"<<endl;
608  cout<<"number of mo from geoparams: "<<NSMod<<endl;
609  cout<<"number of col from geoparams: "<<NCol<<endl;
610  cout<<"number of row from geoparams: "<<NRows<<endl;
611  cout<<"Number of rows: "<<NRows<<", rowsMax: "<<fNMaxRowsAbs<<", Number of collumns:"<<NCol<<", colMax: "<<fNMaxColsAbs<<endl;
612  //Number of rows: 24, rowsMax: 120, Number of collumns:48, colMax: 96
613 
614  //check that! (also for DCal!!!)
615 }
616 
617 //________________________________________________________________________
619 {
620  if(fDebug==1)cout<<"Inside of: AliAnalysisTaskEmcalRun2QA::RetrieveEventObjects()"<<endl;
621 
623 
624  if (!fCentMethod2.IsNull() || !fCentMethod3.IsNull())
625  {
626  if (fBeamType == kAA || fBeamType == kpA )
627  {
628  AliCentrality *aliCent = InputEvent()->GetCentrality();
629  if (aliCent)
630  {
631  if (!fCentMethod2.IsNull())
632  fCent2 = aliCent->GetCentralityPercentile(fCentMethod2);
633  if (!fCentMethod3.IsNull())
634  fCent3 = aliCent->GetCentralityPercentile(fCentMethod3);
635  }
636  }
637  }
638 
639  if (fVZERO)
640  {
641  fV0ATotMult = AliESDUtils::GetCorrV0A(fVZERO->GetMTotV0A(),fVertex[2]);
642  fV0CTotMult = AliESDUtils::GetCorrV0C(fVZERO->GetMTotV0C(),fVertex[2]);
643  }
644 
645  return kTRUE;
646 }
647 
648 //________________________________________________________________________
650 {
651  if(fDebug==1)cout<<"Inside of: AliAnalysisTaskEmcalRun2QA::FillHistograms()"<<endl;
652 
653  // Fill histograms.
654  EventQA_t eventQA;
655 
656  DoClusterLoop();
657  AliDebug(2,Form("%d clusters found in EMCal and %d in DCal", fNTotClusters[0], fNTotClusters[1]));
658  for (Int_t i = 0; i < 2; i++)
659  {
660  eventQA.fMaxCluster[i] = fLeadingCluster[i];
661  }
662 
663  if (fCaloCells)
664  {
665  eventQA.fNCells = DoCellLoop();
666  AliDebug(2,Form("%d cells found in the event", eventQA.fNCells));
667  }
668 
669  eventQA.fCent = fCent;
670  eventQA.fCent2 = fCent2;
671  eventQA.fCent3 = fCent3;
672  eventQA.fV0A = fV0ATotMult;
673  eventQA.fV0C = fV0CTotMult;
674  eventQA.fEP = fEPV0;
675  eventQA.fNClusters[0] = fNTotClusters[0];
676  eventQA.fNClusters[1] = fNTotClusters[1];
677 
678  FillEventQAHisto(eventQA);
679 
680  return kTRUE;
681 }
682 
683 //________________________________________________________________________
685 {
686  if(fDebug==1)cout<<"Inside of: AliAnalysisTaskEmcalRun2QA::FillEventQAHisto()"<<endl;
687  Double_t contents[40]={0};
688 
689  Int_t globalNclusters = eventQA.fNClusters[0] + eventQA.fNClusters[1];
690 
691  AliTLorentzVector globalMaxCluster = eventQA.fMaxCluster[0].E() > eventQA.fMaxCluster[1].E() ? eventQA.fMaxCluster[0] : eventQA.fMaxCluster[1];
692 
693  THnSparse* histEventQA = static_cast<THnSparse*>(fHistManager.FindObject("fHistEventQA"));
694 
695  for (Int_t i = 0; i < histEventQA->GetNdimensions(); i++)
696  {
697  TString title(histEventQA->GetAxis(i)->GetTitle());
698  if (title=="Centrality %")
699  contents[i] = eventQA.fCent;
700  else if (title==Form("Centrality %s %%", fCentMethod2.Data()))
701  contents[i] = eventQA.fCent2;
702  else if (title==Form("Centrality %s %%", fCentMethod3.Data()))
703  contents[i] = eventQA.fCent3;
704  else if (title=="V0A total multiplicity")
705  contents[i] = eventQA.fV0A;
706  else if (title=="V0C total multiplicity")
707  contents[i] = eventQA.fV0C;
708  else if (title=="V0A+V0C total multiplicity")
709  contents[i] = eventQA.fV0A+eventQA.fV0C;
710  else if (title=="#psi_{RP}")
711  contents[i] = eventQA.fEP;
712  else if (title=="No. of clusters")
713  contents[i] = globalNclusters;
714  else if (title=="No. of cells")
715  contents[i] = eventQA.fNCells;
716  else if (title=="#it{p}_{T,track}^{leading} (GeV/c)")
717  contents[i] = eventQA.fMaxTrack.Pt();
718  else if (title=="#eta_{track}^{leading}")
719  contents[i] = eventQA.fMaxTrack.Eta();
720  else if (title=="#phi_{track}^{leading}")
721  contents[i] = eventQA.fMaxTrack.Phi_0_2pi();
722  else if (title=="#it{E}_{cluster}^{leading} (GeV)")
723  contents[i] = globalMaxCluster.E();
724  else if (title=="#eta_{cluster}^{leading}")
725  contents[i] = globalMaxCluster.Eta();
726  else if (title=="#phi_{cluster}^{leading}")
727  contents[i] = globalMaxCluster.Phi();
728  else if (title=="#it{E}_{EMCal cluster}^{leading} (GeV)")
729  contents[i] = eventQA.fMaxCluster[0].E();
730  else if (title=="#eta_{EMCal cluster}^{leading}")
731  contents[i] = eventQA.fMaxCluster[0].Phi_0_2pi();
732  else if (title=="#phi_{EMCal cluster}^{leading}")
733  contents[i] = eventQA.fMaxCluster[0].Eta();
734  else if (title=="#it{E}_{DCal cluster}^{leading} (GeV)")
735  contents[i] = eventQA.fMaxCluster[1].E();
736  else if (title=="#phi_{DCal cluster}^{leading}")
737  contents[i] = eventQA.fMaxCluster[1].Phi_0_2pi();
738  else if (title=="#eta_{DCal cluster}^{leading}")
739  contents[i] = eventQA.fMaxCluster[1].Eta();
740  else
741  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
742  }
743 
744  histEventQA->Fill(contents);
745 }
746 
747 //________________________________________________________________________
749 {
750  // Do cell loop.
751  if (!fCaloCells) return 0;
752 
753  if(fNumOfSuperMod != fGeom->GetNumberOfSuperModules())
754  {
755  cout<<"Will cause major error in histogram arrays - change hard coded value of fNumOfSuperMod!"<<endl;
756  return 0;
757  }
758 
759  TString histname = TString::Format("%s/fHistCellsAbsIdEnergy_%d", fCaloCellsName.Data(), fCentBin);
760 
761  const Int_t ncells = fCaloCells->GetNumberOfCells();
762  Int_t nAccCells = 0;
763  Double_t CellEta,CellPhi;
764  Float_t amp = -1;
765  Int_t absId = -1;
766 
767  //will be set by GetModuleNumberCellIndexesAbsCaloMap
768  Int_t icol = -1;
769  Int_t icolAbs = -1;
770  Int_t irow = -1;
771  Int_t irowAbs = -1;
772  Int_t iRCU = -1;
773 
774  for (Int_t pos = 0; pos < ncells; pos++)
775  {
776  amp = fCaloCells->GetAmplitude(pos);
777  if (amp < fCellEnergyCut) continue;
778  absId = fCaloCells->GetCellNumber(pos);
779  fGeom->EtaPhiFromIndex(absId,CellEta,CellPhi);
780  Int_t smNo = fGeom->GetSuperModuleNumber(absId);
781 
782  //get the cell row and collumn
783  //enum detector { kEMCAL = 0, kPHOS = 1, kCTS = 2, kDCAL = 3, kDCALPHOS = 4 };
784  fCaloUtils->GetModuleNumberCellIndexesAbsCaloMap(absId,0,icol,irow,iRCU,icolAbs,irowAbs);
785  if(icolAbs> fNMaxColsAbs || irowAbs>fNMaxRowsAbs)
786  {
787  cout<<" Problem! wrong calculated number of max col and max rows"<<endl;
788  cout<<"current col: "<<icolAbs<<", max col"<<fNMaxColsAbs<<endl;
789  cout<<"current row: "<<irowAbs<<", max row"<<fNMaxRowsAbs<<endl;
790  }
791 
792  fHistCellIDvsE[smNo] ->Fill(absId,amp);
793  fHistCellIDvsELow[smNo]->Fill(absId,amp);
794  fHistCellEtaPhi[smNo] ->Fill(icolAbs,irowAbs); //icolAbs,irowAbs
795 
796  //all centralities
797  fHistCellIDvsE[fNumOfSuperMod] ->Fill(absId,amp);
798  fHistCellIDvsELow[fNumOfSuperMod]->Fill(absId,amp);
799  fHistCellEtaPhi[fNumOfSuperMod] ->Fill(icolAbs,irowAbs);
800 
801 
802  //IsInDCAL
803  //IsInEMCAL
804 
805  nAccCells++;
806  }
807  return nAccCells;
808 }
809 //________________________________________________________________________
811 {
812  // Do cluster loop.
813  TString histname;
814 
815  memset(fNTotClusters, 0, sizeof(Int_t)*2);
816  for (Int_t i = 0; i < 2; i++) fLeadingCluster[i].SetPxPyPzE(0,0,0,0);
817 
818  for (auto cont_it : fClusterCollArray) {
819  AliClusterContainer* clusters = cont_it.second;
820  // Cluster loop
822  //will be set by GetModuleNumberCellIndexesAbsCaloMap
823  Int_t icol = -1;
824  Int_t icolAbs = -1;
825  Int_t irow = -1;
826  Int_t irowAbs = -1;
827  Int_t iRCU = -1;
828 
829  for (AliClusterIterableMomentumContainer::iterator it = itcont.begin(); it != itcont.end(); it++)
830  {
831  UInt_t rejectionReason = 0;
832 //doesn't work if (!clusters->AcceptCluster(clusters->GetCurrentID(), rejectionReason))
833 //doesn't work if (!clusters->AcceptCluster(it->second->GetID(), rejectionReason))
834  if (!clusters->AcceptCluster(it->second, rejectionReason))
835  {
836  histname = TString::Format("%s/fHistRejectionReason_%d", clusters->GetArrayName().Data(), fCentBin);
837  fHistManager.FillTH2(histname, clusters->GetRejectionReasonBitPosition(rejectionReason), it->first.E());
838  continue;
839  }
840 
841  Float_t pos[3]={0};
842  it->second->GetPosition(pos);
843 
844  //- - - - - - - - - - - - - - - - - - - - - - - - - -
845  Int_t NCells = it->second->GetNCells();
846  Double_t FirstCellEnergy=0;
847  Double_t FirstCellID=0;
848  Double_t CellEta,CellPhi;
849  Double_t ClusterRawEnergy = it->second->E();
850 
851  for (Int_t iCell = 0; iCell < NCells; iCell++)
852  {
853  Int_t cellId = it->second->GetCellAbsId(iCell);
854  Double_t cellFrac = it->second->GetCellAmplitudeFraction(iCell);
855  Double_t Eseed = fCaloCells->GetCellAmplitude(cellId);
856  if(iCell==0)
857  {
858  FirstCellEnergy=Eseed;
859  FirstCellID =cellId;
860  }
861  if(Eseed>FirstCellEnergy)cout<<"Problem: Energy bigger than first."<<iCell<<endl;
862  }
863  //use the id of the first cell (should be the leading cell)
864  Int_t smNo = fGeom->GetSuperModuleNumber(FirstCellID);
865  fGeom->EtaPhiFromIndex(FirstCellID,CellEta,CellPhi);
866  //get the cell row and collumn
867  //enum detector { kEMCAL = 0, kPHOS = 1, kCTS = 2, kDCAL = 3, kDCALPHOS = 4 };
868  fCaloUtils->GetModuleNumberCellIndexesAbsCaloMap(FirstCellID,0,icol,irow,iRCU,icolAbs,irowAbs);
869 
870  //number od cells in cluster vs module (Maybe eta-phi raender weniger etc)
871  fHistClusterIDvsE[smNo] ->Fill(FirstCellID,ClusterRawEnergy);
872  fHistClusterIDvsELow[smNo]->Fill(FirstCellID,ClusterRawEnergy);
873  fHistClusterEtaPhi[smNo] ->Fill(icolAbs,irowAbs);
874 
875  //all centralities/modules
876  fHistClusterIDvsE[fNumOfSuperMod] ->Fill(FirstCellID,ClusterRawEnergy);
877  fHistClusterIDvsELow[fNumOfSuperMod]->Fill(FirstCellID,ClusterRawEnergy);
878  fHistClusterEtaPhi[fNumOfSuperMod] ->Fill(icolAbs,irowAbs);
879 
880  //one dimensional // would be nice to have two dimensional.q
881  fHistNrCellsInCluster[smNo] ->Fill(NCells);
882  //- - - - - - - - - - - - - - - - - - - - - - - - - -
883 
884 
885  histname = TString::Format("%s/fHistClusPosition_%d", clusters->GetArrayName().Data(), fCentBin);
886  fHistManager.FillTH3(histname, pos[0], pos[1], pos[2]);
887 
888  Double_t phi = it->first.Phi_0_2pi();
889 
890  Int_t isDcal = Int_t(phi > fgkEMCalDCalPhiDivide);
891 
892  histname = TString::Format("%s/fHistClusPhiEtaEnergy_%d", clusters->GetArrayName().Data(), fCentBin);
893  fHistManager.FillTH3(histname, it->first.Eta(), it->first.Phi_0_2pi(), it->first.E());
894 
895  if (fLeadingCluster[isDcal].E() < it->first.E()) fLeadingCluster[isDcal] = it->first;
896 
897 
898  histname = TString::Format("%s/fHistClusMCEnergyFraction_%d", clusters->GetArrayName().Data(), fCentBin);
899  if (fHistManager.FindObject(histname))
900  {
901  fHistManager.FillTH1(histname, it->second->GetMCEnergyFraction());
902  }
903 
904  histname = TString::Format("%s/fHistClusEnergy_%d", clusters->GetArrayName().Data(), fCentBin);
905  fHistManager.FillTH1(histname, it->second->E());
906 
907  if (it->second->GetNonLinCorrEnergy() > 0.)
908  {
909  histname = TString::Format("%s/fHistClusNonLinCorrEnergy_%d", clusters->GetArrayName().Data(), fCentBin);
910  fHistManager.FillTH1(histname, it->second->GetNonLinCorrEnergy());
911  }
912 
913  if (it->second->GetHadCorrEnergy() > 0.)
914  {
915  histname = TString::Format("%s/fHistClusHadCorrEnergy_%d", clusters->GetArrayName().Data(), fCentBin);
916  fHistManager.FillTH1(histname, it->second->GetHadCorrEnergy());
917  }
918 
919  fNTotClusters[isDcal]++;
920  }
921  }
922 }
Float_t fPtBinWidth
Histogram pt bin width.
THistManager fHistManager
Histogram manager.
THashList * CreateHistoGroup(const char *groupname)
Create a new group of histograms within a parent group.
AliEMCALGeometry * fGeom
!emcal geometry
EBeamType_t fBeamType
!event beam type
double Double_t
Definition: External.C:58
TList * fOutputList_Cell
Output list for Event QA histograms.
Definition: External.C:236
const char * title
Definition: MakeQAPdf.C:27
AliCalorimeterUtils * fCaloUtils
Cells in cluster per supermodule.
Int_t fNTotClusters[2]
!Total number of accepted clusters in current event (DCal/EMCal)
Double_t fV0CTotMult
!Event V0C total multiplicity
Double_t fEPV0
!event plane V0
TH2 ** fHistClusterEtaPhi
Cluster QA histogram.
bidirectional stl iterator over the EMCAL iterable container
Int_t fDoLeadingObjectPosition
Add axis for leading object position (eta-phi)
Float_t fMaxPt
Histogram pt limit.
TH2 ** fHistCellIDvsE
Output list for Cluster QA histograms.
Int_t fNMaxRowsAbs
Number of rows in one supermodule.
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
Fill a 2D histogram within the container.
void FillTH3(const char *hname, double x, double y, double z, double weight=1., Option_t *opt="")
Fill a 3D histogram within the container.
AliVVZERO * fVZERO
!Event V0 object
TH2 ** fHistCellEtaPhi
Cell QA histogram.
TString fCaloCellsName
name of calo cell collection
Int_t fDoV0QA
Add V0 QA histograms.
Float_t fCellEnergyCut
Energy cell cut.
static Double_t fgkEMCalDCalPhiDivide
phi value used to distinguish between DCal and EMCal
EBeamType_t fForceBeamType
forced beam type
const Int_t nPtBins
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
Create a new TH2 within the container.
TObject * FindObject(const char *name) const
Find an object inside the container.
int Int_t
Definition: External.C:63
Double_t fCent2
!Event centrality with method 2
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
float Float_t
Definition: External.C:68
TH1 ** fHistNrCellsInCluster
Cluster QA histogram.
Base task in the EMCAL framework (lighter version of AliAnalysisTaskEmcal)
Double_t Phi_0_2pi() const
Double_t fVertex[3]
!event vertex
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
Double_t fV0ATotMult
!Event V0A total multiplicity
virtual Bool_t AcceptCluster(Int_t i, UInt_t &rejectionReason) const
Implementation of a task to perform basic QA on EMCal and DCal clusters and cells.
Int_t GetNumberOfSuperModulesUsed() const
Int_t fCentBin
!event centrality bin
Int_t fNMaxColsAbs
Number of columns in one supermodule.
AliTLorentzVector fLeadingCluster[2]
!Leading cluster in current event (EMCal/DCal)
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
void SetNumberOfSuperModulesUsed(Int_t nSM)
Float_t fRtoD
Transformation of rad to deg.
TString fCentMethod2
Centrality method 2.
const AliClusterIterableMomentumContainer all_momentum() const
TH2 ** fHistClusterIDvsELow
Cluster QA histogram.
Definition: External.C:220
std::map< std::string, AliParticleContainer * > fParticleCollArray
particle/track collection array
Double_t fCent3
!Event centrality with method 3
TH2 ** fHistClusterIDvsE
Cell QA histogram.
Int_t fNumOfSuperMod
Number of supermodules.
TH2 ** fHistCellIDvsELow
Cell QA histogram.
TList * fOutputList_Cluster
Output list for Cell QA histograms.
TString fCentMethod3
Centrality method 3.
const Int_t nbins
bool Bool_t
Definition: External.C:53
Class with utils specific to calorimeter clusters/cells.
void AccessGeometry(AliVEvent *inputEvent)
Bool_t fSeparateEMCalDCal
Separate EMCal from DCal in QA plots.
Bool_t fDebug
Prints out some extra information.
std::map< std::string, AliClusterContainer * > fClusterCollArray
cluster collection array
THnSparse * CreateTHnSparse(const char *name, const char *title, int ndim, const int *nbins, const double *min, const double *max, Option_t *opt="")
Create a new THnSparse within the container.
void AccessOADB(AliVEvent *event)
Container structure for EMCAL clusters.
Double_t fCent
!event centrality
Int_t fDoEPQA
Add event plane QA histograms.
Definition: External.C:196
void FillEventQAHisto(const EventQA_t &eventQA)
TH3 * CreateTH3(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, int nbinsz, double zmin, double zmax, Option_t *opt="")
Create a new TH2 within the container.
Bool_t fNeedEmcalGeom
whether or not the task needs the emcal geometry
Int_t GetModuleNumberCellIndexesAbsCaloMap(Int_t absId, Int_t calo, Int_t &icol, Int_t &irow, Int_t &iRCU, Int_t &icolAbs, Int_t &irowAbs) const
Declaration of class AliAnalysisTaskEmcalRun2QA.