Coding violations corrected.
[u/mrichter/AliRoot.git] / JETAN / AliJetAODReader.cxx
CommitLineData
586f2bc3 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// Jet AOD Reader
18// AOD reader for jet analysis
e53baffe 19// This is the reader which must be used if the jet analysis task
20// is executed after the ESD filter task, in order to read its output
21//
586f2bc3 22// Author: Davide Perrino <davide.perrino@cern.ch>
23//-------------------------------------------------------------------------
24
25
26#include <Riostream.h>
27#include <TSystem.h>
28#include <TLorentzVector.h>
29#include <TVector3.h>
ee7de0dd 30#include <TChain.h>
be6e5811 31#include <TFile.h>
32#include <TTask.h>
33#include <TGeoManager.h>
ee7de0dd 34
586f2bc3 35#include "AliJetAODReader.h"
36#include "AliJetAODReaderHeader.h"
37#include "AliAODEvent.h"
38#include "AliAODTrack.h"
0e4398c8 39#include "AliAODMCParticle.h"
be6e5811 40#include "AliJetDummyGeo.h"
41#include "AliJetAODFillUnitArrayTracks.h"
42#include "AliJetAODFillUnitArrayEMCalDigits.h"
e36a3f22 43#include "AliJetHadronCorrectionv1.h"
be6e5811 44#include "AliJetUnitArray.h"
586f2bc3 45
46ClassImp(AliJetAODReader)
47
48AliJetAODReader::AliJetAODReader():
49 AliJetReader(),
586f2bc3 50 fAOD(0x0),
51 fRef(new TRefArray),
52 fDebug(0),
be6e5811 53 fOpt(0),
54 fGeom(0),
55 fHadCorr(0x0),
56 fTpcGrid(0x0),
57 fEmcalGrid(0x0),
58 fGrid0(0),
59 fGrid1(0),
60 fGrid2(0),
61 fGrid3(0),
62 fGrid4(0),
63 fPtCut(0),
be6e5811 64 fApplyElectronCorrection(kFALSE),
be6e5811 65 fApplyMIPCorrection(kTRUE),
66 fApplyFractionHadronicCorrection(kFALSE),
67 fFractionHadronicCorrection(0.3),
68 fNumUnits(0),
69 fMass(0),
70 fSign(0),
71 fNIn(0),
72 fDZ(0),
73 fNeta(0),
74 fNphi(0),
be6e5811 75 fRefArray(0x0),
76 fProcId(kFALSE)
586f2bc3 77{
78 // Constructor
79}
80
81//____________________________________________________________________________
82
83AliJetAODReader::~AliJetAODReader()
84{
85 // Destructor
586f2bc3 86 delete fAOD;
87 delete fRef;
be6e5811 88 delete fTpcGrid;
89 delete fEmcalGrid;
90 if(fDZ)
91 {
92 delete fGrid0;
93 delete fGrid1;
94 delete fGrid2;
95 delete fGrid3;
96 delete fGrid4;
97 }
98
586f2bc3 99}
100
101//____________________________________________________________________________
102
103void AliJetAODReader::OpenInputFiles()
104{
105 // Open the necessary input files
106 // chain for the AODs
107 fChain = new TChain("aodTree");
108
109 // get directory and pattern name from the header
110 const char* dirName=fReaderHeader->GetDirectory();
111 const char* pattern=fReaderHeader->GetPattern();
112
113// // Add files matching patters to the chain
114
115 void *dir = gSystem->OpenDirectory(dirName);
116 const char *name = 0x0;
117 int naod = ((AliJetAODReaderHeader*) fReaderHeader)->GetNaod();
118 int a = 0;
119 while ((name = gSystem->GetDirEntry(dir))){
120 if (a>=naod) continue;
121
122 if (strstr(name,pattern)){
123 char path[256];
124 sprintf(path,"%s/%s/aod.root",dirName,name);
125 fChain->AddFile(path);
126 a++;
127 }
128 }
129
130 gSystem->FreeDirectory(dir);
131
586f2bc3 132 fAOD = 0;
133 fChain->SetBranchAddress("AOD",&fAOD);
134
135 int nMax = fChain->GetEntries();
136
be6e5811 137 printf("\n AliJetAODReader: Total number of events in chain= %d \n",nMax);
586f2bc3 138
139 // set number of events in header
140 if (fReaderHeader->GetLastEvent() == -1)
141 fReaderHeader->SetLastEvent(nMax);
142 else {
143 Int_t nUsr = fReaderHeader->GetLastEvent();
144 fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
145 }
146}
147
148//____________________________________________________________________________
149
150void AliJetAODReader::ConnectTree(TTree* tree, TObject* /*data*/) {
151 // Connect the tree
152 // For AOD reader it's needed only to set the number of events
153 fChain = (TChain*) tree;
154
155 Int_t nMax = fChain->GetEntries();
be6e5811 156 printf("\n AliJetAODReader: Total number of events in chain= %5d \n", nMax);
586f2bc3 157 // set number of events in header
158 if (fReaderHeader->GetLastEvent() == -1)
159 fReaderHeader->SetLastEvent(nMax);
160 else {
161 Int_t nUsr = fReaderHeader->GetLastEvent();
162 fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
163 }
164}
165
0e4398c8 166
167Bool_t AliJetAODReader::AcceptAODMCParticle(AliAODMCParticle *mcP,Short_t flag){
168
169 //
170 // filter for charge and for charged and neutral, no detector
171 // response yet
172 // physical priamries already selected
173
174 Int_t pdg = TMath::Abs(mcP->GetPdgCode());
175
176 // exclude neutrinos anyway
177 if((pdg == 12 || pdg == 14 || pdg == 16)) return kFALSE;
178
179 if(flag==AliJetAODReaderHeader::kAllMC)return kTRUE;
180 if(flag==AliJetAODReaderHeader::kChargedMC){
181 if(mcP->Charge()==0)return kFALSE;
182 return kTRUE;
183 }
184
185 return kTRUE;
186
187}
188
189
190Bool_t AliJetAODReader::FillMomentumArrayMC(){
191
192 //
193 // This routine fetches the MC particles from the AOD
194 // Depending on the flag all particles except neurinos are use
195 // or only charged particles
196 //
197
198 TClonesArray *mcArray = (TClonesArray*)fAOD->FindListObject(AliAODMCParticle::StdBranchName());
199 if(!mcArray){
200 Printf("%s:%d No MC particle branch found",(char*)__FILE__,__LINE__);
201 return kFALSE;
202 }
203
204 Int_t nMC = mcArray->GetEntriesFast();
205
206 // get number of tracks in event (for the loop)
207 if(fDebug>0)printf("AOD MC tracks: %5d \t", nMC);
208
209 // temporary storage of signal and pt cut flag
210 Int_t* sflag = new Int_t[nMC];
211 Int_t* cflag = new Int_t[nMC];
212
213 // get cuts set by user
214 Float_t ptMin = fReaderHeader->GetPtCut();
215 Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
216 Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
217 Int_t mcTrack = 0;
218 Float_t pt, eta;
219 TVector3 p3;
220
221 Short_t readerFlag = ((AliJetAODReaderHeader*)fReaderHeader)->GetReadAODMC();
222
223
224 for (Int_t it = 0; it < nMC; it++) {
225 AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
226 if(!track->IsPhysicalPrimary())continue;
227
228 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
229 eta = p3.Eta();
230 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
231 if(!AcceptAODMCParticle(track,readerFlag))continue;
232 pt = p3.Pt();
233
234
235
236
237 if (mcTrack == 0) new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
238 new ((*fMomentumArray)[mcTrack]) TLorentzVector(p3,p3.Mag());
239 sflag[mcTrack] = 1;
240 cflag[mcTrack] = ( pt > ptMin ) ? 1: 0;
241 fRef->Add(track);
242 mcTrack++;
243 }
244 if(fDebug>0)printf("Used MC tracks: %5d \n", mcTrack);
245 // set the signal flags
246 fSignalFlag.Set(mcTrack,sflag);
247 fCutFlag.Set(mcTrack,cflag);
248
249 delete [] sflag;
250 delete [] cflag;
251
252 return kTRUE;
253}
254
586f2bc3 255//____________________________________________________________________________
256
9e4cc50d 257Bool_t AliJetAODReader::FillMomentumArray()
586f2bc3 258{
259 // Clear momentum array
260 ClearArray();
e53baffe 261 fRef->Clear();
586f2bc3 262 fDebug = fReaderHeader->GetDebug();
0e4398c8 263
264 if(((AliJetAODReaderHeader*)fReaderHeader)->GetReadAODMC()){
265 return FillMomentumArrayMC();
266 }
267
268
586f2bc3 269
270 if (!fAOD) {
271 return kFALSE;
272 }
273
274 // get number of tracks in event (for the loop)
275 Int_t nt = fAOD->GetNTracks();
0e4398c8 276 if(fDebug>0)printf("AOD tracks: %5d \t", nt);
586f2bc3 277
278 // temporary storage of signal and pt cut flag
279 Int_t* sflag = new Int_t[nt];
280 Int_t* cflag = new Int_t[nt];
281
282 // get cuts set by user
283 Float_t ptMin = fReaderHeader->GetPtCut();
284 Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
285 Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
d61f057d 286 UInt_t filterMask = ((AliJetAODReaderHeader*)fReaderHeader)->GetTestFilterMask();
e53baffe 287
586f2bc3 288 //loop over tracks
289 Int_t aodTrack = 0;
290 Float_t pt, eta;
291 TVector3 p3;
d61f057d 292
586f2bc3 293 for (Int_t it = 0; it < nt; it++) {
d61f057d 294 AliAODTrack *track = fAOD->GetTrack(it);
e53baffe 295 UInt_t status = track->GetStatus();
d61f057d 296
297 Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
298 p3.SetXYZ(mom[0],mom[1],mom[2]);
299 pt = p3.Pt();
300 eta = p3.Eta();
e53baffe 301 if (status == 0) continue;
d61f057d 302 if((filterMask>0)&&!(track->TestFilterBit(filterMask)))continue;
d61f057d 303 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
e53baffe 304
6473bac6 305 if (aodTrack == 0){
306 fRef->Delete(); // make sure to delete before placement new...
307 new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
308 }
e53baffe 309 new ((*fMomentumArray)[aodTrack]) TLorentzVector(p3,p3.Mag());
310 sflag[aodTrack] = (TMath::Abs(track->GetLabel()) < 10000) ? 1 : 0;
311 cflag[aodTrack] = ( pt > ptMin ) ? 1: 0;
312 aodTrack++;
d61f057d 313 fRef->Add(track);
586f2bc3 314 }
0e4398c8 315 if(fDebug>0)printf("Used AOD tracks: %5d \n", aodTrack);
586f2bc3 316 // set the signal flags
317 fSignalFlag.Set(aodTrack,sflag);
318 fCutFlag.Set(aodTrack,cflag);
d61f057d 319
680ed75f 320 delete [] sflag;
321 delete [] cflag;
d61f057d 322
586f2bc3 323 return kTRUE;
324}
be6e5811 325
326//__________________________________________________________
327void AliJetAODReader::SetApplyMIPCorrection(Bool_t val)
328{
329 //
330 // Set flag to apply MIP correction fApplyMIPCorrection
331 // - exclusive with fApplyFractionHadronicCorrection
332 //
333
334 fApplyMIPCorrection = val;
335 if(fApplyMIPCorrection == kTRUE)
336 {
337 SetApplyFractionHadronicCorrection(kFALSE);
338 printf("Enabling MIP Correction \n");
339 }
340 else
341 {
342 printf("Disabling MIP Correction \n");
343 }
344}
345
346//__________________________________________________________
347void AliJetAODReader::SetApplyFractionHadronicCorrection(Bool_t val)
348{
349 //
350 // Set flag to apply EMC hadronic correction fApplyFractionHadronicCorrection
351 // - exclusive with fApplyMIPCorrection
352 //
353
354 fApplyFractionHadronicCorrection = val;
355 if(fApplyFractionHadronicCorrection == kTRUE)
356 {
357 SetApplyMIPCorrection(kFALSE);
358 printf("Enabling Fraction Hadronic Correction \n");
359 }
360 else
361 {
362 printf("Disabling Fraction Hadronic Correction \n");
363 }
364}
365
366//__________________________________________________________
367void AliJetAODReader::SetFractionHadronicCorrection(Double_t val)
368{
369 //
370 // Set value to fFractionHadronicCorrection (default is 0.3)
371 // apply EMC hadronic correction fApplyFractionHadronicCorrection
372 // - exclusive with fApplyMIPCorrection
373 //
374
375 fFractionHadronicCorrection = val;
376 if(fFractionHadronicCorrection > 0.0 && fFractionHadronicCorrection <= 1.0)
377 {
378 SetApplyFractionHadronicCorrection(kTRUE);
379 printf("Fraction Hadronic Correction %1.3f \n",fFractionHadronicCorrection);
380 }
381 else
382 {
383 SetApplyFractionHadronicCorrection(kFALSE);
384 }
385}
386
387//____________________________________________________________________________
388void AliJetAODReader::CreateTasks(TChain* tree)
389{
4751efb5 390 //
391 // For reader task initialization
392 //
393
be6e5811 394 fDebug = fReaderHeader->GetDebug();
395 fDZ = fReaderHeader->GetDZ();
396 fTree = tree;
397
398 // Init EMCAL geometry and create UnitArray object
399 SetEMCALGeometry();
400 // cout << "In create task" << endl;
401 InitParameters();
402 InitUnitArray();
403
404 fFillUnitArray = new TTask("fFillUnitArray","Fill unit array jet finder");
405 fFillUAFromTracks = new AliJetAODFillUnitArrayTracks();
406 fFillUAFromTracks->SetReaderHeader(fReaderHeader);
407 fFillUAFromTracks->SetGeom(fGeom);
408 fFillUAFromTracks->SetTPCGrid(fTpcGrid);
409 fFillUAFromTracks->SetEMCalGrid(fEmcalGrid);
410
411 if(fDZ)
412 {
413 fFillUAFromTracks->SetGrid0(fGrid0);
414 fFillUAFromTracks->SetGrid1(fGrid1);
415 fFillUAFromTracks->SetGrid2(fGrid2);
416 fFillUAFromTracks->SetGrid3(fGrid3);
417 fFillUAFromTracks->SetGrid4(fGrid4);
418 }
419 fFillUAFromTracks->SetApplyMIPCorrection(fApplyMIPCorrection);
420 fFillUAFromTracks->SetHadCorrector(fHadCorr);
421 fFillUAFromEMCalDigits = new AliJetAODFillUnitArrayEMCalDigits();
422 fFillUAFromEMCalDigits->SetReaderHeader(fReaderHeader);
423 fFillUAFromEMCalDigits->SetGeom(fGeom);
424 fFillUAFromEMCalDigits->SetTPCGrid(fTpcGrid);
425 fFillUAFromEMCalDigits->SetEMCalGrid(fEmcalGrid);
426 fFillUAFromEMCalDigits->SetApplyFractionHadronicCorrection(fApplyFractionHadronicCorrection);
427 fFillUAFromEMCalDigits->SetFractionHadronicCorrection(fFractionHadronicCorrection);
428 fFillUAFromEMCalDigits->SetApplyElectronCorrection(fApplyElectronCorrection);
e36a3f22 429 // Add the task to global task
be6e5811 430 fFillUnitArray->Add(fFillUAFromTracks);
431 fFillUnitArray->Add(fFillUAFromEMCalDigits);
432 fFillUAFromTracks->SetActive(kFALSE);
433 fFillUAFromEMCalDigits->SetActive(kFALSE);
434
435 cout << "Tasks instantiated at that stage ! " << endl;
436 cout << "You can loop over events now ! " << endl;
437
438}
439
440//____________________________________________________________________________
441Bool_t AliJetAODReader::ExecTasks(Bool_t procid, TRefArray* refArray)
442{
443 //
444 // Main function
445 // Fill the reader part
446 //
447
448 fProcId = procid;
449 fRefArray = refArray;
450//(not used ?) Int_t nEntRef = fRefArray->GetEntries();
451//(not used ?) Int_t nEntUnit = fUnitArray->GetEntries();
452
453 // clear momentum array
454 ClearArray();
e36a3f22 455 fRef->Clear();
be6e5811 456
457 fDebug = fReaderHeader->GetDebug();
458 fOpt = fReaderHeader->GetDetector();
459
460 if(!fAOD) {
461 return kFALSE;
462 }
463
464 // TPC only or Digits+TPC or Clusters+TPC
465 if(fOpt%2==!0 && fOpt!=0){
466 fFillUAFromTracks->SetAOD(fAOD);
467 fFillUAFromTracks->SetActive(kTRUE);
468 fFillUAFromTracks->SetUnitArray(fUnitArray);
469 fFillUAFromTracks->SetRefArray(fRefArray);
e36a3f22 470 fFillUAFromTracks->SetReferences(fRef);
471 fFillUAFromTracks->SetSignalFlag(fSignalFlag);
472 fFillUAFromTracks->SetCutFlag(fCutFlag);
be6e5811 473 fFillUAFromTracks->SetProcId(fProcId);
474 // fFillUAFromTracks->ExecuteTask("tpc"); // => Temporarily changed
475 fFillUAFromTracks->Exec("tpc");
476 if(fOpt==1){
477 fNumCandidate = fFillUAFromTracks->GetMult();
478 fNumCandidateCut = fFillUAFromTracks->GetMultCut();
479 }
e36a3f22 480 fSignalFlag = fFillUAFromTracks->GetSignalFlag();
481 fCutFlag = fFillUAFromTracks->GetCutFlag();
be6e5811 482 }
483
484 // Digits only or Digits+TPC
485 if(fOpt>=2 && fOpt<=3){
486 fFillUAFromEMCalDigits->SetAOD(fAOD);
487 fFillUAFromEMCalDigits->SetActive(kTRUE);
488 fFillUAFromEMCalDigits->SetUnitArray(fUnitArray);
489 fFillUAFromEMCalDigits->SetRefArray(fRefArray);
490 fFillUAFromEMCalDigits->SetProcId(fFillUAFromTracks->GetProcId());
491 fFillUAFromEMCalDigits->SetInitMult(fFillUAFromTracks->GetMult());
492 fFillUAFromEMCalDigits->SetInitMultCut(fFillUAFromTracks->GetMultCut());
493 fFillUAFromEMCalDigits->Exec("digits"); // => Temporarily added
494 fNumCandidate = fFillUAFromEMCalDigits->GetMult();
495 fNumCandidateCut = fFillUAFromEMCalDigits->GetMultCut();
496 }
497
498 // fFillUnitArray->ExecuteTask(); // => Temporarily commented
499
500 return kTRUE;
501}
502
503//____________________________________________________________________________
504Bool_t AliJetAODReader::SetEMCALGeometry()
505{
506 //
507 // Set the EMCal Geometry
508 //
509
510 if (!fTree->GetFile())
511 return kFALSE;
512
513 TString geomFile(fTree->GetFile()->GetName());
514 geomFile.ReplaceAll("AliESDs", "geometry");
515
516 // temporary workaround for PROOF bug #18505
517 geomFile.ReplaceAll("#geometry.root#geometry.root", "#geometry.root");
518 if(fDebug>1) printf("Current geometry file %s \n", geomFile.Data());
519
520 // Define EMCAL geometry to be able to read ESDs
521 fGeom = AliJetDummyGeo::GetInstance();
522 if (fGeom == 0)
523 fGeom = AliJetDummyGeo::GetInstance("EMCAL_COMPLETE","EMCAL");
524
525 // To be setted to run some AliEMCALGeometry functions
526 TGeoManager::Import(geomFile);
527 fGeom->GetTransformationForSM();
528 printf("\n EMCal Geometry set ! \n");
529
530 return kTRUE;
531
532}
533
534//____________________________________________________________________________
535void AliJetAODReader::InitParameters()
536{
537 // Initialise parameters
538 fOpt = fReaderHeader->GetDetector();
be6e5811 539 fNumUnits = fGeom->GetNCells(); // Number of cells in EMCAL
540 if(fDebug>1) printf("\n EMCal parameters initiated ! \n");
541}
542
543//____________________________________________________________________________
544void AliJetAODReader::InitUnitArray()
545{
546 //Initialises unit arrays
547 Int_t nElements = fTpcGrid->GetNEntries();
f23f7908 548 Float_t eta = 0., phi = 0., deltaEta = 0., deltaPhi = 0.;
be6e5811 549 if(fArrayInitialised) fUnitArray->Delete();
550
551 if(fTpcGrid->GetGridType()==0)
552 { // Fill the following quantities :
f23f7908 553 // Good track ID, (Eta,Phi) position ID, eta, phi, energy, px, py, pz, deltaEta, deltaPhi,
be6e5811 554 // detector flag, in/out jet, pt cut, mass, cluster ID)
555 for(Int_t nBin = 1; nBin < nElements+1; nBin++)
556 {
557 // fTpcGrid->GetEtaPhiFromIndex2(nBin,eta,phi);
558 fTpcGrid->GetEtaPhiFromIndex2(nBin,phi,eta);
559 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
f23f7908 560 deltaEta = fTpcGrid->GetDeta();
561 deltaPhi = fTpcGrid->GetDphi();
562 new ((*fUnitArray)[nBin-1]) AliJetUnitArray(nBin-1,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 563 }
564 }
565
566 if(fTpcGrid->GetGridType()==1)
567 {
568 Int_t nGaps = 0;
569 Int_t n0 = 0, n1 = 0, n2 = 0, n3 = 0, n4 = 0;
570
571 if(fDZ)
572 {
573 // Define a grid of cell for the gaps between SM
574 Double_t phimin0 = 0., phimin1 = 0., phimin2 = 0., phimin3 = 0., phimin4 = 0.;
575 Double_t phimax0 = 0., phimax1 = 0., phimax2 = 0., phimax3 = 0., phimax4 = 0.;
576 fGeom->GetPhiBoundariesOfSMGap(0,phimin0,phimax0);
577 fGrid0 = new AliJetGrid(0,95,phimin0,phimax0,-0.7,0.7); // 0.015 x 0.015
578 fGrid0->SetGridType(0);
579 fGrid0->SetMatrixIndexes();
580 fGrid0->SetIndexIJ();
581 n0 = fGrid0->GetNEntries();
582 fGeom->GetPhiBoundariesOfSMGap(1,phimin1,phimax1);
583 fGrid1 = new AliJetGrid(0,95,phimin1,phimax1,-0.7,0.7); // 0.015 x 0.015
584 fGrid1->SetGridType(0);
585 fGrid1->SetMatrixIndexes();
586 fGrid1->SetIndexIJ();
587 n1 = fGrid1->GetNEntries();
588 fGeom->GetPhiBoundariesOfSMGap(2,phimin2,phimax2);
589 fGrid2 = new AliJetGrid(0,95,phimin2,phimax2,-0.7,0.7); // 0.015 x 0.015
590 fGrid2->SetGridType(0);
591 fGrid2->SetMatrixIndexes();
592 fGrid2->SetIndexIJ();
593 n2 = fGrid2->GetNEntries();
594 fGeom->GetPhiBoundariesOfSMGap(3,phimin3,phimax3);
595 fGrid3 = new AliJetGrid(0,95,phimin3,phimax3,-0.7,0.7); // 0.015 x 0.015
596 fGrid3->SetGridType(0);
597 fGrid3->SetMatrixIndexes();
598 fGrid3->SetIndexIJ();
599 n3 = fGrid3->GetNEntries();
600 fGeom->GetPhiBoundariesOfSMGap(4,phimin4,phimax4);
601 fGrid4 = new AliJetGrid(0,95,phimin4,phimax4,-0.7,0.7); // 0.015 x 0.015
602 fGrid4->SetGridType(0);
603 fGrid4->SetMatrixIndexes();
604 fGrid4->SetIndexIJ();
605 n4 = fGrid4->GetNEntries();
606
607 nGaps = n0+n1+n2+n3+n4;
608
609 }
610
611 for(Int_t nBin = 0; nBin < fNumUnits+nElements+nGaps; nBin++)
612 {
613 if(nBin<fNumUnits)
614 {
615 fGeom->EtaPhiFromIndex(nBin, eta, phi); // From EMCal geometry
616 // fEmcalGrid->GetEtaPhiFromIndex2(nBin,phi,eta); // My function from Grid
617 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
f23f7908 618 deltaEta = fEmcalGrid->GetDeta(); // Modify with the exact detector values
619 deltaPhi = fEmcalGrid->GetDphi(); // Modify with the exact detector values
620 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 621 }
622 else {
623 if(nBin>=fNumUnits && nBin<fNumUnits+nElements){
624 fTpcGrid->GetEtaPhiFromIndex2(nBin+1-fNumUnits,phi,eta);
625 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
f23f7908 626 deltaEta = fTpcGrid->GetDeta();
627 deltaPhi = fTpcGrid->GetDphi();
628 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 629 }
630 else {
631 if(fDZ) {
632 if(nBin>=fNumUnits+nElements && nBin<fNumUnits+nElements+nGaps){
633 if(nBin<fNumUnits+nElements+n0)
634 {
f23f7908 635 phi = eta = 0.;
be6e5811 636 fGrid0->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements),phi,eta);
f23f7908 637 deltaEta = fGrid0->GetDeta();
638 deltaPhi = fGrid0->GetDphi();
639 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 640 }
641 else if(nBin>=fNumUnits+nElements+n0 && nBin<fNumUnits+nElements+n0+n1)
642 {
f23f7908 643 phi = eta = 0.;
be6e5811 644 fGrid1->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0),phi,eta);
f23f7908 645 deltaEta = fGrid1->GetDeta();
646 deltaPhi = fGrid1->GetDphi();
647 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 648 }
649 else if(nBin>=fNumUnits+nElements+n0+n1 && nBin<fNumUnits+nElements+n0+n1+n2)
650 {
f23f7908 651 phi = eta = 0.;
be6e5811 652 fGrid2->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1),phi,eta);
f23f7908 653 deltaEta = fGrid2->GetDeta();
654 deltaPhi = fGrid2->GetDphi();
655 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 656 }
657 else if(nBin>=fNumUnits+nElements+n0+n1+n2 && nBin<fNumUnits+nElements+n0+n1+n2+n3)
658 {
f23f7908 659 phi = eta = 0.;
be6e5811 660 fGrid3->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1+n2),phi,eta);
f23f7908 661 deltaEta = fGrid3->GetDeta();
662 deltaPhi = fGrid3->GetDphi();
663 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 664 }
665 else if(nBin>=fNumUnits+nElements+n0+n1+n2+n3 && nBin<fNumUnits+nElements+nGaps)
666 {
f23f7908 667 phi = eta = 0.;
be6e5811 668 fGrid4->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1+n2+n3),phi,eta);
f23f7908 669 deltaEta = fGrid4->GetDeta();
670 deltaPhi = fGrid4->GetDphi();
671 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 672 }
673 }
674 } // end if(fDZ)
675 } // end else 2
676 } // end else 1
677 } // end loop on nBin
678 } // end grid type == 1
679 fArrayInitialised = 1;
680}
0e4398c8 681