]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/AliJetModelBaseTask.cxx
from Redmer
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliJetModelBaseTask.cxx
CommitLineData
762e8424 1// $Id$
2//
3// Jet modelling task.
4//
297edd60 5// Author: S.Aiola, C.Loizides
762e8424 6
ffe32451 7#include "AliJetModelBaseTask.h"
8
762e8424 9#include <TClonesArray.h>
787a3c4f 10#include <TH1I.h>
762e8424 11#include <TLorentzVector.h>
12#include <TRandom3.h>
b16bb001 13#include <TList.h>
762e8424 14
b16bb001 15#include "AliVEvent.h"
e44e8726 16#include "AliAODCaloCluster.h"
b16bb001 17#include "AliESDCaloCluster.h"
18#include "AliVCluster.h"
762e8424 19#include "AliEMCALDigit.h"
ffe32451 20#include "AliEMCALRecPoint.h"
b16bb001 21#include "AliESDCaloCells.h"
22#include "AliAODCaloCells.h"
787a3c4f 23#include "AliAODMCParticle.h"
b16bb001 24#include "AliVCaloCells.h"
ffe32451 25#include "AliPicoTrack.h"
b16bb001 26#include "AliEMCALGeometry.h"
27#include "AliLog.h"
5be3857d 28#include "AliNamedArrayI.h"
762e8424 29
30ClassImp(AliJetModelBaseTask)
31
32//________________________________________________________________________
33AliJetModelBaseTask::AliJetModelBaseTask() :
34 AliAnalysisTaskSE("AliJetModelBaseTask"),
35 fGeomName(),
36 fTracksName(),
37 fOutTracksName(),
38 fCaloName(),
39 fOutCaloName(),
b16bb001 40 fCellsName(),
41 fOutCellsName(),
787a3c4f 42 fMCParticlesName(),
43 fOutMCParticlesName(),
4358e58a 44 fIsMC(kFALSE),
762e8424 45 fSuffix(),
46 fEtaMin(-1),
47 fEtaMax(1),
48 fPhiMin(0),
49 fPhiMax(TMath::Pi() * 2),
50 fPtMin(0),
51 fPtMax(0),
52 fCopyArray(kTRUE),
53 fNClusters(0),
b16bb001 54 fNCells(0),
762e8424 55 fNTracks(0),
7f76e479 56 fMarkMC(99999),
e44e8726 57 fPtSpectrum(0),
b16bb001 58 fQAhistos(kFALSE),
ffe32451 59 fIsInit(0),
762e8424 60 fGeom(0),
61 fClusters(0),
62 fOutClusters(0),
63 fTracks(0),
b16bb001 64 fOutTracks(0),
65 fCaloCells(0),
66 fOutCaloCells(0),
67 fAddedCells(0),
787a3c4f 68 fMCParticles(0),
69 fMCParticlesMap(0),
70 fOutMCParticles(0),
71 fOutMCParticlesMap(0),
72 fMCLabelShift(0),
b16bb001 73 fEsdMode(kFALSE),
74 fOutput(0)
762e8424 75{
76 // Default constructor.
77}
78
79//________________________________________________________________________
b16bb001 80AliJetModelBaseTask::AliJetModelBaseTask(const char *name, Bool_t drawqa) :
762e8424 81 AliAnalysisTaskSE(name),
82 fGeomName(""),
83 fTracksName("PicoTracks"),
84 fOutTracksName("PicoTracksEmbedded"),
85 fCaloName("CaloClustersCorr"),
86 fOutCaloName("CaloClustersCorrEmbedded"),
b16bb001 87 fCellsName(""),
88 fOutCellsName(""),
787a3c4f 89 fMCParticlesName(""),
90 fOutMCParticlesName(""),
4358e58a 91 fIsMC(kFALSE),
762e8424 92 fSuffix("Processed"),
93 fEtaMin(-1),
94 fEtaMax(1),
95 fPhiMin(0),
96 fPhiMax(TMath::Pi() * 2),
97 fPtMin(50),
98 fPtMax(60),
99 fCopyArray(kTRUE),
100 fNClusters(0),
b16bb001 101 fNCells(0),
762e8424 102 fNTracks(1),
7f76e479 103 fMarkMC(99999),
e44e8726 104 fPtSpectrum(0),
b16bb001 105 fQAhistos(drawqa),
ffe32451 106 fIsInit(0),
762e8424 107 fGeom(0),
108 fClusters(0),
109 fOutClusters(0),
110 fTracks(0),
b16bb001 111 fOutTracks(0),
112 fCaloCells(0),
113 fOutCaloCells(0),
114 fAddedCells(0),
787a3c4f 115 fMCParticles(0),
116 fMCParticlesMap(0),
117 fOutMCParticles(0),
118 fOutMCParticlesMap(0),
119 fMCLabelShift(0),
b16bb001 120 fEsdMode(kFALSE),
121 fOutput(0)
762e8424 122{
123 // Standard constructor.
b16bb001 124
125 if (fQAhistos) {
126 DefineOutput(1, TList::Class());
127 }
762e8424 128}
129
130//________________________________________________________________________
131AliJetModelBaseTask::~AliJetModelBaseTask()
132{
133 // Destructor
134}
135
b16bb001 136//________________________________________________________________________
137void AliJetModelBaseTask::UserCreateOutputObjects()
138{
139 // Create user output.
140 if (!fQAhistos)
141 return;
142
143 OpenFile(1);
144 fOutput = new TList();
145 fOutput->SetOwner();
146
147 PostData(1, fOutput);
148}
149
762e8424 150//________________________________________________________________________
ffe32451 151void AliJetModelBaseTask::UserExec(Option_t *)
152{
153 // Execute per event.
154
d63c8c07 155 if (!fIsInit)
156 fIsInit = ExecOnce();
157
158 if (!fIsInit)
159 return;
160
161 if (fCopyArray) {
162 if (fOutTracks)
163 fOutTracks->Delete();
164 if (fOutClusters)
165 fOutClusters->Delete();
787a3c4f 166 if (fOutMCParticles)
167 fOutMCParticles->Delete();
b4ffb4d2 168 }
63fac07f 169
5be3857d 170 // Clear map
171 if (fOutMCParticlesMap)
4358e58a 172 fOutMCParticlesMap->Clear();
d63c8c07 173
cd6431de 174 AliVCaloCells *tempCaloCells = 0;
175
b16bb001 176 if (fCaloCells) {
177 fAddedCells = 0;
cd6431de 178 if (!fCopyArray) {
179 tempCaloCells = fCaloCells;
180 fCaloCells = static_cast<AliVCaloCells*>(tempCaloCells->Clone(Form("%s_old",fCaloCells->GetName())));
181 }
b16bb001 182 }
183
ffe32451 184 Run();
cd6431de 185
fde82e42 186 if (fCaloCells && !fCopyArray) {
cd6431de 187 delete fCaloCells;
188 fCaloCells = tempCaloCells;
189 }
ffe32451 190}
191
b16bb001 192//________________________________________________________________________
193Bool_t AliJetModelBaseTask::ExecOnce()
194{
195 // Init task.
196
197 delete gRandom;
198 gRandom = new TRandom3(0);
199
200 fEsdMode = InputEvent()->InheritsFrom("AliESDEvent");
201
202 if (fPtMax < fPtMin) {
203 AliWarning (Form("PtMax (%f) < PtMin (%f), setting PtMax = PtMin = %f", fPtMax, fPtMin, fPtMin));
204 fPtMax = fPtMin;
205 }
206
207 if (fEtaMax < fEtaMin) {
208 AliWarning (Form("EtaMax (%f) < EtaMin (%f), setting EtaMax = EtaMin = %f", fEtaMax, fEtaMin, fEtaMin));
209 fEtaMax = fEtaMin;
210 }
211
212 if (fPhiMax < fPhiMin) {
213 AliWarning (Form("PhiMax (%f) < PhiMin (%f), setting PhiMax = PhiMin = %f", fPhiMax, fPhiMin, fPhiMin));
214 fPhiMax = fPhiMin;
215 }
216
787a3c4f 217 if (!fCellsName.IsNull()) {
b16bb001 218 fCaloCells = dynamic_cast<AliVCaloCells*>(InputEvent()->FindListObject(fCellsName));
219 if (!fCaloCells) {
787a3c4f 220 AliWarning(Form("%s: Couldn't retrieve calo cells with name %s!", GetName(), fCellsName.Data()));
b16bb001 221 }
787a3c4f 222 else if (!fCaloCells->InheritsFrom("AliVCaloCells")) {
b16bb001 223 AliError(Form("%s: Collection %s does not contain a AliVCaloCells object!", GetName(), fCellsName.Data()));
224 fCaloCells = 0;
225 return kFALSE;
226 }
227
228 if (!fOutCaloCells) {
229 fOutCellsName = fCellsName;
787a3c4f 230 if (fCopyArray)
b16bb001 231 fOutCellsName += fSuffix;
787a3c4f 232 if (fCopyArray || !fCaloCells) {
b16bb001 233 if (fEsdMode)
234 fOutCaloCells = new AliESDCaloCells(fOutCellsName,fOutCellsName);
235 else
236 fOutCaloCells = new AliAODCaloCells(fOutCellsName,fOutCellsName);
237
238 if (InputEvent()->FindListObject(fOutCellsName)) {
239 AliFatal(Form("%s: Collection %s is already present in the event!", GetName(), fOutCellsName.Data()));
240 return kFALSE;
241 }
242 else {
243 InputEvent()->AddObject(fOutCaloCells);
244 }
245 }
246 else {
247 fOutCaloCells = fCaloCells;
248 }
249 }
250 }
251
787a3c4f 252 if (fNTracks > 0 && !fTracksName.IsNull()) {
b16bb001 253 fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
254 if (!fTracks) {
787a3c4f 255 AliWarning(Form("%s: Couldn't retrieve tracks with name %s!", GetName(), fTracksName.Data()));
b16bb001 256 }
787a3c4f 257 else if (!fTracks->GetClass()->GetBaseClass("AliPicoTrack")) {
b16bb001 258 AliError(Form("%s: Collection %s does not contain AliPicoTrack objects!", GetName(), fTracksName.Data()));
259 fTracks = 0;
260 return kFALSE;
261 }
262
263 if (!fOutTracks) {
264 fOutTracksName = fTracksName;
787a3c4f 265 if (fCopyArray)
b16bb001 266 fOutTracksName += fSuffix;
787a3c4f 267 if (fCopyArray || !fTracks) {
268 fOutTracks = new TClonesArray("AliPicoTrack");
b16bb001 269 fOutTracks->SetName(fOutTracksName);
270 if (InputEvent()->FindListObject(fOutTracksName)) {
271 AliFatal(Form("%s: Collection %s is already present in the event!", GetName(), fOutTracksName.Data()));
272 return kFALSE;
273 }
274 else {
275 InputEvent()->AddObject(fOutTracks);
276 }
277 }
278 else {
279 fOutTracks = fTracks;
280 }
281 }
282 }
283
787a3c4f 284 if (fNClusters > 0 && !fCaloName.IsNull()) {
b16bb001 285 fClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
286
287 if (!fClusters) {
787a3c4f 288 AliWarning(Form("%s: Couldn't retrieve clusters with name %s!", GetName(), fCaloName.Data()));
b16bb001 289 }
787a3c4f 290 else if (!fClusters->GetClass()->GetBaseClass("AliVCluster")) {
b16bb001 291 AliError(Form("%s: Collection %s does not contain AliVCluster objects!", GetName(), fCaloName.Data()));
292 fClusters = 0;
293 return kFALSE;
294 }
295
296 if (!fOutClusters) {
297 fOutCaloName = fCaloName;
787a3c4f 298 if (fCopyArray)
b16bb001 299 fOutCaloName += fSuffix;
507f74bc 300 TString className;
301 if (fClusters)
302 className = fClusters->GetClass()->GetName();
303 else if (fEsdMode)
304 className = "AliESDCaloCluster";
305 else
306 className = "AliAODCaloCluster";
307
787a3c4f 308 if (fCopyArray || !fClusters) {
507f74bc 309 fOutClusters = new TClonesArray(className.Data());
b16bb001 310 fOutClusters->SetName(fOutCaloName);
311 if (InputEvent()->FindListObject(fOutCaloName)) {
312 AliFatal(Form("%s: Collection %s is already present in the event!", GetName(), fOutCaloName.Data()));
313 return kFALSE;
314 }
315 else {
316 InputEvent()->AddObject(fOutClusters);
317 }
318 }
319 else {
320 fOutClusters = fClusters;
321 }
322 }
323 }
324
787a3c4f 325 if (!fMCParticlesName.IsNull()) {
326 fMCParticles = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCParticlesName));
327 if (!fMCParticles) {
328 AliWarning(Form("%s: Couldn't retrieve MC particles with name %s!", GetName(), fMCParticlesName.Data()));
329 }
330 else {
331 if (!fMCParticles->GetClass()->GetBaseClass("AliAODMCParticle")) {
332 AliError(Form("%s: Collection %s does not contain AliAODMCParticle objects!", GetName(), fMCParticlesName.Data()));
333 fMCParticles = 0;
334 return kFALSE;
335 }
336
5be3857d 337 fMCParticlesMap = dynamic_cast<AliNamedArrayI*>(InputEvent()->FindListObject(fMCParticlesName + "_Map"));
787a3c4f 338
339 if (!fMCParticlesMap) {
340 AliWarning(Form("%s: Could not retrieve map for MC particles %s! Will assume MC labels consistent with indexes...", GetName(), fMCParticlesName.Data()));
5be3857d 341 fMCParticlesMap = new AliNamedArrayI(fMCParticlesName + "_Map", 9999);
787a3c4f 342 for (Int_t i = 0; i < 9999; i++) {
5be3857d 343 fMCParticlesMap->AddAt(i,i);
787a3c4f 344 }
345 }
346 }
347
348 if (!fOutMCParticles) {
349 fOutMCParticlesName = fMCParticlesName;
350 if (fCopyArray)
351 fOutMCParticlesName += fSuffix;
352 if (fCopyArray || !fMCParticles) {
353
354 fOutMCParticles = new TClonesArray("AliAODMCParticle");
355 fOutMCParticles->SetName(fOutMCParticlesName);
356 if (InputEvent()->FindListObject(fOutMCParticlesName)) {
357 AliFatal(Form("%s: Collection %s is already present in the event!", GetName(), fOutMCParticlesName.Data()));
358 return kFALSE;
359 }
360 else {
361 InputEvent()->AddObject(fOutMCParticles);
362 }
363
5be3857d 364 fOutMCParticlesMap = new AliNamedArrayI(fOutMCParticlesName + "_Map",9999);
787a3c4f 365 if (InputEvent()->FindListObject(fOutMCParticlesName + "_Map")) {
366 AliFatal(Form("%s: Map %s_Map is already present in the event!", GetName(), fOutMCParticlesName.Data()));
367 return kFALSE;
368 }
369 else {
370 InputEvent()->AddObject(fOutMCParticlesMap);
371 }
372 }
373 else {
374 fOutMCParticles = fMCParticles;
375 fOutMCParticlesMap = fMCParticlesMap;
376 }
377 }
378 }
379
2103dc6a 380 if (!fGeom) {
381 if (fGeomName.Length() > 0) {
382 fGeom = AliEMCALGeometry::GetInstance(fGeomName);
383 if (!fGeom) {
384 AliFatal(Form("Could not get geometry with name %s!", fGeomName.Data()));
385 return kFALSE;
386 }
387 } else {
4358e58a 388 fGeom = AliEMCALGeometry::GetInstance();
389 if (!fGeom) {
390 AliFatal("Could not get default geometry!");
391 return kFALSE;
392 }
b16bb001 393 }
cfc2ac24 394 }
b16bb001 395
396 return kTRUE;
397}
398
399//________________________________________________________________________
400Int_t AliJetModelBaseTask::SetNumberOfOutCells(Int_t n)
401{
402 if (fOutCaloCells->GetNumberOfCells() < n) {
403 fOutCaloCells->DeleteContainer();
404 fOutCaloCells->CreateContainer(n);
405 }
406 else {
407 fOutCaloCells->SetNumberOfCells(n);
408 }
409
410 fAddedCells = 0;
411
412 return n;
413}
414
b16bb001 415//________________________________________________________________________
416Int_t AliJetModelBaseTask::AddCell(Double_t e, Double_t eta, Double_t phi)
417{
418 // Add a cell to the event.
419
420 Int_t absId = 0;
421 if (eta < -100 || phi < 0) {
422 GetRandomCell(eta, phi, absId);
423 }
424 else {
425 fGeom->EtaPhiFromIndex(absId, eta, phi);
426 }
427
428 if (absId == -1) {
429 AliWarning(Form("Unable to embed cell in eta = %f, phi = %f!"
430 " Maybe the eta-phi range is not inside the EMCal acceptance (eta = [%f, %f], phi = [%f, %f])",
431 eta, phi, fEtaMin, fEtaMax, fPhiMin, fPhiMax));
432 return 0;
433 }
434
435 if (e < 0) {
436 Double_t pt = GetRandomPt();
437 TLorentzVector nPart;
438 nPart.SetPtEtaPhiM(pt, eta, phi, 0);
439 e = nPart.E();
440 }
441
442 return AddCell(e, absId);
443}
444
445//________________________________________________________________________
787a3c4f 446Int_t AliJetModelBaseTask::AddCell(Double_t e, Int_t absId, Double_t time, Int_t label)
b16bb001 447{
448 // Add a cell to the event.
449
fde82e42 450 if (label <= 0)
787a3c4f 451 label = fMarkMC;
452 else
453 label += fMCLabelShift;
454
fde82e42 455 Short_t pos = -1;
456 if (fCaloCells)
457 pos = fCaloCells->GetCellPosition(absId);
458
4358e58a 459 Double_t efrac = 1;
460 Bool_t increaseOnSuccess = kFALSE;
461
462 if (pos < 0) {
463 increaseOnSuccess = kTRUE;
464 pos = fAddedCells;
465 }
466 else {
467 Short_t cellNumber = -1;
468 Double_t old_e = 0;
469 Double_t old_time = 0;
470 Int_t old_label = 0;
471 Double_t old_efrac = 0;
472 fOutCaloCells->GetCell(pos, cellNumber, old_e, old_time, old_label, old_efrac);
473
fde82e42 474 efrac = e / (old_e + e);
475
476 if (old_label > 0 && e < old_efrac * old_e) {
4358e58a 477 label = old_label;
fde82e42 478 efrac = old_efrac;
4358e58a 479 time = old_time;
480 }
481
4358e58a 482 e += old_e;
483 }
484
485 Bool_t r = fOutCaloCells->SetCell(pos, absId, e, time, label, efrac);
b16bb001 486
487 if (r) {
4358e58a 488 if (increaseOnSuccess)
489 fAddedCells++;
b16bb001 490 return fAddedCells;
491 }
4358e58a 492 else
b16bb001 493 return 0;
b16bb001 494}
495
507f74bc 496//________________________________________________________________________
497AliVCluster* AliJetModelBaseTask::AddCluster(AliVCluster *oc)
498{
499 // Add a cluster to the event.
500
501 const Int_t nClusters = fOutClusters->GetEntriesFast();
502 AliVCluster *dc = static_cast<AliVCluster*>(fOutClusters->New(nClusters));
503 dc->SetType(AliVCluster::kEMCALClusterv1);
504 dc->SetE(oc->E());
505 Float_t pos[3] = {0};
506 oc->GetPosition(pos);
507 dc->SetPosition(pos);
508 dc->SetNCells(oc->GetNCells());
509 dc->SetCellsAbsId(oc->GetCellsAbsId());
510 dc->SetCellsAmplitudeFraction(oc->GetCellsAmplitudeFraction());
511 dc->SetID(oc->GetID());
512 dc->SetDispersion(oc->GetDispersion());
513 dc->SetEmcCpvDistance(-1);
514 dc->SetChi2(-1);
515 dc->SetTOF(oc->GetTOF()); //time-of-flight
516 dc->SetNExMax(oc->GetNExMax()); //number of local maxima
517 dc->SetM02(oc->GetM02());
518 dc->SetM20(oc->GetM20());
519 dc->SetDistanceToBadChannel(oc->GetDistanceToBadChannel());
520
521 //MC
522 UInt_t nlabels = oc->GetNLabels();
523 Int_t *labels = oc->GetLabels();
524
fde82e42 525 if (nlabels != 0 && labels && labels[0] >= 0) {
507f74bc 526 AliESDCaloCluster *esdClus = dynamic_cast<AliESDCaloCluster*>(dc);
527 if (esdClus) {
528 TArrayI parents(nlabels, labels);
529 esdClus->AddLabels(parents);
530 }
531 else {
532 AliAODCaloCluster *aodClus = dynamic_cast<AliAODCaloCluster*>(dc);
533 if (aodClus)
534 aodClus->SetLabel(labels, nlabels);
535 }
536 }
fde82e42 537 else if (fMarkMC != 0) {
538 AliESDCaloCluster *esdClus = dynamic_cast<AliESDCaloCluster*>(dc);
539 if (esdClus) {
540 TArrayI parents(1, &fMarkMC);
541 esdClus->AddLabels(parents);
542 }
543 else {
544 AliAODCaloCluster *aodClus = dynamic_cast<AliAODCaloCluster*>(dc);
545 if (aodClus)
546 aodClus->SetLabel(&fMarkMC,1);
547 }
548 }
507f74bc 549
550 return dc;
551}
b16bb001 552
ffe32451 553//________________________________________________________________________
787a3c4f 554AliVCluster* AliJetModelBaseTask::AddCluster(Double_t e, Double_t eta, Double_t phi, Int_t label)
ffe32451 555{
556 // Add a cluster to the event.
557
558 Int_t absId = 0;
559 if (eta < -100 || phi < 0) {
560 GetRandomCell(eta, phi, absId);
561 }
562 else {
563 fGeom->EtaPhiFromIndex(absId, eta, phi);
564 }
565
566 if (absId == -1) {
567 AliWarning(Form("Unable to embed cluster in eta = %f, phi = %f!"
568 " Maybe the eta-phi range is not inside the EMCal acceptance (eta = [%f, %f], phi = [%f, %f])",
569 eta, phi, fEtaMin, fEtaMax, fPhiMin, fPhiMax));
570 return 0;
571 }
572
573 if (e < 0) {
574 Double_t pt = GetRandomPt();
575 TLorentzVector nPart;
576 nPart.SetPtEtaPhiM(pt, eta, phi, 0);
577 e = nPart.E();
578 }
579
787a3c4f 580 return AddCluster(e, absId, label);
ffe32451 581}
582
583//________________________________________________________________________
787a3c4f 584AliVCluster* AliJetModelBaseTask::AddCluster(Double_t e, Int_t absId, Int_t label)
ffe32451 585{
586 // Add a cluster to the event.
587
588 const Int_t nClusters = fOutClusters->GetEntriesFast();
589
590 TClonesArray digits("AliEMCALDigit", 1);
591
592 AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(digits.New(0));
593 digit->SetId(absId);
594 digit->SetIndexInList(0);
595 digit->SetType(AliEMCALDigit::kHG);
596 digit->SetAmplitude(e);
597
598 AliEMCALRecPoint recPoint("");
599 recPoint.AddDigit(*digit, e, kFALSE);
600 recPoint.EvalGlobalPosition(0, &digits);
601
602 TVector3 gpos;
603 recPoint.GetGlobalPosition(gpos);
604 Float_t g[3];
605 gpos.GetXYZ(g);
606
607 AliVCluster *cluster = static_cast<AliVCluster*>(fOutClusters->New(nClusters));
608 cluster->SetType(AliVCluster::kEMCALClusterv1);
609 cluster->SetE(recPoint.GetEnergy());
610 cluster->SetPosition(g);
611 cluster->SetNCells(1);
612 UShort_t shortAbsId = absId;
613 cluster->SetCellsAbsId(&shortAbsId);
614 Double32_t fract = 1;
615 cluster->SetCellsAmplitudeFraction(&fract);
616 cluster->SetID(nClusters);
617 cluster->SetEmcCpvDistance(-1);
ffe32451 618
f660c2d6 619 //MC label
787a3c4f 620 if (label == 0)
621 label = fMarkMC;
622 else
623 label += fMCLabelShift;
624
f660c2d6 625 if (fEsdMode) {
626 AliESDCaloCluster *esdClus = static_cast<AliESDCaloCluster*>(cluster);
787a3c4f 627 TArrayI parents(1, &label);
f660c2d6 628 esdClus->AddLabels(parents);
629 }
630 else {
631 AliAODCaloCluster *aodClus = static_cast<AliAODCaloCluster*>(cluster);
787a3c4f 632 aodClus->SetLabel(&label, 1);
f660c2d6 633 }
634
ffe32451 635 return cluster;
636}
637
638//________________________________________________________________________
787a3c4f 639AliPicoTrack* AliJetModelBaseTask::AddTrack(Double_t pt, Double_t eta, Double_t phi, Byte_t type, Double_t etaemc, Double_t phiemc, Bool_t ise, Int_t label)
ffe32451 640{
641 // Add a track to the event.
642
643 const Int_t nTracks = fOutTracks->GetEntriesFast();
644
645 if (pt < 0)
646 pt = GetRandomPt();
647 if (eta < -100)
648 eta = GetRandomEta();
649 if (phi < 0)
650 phi = GetRandomPhi();
651
787a3c4f 652 if (label == 0)
653 label = fMarkMC;
654 else
655 label += fMCLabelShift;
656
ffe32451 657 AliPicoTrack *track = new ((*fOutTracks)[nTracks]) AliPicoTrack(pt,
658 eta,
659 phi,
b16bb001 660 1,
787a3c4f 661 label,
b16bb001 662 type,
663 etaemc,
664 phiemc,
665 ise);
666
ffe32451 667 return track;
668}
669
787a3c4f 670//________________________________________________________________________
671AliAODMCParticle* AliJetModelBaseTask::AddMCParticle(AliAODMCParticle *part, Int_t origIndex)
672{
673 const Int_t nPart = fOutMCParticles->GetEntriesFast();
674
675 AliAODMCParticle *aodpart = new ((*fOutMCParticles)[nPart]) AliAODMCParticle(*part);
676
5be3857d 677 fOutMCParticlesMap->AddAt(nPart, origIndex + fMCLabelShift);
5ce8ae64 678 AliDebug(2, Form("Setting bin %d to %d (fMCLabelShift=%d, origIndex=%d)",
679 origIndex + fMCLabelShift, fOutMCParticlesMap->At(origIndex + fMCLabelShift), fMCLabelShift, origIndex));
787a3c4f 680
681 return aodpart;
682}
683
4358e58a 684//________________________________________________________________________
685void AliJetModelBaseTask::CopyCells()
686{
687 if (!fCaloCells)
688 return;
689
fde82e42 690 fAddedCells = 0;
691 fCaloCells->Sort();
4358e58a 692 for (Short_t i = 0; i < fCaloCells->GetNumberOfCells(); i++) {
693 Int_t mclabel = 0;
694 Double_t efrac = 0.;
695 Double_t time = -1;
696 Short_t cellNum = -1;
697 Double_t amp = -1;
698
699 fCaloCells->GetCell(i, cellNum, amp, time, mclabel, efrac);
700
701 if (!fIsMC)
702 mclabel = 0;
703
704 fOutCaloCells->SetCell(i, cellNum, amp, time, mclabel, efrac);
fde82e42 705 fAddedCells++;
4358e58a 706 }
707
4358e58a 708 AliDebug(2, Form("%d cells from the current event", fAddedCells));
709}
710
ffe32451 711//________________________________________________________________________
712void AliJetModelBaseTask::CopyClusters()
713{
714 // Copy all the clusters in the new collection
787a3c4f 715 if (!fClusters)
716 return;
ffe32451 717
ffe32451 718 const Int_t nClusters = fClusters->GetEntriesFast();
56cc432c 719 Int_t nCopiedClusters = 0;
ffe32451 720
b16bb001 721 if (fEsdMode) {
ffe32451 722 for (Int_t i = 0; i < nClusters; ++i) {
723 AliESDCaloCluster *esdcluster = static_cast<AliESDCaloCluster*>(fClusters->At(i));
56cc432c 724 if (!esdcluster || !esdcluster->IsEMCAL())
ffe32451 725 continue;
4358e58a 726 AliESDCaloCluster *clus = new ((*fOutClusters)[nCopiedClusters]) AliESDCaloCluster(*esdcluster);
727 if (!fIsMC) {
728 TArrayI *labels = clus->GetLabelsArray();
729 if (labels)
730 labels->Reset();
731 }
56cc432c 732 nCopiedClusters++;
ffe32451 733 }
734 }
735 else {
736 for (Int_t i = 0; i < nClusters; ++i) {
737 AliAODCaloCluster *aodcluster = static_cast<AliAODCaloCluster*>(fClusters->At(i));
56cc432c 738 if (!aodcluster || !aodcluster->IsEMCAL())
ffe32451 739 continue;
4358e58a 740 AliAODCaloCluster *clus = new ((*fOutClusters)[nCopiedClusters]) AliAODCaloCluster(*aodcluster);
741 if (!fIsMC)
742 clus->SetLabel(0,0);
56cc432c 743 nCopiedClusters++;
ffe32451 744 }
745 }
746}
747
748//________________________________________________________________________
749void AliJetModelBaseTask::CopyTracks()
750{
787a3c4f 751 // Copy all the tracks in the new collection
752
753 if (!fTracks)
754 return;
755
ffe32451 756 const Int_t nTracks = fTracks->GetEntriesFast();
56cc432c 757 Int_t nCopiedTracks = 0;
ffe32451 758 for (Int_t i = 0; i < nTracks; ++i) {
4358e58a 759 AliPicoTrack *picotrack = static_cast<AliPicoTrack*>(fTracks->At(i));
760 if (!picotrack)
ffe32451 761 continue;
4358e58a 762 AliPicoTrack *track = new ((*fOutTracks)[nCopiedTracks]) AliPicoTrack(*picotrack);
763 if (!fIsMC && track->GetLabel() != 0)
764 track->SetLabel(0);
56cc432c 765 nCopiedTracks++;
ffe32451 766 }
767}
768
787a3c4f 769//________________________________________________________________________
770void AliJetModelBaseTask::CopyMCParticles()
771{
772 // Copy all the MC particles in the new collection
773
774 if (!fMCParticles)
775 return;
776
777 const Int_t nPart = fMCParticles->GetEntriesFast();
778 Int_t nCopiedPart = 0;
779 for (Int_t i = 0; i < nPart; ++i) {
780 AliAODMCParticle *part = static_cast<AliAODMCParticle*>(fMCParticles->At(i));
781 if (!part)
782 continue;
783 new ((*fOutMCParticles)[nCopiedPart]) AliAODMCParticle(*part);
784
785 nCopiedPart++;
786 }
787
788 if (!fMCParticlesMap)
789 return;
790
5be3857d 791 for (Int_t i = 0; i < fMCParticlesMap->GetSize(); i++) {
792 fOutMCParticlesMap->AddAt(fMCParticlesMap->At(i), i);
793 if (fMCParticlesMap->At(i) >= 0)
5ce8ae64 794 fMCLabelShift = i;
787a3c4f 795 }
796
5ce8ae64 797 AliDebug(2,Form("MC particles copied. fMCLabelShift=%d",fMCLabelShift));
787a3c4f 798}
799
762e8424 800//________________________________________________________________________
801void AliJetModelBaseTask::GetRandomCell(Double_t &eta, Double_t &phi, Int_t &absId)
802{
acd859a2 803 // Get random cell.
804
762e8424 805 Int_t repeats = 0;
05debff6 806 Double_t rndEta = eta;
807 Double_t rndPhi = phi;
762e8424 808 do {
05debff6 809 if (eta < -100)
2103dc6a 810 rndEta = GetRandomEta(kTRUE);
05debff6 811 if (phi < 0)
2103dc6a 812 rndPhi = GetRandomPhi(kTRUE);
05debff6 813 fGeom->GetAbsCellIdFromEtaPhi(rndEta, rndPhi, absId);
762e8424 814 repeats++;
815 } while (absId == -1 && repeats < 100);
816
817 if (!(absId > -1)) {
818 AliWarning(Form("Could not extract random cluster! Random eta-phi extracted more than 100 times!\n"
819 "eta [%f, %f], phi [%f, %f]\n", fEtaMin, fEtaMax, fPhiMin, fPhiMax));
820 }
05debff6 821 else {
822 eta = rndEta;
823 phi = rndPhi;
824 }
762e8424 825}
826
827//________________________________________________________________________
2103dc6a 828Double_t AliJetModelBaseTask::GetRandomEta(Bool_t emcal)
762e8424 829{
acd859a2 830 // Get random eta.
831
2103dc6a 832 Double_t etamax = fEtaMax;
833 Double_t etamin = fEtaMin;
834
835 if (emcal) {
836 const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
837 const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
838
839 if (etamax > EmcalMaxEta) etamax = EmcalMaxEta;
840 if (etamax < EmcalMinEta) etamax = EmcalMinEta;
841 if (etamin > EmcalMaxEta) etamin = EmcalMaxEta;
842 if (etamin < EmcalMinEta) etamin = EmcalMinEta;
843 }
844
845 return gRandom->Rndm() * (etamax - etamin) + etamin;
762e8424 846}
847
848//________________________________________________________________________
2103dc6a 849Double_t AliJetModelBaseTask::GetRandomPhi(Bool_t emcal)
762e8424 850{
acd859a2 851 // Get random phi.
2103dc6a 852
853 Double_t phimax = fPhiMax;
854 Double_t phimin = fPhiMin;
855
856 if (emcal) {
857 const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
858 const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
859
860 if (phimax > EmcalMaxPhi) phimax = EmcalMaxPhi;
861 if (phimax < EmcalMinPhi) phimax = EmcalMinPhi;
862 if (phimin > EmcalMaxPhi) phimin = EmcalMaxPhi;
863 if (phimin < EmcalMinPhi) phimin = EmcalMinPhi;
864 }
acd859a2 865
2103dc6a 866 return gRandom->Rndm() * (phimax - phimin) + phimin;
762e8424 867}
868
869//________________________________________________________________________
870Double_t AliJetModelBaseTask::GetRandomPt()
871{
acd859a2 872 // Get random pt.
873
e44e8726 874 if (fPtSpectrum)
875 return fPtSpectrum->GetRandom();
876 else
877 return gRandom->Rndm() * (fPtMax - fPtMin) + fPtMin;
762e8424 878}
879
762e8424 880//________________________________________________________________________
881void AliJetModelBaseTask::Run()
882{
acd859a2 883 // Run.
762e8424 884}