AliPhysics  ec707b8 (ec707b8)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
AliAnalysisTaskEMCALPi0CalibSelection.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 
4  * Permission to use, copy, modify and distribute this software and its *
5  * documentation strictly for non-commercial purposes is hereby granted *
6  * without fee, provided that the above copyright notice appears in all *
7  * copies and that both the copyright notice and this permission notice *
8  * appear in the supporting documentation. The authors make no claims *
9  * about the suitability of this software for any purpose. It is *
10  * provided "as is" without express or implied warranty. *
11  **************************************************************************/
12 
13 // Root
14 #include <TRefArray.h>
15 #include <TList.h>
16 #include <TH1F.h>
17 #include <TGeoManager.h>
18 #include <TFile.h>
19 
20 // AliRoot
22 #include "AliAODEvent.h"
23 #include "AliESDEvent.h"
24 #include "AliEMCALGeometry.h"
25 #include "AliVCluster.h"
26 #include "AliVCaloCells.h"
27 #include "AliEMCALRecoUtils.h"
28 #include "AliOADBContainer.h"
29 
33 
36 //______________________________________________________________________________________________
38 AliAnalysisTaskSE(),
39 fEMCALGeo(0x0), fLoadMatrices(0),
40 fEMCALGeoName("EMCAL_COMPLETE12SMV1_DCAL_8SM"),
41 fTriggerName("EMC"),
42 fRecoUtils(new AliEMCALRecoUtils),
43 fOADBFilePath(""), fCalibFilePath(""),
44 fCorrectClusters(kFALSE), fRecalPosition(kTRUE),
45 fCaloClustersArr(0x0), fEMCALCells(0x0),
46 //fCuts(0x0),
47 fOutputContainer(0x0),
48 fVertex(), fFilteredInput(kFALSE),
49 fImportGeometryFromFile(1), fImportGeometryFilePath(""),
50 fEmin(0.5), fEmax(15.),
51 fL0min(0.01), fL0max(0.5),
52 fDTimeCut(100.), fTimeMax(1000000), fTimeMin(-1000000),
53 fAsyCut(1.), fMinNCells(2), fGroupNCells(0),
54 fLogWeight(4.5), fSameSM(kFALSE),
55 fNMaskCellColumns(11), fMaskCellColumns(0x0),
56 fInvMassCutMin(110.), fInvMassCutMax(160.),
57 // Histograms binning
58 fNbins(300),
59 fMinBin(0.), fMaxBin(300.),
60 fNTimeBins(1000), fMinTimeBin(0.), fMaxTimeBin(1000.),
61 // Temporal
62 fMomentum1(), fMomentum2(), fMomentum12(),
63 // Histograms
64 fHmgg(0x0), fHmggDifferentSM(0x0),
65 fHmggMaskFrame(0x0), fHmggDifferentSMMaskFrame(0x0),
66 fHOpeningAngle(0x0), fHOpeningAngleDifferentSM(0x0),
67 fHAsymmetry(0x0), fHAsymmetryDifferentSM(0x0),
68 fhNEvents(0x0),
69 fhClusterTime(0x0), fhClusterPairDiffTime(0x0)
70 {
71  for(Int_t iMod=0; iMod < AliEMCALGeoParams::fgkEMCALModules; iMod++)
72  {
73  for(Int_t iX=0; iX<24; iX++)
74  {
75  for(Int_t iZ=0; iZ<48; iZ++)
76  {
77  fHmpi0[iMod][iZ][iX] = 0 ;
78  }
79  }
80  }
81 
82  fVertex[0]=fVertex[1]=fVertex[2]=-1000;
83 
84  fHTpi0[0]= 0 ;
85  fHTpi0[1]= 0 ;
86  fHTpi0[2]= 0 ;
87  fHTpi0[3]= 0 ;
88 
90  fMaskCellColumns[0] = 6 ; fMaskCellColumns[1] = 7 ; fMaskCellColumns[2] = 8 ;
91  fMaskCellColumns[3] = 35; fMaskCellColumns[4] = 36; fMaskCellColumns[5] = 37;
92  fMaskCellColumns[6] = 12+AliEMCALGeoParams::fgkEMCALCols; fMaskCellColumns[7] = 13+AliEMCALGeoParams::fgkEMCALCols;
93  fMaskCellColumns[8] = 40+AliEMCALGeoParams::fgkEMCALCols; fMaskCellColumns[9] = 41+AliEMCALGeoParams::fgkEMCALCols;
94  fMaskCellColumns[10]= 42+AliEMCALGeoParams::fgkEMCALCols;
95 
96  for(Int_t iSMPair = 0; iSMPair < AliEMCALGeoParams::fgkEMCALModules/2; iSMPair++)
97  {
98  fHmggPairSameSectorSM[iSMPair] = 0;
99  fHmggPairSameSectorSMMaskFrame[iSMPair] = 0;
101  }
102 
103  for(Int_t iSMPair = 0; iSMPair < AliEMCALGeoParams::fgkEMCALModules-2; iSMPair++)
104  {
105  fHmggPairSameSideSM[iSMPair] = 0;
106  fHmggPairSameSideSMMaskFrame[iSMPair] = 0;
107  fhClusterPairDiffTimeSameSide[iSMPair] = 0;
108  }
109 
110  for(Int_t iSM = 0; iSM < AliEMCALGeoParams::fgkEMCALModules; iSM++)
111  {
112  fHmggSM[iSM] = 0;
113  fHmggSM_Zone1[iSM] = 0;
114  fHmggSM_Zone2[iSM] = 0;
115  fHmggSM_Zone3[iSM] = 0;
116  fHmggSM_Zone4[iSM] = 0;
117  fHmggSM_Zone5[iSM] = 0;
118  fHmggSM_Zone6[iSM] = 0;
119  fHmggSM_Zone7[iSM] = 0;
120  fHmggSMMaskFrame[iSM] = 0;
121  fHOpeningAngleSM[iSM] = 0;
122  fHOpeningAnglePairSM[iSM] = 0;
123  fHAsymmetrySM[iSM] = 0;
124  fHAsymmetryPairSM[iSM] = 0;
125  fhTowerDecayPhotonHit[iSM] = 0;
126  fhTowerDecayPhotonEnergy[iSM] = 0;
129  fMatrix[iSM] = 0x0;
130  fhClusterTimeSM[iSM] = 0;
132  }
133 }
134 
140 //______________________________________________________________________________________________
142 AliAnalysisTaskSE(name),
143 fEMCALGeo(0x0), fLoadMatrices(0),
144 fEMCALGeoName("EMCAL_COMPLETE12SMV1_DCAL_8SM"),
145 fTriggerName("EMC"),
146 fRecoUtils(new AliEMCALRecoUtils),
147 fOADBFilePath(""), fCalibFilePath(""),
148 fCorrectClusters(kFALSE), fRecalPosition(kTRUE),
149 fCaloClustersArr(0x0), fEMCALCells(0x0),
150 //fCuts(0x0),
151 fOutputContainer(0x0),
152 fVertex(), fFilteredInput(kFALSE),
153 fImportGeometryFromFile(1), fImportGeometryFilePath(""),
154 fEmin(0.5), fEmax(15.),
155 fL0min(0.01), fL0max(0.5),
156 fDTimeCut(100.), fTimeMax(1000000), fTimeMin(-1000000),
157 fAsyCut(1.), fMinNCells(2), fGroupNCells(0),
158 fLogWeight(4.5), fSameSM(kFALSE),
159 fNMaskCellColumns(11), fMaskCellColumns(0x0),
160 fInvMassCutMin(110.), fInvMassCutMax(160.),
161 // Histograms binning
162 fNbins(300),
163 fMinBin(0.), fMaxBin(300.),
164 fNTimeBins(1000), fMinTimeBin(0.), fMaxTimeBin(1000.),
165 // Temporal
166 fMomentum1(), fMomentum2(), fMomentum12(),
167 // Histograms
168 fHmgg(0x0), fHmggDifferentSM(0x0),
169 fHmggMaskFrame(0x0), fHmggDifferentSMMaskFrame(0x0),
170 fHOpeningAngle(0x0), fHOpeningAngleDifferentSM(0x0),
171 fHAsymmetry(0x0), fHAsymmetryDifferentSM(0x0),
172 fhNEvents(0x0),
173 fhClusterTime(0x0), fhClusterPairDiffTime(0x0)
174 {
175  for(Int_t iMod=0; iMod < AliEMCALGeoParams::fgkEMCALModules; iMod++)
176  {
177  for(Int_t iX=0; iX<24; iX++)
178  {
179  for(Int_t iZ=0; iZ<48; iZ++)
180  {
181  fHmpi0[iMod][iZ][iX] = 0 ;
182  }
183  }
184  }
185 
186  fVertex[0]=fVertex[1]=fVertex[2]=-1000;
187 
188  fHTpi0[0]= 0 ;
189  fHTpi0[1]= 0 ;
190  fHTpi0[2]= 0 ;
191  fHTpi0[3]= 0 ;
192 
193  fMaskCellColumns = new Int_t[fNMaskCellColumns];
194  fMaskCellColumns[0] = 6 ; fMaskCellColumns[1] = 7 ; fMaskCellColumns[2] = 8 ;
195  fMaskCellColumns[3] = 35; fMaskCellColumns[4] = 36; fMaskCellColumns[5] = 37;
196  fMaskCellColumns[6] = 12+AliEMCALGeoParams::fgkEMCALCols; fMaskCellColumns[7] = 13+AliEMCALGeoParams::fgkEMCALCols;
197  fMaskCellColumns[8] = 40+AliEMCALGeoParams::fgkEMCALCols; fMaskCellColumns[9] = 41+AliEMCALGeoParams::fgkEMCALCols;
198  fMaskCellColumns[10]= 42+AliEMCALGeoParams::fgkEMCALCols;
199 
200  for(Int_t iSMPair = 0; iSMPair < AliEMCALGeoParams::fgkEMCALModules/2; iSMPair++)
201  {
202  fHmggPairSameSectorSM[iSMPair] = 0;
203  fHmggPairSameSectorSMMaskFrame[iSMPair] = 0;
205  }
206 
207  for(Int_t iSMPair = 0; iSMPair < AliEMCALGeoParams::fgkEMCALModules-2; iSMPair++)
208  {
209  fHmggPairSameSideSM[iSMPair] = 0;
210  fHmggPairSameSideSMMaskFrame[iSMPair] = 0;
211  fhClusterPairDiffTimeSameSide[iSMPair] = 0;
212  }
213 
214  for(Int_t iSM = 0; iSM < AliEMCALGeoParams::fgkEMCALModules; iSM++)
215  {
216  fHmggSM[iSM] = 0;
217  fHmggSM_Zone1[iSM] = 0;
218  fHmggSM_Zone2[iSM] = 0;
219  fHmggSM_Zone3[iSM] = 0;
220  fHmggSM_Zone4[iSM] = 0;
221  fHmggSM_Zone5[iSM] = 0;
222  fHmggSM_Zone6[iSM] = 0;
223  fHmggSM_Zone7[iSM] = 0;
224  fHmggSMMaskFrame[iSM] = 0;
225  fHOpeningAngleSM[iSM] = 0;
226  fHOpeningAnglePairSM[iSM] = 0;
227  fHAsymmetrySM[iSM] = 0;
228  fHAsymmetryPairSM[iSM] = 0;
229  fhTowerDecayPhotonHit[iSM] = 0;
230  fhTowerDecayPhotonEnergy[iSM] = 0;
233  fMatrix[iSM] = 0x0;
234  fhClusterTimeSM[iSM] = 0;
236  }
237 
238  DefineOutput(1, TList::Class());
239 //DefineOutput(2, TList::Class()); // will contain cuts or local params
240 }
241 
244 //_____________________________________________________________________________
246 {
247  if(fOutputContainer)
248  {
249  fOutputContainer->Delete() ;
250  delete fOutputContainer ;
251  }
252 
253  if(fEMCALGeo) delete fEMCALGeo ;
254  if(fRecoUtils) delete fRecoUtils ;
255  if(fNMaskCellColumns) delete [] fMaskCellColumns;
256 }
257 
261 //____________________________________________________________
263 {
264  if(fRecoUtils->GetParticleType()!=AliEMCALRecoUtils::kPhoton)
265  AliFatal(Form("Wrong particle type for cluster position recalculation! = %d\n", fRecoUtils->GetParticleType()));
266 
267  AliDebug(1,Form("It will use fLogWeight %.3f",fLogWeight));
268 
269  Float_t pos[]={0,0,0};
270 
271  for(Int_t iClu=0; iClu < fCaloClustersArr->GetEntriesFast(); iClu++)
272  {
273  AliVCluster *c1 = (AliVCluster *) fCaloClustersArr->At(iClu);
274 
275  Float_t e1i = c1->E(); // cluster energy before correction
276  if (e1i < fEmin) continue;
277  else if (e1i > fEmax) continue;
278 
279  else if (c1->GetNCells() < fMinNCells) continue;
280 
281  else if (c1->GetM02() < fL0min || c1->GetM02() > fL0max) continue;
282 
283  if(fRecoUtils->ClusterContainsBadChannel(fEMCALGeo, c1->GetCellsAbsId(), c1->GetNCells())) continue;
284 
285  if(DebugLevel() > 2)
286  {
287  AliInfo(Form("Std : i %d, E %f, dispersion %f, m02 %f, m20 %f\n",c1->GetID(),c1->E(),c1->GetDispersion(),c1->GetM02(),c1->GetM20()));
288  c1->GetPosition(pos);
289  AliInfo(Form("Std : i %d, x %f, y %f, z %f\n",c1->GetID(), pos[0], pos[1], pos[2]));
290  }
291 
292  // Correct cluster energy and position if requested, and not corrected previously, by default Off
293  if(fRecoUtils->IsRecalibrationOn())
294  {
295  fRecoUtils->RecalibrateClusterEnergy(fEMCALGeo, c1, fEMCALCells);
296  fRecoUtils->RecalculateClusterShowerShapeParameters(fEMCALGeo, fEMCALCells,c1);
297  fRecoUtils->RecalculateClusterPID(c1);
298  }
299 
300  AliDebug(2,Form("Energy: after recalibration %f",c1->E()));
301 
302  // Recalculate cluster position
303  if ( fRecalPosition ) fRecoUtils->RecalculateClusterPosition(fEMCALGeo, fEMCALCells,c1);
304 
305  // Correct Non-Linearity
306  c1->SetE(fRecoUtils->CorrectClusterEnergyLinearity(c1));
307 
308  AliDebug(2,Form("after linearity correction %f",c1->E()));
309 
310  // In case of MC analysis, to match resolution/calibration in real data
311  //c1->SetE(fRecoUtils->SmearClusterEnergy(c1)); // Not needed anymore
312 
313  AliDebug(2,Form("after smearing %f\n",c1->E()));
314 
315  if(DebugLevel() > 2)
316  {
317  AliInfo(Form("Cor : i %d, E %f, dispersion %f, m02 %f, m20 %f\n",c1->GetID(),c1->E(),c1->GetDispersion(),c1->GetM02(),c1->GetM20()));
318  c1->GetPosition(pos);
319  AliInfo(Form("Cor : i %d, x %f, y %f, z %f\n",c1->GetID(), pos[0], pos[1], pos[2]));
320  }
321  } // cluster loop
322 }
323 
327 //__________________________________________________________
329 {
330  Int_t absId1 = -1;
331  Int_t iSupMod1 = -1;
332  Int_t iphi1 = -1;
333  Int_t ieta1 = -1;
334  Int_t absId2 = -1;
335  Int_t iSupMod2 = -1;
336  Int_t iphi2 = -1;
337  Int_t ieta2 = -1;
338  Bool_t shared = kFALSE;
339 
340  Float_t pos[] = {0,0,0};
341 
342  Int_t bc = InputEvent()->GetBunchCrossNumber();
343  Int_t nSM = (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules();
344 
345  for(Int_t iClu=0; iClu<fCaloClustersArr->GetEntriesFast()-1; iClu++)
346  {
347  AliVCluster *c1 = (AliVCluster *) fCaloClustersArr->At(iClu);
348 
349  // Exclude bad channels
350  if(fRecoUtils->ClusterContainsBadChannel(fEMCALGeo, c1->GetCellsAbsId(), c1->GetNCells())) continue;
351 
352  Float_t e1i = c1->E(); // cluster energy before correction
353 
354  if (e1i < fEmin) continue;
355  else if (e1i > fEmax) continue;
356 
357  else if (!fRecoUtils->IsGoodCluster(c1,fEMCALGeo,fEMCALCells,bc)) continue;
358 
359  else if (c1->GetNCells() < fMinNCells) continue;
360 
361  else if (c1->GetM02() < fL0min || c1->GetM02() > fL0max) continue;
362 
363  if(DebugLevel() > 2)
364  {
365  AliInfo(Form("IMA : i %d, E %f, dispersion %f, m02 %f, m20 %f",c1->GetID(),e1i,c1->GetDispersion(),c1->GetM02(),c1->GetM20()));
366  c1->GetPosition(pos);
367  AliInfo(Form("IMA : i %d, x %f, y %f, z %f",c1->GetID(), pos[0], pos[1], pos[2]));
368  }
369 
370  fRecoUtils->GetMaxEnergyCell(fEMCALGeo, fEMCALCells,c1,absId1,iSupMod1,ieta1,iphi1,shared);
371 
372  c1->GetMomentum(fMomentum1,fVertex);
373 
374  // Check if cluster is in fidutial region, not too close to borders
375  Bool_t in1 = fRecoUtils->CheckCellFiducialRegion(fEMCALGeo, c1, fEMCALCells);
376 
377  // Clusters not facing frame structures
378  Bool_t mask1 = MaskFrameCluster(iSupMod1, ieta1);
379  //if(mask1) printf("Reject eta %d SM %d\n",ieta1, iSupMod1);
380 
381  Double_t time1 = c1->GetTOF()*1.e9;
382 
383  if(time1 > fTimeMax || time1 < fTimeMin) continue;
384 
385  fhClusterTime ->Fill(c1->E(),time1);
386  fhClusterTimeSM[iSupMod1]->Fill(c1->E(),time1);
387 
388  // Combine cluster with other clusters and get the invariant mass
389  for (Int_t jClu=iClu+1; jClu < fCaloClustersArr->GetEntriesFast(); jClu++)
390  {
391  AliAODCaloCluster *c2 = (AliAODCaloCluster *) fCaloClustersArr->At(jClu);
392 
393  Float_t e2i = c2->E();
394  if (e2i < fEmin) continue;
395  else if (e2i > fEmax) continue;
396 
397  else if (!fRecoUtils->IsGoodCluster(c2,fEMCALGeo,fEMCALCells,bc))continue;
398 
399  else if (c2->GetNCells() < fMinNCells) continue;
400 
401  else if (c2->GetM02() < fL0min || c2->GetM02() > fL0max) continue;
402 
403  fRecoUtils->GetMaxEnergyCell(fEMCALGeo, fEMCALCells,c2,absId2,iSupMod2,ieta2,iphi2,shared);
404 
405  c2->GetMomentum(fMomentum2,fVertex);
406 
408  Float_t invmass = fMomentum12.M()*1000;
409 
410  //Asimetry cut
411  Float_t asym = TMath::Abs(fMomentum1.E()-fMomentum2.E())/(fMomentum1.E()+fMomentum2.E());
412 
413  if(asym > fAsyCut) continue;
414 
415  //Time cut
416  Double_t time2 = c2->GetTOF()*1.e9;
417 
418  if(time2 > fTimeMax || time2 < fTimeMin) continue;
419 
420  fhClusterPairDiffTime->Fill(fMomentum12.E(),time1-time2);
421  if(TMath::Abs(time1-time2) > fDTimeCut) continue;
422 
423  if(invmass < fMaxBin && invmass > fMinBin )
424  {
425  //Check if cluster is in fidutial region, not too close to borders
426  Bool_t in2 = fRecoUtils->CheckCellFiducialRegion(fEMCALGeo, c2, fEMCALCells);
427 
428  // Clusters not facing frame structures
429  Bool_t mask2 = MaskFrameCluster(iSupMod2, ieta2);
430  //if(mask2) printf("Reject eta %d SM %d\n",ieta2, iSupMod2);
431 
432  if(in1 && in2)
433  {
434  fHmgg->Fill(invmass,fMomentum12.Pt());
435 
436  if(iSupMod1==iSupMod2)
437  {
438  fHmggSM[iSupMod1]->Fill(invmass,fMomentum12.Pt());
439  fhClusterPairDiffTimeSameSM[iSupMod1]->Fill(fMomentum12.E(),time1-time2);
440 
441  //Is in zone number i
442  Bool_t zone1 = IsInZone1(iSupMod1,ieta1,iphi1);
443  Bool_t zone2 = IsInZone2(iSupMod1,ieta1,iphi1);
444  Bool_t zone3 = IsInZone3(iSupMod1,ieta1,iphi1);
445  Bool_t zone4 = IsInZone4(iSupMod1,ieta1,iphi1);
446  Bool_t zone5 = IsInZone5(iSupMod1,ieta1,iphi1);
447  Bool_t zone6 = IsInZone6(iSupMod1,ieta1,iphi1);
448  Bool_t zone7 = IsInZone7(iSupMod1,ieta1,iphi1);
449 
450 
451  if(zone1) fHmggSM_Zone1[iSupMod1]->Fill(invmass,fMomentum12.Pt());
452  if(zone2) fHmggSM_Zone2[iSupMod1]->Fill(invmass,fMomentum12.Pt());
453  if(zone3) fHmggSM_Zone3[iSupMod1]->Fill(invmass,fMomentum12.Pt());
454  if(zone4) fHmggSM_Zone4[iSupMod1]->Fill(invmass,fMomentum12.Pt());
455  if(zone5) fHmggSM_Zone5[iSupMod1]->Fill(invmass,fMomentum12.Pt());
456  if(zone6) fHmggSM_Zone6[iSupMod1]->Fill(invmass,fMomentum12.Pt());
457  if(zone7) fHmggSM_Zone7[iSupMod1]->Fill(invmass,fMomentum12.Pt());
458 
459  }
460  else
461  fHmggDifferentSM ->Fill(invmass,fMomentum12.Pt());
462 
463  // Same sector
464  Int_t j=0;
465  for(Int_t i = 0; i < nSM/2; i++)
466  {
467  j=2*i;
468  if((iSupMod1==j && iSupMod2==j+1) || (iSupMod1==j+1 && iSupMod2==j))
469  {
470  fHmggPairSameSectorSM[i]->Fill(invmass,fMomentum12.Pt());
471  fhClusterPairDiffTimeSameSector[i]->Fill(fMomentum12.E(),time1-time2);
472  }
473  }
474 
475  // Same side
476  for(Int_t i = 0; i < nSM-2; i++)
477  {
478  if((iSupMod1==i && iSupMod2==i+2) || (iSupMod1==i+2 && iSupMod2==i))
479  {
480  fHmggPairSameSideSM[i]->Fill(invmass,fMomentum12.Pt());
481  fhClusterPairDiffTimeSameSide[i]->Fill(fMomentum12.E(),time1-time2);
482  }
483  }
484 
485 
486  if(!mask1 && !mask2)
487  {
488  fHmggMaskFrame->Fill(invmass,fMomentum12.Pt());
489 
490  if(iSupMod1==iSupMod2) fHmggSMMaskFrame[iSupMod1]->Fill(invmass,fMomentum12.Pt());
491  else fHmggDifferentSMMaskFrame ->Fill(invmass,fMomentum12.Pt());
492 
493  // Same sector
494  j=0;
495  for(Int_t i = 0; i < nSM/2; i++)
496  {
497  j=2*i;
498  if((iSupMod1==j && iSupMod2==j+1) || (iSupMod1==j+1 && iSupMod2==j)) fHmggPairSameSectorSMMaskFrame[i]->Fill(invmass,fMomentum12.Pt());
499  }
500 
501  // Same side
502  for(Int_t i = 0; i < nSM-2; i++)
503  {
504  if((iSupMod1==i && iSupMod2==i+2) || (iSupMod1==i+2 && iSupMod2==i)) fHmggPairSameSideSMMaskFrame[i]->Fill(invmass,fMomentum12.Pt());
505  }
506 
507  }// Pair not facing frame
508 
509  if(invmass > fInvMassCutMin && invmass < fInvMassCutMax) //restrict to clusters really close to pi0 peak
510  {
511 
512  // Check time of cells in both clusters, and fill time histogram
513  for(Int_t icell = 0; icell < c1->GetNCells(); icell++)
514  {
515  Int_t absID = c1->GetCellAbsId(icell);
516  fHTpi0[bc%4]->Fill(absID, fEMCALCells->GetCellTime(absID)*1.e9);
517  }
518 
519  for(Int_t icell = 0; icell < c2->GetNCells(); icell++)
520  {
521  Int_t absID = c2->GetCellAbsId(icell);
522  fHTpi0[bc%4]->Fill(absID, fEMCALCells->GetCellTime(absID)*1.e9);
523  }
524 
525  //Opening angle of 2 photons
526  Float_t opangle = fMomentum1.Angle(fMomentum2.Vect())*TMath::RadToDeg();
527  //printf("*******>>>>>>>> In PEAK pt %f, angle %f \n",fMomentum12.Pt(),opangle);
528 
529 
530  fHOpeningAngle ->Fill(opangle,fMomentum12.Pt());
531  fHAsymmetry ->Fill(asym,fMomentum12.Pt());
532 
533  if(iSupMod1==iSupMod2)
534  {
535  fHOpeningAngleSM[iSupMod1] ->Fill(opangle,fMomentum12.Pt());
536  fHAsymmetrySM[iSupMod1] ->Fill(asym,fMomentum12.Pt());
537  }
538  else
539  {
540  fHOpeningAngleDifferentSM ->Fill(opangle,fMomentum12.Pt());
541  fHAsymmetryDifferentSM ->Fill(asym,fMomentum12.Pt());
542  }
543 
544  if((iSupMod1==0 && iSupMod2==2) || (iSupMod1==2 && iSupMod2==0))
545  {
546  fHOpeningAnglePairSM[0] ->Fill(opangle,fMomentum12.Pt());
547  fHAsymmetryPairSM[0] ->Fill(asym,fMomentum12.Pt());
548 
549  }
550  if((iSupMod1==1 && iSupMod2==3) || (iSupMod1==3 && iSupMod2==1))
551  {
552  fHOpeningAnglePairSM[1] ->Fill(opangle,fMomentum12.Pt());
553  fHAsymmetryPairSM[1] ->Fill(asym,fMomentum12.Pt());
554  }
555 
556  if((iSupMod1==0 && iSupMod2==1) || (iSupMod1==1 && iSupMod2==0))
557  {
558  fHOpeningAnglePairSM[2] ->Fill(opangle,fMomentum12.Pt());
559  fHAsymmetryPairSM[2] ->Fill(asym,fMomentum12.Pt());
560  }
561  if((iSupMod1==2 && iSupMod2==3) || (iSupMod1==3 && iSupMod2==2))
562  {
563  fHOpeningAnglePairSM[3] ->Fill(opangle,fMomentum12.Pt());
564  fHAsymmetryPairSM[3] ->Fill(asym,fMomentum12.Pt());
565  }
566 
567  }// pair in 100 < mass < 160
568 
569  }//in acceptance cuts
570 
571  //In case of filling only channels with second cluster in same SM
572  if(fSameSM && iSupMod1!=iSupMod2) continue;
573 
574  if (fGroupNCells == 0)
575  {
576  fHmpi0[iSupMod1][ieta1][iphi1]->Fill(invmass);
577  fHmpi0[iSupMod2][ieta2][iphi2]->Fill(invmass);
578 
579  if(invmass > fInvMassCutMin && invmass < fInvMassCutMax)//restrict to clusters really close to pi0 peak
580  {
581  fhTowerDecayPhotonHit [iSupMod1]->Fill(ieta1,iphi1);
582  fhTowerDecayPhotonEnergy [iSupMod1]->Fill(ieta1,iphi1,fMomentum1.E());
583  fhTowerDecayPhotonAsymmetry[iSupMod1]->Fill(ieta1,iphi1,asym);
584 
585  fhTowerDecayPhotonHit [iSupMod2]->Fill(ieta2,iphi2);
586  fhTowerDecayPhotonEnergy [iSupMod2]->Fill(ieta2,iphi2,fMomentum2.E());
587  fhTowerDecayPhotonAsymmetry[iSupMod2]->Fill(ieta2,iphi2,asym);
588 
589  if(!mask1)fhTowerDecayPhotonHitMaskFrame[iSupMod1]->Fill(ieta1,iphi1);
590  if(!mask2)fhTowerDecayPhotonHitMaskFrame[iSupMod2]->Fill(ieta2,iphi2);
591 
592  }// pair in mass of pi0
593  }
594  else
595  {
596  //printf("Regroup N %d, eta1 %d, phi1 %d, eta2 %d, phi2 %d \n",fGroupNCells, ieta1, iphi1, ieta2, iphi2);
597  for (Int_t i = -fGroupNCells; i < fGroupNCells+1; i++)
598  {
599  for (Int_t j = -fGroupNCells; j < fGroupNCells+1; j++)
600  {
601  Int_t absId11 = fEMCALGeo->GetAbsCellIdFromCellIndexes(iSupMod1, iphi1+j, ieta1+i);
602  Int_t absId22 = fEMCALGeo->GetAbsCellIdFromCellIndexes(iSupMod2, iphi2+j, ieta2+i);
603 
604  Bool_t ok1 = kFALSE;
605  Bool_t ok2 = kFALSE;
606 
607  for(Int_t icell = 0; icell < c1->GetNCells(); icell++)
608  {
609  if(c1->GetCellsAbsId()[icell] == absId11) ok1=kTRUE;
610  }
611 
612  for(Int_t icell = 0; icell < c2->GetNCells(); icell++)
613  {
614  if(c2->GetCellsAbsId()[icell] == absId22) ok2=kTRUE;
615  }
616 
617  if(ok1 && (ieta1+i >= 0) && (iphi1+j >= 0) && (ieta1+i < 48) && (iphi1+j < 24))
618  {
619  fHmpi0[iSupMod1][ieta1+i][iphi1+j]->Fill(invmass);
620  }
621 
622  if(ok2 && (ieta2+i >= 0) && (iphi2+j >= 0) && (ieta2+i < 48) && (iphi2+j < 24))
623  {
624  fHmpi0[iSupMod2][ieta2+i][iphi2+j]->Fill(invmass);
625  }
626  }// j loop
627  }//i loop
628  }//group cells
629 
630  AliDebug(1,Form("Mass in (SM%d,%d,%d) and (SM%d,%d,%d): %.3f GeV E1_i=%f E1_ii=%f E2_i=%f E2_ii=%f\n",
631  iSupMod1,iphi1,ieta1,iSupMod2,iphi2,ieta2,fMomentum12.M(),e1i,c1->E(),e2i,c2->E()));
632  }
633  }
634  } // end of loop over EMCAL clusters
635 }
636 
637 
642 //______________________________________________________________________
644 {
645  if ( !fRecoUtils->IsRecalibrationOn() || fCalibFilePath == "" ) return ;
646 
647  if(!fEMCALGeo)fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
648 
649  TFile * calibFactorsFile = TFile::Open(fCalibFilePath.Data());
650 
651  if ( !calibFactorsFile ) AliFatal("Cannot recover the calibration factors");
652 
653  for(Int_t ism = 0; ism < fEMCALGeo->GetNumberOfSuperModules(); ism++)
654  {
655  TH2F * histo = (TH2F*) calibFactorsFile->Get(Form("EMCALRecalFactors_SM%d",ism));
656  printf("In AliAnalysisTaskEMCALPi0CalibSelection::InitEnergyCalibrationFactors \n ---Recover calibration factor for : EMCALRecalFactors_SM%d %p\n",ism,histo);
657 
658  if ( histo )
659  fRecoUtils->SetEMCALChannelRecalibrationFactors(ism,histo);
660  else
661  AliWarning(Form("Null histogram with calibration factors for SM%d, 1 will be used for the full SM!",ism));
662  }
663 }
664 
668 //________________________________________________________________
670 {
671  Int_t runnumber = InputEvent()->GetRunNumber() ;
672 
673  //
674  // Load default geo matrices if requested
675  if(fImportGeometryFromFile && !gGeoManager)
676  {
677  if(fImportGeometryFilePath=="") // If not specified, set location depending on run number
678  {
679  // "$ALICE_ROOT/EVE/alice-data/default_geo.root"
680  if (runnumber < 140000) fImportGeometryFilePath = "$ALICE_PHYSICS/OADB/EMCAL/geometry_2010.root";
681  else if(runnumber < 171000) fImportGeometryFilePath = "$ALICE_PHYSICS/OADB/EMCAL/geometry_2011.root";
682  else if(runnumber < 198000) fImportGeometryFilePath = "$ALICE_PHYSICS/OADB/EMCAL/geometry_2012.root"; // 2012-2013
683  else fImportGeometryFilePath = "$ALICE_PHYSICS/OADB/EMCAL/geometry_2015.root"; // >= 2015
684  }
685 
686  AliInfo(Form("Import %s",fImportGeometryFilePath.Data()));
687 
688  TGeoManager::Import(fImportGeometryFilePath) ; // default need file "geometry.root" in local dir!!!!
689  }
690 
691  //
692  if(fLoadMatrices)
693  {
694  AliInfo("Load user defined EMCAL geometry matrices");
695  // OADB if available
696  AliOADBContainer emcGeoMat("AliEMCALgeo");
697 
698  if(fOADBFilePath=="") fOADBFilePath = "$ALICE_PHYSICS/OADB/EMCAL" ;
699 
700  emcGeoMat.InitFromFile(Form("%s/EMCALlocal2master.root",fOADBFilePath.Data()),"AliEMCALgeo");
701 
702  TObjArray *matEMCAL=(TObjArray*)emcGeoMat.GetObject(runnumber,"EmcalMatrices");
703 
704  for(Int_t mod = 0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
705  {
706  if (!fMatrix[mod]) // Get it from OADB
707  {
708  AliDebug(1,Form("EMCAL matrices SM %d, %p",mod,((TGeoHMatrix*) matEMCAL->At(mod))));
709  //((TGeoHMatrix*) matEMCAL->At(mod))->Print();
710 
711  fMatrix[mod] = (TGeoHMatrix*) matEMCAL->At(mod) ;
712  }
713 
714  if(fMatrix[mod])
715  {
716  if(DebugLevel() > 1)
717  fMatrix[mod]->Print();
718 
719  fEMCALGeo->SetMisalMatrix(fMatrix[mod],mod) ;
720  }
721  else if(gGeoManager)
722  {
723  AliWarning(Form("Set matrix for SM %d from gGeoManager",mod));
724  fEMCALGeo->SetMisalMatrix(fEMCALGeo->GetMatrixForSuperModuleFromGeoManager(mod),mod) ;
725  }
726  else
727  {
728  AliError(Form("Alignment matrix for SM %d is not available",mod));
729  }
730  }//SM loop
731  }//Load matrices
732  else if(!gGeoManager)
733  {
734  AliInfo("Get geo matrices from data");
735  //Still not implemented in AOD, just a workaround to be able to work at least with ESDs
736  if(!strcmp(InputEvent()->GetName(),"AliAODEvent"))
737  {
738  AliWarning("Use ideal geometry, values geometry matrix not kept in AODs");
739  }//AOD
740  else
741  {
742  AliDebug(1,"AliAnalysisTaskEMCALClusterize Load Misaligned matrices");
743 
744  for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
745  {
746  if(InputEvent()->GetEMCALMatrix(mod))
747  {
748  if(DebugLevel() > 1)
749  InputEvent()->GetEMCALMatrix(mod)->Print();
750 
751  fEMCALGeo->SetMisalMatrix(InputEvent()->GetEMCALMatrix(mod),mod) ;
752  }
753 
754  }
755  }// ESD
756  }// Load matrices from Data
757  else if(gGeoManager) // Load default matrices
758  {
759  for(Int_t mod = 0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
760  {
761  AliWarning(Form("Set matrix for SM %d from gGeoManager",mod));
762  fEMCALGeo->SetMisalMatrix(fEMCALGeo->GetMatrixForSuperModuleFromGeoManager(mod),mod) ;
763  }
764  } // gGeoManager matrices
765 
766 }
767 
771 //______________________________________________________________________
773 {
774  if(!fRecoUtils->IsRunDepRecalibrationOn()) return;
775 
776  AliOADBContainer *contRFTD=new AliOADBContainer("");
777 
778  contRFTD->InitFromFile(Form("%s/EMCALTemperatureCorrCalib.root",fOADBFilePath.Data()),"AliEMCALRunDepTempCalibCorrections");
779 
780  Int_t runnumber = InputEvent()->GetRunNumber() ;
781 
782  TH1S *htd=(TH1S*)contRFTD->GetObject(runnumber);
783 
784  //If it did not exist for this run, get closes one
785  if (!htd)
786  {
787  AliWarning(Form("No TemperatureCorrCalib Objects for run: %d",runnumber));
788 
789  // let's get the closest runnumber instead then..
790  Int_t lower = 0;
791  Int_t ic = 0;
792  Int_t maxEntry = contRFTD->GetNumberOfEntries();
793 
794  while ( (ic < maxEntry) && (contRFTD->UpperLimit(ic) < runnumber) )
795  {
796  lower = ic;
797  ic++;
798  }
799 
800  Int_t closest = lower;
801  if ( (ic<maxEntry) &&
802  (contRFTD->LowerLimit(ic)-runnumber) < (runnumber - contRFTD->UpperLimit(lower)) )
803  {
804  closest = ic;
805  }
806 
807  AliWarning(Form("TemperatureCorrCalib Objects found closest id %d from run: %d",
808  closest, contRFTD->LowerLimit(closest)));
809 
810  htd = (TH1S*) contRFTD->GetObjectByIndex(closest);
811  }
812 
813  // Fill parameters
814  if(htd)
815  {
816  AliInfo("Recalibrate (Temperature) EMCAL");
817 
818  Int_t nSM = fEMCALGeo->GetNumberOfSuperModules();
819 
820  for (Int_t ism = 0; ism < nSM; ++ism)
821  {
822  for (Int_t icol = 0; icol < 48; ++icol)
823  {
824  for (Int_t irow = 0; irow < 24; ++irow)
825  {
826  Float_t factor = fRecoUtils->GetEMCALChannelRecalibrationFactor(ism,icol,irow);
827 
828  Int_t absID = fEMCALGeo->GetAbsCellIdFromCellIndexes(ism, irow, icol); // original calibration factor
829 
830  AliDebug(3,Form(" ism %d, icol %d, irow %d,absID %d - Calib factor %1.5f - ",ism, icol, irow, absID, factor));
831 
832  factor *= htd->GetBinContent(absID) / 10000. ; // correction dependent on T
833 
834  fRecoUtils->SetEMCALChannelRecalibrationFactor(ism,icol,irow,factor);
835 
836  AliDebug(3,Form("T factor %1.5f - final factor %1.5f",
837  htd->GetBinContent(absID) / 10000.,
838  fRecoUtils->GetEMCALChannelRecalibrationFactor(ism,icol,irow)));
839  } // columns
840  } // rows
841  } // SM loop
842  }
843  else AliInfo("Do NOT recalibrate EMCAL with T variations, no params TH1");
844 
845  delete contRFTD;
846 }
847 
848 
851 //___________________________________________________________________
853 {
854  if(!fEMCALGeo)fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
855  Int_t nSM = (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules();
856 
857  fOutputContainer = new TList();
858  const Int_t buffersize = 255;
859  char hname[buffersize], htitl[buffersize];
860 
861  fhNEvents = new TH1I("hNEvents", "Number of analyzed events" , 1 , 0 , 1 ) ;
863 
864  fHmgg = new TH2F("hmgg","2-cluster invariant mass",fNbins,fMinBin,fMaxBin,100,0,10);
865  fHmgg->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
866  fHmgg->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
867  fOutputContainer->Add(fHmgg);
868 
869  fHmggDifferentSM = new TH2F("hmggDifferentSM","2-cluster invariant mass, different SM",fNbins,fMinBin,fMaxBin,100,0,10);
870  fHmggDifferentSM->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
871  fHmggDifferentSM->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
873 
874  fHOpeningAngle = new TH2F("hopang","2-cluster opening angle",100,0.,50.,100,0,10);
875  fHOpeningAngle->SetXTitle("#alpha_{#gamma #gamma}");
876  fHOpeningAngle->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
878 
879  fHOpeningAngleDifferentSM = new TH2F("hopangDifferentSM","2-cluster opening angle, different SM",100,0,50.,100,0,10);
880  fHOpeningAngleDifferentSM->SetXTitle("#alpha_{#gamma #gamma}");
881  fHOpeningAngleDifferentSM->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
883 
884  fHAsymmetry = new TH2F("hasym","2-cluster opening angle",100,0.,1.,100,0,10);
885  fHAsymmetry->SetXTitle("a");
886  fHAsymmetry->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
888 
889  fHAsymmetryDifferentSM = new TH2F("hasymDifferentSM","2-cluster opening angle, different SM",100,0,1.,100,0,10);
890  fHAsymmetryDifferentSM->SetXTitle("a");
891  fHAsymmetryDifferentSM->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
893 
894  //TString pairname[] = {"A side (0-2)", "C side (1-3)","Row 0 (0-1)", "Row 1 (2-3)"};
895 
896  fHmggMaskFrame = new TH2F("hmggMaskFrame","2-cluster invariant mass, frame masked",fNbins,fMinBin,fMaxBin,100,0,10);
897  fHmggMaskFrame->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
898  fHmggMaskFrame->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
900 
901  fHmggDifferentSMMaskFrame = new TH2F("hmggDifferentSMMaskFrame","2-cluster invariant mass, different SM, frame masked",
902  fNbins,fMinBin,fMaxBin,100,0,10);
903  fHmggDifferentSMMaskFrame->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
904  fHmggDifferentSMMaskFrame->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
906 
907  for(Int_t iSM = 0; iSM < nSM; iSM++)
908  {
909  snprintf(hname, buffersize, "hmgg_SM%d",iSM);
910  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d",iSM);
911  fHmggSM[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
912  fHmggSM[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
913  fHmggSM[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
914  fOutputContainer->Add(fHmggSM[iSM]);
915 
916  snprintf(hname, buffersize, "hmgg_SM%d_MaskFrame",iSM);
917  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d",iSM);
918  fHmggSMMaskFrame[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
919  fHmggSMMaskFrame[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
920  fHmggSMMaskFrame[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
922 
923  snprintf(hname, buffersize, "hmgg_SM%d_Zone1",iSM);
924  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d in zone 1",iSM);
925  fHmggSM_Zone1[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
926  fHmggSM_Zone1[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
927  fHmggSM_Zone1[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
928  fOutputContainer->Add(fHmggSM_Zone1[iSM]);
929 
930  snprintf(hname, buffersize, "hmgg_SM%d_Zone2",iSM);
931  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d in zone 2",iSM);
932  fHmggSM_Zone2[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
933  fHmggSM_Zone2[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
934  fHmggSM_Zone2[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
935  fOutputContainer->Add(fHmggSM_Zone2[iSM]);
936 
937  snprintf(hname, buffersize, "hmgg_SM%d_Zone3",iSM);
938  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d in zone 3",iSM);
939  fHmggSM_Zone3[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
940  fHmggSM_Zone3[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
941  fHmggSM_Zone3[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
942  fOutputContainer->Add(fHmggSM_Zone3[iSM]);
943 
944  snprintf(hname, buffersize, "hmgg_SM%d_Zone4",iSM);
945  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d in zone 4",iSM);
946  fHmggSM_Zone4[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
947  fHmggSM_Zone4[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
948  fHmggSM_Zone4[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
949  fOutputContainer->Add(fHmggSM_Zone4[iSM]);
950 
951  snprintf(hname, buffersize, "hmgg_SM%d_Zone5",iSM);
952  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d in zone 5",iSM);
953  fHmggSM_Zone5[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
954  fHmggSM_Zone5[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
955  fHmggSM_Zone5[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
956  fOutputContainer->Add(fHmggSM_Zone5[iSM]);
957 
958  snprintf(hname, buffersize, "hmgg_SM%d_Zone6",iSM);
959  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d in zone 6",iSM);
960  fHmggSM_Zone6[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
961  fHmggSM_Zone6[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
962  fHmggSM_Zone6[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
963  fOutputContainer->Add(fHmggSM_Zone6[iSM]);
964 
965  snprintf(hname, buffersize, "hmgg_SM%d_Zone7",iSM);
966  snprintf(htitl, buffersize, "Two-gamma inv. mass for super mod %d in zone 7",iSM);
967  fHmggSM_Zone7[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
968  fHmggSM_Zone7[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
969  fHmggSM_Zone7[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
970  fOutputContainer->Add(fHmggSM_Zone7[iSM]);
971 
972  if(iSM < nSM/2)
973  {
974  snprintf(hname,buffersize, "hmgg_PairSameSectorSM%d",iSM);
975  snprintf(htitl,buffersize, "Two-gamma inv. mass for SM pair Sector: %d",iSM);
976  fHmggPairSameSectorSM[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
977  fHmggPairSameSectorSM[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
978  fHmggPairSameSectorSM[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
980 
981  snprintf(hname,buffersize, "hmgg_PairSameSectorSM%d_MaskFrame",iSM);
982  snprintf(htitl,buffersize, "Two-gamma inv. mass for SM pair Sector: %d",iSM);
983  fHmggPairSameSectorSMMaskFrame[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
984  fHmggPairSameSectorSMMaskFrame[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
985  fHmggPairSameSectorSMMaskFrame[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
987 
988  fhClusterPairDiffTimeSameSector[iSM] = new TH2F(Form("hClusterPairDiffTimeSameSector%d",iSM),
989  Form("cluster pair time difference vs E, Sector %d",iSM),
990  100,0,10, 200,-100,100);
991  fhClusterPairDiffTimeSameSector[iSM]->SetXTitle("E_{pair} (GeV)");
992  fhClusterPairDiffTimeSameSector[iSM]->SetYTitle("#Delta t (ns)");
994  }
995 
996  if(iSM < nSM-2)
997  {
998  snprintf(hname,buffersize, "hmgg_PairSameSideSM%d",iSM);
999  snprintf(htitl,buffersize, "Two-gamma inv. mass for SM pair Sector: %d",iSM);
1000  fHmggPairSameSideSM[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
1001  fHmggPairSameSideSM[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
1002  fHmggPairSameSideSM[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
1004 
1005  snprintf(hname,buffersize, "hmgg_PairSameSideSM%d_MaskFrame",iSM);
1006  snprintf(htitl,buffersize, "Two-gamma inv. mass for SM pair Sector: %d",iSM);
1007  fHmggPairSameSideSMMaskFrame[iSM] = new TH2F(hname,htitl,fNbins,fMinBin,fMaxBin,100,0,10);
1008  fHmggPairSameSideSMMaskFrame[iSM]->SetXTitle("m_{#gamma #gamma} (MeV/c^{2})");
1009  fHmggPairSameSideSMMaskFrame[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
1011 
1012  fhClusterPairDiffTimeSameSide[iSM] = new TH2F(Form("hClusterPairDiffTimeSameSide%d",iSM),
1013  Form("cluster pair time difference vs E, Side %d",iSM),
1014  100,0,10, 200,-100,100);
1015  fhClusterPairDiffTimeSameSide[iSM]->SetXTitle("E_{pair} (GeV)");
1016  fhClusterPairDiffTimeSameSide[iSM]->SetYTitle("#Delta t (ns)");
1018  }
1019 
1020  snprintf(hname, buffersize, "hopang_SM%d",iSM);
1021  snprintf(htitl, buffersize, "Opening angle for super mod %d",iSM);
1022  fHOpeningAngleSM[iSM] = new TH2F(hname,htitl,100,0.,50.,100,0,10);
1023  fHOpeningAngleSM[iSM]->SetXTitle("#alpha_{#gamma #gamma} (deg)");
1024  fHOpeningAngleSM[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
1026 
1027  snprintf(hname,buffersize, "hopang_PairSM%d",iSM);
1028  snprintf(htitl,buffersize, "Opening angle for SM pair: %d",iSM);
1029  fHOpeningAnglePairSM[iSM] = new TH2F(hname,htitl,100,0.,50.,100,0,10);
1030  fHOpeningAnglePairSM[iSM]->SetXTitle("#alpha_{#gamma #gamma} (deg)");
1031  fHOpeningAnglePairSM[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
1033 
1034  snprintf(hname, buffersize, "hasym_SM%d",iSM);
1035  snprintf(htitl, buffersize, "Asymmetry for super mod %d",iSM);
1036  fHAsymmetrySM[iSM] = new TH2F(hname,htitl,100,0.,1.,100,0,10);
1037  fHAsymmetrySM[iSM]->SetXTitle("a");
1038  fHAsymmetrySM[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
1039  fOutputContainer->Add(fHAsymmetrySM[iSM]);
1040 
1041  snprintf(hname,buffersize, "hasym_PairSM%d",iSM);
1042  snprintf(htitl,buffersize, "Asymmetry for SM pair: %d",iSM);
1043  fHAsymmetryPairSM[iSM] = new TH2F(hname,htitl,100,0.,1.,100,0,10);
1044  fHAsymmetryPairSM[iSM]->SetXTitle("a");
1045  fHAsymmetryPairSM[iSM]->SetYTitle("p_{T #gamma #gamma} (GeV/c)");
1046  fOutputContainer->Add(fHAsymmetryPairSM[iSM]);
1047 
1048  Int_t colmax = 48;
1049  Int_t rowmax = 24;
1050 
1051  fhTowerDecayPhotonHit[iSM] = new TH2F (Form("hTowerDecPhotonHit_Mod%d",iSM),
1052  Form("Entries in grid of cells in Module %d",iSM),
1053  colmax+2,-1.5,colmax+0.5, rowmax+2,-1.5,rowmax+0.5);
1054  fhTowerDecayPhotonHit[iSM]->SetYTitle("row (phi direction)");
1055  fhTowerDecayPhotonHit[iSM]->SetXTitle("column (eta direction)");
1057 
1058  fhTowerDecayPhotonEnergy[iSM] = new TH2F (Form("hTowerDecPhotonEnergy_Mod%d",iSM),
1059  Form("Accumulated energy in grid of cells in Module %d",iSM),
1060  colmax+2,-1.5,colmax+0.5, rowmax+2,-1.5,rowmax+0.5);
1061  fhTowerDecayPhotonEnergy[iSM]->SetYTitle("row (phi direction)");
1062  fhTowerDecayPhotonEnergy[iSM]->SetXTitle("column (eta direction)");
1064 
1065  fhTowerDecayPhotonAsymmetry[iSM] = new TH2F (Form("hTowerDecPhotonAsymmetry_Mod%d",iSM),
1066  Form("Accumulated asymmetry in grid of cells in Module %d",iSM),
1067  colmax+2,-1.5,colmax+0.5, rowmax+2,-1.5,rowmax+0.5);
1068  fhTowerDecayPhotonAsymmetry[iSM]->SetYTitle("row (phi direction)");
1069  fhTowerDecayPhotonAsymmetry[iSM]->SetXTitle("column (eta direction)");
1071 
1072  fhTowerDecayPhotonHitMaskFrame[iSM] = new TH2F (Form("hTowerDecPhotonHit_Mod%d_MaskFrame",iSM),Form("Entries in grid of cells in Module %d",iSM),
1073  colmax+2,-1.5,colmax+0.5, rowmax+2,-1.5,rowmax+0.5);
1074  fhTowerDecayPhotonHitMaskFrame[iSM]->SetYTitle("row (phi direction)");
1075  fhTowerDecayPhotonHitMaskFrame[iSM]->SetXTitle("column (eta direction)");
1077 
1078  fhClusterTimeSM[iSM] = new TH2F(Form("hClusterTime_SM%d",iSM),"cluster time vs E",100,0,10, 100,0,1000);
1079  fhClusterTimeSM[iSM]->SetXTitle("E (GeV)");
1080  fhClusterTimeSM[iSM]->SetYTitle("t (ns)");
1081  fOutputContainer->Add(fhClusterTimeSM[iSM]);
1082 
1083  fhClusterPairDiffTimeSameSM[iSM] = new TH2F(Form("hClusterPairDiffTimeSameSM%d",iSM),
1084  Form("cluster pair time difference vs E, SM %d",iSM),
1085  100,0,10, 200,-100,100);
1086  fhClusterPairDiffTimeSameSM[iSM]->SetXTitle("E (GeV)");
1087  fhClusterPairDiffTimeSameSM[iSM]->SetYTitle("#Delta t (ns)");
1089  }
1090 
1091  Int_t nchannels = nSM*AliEMCALGeoParams::fgkEMCALRows*AliEMCALGeoParams::fgkEMCALCols;
1092  for(Int_t ibc = 0; ibc < 4; ibc++)
1093  {
1094  fHTpi0[ibc] = new TH2F(Form("hTime_BC%d",ibc),Form("Time of cell clusters under pi0 peak, bunch crossing %d",ibc),
1095  nchannels,0,nchannels, fNTimeBins,fMinTimeBin,fMaxTimeBin);
1096  fOutputContainer->Add(fHTpi0[ibc]);
1097  fHTpi0[ibc]->SetYTitle("time (ns)");
1098  fHTpi0[ibc]->SetXTitle("abs. Id. ");
1099  }
1100 
1101  fhClusterTime = new TH2F("hClusterTime","cluster time vs E",100,0,10, 100,0,1000);
1102  fhClusterTime->SetXTitle("E (GeV)");
1103  fhClusterTime->SetYTitle("t (ns)");
1105 
1106  fhClusterPairDiffTime = new TH2F("hClusterPairDiffTime","cluster pair time difference vs E",100,0,10, 800,-400,400);
1107  fhClusterPairDiffTime->SetXTitle("E_{pair} (GeV)");
1108  fhClusterPairDiffTime->SetYTitle("#Delta t (ns)");
1110 
1111  for(Int_t iMod=0; iMod < nSM; iMod++)
1112  {
1113  for(Int_t iRow=0; iRow < AliEMCALGeoParams::fgkEMCALRows; iRow++)
1114  {
1115  for(Int_t iCol=0; iCol < AliEMCALGeoParams::fgkEMCALCols; iCol++)
1116  {
1117  snprintf(hname,buffersize, "%d_%d_%d",iMod,iCol,iRow);
1118  snprintf(htitl,buffersize, "Two-gamma inv. mass for super mod %d, cell(col,row)=(%d,%d)",iMod,iCol,iRow);
1119  fHmpi0[iMod][iCol][iRow] = new TH1F(hname,htitl,fNbins,fMinBin,fMaxBin);
1120  fHmpi0[iMod][iCol][iRow]->SetXTitle("mass (MeV/c^{2})");
1121  fOutputContainer->Add(fHmpi0[iMod][iCol][iRow]);
1122  }
1123  }
1124  }
1125 
1126  fOutputContainer->SetOwner(kTRUE);
1127 
1128  PostData(1,fOutputContainer);
1129 
1130 // // cuts container, set in terminate but init and post here
1131 //
1132 // fCuts = new TList();
1133 //
1134 // fCuts ->SetOwner(kTRUE);
1135 //
1136 // PostData(2, fCuts);
1137 }
1138 
1144 //______________________________________________________________________________________________________
1146 {
1147  Int_t icol = ieta;
1148  if(iSM%2) icol+=48; // Impair SM, shift index [0-47] to [48-96]
1149 
1151  {
1152  for (Int_t imask = 0; imask < fNMaskCellColumns; imask++)
1153  {
1154  if(icol==fMaskCellColumns[imask]) return kTRUE;
1155  }
1156  }
1157 
1158  return kFALSE;
1159 }
1160 
1161 
1170 //______________________________________________________________________________________________________
1171 Bool_t AliAnalysisTaskEMCALPi0CalibSelection::IsInZone1(Int_t iSupMod, Int_t ieta, Int_t iphi)
1172 {
1173  Int_t icol = ieta;
1174  Int_t irow = iphi;
1175 
1176 // printf("SM = %i\n",iSupMod);
1177 // printf("icol = %i\n",icol);
1178 // printf("irow = %i\n",irow);
1179 
1180  if(iSupMod%2) //Odd SM
1181  {
1182  if((irow >= 2 && irow <= 21) && ((icol >= 42 && icol <= 46) || (icol >= 13 && icol <= 37) || (icol >= 1 && icol <= 9)))
1183  {
1184  if((irow >= 2 && irow <= 3) || (irow >= 20 && irow <= 21))
1185  {
1186 // printf("zone 1\n");
1187  return kTRUE;
1188  }
1189  }
1190  }
1191  else //Even SM
1192  {
1193  if((irow >= 2 && irow <= 21) && ((icol >= 1 && icol <= 5) || (icol >= 10 && icol <= 34) || (icol >= 38 && icol <= 46)))
1194  {
1195  if((irow >= 2 && irow <= 3) || (irow >= 20 && irow <= 21))
1196  {
1197 // printf("zone 1\n");
1198  return kTRUE;
1199  }
1200  }
1201  }
1202 
1203  return kFALSE;
1204 }
1205 
1206 
1215 //______________________________________________________________________________________________________
1216 Bool_t AliAnalysisTaskEMCALPi0CalibSelection::IsInZone2(Int_t iSupMod, Int_t ieta, Int_t iphi)
1217 {
1218  Int_t icol = ieta;
1219  Int_t irow = iphi;
1220 
1221 // printf("SM = %i\n",iSupMod);
1222 // printf("icol = %i\n",icol);
1223 // printf("irow = %i\n",irow);
1224 
1225  if(iSupMod%2) //Odd SM
1226  {
1227  if((irow >= 2 && irow <= 21) && ((icol >= 42 && icol <= 46) || (icol >= 13 && icol <= 37) || (icol >= 1 && icol <= 9)))
1228  {
1229  if((irow >= 2 && irow <= 3) || (irow >= 20 && irow <= 21))
1230  {
1231  return kFALSE;
1232  }
1233  else
1234  {
1235 // printf("zone 2\n");
1236  return kTRUE;
1237  }
1238  }
1239  }
1240  else //Even SM
1241  {
1242  if((irow >= 2 && irow <= 21) && ((icol >= 1 && icol <= 5) || (icol >= 10 && icol <= 34) || (icol >= 38 && icol <= 46)))
1243  {
1244  if((irow >= 2 && irow <= 3) || (irow >= 20 && irow <= 21))
1245  {
1246  return kFALSE;
1247  }
1248  else
1249  {
1250 // printf("zone 2\n");
1251  return kTRUE;
1252  }
1253  }
1254  }
1255 
1256  return kFALSE;
1257 }
1258 
1259 
1268 //______________________________________________________________________________________________________
1269 Bool_t AliAnalysisTaskEMCALPi0CalibSelection::IsInZone3(Int_t iSupMod, Int_t ieta, Int_t iphi)
1270 {
1271  Int_t icol = ieta;
1272  Int_t irow = iphi;
1273 
1274 // printf("SM = %i\n",iSupMod);
1275 // printf("icol = %i\n",icol);
1276 // printf("irow = %i\n",irow);
1277 
1278  if(iSupMod%2) //Odd SM
1279  {
1280  if((irow >= 2 && irow <= 21) && ((icol >= 42 && icol <= 46) || (icol >= 13 && icol <= 37) || (icol >= 1 && icol <= 9)))
1281  {
1282  if((icol >= 1 && icol <= 3) || (icol >= 44 && icol <= 46))
1283  {
1284 // printf("zone 3\n");
1285  return kTRUE;
1286  }
1287  }
1288  }
1289  else //Even SM
1290  {
1291  if((irow >= 2 && irow <= 21) && ((icol >= 1 && icol <= 5) || (icol >= 10 && icol <= 34) || (icol >= 38 && icol <= 46)))
1292  {
1293  if((icol >= 1 && icol <= 3) || (icol >= 44 && icol <= 46))
1294  {
1295 // printf("zone 3\n");
1296  return kTRUE;
1297  }
1298  }
1299  }
1300 
1301  return kFALSE;
1302 }
1303 
1304 
1313 //______________________________________________________________________________________________________
1314 Bool_t AliAnalysisTaskEMCALPi0CalibSelection::IsInZone4(Int_t iSupMod, Int_t ieta, Int_t iphi)
1315 {
1316  Int_t icol = ieta;
1317  Int_t irow = iphi;
1318 
1319 // printf("SM = %i\n",iSupMod);
1320 // printf("icol = %i\n",icol);
1321 // printf("irow = %i\n",irow);
1322 
1323  if(iSupMod%2) //Odd SM
1324  {
1325  if((irow >= 2 && irow <= 21) && ((icol >= 42 && icol <= 46) || (icol >= 13 && icol <= 37) || (icol >= 1 && icol <= 9)))
1326  {
1327  if((icol >= 1 && icol <= 3) || (icol >= 44 && icol <= 46))
1328  {
1329  return kFALSE;
1330  }
1331  else
1332  {
1333 // printf("zone 4\n");
1334  return kTRUE;
1335  }
1336  }
1337  }
1338  else //Even SM
1339  {
1340  if((irow >= 2 && irow <= 21) && ((icol >= 1 && icol <= 5) || (icol >= 10 && icol <= 34) || (icol >= 38 && icol <= 46)))
1341  {
1342  if((icol >= 1 && icol <= 3) || (icol >= 44 && icol <= 46))
1343  {
1344  return kFALSE;
1345  }
1346  else
1347  {
1348 // printf("zone 4\n");
1349  return kTRUE;
1350  }
1351  }
1352  }
1353 
1354  return kFALSE;
1355 }
1356 
1357 
1367 //______________________________________________________________________________________________________
1368 Bool_t AliAnalysisTaskEMCALPi0CalibSelection::IsInZone5(Int_t iSupMod, Int_t ieta, Int_t iphi)
1369 {
1370  Int_t icol = ieta;
1371  Int_t irow = iphi;
1372 
1373  //Center of ellipse
1374  Float_t col0 = 47/2;
1375  Float_t row0 = 23/2;
1376 
1377  //Parameters
1378  Float_t a = 3-col0;
1379  Float_t b = 2-row0;
1380 
1381  if(((icol-col0)*(icol-col0)) / (a*a) + ((irow-row0)*(irow-row0) / (b*b)) > 1)
1382  {
1383  return kTRUE;
1384  }
1385  else
1386  {
1387  return kFALSE;
1388  }
1389 }
1390 
1391 
1401 //______________________________________________________________________________________________________
1402 Bool_t AliAnalysisTaskEMCALPi0CalibSelection::IsInZone6(Int_t iSupMod, Int_t ieta, Int_t iphi)
1403 {
1404  Int_t icol = ieta;
1405  Int_t irow = iphi;
1406 
1407  //Center of ellipse
1408  Float_t col0 = 47/2;
1409  Float_t row0 = 23/2;
1410 
1411  //Paramters
1412  Float_t aLarge = 3-col0;
1413  Float_t bLarge = 2-row0;
1414  Float_t aSmall = 16-col0;
1415  Float_t bSmall = 7-row0;
1416 
1417  if((((icol-col0)*(icol-col0)) / (aLarge*aLarge) + ((irow-row0)*(irow-row0) / (bLarge*bLarge)) < 1) && (((icol-col0)*(icol-col0)) / (aSmall*aSmall) + ((irow-row0)*(irow-row0) / (bSmall*bSmall)) > 1))
1418  {
1419  return kTRUE;
1420  }
1421  else
1422  {
1423  return kFALSE;
1424  }
1425 }
1426 
1427 
1437 //______________________________________________________________________________________________________
1438 Bool_t AliAnalysisTaskEMCALPi0CalibSelection::IsInZone7(Int_t iSupMod, Int_t ieta, Int_t iphi)
1439 {
1440  Int_t icol = ieta;
1441  Int_t irow = iphi;
1442 
1443  //Center of ellipse
1444  Float_t col0 = 47/2;
1445  Float_t row0 = 23/2;
1446 
1447  //Paramters
1448  Float_t a = 16-col0;
1449  Float_t b = 7-row0;
1450 
1451  if(((icol-col0)*(icol-col0)) / (a*a) + ((irow-row0)*(irow-row0) / (b*b)) < 1)
1452  {
1453  return kTRUE;
1454  }
1455  else
1456  {
1457  return kFALSE;
1458  }
1459 
1460 }
1461 
1462 
1468 //__________________________________________________________________________
1470 {
1471  // Event selection
1472 
1473  if(fTriggerName!="")
1474  {
1475  AliESDEvent* esdevent = dynamic_cast<AliESDEvent*> (InputEvent());
1476  AliAODEvent* aodevent = dynamic_cast<AliAODEvent*> (InputEvent());
1477 
1478  TString triggerClass = "";
1479  if (esdevent) triggerClass = esdevent->GetFiredTriggerClasses();
1480  else if(aodevent) triggerClass = aodevent->GetFiredTriggerClasses();
1481 
1482  AliDebug(1,Form("Event %d, FiredClass %s",
1483  (Int_t)Entry(),(((AliESDEvent*)InputEvent())->GetFiredTriggerClasses()).Data()));
1484 
1485  if(!triggerClass.Contains(fTriggerName))
1486  {
1487  AliDebug(1,"Reject event!");
1488  return;
1489  }
1490  else
1491  AliDebug(1,"Accept event!");
1492  }
1493 
1494  // Get the input event
1495 
1496  AliVEvent* event = 0;
1497  if(fFilteredInput) event = AODEvent();
1498  else event = InputEvent();
1499 
1500  if(!event)
1501  {
1502  AliWarning("Input event not available!");
1503  return;
1504  }
1505 
1506  AliDebug(1,Form("<<< %s: Event %d >>>",event->GetName(), (Int_t)Entry()));
1507 
1508  // Get the primary vertex
1509 
1510  event->GetPrimaryVertex()->GetXYZ(fVertex) ;
1511 
1512  AliDebug(1,Form("Vertex: (%.3f,%.3f,%.3f)",fVertex[0],fVertex[1],fVertex[2]));
1513 
1514  //Int_t runNum = aod->GetRunNumber();
1515  //if(DebugLevel() > 1) printf("Run number: %d\n",runNum);
1516 
1517  fhNEvents->Fill(0); //Count the events to be analyzed
1518 
1519  // Acccess once the geometry matrix and temperature corrections and calibration coefficients
1520  if(fhNEvents->GetEntries() == 1)
1521  {
1523 
1525 
1526 // InitEnergyCalibrationFactors();
1527  }
1528 
1529  //Get the list of clusters and cells
1530  fEMCALCells = event->GetEMCALCells();
1531 
1532  fCaloClustersArr = new TRefArray();
1533  event->GetEMCALClusters(fCaloClustersArr);
1534 
1535  AliDebug(1,Form("N CaloClusters: %d - N CaloCells %d",fCaloClustersArr->GetEntriesFast(), fEMCALCells->GetNumberOfCells()));
1536 
1537  // Apply non linearity, new calibration, T calibration to the clusters
1538  if( fCorrectClusters )
1539  CorrectClusters();
1540 
1541  FillHistograms();
1542 
1543  delete fCaloClustersArr;
1544 
1545  PostData(1,fOutputContainer);
1546 }
1547 
1550 //_____________________________________________________
1552 {
1553  printf("Cluster cuts: %2.2f < E < %2.2f GeV; number of cells > %d; Assymetry < %1.2f, pair time diff < %2.2f, %2.2f < t < %2.2f ns\n",
1555 
1556  printf("Group %d cells\n", fGroupNCells) ;
1557 
1558  printf("Cluster maximal cell away from border at least %d cells\n", fRecoUtils->GetNumberOfCellsFromEMCALBorder()) ;
1559 
1560  printf("Histograms: bins %d; energy range: %2.2f < E < %2.2f MeV\n",fNbins,fMinBin,fMaxBin) ;
1561 
1562  printf("Switchs:\n \t Remove Bad Channels? %d; Use filtered input? %d; Correct Clusters? %d, and their position? %d \n \t Mass per channel same SM clusters? %d\n",
1563  fRecoUtils->IsBadChannelsRemovalSwitchedOn(),fFilteredInput,fCorrectClusters, fRecalPosition, fSameSM) ;
1564 
1565  printf("OADB path : %s\n",fOADBFilePath .Data());
1566  printf("Calibration path : %s\n",fCalibFilePath.Data());
1567 
1568  printf("EMCAL Geometry name: < %s >, Load Matrices %d\n",fEMCALGeoName.Data(), fLoadMatrices) ;
1569 
1570  if(fLoadMatrices) { for(Int_t ism = 0; ism < AliEMCALGeoParams::fgkEMCALModules; ism++) if(fMatrix[ism]) fMatrix[ism]->Print() ; }
1571 }
1572 
1576 //_____________________________________________________________________
1578 {
1579  if(n > fNMaskCellColumns)
1580  {
1581  delete [] fMaskCellColumns ;
1582 
1583  fMaskCellColumns = new Int_t[n] ;
1584  }
1585 
1586  fNMaskCellColumns = n ;
1587 }
1588 
1593 //___________________________________________________________________________________
1595 {
1596  if(ipos < fNMaskCellColumns) fMaskCellColumns[ipos] = icol ;
1597  else AliWarning("Mask column not set, position larger than allocated set size first") ;
1598 }
1599 
1602 //______________________________________________________________
1604 {
1605  AliDebug(1,"Not implemented");
1606 // const Int_t buffersize = 255;
1607 // char onePar[buffersize] ;
1608 
1609 // snprintf(onePar,buffersize, "Custer cuts: %2.2f < E < %2.2f GeV; %2.2f < Lambda0_2 < %2.2f GeV; min number of cells %d; Assymetry cut %1.2f, time1-time2 < %2.2f; %2.2f < T < %2.2f ns; %3.1f < Mass < %3.1f",
1610 // fEmin,fEmax, fL0min, fL0max, fMinNCells, fAsyCut, fDTimeCut, fTimeMin, fTimeMax, fInvMassCutMin, fInvMassCutMax) ;
1611 // fCuts->Add(new TObjString(onePar));
1612 // snprintf(onePar,buffersize, "Group %d cells;", fGroupNCells) ;
1613 // fCuts->Add(new TObjString(onePar));
1614 // snprintf(onePar,buffersize, "Cluster maximal cell away from border at least %d cells;", fRecoUtils->GetNumberOfCellsFromEMCALBorder()) ;
1615 // fCuts->Add(new TObjString(onePar));
1616 // snprintf(onePar,buffersize, "Histograms, Mass bins %d; energy range: %2.2f < E < %2.2f GeV;",fNbins,fMinBin,fMaxBin) ;
1617 // fCuts->Add(new TObjString(onePar));
1618 // snprintf(onePar,buffersize, "Histograms, Time bins %d; energy range: %2.2f < E < %2.2f GeV;",fNTimeBins,fMinTimeBin,fMaxTimeBin) ;
1619 // fCuts->Add(new TObjString(onePar));
1620 // snprintf(onePar,buffersize, "Switchs: Remove Bad Channels? %d; Use filtered input? %d; Correct Clusters? %d and their position? %d, Mass per channel same SM clusters? %d ",
1621 // fRecoUtils->IsBadChannelsRemovalSwitchedOn(),fFilteredInput,fCorrectClusters, fRecalPosition, fSameSM) ;
1622 // fCuts->Add(new TObjString(onePar));
1623 // snprintf(onePar,buffersize, "EMCAL Geometry name: < %s >, Load Matrices? %d",fEMCALGeoName.Data(),fLoadMatrices) ;
1624 // fCuts->Add(new TObjString(onePar));
1625 //
1626 // // Post Data
1627 // PostData(2, fCuts);
1628 }
1629 
TH2F * fHmggSM_Zone4[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM in zone 4.
Float_t fDTimeCut
Maximum difference between time of cluster pairs (ns).
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
TH2F * fHmggPairSameSectorSMMaskFrame[AliEMCALGeoParams::fgkEMCALModules/2]
! Two-cluster invariant mass per Pair, mask clusters facing frames.
TGeoHMatrix * fMatrix[AliEMCALGeoParams::fgkEMCALModules]
Bool_t IsInZone3(Int_t iSupMod, Int_t ieta, Int_t iphi)
TH2F * fhClusterPairDiffTimeSameSide[AliEMCALGeoParams::fgkEMCALModules-2]
! Diference in time of clusters same side.
Bool_t IsInZone6(Int_t iSupMod, Int_t ieta, Int_t iphi)
TH2F * fHmggSM_Zone3[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM in zone 3.
Int_t fNTimeBins
N time bins of invariant mass histograms.
TH2F * fhTowerDecayPhotonHit[AliEMCALGeoParams::fgkEMCALModules]
! Cells ordered in column/row for different module, number of times a decay photon hits...
Float_t fInvMassCutMax
Maximum mass cut for clusters to fill time or other histograms.
TH2F * fHAsymmetryPairSM[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster asymmetry vs pt per Pair,with mass close to pi0.
TH1F * fHmpi0[AliEMCALGeoParams::fgkEMCALModules][AliEMCALGeoParams::fgkEMCALCols][AliEMCALGeoParams::fgkEMCALRows]
< Two-cluster invariant mass assigned to each cell.
Bool_t IsInZone4(Int_t iSupMod, Int_t ieta, Int_t iphi)
TString fImportGeometryFilePath
Path fo geometry.root file.
Bool_t fSameSM
Combine clusters in channels on same SM.
Bool_t fLoadMatrices
Matrices set from configuration, not get from geometry.root or from ESDs/AODs.
AliEMCALRecoUtils * fRecoUtils
Access to reconstruction utilities.
AliEMCALGeometry * fEMCALGeo
! EMCAL geometry pointer.
Float_t fMaxTimeBin
Maximum time bins of invariant mass histograms.
Float_t fLogWeight
Logarithmic weight used in cluster recalibration.
Bool_t IsInZone7(Int_t iSupMod, Int_t ieta, Int_t iphi)
TH2F * fHAsymmetry
! Two-cluster asymmetry vs pt of pair, with mass close to pi0.
This task provides the input for the EMCal energy calibration with pi0 invariant mass analysis per ch...
Float_t fMinTimeBin
Minimum time bins of invariant mass histograms.
Bool_t IsInZone5(Int_t iSupMod, Int_t ieta, Int_t iphi)
TH2F * fHOpeningAngleDifferentSM
! Two-cluster opening angle vs pt of pair, each cluster in different SM, with mass close to pi0...
TH2F * fhClusterPairDiffTimeSameSM[AliEMCALGeoParams::fgkEMCALModules]
! Diference in time of clusters same SM.
TH2F * fHOpeningAngle
! Two-cluster opening angle vs pt of pair, with mass close to pi0.
TH2F * fhClusterTimeSM[AliEMCALGeoParams::fgkEMCALModules]
! Timing of clusters vs energy per SM.
TH2F * fHmggSM_Zone5[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM in zone 5.
TH2F * fHOpeningAngleSM[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster opening angle vs pt per SM,with mass close to pi0.
TLorentzVector fMomentum2
Cluster kinematics, temporal.
Float_t fInvMassCutMin
Minimum mass cut for clusters to fill time or other histograms.
Bool_t IsInZone2(Int_t iSupMod, Int_t ieta, Int_t iphi)
TH2F * fHAsymmetrySM[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster asymmetry vs pt per SM,with mass close to pi0.
TString fCalibFilePath
Full path with file with energy calibration factors per channel from previous iteration.
TH2F * fHmggSM[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM.
TH2F * fhTowerDecayPhotonHitMaskFrame[AliEMCALGeoParams::fgkEMCALModules]
! Cells ordered in column/row for different module, number of times a decay photon hits...
TH2F * fHmggSM_Zone6[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM in zone 6.
TH2F * fHmgg
! Two-cluster invariant mass vs pt of pair.
TH2F * fHTpi0[4]
! Time of cell under pi0 mass, for 4 bunch crossings.
AliAnalysisTaskEMCALPi0CalibSelection()
Default constructor. Arrays initialization is done here.
Float_t fMinBin
Minimum mass bins of invariant mass histograms.
Bool_t fRecalPosition
Switch on/off cluster position calculation, in case alignment matrices are not available.
TString fTriggerName
Trigger name must contain this name.
TH2F * fHmggDifferentSMMaskFrame
! Two-cluster invariant mass vs pt of pair, each cluster in different SM,mask clusters facing frames...
TH2F * fhClusterPairDiffTime
! Diference in time of clusters.
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Int_t fNbins
N mass bins of invariant mass histograms.
TH2F * fHmggSM_Zone2[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM in zone 2.
TH2F * fhClusterTime
! Timing of clusters vs energy.
TLorentzVector fMomentum12
Cluster pair kinematics, temporal.
TH2F * fhTowerDecayPhotonAsymmetry[AliEMCALGeoParams::fgkEMCALModules]
! Cells ordered in column/row for different module, accumulated asymmetry in the tower by decay photo...
TH2F * fHmggPairSameSideSMMaskFrame[AliEMCALGeoParams::fgkEMCALModules-2]
! Two-cluster invariant mass per Pair, mask clusters facing frames.
Bool_t fCorrectClusters
Correct clusters energy, position etc.
void UserCreateOutputObjects()
Create output container, init geometry.
TH2F * fHmggSMMaskFrame[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM, mask clusters facing frames.
void Terminate(Option_t *opt)
Create cuts/param objects and publish to slot. Comment out for the moment.
TH2F * fhTowerDecayPhotonEnergy[AliEMCALGeoParams::fgkEMCALModules]
! Cells ordered in column/row for different module, accumulated energy in the tower by decay photons...
TH2F * fHmggPairSameSideSM[AliEMCALGeoParams::fgkEMCALModules-2]
! Two-cluster invariant mass per Pair.
Float_t fMaxBin
Maximum mass bins of invariant mass histograms.
TH2F * fHmggDifferentSM
! Two-cluster invariant mass vs pt of pair, each cluster in different SM.
TH2F * fhClusterPairDiffTimeSameSector[AliEMCALGeoParams::fgkEMCALModules/2]
! Diference in time of clusters same sector.
TH1I * fhNEvents
! Number of events counter histogram.
TLorentzVector fMomentum1
Cluster kinematics, temporal.
TH2F * fHmggMaskFrame
! Two-cluster invariant mass vs pt of pair, mask clusters facing frames.
Bool_t fImportGeometryFromFile
Import geometry settings in geometry.root file.
Bool_t IsInZone1(Int_t iSupMod, Int_t ieta, Int_t iphi)
TString fOADBFilePath
Default path $ALICE_PHYSICS/OADB/EMCAL, if needed change.
TH2F * fHOpeningAnglePairSM[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster opening angle vs pt per Pair,with mass close to pi0.
Bool_t fFilteredInput
Read input produced with filter.
TH2F * fHAsymmetryDifferentSM
! Two-cluster asymmetry vs pt of pair, each cluster in different SM, with mass close to pi0...
TH2F * fHmggSM_Zone1[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM in zone 1.
TH2F * fHmggSM_Zone7[AliEMCALGeoParams::fgkEMCALModules]
! Two-cluster invariant mass per SM in zone 7.
TH2F * fHmggPairSameSectorSM[AliEMCALGeoParams::fgkEMCALModules/2]
! Two-cluster invariant mass per Pair.