charge selection on trigger particle
[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliAnalysisTaskEMCALClusterizeFast.cxx
CommitLineData
2f7259cf 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
cd231d42 16//_________________________________________________________________________
17// This analysis provides a new list of clusters to be used in other analysis
18//
19// Author: Constantin Loizides, Salvatore Aiola
20// Adapted from analysis class from Deepa Thomas
21//
22// $Id$
23//_________________________________________________________________________
2f7259cf 24
25// --- Root ---
26#include <TClonesArray.h>
27#include <TGeoManager.h>
28#include <TObjArray.h>
29#include <TString.h>
30#include <TTree.h>
7a7d5a64 31
32// --- AliRoot ---
2f7259cf 33#include "AliAODCaloCluster.h"
34#include "AliAODEvent.h"
35#include "AliAnalysisManager.h"
36#include "AliCDBEntry.h"
37#include "AliCDBManager.h"
38#include "AliCaloCalibPedestal.h"
39#include "AliEMCALAfterBurnerUF.h"
40#include "AliEMCALCalibData.h"
41#include "AliEMCALClusterizerNxN.h"
42#include "AliEMCALClusterizerv1.h"
3786256f 43#include "AliEMCALClusterizerv2.h"
7a7d5a64 44#include "AliEMCALClusterizerFixedWindow.h"
2f7259cf 45#include "AliEMCALDigit.h"
46#include "AliEMCALGeometry.h"
47#include "AliEMCALRecParam.h"
48#include "AliEMCALRecPoint.h"
49#include "AliEMCALRecoUtils.h"
50#include "AliESDEvent.h"
6b7ed7bd 51#include "AliInputEventHandler.h"
2f7259cf 52#include "AliLog.h"
53
54#include "AliAnalysisTaskEMCALClusterizeFast.h"
55
56ClassImp(AliAnalysisTaskEMCALClusterizeFast)
57
58//________________________________________________________________________
2431b20f 59AliAnalysisTaskEMCALClusterizeFast::AliAnalysisTaskEMCALClusterizeFast()
60 : AliAnalysisTaskSE(),
349bf69e 61 fRun(-1),
2431b20f 62 fDigitsArr(0),
63 fClusterArr(0),
64 fRecParam(0),
65 fClusterizer(0),
66 fUnfolder(0),
67 fJustUnfold(kFALSE),
68 fGeomName(),
69 fGeomMatrixSet(kFALSE),
70 fLoadGeomMatrices(kFALSE),
71 fOCDBpath(),
72 fCalibData(0),
73 fPedestalData(0),
74 fOutputAODBranch(0),
75 fOutputAODBrName(),
76 fRecoUtils(0),
77 fLoadCalib(0),
4cf9204b 78 fLoadPed(0),
e50a444f 79 fAttachClusters(1),
95f69c66 80 fRecalibOnly(0),
7a7d5a64 81 fSubBackground(0),
82 fCreatePattern(0),
e50a444f 83 fOverwrite(0),
7a7d5a64 84 fNewClusterArrayName("newCaloClusters"),
85 fNPhi(4),
86 fNEta(4),
3c56da8f 87 fShiftPhi(2),
88 fShiftEta(2),
89 fTRUShift(0),
0408f54f 90 fClusterizeFastORs(0),
be268574 91 fTrackName(),
92 fCutL0Times(kTRUE)
2431b20f 93{
94 // Constructor
c6d60de7 95
96 for(Int_t i = 0; i < 12; ++i)
97 fGeomMatrix[i] = 0;
2431b20f 98}
99
100//________________________________________________________________________
2f7259cf 101AliAnalysisTaskEMCALClusterizeFast::AliAnalysisTaskEMCALClusterizeFast(const char *name)
102 : AliAnalysisTaskSE(name),
349bf69e 103 fRun(-1),
2f7259cf 104 fDigitsArr(0),
105 fClusterArr(0),
106 fRecParam(new AliEMCALRecParam),
107 fClusterizer(0),
108 fUnfolder(0),
109 fJustUnfold(kFALSE),
110 fGeomName("EMCAL_FIRSTYEARV1"),
111 fGeomMatrixSet(kFALSE),
112 fLoadGeomMatrices(kFALSE),
113 fOCDBpath(),
114 fCalibData(0),
115 fPedestalData(0),
116 fOutputAODBranch(0),
117 fOutputAODBrName(),
10d33986 118 fRecoUtils(0),
119 fLoadCalib(0),
4cf9204b 120 fLoadPed(0),
e50a444f 121 fAttachClusters(1),
95f69c66 122 fRecalibOnly(0),
7a7d5a64 123 fSubBackground(0),
124 fCreatePattern(0),
e50a444f 125 fOverwrite(0),
7a7d5a64 126 fNewClusterArrayName("newCaloClusters"),
127 fNPhi(4),
128 fNEta(4),
3c56da8f 129 fShiftPhi(2),
130 fShiftEta(2),
131 fTRUShift(0),
0408f54f 132 fClusterizeFastORs(0),
be268574 133 fTrackName(),
134 fCutL0Times(kTRUE)
2f7259cf 135{
136 // Constructor
137
0a7eb283 138 fBranchNames = "ESD:AliESDHeader.,AliESDRun.,EMCALCells.,EMCALTrigger. AOD:header,emcalCells";
95f69c66 139 for(Int_t i = 0; i < 12; ++i)
2f7259cf 140 fGeomMatrix[i] = 0;
141}
142
143//________________________________________________________________________
144AliAnalysisTaskEMCALClusterizeFast::~AliAnalysisTaskEMCALClusterizeFast()
145{
146 // Destructor.
147
2f7259cf 148 delete fClusterizer;
149 delete fUnfolder;
150 delete fRecoUtils;
151}
152
153//-------------------------------------------------------------------
154void AliAnalysisTaskEMCALClusterizeFast::UserCreateOutputObjects()
155{
156 // Create output objects.
157
158 if (!fOutputAODBrName.IsNull()) {
159 fOutputAODBranch = new TClonesArray("AliAODCaloCluster", 0);
160 fOutputAODBranch->SetName(fOutputAODBrName);
161 AddAODBranch("TClonesArray", &fOutputAODBranch);
162 AliInfo(Form("Created Branch: %s",fOutputAODBrName.Data()));
163 }
164}
165
166//________________________________________________________________________
167void AliAnalysisTaskEMCALClusterizeFast::UserExec(Option_t *)
168{
169 // Main loop, called for each event
170
171 // remove the contents of output list set in the previous event
172 if (fOutputAODBranch)
173 fOutputAODBranch->Clear("C");
174
175 AliESDEvent *esdevent = dynamic_cast<AliESDEvent*>(InputEvent());
176 AliAODEvent *aodevent = dynamic_cast<AliAODEvent*>(InputEvent());
177
178 if (!esdevent&&!aodevent) {
179 Error("UserExec","Event not available");
180 return;
181 }
182
183 LoadBranches();
6b7ed7bd 184
185 UInt_t offtrigger = 0;
186 if (esdevent) {
9809ed8c 187 UInt_t mask1 = esdevent->GetESDRun()->GetDetectorsInDAQ();
188 UInt_t mask2 = esdevent->GetESDRun()->GetDetectorsInReco();
189 Bool_t desc1 = (mask1 >> 18) & 0x1;
190 Bool_t desc2 = (mask2 >> 18) & 0x1;
191 if (desc1==0 || desc2==0) { //AliDAQ::OfflineModuleName(180=="EMCAL"
192 AliError(Form("EMCAL not in DAQ/RECO: %u (%u)/%u (%u)",
193 mask1, esdevent->GetESDRun()->GetDetectorsInReco(),
194 mask2, esdevent->GetESDRun()->GetDetectorsInDAQ()));
6b7ed7bd 195 return;
196 }
197 AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
198 offtrigger = ((AliInputEventHandler*)(am->GetInputEventHandler()))->IsEventSelected();
199 } else if (aodevent) {
200 offtrigger = aodevent->GetHeader()->GetOfflineTrigger();
201 }
179c3e32 202
2ecad9b5 203 if (esdevent) {
179c3e32 204 AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
205 Bool_t mcmode = 0;
206 if (am->GetMCtruthEventHandler())
207 mcmode = 1;
208 if (!mcmode) {
179c3e32 209 if (offtrigger & AliVEvent::kFastOnly) {
210 AliWarning(Form("EMCAL not in fast only partition"));
211 return;
212 }
213 }
6b7ed7bd 214 }
2f7259cf 215
216 Init();
217
218 if (fJustUnfold) {
219 AliWarning("Unfolding not implemented");
84ccad86 220 return;
221 }
222
95f69c66 223 FillDigitsArray();
84ccad86 224
95f69c66 225 if (fRecalibOnly) {
226 UpdateCells();
227 return; // not requested to run clusterizer
84ccad86 228 }
7099ec08 229
95f69c66 230 Clusterize();
231 UpdateCells();
232 UpdateClusters();
2431b20f 233
95f69c66 234 if (fOutputAODBranch)
235 RecPoints2Clusters(fOutputAODBranch);
2431b20f 236}
237
238//________________________________________________________________________
95f69c66 239void AliAnalysisTaskEMCALClusterizeFast::Clusterize()
2431b20f 240{
95f69c66 241 // Clusterize
7099ec08 242
95f69c66 243 if (fSubBackground) {
244 fClusterizer->SetInputCalibrated(kTRUE);
245 fClusterizer->SetCalibrationParameters(0);
2f7259cf 246 }
d3b5a5ec 247
95f69c66 248 fClusterizer->Digits2Clusters("");
249 if (fSubBackground) {
4889ed8b 250 if (fCalibData) {
251 fClusterizer->SetInputCalibrated(kFALSE);
252 fClusterizer->SetCalibrationParameters(fCalibData);
253 }
4cf9204b 254 }
4cf9204b 255}
256
257//________________________________________________________________________
95f69c66 258void AliAnalysisTaskEMCALClusterizeFast::FillDigitsArray()
2f7259cf 259{
6265e793 260 // Fill digits array
55a8e3e4 261
6265e793 262 fDigitsArr->Clear("C");
263
5133b922 264 if (fCreatePattern) { // Fill digits from a pattern
7a7d5a64 265 AliEMCALGeometry *fGeom = AliEMCALGeometry::GetInstance(fGeomName);
7a7d5a64 266 Int_t maxd = fGeom->GetNCells() / 4;
6265e793 267 for (Int_t idigit = 0; idigit < maxd; idigit++){
7a7d5a64 268 if (idigit % 24 == 12) idigit += 12;
269 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
270 digit->SetId(idigit * 4);
271 digit->SetTime(600);
272 digit->SetTimeR(600);
273 digit->SetIndexInList(idigit);
274 digit->SetType(AliEMCALDigit::kHG);
275 digit->SetAmplitude(0.1);
b0e354fe 276 }
06ca320f 277
5133b922 278 } else if (fClusterizeFastORs) { // Fill digits from FastORs
6265e793 279
280 AliEMCALGeometry *fGeom = AliEMCALGeometry::GetInstance(fGeomName);
6265e793 281
35799643 282 AliVCaloTrigger *triggers = InputEvent()->GetCaloTrigger("EMCAL");
6265e793 283
284 if (!triggers || !(triggers->GetEntries() > 0))
285 return;
286
287 Int_t idigit = 0;
288 triggers->Reset();
289
5133b922 290 while ((triggers->Next())) {
6265e793 291 Float_t triggerAmplitude = 0;
292 triggers->GetAmplitude(triggerAmplitude);
293 if (triggerAmplitude <= 0)
294 continue;
295
296 Int_t triggerTime = 0;
297 Int_t ntimes = 0;
298 triggers->GetNL0Times(ntimes);
be268574 299 if (!(ntimes > 0) && fCutL0Times)
300 continue;
301
302 Int_t trgtimes[25];
303 triggers->GetL0Times(trgtimes);
304 triggerTime = trgtimes[0];
305
6265e793 306 Int_t triggerCol = 0, triggerRow = 0;
307 triggers->GetPosition(triggerCol, triggerRow);
308
309 Int_t find = -1;
310 fGeom->GetAbsFastORIndexFromPositionInEMCAL(triggerCol, triggerRow, find);
311
312 if (find<0)
313 continue;
314
315 Int_t cidx[4] = {-1};
316 Bool_t ret = fGeom->GetCellIndexFromFastORIndex(find, cidx);
317
318 if (!ret)
319 continue;
320
5133b922 321 for (Int_t idxpos = 0; idxpos < 4; idxpos++) {
6265e793 322 Int_t triggerNumber = cidx[idxpos];
323 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
324 digit->SetId(triggerNumber);
325 digit->SetTime(triggerTime);
326 digit->SetTimeR(triggerTime);
327 digit->SetIndexInList(idigit);
328 digit->SetType(AliEMCALDigit::kHG);
329 digit->SetAmplitude(triggerAmplitude);
330 idigit++;
331 }
332 }
06ca320f 333
334 } else { // Fill digits from cells.
7a7d5a64 335
7a7d5a64 336 AliVCaloCells *cells = InputEvent()->GetEMCALCells();
55a8e3e4 337 Double_t avgE = 0; // for background subtraction
338 const Int_t ncells = cells->GetNumberOfCells();
7a7d5a64 339 for (Int_t icell = 0, idigit = 0; icell < ncells; ++icell) {
77e93dc2 340 Double_t cellAmplitude=0, cellTime=0, cellEFrac = 0;
341 Short_t cellNumber=0, cellMCLabel=-1;
342 if (cells->GetCell(icell, cellNumber, cellAmplitude, cellTime, cellMCLabel, cellEFrac) != kTRUE)
7a7d5a64 343 break;
344 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->New(idigit));
345 digit->SetId(cellNumber);
346 digit->SetTime(cellTime);
347 digit->SetTimeR(cellTime);
348 digit->SetIndexInList(idigit);
349 digit->SetType(AliEMCALDigit::kHG);
350 if (fRecalibOnly||fSubBackground) {
351 Float_t energy = cellAmplitude;
55a8e3e4 352 Float_t time = cellTime;
7a7d5a64 353 fClusterizer->Calibrate(energy,time,cellNumber);
95f69c66 354 digit->SetAmplitude(energy);
7a7d5a64 355 avgE += energy;
356 } else {
357 digit->SetAmplitude(cellAmplitude);
95f69c66 358 }
7a7d5a64 359 idigit++;
b0e354fe 360 }
7a7d5a64 361
55a8e3e4 362 fDigitsArr->Sort();
363
7a7d5a64 364 if (fSubBackground) {
365 avgE /= AliEMCALGeometry::GetInstance(fGeomName)->GetNumberOfSuperModules()*48*24;
366 Int_t ndigis = fDigitsArr->GetEntries();
367 for (Int_t i = 0; i < ndigis; ++i) {
368 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(i));
369 Double_t energy = digit->GetAmplitude() - avgE;
370 if (energy<=0.001) {
371 digit->SetAmplitude(0);
372 } else {
373 digit->SetAmplitude(energy);
374 }
375 }
376 }
06ca320f 377 }
2f7259cf 378}
379
380//________________________________________________________________________________________
95f69c66 381void AliAnalysisTaskEMCALClusterizeFast::RecPoints2Clusters(TClonesArray *clus)
4cf9204b 382{
383 // Cluster energy, global position, cells and their amplitude fractions are restored.
384
95f69c66 385 Bool_t esdobjects = 0;
386 if (strcmp(clus->GetClass()->GetName(),"AliESDCaloCluster")==0)
387 esdobjects = 1;
388
b98487bc 389 AliVCaloCells *cells = InputEvent()->GetEMCALCells();
390 AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance(fGeomName);
e278a6e9 391 AliVEvent *event = InputEvent();
392 if (!event) {
393 AliError("Cannot get the event");
31f5b395 394 return;
395 }
0408f54f 396
397 // tracks array for track/cluster matching
398 TClonesArray *tarr = 0;
399 if (!fTrackName.IsNull()) {
e278a6e9 400 tarr = dynamic_cast<TClonesArray*>(event->FindListObject(fTrackName));
0408f54f 401 if (!tarr) {
402 AliError(Form("Cannot get tracks named %s", fTrackName.Data()));
403 }
404 }
e93ec1f7 405
eb33b7f8 406 const Int_t Ncls = fClusterArr->GetEntries();
0408f54f 407 AliDebug(1, Form("total no of clusters %d", Ncls));
4cf9204b 408 for(Int_t i=0, nout=clus->GetEntries(); i < Ncls; ++i) {
409 AliEMCALRecPoint *recpoint = static_cast<AliEMCALRecPoint*>(fClusterArr->At(i));
410 Int_t ncells_true = 0;
411 const Int_t ncells = recpoint->GetMultiplicity();
412 UShort_t absIds[ncells];
413 Double32_t ratios[ncells];
414 Int_t *dlist = recpoint->GetDigitsList();
415 Float_t *elist = recpoint->GetEnergiesList();
416 for (Int_t c = 0; c < ncells; ++c) {
417 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(dlist[c]));
418 absIds[ncells_true] = digit->GetId();
419 ratios[ncells_true] = elist[c]/digit->GetAmplitude();
84ccad86 420 ++ncells_true;
4cf9204b 421 }
422
423 if (ncells_true < 1) {
424 AliWarning("Skipping cluster with no cells");
425 continue;
426 }
427
428 // calculate new cluster position
429 TVector3 gpos;
4cf9204b 430 recpoint->GetGlobalPosition(gpos);
84ccad86 431 Float_t g[3];
4cf9204b 432 gpos.GetXYZ(g);
433
e93ec1f7 434 AliDebug(1, Form("energy %f", recpoint->GetEnergy()));
435
95f69c66 436 AliVCluster *c = static_cast<AliVCluster*>(clus->New(nout++));
4cf9204b 437 c->SetType(AliVCluster::kEMCALClusterv1);
438 c->SetE(recpoint->GetEnergy());
439 c->SetPosition(g);
440 c->SetNCells(ncells_true);
eb33b7f8 441 if (esdobjects) {
442 AliESDCaloCluster *cesd = static_cast<AliESDCaloCluster*>(c);
443 cesd->SetCellsAbsId(absIds);
444 cesd->SetCellsAmplitudeFraction(ratios);
445 cesd->SetID(recpoint->GetUniqueID());
446 } else {
447 AliAODCaloCluster *caod = static_cast<AliAODCaloCluster*>(c);
448 caod->SetCellsAbsId(absIds);
449 caod->SetCellsAmplitudeFraction(ratios);
450 }
4cf9204b 451 c->SetDispersion(recpoint->GetDispersion());
0408f54f 452 c->SetEmcCpvDistance(-1);
453 c->SetChi2(-1);
4cf9204b 454 c->SetTOF(recpoint->GetTime()) ; //time-of-flight
455 c->SetNExMax(recpoint->GetNExMax()); //number of local maxima
456 Float_t elipAxis[2];
457 recpoint->GetElipsAxis(elipAxis);
e93ec1f7 458 c->SetM02(elipAxis[0]*elipAxis[0]);
459 c->SetM20(elipAxis[1]*elipAxis[1]);
460 if (fPedestalData) {
b98487bc 461 c->SetDistanceToBadChannel(recpoint->GetDistanceToBadTower());
eb33b7f8 462 } else {
e93ec1f7 463 if (fRecoUtils && fRecoUtils->IsBadChannelsRemovalSwitchedOn()) {
464 fRecoUtils->RecalculateClusterDistanceToBadChannel(geom, cells, c);
465 }
b98487bc 466 }
84ccad86 467
0408f54f 468 if (tarr) {
469 Double_t dEtaMin = 1e9;
470 Double_t dPhiMin = 1e9;
471 Int_t imin = -1;
472 Double_t ceta = gpos.Eta();
473 Double_t cphi = gpos.Phi();
474 const Int_t ntrks = tarr->GetEntries();
475 for(Int_t t = 0; t<ntrks; ++t) {
476 AliVTrack *track = static_cast<AliVTrack*>(tarr->At(t));
477 if (!track)
478 continue;
479 const AliExternalTrackParam *outp = track->GetOuterParam();
480 if (!outp)
481 continue;
482 Double_t trkPos[3] = {0.,0.,0.};
483 if (!outp->GetXYZ(trkPos))
484 continue;
485 TVector3 vec(trkPos[0],trkPos[1],trkPos[2]);
486 Double_t veta = vec.Eta();
487 Double_t vphi = vec.Phi();
488 if(vphi<0)
489 vphi += 2*TMath::Pi();
490 if (TMath::Abs(veta)>0.75 || (vphi<70*TMath::DegToRad()) || (vphi>190*TMath::DegToRad()))
491 continue;
492 Double_t dR = vec.DeltaR(gpos);
493 if(dR > 25)
494 continue;
495 Float_t tmpEta=0, tmpPhi=0;
496 if (0) {
497 AliExternalTrackParam trkParTemp(*outp); // retrieve the starting point every time before the extrapolation
498 Bool_t ret = fRecoUtils->ExtrapolateTrackToCluster(&trkParTemp, c, fRecoUtils->GetMass(), fRecoUtils->GetStep(), tmpEta, tmpPhi);
499 if (!ret)
500 continue;
501 } else {
502 tmpEta = ceta - veta;
503 tmpPhi = cphi - vphi;
504 }
505 if (TMath::Abs(tmpEta)<TMath::Abs(dEtaMin) && TMath::Abs(tmpPhi)<TMath::Abs(dPhiMin)) {
506 dEtaMin = tmpEta;
507 dPhiMin = tmpPhi;
508 imin = t;
509 }
95f69c66 510 }
eb33b7f8 511 c->SetEmcCpvDistance(imin);
512 c->SetTrackDistance(dPhiMin, dEtaMin);
2f7259cf 513 }
514 }
515}
516
95f69c66 517//________________________________________________________________________
518void AliAnalysisTaskEMCALClusterizeFast::UpdateCells()
519{
520 // Update cells in case re-calibration was done.
521
522 if (!fCalibData&&!fSubBackground)
523 return;
524
525 AliVCaloCells *cells = InputEvent()->GetEMCALCells();
55a8e3e4 526 const Int_t ncells = cells->GetNumberOfCells();
527 const Int_t ndigis = fDigitsArr->GetEntries();
95f69c66 528 if (ncells!=ndigis) {
529 cells->DeleteContainer();
530 cells->CreateContainer(ndigis);
531 }
532 for (Int_t idigit = 0; idigit < ndigis; ++idigit) {
533 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(fDigitsArr->At(idigit));
534 Double_t cellAmplitude = digit->GetCalibAmp();
55a8e3e4 535 Short_t cellNumber = digit->GetId();
536 Double_t cellTime = digit->GetTime();
95f69c66 537 cells->SetCell(idigit, cellNumber, cellAmplitude, cellTime);
538 }
539}
540
541//________________________________________________________________________
542void AliAnalysisTaskEMCALClusterizeFast::UpdateClusters()
543{
544 // Update cells in case re-calibration was done.
7a7d5a64 545
95f69c66 546 if (!fAttachClusters)
547 return;
7a7d5a64 548
06ca320f 549 TClonesArray *clus = 0;
7a7d5a64 550
e50a444f 551 if (fOverwrite) {
7a7d5a64 552 clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("caloClusters"));
553 if (!clus)
554 clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("CaloClusters"));
06ca320f 555 if (!clus)
7a7d5a64 556 return;
557
55a8e3e4 558 const Int_t nents = clus->GetEntries();
7a7d5a64 559 for (Int_t i=0;i<nents;++i) {
560 AliVCluster *c = static_cast<AliVCluster*>(clus->At(i));
561 if (!c)
562 continue;
563 if (c->IsEMCAL())
564 delete clus->RemoveAt(i);
565 }
06ca320f 566 } else {
7a7d5a64 567 clus = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fNewClusterArrayName));
06ca320f 568 if (!clus) {
7a7d5a64 569 clus = new TClonesArray("AliESDCaloCluster");
570 clus->SetName(fNewClusterArrayName);
571 InputEvent()->AddObject(clus);
06ca320f 572 } else {
0a7eb283 573 clus->Delete();
7a7d5a64 574 }
575 }
576
95f69c66 577 RecPoints2Clusters(clus);
578}
579
2f7259cf 580//________________________________________________________________________________________
581void AliAnalysisTaskEMCALClusterizeFast::Init()
582{
583 //Select clusterization/unfolding algorithm and set all the needed parameters
01bc1ce8 584
2f7259cf 585 AliVEvent * event = InputEvent();
586 if (!event) {
587 AliWarning("Event not available!!!");
588 return;
589 }
590
591 if (event->GetRunNumber()==fRun)
592 return;
593 fRun = event->GetRunNumber();
594
595 if (fJustUnfold){
596 // init the unfolding afterburner
597 delete fUnfolder;
1759f743 598 fUnfolder = new AliEMCALAfterBurnerUF(fRecParam->GetW0(),fRecParam->GetLocMaxCut(),fRecParam->GetMinECut());
2f7259cf 599 return;
600 }
601
2f7259cf 602 AliEMCALGeometry *geometry = AliEMCALGeometry::GetInstance(fGeomName);
603 if (!geometry) {
604 AliFatal("Geometry not available!!!");
605 return;
606 }
607
608 if (!fGeomMatrixSet) {
609 if (fLoadGeomMatrices) {
e51bc3f5 610 for(Int_t mod=0; mod < geometry->GetNumberOfSuperModules(); ++mod) {
06ca320f 611 if (fGeomMatrix[mod]){
612 if (DebugLevel() > 2)
2f7259cf 613 fGeomMatrix[mod]->Print();
614 geometry->SetMisalMatrix(fGeomMatrix[mod],mod);
615 }
616 }
29b496d5 617 } else { // get matrix from file (work around bug in aliroot)
2f7259cf 618 for(Int_t mod=0; mod < geometry->GetEMCGeometry()->GetNumberOfSuperModules(); ++mod) {
29b496d5 619 const TGeoHMatrix *gm = 0;
441c4d3c 620 AliESDEvent *esdevent = dynamic_cast<AliESDEvent*>(event);
29b496d5 621 if (esdevent) {
622 gm = esdevent->GetEMCALMatrix(mod);
623 } else {
624 AliAODHeader *aodheader = dynamic_cast<AliAODHeader*>(event->GetHeader());
625 if (aodheader) {
626 gm = aodheader->GetEMCALMatrix(mod);
627 }
628 }
629 if (gm) {
06ca320f 630 if (DebugLevel() > 2)
29b496d5 631 gm->Print();
632 geometry->SetMisalMatrix(gm,mod);
2f7259cf 633 }
634 }
635 }
636 fGeomMatrixSet=kTRUE;
637 }
638
639 // setup digit array if needed
640 if (!fDigitsArr) {
641 fDigitsArr = new TClonesArray("AliEMCALDigit", 1000);
642 fDigitsArr->SetOwner(1);
643 }
644
645 // then setup clusterizer
646 delete fClusterizer;
647 if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv1)
648 fClusterizer = new AliEMCALClusterizerv1(geometry);
06ca320f 649 else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerNxN) {
2f7259cf 650 AliEMCALClusterizerNxN *clusterizer = new AliEMCALClusterizerNxN(geometry);
3c56da8f 651 clusterizer->SetNRowDiff(fRecParam->GetNRowDiff());
652 clusterizer->SetNColDiff(fRecParam->GetNColDiff());
653 fClusterizer = clusterizer;
06ca320f 654 } else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerv2)
3c56da8f 655 fClusterizer = new AliEMCALClusterizerv2(geometry);
06ca320f 656 else if (fRecParam->GetClusterizerFlag() == AliEMCALRecParam::kClusterizerFW){
7a7d5a64 657 AliEMCALClusterizerFixedWindow *clusterizer = new AliEMCALClusterizerFixedWindow(geometry);
0fda97fe 658 clusterizer->SetNphi(fNPhi);
659 clusterizer->SetNeta(fNEta);
660 clusterizer->SetShiftPhi(fShiftPhi);
661 clusterizer->SetShiftEta(fShiftEta);
3c56da8f 662 clusterizer->SetTRUshift(fTRUShift);
7a7d5a64 663 fClusterizer = clusterizer;
664 }
665 else{
2f7259cf 666 AliFatal(Form("Clusterizer < %d > not available", fRecParam->GetClusterizerFlag()));
667 }
668 fClusterizer->InitParameters(fRecParam);
01bc1ce8 669
670 if ((!fCalibData&&fLoadCalib) || (!fPedestalData&&fLoadPed)) {
671 AliCDBManager *cdb = AliCDBManager::Instance();
672 if (!cdb->IsDefaultStorageSet() && !fOCDBpath.IsNull())
673 cdb->SetDefaultStorage(fOCDBpath);
674 if (fRun!=cdb->GetRun())
675 cdb->SetRun(fRun);
676 }
349bf69e 677 if (!fCalibData&&fLoadCalib&&fRun>0) {
2f7259cf 678 AliCDBEntry *entry = static_cast<AliCDBEntry*>(AliCDBManager::Instance()->Get("EMCAL/Calib/Data"));
679 if (entry)
680 fCalibData = static_cast<AliEMCALCalibData*>(entry->GetObject());
681 if (!fCalibData)
682 AliFatal("Calibration parameters not found in CDB!");
683 }
349bf69e 684 if (!fPedestalData&&fLoadPed&&fRun>0) {
2f7259cf 685 AliCDBEntry *entry = static_cast<AliCDBEntry*>(AliCDBManager::Instance()->Get("EMCAL/Calib/Pedestals"));
686 if (entry)
687 fPedestalData = static_cast<AliCaloCalibPedestal*>(entry->GetObject());
688 }
10d33986 689 if (fCalibData) {
690 fClusterizer->SetInputCalibrated(kFALSE);
691 fClusterizer->SetCalibrationParameters(fCalibData);
10d33986 692 } else {
693 fClusterizer->SetInputCalibrated(kTRUE);
694 }
e93ec1f7 695 fClusterizer->SetCaloCalibPedestal(fPedestalData);
10d33986 696 fClusterizer->SetJustClusters(kTRUE);
2f7259cf 697 fClusterizer->SetDigitsArr(fDigitsArr);
698 fClusterizer->SetOutput(0);
699 fClusterArr = const_cast<TObjArray *>(fClusterizer->GetRecPoints());
700}