PWG4/UserTasks/EmcalTasks/AliAnalysisTaskEMCALClusterizeFast.cxx

   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 /* $Id$ */
  17
  18 // --- Root ---
  19 #include <TClonesArray.h>
  20 #include <TGeoManager.h>
  21 #include <TObjArray.h>
  22 #include <TString.h>
  23 #include <TTree.h>
  24
  25 // --- AliRoot ---
  26 #include "AliAODCaloCluster.h"
  27 #include "AliAODEvent.h"
  28 #include "AliAnalysisManager.h"
  29 #include "AliCDBEntry.h"
  30 #include "AliCDBManager.h"
  31 #include "AliCaloCalibPedestal.h"
  32 #include "AliEMCALAfterBurnerUF.h"
  33 #include "AliEMCALCalibData.h"
  34 #include "AliEMCALClusterizerNxN.h"
  35 #include "AliEMCALClusterizerv1.h"
  36 #include "AliEMCALClusterizerv2.h"
  37 #include "AliEMCALClusterizerFixedWindow.h"
  38 #include "AliEMCALDigit.h"
  39 #include "AliEMCALGeometry.h"
  40 #include "AliEMCALRecParam.h"
  41 #include "AliEMCALRecPoint.h"
  42 #include "AliEMCALRecoUtils.h"
  43 #include "AliESDEvent.h"
  44 #include "AliInputEventHandler.h"
  45 #include "AliLog.h"
  46
  47 #include "AliAnalysisTaskEMCALClusterizeFast.h"
  48
  49 ClassImp(AliAnalysisTaskEMCALClusterizeFast)
  50
  51 //________________________________________________________________________
  52 AliAnalysisTaskEMCALClusterizeFast::AliAnalysisTaskEMCALClusterizeFast()
  53   : AliAnalysisTaskSE(),
  54     fRun(-1),
  55     fDigitsArr(0),
  56     fClusterArr(0),
  57     fRecParam(0),
  58     fClusterizer(0),
  59     fUnfolder(0),
  60     fJustUnfold(kFALSE),
  61     fGeomName(),
  62     fGeomMatrixSet(kFALSE),
  63     fLoadGeomMatrices(kFALSE),
  64     fOCDBpath(),
  65     fCalibData(0),
  66     fPedestalData(0),
  67     fOutputAODBranch(0),
  68     fOutputAODBrName(),
  69     fRecoUtils(0),
  70     fLoadCalib(0),
  71     fLoadPed(0),
  72     fAttachClusters(1),
  73     fRecalibOnly(0),
  74     fSubBackground(0),
  75     fCreatePattern(0),
  76     fOverwrite(0),
  77     fNewClusterArrayName("newCaloClusters"),
  78     fNPhi(4),
  79     fNEta(4),
  80     fShiftPhi(2),
  81     fShiftEta(2),
  82     fTRUShift(0),
  83     fClusterizeFastORs(0)
  84 {
  85   // Constructor
  86 }
  87
  88 //________________________________________________________________________
  89 AliAnalysisTaskEMCALClusterizeFast::AliAnalysisTaskEMCALClusterizeFast(const char *name)
  90   : AliAnalysisTaskSE(name),
  91     fRun(-1),
  92     fDigitsArr(0),
  93     fClusterArr(0),
  94     fRecParam(new AliEMCALRecParam),
  95     fClusterizer(0),
  96     fUnfolder(0),
  97     fJustUnfold(kFALSE),
  98     fGeomName("EMCAL_FIRSTYEARV1"),
  99     fGeomMatrixSet(kFALSE),
 100     fLoadGeomMatrices(kFALSE),
 101     fOCDBpath(),
 102     fCalibData(0),
 103     fPedestalData(0),
 104     fOutputAODBranch(0),
 105     fOutputAODBrName(),
 106     fRecoUtils(0),
 107     fLoadCalib(0),
 108     fLoadPed(0),
 109     fAttachClusters(1),
 110     fRecalibOnly(0),
 111     fSubBackground(0),
 112     fCreatePattern(0),
 113     fOverwrite(0),
 114     fNewClusterArrayName("newCaloClusters"),
 115     fNPhi(4),
 116     fNEta(4),
 117     fShiftPhi(2),
 118     fShiftEta(2),
 119     fTRUShift(0),
 120     fClusterizeFastORs(0)
 121 {
 122   // Constructor
 123
 124   fBranchNames     = "ESD:AliESDHeader.,AliESDRun.,EMCALCells.,EMCALTrigger. AOD:header,emcalCells";
 125   for(Int_t i = 0; i < 12; ++i)
 126     fGeomMatrix[i] = 0;
 127 }
 128
 129 //________________________________________________________________________
 130 AliAnalysisTaskEMCALClusterizeFast::~AliAnalysisTaskEMCALClusterizeFast()
 131 {
 132   // Destructor.
 133
 134   delete fDigitsArr;
 135   delete fClusterizer;
 136   delete fUnfolder;
 137   delete fRecoUtils;
 138 }
 139
 140 //-------------------------------------------------------------------
 141 void AliAnalysisTaskEMCALClusterizeFast::UserCreateOutputObjects()
 142 {
 143   // Create output objects.
 144
 145   if (!fOutputAODBrName.IsNull()) {
 146     fOutputAODBranch = new TClonesArray("AliAODCaloCluster", 0);
 147     fOutputAODBranch->SetName(fOutputAODBrName);
 148     AddAODBranch("TClonesArray", &fOutputAODBranch);
 149     AliInfo(Form("Created Branch: %s",fOutputAODBrName.Data()));
 150   }
 151 }
 152
 153 //________________________________________________________________________
 154 void AliAnalysisTaskEMCALClusterizeFast::UserExec(Option_t *)
 155 {
 156   // Main loop, called for each event
 157
 158   // remove the contents of output list set in the previous event
 159   if (fOutputAODBranch)
 160     fOutputAODBranch->Clear("C");
 161
 162   AliESDEvent *esdevent = dynamic_cast<AliESDEvent*>(InputEvent());
 163   AliAODEvent *aodevent = dynamic_cast<AliAODEvent*>(InputEvent());
 164
 165   if (!esdevent&&!aodevent) {
 166     Error("UserExec","Event not available");
 167     return;
 168   }
 169
 170   LoadBranches();
 171
 172   UInt_t offtrigger = 0;
 173   if (esdevent) {
 174     UInt_t mask1 = esdevent->GetESDRun()->GetDetectorsInDAQ();
 175     UInt_t mask2 = esdevent->GetESDRun()->GetDetectorsInReco();
 176     Bool_t desc1 = (mask1 >> 18) & 0x1;
 177     Bool_t desc2 = (mask2 >> 18) & 0x1;
 178     if (desc1==0 || desc2==0) { //AliDAQ::OfflineModuleName(180=="EMCAL"
 179       AliError(Form("EMCAL not in DAQ/RECO: %u (%u)/%u (%u)",
 180                     mask1, esdevent->GetESDRun()->GetDetectorsInReco(),
 181                     mask2, esdevent->GetESDRun()->GetDetectorsInDAQ()));
 182       return;
 183     }
 184     AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
 185     offtrigger = ((AliInputEventHandler*)(am->GetInputEventHandler()))->IsEventSelected();
 186   } else if (aodevent) {
 187     offtrigger =  aodevent->GetHeader()->GetOfflineTrigger();
 188   }
 189   if (offtrigger & AliVEvent::kFastOnly) {
 190     AliWarning(Form("EMCAL not in fast only partition"));
 191     return;
 192   }
 193
 194   Init();
 195
 196   if (fJustUnfold) {
 197     AliWarning("Unfolding not implemented");
 198     return;
 199   }
 200
 201   FillDigitsArray();
 202
 203   if (fRecalibOnly) {
 204     UpdateCells();
 205     return; // not requested to run clusterizer
 206   }
 207
 208   Clusterize();
 209   UpdateCells();
 210   UpdateClusters();
 211
 212   if (fOutputAODBranch)
 213     RecPoints2Clusters(fOutputAODBranch);
 214 }
 215
 216 //________________________________________________________________________
 217 void AliAnalysisTaskEMCALClusterizeFast::Clusterize()
 218 {
 219   // Clusterize
 220
 221   if (fSubBackground) {
 222     fClusterizer->SetInputCalibrated(kTRUE);
 223     fClusterizer->SetCalibrationParameters(0);
 224   }
 225
 226   fClusterizer->Digits2Clusters("");
 227   if (fSubBackground) {
 228     if (fCalibData) {
 229       fClusterizer->SetInputCalibrated(kFALSE);
 230       fClusterizer->SetCalibrationParameters(fCalibData);
 231     }
 232   }
 233 }
 234
 235 //________________________________________________________________________
 236 void AliAnalysisTaskEMCALClusterizeFast::FillDigitsArray()
 237 {
 238   // Fill digits array
 239
 240   fDigitsArr->Clear("C");
 241
 242   if (fCreatePattern) {  // Fill digits from a pattern
 243     AliEMCALGeometry *fGeom = AliEMCALGeometry::GetInstance(fGeomName);
 244     Int_t maxd = fGeom->GetNCells() / 4;
 245     for (Int_t idigit = 0; idigit < maxd; idigit++){
 246       if (idigit % 24 == 12) idigit += 12;
 247       AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
 248       digit->SetId(idigit * 4);
 249       digit->SetTime(600);
 250       digit->SetTimeR(600);
 251       digit->SetIndexInList(idigit);
 252       digit->SetType(AliEMCALDigit::kHG);
 253       digit->SetAmplitude(0.1);
 254     }
 255
 256   } else if (fClusterizeFastORs) { // Fill digits from FastORs
 257
 258     AliEMCALGeometry *fGeom = AliEMCALGeometry::GetInstance(fGeomName);
 259
 260     AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
 261     if (!esd){
 262       AliError("Cannot get the ESD event");
 263       return;
 264     }
 265
 266     AliESDCaloTrigger *triggers = esd->GetCaloTrigger("EMCAL");
 267
 268     if (!triggers || !(triggers->GetEntries() > 0))
 269       return;
 270
 271     Int_t idigit = 0;
 272     triggers->Reset();
 273
 274     while ((triggers->Next())) {
 275       Float_t triggerAmplitude = 0;
 276       triggers->GetAmplitude(triggerAmplitude);
 277       if (triggerAmplitude <= 0)
 278         continue;
 279
 280       Int_t triggerTime = 0;
 281       Int_t ntimes = 0;
 282       triggers->GetNL0Times(ntimes);
 283       if (ntimes > 0) {
 284         Int_t trgtimes[25];
 285         triggers->GetL0Times(trgtimes);
 286         triggerTime = trgtimes[0];
 287       }
 288
 289       Int_t triggerCol = 0, triggerRow = 0;
 290       triggers->GetPosition(triggerCol, triggerRow);
 291
 292       Int_t find = -1;
 293       fGeom->GetAbsFastORIndexFromPositionInEMCAL(triggerCol, triggerRow, find);
 294
 295       if (find<0)
 296         continue;
 297
 298       Int_t cidx[4] = {-1};
 299       Bool_t ret = fGeom->GetCellIndexFromFastORIndex(find, cidx);
 300
 301       if (!ret)
 302         continue;
 303
 304       for (Int_t idxpos = 0; idxpos < 4; idxpos++) {
 305         Int_t triggerNumber = cidx[idxpos];
 306         AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
 307         digit->SetId(triggerNumber);
 308         digit->SetTime(triggerTime);
 309         digit->SetTimeR(triggerTime);
 310         digit->SetIndexInList(idigit);
 311         digit->SetType(AliEMCALDigit::kHG);
 312         digit->SetAmplitude(triggerAmplitude);
 313         idigit++;
 314       }
 315     }
 316
 317   } else { // Fill digits from cells.
 318
 319     AliVCaloCells *cells = InputEvent()->GetEMCALCells();
 320     Double_t avgE = 0; // for background subtraction
 321     Int_t ncells = cells->GetNumberOfCells();
 322     for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell) {
 323       Double_t cellAmplitude=0, cellTime=0;
 324       Short_t cellNumber=0;
 325       if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime) != kTRUE)
 326         break;
 327       AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
 328       digit->SetId(cellNumber);
 329       digit->SetTime(cellTime);
 330       digit->SetTimeR(cellTime);
 331       digit->SetIndexInList(idigit);
 332       digit->SetType(AliEMCALDigit::kHG);
 333       if (fRecalibOnly||fSubBackground) {
 334         Float_t energy = cellAmplitude;
 335         Float_t time    = cellTime;
 336         fClusterizer->Calibrate(energy,time,cellNumber);
 337         digit->SetAmplitude(energy);
 338         avgE += energy;
 339       } else {
 340         digit->SetAmplitude(cellAmplitude);
 341       }
 342       idigit++;
 343     }
 344
 345     if (fSubBackground) {
 346       avgE /= AliEMCALGeometry::GetInstance(fGeomName)->GetNumberOfSuperModules()*48*24;
 347       Int_t ndigis = fDigitsArr->GetEntries();
 348       for (Int_t i = 0; i < ndigis; ++i) {
 349         AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(i));
 350         Double_t energy = digit->GetAmplitude() - avgE;
 351         if (energy<=0.001) {
 352           digit->SetAmplitude(0);
 353         } else {
 354           digit->SetAmplitude(energy);
 355         }
 356       }
 357     }
 358   }
 359 }
 360
 361 //________________________________________________________________________________________
 362 void AliAnalysisTaskEMCALClusterizeFast::RecPoints2Clusters(TClonesArray *clus)
 363 {
 364   // Cluster energy, global position, cells and their amplitude fractions are restored.
 365
 366   Bool_t esdobjects = 0;
 367   if (strcmp(clus->GetClass()->GetName(),"AliESDCaloCluster")==0)
 368     esdobjects = 1;
 369
 370   AliVCaloCells *cells = InputEvent()->GetEMCALCells();
 371   AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance(fGeomName);
 372
 373   AliDebug(1, Form("total no of clusters %d", fClusterArr->GetEntriesFast()));
 374
 375   Int_t Ncls = fClusterArr->GetEntriesFast();
 376   for(Int_t i=0, nout=clus->GetEntries(); i < Ncls; ++i) {
 377     AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
 378     Int_t ncells_true = 0;
 379     const Int_t ncells = recpoint->GetMultiplicity();
 380     UShort_t   absIds[ncells];
 381     Double32_t ratios[ncells];
 382     Int_t *dlist = recpoint->GetDigitsList();
 383     Float_t *elist = recpoint->GetEnergiesList();
 384     for (Int_t c = 0; c < ncells; ++c) {
 385       AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
 386       absIds[ncells_true] = digit->GetId();
 387       ratios[ncells_true] = elist[c]/digit->GetAmplitude();
 388       ++ncells_true;
 389     }
 390
 391     if (ncells_true < 1) {
 392       AliWarning("Skipping cluster with no cells");
 393       continue;
 394     }
 395
 396     // calculate new cluster position
 397     TVector3 gpos;
 398     recpoint->GetGlobalPosition(gpos);
 399     Float_t g[3];
 400     gpos.GetXYZ(g);
 401
 402     AliDebug(1, Form("energy %f", recpoint->GetEnergy()));
 403
 404     AliVCluster *c = static_cast<AliVCluster*>(clus->New(nout++));
 405     c->SetType(AliVCluster::kEMCALClusterv1);
 406     c->SetE(recpoint->GetEnergy());
 407     c->SetPosition(g);
 408     c->SetNCells(ncells_true);
 409     c->SetDispersion(recpoint->GetDispersion());
 410     c->SetEmcCpvDistance(-1);            //not yet implemented
 411     c->SetChi2(-1);                      //not yet implemented
 412     c->SetTOF(recpoint->GetTime()) ;     //time-of-flight
 413     c->SetNExMax(recpoint->GetNExMax()); //number of local maxima
 414     Float_t elipAxis[2];
 415     recpoint->GetElipsAxis(elipAxis);
 416     c->SetM02(elipAxis[0]*elipAxis[0]);
 417     c->SetM20(elipAxis[1]*elipAxis[1]);
 418     if (fPedestalData) {
 419       c->SetDistanceToBadChannel(recpoint->GetDistanceToBadTower());
 420     }
 421     else {
 422       if (fRecoUtils && fRecoUtils->IsBadChannelsRemovalSwitchedOn()) {
 423         fRecoUtils->RecalculateClusterDistanceToBadChannel(geom, cells, c);
 424       }
 425     }
 426
 427     if (esdobjects) {
 428       AliESDCaloCluster *cesd = static_cast<AliESDCaloCluster*>(c);
 429       cesd->SetCellsAbsId(absIds);
 430       cesd->SetCellsAmplitudeFraction(ratios);
 431       cesd->SetID(recpoint->GetUniqueID());
 432     } else {
 433       AliAODCaloCluster *caod = static_cast<AliAODCaloCluster*>(c);
 434       caod->SetCellsAbsId(absIds);
 435       caod->SetCellsAmplitudeFraction(ratios);
 436     }
 437   }
 438
 439   AliESDEvent *esdevent = dynamic_cast<AliESDEvent*>(InputEvent());
 440   if (!esdevent)
 441     return;
 442   if (!fRecoUtils)
 443     return;
 444
 445   AliAnalysisManager::GetAnalysisManager()->LoadBranch("Tracks");
 446   fRecoUtils->FindMatches(esdevent,clus);
 447
 448   if (!esdobjects) {
 449     Int_t Nclus = clus->GetEntries();
 450     for(Int_t i=0; i < Nclus; ++i) {
 451       AliAODCaloCluster *c = static_cast<AliAODCaloCluster*>(clus->At(i));
 452       Int_t trackIndex = fRecoUtils->GetMatchedTrackIndex(i);
 453       if (trackIndex >= 0) {
 454         Float_t dR, dZ;
 455         fRecoUtils->GetMatchedResiduals(i,dR, dZ);
 456         c->AddTrackMatched(0x0);   //esdevent->GetTrack(trackIndex));
 457         c->SetTrackDistance(dR,dZ); // not implemented
 458         c->SetEmcCpvDistance(dR);
 459         c->SetChi2(dZ);
 460         if (DebugLevel() > 1)
 461           AliInfo(Form("Matched Track index %d to new cluster %d\n",trackIndex,i));
 462       }
 463     }
 464   } else {
 465     Int_t Nclus = clus->GetEntries();
 466     for(Int_t i=0; i < Nclus; ++i) {
 467       AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->At(i));
 468       Int_t trackIndex = fRecoUtils->GetMatchedTrackIndex(i);
 469       if (trackIndex >= 0) {
 470         Float_t dR, dZ;
 471         fRecoUtils->GetMatchedResiduals(i,dR, dZ);
 472         c->SetTrackDistance(dR,dZ);
 473         c->SetEmcCpvDistance(dR); //to be consistent with AODs
 474         c->SetChi2(dZ);           //to be consistent with AODs
 475         TArrayI tm(1,&trackIndex);
 476         c->AddTracksMatched(tm);
 477         if (DebugLevel() > 1)
 478           AliInfo(Form("Matched Track index %d to new cluster %d\n",trackIndex,i));
 479       }
 480     }
 481   }
 482 }
 483
 484 //________________________________________________________________________
 485 void AliAnalysisTaskEMCALClusterizeFast::UpdateCells()
 486 {
 487   // Update cells in case re-calibration was done.
 488
 489   if (!fCalibData&&!fSubBackground)
 490     return;
 491
 492   AliVCaloCells *cells = InputEvent()->GetEMCALCells();
 493   Int_t ncells = cells->GetNumberOfCells();
 494   Int_t ndigis = fDigitsArr->GetEntries();
 495
 496   if (ncells!=ndigis) {
 497     cells->DeleteContainer();
 498     cells->CreateContainer(ndigis);
 499   }
 500   for (Int_t idigit = 0; idigit < ndigis; ++idigit) {
 501     AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(idigit));
 502     Double_t cellAmplitude = digit->GetCalibAmp();
 503     Short_t cellNumber = digit->GetId();
 504     Double_t cellTime = digit->GetTime();
 505     cells->SetCell(idigit, cellNumber, cellAmplitude, cellTime);
 506   }
 507 }
 508
 509 //________________________________________________________________________
 510 void AliAnalysisTaskEMCALClusterizeFast::UpdateClusters()
 511 {
 512   // Update cells in case re-calibration was done.
 513
 514   if (!fAttachClusters)
 515     return;
 516
 517   TClonesArray *clus = 0;
 518
 519   if (fOverwrite) {
 520     clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("caloClusters"));
 521     if (!clus)
 522       clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("CaloClusters"));
 523     if (!clus)
 524       return;
 525
 526     Int_t nents = clus->GetEntries();
 527     for (Int_t i=0;i<nents;++i) {
 528       AliVCluster *c = static_cast<AliVCluster*>(clus->At(i));
 529       if (!c)
 530         continue;
 531       if (c->IsEMCAL())
 532         delete clus->RemoveAt(i);
 533     }
 534   } else {
 535     clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fNewClusterArrayName));
 536     if (!clus) {
 537       clus = new TClonesArray("AliESDCaloCluster");
 538       clus->SetName(fNewClusterArrayName);
 539       InputEvent()->AddObject(clus);
 540     } else {
 541       clus->Delete();
 542     }
 543   }
 544
 545   RecPoints2Clusters(clus);
 546 }
 547
 548 //________________________________________________________________________________________
 549 void AliAnalysisTaskEMCALClusterizeFast::Init()
 550 {
 551   //Select clusterization/unfolding algorithm and set all the needed parameters
 552
 553   AliVEvent * event = InputEvent();
 554   if (!event) {
 555     AliWarning("Event not available!!!");
 556     return;
 557   }
 558
 559   if (event->GetRunNumber()==fRun)
 560     return;
 561   fRun = event->GetRunNumber();
 562
 563   if (fJustUnfold){
 564     // init the unfolding afterburner
 565     delete fUnfolder;
 566     fUnfolder = new AliEMCALAfterBurnerUF(fRecParam->GetW0(),fRecParam->GetLocMaxCut());
 567     return;
 568   }
 569
 570   AliEMCALGeometry *geometry = AliEMCALGeometry::GetInstance(fGeomName);
 571   if (!geometry) {
 572     AliFatal("Geometry not available!!!");
 573     return;
 574   }
 575
 576   if (!fGeomMatrixSet) {
 577     if (fLoadGeomMatrices) {
 578       for(Int_t mod=0; mod < geometry->GetNumberOfSuperModules(); ++mod) {
 579         if (fGeomMatrix[mod]){
 580           if (DebugLevel() > 2)
 581             fGeomMatrix[mod]->Print();
 582           geometry->SetMisalMatrix(fGeomMatrix[mod],mod);
 583         }
 584       }
 585     } else { // get matrix from file (work around bug in aliroot)
 586       for(Int_t mod=0; mod < geometry->GetEMCGeometry()->GetNumberOfSuperModules(); ++mod) {
 587         const TGeoHMatrix *gm = 0;
 588         AliESDEvent *esdevent = dynamic_cast<AliESDEvent*>(event);
 589         if (esdevent) {
 590           gm = esdevent->GetEMCALMatrix(mod);
 591         } else {
 592           AliAODHeader *aodheader = dynamic_cast<AliAODHeader*>(event->GetHeader());
 593           if (aodheader) {
 594             gm = aodheader->GetEMCALMatrix(mod);
 595           }
 596         }
 597         if (gm) {
 598           if (DebugLevel() > 2)
 599             gm->Print();
 600           geometry->SetMisalMatrix(gm,mod);
 601         }
 602       }
 603     }
 604     fGeomMatrixSet=kTRUE;
 605   }
 606
 607   // setup digit array if needed
 608   if (!fDigitsArr) {
 609     fDigitsArr = new TClonesArray("AliEMCALDigit", 1000);
 610     fDigitsArr->SetOwner(1);
 611   }
 612
 613   // then setup clusterizer
 614   delete fClusterizer;
 615   if     (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
 616     fClusterizer = new AliEMCALClusterizerv1(geometry);
 617   else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN) {
 618    AliEMCALClusterizerNxN *clusterizer = new AliEMCALClusterizerNxN(geometry);
 619    clusterizer->SetNRowDiff(fRecParam->GetNRowDiff());
 620    clusterizer->SetNColDiff(fRecParam->GetNColDiff());
 621     fClusterizer = clusterizer;
 622   } else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
 623     fClusterizer = new AliEMCALClusterizerv2(geometry);
 624   else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerFW){
 625     AliEMCALClusterizerFixedWindow *clusterizer = new AliEMCALClusterizerFixedWindow(geometry);
 626     clusterizer->SetNphi(fNPhi);
 627     clusterizer->SetNeta(fNEta);
 628     clusterizer->SetShiftPhi(fShiftPhi);
 629     clusterizer->SetShiftEta(fShiftEta);
 630     clusterizer->SetTRUshift(fTRUShift);
 631     fClusterizer = clusterizer;
 632   }
 633   else{
 634     AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
 635   }
 636   fClusterizer->InitParameters(fRecParam);
 637
 638   if ((!fCalibData&&fLoadCalib) || (!fPedestalData&&fLoadPed)) {
 639     AliCDBManager *cdb = AliCDBManager::Instance();
 640     if (!cdb->IsDefaultStorageSet() && !fOCDBpath.IsNull())
 641       cdb->SetDefaultStorage(fOCDBpath);
 642     if (fRun!=cdb->GetRun())
 643       cdb->SetRun(fRun);
 644   }
 645   if (!fCalibData&&fLoadCalib&&fRun>0) {
 646     AliCDBEntry *entry = static_cast<AliCDBEntry*>(AliCDBManager::Instance()->Get("EMCAL/Calib/Data"));
 647     if (entry)
 648       fCalibData =  static_cast<AliEMCALCalibData*>(entry->GetObject());
 649     if (!fCalibData)
 650       AliFatal("Calibration parameters not found in CDB!");
 651   }
 652   if (!fPedestalData&&fLoadPed&&fRun>0) {
 653     AliCDBEntry *entry = static_cast<AliCDBEntry*>(AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals"));
 654     if (entry)
 655       fPedestalData =  static_cast<AliCaloCalibPedestal*>(entry->GetObject());
 656   }
 657   if (fCalibData) {
 658     fClusterizer->SetInputCalibrated(kFALSE);
 659     fClusterizer->SetCalibrationParameters(fCalibData);
 660   } else {
 661     fClusterizer->SetInputCalibrated(kTRUE);
 662   }
 663   fClusterizer->SetCaloCalibPedestal(fPedestalData);
 664   fClusterizer->SetJustClusters(kTRUE);
 665   fClusterizer->SetDigitsArr(fDigitsArr);
 666   fClusterizer->SetOutput(0);
 667   fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
 668 }