]> 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
Correct setting of non linearity correction in Tender
[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();
292 fReCalibCluster= kFALSE;
293 fReClusterize = kTRUE;
294 }
295
296 // Reclusterize
acf53135 297 if(fReClusterize)
8b775c10 298 {
299 FillDigitsArray();
300 Clusterize();
301 UpdateClusters();
302 }
303
766cc9de 304 // Store good clusters
305 TClonesArray *clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
306 if (!clusArr)
307 clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
308 if (!clusArr) {
309 AliWarning(Form("No cluster array, number of cells in event = %d, returning", cells->GetNumberOfCells()));
310 return;
311 }
8b775c10 312
766cc9de 313 Int_t nclusters = clusArr->GetEntriesFast();
acf53135 314
766cc9de 315 for (Int_t icluster=0; icluster < nclusters; ++icluster) {
316 AliVCluster *clust = static_cast<AliVCluster*>(clusArr->At(icluster));
317 if (!clust)
318 continue;
319 if (!clust->IsEMCAL())
320 continue;
766cc9de 321 if (fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo, clust->GetCellsAbsId(), clust->GetNCells())) {
322 delete clusArr->RemoveAt(icluster);
323 continue; //todo is it really needed to remove it? Or should we flag it?
324 }
6b7df55d 325
766cc9de 326 if (fFiducial){
327 if (!fEMCALRecoUtils->CheckCellFiducialRegion(fEMCALGeo, clust, cells)){
328 delete clusArr->RemoveAt(icluster);
329 continue; // todo it would be nice to store the distance
330 }
331 }
6b7df55d 332
766cc9de 333 if(fRecalDistToBadChannels)
334 fEMCALRecoUtils->RecalculateClusterDistanceToBadChannel(fEMCALGeo, cells, clust);
335 if(fReCalibCluster)
336 fEMCALRecoUtils->RecalibrateClusterEnergy(fEMCALGeo, clust, cells);
337 if(fRecalClusPos)
338 fEMCALRecoUtils->RecalculateClusterPosition(fEMCALGeo, cells, clust);
57131575 339
340 Float_t correctedEnergy = fEMCALRecoUtils->CorrectClusterEnergyLinearity(clust);
341 clust->SetE(correctedEnergy);
342
766cc9de 343 }
57131575 344
345
766cc9de 346 clusArr->Compress();
6b7df55d 347
766cc9de 348 // Track matching
349 if (!TGeoGlobalMagField::Instance()->GetField()) {
350 const AliESDRun *erun = event->GetESDRun();
351 AliMagF *field = AliMagF::CreateFieldMap(erun->GetCurrentL3(),
352 erun->GetCurrentDip(),
353 AliMagF::kConvLHC,
354 kFALSE,
355 erun->GetBeamEnergy(),
356 erun->GetBeamType());
357 TGeoGlobalMagField::Instance()->SetField(field);
358 }
a6a717c7 359
766cc9de 360 fEMCALRecoUtils->FindMatches(event,0x0,fEMCALGeo);
8b775c10 361
57131575 362 fEMCALRecoUtils->SetClusterMatchedToTrack(event);
363 fEMCALRecoUtils->SetTracksMatchedToCluster(event);
6b7df55d 364
6d67b5a7 365}
366
367//_____________________________________________________
368Bool_t AliEMCALTenderSupply::InitMisalignMatrix()
369{
766cc9de 370 // Initialising misalignment matrices
6b7df55d 371
766cc9de 372 AliESDEvent *event = fTender->GetEvent();
373 if (!event)
374 return kFALSE;
6b7df55d 375
766cc9de 376 if (fGeomMatrixSet) {
377 AliInfo("Misalignment matrix already set");
378 return kTRUE;
379 }
6b7df55d 380
766cc9de 381 if (fDebugLevel>0)
382 AliInfo("Initialising misalignment matrix");
6b7df55d 383
766cc9de 384 fEMCALRecoUtils->SetPositionAlgorithm(AliEMCALRecoUtils::kPosTowerGlobal);
6b7df55d 385
766cc9de 386 if (fLoadGeomMatrices) {
387 for(Int_t mod=0; mod < fEMCALGeo->GetNumberOfSuperModules(); ++mod) {
388 if (fEMCALMatrix[mod]){
389 if(fDebugLevel > 2)
390 fEMCALMatrix[mod]->Print();
391 fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
392 }
393 }
394 fGeomMatrixSet = kTRUE;
395 return kTRUE;
396 }
6b7df55d 397
766cc9de 398 Int_t runGM = event->GetRunNumber();
6b7df55d 399 TObjArray *mobj = 0;
acf53135 400
401 if(fMisalignSurvey == kdefault){ //take default alignment corresponding to run no
57131575 402 printf("***DEFAULT MATRICES***\n!");
acf53135 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;
694 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
695 digit->SetId(cellNumber);
696 digit->SetTime(cellTime);
697 digit->SetTimeR(cellTime);
698 digit->SetIndexInList(idigit);
699 digit->SetType(AliEMCALDigit::kHG);
700 digit->SetAmplitude(cellAmplitude);
701 idigit++;
702 }
b20bc239 703}
704
705//_____________________________________________________
706void AliEMCALTenderSupply::Clusterize()
707{
766cc9de 708 // Clusterize.
6b7df55d 709
766cc9de 710 fClusterizer->Digits2Clusters("");
b20bc239 711}
712
713//_____________________________________________________
714void AliEMCALTenderSupply::UpdateClusters()
715{
766cc9de 716 // Update ESD cluster list.
6b7df55d 717
766cc9de 718 AliESDEvent *event = fTender->GetEvent();
719 if (!event)
720 return;
6b7df55d 721
766cc9de 722 TClonesArray *clus = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
723 if (!clus)
724 clus = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
725 if (!clus) {
6b7df55d 726 AliError(" Null pointer to calo clusters array, returning");
766cc9de 727 return;
728 }
6b7df55d 729
766cc9de 730 Int_t nents = clus->GetEntriesFast();
731 for (Int_t i=0; i < nents; ++i) {
732 AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->At(i));
733 if (!c)
734 continue;
735 if (c->IsEMCAL())
736 delete clus->RemoveAt(i);
737 }
738 clus->Compress();
739 RecPoints2Clusters(clus);
b20bc239 740}
741
742//_____________________________________________________
743void AliEMCALTenderSupply::RecPoints2Clusters(TClonesArray *clus)
744{
766cc9de 745 // Convert AliEMCALRecoPoints to AliESDCaloClusters.
746 // Cluster energy, global position, cells and their amplitude fractions are restored.
747
748 AliESDEvent *event = fTender->GetEvent();
749 if (!event)
750 return;
6b7df55d 751
766cc9de 752 Int_t ncls = fClusterArr->GetEntriesFast();
753 for(Int_t i=0, nout=clus->GetEntriesFast(); i < ncls; ++i) {
754 AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
755
756 Int_t ncells_true = 0;
757 const Int_t ncells = recpoint->GetMultiplicity();
758 UShort_t absIds[ncells];
759 Double32_t ratios[ncells];
760 Int_t *dlist = recpoint->GetDigitsList();
761 Float_t *elist = recpoint->GetEnergiesList();
762 for (Int_t c = 0; c < ncells; ++c) {
763 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
764 absIds[ncells_true] = digit->GetId();
765 ratios[ncells_true] = elist[c]/digit->GetAmplitude();
766 if (ratios[ncells_true] < 0.001)
767 continue;
768 ++ncells_true;
769 }
770
771 if (ncells_true < 1) {
772 AliWarning("Skipping cluster with no cells");
773 continue;
774 }
775
776 // calculate new cluster position
777 TVector3 gpos;
778 recpoint->GetGlobalPosition(gpos);
779 Float_t g[3];
780 gpos.GetXYZ(g);
781
782 AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->New(nout++));
783 c->SetType(AliVCluster::kEMCALClusterv1);
784 c->SetE(recpoint->GetEnergy());
785 c->SetPosition(g);
786 c->SetNCells(ncells_true);
787 c->SetDispersion(recpoint->GetDispersion());
788 c->SetEmcCpvDistance(-1); //not yet implemented
789 c->SetChi2(-1); //not yet implemented
790 c->SetTOF(recpoint->GetTime()) ; //time-of-flight
791 c->SetNExMax(recpoint->GetNExMax()); //number of local maxima
792 Float_t elipAxis[2];
793 recpoint->GetElipsAxis(elipAxis);
794 c->SetM02(elipAxis[0]*elipAxis[0]) ;
795 c->SetM20(elipAxis[1]*elipAxis[1]) ;
796 AliESDCaloCluster *cesd = static_cast<AliESDCaloCluster*>(c);
797 cesd->SetCellsAbsId(absIds);
798 cesd->SetCellsAmplitudeFraction(ratios);
799 }
800}
b20bc239 801
766cc9de 802//_____________________________________________________
803void AliEMCALTenderSupply::GetPass()
804{
805 // Get passx from filename.
6b7df55d 806
766cc9de 807 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
808 fInputTree = mgr->GetTree();
6b7df55d 809
766cc9de 810 if (!fInputTree) {
6b7df55d 811 AliError("Pointer to tree = 0, returning");
766cc9de 812 return;
813 }
6b7df55d 814
766cc9de 815 fInputFile = fInputTree->GetCurrentFile();
816 if (!fInputFile) {
6b7df55d 817 AliError("Null pointer input file, returning");
766cc9de 818 return;
819 }
6b7df55d 820
766cc9de 821 TString fname(fInputFile->GetName());
822 if (fname.Contains("pass1")) fFilepass = TString("pass1");
823 else if(fname.Contains("pass2")) fFilepass = TString("pass2");
824 else if(fname.Contains("pass3")) fFilepass = TString("pass3");
825 else {
826 AliError(Form("Pass number string not found: %s", fname.Data()));
827 return;
828 }
b20bc239 829}
830