JETAN/AliJetFillUnitArrayTracks.cxx

   1
   2 /**************************************************************************
   3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   4  *                                                                        *
   5  * Author: The ALICE Off-line Project.                                    *
   6  * Contributors are mentioned in the code where appropriate.              *
   7  *                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is          *
  14  * provided "as is" without express or implied warranty.                  *
  15  **************************************************************************/
  16
  17 //======================================================================
  18 // *** November 2006
  19 // Author: magali.estienne@ires.in2p3.fr
  20 // 1) Define 2 grids and take (eta,phi) from the grid or use the grid for the TPC and
  21 // EtaPhiFromIndex and TowerIndexFromEtaPhi for the particles in EMCAL acceptance
  22 //     2 options are possible : fGrid==0, work only with full TPC acceptance (for the moment)
  23 //                              fGrid==1, work with a part of the TPC acceptance
  24 // 2) Try to implement 2 full grids for TPC and EMCal separately and to merge them
  25 // 3) Need to include Dead-zone -> Wait for exact positions in the new detector geometry
  26 // Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
  27 //======================================================================
  28 // ***September 2006
  29 // TTask : Fill Unit Array for the Tracks information
  30 // Called by ESD reader for jet analysis
  31 // Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
  32 //======================================================================
  33 // *** July 2006
  34 // 1) When the tracks are in the EMCal acceptance, the functions EtaPhiFromIndex
  35 // and TowerIndexFromEtaPhi in the AliEMCALGeometry class are used to extract the
  36 // index or the eta, phi position of a grid.
  37 // 2) Define a grid for TPC
  38 // Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
  39 //======================================================================
  40 // ***July 2006
  41 // Fill Unit Array class
  42 // Class used by AliJetESDReader to fill a UnitArray from the information extracted
  43 // from the particle tracks
  44 // Author: magali.estienne@ires.in2p3.fr
  45 //======================================================================
  46
  47
  48 // --- Standard library ---
  49 #include <Riostream.h>
  50
  51 // --- ROOT system ---
  52 #include <TSystem.h>
  53 #include <TLorentzVector.h>
  54 #include <TVector3.h>
  55 //#include "Math/Vector3D.h"
  56 //#include "Math/Vector3Dfwd.h"
  57 #include "TTask.h"
  58 #include <TGeoManager.h>
  59 #include <TMatrixD.h>
  60 #include <TArrayD.h>
  61
  62 // --- AliRoot header files ---
  63 #include "AliJetFinder.h"
  64 #include "AliJetReaderHeader.h"
  65 #include "AliJetReader.h"
  66 #include "AliJetESDReader.h"
  67 #include "AliJetESDReaderHeader.h"
  68 #include "AliESD.h"
  69 #include "AliEMCALGeometry.h"
  70 #include "AliJetUnitArray.h"
  71 #include "AliJetFillUnitArrayTracks.h"
  72 #include "AliJetGrid.h"
  73
  74 ClassImp(AliJetFillUnitArrayTracks)
  75
  76 AliEMCALGeometry *AliJetFillUnitArrayTracks::fGeom=0;
  77
  78 //_____________________________________________________________________________
  79 AliJetFillUnitArrayTracks::AliJetFillUnitArrayTracks()
  80   : TTask("AliJetFillUnitArrayTracks","Fill Unit Array with tpc/its and emcal information")
  81 {
  82   // constructor
  83    fNIn = 0;
  84   fOpt = 0;
  85   fDebug = 0;
  86   fNphi = 0;
  87   fNeta = 0;
  88   fPhiMin = 0;
  89   fPhiMax = 0;
  90   fEtaMin = 0;
  91   fEtaMax = 0;
  92   fEtaBinInTPCAcc = 0;
  93   fPhiBinInTPCAcc = 0;
  94   fEtaBinInEMCalAcc = 0;
  95   fPhiBinInEMCalAcc = 0;
  96   fNbinPhi = 0;
  97 }
  98
  99 //____________________________________________________________________________
 100 void AliJetFillUnitArrayTracks::SetEMCALGeometry()
 101 {
 102   // Set EMCAL geometry information
 103   fGeom = AliEMCALGeometry::GetInstance();
 104   if (fGeom == 0)
 105     fGeom = AliEMCALGeometry::GetInstance("SHISH_77_TRD1_2X2_FINAL_110DEG","EMCAL");
 106   if(fDebug>1) printf("\n EMCAL Geometry setted ! \n");
 107 }
 108
 109 //____________________________________________________________________________
 110 void AliJetFillUnitArrayTracks::InitParameters()
 111 {
 112   fHCorrection    = 0;     // For hadron correction
 113   fHadCorr        = 0;     // For hadron correction
 114   fNumUnits = fGeom->GetNCells();      // Number of towers in EMCAL
 115   fDebug = fReaderHeader->GetDebug();
 116
 117   fEtaMinCal = fGeom->GetArm1EtaMin();
 118   fEtaMaxCal = fGeom->GetArm1EtaMax();
 119   fPhiMinCal = fGeom->GetArm1PhiMin();
 120   fPhiMaxCal = fGeom->GetArm1PhiMax()-10.; // A verifier quelle doit etre la derniere valeur !
 121
 122   if(fDebug>30){
 123     cout << "fEtaMinCal : " << fEtaMinCal << endl;
 124     cout << "fEtaMaxCal : " << fEtaMaxCal << endl;
 125     cout << "fPhiMinCal : " << fPhiMinCal << endl;
 126     cout << "fPhiMaxCal : " << fPhiMaxCal << endl;
 127   }
 128
 129   fTPCGrid->GetAccParam(fNphi,fNeta,fPhiMin,
 130                         fPhiMax,fEtaMin,fEtaMax);
 131   fTPCGrid->GetBinParam(fPhiBinInTPCAcc,fEtaBinInTPCAcc,
 132                         fPhiBinInEMCalAcc,fEtaBinInEMCalAcc,fNbinPhi);
 133
 134   fEta   = fTPCGrid->GetArrayEta();
 135   fPhi   = fTPCGrid->GetArrayPhi();
 136   fIndex = fTPCGrid->GetIndexObject();
 137
 138   if(fDebug>20){
 139     for(Int_t i=0; i<fNphi+1; i++) cout << "phi[" << i << "] : " << (*fPhi)[i] << endl;
 140     for(Int_t i=0; i<fNeta+1; i++) cout << "eta[" << i << "] : " << (*fEta)[i] << endl;
 141
 142     for(Int_t i=0; i<fNphi+1; i++)
 143       for(Int_t j=0; j<fNeta+1; j++) {cout << "fIndex[" << i << "," << j << "] : " <<
 144           (*fIndex)(i,j) << endl; }
 145   }
 146   if(fDebug>1) printf("\n Parameters initiated ! \n");
 147 }
 148
 149 //_____________________________________________________________________________
 150 AliJetFillUnitArrayTracks::~AliJetFillUnitArrayTracks()
 151 {
 152   // destructor
 153 }
 154
 155 //_____________________________________________________________________________
 156 void AliJetFillUnitArrayTracks::Exec(Option_t* option)
 157 {
 158   //
 159   // Main method.
 160   // Explain
 161
 162   fDebug = fReaderHeader->GetDebug();
 163   if(fDebug>1) printf("In AliJetFillUnitArrayTracks !");
 164   if(fDebug>3) printf("\nfGeom->GetEntries() = %d\n", fGeom->GetNCells());
 165   // Set EMCal Geometry
 166   SetEMCALGeometry();
 167   // Set parameters
 168   InitParameters();
 169
 170   TClonesArray *lvArray = fMomentumArray;     // Correct checked !
 171   Int_t nInT =  lvArray->GetEntries();        // Correct checked !
 172   Float_t ptMin = fReaderHeader->GetPtCut();  // Correct checked !
 173
 174   // sflag -> not yet implemented !!!!
 175
 176   if(fDebug>3) cout << "nInT : " << nInT << endl;
 177
 178   if (nInT == 0) return;
 179
 180   // local arrays for input
 181   Float_t* enT  = new Float_t[nInT];
 182   Float_t* ptT  = new Float_t[nInT];
 183   Float_t* etaT = new Float_t[nInT];
 184   Float_t* phiT = new Float_t[nInT];
 185   Float_t* thetaT = new Float_t[nInT];
 186
 187   Int_t trackInEmcalAcc = 0;
 188   Int_t trackInTpcAcc = 0;
 189   Int_t trackInTpcAccOnly = 0;
 190
 191   Int_t nElements = fTPCGrid->GetNEntries();
 192   Int_t nElements2 = fEMCalGrid->GetNEntries();
 193   fGrid = fTPCGrid->GetGridType();
 194
 195   if(fDebug>3){
 196     cout << "nElements : " << nElements << endl;
 197     cout << "nElements2 : " << nElements2 << endl;
 198     cout << "fNumUnits : " << fNumUnits << endl;
 199     cout << "sum : " << fNumUnits+nElements << endl;
 200   }
 201
 202   // Set energy exactly to zero
 203   if(fGrid==0)
 204     for(Int_t k=0; k<nElements; k++)
 205       fUnitArray[k].SetUnitEnergy(0.);
 206
 207   if(fGrid==1)
 208     for(Int_t k=0; k<fNumUnits+nElements; k++)
 209       fUnitArray[k].SetUnitEnergy(0.);
 210
 211   // load input vectors
 212   for (Int_t i = 0; i < nInT; i++)
 213     {
 214       TLorentzVector *lv = (TLorentzVector*) lvArray->At(i);
 215       enT[i]  = lv->Energy();
 216       ptT[i]  = lv->Pt();
 217       etaT[i] = lv->Eta();
 218       phiT[i] = ((lv->Phi() < 0) ? (lv->Phi()) + 2. * TMath::Pi() : lv->Phi());
 219       thetaT[i] = 2.0*TMath::ATan(TMath::Exp(-etaT[i]));
 220
 221       if(fDebug>20){
 222         cout << "enT[" << i << "] : " <<  enT[i] << endl;
 223         cout << "ptT[" << i << "] : " <<  ptT[i] << endl;
 224         cout << "etaT[" << i << "] : " <<  etaT[i] << endl;
 225         cout << "phiT[" << i << "] : " <<  phiT[i] << endl;
 226         cout << "thetaT[" << i << "] : " <<  thetaT[i] << endl;
 227         cout << "fEtaMinCal : " << fEtaMinCal << ", fEtaMaxCal : " << fEtaMaxCal << endl;
 228         cout << "fPhiMinCal : " << fPhiMinCal << ", fPhiMaxCal : " << fPhiMaxCal << endl;
 229         cout << "fEtaMin : " << fEtaMin << ", fEtaMax : " << fEtaMax << endl;
 230         cout << "fPhiMin : " << fPhiMin << ", fPhiMax : " << fPhiMax << endl;
 231       }
 232
 233       if(fGrid==0)
 234         {
 235           // For the moment, only TPC filled from its grid in its total acceptance
 236           if(fDebug>2)
 237             cout << "In total TPC acceptance +++++++++++++++++++++++++++++++++++++++++++" << endl;
 238
 239             trackInTpcAccOnly += 1;
 240
 241             Int_t idTPC = fTPCGrid->GetIndex(phiT[i],etaT[i]);
 242
 243             Float_t unitEnergy = 0.;
 244             unitEnergy = fUnitArray[idTPC].GetUnitEnergy();
 245
 246             if(unitEnergy > 0. && fDebug >10){
 247               cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 248               cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 249               cout << "TPC id : " << idTPC << endl;
 250               cout << "etaT[" << i << "] : " << etaT[i] << endl;
 251               cout << "phiT[" << i << "] : " << phiT[i] << endl;
 252               cout << "unitEnergy in TPC acceptance : " << unitEnergy << endl;
 253               cout << "fUnitArray[idTPC].GetUnitEnergy(): " <<
 254                 fUnitArray[idTPC].GetUnitEnergy() << endl;
 255               cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 256               cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 257             }
 258
 259             // Fill energy in TPC acceptance
 260             fUnitArray[idTPC].SetUnitEnergy(unitEnergy + ptT[i]);
 261
 262             if(fDebug > 10){
 263               cout << "ptT[" << i << "] : " << ptT[i] << endl;
 264               cout << "unitEnergy in TPC acceptance after : " <<
 265                 fUnitArray[idTPC].GetUnitEnergy() << endl;
 266             }
 267
 268             // Pt cut flag
 269             if(fUnitArray[idTPC].GetUnitEnergy()<ptMin)
 270               fUnitArray[idTPC].SetUnitCutFlag(kPtSmaller);
 271             else fUnitArray[idTPC].SetUnitCutFlag(kPtHigher);
 272             // Detector flag
 273             fUnitArray[idTPC].SetUnitDetectorFlag(kTpc);
 274         }
 275
 276       if(fGrid==1)
 277         {
 278           // Fill track information in EMCAL acceptance
 279           if((etaT[i] >= fEtaMin && etaT[i] <= fEtaMax) &&
 280              (phiT[i] >= fPhiMin && phiT[i] <= fPhiMax))// &&
 281             //   GetCutFlag(i) == 1)
 282             //   ptT[i] > ptMin)
 283             {
 284               trackInEmcalAcc += 1;
 285
 286               if(fDebug>20){
 287                 cout << "before : " << endl;
 288                 cout << "etaT[i] : " << etaT[i] << endl;
 289                 cout << "phiT[i] : " << phiT[i]*180/TMath::Pi() << endl;
 290               }
 291
 292               // This function should be modified soon
 293 //            Int_t towerID = fGeom->TowerIndexFromEtaPhi(etaT[i],180.0/TMath::Pi()*phiT[i]); // Obsolete
 294               Int_t towerID = fGeom->TowerIndexFromEtaPhi2(etaT[i],180.0/TMath::Pi()*phiT[i]); // Mine modified
 295 //            Int_t towerID = fEMCalGrid->GetIndexFromEtaPhi(phiT[i],etaT[i]);  // Using an EMCal grid -> calculated
 296 //            Int_t towerID = fEMCalGrid->GetIndex(phiT[i],etaT[i]);  // Using an EMCal grid -> tabulated (faster)
 297
 298               Float_t unitEnergy = fUnitArray[towerID].GetUnitEnergy();
 299
 300               if(fDebug>20) {
 301                 cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 302                 cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 303                 cout << "after : " << endl;
 304                 cout << "towerID : " << towerID << endl;
 305                 cout << "etaT[i] : " << etaT[i] << endl;
 306                 cout << "phiT[i](rad) : " << phiT[i] << endl;
 307                 cout << "phiT[i] : " << phiT[i]*180/TMath::Pi() << endl;
 308                 cout << "unitEnergy in emcal acceptance : " << unitEnergy << endl;
 309                 cout << "fHadCorr : " << fHadCorr << endl;
 310                 cout << "fHCorrection : " << fHCorrection << endl;
 311                 cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 312                 cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 313               }
 314
 315               //OLD WAY:   //Do Hadron Correction
 316               if (fHCorrection != 0)
 317                 {
 318                   Float_t   hCEnergy = fHadCorr->GetEnergy(enT[i], (Double_t)etaT[i]);
 319                   unitEnergy -= hCEnergy*TMath::Sin(thetaT[i]);
 320                   if(fDebug>20){
 321                     cout << "Inside loop for hadron correction ==========================" << endl;
 322                     cout << "enT[" << i << "] : " << enT[i] << endl;
 323                     cout << "ptT[" << i << "] : " << ptT[i] << endl;
 324                     cout << "etaT[" << i << "] : " << etaT[i] << endl;
 325                     cout << "Energy correction : " << hCEnergy << endl;
 326                     cout << "unitEnergy : " << unitEnergy << endl;
 327                     cout << "fUnitArray[towerID].GetUnitEnergy() : " <<
 328                       fUnitArray[towerID].GetUnitEnergy() << endl;
 329                   }
 330                 } //end Hadron Correction loop
 331
 332               fUnitArray[towerID].SetUnitEnergy(unitEnergy + ptT[i]);
 333
 334               // Put a pt cut flag
 335               if(fUnitArray[towerID].GetUnitEnergy()<ptMin)
 336                 fUnitArray[towerID].SetUnitCutFlag(kPtSmaller);
 337               else fUnitArray[towerID].SetUnitCutFlag(kPtHigher);
 338               // Detector flag
 339               fUnitArray[towerID].SetUnitDetectorFlag(kTpc);
 340
 341               if(fDebug>10){
 342                 cout << "After pT filled ===============================================" << endl;
 343                 cout << "PtCut : " << ptMin << endl;
 344                 cout << "ptT[" << i << "] : " << ptT[i] << endl;
 345                 cout << "unitEnergy : " << unitEnergy << endl;
 346                 cout << "fUnitArray[towerID].GetUnitEnergy() : " << fUnitArray[towerID].GetUnitEnergy() << endl;
 347                 cout << "fUnitArray[towerID].GetUnitCutFlag() : " << fUnitArray[towerID].GetUnitCutFlag() << endl;
 348                 cout << "fUnitArray[towerID].GetUnitEta() : " << fUnitArray[towerID].GetUnitEta() << endl;
 349                 cout << "fUnitArray[towerID].GetUnitPhi() : " << fUnitArray[towerID].GetUnitPhi() << endl;
 350               }
 351             } // end loop on EMCal acceptance cut + tracks quality
 352           else{
 353             // Outside EMCal acceptance
 354             if(fDebug>2)
 355               cout << "Outside EMCal acceptance +++++++++++++++++++++++++++++++++++++++++++" << endl;
 356
 357             trackInTpcAcc += 1;
 358
 359             //  Int_t idTPC = GetIndexFromEtaPhi(etaT[i],phiT[i]);
 360             Int_t idTPC = fTPCGrid->GetIndex(phiT[i],etaT[i]);
 361
 362             Float_t unitEnergy2 = 0.;
 363             unitEnergy2 = fUnitArray[fNumUnits-1+idTPC].GetUnitEnergy(); // check if fNumUnits or fNumUnits-1
 364
 365             if(unitEnergy2 > 0. && fDebug >10){
 366               cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 367               cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 368               cout << "TPC id : " << idTPC << endl;
 369               cout << "Total id : " << fNumUnits-1+idTPC << endl;
 370               cout << "etaT[" << i << "] : " << etaT[i] << endl;
 371               cout << "phiT[" << i << "] : " << phiT[i] << endl;
 372               cout << "unitEnergy outside emcal acceptance : " << unitEnergy2 << endl;
 373               cout << "fUnitArray[fNumUnits-1+idTPC].GetUnitEnergy(): " <<
 374                 fUnitArray[fNumUnits-1+idTPC].GetUnitEnergy() << endl;
 375               cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 376               cout << "§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§§" << endl;
 377             }
 378
 379             // Fill energy outside emcal acceptance
 380             fUnitArray[fNumUnits-1+idTPC].SetUnitEnergy(unitEnergy2 + ptT[i]);
 381
 382             // Pt cut flag
 383             if(fUnitArray[fNumUnits-1+idTPC].GetUnitEnergy()<ptMin)
 384               fUnitArray[fNumUnits-1+idTPC].SetUnitCutFlag(kPtSmaller);
 385             else fUnitArray[fNumUnits-1+idTPC].SetUnitCutFlag(kPtHigher);
 386             // Detector flag
 387             fUnitArray[fNumUnits-1+idTPC].SetUnitDetectorFlag(kTpc);
 388           }
 389         } // end fGrid==1
 390     } // end loop on entries (tpc tracks)
 391
 392   if(fDebug>3){
 393     printf("Number of tracks in EMCAL acceptance: %d\n", trackInEmcalAcc);
 394     printf("Number of tracks in TPC acceptance: %d\n", trackInTpcAcc);
 395   }
 396
 397   if(fGrid==0) fNIn = trackInTpcAccOnly;
 398   if(fGrid==1) fNIn = trackInEmcalAcc+trackInTpcAcc;
 399
 400   fOpt = fReaderHeader->GetDetector();
 401
 402   if(fDebug>1) printf("fOpt in FillUnitArrayFromTPCTracks : %d\n", fOpt);
 403
 404       if(fOpt==1 || option=="tpc") // if only TPC
 405         { //Set all unit flags, Eta, Phi
 406
 407           if(fGrid==0)
 408             {
 409               if(fDebug>1) printf("In if(fOpt==1 || option==tpc)\n");
 410               for(Int_t i=0; i<nElements; i++)
 411                 {
 412                   Float_t eta = 0.;
 413                   Float_t phi = 0.;
 414                   if (fDebug>10) Info("FillUnitArray","Setting all units outside jets \n");
 415                   fUnitArray[i].SetUnitFlag(kOutJet); // returns 0, 1, 2...
 416                   fUnitArray[i].SetUnitEntries(nElements);
 417                   fUnitArray[i].SetUnitID(i);
 418                   fTPCGrid->GetEtaPhiFromIndex2(fUnitArray[i].GetUnitID(),phi,eta);
 419                   fUnitArray[i].SetUnitEta(eta);
 420                   fUnitArray[i].SetUnitPhi(phi);
 421                   fUnitArray[i].SetUnitDeta(fTPCGrid->GetDeta());
 422                   fUnitArray[i].SetUnitDphi(fTPCGrid->GetDphi());
 423                 }
 424             }
 425
 426           if(fGrid==1)
 427             {
 428               if(fDebug>1) printf("In if(fOpt==1 || option==tpc)\n");
 429               for(Int_t i=0; i<fNumUnits+nElements; i++)
 430                 {
 431                   Float_t eta = 0.;
 432                   Float_t phi = 0.;
 433                   if (fDebug>10) Info("FillUnitArray","Setting all units outside jets \n");
 434                   //Set all units to be outside a jet initially
 435                   fUnitArray[i].SetUnitFlag(kOutJet); // returns 0, 1, 2...
 436                   fUnitArray[i].SetUnitEntries(fNumUnits+nElements);
 437                   fUnitArray[i].SetUnitID(i);
 438
 439                   if(fUnitArray[i].GetUnitID()<fNumUnits)
 440                     {
 441                       // fGeom->EtaPhiFromIndex2(fUnitArray[i].GetUnitID(), eta, phi); // My function in HEADPythia63
 442                       fGeom->EtaPhiFromIndex(fUnitArray[i].GetUnitID(), eta, phi); // From EMCal geometry
 443                       phi = ((phi < 0) ? phi + 2. * TMath::Pi() : phi);
 444                       // fEMCalGrid->GetEtaPhiFromIndex2(i,phi,eta); // My function from Grid
 445                       // fEMCalGrid->GetEtaPhiFromIndex2(fUnitArray[i].GetUnitID(),phi,eta); // My function from Grid
 446                       fUnitArray[i].SetUnitEta(eta);
 447                       //fUnitArray[i].SetUnitPhi(phi*TMath::Pi()/180.0);
 448                       fUnitArray[i].SetUnitPhi(phi);
 449                       fUnitArray[i].SetUnitDeta(fEMCalGrid->GetDeta());
 450                       fUnitArray[i].SetUnitDphi(fEMCalGrid->GetDphi());
 451                     }
 452                   else {
 453                     fTPCGrid->GetEtaPhiFromIndex2(fUnitArray[i].GetUnitID()+1-fNumUnits,phi,eta);
 454                     fUnitArray[i].SetUnitEta(eta);
 455                     fUnitArray[i].SetUnitPhi(phi);
 456                     fUnitArray[i].SetUnitDeta(fTPCGrid->GetDeta());
 457                     fUnitArray[i].SetUnitDphi(fTPCGrid->GetDphi());
 458
 459                     if(fDebug>10)
 460                       {
 461                         if(fUnitArray[i].GetUnitEnergy()!=0.){
 462                           cout << "(fUnitArray[" << i << "].GetUnitID()+1-fNumUnits : " <<
 463                             fUnitArray[i].GetUnitID()+1-fNumUnits << endl;
 464                           cout << "(fUnitArray[" << i << "].GetUnitEnergy() : " <<
 465                             fUnitArray[i].GetUnitEnergy() << endl;
 466                           cout << "(fUnitArray[" << i << "].GetUnitEta() : " <<
 467                             fUnitArray[i].GetUnitEta() << endl;
 468                           cout << "(fUnitArray[" << i << "].GetUnitPhi() : " <<
 469                             fUnitArray[i].GetUnitPhi() << endl;
 470                         }
 471                       }
 472                   }
 473                   fUnitArray[i].SetUnitClusterID(0);
 474                 }//end loop over all units in array (same as all towers in EMCAL)
 475             }
 476
 477           if(fDebug>20)
 478             {
 479               for(Int_t i=0; i<fNumUnits+nElements; i++)
 480                 {
 481                   if(fUnitArray[i].GetUnitEnergy()!=0.){
 482                     cout << "######################################################### " << endl;
 483                     cout << "Final UnitArray filled with energy != 0" << i << endl;
 484                     cout << "Pointeur UnitArray : " << fUnitArray << " ID : " <<
 485                       fUnitArray[i].GetUnitID() << " Energy : " << fUnitArray[i].GetUnitEnergy() <<
 486                       " Eta : " << fUnitArray[i].GetUnitEta() << " Phi : " << fUnitArray[i].GetUnitPhi() << endl;
 487                   }
 488                 }
 489             }
 490
 491         } // end  if(fOpt==1 || option=="tpc")
 492
 493       delete[] enT;
 494       delete[] ptT;
 495       delete[] etaT;
 496       delete[] phiT;
 497       delete[] thetaT;
 498
 499 }
 500
 501 //__________________________________________________________
 502 void AliJetFillUnitArrayTracks::GetEtaPhiFromIndex(Int_t index, Float_t &eta, Float_t &phi)
 503 {
 504   for(Int_t j=0; j<fNphi+1; j++) {
 505     for(Int_t i=0; i<fNeta+1; i++) {
 506
 507       // TPC grid only
 508       //-------------------------------------
 509       if(fGrid==0) {
 510         if(j*(fNeta+1)+i == index) {
 511           eta = fEta2->At(i);
 512           phi = fPhi2->At(j);
 513         }
 514       }
 515
 516       // TPC-EMCAL grid
 517       //-------------------------------------
 518       Int_t ii = 0;
 519       if(i==0) ii = 0;
 520       if(i>0 && i<(fEtaBinInTPCAcc-fEtaBinInEMCalAcc)/2) ii = i;
 521       if(i>=(fEtaBinInTPCAcc+fEtaBinInEMCalAcc)/2 && i<fNeta+1) ii = i-fEtaBinInEMCalAcc;
 522
 523       if(fGrid==1) {
 524         if(j<(fNbinPhi+1) && j*(fNeta+1)+i == index) {
 525           eta = fEta2->At(i);
 526           phi = fPhi2->At(j);
 527         }
 528
 529         if((j>=(fNbinPhi+1) && j<(fNbinPhi+1+fPhiBinInEMCalAcc)) &&
 530            ((fNbinPhi+1)*(fNeta+1) + (j-fNbinPhi-1)*(fEtaBinInTPCAcc-fEtaBinInEMCalAcc) + ii)== index ) {
 531           if(ii==0) {Int_t ind = 0; eta = fEta2->At(ind);}
 532           else eta = fEta2->At(i);
 533           phi = fPhi2->At(j);
 534         }
 535
 536         if(j>=(fNbinPhi+1+fPhiBinInEMCalAcc) && ((fNbinPhi+1)*(fNeta+1)+fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))+(j-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)+i == index)) {
 537           eta = fEta2->At(i);
 538           phi = fPhi2->At(j);
 539         }
 540       }
 541     }
 542   }
 543 }
 544
 545
 546
 547
 548
 549
 550
 551
 552
 553
 554
 555