1 /**************************************************************************
2 * Copyright(c) 2001-2002, 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 /* $Id: AliJetGrid.cxx,v 1.0 07/09/2006 */
18 //=========================================================================
19 // *** 07 September 2006
20 // This class allows to fill a really general grid in (eta,phi)
21 // Two types of grid can be setted :
22 // if fGrid==0 -> the complete acceptance of a rectangular detector acceptance is filled
23 // if fGrid==1 -> the previous grid - an other rectangular detector acceptance is filled
24 // A (eta,phi(rad)) position can be extracted from an index value without looping over all entries
25 // An index position can be extracted from a (eta,phi(rad)) position without looping over all entries
28 // > AliJetGrid *grid = new AliJetGrid((NbinPhi-1),(NbinEta-1),PhiMin,PhiMax,EtaMin,EtaMax);
29 // > grid->SetGridType(...); // 0 or 1 for the moment
30 // > grid->SetMatrixIndexes();
31 // > grid->SetIndexIJ();
33 // Author : magali.estienne@subatech.in2p3.fr
34 //=========================================================================
37 #include <Riostream.h>
44 #include "AliJetGrid.h"
48 //__________________________________________________________
49 AliJetGrid::AliJetGrid():
74 // Default constructor
77 //__________________________________________________________
78 AliJetGrid::AliJetGrid(Int_t nphi,Int_t neta,Double_t phiMin,Double_t phiMax,Double_t etaMin,Double_t etaMax):
103 // Standard constructor
104 fPhi = new TArrayD(fNphi+1);
105 fEta = new TArrayD(fNeta+1);
106 fIndexI = new TArrayI((fNeta+1)*(fNphi+1)+1);
107 fIndexJ = new TArrayI((fNeta+1)*(fNphi+1)+1);
109 for(Int_t i=0; i<fNphi+1; i++) {
110 if(fNphi!=0) (*fPhi)[i] = (phiMax-phiMin)/fNphi*i+phiMin;
111 else (*fPhi)[i] = phiMin+(phiMax-phiMin)/2;
113 for(Int_t i=0; i<fNeta+1; i++) {
114 if(fNeta!=0) (*fEta)[i] = (etaMax-etaMin)/fNeta*i+etaMin;
115 else (*fEta)[i] = etaMin+(etaMax-etaMin)/2;
119 for(Int_t i=0; i<(fNphi+1); i++) cout << (*fPhi)[i] << endl;
120 for(Int_t i=0; i<(fNeta+1); i++) cout << (*fEta)[i] << endl;
123 fIndex = new TMatrixD(fNphi+1,fNeta+1);
127 //__________________________________________________________
128 AliJetGrid::AliJetGrid(const AliJetGrid& grid) :
136 fIndexI(grid.fIndexI),
137 fIndexJ(grid.fIndexJ),
138 fPhiMin(grid.fPhiMin),
139 fPhiMax(grid.fPhiMax),
140 fEtaMin(grid.fEtaMin),
141 fEtaMax(grid.fEtaMax),
142 fEtaBinInTPCAcc(grid.fEtaBinInTPCAcc),
143 fPhiBinInTPCAcc(grid.fPhiBinInTPCAcc),
144 fEtaBinInEMCalAcc(grid.fEtaBinInEMCalAcc),
145 fPhiBinInEMCalAcc(grid.fPhiBinInEMCalAcc),
146 fNbinEta(grid.fNbinEta),
147 fNbinPhi(grid.fNbinPhi),
148 fMaxPhi(grid.fMaxPhi),
149 fMinPhi(grid.fMinPhi),
150 fMaxEta(grid.fMaxEta),
151 fMinEta(grid.fMinEta),
157 fPhi = new TArrayD(fNphi+1);
158 for(Int_t i=0; i<fNphi+1; i++) (*fPhi)[i] = grid.fPhi->At(i);
159 fEta = new TArrayD(fNeta+1);
160 for(Int_t i=0; i<fNeta+1; i++) (*fEta)[i] = grid.fEta->At(i);
162 fIndex = new TMatrixD(fNphi+1,fNeta+1);
163 for(Int_t i=0; i<fNphi+1; i++) {
164 for(Int_t j=0; j<fNeta+1; j++) (*fIndex)(i,j)=(*grid.fIndex)(i,j);
169 AliJetGrid& AliJetGrid::operator=(const AliJetGrid& other)
172 if (this != &other) {
179 fIndexI = other.fIndexI;
180 fIndexJ = other.fIndexJ;
181 fPhiMin = other.fPhiMin;
182 fPhiMax = other.fPhiMax;
183 fEtaMin = other.fEtaMin;
184 fEtaMax = other.fEtaMax;
185 fEtaBinInTPCAcc = other.fEtaBinInTPCAcc;
186 fPhiBinInTPCAcc = other.fPhiBinInTPCAcc;
187 fEtaBinInEMCalAcc = other.fEtaBinInEMCalAcc;
188 fPhiBinInEMCalAcc = other.fPhiBinInEMCalAcc;
189 fNbinEta = other.fNbinEta;
190 fNbinPhi = other.fNbinPhi;
191 fMaxPhi = other.fMaxPhi;
192 fMinPhi = other.fMinPhi;
193 fMaxEta = other.fMaxEta;
194 fMinEta = other.fMinEta;
195 fDebug = other.fDebug;
196 fPhi = new TArrayD(fNphi+1);
197 for(Int_t i=0; i<fNphi+1; i++) (*fPhi)[i] = other.fPhi->At(i);
198 fEta = new TArrayD(fNeta+1);
199 for(Int_t i=0; i<fNeta+1; i++) (*fEta)[i] = other.fEta->At(i);
201 fIndex = new TMatrixD(fNphi+1,fNeta+1);
202 for(Int_t i=0; i<fNphi+1; i++) {
203 for(Int_t j=0; j<fNeta+1; j++) (*fIndex)(i,j)=(*other.fIndex)(i,j);
210 //__________________________________________________________
211 AliJetGrid::~AliJetGrid() {
221 //__________________________________________________________
222 void AliJetGrid::InitParams(Double_t phiMinCal,Double_t phiMaxCal,Double_t etaMinCal,Double_t etaMaxCal)
224 // To set initial parameters
226 fPhiMin = phiMinCal; // rad
227 fPhiMax = phiMaxCal; // rad
230 fNbinPhi = static_cast<int>(fPhiMin/(2*TMath::Pi()/fNphi));
232 // Define some binning numbers
234 for(Int_t i=0; i<fNphi+1; i++) fPhiBinInTPCAcc++;
235 fPhiBinInEMCalAcc = 0;
237 for(Int_t i=0; i<fNeta+1; i++) fEtaBinInTPCAcc++;
238 fEtaBinInEMCalAcc = 0;
242 for(Int_t i=0; i<fNphi+1; i++) {
244 if(fPhi->At(i) >= fPhiMin &&
245 fPhi->At(i) <= fPhiMax)
248 for(Int_t i=0; i<fNeta+1; i++) {
250 if((fEta->At(i) >= fEtaMin) &&
251 (fEta->At(i) <= fEtaMax))
258 //__________________________________________________________
259 TArrayD* AliJetGrid::GetArrayEta()
261 // Returns an array with the eta points
267 //__________________________________________________________
268 TArrayD* AliJetGrid::GetArrayPhi()
270 // Returns an array with the phi points
276 //__________________________________________________________
277 TMatrixD* AliJetGrid::GetIndexObject()
279 // Returns a pointer to the matrix
285 //__________________________________________________________
286 void AliJetGrid::GetAccParam(Int_t &nphi, Int_t &neta, Float_t &minphi, Float_t &maxphi,
287 Float_t &mineta, Float_t &maxeta) const
289 // Returns EMCAL acceptance initially setted
299 //__________________________________________________________
300 void AliJetGrid::GetBinParam(Int_t &phibintpc, Int_t &etabintpc,
301 Int_t &phibinemc, Int_t &etabinemc, Int_t &nbinphi) const
303 // Returns number of bins in TPC and EMCAL geometry
305 etabintpc = fEtaBinInTPCAcc;
306 phibintpc = fPhiBinInTPCAcc;
307 etabinemc = fEtaBinInEMCalAcc;
308 phibinemc = fPhiBinInEMCalAcc;
312 //__________________________________________________________
313 Int_t AliJetGrid::GetIndexFromEtaPhi(Double_t phi,Double_t eta) const
315 // Tells the index value of a corresponding (eta,phi) real position
316 // Loop over all entries -> takes time.
317 // Used one time at the begining to fill the grids
319 /* this is how bins are numbered
321 in all TPC or in TPC - EMCAL
324 ... ... . 10 | | | 11
325 ---+---+--- ---------------
326 ^ 6 | 7 | 8 or 8 | | | 9
327 | ---+---+--- ---------------
328 3 | 4 | 5 4 | 5 | 6 | 7
329 phi ---+---+--- ---------------
330 0 | 1 | 2 0 | 1 | 2 | 3
336 Int_t etaBin=0,phiBin=0,absID=0;
337 Int_t etaBin2=0,etaBin3=0;
339 // Fill all the grid in eta/phi (all TPC acceptance)
340 //-----------------------------------------------------
342 if(eta <= fEta->At(0)) {
344 } else if(eta >= fEta->At(fNeta)) {
347 for(Int_t i=0; i<fNeta+1; i++) {
348 if(eta < fEta->At(i)) {
354 if(phi <= fPhi->At(0)) {
356 } else if(phi >= fPhi->At(fNphi)) {
359 for(Int_t i=0; i<fNphi+1; i++) {
360 if(phi < fPhi->At(i)) {
367 // Calculate absolute id
368 absID = phiBin*(fNeta+1) + etaBin;
372 // Fill the grid but do not count id in EMCal acceptance
373 //------------------------------------------------------
374 if((eta >= fEtaMin && eta <= fEtaMax) &&
375 (phi >= fPhiMin && phi <= fPhiMax)){
376 etaBin = etaBin2 = etaBin3 = -1;
382 if(eta <= fEta->At(0)) {
384 } else if(eta >= fEta->At(fNeta)) {
387 for(Int_t i=0; i<fNeta+1; i++) {
388 if(eta < fEta->At(i)) {
396 if((phi>=fPhiMin) && (phi<=fPhiMax)) {
397 if(eta <= fEta->At(0)) {
399 } else if(eta >= fEta->At(fNeta)) {
400 etaBin2 = fNeta-fEtaBinInEMCalAcc;
402 for(Int_t i=0; i<fNeta+1; i++) {
403 if(eta < fEta->At(i) && ((fEta->At(0)<eta && eta<fEtaMin) ||
404 (fEtaMax<eta && eta<fEta->At(fNeta)))) {
405 // cout << "i : " << i << endl;
408 else etaBin2 = i-1-fEtaBinInEMCalAcc;
415 if(eta <= fEta->At(0)) {
417 } else if(eta >= fEta->At(fNeta)) {
420 for(Int_t i=0; i<fNeta+1; i++) {
421 if(eta < fEta->At(i)) {
430 if(phi <= fPhi->At(0)) {
432 } else if(phi >= fPhi->At(fNphi)) {
435 for(Int_t i=0; i<fNphi+1; i++) {
436 if(phi < fPhi->At(i)) {
444 //-------------------------------------------------------
446 absID = phiBin*(fNeta+1) + etaBin;
447 if(phi>=fPhiMin && phi<=fPhiMax) {
448 if(eta >= fEtaMin && eta <= fEtaMax) absID = -1;
450 absID = (fNbinPhi+1)*(fNeta+1)
451 + (phiBin-fNbinPhi-1)*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
456 absID = (fNbinPhi+1)*(fNeta+1)
457 + fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
458 + (phiBin-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)
467 //__________________________________________________________
468 void AliJetGrid::GetEtaPhiFromIndex(Int_t index, Float_t &eta, Float_t &phi)
470 // Get (eta,phi) position for a given index BUT loop over all entries (takes time)
472 for(Int_t j=0; j<fNphi+1; j++) {
473 for(Int_t i=0; i<fNeta+1; i++) {
476 //-------------------------------------
478 if(j*(fNeta+1)+i == index) {
485 //-------------------------------------
488 if(i>0 && i<(fEtaBinInTPCAcc-fEtaBinInEMCalAcc)/2) ii = i;
489 if(i>=(fEtaBinInTPCAcc+fEtaBinInEMCalAcc)/2 && i<fNeta+1) ii = i-fEtaBinInEMCalAcc;
492 if(j<(fNbinPhi+1) && j*(fNeta+1)+i == index) {
497 if((j>=(fNbinPhi+1) && j<(fNbinPhi+1+fPhiBinInEMCalAcc)) &&
498 ((fNbinPhi+1)*(fNeta+1) + (j-fNbinPhi-1)*(fEtaBinInTPCAcc-fEtaBinInEMCalAcc) + ii)== index ) {
499 if(ii==0) {Int_t ind = 0; eta = fEta->At(ind);}
500 else eta = fEta->At(i);
504 if(j>=(fNbinPhi+1+fPhiBinInEMCalAcc) &&
505 ((fNbinPhi+1)*(fNeta+1)+fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
506 +(j-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)+i == index)) {
515 //__________________________________________________________
516 Int_t AliJetGrid::GetIndex(Double_t phi, Double_t eta)
518 // Get index value for a (eta,phi) position - Direct value
520 Int_t ieta = GetIndexJFromEta(eta);
521 Int_t iphi = GetIndexIFromPhi(phi);
523 Int_t index = GetMatrixIndex(iphi,ieta);
526 cout << "(phi,eta) : " << phi << ", " << eta << endl;
527 cout << "index : " << index << endl;
533 //__________________________________________________________
534 Int_t AliJetGrid::GetIndexJFromEta(Double_t eta)
540 Double_t temp = (eta+fMaxEta)/(fMaxEta-fMinEta)*fNeta;
543 cout << "eta : " << eta << endl;
544 cout << "fMaxEta : " << fMaxEta << endl;
545 cout << "fMinEta : " << fMinEta << endl;
546 cout << "fNeta : " << fNeta << endl;
547 cout << "temp eta before cast : " << temp << endl;
549 idEta = static_cast<Int_t>(temp+0.5);
550 if(fDebug>20) cout << "temp eta after cast : " << idEta << endl;
553 //__________________________________________________________
554 Int_t AliJetGrid::GetIndexIFromPhi(Double_t phi)
561 if(fMinPhi==0) temp = phi/(fMaxPhi-fMinPhi)*fNphi;
562 else temp = (phi-fMinPhi)/(fMaxPhi-fMinPhi)*fNphi;
566 cout << "phi : " << phi << endl;
567 cout << "fMaxPhi : " << fMaxPhi << endl;
568 cout << "fMinPhi : " << fMinPhi << endl;
569 cout << "fNphi : " << fNphi << endl;
570 cout << "temp phi before cast : " << temp << endl;
572 idPhi = static_cast<Int_t>(temp+0.5);
573 if(fDebug>20) cout << "temp phi after cast : " << idPhi << endl;
579 //__________________________________________________________
580 void AliJetGrid::SetMatrixIndex(Int_t i,Double_t par)
582 // Allows to set parameters using only one index (if fGrid==0) !!
585 Int_t iphi = (Int_t)i/fNeta;
586 Int_t ieta = i-iphi*fNeta;
587 SetMatrixIndex(iphi,ieta,par);
592 //__________________________________________________________
593 void AliJetGrid::SetMatrixIndexes()
595 // Fill the final matrix object with the corresponding index in eta/phi
597 for(Int_t i=0; i<fNphi+1; i++){
598 for(Int_t j=0; j<fNeta+1; j++){
599 (*fIndex)(i,j) = GetIndexFromEtaPhi(fPhi->At(i),fEta->At(j))+1;
601 cout << "(*fIndex)(" << i << "," << j << ") : " << (*fIndex)(i,j) <<
602 ", phi : " << fPhi->At(i) << ", eta : " << fEta->At(j) << endl;
606 printf("##########################################################\n");
607 printf("TMatrix object filled !\n");
608 printf("Size of the object : phi x eta = (fNphi+1) x (fNeta+1) = %d\n",(fNphi+1)*(fNeta+1));
609 printf("##########################################################\n");
612 //__________________________________________________________
613 void AliJetGrid::SetIndexIJ()
615 // Set the grid index
617 for(Int_t i=0; i<fNphi+1; i++){
618 for(Int_t j=0; j<fNeta+1; j++){
619 Int_t id = static_cast<Int_t>((*fIndex)(i,j));
629 printf("##########################################################\n");
630 printf(" In SetIndexIJ - Grid indexes setted !\n");
631 printf("##########################################################\n");
634 //__________________________________________________________
635 void AliJetGrid::GetIJFromIndex(Int_t index, Int_t& i, Int_t& j) const
637 // Returns i position id of eta and j position id of phi for a given grid index
638 i = (*fIndexI)[index];
639 j = (*fIndexJ)[index];
642 //__________________________________________________________
643 void AliJetGrid::GetEtaPhiFromIndex2(Int_t index, Float_t &phi, Float_t &eta)
645 // Returns eta, phi values for a given grid index
647 phi = fPhi->At((*fIndexI)[index]);
648 eta = fEta->At((*fIndexJ)[index]);
651 //__________________________________________________________
652 Int_t AliJetGrid::GetNEntries()
654 // Returns the number of entries of the grid
657 cout << "fMaxPhi : " << fMaxPhi << endl;
658 cout << "fMaxEta : " << fMaxEta << endl;
661 Int_t indexNum = GetIndex(fMaxPhi,fMaxEta);
662 if(fDebug>20) cout << "indexNum : " << indexNum << endl;
667 //__________________________________________________________
668 Int_t AliJetGrid::GetNEntries2()
670 // Returns the number of entries of the grid
672 Int_t indexNum = GetIndex(fMaxPhi-1.,fMaxEta-0.5);
673 if(fDebug>20) cout << "indexNum : " << indexNum << endl;