ANALYSIS/TenderSupplies/AliEMCALTenderSupply.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
  17 ///////////////////////////////////////////////////////////////////////////////
  18 //                                                                           //
  19 //  EMCAL tender, apply corrections to EMCAL clusters                        //
  20 //  and do track matching.                                                   //
  21 //  Author: Deepa Thomas (Utrecht University)                                //
  22 //                                                                           //
  23 ///////////////////////////////////////////////////////////////////////////////
  24
  25 #include <TROOT.h>
  26 #include <TFile.h>
  27 #include <TTree.h>
  28 #include <TInterpreter.h>
  29 #include <TObjArray.h>
  30 #include <TClonesArray.h>
  31 #include <TGeoGlobalMagField.h>
  32 #include <AliLog.h>
  33 #include <AliESDEvent.h>
  34 #include <AliAnalysisManager.h>
  35 #include <AliTender.h>
  36 #include "AliOADBContainer.h"
  37 #include "AliCDBManager.h"
  38 #include "AliCDBStorage.h"
  39 #include "AliCDBEntry.h"
  40 #include "AliMagF.h"
  41 #include "AliESDCaloCluster.h"
  42 #include "AliEMCALTenderSupply.h"
  43 #include "AliEMCALGeometry.h"
  44 #include "AliEMCALRecoUtils.h"
  45 #include "AliEMCALClusterizer.h"
  46 #include "AliEMCALRecParam.h"
  47 #include "AliEMCALRecPoint.h"
  48 #include "AliEMCALAfterBurnerUF.h"
  49 #include "AliEMCALClusterizerNxN.h"
  50 #include "AliEMCALClusterizerv1.h"
  51 #include "AliEMCALClusterizerv2.h"
  52 #include "AliEMCALDigit.h"
  53 #include "AliEMCALRecParam.h"
  54
  55 ClassImp(AliEMCALTenderSupply)
  56
  57 AliEMCALTenderSupply::AliEMCALTenderSupply() :
  58 AliTenderSupply()
  59 ,fEMCALGeo(0x0)
  60 ,fEMCALGeoName("EMCAL_COMPLETEV1")
  61 ,fEMCALRecoUtils(0)
  62 ,fConfigName("")
  63 ,fDebugLevel(0)
  64 ,fNonLinearFunc(-1)
  65 ,fNonLinearThreshold(-1)
  66 ,fReCalibCluster(kFALSE)
  67 ,fUpdateCell(kFALSE)
  68 ,fCalibrateEnergy(kFALSE)
  69 ,fCalibrateTime(kFALSE)
  70 ,fDoNonLinearity(kFALSE)
  71 ,fBadCellRemove(kFALSE)
  72 ,fRejectExoticCells(kFALSE)
  73 ,fRejectExoticClusters(kFALSE)
  74 ,fClusterBadChannelCheck(kFALSE)
  75 ,fRecalClusPos(kFALSE)
  76 ,fFiducial(kFALSE)
  77 ,fNCellsFromEMCALBorder(-1)
  78 ,fRecalDistToBadChannels(kFALSE)
  79 ,fRecalShowerShape(kFALSE)
  80 ,fInputTree(0)
  81 ,fInputFile(0)
  82 ,fFilepass(0)
  83 ,fMass(-1)
  84 ,fStep(-1)
  85 ,fCutEtaPhiSum(kTRUE)
  86 ,fCutEtaPhiSeparate(kFALSE)
  87 ,fRcut(-1)
  88 ,fEtacut(-1)
  89 ,fPhicut(-1)
  90 ,fBasePath("")
  91 ,fReClusterize(kFALSE)
  92 ,fClusterizer(0)
  93 ,fGeomMatrixSet(kFALSE)
  94 ,fLoadGeomMatrices(kFALSE)
  95 ,fRecParam(0x0)
  96 ,fDoTrackMatch(kFALSE)
  97 ,fDoUpdateOnly(kFALSE)
  98 ,fUnfolder(0)
  99 ,fDigitsArr(0)
 100 ,fClusterArr(0)
 101 ,fMisalignSurvey(kdefault)
 102 ,fExoticCellFraction(-1)
 103 ,fExoticCellDiffTime(-1)
 104 ,fExoticCellMinAmplitude(-1)
 105 {
 106   // Default constructor.
 107   for(Int_t i = 0; i < 10; i++) fEMCALMatrix[i] = 0 ;
 108   fEMCALRecoUtils  = new AliEMCALRecoUtils;
 109 }
 110
 111 //_____________________________________________________
 112 AliEMCALTenderSupply::AliEMCALTenderSupply(const char *name, const AliTender *tender) :
 113 AliTenderSupply(name,tender)
 114 ,fEMCALGeo(0x0)
 115 ,fEMCALGeoName("EMCAL_COMPLETEV1")
 116 ,fEMCALRecoUtils(0)
 117 ,fConfigName("")
 118 ,fDebugLevel(0)
 119 ,fNonLinearFunc(-1)
 120 ,fNonLinearThreshold(-1)
 121 ,fReCalibCluster(kFALSE)
 122 ,fUpdateCell(kFALSE)
 123 ,fCalibrateEnergy(kFALSE)
 124 ,fCalibrateTime(kFALSE)
 125 ,fDoNonLinearity(kFALSE)
 126 ,fBadCellRemove(kFALSE)
 127 ,fRejectExoticCells(kFALSE)
 128 ,fRejectExoticClusters(kFALSE)
 129 ,fClusterBadChannelCheck(kFALSE)
 130 ,fRecalClusPos(kFALSE)
 131 ,fFiducial(kFALSE)
 132 ,fNCellsFromEMCALBorder(-1)
 133 ,fRecalDistToBadChannels(kFALSE)
 134 ,fRecalShowerShape(kFALSE)
 135 ,fInputTree(0)
 136 ,fInputFile(0)
 137 ,fFilepass(0)
 138 ,fMass(-1)
 139 ,fStep(-1)
 140 ,fCutEtaPhiSum(kTRUE)
 141 ,fCutEtaPhiSeparate(kFALSE)
 142 ,fRcut(-1)
 143 ,fEtacut(-1)
 144 ,fPhicut(-1)
 145 ,fBasePath("")
 146 ,fReClusterize(kFALSE)
 147 ,fClusterizer(0)
 148 ,fGeomMatrixSet(kFALSE)
 149 ,fLoadGeomMatrices(kFALSE)
 150 ,fRecParam(0x0)
 151 ,fDoTrackMatch(kFALSE)
 152 ,fDoUpdateOnly(kFALSE)
 153 ,fUnfolder(0)
 154 ,fDigitsArr(0)
 155 ,fClusterArr(0)
 156 ,fMisalignSurvey(kdefault)
 157 ,fExoticCellFraction(-1)
 158 ,fExoticCellDiffTime(-1)
 159 ,fExoticCellMinAmplitude(-1)
 160 {
 161   // Named constructor
 162
 163   for(Int_t i = 0; i < 10; i++) fEMCALMatrix[i] = 0 ;
 164   fEMCALRecoUtils  = new AliEMCALRecoUtils;
 165 }
 166
 167 //_____________________________________________________
 168 AliEMCALTenderSupply::~AliEMCALTenderSupply()
 169 {
 170   //Destructor
 171
 172   delete fEMCALRecoUtils;
 173   delete fRecParam;
 174   delete fUnfolder;
 175   if (!fClusterizer) {
 176     fDigitsArr->Clear("C");
 177     delete fDigitsArr;
 178   } else {
 179     delete fClusterizer;
 180     fDigitsArr = 0;
 181   }
 182 }
 183
 184 //_____________________________________________________
 185 void AliEMCALTenderSupply::Init()
 186 {
 187   // Initialise EMCAL tender.
 188
 189   if (fDebugLevel>0)
 190     AliWarning("Init EMCAL Tender supply");
 191
 192   if (fConfigName.Length()>0 && gROOT->LoadMacro(fConfigName) >=0) {
 193     AliDebug(1, Form("Loading settings from macro %s", fConfigName.Data()));
 194     AliEMCALTenderSupply *tender = (AliEMCALTenderSupply*)gInterpreter->ProcessLine("ConfigEMCALTenderSupply()");
 195     fDebugLevel             = tender->fDebugLevel;
 196     fEMCALGeoName           = tender->fEMCALGeoName;
 197     fEMCALRecoUtils         = tender->fEMCALRecoUtils;
 198     fConfigName             = tender->fConfigName;
 199     fNonLinearFunc          = tender->fNonLinearFunc;
 200     fNonLinearThreshold     = tender->fNonLinearThreshold;
 201     fReCalibCluster         = tender->fReCalibCluster;
 202     fUpdateCell             = tender->fUpdateCell;
 203     fRecalClusPos           = tender->fRecalClusPos;
 204     fCalibrateEnergy        = tender->fCalibrateEnergy;
 205     fCalibrateTime          = tender->fCalibrateTime;
 206     fFiducial               = tender->fFiducial;
 207     fNCellsFromEMCALBorder  = tender->fNCellsFromEMCALBorder;
 208     fRecalDistToBadChannels = tender->fRecalDistToBadChannels;
 209     fRecalShowerShape       = tender->fRecalShowerShape;
 210     fClusterBadChannelCheck = tender->fClusterBadChannelCheck;
 211     fBadCellRemove          = tender->fBadCellRemove;
 212     fRejectExoticCells      = tender->fRejectExoticCells;
 213     fRejectExoticClusters   = tender->fRejectExoticClusters;
 214     fMass                   = tender->fMass;
 215     fStep                   = tender->fStep;
 216     fCutEtaPhiSum           = tender->fCutEtaPhiSum;
 217     fCutEtaPhiSeparate      = tender->fCutEtaPhiSeparate;
 218     fRcut                   = tender->fRcut;
 219     fEtacut                 = tender->fEtacut;
 220     fPhicut                 = tender->fPhicut;
 221     fReClusterize           = tender->fReClusterize;
 222     fLoadGeomMatrices       = tender->fLoadGeomMatrices;
 223     fRecParam               = tender->fRecParam;
 224     fDoNonLinearity         = tender->fDoNonLinearity;
 225     fDoTrackMatch           = tender->fDoTrackMatch;
 226     fDoUpdateOnly           = tender->fDoUpdateOnly;
 227     fMisalignSurvey         = tender->fMisalignSurvey;
 228     fExoticCellFraction     = tender->fExoticCellFraction;
 229     fExoticCellDiffTime     = tender->fExoticCellDiffTime;
 230     fExoticCellMinAmplitude = tender->fExoticCellMinAmplitude;
 231
 232     for(Int_t i = 0; i < 10; i++)
 233       fEMCALMatrix[i] = tender->fEMCALMatrix[i] ;
 234   }
 235
 236   if (fDebugLevel>0){
 237     AliInfo( "Emcal Tender settings: ======================================" );
 238     AliInfo( "------------ Switches --------------------------" );
 239     AliInfo( Form( "BadCellRemove : %d", fBadCellRemove ));
 240     AliInfo( Form( "ExoticCellRemove : %d", fRejectExoticCells ));
 241     AliInfo( Form( "CalibrateEnergy : %d", fCalibrateEnergy ));
 242     AliInfo( Form( "CalibrateTime : %d", fCalibrateTime ));
 243     AliInfo( Form( "UpdateCell : %d", fUpdateCell ));
 244     AliInfo( Form( "DoUpdateOnly : %d", fDoUpdateOnly ));
 245     AliInfo( Form( "Reclustering : %d", fReClusterize ));
 246     AliInfo( Form( "ClusterBadChannelCheck : %d", fClusterBadChannelCheck ));
 247     AliInfo( Form( "ClusterExoticChannelCheck : %d", fRejectExoticClusters ));
 248     AliInfo( Form( "CellFiducialRegion : %d", fFiducial ));
 249     AliInfo( Form( "ReCalibrateCluster : %d", fReCalibCluster ));
 250     AliInfo( Form( "RecalculateClusPos : %d", fRecalClusPos ));
 251     AliInfo( Form( "RecalShowerShape : %d", fRecalShowerShape ));
 252     AliInfo( Form( "NonLinearityCorrection : %d", fDoNonLinearity ));
 253     AliInfo( Form( "RecalDistBadChannel : %d", fRecalDistToBadChannels ));
 254     AliInfo( Form( "TrackMatch : %d", fDoTrackMatch ));
 255     AliInfo( "------------ Variables -------------------------" );
 256     AliInfo( Form( "DebugLevel : %d", fDebugLevel ));
 257     AliInfo( Form( "BasePath : %s", fBasePath.Data() ));
 258     AliInfo( Form( "ConfigFileName : %s", fConfigName.Data() ));
 259     AliInfo( Form( "EMCALGeometryName : %s", fEMCALGeoName.Data() ));
 260     AliInfo( Form( "NonLinearityFunction : %d", fNonLinearFunc ));
 261     AliInfo( Form( "NonLinearityThreshold : %d", fNonLinearThreshold ));
 262     AliInfo( Form( "MisalignmentMatrixSurvey : %d", fMisalignSurvey ));
 263     AliInfo( Form( "NumberOfCellsFromEMCALBorder : %d", fNCellsFromEMCALBorder ));
 264     AliInfo( Form( "RCut : %f", fRcut ));
 265     AliInfo( Form( "Mass : %f", fMass ));
 266     AliInfo( Form( "Step : %f", fStep ));
 267     AliInfo( Form( "EtaCut : %f", fEtacut ));
 268     AliInfo( Form( "PhiCut : %f", fPhicut ));
 269     AliInfo( Form( "ExoticCellFraction : %f", fExoticCellFraction ));
 270     AliInfo( Form( "ExoticCellDiffTime : %f", fExoticCellDiffTime ));
 271     AliInfo( Form( "ExoticCellMinAmplitude : %f", fExoticCellMinAmplitude ));
 272     AliInfo( "=============================================================" );
 273   }
 274
 275   // Init geometry
 276   fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
 277
 278   // digits array
 279   fDigitsArr       = new TClonesArray("AliEMCALDigit",1000);
 280
 281   // Initialising non-linearity parameters
 282   if( fNonLinearThreshold != -1 )
 283     fEMCALRecoUtils->SetNonLinearityThreshold(fNonLinearThreshold);
 284   if( fNonLinearFunc != -1 )
 285     fEMCALRecoUtils->SetNonLinearityFunction(fNonLinearFunc);
 286
 287   // missalignment function
 288   fEMCALRecoUtils->SetPositionAlgorithm(AliEMCALRecoUtils::kPosTowerGlobal);
 289
 290   // fiducial cut
 291   // do not do the eta0 fiducial cut
 292   if( fNCellsFromEMCALBorder != -1 )
 293     fEMCALRecoUtils->SetNumberOfCellsFromEMCALBorder(fNCellsFromEMCALBorder);
 294   fEMCALRecoUtils->SwitchOnNoFiducialBorderInEMCALEta0();
 295
 296   // exotic cell rejection
 297   if( fExoticCellFraction != -1 )
 298     fEMCALRecoUtils->SetExoticCellFractionCut( fExoticCellFraction );
 299   if( fExoticCellDiffTime != -1 )
 300     fEMCALRecoUtils->SetExoticCellDiffTimeCut( fExoticCellDiffTime );
 301   if( fExoticCellMinAmplitude != -1 )
 302     fEMCALRecoUtils->SetExoticCellMinAmplitudeCut( fExoticCellMinAmplitude );
 303
 304   // Setting track matching parameters ... mass, step size and residual cut
 305   if( fMass != -1 )
 306     fEMCALRecoUtils->SetMass(fMass);
 307   if( fStep != -1 )
 308     fEMCALRecoUtils->SetStep(fStep);
 309
 310   // spatial cut based on separate eta/phi or common processing
 311   if(fCutEtaPhiSum){
 312     fEMCALRecoUtils->SwitchOnCutEtaPhiSum();
 313     if( fRcut != -1 )
 314       fEMCALRecoUtils->SetCutR(fRcut);
 315   } else if (fCutEtaPhiSeparate) {
 316     fEMCALRecoUtils->SwitchOnCutEtaPhiSeparate();
 317     if( fEtacut != -1 )
 318       fEMCALRecoUtils->SetCutEta(fEtacut);
 319     if( fPhicut != -1 )
 320       fEMCALRecoUtils->SetCutPhi(fPhicut);
 321   }
 322 }
 323
 324 //_____________________________________________________
 325 void AliEMCALTenderSupply::ProcessEvent()
 326 {
 327   // Event loop.
 328
 329   AliESDEvent *event = fTender->GetEvent();
 330   if (!event) {
 331     AliError("ESD event ptr = 0, returning");
 332     return;
 333   }
 334
 335   // Initialising parameters once per run number
 336   if (fTender->RunChanged()){
 337
 338     // initiate default reco params if not supplied by user
 339     if (!fRecParam)
 340       InitRecParam();
 341
 342     // get pass
 343     GetPass();
 344
 345     // define what recalib parameters are needed for various switches
 346     // this is based on implementation in AliEMCALRecoUtils
 347     Bool_t needBadChannels = fBadCellRemove   | fClusterBadChannelCheck | fRecalDistToBadChannels | fReClusterize;
 348     Bool_t needRecalib     = fCalibrateEnergy | fReClusterize;
 349     Bool_t needTimecalib   = fCalibrateTime   | fReClusterize;
 350     Bool_t needMisalign    = fRecalClusPos    | fReClusterize;
 351     Bool_t needClusterizer = fReClusterize;
 352
 353     // Init bad channels
 354     if( needBadChannels )
 355     {
 356       Int_t fInitBC = InitBadChannels();
 357       if (fInitBC==0)
 358       {
 359         AliError("InitBadChannels returned false, returning");
 360         return;
 361       }
 362       if (fInitBC==1)
 363       {
 364         AliWarning("InitBadChannels OK");
 365       }
 366       if (fInitBC>1)
 367       {
 368         AliWarning(Form("No external hot channel set: %d - %s", event->GetRunNumber(), fFilepass.Data()));
 369       }
 370     }
 371
 372     // init recalibration factors
 373     if( needRecalib )
 374     {
 375       Int_t fInitRecalib = InitRecalib();
 376       if (fInitRecalib==0)
 377       {
 378         AliError("InitRecalib returned false, returning");
 379         return;
 380       }
 381       if (fInitRecalib==1)
 382       {
 383         AliWarning("InitRecalib OK");
 384       }
 385       if (fInitRecalib>1)
 386       {
 387         AliWarning(Form("No recalibration available: %d - %s", event->GetRunNumber(), fFilepass.Data()));
 388         fReCalibCluster = kFALSE;
 389       }
 390     }
 391
 392     // init time calibration
 393     if( needTimecalib )
 394     {
 395       Int_t initTC = InitTimeCalibration();
 396       if ( !initTC )
 397       {
 398         AliError("InitTimeCalibration returned false, returning");
 399         return;
 400       }
 401       if (initTC==1)
 402       {
 403         AliWarning("InitTimeCalib OK");
 404       }
 405       if( initTC > 1 )
 406       {
 407         AliWarning(Form("No external time calibration set: %d - %s", event->GetRunNumber(), fFilepass.Data()));
 408       }
 409     }
 410
 411     // init misalignment matrix
 412     if( needMisalign )
 413     {
 414       if (!InitMisalignMatrix()) {
 415
 416         AliError("InitMisalignmentMatrix returned false, returning");
 417         return;
 418       }
 419       else
 420         AliWarning("InitMisalignMatrix OK");
 421     }
 422
 423     // init clusterizer
 424     if( needClusterizer )
 425     {
 426       if (!InitClusterization())
 427       {
 428         AliError("InitClusterization returned false, returning");
 429         return;
 430       }
 431       else
 432         AliWarning("InitClusterization OK");
 433     }
 434
 435     if(fDebugLevel>1)
 436       fEMCALRecoUtils->Print("");
 437   }
 438
 439   // disable implied switches -------------------------------------------------
 440   // AliEMCALRecoUtils or clusterizer functions alredy include some
 441   // recalibration so based on those implied callibration te switches
 442   // here are disabled to avoid duplication
 443
 444   // clusterizer does cluster energy recalibration, position recomputation
 445   // and shower shape
 446   if( fReClusterize )
 447   {
 448     fReCalibCluster   = kFALSE;
 449     fRecalClusPos     = kFALSE;
 450     fRecalShowerShape = kFALSE;
 451   }
 452
 453   // CONFIGURE THE RECO UTILS -------------------------------------------------
 454   // configure the reco utils
 455   // this option does energy recalibration
 456   if( fCalibrateEnergy )
 457     fEMCALRecoUtils->SwitchOnRecalibration();
 458   else
 459     fEMCALRecoUtils->SwitchOffRecalibration();
 460
 461   // allows time calibration
 462   if( fCalibrateTime )
 463     fEMCALRecoUtils->SwitchOnTimeRecalibration();
 464   else
 465     fEMCALRecoUtils->SwitchOffTimeRecalibration();
 466
 467   // allows to zero bad cells
 468   if( fBadCellRemove )
 469     fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
 470   else
 471     fEMCALRecoUtils->SwitchOffBadChannelsRemoval();
 472
 473   // distance to bad channel recalibration
 474   if( fRecalDistToBadChannels )
 475     fEMCALRecoUtils->SwitchOnDistToBadChannelRecalculation();
 476   else
 477     fEMCALRecoUtils->SwitchOffDistToBadChannelRecalculation();
 478
 479   // exclude exotic cells
 480   if( fRejectExoticCells )
 481     fEMCALRecoUtils->SwitchOnRejectExoticCell();
 482   else
 483     fEMCALRecoUtils->SwitchOffRejectExoticCell();
 484
 485   // exclude clusters with exotic cells
 486   if( fRejectExoticClusters )
 487     fEMCALRecoUtils->SwitchOnRejectExoticCluster();
 488   else
 489     fEMCALRecoUtils->SwitchOffRejectExoticCluster();
 490
 491   // TODO: not implemented switches
 492   // SwitchOnClusterEnergySmearing
 493   // SwitchOnRunDepCorrection
 494
 495   // START PROCESSING ---------------------------------------------------------
 496   // Test if cells present
 497   AliESDCaloCells *cells= event->GetEMCALCells();
 498   if (cells->GetNumberOfCells()<=0)
 499   {
 500     if(fDebugLevel>1)
 501       AliWarning(Form("Number of EMCAL cells = %d, returning", cells->GetNumberOfCells()));
 502     return;
 503   }
 504
 505   if (fDebugLevel>2)
 506     AliInfo(Form("Re-calibrate cluster %d\n",fReCalibCluster));
 507
 508   // mark the cells not recalibrated in case of selected
 509   // time, energy recalibration or bad channel removal
 510   if( fCalibrateEnergy || fCalibrateTime || fBadCellRemove )
 511     fEMCALRecoUtils->ResetCellsCalibrated();
 512
 513   // CELL RECALIBRATION -------------------------------------------------------
 514   // cell objects will be updated
 515   // the cell calibrations are also processed locally any time those are needed
 516   // in case that the cell objects are not to be updated here for later use
 517   if( fUpdateCell )
 518   {
 519     // do the update
 520     // includes exotic cell check in the UpdateCells function - is not provided
 521     // by the reco utils
 522     UpdateCells();
 523
 524     // switch off recalibrations so those are not done multiple times
 525     // this is just for safety, the recalibrated flag of cell object
 526     // should not allow for farther processing anyways
 527     fEMCALRecoUtils->SwitchOffRecalibration();
 528     fEMCALRecoUtils->SwitchOffTimeRecalibration();
 529
 530     if (fDoUpdateOnly)
 531       return;
 532   }
 533
 534   // RECLUSTERIZATION ---------------------------------------------------------
 535   if (fReClusterize)
 536   {
 537     FillDigitsArray();
 538     Clusterize();
 539     UpdateClusters();
 540   }
 541
 542   // Store good clusters
 543   TClonesArray *clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
 544   if (!clusArr)
 545     clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
 546   if (!clusArr) {
 547     AliWarning(Form("No cluster array, number of cells in event = %d, returning", cells->GetNumberOfCells()));
 548     return;
 549   }
 550
 551   // PROCESS SINGLE CLUSTER RECALIBRATION -------------------------------------
 552   // now go through clusters one by one and process separate correction
 553   // as those were defined or not
 554   Int_t nclusters = clusArr->GetEntriesFast();
 555   for (Int_t icluster=0; icluster < nclusters; ++icluster)
 556   {
 557     AliVCluster *clust = static_cast<AliVCluster*>(clusArr->At(icluster));
 558     if (!clust)
 559       continue;
 560     if  (!clust->IsEMCAL())
 561       continue;
 562
 563     // REMOVE CLUSTERS WITH BAD CELLS -----------------------------
 564     if( fClusterBadChannelCheck )
 565     {
 566       // careful, the the ClusterContainsBadChannel is dependent on
 567       // SwitchOnBadChannelsRemoval, switching it ON automatically
 568       // and returning to original value after processing
 569       Bool_t badRemoval = fEMCALRecoUtils->IsBadChannelsRemovalSwitchedOn();
 570       fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
 571
 572       Bool_t badResult = fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo, clust->GetCellsAbsId(), clust->GetNCells());
 573
 574       // switch the bad channels removal back
 575       if( ! badRemoval )
 576         fEMCALRecoUtils->SwitchOffBadChannelsRemoval();
 577
 578       if( badResult )
 579       {
 580         delete clusArr->RemoveAt(icluster);
 581         continue; //TODO is it really needed to remove it? Or should we flag it?
 582       }
 583     }
 584
 585     // REMOVE EXOTIC CLUSTERS -------------------------------------
 586     // does process local cell recalibration energy and time without replacing
 587     // the global cell values, in case of no cell recalib done yet
 588     if( fRejectExoticClusters )
 589     {
 590       // careful, the IsExoticCluster is dependent on
 591       // SwitchOnRejectExoticCell, switching it ON automatically
 592       // and returning to original value after processing
 593       Bool_t exRemoval = fEMCALRecoUtils->IsRejectExoticCell();
 594       fEMCALRecoUtils->SwitchOnRejectExoticCell();
 595
 596       // get bunch crossing
 597       Int_t bunchCrossNo = fTender->GetEvent()->GetBunchCrossNumber();
 598
 599       Bool_t exResult = fEMCALRecoUtils->IsExoticCluster(clust, cells, bunchCrossNo );
 600
 601       // switch the exotic channels removal back
 602       if( ! exRemoval )
 603         fEMCALRecoUtils->SwitchOffRejectExoticCell();
 604
 605       if( exResult )
 606       {
 607         delete clusArr->RemoveAt(icluster);
 608         continue; //TODO is it really needed to remove it? Or should we flag it?
 609       }
 610     }
 611
 612     // FIDUCIAL CUT -----------------------------------------------
 613     if (fFiducial)
 614     {
 615       // depends on SetNumberOfCellsFromEMCALBorder
 616       // SwitchOnNoFiducialBorderInEMCALEta0
 617       if (!fEMCALRecoUtils->CheckCellFiducialRegion(fEMCALGeo, clust, cells)){
 618         delete clusArr->RemoveAt(icluster);
 619         continue; //TODO it would be nice to store the distance
 620       }
 621     }
 622
 623     // CLUSTER ENERGY ---------------------------------------------
 624     // does process local cell recalibration energy and time without replacing
 625     // the global cell values, in case of no cell recalib done yet
 626     if( fReCalibCluster )
 627       fEMCALRecoUtils->RecalibrateClusterEnergy(fEMCALGeo, clust, cells);
 628
 629     // CLUSTER POSITION -------------------------------------------
 630     // does process local cell energy recalibration, if enabled and cells
 631     // not calibratied yet
 632     if( fRecalClusPos )
 633       fEMCALRecoUtils->RecalculateClusterPosition(fEMCALGeo, cells, clust);
 634
 635     // SHOWER SHAPE -----------------------------------------------
 636     if( fRecalShowerShape )
 637       fEMCALRecoUtils->RecalculateClusterShowerShapeParameters(fEMCALGeo, cells, clust);
 638
 639     // NONLINEARITY -----------------------------------------------
 640     if( fDoNonLinearity )
 641     {
 642       Float_t correctedEnergy = fEMCALRecoUtils->CorrectClusterEnergyLinearity(clust);
 643       clust->SetE(correctedEnergy);
 644     }
 645
 646     // DISTANCE TO BAD CHANNELS -----------------------------------
 647     if( fRecalDistToBadChannels )
 648       fEMCALRecoUtils->RecalculateClusterDistanceToBadChannel(fEMCALGeo, cells, clust);
 649   }
 650
 651   clusArr->Compress();
 652
 653   if (!fDoTrackMatch)
 654     return;
 655
 656   // TRACK MATCHING -----------------------------------------------------------
 657   if (!TGeoGlobalMagField::Instance()->GetField())
 658   {
 659     event->InitMagneticField();
 660   }
 661
 662   fEMCALRecoUtils->FindMatches(event,0x0,fEMCALGeo);
 663   fEMCALRecoUtils->SetClusterMatchedToTrack(event);
 664   fEMCALRecoUtils->SetTracksMatchedToCluster(event);
 665 }
 666
 667 //_____________________________________________________
 668 Bool_t AliEMCALTenderSupply::InitMisalignMatrix()
 669 {
 670   // Initialising misalignment matrices
 671
 672   AliESDEvent *event = fTender->GetEvent();
 673   if (!event)
 674     return kFALSE;
 675
 676   if (fGeomMatrixSet)
 677   {
 678     AliInfo("Misalignment matrix already set");
 679     return kTRUE;
 680   }
 681
 682   if (fDebugLevel>0)
 683     AliInfo("Initialising misalignment matrix");
 684
 685   if (fLoadGeomMatrices) {
 686     for(Int_t mod=0; mod < fEMCALGeo->GetNumberOfSuperModules(); ++mod)
 687     {
 688       if (fEMCALMatrix[mod]){
 689         if(fDebugLevel > 2)
 690           fEMCALMatrix[mod]->Print();
 691         fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
 692       }
 693     }
 694     fGeomMatrixSet = kTRUE;
 695     return kTRUE;
 696   }
 697
 698   Int_t runGM = event->GetRunNumber();
 699   TObjArray *mobj = 0;
 700
 701  if(fMisalignSurvey == kdefault)
 702  { //take default alignment corresponding to run no
 703     AliOADBContainer emcalgeoCont(Form("emcal"));
 704     emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 705     mobj=(TObjArray*)emcalgeoCont.GetObject(runGM,"EmcalMatrices");
 706  }
 707
 708  if(fMisalignSurvey == kSurveybyS)
 709  { //take alignment at sector level
 710   if (runGM <= 140000) { //2010 data
 711     AliOADBContainer emcalgeoCont(Form("emcal2010"));
 712     emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 713     mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
 714
 715   } else if (runGM>140000)
 716   { // 2011 LHC11a pass1 data
 717     AliOADBContainer emcalgeoCont(Form("emcal2011"));
 718     emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 719     mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byS");
 720   }
 721  }
 722
 723  if(fMisalignSurvey == kSurveybyM)
 724  { //take alignment at module level
 725   if (runGM <= 140000) { //2010 data
 726     AliOADBContainer emcalgeoCont(Form("emcal2010"));
 727     emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 728     mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
 729
 730   } else if (runGM>140000)
 731   { // 2011 LHC11a pass1 data
 732     AliOADBContainer emcalgeoCont(Form("emcal2011"));
 733     emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 734     mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byM");
 735   }
 736  }
 737
 738   if(!mobj)
 739   {
 740     AliFatal("Geometry matrix array not found");
 741     return kFALSE;
 742   }
 743
 744  for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
 745  {
 746    fEMCALMatrix[mod] = (TGeoHMatrix*) mobj->At(mod);
 747    fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
 748    fEMCALMatrix[mod]->Print();
 749  }
 750
 751   return kTRUE;
 752 }
 753
 754 //_____________________________________________________
 755 Int_t AliEMCALTenderSupply::InitBadChannels()
 756 {
 757   // Initialising bad channel maps
 758   AliESDEvent *event = fTender->GetEvent();
 759   if (!event)
 760     return 0;
 761
 762   if (fDebugLevel>0)
 763     AliInfo("Initialising Bad channel map");
 764
 765   // init default maps first
 766   if( !fEMCALRecoUtils->GetEMCALBadChannelStatusMapArray() )
 767     fEMCALRecoUtils->InitEMCALBadChannelStatusMap() ;
 768
 769   Int_t runBC = event->GetRunNumber();
 770
 771   AliOADBContainer *contBC = new AliOADBContainer("");
 772   if (fBasePath!="")
 773   { //if fBasePath specified in the ->SetBasePath()
 774     if (fDebugLevel>0) AliInfo(Form("Loading Bad Channels OADB from given path %s",fBasePath.Data()));
 775
 776     TFile *fbad=new TFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"read");
 777     if (!fbad || fbad->IsZombie())
 778     {
 779       AliFatal(Form("EMCALBadChannels.root was not found in the path provided: %s",fBasePath.Data()));
 780       return 0;
 781     }
 782
 783     if (fbad) delete fbad;
 784
 785     contBC->InitFromFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"AliEMCALBadChannels");
 786   }
 787   else
 788   { // Else choose the one in the $ALICE_ROOT directory
 789     if (fDebugLevel>0) AliInfo("Loading Bad Channels OADB from $ALICE_ROOT/OADB/EMCAL");
 790
 791     TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","read");
 792     if (!fbad || fbad->IsZombie())
 793     {
 794       AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root was not found");
 795       return 0;
 796     }
 797
 798     if (fbad) delete fbad;
 799
 800     contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","AliEMCALBadChannels");
 801   }
 802
 803   TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
 804   if (!arrayBC)
 805   {
 806     AliError(Form("No external hot channel set for run number: %d", runBC));
 807     return 2;
 808   }
 809
 810   TObjArray *arrayBCpass=(TObjArray*)arrayBC->FindObject(fFilepass); // Here, it looks for a specific pass
 811   if (!arrayBCpass)
 812   {
 813     AliError(Form("No external hot channel set for: %d -%s", runBC,fFilepass.Data()));
 814     return 2;
 815   }
 816
 817   if (fDebugLevel>0) arrayBCpass->Print();
 818
 819   Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
 820   for (Int_t i=0; i<sms; ++i)
 821   {
 822     TH2I *h = fEMCALRecoUtils->GetEMCALChannelStatusMap(i);
 823     if (h)
 824       delete h;
 825     h=(TH2I*)arrayBCpass->FindObject(Form("EMCALBadChannelMap_Mod%d",i));
 826
 827     if (!h)
 828     {
 829       AliError(Form("Can not get EMCALBadChannelMap_Mod%d",i));
 830       continue;
 831     }
 832     h->SetDirectory(0);
 833     fEMCALRecoUtils->SetEMCALChannelStatusMap(i,h);
 834   }
 835   return 1;
 836 }
 837
 838 //_____________________________________________________
 839 Int_t AliEMCALTenderSupply::InitRecalib()
 840 {
 841   // Initialising recalibration factors.
 842
 843   AliESDEvent *event = fTender->GetEvent();
 844   if (!event)
 845     return 0;
 846
 847   if (fDebugLevel>0)
 848     AliInfo("Initialising recalibration factors");
 849
 850   // init default maps first
 851   if( !fEMCALRecoUtils->GetEMCALRecalibrationFactorsArray() )
 852     fEMCALRecoUtils->InitEMCALRecalibrationFactors() ;
 853
 854   Int_t runRC = event->GetRunNumber();
 855
 856   AliOADBContainer *contRF=new AliOADBContainer("");
 857   if (fBasePath!="")
 858   { //if fBasePath specified in the ->SetBasePath()
 859     if (fDebugLevel>0)  AliInfo(Form("Loading Recalib OADB from given path %s",fBasePath.Data()));
 860
 861     TFile *fRecalib= new TFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"read");
 862     if (!fRecalib || fRecalib->IsZombie())
 863     {
 864       AliFatal(Form("EMCALRecalib.root not found in %s",fBasePath.Data()));
 865       return 0;
 866     }
 867
 868     if (fRecalib) delete fRecalib;
 869
 870     contRF->InitFromFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"AliEMCALRecalib");
 871   }
 872     else
 873     { // Else choose the one in the $ALICE_ROOT directory
 874       if (fDebugLevel>0)  AliInfo("Loading Recalib OADB from $ALICE_ROOT/OADB/EMCAL");
 875
 876       TFile *fRecalib= new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","read");
 877       if (!fRecalib || fRecalib->IsZombie())
 878       {
 879         AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root was not found");
 880         return 0;
 881       }
 882
 883       if (fRecalib) delete fRecalib;
 884
 885       contRF->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","AliEMCALRecalib");
 886     }
 887
 888   TObjArray *recal=(TObjArray*)contRF->GetObject(runRC); //GetObject(int runnumber)
 889   if (!recal)
 890   {
 891     AliError(Form("No Objects for run: %d",runRC));
 892     return 2;
 893   }
 894
 895   TObjArray *recalpass=(TObjArray*)recal->FindObject(fFilepass);
 896   if (!recalpass)
 897   {
 898     AliError(Form("No Objects for run: %d - %s",runRC,fFilepass.Data()));
 899     return 2;
 900   }
 901
 902   TObjArray *recalib=(TObjArray*)recalpass->FindObject("Recalib");
 903   if (!recalib)
 904   {
 905     AliError(Form("No Recalib histos found for  %d - %s",runRC,fFilepass.Data()));
 906     return 2;
 907   }
 908
 909   if (fDebugLevel>0) recalib->Print();
 910
 911   Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
 912   for (Int_t i=0; i<sms; ++i)
 913   {
 914     TH2F *h = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactors(i);
 915     if (h)
 916       delete h;
 917     h = (TH2F*)recalib->FindObject(Form("EMCALRecalFactors_SM%d",i));
 918     if (!h)
 919     {
 920       AliError(Form("Could not load EMCALRecalFactors_SM%d",i));
 921       continue;
 922     }
 923     h->SetDirectory(0);
 924     fEMCALRecoUtils->SetEMCALChannelRecalibrationFactors(i,h);
 925   }
 926   return 1;
 927 }
 928
 929 //_____________________________________________________
 930 Int_t AliEMCALTenderSupply::InitTimeCalibration()
 931 {
 932   // Initialising bad channel maps
 933   AliESDEvent *event = fTender->GetEvent();
 934   if (!event)
 935     return 0;
 936
 937   if (fDebugLevel>0)
 938     AliInfo("Initialising time calibration map");
 939
 940   // init default maps first
 941   if( !fEMCALRecoUtils->GetEMCALTimeRecalibrationFactorsArray() )
 942     fEMCALRecoUtils->InitEMCALTimeRecalibrationFactors() ;
 943
 944   Int_t runBC = event->GetRunNumber();
 945
 946   AliOADBContainer *contBC = new AliOADBContainer("");
 947   if (fBasePath!="")
 948   { //if fBasePath specified in the ->SetBasePath()
 949     if (fDebugLevel>0) AliInfo(Form("Loading time calibration OADB from given path %s",fBasePath.Data()));
 950
 951     TFile *fbad=new TFile(Form("%s/EMCALTimeCalib.root",fBasePath.Data()),"read");
 952     if (!fbad || fbad->IsZombie())
 953     {
 954       AliFatal(Form("EMCALTimeCalib.root was not found in the path provided: %s",fBasePath.Data()));
 955       return 0;
 956     }
 957
 958     if (fbad) delete fbad;
 959
 960     contBC->InitFromFile(Form("%s/EMCALTimeCalib.root",fBasePath.Data()),"AliEMCALTimeCalib");
 961   }
 962   else
 963   { // Else choose the one in the $ALICE_ROOT directory
 964     if (fDebugLevel>0) AliInfo("Loading time calibration OADB from $ALICE_ROOT/OADB/EMCAL");
 965
 966     TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root","read");
 967     if (!fbad || fbad->IsZombie())
 968     {
 969       AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root was not found");
 970       return 0;
 971     }
 972
 973     if (fbad) delete fbad;
 974
 975     contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root","AliEMCALTimeCalib");
 976   }
 977
 978   TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
 979   if (!arrayBC)
 980   {
 981     AliError(Form("No external time calibration set for run number: %d", runBC));
 982     return 2;
 983   }
 984
 985   TObjArray *arrayBCpass=(TObjArray*)arrayBC->FindObject(fFilepass); // Here, it looks for a specific pass
 986   if (!arrayBCpass)
 987   {
 988     AliError(Form("No external time calibration set for: %d -%s", runBC,fFilepass.Data()));
 989     return 2;
 990   }
 991
 992   if (fDebugLevel>0) arrayBCpass->Print();
 993
 994   for( Int_t i = 0; i < 4; i++ )
 995   {
 996     TH1F *h = fEMCALRecoUtils->GetEMCALChannelTimeRecalibrationFactors( i );
 997     if( h )
 998       delete h;
 999
1000     h = (TH1F*)arrayBCpass->FindObject(Form("hAllTimeAvBC%d",i));
1001
1002     if (!h)
1003     {
1004       AliError(Form("Can not get hAllTimeAvBC%d",i));
1005       continue;
1006     }
1007     h->SetDirectory(0);
1008     fEMCALRecoUtils->SetEMCALChannelTimeRecalibrationFactors(i,h);
1009   }
1010   return 1;
1011 }
1012
1013 //_____________________________________________________
1014 void AliEMCALTenderSupply::UpdateCells()
1015 {
1016   //Remove bad cells from the cell list
1017   //Recalibrate energy and time cells
1018   //This is required for later reclusterization
1019
1020   AliESDCaloCells *cells = fTender->GetEvent()->GetEMCALCells();
1021   Int_t bunchCrossNo = fTender->GetEvent()->GetBunchCrossNumber();
1022
1023   fEMCALRecoUtils->RecalibrateCells(cells, bunchCrossNo);
1024
1025   // remove exotic cells - loop through cells and zero the exotic ones
1026   // just like with bad cell rejection in reco utils (inside RecalibrateCells)
1027   if( fRejectExoticCells )
1028   {
1029     Int_t    absId  =-1;
1030     Bool_t   isExot = kFALSE;
1031
1032     // loop through cells
1033     Int_t nEMcell  = cells->GetNumberOfCells() ;
1034     for (Int_t iCell = 0; iCell < nEMcell; iCell++)
1035     {
1036       absId  = cells->GetCellNumber(iCell);
1037
1038       isExot = fEMCALRecoUtils->IsExoticCell( absId, cells, bunchCrossNo );
1039       // zero if exotic
1040       if( isExot )
1041         cells->SetCell( iCell, absId, 0.0, 0.0 );
1042     } // cell loop
1043   } // reject exotic cells
1044
1045   cells->Sort();
1046 }
1047
1048 //_____________________________________________________
1049 void AliEMCALTenderSupply::InitRecParam()
1050 {
1051   // Initalize the reconstructor parameters
1052   // Depending on the data type, different type of clusterizer
1053
1054   if (fDebugLevel>0)
1055     AliInfo("Initialize the recParam");
1056
1057   fRecParam = new AliEMCALRecParam;
1058   const AliESDRun *run = fTender->GetEvent()->GetESDRun();
1059   TString bt(run->GetBeamType());
1060   if (bt=="A-A")
1061   {
1062     fRecParam->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv2);
1063     fRecParam->SetClusteringThreshold(0.100);
1064     fRecParam->SetMinECut(0.050);
1065   }
1066   else
1067   {
1068     fRecParam->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv1);
1069     fRecParam->SetClusteringThreshold(0.100);
1070     fRecParam->SetMinECut(0.050);
1071   }
1072 }
1073
1074 //_____________________________________________________
1075 Bool_t AliEMCALTenderSupply::InitClusterization()
1076 {
1077   // Initialing clusterization/unfolding algorithm and set all the needed parameters.
1078
1079   AliESDEvent *event = fTender->GetEvent();
1080   if (!event)
1081     return kFALSE;
1082
1083   if (fDebugLevel>0)
1084     AliInfo(Form("Initialising reclustering parameters: Clusterizer type: %d",fRecParam->GetClusterizerFlag()));
1085
1086   //---setup clusterizer
1087   delete fClusterizer;
1088   if     (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
1089     fClusterizer = new AliEMCALClusterizerv1 (fEMCALGeo);
1090   else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
1091     fClusterizer = new AliEMCALClusterizerv2(fEMCALGeo);
1092   else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN)
1093   {
1094     AliEMCALClusterizerNxN *clusterizer = new AliEMCALClusterizerNxN(fEMCALGeo);
1095     clusterizer->SetNRowDiff(fRecParam->GetNRowDiff());
1096     clusterizer->SetNColDiff(fRecParam->GetNColDiff());
1097     fClusterizer = clusterizer;
1098   }
1099   else
1100   {
1101     AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
1102     return kFALSE;
1103   }
1104
1105   // Set the clustering parameters
1106   fClusterizer->SetECAClusteringThreshold( fRecParam->GetClusteringThreshold() );
1107   fClusterizer->SetECALogWeight          ( fRecParam->GetW0()                  );
1108   fClusterizer->SetMinECut               ( fRecParam->GetMinECut()             );
1109   fClusterizer->SetUnfolding             ( fRecParam->GetUnfold()              );
1110   fClusterizer->SetECALocalMaxCut        ( fRecParam->GetLocMaxCut()           );
1111   fClusterizer->SetTimeCut               ( fRecParam->GetTimeCut()             );
1112   fClusterizer->SetTimeMin               ( fRecParam->GetTimeMin()             );
1113   fClusterizer->SetTimeMax               ( fRecParam->GetTimeMax()             );
1114   fClusterizer->SetInputCalibrated       ( kTRUE                               );
1115   fClusterizer->SetJustClusters          ( kTRUE                               );
1116
1117   // In case of unfolding after clusterization is requested, set the corresponding parameters
1118   if (fRecParam->GetUnfold())
1119   {
1120     for (Int_t i = 0; i < 8; ++i)
1121     {
1122       fClusterizer->SetSSPars(i, fRecParam->GetSSPars(i));
1123     }
1124     for (Int_t i = 0; i < 3; ++i)
1125     {
1126       fClusterizer->SetPar5  (i, fRecParam->GetPar5(i));
1127       fClusterizer->SetPar6  (i, fRecParam->GetPar6(i));
1128     }
1129     fClusterizer->InitClusterUnfolding();
1130   }
1131
1132   fClusterizer->SetDigitsArr(fDigitsArr);
1133   fClusterizer->SetOutput(0);
1134   fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
1135   return kTRUE;
1136 }
1137
1138 //_____________________________________________________
1139 void AliEMCALTenderSupply::FillDigitsArray()
1140 {
1141   // Fill digits from cells to a TClonesArray.
1142
1143   AliESDEvent *event = fTender->GetEvent();
1144   if (!event)
1145     return;
1146
1147   fDigitsArr->Clear("C");
1148   AliESDCaloCells *cells = event->GetEMCALCells();
1149   Int_t ncells = cells->GetNumberOfCells();
1150   for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell)
1151   {
1152     Double_t cellAmplitude=0, cellTime=0;
1153     Short_t  cellNumber=0;
1154
1155     if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime) != kTRUE)
1156       break;
1157
1158     // Do not add if already too low (some cells set to 0 if bad channels)
1159     if (cellAmplitude < fRecParam->GetMinECut() )
1160       continue;
1161
1162     // If requested, do not include exotic cells
1163    if (fEMCALRecoUtils->IsExoticCell(cellNumber,cells,event->GetBunchCrossNumber()))
1164       continue;
1165
1166     AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
1167     digit->SetId(cellNumber);
1168     digit->SetTime(cellTime);
1169     digit->SetTimeR(cellTime);
1170     digit->SetIndexInList(idigit);
1171     digit->SetType(AliEMCALDigit::kHG);
1172     digit->SetAmplitude(cellAmplitude);
1173     idigit++;
1174   }
1175 }
1176
1177 //_____________________________________________________
1178 void AliEMCALTenderSupply::Clusterize()
1179 {
1180   // Clusterize.
1181
1182   fClusterizer->Digits2Clusters("");
1183 }
1184
1185 //_____________________________________________________
1186 void AliEMCALTenderSupply::UpdateClusters()
1187 {
1188   // Update ESD cluster list.
1189
1190   AliESDEvent *event = fTender->GetEvent();
1191   if (!event)
1192     return;
1193
1194   TClonesArray *clus = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
1195   if (!clus)
1196     clus = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
1197   if (!clus)
1198   {
1199     AliError(" Null pointer to calo clusters array, returning");
1200     return;
1201   }
1202
1203   Int_t nents = clus->GetEntriesFast();
1204   for (Int_t i=0; i < nents; ++i)
1205   {
1206     AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->At(i));
1207     if (!c)
1208       continue;
1209     if (c->IsEMCAL())
1210       delete clus->RemoveAt(i);
1211   }
1212
1213   clus->Compress();
1214
1215   RecPoints2Clusters(clus);
1216
1217 }
1218
1219 //_____________________________________________________
1220 void AliEMCALTenderSupply::RecPoints2Clusters(TClonesArray *clus)
1221 {
1222   // Convert AliEMCALRecoPoints to AliESDCaloClusters.
1223   // Cluster energy, global position, cells and their amplitude fractions are restored.
1224
1225   AliESDEvent *event = fTender->GetEvent();
1226   if (!event)
1227     return;
1228
1229   Int_t ncls = fClusterArr->GetEntriesFast();
1230   for(Int_t i=0, nout=clus->GetEntriesFast(); i < ncls; ++i)
1231   {
1232     AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
1233
1234     Int_t ncellsTrue = 0;
1235     const Int_t ncells = recpoint->GetMultiplicity();
1236     UShort_t   absIds[ncells];
1237     Double32_t ratios[ncells];
1238     Int_t *dlist = recpoint->GetDigitsList();
1239     Float_t *elist = recpoint->GetEnergiesList();
1240     for (Int_t c = 0; c < ncells; ++c)
1241     {
1242       AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
1243       absIds[ncellsTrue] = digit->GetId();
1244       ratios[ncellsTrue] = elist[c]/digit->GetAmplitude();
1245       if (ratios[ncellsTrue] < 0.001)
1246         continue;
1247       ++ncellsTrue;
1248     }
1249
1250     if (ncellsTrue < 1)
1251     {
1252       AliWarning("Skipping cluster with no cells");
1253       continue;
1254     }
1255
1256     // calculate new cluster position
1257     TVector3 gpos;
1258     recpoint->GetGlobalPosition(gpos);
1259     Float_t g[3];
1260     gpos.GetXYZ(g);
1261
1262     AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->New(nout++));
1263     c->SetID(nout-1);
1264     c->SetType(AliVCluster::kEMCALClusterv1);
1265     c->SetE(recpoint->GetEnergy());
1266     c->SetPosition(g);
1267     c->SetNCells(ncellsTrue);
1268     c->SetDispersion(recpoint->GetDispersion());
1269     c->SetEmcCpvDistance(-1);            //not yet implemented
1270     c->SetChi2(-1);                      //not yet implemented
1271     c->SetTOF(recpoint->GetTime()) ;     //time-of-flight
1272     c->SetNExMax(recpoint->GetNExMax()); //number of local maxima
1273     Float_t elipAxis[2];
1274     recpoint->GetElipsAxis(elipAxis);
1275     c->SetM02(elipAxis[0]*elipAxis[0]) ;
1276     c->SetM20(elipAxis[1]*elipAxis[1]) ;
1277     AliESDCaloCluster *cesd = static_cast<AliESDCaloCluster*>(c);
1278     cesd->SetCellsAbsId(absIds);
1279     cesd->SetCellsAmplitudeFraction(ratios);
1280   }
1281 }
1282
1283 //_____________________________________________________
1284 void AliEMCALTenderSupply::GetPass()
1285 {
1286   // Get passx from filename.
1287
1288   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1289   fInputTree = mgr->GetTree();
1290
1291   if (!fInputTree)
1292   {
1293     AliError("Pointer to tree = 0, returning");
1294     return;
1295   }
1296
1297   fInputFile = fInputTree->GetCurrentFile();
1298   if (!fInputFile) {
1299     AliError("Null pointer input file, returning");
1300     return;
1301   }
1302
1303   TString fname(fInputFile->GetName());
1304   if      (fname.Contains("pass1")) fFilepass = TString("pass1");
1305   else if (fname.Contains("pass2")) fFilepass = TString("pass2");
1306   else if (fname.Contains("pass3")) fFilepass = TString("pass3");
1307   else
1308   {
1309     AliError(Form("Pass number string not found: %s", fname.Data()));
1310     return;
1311   }
1312 }
1313