]> git.uio.no Git - u/mrichter/AliRoot.git/blame - JETAN/AliJetAODReader.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
New class to extract pt distribution of given particle from pt distribution of a...
[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
921b3f4d 305 if (aodTrack == 0) new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
306
e53baffe 307 new ((*fMomentumArray)[aodTrack]) TLorentzVector(p3,p3.Mag());
308 sflag[aodTrack] = (TMath::Abs(track->GetLabel()) < 10000) ? 1 : 0;
309 cflag[aodTrack] = ( pt > ptMin ) ? 1: 0;
310 aodTrack++;
d61f057d 311 fRef->Add(track);
586f2bc3 312 }
0e4398c8 313 if(fDebug>0)printf("Used AOD tracks: %5d \n", aodTrack);
586f2bc3 314 // set the signal flags
315 fSignalFlag.Set(aodTrack,sflag);
316 fCutFlag.Set(aodTrack,cflag);
d61f057d 317
680ed75f 318 delete [] sflag;
319 delete [] cflag;
d61f057d 320
586f2bc3 321 return kTRUE;
322}
be6e5811 323
324//__________________________________________________________
325void AliJetAODReader::SetApplyMIPCorrection(Bool_t val)
326{
327 //
328 // Set flag to apply MIP correction fApplyMIPCorrection
329 // - exclusive with fApplyFractionHadronicCorrection
330 //
331
332 fApplyMIPCorrection = val;
333 if(fApplyMIPCorrection == kTRUE)
334 {
335 SetApplyFractionHadronicCorrection(kFALSE);
336 printf("Enabling MIP Correction \n");
337 }
338 else
339 {
340 printf("Disabling MIP Correction \n");
341 }
342}
343
344//__________________________________________________________
345void AliJetAODReader::SetApplyFractionHadronicCorrection(Bool_t val)
346{
347 //
348 // Set flag to apply EMC hadronic correction fApplyFractionHadronicCorrection
349 // - exclusive with fApplyMIPCorrection
350 //
351
352 fApplyFractionHadronicCorrection = val;
353 if(fApplyFractionHadronicCorrection == kTRUE)
354 {
355 SetApplyMIPCorrection(kFALSE);
356 printf("Enabling Fraction Hadronic Correction \n");
357 }
358 else
359 {
360 printf("Disabling Fraction Hadronic Correction \n");
361 }
362}
363
364//__________________________________________________________
365void AliJetAODReader::SetFractionHadronicCorrection(Double_t val)
366{
367 //
368 // Set value to fFractionHadronicCorrection (default is 0.3)
369 // apply EMC hadronic correction fApplyFractionHadronicCorrection
370 // - exclusive with fApplyMIPCorrection
371 //
372
373 fFractionHadronicCorrection = val;
374 if(fFractionHadronicCorrection > 0.0 && fFractionHadronicCorrection <= 1.0)
375 {
376 SetApplyFractionHadronicCorrection(kTRUE);
377 printf("Fraction Hadronic Correction %1.3f \n",fFractionHadronicCorrection);
378 }
379 else
380 {
381 SetApplyFractionHadronicCorrection(kFALSE);
382 }
383}
384
385//____________________________________________________________________________
386void AliJetAODReader::CreateTasks(TChain* tree)
387{
4751efb5 388 //
389 // For reader task initialization
390 //
391
be6e5811 392 fDebug = fReaderHeader->GetDebug();
393 fDZ = fReaderHeader->GetDZ();
394 fTree = tree;
395
396 // Init EMCAL geometry and create UnitArray object
397 SetEMCALGeometry();
398 // cout << "In create task" << endl;
399 InitParameters();
400 InitUnitArray();
401
402 fFillUnitArray = new TTask("fFillUnitArray","Fill unit array jet finder");
403 fFillUAFromTracks = new AliJetAODFillUnitArrayTracks();
404 fFillUAFromTracks->SetReaderHeader(fReaderHeader);
405 fFillUAFromTracks->SetGeom(fGeom);
406 fFillUAFromTracks->SetTPCGrid(fTpcGrid);
407 fFillUAFromTracks->SetEMCalGrid(fEmcalGrid);
408
409 if(fDZ)
410 {
411 fFillUAFromTracks->SetGrid0(fGrid0);
412 fFillUAFromTracks->SetGrid1(fGrid1);
413 fFillUAFromTracks->SetGrid2(fGrid2);
414 fFillUAFromTracks->SetGrid3(fGrid3);
415 fFillUAFromTracks->SetGrid4(fGrid4);
416 }
417 fFillUAFromTracks->SetApplyMIPCorrection(fApplyMIPCorrection);
418 fFillUAFromTracks->SetHadCorrector(fHadCorr);
419 fFillUAFromEMCalDigits = new AliJetAODFillUnitArrayEMCalDigits();
420 fFillUAFromEMCalDigits->SetReaderHeader(fReaderHeader);
421 fFillUAFromEMCalDigits->SetGeom(fGeom);
422 fFillUAFromEMCalDigits->SetTPCGrid(fTpcGrid);
423 fFillUAFromEMCalDigits->SetEMCalGrid(fEmcalGrid);
424 fFillUAFromEMCalDigits->SetApplyFractionHadronicCorrection(fApplyFractionHadronicCorrection);
425 fFillUAFromEMCalDigits->SetFractionHadronicCorrection(fFractionHadronicCorrection);
426 fFillUAFromEMCalDigits->SetApplyElectronCorrection(fApplyElectronCorrection);
e36a3f22 427 // Add the task to global task
be6e5811 428 fFillUnitArray->Add(fFillUAFromTracks);
429 fFillUnitArray->Add(fFillUAFromEMCalDigits);
430 fFillUAFromTracks->SetActive(kFALSE);
431 fFillUAFromEMCalDigits->SetActive(kFALSE);
432
433 cout << "Tasks instantiated at that stage ! " << endl;
434 cout << "You can loop over events now ! " << endl;
435
436}
437
438//____________________________________________________________________________
439Bool_t AliJetAODReader::ExecTasks(Bool_t procid, TRefArray* refArray)
440{
441 //
442 // Main function
443 // Fill the reader part
444 //
445
446 fProcId = procid;
447 fRefArray = refArray;
448//(not used ?) Int_t nEntRef = fRefArray->GetEntries();
449//(not used ?) Int_t nEntUnit = fUnitArray->GetEntries();
450
451 // clear momentum array
452 ClearArray();
e36a3f22 453 fRef->Clear();
be6e5811 454
455 fDebug = fReaderHeader->GetDebug();
456 fOpt = fReaderHeader->GetDetector();
457
458 if(!fAOD) {
459 return kFALSE;
460 }
461
462 // TPC only or Digits+TPC or Clusters+TPC
463 if(fOpt%2==!0 && fOpt!=0){
464 fFillUAFromTracks->SetAOD(fAOD);
465 fFillUAFromTracks->SetActive(kTRUE);
466 fFillUAFromTracks->SetUnitArray(fUnitArray);
467 fFillUAFromTracks->SetRefArray(fRefArray);
e36a3f22 468 fFillUAFromTracks->SetReferences(fRef);
469 fFillUAFromTracks->SetSignalFlag(fSignalFlag);
470 fFillUAFromTracks->SetCutFlag(fCutFlag);
be6e5811 471 fFillUAFromTracks->SetProcId(fProcId);
472 // fFillUAFromTracks->ExecuteTask("tpc"); // => Temporarily changed
473 fFillUAFromTracks->Exec("tpc");
474 if(fOpt==1){
475 fNumCandidate = fFillUAFromTracks->GetMult();
476 fNumCandidateCut = fFillUAFromTracks->GetMultCut();
477 }
e36a3f22 478 fSignalFlag = fFillUAFromTracks->GetSignalFlag();
479 fCutFlag = fFillUAFromTracks->GetCutFlag();
be6e5811 480 }
481
482 // Digits only or Digits+TPC
483 if(fOpt>=2 && fOpt<=3){
484 fFillUAFromEMCalDigits->SetAOD(fAOD);
485 fFillUAFromEMCalDigits->SetActive(kTRUE);
486 fFillUAFromEMCalDigits->SetUnitArray(fUnitArray);
487 fFillUAFromEMCalDigits->SetRefArray(fRefArray);
488 fFillUAFromEMCalDigits->SetProcId(fFillUAFromTracks->GetProcId());
489 fFillUAFromEMCalDigits->SetInitMult(fFillUAFromTracks->GetMult());
490 fFillUAFromEMCalDigits->SetInitMultCut(fFillUAFromTracks->GetMultCut());
491 fFillUAFromEMCalDigits->Exec("digits"); // => Temporarily added
492 fNumCandidate = fFillUAFromEMCalDigits->GetMult();
493 fNumCandidateCut = fFillUAFromEMCalDigits->GetMultCut();
494 }
495
496 // fFillUnitArray->ExecuteTask(); // => Temporarily commented
497
498 return kTRUE;
499}
500
501//____________________________________________________________________________
502Bool_t AliJetAODReader::SetEMCALGeometry()
503{
504 //
505 // Set the EMCal Geometry
506 //
507
508 if (!fTree->GetFile())
509 return kFALSE;
510
511 TString geomFile(fTree->GetFile()->GetName());
512 geomFile.ReplaceAll("AliESDs", "geometry");
513
514 // temporary workaround for PROOF bug #18505
515 geomFile.ReplaceAll("#geometry.root#geometry.root", "#geometry.root");
516 if(fDebug>1) printf("Current geometry file %s \n", geomFile.Data());
517
518 // Define EMCAL geometry to be able to read ESDs
519 fGeom = AliJetDummyGeo::GetInstance();
520 if (fGeom == 0)
521 fGeom = AliJetDummyGeo::GetInstance("EMCAL_COMPLETE","EMCAL");
522
523 // To be setted to run some AliEMCALGeometry functions
524 TGeoManager::Import(geomFile);
525 fGeom->GetTransformationForSM();
526 printf("\n EMCal Geometry set ! \n");
527
528 return kTRUE;
529
530}
531
532//____________________________________________________________________________
533void AliJetAODReader::InitParameters()
534{
535 // Initialise parameters
536 fOpt = fReaderHeader->GetDetector();
be6e5811 537 fNumUnits = fGeom->GetNCells(); // Number of cells in EMCAL
538 if(fDebug>1) printf("\n EMCal parameters initiated ! \n");
539}
540
541//____________________________________________________________________________
542void AliJetAODReader::InitUnitArray()
543{
544 //Initialises unit arrays
545 Int_t nElements = fTpcGrid->GetNEntries();
f23f7908 546 Float_t eta = 0., phi = 0., deltaEta = 0., deltaPhi = 0.;
be6e5811 547 if(fArrayInitialised) fUnitArray->Delete();
548
549 if(fTpcGrid->GetGridType()==0)
550 { // Fill the following quantities :
f23f7908 551 // Good track ID, (Eta,Phi) position ID, eta, phi, energy, px, py, pz, deltaEta, deltaPhi,
be6e5811 552 // detector flag, in/out jet, pt cut, mass, cluster ID)
553 for(Int_t nBin = 1; nBin < nElements+1; nBin++)
554 {
555 // fTpcGrid->GetEtaPhiFromIndex2(nBin,eta,phi);
556 fTpcGrid->GetEtaPhiFromIndex2(nBin,phi,eta);
557 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
f23f7908 558 deltaEta = fTpcGrid->GetDeta();
559 deltaPhi = fTpcGrid->GetDphi();
560 new ((*fUnitArray)[nBin-1]) AliJetUnitArray(nBin-1,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 561 }
562 }
563
564 if(fTpcGrid->GetGridType()==1)
565 {
566 Int_t nGaps = 0;
567 Int_t n0 = 0, n1 = 0, n2 = 0, n3 = 0, n4 = 0;
568
569 if(fDZ)
570 {
571 // Define a grid of cell for the gaps between SM
572 Double_t phimin0 = 0., phimin1 = 0., phimin2 = 0., phimin3 = 0., phimin4 = 0.;
573 Double_t phimax0 = 0., phimax1 = 0., phimax2 = 0., phimax3 = 0., phimax4 = 0.;
574 fGeom->GetPhiBoundariesOfSMGap(0,phimin0,phimax0);
575 fGrid0 = new AliJetGrid(0,95,phimin0,phimax0,-0.7,0.7); // 0.015 x 0.015
576 fGrid0->SetGridType(0);
577 fGrid0->SetMatrixIndexes();
578 fGrid0->SetIndexIJ();
579 n0 = fGrid0->GetNEntries();
580 fGeom->GetPhiBoundariesOfSMGap(1,phimin1,phimax1);
581 fGrid1 = new AliJetGrid(0,95,phimin1,phimax1,-0.7,0.7); // 0.015 x 0.015
582 fGrid1->SetGridType(0);
583 fGrid1->SetMatrixIndexes();
584 fGrid1->SetIndexIJ();
585 n1 = fGrid1->GetNEntries();
586 fGeom->GetPhiBoundariesOfSMGap(2,phimin2,phimax2);
587 fGrid2 = new AliJetGrid(0,95,phimin2,phimax2,-0.7,0.7); // 0.015 x 0.015
588 fGrid2->SetGridType(0);
589 fGrid2->SetMatrixIndexes();
590 fGrid2->SetIndexIJ();
591 n2 = fGrid2->GetNEntries();
592 fGeom->GetPhiBoundariesOfSMGap(3,phimin3,phimax3);
593 fGrid3 = new AliJetGrid(0,95,phimin3,phimax3,-0.7,0.7); // 0.015 x 0.015
594 fGrid3->SetGridType(0);
595 fGrid3->SetMatrixIndexes();
596 fGrid3->SetIndexIJ();
597 n3 = fGrid3->GetNEntries();
598 fGeom->GetPhiBoundariesOfSMGap(4,phimin4,phimax4);
599 fGrid4 = new AliJetGrid(0,95,phimin4,phimax4,-0.7,0.7); // 0.015 x 0.015
600 fGrid4->SetGridType(0);
601 fGrid4->SetMatrixIndexes();
602 fGrid4->SetIndexIJ();
603 n4 = fGrid4->GetNEntries();
604
605 nGaps = n0+n1+n2+n3+n4;
606
607 }
608
609 for(Int_t nBin = 0; nBin < fNumUnits+nElements+nGaps; nBin++)
610 {
611 if(nBin<fNumUnits)
612 {
613 fGeom->EtaPhiFromIndex(nBin, eta, phi); // From EMCal geometry
614 // fEmcalGrid->GetEtaPhiFromIndex2(nBin,phi,eta); // My function from Grid
615 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
f23f7908 616 deltaEta = fEmcalGrid->GetDeta(); // Modify with the exact detector values
617 deltaPhi = fEmcalGrid->GetDphi(); // Modify with the exact detector values
618 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 619 }
620 else {
621 if(nBin>=fNumUnits && nBin<fNumUnits+nElements){
622 fTpcGrid->GetEtaPhiFromIndex2(nBin+1-fNumUnits,phi,eta);
623 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
f23f7908 624 deltaEta = fTpcGrid->GetDeta();
625 deltaPhi = fTpcGrid->GetDphi();
626 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 627 }
628 else {
629 if(fDZ) {
630 if(nBin>=fNumUnits+nElements && nBin<fNumUnits+nElements+nGaps){
631 if(nBin<fNumUnits+nElements+n0)
632 {
f23f7908 633 phi = eta = 0.;
be6e5811 634 fGrid0->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements),phi,eta);
f23f7908 635 deltaEta = fGrid0->GetDeta();
636 deltaPhi = fGrid0->GetDphi();
637 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 638 }
639 else if(nBin>=fNumUnits+nElements+n0 && nBin<fNumUnits+nElements+n0+n1)
640 {
f23f7908 641 phi = eta = 0.;
be6e5811 642 fGrid1->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0),phi,eta);
f23f7908 643 deltaEta = fGrid1->GetDeta();
644 deltaPhi = fGrid1->GetDphi();
645 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 646 }
647 else if(nBin>=fNumUnits+nElements+n0+n1 && nBin<fNumUnits+nElements+n0+n1+n2)
648 {
f23f7908 649 phi = eta = 0.;
be6e5811 650 fGrid2->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1),phi,eta);
f23f7908 651 deltaEta = fGrid2->GetDeta();
652 deltaPhi = fGrid2->GetDphi();
653 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 654 }
655 else if(nBin>=fNumUnits+nElements+n0+n1+n2 && nBin<fNumUnits+nElements+n0+n1+n2+n3)
656 {
f23f7908 657 phi = eta = 0.;
be6e5811 658 fGrid3->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1+n2),phi,eta);
f23f7908 659 deltaEta = fGrid3->GetDeta();
660 deltaPhi = fGrid3->GetDphi();
661 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 662 }
663 else if(nBin>=fNumUnits+nElements+n0+n1+n2+n3 && nBin<fNumUnits+nElements+nGaps)
664 {
f23f7908 665 phi = eta = 0.;
be6e5811 666 fGrid4->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1+n2+n3),phi,eta);
f23f7908 667 deltaEta = fGrid4->GetDeta();
668 deltaPhi = fGrid4->GetDphi();
669 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
be6e5811 670 }
671 }
672 } // end if(fDZ)
673 } // end else 2
674 } // end else 1
675 } // end loop on nBin
676 } // end grid type == 1
677 fArrayInitialised = 1;
678}
0e4398c8 679