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 //  EMCAL tender, apply corrections to EMCAL clusters and do track matching. //
  19 //                                                                           //
  20 //  Author: Deepa Thomas (Utrecht University)                                //
  21 //  Later mods/rewrite: Jiri Kral (University of Jyvaskyla)                  //
  22 //  S. Aiola, C. Loizides : Make it work for AODs                            //
  23 //                                                                           //
  24 ///////////////////////////////////////////////////////////////////////////////
  25
  26 #include <TClonesArray.h>
  27 #include <TFile.h>
  28 #include <TGeoGlobalMagField.h>
  29 #include <TInterpreter.h>
  30 #include <TObjArray.h>
  31 #include <TROOT.h>
  32 #include <TTree.h>
  33 #include "AliAODEvent.h"
  34 #include "AliAODMCParticle.h"
  35 #include "AliAnalysisManager.h"
  36 #include "AliEMCALAfterBurnerUF.h"
  37 #include "AliEMCALClusterizer.h"
  38 #include "AliEMCALClusterizerNxN.h"
  39 #include "AliEMCALClusterizerv1.h"
  40 #include "AliEMCALClusterizerv2.h"
  41 #include "AliEMCALDigit.h"
  42 #include "AliEMCALGeometry.h"
  43 #include "AliEMCALRecParam.h"
  44 #include "AliEMCALRecParam.h"
  45 #include "AliEMCALRecPoint.h"
  46 #include "AliEMCALRecoUtils.h"
  47 #include "AliESDCaloCluster.h"
  48 #include "AliESDEvent.h"
  49 #include "AliLog.h"
  50 #include "AliMagF.h"
  51 #include "AliOADBContainer.h"
  52 #include "AliTender.h"
  53 #include "AliEMCALTenderSupply.h"
  54
  55 ClassImp(AliEMCALTenderSupply)
  56
  57 AliEMCALTenderSupply::AliEMCALTenderSupply() :
  58   AliTenderSupply()
  59   ,fTask(0)
  60   ,fRun(0)
  61   ,fEMCALGeo(0x0)
  62   ,fEMCALGeoName("")
  63   ,fEMCALRecoUtils(0)
  64   ,fConfigName("")
  65   ,fDebugLevel(0)
  66   ,fNonLinearFunc(-1)
  67   ,fNonLinearThreshold(-1)
  68   ,fReCalibCluster(kFALSE)
  69   ,fUpdateCell(kFALSE)
  70   ,fCalibrateEnergy(kFALSE)
  71   ,fCalibrateTime(kFALSE)
  72   ,fCalibrateTimeParamAvailable(kFALSE)
  73   ,fDoNonLinearity(kFALSE)
  74   ,fBadCellRemove(kFALSE)
  75   ,fRejectExoticCells(kFALSE)
  76   ,fRejectExoticClusters(kFALSE)
  77   ,fClusterBadChannelCheck(kFALSE)
  78   ,fRecalClusPos(kFALSE)
  79   ,fFiducial(kFALSE)
  80   ,fNCellsFromEMCALBorder(-1)
  81   ,fRecalDistToBadChannels(kFALSE)
  82   ,fRecalShowerShape(kFALSE)
  83   ,fInputTree(0)
  84   ,fInputFile(0)
  85   ,fGetPassFromFileName(kTRUE)
  86   ,fFilepass(0)
  87   ,fMass(-1)
  88   ,fStep(-1)
  89   ,fCutEtaPhiSum(kTRUE)
  90   ,fCutEtaPhiSeparate(kFALSE)
  91   ,fRcut(-1)
  92   ,fEtacut(-1)
  93   ,fPhicut(-1)
  94   ,fBasePath("")
  95   ,fReClusterize(kFALSE)
  96   ,fClusterizer(0)
  97   ,fGeomMatrixSet(kFALSE)
  98   ,fLoadGeomMatrices(kFALSE)
  99   ,fRecParam(0x0)
 100   ,fDoTrackMatch(kFALSE)
 101   ,fDoUpdateOnly(kFALSE)
 102   ,fUnfolder(0)
 103   ,fDigitsArr(0)
 104   ,fClusterArr(0)
 105   ,fMisalignSurvey(kdefault)
 106   ,fExoticCellFraction(-1)
 107   ,fExoticCellDiffTime(-1)
 108   ,fExoticCellMinAmplitude(-1)
 109   ,fSetCellMCLabelFromCluster(0)
 110   ,fTempClusterArr(0)
 111   ,fRemapMCLabelForAODs(0)
 112 {
 113   // Default constructor.
 114
 115   for(Int_t i = 0; i < 12;    i++) fEMCALMatrix[i]      = 0 ;
 116   for(Int_t j = 0; j < 12672; j++) fOrgClusterCellId[j] = -1;
 117 }
 118
 119 //_____________________________________________________
 120 AliEMCALTenderSupply::AliEMCALTenderSupply(const char *name, const AliTender *tender) :
 121   AliTenderSupply(name,tender)
 122   ,fTask(0)
 123   ,fRun(0)
 124   ,fEMCALGeo(0x0)
 125   ,fEMCALGeoName("")
 126   ,fEMCALRecoUtils(0)
 127   ,fConfigName("")
 128   ,fDebugLevel(0)
 129   ,fNonLinearFunc(-1)
 130   ,fNonLinearThreshold(-1)
 131   ,fReCalibCluster(kFALSE)
 132   ,fUpdateCell(kFALSE)
 133   ,fCalibrateEnergy(kFALSE)
 134   ,fCalibrateTime(kFALSE)
 135   ,fCalibrateTimeParamAvailable(kFALSE)
 136   ,fDoNonLinearity(kFALSE)
 137   ,fBadCellRemove(kFALSE)
 138   ,fRejectExoticCells(kFALSE)
 139   ,fRejectExoticClusters(kFALSE)
 140   ,fClusterBadChannelCheck(kFALSE)
 141   ,fRecalClusPos(kFALSE)
 142   ,fFiducial(kFALSE)
 143   ,fNCellsFromEMCALBorder(-1)
 144   ,fRecalDistToBadChannels(kFALSE)
 145   ,fRecalShowerShape(kFALSE)
 146   ,fInputTree(0)
 147   ,fInputFile(0)
 148   ,fGetPassFromFileName(kTRUE)
 149   ,fFilepass("")
 150   ,fMass(-1)
 151   ,fStep(-1)
 152   ,fCutEtaPhiSum(kTRUE)
 153   ,fCutEtaPhiSeparate(kFALSE)
 154   ,fRcut(-1)
 155   ,fEtacut(-1)
 156   ,fPhicut(-1)
 157   ,fBasePath("")
 158   ,fReClusterize(kFALSE)
 159   ,fClusterizer(0)
 160   ,fGeomMatrixSet(kFALSE)
 161   ,fLoadGeomMatrices(kFALSE)
 162   ,fRecParam(0x0)
 163   ,fDoTrackMatch(kFALSE)
 164   ,fDoUpdateOnly(kFALSE)
 165   ,fUnfolder(0)
 166   ,fDigitsArr(0)
 167   ,fClusterArr(0)
 168   ,fMisalignSurvey(kdefault)
 169   ,fExoticCellFraction(-1)
 170   ,fExoticCellDiffTime(-1)
 171   ,fExoticCellMinAmplitude(-1)
 172   ,fSetCellMCLabelFromCluster(0)
 173   ,fTempClusterArr(0)
 174   ,fRemapMCLabelForAODs(0)
 175 {
 176   // Named constructor
 177
 178   for(Int_t i = 0; i < 12;    i++) fEMCALMatrix[i]      = 0 ;
 179   for(Int_t j = 0; j < 12672; j++) fOrgClusterCellId[j] = -1;
 180 }
 181
 182 //_____________________________________________________
 183 AliEMCALTenderSupply::AliEMCALTenderSupply(const char *name, AliAnalysisTaskSE *task) :
 184   AliTenderSupply(name)
 185   ,fTask(task)
 186   ,fRun(0)
 187   ,fEMCALGeo(0x0)
 188   ,fEMCALGeoName("")
 189   ,fEMCALRecoUtils(0)
 190   ,fConfigName("")
 191   ,fDebugLevel(0)
 192   ,fNonLinearFunc(-1)
 193   ,fNonLinearThreshold(-1)
 194   ,fReCalibCluster(kFALSE)
 195   ,fUpdateCell(kFALSE)
 196   ,fCalibrateEnergy(kFALSE)
 197   ,fCalibrateTime(kFALSE)
 198   ,fCalibrateTimeParamAvailable(kFALSE)
 199   ,fDoNonLinearity(kFALSE)
 200   ,fBadCellRemove(kFALSE)
 201   ,fRejectExoticCells(kFALSE)
 202   ,fRejectExoticClusters(kFALSE)
 203   ,fClusterBadChannelCheck(kFALSE)
 204   ,fRecalClusPos(kFALSE)
 205   ,fFiducial(kFALSE)
 206   ,fNCellsFromEMCALBorder(-1)
 207   ,fRecalDistToBadChannels(kFALSE)
 208   ,fRecalShowerShape(kFALSE)
 209   ,fInputTree(0)
 210   ,fInputFile(0)
 211   ,fGetPassFromFileName(kTRUE)
 212   ,fFilepass("")
 213   ,fMass(-1)
 214   ,fStep(-1)
 215   ,fCutEtaPhiSum(kTRUE)
 216   ,fCutEtaPhiSeparate(kFALSE)
 217   ,fRcut(-1)
 218   ,fEtacut(-1)
 219   ,fPhicut(-1)
 220   ,fBasePath("")
 221   ,fReClusterize(kFALSE)
 222   ,fClusterizer(0)
 223   ,fGeomMatrixSet(kFALSE)
 224   ,fLoadGeomMatrices(kFALSE)
 225   ,fRecParam(0x0)
 226   ,fDoTrackMatch(kFALSE)
 227   ,fDoUpdateOnly(kFALSE)
 228   ,fUnfolder(0)
 229   ,fDigitsArr(0)
 230   ,fClusterArr(0)
 231   ,fMisalignSurvey(kdefault)
 232   ,fExoticCellFraction(-1)
 233   ,fExoticCellDiffTime(-1)
 234   ,fExoticCellMinAmplitude(-1)
 235   ,fSetCellMCLabelFromCluster(0)
 236   ,fTempClusterArr(0)
 237   ,fRemapMCLabelForAODs(0)
 238 {
 239   // Named constructor.
 240
 241   for(Int_t i = 0; i < 12;    i++) fEMCALMatrix[i]      = 0 ;
 242   for(Int_t j = 0; j < 12672; j++) fOrgClusterCellId[j] = -1;
 243 }
 244
 245 //_____________________________________________________
 246 AliEMCALTenderSupply::~AliEMCALTenderSupply()
 247 {
 248   //Destructor
 249
 250   if (!AliAnalysisManager::GetAnalysisManager())  return;
 251
 252   if (!AliAnalysisManager::GetAnalysisManager()->IsProofMode())
 253   {
 254     delete fEMCALRecoUtils;
 255     delete fRecParam;
 256     delete fUnfolder;
 257
 258     if (!fClusterizer)
 259     {
 260       if (fDigitsArr)
 261       {
 262         fDigitsArr->Clear("C");
 263         delete fDigitsArr;
 264       }
 265     }
 266     else
 267     {
 268       delete fClusterizer;
 269       fDigitsArr = 0;
 270     }
 271   }
 272 }
 273
 274 //_____________________________________________________
 275 void AliEMCALTenderSupply::SetDefaults()
 276 {
 277   // Set default settings.
 278
 279   SwitchOnReclustering();
 280   SwitchOnTrackMatch();
 281 }
 282
 283 //_____________________________________________________
 284 Bool_t AliEMCALTenderSupply::RunChanged() const
 285 {
 286   // Get run number.
 287
 288   return (fTender && fTender->RunChanged()) || (fTask && fRun != fTask->InputEvent()->GetRunNumber());
 289 }
 290
 291 //_____________________________________________________
 292 void AliEMCALTenderSupply::Init()
 293 {
 294   // Initialise EMCAL tender.
 295
 296   if (fDebugLevel>0)
 297     AliWarning("Init EMCAL Tender supply");
 298
 299   if (fConfigName.Length()>0 && gROOT->LoadMacro(fConfigName) >=0) {
 300     AliDebug(1, Form("Loading settings from macro %s", fConfigName.Data()));
 301     AliEMCALTenderSupply *tender = (AliEMCALTenderSupply*)gInterpreter->ProcessLine("ConfigEMCALTenderSupply()");
 302     fDebugLevel             = tender->fDebugLevel;
 303     fEMCALGeoName           = tender->fEMCALGeoName;
 304     fEMCALRecoUtils         = tender->fEMCALRecoUtils;
 305     fConfigName             = tender->fConfigName;
 306     fNonLinearFunc          = tender->fNonLinearFunc;
 307     fNonLinearThreshold     = tender->fNonLinearThreshold;
 308     fReCalibCluster         = tender->fReCalibCluster;
 309     fUpdateCell             = tender->fUpdateCell;
 310     fRecalClusPos           = tender->fRecalClusPos;
 311     fCalibrateEnergy        = tender->fCalibrateEnergy;
 312     fCalibrateTime          = tender->fCalibrateTime;
 313     fCalibrateTimeParamAvailable = tender->fCalibrateTimeParamAvailable;
 314     fFiducial               = tender->fFiducial;
 315     fNCellsFromEMCALBorder  = tender->fNCellsFromEMCALBorder;
 316     fRecalDistToBadChannels = tender->fRecalDistToBadChannels;
 317     fRecalShowerShape       = tender->fRecalShowerShape;
 318     fClusterBadChannelCheck = tender->fClusterBadChannelCheck;
 319     fBadCellRemove          = tender->fBadCellRemove;
 320     fRejectExoticCells      = tender->fRejectExoticCells;
 321     fRejectExoticClusters   = tender->fRejectExoticClusters;
 322     fMass                   = tender->fMass;
 323     fStep                   = tender->fStep;
 324     fCutEtaPhiSum           = tender->fCutEtaPhiSum;
 325     fCutEtaPhiSeparate      = tender->fCutEtaPhiSeparate;
 326     fRcut                   = tender->fRcut;
 327     fEtacut                 = tender->fEtacut;
 328     fPhicut                 = tender->fPhicut;
 329     fReClusterize           = tender->fReClusterize;
 330     fLoadGeomMatrices       = tender->fLoadGeomMatrices;
 331     fRecParam               = tender->fRecParam;
 332     fDoNonLinearity         = tender->fDoNonLinearity;
 333     fDoTrackMatch           = tender->fDoTrackMatch;
 334     fDoUpdateOnly           = tender->fDoUpdateOnly;
 335     fMisalignSurvey         = tender->fMisalignSurvey;
 336     fExoticCellFraction     = tender->fExoticCellFraction;
 337     fExoticCellDiffTime     = tender->fExoticCellDiffTime;
 338     fExoticCellMinAmplitude = tender->fExoticCellMinAmplitude;
 339
 340     for(Int_t i = 0; i < 12; i++)
 341       fEMCALMatrix[i] = tender->fEMCALMatrix[i] ;
 342   }
 343
 344   if (fDebugLevel>0){
 345     AliInfo("Emcal Tender settings: ======================================");
 346     AliInfo("------------ Switches --------------------------");
 347     AliInfo(Form("BadCellRemove : %d", fBadCellRemove));
 348     AliInfo(Form("ExoticCellRemove : %d", fRejectExoticCells));
 349     AliInfo(Form("CalibrateEnergy : %d", fCalibrateEnergy));
 350     AliInfo(Form("CalibrateTime : %d", fCalibrateTime));
 351     AliInfo(Form("UpdateCell : %d", fUpdateCell));
 352     AliInfo(Form("DoUpdateOnly : %d", fDoUpdateOnly));
 353     AliInfo(Form("Reclustering : %d", fReClusterize));
 354     AliInfo(Form("ClusterBadChannelCheck : %d", fClusterBadChannelCheck));
 355     AliInfo(Form("ClusterExoticChannelCheck : %d", fRejectExoticClusters));
 356     AliInfo(Form("CellFiducialRegion : %d", fFiducial));
 357     AliInfo(Form("ReCalibrateCluster : %d", fReCalibCluster));
 358     AliInfo(Form("RecalculateClusPos : %d", fRecalClusPos));
 359     AliInfo(Form("RecalShowerShape : %d", fRecalShowerShape));
 360     AliInfo(Form("NonLinearityCorrection : %d", fDoNonLinearity));
 361     AliInfo(Form("RecalDistBadChannel : %d", fRecalDistToBadChannels));
 362     AliInfo(Form("TrackMatch : %d", fDoTrackMatch));
 363     AliInfo("------------ Variables -------------------------");
 364     AliInfo(Form("DebugLevel : %d", fDebugLevel));
 365     AliInfo(Form("BasePath : %s", fBasePath.Data()));
 366     AliInfo(Form("ConfigFileName : %s", fConfigName.Data()));
 367     AliInfo(Form("EMCALGeometryName : %s", fEMCALGeoName.Data()));
 368     AliInfo(Form("NonLinearityFunction : %d", fNonLinearFunc));
 369     AliInfo(Form("NonLinearityThreshold : %d", fNonLinearThreshold));
 370     AliInfo(Form("MisalignmentMatrixSurvey : %d", fMisalignSurvey));
 371     AliInfo(Form("NumberOfCellsFromEMCALBorder : %d", fNCellsFromEMCALBorder));
 372     AliInfo(Form("RCut : %f", fRcut));
 373     AliInfo(Form("Mass : %f", fMass));
 374     AliInfo(Form("Step : %f", fStep));
 375     AliInfo(Form("EtaCut : %f", fEtacut));
 376     AliInfo(Form("PhiCut : %f", fPhicut));
 377     AliInfo(Form("ExoticCellFraction : %f", fExoticCellFraction));
 378     AliInfo(Form("ExoticCellDiffTime : %f", fExoticCellDiffTime));
 379     AliInfo(Form("ExoticCellMinAmplitude : %f", fExoticCellMinAmplitude));
 380     AliInfo("=============================================================");
 381   }
 382
 383   // init reco utils
 384
 385   if (!fEMCALRecoUtils)
 386     fEMCALRecoUtils  = new AliEMCALRecoUtils;
 387
 388   // init geometry if requested
 389   if (fEMCALGeoName.Length()>0)
 390     fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
 391
 392   // digits array
 393   fDigitsArr       = new TClonesArray("AliEMCALDigit",1000);
 394
 395   // initialising non-linearity parameters
 396   if (fNonLinearThreshold != -1)
 397     fEMCALRecoUtils->SetNonLinearityThreshold(fNonLinearThreshold);
 398   if (fNonLinearFunc != -1)
 399     fEMCALRecoUtils->SetNonLinearityFunction(fNonLinearFunc);
 400
 401   // missalignment function
 402   fEMCALRecoUtils->SetPositionAlgorithm(AliEMCALRecoUtils::kPosTowerGlobal);
 403
 404   // fiducial cut
 405   // do not do the eta0 fiducial cut
 406   if (fNCellsFromEMCALBorder != -1)
 407     fEMCALRecoUtils->SetNumberOfCellsFromEMCALBorder(fNCellsFromEMCALBorder);
 408   fEMCALRecoUtils->SwitchOnNoFiducialBorderInEMCALEta0();
 409
 410   // exotic cell rejection
 411   if (fExoticCellFraction != -1)
 412     fEMCALRecoUtils->SetExoticCellFractionCut(fExoticCellFraction);
 413   if (fExoticCellDiffTime != -1)
 414     fEMCALRecoUtils->SetExoticCellDiffTimeCut(fExoticCellDiffTime);
 415   if (fExoticCellMinAmplitude != -1)
 416     fEMCALRecoUtils->SetExoticCellMinAmplitudeCut(fExoticCellMinAmplitude);
 417
 418   // setting track matching parameters ... mass, step size and residual cut
 419   if (fMass != -1)
 420     fEMCALRecoUtils->SetMass(fMass);
 421   if (fStep != -1)
 422     fEMCALRecoUtils->SetStep(fStep);
 423
 424   // spatial cut based on separate eta/phi or common processing
 425   if (fCutEtaPhiSum) {
 426     fEMCALRecoUtils->SwitchOnCutEtaPhiSum();
 427     if (fRcut != -1)
 428       fEMCALRecoUtils->SetCutR(fRcut);
 429   } else if (fCutEtaPhiSeparate) {
 430     fEMCALRecoUtils->SwitchOnCutEtaPhiSeparate();
 431     if (fEtacut != -1)
 432       fEMCALRecoUtils->SetCutEta(fEtacut);
 433     if (fPhicut != -1)
 434       fEMCALRecoUtils->SetCutPhi(fPhicut);
 435   }
 436 }
 437
 438 //_____________________________________________________
 439 AliVEvent* AliEMCALTenderSupply::GetEvent()
 440 {
 441   // Return the event pointer.
 442
 443   if (fTender) {
 444     return fTender->GetEvent();
 445   } else if (fTask) {
 446     return fTask->InputEvent();
 447   }
 448
 449   return 0;
 450 }
 451
 452 //_____________________________________________________
 453 void AliEMCALTenderSupply::ProcessEvent()
 454 {
 455   // Event loop.
 456
 457   AliVEvent *event = GetEvent();
 458
 459   if (!event) {
 460     AliError("Event ptr = 0, returning");
 461     return;
 462   }
 463
 464   // Initialising parameters once per run number
 465   if (RunChanged()) {
 466
 467     fRun = event->GetRunNumber();
 468     AliWarning(Form("Run changed, initializing parameters for %d", fRun));
 469     if (dynamic_cast<AliAODEvent*>(event)) {
 470       AliWarning("=============================================================");
 471       AliWarning("===  Running on AOD is not equivalent to running on ESD!  ===");
 472       AliWarning("=============================================================");
 473     }
 474
 475     // init geometry if not already done
 476     if (fEMCALGeoName.Length()==0) {
 477       fEMCALGeoName = "EMCAL_FIRSTYEARV1";
 478       if (fRun>139517) {
 479         fEMCALGeoName = "EMCAL_COMPLETEV1";
 480       }
 481       if (fRun>170593) {
 482         fEMCALGeoName = "EMCAL_COMPLETE12SMV1";
 483       }
 484       fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName);
 485       if (!fEMCALGeo) {
 486         AliFatal(Form("Can not create geometry: %s", fEMCALGeoName.Data()));
 487         return;
 488       }
 489     }
 490
 491     // get pass
 492     if (fGetPassFromFileName)
 493       GetPass();
 494
 495     // define what recalib parameters are needed for various switches
 496     // this is based on implementation in AliEMCALRecoUtils
 497     Bool_t needRecoParam   = fReClusterize;
 498     Bool_t needBadChannels = fBadCellRemove   | fClusterBadChannelCheck | fRecalDistToBadChannels | fReClusterize;
 499     Bool_t needRecalib     = fCalibrateEnergy | fReClusterize;
 500     Bool_t needTimecalib   = fCalibrateTime   | fReClusterize;
 501     Bool_t needMisalign    = fRecalClusPos    | fReClusterize;
 502     Bool_t needClusterizer = fReClusterize;
 503
 504     // init bad channels
 505     if (needBadChannels) {
 506       Int_t fInitBC = InitBadChannels();
 507       if (fInitBC==0)
 508         AliError("InitBadChannels returned false, returning");
 509       if (fInitBC==1)
 510         AliWarning("InitBadChannels OK");
 511       if (fInitBC>1)
 512         AliWarning(Form("No external hot channel set: %d - %s", event->GetRunNumber(), fFilepass.Data()));
 513     }
 514
 515     // init recalibration factors
 516     if (needRecalib) {
 517       Int_t fInitRecalib = InitRecalib();
 518       if (fInitRecalib==0)
 519         AliError("InitRecalib returned false, returning");
 520       if (fInitRecalib==1)
 521         AliWarning("InitRecalib OK");
 522       if (fInitRecalib>1)
 523         AliWarning(Form("No recalibration available: %d - %s", event->GetRunNumber(), fFilepass.Data()));
 524
 525       Int_t fInitRunDepRecalib = InitRunDepRecalib();
 526       if (fInitRunDepRecalib==0)
 527         AliError("InitrunDepRecalib returned false, returning");
 528       if (fInitRunDepRecalib==1)
 529         AliWarning("InitRecalib OK");
 530       if (fInitRunDepRecalib>1)
 531         AliWarning(Form("No Temperature recalibration available: %d - %s", event->GetRunNumber(), fFilepass.Data()));
 532     }
 533
 534     // init time calibration
 535     if (needTimecalib) {
 536       Int_t initTC = InitTimeCalibration();
 537       if (!initTC)
 538         AliError("InitTimeCalibration returned false, returning");
 539       if (initTC==1) {
 540         fCalibrateTimeParamAvailable = kTRUE;
 541         AliWarning("InitTimeCalib OK");
 542       }
 543       if (initTC > 1)
 544         AliWarning(Form("No external time calibration available: %d - %s", event->GetRunNumber(), fFilepass.Data()));
 545     }
 546
 547     // init misalignment matrix
 548     if (needMisalign) {
 549       if (!InitMisalignMatrix())
 550         AliError("InitMisalignmentMatrix returned false, returning");
 551       else
 552         AliWarning("InitMisalignMatrix OK");
 553     }
 554
 555     // initiate reco params with some defaults
 556     // will not overwrite, if those have been provided by user
 557     if (needRecoParam) {
 558       Int_t initRC = InitRecParam();
 559
 560       if (initRC == 0)
 561         AliInfo("Defaults reco params loaded.");
 562       if (initRC > 1)
 563         AliWarning("User defined reco params.");
 564     }
 565
 566     // init clusterizer
 567     if (needClusterizer) {
 568       if (!InitClusterization())
 569         AliError("InitClusterization returned false, returning");
 570       else
 571         AliWarning("InitClusterization OK");
 572     }
 573
 574     if (fDebugLevel>1)
 575       fEMCALRecoUtils->Print("");
 576   }
 577
 578   // disable implied switches -------------------------------------------------
 579   // AliEMCALRecoUtils or clusterizer functions alredy include some
 580   // recalibration so based on those implied callibration te switches
 581   // here are disabled to avoid duplication
 582
 583   // clusterizer does cluster energy recalibration, position recomputation
 584   // and shower shape
 585   if (fReClusterize) {
 586     fReCalibCluster   = kFALSE;
 587     fRecalClusPos     = kFALSE;
 588     fRecalShowerShape = kFALSE;
 589   }
 590
 591   // CONFIGURE THE RECO UTILS -------------------------------------------------
 592   // configure the reco utils
 593   // this option does energy recalibration
 594   if (fCalibrateEnergy)
 595     fEMCALRecoUtils->SwitchOnRecalibration();
 596   else
 597     fEMCALRecoUtils->SwitchOffRecalibration();
 598
 599   // allows time calibration
 600   if (fCalibrateTime)
 601     fEMCALRecoUtils->SwitchOnTimeRecalibration();
 602   else
 603     fEMCALRecoUtils->SwitchOffTimeRecalibration();
 604
 605   // allows to zero bad cells
 606   if (fBadCellRemove)
 607     fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
 608   else
 609     fEMCALRecoUtils->SwitchOffBadChannelsRemoval();
 610
 611   // distance to bad channel recalibration
 612   if (fRecalDistToBadChannels)
 613     fEMCALRecoUtils->SwitchOnDistToBadChannelRecalculation();
 614   else
 615     fEMCALRecoUtils->SwitchOffDistToBadChannelRecalculation();
 616
 617   // exclude exotic cells
 618   if (fRejectExoticCells)
 619     fEMCALRecoUtils->SwitchOnRejectExoticCell();
 620   else
 621     fEMCALRecoUtils->SwitchOffRejectExoticCell();
 622
 623   // exclude clusters with exotic cells
 624   if (fRejectExoticClusters)
 625     fEMCALRecoUtils->SwitchOnRejectExoticCluster();
 626   else
 627     fEMCALRecoUtils->SwitchOffRejectExoticCluster();
 628
 629   // TODO: not implemented switches
 630   // SwitchOnClusterEnergySmearing
 631   // SwitchOnRunDepCorrection
 632
 633   // START PROCESSING ---------------------------------------------------------
 634   // Test if cells present
 635   AliVCaloCells *cells= event->GetEMCALCells();
 636   if (cells->GetNumberOfCells()<=0)
 637   {
 638     if (fDebugLevel>1)
 639       AliWarning(Form("Number of EMCAL cells = %d, returning", cells->GetNumberOfCells()));
 640     return;
 641   }
 642
 643   if (fDebugLevel>2)
 644     AliInfo(Form("Re-calibrate cluster %d\n",fReCalibCluster));
 645
 646   // mark the cells not recalibrated in case of selected
 647   // time, energy recalibration or bad channel removal
 648   if (fCalibrateEnergy || fCalibrateTime || fBadCellRemove)
 649     fEMCALRecoUtils->ResetCellsCalibrated();
 650
 651  // CELL RECALIBRATION -------------------------------------------------------
 652   // cell objects will be updated
 653   // the cell calibrations are also processed locally any time those are needed
 654   // in case that the cell objects are not to be updated here for later use
 655   if (fUpdateCell) {
 656     // do the update
 657     // includes exotic cell check in the UpdateCells function - is not provided
 658     // by the reco utils
 659     UpdateCells();
 660
 661     // switch off recalibrations so those are not done multiple times
 662     // this is just for safety, the recalibrated flag of cell object
 663     // should not allow for farther processing anyways
 664     fEMCALRecoUtils->SwitchOffRecalibration();
 665     fEMCALRecoUtils->SwitchOffTimeRecalibration();
 666
 667     if (fDoUpdateOnly)
 668       return;
 669   }
 670
 671   // RECLUSTERIZATION ---------------------------------------------------------
 672   if (fReClusterize)
 673   {
 674     FillDigitsArray();
 675     Clusterize();
 676     UpdateClusters();
 677   }
 678
 679   // Store good clusters
 680   TClonesArray *clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
 681   if (!clusArr)
 682     clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
 683   if (!clusArr) {
 684     AliWarning(Form("No cluster array, number of cells in event = %d, returning", cells->GetNumberOfCells()));
 685     return;
 686   }
 687
 688   // PROCESS SINGLE CLUSTER RECALIBRATION -------------------------------------
 689   // now go through clusters one by one and process separate correction
 690   // as those were defined or not
 691   Int_t nclusters = clusArr->GetEntriesFast();
 692   for (Int_t icluster=0; icluster < nclusters; ++icluster)
 693   {
 694     AliVCluster *clust = static_cast<AliVCluster*>(clusArr->At(icluster));
 695     if (!clust)
 696       continue;
 697     if  (!clust->IsEMCAL())
 698       continue;
 699
 700     // REMOVE CLUSTERS WITH BAD CELLS -----------------------------
 701     if (fClusterBadChannelCheck) {
 702       // careful, the the ClusterContainsBadChannel is dependent on
 703       // SwitchOnBadChannelsRemoval, switching it ON automatically
 704       // and returning to original value after processing
 705       Bool_t badRemoval = fEMCALRecoUtils->IsBadChannelsRemovalSwitchedOn();
 706       fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
 707
 708       Bool_t badResult = fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo, clust->GetCellsAbsId(), clust->GetNCells());
 709
 710       // switch the bad channels removal back
 711       if (!badRemoval)
 712         fEMCALRecoUtils->SwitchOffBadChannelsRemoval();
 713
 714       if (badResult)
 715       {
 716         delete clusArr->RemoveAt(icluster);
 717         continue; //TODO is it really needed to remove it? Or should we flag it?
 718       }
 719     }
 720
 721     // REMOVE EXOTIC CLUSTERS -------------------------------------
 722     // does process local cell recalibration energy and time without replacing
 723     // the global cell values, in case of no cell recalib done yet
 724     if (fRejectExoticClusters)
 725     {
 726       // careful, the IsExoticCluster is dependent on
 727       // SwitchOnRejectExoticCell, switching it ON automatically
 728       // and returning to original value after processing
 729       Bool_t exRemoval = fEMCALRecoUtils->IsRejectExoticCell();
 730       fEMCALRecoUtils->SwitchOnRejectExoticCell();
 731
 732       // get bunch crossing
 733       Int_t bunchCrossNo = event->GetBunchCrossNumber();
 734
 735       Bool_t exResult = fEMCALRecoUtils->IsExoticCluster(clust, cells, bunchCrossNo);
 736
 737       // switch the exotic channels removal back
 738       if (!exRemoval)
 739         fEMCALRecoUtils->SwitchOffRejectExoticCell();
 740
 741       if (exResult) {
 742         delete clusArr->RemoveAt(icluster);
 743         continue; //TODO is it really needed to remove it? Or should we flag it?
 744       }
 745     }
 746
 747     // FIDUCIAL CUT -----------------------------------------------
 748     if (fFiducial) {
 749       // depends on SetNumberOfCellsFromEMCALBorder
 750       // SwitchOnNoFiducialBorderInEMCALEta0
 751       if (!fEMCALRecoUtils->CheckCellFiducialRegion(fEMCALGeo, clust, cells)){
 752         delete clusArr->RemoveAt(icluster);
 753         continue; //TODO it would be nice to store the distance
 754       }
 755     }
 756
 757     // CLUSTER ENERGY ---------------------------------------------
 758     // does process local cell recalibration energy and time without replacing
 759     // the global cell values, in case of no cell recalib done yet
 760     if (fReCalibCluster) {
 761       fEMCALRecoUtils->RecalibrateClusterEnergy(fEMCALGeo, clust, cells);
 762       if (clust->E() < 1e-9) {
 763         delete clusArr->RemoveAt(icluster);
 764         continue;
 765       }
 766     }
 767
 768     // CLUSTER POSITION -------------------------------------------
 769     // does process local cell energy recalibration, if enabled and cells
 770     // not calibrated yet
 771     if (fRecalClusPos)
 772       fEMCALRecoUtils->RecalculateClusterPosition(fEMCALGeo, cells, clust);
 773
 774     // SHOWER SHAPE -----------------------------------------------
 775     if (fRecalShowerShape)
 776       fEMCALRecoUtils->RecalculateClusterShowerShapeParameters(fEMCALGeo, cells, clust);
 777
 778     // NONLINEARITY -----------------------------------------------
 779     if (fDoNonLinearity) {
 780       Float_t correctedEnergy = fEMCALRecoUtils->CorrectClusterEnergyLinearity(clust);
 781       clust->SetE(correctedEnergy);
 782     }
 783
 784     // DISTANCE TO BAD CHANNELS -----------------------------------
 785     if (fRecalDistToBadChannels)
 786       fEMCALRecoUtils->RecalculateClusterDistanceToBadChannel(fEMCALGeo, cells, clust);
 787   }
 788
 789   clusArr->Compress();
 790
 791   if (!fDoTrackMatch)
 792     return;
 793
 794   // TRACK MATCHING -----------------------------------------------------------
 795   if (!TGeoGlobalMagField::Instance()->GetField()) {
 796     event->InitMagneticField();
 797   }
 798
 799   fEMCALRecoUtils->FindMatches(event,0x0,fEMCALGeo);
 800   fEMCALRecoUtils->SetClusterMatchedToTrack(event);
 801   fEMCALRecoUtils->SetTracksMatchedToCluster(event);
 802 }
 803
 804 //_____________________________________________________
 805 Bool_t AliEMCALTenderSupply::InitMisalignMatrix()
 806 {
 807   // Initialising misalignment matrices
 808
 809   AliVEvent *event = GetEvent();
 810
 811   if (!event)
 812     return kFALSE;
 813
 814   if (fGeomMatrixSet)
 815   {
 816     AliInfo("Misalignment matrix already set");
 817     return kTRUE;
 818   }
 819
 820   if (fDebugLevel>0)
 821     AliInfo("Initialising misalignment matrix");
 822
 823   if (fLoadGeomMatrices) {
 824     for(Int_t mod=0; mod < fEMCALGeo->GetNumberOfSuperModules(); ++mod)
 825     {
 826       if (fEMCALMatrix[mod]){
 827         if (fDebugLevel > 2)
 828           fEMCALMatrix[mod]->Print();
 829         fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
 830       }
 831     }
 832     fGeomMatrixSet = kTRUE;
 833     return kTRUE;
 834   }
 835
 836   Int_t runGM = event->GetRunNumber();
 837   TObjArray *mobj = 0;
 838
 839   if (fMisalignSurvey == kdefault)
 840   { //take default alignment corresponding to run no
 841     AliOADBContainer emcalgeoCont(Form("emcal"));
 842     emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 843     mobj=(TObjArray*)emcalgeoCont.GetObject(runGM,"EmcalMatrices");
 844   }
 845
 846   if (fMisalignSurvey == kSurveybyS)
 847   { //take alignment at sector level
 848     if (runGM <= 140000) { //2010 data
 849       AliOADBContainer emcalgeoCont(Form("emcal2010"));
 850       emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 851       mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
 852     }
 853     else if (runGM>140000)
 854     { // 2011 LHC11a pass1 data
 855       AliOADBContainer emcalgeoCont(Form("emcal2011"));
 856       emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 857       mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byS");
 858     }
 859   }
 860
 861   if (fMisalignSurvey == kSurveybyM)
 862   { //take alignment at module level
 863     if (runGM <= 140000) { //2010 data
 864       AliOADBContainer emcalgeoCont(Form("emcal2010"));
 865       emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 866       mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
 867     }
 868     else if (runGM>140000)
 869     { // 2011 LHC11a pass1 data
 870       AliOADBContainer emcalgeoCont(Form("emcal2011"));
 871       emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
 872       mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byM");
 873     }
 874   }
 875
 876   if (!mobj) {
 877     AliFatal("Geometry matrix array not found");
 878     return kFALSE;
 879   }
 880
 881   for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
 882   {
 883     fEMCALMatrix[mod] = (TGeoHMatrix*) mobj->At(mod);
 884     fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
 885     fEMCALMatrix[mod]->Print();
 886   }
 887
 888   return kTRUE;
 889 }
 890
 891 //_____________________________________________________
 892 Int_t AliEMCALTenderSupply::InitBadChannels()
 893 {
 894   // Initialising bad channel maps
 895
 896   AliVEvent *event = GetEvent();
 897
 898   if (!event)
 899     return 0;
 900
 901   if (fDebugLevel>0)
 902     AliInfo("Initialising Bad channel map");
 903
 904   // init default maps first
 905   if (!fEMCALRecoUtils->GetEMCALBadChannelStatusMapArray())
 906     fEMCALRecoUtils->InitEMCALBadChannelStatusMap() ;
 907
 908   Int_t runBC = event->GetRunNumber();
 909
 910   AliOADBContainer *contBC = new AliOADBContainer("");
 911   if (fBasePath!="")
 912   { //if fBasePath specified in the ->SetBasePath()
 913     if (fDebugLevel>0) AliInfo(Form("Loading Bad Channels OADB from given path %s",fBasePath.Data()));
 914
 915     TFile *fbad=new TFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"read");
 916     if (!fbad || fbad->IsZombie())
 917     {
 918       AliFatal(Form("EMCALBadChannels.root was not found in the path provided: %s",fBasePath.Data()));
 919       return 0;
 920     }
 921
 922     if (fbad) delete fbad;
 923
 924     contBC->InitFromFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"AliEMCALBadChannels");
 925   }
 926   else
 927   { // Else choose the one in the $ALICE_ROOT directory
 928     if (fDebugLevel>0) AliInfo("Loading Bad Channels OADB from $ALICE_ROOT/OADB/EMCAL");
 929
 930     TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","read");
 931     if (!fbad || fbad->IsZombie())
 932     {
 933       AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root was not found");
 934       return 0;
 935     }
 936
 937     if (fbad) delete fbad;
 938
 939     contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","AliEMCALBadChannels");
 940   }
 941
 942   TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
 943   if (!arrayBC)
 944   {
 945     AliError(Form("No external hot channel set for run number: %d", runBC));
 946     return 2;
 947   }
 948
 949   Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
 950   for (Int_t i=0; i<sms; ++i)
 951   {
 952     TH2I *h = fEMCALRecoUtils->GetEMCALChannelStatusMap(i);
 953     if (h)
 954       delete h;
 955     h=(TH2I*)arrayBC->FindObject(Form("EMCALBadChannelMap_Mod%d",i));
 956
 957     if (!h)
 958     {
 959       AliError(Form("Can not get EMCALBadChannelMap_Mod%d",i));
 960       continue;
 961     }
 962     h->SetDirectory(0);
 963     fEMCALRecoUtils->SetEMCALChannelStatusMap(i,h);
 964   }
 965   return 1;
 966 }
 967
 968 //_____________________________________________________
 969 Int_t AliEMCALTenderSupply::InitRecalib()
 970 {
 971   // Initialising recalibration factors.
 972
 973   AliVEvent *event = GetEvent();
 974
 975   if (!event)
 976     return 0;
 977
 978   if (fDebugLevel>0)
 979     AliInfo("Initialising recalibration factors");
 980
 981   // init default maps first
 982   if (!fEMCALRecoUtils->GetEMCALRecalibrationFactorsArray())
 983     fEMCALRecoUtils->InitEMCALRecalibrationFactors() ;
 984
 985   Int_t runRC = event->GetRunNumber();
 986
 987   AliOADBContainer *contRF=new AliOADBContainer("");
 988   if (fBasePath!="")
 989   { //if fBasePath specified in the ->SetBasePath()
 990     if (fDebugLevel>0)  AliInfo(Form("Loading Recalib OADB from given path %s",fBasePath.Data()));
 991
 992     TFile *fRecalib= new TFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"read");
 993     if (!fRecalib || fRecalib->IsZombie())
 994     {
 995       AliFatal(Form("EMCALRecalib.root not found in %s",fBasePath.Data()));
 996       return 0;
 997     }
 998
 999     if (fRecalib) delete fRecalib;
1000
1001     contRF->InitFromFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"AliEMCALRecalib");
1002   }
1003   else
1004   { // Else choose the one in the $ALICE_ROOT directory
1005     if (fDebugLevel>0)  AliInfo("Loading Recalib OADB from $ALICE_ROOT/OADB/EMCAL");
1006
1007     TFile *fRecalib= new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","read");
1008     if (!fRecalib || fRecalib->IsZombie())
1009     {
1010       AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root was not found");
1011       return 0;
1012     }
1013
1014     if (fRecalib) delete fRecalib;
1015
1016     contRF->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","AliEMCALRecalib");
1017   }
1018
1019   TObjArray *recal=(TObjArray*)contRF->GetObject(runRC);
1020   if (!recal)
1021   {
1022     AliError(Form("No Objects for run: %d",runRC));
1023     return 2;
1024   }
1025
1026   TObjArray *recalpass=(TObjArray*)recal->FindObject(fFilepass);
1027   if (!recalpass)
1028   {
1029     AliError(Form("No Objects for run: %d - %s",runRC,fFilepass.Data()));
1030     return 2;
1031   }
1032
1033   TObjArray *recalib=(TObjArray*)recalpass->FindObject("Recalib");
1034   if (!recalib)
1035   {
1036     AliError(Form("No Recalib histos found for  %d - %s",runRC,fFilepass.Data()));
1037     return 2;
1038   }
1039
1040   if (fDebugLevel>0) recalib->Print();
1041
1042   Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
1043   for (Int_t i=0; i<sms; ++i)
1044   {
1045     TH2F *h = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactors(i);
1046     if (h)
1047       delete h;
1048     h = (TH2F*)recalib->FindObject(Form("EMCALRecalFactors_SM%d",i));
1049     if (!h)
1050     {
1051       AliError(Form("Could not load EMCALRecalFactors_SM%d",i));
1052       continue;
1053     }
1054     h->SetDirectory(0);
1055     fEMCALRecoUtils->SetEMCALChannelRecalibrationFactors(i,h);
1056   }
1057   return 1;
1058 }
1059
1060 //_____________________________________________________
1061 Int_t AliEMCALTenderSupply::InitRunDepRecalib()
1062 {
1063   // Initialising recalibration factors.
1064
1065   AliVEvent *event = GetEvent();
1066
1067   if (!event)
1068     return 0;
1069
1070   if (fDebugLevel>0)
1071     AliInfo("Initialising recalibration factors");
1072
1073   // init default maps first
1074   if (!fEMCALRecoUtils->GetEMCALRecalibrationFactorsArray())
1075     fEMCALRecoUtils->InitEMCALRecalibrationFactors() ;
1076
1077   Int_t runRC = event->GetRunNumber();
1078
1079   AliOADBContainer *contRF=new AliOADBContainer("");
1080   if (fBasePath!="")
1081   { //if fBasePath specified in the ->SetBasePath()
1082     if (fDebugLevel>0)  AliInfo(Form("Loading Recalib OADB from given path %s",fBasePath.Data()));
1083
1084     TFile *fRunDepRecalib= new TFile(Form("%s/EMCALTemperatureCorrCalib.root",fBasePath.Data()),"read");
1085     if (!fRunDepRecalib || fRunDepRecalib->IsZombie())
1086     {
1087       AliFatal(Form("EMCALTemperatureCorrCalib.root not found in %s",fBasePath.Data()));
1088       return 0;
1089     }
1090
1091     if (fRunDepRecalib) delete fRunDepRecalib;
1092
1093     contRF->InitFromFile(Form("%s/EMCALTemperatureCorrCalib.root",fBasePath.Data()),"AliEMCALRunDepTempCalibCorrections");
1094   }
1095   else
1096   { // Else choose the one in the $ALICE_ROOT directory
1097     if (fDebugLevel>0)  AliInfo("Loading Recalib OADB from $ALICE_ROOT/OADB/EMCAL");
1098
1099     TFile *fRunDepRecalib= new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALTemperatureCorrCalib.root","read");
1100     if (!fRunDepRecalib || fRunDepRecalib->IsZombie())
1101     {
1102       AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALTemperatureCorrCalib.root was not found");
1103       return 0;
1104     }
1105
1106     if (fRunDepRecalib) delete fRunDepRecalib;
1107
1108     contRF->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALTemperatureCorrCalib.root","AliEMCALRunDepTempCalibCorrections");
1109   }
1110
1111   TH1S *rundeprecal=(TH1S*)contRF->GetObject(runRC);
1112   if (!rundeprecal)
1113   {
1114     AliWarning(Form("No TemperatureCorrCalib Objects for run: %d",runRC));
1115     // let's get the closest runnumber instead then..
1116     Int_t lower = 0;
1117     Int_t ic = 0;
1118     Int_t maxEntry = contRF->GetNumberOfEntries();
1119
1120     while ((ic < maxEntry) && (contRF->UpperLimit(ic) < runRC)) {
1121       lower = ic;
1122       ic++;
1123     }
1124
1125     Int_t closest = lower;
1126     if ((ic<maxEntry) &&
1127          (contRF->LowerLimit(ic)-runRC) < (runRC - contRF->UpperLimit(lower))) {
1128          closest = ic;
1129     }
1130
1131     AliWarning(Form("TemperatureCorrCalib Objects found closest id %d from run: %d", closest, contRF->LowerLimit(closest)));
1132     rundeprecal = (TH1S*) contRF->GetObjectByIndex(closest);
1133   }
1134
1135   if (fDebugLevel>0) rundeprecal->Print();
1136
1137   Int_t nSM = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
1138
1139   for (Int_t ism=0; ism<nSM; ++ism)
1140   {
1141     for (Int_t icol=0; icol<48; ++icol)
1142     {
1143       for (Int_t irow=0; irow<24; ++irow)
1144       {
1145         Float_t factor = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactor(ism,icol,irow);
1146
1147         Int_t absID = fEMCALGeo->GetAbsCellIdFromCellIndexes(ism, irow, icol); // original calibration factor
1148         factor *= rundeprecal->GetBinContent(absID) / 10000. ; // correction dependent on T
1149
1150         fEMCALRecoUtils->SetEMCALChannelRecalibrationFactor(ism,icol,irow,factor);
1151       } // columns
1152     } // rows
1153   } // SM loop
1154
1155   return 1;
1156 }
1157
1158
1159 //_____________________________________________________
1160 Int_t AliEMCALTenderSupply::InitTimeCalibration()
1161 {
1162   // Initialising bad channel maps
1163   AliVEvent *event = GetEvent();
1164
1165   if (!event)
1166     return 0;
1167
1168   if (fDebugLevel>0)
1169     AliInfo("Initialising time calibration map");
1170
1171   // init default maps first
1172   if (!fEMCALRecoUtils->GetEMCALTimeRecalibrationFactorsArray())
1173     fEMCALRecoUtils->InitEMCALTimeRecalibrationFactors() ;
1174
1175   Int_t runBC = event->GetRunNumber();
1176
1177   AliOADBContainer *contBC = new AliOADBContainer("");
1178   if (fBasePath!="")
1179   { //if fBasePath specified in the ->SetBasePath()
1180     if (fDebugLevel>0) AliInfo(Form("Loading time calibration OADB from given path %s",fBasePath.Data()));
1181
1182     TFile *fbad=new TFile(Form("%s/EMCALTimeCalib.root",fBasePath.Data()),"read");
1183     if (!fbad || fbad->IsZombie())
1184     {
1185       AliFatal(Form("EMCALTimeCalib.root was not found in the path provided: %s",fBasePath.Data()));
1186       return 0;
1187     }
1188
1189     if (fbad) delete fbad;
1190
1191     contBC->InitFromFile(Form("%s/EMCALTimeCalib.root",fBasePath.Data()),"AliEMCALTimeCalib");
1192   }
1193   else
1194   { // Else choose the one in the $ALICE_ROOT directory
1195     if (fDebugLevel>0) AliInfo("Loading time calibration OADB from $ALICE_ROOT/OADB/EMCAL");
1196
1197     TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root","read");
1198     if (!fbad || fbad->IsZombie())
1199     {
1200       AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root was not found");
1201       return 0;
1202     }
1203
1204     if (fbad) delete fbad;
1205
1206     contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root","AliEMCALTimeCalib");
1207   }
1208
1209   TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
1210   if (!arrayBC)
1211   {
1212     AliError(Form("No external time calibration set for run number: %d", runBC));
1213     return 2;
1214   }
1215
1216   // Here, it looks for a specific pass
1217   TString pass = fFilepass;
1218   if (fFilepass=="calo_spc") pass ="pass1";
1219   TObjArray *arrayBCpass=(TObjArray*)arrayBC->FindObject(pass);
1220   if (!arrayBCpass)
1221   {
1222     AliError(Form("No external time calibration set for: %d -%s", runBC,pass.Data()));
1223     return 2;
1224   }
1225
1226   if (fDebugLevel>0) arrayBCpass->Print();
1227
1228   for(Int_t i = 0; i < 4; i++)
1229   {
1230     TH1F *h = fEMCALRecoUtils->GetEMCALChannelTimeRecalibrationFactors(i);
1231     if (h)
1232       delete h;
1233
1234     h = (TH1F*)arrayBCpass->FindObject(Form("hAllTimeAvBC%d",i));
1235
1236     if (!h)
1237     {
1238       AliError(Form("Can not get hAllTimeAvBC%d",i));
1239       continue;
1240     }
1241     h->SetDirectory(0);
1242     fEMCALRecoUtils->SetEMCALChannelTimeRecalibrationFactors(i,h);
1243   }
1244   return 1;
1245 }
1246
1247 //_____________________________________________________
1248 void AliEMCALTenderSupply::UpdateCells()
1249 {
1250   //Remove bad cells from the cell list
1251   //Recalibrate energy and time cells
1252   //This is required for later reclusterization
1253
1254   AliVEvent *event = GetEvent();
1255
1256   if (!event) return ;
1257
1258   AliVCaloCells *cells = event->GetEMCALCells();
1259   Int_t bunchCrossNo = event->GetBunchCrossNumber();
1260
1261   fEMCALRecoUtils->RecalibrateCells(cells, bunchCrossNo);
1262
1263   // remove exotic cells - loop through cells and zero the exotic ones
1264   // just like with bad cell rejection in reco utils (inside RecalibrateCells)
1265   if (fRejectExoticCells)
1266   {
1267     Short_t  absId  =-1;
1268     Double_t ecell = 0;
1269     Double_t tcell = 0;
1270     Double_t efrac = 0;
1271     Int_t  mclabel = -1;
1272     Bool_t   isExot = kFALSE;
1273
1274     // loop through cells
1275     Int_t nEMcell  = cells->GetNumberOfCells() ;
1276     for (Int_t iCell = 0; iCell < nEMcell; iCell++)
1277     {
1278       cells->GetCell(iCell, absId, ecell, tcell, mclabel, efrac);
1279
1280       isExot = fEMCALRecoUtils->IsExoticCell(absId, cells, bunchCrossNo);
1281       // zero if exotic
1282       if (isExot)
1283         cells->SetCell(iCell, absId, 0.0, -1.0, mclabel, efrac);
1284     } // cell loop
1285   } // reject exotic cells
1286
1287   cells->Sort();
1288 }
1289
1290 //_____________________________________________________
1291 TString AliEMCALTenderSupply::GetBeamType()
1292 {
1293   // Get beam type : pp-AA-pA
1294   // ESDs have it directly, AODs get it from hardcoded run number ranges
1295
1296   AliVEvent *event = GetEvent();
1297
1298   if (!event) {
1299     AliError("Couldn't retrieve event!");
1300     return "";
1301   }
1302
1303   TString beamType;
1304
1305   AliESDEvent *esd = dynamic_cast<AliESDEvent*>(event);
1306   if (esd) {
1307     const AliESDRun *run = esd->GetESDRun();
1308     beamType = run->GetBeamType();
1309   }
1310   else
1311   {
1312     Int_t runNumber = event->GetRunNumber();
1313     if ((runNumber >= 136851 && runNumber <= 139517)  // LHC10h
1314   || (runNumber >= 166529 && runNumber <= 170593))  // LHC11h
1315     {
1316       beamType = "A-A";
1317     }
1318     else
1319     {
1320       beamType = "p-p";
1321     }
1322   }
1323
1324   return beamType;
1325 }
1326
1327 //_____________________________________________________
1328 Int_t AliEMCALTenderSupply::InitRecParam()
1329 {
1330   // Init reco params if not yet exist (probably shipped by the user already)
1331
1332   if (fRecParam != 0)
1333     return 2;
1334
1335   TString beamType = GetBeamType();
1336
1337   // set some default reco params
1338   fRecParam = new AliEMCALRecParam();
1339   fRecParam->SetClusteringThreshold(0.100);
1340   fRecParam->SetMinECut(0.050);
1341
1342   if (!fCalibrateTimeParamAvailable) {
1343     fRecParam->SetTimeCut(250*1.e-9);
1344     fRecParam->SetTimeMin(425*1.e-9);
1345     fRecParam->SetTimeMax(825*1.e-9);
1346   } else {
1347     fRecParam->SetTimeCut(100*1.e-9);
1348     fRecParam->SetTimeMin(-50*1.e-9);
1349     fRecParam->SetTimeMax(50*1.e-9);
1350   }
1351
1352   if (beamType == "A-A") {
1353     fRecParam->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv2);
1354   } else {
1355     fRecParam->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv1);
1356   }
1357
1358   return 0;
1359 }
1360
1361 //_____________________________________________________
1362 Bool_t AliEMCALTenderSupply::InitClusterization()
1363 {
1364   // Initialing clusterization/unfolding algorithm and set all the needed parameters.
1365
1366   AliVEvent *event = GetEvent();
1367
1368   if (!event)
1369     return kFALSE;
1370
1371   if (fDebugLevel>0)
1372     AliInfo(Form("Initialising reclustering parameters: Clusterizer type: %d",fRecParam->GetClusterizerFlag()));
1373
1374   //---setup clusterizer
1375   if (fClusterizer) {
1376     // avoid deleting fDigitsArr
1377     fClusterizer->SetDigitsArr(0);
1378     delete fClusterizer;
1379     fClusterizer = 0;
1380   }
1381   if     (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
1382     fClusterizer = new AliEMCALClusterizerv1 (fEMCALGeo);
1383   else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
1384     fClusterizer = new AliEMCALClusterizerv2(fEMCALGeo);
1385   else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN)
1386   {
1387     AliEMCALClusterizerNxN *clusterizer = new AliEMCALClusterizerNxN(fEMCALGeo);
1388     clusterizer->SetNRowDiff(fRecParam->GetNRowDiff());
1389     clusterizer->SetNColDiff(fRecParam->GetNColDiff());
1390     fClusterizer = clusterizer;
1391   }
1392   else
1393   {
1394     AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
1395     return kFALSE;
1396   }
1397
1398   // Set the clustering parameters
1399   fClusterizer->SetECAClusteringThreshold(fRecParam->GetClusteringThreshold());
1400   fClusterizer->SetECALogWeight          (fRecParam->GetW0()                 );
1401   fClusterizer->SetMinECut               (fRecParam->GetMinECut()            );
1402   fClusterizer->SetUnfolding             (fRecParam->GetUnfold()             );
1403   fClusterizer->SetECALocalMaxCut        (fRecParam->GetLocMaxCut()          );
1404   fClusterizer->SetTimeCut               (fRecParam->GetTimeCut()            );
1405   fClusterizer->SetTimeMin               (fRecParam->GetTimeMin()            );
1406   fClusterizer->SetTimeMax               (fRecParam->GetTimeMax()            );
1407   fClusterizer->SetInputCalibrated       (kTRUE                              );
1408   fClusterizer->SetJustClusters          (kTRUE                              );
1409
1410   // In case of unfolding after clusterization is requested, set the corresponding parameters
1411   if (fRecParam->GetUnfold())
1412   {
1413     for (Int_t i = 0; i < 8; ++i)
1414     {
1415       fClusterizer->SetSSPars(i, fRecParam->GetSSPars(i));
1416     }
1417     for (Int_t i = 0; i < 3; ++i)
1418     {
1419       fClusterizer->SetPar5  (i, fRecParam->GetPar5(i));
1420       fClusterizer->SetPar6  (i, fRecParam->GetPar6(i));
1421     }
1422     fClusterizer->InitClusterUnfolding();
1423   }
1424
1425   fClusterizer->SetDigitsArr(fDigitsArr);
1426   fClusterizer->SetOutput(0);
1427   fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
1428   return kTRUE;
1429 }
1430
1431 //_____________________________________________________
1432 void AliEMCALTenderSupply::FillDigitsArray()
1433 {
1434   // Fill digits from cells to a TClonesArray.
1435
1436   AliVEvent *event = GetEvent();
1437
1438  if (!event)
1439     return;
1440
1441   // In case of MC productions done before aliroot tag v5-02-Rev09
1442   // assing the cluster label to all the cells belonging to this cluster
1443   // very rough
1444   Int_t cellLabels[12672];
1445   if (fSetCellMCLabelFromCluster)
1446   {
1447     for (Int_t i = 0; i < 12672; i++)
1448     {
1449       cellLabels       [i] = 0 ;
1450       fOrgClusterCellId[i] =-1 ;
1451     }
1452
1453     Int_t nClusters = event->GetNumberOfCaloClusters();
1454     for (Int_t i = 0; i < nClusters; i++)
1455     {
1456       AliVCluster *clus =  event->GetCaloCluster(i);
1457
1458       if (!clus) continue;
1459
1460       if (!clus->IsEMCAL()) continue ;
1461
1462       Int_t      label = clus->GetLabel();
1463       UShort_t * index = clus->GetCellsAbsId() ;
1464
1465       for(Int_t icell=0; icell < clus->GetNCells(); icell++)
1466       {
1467         cellLabels       [index[icell]] = label;
1468         fOrgClusterCellId[index[icell]] = i ; // index of the original cluster
1469       }
1470     }// cluster loop
1471   }
1472
1473   fDigitsArr->Clear("C");
1474   AliVCaloCells *cells = event->GetEMCALCells();
1475   Int_t ncells = cells->GetNumberOfCells();
1476   for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell)
1477   {
1478     Double_t cellAmplitude=0, cellTime=0, efrac = 0;
1479     Short_t  cellNumber=0;
1480     Int_t mcLabel=-1;
1481
1482     if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime, mcLabel, efrac) != kTRUE)
1483       break;
1484
1485     // Do not add if energy already too low (some cells set to 0 if bad channels)
1486     if (cellAmplitude < fRecParam->GetMinECut())
1487       continue;
1488
1489     // If requested, do not include exotic cells
1490    if (fEMCALRecoUtils->IsExoticCell(cellNumber,cells,event->GetBunchCrossNumber()))
1491       continue;
1492
1493     if      (fSetCellMCLabelFromCluster) mcLabel = cellLabels[cellNumber];
1494     else if (fRemapMCLabelForAODs     ) RemapMCLabelForAODs(mcLabel);
1495
1496     if (mcLabel > 0 && efrac < 1e-6) efrac = 1;
1497
1498     new((*fDigitsArr)[idigit]) AliEMCALDigit(mcLabel, mcLabel, cellNumber,
1499                                              (Float_t)cellAmplitude, (Float_t)cellTime,
1500                                              AliEMCALDigit::kHG,idigit, 0, 0, efrac*cellAmplitude);
1501     idigit++;
1502   }
1503 }
1504
1505 //_____________________________________________________
1506 void AliEMCALTenderSupply::Clusterize()
1507 {
1508   // Clusterize.
1509
1510   fClusterizer->Digits2Clusters("");
1511 }
1512
1513 //_____________________________________________________
1514 void AliEMCALTenderSupply::UpdateClusters()
1515 {
1516   // Update ESD cluster list.
1517
1518   AliVEvent *event = GetEvent();
1519
1520   if (!event)
1521     return;
1522
1523   TClonesArray *clus = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
1524   if (!clus)
1525     clus = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
1526   if (!clus)
1527   {
1528     AliError(" Null pointer to calo clusters array, returning");
1529     return;
1530   }
1531
1532   // Before destroying the orignal list, assign to the rec points the MC labels
1533   // of the original clusters, if requested
1534   if (fSetCellMCLabelFromCluster == 2)
1535     SetClustersMCLabelFromOriginalClusters() ;
1536
1537   Int_t nents = clus->GetEntriesFast();
1538   for (Int_t i=0; i < nents; ++i)
1539   {
1540     AliVCluster *c = dynamic_cast<AliVCluster*>(clus->At(i));
1541     if (!c)
1542       continue;
1543     if (c->IsEMCAL())
1544     {
1545       delete clus->RemoveAt(i);
1546     }
1547   }
1548
1549   clus->Compress();
1550
1551   RecPoints2Clusters(clus);
1552 }
1553
1554 //_____________________________________________________
1555 void AliEMCALTenderSupply::RecPoints2Clusters(TClonesArray *clus)
1556 {
1557   // Convert AliEMCALRecoPoints to AliESDCaloClusters/AliAODCaloClusters.
1558   // Cluster energy, global position, cells and their amplitude fractions are restored.
1559
1560   AliVEvent *event = GetEvent();
1561
1562   if (!event)
1563     return;
1564
1565   Int_t ncls = fClusterArr->GetEntriesFast();
1566   for(Int_t i=0, nout=clus->GetEntriesFast(); i < ncls; ++i)
1567   {
1568     AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
1569
1570     Int_t ncellsTrue = 0;
1571     const Int_t ncells = recpoint->GetMultiplicity();
1572     UShort_t   absIds[ncells];
1573     Double32_t ratios[ncells];
1574     Int_t *dlist = recpoint->GetDigitsList();
1575     Float_t *elist = recpoint->GetEnergiesList();
1576     for (Int_t c = 0; c < ncells; ++c)
1577     {
1578       AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
1579       absIds[ncellsTrue] = digit->GetId();
1580       ratios[ncellsTrue] = elist[c]/digit->GetAmplitude();
1581       if (ratios[ncellsTrue] < 0.001)
1582         continue;
1583       ++ncellsTrue;
1584     }
1585
1586     if (ncellsTrue < 1)
1587     {
1588       AliWarning("Skipping cluster with no cells");
1589       continue;
1590     }
1591
1592     // calculate new cluster position
1593     TVector3 gpos;
1594     recpoint->GetGlobalPosition(gpos);
1595     Float_t g[3];
1596     gpos.GetXYZ(g);
1597
1598     AliVCluster *c = static_cast<AliVCluster*>(clus->New(nout++));
1599     c->SetID(nout-1);
1600     c->SetType(AliVCluster::kEMCALClusterv1);
1601     c->SetE(recpoint->GetEnergy());
1602     c->SetPosition(g);
1603     c->SetNCells(ncellsTrue);
1604     c->SetDispersion(recpoint->GetDispersion());
1605     c->SetEmcCpvDistance(-1);            //not yet implemented
1606     c->SetChi2(-1);                      //not yet implemented
1607     c->SetTOF(recpoint->GetTime()) ;     //time-of-flight
1608     c->SetNExMax(recpoint->GetNExMax()); //number of local maxima
1609     Float_t elipAxis[2];
1610     recpoint->GetElipsAxis(elipAxis);
1611     c->SetM02(elipAxis[0]*elipAxis[0]) ;
1612     c->SetM20(elipAxis[1]*elipAxis[1]) ;
1613     c->SetCellsAbsId(absIds);
1614     c->SetCellsAmplitudeFraction(ratios);
1615
1616     //MC labels
1617     Int_t  parentMult = 0;
1618     Int_t *parentList = recpoint->GetParents(parentMult);
1619     if (parentMult > 0) c->SetLabel(parentList, parentMult);
1620
1621   }
1622 }
1623
1624 //___________________________________________________________
1625 void AliEMCALTenderSupply::RemapMCLabelForAODs(Int_t & label)
1626 {
1627   // MC label for Cells not remapped after ESD filtering, do it here.
1628
1629   if (label < 0) return;
1630
1631   AliAODEvent  * evt = dynamic_cast<AliAODEvent*> (GetEvent()) ;
1632   if (!evt) return ;
1633
1634   TClonesArray * arr = dynamic_cast<TClonesArray*>(evt->FindListObject("mcparticles")) ;
1635   if (!arr) return ;
1636
1637   if (label < arr->GetEntriesFast())
1638   {
1639     AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(label));
1640     if (!particle) return ;
1641
1642     if (label == particle->Label()) return ; // label already OK
1643   }
1644
1645   // loop on the particles list and check if there is one with the same label
1646   for (Int_t ind = 0; ind < arr->GetEntriesFast(); ind++)
1647   {
1648     AliAODMCParticle * particle = dynamic_cast<AliAODMCParticle *>(arr->At(ind));
1649     if (!particle) continue ;
1650
1651     if (label == particle->Label())
1652     {
1653       label = ind;
1654       return;
1655     }
1656   }
1657
1658   label = -1;
1659 }
1660
1661 //_____________________________________________________________________________________________
1662 void AliEMCALTenderSupply::SetClustersMCLabelFromOriginalClusters()
1663 {
1664   // Get the original clusters that contribute to the new rec point cluster,
1665   // assign the labels of such clusters to the new rec point cluster.
1666   // Only approximatedly valid  when output are V1 clusters, or input/output clusterizer
1667   // are the same handle with care
1668   // Copy from same method in AliAnalysisTaskEMCALClusterize, but here modify the recpoint and
1669   // not the output calocluster
1670
1671   Int_t ncls = fClusterArr->GetEntriesFast();
1672   for(Int_t irp=0; irp < ncls; ++irp)
1673   {
1674     TArrayI clArray(300) ; //Weird if more than a few clusters are in the origin ...
1675     clArray.Reset();
1676     Int_t nClu = 0;
1677     Int_t nLabTotOrg = 0;
1678     Float_t emax = -1;
1679     Int_t idMax = -1;
1680
1681     AliEMCALRecPoint *clus = static_cast<AliEMCALRecPoint*>(fClusterArr->At(irp));
1682
1683     //Find the clusters that originally had the cells
1684     const Int_t ncells = clus->GetMultiplicity();
1685     Int_t *digList     = clus->GetDigitsList();
1686
1687     for (Int_t iLoopCell = 0 ; iLoopCell < ncells ; iLoopCell++)
1688     {
1689       AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(digList[iLoopCell]));
1690       Int_t idCell = digit->GetId();
1691
1692       if (idCell>=0)
1693       {
1694         Int_t idCluster = fOrgClusterCellId[idCell];
1695         Bool_t set = kTRUE;
1696         for (Int_t icl =0; icl < nClu; icl++)
1697         {
1698           if (((Int_t)clArray.GetAt(icl))==-1) continue;
1699           if (idCluster == ((Int_t)clArray.GetAt(icl))) set = kFALSE;
1700         }
1701         if (set && idCluster >= 0)
1702         {
1703           clArray.SetAt(idCluster,nClu++);
1704           nLabTotOrg+=(GetEvent()->GetCaloCluster(idCluster))->GetNLabels();
1705
1706           //Search highest E cluster
1707           AliVCluster * clOrg = GetEvent()->GetCaloCluster(idCluster);
1708           if (emax < clOrg->E())
1709           {
1710             emax  = clOrg->E();
1711             idMax = idCluster;
1712           }
1713         }
1714       }
1715     }// cell loop
1716
1717     // Put the first in the list the cluster with highest energy
1718     if (idMax != ((Int_t)clArray.GetAt(0))) // Max not at first position
1719     {
1720       Int_t maxIndex = -1;
1721       Int_t firstCluster = ((Int_t)clArray.GetAt(0));
1722       for (Int_t iLoopCluster = 0 ; iLoopCluster < nClu ; iLoopCluster++)
1723       {
1724         if (idMax == ((Int_t)clArray.GetAt(iLoopCluster))) maxIndex = iLoopCluster;
1725       }
1726
1727       if (firstCluster >=0 && idMax >=0)
1728       {
1729         clArray.SetAt(idMax,0);
1730         clArray.SetAt(firstCluster,maxIndex);
1731       }
1732     }
1733
1734     // Get the labels list in the original clusters, assign all to the new cluster
1735     TArrayI clMCArray(nLabTotOrg) ;
1736     clMCArray.Reset();
1737
1738     Int_t nLabTot = 0;
1739     for (Int_t iLoopCluster = 0 ; iLoopCluster < nClu ; iLoopCluster++)
1740     {
1741       Int_t idCluster = (Int_t) clArray.GetAt(iLoopCluster);
1742       AliVCluster * clOrg = GetEvent()->GetCaloCluster(idCluster);
1743       Int_t nLab = clOrg->GetNLabels();
1744
1745       for (Int_t iLab = 0 ; iLab < nLab ; iLab++)
1746       {
1747         Int_t lab = clOrg->GetLabelAt(iLab) ;
1748         if (lab>=0)
1749         {
1750           Bool_t set = kTRUE;
1751           for(Int_t iLabTot =0; iLabTot < nLabTot; iLabTot++)
1752           {
1753             if (lab == ((Int_t)clMCArray.GetAt(iLabTot))) set = kFALSE;
1754           }
1755           if (set) clMCArray.SetAt(lab,nLabTot++);
1756         }
1757       }
1758     }// cluster loop
1759
1760     // Set the final list of labels to rec point
1761
1762     Int_t *labels = new Int_t[nLabTot];
1763     for(Int_t il = 0; il < nLabTot; il++) labels[il] = clMCArray.GetArray()[il];
1764     clus->SetParents(nLabTot,labels);
1765
1766   } // rec point array
1767 }
1768
1769 //_____________________________________________________
1770 void AliEMCALTenderSupply::GetPass()
1771 {
1772   // Get passx from filename
1773
1774   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1775   fInputTree = mgr->GetTree();
1776
1777   if (!fInputTree)
1778   {
1779     AliError("Pointer to tree = 0, returning");
1780     return;
1781   }
1782
1783   fInputFile = fInputTree->GetCurrentFile();
1784   if (!fInputFile) {
1785     AliError("Null pointer input file, returning");
1786     return;
1787   }
1788
1789   TString fname(fInputFile->GetName());
1790   if      (fname.Contains("pass1")) fFilepass = TString("pass1");
1791   else if (fname.Contains("pass2")) fFilepass = TString("pass2");
1792   else if (fname.Contains("pass3")) fFilepass = TString("pass3");
1793   else if (fname.Contains("pass4")) fFilepass = TString("pass4");
1794   else if (fname.Contains("pass5")) fFilepass = TString("pass5");
1795   else if (fname.Contains("LHC11c") &&
1796            fname.Contains("spc_calo")) fFilepass = TString("spc_calo");
1797   else if (fname.Contains("calo") || fname.Contains("high_lumi"))
1798
1799   {
1800     //printf("AliEMCALTenderSupply::GetPass() - Path contains <calo> or <high-lumi>, set as <pass1>\n");
1801     fFilepass = TString("pass1");
1802   }
1803   else
1804   {
1805     AliError(Form("Pass number string not found: %s", fname.Data()));
1806     return;
1807   }
1808 }