]> git.uio.no Git - u/mrichter/AliRoot.git/blame - JETAN/AliJetGrid.cxx
Compact plots in a single canvas
[u/mrichter/AliRoot.git] / JETAN / AliJetGrid.cxx
CommitLineData
4a01bb2c 1/**************************************************************************
2 * Copyright(c) 2001-2002, 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/* $Id: AliJetGrid.cxx,v 1.0 07/09/2006 */
17
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
26//
27// How to run it :
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();
32//
33// Author : magali.estienne@ires.in2p3.fr
34//=========================================================================
35
36// Standard headers
37#include <Riostream.h>
38// Root headers
39#include <TMatrixD.h>
40#include <TArrayD.h>
41#include <TArrayI.h>
42// AliRoot headers
43#include "AliJetGrid.h"
44
45ClassImp(AliJetGrid)
46
47//__________________________________________________________
48AliJetGrid::AliJetGrid() {
49
50 // Default constructor
51
52 fNphi = 0;
53 fNeta = 0;
54 fPhi = 0;
55 fEta = 0;
56 fIndex = 0x0;
57 fPhiMin = 0.; // acceptance removed inside
58 fPhiMax = 0.; // acceptance removed inside
59 fEtaMin = 0.; // acceptance removed inside
60 fEtaMax = 0.; // acceptance removed inside
61 fMinPhi = 0.; // total acceptance
62 fMaxPhi = 0.; // total acceptance
63 fMinEta = 0.; // total acceptance
64 fMaxEta = 0.; // total acceptance
65 fEtaBinInTPCAcc = 0;
66 fPhiBinInTPCAcc = 0;
67 fEtaBinInEMCalAcc = 0;
68 fPhiBinInEMCalAcc = 0;
69 fNbinPhi = 0;
70
71 fDebug = 1;
72}
73
74//__________________________________________________________
75AliJetGrid::AliJetGrid(Int_t nphi,Int_t neta,Double_t phiMin,Double_t phiMax,Double_t etaMin,Double_t etaMax) {
76
77 // Standard constructor
78
79 fNphi = nphi;
80 fNeta = neta;
81 fPhiMin = 0.; // rad - acceptance removed inside
82 fPhiMax = 0.; // rad - acceptance removed inside
83 fEtaMin = 0.; // acceptance removed inside
84 fEtaMax = 0.; // acceptance removed inside
85 fNbinPhi = 0;
86
87 fMaxPhi = phiMax; // total acceptance
88 fMinPhi = phiMin; // total acceptance
89 fMaxEta = etaMax; // total acceptance
90 fMinEta = etaMin; // total acceptance
91
92 fPhi = new TArrayD(fNphi+1);
93 fEta = new TArrayD(fNeta+1);
94 fIndexI = new TArrayI((fNeta+1)*(fNphi+1)+1);
95 fIndexJ = new TArrayI((fNeta+1)*(fNphi+1)+1);
96
97 for(Int_t i=0; i<fNphi+1; i++) (*fPhi)[i] = (phiMax-phiMin)/fNphi*i+phiMin;
98 for(Int_t i=0; i<fNeta+1; i++) (*fEta)[i] = (etaMax-etaMin)/fNeta*i+etaMin;
99
100 if(fDebug > 3){
101 for(Int_t i=0; i<(fNphi+1); i++) cout << (*fPhi)[i] << endl;
102 for(Int_t i=0; i<(fNeta+1); i++) cout << (*fEta)[i] << endl;
103 }
104
105 fIndex = new TMatrixD(fNphi+1,fNeta+1);
106
107}
108
109//__________________________________________________________
110AliJetGrid::AliJetGrid(const AliJetGrid& grid):TNamed(grid) {
111
112 // Copy constructor
113
114 fNphi = grid.fNphi;
115 fNeta = grid.fNeta;
116 fPhiMin = grid.fPhiMin;
117 fPhiMax = grid.fPhiMax;
118 fEtaMin = grid.fEtaMin;
119 fEtaMax = grid.fEtaMax;
120 fEtaBinInTPCAcc = grid.fEtaBinInTPCAcc;
121 fPhiBinInTPCAcc = grid.fPhiBinInTPCAcc;
122 fEtaBinInEMCalAcc = grid.fEtaBinInEMCalAcc;
123 fPhiBinInEMCalAcc = grid.fPhiBinInEMCalAcc;
124 fNbinPhi = grid.fNbinPhi;
125 fMaxPhi = grid.fMaxPhi;
126 fMinPhi = grid.fMinPhi;
127 fMaxEta = grid.fMaxEta;
128 fMinEta = grid.fMinEta;
129
130 fPhi = new TArrayD(fNphi+1);
131 for(Int_t i=0; i<fNphi+1; i++) (*fPhi)[i] = grid.fPhi->At(i);
132 fEta = new TArrayD(fNeta+1);
133 for(Int_t i=0; i<fNeta+1; i++) (*fEta)[i] = grid.fEta->At(i);
134
135 fIndex = new TMatrixD(fNphi+1,fNeta+1);
136 for(Int_t i=0; i<fNphi+1; i++) {
137 for(Int_t j=0; j<fNeta+1; j++) (*fIndex)(i,j)=(*grid.fIndex)(i,j);
138 }
139}
140
141//__________________________________________________________
142AliJetGrid::~AliJetGrid() {
143
144 // Destructor
145 delete fPhi;
146 delete fEta;
147 delete fIndexI;
148 delete fIndexJ;
149 delete fIndex;
150}
151
152//__________________________________________________________
153void AliJetGrid::InitParams(Double_t phiMinCal,Double_t phiMaxCal,Double_t etaMinCal,Double_t etaMaxCal)
154{ // To set initial parameters
155
156 fPhiMin = phiMinCal; // rad
157 fPhiMax = phiMaxCal; // rad
158 fEtaMin = etaMinCal;
159 fEtaMax = etaMaxCal;
160 fNbinPhi = static_cast<int>(fPhiMin/(2*TMath::Pi()/fNphi));
161
162 // Define some binning numbers
163 if(fGrid==0){
164 for(Int_t i=0; i<fNphi+1; i++) fPhiBinInTPCAcc++;
165 fPhiBinInEMCalAcc = 0;
166
167 for(Int_t i=0; i<fNeta+1; i++) fEtaBinInTPCAcc++;
168 fEtaBinInEMCalAcc = 0;
169 }
170
171 if(fGrid==1){
172 for(Int_t i=0; i<fNphi+1; i++) {
173 fPhiBinInTPCAcc++;
174 if(fPhi->At(i) >= fPhiMin &&
175 fPhi->At(i) <= fPhiMax)
176 fPhiBinInEMCalAcc++;
177 }
178 for(Int_t i=0; i<fNeta+1; i++) {
179 fEtaBinInTPCAcc++;
180 if((fEta->At(i) >= fEtaMin) &&
181 (fEta->At(i) <= fEtaMax))
182 fEtaBinInEMCalAcc++;
183 }
184 }
185
186}
187
188//__________________________________________________________
189TArrayD* AliJetGrid::GetArrayEta()
190{ // Returns an array with the eta points
191
192 return fEta;
193
194}
195
196//__________________________________________________________
197TArrayD* AliJetGrid::GetArrayPhi()
198{ // Returns an array with the phi points
199
200 return fPhi;
201
202}
203
204//__________________________________________________________
205TMatrixD* AliJetGrid::GetIndexObject()
206{ // Returns a pointer to the matrix
207
208 return fIndex;
209
210}
211
212//__________________________________________________________
213void AliJetGrid::GetAccParam(Int_t &nphi, Int_t &neta, Float_t &minphi, Float_t &maxphi,
214 Float_t &mineta, Float_t &maxeta)
215{ // Returns EMCAL acceptance initially setted
216
217 nphi = fNphi;
218 neta = fNeta;
219 minphi = fPhiMin;
220 maxphi = fPhiMax;
221 mineta = fEtaMin;
222 maxeta = fEtaMax;
223}
224
225//__________________________________________________________
226void AliJetGrid::GetBinParam(Int_t &phibintpc, Int_t &etabintpc,
227 Int_t &phibinemc, Int_t &etabinemc, Int_t &nbinphi)
228{ // Returns number of bins in TPC and EMCAL geometry
229
230 etabintpc = fEtaBinInTPCAcc;
231 phibintpc = fPhiBinInTPCAcc;
232 etabinemc = fEtaBinInEMCalAcc;
233 phibinemc = fPhiBinInEMCalAcc;
234 nbinphi = fNbinPhi;
235}
236
237//__________________________________________________________
238Int_t AliJetGrid::GetIndexFromEtaPhi(Double_t phi,Double_t eta) const
239{ // Tells the index value of a corresponding (eta,phi) real position
240 // Loop over all entries -> takes time.
241 // Used one time at the begining to fill the grids
242
243 /* this is how bins are numbered
244
245 in all TPC or in TPC - EMCAL
246 ... ... ... ..
247 ---------------
248 ... ... . 10 | | | 11
249 ---+---+--- ---------------
250 ^ 6 | 7 | 8 or 8 | | | 9
251 | ---+---+--- ---------------
252 3 | 4 | 5 4 | 5 | 6 | 7
253 phi ---+---+--- ---------------
254 0 | 1 | 2 0 | 1 | 2 | 3
255
256
257 eta ->
258 */
259
260 Int_t etaBin=0,phiBin=0,absID=0;
261 Int_t etaBin2=0,etaBin3=0;
262
263 // Fill all the grid in eta/phi (all TPC acceptance)
264 //-----------------------------------------------------
265 if(fGrid==0){
266 if(eta <= fEta->At(0)) {
267 etaBin = 0;
268 } else if(eta >= fEta->At(fNeta)) {
269 etaBin = fNeta;
270 } else {
271 for(Int_t i=0; i<fNeta+1; i++) {
272 if(eta < fEta->At(i)) {
273 etaBin = i-1;
274 break;
275 }
276 }
277 }
278 if(phi <= fPhi->At(0)) {
279 phiBin = 0;
280 } else if(phi >= fPhi->At(fNphi)) {
281 phiBin = fNphi;
282 } else {
283 for(Int_t i=0; i<fNphi+1; i++) {
284 if(phi < fPhi->At(i)) {
285 phiBin = i-1;
286 break;
287 }
288 }
289 }
290
291 // Calculate absolute id
292 absID = phiBin*(fNeta+1) + etaBin;
293
294 }
295
296 // Fill the grid but do not count id in EMCal acceptance
297 //------------------------------------------------------
298 if((eta >= fEtaMin && eta <= fEtaMax) &&
299 (phi >= fPhiMin && phi <= fPhiMax)){
300 etaBin = etaBin2 = etaBin3 = -1;
301 phiBin = -1;
302 }
303
304 if(fGrid == 1){
305 if(phi<fPhiMin) {
306 if(eta <= fEta->At(0)) {
307 etaBin = 0;
308 } else if(eta >= fEta->At(fNeta)) {
309 etaBin = fNeta;
310 } else {
311 for(Int_t i=0; i<fNeta+1; i++) {
312 if(eta < fEta->At(i)) {
313 etaBin = i-1;
314 break;
315 }
316 }
317 }
318 }
319 else {
320 if((phi>=fPhiMin) && (phi<=fPhiMax)) {
321 if(eta <= fEta->At(0)) {
322 etaBin2 = 0;
323 } else if(eta >= fEta->At(fNeta)) {
324 etaBin2 = fNeta-fEtaBinInEMCalAcc;
325 } else {
326 for(Int_t i=0; i<fNeta+1; i++) {
327 if(eta < fEta->At(i) && ((fEta->At(0)<eta && eta<fEtaMin) ||
328 (fEtaMax<eta && eta<fEta->At(fNeta)))) {
329 // cout << "i : " << i << endl;
330 if(eta<fEtaMin)
331 etaBin2 = i-1;
332 else etaBin2 = i-1-fEtaBinInEMCalAcc;
333 break;
334 }
335 }
336 }
337 }
338 else {
339 if(eta <= fEta->At(0)) {
340 etaBin3 = 0;
341 } else if(eta >= fEta->At(fNeta)) {
342 etaBin3 = fNeta;
343 } else {
344 for(Int_t i=0; i<fNeta+1; i++) {
345 if(eta < fEta->At(i)) {
346 etaBin3 = i-1;
347 break;
348 }
349 }
350 }
351 }
352 }
353
354 if(phi <= fPhi->At(0)) {
355 phiBin = 0;
356 } else if(phi >= fPhi->At(fNphi)) {
357 phiBin = fNphi;
358 } else {
359 for(Int_t i=0; i<fNphi+1; i++) {
360 if(phi < fPhi->At(i)) {
361 phiBin = i-1;
362 break;
363 }
364 }
365 }
366
367 // Calculate absID
368 //-------------------------------------------------------
369 if(phi<fPhiMin)
370 absID = phiBin*(fNeta+1) + etaBin;
371 if(phi>=fPhiMin && phi<=fPhiMax) {
372 if(eta >= fEtaMin && eta <= fEtaMax) absID = -1;
373 else{
374 absID = (fNbinPhi+1)*(fNeta+1)
375 + (phiBin-fNbinPhi-1)*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
376 + etaBin2;
377 }
378 }
379 if(phi>fPhiMax)
380 absID = (fNbinPhi+1)*(fNeta+1)
381 + fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
382 + (phiBin-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)
383 + etaBin3;
384
385 } // END OPTION==1
386
387 return absID;
388
389}
390
391//__________________________________________________________
392void AliJetGrid::GetEtaPhiFromIndex(Int_t index, Float_t &eta, Float_t &phi)
393{ // Get (eta,phi) position for a given index BUT loop over all entries (takes time)
394
395 for(Int_t j=0; j<fNphi+1; j++) {
396 for(Int_t i=0; i<fNeta+1; i++) {
397
398 // TPC grid only
399 //-------------------------------------
400 if(fGrid==0) {
401 if(j*(fNeta+1)+i == index) {
402 eta = fEta->At(i);
403 phi = fPhi->At(j);
404 }
405 }
406
407 // TPC-EMCAL grid
408 //-------------------------------------
409 Int_t ii = 0;
410 if(i==0) ii = 0;
411 if(i>0 && i<(fEtaBinInTPCAcc-fEtaBinInEMCalAcc)/2) ii = i;
412 if(i>=(fEtaBinInTPCAcc+fEtaBinInEMCalAcc)/2 && i<fNeta+1) ii = i-fEtaBinInEMCalAcc;
413
414 if(fGrid==1) {
415 if(j<(fNbinPhi+1) && j*(fNeta+1)+i == index) {
416 eta = fEta->At(i);
417 phi = fPhi->At(j);
418 }
419
420 if((j>=(fNbinPhi+1) && j<(fNbinPhi+1+fPhiBinInEMCalAcc)) &&
421 ((fNbinPhi+1)*(fNeta+1) + (j-fNbinPhi-1)*(fEtaBinInTPCAcc-fEtaBinInEMCalAcc) + ii)== index ) {
422 if(ii==0) {Int_t ind = 0; eta = fEta->At(ind);}
423 else eta = fEta->At(i);
424 phi = fPhi->At(j);
425 }
426
427 if(j>=(fNbinPhi+1+fPhiBinInEMCalAcc) &&
428 ((fNbinPhi+1)*(fNeta+1)+fPhiBinInEMCalAcc*((fEtaBinInTPCAcc-fEtaBinInEMCalAcc))
429 +(j-(fNbinPhi+1+fPhiBinInEMCalAcc))*(fNeta+1)+i == index)) {
430 eta = fEta->At(i);
431 phi = fPhi->At(j);
432 }
433 }
434 }
435 }
436}
437
438//__________________________________________________________
439Int_t AliJetGrid::GetIndex(Double_t phi, Double_t eta)
440{ // Get index value for a (eta,phi) position - Direct value
441
442 Int_t ieta = GetIndexJFromEta(eta);
443 Int_t iphi = GetIndexIFromPhi(phi);
444
445 Int_t index = GetMatrixIndex(iphi,ieta);
446
447 if(fDebug>10){
448 cout << "(phi,eta) : " << phi << ", " << eta << endl;
449 cout << "index : " << index << endl;
450 }
451 return index;
452
453}
454
455//__________________________________________________________
456Int_t AliJetGrid::GetIndexJFromEta(Double_t eta)
457{ // Get eta id
458 // Eta discretized
459
460 Int_t idEta =0;
461 Double_t temp = (eta+fMaxEta)/(fMaxEta-fMinEta)*fNeta;
462 if(fDebug>20)
463 {
464 cout << "eta : " << eta << endl;
465 cout << "fMaxEta : " << fMaxEta << endl;
466 cout << "fMinEta : " << fMinEta << endl;
467 cout << "fNeta : " << fNeta << endl;
468 cout << "temp eta before cast : " << temp << endl;
469 }
470 idEta = static_cast<Int_t>(temp+0.5);
471 if(fDebug>20) cout << "temp eta after cast : " << idEta << endl;
472 return idEta;
473}
474//__________________________________________________________
475Int_t AliJetGrid::GetIndexIFromPhi(Double_t phi)
476{ // Get phi id
477 // Phi discretized
478
479 Int_t idPhi = 0;
480 Double_t temp = 0.;
481 if(fMinPhi==0) temp = phi/(fMaxPhi-fMinPhi)*fNphi;
482 else temp = (phi-fMinPhi)/(fMaxPhi-fMinPhi)*fNphi;
483
484 if(fDebug>20)
485 {
486 cout << "phi : " << phi << endl;
487 cout << "fMaxPhi : " << fMaxPhi << endl;
488 cout << "fMinPhi : " << fMinPhi << endl;
489 cout << "fNphi : " << fNphi << endl;
490 cout << "temp phi before cast : " << temp << endl;
491 }
492 idPhi = static_cast<Int_t>(temp+0.5);
493 if(fDebug>20) cout << "temp phi after cast : " << idPhi << endl;
494 return idPhi;
495
496
497}
498
499//__________________________________________________________
500void AliJetGrid::SetMatrixIndex(Int_t i,Double_t par)
501{ // Allows to set parameters using only one index (if fGrid==0) !!
502 // Not used !
503
504 Int_t iphi = (Int_t)i/fNeta;
505 Int_t ieta = i-iphi*fNeta;
506 SetMatrixIndex(iphi,ieta,par);
507
508 return;
509}
510
511//__________________________________________________________
512void AliJetGrid::SetMatrixIndexes()
513{ // Fill the final matrix object with the corresponding index in eta/phi
514
515 for(Int_t i=0; i<fNphi+1; i++){
516 for(Int_t j=0; j<fNeta+1; j++){
517 (*fIndex)(i,j) = GetIndexFromEtaPhi(fPhi->At(i),fEta->At(j))+1;
518 if(fDebug>2){
519 cout << "(*fIndex)(" << i << "," << j << ") : " << (*fIndex)(i,j) <<
520 ", phi : " << fPhi->At(i) << ", eta : " << fEta->At(j) << endl;
521 }
522 }
523 }
524 printf("##########################################################\n");
525 printf("TMatrix object filled !\n");
526 printf("Size of the object : phi x eta = (fNphi+1) x (fNeta+1) = %d\n",(fNphi+1)*(fNeta+1));
527 printf("##########################################################\n");
528}
529
530//__________________________________________________________
531void AliJetGrid::SetIndexIJ()
532{ //
533
534 for(Int_t i=0; i<fNphi+1; i++){
535 for(Int_t j=0; j<fNeta+1; j++){
536 Int_t id = static_cast<Int_t>((*fIndex)(i,j));
537
538 if(id!=-1)
539 {
540 (*fIndexI)[id] = i;
541 (*fIndexJ)[id] = j;
542 }
543 }
544 }
545
546 printf("##########################################################\n");
547 printf(" In SetIndexIJ - Grid indexes setted !\n");
548 printf("##########################################################\n");
549}
550
551//__________________________________________________________
552void AliJetGrid::GetIJFromIndex(Int_t index, Int_t i, Int_t j)
553{ // returns i position id of eta and j position id of phi for a given grid index
554 i = (*fIndexI)[index];
555 j = (*fIndexJ)[index];
556}
557
558//__________________________________________________________
559void AliJetGrid::GetEtaPhiFromIndex2(Int_t index, Float_t &phi, Float_t &eta)
560{ // returns eta, phi values for a given grid index
561
562 phi = fPhi->At((*fIndexI)[index]);
563 eta = fEta->At((*fIndexJ)[index]);
564}
565
566//__________________________________________________________
567Int_t AliJetGrid::GetNEntries()
568{ // Returns the number of entries of the grid
569
570 if(fDebug>20){
571 cout << "fMaxPhi : " << fMaxPhi << endl;
572 cout << "fMaxEta : " << fMaxEta << endl;
573 }
574
575 Int_t indexNum = GetIndex(fMaxPhi,fMaxEta);
576 if(fDebug>20) cout << "indexNum : " << indexNum << endl;
577 return indexNum;
578
579}
580
581//__________________________________________________________
582Int_t AliJetGrid::GetNEntries2()
583{ // Returns the number of entries of the grid
584
585 Int_t indexNum = GetIndex(fMaxPhi-1.,fMaxEta-0.5);
586 if(fDebug>20) cout << "indexNum : " << indexNum << endl;
587 return indexNum;
588
589}
590
591
592
593
594
595