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@ires.in2p3.fr
34 //=========================================================================
37 #include <Riostream.h>
44 #include "AliJetGrid.h"
49 //__________________________________________________________
50 AliJetGrid::AliJetGrid():
75 // Default constructor
78 //__________________________________________________________
79 AliJetGrid::AliJetGrid(Int_t nphi,Int_t neta,Double_t phiMin,Double_t phiMax,Double_t etaMin,Double_t etaMax):
104 // Standard constructor
105 fPhi = new TArrayD(fNphi+1);
106 fEta = new TArrayD(fNeta+1);
107 fIndexI = new TArrayI((fNeta+1)*(fNphi+1)+1);
108 fIndexJ = new TArrayI((fNeta+1)*(fNphi+1)+1);
110 for(Int_t i=0; i<fNphi+1; i++) (*fPhi)[i] = (phiMax-phiMin)/fNphi*i+phiMin;
111 for(Int_t i=0; i<fNeta+1; i++) (*fEta)[i] = (etaMax-etaMin)/fNeta*i+etaMin;
114 for(Int_t i=0; i<(fNphi+1); i++) cout << (*fPhi)[i] << endl;
115 for(Int_t i=0; i<(fNeta+1); i++) cout << (*fEta)[i] << endl;
118 fIndex = new TMatrixD(fNphi+1,fNeta+1);
122 //__________________________________________________________
123 AliJetGrid::AliJetGrid(const AliJetGrid& grid):TNamed(grid) {
129 fPhiMin = grid.fPhiMin;
130 fPhiMax = grid.fPhiMax;
131 fEtaMin = grid.fEtaMin;
132 fEtaMax = grid.fEtaMax;
133 fEtaBinInTPCAcc = grid.fEtaBinInTPCAcc;
134 fPhiBinInTPCAcc = grid.fPhiBinInTPCAcc;
135 fEtaBinInEMCalAcc = grid.fEtaBinInEMCalAcc;
136 fPhiBinInEMCalAcc = grid.fPhiBinInEMCalAcc;
137 fNbinPhi = grid.fNbinPhi;
138 fMaxPhi = grid.fMaxPhi;
139 fMinPhi = grid.fMinPhi;
140 fMaxEta = grid.fMaxEta;
141 fMinEta = grid.fMinEta;
143 fPhi = new TArrayD(fNphi+1);
144 for(Int_t i=0; i<fNphi+1; i++) (*fPhi)[i] = grid.fPhi->At(i);
145 fEta = new TArrayD(fNeta+1);
146 for(Int_t i=0; i<fNeta+1; i++) (*fEta)[i] = grid.fEta->At(i);
148 fIndex = new TMatrixD(fNphi+1,fNeta+1);
149 for(Int_t i=0; i<fNphi+1; i++) {
150 for(Int_t j=0; j<fNeta+1; j++) (*fIndex)(i,j)=(*grid.fIndex)(i,j);
154 //__________________________________________________________
155 AliJetGrid::~AliJetGrid() {
165 //__________________________________________________________
166 void AliJetGrid::InitParams(Double_t phiMinCal,Double_t phiMaxCal,Double_t etaMinCal,Double_t etaMaxCal)
167 { // To set initial parameters
169 fPhiMin = phiMinCal; // rad
170 fPhiMax = phiMaxCal; // rad
173 fNbinPhi = static_cast<int>(fPhiMin/(2*TMath::Pi()/fNphi));
175 // Define some binning numbers
177 for(Int_t i=0; i<fNphi+1; i++) fPhiBinInTPCAcc++;
178 fPhiBinInEMCalAcc = 0;
180 for(Int_t i=0; i<fNeta+1; i++) fEtaBinInTPCAcc++;
181 fEtaBinInEMCalAcc = 0;
185 for(Int_t i=0; i<fNphi+1; i++) {
187 if(fPhi->At(i) >= fPhiMin &&
188 fPhi->At(i) <= fPhiMax)
191 for(Int_t i=0; i<fNeta+1; i++) {
193 if((fEta->At(i) >= fEtaMin) &&
194 (fEta->At(i) <= fEtaMax))
201 //__________________________________________________________
202 TArrayD* AliJetGrid::GetArrayEta()
203 { // Returns an array with the eta points
209 //__________________________________________________________
210 TArrayD* AliJetGrid::GetArrayPhi()
211 { // Returns an array with the phi points
217 //__________________________________________________________
218 TMatrixD* AliJetGrid::GetIndexObject()
219 { // Returns a pointer to the matrix
225 //__________________________________________________________
226 void AliJetGrid::GetAccParam(Int_t &nphi, Int_t &neta, Float_t &minphi, Float_t &maxphi,
227 Float_t &mineta, Float_t &maxeta)
228 { // Returns EMCAL acceptance initially setted
238 //__________________________________________________________
239 void AliJetGrid::GetBinParam(Int_t &phibintpc, Int_t &etabintpc,
240 Int_t &phibinemc, Int_t &etabinemc, Int_t &nbinphi)
241 { // Returns number of bins in TPC and EMCAL geometry
243 etabintpc = fEtaBinInTPCAcc;
244 phibintpc = fPhiBinInTPCAcc;
245 etabinemc = fEtaBinInEMCalAcc;
246 phibinemc = fPhiBinInEMCalAcc;
250 //__________________________________________________________
251 Int_t AliJetGrid::GetIndexFromEtaPhi(Double_t phi,Double_t eta) const
252 { // Tells the index value of a corresponding (eta,phi) real position
253 // Loop over all entries -> takes time.
254 // Used one time at the begining to fill the grids
256 /* this is how bins are numbered
258 in all TPC or in TPC - EMCAL
261 ... ... . 10 | | | 11
262 ---+---+--- ---------------
263 ^ 6 | 7 | 8 or 8 | | | 9
264 | ---+---+--- ---------------
265 3 | 4 | 5 4 | 5 | 6 | 7
266 phi ---+---+--- ---------------
267 0 | 1 | 2 0 | 1 | 2 | 3
273 Int_t etaBin=0,phiBin=0,absID=0;
274 Int_t etaBin2=0,etaBin3=0;
276 // Fill all the grid in eta/phi (all TPC acceptance)
277 //-----------------------------------------------------
279 if(eta <= fEta->At(0)) {
281 } else if(eta >= fEta->At(fNeta)) {
284 for(Int_t i=0; i<fNeta+1; i++) {
285 if(eta < fEta->At(i)) {
291 if(phi <= fPhi->At(0)) {
293 } else if(phi >= fPhi->At(fNphi)) {
296 for(Int_t i=0; i<fNphi+1; i++) {
297 if(phi < fPhi->At(i)) {
304 // Calculate absolute id
305 absID = phiBin*(fNeta+1) + etaBin;
309 // Fill the grid but do not count id in EMCal acceptance
310 //------------------------------------------------------
311 if((eta >= fEtaMin && eta <= fEtaMax) &&
312 (phi >= fPhiMin && phi <= fPhiMax)){
313 etaBin = etaBin2 = etaBin3 = -1;
319 if(eta <= fEta->At(0)) {
321 } else if(eta >= fEta->At(fNeta)) {
324 for(Int_t i=0; i<fNeta+1; i++) {
325 if(eta < fEta->At(i)) {
333 if((phi>=fPhiMin) && (phi<=fPhiMax)) {
334 if(eta <= fEta->At(0)) {
336 } else if(eta >= fEta->At(fNeta)) {
337 etaBin2 = fNeta-fEtaBinInEMCalAcc;
339 for(Int_t i=0; i<fNeta+1; i++) {
340 if(eta < fEta->At(i) && ((fEta->At(0)<eta && eta<fEtaMin) ||
341 (fEtaMax<eta && eta<fEta->At(fNeta)))) {
342 // cout << "i : " << i << endl;
345 else etaBin2 = i-1-fEtaBinInEMCalAcc;
352 if(eta <= fEta->At(0)) {
354 } else if(eta >= fEta->At(fNeta)) {
357 for(Int_t i=0; i<fNeta+1; i++) {
358 if(eta < fEta->At(i)) {
367 if(phi <= fPhi->At(0)) {
369 } else if(phi >= fPhi->At(fNphi)) {
372 for(Int_t i=0; i<fNphi+1; i++) {
373 if(phi < fPhi->At(i)) {
381 //-------------------------------------------------------
383 absID = phiBin*(fNeta+1) + etaBin;
384 if(phi>=fPhiMin && phi<=fPhiMax) {
385 if(eta >= fEtaMin && eta <= fEtaMax) absID = -1;
387 absID = (fNbinPhi+1)*(fNeta+1)
388 + (phiBin-fNbinPhi-1)*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
393 absID = (fNbinPhi+1)*(fNeta+1)
394 + fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
395 + (phiBin-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)
404 //__________________________________________________________
405 void AliJetGrid::GetEtaPhiFromIndex(Int_t index, Float_t &eta, Float_t &phi)
406 { // Get (eta,phi) position for a given index BUT loop over all entries (takes time)
408 for(Int_t j=0; j<fNphi+1; j++) {
409 for(Int_t i=0; i<fNeta+1; i++) {
412 //-------------------------------------
414 if(j*(fNeta+1)+i == index) {
421 //-------------------------------------
424 if(i>0 && i<(fEtaBinInTPCAcc-fEtaBinInEMCalAcc)/2) ii = i;
425 if(i>=(fEtaBinInTPCAcc+fEtaBinInEMCalAcc)/2 && i<fNeta+1) ii = i-fEtaBinInEMCalAcc;
428 if(j<(fNbinPhi+1) && j*(fNeta+1)+i == index) {
433 if((j>=(fNbinPhi+1) && j<(fNbinPhi+1+fPhiBinInEMCalAcc)) &&
434 ((fNbinPhi+1)*(fNeta+1) + (j-fNbinPhi-1)*(fEtaBinInTPCAcc-fEtaBinInEMCalAcc) + ii)== index ) {
435 if(ii==0) {Int_t ind = 0; eta = fEta->At(ind);}
436 else eta = fEta->At(i);
440 if(j>=(fNbinPhi+1+fPhiBinInEMCalAcc) &&
441 ((fNbinPhi+1)*(fNeta+1)+fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
442 +(j-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)+i == index)) {
451 //__________________________________________________________
452 Int_t AliJetGrid::GetIndex(Double_t phi, Double_t eta)
453 { // Get index value for a (eta,phi) position - Direct value
455 Int_t ieta = GetIndexJFromEta(eta);
456 Int_t iphi = GetIndexIFromPhi(phi);
458 Int_t index = GetMatrixIndex(iphi,ieta);
461 cout << "(phi,eta) : " << phi << ", " << eta << endl;
462 cout << "index : " << index << endl;
468 //__________________________________________________________
469 Int_t AliJetGrid::GetIndexJFromEta(Double_t eta)
474 Double_t temp = (eta+fMaxEta)/(fMaxEta-fMinEta)*fNeta;
477 cout << "eta : " << eta << endl;
478 cout << "fMaxEta : " << fMaxEta << endl;
479 cout << "fMinEta : " << fMinEta << endl;
480 cout << "fNeta : " << fNeta << endl;
481 cout << "temp eta before cast : " << temp << endl;
483 idEta = static_cast<Int_t>(temp+0.5);
484 if(fDebug>20) cout << "temp eta after cast : " << idEta << endl;
487 //__________________________________________________________
488 Int_t AliJetGrid::GetIndexIFromPhi(Double_t phi)
494 if(fMinPhi==0) temp = phi/(fMaxPhi-fMinPhi)*fNphi;
495 else temp = (phi-fMinPhi)/(fMaxPhi-fMinPhi)*fNphi;
499 cout << "phi : " << phi << endl;
500 cout << "fMaxPhi : " << fMaxPhi << endl;
501 cout << "fMinPhi : " << fMinPhi << endl;
502 cout << "fNphi : " << fNphi << endl;
503 cout << "temp phi before cast : " << temp << endl;
505 idPhi = static_cast<Int_t>(temp+0.5);
506 if(fDebug>20) cout << "temp phi after cast : " << idPhi << endl;
512 //__________________________________________________________
513 void AliJetGrid::SetMatrixIndex(Int_t i,Double_t par)
514 { // Allows to set parameters using only one index (if fGrid==0) !!
517 Int_t iphi = (Int_t)i/fNeta;
518 Int_t ieta = i-iphi*fNeta;
519 SetMatrixIndex(iphi,ieta,par);
524 //__________________________________________________________
525 void AliJetGrid::SetMatrixIndexes()
526 { // Fill the final matrix object with the corresponding index in eta/phi
528 for(Int_t i=0; i<fNphi+1; i++){
529 for(Int_t j=0; j<fNeta+1; j++){
530 (*fIndex)(i,j) = GetIndexFromEtaPhi(fPhi->At(i),fEta->At(j))+1;
532 cout << "(*fIndex)(" << i << "," << j << ") : " << (*fIndex)(i,j) <<
533 ", phi : " << fPhi->At(i) << ", eta : " << fEta->At(j) << endl;
537 printf("##########################################################\n");
538 printf("TMatrix object filled !\n");
539 printf("Size of the object : phi x eta = (fNphi+1) x (fNeta+1) = %d\n",(fNphi+1)*(fNeta+1));
540 printf("##########################################################\n");
543 //__________________________________________________________
544 void AliJetGrid::SetIndexIJ()
547 for(Int_t i=0; i<fNphi+1; i++){
548 for(Int_t j=0; j<fNeta+1; j++){
549 Int_t id = static_cast<Int_t>((*fIndex)(i,j));
559 printf("##########################################################\n");
560 printf(" In SetIndexIJ - Grid indexes setted !\n");
561 printf("##########################################################\n");
564 //__________________________________________________________
565 void AliJetGrid::GetIJFromIndex(Int_t index, Int_t i, Int_t j)
566 { // returns i position id of eta and j position id of phi for a given grid index
567 i = (*fIndexI)[index];
568 j = (*fIndexJ)[index];
571 //__________________________________________________________
572 void AliJetGrid::GetEtaPhiFromIndex2(Int_t index, Float_t &phi, Float_t &eta)
573 { // returns eta, phi values for a given grid index
575 phi = fPhi->At((*fIndexI)[index]);
576 eta = fEta->At((*fIndexJ)[index]);
579 //__________________________________________________________
580 Int_t AliJetGrid::GetNEntries()
581 { // Returns the number of entries of the grid
584 cout << "fMaxPhi : " << fMaxPhi << endl;
585 cout << "fMaxEta : " << fMaxEta << endl;
588 Int_t indexNum = GetIndex(fMaxPhi,fMaxEta);
589 if(fDebug>20) cout << "indexNum : " << indexNum << endl;
594 //__________________________________________________________
595 Int_t AliJetGrid::GetNEntries2()
596 { // Returns the number of entries of the grid
598 Int_t indexNum = GetIndex(fMaxPhi-1.,fMaxEta-0.5);
599 if(fDebug>20) cout << "indexNum : " << indexNum << endl;