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