AliPhysics  32b88a8 (32b88a8)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliEmcalCorrectionClusterizer.cxx
Go to the documentation of this file.
1 // AliEmcalCorrectionClusterizer
2 //
3 /**************************************************************************
4  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5  * *
6  * Author: The ALICE Off-line Project. *
7  * Contributors are mentioned in the code where appropriate. *
8  * *
9  * Permission to use, copy, modify and distribute this software and its *
10  * documentation strictly for non-commercial purposes is hereby granted *
11  * without fee, provided that the above copyright notice appears in all *
12  * copies and that both the copyright notice and this permission notice *
13  * appear in the supporting documentation. The authors make no claims *
14  * about the suitability of this software for any purpose. It is *
15  * provided "as is" without express or implied warranty. *
16  **************************************************************************/
17 
18 // --- Root ---
19 #include <TGeoManager.h>
20 #include <TObjArray.h>
21 #include <TString.h>
22 #include <TTree.h>
23 #include <TArrayI.h>
24 
25 // --- AliRoot ---
26 #include "AliCDBEntry.h"
27 #include "AliCDBManager.h"
28 #include "AliCaloCalibPedestal.h"
29 #include "AliEMCALAfterBurnerUF.h"
30 #include "AliEMCALCalibData.h"
31 #include "AliEMCALClusterizerNxN.h"
32 #include "AliEMCALClusterizerv1.h"
33 #include "AliEMCALClusterizerv2.h"
34 #include "AliEMCALClusterizerFixedWindow.h"
35 #include "AliEMCALDigit.h"
36 #include "AliEMCALGeometry.h"
37 #include "AliEMCALRecParam.h"
38 #include "AliEMCALRecPoint.h"
39 #include "AliInputEventHandler.h"
40 
42 
46 
47 // Actually registers the class with the base class
49 
50 //________________________________________________________________________
52  AliEmcalCorrectionComponent("AliEmcalCorrectionClusterizer"),
53  fDigitsArr(0),
54  fClusterArr(0),
55  fRecParam(new AliEMCALRecParam),
56  fClusterizer(0),
57  fUnfolder(0),
58  fJustUnfold(kFALSE),
59  fGeomName(),
60  fGeomMatrixSet(kFALSE),
61  fLoadGeomMatrices(kFALSE),
62  fOCDBpath(),
63  fCalibData(0),
64  fPedestalData(0),
65  fLoadCalib(kFALSE),
66  fLoadPed(kFALSE),
67  fSubBackground(kFALSE),
68  fNPhi(4),
69  fNEta(4),
70  fShiftPhi(2),
71  fShiftEta(2),
72  fTRUShift(0),
73  fInputCellType(kFEEData),
74  fSetCellMCLabelFromCluster(0),
75  fSetCellMCLabelFromEdepFrac(0),
76  fCaloClusters(0),
77  fEsd(0),
78  fAod(0),
79  fRecalDistToBadChannels(kFALSE),
80  fRecalShowerShape(kFALSE)
81 {
82  // Default constructor
83  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
84 
85  for(Int_t i = 0; i < AliEMCALGeoParams::fgkEMCALModules; i++) fGeomMatrix[i] = 0 ;
86  for(Int_t j = 0; j < fgkTotalCellNumber; j++)
87  { fOrgClusterCellId[j] =-1; fCellLabels[j] =-1 ; }
88 }
89 
90 //________________________________________________________________________
92 {
93  // Destructor
94 
95  delete fClusterizer;
96  delete fUnfolder;
97  delete fRecParam;
98 }
99 
100 //________________________________________________________________________
102 {
103  // Initialization
104  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
106  // Do base class initializations and if it fails -> bail out
107  //AliAnalysisTaskEmcal::ExecOnce();
108  //if (!fInitialized) return;
109 
110  std::string clusterizerTypeStr = "";
111  GetProperty("clusterizer", clusterizerTypeStr);
112  UInt_t clusterizerType = clusterizerTypeMap.at(clusterizerTypeStr);
113  Double_t cellE = 0.05;
114  GetProperty("cellE", cellE);
115  Double_t seedE = 0.1;
116  GetProperty("seedE", seedE);
117  Float_t timeMin = -1; //minimum time of physical signal in a cell/digit (s) (in run macro, -50e-6)
118  GetProperty("cellTimeMin", timeMin);
119  Float_t timeMax = +1; //maximum time of physical signal in a cell/digit (s) (in run macro, 50e-6)
120  GetProperty("cellTimeMax", timeMax);
121  Float_t timeCut = 1; //maximum time difference between the digits inside EMC cluster (s) (in run macro, 1e-6)
122  GetProperty("clusterTimeLength", timeCut);
123  Float_t w0 = 4.5;
124  GetProperty("w0", w0);
125  GetProperty("recalDistToBadChannels", fRecalDistToBadChannels);
126  GetProperty("recalShowerShape", fRecalShowerShape);
127 
129 
130  fRecParam->SetClusterizerFlag(clusterizerType);
131  fRecParam->SetMinECut(cellE);
132  fRecParam->SetClusteringThreshold(seedE);
133  fRecParam->SetW0(w0);
134  fRecParam->SetTimeMin(timeMin);
135  fRecParam->SetTimeMax(timeMax);
136  fRecParam->SetTimeCut(timeCut);
137 
138  if (clusterizerType == AliEMCALRecParam::kClusterizerNxN)
139  fRecParam->SetNxM(1,1); // -> (1,1) means 3x3!
140 
142  Printf("inputCellType: %d",fInputCellType);
143 
144  return kTRUE;
145 }
146 
147 //________________________________________________________________________
149 {
150  // Run
151  AliDebug(3, Form("%s", __PRETTY_FUNCTION__));
153 
154  // Main loop, called for each event
155 
156  fEsd = dynamic_cast<AliESDEvent*>(fEvent);
157  fAod = dynamic_cast<AliAODEvent*>(fEvent);
158 
159  fCaloClusters = fClusCont->GetArray();
160 
161  UInt_t offtrigger = 0;
162  if (fEsd) {
163  UInt_t mask1 = fEsd->GetESDRun()->GetDetectorsInDAQ();
164  UInt_t mask2 = fEsd->GetESDRun()->GetDetectorsInReco();
165  Bool_t desc1 = (mask1 >> 18) & 0x1;
166  Bool_t desc2 = (mask2 >> 18) & 0x1;
167  if (desc1==0 || desc2==0) { //AliDAQ::OfflineModuleName(180=="EMCAL"
168  AliError(Form("EMCAL not in DAQ/RECO: %u (%u)/%u (%u)",
169  mask1, fEsd->GetESDRun()->GetDetectorsInReco(),
170  mask2, fEsd->GetESDRun()->GetDetectorsInDAQ()));
171  return kFALSE;
172  }
173  AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
174  offtrigger = ((AliInputEventHandler*)(am->GetInputEventHandler()))->IsEventSelected();
175  } else {
176  offtrigger = ((AliVAODHeader*)fAod->GetHeader())->GetOfflineTrigger();
177  }
178 
179  if (!fMCEvent) {
180  if (offtrigger & AliVEvent::kFastOnly) {
181  AliError(Form("EMCAL not in fast only partition"));
182  return kFALSE;
183  }
184  }
185 
186  Init();
187 
188  if (fJustUnfold) {
189  AliWarning("Unfolding not implemented");
190  return kTRUE;
191  }
192 
193  FillDigitsArray();
194 
195  Clusterize();
196 
197  UpdateClusters();
198 
200 
201  return kTRUE;
202 }
203 
204 //________________________________________________________________________
206 {
207  // Clusterize
208 
209  if (fSubBackground) {
210  fClusterizer->SetInputCalibrated(kTRUE);
211  fClusterizer->SetCalibrationParameters(0);
212  }
213 
214  fClusterizer->Digits2Clusters("");
215 
216  if (fSubBackground) {
217  if (fCalibData) {
218  fClusterizer->SetInputCalibrated(kFALSE);
219  fClusterizer->SetCalibrationParameters(fCalibData);
220  }
221  }
222 }
223 
224 //________________________________________________________________________
226 {
227  // Fill digits array
228 
229  fDigitsArr->Clear("C");
230  switch (fInputCellType) {
231 
232  case kFEEData :
233  case kFEEDataMCOnly :
234  case kFEEDataExcludeMC :
235  {
236  // In case of MC productions done before aliroot tag v5-02-Rev09
237  // passing the cluster label to all the cells belonging to this cluster
238  // very rough
239  // Copied and simplified from AliEMCALTenderSupply
241  {
242  for (Int_t i = 0; i < fgkTotalCellNumber; i++)
243  {
244  fCellLabels [i] =-1 ;
245  fOrgClusterCellId[i] =-1 ;
246  }
247 
248  Int_t nClusters = fEvent->GetNumberOfCaloClusters();
249  for (Int_t i = 0; i < nClusters; i++)
250  {
251  AliVCluster *clus = fEvent->GetCaloCluster(i);
252 
253  if (!clus) continue;
254 
255  if (!clus->IsEMCAL()) continue ;
256 
257  Int_t label = clus->GetLabel();
258  UShort_t * index = clus->GetCellsAbsId() ;
259 
260  for(Int_t icell=0; icell < clus->GetNCells(); icell++)
261  {
263  fCellLabels[index[icell]] = label;
264 
265  fOrgClusterCellId[index[icell]] = i ; // index of the original cluster
266  } // cell in cluster loop
267  } // cluster loop
268  }
269 
270  Double_t avgE = 0; // for background subtraction
271  const Int_t ncells = fCaloCells->GetNumberOfCells();
272  for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell)
273  {
274  Double_t cellAmplitude=0, cellTime=0, cellEFrac = 0;
275  Short_t cellNumber=0;
276  Int_t cellMCLabel=-1;
277  if (fCaloCells->GetCell(icell, cellNumber, cellAmplitude, cellTime, cellMCLabel, cellEFrac) != kTRUE)
278  break;
279 
280  if (fSetCellMCLabelFromCluster) cellMCLabel = fCellLabels[cellNumber];
281 
282  if (cellMCLabel > 0 && cellEFrac < 1e-6)
283  cellEFrac = 1;
284 
285  if (cellAmplitude < 1e-6 || cellNumber < 0)
286  continue;
287 
289  if (cellMCLabel <= 0)
290  continue;
291  else {
292  cellAmplitude *= cellEFrac;
293  cellEFrac = 1;
294  }
295  }
296  else if (fInputCellType == kFEEDataExcludeMC) {
297  if (cellMCLabel > 0)
298  continue;
299  else {
300  cellAmplitude *= 1 - cellEFrac;
301  cellEFrac = 0;
302  }
303  }
304 
305  AliEMCALDigit *digit = new((*fDigitsArr)[idigit]) AliEMCALDigit(cellMCLabel, cellMCLabel, cellNumber,
306  (Float_t)cellAmplitude, (Float_t)cellTime,
307  AliEMCALDigit::kHG,idigit, 0, 0, cellEFrac*cellAmplitude);
308 
309  if (fSubBackground)
310  {
311  Float_t energy = cellAmplitude;
312  Float_t time = cellTime;
313  fClusterizer->Calibrate(energy,time,cellNumber);
314  digit->SetAmplitude(energy);
315  avgE += energy;
316  }
317 
318  // New way to set the cell MC labels,
319  // valid only for MC productions with aliroot > v5-07-21
320  if(fSetCellMCLabelFromEdepFrac && fOrgClusterCellId[cellNumber] >= 0) // index can be negative if noisy cell that did not form cluster
321  {
322  fCellLabels[cellNumber] = idigit;
323 
324  AliVCluster *clus = 0;
325  Int_t iclus = fOrgClusterCellId[cellNumber];
326 
327  if(iclus < 0)
328  {
329  AliInfo("Negative original cluster index, skip \n");
330  continue;
331  }
332 
333  clus = fEvent->GetCaloCluster(iclus);
334 
335  for(Int_t icluscell=0; icluscell < clus->GetNCells(); icluscell++ )
336  {
337  if(cellNumber != clus->GetCellAbsId(icluscell)) continue ;
338 
339  // Get the energy deposition fraction.
340  Float_t eDepFrac[4];
341  clus->GetCellMCEdepFractionArray(icluscell,eDepFrac);
342 
343  // Select the MC label contributing, only if enough energy deposition
344  TArrayI labeArr(0);
345  TArrayF eDepArr(0);
346  Int_t nLabels = 0;
347  for(Int_t imc = 0; imc < 4; imc++)
348  {
349  if(eDepFrac[imc] > 0 && clus->GetNLabels() > imc)
350  {
351  nLabels++;
352 
353  labeArr.Set(nLabels);
354  labeArr.AddAt(clus->GetLabelAt(imc), nLabels-1);
355 
356  eDepArr.Set(nLabels);
357  eDepArr.AddAt(eDepFrac[imc]*cellAmplitude, nLabels-1);
358  // use as deposited energy a fraction of the simulated energy (smeared and with noise)
359  }
360  }
361 
362  if(nLabels > 0)
363  {
364  digit->SetListOfParents(nLabels,labeArr.GetArray(),eDepArr.GetArray());
365  }
366  }
367  }
368 
369 
370  idigit++;
371  }
372 
373  if (fSubBackground) {
374  avgE /= fGeom->GetNumberOfSuperModules()*48*24;
375  Int_t ndigis = fDigitsArr->GetEntries();
376  for (Int_t i = 0; i < ndigis; ++i) {
377  AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(i));
378  Double_t energy = digit->GetAmplitude() - avgE;
379  if (energy<=0.001) {
380  digit->SetAmplitude(0);
381  } else {
382  digit->SetAmplitude(energy);
383  }
384  }
385  }
386  }
387  break;
388 
389  case kPattern :
390  {
391  // Fill digits from a pattern
392  Int_t maxd = fGeom->GetNCells() / 4;
393  for (Int_t idigit = 0; idigit < maxd; idigit++){
394  if (idigit % 24 == 12) idigit += 12;
395  AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
396  digit->SetId(idigit * 4);
397  digit->SetTime(600);
398  digit->SetTimeR(600);
399  digit->SetIndexInList(idigit);
400  digit->SetType(AliEMCALDigit::kHG);
401  digit->SetAmplitude(0.1);
402  }
403  }
404  break;
405 
406  case kL0FastORs :
407  case kL0FastORsTC :
408  case kL1FastORs :
409  {
410  // Fill digits from FastORs
411 
412  AliVCaloTrigger *triggers = fEvent->GetCaloTrigger("EMCAL");
413 
414  if (!triggers || !(triggers->GetEntries() > 0))
415  return;
416 
417  Int_t idigit = 0;
418  triggers->Reset();
419 
420  while ((triggers->Next())) {
421  Float_t L0Amplitude = 0;
422  triggers->GetAmplitude(L0Amplitude);
423 
424  if (L0Amplitude <= 0 && fInputCellType != kL1FastORs)
425  continue;
426 
427  Int_t L1Amplitude = 0;
428  triggers->GetL1TimeSum(L1Amplitude);
429 
430  if (L1Amplitude <= 0 && fInputCellType == kL1FastORs)
431  continue;
432 
433  Int_t triggerTime = 0;
434  Int_t ntimes = 0;
435  triggers->GetNL0Times(ntimes);
436 
437  if (ntimes < 1 && fInputCellType == kL0FastORsTC)
438  continue;
439 
440  if (ntimes > 0) {
441  Int_t trgtimes[25];
442  triggers->GetL0Times(trgtimes);
443  triggerTime = trgtimes[0];
444  }
445 
446  Int_t triggerCol = 0, triggerRow = 0;
447  triggers->GetPosition(triggerCol, triggerRow);
448 
449  Int_t find = -1;
450  fGeom->GetAbsFastORIndexFromPositionInEMCAL(triggerCol, triggerRow, find);
451 
452  if (find < 0)
453  continue;
454 
455  Int_t cidx[4] = {-1};
456  Bool_t ret = fGeom->GetCellIndexFromFastORIndex(find, cidx);
457 
458  if (!ret)
459  continue;
460 
461  Float_t triggerAmplitude = 0;
462 
463  if (fInputCellType == kL1FastORs) {
464  triggerAmplitude = 0.25 * L1Amplitude; // it will add 4 cells for 1 amplitude
465  }
466  else {
467  triggerAmplitude = L0Amplitude; // 10 bit truncated, so it is already divided by 4
468  }
469 
470  for (Int_t idxpos = 0; idxpos < 4; idxpos++) {
471  Int_t triggerNumber = cidx[idxpos];
472  AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
473  digit->SetId(triggerNumber);
474  digit->SetTime(triggerTime);
475  digit->SetTimeR(triggerTime);
476  digit->SetIndexInList(idigit);
477  digit->SetType(AliEMCALDigit::kHG);
478  digit->SetAmplitude(triggerAmplitude);
479  idigit++;
480  }
481  }
482  }
483  break;
484  }
485 }
486 
487 //________________________________________________________________________________________
489 {
490  // Cluster energy, global position, cells and their amplitude fractions are restored.
491 
492  const Int_t Ncls = fClusterArr->GetEntries();
493  AliDebug(1, Form("total no of clusters %d", Ncls));
494 
495  for(Int_t i=0, nout=clus->GetEntries(); i < Ncls; ++i)
496  {
497  AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
498 
499  Int_t ncellsTrue = 0;
500  const Int_t ncells = recpoint->GetMultiplicity();
501  UShort_t absIds[ncells];
502  Double32_t ratios[ncells];
503  Int_t *dlist = recpoint->GetDigitsList();
504  Float_t *elist = recpoint->GetEnergiesList();
505  Double_t mcEnergy = 0;
506 
507  for (Int_t c = 0; c < ncells; ++c)
508  {
509  AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
510  absIds[ncellsTrue] = digit->GetId();
511  ratios[ncellsTrue] = elist[c]/digit->GetAmplitude();
512 
513  if (digit->GetIparent(1) > 0)
514  mcEnergy += digit->GetDEParent(1)/recpoint->GetEnergy();
515 
516  ++ncellsTrue;
517  }
518 
519  if (ncellsTrue < 1)
520  {
521  AliWarning("Skipping cluster with no cells");
522  continue;
523  }
524 
525  // calculate new cluster position
526  TVector3 gpos;
527  recpoint->GetGlobalPosition(gpos);
528  Float_t g[3];
529  gpos.GetXYZ(g);
530 
531  AliDebug(1, Form("energy %f", recpoint->GetEnergy()));
532 
533  AliVCluster *c = static_cast<AliVCluster*>(clus->New(nout++));
534  c->SetType(AliVCluster::kEMCALClusterv1);
535  c->SetE(recpoint->GetEnergy());
536  c->SetPosition(g);
537  c->SetNCells(ncellsTrue);
538  c->SetCellsAbsId(absIds);
539  c->SetCellsAmplitudeFraction(ratios);
540  c->SetID(recpoint->GetUniqueID());
541  c->SetDispersion(recpoint->GetDispersion());
542  c->SetEmcCpvDistance(-1);
543  c->SetChi2(-1);
544  c->SetTOF(recpoint->GetTime()) ; //time-of-flight
545  c->SetNExMax(recpoint->GetNExMax()); //number of local maxima
546  Float_t elipAxis[2];
547  recpoint->GetElipsAxis(elipAxis);
548  c->SetM02(elipAxis[0]*elipAxis[0]);
549  c->SetM20(elipAxis[1]*elipAxis[1]);
550  c->SetMCEnergyFraction(mcEnergy);
551 
552  //
553  // MC labels
554  //
555  Int_t parentMult = 0;
556  Int_t *parentList = recpoint->GetParents(parentMult);
557  Float_t *parentListDE = recpoint->GetParentsDE(); // deposited energy
558 
559  c->SetLabel(parentList, parentMult);
560  c->SetClusterMCEdepFractionFromEdepArray(parentListDE);
561 
562  //
563  // Set the cell energy deposition fraction map:
564  //
565  if( parentMult > 0 && fSetCellMCLabelFromEdepFrac )
566  {
567  UInt_t * mcEdepFracPerCell = new UInt_t[ncellsTrue];
568 
569  // Get the digit that originated this cell cluster
570  //AliVCaloCells* cells = InputEvent()->GetEMCALCells();
571 
572  for(Int_t icell = 0; icell < ncellsTrue ; icell++)
573  {
574  Int_t idigit = fCellLabels[absIds[icell]];
575 
576  const AliEMCALDigit * dig = (const AliEMCALDigit*)fDigitsArr->At(idigit);
577 
578  // Find the 4 MC labels that contributed to the cluster and their
579  // deposited energy in the current digit
580 
581  mcEdepFracPerCell[icell] = 0; // init
582 
583  Int_t nparents = dig->GetNiparent();
584  if ( nparents > 0 )
585  {
586  Int_t digLabel =-1 ;
587  Float_t edep = 0 ;
588  Float_t edepTot = 0 ;
589  Float_t mcEDepFrac[4] = {0,0,0,0};
590 
591  // all parents in digit
592  for ( Int_t jndex = 0 ; jndex < nparents ; jndex++ )
593  {
594  digLabel = dig->GetIparent (jndex+1);
595  edep = dig->GetDEParent(jndex+1);
596  edepTot += edep;
597 
598  if ( digLabel == parentList[0] ) mcEDepFrac[0] = edep;
599  else if ( digLabel == parentList[1] ) mcEDepFrac[1] = edep;
600  else if ( digLabel == parentList[2] ) mcEDepFrac[2] = edep;
601  else if ( digLabel == parentList[3] ) mcEDepFrac[3] = edep;
602  } // all prarents in digit
603 
604  // Divide energy deposit by total deposited energy
605  // Do this only when deposited energy is significant, use 10 MeV although 50 MeV should be expected
606  if(edepTot > 0.01)
607  {
608  mcEdepFracPerCell[icell] = c->PackMCEdepFraction(mcEDepFrac);
609  }
610  } // at least one parent label in digit
611  } // cell in cluster loop
612 
613  c->SetCellsMCEdepFractionMap(mcEdepFracPerCell);
614 
615  delete [] mcEdepFracPerCell;
616 
617  } // at least one parent in cluster, do the cell primary packing
618  }
619 }
620 
621 //________________________________________________________________________
623 {
624  // Update cells in case re-calibration was done.
625 
626  const Int_t nents = fCaloClusters->GetEntries();
627  for (Int_t i=0;i<nents;++i) {
628  AliVCluster *c = static_cast<AliVCluster*>(fCaloClusters->At(i));
629  if (!c)
630  continue;
631  if (c->IsEMCAL())
632  delete fCaloClusters->RemoveAt(i);
633  }
634 
635  fCaloClusters->Compress();
636 
638 }
639 
640 //________________________________________________________________________________________
642 {
643  // Go through clusters one by one and process separate correction
644 
645  Int_t nclusters = fCaloClusters->GetEntriesFast();
646  for (Int_t icluster=0; icluster < nclusters; ++icluster) {
647  AliVCluster *clust = static_cast<AliVCluster*>(fCaloClusters->At(icluster));
648  if (!clust)
649  continue;
650 
651  // SHOWER SHAPE -----------------------------------------------
652  if (fRecalShowerShape)
653  fRecoUtils->RecalculateClusterShowerShapeParameters(fGeom, fCaloCells, clust);
654 
655  // DISTANCE TO BAD CHANNELS -----------------------------------
657  fRecoUtils->RecalculateClusterDistanceToBadChannel(fGeom, fCaloCells, clust);
658  }
659 
660  fCaloClusters->Compress();
661 }
662 
663 //________________________________________________________________________________________
665 {
666  // Select clusterization/unfolding algorithm and set all the needed parameters.
667 
668  if (fEvent->GetRunNumber()==fRun)
669  return;
670  fRun = fEvent->GetRunNumber();
671 
672  if (fJustUnfold){
673  // init the unfolding afterburner
674  delete fUnfolder;
675  fUnfolder = new AliEMCALAfterBurnerUF(fRecParam->GetW0(),fRecParam->GetLocMaxCut(),fRecParam->GetMinECut());
676  return;
677  }
678 
679  if (fGeomName.Length()>0)
680  fGeom = AliEMCALGeometry::GetInstance(fGeomName);
681  else
682  fGeom = AliEMCALGeometry::GetInstanceFromRunNumber(fRun);
683  if (!fGeom) {
684  AliFatal("Geometry not available!!!");
685  return;
686  }
687 
688  if (!fGeomMatrixSet) {
689  if (fLoadGeomMatrices) {
690  for(Int_t mod=0; mod < fGeom->GetNumberOfSuperModules(); ++mod) {
691  if (fGeomMatrix[mod]){
692  // AliDebug(3, fGeomMatrix[mod]); need to "print"
693  fGeom->SetMisalMatrix(fGeomMatrix[mod],mod);
694  }
695  }
696  } else { // get matrix from file (work around bug in aliroot)
697  for(Int_t mod=0; mod < fGeom->GetEMCGeometry()->GetNumberOfSuperModules(); ++mod) {
698  const TGeoHMatrix *gm = 0;
699  if (fEsd) {
700  gm = fEsd->GetEMCALMatrix(mod);
701  } else {
702  AliAODHeader *aodheader = dynamic_cast<AliAODHeader*>(fAod->GetHeader());
703  if(!aodheader) AliFatal("Not a standard AOD");
704  if (aodheader) {
705  gm = aodheader->GetEMCALMatrix(mod);
706  }
707  }
708  if (gm) {
709  //AliDebug(3, gm); need to "print"
710  fGeom->SetMisalMatrix(gm,mod);
711  }
712  }
713  }
714  fGeomMatrixSet=kTRUE;
715  }
716 
717  // setup digit array if needed
718  if (!fDigitsArr) {
719  fDigitsArr = new TClonesArray("AliEMCALDigit", 1000);
720  fDigitsArr->SetOwner(1);
721  }
722 
723  // then setup clusterizer
724  if (fClusterizer) {
725  // avoid to delete digits array
726  fClusterizer->SetDigitsArr(0);
727  delete fClusterizer;
728  }
729  if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
730  fClusterizer = new AliEMCALClusterizerv1(fGeom);
731  else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN) {
732  AliEMCALClusterizerNxN *clusterizer = new AliEMCALClusterizerNxN(fGeom);
733  clusterizer->SetNRowDiff(fRecParam->GetNRowDiff()); //MV: already done in AliEMCALClusterizer::InitParameters
734  clusterizer->SetNColDiff(fRecParam->GetNColDiff()); //MV: already done in AliEMCALClusterizer::InitParameters
735  fClusterizer = clusterizer;
736  }
737  else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
738  fClusterizer = new AliEMCALClusterizerv2(fGeom);
739  else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerFW) {
740  AliEMCALClusterizerFixedWindow *clusterizer = new AliEMCALClusterizerFixedWindow(fGeom);
741  clusterizer->SetNphi(fNPhi);
742  clusterizer->SetNeta(fNEta);
743  clusterizer->SetShiftPhi(fShiftPhi);
744  clusterizer->SetShiftEta(fShiftEta);
745  clusterizer->SetTRUshift(fTRUShift);
746  fClusterizer = clusterizer;
747  }
748  else {
749  AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
750  }
751  fClusterizer->InitParameters(fRecParam);
752 
753  if ((!fCalibData&&fLoadCalib) || (!fPedestalData&&fLoadPed)) {
754  AliCDBManager *cdb = AliCDBManager::Instance();
755  if (!cdb->IsDefaultStorageSet() && !fOCDBpath.IsNull())
756  cdb->SetDefaultStorage(fOCDBpath);
757  if (fRun!=cdb->GetRun())
758  cdb->SetRun(fRun);
759  }
760  if (!fCalibData&&fLoadCalib&&fRun>0) {
761  AliCDBEntry *entry = static_cast<AliCDBEntry*>(AliCDBManager::Instance()->Get("EMCAL/Calib/Data"));
762  if (entry)
763  fCalibData = static_cast<AliEMCALCalibData*>(entry->GetObject());
764  if (!fCalibData)
765  AliFatal("Calibration parameters not found in CDB!");
766  }
767  if (!fPedestalData&&fLoadPed&&fRun>0) {
768  AliCDBEntry *entry = static_cast<AliCDBEntry*>(AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals"));
769  if (entry)
770  fPedestalData = static_cast<AliCaloCalibPedestal*>(entry->GetObject());
771  }
772  if (fCalibData) {
773  fClusterizer->SetInputCalibrated(kFALSE);
774  fClusterizer->SetCalibrationParameters(fCalibData);
775  } else {
776  fClusterizer->SetInputCalibrated(kTRUE);
777  }
778  fClusterizer->SetCaloCalibPedestal(fPedestalData);
779  fClusterizer->SetJustClusters(kTRUE);
780  fClusterizer->SetDigitsArr(fDigitsArr);
781  fClusterizer->SetOutput(0);
782  fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
783 
784 }
void AddContainer(inputObjectType type)
AliEMCALGeometry * fGeom
!geometry object
double Double_t
Definition: External.C:58
AliEMCALRecParam * fRecParam
recpoints array
std::map< std::string, AliEMCALRecParam::AliEMCALClusterizerFlag > clusterizerTypeMap
AliEMCALAfterBurnerUF * fUnfolder
clusterizer
TCanvas * c
Definition: TestFitELoss.C:172
AliVCaloCells * fCaloCells
! pointer to calo cells
AliEMCALRecoUtils * fRecoUtils
! pointer to reco utils
AliClusterContainer * fClusCont
! pointer to the cluster container
TGeoHMatrix * fGeomMatrix[AliEMCALGeoParams::fgkEMCALModules]
AliESDEvent * fEsd
calo clusters array
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
void GetProperty(std::string propertyName, T &property, bool requiredProperty=true, std::string correctionName="")
Retrieve property.
short Short_t
Definition: External.C:23
energy
static RegisterCorrectionComponent< AliEmcalCorrectionClusterizer > reg
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
Int_t fOrgClusterCellId[fgkTotalCellNumber]
unsigned short UShort_t
Definition: External.C:28
bool Bool_t
Definition: External.C:53