AliPhysics  5364b50 (5364b50)
AliAnalysisTaskEMCALClusterize.cxx
Go to the documentation of this file.
1  /**************************************************************************
2  * Copyright(c) 1998-1999, 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 // --- Root ---
17 #include <TString.h>
18 #include <TRefArray.h>
19 #include <TClonesArray.h>
20 #include <TTree.h>
21 #include <TGeoManager.h>
22 #include <TROOT.h>
23 #include <TInterpreter.h>
24 #include <TFile.h>
25 
26 // --- AliRoot Analysis Steering
27 #include "AliAnalysisTask.h"
28 #include "AliAnalysisManager.h"
29 #include "AliESDEvent.h"
30 #include "AliGeomManager.h"
31 #include "AliVCaloCells.h"
32 #include "AliAODCaloCluster.h"
33 #include "AliCDBManager.h"
34 #include "AliCDBStorage.h"
35 #include "AliCDBEntry.h"
36 #include "AliLog.h"
37 #include "AliVEventHandler.h"
38 #include "AliAODInputHandler.h"
39 #include "AliOADBContainer.h"
40 #include "AliAODMCParticle.h"
41 #include "AliCentrality.h"
42 #include "AliMultSelection.h"
43 #include "AliDataFile.h"
44 
45 // --- EMCAL
46 #include "AliEMCALAfterBurnerUF.h"
47 #include "AliEMCALGeometry.h"
48 #include "AliEMCALClusterizerNxN.h"
49 #include "AliEMCALClusterizerv1.h"
50 #include "AliEMCALClusterizerv2.h"
51 #include "AliEMCALRecPoint.h"
52 #include "AliEMCALDigit.h"
53 
55 
59 
60 //______________________________________________________________________________
61 // Constructor.
62 //______________________________________________________________________________
64 : AliAnalysisTaskSE(name)
65 , fEvent(0)
66 , fGeom(0), fGeomName("")
67 , fGeomMatrixSet(kFALSE), fLoadGeomMatrices(kFALSE)
68 , fOCDBpath(""), fAccessOCDB(kFALSE)
69 , fDigitsArr(0), fClusterArr(0)
70 , fCaloClusterArr(0), fCaloCells(0)
71 , fRecParam(0), fClusterizer(0)
72 , fUnfolder(0), fJustUnfold(kFALSE)
73 , fOutputAODBranch(0), fOutputAODBranchName("")
74 , fOutputAODCells (0), fOutputAODCellsName (""), fInputCaloCellsName ("")
75 , fOutputAODBranchSet(0)
76 , fFillAODFile(kFALSE), fFillAODHeader(0)
77 , fFillAODCaloCells(0), fRun(-1)
78 , fRecoUtils(0), fConfigName("")
79 , fOrgClusterCellId()
80 , fCellLabels(), fCellSecondLabels(), fCellTime()
81 , fCellMatchdEta(), fCellMatchdPhi()
82 , fRecalibrateWithClusterTime(0)
83 , fMaxEvent(0), fDoTrackMatching(kFALSE), fUpdateCell(0)
84 , fSelectCell(kFALSE), fSelectCellMinE(0), fSelectCellMinFrac(0)
85 , fRejectBelowThreshold(kFALSE)
86 , fRemoveLEDEvents(kTRUE),fRemoveExoticEvents(kFALSE)
87 , fImportGeometryFromFile(kTRUE), fImportGeometryFilePath("")
88 , fOADBSet(kFALSE), fAccessOADB(kTRUE), fOADBFilePath("")
89 , fConstantTimeShift(0)
90 , fCentralityClass(""), fUseAliCentrality(0), fSelectEMCALEvent(0)
91 , fEMCALEnergyCut(0.), fEMCALNcellsCut (0)
92 , fSetCellMCLabelFromCluster(0)
93 , fSetCellMCLabelFromEdepFrac(0)
94 , fRemapMCLabelForAODs(0)
95 , fInputFromFilter(0)
96 , fTCardCorrEmulation(0), fTCardCorrClusEnerConserv(0)
97 , fRandom(0), fRandomizeTCard(1)
98 , fTCardCorrMinAmp(0.01), fTCardCorrMaxInduced(100)
99 , fPrintOnce(0)
100 
101 {
102  for(Int_t i = 0; i < 22; i++)
103  {
104  fGeomMatrix[i] = 0;
105  fTCardCorrInduceEnerProb [i] = 0;
108 
109  for(Int_t j = 0; j < 4 ; j++)
110  {
111  fTCardCorrInduceEner [j][i] = 0 ;
112  fTCardCorrInduceEnerFrac [j][i] = 0 ;
113  fTCardCorrInduceEnerFracP1 [j][i] = 0 ;
114  fTCardCorrInduceEnerFracWidth[j][i] = 0 ;
115  }
116  }
117 
118  ResetArrays();
119 
120  fCentralityBin[0] = fCentralityBin[1]=-1;
121 }
122 
123 //______________________________________________________________
125 //______________________________________________________________
127 : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisTaskEMCALClusterize")
128 , fEvent(0)
129 , fGeom(0), fGeomName("")
130 , fGeomMatrixSet(kFALSE), fLoadGeomMatrices(kFALSE)
131 , fOCDBpath(""), fAccessOCDB(kFALSE)
132 , fDigitsArr(0), fClusterArr(0)
134 , fRecParam(0), fClusterizer(0)
135 , fUnfolder(0), fJustUnfold(kFALSE)
139 , fFillAODFile(kFALSE), fFillAODHeader(0)
140 , fFillAODCaloCells(0), fRun(-1)
141 , fRecoUtils(0), fConfigName("")
146 , fMaxEvent(0), fDoTrackMatching(kFALSE), fUpdateCell(0)
148 , fRejectBelowThreshold(kFALSE)
149 , fRemoveLEDEvents(kTRUE), fRemoveExoticEvents(kFALSE)
151 , fOADBSet(kFALSE), fAccessOADB(kTRUE), fOADBFilePath("")
158 , fInputFromFilter(0)
160 , fRandom(0), fRandomizeTCard(1)
162 , fPrintOnce(0)
163 {
164  for(Int_t i = 0; i < 22; i++)
165  {
166  fGeomMatrix[i] = 0;
167  fTCardCorrInduceEnerProb [i] = 0;
170 
171  for(Int_t j = 0; j < 4 ; j++)
172  {
173  fTCardCorrInduceEner [j][i] = 0 ;
174  fTCardCorrInduceEnerFrac [j][i] = 0 ;
175  fTCardCorrInduceEnerFracP1 [j][i] = 0 ;
176  fTCardCorrInduceEnerFracWidth[j][i] = 0 ;
177  }
178  }
179 
180  ResetArrays();
181 
182  fCentralityBin[0] = fCentralityBin[1]=-1;
183 }
184 
185 //_______________________________________________________________
187 //_______________________________________________________________
189 {
190  if (fDigitsArr)
191  {
192  fDigitsArr->Clear("C");
193  delete fDigitsArr;
194  }
195 
196  if (fClusterArr)
197  {
198  fClusterArr->Delete();
199  delete fClusterArr;
200  }
201 
202  if (fCaloClusterArr)
203  {
204  fCaloClusterArr->Delete();
205  delete fCaloClusterArr;
206  }
207 
208  if(fClusterizer) delete fClusterizer;
209  if(fUnfolder) delete fUnfolder;
210  if(fRecoUtils) delete fRecoUtils;
211 }
212 
213 //_______________________________________________________
221 //_______________________________________________________________________________
223 {
224  if ( absID < 0 || absID >= 24*48*fGeom->GetNumberOfSuperModules() )
225  return kFALSE;
226 
227  Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
228  if (!fGeom->GetCellIndex(absID,imod,iTower,iIphi,iIeta))
229  return kFALSE;
230 
231  fGeom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,iphi,ieta);
232 
233  // Do not include bad channels found in analysis,
234  if ( fRecoUtils->IsBadChannelsRemovalSwitchedOn() &&
235  fRecoUtils->GetEMCALChannelStatus(imod, ieta, iphi) )
236  return kFALSE;
237 
238  return kTRUE;
239 }
240 
241 //_______________________________________________________
245 //_______________________________________________________
247 {
248  if(!fSelectEMCALEvent) return kTRUE; // accept
249 
250  if(fEMCALEnergyCut <= 0) return kTRUE; // accept
251 
252  Int_t nCluster = fEvent -> GetNumberOfCaloClusters();
253  Int_t bc = fEvent -> GetBunchCrossNumber();
254 
255  for(Int_t icalo = 0; icalo < nCluster; icalo++)
256  {
257  AliVCluster *clus = (AliVCluster*) (fEvent->GetCaloCluster(icalo));
258 
259  if( ( clus->IsEMCAL() ) && ( clus->GetNCells() > fEMCALNcellsCut ) && ( clus->E() > fEMCALEnergyCut ) &&
260  fRecoUtils->IsGoodCluster(clus,fGeom,fCaloCells,bc))
261  {
262 
263  AliDebug(1, Form("Accept : E %2.2f > %2.2f, nCells %d > %d",
264  clus->E(), fEMCALEnergyCut, clus->GetNCells(), fEMCALNcellsCut));
265 
266  return kTRUE;
267  }
268 
269  }// loop
270 
271  AliDebug(1,"Reject");
272 
273  return kFALSE;
274 }
275 
276 //_______________________________________________
279 //_______________________________________________
281 {
282  // Set it only once, unless run changed
283  if ( fOADBSet ) return ;
284 
285  TString pass = GetPass();
286 
287  AliInfo(Form("Get AODB parameters from EMCAL in %s for run %d, and <%s>",fOADBFilePath.Data(),fRun,pass.Data()));
288 
289  Int_t nSM = fGeom->GetNumberOfSuperModules();
290 
291  // Bad map
292  if(fRecoUtils->IsBadChannelsRemovalSwitchedOn())
293  {
294  AliOADBContainer *contBC=new AliOADBContainer("");
295  if(fOADBFilePath!="")
296  contBC->InitFromFile(Form("%s/EMCALBadChannels.root",fOADBFilePath.Data()),"AliEMCALBadChannels");
297  else
298  contBC->InitFromFile(AliDataFile::GetFileNameOADB("EMCAL/EMCALBadChannels.root").data(),"AliEMCALBadChannels");
299 
300  TObjArray *arrayBC=(TObjArray*)contBC->GetObject(fRun);
301 
302  if(arrayBC)
303  {
304  fRecoUtils->SwitchOnDistToBadChannelRecalculation();
305  AliInfo("Remove EMCAL bad cells");
306 
307  for (Int_t i=0; i < nSM; ++i)
308  {
309  TH2I *hbm = fRecoUtils->GetEMCALChannelStatusMap(i);
310 
311  if (hbm)
312  delete hbm;
313 
314  hbm=(TH2I*)arrayBC->FindObject(Form("EMCALBadChannelMap_Mod%d",i));
315 
316  if (!hbm)
317  {
318  AliError(Form("Can not get EMCALBadChannelMap_Mod%d",i));
319  continue;
320  }
321 
322  hbm->SetDirectory(0);
323  fRecoUtils->SetEMCALChannelStatusMap(i,hbm);
324 
325  } // loop
326  } else AliInfo("Do NOT remove EMCAL bad channels"); // run array
327 
328  delete contBC;
329  } // Remove bad
330 
331  // Energy Recalibration
332  if(fRecoUtils->IsRecalibrationOn())
333  {
334  AliOADBContainer *contRF=new AliOADBContainer("");
335 
336  if(fOADBFilePath!="")
337  contRF->InitFromFile(Form("%s/EMCALRecalib.root",fOADBFilePath.Data()),"AliEMCALRecalib");
338  else
339  contRF->InitFromFile(AliDataFile::GetFileNameOADB("EMCAL/EMCALRecalib.root").data(),"AliEMCALRecalib");
340 
341  TObjArray *recal=(TObjArray*)contRF->GetObject(fRun);
342 
343  if(recal)
344  {
345  TObjArray *recalpass=(TObjArray*)recal->FindObject(pass);
346 
347  if(recalpass)
348  {
349  TObjArray *recalib=(TObjArray*)recalpass->FindObject("Recalib");
350 
351  if(recalib)
352  {
353  AliInfo("Recalibrate EMCAL");
354  for (Int_t i=0; i<nSM; ++i)
355  {
356  TH2F *h = fRecoUtils->GetEMCALChannelRecalibrationFactors(i);
357 
358  if (h)
359  delete h;
360 
361  h = (TH2F*)recalib->FindObject(Form("EMCALRecalFactors_SM%d",i));
362 
363  if (!h)
364  {
365  AliError(Form("Could not load EMCALRecalFactors_SM%d",i));
366  continue;
367  }
368 
369  h->SetDirectory(0);
370 
371  fRecoUtils->SetEMCALChannelRecalibrationFactors(i,h);
372  } // SM loop
373  } else AliInfo("Do NOT recalibrate EMCAL, no params object array"); // array ok
374  } else AliInfo("Do NOT recalibrate EMCAL, no params for pass"); // array pass ok
375  } else AliInfo("Do NOT recalibrate EMCAL, no params for run"); // run number array ok
376 
377  delete contRF;
378  } // Recalibration on
379 
380  // Energy Recalibration, apply on top of previous calibration factors
381  if(fRecoUtils->IsRunDepRecalibrationOn())
382  {
383  AliOADBContainer *contRFTD=new AliOADBContainer("");
384 
385  if(fOADBFilePath!="")
386  contRFTD->InitFromFile(Form("%s/EMCALTemperatureCorrCalib.root",fOADBFilePath.Data()),"AliEMCALRunDepTempCalibCorrections");
387  else
388  contRFTD->InitFromFile(AliDataFile::GetFileNameOADB("EMCAL/EMCALTemperatureCorrCalib.root").data(),"AliEMCALRunDepTempCalibCorrections");
389 
390  TH1S *htd=(TH1S*)contRFTD->GetObject(fRun);
391 
392  //If it did not exist for this run, get closes one
393  if (!htd)
394  {
395  AliWarning(Form("No TemperatureCorrCalib Objects for run: %d",fRun));
396  // let's get the closest fRun instead then..
397  Int_t lower = 0;
398  Int_t ic = 0;
399  Int_t maxEntry = contRFTD->GetNumberOfEntries();
400 
401  while ( (ic < maxEntry) && (contRFTD->UpperLimit(ic) < fRun) ) {
402  lower = ic;
403  ic++;
404  }
405 
406  Int_t closest = lower;
407  if ( (ic<maxEntry) &&
408  (contRFTD->LowerLimit(ic)-fRun) < (fRun - contRFTD->UpperLimit(lower)) ) {
409  closest = ic;
410  }
411 
412  AliWarning(Form("TemperatureCorrCalib Objects found closest id %d from run: %d", closest, contRFTD->LowerLimit(closest)));
413  htd = (TH1S*) contRFTD->GetObjectByIndex(closest);
414  }
415 
416  if(htd)
417  {
418  AliInfo("Recalibrate (Temperature) EMCAL");
419 
420  for (Int_t ism=0; ism<nSM; ++ism)
421  {
422  for (Int_t icol=0; icol<48; ++icol)
423  {
424  for (Int_t irow=0; irow<24; ++irow)
425  {
426  Float_t factor = fRecoUtils->GetEMCALChannelRecalibrationFactor(ism,icol,irow);
427 
428  Int_t absID = fGeom->GetAbsCellIdFromCellIndexes(ism, irow, icol); // original calibration factor
429  factor *= htd->GetBinContent(absID) / 10000. ; // correction dependent on T
430  //printf("\t ism %d, icol %d, irow %d,absID %d, corrA %2.3f, corrB %2.3f, corrAB %2.3f\n",ism, icol, irow, absID,
431  // GetEMCALChannelRecalibrationFactor(ism,icol,irow) , htd->GetBinContent(absID) / 10000., factor);
432  fRecoUtils->SetEMCALChannelRecalibrationFactor(ism,icol,irow,factor);
433  } // columns
434  } // rows
435  } // SM loop
436  } else AliInfo("Do NOT recalibrate EMCAL with T variations, no params TH1");
437 
438  delete contRFTD;
439  } // Run by Run T calibration
440 
441  // Time Recalibration
442  if(fRecoUtils->IsTimeRecalibrationOn())
443  {
444  AliOADBContainer *contTRF=new AliOADBContainer("");
445 
446  if(fOADBFilePath!="")
447  contTRF->InitFromFile(Form("%s/EMCALTimeCalib.root",fOADBFilePath.Data()),"AliEMCALTimeCalib");
448  else
449  contTRF->InitFromFile(AliDataFile::GetFileNameOADB("EMCAL/EMCALTimeCalib.root").data(),"AliEMCALTimeCalib");
450 
451  TObjArray *trecal=(TObjArray*)contTRF->GetObject(fRun);
452 
453  if(trecal)
454  {
455  // pass number should be pass1 except on Run1 and special cases
456  TString passM = pass;
457  if ( pass=="spc_calo" ) passM = "pass3";
458  if ( fRun > 209121 ) passM = "pass1"; // run2 periods
459  if ( pass == "muon_calo_pass1" && fRun > 209121 && fRun < 244284 )
460  passM = "pass0";//period LHC15a-m
461 
462  TObjArray *trecalpass=(TObjArray*)trecal->FindObject(passM);
463 
464  if(trecalpass)
465  {
466  AliInfo("Time Recalibrate EMCAL");
467  for (Int_t ibc = 0; ibc < 4; ++ibc)
468  {
469  TH1F *h = fRecoUtils->GetEMCALChannelTimeRecalibrationFactors(ibc);
470 
471  if (h)
472  delete h;
473 
474  h = (TH1F*)trecalpass->FindObject(Form("hAllTimeAvBC%d",ibc));
475 
476  if (!h)
477  {
478  AliError(Form("Could not load hAllTimeAvBC%d",ibc));
479  continue;
480  }
481 
482  h->SetDirectory(0);
483 
484  fRecoUtils->SetEMCALChannelTimeRecalibrationFactors(ibc,h);
485  } // bunch crossing loop
486  } else AliInfo("Do NOT recalibrate time EMCAL, no params for pass"); // array pass ok
487  } else AliInfo("Do NOT recalibrate time EMCAL, no params for run"); // run number array ok
488 
489  delete contTRF;
490  } // Time recalibration on
491 
492  // L1 Phase Time Recalibration
493  if( fRecoUtils->IsL1PhaseInTimeRecalibrationOn() )
494  {
495  // init default maps first
496  if (!fRecoUtils->GetEMCALL1PhaseInTimeRecalibrationArray())
497  fRecoUtils->InitEMCALL1PhaseInTimeRecalibration() ;
498 
499  AliOADBContainer *contBC = new AliOADBContainer("");
500 
501  TFile *timeFile;
502  if(fOADBFilePath!="")
503  timeFile = new TFile(Form("%s/EMCALTimeL1PhaseCalib.root",fOADBFilePath.Data()),"read");
504  else
505  timeFile = new TFile(AliDataFile::GetFileNameOADB("EMCAL/EMCALTimeL1PhaseCalib.root").data(),"read");
506 
507  if (!timeFile || timeFile->IsZombie())
508  {
509  AliFatal(Form("EMCALTimeL1PhaseCalib.root was not found in the path provided: %s",fOADBFilePath.Data()));
510  return ;
511  }
512 
513  if (timeFile) delete timeFile;
514 
515  if(fOADBFilePath!="")
516  contBC->InitFromFile(Form("%s/EMCALTimeL1PhaseCalib.root",fOADBFilePath.Data()),"AliEMCALTimeL1PhaseCalib");
517  else
518  contBC->InitFromFile(AliDataFile::GetFileNameOADB("EMCAL/EMCALTimeL1PhaseCalib.root").data(),"AliEMCALTimeL1PhaseCalib");
519 
520  TObjArray *arrayBC=(TObjArray*)contBC->GetObject(fRun);
521  if (!arrayBC)
522  {
523  AliError(Form("No external L1 phase in time calibration set for run number: %d", fRun));
524  fRecoUtils->SwitchOffL1PhaseInTimeRecalibration();
525  }
526  else
527  {
528  // Only 1 L1 phase correction possible, except special cases
529  TString pass2 = "pass1";
530 
531  if ( pass=="muon_calo_pass1" && fRun > 209121 && fRun < 244284 )
532  pass2 = "pass0"; // period LHC15a-m
533 
534  TObjArray *arrayBCpass=(TObjArray*)arrayBC->FindObject(pass2);
535  if (!arrayBCpass)
536  {
537  AliError(Form("No external L1 phase in time calibration set for: %d -%s", fRun,pass2.Data()));
538  fRecoUtils->SwitchOffL1PhaseInTimeRecalibration();
539  }
540  else AliInfo("Recalibrate L1 Phase time");
541 
542  if(arrayBCpass)
543  {
544  if ( DebugLevel()>0 ) arrayBCpass->Print();
545 
546  TH1C *h = fRecoUtils->GetEMCALL1PhaseInTimeRecalibrationForAllSM();
547  if (h) delete h;
548 
549  h = (TH1C*)arrayBCpass->FindObject(Form("h%d",fRun));
550 
551  if (!h)
552  {
553  AliFatal(Form("There is no calibration histogram h%d for this run",fRun));
554  return;
555  }
556 
557  h->SetDirectory(0);
558  fRecoUtils->SetEMCALL1PhaseInTimeRecalibrationForAllSM(h);
559  }
560  }
561 
562  delete contBC;
563  } // L1 Phase Time Recalibration
564 
565  // Parameters already set once, so do not it again, unless run changes
566  fOADBSet = kTRUE;
567 }
568 
569 //_________________________________________________
572 //_________________________________________________
574 {
575  // Set once per run
576  if ( fOADBSet ) return;
577 
578  AliDebug(1,"Begin");
579 
580  AliCDBManager *cdb = AliCDBManager::Instance();
581 
582  if (fOCDBpath.Length())
583  {
584  cdb->SetDefaultStorage(fOCDBpath.Data());
585  AliInfo(Form("Default storage %s",fOCDBpath.Data()));
586  }
587 
588  cdb->SetRun(fRun);
589 
590  //
591  // EMCAL from RAW OCDB
592  if (fOCDBpath.Contains("alien:"))
593  {
594  cdb->SetSpecificStorage("EMCAL/Calib/Data","raw://");
595  cdb->SetSpecificStorage("EMCAL/Calib/Time","raw://");
596  cdb->SetSpecificStorage("EMCAL/Calib/Pedestals","raw://");
597  }
598 
599  TString path = cdb->GetDefaultStorage()->GetBaseFolder();
600 
601  fOADBSet = kTRUE;
602 
603  return ;
604 }
605 
606 //_____________________________________________________
610 //_____________________________________________________
612 {
613  // Get the number of stored digits, to assign the index of the new ones
614  Int_t idigit = fDigitsArr->GetEntriesFast();
615 
616  Double_t fixTime = 615.*1e-9;
617 
618  for(Int_t j = 0; j < fgkNEMCalCells; j++)
619  {
620  // Newly created?
621  if ( !fTCardCorrCellsNew[j] ) continue;
622 
623  // Accept only if at least 10 MeV
624  if ( fTCardCorrCellsEner[j] < 0.01 ) continue;
625 
626  // Check if it was not masked
627  Float_t amp = 0;
628  Double_t time = 0;
629  if ( !fRecoUtils->AcceptCalibrateCell(j,0,amp,time,fCaloCells) ) continue;
630 
631 // printf("add new digit absId %d, accept? %d, digit %d, induced amp %2.2f\n",
632 // j,fRecoUtils->AcceptCalibrateCell(j,0,amp,time,fCaloCells),idigit,fTCardCorrCellsEner[j]);
633 
634  // Now add the cell to the digits list
635  //AliEMCALDigit* digit =
636  new((*fDigitsArr)[idigit]) AliEMCALDigit( -2, -2, j, fTCardCorrCellsEner[j], fixTime,AliEMCALDigit::kHG,idigit, 0, 0, 0);
637 
638  //digit->SetListOfParents(0,0x0,0x0); // not needed
639 
640  idigit++;
641  }// loop on all possible cells
642 }
643 
644 
645 //_____________________________________________________
650 //_____________________________________________________
652 {
653  fEvent = 0x0;
654 
655  AliAODInputHandler* aodIH = dynamic_cast<AliAODInputHandler*>((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
656  Int_t eventN = Entry();
657  if(aodIH) eventN = aodIH->GetReadEntry();
658 
659  if (eventN > fMaxEvent)
660  return ;
661 
662  //printf("Clusterizer --- Event %d-- \n",eventN);
663 
664  // Check if input event are embedded events
665  // If so, take output event
666  if (aodIH && aodIH->GetMergeEvents())
667  {
668  fEvent = AODEvent();
669 
670  if(!aodIH->GetMergeEMCALCells())
671  AliFatal("Events merged but not EMCAL cells, check analysis settings!");
672 
673  AliDebug(1,"Use embedded events");
674 
675  AliDebug(1,Form("\t InputEvent N Clusters %d, N Cells %d",
676  InputEvent()->GetNumberOfCaloClusters(),InputEvent()->GetEMCALCells()->GetNumberOfCells()));
677 
678  AliDebug(1,Form("\t MergedEvent N Clusters %d, N Cells %d",
679  aodIH->GetEventToMerge()->GetNumberOfCaloClusters(), aodIH->GetEventToMerge()->GetEMCALCells()->GetNumberOfCells()));
680 
681 // if(DebugLevel() > 1)
682 // {
683 // for (Int_t icl=0; icl < aodIH->GetEventToMerge()->GetNumberOfCaloClusters(); icl++)
684 // {
685 // AliAODCaloCluster *sigCluster = aodIH->GetEventToMerge()->GetCaloCluster(icl);
686 // if(sigCluster->IsEMCAL()) AliInfo(Form("\t \t Signal cluster: i %d, E %f",icl,sigCluster->E()));
687 // }
688 // }
689 
690  AliDebug(1,Form("\t OutputEvent N Clusters %d, N Cells %d",
691  AODEvent()->GetNumberOfCaloClusters(), AODEvent()->GetEMCALCells()->GetNumberOfCells()));
692  }
693  else if(fInputFromFilter)
694  {
695  //printf("Get Input From Filtered AOD\n");
696  fEvent = AODEvent();
697  }
698  else
699  {
700  fEvent = InputEvent();
703  }
704 
705  if (!fEvent)
706  {
707  AliError("Event not available");
708  return ;
709  }
710 
711  //Recover the pointer to CaloCells container
712  if ( fInputCaloCellsName.Length() == 0 )
713  fCaloCells = fEvent->GetEMCALCells();
714  else
715  {
716  fCaloCells = (AliVCaloCells*) fEvent->FindListObject(fInputCaloCellsName);
717  if ( !fCaloCells )
718  AliWarning(Form("CaloCells branch <%s> not found use STD!",fInputCaloCellsName.Data()));
719  else
720  fCaloCells = fEvent->GetEMCALCells();
721  }
722 
723  //Process events if there is a high energy cluster
724  if(!AcceptEventEMCAL()) { fEvent = 0x0 ; return ; }
725 
726  //-------------------------------------------------------------------------------------
727  // Reject events if LED was firing, use only for LHC11a data
728  // Reject event if triggered by exotic cell and remove exotic cells if not triggered
729  //-------------------------------------------------------------------------------------
730 
731  if( IsLEDEvent( fRun ) ) { fEvent = 0x0 ; return ; }
732 
733  if( IsExoticEvent() ) { fEvent = 0x0 ; return ; }
734 
735  //-------------------------------------------------------------------------------------
736  // Set the cluster array in the event (output or input)
737  //-------------------------------------------------------------------------------------
738 
739  if ( fFillAODFile )
740  {
741  //Magic line to write events to AOD filem put after event rejection
742  AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
743  }
744  else if( !fOutputAODBranchSet )
745  {
746  // Create array of clusters/cells and put it in the input event, if output AOD not selected, only once
747  InputEvent()->AddObject(fOutputAODBranch);
748  AliInfo(Form("Add AOD clusters branch <%s> to input event",fOutputAODBranchName.Data()));
749 
750  if ( fOutputAODBranchName.Length() > 0 )
751  {
752  InputEvent()->AddObject(fOutputAODCells);
753  AliInfo(Form("Add AOD cells branch <%s> to input event",fOutputAODCellsName.Data()));
754  }
755 
756  fOutputAODBranchSet = kTRUE;
757  }
758 }
759 
760 //____________________________________________________
765 //____________________________________________________
767 {
768  Int_t nClusters = fEvent->GetNumberOfCaloClusters();
769  Int_t nClustersOrg = 0;
770 
771  AliAODInputHandler* aodIH = dynamic_cast<AliAODInputHandler*>((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
772  if(aodIH && aodIH->GetEventToMerge()) //Embedding
773  nClusters = aodIH->GetEventToMerge()->GetNumberOfCaloClusters(); //Get clusters directly from embedded signal
774 
775  ResetArrays();
776 
777  // Loop on original clusters, get MC labels, cluster time (OLD AODs),
778  // or track matching residuals (if matching is not requested)
780  {
781  for (Int_t i = 0; i < nClusters; i++)
782  {
783  AliVCluster *clus = 0;
784  if(aodIH && aodIH->GetEventToMerge()) //Embedding
785  clus = aodIH->GetEventToMerge()->GetCaloCluster(i); //Get clusters directly from embedded signal
786  else
787  clus = fEvent->GetCaloCluster(i);
788 
789  if(!clus) return;
790 
791  nClustersOrg++;
792 
793  if(!clus->IsEMCAL()) continue;
794 
795  Int_t label = clus->GetLabel();
796  Int_t label2 = -1 ;
797  if (clus->GetNLabels() >=2 ) label2 = clus->GetLabelAt(1) ;
798 
799  //printf("Org cluster %d) ID = %d, E %2.2f, Time %3.0f, N Cells %d, N MC labels %d, main MC label %d, all MC labels:\n",
800  // i, clus->GetID(), clus->E(), clus->GetTOF()*1e9, clus->GetNCells(), clus->GetNLabels(), label);
801  //
802  //for(Int_t imc = 0; imc < clus->GetNLabels(); imc++)
803  // printf("%d) Label %d, E dep frac %0.2f; ",
804  // imc, clus->GetLabelAt(imc),clus->GetClusterMCEdepFraction(imc));
805  //if(clus->GetNLabels() > 0) printf("\n");
806 
807  UShort_t * index = clus->GetCellsAbsId() ;
808  for(Int_t icell=0; icell < clus->GetNCells(); icell++ )
809  {
810  fOrgClusterCellId[index[icell]] = i;
811  fCellTime[index[icell]] = clus->GetTOF();
812  fCellMatchdEta[index[icell]] = clus->GetTrackDz();
813  fCellMatchdPhi[index[icell]] = clus->GetTrackDx();
814 
816  {
817  fCellLabels[index[icell]] = label;
818  fCellSecondLabels[index[icell]] = label2;
819  }
820  }
821  }
822  }
823 
824  // Do here induced cell energy assignation by T-Card correlation emulation, ONLY MC
826 
827  //
828  // Transform CaloCells into Digits
829  //
830  Int_t idigit = 0;
831  Int_t id = -1;
832  Float_t amp = -1;
833  Double_t time = -1;
834 
835  Int_t bc = InputEvent()->GetBunchCrossNumber();
836 
837  for (Int_t icell = 0; icell < fCaloCells->GetNumberOfCells(); icell++)
838  {
839  // Get cell values, recalibrate and not include bad channels found in analysis, nor cells with too low energy, nor exotic cell
840  id = fCaloCells->GetCellNumber(icell);
841  Bool_t accept = fRecoUtils->AcceptCalibrateCell(id,bc,amp,time,fCaloCells);
842  time-=fConstantTimeShift*1e-9; // only in case of simulations done before 2015
843 
844  // Do not include cells with too low energy, nor exotic cell
845  if( amp < fRecParam->GetMinECut() ||
846  time > fRecParam->GetTimeMax() ||
847  time < fRecParam->GetTimeMin() ) accept = kFALSE;
848 
849  // In case of old AOD analysis cell time is -1 s, approximate replacing by time of the cluster the digit belongs.
851  {
852  time = fCellTime[id];
853  //printf("cell %d time org %f - ",id, time*1.e9);
854  fRecoUtils->RecalibrateCellTime(id,bc,time);
855  //printf("recal %f\n",time*1.e9);
856  }
857 
858  //Exotic?
859  if (accept && fRecoUtils->IsExoticCell(id,fCaloCells,bc))
860  accept = kFALSE;
861 
862  if( !accept )
863  {
864  AliDebug(2,Form("Remove channel absId %d, index %d of %d, amp %f, time %f",
865  id,icell, fCaloCells->GetNumberOfCells(), amp, time*1.e9));
866  continue;
867  }
868 
869  //
870  // MC
871  //
872 
873  // Old way
874  Int_t mcLabel = fCaloCells->GetMCLabel(icell);
875  Float_t eDep = amp;
876 
877  //printf("--- Selected cell %d--- amp %f\n",id,amp);
878 
879  // New way
880  TArrayI labeArr(0);
881  TArrayF eDepArr(0);
882  Int_t nLabels = 0;
883 
885  {
886  // Old way to recover/set the cell MC label
887  // Only possibility for old Run1 productions
888 
889  if ( fSetCellMCLabelFromCluster == 1 ) mcLabel = fCellLabels[id]; // Older aliroot MC productions
890 
891  else if( fSetCellMCLabelFromCluster == 0 &&
893 
894  else mcLabel = -1; // found later
895 
896  // Last parameter of the digit object should be MC deposited energy,
897  // since it is not available in aliroot before year 2016, add just the cell amplitude so that
898  // we give more weight to the MC label of the cell with highest energy in the cluster
899 
900  Float_t efrac = fCaloCells->GetEFraction(icell);
901 
902  // When checking the MC of digits, give weight to cells with embedded signal
903  if (mcLabel > 0 && efrac < 1.e-6) efrac = 1;
904 
905  eDep *= efrac ;
906  }
907  else // fSetCellMCLabelFromEdepFrac = true
908  {
909  // New way, valid only for MC productions with aliroot > v5-07-21
910  mcLabel = -1;
911 
912  // Map the digit to cell index for later to calculate the cell MC energy deposition map
913  fCellLabels[id] = idigit;
914  //printf("\t absId %d, idigit %d\n",id,idigit);
915 
916  if(fOrgClusterCellId[id] < 0) continue; // index can be negative if noisy cell that did not form cluster
917 
918  AliVCluster *clus = 0;
919  Int_t iclus = fOrgClusterCellId[id];
920 
921  if(iclus < 0)
922  {
923  AliInfo("Negative original cluster index, skip \n");
924  continue;
925  }
926 
927  if(aodIH && aodIH->GetEventToMerge()) //Embedding
928  clus = aodIH->GetEventToMerge()->GetCaloCluster(iclus); //Get clusters directly from embedded signal
929  else
930  clus = fEvent->GetCaloCluster(iclus);
931 
932  fRecoUtils->RecalculateCellLabelsRemoveAddedGenerator(id, clus, MCEvent(), amp, labeArr, eDepArr);
933  nLabels = labeArr.GetSize();
934 
935  } // cell MC label, new
936 
937  // Apply here the found induced energies
938  if ( fTCardCorrEmulation )
939  {
940 // if( TMath::Abs(fTCardCorrCellsEner[id]) > 0.001 )
941 // printf("add energy to digit %d, absId %d: amp %2.2f + %2.2f\n",idigit,id,amp,fTCardCorrCellsEner[id]);
942  amp+=fTCardCorrCellsEner[id];
943  }
944 
945  //
946  // Create the digit
947  //
948  if(amp <= 0.01) continue ; // accept if > 10 MeV
949 
950  //if(amp > 1.5)printf("*** Add digit *** amp %f, nlabels %d, label %d, found %d, edep fract tot %f, ncluslabels %d\n",amp, nLabels,mcLabel,found,edepTotFrac,nClusLabels);
951 
952  AliEMCALDigit* digit = new((*fDigitsArr)[idigit]) AliEMCALDigit( mcLabel, mcLabel, id, amp, time,AliEMCALDigit::kHG,idigit, 0, 0, eDep);
953 
954  if(nLabels > 0)
955  digit->SetListOfParents(nLabels,labeArr.GetArray(),eDepArr.GetArray());
956 
957  idigit++;
958 
959  }
960 
961  fDigitsArr->Sort();
962 
963  // Call after Sort, if not it screws up digits index order in clusterization
965 
966 
967  //-------------------------------------------------------------------------------------
968  // Do the clusterization
969  //-------------------------------------------------------------------------------------
970 
971  fClusterizer->Digits2Clusters("");
972 
973  //-------------------------------------------------------------------------------------
974  // Transform the recpoints into AliVClusters
975  //-------------------------------------------------------------------------------------
976 
978 
979  if(!fCaloClusterArr)
980  {
981  AliWarning(Form("No array with CaloClusters, input RecPoints entries %d",fClusterArr->GetEntriesFast()));
982  return;
983  }
984 
985  AliDebug(1,Form("N clusters: before recluster %d, after recluster %d, recpoints %d",
986  nClustersOrg, fCaloClusterArr->GetEntriesFast(),fClusterArr->GetEntriesFast()));
987 
988 // if(fCaloClusterArr->GetEntriesFast() != fClusterArr->GetEntriesFast())
989 // {
990 // AliInfo("\t Some RecRoints not transformed into CaloClusters (clusterizer %d, unfold %d): Input entries %d - Output entries %d - %d (not fast)\n",
991 // fRecParam->GetClusterizerFlag(),fRecParam->GetUnfold(),
992 // fClusterArr->GetEntriesFast(), fCaloClusterArr->GetEntriesFast(), fCaloClusterArr->GetEntries());
993 // }
994 }
995 
996 
997 //_____________________________________________________
999 //_____________________________________________________
1001 {
1002  Double_t cellAmplitude = 0;
1003  Double_t cellTime = 0;
1004  Short_t cellNumber = 0;
1005  Int_t cellMCLabel = 0;
1006  Double_t cellEFrac = 0;
1007  Int_t nClustersOrg = 0;
1008 
1009  // Fill the array with the EMCAL clusters, copy them
1010  for (Int_t i = 0; i < fEvent->GetNumberOfCaloClusters(); i++)
1011  {
1012  AliVCluster *clus = fEvent->GetCaloCluster(i);
1013  if(clus->IsEMCAL())
1014  {
1015  //recalibrate/remove bad channels/etc if requested
1016  if(fRecoUtils->ClusterContainsBadChannel(fGeom,clus->GetCellsAbsId(), clus->GetNCells()))
1017  {
1018  continue;
1019  }
1020 
1021  if(fRecoUtils->IsRecalibrationOn())
1022  {
1023  //Calibrate cluster
1024  fRecoUtils->RecalibrateClusterEnergy(fGeom, clus, fCaloCells);
1025 
1026  //CalibrateCells
1027  for (Int_t icell = 0; icell < fCaloCells->GetNumberOfCells(); icell++)
1028  {
1029  if (fCaloCells->GetCell(icell, cellNumber, cellAmplitude, cellTime, cellMCLabel, cellEFrac) != kTRUE)
1030  break;
1031 
1032  Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
1033  fGeom->GetCellIndex(cellNumber,imod,iTower,iIphi,iIeta);
1034  fGeom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,iphi,ieta);
1035 
1036  //Do not include bad channels found in analysis?
1037  if( fRecoUtils->IsBadChannelsRemovalSwitchedOn() &&
1038  fRecoUtils->GetEMCALChannelStatus(imod, ieta, iphi))
1039  continue;
1040 
1041  fCaloCells->SetCell(icell, cellNumber, cellAmplitude*fRecoUtils->GetEMCALChannelRecalibrationFactor(imod,ieta,iphi),cellTime);
1042 
1043  }// cells loop
1044  }// recalibrate
1045 
1046  //Cast to ESD or AOD, needed to create the cluster array
1047  AliESDCaloCluster * esdCluster = dynamic_cast<AliESDCaloCluster*> (clus);
1048  AliAODCaloCluster * aodCluster = dynamic_cast<AliAODCaloCluster*> (clus);
1049 
1050  if (esdCluster)
1051  {
1052  fCaloClusterArr->Add( new AliESDCaloCluster(*esdCluster) );
1053  }//ESD
1054  else if(aodCluster)
1055  {
1056  fCaloClusterArr->Add( new AliAODCaloCluster(*aodCluster) );
1057  }//AOD
1058  else
1059  AliWarning("Wrong CaloCluster type?");
1060 
1061  nClustersOrg++;
1062  }
1063  }
1064 
1065  //Do the unfolding
1066  fUnfolder->UnfoldClusters(fCaloClusterArr, fCaloCells);
1067 
1068  //CLEAN-UP
1069  fUnfolder->Clear();
1070 }
1071 
1072 //_______________________________________________________________
1085 //_______________________________________________________________
1087 (Bool_t bMC , Bool_t bExotic, Bool_t bNonLin,
1088  Bool_t bRecalE, Bool_t bBad , Bool_t bRecalT, Int_t debug)
1089 {
1090  if ( debug > 0 ) printf("ConfigureEMCALRecoUtils() - **** Start ***\n");
1091 
1092  // Init
1093  if(!fRecoUtils) fRecoUtils = new AliEMCALRecoUtils ;
1094 
1095  // Exotic cells removal
1096 
1097  if(bExotic)
1098  {
1099  if ( debug > 0 ) printf("ConfigureEMCALRecoUtils() - Remove exotics in EMCAL\n");
1100  fRecoUtils->SwitchOnRejectExoticCell() ;
1101  fRecoUtils->SwitchOnRejectExoticCluster();
1102 
1103  // fRecoUtils->SetExoticCellDiffTimeCut(50); // If |t cell max - t cell in cross| > 50 do not add its energy, avoid
1104  fRecoUtils->SetExoticCellFractionCut(0.97); // 1-Ecross/Ecell > 0.97 -> out
1105  fRecoUtils->SetExoticCellMinAmplitudeCut(4.); // 4 GeV
1106  }
1107 
1108  // Recalibration factors
1109 
1110  if(bRecalE && ! bMC)
1111  {
1112  if ( debug > 0 ) printf("ConfigureEMCALRecoUtils() - Switch on energy recalibration in EMCAL\n");
1113  fRecoUtils->SwitchOnRecalibration();
1114  fRecoUtils->SwitchOnRunDepCorrection();
1115  }
1116 
1117  // Remove EMCAL hot channels
1118 
1119  if(bBad)
1120  {
1121  if ( debug > 0 ) printf("ConfigureEMCALRecoUtils() - Switch on bad channels removal in EMCAL\n");
1122  fRecoUtils->SwitchOnBadChannelsRemoval();
1123  fRecoUtils->SwitchOnDistToBadChannelRecalculation();
1124  }
1125 
1126  // *** Time recalibration settings ***
1127 
1128  if(bRecalT && ! bMC)
1129  {
1130  if ( debug > 0 ) printf("ConfigureEMCALRecoUtils() - Switch on time recalibration in EMCAL\n");
1131  fRecoUtils->SwitchOnTimeRecalibration();
1132  fRecoUtils->SwitchOnL1PhaseInTimeRecalibration() ;
1133  }
1134 
1135  // Recalculate position with method
1136 
1137  fRecoUtils->SetPositionAlgorithm(AliEMCALRecoUtils::kPosTowerGlobal);
1138 
1139  // Non linearity
1140 
1141  if( bNonLin )
1142  {
1143  if(!bMC)
1144  {
1145  if ( debug > 0 ) printf("ConfigureEMCALRecoUtils() xxx SET Non linearity correction kBeamTestCorrected xxx\n");
1146  fRecoUtils->SetNonLinearityFunction(AliEMCALRecoUtils::kBeamTestCorrectedv3);
1147  }
1148  else
1149  {
1150  if ( debug > 0 ) printf("ConfigureEMCALRecoUtils() xxx SET Non linearity correction kPi0MCv3 xxx\n");
1151  fRecoUtils->SetNonLinearityFunction(AliEMCALRecoUtils::kPi0MCv3);
1152  }
1153  }
1154  else
1155  {
1156  if ( debug > 0 ) printf("ConfigureEMCALRecoUtils() xxx DON'T SET Non linearity correction xxx\n");
1157  fRecoUtils->SetNonLinearityFunction(AliEMCALRecoUtils::kNoCorrection);
1158  }
1159 
1160 }
1161 
1162 //_____________________________________________________
1164 //_____________________________________________________
1166 {
1167  AliVCaloCells &eventEMcells = *(fEvent->GetEMCALCells());
1168  Int_t nEMcell = eventEMcells.GetNumberOfCells() ;
1169 
1170  AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());
1171  aodEMcells.CreateContainer(nEMcell);
1172  aodEMcells.SetType(AliVCaloCells::kEMCALCell);
1173  Double_t calibFactor = 1.;
1174  for (Int_t iCell = 0; iCell < nEMcell; iCell++)
1175  {
1176  Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
1177  fGeom->GetCellIndex(eventEMcells.GetCellNumber(iCell),imod,iTower,iIphi,iIeta);
1178  fGeom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,iphi,ieta);
1179 
1180  if(fRecoUtils->IsRecalibrationOn())
1181  {
1182  calibFactor = fRecoUtils->GetEMCALChannelRecalibrationFactor(imod,ieta,iphi);
1183  }
1184 
1185  if(!fRecoUtils->GetEMCALChannelStatus(imod, ieta, iphi))
1186  { //Channel is not declared as bad
1187  aodEMcells.SetCell(iCell,eventEMcells.GetCellNumber(iCell),eventEMcells.GetAmplitude(iCell)*calibFactor,
1188  eventEMcells.GetTime(iCell),eventEMcells.GetMCLabel(iCell),eventEMcells.GetEFraction(iCell));
1189  }
1190  else
1191  {
1192  aodEMcells.SetCell(iCell,eventEMcells.GetCellNumber(iCell),0,-1,-1,0);
1193  }
1194  }
1195 
1196  aodEMcells.Sort();
1197 }
1198 
1199 //__________________________________________________
1201 //__________________________________________________
1203 {
1204  AliESDEvent* esdevent = dynamic_cast<AliESDEvent*> (fEvent);
1205  AliAODEvent* aodevent = dynamic_cast<AliAODEvent*> (fEvent);
1206 
1207  Double_t pos[3] ;
1208  Double_t covVtx[6];
1209  for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
1210 
1211  AliAODHeader* header = dynamic_cast<AliAODHeader*>(AODEvent()->GetHeader());
1212  if(!header) AliFatal("Not a standard AOD");
1213  header->SetRunNumber(fRun);
1214 
1215  if(esdevent)
1216  {
1217  TTree* tree = fInputHandler->GetTree();
1218  if (tree)
1219  {
1220  TFile* file = tree->GetCurrentFile();
1221  if (file) header->SetESDFileName(file->GetName());
1222  }
1223  }
1224  else if (aodevent) {
1225  AliAODHeader * aodheader = dynamic_cast<AliAODHeader*>(aodevent->GetHeader());
1226  if(!aodheader) AliFatal("Not a standard AOD");
1227  header->SetESDFileName(aodheader->GetESDFileName());
1228  }
1229 
1230  header->SetBunchCrossNumber(fEvent->GetBunchCrossNumber());
1231  header->SetOrbitNumber(fEvent->GetOrbitNumber());
1232  header->SetPeriodNumber(fEvent->GetPeriodNumber());
1233  header->SetEventType(fEvent->GetEventType());
1234 
1235  // Centrality
1236  if(fUseAliCentrality)
1237  {
1238  if(fEvent->GetCentrality())
1239  header->SetCentrality(new AliCentrality(*(fEvent->GetCentrality())));
1240  else
1241  header->SetCentrality(0);
1242  }
1243 
1244  // Trigger
1245  header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection
1246 
1247  header->SetFiredTriggerClasses(fEvent->GetFiredTriggerClasses());
1248 
1249  header->SetTriggerMask(fEvent->GetTriggerMask());
1250  header->SetTriggerCluster(fEvent->GetTriggerCluster());
1251 
1252  if (esdevent)
1253  {
1254  header->SetL0TriggerInputs(esdevent->GetHeader()->GetL0TriggerInputs());
1255  header->SetL1TriggerInputs(esdevent->GetHeader()->GetL1TriggerInputs());
1256  header->SetL2TriggerInputs(esdevent->GetHeader()->GetL2TriggerInputs());
1257  }
1258  else if (aodevent)
1259  {
1260  header->SetL0TriggerInputs(aodevent->GetHeader()->GetL0TriggerInputs());
1261  header->SetL1TriggerInputs(aodevent->GetHeader()->GetL1TriggerInputs());
1262  header->SetL2TriggerInputs(aodevent->GetHeader()->GetL2TriggerInputs());
1263  }
1264 
1265  header->SetMagneticField(fEvent->GetMagneticField());
1266  //header->SetMuonMagFieldScale(esdevent->GetCurrentDip()/6000.);
1267 
1268  header->SetZDCN1Energy(fEvent->GetZDCN1Energy());
1269  header->SetZDCP1Energy(fEvent->GetZDCP1Energy());
1270  header->SetZDCN2Energy(fEvent->GetZDCN2Energy());
1271  header->SetZDCP2Energy(fEvent->GetZDCP2Energy());
1272  header->SetZDCEMEnergy(fEvent->GetZDCEMEnergy(0),fEvent->GetZDCEMEnergy(1));
1273 
1274  Float_t diamxy[2]={(Float_t)fEvent->GetDiamondX(),(Float_t)fEvent->GetDiamondY()};
1275  Float_t diamcov[3];
1276  fEvent->GetDiamondCovXY(diamcov);
1277  header->SetDiamond(diamxy,diamcov);
1278  if (esdevent) header->SetDiamondZ(esdevent->GetDiamondZ(),esdevent->GetSigma2DiamondZ());
1279  else if (aodevent) header->SetDiamondZ(aodevent->GetDiamondZ(),aodevent->GetSigma2DiamondZ());
1280 
1281  //
1282  Int_t nVertices = 1 ;/* = prim. vtx*/;
1283  Int_t nCaloClus = fEvent->GetNumberOfCaloClusters();
1284 
1285  AODEvent()->ResetStd(0, nVertices, 0, 0, 0, nCaloClus, 0, 0);
1286 
1287  // Access to the AOD container of vertices
1288  TClonesArray &vertices = *(AODEvent()->GetVertices());
1289  Int_t jVertices=0;
1290 
1291  // Add primary vertex. The primary tracks will be defined
1292  // after the loops on the composite objects (V0, cascades, kinks)
1293  fEvent->GetPrimaryVertex()->GetXYZ(pos);
1294  Float_t chi = 0;
1295  if (esdevent)
1296  {
1297  esdevent->GetPrimaryVertex()->GetCovMatrix(covVtx);
1298  chi = esdevent->GetPrimaryVertex()->GetChi2toNDF();
1299  }
1300  else if (aodevent)
1301  {
1302  aodevent->GetPrimaryVertex()->GetCovMatrix(covVtx);
1303  chi = aodevent->GetPrimaryVertex()->GetChi2perNDF();//Different from ESD?
1304  }
1305 
1306  AliAODVertex * primary = new(vertices[jVertices++])
1307  AliAODVertex(pos, covVtx, chi, NULL, -1, AliAODVertex::kPrimary);
1308  primary->SetName(fEvent->GetPrimaryVertex()->GetName());
1309  primary->SetTitle(fEvent->GetPrimaryVertex()->GetTitle());
1310 }
1311 
1312 //___________________________________________________________
1316 //___________________________________________________________
1318 {
1319  // First recalculate track-matching for the new clusters
1320  if(fDoTrackMatching)
1321  {
1322  fRecoUtils->FindMatches(fEvent,fCaloClusterArr,fGeom,MCEvent());
1323  }
1324 
1325  // Put the new clusters in the AOD list
1326 
1327  Int_t kNumberOfCaloClusters = fCaloClusterArr->GetEntriesFast();
1328 
1329  for(Int_t i = 0; i < kNumberOfCaloClusters; i++)
1330  {
1331  AliAODCaloCluster *newCluster = (AliAODCaloCluster *) fCaloClusterArr->At(i);
1332 
1333  newCluster->SetID(i);
1334 
1335  // Correct cluster energy non linearity
1336  newCluster->SetE(fRecoUtils->CorrectClusterEnergyLinearity(newCluster));
1337 
1338  //Add matched track
1339  if(fDoTrackMatching)
1340  {
1341  Int_t trackIndex = fRecoUtils->GetMatchedTrackIndex(i);
1342  if(trackIndex >= 0)
1343  {
1344  newCluster->AddTrackMatched(fEvent->GetTrack(trackIndex));
1345  AliDebug(2,Form("Matched Track index %d to new cluster %d",trackIndex,i));
1346  }
1347 
1348  Float_t dR = 999., dZ = 999.;
1349  fRecoUtils->GetMatchedResiduals(newCluster->GetID(),dZ,dR);
1350  newCluster->SetTrackDistance(dR,dZ);
1351  }
1352  else
1353  {// Assign previously assigned matched track in reco, very very rough
1354  Int_t absId0 = newCluster->GetCellsAbsId()[0]; // Assign match of first cell in cluster
1355  newCluster->SetTrackDistance(fCellMatchdPhi[absId0],fCellMatchdEta[absId0]);
1356  }
1357 
1358  //printf("New cluster E %f, Time %e, Id = ", newCluster->E(), newCluster->GetTOF() );
1359  //for(Int_t icell=0; icell < newCluster->GetNCells(); icell++ ) printf(" %d,", newCluster->GetCellsAbsId() [icell] );
1360  //printf("\n");
1361  //
1362  //printf("New cluster %d) ID = %d, E %2.2f, Time %3.0f, N Cells %d, N MC labels %d, main MC label %d, all MC labels:\n",
1363  // i, newCluster->GetID(), newCluster->E(), newCluster->GetTOF()*1e9, newCluster->GetNCells(), newCluster->GetNLabels(), newCluster->GetLabel());
1364  //
1365  //for(Int_t imc = 0; imc < newCluster->GetNLabels(); imc++)
1366  // printf("%d) Label %d, E dep frac %0.2f; ",
1367  // imc,newCluster->GetLabelAt(imc),newCluster->GetClusterMCEdepFraction(imc));
1368  //if(newCluster->GetNLabels() > 0) printf("\n");
1369 
1370  // Calculate distance to bad channel for new cluster. Make sure you give the list of bad channels.
1371  fRecoUtils->RecalculateClusterDistanceToBadChannel(fGeom, fCaloCells, newCluster);
1372 
1373  new((*fOutputAODBranch)[i]) AliAODCaloCluster(*newCluster);
1374 
1375  AliDebug(2,Form("New cluster %d of %d, energy %f, mc label %d",
1376  newCluster->GetID(), kNumberOfCaloClusters, newCluster->E(), newCluster->GetLabel()));
1377 
1378  } // cluster loop
1379 
1380  fOutputAODBranch->Expand(kNumberOfCaloClusters); // resize TObjArray to 'remove' slots
1381 }
1382 
1383 
1384 //________________________________________________________________
1386 //________________________________________________________________
1388 {
1389  if(fUseAliCentrality)
1390  {
1391  if(GetCentrality()) return GetCentrality()->GetCentralityPercentile(fCentralityClass) ;
1392  else return -1. ;
1393  }
1394  else
1395  {
1396  if(GetMultSelCen()) return GetMultSelCen()->GetMultiplicityPercentile(fCentralityClass,kTRUE) ;
1397  else return -1. ;
1398  }
1399 
1400  return -1.;
1401 }
1402 
1403 //_______________________________________________
1405 //_______________________________________________
1407 {
1408  if (!AliAnalysisManager::GetAnalysisManager()->GetTree())
1409  {
1410  AliError("AliAnalysisTaskEMCALClusterize::GetPass() - Pointer to tree = 0, returning null");
1411  return TString("");
1412  }
1413 
1414  if (!AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile())
1415  {
1416  AliError("AliAnalysisTaskEMCALClusterize::GetPass() - Null pointer input file, returning null");
1417  return TString("");
1418  }
1419 
1420  TString pass(AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName());
1421  if (pass.Contains("ass1")) return TString("pass1");
1422  else if (pass.Contains("ass2")) return TString("pass2");
1423  else if (pass.Contains("ass3")) return TString("pass3");
1424  else if (pass.Contains("ass4")) return TString("pass4");
1425  else if (pass.Contains("ass5")) return TString("pass5");
1426  else if (pass.Contains("LHC11c") && pass.Contains("spc_calo") ) return TString("spc_calo");
1427  else if (pass.Contains("calo") || pass.Contains("high_lumi"))
1428  {
1429  AliInfo("Path contains <calo> or <high-lumi>, set as <pass1>");
1430  return TString("pass1");
1431  }
1432  else if (pass.Contains("LHC14a1a"))
1433  {
1434  AliInfo("Enable EMCal energy calibration for this MC production!!");
1435 
1436  fRecoUtils->SwitchOnRecalibration();
1437 
1438  return TString("LHC14a1a");
1439  }
1440 
1441  // No condition fullfilled, give a default value
1442  AliInfo("Pass number string not found");
1443 
1444  return TString("");
1445 }
1446 
1447 //_________________________________________
1450 //_________________________________________
1452 {
1453  if(fDebug >=0) (AliAnalysisManager::GetAnalysisManager())->AddClassDebug(this->ClassName(),fDebug);
1454 
1455  fOADBSet = kFALSE;
1456 
1457  fBranchNames = "ESD:AliESDHeader.,EMCALCells.";
1458 
1459  if(!fRecParam) fRecParam = new AliEMCALRecParam;
1460  if(!fRecoUtils) fRecoUtils = new AliEMCALRecoUtils();
1461 
1462  if(fMaxEvent <= 0) fMaxEvent = 1000000000;
1463  if(fSelectCellMinE <= 0) fSelectCellMinE = 0.005;
1464  if(fSelectCellMinFrac <= 0) fSelectCellMinFrac = 0.001;
1465  fRejectBelowThreshold = kFALSE;
1466 
1467  // Centrality
1468  if(fCentralityClass == "") fCentralityClass = "V0M";
1469 
1470  if (fOCDBpath == "") fOCDBpath = "raw://" ;
1471  if (fOutputAODBranchName == "") fOutputAODBranchName = "newEMCALClusters" ;
1472 
1473  if(gROOT->LoadMacro(fConfigName) >=0)
1474  {
1475  AliInfo(Form("Configure analysis with %s",fConfigName.Data()));
1476  AliAnalysisTaskEMCALClusterize *clus = (AliAnalysisTaskEMCALClusterize*)gInterpreter->ProcessLine("ConfigEMCALClusterize()");
1477  fGeomName = clus->fGeomName;
1479  fOCDBpath = clus->fOCDBpath;
1480  fAccessOCDB = clus->fAccessOCDB;
1481  fRecParam = clus->fRecParam;
1482  fJustUnfold = clus->fJustUnfold;
1483  fFillAODFile = clus->fFillAODFile;
1484  fRecoUtils = clus->fRecoUtils;
1485  fConfigName = clus->fConfigName;
1486  fMaxEvent = clus->fMaxEvent;
1488  fUpdateCell = clus->fUpdateCell;
1490  for(Int_t i = 0; i < 22; i++) fGeomMatrix[i] = clus->fGeomMatrix[i] ;
1492  fCentralityBin[0] = clus->fCentralityBin[0];
1493  fCentralityBin[1] = clus->fCentralityBin[1];
1495  }
1496 }
1497 
1498 //_______________________________________________________
1501 //_______________________________________________________
1503 {
1504  if (fJustUnfold)
1505  {
1506  // init the unfolding afterburner
1507  delete fUnfolder;
1508  fUnfolder = new AliEMCALAfterBurnerUF(fRecParam->GetW0(),fRecParam->GetLocMaxCut(),fRecParam->GetMinECut());
1509  return;
1510  }
1511 
1512  //First init the clusterizer
1513  delete fClusterizer;
1514  if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
1515  fClusterizer = new AliEMCALClusterizerv1 (fGeom);
1516  else if(fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
1517  fClusterizer = new AliEMCALClusterizerv2(fGeom);
1518  else if(fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN)
1519  {
1520  fClusterizer = new AliEMCALClusterizerNxN(fGeom);
1521  fClusterizer->SetNRowDiff(fRecParam->GetNRowDiff());
1522  fClusterizer->SetNColDiff(fRecParam->GetNColDiff());
1523  }
1524  else
1525  {
1526  AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
1527  }
1528 
1529  //Now set the parameters
1530  fClusterizer->SetECAClusteringThreshold( fRecParam->GetClusteringThreshold() );
1531  fClusterizer->SetECALogWeight ( fRecParam->GetW0() );
1532  fClusterizer->SetMinECut ( fRecParam->GetMinECut() );
1533  fClusterizer->SetUnfolding ( fRecParam->GetUnfold() );
1534  fClusterizer->SetECALocalMaxCut ( fRecParam->GetLocMaxCut() );
1535  fClusterizer->SetTimeCut ( fRecParam->GetTimeCut() );
1536  fClusterizer->SetTimeMin ( fRecParam->GetTimeMin() );
1537  fClusterizer->SetTimeMax ( fRecParam->GetTimeMax() );
1538  fClusterizer->SetTimeCalibration ( fRecParam->IsTimeCalibrationOn() );
1539  fClusterizer->SetInputCalibrated ( kTRUE );
1540  fClusterizer->SetJustClusters ( kTRUE );
1541 
1542  // Initialize the cluster rec points and digits arrays and get them.
1543  fClusterizer->SetOutput(0);
1544  fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
1545  fDigitsArr = fClusterizer->GetDigits();
1546 
1547  //In case of unfolding after clusterization is requested, set the corresponding parameters
1548  if(fRecParam->GetUnfold())
1549  {
1550  Int_t i=0;
1551  for (i = 0; i < 8; i++)
1552  {
1553  fClusterizer->SetSSPars(i, fRecParam->GetSSPars(i));
1554  }//end of loop over parameters
1555 
1556  for (i = 0; i < 3; i++)
1557  {
1558  fClusterizer->SetPar5 (i, fRecParam->GetPar5(i));
1559  fClusterizer->SetPar6 (i, fRecParam->GetPar6(i));
1560  }//end of loop over parameters
1561 
1562  fClusterizer->SetRejectBelowThreshold(fRejectBelowThreshold);//here we set option of unfolding: split or reject energy
1563  fClusterizer->InitClusterUnfolding();
1564 
1565  }// to unfold
1566 }
1567 
1568 //________________________________________________________________
1572 //________________________________________________________________
1574 {
1575  if(fGeomMatrixSet) return;
1576 
1577  if (!fGeom)
1578  {
1579  if(fGeomName=="")
1580  {
1581  fGeom = AliEMCALGeometry::GetInstanceFromRunNumber(fRun);
1582  AliInfo(Form("Get EMCAL geometry name <%s> for run %d",fGeom->GetName(),fRun));
1583  }
1584  else
1585  {
1586  fGeom = AliEMCALGeometry::GetInstance(fGeomName);
1587  AliInfo(Form("Set EMCAL geometry name to <%s>",fGeomName.Data()));
1588  }
1589 
1590  // Init geometry, I do not like much to do it like this ...
1591  if(fImportGeometryFromFile && !gGeoManager)
1592  {
1593  if(fImportGeometryFilePath=="") // If not specified, set location depending on run number
1594  {
1595  // "$ALICE_ROOT/EVE/alice-data/default_geo.root"
1596  if (fRun < 140000) fImportGeometryFilePath = AliDataFile::GetFileNameOADB("EMCAL/geometry_2010.root").data();
1597  else if (fRun < 171000) fImportGeometryFilePath = AliDataFile::GetFileNameOADB("EMCAL/geometry_2011.root").data();
1598  else if (fRun < 198000) fImportGeometryFilePath = AliDataFile::GetFileNameOADB("EMCAL/geometry_2012.root").data(); // 2012-2013
1599  else fImportGeometryFilePath = AliDataFile::GetFileNameOADB("EMCAL/geometry_2015.root").data(); // >=2015
1600  }
1601 
1602  AliInfo(Form("Import %s",fImportGeometryFilePath.Data()));
1603 
1604  TGeoManager::Import(fImportGeometryFilePath) ;
1605  }
1606 
1607  AliDebug(1,Form("Init for run=%d",fRun));
1608  if (!gGeoManager) AliDebug(1,"Careful!, gGeoManager not loaded, load misalign matrices");
1609  } // geometry pointer did not exist before
1610 
1611  if(fLoadGeomMatrices)
1612  {
1613  AliInfo("Load user defined EMCAL geometry matrices");
1614 
1615  // OADB if available
1616  AliOADBContainer emcGeoMat("AliEMCALgeo");
1617  if(fOADBFilePath!="")
1618  emcGeoMat.InitFromFile(Form("%s/EMCALlocal2master.root",fOADBFilePath.Data()),"AliEMCALgeo");
1619  else
1620  emcGeoMat.InitFromFile(AliDataFile::GetFileNameOADB("EMCAL/EMCALlocal2master.root").data(),"AliEMCALgeo");
1621 
1622  TObjArray *matEMCAL=(TObjArray*)emcGeoMat.GetObject(fRun,"EmcalMatrices");
1623 
1624  for(Int_t mod=0; mod < (fGeom->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
1625  {
1626 
1627  if (!fGeomMatrix[mod]) // Get it from OADB
1628  {
1629  AliDebug(2,Form("EMCAL matrices SM %d, %p",mod,((TGeoHMatrix*) matEMCAL->At(mod))));
1630  //((TGeoHMatrix*) matEMCAL->At(mod))->Print();
1631 
1632  fGeomMatrix[mod] = (TGeoHMatrix*) matEMCAL->At(mod) ;
1633  }
1634 
1635  if(fGeomMatrix[mod])
1636  {
1637  if(DebugLevel() > 1)
1638  fGeomMatrix[mod]->Print();
1639 
1640  fGeom->SetMisalMatrix(fGeomMatrix[mod],mod) ;
1641  }
1642  else if(gGeoManager)
1643  {
1644  AliWarning(Form("Set matrix for SM %d from gGeoManager",mod));
1645  fGeom->SetMisalMatrix(fGeom->GetMatrixForSuperModuleFromGeoManager(mod),mod) ;
1646  }
1647  else
1648  {
1649  AliError(Form("Alignment atrix for SM %d is not available",mod));
1650  }
1651 
1652  fGeomMatrixSet=kTRUE;
1653 
1654  }//SM loop
1655  }//Load matrices
1656  else if(!gGeoManager)
1657  {
1658  AliInfo("Get geo matrices from data");
1659  //Still not implemented in AOD, just a workaround to be able to work at least with ESDs
1660  if(!strcmp(fEvent->GetName(),"AliAODEvent"))
1661  {
1662  AliWarning("Use ideal geometry, values geometry matrix not kept in AODs");
1663  }//AOD
1664  else
1665  {
1666  AliDebug(1,"Load Misaligned matrices");
1667 
1668  AliESDEvent* esd = dynamic_cast<AliESDEvent*>(fEvent) ;
1669 
1670  if(!esd)
1671  {
1672  AliError("This event does not contain ESDs?");
1673  if(fFillAODFile) AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kFALSE);
1674  return;
1675  }
1676 
1677  for(Int_t mod=0; mod < (fGeom->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
1678  {
1679  if(DebugLevel() > 1)
1680  esd->GetEMCALMatrix(mod)->Print();
1681 
1682  if(esd->GetEMCALMatrix(mod)) fGeom->SetMisalMatrix(esd->GetEMCALMatrix(mod),mod) ;
1683 
1684  }
1685 
1686  fGeomMatrixSet=kTRUE;
1687 
1688  } // ESD
1689  } // Load matrices from Data
1690 }
1691 
1692 //____________________________________________________
1697 //____________________________________________________
1699 {
1700  if(!fRemoveExoticEvents) return kFALSE;
1701 
1702  // Loop on cells
1703 
1704  Float_t totCellE = 0;
1705  Int_t bc = InputEvent()->GetBunchCrossNumber();
1706 
1707  for(Int_t icell = 0; icell < fCaloCells->GetNumberOfCells(); icell++)
1708  {
1709  Float_t ecell = 0 ;
1710  Double_t tcell = 0 ;
1711 
1712  Int_t absID = fCaloCells->GetCellNumber(icell);
1713  Bool_t accept = fRecoUtils->AcceptCalibrateCell(absID,bc,ecell,tcell,fCaloCells);
1714  tcell-=fConstantTimeShift*1e-9;// Only for MC simulations done before 2015
1715 
1716  if(accept && !fRecoUtils->IsExoticCell(absID,fCaloCells,bc)) totCellE += ecell;
1717  }
1718 
1719  // TString triggerclasses = event->GetFiredTriggerClasses();
1720  // printf("AliAnalysisTaskEMCALClusterize - reject event %d with cluster - reject event with ncells in SM3 %d and SM4 %d\n",(Int_t)Entry(),ncellsSM3, ncellsSM4);
1721  // if(fFillAODFile) AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kFALSE);;
1722  // return;
1723  //
1724 
1725  //printf("TotE cell %f\n",totCellE);
1726  if(totCellE < 1) return kTRUE;
1727 
1728  return kFALSE;
1729 }
1730 
1731 //________________________________________________________________
1733 //________________________________________________________________
1735 {
1736  if(!fRemoveLEDEvents) return kFALSE;
1737 
1738  //check events of LHC11a period
1739  if(run < 146858 || run > 146860) return kFALSE ;
1740 
1741  // Count number of cells with energy larger than 0.1 in SM3, cut on this number
1742  Int_t ncellsSM3 = 0;
1743  for(Int_t icell = 0; icell < fCaloCells->GetNumberOfCells(); icell++)
1744  {
1745  if ( fCaloCells->GetAmplitude (icell) > 0.1 &&
1746  fCaloCells->GetCellNumber(icell)/(24*48)==3 ) ncellsSM3++;
1747  }
1748 
1749  TString triggerclasses = fEvent->GetFiredTriggerClasses();
1750 
1751  Int_t ncellcut = 21;
1752  if(triggerclasses.Contains("EMC")) ncellcut = 35;
1753 
1754  if( ncellsSM3 >= ncellcut)
1755  {
1756  AliInfo(Form("Reject event %d with ncells in SM3 %d",(Int_t)Entry(),ncellsSM3));
1757  if(fFillAODFile) AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kFALSE);;
1758  return kTRUE;
1759  }
1760 
1761  return kFALSE;
1762 }
1763 
1764 
1765 //_______________________________________________________
1768 //_______________________________________________________
1770 {
1771  Int_t id = -1;
1772  Float_t amp = -1;
1773 
1774  // Loop on all cells with signal
1775  for (Int_t icell = 0; icell < fCaloCells->GetNumberOfCells(); icell++)
1776  {
1777  id = fCaloCells->GetCellNumber(icell);
1778  amp = fCaloCells->GetAmplitude (icell); // fCaloCells->GetCellAmplitude(id);
1779 
1780  if ( amp <= fTCardCorrMinAmp ) continue ;
1781 
1782  //
1783  // First get the SM, col-row of this tower
1784  Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
1785  fGeom->GetCellIndex(id,imod,iTower,iIphi,iIeta);
1786  fGeom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi, iIeta,iphi,ieta);
1787 
1788  //
1789  // Determine randomly if we want to create a correlation for this cell,
1790  // depending the SM number of the cell
1791  Float_t rand = fRandom.Uniform(0, 1);
1792 
1793  if ( rand > fTCardCorrInduceEnerProb[imod] ) continue;
1794 
1795  AliDebug(1,Form("Reference cell absId %d, iEta %d, iPhi %d, amp %2.3f",id,ieta,iphi,amp));
1796 
1797  //
1798  // Get the absId of the cells in the cross and same T-Card
1799  Int_t absIDup = -1;
1800  Int_t absIDdo = -1;
1801  Int_t absIDlr = -1;
1802  Int_t absIDuplr = -1;
1803  Int_t absIDdolr = -1;
1804 
1805  Int_t absIDup2 = -1;
1806  Int_t absIDup2lr = -1;
1807  Int_t absIDdo2 = -1;
1808  Int_t absIDdo2lr = -1;
1809 
1810  // Only 2 columns in the T-Card, +1 for even and -1 for odd with respect reference cell
1811  Int_t colShift = 0;
1812  if ( (ieta%2) && ieta <= AliEMCALGeoParams::fgkEMCALCols-1 ) colShift = -1;
1813  if ( !(ieta%2) && ieta >= 0 ) colShift = +1;
1814 
1815  absIDlr = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta+colShift);
1816 
1817  // Check up / down cells from reference cell not out of SM and in same T-Card
1818  if ( iphi < AliEMCALGeoParams::fgkEMCALRows-1 )
1819  {
1820  absIDup = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi+1, ieta);
1821  absIDuplr = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi+1, ieta+colShift);
1822  }
1823 
1824  if ( iphi > 0 )
1825  {
1826  absIDdo = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi-1, ieta);
1827  absIDdolr = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi-1, ieta+colShift);
1828  }
1829 
1830  // Check 2 up / 2 down cells from reference cell not out of SM
1831  if ( iphi < AliEMCALGeoParams::fgkEMCALRows-2 )
1832  {
1833  absIDup2 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi+2, ieta);
1834  absIDup2lr = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi+2, ieta+colShift);
1835  }
1836 
1837  if ( iphi > 1 )
1838  {
1839  absIDdo2 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi-2, ieta);
1840  absIDdo2lr = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi-2, ieta+colShift);
1841  }
1842 
1843  // In same T-Card?
1844  if ( TMath::FloorNint(iphi/8) != TMath::FloorNint((iphi+1)/8) ) { absIDup = -1 ; absIDuplr = -1 ; }
1845  if ( TMath::FloorNint(iphi/8) != TMath::FloorNint((iphi-1)/8) ) { absIDdo = -1 ; absIDdolr = -1 ; }
1846  if ( TMath::FloorNint(iphi/8) != TMath::FloorNint((iphi+2)/8) ) { absIDup2 = -1 ; absIDup2lr = -1 ; }
1847  if ( TMath::FloorNint(iphi/8) != TMath::FloorNint((iphi-2)/8) ) { absIDdo2 = -1 ; absIDdo2lr = -1 ; }
1848 
1849  //
1850  // Check if they are not declared bad or exist
1851  Bool_t okup = AcceptCell(absIDup );
1852  Bool_t okdo = AcceptCell(absIDdo );
1853  Bool_t oklr = AcceptCell(absIDlr );
1854  Bool_t okuplr = AcceptCell(absIDuplr );
1855  Bool_t okdolr = AcceptCell(absIDdolr );
1856  Bool_t okup2 = AcceptCell(absIDup2 );
1857  Bool_t okdo2 = AcceptCell(absIDdo2 );
1858  Bool_t okup2lr= AcceptCell(absIDup2lr);
1859  Bool_t okdo2lr= AcceptCell(absIDdo2lr);
1860 
1861  AliDebug(1,Form("Same T-Card cells:\n \t up %d (%d), down %d (%d), left-right %d (%d), up-lr %d (%d), down-lr %d (%d)\n"
1862  "\t up2 %d (%d), down2 %d (%d), up2-lr %d (%d), down2-lr %d (%d)",
1863  absIDup ,okup ,absIDdo ,okdo ,absIDlr,oklr,absIDuplr ,okuplr ,absIDdolr ,okdolr ,
1864  absIDup2,okup2,absIDdo2,okdo2, absIDup2lr,okup2lr,absIDdo2lr,okdo2lr));
1865 
1866  //
1867  // Generate some energy for the nearby cells in same TCard , depending on this cell energy
1868  // Check if originally the tower had no or little energy, in which case tag it as new
1869  Float_t fracupdown = fTCardCorrInduceEnerFrac[0][imod]+amp*fTCardCorrInduceEnerFracP1[0][imod];
1870  Float_t fracupdownleri = fTCardCorrInduceEnerFrac[1][imod]+amp*fTCardCorrInduceEnerFracP1[1][imod];
1871  Float_t fracleri = fTCardCorrInduceEnerFrac[2][imod]+amp*fTCardCorrInduceEnerFracP1[2][imod];
1872  Float_t frac2nd = fTCardCorrInduceEnerFrac[3][imod]+amp*fTCardCorrInduceEnerFracP1[3][imod];
1873 
1874  AliDebug(1,Form("Fraction for SM %d (min %2.3f, max %2.3f):\n"
1875  "\t up-down : c %2.3e, p1 %2.3e, p2 %2.4e, sig %2.3e, fraction %2.3f\n"
1876  "\t up-down-lr: c %2.3e, p1 %2.3e, p2 %2.4e, sig %2.3e, fraction %2.3f\n"
1877  "\t left-right: c %2.3e, p1 %2.3e, p2 %2.4e, sig %2.3e, fraction %2.3f\n"
1878  "\t 2nd row : c %2.3e, p1 %2.3e, p2 %2.4e, sig %2.3e, fraction %2.3f",
1880  fTCardCorrInduceEner[0][imod],fTCardCorrInduceEnerFrac[0][imod],fTCardCorrInduceEnerFracP1[0][imod],fTCardCorrInduceEnerFracWidth[0][imod],fracupdown,
1881  fTCardCorrInduceEner[1][imod],fTCardCorrInduceEnerFrac[1][imod],fTCardCorrInduceEnerFracP1[1][imod],fTCardCorrInduceEnerFracWidth[1][imod],fracupdownleri,
1882  fTCardCorrInduceEner[2][imod],fTCardCorrInduceEnerFrac[2][imod],fTCardCorrInduceEnerFracP1[2][imod],fTCardCorrInduceEnerFracWidth[2][imod],fracleri,
1883  fTCardCorrInduceEner[3][imod],fTCardCorrInduceEnerFrac[3][imod],fTCardCorrInduceEnerFracP1[3][imod],fTCardCorrInduceEnerFracWidth[3][imod],frac2nd));
1884 
1885  if( fracupdown < fTCardCorrInduceEnerFracMin[imod] ) fracupdown = fTCardCorrInduceEnerFracMin[imod];
1886  if( fracupdown > fTCardCorrInduceEnerFracMax[imod] ) fracupdown = fTCardCorrInduceEnerFracMax[imod];
1887  if( fracupdownleri < fTCardCorrInduceEnerFracMin[imod] ) fracupdownleri = fTCardCorrInduceEnerFracMin[imod];
1888  if( fracupdownleri > fTCardCorrInduceEnerFracMax[imod] ) fracupdownleri = fTCardCorrInduceEnerFracMax[imod];
1889  if( fracleri < fTCardCorrInduceEnerFracMin[imod] ) fracleri = fTCardCorrInduceEnerFracMin[imod];
1890  if( fracleri > fTCardCorrInduceEnerFracMax[imod] ) fracleri = fTCardCorrInduceEnerFracMax[imod];
1891  if( frac2nd < fTCardCorrInduceEnerFracMin[imod] ) frac2nd = fTCardCorrInduceEnerFracMin[imod];
1892  if( frac2nd > fTCardCorrInduceEnerFracMax[imod] ) frac2nd = fTCardCorrInduceEnerFracMax[imod];
1893 
1894  // Randomize the induced fraction, if requested
1895  if(fRandomizeTCard)
1896  {
1897  fracupdown = fRandom.Gaus(fracupdown ,fTCardCorrInduceEnerFracWidth[0][imod]);
1898  fracupdownleri = fRandom.Gaus(fracupdownleri,fTCardCorrInduceEnerFracWidth[1][imod]);
1899  fracleri = fRandom.Gaus(fracleri ,fTCardCorrInduceEnerFracWidth[2][imod]);
1900  frac2nd = fRandom.Gaus(frac2nd ,fTCardCorrInduceEnerFracWidth[3][imod]);
1901 
1902  AliDebug(1,Form("Randomized fraction: up-down %2.3f; up-down-left-right %2.3f; left-right %2.3f; 2nd row %2.3f",
1903  fracupdown,fracupdownleri,fracleri,frac2nd));
1904  }
1905 
1906  // Calculate induced energy
1907  Float_t indEupdown = fTCardCorrInduceEner[0][imod]+amp*fracupdown;
1908  Float_t indEupdownleri = fTCardCorrInduceEner[1][imod]+amp*fracupdownleri;
1909  Float_t indEleri = fTCardCorrInduceEner[2][imod]+amp*fracleri;
1910  Float_t indE2nd = fTCardCorrInduceEner[3][imod]+amp*frac2nd;
1911 
1912  AliDebug(1,Form("Induced energy: up-down %2.3f; up-down-left-right %2.3f; left-right %2.3f; 2nd row %2.3f",
1913  indEupdown,indEupdownleri,indEleri,indE2nd));
1914 
1915  // Check if we induce too much energy, in such case use a constant value
1916  if ( fTCardCorrMaxInduced < indE2nd ) indE2nd = fTCardCorrMaxInduced;
1917  if ( fTCardCorrMaxInduced < indEupdownleri ) indEupdownleri = fTCardCorrMaxInduced;
1918  if ( fTCardCorrMaxInduced < indEupdown ) indEupdown = fTCardCorrMaxInduced;
1919  if ( fTCardCorrMaxInduced < indEleri ) indEleri = fTCardCorrMaxInduced;
1920 
1921  AliDebug(1,Form("Induced energy, saturated?: up-down %2.3f; up-down-left-right %2.3f; left-right %2.3f; 2nd row %2.3f",
1922  indEupdown,indEupdownleri,indEleri,indE2nd));
1923 
1924  //
1925  // Add the induced energy, check if cell existed
1926  if ( okup )
1927  {
1928  fTCardCorrCellsEner[absIDup] += indEupdown;
1929 
1930  if ( fCaloCells->GetCellAmplitude(absIDup) < 0.01 ) fTCardCorrCellsNew[absIDup] = kTRUE;
1931  }
1932 
1933  if ( okdo )
1934  {
1935  fTCardCorrCellsEner[absIDdo] += indEupdown;
1936 
1937  if ( fCaloCells->GetCellAmplitude(absIDdo) < 0.01 ) fTCardCorrCellsNew[absIDdo] = kTRUE;
1938  }
1939 
1940  if ( oklr )
1941  {
1942  fTCardCorrCellsEner[absIDlr] += indEleri;
1943 
1944  if ( fCaloCells->GetCellAmplitude(absIDlr) < 0.01 ) fTCardCorrCellsNew[absIDlr] = kTRUE;
1945  }
1946 
1947  if ( okuplr )
1948  {
1949  fTCardCorrCellsEner[absIDuplr] += indEupdownleri;
1950 
1951  if ( fCaloCells->GetCellAmplitude(absIDuplr ) < 0.01 ) fTCardCorrCellsNew[absIDuplr] = kTRUE;
1952  }
1953 
1954  if ( okdolr )
1955  {
1956  fTCardCorrCellsEner[absIDdolr] += indEupdownleri;
1957 
1958  if ( fCaloCells->GetCellAmplitude(absIDdolr ) < 0.01 ) fTCardCorrCellsNew[absIDdolr] = kTRUE;
1959  }
1960 
1961  if ( okup2 )
1962  {
1963  fTCardCorrCellsEner[absIDup2] += indE2nd;
1964 
1965  if ( fCaloCells->GetCellAmplitude(absIDup2) < 0.01 ) fTCardCorrCellsNew[absIDup2] = kTRUE;
1966  }
1967 
1968  if ( okup2lr )
1969  {
1970  fTCardCorrCellsEner[absIDup2lr] += indE2nd;
1971 
1972  if ( fCaloCells->GetCellAmplitude(absIDup2lr) < 0.01 ) fTCardCorrCellsNew[absIDup2lr] = kTRUE;
1973  }
1974 
1975  if ( okdo2 )
1976  {
1977  fTCardCorrCellsEner[absIDdo2] += indE2nd;
1978 
1979  if ( fCaloCells->GetCellAmplitude(absIDdo2) < 0.01 ) fTCardCorrCellsNew[absIDdo2] = kTRUE;
1980  }
1981 
1982  if ( okdo2lr )
1983  {
1984  fTCardCorrCellsEner[absIDdo2lr] += indE2nd;
1985 
1986  if ( fCaloCells->GetCellAmplitude(absIDdo2lr) < 0.01 ) fTCardCorrCellsNew[absIDdo2lr] = kTRUE;
1987  }
1988 
1989  //
1990  // Subtract the added energy to main cell, if energy conservation is requested
1992  {
1993  if ( oklr ) fTCardCorrCellsEner[id] -= indEleri;
1994  if ( okuplr ) fTCardCorrCellsEner[id] -= indEupdownleri;
1995  if ( okdolr ) fTCardCorrCellsEner[id] -= indEupdownleri;
1996  if ( okup ) fTCardCorrCellsEner[id] -= indEupdown;
1997  if ( okdo ) fTCardCorrCellsEner[id] -= indEupdown;
1998  if ( okup2 ) fTCardCorrCellsEner[id] -= indE2nd;
1999  if ( okup2lr ) fTCardCorrCellsEner[id] -= indE2nd;
2000  if ( okdo2 ) fTCardCorrCellsEner[id] -= indE2nd;
2001  if ( okdo2lr ) fTCardCorrCellsEner[id] -= indE2nd;
2002  } // conserve energy
2003 
2004  } // cell loop
2005 
2006 }
2007 
2008 //_______________________________________________________
2010 //_______________________________________________________
2012 {
2013  AliInfo(Form("Geometry: name <%s>, matrix set <%d>, load matrix <%d>, import geo <%d> from path <%s>",
2015 
2016  if ( fAccessOCDB ) AliInfo(Form("OCDB path name <%s>", fOCDBpath.Data()));
2017  if ( fAccessOADB ) AliInfo(Form("OADB path name <%s>", fOADBFilePath.Data()));
2018 
2019  if ( fInputCaloCellsName.Length() > 0 )
2020  AliInfo(Form("Input CaloCells <%s>", fInputCaloCellsName.Data()));
2021 
2022  AliInfo(Form("Just Unfold clusters <%d>, new clusters list name <%s>, new cells name <%s>",
2024 
2025  if ( fFillAODFile ) AliInfo(Form("Fill new AOD file with: header <%d>, cells <%d>",fFillAODHeader,fFillAODCaloCells));
2026 
2027  AliInfo(Form("Use cell time for cluster <%d>, Apply constant time shift <%2.2f>, Do track-matching <%d>, Update cells <%d>, Input from ESD filter <%d>",
2029 
2030  AliInfo(Form("Reject events: larger than <%d>, LED <%d>, exotics <%d>", fMaxEvent, fRemoveLEDEvents, fRemoveExoticEvents));
2031 
2032  if (fCentralityBin[0] != -1 && fCentralityBin[1] != -1 )
2033  AliInfo(Form("Centrality bin [%2.2f,%2.2f], class <%s>, use AliCentrality? <%d>",
2035 
2036  if ( fSelectEMCALEvent )
2037  AliInfo(Form("Select events with signal in EMCal: E min <%2.2f>, n cell min <%d>", fEMCALEnergyCut, fEMCALNcellsCut));
2038 
2039  AliInfo(Form("MC label from cluster <%d>, Use EdepFrac <%d>, remap AODs <%d>",
2041 }
2042 
2043 //_______________________________________________________
2045 //_______________________________________________________
2047 {
2048  if(!fTCardCorrEmulation)
2049  {
2050  AliInfo("T-Card emulation not activated");
2051  return;
2052  }
2053 
2054  AliInfo(Form("T-Card emulation activated, energy conservation <%d>, randomize E <%d>, induced energy parameters:",
2056 
2057  AliInfo(Form("T-Card emulation super-modules fraction: Min cell E %2.2f Max induced E %2.2f",
2059 
2060  for(Int_t ism = 0; ism < 22; ism++)
2061  {
2062  printf("\t sm %d, fraction %2.3f, E frac abs min %2.3e max %2.3e \n",
2064 
2065  for(Int_t icell = 0; icell < 4; icell++)
2066  {
2067  printf("\t \t cell type %d, c %2.4e, p0 %2.4e, p1 %2.4e, sigma %2.4e \n",
2068  icell,fTCardCorrInduceEner[icell][ism],fTCardCorrInduceEnerFrac[icell][ism],
2070  }
2071  }
2072 }
2073 
2074 //_______________________________________________________
2078 //_______________________________________________________
2080 {
2081  Int_t j = 0;
2082  for(Int_t i = 0; i < fClusterArr->GetEntriesFast(); i++)
2083  {
2084  AliEMCALRecPoint *recPoint = (AliEMCALRecPoint*) fClusterArr->At(i);
2085 
2086  const Int_t ncells = recPoint->GetMultiplicity();
2087  Int_t ncellsTrue = 0;
2088 
2089  if(recPoint->GetEnergy() < fRecParam->GetClusteringThreshold()) continue;
2090 
2091  // cells and their amplitude fractions
2092  UShort_t absIds[ncells];
2093  Double32_t ratios[ncells];
2094 
2095  //For later check embedding
2096  AliAODInputHandler* aodIH = dynamic_cast<AliAODInputHandler*>((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
2097 
2098  Float_t clusterE = 0;
2099  for (Int_t c = 0; c < ncells; c++)
2100  {
2101  AliEMCALDigit *digit = (AliEMCALDigit*) fDigitsArr->At(recPoint->GetDigitsList()[c]);
2102 
2103  absIds[ncellsTrue] = digit->GetId();
2104  ratios[ncellsTrue] = recPoint->GetEnergiesList()[c]/digit->GetAmplitude();
2105 
2106  if ( !fRecParam->GetUnfold() && (ratios[ncellsTrue] > 1 || ratios[ncellsTrue] < 1) )
2107  AliWarning(Form("recpoint cell E %2.3f but digit E %2.3f and no unfolding", recPoint->GetEnergiesList()[c], digit->GetAmplitude()));
2108 
2109  // In case of unfolding, remove digits with unfolded energy too low
2110  if(fSelectCell)
2111  {
2112  if (recPoint->GetEnergiesList()[c] < fSelectCellMinE || ratios[ncellsTrue] < fSelectCellMinFrac)
2113  {
2114 
2115  AliDebug(2,Form("Too small energy in cell of cluster: cluster cell %f, digit %f",
2116  recPoint->GetEnergiesList()[c],digit->GetAmplitude()));
2117  continue;
2118  } // if cuts
2119  }// Select cells
2120 
2121  //Recalculate cluster energy and number of cluster cells in case any of the cells was rejected
2122  clusterE +=recPoint->GetEnergiesList()[c];
2123 
2124  // In case of embedding, fill ratio with amount of signal,
2125  if (aodIH && aodIH->GetMergeEvents())
2126  {
2127  //AliVCaloCells* inEMCALCells = InputEvent()->GetEMCALCells();
2128  AliVCaloCells* meEMCALCells = aodIH->GetEventToMerge()->GetEMCALCells();
2129  AliVCaloCells* ouEMCALCells = AODEvent()->GetEMCALCells();
2130 
2131  Float_t sigAmplitude = meEMCALCells->GetCellAmplitude(absIds[ncellsTrue]);
2132  //Float_t bkgAmplitude = inEMCALCells->GetCellAmplitude(absIds[ncellsTrue]);
2133  Float_t sumAmplitude = ouEMCALCells->GetCellAmplitude(absIds[ncellsTrue]);
2134  //printf("\t AbsID %d, amplitude : bkg %f, sigAmplitude %f, summed %f - %f\n",absIds[ncellsTrue], bkgAmplitude, sigAmplitude, sumAmplitude, digit->GetAmplitude());
2135 
2136  if(sumAmplitude > 0) ratios[ncellsTrue] = sigAmplitude/sumAmplitude;
2137  //printf("\t \t ratio %f\n",ratios[ncellsTrue]);
2138 
2139  }//Embedding
2140 
2141  ncellsTrue++;
2142 
2143  }// cluster cell loop
2144 
2145  if (ncellsTrue < 1)
2146  {
2147  AliDebug(2,Form("Skipping cluster with no cells avobe threshold E = %f, ncells %d",
2148  recPoint->GetEnergy(), ncells));
2149  continue;
2150  }
2151 
2152  //if(ncellsTrue != ncells) printf("Old E %f, ncells %d; New E %f, ncells %d\n",recPoint->GetEnergy(),ncells,clusterE,ncellsTrue);
2153 
2154  if(clusterE < fRecParam->GetClusteringThreshold())
2155  {
2156  AliDebug(2,Form("Remove cluster with energy below seed threshold %f",clusterE));
2157  continue;
2158  }
2159 
2160  TVector3 gpos;
2161  Float_t g[3];
2162 
2163  // calculate new cluster position
2164 
2165  recPoint->EvalGlobalPosition(fRecParam->GetW0(), fDigitsArr);
2166  recPoint->GetGlobalPosition(gpos);
2167  gpos.GetXYZ(g);
2168 
2169  // create a new cluster
2170 
2171  (*fCaloClusterArr)[j] = new AliAODCaloCluster() ;
2172  AliAODCaloCluster *clus = dynamic_cast<AliAODCaloCluster *>( fCaloClusterArr->At(j) ) ;
2173  j++;
2174  clus->SetType(AliVCluster::kEMCALClusterv1);
2175  clus->SetE(clusterE);
2176  clus->SetPosition(g);
2177  clus->SetNCells(ncellsTrue);
2178  clus->SetCellsAbsId(absIds);
2179  clus->SetCellsAmplitudeFraction(ratios);
2180  clus->SetChi2(-1); //not yet implemented
2181  clus->SetTOF(recPoint->GetTime()) ; //time-of-flight
2182  clus->SetNExMax(recPoint->GetNExMax()); //number of local maxima
2183  clus->SetDistanceToBadChannel(recPoint->GetDistanceToBadTower());
2184 
2185  if(ncells == ncellsTrue)
2186  {
2187  Float_t elipAxis[2];
2188  recPoint->GetElipsAxis(elipAxis);
2189  clus->SetM02(elipAxis[0]*elipAxis[0]) ;
2190  clus->SetM20(elipAxis[1]*elipAxis[1]) ;
2191  clus->SetDispersion(recPoint->GetDispersion());
2192  }
2193  else if(fSelectCell)
2194  {
2195  // In case some cells rejected, in unfolding case, recalculate
2196  // shower shape parameters and position
2197  AliDebug(2,Form("Cells removed from cluster (ncells %d, ncellsTrue %d), recalculate Shower Shape",ncells,ncellsTrue));
2198 
2199  AliVCaloCells* cells = 0x0;
2200  if (aodIH && aodIH->GetMergeEvents()) cells = AODEvent() ->GetEMCALCells();
2201  else cells = InputEvent()->GetEMCALCells();
2202 
2203  fRecoUtils->RecalculateClusterShowerShapeParameters(fGeom,cells,clus);
2204  fRecoUtils->RecalculateClusterPID(clus);
2205  fRecoUtils->RecalculateClusterPosition(fGeom,cells,clus);
2206  }
2207 
2208  // MC
2209 
2212  else
2213  {
2214  //
2215  // Normal case, trust what the clusterizer has found
2216  //
2217  Int_t parentMult = 0;
2218  Int_t *parentList = recPoint->GetParents(parentMult);
2219  Float_t *parentListDE = recPoint->GetParentsDE(); // deposited energy
2220 
2221  clus->SetLabel(parentList, parentMult);
2223  clus->SetClusterMCEdepFractionFromEdepArray(parentListDE);
2224 
2225  //
2226  // Set the cell energy deposition fraction map:
2227  //
2228  if( parentMult > 0 && fSetCellMCLabelFromEdepFrac )
2229  {
2230  UInt_t * mcEdepFracPerCell = new UInt_t[ncellsTrue];
2231 
2232  // Get the digit that originated this cell cluster
2233 // AliVCaloCells* cells = 0x0;
2234 // if (aodIH && aodIH->GetMergeEvents()) cells = AODEvent() ->GetEMCALCells();
2235 // else cells = InputEvent()->GetEMCALCells();
2236 
2237  for(Int_t icell = 0; icell < ncellsTrue ; icell++)
2238  {
2239  Int_t idigit = fCellLabels[absIds[icell]];
2240 
2241  const AliEMCALDigit * dig = (const AliEMCALDigit*)fDigitsArr->At(idigit);
2242 
2243  // Find the 4 MC labels that contributed to the cluster and their
2244  // deposited energy in the current digit
2245 
2246  mcEdepFracPerCell[icell] = 0; // init
2247 
2248  Int_t nparents = dig->GetNiparent();
2249  if ( nparents > 0 )
2250  {
2251  Int_t digLabel =-1 ;
2252  Float_t edep = 0 ;
2253  Float_t edepTot = 0 ;
2254  Float_t mcEDepFrac[4] = {0,0,0,0};
2255 
2256  // all parents in digit
2257  for ( Int_t jndex = 0 ; jndex < nparents ; jndex++ )
2258  {
2259  digLabel = dig->GetIparent (jndex+1);
2260  edep = dig->GetDEParent(jndex+1);
2261  edepTot += edep;
2262 
2263  if ( digLabel == parentList[0] ) mcEDepFrac[0] = edep;
2264  else if ( digLabel == parentList[1] ) mcEDepFrac[1] = edep;
2265  else if ( digLabel == parentList[2] ) mcEDepFrac[2] = edep;
2266  else if ( digLabel == parentList[3] ) mcEDepFrac[3] = edep;
2267  } // all prarents in digit
2268 
2269  // Divide energy deposit by total deposited energy
2270  // Do this only when deposited energy is significant, use 10 MeV although 50 MeV should be expected
2271  if(edepTot > 0.01)
2272  {
2273  mcEdepFracPerCell[icell] = clus->PackMCEdepFraction(mcEDepFrac);
2274  }
2275  } // at least one parent label in digit
2276  } // cell in cluster loop
2277 
2278  clus->SetCellsMCEdepFractionMap(mcEdepFracPerCell);
2279 
2280  delete [] mcEdepFracPerCell;
2281 
2282  } // at least one parent in cluster, do the cell primary packing
2283  }
2284  } // recPoints loop
2285 }
2286 
2287 //_____________________________________________________________________
2290 //_____________________________________________________________________
2292 {
2293  if(label < 0) return ;
2294 
2295  AliAODEvent * evt = dynamic_cast<AliAODEvent*> (fEvent) ;
2296  if(!evt) return ;
2297 
2298  TClonesArray * arr = dynamic_cast<TClonesArray*>(evt->FindListObject("mcparticles")) ;
2299  if(!arr) return ;
2300 
2301  if(label < arr->GetEntriesFast())
2302  {
2303  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(label));
2304  if(!particle) return ;
2305 
2306  if(label == particle->Label()) return ; // label already OK
2307  //else printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - Label %d - AOD stack %d \n",label, particle->Label());
2308  }
2309  //else printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - Label %d > AOD labels %d \n",label, arr->GetEntriesFast());
2310 
2311  // loop on the particles list and check if there is one with the same label
2312  for(Int_t ind = 0; ind < arr->GetEntriesFast(); ind++ )
2313  {
2314  AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
2315  if(!particle) continue ;
2316 
2317  if(label == particle->Label())
2318  {
2319  label = ind;
2320  //printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - New Label Index %d \n",label);
2321  return;
2322  }
2323  }
2324 
2325  label = -1;
2326 
2327  //printf("AliAnalysisTaskEMCALClusterize::RemapMCLabelForAODs() - Label not found set to -1 \n");
2328 }
2329 
2330 //________________________________________________
2332 //________________________________________________
2334 {
2335  for(Int_t j = 0; j < fgkNEMCalCells; j++)
2336  {
2337  fOrgClusterCellId[j] = -1;
2338  fCellLabels[j] = -1;
2339  fCellSecondLabels[j] = -1;
2340  fCellTime[j] = 0.;
2341  fCellMatchdEta[j] = -999;
2342  fCellMatchdPhi[j] = -999;
2343 
2344  fTCardCorrCellsEner[j] = 0.;
2345  fTCardCorrCellsNew [j] = kFALSE;
2346  }
2347 }
2348 
2349 //_____________________________________________________________________________________________________
2353 //_____________________________________________________________________________________________________
2355  AliAODCaloCluster * clus)
2356 {
2357  Int_t parentMult = 0;
2358  Int_t *parentList = recPoint->GetParents(parentMult);
2359  clus->SetLabel(parentList, parentMult);
2360 
2361  //Write the second major contributor to each MC cluster.
2362  Int_t iNewLabel ;
2363  for ( Int_t iLoopCell = 0 ; iLoopCell < clus->GetNCells() ; iLoopCell++ )
2364  {
2365 
2366  Int_t idCell = clus->GetCellAbsId(iLoopCell) ;
2367  if(idCell>=0)
2368  {
2369  iNewLabel = 1 ; //iNewLabel makes sure we don't write twice the same label.
2370  for ( UInt_t iLoopLabels = 0 ; iLoopLabels < clus->GetNLabels() ; iLoopLabels++ )
2371  {
2372  if ( fCellSecondLabels[idCell] == -1 ) iNewLabel = 0; // -1 is never a good second label.
2373  if ( fCellSecondLabels[idCell] == clus->GetLabelAt(iLoopLabels) ) iNewLabel = 0;
2374  }
2375  if (iNewLabel == 1)
2376  {
2377  Int_t * newLabelArray = new Int_t[clus->GetNLabels()+1] ;
2378  for ( UInt_t iLoopNewLabels = 0 ; iLoopNewLabels < clus->GetNLabels() ; iLoopNewLabels++ )
2379  {
2380  newLabelArray[iLoopNewLabels] = clus->GetLabelAt(iLoopNewLabels) ;
2381  }
2382 
2383  newLabelArray[clus->GetNLabels()] = fCellSecondLabels[idCell] ;
2384  clus->SetLabel(newLabelArray,clus->GetNLabels()+1) ;
2385  delete [] newLabelArray;
2386  }
2387  } // positive cell number
2388  }
2389 }
2390 
2391 //___________________________________________________________________________________________________
2395 //___________________________________________________________________________________________________
2397 {
2398  TArrayI clArray(300) ; //Weird if more than a few clusters are in the origin ...
2399  clArray.Reset();
2400  Int_t nClu = 0;
2401  Int_t nLabTotOrg = 0;
2402  Float_t emax = -1;
2403  Int_t idMax = -1;
2404 
2405  AliVEvent * event = fEvent;
2406 
2407  // In case of embedding the MC signal is in the added event
2408  AliAODInputHandler* aodIH = dynamic_cast<AliAODInputHandler*>((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
2409  if(aodIH && aodIH->GetEventToMerge()) //Embedding
2410  event = aodIH->GetEventToMerge(); //Get clusters directly from embedded signal
2411 
2412 
2413  // Find the clusters that originally had the cells
2414  for ( Int_t iLoopCell = 0 ; iLoopCell < clus->GetNCells() ; iLoopCell++ )
2415  {
2416  Int_t idCell = clus->GetCellAbsId(iLoopCell) ;
2417 
2418  if(idCell>=0)
2419  {
2420  Int_t idCluster = fOrgClusterCellId[idCell];
2421 
2422  Bool_t set = kTRUE;
2423  for(Int_t icl =0; icl < nClu; icl++)
2424  {
2425  if(((Int_t)clArray.GetAt(icl))==-1 ) continue;
2426  if( idCluster == ((Int_t)clArray.GetAt(icl)) ) set = kFALSE;
2427  // printf("\t \t icell %d Cluster in array %d, IdCluster %d, in array %d, set %d\n",
2428  // iLoopCell, icl, idCluster,((Int_t)clArray.GetAt(icl)),set);
2429  }
2430  if( set && idCluster >= 0)
2431  {
2432  clArray.SetAt(idCluster,nClu++);
2433  //printf("******** idCluster %d \n",idCluster);
2434  nLabTotOrg+=(event->GetCaloCluster(idCluster))->GetNLabels();
2435 
2436  //printf("Cluster in array %d, IdCluster %d\n",nClu-1, idCluster);
2437 
2438  //Search highest E cluster
2439  AliVCluster * clOrg = event->GetCaloCluster(idCluster);
2440  //printf("\t E %f\n",clOrg->E());
2441  if(emax < clOrg->E())
2442  {
2443  emax = clOrg->E();
2444  idMax = idCluster;
2445  }
2446  }
2447  }
2448  }// cell loop
2449 
2450  if(nClu==0 || nLabTotOrg == 0)
2451  {
2452  //if(clus->E() > 0.25) printf("AliAnalysisTaskEMCALClusterize::SetClustersMCLabelFromOriginalClusters() - Check: N org clusters %d, n tot labels %d, cluster E %f, n cells %d\n",nClu,nLabTotOrg,clus->E(), clus->GetNCells());
2453  //for(Int_t icell = 0; icell < clus->GetNCells(); icell++) printf("\t cell %d",clus->GetCellsAbsId()[icell]);
2454  //printf("\n");
2455  }
2456 
2457  // Put the first in the list the cluster with highest energy
2458  if(idMax != ((Int_t)clArray.GetAt(0))) // Max not at first position
2459  {
2460  Int_t maxIndex = -1;
2461  Int_t firstCluster = ((Int_t)clArray.GetAt(0));
2462  for ( Int_t iLoopCluster = 0 ; iLoopCluster < nClu ; iLoopCluster++ )
2463  {
2464  if(idMax == ((Int_t)clArray.GetAt(iLoopCluster))) maxIndex = iLoopCluster;
2465  }
2466 
2467  if(firstCluster >=0 && idMax >=0)
2468  {
2469  clArray.SetAt(idMax,0);
2470  clArray.SetAt(firstCluster,maxIndex);
2471  }
2472  }
2473 
2474  // Get the labels list in the original clusters, assign all to the new cluster
2475  TArrayI clMCArray(nLabTotOrg) ;
2476  clMCArray.Reset();
2477 
2478  Int_t nLabTot = 0;
2479  for ( Int_t iLoopCluster = 0 ; iLoopCluster < nClu ; iLoopCluster++ )
2480  {
2481  Int_t idCluster = (Int_t) clArray.GetAt(iLoopCluster);
2482  //printf("New Cluster in Array %d, idCluster %d \n",iLoopCluster,idCluster);
2483  AliVCluster * clOrg = event->GetCaloCluster(idCluster);
2484  Int_t nLab = clOrg->GetNLabels();
2485 
2486  for ( Int_t iLab = 0 ; iLab < nLab ; iLab++ )
2487  {
2488  Int_t lab = clOrg->GetLabelAt(iLab) ;
2489  if(lab>=0)
2490  {
2491  Bool_t set = kTRUE;
2492  //printf("\t \t Set Label %d \n", lab);
2493  for(Int_t iLabTot =0; iLabTot < nLabTot; iLabTot++)
2494  {
2495  if( lab == ((Int_t)clMCArray.GetAt(iLabTot)) ) set = kFALSE;
2496  //printf("iLoopCluster %d, Label ID in Org Cluster %d,label %d Label ID in array %d, label in array %d, set %d\n",
2497  // iLoopCluster, iLab, lab, iLabTot, ((Int_t)clMCArray.GetAt(iLabTot)),set);
2498  }
2499  if( set ) clMCArray.SetAt(lab,nLabTot++);
2500  }
2501  }
2502  }// cluster loop
2503 
2504  // Set the final list of labels
2505 
2506  clus->SetLabel(clMCArray.GetArray(), nLabTot);
2507 
2508 // printf("Final list of labels for new cluster : \n");
2509 // for(Int_t ice = 0; ice < clus->GetNCells() ; ice++)
2510 // {
2511 // printf("\t Cell %d ",clus->GetCellsAbsId()[ice]);
2512 // Int_t label = InputEvent()->GetEMCALCells()->GetCellMCLabel(clus->GetCellsAbsId()[ice]);
2513 // printf(" org %d ",label);
2514 // RemapMCLabelForAODs(label);
2515 // printf(" new %d \n",label);
2516 // }
2517 // for(Int_t imc = 0; imc < clus->GetNLabels(); imc++) printf("\t Label %d\n",clus->GetLabelAt(imc));
2518 }
2519 
2520 //____________________________________________________________
2521 // Init geometry, create list of output clusters and cells.
2522 //____________________________________________________________
2524 {
2525  // Clusters
2526  fOutputAODBranch = new TClonesArray("AliAODCaloCluster", 0);
2527 
2528  if(fOutputAODBranchName.Length()==0)
2529  {
2530  fOutputAODBranchName = "newEMCALClustersArray";
2531  AliInfo("Cluster branch name not set, set it to newEMCALClustersArray");
2532  }
2533 
2535 
2536  if( fFillAODFile )
2537  {
2538  //fOutputAODBranch = new TClonesArray("AliAODCaloCluster", 0);
2539 
2540  //fOutputAODBranch->SetOwner(kFALSE);
2541 
2542  AddAODBranch("TClonesArray", &fOutputAODBranch);
2543  }
2544 
2545  // Cells
2546  if ( fOutputAODBranchName.Length() > 0 )
2547  {
2548  fOutputAODCells = new AliAODCaloCells(fOutputAODCellsName,fOutputAODCellsName,AliAODCaloCells::kEMCALCell);
2549 
2550  if( fFillAODFile ) AddAODBranch("AliAODCaloCells", &fOutputAODCells);
2551  }
2552 }
2553 
2554 //_______________________________________________________________________
2557 //________________________________________________________________________
2559 {
2560  // Update cells only in case re-calibration was done
2561  // or bad map applied or additional T-Card cells added.
2562  if(!fRecoUtils->IsBadChannelsRemovalSwitchedOn() &&
2563  !fRecoUtils->IsRecalibrationOn() &&
2564  !fRecoUtils->IsRunDepRecalibrationOn() &&
2565  !fRecoUtils->IsTimeRecalibrationOn() &&
2566  !fRecoUtils->IsL1PhaseInTimeRecalibrationOn() &&
2567  !fTCardCorrEmulation ) return;
2568 
2569  const Int_t ncells = fCaloCells->GetNumberOfCells();
2570  const Int_t ndigis = fDigitsArr->GetEntries();
2571 
2572  if ( fOutputAODCellsName.Length() > 0 )
2573  {
2574  fOutputAODCells->DeleteContainer();
2575  fOutputAODCells->CreateContainer(ndigis);
2576  }
2577  else if ( ncells != ndigis ) // update case
2578  {
2579  fCaloCells->DeleteContainer();
2580  fCaloCells->CreateContainer(ndigis);
2581  }
2582 
2583  for (Int_t idigit = 0; idigit < ndigis; ++idigit)
2584  {
2585  AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(idigit));
2586 
2587  Double_t cellAmplitude = digit->GetAmplitude();
2588  Short_t cellNumber = digit->GetId();
2589  Double_t cellTime = digit->GetTime();
2590 
2591  Bool_t highGain = kFALSE;
2592  if( digit->GetType() == AliEMCALDigit::kHG ) highGain = kTRUE;
2593 
2594  // Only for MC
2595  // Get the label of the primary particle that generated the cell
2596  // Assign the particle that deposited more energy
2597  Int_t nparents = digit->GetNiparent();
2598  Int_t cellMcEDepFrac =-1 ;
2599  Float_t cellMcLabel =-1.;
2600  if ( nparents > 0 )
2601  {
2602  for ( Int_t jndex = 0 ; jndex < nparents ; jndex++ )
2603  {
2604  if(cellMcEDepFrac >= digit->GetDEParent(jndex+1)) continue ;
2605 
2606  cellMcLabel = digit->GetIparent (jndex+1);
2607  cellMcEDepFrac= digit->GetDEParent(jndex+1);
2608  } // all primaries in digit
2609  } // select primary label
2610 
2611  if ( cellMcEDepFrac < 0 ) cellMcEDepFrac = 0.;
2612 
2613  if ( fUpdateCell )
2614  fCaloCells ->SetCell(idigit, cellNumber, cellAmplitude, cellTime, cellMcLabel, cellMcEDepFrac, highGain);
2615  else
2616  fOutputAODCells->SetCell(idigit, cellNumber, cellAmplitude, cellTime, cellMcLabel, cellMcEDepFrac, highGain);
2617  }
2618 }
2619 
2620 //_______________________________________________________
2631 //_______________________________________________________
2633 {
2634  //-------
2635  // Step 1
2636 
2637  // Remove the contents of AOD branch output list set in the previous event
2638  fOutputAODBranch->Clear("C");
2639 
2640  LoadBranches();
2641 
2642  // Check if there is a centrality value, PbPb analysis, and if a centrality bin selection is requested
2643  // If we need a centrality bin, we select only those events in the corresponding bin.
2644  if( fCentralityBin[0] >= 0 && fCentralityBin[1] >= 0 )
2645  {
2646  Float_t cen = GetEventCentrality();
2647  if(cen > fCentralityBin[1] || cen < fCentralityBin[0]) return ; //reject events out of bin.
2648  }
2649 
2650  // Intermediate array with new clusters : init the array only once or clear from previous event
2651  if(!fCaloClusterArr) fCaloClusterArr = new TObjArray(10000);
2652  else fCaloClusterArr->Delete();//Clear("C"); it leaks?
2653 
2654 
2655  // In case of analysis in multiple runs, check the OADB again
2656  if ( InputEvent()->GetRunNumber() != fRun )
2657  {
2658  fRun = InputEvent()->GetRunNumber();
2659 
2660  fOADBSet = kFALSE; // recover the OADB for this run
2661 
2662  AliInfo(Form("Set run to %d",fRun));
2663  }
2664 
2665  InitGeometry(); // only once, must be done before OADB, geo OADB accessed here
2666 
2667  // Get the event, do some checks and settings
2668  CheckAndGetEvent() ;
2669 
2670  if (!fEvent)
2671  {
2672  AliDebug(1,Form("Skip Event %d", (Int_t) Entry()));
2673  return ;
2674  }
2675 
2676  // Init pointers, geometry, clusterizer, ocdb, aodb
2677 
2678  if(fAccessOCDB) AccessOCDB();
2679  if(fAccessOADB) AccessOADB(); // only once
2680 
2682 
2683  // Print once the analysis parameters
2684  if ( fDebug > 0 || !fPrintOnce )
2685  {
2686  fRecParam->Print("reco");
2687 
2688  PrintParam();
2689 
2690  PrintTCardParam();
2691 
2692  fPrintOnce = kTRUE;
2693  }
2694 
2695  //-------
2696  // Step 2
2697 
2698  // Make clusters
2700  else ClusterizeCells() ;
2701 
2702  //-------
2703  // Step 3
2704 
2706 
2707  if ( fUpdateCell || fOutputAODCellsName.Length() > 0 )
2708  UpdateCells();
2709 
2710 }
2711 
2712 
TString fOutputAODCellsName
New of output cells AOD branch name.
void FillAODCaloCells()
Put calo cells in standard branch.
Float_t fTCardCorrMinAmp
Minimum cell energy to induce signal on adjacent cells.
AliVCaloCells * fCaloCells
! CaloCells container
Bool_t fTCardCorrEmulation
Activate T-Card cells energy correlation.
Float_t fTCardCorrInduceEnerFracP1[4][22]
Induced energy loss gauss fraction param1 on 0-same row, diff col, 1-up/down cells left/right col 2-l...
double Double_t
Definition: External.C:58
virtual void ResetArrays()
Reset arrays containing information for all possible cells.
Definition: External.C:236
const int debug
Definition: scanAll.C:15
TString GetPass()
Get or guess pass number/string from path of filename.
Bool_t fFillAODHeader
Copy header to standard branch.
TClonesArray * fDigitsArr
! Digits array
Bool_t fTCardCorrClusEnerConserv
When making correlation, subtract from the reference cell the induced energy on the neighbour cells...
AliEMCALClusterizer * fClusterizer
! EMCAL clusterizer
Bool_t fJustUnfold
Just unfold, do not recluster.
Bool_t fImportGeometryFromFile
Import geometry settings in geometry.root file.
Bool_t fSelectCell
Reject cells from cluster if energy is too low and recalculate position/energy and other...
TString fConfigName
Name of analysis configuration file.
Float_t fTCardCorrInduceEnerFracMax[22]
In case fTCardCorrInduceEnerFracP1 is non null, restrict the maximum fraction of induced energy per S...
TCanvas * c
Definition: TestFitELoss.C:172
void PrintTCardParam()
Print parameters for T-Card correlation emulation.
Int_t fCellSecondLabels[fgkNEMCalCells]
Array with Second MC label to be passed to digit.
void SetClustersMCLabelFrom2SelectedLabels(AliEMCALRecPoint *recPoint, AliAODCaloCluster *clus)
Float_t fCellMatchdEta[fgkNEMCalCells]
Array with cluster-track dPhi.
TClonesArray * fOutputAODBranch
! AOD Branch with output clusters
TString fOADBFilePath
Default path $ALICE_PHYSICS/OADB/EMCAL, if needed change.
Float_t fTCardCorrCellsEner[fgkNEMCalCells]
Array with induced cell energy in T-Card neighbour cells.
TString fImportGeometryFilePath
path fo geometry.root file
Float_t fTCardCorrInduceEnerFrac[4][22]
Induced energy loss gauss fraction param0 on 0-same row, diff col, 1-up/down cells left/right col 2-l...
void ConfigureEMCALRecoUtils(Bool_t bMC=kFALSE, Bool_t bExotic=kTRUE, Bool_t bNonLin=kFALSE, Bool_t bRecalE=kTRUE, Bool_t bBad=kTRUE, Bool_t bRecalT=kTRUE, Int_t debug=-1)
Bool_t fPrintOnce
Print once analysis parameters.
Bool_t fSelectEMCALEvent
Process the event if there is some high energy cluster.
Float_t fTCardCorrMaxInduced
Maximum induced energy signal on adjacent cells.
Float_t GetEventCentrality() const
Get centrality/multiplicity percentile.
Bool_t fInputFromFilter
Get the input from AODs from the filter.
TString fOutputAODBranchName
New of output clusters AOD branch.
AliEMCALGeometry * fGeom
EMCAL geometry.
void ClusterUnfolding()
Take the event clusters and unfold them.
int Int_t
Definition: External.C:63
Bool_t fAccessOADB
Get calibration from OADB for EMCAL.
Bool_t fRemoveLEDEvents
Remove LED events, use only for LHC11a.
unsigned int UInt_t
Definition: External.C:33
AliEMCALRecoUtils * fRecoUtils
Access to factorized reconstruction algorithms.
float Float_t
Definition: External.C:68
Bool_t fRejectBelowThreshold
split (false-default) or reject (true) cell energy below threshold after UF
Float_t fConstantTimeShift
Apply a 600 ns time shift in case of simulation, shift in ns.
void SetClustersMCLabelFromOriginalClusters(AliAODCaloCluster *clus)
Bool_t fUpdateCell
On/Off the upate of the CaloCells container.
Float_t fSelectCellMinE
Min energy cell threshold, after unfolding.
Bool_t fOADBSet
AODB parameters already set.
AliAODCaloCells * fOutputAODCells
! AOD Branch with output cells
Float_t fTCardCorrInduceEner[4][22]
Induced energy loss gauss constant on 0-same row, diff col, 1-up/down cells left/right col 2-left/rig...
Int_t fCellLabels[fgkNEMCalCells]
Array with MC label to be passed to digit.
Bool_t fRemoveExoticEvents
Remove exotic events.
Bool_t fDoTrackMatching
On/Off the matching recalculation to speed up analysis in PbPb.
Float_t fCellMatchdPhi[fgkNEMCalCells]
Array with cluster-track dEta.
Float_t fTCardCorrInduceEnerProb[22]
Probability to induce energy loss per SM.
TString fInputCaloCellsName
Input cells branch name, if different from default branch.
TGeoHMatrix * fGeomMatrix[22]
Geometry matrices with alignments.
Bool_t IsLEDEvent(const Int_t run)
Check if event is LED, is so remove it. Affected LHC11a runs.
AliEMCALRecParam * fRecParam
Reconstruction parameters container.
TString fGeomName
Name of geometry to use.
Bool_t fOutputAODBranchSet
Set the AOD clusters branch in the input event once.
void PrintParam()
Print clusterization task parameters.
short Short_t
Definition: External.C:23
TString fOCDBpath
Path with OCDB location.
AliEMCALAfterBurnerUF * fUnfolder
! Unfolding procedure
Bool_t fGeomMatrixSet
Set geometry matrices only once, for the first event.
Float_t fCentralityBin[2]
Minimum and maximum value of the centrality for the analysis.
Bool_t fFillAODCaloCells
Copy calocells to standard branch.
TObjArray * fCaloClusterArr
! CaloClusters array
Double_t fCellTime[fgkNEMCalCells]
Array with cluster time to be passed to digit in case of AODs.
TString fCentralityClass
Name of selected centrality class.
TFile * file
TList with histograms for a given trigger.
Bool_t fAccessOCDB
Need to access info from OCDB (not really)
TObjArray * fClusterArr
! Recpoints array
unsigned short UShort_t
Definition: External.C:28
Bool_t fRecalibrateWithClusterTime
Use fCellTime to store time of cells in cluster.
Int_t GetRunNumber(TString)
Definition: PlotMuonQA.C:2235
const char Option_t
Definition: External.C:48
bool Bool_t
Definition: External.C:53
void FillAODHeader()
Put event header information in standard AOD branch.
static const Int_t fgkNEMCalCells
Total number of cells in the calorimeter, 10*48*24 (EMCal) + 4*48*8 (EMCal/DCal 1/3) + 6*32*24 (DCal)...
Float_t fTCardCorrInduceEnerFracMin[22]
In case fTCardCorrInduceEnerFracP1 is non null, restrict the minimum fraction of induced energy per S...
Int_t fEMCALNcellsCut
At least an EMCAL cluster with fNCellsCut cells over fEnergyCut.
Reclusterize EMCal clusters, put them in a new branch for other following analysis.
Bool_t fTCardCorrCellsNew[fgkNEMCalCells]
Array with induced cell energy in T-Card neighbour cells, that before had no signal.
Bool_t fLoadGeomMatrices
Matrices set from configuration, not get from geometry.root or from ESDs/AODs.
Float_t fEMCALEnergyCut
At least an EMCAL cluster with this energy in the event.
Bool_t fRandomizeTCard
Use random induced energy.
Bool_t fUseAliCentrality
Use the centrality estimator from AliCentrality or AliMultSelection.
Float_t fTCardCorrInduceEnerFracWidth[4][22]
Induced energy loss gauss witdth on 0-same row, diff col, 1-up/down cells left/right col 2-left/righ ...
Int_t fOrgClusterCellId[fgkNEMCalCells]
Array ID of cluster to wich the cell belongs in unmodified clusters.
Float_t fSelectCellMinFrac
Min fraction of cell energy after unfolding cut.
Bool_t fRemapMCLabelForAODs
Remap AOD cells MC label. Needed in old AOD productions.