]> git.uio.no Git - u/mrichter/AliRoot.git/blame - JETAN/AliJetAODFillUnitArrayTracks.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Adding new macro, loadFromOCDB.C, for loading magnetic field, mapping and
[u/mrichter/AliRoot.git] / JETAN / AliJetAODFillUnitArrayTracks.cxx
CommitLineData
be6e5811 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// ***July 2006
19// Fill Unit Array class
20// Class used by AliJetESDReader to fill a UnitArray from the information extracted
21// from the particle tracks
22// Author: magali.estienne@ires.in2p3.fr
23//======================================================================
24
25
26// --- ROOT system ---
27#include <TVector3.h>
28#include <TProcessID.h>
29
30#include "AliJetUnitArray.h"
e36a3f22 31#include "AliJetHadronCorrectionv1.h"
be6e5811 32#include "AliJetAODFillUnitArrayTracks.h"
33
34// --- ROOT system ---
35class TSystem;
36class TLorentzVector;
37class TGeoManager;
38
39// --- AliRoot header files ---
40class AliJetFinder;
41
42ClassImp(AliJetAODFillUnitArrayTracks)
43
44//_____________________________________________________________________________
45AliJetAODFillUnitArrayTracks::AliJetAODFillUnitArrayTracks()
46 : AliJetFillUnitArray(),
47 fNumUnits(0),
48 fHadCorr(0),
49 fApplyMIPCorrection(kTRUE),
50 fAOD(0x0),
51 fGrid0(0x0),
52 fGrid1(0x0),
53 fGrid2(0x0),
54 fGrid3(0x0),
55 fGrid4(0x0)
56{
57 // constructor
58}
59
60//_____________________________________________________________________________
61AliJetAODFillUnitArrayTracks::AliJetAODFillUnitArrayTracks(AliAODEvent* aod)
62 : AliJetFillUnitArray(),
63 fNumUnits(0),
64 fHadCorr(0),
65 fApplyMIPCorrection(kTRUE),
66 fAOD(aod),
67 fGrid0(0x0),
68 fGrid1(0x0),
69 fGrid2(0x0),
70 fGrid3(0x0),
71 fGrid4(0x0)
72{
73 // constructor
74}
75
76//_____________________________________________________________________________
77AliJetAODFillUnitArrayTracks::AliJetAODFillUnitArrayTracks(const AliJetAODFillUnitArrayTracks &det)
78 : AliJetFillUnitArray(det),
79 fNumUnits(det.fNumUnits),
80 fHadCorr(det.fHadCorr),
81 fApplyMIPCorrection(det.fApplyMIPCorrection),
82 fAOD(det.fAOD),
83 fGrid0(det.fGrid0),
84 fGrid1(det.fGrid1),
85 fGrid2(det.fGrid2),
86 fGrid3(det.fGrid3),
87 fGrid4(det.fGrid4)
88{
89 // Copy constructor
90}
91
92//_____________________________________________________________________________
93AliJetAODFillUnitArrayTracks& AliJetAODFillUnitArrayTracks::operator=(const AliJetAODFillUnitArrayTracks& other)
94{
95 // Assignment
96
97 fNumUnits = other.fNumUnits;
98 fHadCorr = other.fHadCorr;
99 fApplyMIPCorrection = other.fApplyMIPCorrection;
100 fAOD = other.fAOD;
101 fGrid0 = other.fGrid0;
102 fGrid1 = other.fGrid1;
103 fGrid2 = other.fGrid2;
104 fGrid3 = other.fGrid3;
105 fGrid4 = other.fGrid4;
106
107 return (*this);
108}
109
110//____________________________________________________________________________
111void AliJetAODFillUnitArrayTracks::InitParameters()
112{
e36a3f22 113 // fHadCorr = 0; // For hadron correction
be6e5811 114 fNumUnits = fGeom->GetNCells(); // Number of towers in EMCAL
115
116 fTPCGrid->GetAccParam(fNphi,fNeta,fPhiMin,
117 fPhiMax,fEtaMin,fEtaMax);
118 fTPCGrid->GetBinParam(fPhiBinInTPCAcc,fEtaBinInTPCAcc,
119 fPhiBinInEMCalAcc,fEtaBinInEMCalAcc,fNbinPhi);
120
121 fIndex = fTPCGrid->GetIndexObject();
122
123 if(fDebug>20){
124 for(Int_t i=0; i<fNphi+1; i++)
125 for(Int_t j=0; j<fNeta+1; j++) {cout << "fIndex[" << i << "," << j << "] : " <<
126 (*fIndex)(i,j) << endl; }
127 }
128 if(fDebug>1) printf("\n Parameters initiated ! \n");
129}
130
131//_____________________________________________________________________________
132AliJetAODFillUnitArrayTracks::~AliJetAODFillUnitArrayTracks()
133{
134 // destructor
135}
136
137//_____________________________________________________________________________
138void AliJetAODFillUnitArrayTracks::Exec(Option_t* const /*option*/)
139//void AliJetAODFillUnitArrayTracks::Exec(Option_t* option)
140{
141 //
142 // Main method.
143 // Fill the unit array with the charged particle information in AOD
144 //
145
146 fDebug = fReaderHeader->GetDebug();
147
148 // Set parameters
149 InitParameters();
e36a3f22 150 fRef->Clear();
be6e5811 151
152 // get number of tracks in event (for the loop)
153 Int_t goodTrack = 0;
154 Int_t nt = 0;
155//(not used ?) Int_t nt2 = 0;
156 Int_t nmax = 0;
157 Float_t pt,eta,phi;
158 // Int_t sflag = 0;
159 TVector3 p3;
160
161 nt = fAOD->GetNumberOfTracks();
e36a3f22 162 if(fDebug>1) cout << "Number of Tracks in AOD : " << nt << endl;
163
be6e5811 164 // temporary storage of signal and pt cut flag
165 Int_t* sflag = new Int_t[nt];
166 Int_t* cflag = new Int_t[nt];
167
168 // get cuts set by user
169 Float_t ptMin = fReaderHeader->GetPtCut();
170 Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
171 Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
172 fOpt = fReaderHeader->GetDetector();
173 fDZ = fReaderHeader->GetDZ();
174 UInt_t filterMask = ((AliJetAODReaderHeader*)fReaderHeader)->GetTestFilterMask();
175
176 Int_t nTracksEmcal = 0;
177 Int_t nTracksEmcalDZ = 0;
178 Int_t nTracksTpc = 0;
179 Int_t nTracksTpcOnly = 0;
180 Int_t nTracksEmcalCut = 0;
181 Int_t nTracksEmcalDZCut = 0;
182 Int_t nTracksTpcCut = 0;
183 Int_t nTracksTpcOnlyCut = 0;
184
be6e5811 185 fGrid = fTPCGrid->GetGridType();
186
be6e5811 187 //loop over tracks
188 nmax = nt;
189 for (Int_t it = 0; it < nmax; it++) {
190 AliAODTrack *track;
191 track = fAOD->GetTrack(it);
192 UInt_t status = track->GetStatus();
193
194 Double_t mom[3];
195 track->GetPxPyPz(mom);
196
197 p3.SetXYZ(mom[0],mom[1],mom[2]);
198 pt = p3.Pt();
199
be6e5811 200 eta = p3.Eta();
201 phi = ( (p3.Phi()) < 0) ? (p3.Phi()) + 2. * TMath::Pi() : (p3.Phi());
202 if (status == 0) continue;
203 if((filterMask>0)&&!(track->TestFilterBit(filterMask)))continue;
204
205 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
206
e36a3f22 207 if (goodTrack == 0) new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
208
be6e5811 209 sflag[goodTrack]=0;
210 if (TMath::Abs(track->GetLabel()) < 10000) sflag[goodTrack]=1;
211 cflag[goodTrack]=0;
212 if (pt > ptMin) cflag[goodTrack]=1; // pt cut
e36a3f22 213 fRef->Add(track);
be6e5811 214
215 if(fGrid==0)
216 {
217 // Only TPC filled from its grid in its total acceptance
218
219 Int_t idTPC = fTPCGrid->GetIndex(phi,eta);
220 Bool_t ok = kFALSE;
221 Bool_t ref = kFALSE;
222
223 AliJetUnitArray *uArray = (AliJetUnitArray*)fUnitArray->At(idTPC-1);
224 TRefArray *reference = uArray->GetUnitTrackRef();
a9c14d46 225 if (reference->GetEntries() == 0) {
226 new(reference) TRefArray(TProcessID::GetProcessWithUID(track));
227 }
228
be6e5811 229 reference->Add(track);
230
231 Float_t unitEnergy = 0.;
232 unitEnergy = uArray->GetUnitEnergy();
233 // nTracksTpcOnly is necessary to count the number of candidate cells
234 // but it doesn't give the real multiplicity -> it will be extracted
235 // in the jet finder through the reference to tracks
236 if(unitEnergy==0.){
237 nTracksTpcOnly++;
238 ok = kTRUE;
239 ref = kTRUE;
240 }
241
242 // Fill energy in TPC acceptance
243 uArray->SetUnitEnergy(unitEnergy + pt);
244
245 // Pt cut flag
246 if(uArray->GetUnitEnergy()<ptMin){
247 uArray->SetUnitCutFlag(kPtSmaller);
248 }
249 else {
250 uArray->SetUnitCutFlag(kPtHigher);
251 if(ok) nTracksTpcOnlyCut++;
252 }
253
254 // Signal flag
255 if(sflag[goodTrack] == 1) {
256 uArray->SetUnitSignalFlag(kGood);
257 uArray->SetUnitSignalFlagC(kFALSE,kGood);
258 } else uArray->SetUnitSignalFlagC(kFALSE,kBad);
259
260 if(uArray->GetUnitEnergy()>0 && ref){
261 if(!fProcId) {
262 fProcId = kTRUE;
263 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray));
264 }
265 fRefArray->Add(uArray);
266 }
267 }
268
269 if(fGrid==1)
270 {
271
272 Int_t nElements = fTPCGrid->GetNEntries();
273 // Fill track information in EMCAL acceptance
274 if((eta >= fEtaMin && eta <= fEtaMax) &&
275 (phi >= fPhiMin && phi <= fPhiMax))// &&
276 {
277 // Include dead-zones
278 if(fDZ)
279 {
280 Double_t phimin0 = 0., phimin1 = 0., phimin2 = 0., phimin3 = 0., phimin4 = 0.;
281 Double_t phimax0 = 0., phimax1 = 0., phimax2 = 0., phimax3 = 0., phimax4 = 0.;
282 fGeom->GetPhiBoundariesOfSMGap(0,phimin0,phimax0);
283 fGeom->GetPhiBoundariesOfSMGap(1,phimin1,phimax1);
284 fGeom->GetPhiBoundariesOfSMGap(2,phimin2,phimax2);
285 fGeom->GetPhiBoundariesOfSMGap(3,phimin3,phimax3);
286 fGeom->GetPhiBoundariesOfSMGap(4,phimin4,phimax4);
287 Int_t n0 = fGrid0->GetNEntries();
288 Int_t n1 = fGrid1->GetNEntries();
289 Int_t n2 = fGrid2->GetNEntries();
290 Int_t n3 = fGrid3->GetNEntries();
be6e5811 291
be6e5811 292 if(phi >= phimin0 && phi <= phimax0){
293 Int_t id0 = fGrid0->GetIndex(phi,eta)-1;
294 AliJetUnitArray *uArray0 = (AliJetUnitArray*)fUnitArray->At(id0+fNumUnits+nElements);
295 TRefArray *reference0 = uArray0->GetUnitTrackRef();
a9c14d46 296 if (reference0->GetEntries() == 0) {
297 new(reference0) TRefArray(TProcessID::GetProcessWithUID(track));
298 }
299
be6e5811 300 reference0->Add(track);
301 uArray0->SetUnitTrackID(it);
302
303 Float_t uEnergy0 = uArray0->GetUnitEnergy();
304 Bool_t ok0 = kFALSE;
305 Bool_t ref0 = kFALSE;
306 if(uEnergy0==0.){
307 nTracksEmcalDZ++;
308 ok0 = kTRUE;
309 ref0 = kTRUE;
310 }
311 uArray0->SetUnitEnergy(uEnergy0+pt);
312
313 if(uArray0->GetUnitEnergy()<ptMin)
314 uArray0->SetUnitCutFlag(kPtSmaller);
315 else {
316 uArray0->SetUnitCutFlag(kPtHigher);
317 if(ok0) nTracksEmcalDZCut++;
318 }
319 if(sflag[goodTrack] == 1) {
320 uArray0->SetUnitSignalFlag(kGood);
321 uArray0->SetUnitSignalFlagC(kFALSE,kGood);
322 } else uArray0->SetUnitSignalFlagC(kFALSE,kBad);
323
324 if(uArray0->GetUnitEnergy()>0 && ref0){
325 if(!fProcId){
326 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray0));
327 fProcId = kTRUE;
328 }
329 fRefArray->Add(uArray0);
330 }
331 }
332
333 if(phi >= phimin1 && phi <= phimax1){
334 Int_t id1 = fGrid1->GetIndex(phi,eta)-1+n0;
335 AliJetUnitArray *uArray1 = (AliJetUnitArray*)fUnitArray->At(id1+fNumUnits+nElements);
336 TRefArray *reference1 = uArray1->GetUnitTrackRef();
a9c14d46 337 if (reference1->GetEntries() == 0) {
338 new(reference1) TRefArray(TProcessID::GetProcessWithUID(track));
339 }
340
be6e5811 341 reference1->Add(track);
342
343 Float_t uEnergy1 = uArray1->GetUnitEnergy();
344 Bool_t ok1 = kFALSE;
345 Bool_t ref1 = kFALSE;
346 if(uEnergy1==0.){
347 nTracksEmcalDZ++;
348 ok1 = kTRUE;
349 ref1 = kTRUE;
350 }
351 uArray1->SetUnitEnergy(uEnergy1+pt);
352
353 if(uArray1->GetUnitEnergy()<ptMin)
354 uArray1->SetUnitCutFlag(kPtSmaller);
355 else {
356 uArray1->SetUnitCutFlag(kPtHigher);
357 if(ok1) nTracksEmcalDZCut++;
358 }
359 if(sflag[goodTrack] == 1) {
360 uArray1->SetUnitSignalFlag(kGood);
361 uArray1->SetUnitSignalFlagC(kFALSE,kGood);
362 } else uArray1->SetUnitSignalFlagC(kFALSE,kBad);
363
364 if(uArray1->GetUnitEnergy()>0 && ref1){
365 if(!fProcId){
366 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray1));
367 fProcId = kTRUE;
368 }
369 fRefArray->Add(uArray1);
370 }
371 }
372
373 if(phi >= phimin2 && phi <= phimax2){
374 Int_t id2 = fGrid2->GetIndex(phi,eta)-1+n0+n1;
375 AliJetUnitArray *uArray2 = (AliJetUnitArray*)fUnitArray->At(id2+fNumUnits+nElements);
376 TRefArray *reference2 = uArray2->GetUnitTrackRef();
a9c14d46 377 if (reference2->GetEntries() == 0) {
378 new(reference2) TRefArray(TProcessID::GetProcessWithUID(track));
379 }
380
be6e5811 381 reference2->Add(track);
382
383 Float_t uEnergy2 = uArray2->GetUnitEnergy();
384 Bool_t ok2 = kFALSE;
385 Bool_t ref2 = kFALSE;
386 if(uEnergy2==0.){
387 nTracksEmcalDZ++;
388 ok2 = kTRUE;
389 ref2 = kTRUE;
390 }
391 uArray2->SetUnitEnergy(uEnergy2+pt);
392
393 if(uArray2->GetUnitEnergy()<ptMin)
394 uArray2->SetUnitCutFlag(kPtSmaller);
395 else {
396 uArray2->SetUnitCutFlag(kPtHigher);
397 if(ok2) nTracksEmcalDZCut++;
398 }
399 if(sflag[goodTrack] == 1) {
400 uArray2->SetUnitSignalFlag(kGood);
401 uArray2->SetUnitSignalFlagC(kFALSE,kGood);
402 } else uArray2->SetUnitSignalFlagC(kFALSE,kBad);
403
404 if(uArray2->GetUnitEnergy()>0 && ref2){
405 if(!fProcId){
406 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray2));
407 fProcId = kTRUE;
408 }
409 fRefArray->Add(uArray2);
410 }
411 }
412
413 if(phi >= phimin3 && phi <= phimax3){
414 Int_t id3 = fGrid3->GetIndex(phi,eta)-1+n0+n1+n2;
415 AliJetUnitArray *uArray3 = (AliJetUnitArray*)fUnitArray->At(id3+fNumUnits+nElements);
416 TRefArray *reference3 = uArray3->GetUnitTrackRef();
a9c14d46 417 if (reference3->GetEntries() == 0) {
418 new(reference3) TRefArray(TProcessID::GetProcessWithUID(track));
419 }
420
be6e5811 421 reference3->Add(track);
422
423 Float_t uEnergy3 = uArray3->GetUnitEnergy();
424 Bool_t ok3 = kFALSE;
425 Bool_t ref3 = kFALSE;
426 if(uEnergy3==0.){
427 nTracksEmcalDZ++;
428 ok3 = kTRUE;
429 ref3 = kTRUE;
430 }
431 uArray3->SetUnitEnergy(uEnergy3+pt);
432
433 if(uArray3->GetUnitEnergy()<ptMin)
434 uArray3->SetUnitCutFlag(kPtSmaller);
435 else {
436 uArray3->SetUnitCutFlag(kPtHigher);
437 if(ok3) nTracksEmcalDZCut++;
438 }
439 if(sflag[goodTrack] == 1) {
440 uArray3->SetUnitSignalFlag(kGood);
441 uArray3->SetUnitSignalFlagC(kFALSE,kGood);
442 } else uArray3->SetUnitSignalFlagC(kFALSE,kBad);
443
444 if(uArray3->GetUnitEnergy()>0 && ref3){
445 if(!fProcId){
446 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray3));
447 fProcId = kTRUE;
448 }
449 fRefArray->Add(uArray3);
450 }
451 }
452
453 if(phi >= phimin4 && phi <= phimax4){
454 Int_t id4 = fGrid4->GetIndex(phi,eta)-1+n0+n1+n2+n3;
455 AliJetUnitArray *uArray4 = (AliJetUnitArray*)fUnitArray->At(id4+fNumUnits+nElements);
456 TRefArray *reference4 = uArray4->GetUnitTrackRef();
a9c14d46 457 if (reference4->GetEntries() == 0) {
458 new(reference4) TRefArray(TProcessID::GetProcessWithUID(track));
459 }
460
be6e5811 461 reference4->Add(track);
462
463 Float_t uEnergy4 = uArray4->GetUnitEnergy();
464 Bool_t ok4 = kFALSE;
465 Bool_t ref4 = kFALSE;
466 if(uEnergy4==0.){
467 nTracksEmcalDZ++;
468 ok4 = kTRUE;
469 ref4 = kTRUE;
470 }
471 uArray4->SetUnitEnergy(uEnergy4+pt);
472
473 if(uArray4->GetUnitEnergy()<ptMin)
474 uArray4->SetUnitCutFlag(kPtSmaller);
475 else {
476 uArray4->SetUnitCutFlag(kPtHigher);
477 if(ok4) nTracksEmcalDZCut++;
478 }
479 if(sflag[goodTrack] == 1) {
480 uArray4->SetUnitSignalFlag(kGood);
481 uArray4->SetUnitSignalFlagC(kFALSE,kGood);
482 } else uArray4->SetUnitSignalFlagC(kFALSE,kBad);
483
484 if(uArray4->GetUnitEnergy()>0 && ref4){
485 if(!fProcId){
486 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray4));
487 fProcId = kTRUE;
488 }
489 fRefArray->Add(uArray4);
490 }
491 }
492 } // end fDZ
493
494 Int_t towerID = 0;
495 fGeom->GetAbsCellIdFromEtaPhi(eta,phi,towerID);
496 if(towerID==-1) continue;
497
498 AliJetUnitArray *uArray = (AliJetUnitArray*)fUnitArray->At(towerID);
499 TRefArray *reference = uArray->GetUnitTrackRef();
a9c14d46 500 if (reference->GetEntries() == 0) {
501 new(reference) TRefArray(TProcessID::GetProcessWithUID(track));
502 }
503
be6e5811 504 reference->Add(track);
505
506 Float_t unitEnergy = uArray->GetUnitEnergy();
507 Bool_t ok = kFALSE;
508 Bool_t ref = kFALSE;
509 if(unitEnergy==0.){
510 nTracksEmcal++;
511 ok=kTRUE;
512 ref=kTRUE;
513 }
514
515 // Do Hadron Correction
e36a3f22 516 // For the moment I apply MIP correction if p >= 0.5 GeV/c
517 // and only for tracks inside EMCal acceptance
518 // End of if(fGrid==1) -> hadron correction for all tracks
519 if (fApplyMIPCorrection != 0 && p3.Mag() >= 0.5)
be6e5811 520 {
e36a3f22 521 ((AliJetHadronCorrectionv1*)fHadCorr)->SetGeometry("EMCAL_COMPLETE",1.);
522
523 Double_t etaOut = 0.;
524 Double_t phiOut = 0.;
525 ((AliJetHadronCorrectionv1*)fHadCorr)->TrackPositionEMCal(track,etaOut,phiOut);
526
527 // One has to make sure that the track hits the calorimeter
528 // We can then correct on average for these particles
529 if((etaOut >= fEtaMin && etaOut <= fEtaMax) &&
530 (phiOut >= fPhiMin && phiOut <= fPhiMax))// &&
531 {
532 Double_t hCEnergy = (Double_t)fHadCorr->GetEnergy(p3.Mag(), (Double_t)eta,0);
533 unitEnergy -= hCEnergy*TMath::Sin(2.0*TMath::ATan(TMath::Exp(-eta)));
534 }
be6e5811 535 } //end Hadron Correction loop
536
e36a3f22 537 if((unitEnergy + pt) > 0.) uArray->SetUnitEnergy(unitEnergy + pt);
538 else uArray->SetUnitEnergy(0.);
be6e5811 539
540 // Put a pt cut flag
541 if(uArray->GetUnitEnergy()<ptMin){
542 uArray->SetUnitCutFlag(kPtSmaller);
543 }
544 else {
545 uArray->SetUnitCutFlag(kPtHigher);
546 if(ok) nTracksEmcalCut++;
547 }
548
549 // Signal flag
550 if(sflag[goodTrack] == 1) {
551 uArray->SetUnitSignalFlag(kGood);
552 uArray->SetUnitSignalFlagC(kFALSE,kGood);
553 } else uArray->SetUnitSignalFlagC(kFALSE,kBad);
554
555 if(uArray->GetUnitEnergy()>0 && ref){
556 if(!fProcId){
557 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray));
558 fProcId = kTRUE;
559 }
560 fRefArray->Add(uArray);
561 }
562 } // end loop on EMCal acceptance cut + tracks quality
563 else{
564 // Outside EMCal acceptance
565 // Int_t idTPC = GetIndexFromEtaPhi(etaT[i],phiT[i]);
566 Int_t idTPC = fTPCGrid->GetIndex(phi,eta);
567
568 AliJetUnitArray *uArray = (AliJetUnitArray*)fUnitArray->At(fNumUnits-1+idTPC);
569 TRefArray *reference = uArray->GetUnitTrackRef();
a9c14d46 570 if (reference->GetEntries() == 0) {
571 new(reference) TRefArray(TProcessID::GetProcessWithUID(track));
572 }
be6e5811 573 reference->Add(track);
574
575 Float_t unitEnergy2 = uArray->GetUnitEnergy(); // check if fNumUnits or fNumUnits-1
576 Bool_t okout = kFALSE;
577 Bool_t refout = kFALSE;
578 if(unitEnergy2==0.){
579 nTracksTpc++;
580 okout=kTRUE;
581 refout=kTRUE;
582 }
583 // Fill energy outside emcal acceptance
584 uArray->SetUnitEnergy(unitEnergy2 + pt);
585
586 // Pt cut flag
587 if(uArray->GetUnitEnergy()<ptMin){
588 uArray->SetUnitCutFlag(kPtSmaller);
589 }
590 else {
591 uArray->SetUnitCutFlag(kPtHigher);
592 if(okout) nTracksTpcCut++;
593 }
594 // Signal flag
595 if(sflag[goodTrack] == 1) {
596 uArray->SetUnitSignalFlag(kGood);
597 uArray->SetUnitSignalFlagC(kFALSE,kGood);
598 } else uArray->SetUnitSignalFlagC(kFALSE,kBad);
599
600 if(uArray->GetUnitEnergy()>0 && refout){
601 if(!fProcId){
602 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray));
603 fProcId = kTRUE;
604 }
605 fRefArray->Add(uArray);
606 }
607 }
e36a3f22 608
609 /*
610 // Do hadron correction
611 // In principle, the best thing to do as particles which
612 // eta,phi (at vertex) correspond to EMCal dead zone can deposit energy
613 // in the calorimeter as well as particles with eta,phi (at vertex)
614 // outside the EMCal acceptance
615 cout << "Hadronic correction for all tracks goes here" << endl;
616
617 // For the moment I apply MIP correction if p >= 0.5 GeV/c
618 if (fApplyMIPCorrection != 0 && p3.Mag() >= 0.5)
619 {
620 ((AliJetHadronCorrectionv1*)fHadCorr)->SetGeometry("EMCAL_COMPLETE",1.);
621
622 Double_t etaOut = 0.;
623 Double_t phiOut = 0.;
624 Int_t towerID2 = 0;
625 Float_t unitEnergy = 0.;
626 ((AliJetHadronCorrectionv1*)fHadCorr)->TrackPositionEMCal(track,etaOut,phiOut);
627
628 // One has to make sure that the track hits the calorimeter
629 // We can then correct on average for these particles
630 if((etaOut >= fEtaMin && etaOut <= fEtaMax) &&
631 (phiOut >= fPhiMin && phiOut <= fPhiMax))// &&
632 {
633 Int_t towerID = 0;
634 Int_t towerID2 = 0;
635 fGeom->GetAbsCellIdFromEtaPhi(etaOut,phiOut,towerID2);
636 if(towerID2==-1) continue;
637
638 AliJetUnitArray *uArray = (AliJetUnitArray*)fUnitArray->At(towerID2);
639 unitEnergy = uArray->GetUnitEnergy();
640 Bool_t ok = kFALSE;
641 Bool_t ref = kFALSE;
642 if(unitEnergy==0.){
643 nTracksEmcal++;
644 ok=kTRUE;
645 ref=kTRUE;
646 }
647
648 Double_t hCEnergy = (Double_t)fHadCorr->GetEnergy(p3.Mag(), (Double_t)eta,0);
649 unitEnergy -= hCEnergy*TMath::Sin(2.0*TMath::ATan(TMath::Exp(-eta)));
650
651 // The energy in this unitarray can be negative
652 // If the Digits are then filled, make sure the energy sum
653 // is positive. If not, put the value to zero in AliJetAODUnitArrayEMCalDigit.
654 // If they are not filled, put this value to zero
655 uArray->SetUnitEnergy(unitEnergy);
656
657 if(uArray->GetUnitEnergy()!=0. && ref){
658 if(!fProcId){
659 new(fRefArray) TRefArray(TProcessID::GetProcessWithUID(uArray));
660 fProcId = kTRUE;
661 }
662 fRefArray->Add(uArray);
663 }
664 }
665 else cout << "Track out of EMCal acceptance" << endl;
666
667 } //end Hadron Correction loop
668 */
669
be6e5811 670 } // end fGrid==1
671
672 goodTrack++;
673
674 } // end loop on entries (tpc tracks)
675
676
677 if(fDebug>0)
678 {
679 cout << "End of Tracks %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%" << endl;
680 cout << "goodTracks: " << goodTrack << endl;
681 }
682
e36a3f22 683 fSignalFlag.Set(goodTrack,sflag);
684 fCutFlag.Set(goodTrack,cflag);
685
be6e5811 686 delete[] sflag;
687 delete[] cflag;
688
689 // } // end loop on entries (tpc tracks)
690
691 if(fGrid==0) {
692 fNTracks = nTracksTpcOnly;
693 fNTracksCut = nTracksTpcOnlyCut;
694 if(fDebug>10){
695 cout << "fNTracks : " << fNTracks << endl;
696 cout << "fNTracksCut : " << fNTracksCut << endl;
697 }
698 }
699 if(fGrid==1) {
700 fNTracks = nTracksEmcal+nTracksEmcalDZ+nTracksTpc;
701 fNTracksCut = nTracksEmcalCut+nTracksEmcalDZCut+nTracksTpcCut;
702 if(fDebug>10){
703 cout << "fNTracks : " << fNTracks << endl;
704 cout << "fNTracksCut : " << fNTracksCut << endl;
705 }
706 }
707
708}
709
710
711
712
713
714
715
716
717
718
719
720
721