]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ANALYSIS/TenderSupplies/AliEMCALTenderSupply.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Set ID number for the new cluster in case of reclusterization, in case of reclusteriz...
[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
16
17///////////////////////////////////////////////////////////////////////////////
18// //
19// EMCAL tender, apply corrections to EMCAl clusters //
766cc9de 20// and do track matching. //
21// Author: Deepa Thomas (Utrecht University) //
6d67b5a7 22// //
23///////////////////////////////////////////////////////////////////////////////
24
766cc9de 25#include <TROOT.h>
26#include <TFile.h>
27#include <TTree.h>
28#include <TInterpreter.h>
29#include <TObjArray.h>
30#include <TClonesArray.h>
31#include <TGeoGlobalMagField.h>
6d67b5a7 32#include <AliLog.h>
33#include <AliESDEvent.h>
34#include <AliAnalysisManager.h>
35#include <AliTender.h>
36#include "AliOADBContainer.h"
b20bc239 37#include "AliCDBManager.h"
38#include "AliCDBStorage.h"
39#include "AliCDBEntry.h"
6d67b5a7 40#include "AliMagF.h"
6d67b5a7 41#include "AliESDCaloCluster.h"
42#include "AliEMCALTenderSupply.h"
6d67b5a7 43#include "AliEMCALGeometry.h"
44#include "AliEMCALRecoUtils.h"
b20bc239 45#include "AliEMCALClusterizer.h"
46#include "AliEMCALRecParam.h"
47#include "AliEMCALRecPoint.h"
48#include "AliEMCALAfterBurnerUF.h"
49#include "AliEMCALClusterizerNxN.h"
50#include "AliEMCALClusterizerv1.h"
51#include "AliEMCALClusterizerv2.h"
52#include "AliEMCALDigit.h"
6d67b5a7 53
3b7d451e 54ClassImp(AliEMCALTenderSupply)
55
6d67b5a7 56AliEMCALTenderSupply::AliEMCALTenderSupply() :
6b7df55d 57AliTenderSupply()
58,fEMCALGeo(0x0)
59,fEMCALGeoName("EMCAL_COMPLETEV1")
60,fEMCALRecoUtils(0)
61,fConfigName("")
62,fDebugLevel(0)
63,fNonLinearFunc(AliEMCALRecoUtils::kNoCorrection)
64,fNonLinearThreshold(30)
65,fReCalibCluster(kFALSE)
acf53135 66,fUpdateCell(kFALSE)
6b7df55d 67,fRecalClusPos(kFALSE)
68,fFiducial(kFALSE)
69,fNCellsFromEMCALBorder(1)
70,fRecalDistToBadChannels(kFALSE)
71,fInputTree(0)
72,fInputFile(0)
73,fFilepass(0)
74,fMass(0.139)
75,fStep(50)
76,fCutEtaPhiSum(kTRUE)
77,fCutEtaPhiSeparate(kFALSE)
78,fRcut(0.05)
79,fEtacut(0.025)
80,fPhicut(0.05)
acf53135 81,fBasePath("")
6b7df55d 82,fReClusterize(kFALSE)
83,fClusterizer(0)
84,fGeomMatrixSet(kFALSE)
85,fLoadGeomMatrices(kFALSE)
acf53135 86,fRecParam(0x0)
87,fRecParamSet(kFALSE)
6b7df55d 88,fOCDBpath(" ")
89,fUnfolder(0)
90,fDigitsArr(0)
91,fClusterArr(0)
acf53135 92,fMisalignSurvey(kdefault)
6d67b5a7 93{
766cc9de 94 // Default constructor.
6b7df55d 95 for(Int_t i = 0; i < 10; i++) fEMCALMatrix[i] = 0 ;
6d67b5a7 96}
97
98//_____________________________________________________
99AliEMCALTenderSupply::AliEMCALTenderSupply(const char *name, const AliTender *tender) :
6b7df55d 100AliTenderSupply(name,tender)
101,fEMCALGeo(0x0)
102,fEMCALGeoName("EMCAL_COMPLETEV1")
103,fEMCALRecoUtils(0)
104,fConfigName("")
105,fDebugLevel(0)
106,fNonLinearFunc(AliEMCALRecoUtils::kNoCorrection)
107,fNonLinearThreshold(30)
108,fReCalibCluster(kFALSE)
acf53135 109,fUpdateCell(kFALSE)
6b7df55d 110,fRecalClusPos(kFALSE)
111,fFiducial(kFALSE)
112,fNCellsFromEMCALBorder(1)
113,fRecalDistToBadChannels(kFALSE)
114,fInputTree(0)
115,fInputFile(0)
116,fFilepass(0)
117,fMass(0.139)
118,fStep(50)
119,fCutEtaPhiSum(kTRUE)
120,fCutEtaPhiSeparate(kFALSE)
121,fRcut(0.05)
122,fEtacut(0.025)
123,fPhicut(0.05)
acf53135 124,fBasePath("")
6b7df55d 125,fReClusterize(kFALSE)
126,fClusterizer(0)
127,fGeomMatrixSet(kFALSE)
128,fLoadGeomMatrices(kFALSE)
acf53135 129,fRecParam(0x0)
130,fRecParamSet(kFALSE)
6b7df55d 131,fOCDBpath(" ")
132,fUnfolder(0)
133,fDigitsArr(0)
134,fClusterArr(0)
acf53135 135,fMisalignSurvey(kdefault)
6d67b5a7 136{
766cc9de 137 // Named constructor
6b7df55d 138
766cc9de 139 for(Int_t i = 0; i < 10; i++) fEMCALMatrix[i] = 0 ;
6b7df55d 140
766cc9de 141 fEMCALRecoUtils = new AliEMCALRecoUtils();
acf53135 142 fRecParam = new AliEMCALRecParam;
766cc9de 143 fDigitsArr = new TClonesArray("AliEMCALDigit",1000);
6d67b5a7 144}
145
146//_____________________________________________________
147AliEMCALTenderSupply::~AliEMCALTenderSupply()
148{
766cc9de 149 //Destructor
6b7df55d 150
766cc9de 151 delete fEMCALRecoUtils;
8b775c10 152 delete fRecParam;
766cc9de 153 delete fClusterizer;
154 delete fUnfolder;
155 if (fDigitsArr){
156 fDigitsArr->Clear("C");
157 delete fDigitsArr;
158 }
6d67b5a7 159}
160
161//_____________________________________________________
162void AliEMCALTenderSupply::Init()
163{
766cc9de 164 // Initialise EMCAL tender.
acf53135 165
766cc9de 166 if (fDebugLevel>0)
167 AliInfo("Init EMCAL Tender supply");
6b7df55d 168
766cc9de 169 if(fConfigName.Length()>0 && gROOT->LoadMacro(fConfigName) >=0) {
170 AliDebug(1, Form("Loading settings from macro %s", fConfigName.Data()));
171 AliEMCALTenderSupply *tender = (AliEMCALTenderSupply*)gInterpreter->ProcessLine("ConfigEMCALTenderSupply()");
172 fDebugLevel = tender->fDebugLevel;
173 fEMCALGeoName = tender->fEMCALGeoName;
766cc9de 174 fEMCALRecoUtils = tender->fEMCALRecoUtils;
175 fConfigName = tender->fConfigName;
176 fNonLinearFunc = tender->fNonLinearFunc;
177 fNonLinearThreshold = tender->fNonLinearThreshold;
178 fReCalibCluster = tender->fReCalibCluster;
766cc9de 179 fRecalClusPos = tender->fRecalClusPos;
180 fFiducial = tender->fFiducial;
181 fNCellsFromEMCALBorder = tender->fNCellsFromEMCALBorder;
182 fRecalDistToBadChannels = tender->fRecalDistToBadChannels;
183 fMass = tender->fMass;
184 fStep = tender->fStep;
185 fRcut = tender->fRcut;
186 fEtacut = tender->fEtacut;
187 fPhicut = tender->fPhicut;
188 fReClusterize = tender->fReClusterize;
189 fLoadGeomMatrices = tender->fLoadGeomMatrices;
190 fRecParam = tender->fRecParam;
191 fOCDBpath = tender->fOCDBpath;
192 for(Int_t i = 0; i < 10; i++)
193 fEMCALMatrix[i] = tender->fEMCALMatrix[i] ;
194 }
6b7df55d 195
766cc9de 196 // Init goemetry
197 fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
6b7df55d 198
766cc9de 199 // Initialising Non linearity parameters
200 fEMCALRecoUtils->SetNonLinearityThreshold(fNonLinearThreshold);
201 fEMCALRecoUtils->SetNonLinearityFunction(fNonLinearFunc);
6b7df55d 202
766cc9de 203 // Setting mass, step size and residual cut
204 fEMCALRecoUtils->SetMass(fMass);
205 fEMCALRecoUtils->SetStep(fStep);
206 if(fCutEtaPhiSum){
207 fEMCALRecoUtils->SwitchOnCutEtaPhiSum();
208 fEMCALRecoUtils->SetCutR(fRcut);
209 } else if (fCutEtaPhiSeparate) {
210 fEMCALRecoUtils->SwitchOnCutEtaPhiSeparate();
211 fEMCALRecoUtils->SetCutEta(fEtacut);
212 fEMCALRecoUtils->SetCutPhi(fPhicut);
213 }
6d67b5a7 214}
215
216//_____________________________________________________
217void AliEMCALTenderSupply::ProcessEvent()
218{
766cc9de 219 // Event loop.
6b7df55d 220
766cc9de 221 AliESDEvent *event = fTender->GetEvent();
222 if (!event) {
223 AliError("ESD event ptr = 0, returning");
224 return;
225 }
8b775c10 226
766cc9de 227 // Initialising parameters once per run number
228 if(fTender->RunChanged()){
229 GetPass();
acf53135 230 // Init bad channels
231 Int_t fInitBC=InitBadChannels();
232
233 if (fInitBC==0)
8b775c10 234 {
235 AliError("InitBadChannels returned false, returning");
236 return;
237 }
acf53135 238 if(fInitBC>1)
8b775c10 239 {
240 AliInfo(Form("No external hot channel set: %d - %s", event->GetRunNumber(), fFilepass.Data()));
241 }
242
243 if (fReCalibCluster || fUpdateCell) {
acf53135 244 Int_t fInitRecalib=InitRecalib();
245 if (fInitRecalib==0)
8b775c10 246 {
247 AliError("InitRecalib returned false, returning");
248 return;
249 }
acf53135 250 if(fInitRecalib >1)
8b775c10 251 {
252 AliInfo(Form("No recalibration available: %d - %s", event->GetRunNumber(), fFilepass.Data()));
253 fUpdateCell = kFALSE;
254 fReCalibCluster = kFALSE;
255 }
766cc9de 256 }
8b775c10 257
acf53135 258 if (fRecalClusPos || fReClusterize || fUpdateCell) {
766cc9de 259 if (!InitMisalignMatrix()) {
260 AliError("InitMisalignmentMatrix returned false, returning");
261 return;
262 }
263 }
8b775c10 264
acf53135 265 if (fReClusterize || fUpdateCell) {
766cc9de 266 if (!InitClusterization()) {
267 AliError("InitClusterization returned false, returning");
268 return;
269 }
270 }
6b7df55d 271
766cc9de 272 if(fDebugLevel>1)
273 fEMCALRecoUtils->Print("");
8b775c10 274
766cc9de 275 }
8b775c10 276
766cc9de 277 // Test if cells present
278 AliESDCaloCells *cells= event->GetEMCALCells();
279 if (cells->GetNumberOfCells()<=0) {
8b775c10 280 if(fDebugLevel>1)
281 AliWarning(Form("Number of EMCAL cells = %d, returning", cells->GetNumberOfCells()));
766cc9de 282 return;
283 }
6b7df55d 284
acf53135 285 if(fDebugLevel>2)
8b775c10 286 AliInfo(Form("Re-calibrate cluster %d\n",fReCalibCluster));
287
288 // Recalibrate cells
acf53135 289 if(fUpdateCell)
8b775c10 290 {
291 UpdateCells();
8b775c10 292 fReClusterize = kTRUE;
293 }
294
295 // Reclusterize
acf53135 296 if(fReClusterize)
8b775c10 297 {
94f7c836 298 fRecalClusPos = kFALSE; // Done in reclusterization, do nothing
299 fReCalibCluster= kFALSE; // Done in reclusterization, do nothing
8b775c10 300 FillDigitsArray();
301 Clusterize();
302 UpdateClusters();
303 }
304
766cc9de 305 // Store good clusters
306 TClonesArray *clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
307 if (!clusArr)
308 clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
309 if (!clusArr) {
310 AliWarning(Form("No cluster array, number of cells in event = %d, returning", cells->GetNumberOfCells()));
311 return;
312 }
8b775c10 313
766cc9de 314 Int_t nclusters = clusArr->GetEntriesFast();
acf53135 315
766cc9de 316 for (Int_t icluster=0; icluster < nclusters; ++icluster) {
317 AliVCluster *clust = static_cast<AliVCluster*>(clusArr->At(icluster));
318 if (!clust)
319 continue;
320 if (!clust->IsEMCAL())
321 continue;
766cc9de 322 if (fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo, clust->GetCellsAbsId(), clust->GetNCells())) {
323 delete clusArr->RemoveAt(icluster);
324 continue; //todo is it really needed to remove it? Or should we flag it?
325 }
6b7df55d 326
766cc9de 327 if (fFiducial){
328 if (!fEMCALRecoUtils->CheckCellFiducialRegion(fEMCALGeo, clust, cells)){
329 delete clusArr->RemoveAt(icluster);
330 continue; // todo it would be nice to store the distance
331 }
332 }
6b7df55d 333
766cc9de 334 if(fRecalDistToBadChannels)
335 fEMCALRecoUtils->RecalculateClusterDistanceToBadChannel(fEMCALGeo, cells, clust);
336 if(fReCalibCluster)
337 fEMCALRecoUtils->RecalibrateClusterEnergy(fEMCALGeo, clust, cells);
338 if(fRecalClusPos)
339 fEMCALRecoUtils->RecalculateClusterPosition(fEMCALGeo, cells, clust);
57131575 340
341 Float_t correctedEnergy = fEMCALRecoUtils->CorrectClusterEnergyLinearity(clust);
342 clust->SetE(correctedEnergy);
343
766cc9de 344 }
57131575 345
346
766cc9de 347 clusArr->Compress();
6b7df55d 348
766cc9de 349 // Track matching
350 if (!TGeoGlobalMagField::Instance()->GetField()) {
351 const AliESDRun *erun = event->GetESDRun();
352 AliMagF *field = AliMagF::CreateFieldMap(erun->GetCurrentL3(),
353 erun->GetCurrentDip(),
354 AliMagF::kConvLHC,
355 kFALSE,
356 erun->GetBeamEnergy(),
357 erun->GetBeamType());
358 TGeoGlobalMagField::Instance()->SetField(field);
359 }
a6a717c7 360
766cc9de 361 fEMCALRecoUtils->FindMatches(event,0x0,fEMCALGeo);
8b775c10 362
57131575 363 fEMCALRecoUtils->SetClusterMatchedToTrack(event);
364 fEMCALRecoUtils->SetTracksMatchedToCluster(event);
6b7df55d 365
6d67b5a7 366}
367
368//_____________________________________________________
369Bool_t AliEMCALTenderSupply::InitMisalignMatrix()
370{
766cc9de 371 // Initialising misalignment matrices
6b7df55d 372
766cc9de 373 AliESDEvent *event = fTender->GetEvent();
374 if (!event)
375 return kFALSE;
6b7df55d 376
766cc9de 377 if (fGeomMatrixSet) {
378 AliInfo("Misalignment matrix already set");
379 return kTRUE;
380 }
6b7df55d 381
766cc9de 382 if (fDebugLevel>0)
383 AliInfo("Initialising misalignment matrix");
6b7df55d 384
766cc9de 385 fEMCALRecoUtils->SetPositionAlgorithm(AliEMCALRecoUtils::kPosTowerGlobal);
6b7df55d 386
766cc9de 387 if (fLoadGeomMatrices) {
388 for(Int_t mod=0; mod < fEMCALGeo->GetNumberOfSuperModules(); ++mod) {
389 if (fEMCALMatrix[mod]){
390 if(fDebugLevel > 2)
391 fEMCALMatrix[mod]->Print();
392 fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
393 }
394 }
395 fGeomMatrixSet = kTRUE;
396 return kTRUE;
397 }
6b7df55d 398
766cc9de 399 Int_t runGM = event->GetRunNumber();
6b7df55d 400 TObjArray *mobj = 0;
acf53135 401
402 if(fMisalignSurvey == kdefault){ //take default alignment corresponding to run no
403 AliOADBContainer emcalgeoCont(Form("emcal"));
404 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
405 mobj=(TObjArray*)emcalgeoCont.GetObject(runGM,"EmcalMatrices");
406 }
407
408 if(fMisalignSurvey == kSurveybyS){ //take alignment at sector level
766cc9de 409 if (runGM <= 140000) { //2010 data
6b7df55d 410 AliOADBContainer emcalgeoCont(Form("emcal2010"));
411 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
412 mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
413
acf53135 414 } else if (runGM>140000) { // 2011 LHC11a pass1 data
766cc9de 415 AliOADBContainer emcalgeoCont(Form("emcal2011"));
6b7df55d 416 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
acf53135 417 mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byS");
418 }
419 }
420
57131575 421 if(fMisalignSurvey == kSurveybyM){ //take alignment at module level
acf53135 422 if (runGM <= 140000) { //2010 data
423 AliOADBContainer emcalgeoCont(Form("emcal2010"));
424 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
425 mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
426
427 } else if (runGM>140000) { // 2011 LHC11a pass1 data
57131575 428
acf53135 429 AliOADBContainer emcalgeoCont(Form("emcal2011"));
430 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
431 mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byM");
766cc9de 432 }
acf53135 433 }
434
435 for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
436 fEMCALMatrix[mod] = (TGeoHMatrix*) mobj->At(mod);
437 fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
438 fEMCALMatrix[mod]->Print();
439 }
766cc9de 440 return kTRUE;
6d67b5a7 441}
442
443//_____________________________________________________
acf53135 444Int_t AliEMCALTenderSupply::InitBadChannels()
6d67b5a7 445{
766cc9de 446 // Initialising bad channel maps
766cc9de 447 AliESDEvent *event = fTender->GetEvent();
448 if (!event)
acf53135 449 return 0;
6b7df55d 450
766cc9de 451 if (fDebugLevel>0)
452 AliInfo("Initialising Bad channel map");
6b7df55d 453
766cc9de 454 if (fFiducial){
455 fEMCALRecoUtils->SetNumberOfCellsFromEMCALBorder(fNCellsFromEMCALBorder);
456 fEMCALRecoUtils->SwitchOnNoFiducialBorderInEMCALEta0();
457 }
6b7df55d 458
766cc9de 459 fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
460 if (fRecalDistToBadChannels)
461 fEMCALRecoUtils->SwitchOnDistToBadChannelRecalculation();
6b7df55d 462
766cc9de 463 Int_t runBC = event->GetRunNumber();
6b7df55d 464
acf53135 465 AliOADBContainer *contBC=new AliOADBContainer("");
466 if(fBasePath!=""){ //if fBasePath specified in the ->SetBasePath()
467 if (fDebugLevel>0) AliInfo(Form("Loading Bad Channels OADB from given path %s",fBasePath.Data()));
468 TFile *fbad=new TFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"read");
469 if(fbad->IsZombie()){
470 AliFatal(Form("EMCALBadChannels.root was not found in the path provided: %s, aborting",fBasePath.Data()));
471 return 0;
472 }
473 if(fbad) delete fbad;
474 contBC->InitFromFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"AliEMCALBadChannels");
475 }
476 else { // Else choose the one in the $ALICE_ROOT directory
477 if (fDebugLevel>0) AliInfo("Loading Bad Channels OADB from $ALICE_ROOT/OADB/EMCAL");
478 TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","read");
479 if(fbad->IsZombie()){
480 AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root was not found, aborting");
481 return 0;
482 }
483 if(fbad) delete fbad;
484 contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","AliEMCALBadChannels");
766cc9de 485 }
acf53135 486
487 TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
488 if(!arrayBC){
489 AliError(Form("No external hot channel set for run number: %d", runBC));
490 return 2;
766cc9de 491 }
6b7df55d 492
acf53135 493 TObjArray *arrayBCpass=(TObjArray*)arrayBC->FindObject(fFilepass); // Here, it looks for a specific pass
494 if(!arrayBCpass){
495 AliError(Form("No external hot channel set for: %d -%s", runBC,fFilepass.Data()));
496 return 2;
497 }
498
499 if(fDebugLevel>0) arrayBCpass->Print();
500
501 Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
502 for (Int_t i=0; i<sms; ++i) {
503 TH2I *h = fEMCALRecoUtils->GetEMCALChannelStatusMap(i);
504 if (h)
505 delete h;
506 h=(TH2I*)arrayBCpass->FindObject(Form("EMCALBadChannelMap_Mod%d",i));
507
508 if (!h) {
509 AliError(Form("Can not get EMCALBadChannelMap_Mod%d",i));
510 continue;
511 }
512 h->SetDirectory(0);
513 fEMCALRecoUtils->SetEMCALChannelStatusMap(i,h);
766cc9de 514 }
6b7df55d 515
acf53135 516 return 1;
6d67b5a7 517}
518
acf53135 519
6d67b5a7 520//_____________________________________________________
acf53135 521Int_t AliEMCALTenderSupply::InitRecalib()
6d67b5a7 522{
766cc9de 523 // Initialising recalibration factors.
6b7df55d 524
766cc9de 525 AliESDEvent *event = fTender->GetEvent();
526 if (!event)
acf53135 527 return 0;
6b7df55d 528
766cc9de 529 if (fDebugLevel>0)
530 AliInfo("Initialising recalibration factors");
6b7df55d 531
766cc9de 532 fEMCALRecoUtils->SwitchOnRecalibration();
533 Int_t runRC = event->GetRunNumber();
acf53135 534
535 AliOADBContainer *contRF=new AliOADBContainer("");
536 if(fBasePath!="") { //if fBasePath specified in the ->SetBasePath()
537 if (fDebugLevel>0) AliInfo(Form("Loading Recalib OADB from given path %s",fBasePath.Data()));
538 TFile *fRecalib= new TFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"read");
766cc9de 539 if (fRecalib->IsZombie()) {
acf53135 540 AliFatal(Form("EMCALRecalib.root not found in %s, aborting",fBasePath.Data()));
541 return 0;
766cc9de 542 }
acf53135 543 if(fRecalib) delete fRecalib;
544 contRF->InitFromFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"AliEMCALRecalib");
545 }
546 else{ // Else choose the one in the $ALICE_ROOT directory
547 if (fDebugLevel>0) AliInfo("Loading Recalib OADB from $ALICE_ROOT/OADB/EMCAL");
548 TFile *fRecalib= new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","read");
549 if (fRecalib->IsZombie()) {
550 AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root was not found, aborting");
551 return 0;
766cc9de 552 }
acf53135 553 if(fRecalib) delete fRecalib;
554 contRF->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","AliEMCALRecalib");
766cc9de 555 }
acf53135 556
557 TObjArray *recal=(TObjArray*)contRF->GetObject(runRC); //GetObject(int runnumber)
558 if(!recal){
559 AliError(Form("No Objects for run: %d",runRC));
560 return 2;
561 }
562
563 TObjArray *recalpass=(TObjArray*)recal->FindObject(fFilepass);
564 if(!recalpass){
565 AliError(Form("No Objects for run: %d - %s",runRC,fFilepass.Data()));
566 return 2;
567 }
568
569 TObjArray *recalib=(TObjArray*)recalpass->FindObject("Recalib");
570 if(!recalib){
571 AliError(Form("No Recalib histos found for %d - %s",runRC,fFilepass.Data()));
572 return 2;
573 }
574
575 if(fDebugLevel>0) recalib->Print();
576
577 Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
578 for (Int_t i=0; i<sms; ++i) {
579 TH2F *h = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactors(i);
580 if (h)
581 delete h;
582 h = (TH2F*)recalib->FindObject(Form("EMCALRecalFactors_SM%d",i));
583 if (!h) {
584 AliError(Form("Could not load EMCALRecalFactors_SM%d",i));
585 continue;
766cc9de 586 }
acf53135 587 h->SetDirectory(0);
588 fEMCALRecoUtils->SetEMCALChannelRecalibrationFactors(i,h);
766cc9de 589 }
acf53135 590 return 1;
c958a2f7 591}
592
acf53135 593//_____________________________________________________
594void AliEMCALTenderSupply::UpdateCells()
595{
8b775c10 596 //Remove bad cells from the cell list
597 //Recalibrate energy and time cells
acf53135 598 //This is required for later reclusterization
599
600 AliESDCaloCells *cells = fTender->GetEvent()->GetEMCALCells();
8b775c10 601 Int_t bunchCrossNo = fTender->GetEvent()->GetBunchCrossNumber();
acf53135 602
8b775c10 603 fEMCALRecoUtils->SwitchOnRecalibration();
604 fEMCALRecoUtils->SwitchOnTimeRecalibration();
605
606 fEMCALRecoUtils->RecalibrateCells(cells, bunchCrossNo);
607
608}
acf53135 609
610//_____________________________________________________
611void AliEMCALTenderSupply::InitRecParam()
612{
613 if (fDebugLevel>0)
614 AliInfo("Initialize the recParam");
615 fRecParam = new AliEMCALRecParam;
616}
617
618
b20bc239 619//_____________________________________________________
766cc9de 620Bool_t AliEMCALTenderSupply::InitClusterization()
b20bc239 621{
766cc9de 622 // Initialing clusterization/unfolding algorithm and set all the needed parameters.
6b7df55d 623
766cc9de 624 AliESDEvent *event=fTender->GetEvent();
625 if (!event)
626 return kFALSE;
6b7df55d 627
766cc9de 628 if (fDebugLevel>0)
629 AliInfo(Form("Initialising reclustering parameters: Clusterizer-%d",fRecParam->GetClusterizerFlag()));
6b7df55d 630
766cc9de 631 //---setup clusterizer
632 delete fClusterizer;
633 if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
634 fClusterizer = new AliEMCALClusterizerv1 (fEMCALGeo);
635 else if(fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
636 fClusterizer = new AliEMCALClusterizerv2(fEMCALGeo);
637 else if(fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN) {
638 AliEMCALClusterizerNxN *clusterizer = new AliEMCALClusterizerNxN(fEMCALGeo);
639 clusterizer->SetNRowDiff(fRecParam->GetNRowDiff());
640 clusterizer->SetNColDiff(fRecParam->GetNColDiff());
641 fClusterizer = clusterizer;
642 } else {
643 AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
644 return kFALSE;
645 }
6b7df55d 646
766cc9de 647 // Set the clustering parameters
648 fClusterizer->SetECAClusteringThreshold( fRecParam->GetClusteringThreshold() );
649 fClusterizer->SetECALogWeight ( fRecParam->GetW0() );
650 fClusterizer->SetMinECut ( fRecParam->GetMinECut() );
651 fClusterizer->SetUnfolding ( fRecParam->GetUnfold() );
652 fClusterizer->SetECALocalMaxCut ( fRecParam->GetLocMaxCut() );
653 fClusterizer->SetTimeCut ( fRecParam->GetTimeCut() );
654 fClusterizer->SetTimeMin ( fRecParam->GetTimeMin() );
655 fClusterizer->SetTimeMax ( fRecParam->GetTimeMax() );
656 fClusterizer->SetInputCalibrated ( kTRUE );
657 fClusterizer->SetJustClusters ( kTRUE );
6b7df55d 658
766cc9de 659 // In case of unfolding after clusterization is requested, set the corresponding parameters
660 if (fRecParam->GetUnfold()) {
661 for (Int_t i = 0; i < 8; ++i) {
662 fClusterizer->SetSSPars(i, fRecParam->GetSSPars(i));
663 }
664 for (Int_t i = 0; i < 3; ++i) {
665 fClusterizer->SetPar5 (i, fRecParam->GetPar5(i));
666 fClusterizer->SetPar6 (i, fRecParam->GetPar6(i));
667 }
668 fClusterizer->InitClusterUnfolding();
669 }
6b7df55d 670
766cc9de 671 fClusterizer->SetDigitsArr(fDigitsArr);
672 fClusterizer->SetOutput(0);
673 fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
674 return kTRUE;
b20bc239 675}
766cc9de 676
b20bc239 677//_____________________________________________________
678void AliEMCALTenderSupply::FillDigitsArray()
679{
766cc9de 680 // Fill digits from cells to a TClonesArray.
6b7df55d 681
766cc9de 682 AliESDEvent *event = fTender->GetEvent();
683 if (!event)
684 return;
6b7df55d 685
766cc9de 686 fDigitsArr->Clear("C");
687 AliESDCaloCells *cells = event->GetEMCALCells();
688 Int_t ncells = cells->GetNumberOfCells();
689 for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell) {
690 Double_t cellAmplitude=0, cellTime=0;
691 Short_t cellNumber=0;
692 if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime) != kTRUE)
693 break;
94f7c836 694
695 // Do not add if already too low (some cells set to 0 if bad channels)
696 if(cellAmplitude < fRecParam->GetMinECut() )
697 continue;
698
699 // If requested, do not include exotic clusters
700 if(fEMCALRecoUtils->IsExoticCell(cellNumber,cells,event->GetBunchCrossNumber()))
701 continue;
702
766cc9de 703 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
704 digit->SetId(cellNumber);
705 digit->SetTime(cellTime);
706 digit->SetTimeR(cellTime);
707 digit->SetIndexInList(idigit);
708 digit->SetType(AliEMCALDigit::kHG);
709 digit->SetAmplitude(cellAmplitude);
710 idigit++;
711 }
b20bc239 712}
713
714//_____________________________________________________
715void AliEMCALTenderSupply::Clusterize()
716{
766cc9de 717 // Clusterize.
6b7df55d 718
766cc9de 719 fClusterizer->Digits2Clusters("");
b20bc239 720}
721
722//_____________________________________________________
723void AliEMCALTenderSupply::UpdateClusters()
724{
766cc9de 725 // Update ESD cluster list.
6b7df55d 726
766cc9de 727 AliESDEvent *event = fTender->GetEvent();
728 if (!event)
729 return;
6b7df55d 730
766cc9de 731 TClonesArray *clus = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
732 if (!clus)
733 clus = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
734 if (!clus) {
6b7df55d 735 AliError(" Null pointer to calo clusters array, returning");
766cc9de 736 return;
737 }
6b7df55d 738
766cc9de 739 Int_t nents = clus->GetEntriesFast();
740 for (Int_t i=0; i < nents; ++i) {
741 AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->At(i));
742 if (!c)
743 continue;
744 if (c->IsEMCAL())
745 delete clus->RemoveAt(i);
746 }
747 clus->Compress();
748 RecPoints2Clusters(clus);
b20bc239 749}
750
751//_____________________________________________________
752void AliEMCALTenderSupply::RecPoints2Clusters(TClonesArray *clus)
753{
766cc9de 754 // Convert AliEMCALRecoPoints to AliESDCaloClusters.
755 // Cluster energy, global position, cells and their amplitude fractions are restored.
756
757 AliESDEvent *event = fTender->GetEvent();
758 if (!event)
759 return;
94f7c836 760
766cc9de 761 Int_t ncls = fClusterArr->GetEntriesFast();
762 for(Int_t i=0, nout=clus->GetEntriesFast(); i < ncls; ++i) {
763 AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
764
765 Int_t ncells_true = 0;
766 const Int_t ncells = recpoint->GetMultiplicity();
767 UShort_t absIds[ncells];
768 Double32_t ratios[ncells];
769 Int_t *dlist = recpoint->GetDigitsList();
770 Float_t *elist = recpoint->GetEnergiesList();
771 for (Int_t c = 0; c < ncells; ++c) {
772 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
773 absIds[ncells_true] = digit->GetId();
774 ratios[ncells_true] = elist[c]/digit->GetAmplitude();
775 if (ratios[ncells_true] < 0.001)
776 continue;
777 ++ncells_true;
778 }
779
780 if (ncells_true < 1) {
781 AliWarning("Skipping cluster with no cells");
782 continue;
783 }
784
785 // calculate new cluster position
786 TVector3 gpos;
787 recpoint->GetGlobalPosition(gpos);
788 Float_t g[3];
789 gpos.GetXYZ(g);
790
791 AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->New(nout++));
94f7c836 792 c->SetID(nout-1);
766cc9de 793 c->SetType(AliVCluster::kEMCALClusterv1);
794 c->SetE(recpoint->GetEnergy());
795 c->SetPosition(g);
796 c->SetNCells(ncells_true);
797 c->SetDispersion(recpoint->GetDispersion());
798 c->SetEmcCpvDistance(-1); //not yet implemented
799 c->SetChi2(-1); //not yet implemented
800 c->SetTOF(recpoint->GetTime()) ; //time-of-flight
801 c->SetNExMax(recpoint->GetNExMax()); //number of local maxima
802 Float_t elipAxis[2];
803 recpoint->GetElipsAxis(elipAxis);
804 c->SetM02(elipAxis[0]*elipAxis[0]) ;
805 c->SetM20(elipAxis[1]*elipAxis[1]) ;
806 AliESDCaloCluster *cesd = static_cast<AliESDCaloCluster*>(c);
807 cesd->SetCellsAbsId(absIds);
808 cesd->SetCellsAmplitudeFraction(ratios);
94f7c836 809
766cc9de 810 }
811}
b20bc239 812
766cc9de 813//_____________________________________________________
814void AliEMCALTenderSupply::GetPass()
815{
816 // Get passx from filename.
6b7df55d 817
766cc9de 818 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
819 fInputTree = mgr->GetTree();
6b7df55d 820
766cc9de 821 if (!fInputTree) {
6b7df55d 822 AliError("Pointer to tree = 0, returning");
766cc9de 823 return;
824 }
6b7df55d 825
766cc9de 826 fInputFile = fInputTree->GetCurrentFile();
827 if (!fInputFile) {
6b7df55d 828 AliError("Null pointer input file, returning");
766cc9de 829 return;
830 }
6b7df55d 831
766cc9de 832 TString fname(fInputFile->GetName());
833 if (fname.Contains("pass1")) fFilepass = TString("pass1");
834 else if(fname.Contains("pass2")) fFilepass = TString("pass2");
835 else if(fname.Contains("pass3")) fFilepass = TString("pass3");
836 else {
837 AliError(Form("Pass number string not found: %s", fname.Data()));
838 return;
839 }
b20bc239 840}
841