1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //-------------------------------------------------------------------------
18 // AOD reader for jet analysis
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
22 // Author: Davide Perrino <davide.perrino@cern.ch>
25 // implemented standard geometry (AliEMCALGeometry) instead of dummy one (AliJetDummyGeo)
26 // moved geometry definition in AliJetReader
27 // marco.bregant@subatech.in2p3.fr
28 //-------------------------------------------------------------------------
31 #include <Riostream.h>
33 #include <TLorentzVector.h>
38 #include <TGeoManager.h>
39 #include <TGeoMatrix.h>
41 #include "AliJetAODReader.h"
42 #include "AliJetAODReaderHeader.h"
43 #include "AliAODEvent.h"
44 #include "AliAODTrack.h"
45 #include "AliAODMCParticle.h"
46 #include "AliEMCALGeometry.h"
47 #include "AliJetAODFillUnitArrayTracks.h"
48 #include "AliJetAODFillUnitArrayEMCalDigits.h"
49 #include "AliJetHadronCorrectionv1.h"
50 #include "AliJetUnitArray.h"
51 #include "AliOADBContainer.h"
53 ClassImp(AliJetAODReader)
55 AliJetAODReader::AliJetAODReader():
69 fApplyElectronCorrection(kFALSE),
70 fApplyMIPCorrection(kTRUE),
71 fApplyFractionHadronicCorrection(kFALSE),
72 fFractionHadronicCorrection(0.3),
86 //____________________________________________________________________________
88 AliJetAODReader::~AliJetAODReader()
106 //____________________________________________________________________________
108 void AliJetAODReader::OpenInputFiles()
110 // Open the necessary input files
111 // chain for the AODs
112 fChain = new TChain("aodTree");
114 // get directory and pattern name from the header
115 const char* dirName=fReaderHeader->GetDirectory();
116 const char* pattern=fReaderHeader->GetPattern();
118 // // Add files matching patters to the chain
120 void *dir = gSystem->OpenDirectory(dirName);
121 const char *name = 0x0;
122 int naod = ((AliJetAODReaderHeader*) fReaderHeader)->GetNaod();
124 while ((name = gSystem->GetDirEntry(dir))){
125 if (a>=naod) continue;
127 if (strstr(name,pattern)){
128 fChain->AddFile(Form("%s/%s/aod.root",dirName,name));
133 gSystem->FreeDirectory(dir);
136 fChain->SetBranchAddress("AOD",&fAOD);
138 int nMax = fChain->GetEntries();
140 printf("\n AliJetAODReader: Total number of events in chain= %d \n",nMax);
142 // set number of events in header
143 if (fReaderHeader->GetLastEvent() == -1)
144 fReaderHeader->SetLastEvent(nMax);
146 Int_t nUsr = fReaderHeader->GetLastEvent();
147 fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
151 //____________________________________________________________________________
153 void AliJetAODReader::ConnectTree(TTree* tree, TObject* /*data*/) {
155 // For AOD reader it's needed only to set the number of events
156 fChain = (TChain*) tree;
158 Int_t nMax = fChain->GetEntries();
159 printf("\n AliJetAODReader: Total number of events in chain= %5d \n", nMax);
160 // set number of events in header
161 if (fReaderHeader->GetLastEvent() == -1)
162 fReaderHeader->SetLastEvent(nMax);
164 Int_t nUsr = fReaderHeader->GetLastEvent();
165 fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
170 Bool_t AliJetAODReader::AcceptAODMCParticle(AliAODMCParticle *mcP,Short_t flag){
173 // filter for charge and for charged and neutral, no detector
175 // physical priamries already selected
177 Int_t pdg = TMath::Abs(mcP->GetPdgCode());
179 // exclude neutrinos anyway
180 if((pdg == 12 || pdg == 14 || pdg == 16)) return kFALSE;
182 if(flag==AliJetAODReaderHeader::kAllMC)return kTRUE;
183 if(flag==AliJetAODReaderHeader::kChargedMC){
184 if(mcP->Charge()==0)return kFALSE;
193 Bool_t AliJetAODReader::FillMomentumArrayMC(){
196 // This routine fetches the MC particles from the AOD
197 // Depending on the flag all particles except neurinos are use
198 // or only charged particles
201 TClonesArray *mcArray = (TClonesArray*)fAOD->FindListObject(AliAODMCParticle::StdBranchName());
203 Printf("%s:%d No MC particle branch found",(char*)__FILE__,__LINE__);
207 Int_t nMC = mcArray->GetEntriesFast();
209 // get number of tracks in event (for the loop)
210 if(fDebug>0)printf("AOD MC tracks: %5d \t", nMC);
212 // temporary storage of signal and pt cut flag
213 Int_t* sflag = new Int_t[nMC];
214 Int_t* cflag = new Int_t[nMC];
216 // get cuts set by user
217 Float_t ptMin = fReaderHeader->GetPtCut();
218 Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
219 Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
224 Short_t readerFlag = ((AliJetAODReaderHeader*)fReaderHeader)->GetReadAODMC();
227 for (Int_t it = 0; it < nMC; it++) {
228 AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
229 if(!track->IsPhysicalPrimary())continue;
231 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
233 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
234 if(!AcceptAODMCParticle(track,readerFlag))continue;
241 fRef->Delete(); // make sure to delete before placement new...
242 new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
244 new ((*fMomentumArray)[mcTrack]) TLorentzVector(p3,p3.Mag());
246 cflag[mcTrack] = ( pt > ptMin ) ? 1: 0;
250 if(fDebug>0)printf("Used MC tracks: %5d \n", mcTrack);
251 // set the signal flags
252 fSignalFlag.Set(mcTrack,sflag);
253 fCutFlag.Set(mcTrack,cflag);
261 //____________________________________________________________________________
263 Bool_t AliJetAODReader::FillMomentumArray()
265 // Clear momentum array
268 fDebug = fReaderHeader->GetDebug();
274 // get cuts set by user
275 Float_t ptMin = fReaderHeader->GetPtCut();
276 Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
277 Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
278 UInt_t filterMask = ((AliJetAODReaderHeader*)fReaderHeader)->GetTestFilterMask();
280 // ----- number of tracks -----
281 Int_t nTracksStd = 0;
282 Short_t mcReaderFlag = ((AliJetAODReaderHeader*)fReaderHeader)->GetReadAODMC();
283 TClonesArray *mcArray = 0x0;
284 // check if we have to read from MC branch
285 if (((AliJetAODReaderHeader*)fReaderHeader)->GetReadStdBranch()) {
287 mcArray = (TClonesArray*)fAOD->FindListObject(AliAODMCParticle::StdBranchName());
289 Printf("%s:%d No MC particle branch found",(char*)__FILE__,__LINE__);
291 nTracksStd = mcArray->GetEntriesFast();
294 nTracksStd = fAOD->GetNTracks();
295 // printf("no. of standard tracks: %i\n", nTracksStd);
298 Int_t nTracksNonStd = 0;
299 TClonesArray *nonStdTracks = 0x0;
300 if (((AliJetAODReaderHeader*)fReaderHeader)->GetReadNonStdBranch()) {
302 (TClonesArray*) fAOD->FindListObject(((AliJetAODReaderHeader*)fReaderHeader)->GetNonStdBranch());
304 nTracksNonStd = nonStdTracks->GetEntries();
305 // printf("no. of non-standard tracks: %i\n", nTracksNonStd);
307 Int_t nTracks = nTracksStd + nTracksNonStd;
309 // temporary storage of signal and pt cut flag
310 Int_t* sflag = new Int_t[nTracks];
311 Int_t* cflag = new Int_t[nTracks];
318 // ----- looping over standard branch -----
320 for (Int_t iTrack = 0; iTrack < nTracksStd; iTrack++) {
321 AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(iTrack);
322 if(!track->IsPhysicalPrimary())continue;
324 p3.SetXYZ(track->Px(),track->Py(),track->Pz());
326 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
327 if(!AcceptAODMCParticle(track,mcReaderFlag))continue;
331 fRef->Delete(); // make sure to delete before placement new...
332 new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
334 new ((*fMomentumArray)[aodTrack]) TLorentzVector(p3,p3.Mag());
336 cflag[aodTrack] = ( pt > ptMin ) ? 1: 0;
342 for (Int_t iTrack = 0; iTrack < nTracksStd; iTrack++) {
343 AliAODTrack *track = fAOD->GetTrack(iTrack);
344 UInt_t status = track->GetStatus();
346 Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
347 p3.SetXYZ(mom[0],mom[1],mom[2]);
350 if (status == 0) continue;
351 if((filterMask>0)&&!(track->TestFilterBit(filterMask)))continue;
352 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
355 fRef->Delete(); // make sure to delete before placement new...
356 new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
358 new ((*fMomentumArray)[aodTrack]) TLorentzVector(p3,p3.Mag());
359 sflag[aodTrack] = (TMath::Abs(track->GetLabel()) < 10000) ? 1 : 0;
360 cflag[aodTrack] = ( pt > ptMin ) ? 1: 0;
366 // ----- reading of non-standard branch -----
367 for (Int_t iTrack = 0; iTrack < nTracksNonStd; iTrack++) {
368 AliVParticle *track = dynamic_cast<AliVParticle*> ((*nonStdTracks)[iTrack]);
372 Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
373 p3.SetXYZ(mom[0],mom[1],mom[2]);
377 // which cuts to apply if not AOD track (e.g. MC) ???
378 AliAODTrack *trackAOD = dynamic_cast<AliAODTrack*> (track);
380 if (trackAOD->GetStatus() == 0)
382 if ((filterMask > 0) && !(trackAOD->TestFilterBit(filterMask)))
385 if ((eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
388 fRef->Delete(); // make sure to delete before placement new...
389 new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
391 new ((*fMomentumArray)[aodTrack]) TLorentzVector(p3,p3.Mag());
392 sflag[aodTrack] = (TMath::Abs(track->GetLabel()) < 10000) ? 1 : 0;
393 cflag[aodTrack] = ( pt > ptMin ) ? 1 : 0;
396 // printf("added non-standard track\n");
399 // set the signal flags
400 fSignalFlag.Set(aodTrack,sflag);
401 fCutFlag.Set(aodTrack,cflag);
409 //__________________________________________________________
410 void AliJetAODReader::SetApplyMIPCorrection(Bool_t val)
413 // Set flag to apply MIP correction fApplyMIPCorrection
414 // - exclusive with fApplyFractionHadronicCorrection
417 fApplyMIPCorrection = val;
418 if(fApplyMIPCorrection == kTRUE)
420 SetApplyFractionHadronicCorrection(kFALSE);
421 printf("Enabling MIP Correction \n");
425 printf("Disabling MIP Correction \n");
429 //__________________________________________________________
430 void AliJetAODReader::SetApplyFractionHadronicCorrection(Bool_t val)
433 // Set flag to apply EMC hadronic correction fApplyFractionHadronicCorrection
434 // - exclusive with fApplyMIPCorrection
437 fApplyFractionHadronicCorrection = val;
438 if(fApplyFractionHadronicCorrection == kTRUE)
440 SetApplyMIPCorrection(kFALSE);
441 printf("Enabling Fraction Hadronic Correction \n");
445 printf("Disabling Fraction Hadronic Correction \n");
449 //__________________________________________________________
450 void AliJetAODReader::SetFractionHadronicCorrection(Double_t val)
453 // Set value to fFractionHadronicCorrection (default is 0.3)
454 // apply EMC hadronic correction fApplyFractionHadronicCorrection
455 // - exclusive with fApplyMIPCorrection
458 fFractionHadronicCorrection = val;
459 if(fFractionHadronicCorrection > 0.0 && fFractionHadronicCorrection <= 1.0)
461 SetApplyFractionHadronicCorrection(kTRUE);
462 printf("Fraction Hadronic Correction %1.3f \n",fFractionHadronicCorrection);
466 SetApplyFractionHadronicCorrection(kFALSE);
470 //____________________________________________________________________________
471 void AliJetAODReader::CreateTasks(TChain* tree)
474 // For reader task initialization
477 fDebug = fReaderHeader->GetDebug();
478 fDZ = fReaderHeader->GetDZ();
481 // Init EMCAL geometry, if needed
484 else Info(" SetEMCALGeometry","was already done.. it's called just once !!");
489 fFillUnitArray = new TTask("fFillUnitArray","Fill unit array jet finder");
490 fFillUAFromTracks = new AliJetAODFillUnitArrayTracks();
491 fFillUAFromTracks->SetReaderHeader(fReaderHeader);
492 fFillUAFromTracks->SetGeom(fGeom);
493 fFillUAFromTracks->SetTPCGrid(fTpcGrid);
494 fFillUAFromTracks->SetEMCalGrid(fEmcalGrid);
498 fFillUAFromTracks->SetGrid0(fGrid0);
499 fFillUAFromTracks->SetGrid1(fGrid1);
500 fFillUAFromTracks->SetGrid2(fGrid2);
501 fFillUAFromTracks->SetGrid3(fGrid3);
502 fFillUAFromTracks->SetGrid4(fGrid4);
504 fFillUAFromTracks->SetApplyMIPCorrection(fApplyMIPCorrection);
505 fFillUAFromTracks->SetHadCorrector(fHadCorr);
506 fFillUAFromEMCalDigits = new AliJetAODFillUnitArrayEMCalDigits();
507 fFillUAFromEMCalDigits->SetReaderHeader(fReaderHeader);
508 fFillUAFromEMCalDigits->SetGeom(fGeom);
509 fFillUAFromEMCalDigits->SetTPCGrid(fTpcGrid);
510 fFillUAFromEMCalDigits->SetEMCalGrid(fEmcalGrid);
511 fFillUAFromEMCalDigits->SetApplyFractionHadronicCorrection(fApplyFractionHadronicCorrection);
512 fFillUAFromEMCalDigits->SetFractionHadronicCorrection(fFractionHadronicCorrection);
513 fFillUAFromEMCalDigits->SetApplyElectronCorrection(fApplyElectronCorrection);
514 // Add the task to global task
515 fFillUnitArray->Add(fFillUAFromTracks);
516 fFillUnitArray->Add(fFillUAFromEMCalDigits);
517 fFillUAFromTracks->SetActive(kFALSE);
518 fFillUAFromEMCalDigits->SetActive(kFALSE);
520 cout << "Tasks instantiated at that stage ! " << endl;
521 cout << "You can loop over events now ! " << endl;
525 //____________________________________________________________________________
526 Bool_t AliJetAODReader::ExecTasks(Bool_t procid, TRefArray* refArray)
530 // Fill the reader part
534 fRefArray = refArray;
535 //(not used ?) Int_t nEntRef = fRefArray->GetEntries();
536 //(not used ?) Int_t nEntUnit = fUnitArray->GetEntries();
538 // clear momentum array
542 fDebug = fReaderHeader->GetDebug();
543 fOpt = fReaderHeader->GetDetector();
549 // TPC only or Digits+TPC or Clusters+TPC
550 if(fOpt%2==!0 && fOpt!=0){
551 fFillUAFromTracks->SetAOD(fAOD);
552 fFillUAFromTracks->SetActive(kTRUE);
553 fFillUAFromTracks->SetUnitArray(fUnitArray);
554 fFillUAFromTracks->SetRefArray(fRefArray);
555 fFillUAFromTracks->SetReferences(fRef);
556 fFillUAFromTracks->SetSignalFlag(fSignalFlag);
557 fFillUAFromTracks->SetCutFlag(fCutFlag);
558 fFillUAFromTracks->SetProcId(fProcId);
559 // fFillUAFromTracks->ExecuteTask("tpc"); // => Temporarily changed
560 fFillUAFromTracks->Exec("tpc");
562 fNumCandidate = fFillUAFromTracks->GetMult();
563 fNumCandidateCut = fFillUAFromTracks->GetMultCut();
565 fSignalFlag = fFillUAFromTracks->GetSignalFlag();
566 fCutFlag = fFillUAFromTracks->GetCutFlag();
569 // Digits only or Digits+TPC
570 if(fOpt>=2 && fOpt<=3){
571 fFillUAFromEMCalDigits->SetAOD(fAOD);
572 fFillUAFromEMCalDigits->SetActive(kTRUE);
573 fFillUAFromEMCalDigits->SetUnitArray(fUnitArray);
574 fFillUAFromEMCalDigits->SetRefArray(fRefArray);
575 fFillUAFromEMCalDigits->SetProcId(fFillUAFromTracks->GetProcId());
576 fFillUAFromEMCalDigits->SetInitMult(fFillUAFromTracks->GetMult());
577 fFillUAFromEMCalDigits->SetInitMultCut(fFillUAFromTracks->GetMultCut());
578 fFillUAFromEMCalDigits->Exec("digits"); // => Temporarily added
579 fNumCandidate = fFillUAFromEMCalDigits->GetMult();
580 fNumCandidateCut = fFillUAFromEMCalDigits->GetMultCut();
583 // fFillUnitArray->ExecuteTask(); // => Temporarily commented
590 //____________________________________________________________________________
591 void AliJetAODReader::InitParameters()
593 // Initialise parameters
594 fOpt = fReaderHeader->GetDetector();
595 fNumUnits = fGeom->GetEMCGeometry()->GetNCells(); // Number of cells in EMCAL
596 if(fDebug>1) printf("\n EMCal parameters initiated ! \n");
599 //____________________________________________________________________________
600 void AliJetAODReader::InitUnitArray()
602 //Initialises unit arrays
603 Int_t nElements = fTpcGrid->GetNEntries();
604 Float_t eta = 0., phi = 0., deltaEta = 0., deltaPhi = 0.;
605 if(fArrayInitialised) fUnitArray->Delete();
607 if(fTpcGrid->GetGridType()==0)
608 { // Fill the following quantities :
609 // Good track ID, (Eta,Phi) position ID, eta, phi, energy, px, py, pz, deltaEta, deltaPhi,
610 // detector flag, in/out jet, pt cut, mass, cluster ID)
611 for(Int_t nBin = 1; nBin < nElements+1; nBin++)
613 // fTpcGrid->GetEtaPhiFromIndex2(nBin,eta,phi);
614 fTpcGrid->GetEtaPhiFromIndex2(nBin,phi,eta);
615 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
616 deltaEta = fTpcGrid->GetDeta();
617 deltaPhi = fTpcGrid->GetDphi();
618 new ((*fUnitArray)[nBin-1]) AliJetUnitArray(nBin-1,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
622 if(fTpcGrid->GetGridType()==1)
625 Int_t n0 = 0, n1 = 0, n2 = 0, n3 = 0, n4 = 0;
629 // Define a grid of cell for the gaps between SM
630 Double_t phimin0 = 0., phimin1 = 0., phimin2 = 0., phimin3 = 0., phimin4 = 0.;
631 Double_t phimax0 = 0., phimax1 = 0., phimax2 = 0., phimax3 = 0., phimax4 = 0.;
632 fGeom->GetEMCGeometry()->GetPhiBoundariesOfSMGap(0,phimin0,phimax0);
633 fGrid0 = new AliJetGrid(0,95,phimin0,phimax0,-0.7,0.7); // 0.015 x 0.015
634 fGrid0->SetGridType(0);
635 fGrid0->SetMatrixIndexes();
636 fGrid0->SetIndexIJ();
637 n0 = fGrid0->GetNEntries();
638 fGeom->GetEMCGeometry()->GetPhiBoundariesOfSMGap(1,phimin1,phimax1);
639 fGrid1 = new AliJetGrid(0,95,phimin1,phimax1,-0.7,0.7); // 0.015 x 0.015
640 fGrid1->SetGridType(0);
641 fGrid1->SetMatrixIndexes();
642 fGrid1->SetIndexIJ();
643 n1 = fGrid1->GetNEntries();
644 fGeom->GetEMCGeometry()->GetPhiBoundariesOfSMGap(2,phimin2,phimax2);
645 fGrid2 = new AliJetGrid(0,95,phimin2,phimax2,-0.7,0.7); // 0.015 x 0.015
646 fGrid2->SetGridType(0);
647 fGrid2->SetMatrixIndexes();
648 fGrid2->SetIndexIJ();
649 n2 = fGrid2->GetNEntries();
650 fGeom->GetEMCGeometry()->GetPhiBoundariesOfSMGap(3,phimin3,phimax3);
651 fGrid3 = new AliJetGrid(0,95,phimin3,phimax3,-0.7,0.7); // 0.015 x 0.015
652 fGrid3->SetGridType(0);
653 fGrid3->SetMatrixIndexes();
654 fGrid3->SetIndexIJ();
655 n3 = fGrid3->GetNEntries();
656 fGeom->GetEMCGeometry()->GetPhiBoundariesOfSMGap(4,phimin4,phimax4);
657 fGrid4 = new AliJetGrid(0,95,phimin4,phimax4,-0.7,0.7); // 0.015 x 0.015
658 fGrid4->SetGridType(0);
659 fGrid4->SetMatrixIndexes();
660 fGrid4->SetIndexIJ();
661 n4 = fGrid4->GetNEntries();
663 nGaps = n0+n1+n2+n3+n4;
667 for(Int_t nBin = 0; nBin < fNumUnits+nElements+nGaps; nBin++)
671 fGeom->EtaPhiFromIndex(nBin, eta, phi); // From EMCal geometry
672 // fEmcalGrid->GetEtaPhiFromIndex2(nBin,phi,eta); // My function from Grid
673 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
674 deltaEta = fEmcalGrid->GetDeta(); // Modify with the exact detector values
675 deltaPhi = fEmcalGrid->GetDphi(); // Modify with the exact detector values
676 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
679 if(nBin>=fNumUnits && nBin<fNumUnits+nElements){
680 fTpcGrid->GetEtaPhiFromIndex2(nBin+1-fNumUnits,phi,eta);
681 phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
682 deltaEta = fTpcGrid->GetDeta();
683 deltaPhi = fTpcGrid->GetDphi();
684 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
688 if(nBin>=fNumUnits+nElements && nBin<fNumUnits+nElements+nGaps){
689 if(nBin<fNumUnits+nElements+n0)
692 fGrid0->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements),phi,eta);
693 deltaEta = fGrid0->GetDeta();
694 deltaPhi = fGrid0->GetDphi();
695 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
697 else if(nBin>=fNumUnits+nElements+n0 && nBin<fNumUnits+nElements+n0+n1)
700 fGrid1->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0),phi,eta);
701 deltaEta = fGrid1->GetDeta();
702 deltaPhi = fGrid1->GetDphi();
703 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
705 else if(nBin>=fNumUnits+nElements+n0+n1 && nBin<fNumUnits+nElements+n0+n1+n2)
708 fGrid2->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1),phi,eta);
709 deltaEta = fGrid2->GetDeta();
710 deltaPhi = fGrid2->GetDphi();
711 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
713 else if(nBin>=fNumUnits+nElements+n0+n1+n2 && nBin<fNumUnits+nElements+n0+n1+n2+n3)
716 fGrid3->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1+n2),phi,eta);
717 deltaEta = fGrid3->GetDeta();
718 deltaPhi = fGrid3->GetDphi();
719 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
721 else if(nBin>=fNumUnits+nElements+n0+n1+n2+n3 && nBin<fNumUnits+nElements+nGaps)
724 fGrid4->GetEtaPhiFromIndex2(nBin+1-(fNumUnits+nElements+n0+n1+n2+n3),phi,eta);
725 deltaEta = fGrid4->GetDeta();
726 deltaPhi = fGrid4->GetDphi();
727 new ((*fUnitArray)[nBin]) AliJetUnitArray(nBin,0,eta,phi,0.,deltaEta,deltaPhi,kTpc,kOutJet,kPtSmaller,kPtSmaller,kBad,0.,-1);
733 } // end loop on nBin
734 } // end grid type == 1
735 fArrayInitialised = 1;