1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
17 ///////////////////////////////////////////////////////////////////////////////
19 // EMCAL tender, apply corrections to EMCAL clusters //
20 // and do track matching. //
21 // Author: Deepa Thomas (Utrecht University) //
22 // Later mods/rewrite: Jiri Kral (University of Jyvaskyla) //
24 ///////////////////////////////////////////////////////////////////////////////
29 #include <TInterpreter.h>
30 #include <TObjArray.h>
31 #include <TClonesArray.h>
32 #include <TGeoGlobalMagField.h>
34 #include <AliESDEvent.h>
35 #include <AliAnalysisManager.h>
36 #include <AliTender.h>
37 #include "AliOADBContainer.h"
38 #include "AliCDBManager.h"
39 #include "AliCDBStorage.h"
40 #include "AliCDBEntry.h"
42 #include "AliESDCaloCluster.h"
43 #include "AliEMCALTenderSupply.h"
44 #include "AliEMCALGeometry.h"
45 #include "AliEMCALRecoUtils.h"
46 #include "AliEMCALClusterizer.h"
47 #include "AliEMCALRecParam.h"
48 #include "AliEMCALRecPoint.h"
49 #include "AliEMCALAfterBurnerUF.h"
50 #include "AliEMCALClusterizerNxN.h"
51 #include "AliEMCALClusterizerv1.h"
52 #include "AliEMCALClusterizerv2.h"
53 #include "AliEMCALDigit.h"
54 #include "AliEMCALRecParam.h"
56 ClassImp(AliEMCALTenderSupply)
58 AliEMCALTenderSupply::AliEMCALTenderSupply() :
61 ,fEMCALGeoName("EMCAL_COMPLETEV1")
66 ,fNonLinearThreshold(-1)
67 ,fReCalibCluster(kFALSE)
69 ,fCalibrateEnergy(kFALSE)
70 ,fCalibrateTime(kFALSE)
71 ,fDoNonLinearity(kFALSE)
72 ,fBadCellRemove(kFALSE)
73 ,fRejectExoticCells(kFALSE)
74 ,fRejectExoticClusters(kFALSE)
75 ,fClusterBadChannelCheck(kFALSE)
76 ,fRecalClusPos(kFALSE)
78 ,fNCellsFromEMCALBorder(-1)
79 ,fRecalDistToBadChannels(kFALSE)
80 ,fRecalShowerShape(kFALSE)
87 ,fCutEtaPhiSeparate(kFALSE)
92 ,fReClusterize(kFALSE)
94 ,fGeomMatrixSet(kFALSE)
95 ,fLoadGeomMatrices(kFALSE)
97 ,fDoTrackMatch(kFALSE)
98 ,fDoUpdateOnly(kFALSE)
102 ,fMisalignSurvey(kdefault)
103 ,fExoticCellFraction(-1)
104 ,fExoticCellDiffTime(-1)
105 ,fExoticCellMinAmplitude(-1)
106 ,fRecoParamsOCDBLoaded(kFALSE)
108 // Default constructor.
109 for(Int_t i = 0; i < 10; i++) fEMCALMatrix[i] = 0 ;
110 fEMCALRecoUtils = new AliEMCALRecoUtils;
113 //_____________________________________________________
114 AliEMCALTenderSupply::AliEMCALTenderSupply(const char *name, const AliTender *tender) :
115 AliTenderSupply(name,tender)
117 ,fEMCALGeoName("EMCAL_COMPLETEV1")
122 ,fNonLinearThreshold(-1)
123 ,fReCalibCluster(kFALSE)
125 ,fCalibrateEnergy(kFALSE)
126 ,fCalibrateTime(kFALSE)
127 ,fDoNonLinearity(kFALSE)
128 ,fBadCellRemove(kFALSE)
129 ,fRejectExoticCells(kFALSE)
130 ,fRejectExoticClusters(kFALSE)
131 ,fClusterBadChannelCheck(kFALSE)
132 ,fRecalClusPos(kFALSE)
134 ,fNCellsFromEMCALBorder(-1)
135 ,fRecalDistToBadChannels(kFALSE)
136 ,fRecalShowerShape(kFALSE)
142 ,fCutEtaPhiSum(kTRUE)
143 ,fCutEtaPhiSeparate(kFALSE)
148 ,fReClusterize(kFALSE)
150 ,fGeomMatrixSet(kFALSE)
151 ,fLoadGeomMatrices(kFALSE)
153 ,fDoTrackMatch(kFALSE)
154 ,fDoUpdateOnly(kFALSE)
158 ,fMisalignSurvey(kdefault)
159 ,fExoticCellFraction(-1)
160 ,fExoticCellDiffTime(-1)
161 ,fExoticCellMinAmplitude(-1)
162 ,fRecoParamsOCDBLoaded(kFALSE)
166 for(Int_t i = 0; i < 10; i++) fEMCALMatrix[i] = 0 ;
167 fEMCALRecoUtils = new AliEMCALRecoUtils;
170 //_____________________________________________________
171 AliEMCALTenderSupply::~AliEMCALTenderSupply()
175 delete fEMCALRecoUtils;
179 fDigitsArr->Clear("C");
187 //_____________________________________________________
188 void AliEMCALTenderSupply::Init()
190 // Initialise EMCAL tender.
193 AliWarning("Init EMCAL Tender supply");
195 if (fConfigName.Length()>0 && gROOT->LoadMacro(fConfigName) >=0) {
196 AliDebug(1, Form("Loading settings from macro %s", fConfigName.Data()));
197 AliEMCALTenderSupply *tender = (AliEMCALTenderSupply*)gInterpreter->ProcessLine("ConfigEMCALTenderSupply()");
198 fDebugLevel = tender->fDebugLevel;
199 fEMCALGeoName = tender->fEMCALGeoName;
200 fEMCALRecoUtils = tender->fEMCALRecoUtils;
201 fConfigName = tender->fConfigName;
202 fNonLinearFunc = tender->fNonLinearFunc;
203 fNonLinearThreshold = tender->fNonLinearThreshold;
204 fReCalibCluster = tender->fReCalibCluster;
205 fUpdateCell = tender->fUpdateCell;
206 fRecalClusPos = tender->fRecalClusPos;
207 fCalibrateEnergy = tender->fCalibrateEnergy;
208 fCalibrateTime = tender->fCalibrateTime;
209 fFiducial = tender->fFiducial;
210 fNCellsFromEMCALBorder = tender->fNCellsFromEMCALBorder;
211 fRecalDistToBadChannels = tender->fRecalDistToBadChannels;
212 fRecalShowerShape = tender->fRecalShowerShape;
213 fClusterBadChannelCheck = tender->fClusterBadChannelCheck;
214 fBadCellRemove = tender->fBadCellRemove;
215 fRejectExoticCells = tender->fRejectExoticCells;
216 fRejectExoticClusters = tender->fRejectExoticClusters;
217 fMass = tender->fMass;
218 fStep = tender->fStep;
219 fCutEtaPhiSum = tender->fCutEtaPhiSum;
220 fCutEtaPhiSeparate = tender->fCutEtaPhiSeparate;
221 fRcut = tender->fRcut;
222 fEtacut = tender->fEtacut;
223 fPhicut = tender->fPhicut;
224 fReClusterize = tender->fReClusterize;
225 fLoadGeomMatrices = tender->fLoadGeomMatrices;
226 fRecParam = tender->fRecParam;
227 fDoNonLinearity = tender->fDoNonLinearity;
228 fDoTrackMatch = tender->fDoTrackMatch;
229 fDoUpdateOnly = tender->fDoUpdateOnly;
230 fMisalignSurvey = tender->fMisalignSurvey;
231 fExoticCellFraction = tender->fExoticCellFraction;
232 fExoticCellDiffTime = tender->fExoticCellDiffTime;
233 fExoticCellMinAmplitude = tender->fExoticCellMinAmplitude;
234 fRecoParamsOCDBLoaded = tender->fRecoParamsOCDBLoaded;
236 for(Int_t i = 0; i < 10; i++)
237 fEMCALMatrix[i] = tender->fEMCALMatrix[i] ;
241 AliInfo( "Emcal Tender settings: ======================================" );
242 AliInfo( "------------ Switches --------------------------" );
243 AliInfo( Form( "BadCellRemove : %d", fBadCellRemove ));
244 AliInfo( Form( "ExoticCellRemove : %d", fRejectExoticCells ));
245 AliInfo( Form( "CalibrateEnergy : %d", fCalibrateEnergy ));
246 AliInfo( Form( "CalibrateTime : %d", fCalibrateTime ));
247 AliInfo( Form( "UpdateCell : %d", fUpdateCell ));
248 AliInfo( Form( "DoUpdateOnly : %d", fDoUpdateOnly ));
249 AliInfo( Form( "Reclustering : %d", fReClusterize ));
250 AliInfo( Form( "ClusterBadChannelCheck : %d", fClusterBadChannelCheck ));
251 AliInfo( Form( "ClusterExoticChannelCheck : %d", fRejectExoticClusters ));
252 AliInfo( Form( "CellFiducialRegion : %d", fFiducial ));
253 AliInfo( Form( "ReCalibrateCluster : %d", fReCalibCluster ));
254 AliInfo( Form( "RecalculateClusPos : %d", fRecalClusPos ));
255 AliInfo( Form( "RecalShowerShape : %d", fRecalShowerShape ));
256 AliInfo( Form( "NonLinearityCorrection : %d", fDoNonLinearity ));
257 AliInfo( Form( "RecalDistBadChannel : %d", fRecalDistToBadChannels ));
258 AliInfo( Form( "TrackMatch : %d", fDoTrackMatch ));
259 AliInfo( "------------ Variables -------------------------" );
260 AliInfo( Form( "DebugLevel : %d", fDebugLevel ));
261 AliInfo( Form( "BasePath : %s", fBasePath.Data() ));
262 AliInfo( Form( "ConfigFileName : %s", fConfigName.Data() ));
263 AliInfo( Form( "EMCALGeometryName : %s", fEMCALGeoName.Data() ));
264 AliInfo( Form( "NonLinearityFunction : %d", fNonLinearFunc ));
265 AliInfo( Form( "NonLinearityThreshold : %d", fNonLinearThreshold ));
266 AliInfo( Form( "MisalignmentMatrixSurvey : %d", fMisalignSurvey ));
267 AliInfo( Form( "NumberOfCellsFromEMCALBorder : %d", fNCellsFromEMCALBorder ));
268 AliInfo( Form( "RCut : %f", fRcut ));
269 AliInfo( Form( "Mass : %f", fMass ));
270 AliInfo( Form( "Step : %f", fStep ));
271 AliInfo( Form( "EtaCut : %f", fEtacut ));
272 AliInfo( Form( "PhiCut : %f", fPhicut ));
273 AliInfo( Form( "ExoticCellFraction : %f", fExoticCellFraction ));
274 AliInfo( Form( "ExoticCellDiffTime : %f", fExoticCellDiffTime ));
275 AliInfo( Form( "ExoticCellMinAmplitude : %f", fExoticCellMinAmplitude ));
276 AliInfo( "=============================================================" );
280 fEMCALGeo = AliEMCALGeometry::GetInstance(fEMCALGeoName) ;
283 fDigitsArr = new TClonesArray("AliEMCALDigit",1000);
285 // Initialising non-linearity parameters
286 if( fNonLinearThreshold != -1 )
287 fEMCALRecoUtils->SetNonLinearityThreshold(fNonLinearThreshold);
288 if( fNonLinearFunc != -1 )
289 fEMCALRecoUtils->SetNonLinearityFunction(fNonLinearFunc);
291 // missalignment function
292 fEMCALRecoUtils->SetPositionAlgorithm(AliEMCALRecoUtils::kPosTowerGlobal);
295 // do not do the eta0 fiducial cut
296 if( fNCellsFromEMCALBorder != -1 )
297 fEMCALRecoUtils->SetNumberOfCellsFromEMCALBorder(fNCellsFromEMCALBorder);
298 fEMCALRecoUtils->SwitchOnNoFiducialBorderInEMCALEta0();
300 // exotic cell rejection
301 if( fExoticCellFraction != -1 )
302 fEMCALRecoUtils->SetExoticCellFractionCut( fExoticCellFraction );
303 if( fExoticCellDiffTime != -1 )
304 fEMCALRecoUtils->SetExoticCellDiffTimeCut( fExoticCellDiffTime );
305 if( fExoticCellMinAmplitude != -1 )
306 fEMCALRecoUtils->SetExoticCellMinAmplitudeCut( fExoticCellMinAmplitude );
308 // Setting track matching parameters ... mass, step size and residual cut
310 fEMCALRecoUtils->SetMass(fMass);
312 fEMCALRecoUtils->SetStep(fStep);
314 // spatial cut based on separate eta/phi or common processing
316 fEMCALRecoUtils->SwitchOnCutEtaPhiSum();
318 fEMCALRecoUtils->SetCutR(fRcut);
319 } else if (fCutEtaPhiSeparate) {
320 fEMCALRecoUtils->SwitchOnCutEtaPhiSeparate();
322 fEMCALRecoUtils->SetCutEta(fEtacut);
324 fEMCALRecoUtils->SetCutPhi(fPhicut);
328 //_____________________________________________________
329 void AliEMCALTenderSupply::ProcessEvent()
333 AliESDEvent *event = fTender->GetEvent();
335 AliError("ESD event ptr = 0, returning");
339 // Initialising parameters once per run number
340 if (fTender->RunChanged()){
342 AliWarning( "Run changed, initializing parameters" );
347 // define what recalib parameters are needed for various switches
348 // this is based on implementation in AliEMCALRecoUtils
349 Bool_t needRecoParam = fReClusterize;
350 Bool_t needBadChannels = fBadCellRemove | fClusterBadChannelCheck | fRecalDistToBadChannels | fReClusterize;
351 Bool_t needRecalib = fCalibrateEnergy | fReClusterize;
352 Bool_t needTimecalib = fCalibrateTime | fReClusterize;
353 Bool_t needMisalign = fRecalClusPos | fReClusterize;
354 Bool_t needClusterizer = fReClusterize;
356 // initiate reco params from OCDB or load some defaults on OCDB failure
357 // will not overwrive, if those have been provided by user
359 Int_t initRC = InitRecParam();
362 AliError("Reco params load from OCDB failed! Defaults loaded.");
364 AliWarning("Reco params loaded from OCDB.");
366 AliWarning("Reco params not loaded from OCDB (user defined or previously failed to load).");
370 if( needBadChannels ){
371 Int_t fInitBC = InitBadChannels();
373 AliError("InitBadChannels returned false, returning");
375 AliWarning("InitBadChannels OK");
377 AliWarning(Form("No external hot channel set: %d - %s", event->GetRunNumber(), fFilepass.Data()));
380 // init recalibration factors
382 Int_t fInitRecalib = InitRecalib();
384 AliError("InitRecalib returned false, returning");
386 AliWarning("InitRecalib OK");
388 AliWarning(Form("No recalibration available: %d - %s", event->GetRunNumber(), fFilepass.Data()));
389 fReCalibCluster = kFALSE;
392 // init time calibration
394 Int_t initTC = InitTimeCalibration();
396 AliError("InitTimeCalibration returned false, returning");
398 AliWarning("InitTimeCalib OK");
400 AliWarning(Form("No external time calibration set: %d - %s", event->GetRunNumber(), fFilepass.Data()));
403 // init misalignment matrix
405 if (!InitMisalignMatrix())
406 AliError("InitMisalignmentMatrix returned false, returning");
408 AliWarning("InitMisalignMatrix OK");
412 if( needClusterizer ) {
413 if (!InitClusterization())
414 AliError("InitClusterization returned false, returning");
416 AliWarning("InitClusterization OK");
420 fEMCALRecoUtils->Print("");
423 // disable implied switches -------------------------------------------------
424 // AliEMCALRecoUtils or clusterizer functions alredy include some
425 // recalibration so based on those implied callibration te switches
426 // here are disabled to avoid duplication
428 // clusterizer does cluster energy recalibration, position recomputation
431 fReCalibCluster = kFALSE;
432 fRecalClusPos = kFALSE;
433 fRecalShowerShape = kFALSE;
436 // CONFIGURE THE RECO UTILS -------------------------------------------------
437 // configure the reco utils
438 // this option does energy recalibration
439 if( fCalibrateEnergy )
440 fEMCALRecoUtils->SwitchOnRecalibration();
442 fEMCALRecoUtils->SwitchOffRecalibration();
444 // allows time calibration
446 fEMCALRecoUtils->SwitchOnTimeRecalibration();
448 fEMCALRecoUtils->SwitchOffTimeRecalibration();
450 // allows to zero bad cells
452 fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
454 fEMCALRecoUtils->SwitchOffBadChannelsRemoval();
456 // distance to bad channel recalibration
457 if( fRecalDistToBadChannels )
458 fEMCALRecoUtils->SwitchOnDistToBadChannelRecalculation();
460 fEMCALRecoUtils->SwitchOffDistToBadChannelRecalculation();
462 // exclude exotic cells
463 if( fRejectExoticCells )
464 fEMCALRecoUtils->SwitchOnRejectExoticCell();
466 fEMCALRecoUtils->SwitchOffRejectExoticCell();
468 // exclude clusters with exotic cells
469 if( fRejectExoticClusters )
470 fEMCALRecoUtils->SwitchOnRejectExoticCluster();
472 fEMCALRecoUtils->SwitchOffRejectExoticCluster();
474 // TODO: not implemented switches
475 // SwitchOnClusterEnergySmearing
476 // SwitchOnRunDepCorrection
478 // START PROCESSING ---------------------------------------------------------
479 // Test if cells present
480 AliESDCaloCells *cells= event->GetEMCALCells();
481 if (cells->GetNumberOfCells()<=0)
484 AliWarning(Form("Number of EMCAL cells = %d, returning", cells->GetNumberOfCells()));
489 AliInfo(Form("Re-calibrate cluster %d\n",fReCalibCluster));
491 // mark the cells not recalibrated in case of selected
492 // time, energy recalibration or bad channel removal
493 if( fCalibrateEnergy || fCalibrateTime || fBadCellRemove )
494 fEMCALRecoUtils->ResetCellsCalibrated();
496 // CELL RECALIBRATION -------------------------------------------------------
497 // cell objects will be updated
498 // the cell calibrations are also processed locally any time those are needed
499 // in case that the cell objects are not to be updated here for later use
503 // includes exotic cell check in the UpdateCells function - is not provided
507 // switch off recalibrations so those are not done multiple times
508 // this is just for safety, the recalibrated flag of cell object
509 // should not allow for farther processing anyways
510 fEMCALRecoUtils->SwitchOffRecalibration();
511 fEMCALRecoUtils->SwitchOffTimeRecalibration();
517 // RECLUSTERIZATION ---------------------------------------------------------
525 // Store good clusters
526 TClonesArray *clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
528 clusArr = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
530 AliWarning(Form("No cluster array, number of cells in event = %d, returning", cells->GetNumberOfCells()));
534 // PROCESS SINGLE CLUSTER RECALIBRATION -------------------------------------
535 // now go through clusters one by one and process separate correction
536 // as those were defined or not
537 Int_t nclusters = clusArr->GetEntriesFast();
538 for (Int_t icluster=0; icluster < nclusters; ++icluster)
540 AliVCluster *clust = static_cast<AliVCluster*>(clusArr->At(icluster));
543 if (!clust->IsEMCAL())
546 // REMOVE CLUSTERS WITH BAD CELLS -----------------------------
547 if( fClusterBadChannelCheck )
549 // careful, the the ClusterContainsBadChannel is dependent on
550 // SwitchOnBadChannelsRemoval, switching it ON automatically
551 // and returning to original value after processing
552 Bool_t badRemoval = fEMCALRecoUtils->IsBadChannelsRemovalSwitchedOn();
553 fEMCALRecoUtils->SwitchOnBadChannelsRemoval();
555 Bool_t badResult = fEMCALRecoUtils->ClusterContainsBadChannel(fEMCALGeo, clust->GetCellsAbsId(), clust->GetNCells());
557 // switch the bad channels removal back
559 fEMCALRecoUtils->SwitchOffBadChannelsRemoval();
563 delete clusArr->RemoveAt(icluster);
564 continue; //TODO is it really needed to remove it? Or should we flag it?
568 // REMOVE EXOTIC CLUSTERS -------------------------------------
569 // does process local cell recalibration energy and time without replacing
570 // the global cell values, in case of no cell recalib done yet
571 if( fRejectExoticClusters )
573 // careful, the IsExoticCluster is dependent on
574 // SwitchOnRejectExoticCell, switching it ON automatically
575 // and returning to original value after processing
576 Bool_t exRemoval = fEMCALRecoUtils->IsRejectExoticCell();
577 fEMCALRecoUtils->SwitchOnRejectExoticCell();
579 // get bunch crossing
580 Int_t bunchCrossNo = fTender->GetEvent()->GetBunchCrossNumber();
582 Bool_t exResult = fEMCALRecoUtils->IsExoticCluster(clust, cells, bunchCrossNo );
584 // switch the exotic channels removal back
586 fEMCALRecoUtils->SwitchOffRejectExoticCell();
590 delete clusArr->RemoveAt(icluster);
591 continue; //TODO is it really needed to remove it? Or should we flag it?
595 // FIDUCIAL CUT -----------------------------------------------
598 // depends on SetNumberOfCellsFromEMCALBorder
599 // SwitchOnNoFiducialBorderInEMCALEta0
600 if (!fEMCALRecoUtils->CheckCellFiducialRegion(fEMCALGeo, clust, cells)){
601 delete clusArr->RemoveAt(icluster);
602 continue; //TODO it would be nice to store the distance
606 // CLUSTER ENERGY ---------------------------------------------
607 // does process local cell recalibration energy and time without replacing
608 // the global cell values, in case of no cell recalib done yet
609 if( fReCalibCluster )
610 fEMCALRecoUtils->RecalibrateClusterEnergy(fEMCALGeo, clust, cells);
612 // CLUSTER POSITION -------------------------------------------
613 // does process local cell energy recalibration, if enabled and cells
614 // not calibratied yet
616 fEMCALRecoUtils->RecalculateClusterPosition(fEMCALGeo, cells, clust);
618 // SHOWER SHAPE -----------------------------------------------
619 if( fRecalShowerShape )
620 fEMCALRecoUtils->RecalculateClusterShowerShapeParameters(fEMCALGeo, cells, clust);
622 // NONLINEARITY -----------------------------------------------
623 if( fDoNonLinearity )
625 Float_t correctedEnergy = fEMCALRecoUtils->CorrectClusterEnergyLinearity(clust);
626 clust->SetE(correctedEnergy);
629 // DISTANCE TO BAD CHANNELS -----------------------------------
630 if( fRecalDistToBadChannels )
631 fEMCALRecoUtils->RecalculateClusterDistanceToBadChannel(fEMCALGeo, cells, clust);
639 // TRACK MATCHING -----------------------------------------------------------
640 if (!TGeoGlobalMagField::Instance()->GetField())
642 event->InitMagneticField();
645 fEMCALRecoUtils->FindMatches(event,0x0,fEMCALGeo);
646 fEMCALRecoUtils->SetClusterMatchedToTrack(event);
647 fEMCALRecoUtils->SetTracksMatchedToCluster(event);
650 //_____________________________________________________
651 Bool_t AliEMCALTenderSupply::InitMisalignMatrix()
653 // Initialising misalignment matrices
655 AliESDEvent *event = fTender->GetEvent();
661 AliInfo("Misalignment matrix already set");
666 AliInfo("Initialising misalignment matrix");
668 if (fLoadGeomMatrices) {
669 for(Int_t mod=0; mod < fEMCALGeo->GetNumberOfSuperModules(); ++mod)
671 if (fEMCALMatrix[mod]){
673 fEMCALMatrix[mod]->Print();
674 fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
677 fGeomMatrixSet = kTRUE;
681 Int_t runGM = event->GetRunNumber();
684 if(fMisalignSurvey == kdefault)
685 { //take default alignment corresponding to run no
686 AliOADBContainer emcalgeoCont(Form("emcal"));
687 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
688 mobj=(TObjArray*)emcalgeoCont.GetObject(runGM,"EmcalMatrices");
691 if(fMisalignSurvey == kSurveybyS)
692 { //take alignment at sector level
693 if (runGM <= 140000) { //2010 data
694 AliOADBContainer emcalgeoCont(Form("emcal2010"));
695 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
696 mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
698 } else if (runGM>140000)
699 { // 2011 LHC11a pass1 data
700 AliOADBContainer emcalgeoCont(Form("emcal2011"));
701 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
702 mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byS");
706 if(fMisalignSurvey == kSurveybyM)
707 { //take alignment at module level
708 if (runGM <= 140000) { //2010 data
709 AliOADBContainer emcalgeoCont(Form("emcal2010"));
710 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
711 mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey10");
713 } else if (runGM>140000)
714 { // 2011 LHC11a pass1 data
715 AliOADBContainer emcalgeoCont(Form("emcal2011"));
716 emcalgeoCont.InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALlocal2master.root",Form("AliEMCALgeo"));
717 mobj=(TObjArray*)emcalgeoCont.GetObject(100,"survey11byM");
723 AliFatal("Geometry matrix array not found");
727 for(Int_t mod=0; mod < (fEMCALGeo->GetEMCGeometry())->GetNumberOfSuperModules(); mod++)
729 fEMCALMatrix[mod] = (TGeoHMatrix*) mobj->At(mod);
730 fEMCALGeo->SetMisalMatrix(fEMCALMatrix[mod],mod);
731 fEMCALMatrix[mod]->Print();
737 //_____________________________________________________
738 Int_t AliEMCALTenderSupply::InitBadChannels()
740 // Initialising bad channel maps
741 AliESDEvent *event = fTender->GetEvent();
746 AliInfo("Initialising Bad channel map");
748 // init default maps first
749 if( !fEMCALRecoUtils->GetEMCALBadChannelStatusMapArray() )
750 fEMCALRecoUtils->InitEMCALBadChannelStatusMap() ;
752 Int_t runBC = event->GetRunNumber();
754 AliOADBContainer *contBC = new AliOADBContainer("");
756 { //if fBasePath specified in the ->SetBasePath()
757 if (fDebugLevel>0) AliInfo(Form("Loading Bad Channels OADB from given path %s",fBasePath.Data()));
759 TFile *fbad=new TFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"read");
760 if (!fbad || fbad->IsZombie())
762 AliFatal(Form("EMCALBadChannels.root was not found in the path provided: %s",fBasePath.Data()));
766 if (fbad) delete fbad;
768 contBC->InitFromFile(Form("%s/EMCALBadChannels.root",fBasePath.Data()),"AliEMCALBadChannels");
771 { // Else choose the one in the $ALICE_ROOT directory
772 if (fDebugLevel>0) AliInfo("Loading Bad Channels OADB from $ALICE_ROOT/OADB/EMCAL");
774 TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","read");
775 if (!fbad || fbad->IsZombie())
777 AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root was not found");
781 if (fbad) delete fbad;
783 contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALBadChannels.root","AliEMCALBadChannels");
786 TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
789 AliError(Form("No external hot channel set for run number: %d", runBC));
793 TObjArray *arrayBCpass=(TObjArray*)arrayBC->FindObject(fFilepass); // Here, it looks for a specific pass
796 AliError(Form("No external hot channel set for: %d -%s", runBC,fFilepass.Data()));
800 if (fDebugLevel>0) arrayBCpass->Print();
802 Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
803 for (Int_t i=0; i<sms; ++i)
805 TH2I *h = fEMCALRecoUtils->GetEMCALChannelStatusMap(i);
808 h=(TH2I*)arrayBCpass->FindObject(Form("EMCALBadChannelMap_Mod%d",i));
812 AliError(Form("Can not get EMCALBadChannelMap_Mod%d",i));
816 fEMCALRecoUtils->SetEMCALChannelStatusMap(i,h);
821 //_____________________________________________________
822 Int_t AliEMCALTenderSupply::InitRecalib()
824 // Initialising recalibration factors.
826 AliESDEvent *event = fTender->GetEvent();
831 AliInfo("Initialising recalibration factors");
833 // init default maps first
834 if( !fEMCALRecoUtils->GetEMCALRecalibrationFactorsArray() )
835 fEMCALRecoUtils->InitEMCALRecalibrationFactors() ;
837 Int_t runRC = event->GetRunNumber();
839 AliOADBContainer *contRF=new AliOADBContainer("");
841 { //if fBasePath specified in the ->SetBasePath()
842 if (fDebugLevel>0) AliInfo(Form("Loading Recalib OADB from given path %s",fBasePath.Data()));
844 TFile *fRecalib= new TFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"read");
845 if (!fRecalib || fRecalib->IsZombie())
847 AliFatal(Form("EMCALRecalib.root not found in %s",fBasePath.Data()));
851 if (fRecalib) delete fRecalib;
853 contRF->InitFromFile(Form("%s/EMCALRecalib.root",fBasePath.Data()),"AliEMCALRecalib");
856 { // Else choose the one in the $ALICE_ROOT directory
857 if (fDebugLevel>0) AliInfo("Loading Recalib OADB from $ALICE_ROOT/OADB/EMCAL");
859 TFile *fRecalib= new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","read");
860 if (!fRecalib || fRecalib->IsZombie())
862 AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root was not found");
866 if (fRecalib) delete fRecalib;
868 contRF->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALRecalib.root","AliEMCALRecalib");
871 TObjArray *recal=(TObjArray*)contRF->GetObject(runRC); //GetObject(int runnumber)
874 AliError(Form("No Objects for run: %d",runRC));
878 TObjArray *recalpass=(TObjArray*)recal->FindObject(fFilepass);
881 AliError(Form("No Objects for run: %d - %s",runRC,fFilepass.Data()));
885 TObjArray *recalib=(TObjArray*)recalpass->FindObject("Recalib");
888 AliError(Form("No Recalib histos found for %d - %s",runRC,fFilepass.Data()));
892 if (fDebugLevel>0) recalib->Print();
894 Int_t sms = fEMCALGeo->GetEMCGeometry()->GetNumberOfSuperModules();
895 for (Int_t i=0; i<sms; ++i)
897 TH2F *h = fEMCALRecoUtils->GetEMCALChannelRecalibrationFactors(i);
900 h = (TH2F*)recalib->FindObject(Form("EMCALRecalFactors_SM%d",i));
903 AliError(Form("Could not load EMCALRecalFactors_SM%d",i));
907 fEMCALRecoUtils->SetEMCALChannelRecalibrationFactors(i,h);
912 //_____________________________________________________
913 Int_t AliEMCALTenderSupply::InitTimeCalibration()
915 // Initialising bad channel maps
916 AliESDEvent *event = fTender->GetEvent();
921 AliInfo("Initialising time calibration map");
923 // init default maps first
924 if( !fEMCALRecoUtils->GetEMCALTimeRecalibrationFactorsArray() )
925 fEMCALRecoUtils->InitEMCALTimeRecalibrationFactors() ;
927 Int_t runBC = event->GetRunNumber();
929 AliOADBContainer *contBC = new AliOADBContainer("");
931 { //if fBasePath specified in the ->SetBasePath()
932 if (fDebugLevel>0) AliInfo(Form("Loading time calibration OADB from given path %s",fBasePath.Data()));
934 TFile *fbad=new TFile(Form("%s/EMCALTimeCalib.root",fBasePath.Data()),"read");
935 if (!fbad || fbad->IsZombie())
937 AliFatal(Form("EMCALTimeCalib.root was not found in the path provided: %s",fBasePath.Data()));
941 if (fbad) delete fbad;
943 contBC->InitFromFile(Form("%s/EMCALTimeCalib.root",fBasePath.Data()),"AliEMCALTimeCalib");
946 { // Else choose the one in the $ALICE_ROOT directory
947 if (fDebugLevel>0) AliInfo("Loading time calibration OADB from $ALICE_ROOT/OADB/EMCAL");
949 TFile *fbad=new TFile("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root","read");
950 if (!fbad || fbad->IsZombie())
952 AliFatal("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root was not found");
956 if (fbad) delete fbad;
958 contBC->InitFromFile("$ALICE_ROOT/OADB/EMCAL/EMCALTimeCalib.root","AliEMCALTimeCalib");
961 TObjArray *arrayBC=(TObjArray*)contBC->GetObject(runBC);
964 AliError(Form("No external time calibration set for run number: %d", runBC));
968 TObjArray *arrayBCpass=(TObjArray*)arrayBC->FindObject(fFilepass); // Here, it looks for a specific pass
971 AliError(Form("No external time calibration set for: %d -%s", runBC,fFilepass.Data()));
975 if (fDebugLevel>0) arrayBCpass->Print();
977 for( Int_t i = 0; i < 4; i++ )
979 TH1F *h = fEMCALRecoUtils->GetEMCALChannelTimeRecalibrationFactors( i );
983 h = (TH1F*)arrayBCpass->FindObject(Form("hAllTimeAvBC%d",i));
987 AliError(Form("Can not get hAllTimeAvBC%d",i));
991 fEMCALRecoUtils->SetEMCALChannelTimeRecalibrationFactors(i,h);
996 //_____________________________________________________
997 void AliEMCALTenderSupply::UpdateCells()
999 //Remove bad cells from the cell list
1000 //Recalibrate energy and time cells
1001 //This is required for later reclusterization
1003 AliESDCaloCells *cells = fTender->GetEvent()->GetEMCALCells();
1004 Int_t bunchCrossNo = fTender->GetEvent()->GetBunchCrossNumber();
1006 fEMCALRecoUtils->RecalibrateCells(cells, bunchCrossNo);
1008 // remove exotic cells - loop through cells and zero the exotic ones
1009 // just like with bad cell rejection in reco utils (inside RecalibrateCells)
1010 if( fRejectExoticCells )
1013 Bool_t isExot = kFALSE;
1015 // loop through cells
1016 Int_t nEMcell = cells->GetNumberOfCells() ;
1017 for (Int_t iCell = 0; iCell < nEMcell; iCell++)
1019 absId = cells->GetCellNumber(iCell);
1021 isExot = fEMCALRecoUtils->IsExoticCell( absId, cells, bunchCrossNo );
1024 cells->SetCell( iCell, absId, 0.0, 0.0 );
1026 } // reject exotic cells
1031 //_____________________________________________________
1032 Int_t AliEMCALTenderSupply::InitRecParam()
1034 // Initalize the reconstructor parameters from OCDB
1035 // load some default on OCDB failure
1038 AliCDBManager * man = 0x0 ;
1039 TObjArray * arr = 0x0 ;
1040 AliEMCALRecParam * pars = 0x0 ;
1041 AliCDBEntry * entry = 0x0;
1042 const AliESDRun * run = 0x0 ;
1045 // clean the previous reco params, if those came from OCDB
1046 // we do not want to erase user provided params, do we
1047 if( fRecoParamsOCDBLoaded ){
1048 if( fRecParam != 0 ){
1052 // zero the OCDB loaded flag
1053 fRecoParamsOCDBLoaded = kFALSE;
1056 // exit if reco params exist (probably shipped by the user already)
1057 if( fRecParam != 0 )
1061 AliInfo("Initialize the recParam");
1064 run = fTender->GetEvent()->GetESDRun();
1065 beamType = run->GetBeamType();
1066 runNum = fTender->GetEvent()->GetRunNumber();
1068 // OCDB manager should already exist
1069 // and have a default storage defined (done by AliTender)
1070 man = AliCDBManager::Instance();
1072 // load the file data
1074 entry = man->Get("EMCAL/Calib/RecoParam", runNum);
1077 arr = (TObjArray*)(entry->GetObject());
1080 // load given parameters based on beam type
1081 if( beamType == "A-A" ){
1082 if( fDebugLevel > 0 )
1083 AliInfo( "Initializing A-A reco params." );
1084 pars = (AliEMCALRecParam*)arr->FindObject( "High Flux - Pb+Pb" );
1087 if( fDebugLevel > 0 )
1088 AliInfo( "Initializing p-p reco params." );
1089 pars = (AliEMCALRecParam*)arr->FindObject( "Low Flux - p+p" );
1092 // set the parameters, if found
1094 if( fDebugLevel > 0 )
1095 AliInfo( "OCDB reco params set." );
1098 fRecoParamsOCDBLoaded = kTRUE;
1105 // set some defaults if OCDB did not succede
1106 if( !fRecoParamsOCDBLoaded ){
1107 fRecParam = new AliEMCALRecParam();
1109 if ( beamType == "A-A"){
1110 fRecParam->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv2);
1111 fRecParam->SetClusteringThreshold(0.100);
1112 fRecParam->SetMinECut(0.050);
1116 fRecParam->SetClusterizerFlag(AliEMCALRecParam::kClusterizerv1);
1117 fRecParam->SetClusteringThreshold(0.100);
1118 fRecParam->SetMinECut(0.050);
1122 if( fRecoParamsOCDBLoaded )
1128 //_____________________________________________________
1129 Bool_t AliEMCALTenderSupply::InitClusterization()
1131 // Initialing clusterization/unfolding algorithm and set all the needed parameters.
1133 AliESDEvent *event = fTender->GetEvent();
1138 AliInfo(Form("Initialising reclustering parameters: Clusterizer type: %d",fRecParam->GetClusterizerFlag()));
1140 //---setup clusterizer
1141 delete fClusterizer;
1142 if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
1143 fClusterizer = new AliEMCALClusterizerv1 (fEMCALGeo);
1144 else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
1145 fClusterizer = new AliEMCALClusterizerv2(fEMCALGeo);
1146 else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN)
1148 AliEMCALClusterizerNxN *clusterizer = new AliEMCALClusterizerNxN(fEMCALGeo);
1149 clusterizer->SetNRowDiff(fRecParam->GetNRowDiff());
1150 clusterizer->SetNColDiff(fRecParam->GetNColDiff());
1151 fClusterizer = clusterizer;
1155 AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
1159 // Set the clustering parameters
1160 fClusterizer->SetECAClusteringThreshold( fRecParam->GetClusteringThreshold() );
1161 fClusterizer->SetECALogWeight ( fRecParam->GetW0() );
1162 fClusterizer->SetMinECut ( fRecParam->GetMinECut() );
1163 fClusterizer->SetUnfolding ( fRecParam->GetUnfold() );
1164 fClusterizer->SetECALocalMaxCut ( fRecParam->GetLocMaxCut() );
1165 fClusterizer->SetTimeCut ( fRecParam->GetTimeCut() );
1166 fClusterizer->SetTimeMin ( fRecParam->GetTimeMin() );
1167 fClusterizer->SetTimeMax ( fRecParam->GetTimeMax() );
1168 fClusterizer->SetInputCalibrated ( kTRUE );
1169 fClusterizer->SetJustClusters ( kTRUE );
1171 // In case of unfolding after clusterization is requested, set the corresponding parameters
1172 if (fRecParam->GetUnfold())
1174 for (Int_t i = 0; i < 8; ++i)
1176 fClusterizer->SetSSPars(i, fRecParam->GetSSPars(i));
1178 for (Int_t i = 0; i < 3; ++i)
1180 fClusterizer->SetPar5 (i, fRecParam->GetPar5(i));
1181 fClusterizer->SetPar6 (i, fRecParam->GetPar6(i));
1183 fClusterizer->InitClusterUnfolding();
1186 fClusterizer->SetDigitsArr(fDigitsArr);
1187 fClusterizer->SetOutput(0);
1188 fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
1192 //_____________________________________________________
1193 void AliEMCALTenderSupply::FillDigitsArray()
1195 // Fill digits from cells to a TClonesArray.
1197 AliESDEvent *event = fTender->GetEvent();
1201 fDigitsArr->Clear("C");
1202 AliESDCaloCells *cells = event->GetEMCALCells();
1203 Int_t ncells = cells->GetNumberOfCells();
1204 for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell)
1206 Double_t cellAmplitude=0, cellTime=0;
1207 Short_t cellNumber=0;
1209 if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime) != kTRUE)
1212 // Do not add if already too low (some cells set to 0 if bad channels)
1213 if (cellAmplitude < fRecParam->GetMinECut() )
1216 // If requested, do not include exotic cells
1217 if (fEMCALRecoUtils->IsExoticCell(cellNumber,cells,event->GetBunchCrossNumber()))
1220 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
1221 digit->SetId(cellNumber);
1222 digit->SetTime(cellTime);
1223 digit->SetTimeR(cellTime);
1224 digit->SetIndexInList(idigit);
1225 digit->SetType(AliEMCALDigit::kHG);
1226 digit->SetAmplitude(cellAmplitude);
1231 //_____________________________________________________
1232 void AliEMCALTenderSupply::Clusterize()
1236 fClusterizer->Digits2Clusters("");
1239 //_____________________________________________________
1240 void AliEMCALTenderSupply::UpdateClusters()
1242 // Update ESD cluster list.
1244 AliESDEvent *event = fTender->GetEvent();
1248 TClonesArray *clus = dynamic_cast<TClonesArray*>(event->FindListObject("caloClusters"));
1250 clus = dynamic_cast<TClonesArray*>(event->FindListObject("CaloClusters"));
1253 AliError(" Null pointer to calo clusters array, returning");
1257 Int_t nents = clus->GetEntriesFast();
1258 for (Int_t i=0; i < nents; ++i)
1260 AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->At(i));
1264 delete clus->RemoveAt(i);
1269 RecPoints2Clusters(clus);
1273 //_____________________________________________________
1274 void AliEMCALTenderSupply::RecPoints2Clusters(TClonesArray *clus)
1276 // Convert AliEMCALRecoPoints to AliESDCaloClusters.
1277 // Cluster energy, global position, cells and their amplitude fractions are restored.
1279 AliESDEvent *event = fTender->GetEvent();
1283 Int_t ncls = fClusterArr->GetEntriesFast();
1284 for(Int_t i=0, nout=clus->GetEntriesFast(); i < ncls; ++i)
1286 AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
1288 Int_t ncellsTrue = 0;
1289 const Int_t ncells = recpoint->GetMultiplicity();
1290 UShort_t absIds[ncells];
1291 Double32_t ratios[ncells];
1292 Int_t *dlist = recpoint->GetDigitsList();
1293 Float_t *elist = recpoint->GetEnergiesList();
1294 for (Int_t c = 0; c < ncells; ++c)
1296 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
1297 absIds[ncellsTrue] = digit->GetId();
1298 ratios[ncellsTrue] = elist[c]/digit->GetAmplitude();
1299 if (ratios[ncellsTrue] < 0.001)
1306 AliWarning("Skipping cluster with no cells");
1310 // calculate new cluster position
1312 recpoint->GetGlobalPosition(gpos);
1316 AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(clus->New(nout++));
1318 c->SetType(AliVCluster::kEMCALClusterv1);
1319 c->SetE(recpoint->GetEnergy());
1321 c->SetNCells(ncellsTrue);
1322 c->SetDispersion(recpoint->GetDispersion());
1323 c->SetEmcCpvDistance(-1); //not yet implemented
1324 c->SetChi2(-1); //not yet implemented
1325 c->SetTOF(recpoint->GetTime()) ; //time-of-flight
1326 c->SetNExMax(recpoint->GetNExMax()); //number of local maxima
1327 Float_t elipAxis[2];
1328 recpoint->GetElipsAxis(elipAxis);
1329 c->SetM02(elipAxis[0]*elipAxis[0]) ;
1330 c->SetM20(elipAxis[1]*elipAxis[1]) ;
1331 AliESDCaloCluster *cesd = static_cast<AliESDCaloCluster*>(c);
1332 cesd->SetCellsAbsId(absIds);
1333 cesd->SetCellsAmplitudeFraction(ratios);
1337 //_____________________________________________________
1338 void AliEMCALTenderSupply::GetPass()
1340 // Get passx from filename.
1342 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1343 fInputTree = mgr->GetTree();
1347 AliError("Pointer to tree = 0, returning");
1351 fInputFile = fInputTree->GetCurrentFile();
1353 AliError("Null pointer input file, returning");
1357 TString fname(fInputFile->GetName());
1358 if (fname.Contains("pass1")) fFilepass = TString("pass1");
1359 else if (fname.Contains("pass2")) fFilepass = TString("pass2");
1360 else if (fname.Contains("pass3")) fFilepass = TString("pass3");
1363 AliError(Form("Pass number string not found: %s", fname.Data()));