Logo ROOT   6.18/05
Reference Guide
TEveCaloLegoGL.cxx
Go to the documentation of this file.
1// @(#)root/eve:$Id$
2// Author: Alja Mrak-Tadel 2007
3
4/*************************************************************************
5 * Copyright (C) 1995-2007, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include "TAxis.h"
13#include "TH2.h"
14#include "THLimitsFinder.h"
15
16#include "TGLViewer.h"
17#include "TGLIncludes.h"
18#include "TGLPhysicalShape.h"
19#include "TGLRnrCtx.h"
20#include "TGLSelectRecord.h"
21#include "TGLScene.h"
22#include "TGLCamera.h"
23#include "TGLUtil.h"
24#include "TColor.h"
25#include "TROOT.h"
26
27
28#include "TEveCaloLegoGL.h"
29#include "TEveCalo.h"
30#include "TEveManager.h"
31#include "TEveRGBAPalette.h"
32
33#include <algorithm>
34
35/** \class TEveCaloLegoGL
36\ingroup TEve
37OpenGL renderer class for TEveCaloLego.
38*/
39
41
42////////////////////////////////////////////////////////////////////////////////
43/// Constructor.
44
46 TGLObject(),
47
48 fGridColor(-1),
49 fFontColor(-1),
50
51 fEtaAxis(0),
52 fPhiAxis(0),
53 fZAxis(0),
54 fM(0),
55 fDLCacheOK(kFALSE),
56 fMaxVal(0),
57 fValToPixel(0),
58 fCurrentPixelsPerBin(0),
59 fCells3D(kTRUE),
60 fBinStep(-1)
61{
63
64 fEtaAxis = new TAxis();
65 fPhiAxis = new TAxis();
66 fZAxis = new TAxis();
67
69}
70
71////////////////////////////////////////////////////////////////////////////////
72/// Destructor.
73
75{
77
78 delete fEtaAxis;
79 delete fPhiAxis;
80 delete fZAxis;
81}
82
83////////////////////////////////////////////////////////////////////////////////
84/// Set model object.
85
87{
88 fM = SetModelDynCast<TEveCaloLego>(obj);
89 return kTRUE;
90}
91
92////////////////////////////////////////////////////////////////////////////////
93/// Set bounding box.
94
96{
98}
99
100////////////////////////////////////////////////////////////////////////////////
101/// Drop all display-list definitions.
102
104{
106 for (SliceDLMap_i i = fDLMap.begin(); i != fDLMap.end(); ++i)
107 i->second = 0;
108
110}
111
112////////////////////////////////////////////////////////////////////////////////
113/// Unregister all display-lists.
114
116{
117 // all lego cells
119 if (! fDLMap.empty()) {
120 for (SliceDLMap_i i = fDLMap.begin(); i != fDLMap.end(); ++i) {
121 if (i->second) {
122 PurgeDLRange(i->second, 1);
123 i->second = 0;
124 }
125 }
126 }
128}
129
130////////////////////////////////////////////////////////////////////////////////
131/// Draw an axis-aligned box using quads.
132/// ~~~ {.cpp}
133/// z
134/// |
135/// |
136/// |________y
137/// / 6-------7
138/// / /| /|
139/// x 5-------4 |
140/// | 2-----|-3
141/// |/ |/
142/// 1-------0
143/// ~~~
144
146 Float_t xw, Float_t yw, Float_t h) const
147{
148 Float_t x2 = x1 + xw;
149 Float_t y2 = y1 + yw;
150 Float_t z2 = z1 + h;
151
152 if (x1 < fM->GetEtaMin()) x1 = fM->GetEtaMin();
153 if (x2 > fM->GetEtaMax()) x2 = fM->GetEtaMax();
154
155 if (y1 < fM->GetPhiMin()) y1 = fM->GetPhiMin();
156 if (y2 > fM->GetPhiMax()) y2 = fM->GetPhiMax();
157
158 glBegin(GL_QUADS);
159 {
160 // bottom 0123
161 glNormal3f(0, 0, -1);
162 glVertex3f(x2, y2, z1);
163 glVertex3f(x2, y1, z1);
164 glVertex3f(x1, y1, z1);
165 glVertex3f(x1, y2, z1);
166 // top 4765
167 glNormal3f(0, 0, 1);
168 glVertex3f(x2, y2, z2);
169 glVertex3f(x1, y2, z2);
170 glVertex3f(x1, y1, z2);
171 glVertex3f(x2, y1, z2);
172
173 // back 0451
174 glNormal3f(1, 0, 0);
175 glVertex3f(x2, y2, z1);
176 glVertex3f(x2, y2, z2);
177 glVertex3f(x2, y1, z2);
178 glVertex3f(x2, y1, z1);
179 // front 3267
180 glNormal3f(-1, 0, 0);
181 glVertex3f(x1, y2, z1);
182 glVertex3f(x1, y1, z1);
183 glVertex3f(x1, y1, z2);
184 glVertex3f(x1, y2, z2);
185
186 // left 0374
187 glNormal3f(0, 1, 0);
188 glVertex3f(x2, y2, z1);
189 glVertex3f(x1, y2, z1);
190 glVertex3f(x1, y2, z2);
191 glVertex3f(x2, y2, z2);
192 // right 1562
193 glNormal3f(0, -1, 0);
194 glVertex3f(x2, y1, z1);
195 glVertex3f(x2, y1, z2);
196 glVertex3f(x1, y1, z2);
197 glVertex3f(x1, y1, z1);
198 }
199 glEnd();
200}
201
202////////////////////////////////////////////////////////////////////////////////
203/// Create display-list that draws histogram bars for non-rebinned data.
204/// It is used for filled and outline passes.
205
207{
209 Int_t prevTower = 0;
210 Float_t offset = 0;
211
212 // ids in eta phi rng
213 Int_t nSlices = fM->fData->GetNSlices();
214 for (Int_t s = 0; s < nSlices; ++s)
215 {
216 if (dlMap.empty() || dlMap[s] == 0)
217 dlMap[s] = glGenLists(1);
218
219 glNewList(dlMap[s], GL_COMPILE);
220
221 for (UInt_t i = 0; i < cellList.size(); ++i)
222 {
223 if (cellList[i].fSlice > s) continue;
224 if (cellList[i].fTower != prevTower) {
225 offset = 0;
226 prevTower = cellList[i].fTower;
227 }
228
229 fM->fData->GetCellData(cellList[i], cellData);
230 if (s == cellList[i].fSlice)
231 {
232 if (selection) glLoadName(i);
233
234 WrapTwoPi(cellData.fPhiMin, cellData.fPhiMax);
235 MakeQuad(cellData.EtaMin(), cellData.PhiMin(), offset,
236 cellData.EtaDelta(), cellData.PhiDelta(), cellData.Value(fM->fPlotEt));
237 }
238 offset += cellData.Value(fM->fPlotEt);
239 }
240 glEndList();
241 }
242}
243
244////////////////////////////////////////////////////////////////////////////////
245/// Create display-list that draws histogram bars for rebinned data.
246/// It is used for filled and outline passes.
247
249{
250 Int_t nSlices = fM->fData->GetNSlices();
251 Float_t *vals;
252 Float_t offset;
253 Float_t y0, y1;
254
255 for (Int_t s = 0; s < nSlices; ++s)
256 {
257 if (dlMap.empty() || dlMap[s] == 0)
258 dlMap[s] = glGenLists(1);
259
260 glNewList(dlMap[s], GL_COMPILE);
261
262 for (Int_t i = 1; i <= fEtaAxis->GetNbins(); ++i)
263 {
264 for (Int_t j = 1; j <= fPhiAxis->GetNbins(); ++j)
265 {
266 const Int_t bin = (i)+(j)*(fEtaAxis->GetNbins()+2);
267
268 if (rebinData.fBinData[bin] !=-1)
269 {
270 vals = rebinData.GetSliceVals(bin);
271 offset =0;
272 for (Int_t t = 0; t < s; ++t)
273 offset += vals[t];
274
275 y0 = fPhiAxis->GetBinLowEdge(j);
276 y1 = fPhiAxis->GetBinUpEdge(j);
277 WrapTwoPi(y0, y1);
278
279 if (selection) glLoadName(bin);
280
281 MakeQuad(fEtaAxis->GetBinLowEdge(i), y0, offset,
282 fEtaAxis->GetBinWidth(i), y1-y0, vals[s]);
283 }
284 }
285 }
286 glEndList();
287 }
288}
289
290////////////////////////////////////////////////////////////////////////////////
291/// Set the axis 3D title position.
292
294{
295 const GLdouble *pm = rnrCtx.RefCamera().RefLastNoPickProjM().CArr();
296 GLdouble mm[16];
297 GLint vp[4];
298 glGetDoublev(GL_MODELVIEW_MATRIX, mm);
299 glGetIntegerv(GL_VIEWPORT, vp);
300 GLdouble projX[4], projY[4], projZ[4];
301
302 GLdouble cornerX[4];
303 GLdouble cornerY[4];
304 cornerX[0] = x0; cornerY[0] = y0;
305 cornerX[1] = x1; cornerY[1] = y0;
306 cornerX[2] = x1; cornerY[2] = y1;
307 cornerX[3] = x0; cornerY[3] = y1;
308
309 gluProject(cornerX[0], cornerY[0], 0, mm, pm, vp, &projX[0], &projY[0], &projZ[0]);
310 gluProject(cornerX[1], cornerY[1], 0, mm, pm, vp, &projX[1], &projY[1], &projZ[1]);
311 gluProject(cornerX[2], cornerY[2], 0, mm, pm, vp, &projX[2], &projY[2], &projZ[2]);
312 gluProject(cornerX[3], cornerY[3], 0, mm, pm, vp, &projX[3], &projY[3], &projZ[3]);
313
314
315 // Z axis location (left most corner)
316 //
317 Int_t idxLeft = 0;
318 Float_t xt = projX[0];
319 for (Int_t i = 1; i < 4; ++i) {
320 if (projX[i] < xt) {
321 xt = projX[i];
322 idxLeft = i;
323 }
324 }
325 fZAxisTitlePos.Set(cornerX[idxLeft], cornerY[idxLeft], 1.05 * fMaxVal);
326
327
328 // XY axis location (closest to eye) first
329 //
330 Float_t zt = 1.f;
331 Float_t zMin = 0.f;
332 Int_t idxFront = 0;
333 for (Int_t i = 0; i < 4; ++i) {
334 if (projZ[i] < zt) {
335 zt = projZ[i];
336 idxFront = i;
337 }
338 if (projZ[i] > zMin) zMin = projZ[i];
339 }
340
341
342 Int_t xyIdx = idxFront;
343 if (zMin - zt < 1e-2) xyIdx = 0; // avoid flipping in front view
344
345
346 switch (xyIdx) {
347 case 0:
349 fXAxisTitlePos.fY = y0;
350 fYAxisTitlePos.fX = x0;
351 fYAxisTitlePos.fY = y1;
352 break;
353 case 1:
354 fXAxisTitlePos.fX = x0;
355 fXAxisTitlePos.fY = y0;
357 fYAxisTitlePos.fY = y1;
358 break;
359 case 2:
360 fXAxisTitlePos.fX = x0;
361 fXAxisTitlePos.fY = y1;
363 fYAxisTitlePos.fY = y0;
364 break;
365 case 3:
367 fXAxisTitlePos.fY = y1;
368 fYAxisTitlePos.fX = x0;
369 fYAxisTitlePos.fY = y0;
370 break;
371 }
372
373 // move title 5% over the axis length
374 Float_t off = 0.05;
375 Float_t tOffX = (x1-x0) * off; if (fYAxisTitlePos.fX > x0) tOffX = -tOffX;
376 Float_t tOffY = (y1-y0) * off; if (fXAxisTitlePos.fY > y0) tOffY = -tOffY;
377 fXAxisTitlePos.fX += tOffX;
378 fYAxisTitlePos.fY += tOffY;
379
380
381 // frame box
382 //
383 if (fM->fBoxMode)
384 {
385 // get corner closest to eye excluding left corner
386 Double_t zm = 1.f;
387 Int_t idxDepthT = 0;
388 for (Int_t i = 0; i < 4; ++i)
389 {
390 if (projZ[i] < zm && projZ[i] >= zt && i != idxFront )
391 {
392 zm = projZ[i];
393 idxDepthT = i;
394 }
395 }
396 if (idxFront == idxLeft) idxFront =idxDepthT;
397
398 switch (idxFront)
399 {
400 case 0:
401 fBackPlaneXConst[0].Set(x1, y0, 0); fBackPlaneXConst[1].Set(x1, y1, 0);
402 fBackPlaneYConst[0].Set(x0, y1, 0); fBackPlaneYConst[1].Set(x1, y1, 0);
403 break;
404 case 1:
405 fBackPlaneXConst[0].Set(x0, y0, 0); fBackPlaneXConst[1].Set(x0, y1, 0);
406 fBackPlaneYConst[0].Set(x0, y1, 0); fBackPlaneYConst[1].Set(x1, y1, 0);
407 break;
408 case 2:
409 fBackPlaneXConst[0].Set(x0, y0, 0); fBackPlaneXConst[1].Set(x0, y1, 0);
410 fBackPlaneYConst[0].Set(x0, y0, 0); fBackPlaneYConst[1].Set(x1, y0, 0);
411 break;
412 case 3:
413 fBackPlaneXConst[0].Set(x1, y0, 0); fBackPlaneXConst[1].Set(x1, y1, 0);
414 fBackPlaneYConst[0].Set(x0, y0, 0); fBackPlaneYConst[1].Set(x1, y0, 0);
415 break;
416 }
417 }
418}
419
420////////////////////////////////////////////////////////////////////////////////
421/// Draw z-axis and z-box at the appropriate grid corner-point including
422/// tick-marks and labels.
423
425{
426 // set font size first depending on size of projected axis
427
429 GLdouble pm[16];
430 glGetDoublev(GL_MODELVIEW_MATRIX, mm.Arr());
431 glGetDoublev(GL_PROJECTION_MATRIX, pm);
432 Int_t* vp = rnrCtx.RefCamera().RefViewport().CArr();
433 GLdouble dn[3];
434 GLdouble up[3];
435 gluProject(fXAxisTitlePos.fX, fXAxisTitlePos.fY, fXAxisTitlePos.fZ, mm.Arr(), pm, vp, &up[0], &up[1], &up[2]);
436 gluProject(fYAxisTitlePos.fX, fYAxisTitlePos.fY, fYAxisTitlePos.fZ, mm.Arr(), pm, vp, &dn[0], &dn[1], &dn[2]);
437 Float_t len = TMath::Sqrt((up[0] - dn[0]) * (up[0] - dn[0])
438 + (up[1] - dn[1]) * (up[1] - dn[1])
439 + (up[2] - dn[2]) * (up[2] - dn[2]));
440 len = TMath::Min(len, rnrCtx.RefCamera().RefViewport().Diagonal()*0.7f);
441 len /= TMath::Sqrt2();
442
444 fAxisPainter.RefTMOff(0) = rnrCtx.RefCamera().ViewportDeltaToWorld(worldRef, -len, 0, &mm);
447
448 Float_t tickLength = TMath::Max(fM->GetData()->GetEtaBins()->GetTickLength(), 0.02f);
449 Float_t labelOffset = TMath::Max(fM->GetData()->GetEtaBins()->GetLabelOffset(), 0.02f);
450
451 // Z axis
452 //
453 if (fM->fData->Empty() == kFALSE)
454 {
455 Int_t ondiv;
456 Double_t omin=0, omax=0, bw1;
457 THLimitsFinder::Optimize(0, fMaxVal, fM->fNZSteps, omin, omax, ondiv, bw1);
460 // check z axis title does not overalp with label
462 fZAxisTitlePos.fZ = omax + zto.Z();
463
464
468 fZAxis->SetNdivisions(fM->fNZSteps*100 + 10);
470 fZAxis->SetTitle(fM->GetPlotEt() ? "Et[GeV]" : "E[GeV]");
471
473 fAxisPainter.RefDir().Set(0., 0., 1.);
475 glPushMatrix();
476 glTranslatef(fZAxisTitlePos.fX, fZAxisTitlePos.fY, 0);
477
478 // tickmark vector = 10 pixels left
480 fZAxis->SetLabelOffset(labelOffset);
481 fZAxis->SetTickLength(tickLength);
483 glPopMatrix();
484
485 // draw box frame
486 //
487 if (fM->fBoxMode) {
488
489 glPushAttrib(GL_ENABLE_BIT | GL_LINE_BIT);
490
491 // box verticals
493 glBegin(GL_LINES);
495
496 glVertex3f(fBackPlaneXConst[0].fX ,fBackPlaneXConst[0].fY ,0);
497 glVertex3f(fBackPlaneXConst[0].fX ,fBackPlaneXConst[0].fY ,fMaxVal);
498 glVertex3f(fBackPlaneXConst[1].fX ,fBackPlaneXConst[1].fY ,0);
499 glVertex3f(fBackPlaneXConst[1].fX ,fBackPlaneXConst[1].fY ,fMaxVal);
500
501
502 glVertex3f(fBackPlaneYConst[0].fX ,fBackPlaneYConst[0].fY ,0);
503 glVertex3f(fBackPlaneYConst[0].fX ,fBackPlaneYConst[0].fY ,fMaxVal);
504 glVertex3f(fBackPlaneYConst[1].fX ,fBackPlaneYConst[1].fY ,0);
505 glVertex3f(fBackPlaneYConst[1].fX ,fBackPlaneYConst[1].fY ,fMaxVal);
506
507 // box top
508 glVertex3f(fBackPlaneXConst[0].fX ,fBackPlaneXConst[0].fY ,fMaxVal);
509 glVertex3f(fBackPlaneXConst[1].fX ,fBackPlaneXConst[1].fY ,fMaxVal);
510 glVertex3f(fBackPlaneYConst[0].fX ,fBackPlaneYConst[0].fY ,fMaxVal);
511 glVertex3f(fBackPlaneYConst[1].fX ,fBackPlaneYConst[1].fY ,fMaxVal);
512
513 glEnd();
514
515 // box horizontals stippled
516 glEnable(GL_LINE_STIPPLE);
517 glLineStipple(1, 0x5555);
518 glBegin(GL_LINES);
519 Float_t hz = bw1;
520 for (Int_t i = 1; i <= ondiv; ++i, hz += bw1) {
521 glVertex3f(fBackPlaneXConst[0].fX ,fBackPlaneXConst[0].fY ,hz);
522 glVertex3f(fBackPlaneXConst[1].fX ,fBackPlaneXConst[1].fY ,hz);
523 glVertex3f(fBackPlaneYConst[0].fX ,fBackPlaneYConst[0].fY ,hz);
524 glVertex3f(fBackPlaneYConst[1].fX ,fBackPlaneYConst[1].fY ,hz);
525 }
526 glEnd();
527
528 glPopAttrib();
529 }
530 }
531
532 // XY Axis
533 //
534
535 Float_t yOff = fM->GetPhiRng();
536 if (fXAxisTitlePos.fY < fM->GetPhiMax()) yOff = -yOff;
537
538 Float_t xOff = fM->GetEtaRng();
539 if (fYAxisTitlePos.fX < fM->GetEtaMax()) xOff = -xOff;
540
541 TAxis ax;
546 ax.SetLabelOffset(labelOffset);
547 ax.SetTickLength(tickLength);
549 fAxisPainter.RefTMOff(1).Set(0, 0, -fMaxVal);
551
552 // eta
553 glPushMatrix();
554 fAxisPainter.RefDir().Set(1, 0, 0);
555 fAxisPainter.RefTMOff(0).Set(0, yOff, 0);
556 glTranslatef(0, fXAxisTitlePos.fY, 0);
557
559 ax.SetLimits(fM->GetEtaMin(), fM->GetEtaMax());
562 fAxisPainter.PaintAxis(rnrCtx, &ax);
563 glPopMatrix();
564
565 // phi
566 fAxisPainter.RefDir().Set(0, 1, 0);
567 fAxisPainter.RefTMOff(0).Set(xOff, 0, 0);
569 ax.SetLimits(fM->GetPhiMin(), fM->GetPhiMax());
571 glPushMatrix();
572 glTranslatef(fYAxisTitlePos.fX, 0, 0);
574 fAxisPainter.PaintAxis(rnrCtx, &ax);
575 glPopMatrix();
576
577} // DrawAxis3D
578
579////////////////////////////////////////////////////////////////////////////////
580/// Get scale for matrix.
581
583{
584 Double_t em, eM, pm, pM;
585 fM->fData->GetEtaLimits(em, eM);
586 fM->fData->GetPhiLimits(pm, pM);
587 Double_t unit = ((eM - em) < (pM - pm)) ? (eM - em) : (pM - pm);
588 sx = (eM - em) / (fM->GetEtaRng() * unit);
589 sy = (pM - pm) / (fM->GetPhiRng() * unit);
590
591 sz = 1;
592 if (fM->fScaleAbs)
593 {
594 sz = fM->GetMaxTowerH() / fM->fMaxValAbs;
595 }
596 else if (!fM->fData->Empty())
597 {
598 sz = fM->GetMaxTowerH() / fMaxVal;
599 }
600}
601
602////////////////////////////////////////////////////////////////////////////////
603/// Draw XY axis.
604
606{
607 if (fM->GetData()->Empty())
609
610 TGLCamera& cam = rnrCtx.RefCamera();
611
612 TAxis ax;
620
621 // set fonts
623
624 // get projected length of diagonal to determine
626 GLdouble pm[16];
627 GLint vp[4];
628 glGetDoublev(GL_MODELVIEW_MATRIX, mm.Arr());
629 glGetDoublev(GL_PROJECTION_MATRIX, pm);
630 glGetIntegerv(GL_VIEWPORT, vp);
631
632 GLdouble dn[3];
633 GLdouble up[3];
634 gluProject(fM->GetEtaMin(), fM->GetPhiMin(), 0, mm.Arr(), pm, vp, &dn[0], &dn[1], &dn[2]);
635 gluProject(fM->GetEtaMax(), fM->GetPhiMax(), 0, mm.Arr(), pm, vp, &up[0], &up[1], &up[2]);
636 Double_t len = TMath::Sqrt((up[0] - dn[0]) * (up[0] - dn[0])
637 + (up[1] - dn[1]) * (up[1] - dn[1])
638 + (up[2] - dn[2]) * (up[2] - dn[2]));
639
640 // lock upper limit to of relative font size relative to viewport diagonal
641 Double_t vpLimit = cam.RefViewport().Diagonal()*0.5/TMath::Sqrt2();
642 len = TMath::Min(len, vpLimit);
643
644 // eta
648 ax.SetLimits(fM->GetEtaMin(), fM->GetEtaMax());
651 fAxisPainter.RefDir().Set(1, 0, 0);
652
654 fAxisPainter.RefTMOff(0).Set(0, -TMath::Min(fM->GetPhiRng(), tmOffFrustX), 0);
656
657 glPushMatrix();
658 glTranslatef(0, fM->GetPhiMin(), 0);
659 fAxisPainter.PaintAxis(rnrCtx, &ax);
660 glPopMatrix();
661
662 // phi
664 ax.SetLimits(fM->GetPhiMin(), fM->GetPhiMax());
667 fAxisPainter.RefDir().Set(0, 1, 0);
669 fAxisPainter.RefTMOff(0).Set(-TMath::Min(fM->GetEtaRng(), tmOffFrustY), 0, 0);
671
672 glPushMatrix();
673 glTranslatef(fM->GetEtaMin(), 0, 0);
674 fAxisPainter.PaintAxis(rnrCtx, &ax);
675 glPopMatrix();
676
678}
679
680////////////////////////////////////////////////////////////////////////////////
681/// Calculate view-dependent grid density.
682
684{
685 TGLCamera &camera = rnrCtx.RefCamera();
690 Float_t frustD = TMath::Hypot(r-l, t-b);
691
692 GLint vp[4]; glGetIntegerv(GL_VIEWPORT, vp);
693 Float_t viewportD = TMath::Sqrt((vp[1] - vp[0]) * (vp[1] - vp[0]) + (vp[3] - vp[1]) * (vp[3] - vp[1]));
694 Float_t deltaToViewport = viewportD/frustD;
695
696 // average bin width
697 GLdouble em, eM, pm, pM;
698 fM->GetData()->GetEtaLimits(pm, pM);
699 fM->GetData()->GetPhiLimits(em, eM);
704
705 Float_t averageBinWidth = TMath::Hypot(eM - em, pM - pm)/TMath::Sqrt((i0 - i1) * (i0 - i1) + (j0 - j1) * (j0 - j1));
706 Float_t ppb = deltaToViewport*averageBinWidth;
707
708 Int_t ngroup = 1;
709 if (fM->fAutoRebin && fM->fPixelsPerBin > ppb)
710 {
711 // limit rebin realtive to number of axis bins
712 Int_t maxGroup = TMath::Min(fM->fData->GetEtaBins()->GetNbins(), fM->fData->GetPhiBins()->GetNbins())/4;
713 if (maxGroup > 1) {
714 ngroup = TMath::Nint(fM->fPixelsPerBin*0.5/ppb); // symetrical rebin factor 2
715 if (ngroup > maxGroup) ngroup = maxGroup;
716 }
717 }
719
720 return ngroup;
721}
722
723////////////////////////////////////////////////////////////////////////////////
724/// Rebin eta, phi axis.
725
726void TEveCaloLegoGL::RebinAxis(TAxis *orig, TAxis *curr) const
727{
728 Double_t center = 0.5 * (orig->GetXmin() + orig->GetXmax());
729 Int_t idx0 = orig->FindBin(center);
730 Double_t bc = orig->GetBinCenter(idx0);
731 if (bc > center) --idx0;
732
733 Int_t nbR = TMath::FloorNint(idx0/fBinStep) + TMath::FloorNint((orig->GetNbins() - idx0)/fBinStep);
734 Int_t off = idx0 - TMath::FloorNint(idx0/fBinStep)*fBinStep;
735 std::vector<Double_t> bins(nbR + 1);
736 for (Int_t i = 0; i <= nbR; ++i)
737 {
738 bins[i] = orig->GetBinUpEdge(off + i*fBinStep);
739 }
740 curr->Set(nbR, &bins[0]);
741}
742
743////////////////////////////////////////////////////////////////////////////////
744/// Draw basic histogram components: x-y grid
745
747{
748 Float_t eta0 = fM->fEtaMin;
749 Float_t eta1 = fM->fEtaMax;
750 Float_t phi0 = fM->GetPhiMin();
751 Float_t phi1 = fM->GetPhiMax();
752
753 // XY grid
754 //
757 glBegin(GL_LINES);
758 glVertex2f(eta0, phi0);
759 glVertex2f(eta0, phi1);
760 glVertex2f(eta1, phi0);
761 glVertex2f(eta1, phi1);
762
763 glVertex2f(eta0, phi0);
764 glVertex2f(eta1, phi0);
765 glVertex2f(eta0, phi1);
766 glVertex2f(eta1, phi1);
767
768 // eta grid
769 Float_t val;
770 Int_t neb = fEtaAxis->GetNbins();
771 for (Int_t i = 0; i<= neb; i++)
772 {
773 val = fEtaAxis->GetBinUpEdge(i);
774 if (val > eta0 && val < eta1 )
775 {
776 glVertex2f(val, phi0);
777 glVertex2f(val, phi1);
778 }
779 }
780
781 // phi grid
782 Int_t npb = fPhiAxis->GetNbins();
783 for (Int_t i = 1; i <= npb; i++) {
784 val = fPhiAxis->GetBinUpEdge(i);
785 if (val > phi0 && val < phi1)
786 {
787 glVertex2f(eta0, val);
788 glVertex2f(eta1, val);
789 }
790 }
791
792 glEnd();
793
794 // XYZ axes
795 //
796 glPushAttrib(GL_ENABLE_BIT | GL_LINE_BIT | GL_POLYGON_BIT);
798 if (fCells3D)
799 {
800 SetAxis3DTitlePos(rnrCtx, eta0, eta1, phi0, phi1);
801 DrawAxis3D(rnrCtx);
802 }
803 else
804 {
805 DrawAxis2D(rnrCtx);
806 }
807 glPopAttrib();
808}
809
810////////////////////////////////////////////////////////////////////////////////
811/// Render the calo lego-plot with OpenGL.
812
814{
815 // quads
816 {
817 for (SliceDLMap_i i = fDLMap.begin(); i != fDLMap.end(); ++i) {
819 glLoadName(i->first);
820 glPushName(0);
821 glCallList(i->second);
822 glPopName();
823 }
824 }
825 // outlines
826 {
827 if (rnrCtx.SceneStyle() == TGLRnrCtx::kFill) {
828 glPolygonMode(GL_FRONT_AND_BACK, GL_LINE);
829 glDisable(GL_POLYGON_OFFSET_FILL);
831 for (SliceDLMap_i i = fDLMap.begin(); i != fDLMap.end(); ++i)
832 glCallList(i->second);
833 }
834 }
835}
836
837////////////////////////////////////////////////////////////////////////////////
838/// Prepare cells 2D data non-rebinned for drawing.
839
841{
842 Int_t max_energy_slice, cellID=0;
843 Float_t sum, max_energy;
844
845 TEveCaloData::vCellId_t::iterator currentCell = cellList.begin();
846 TEveCaloData::vCellId_t::iterator nextCell = currentCell;
847 ++nextCell;
848
849 while (true)
850 {
851 TEveCaloData::CellData_t currentCellData;
852 TEveCaloData::CellData_t nextCellData;
853
854 fM->fData->GetCellData(*currentCell, currentCellData);
855 sum = max_energy = currentCellData.Value(fM->fPlotEt);
856 max_energy_slice = currentCell->fSlice;
857 while (nextCell != cellList.end() && currentCell->fTower == nextCell->fTower)
858 {
859 fM->fData->GetCellData(*nextCell, nextCellData);
860 Float_t energy = nextCellData.Value(fM->fPlotEt);
861 sum += energy;
862 if (energy > max_energy)
863 {
864 max_energy = energy;
865 max_energy_slice = nextCell->fSlice;
866 }
867 ++nextCell;
868 ++cellID;
869 }
870
871 WrapTwoPi(currentCellData.fPhiMin, currentCellData.fPhiMax);
872 cells2D.push_back(Cell2D_t(cellID, sum, max_energy_slice));
873 cells2D.back().SetGeom(currentCellData.fEtaMin, currentCellData.fEtaMax,
874 currentCellData.fPhiMin, currentCellData.fPhiMax);
875
876 if (nextCell == cellList.end())
877 break;
878
879 currentCell = nextCell;
880 ++nextCell;
881 ++cellID;
882 }
883}
884
885////////////////////////////////////////////////////////////////////////////////
886/// Prepare cells 2D rebinned data for drawing.
887
889{
890 const Int_t nEta = fEtaAxis->GetNbins();
891 const Int_t nPhi = fPhiAxis->GetNbins();
892 std::vector<Float_t> vec;
893 vec.assign((nEta + 2)*(nPhi + 2), 0.f);
894 std::vector<Float_t> max_e;
895 std::vector<Int_t> max_e_slice;
896 max_e.assign((nEta + 2) * (nPhi + 2), 0.f);
897 max_e_slice.assign((nEta + 2) * (nPhi + 2), -1);
898
899 for (UInt_t bin = 0; bin < rebinData.fBinData.size(); ++bin) {
900 Float_t ssum = 0;
901 if (rebinData.fBinData[bin] != -1) {
902 Float_t *val = rebinData.GetSliceVals(bin);
903 for (Int_t s = 0; s < rebinData.fNSlices; ++s) {
904 ssum += val[s];
905 if (val[s] > max_e[bin]) {
906 max_e[bin] = val[s];
907 max_e_slice[bin] = s;
908 }
909 }
910 }
911 vec[bin] = ssum;
912 }
913
914 // smallest threshold
915 Float_t threshold = fM->GetDataSliceThreshold(0);
916 for (Int_t s = 1; s < fM->fData->GetNSlices(); ++s) {
917 if (threshold > fM->GetDataSliceThreshold(s))
918 threshold = fM->GetDataSliceThreshold(s);
919 }
920
921 // write cells
922 for (Int_t i = 1; i <= fEtaAxis->GetNbins(); ++i) {
923 for (Int_t j = 1; j <= fPhiAxis->GetNbins(); ++j) {
924 const Int_t bin = j * (nEta + 2) + i;
925 if (vec[bin] > threshold && rebinData.fBinData[bin] != -1) {
926 cells2D.push_back(Cell2D_t(bin, vec[bin], max_e_slice[bin]));
927 cells2D.back().SetGeom(fEtaAxis->GetBinLowEdge(i), fEtaAxis->GetBinUpEdge(i),
929 }
930 }
931 }
932}
933
934////////////////////////////////////////////////////////////////////////////////
935/// Draw cells in top view.
936
938{
939 Float_t bws = -1; //smallest bin
940 Float_t logMax = -1;
941
943
945 {
946 fM->AssertPalette();
947 UChar_t col[4];
948
949 for (vCell2D_i i = cells2D.begin(); i != cells2D.end(); ++i)
950 {
951 if (rnrCtx.SecSelection()) glLoadName(i->fId);
952 glBegin(GL_POLYGON);
953 fM->fPalette->ColorFromValue(TMath::FloorNint(i->fSumVal), col);
954 col[3] = fM->GetData()->GetSliceTransparency(i->fMaxSlice);
956 Float_t z = fM->GetHasFixedHeightIn2DMode() ? baseOffset : i->fSumVal;
957 glVertex3f(i->fX0, i->fY0, z);
958 glVertex3f(i->fX1, i->fY0, z);
959 glVertex3f(i->fX1, i->fY1, z);
960 glVertex3f(i->fX0, i->fY1, z);
961 glEnd();
962 }
963 }
964 else
965 {
966 Float_t maxv = 0;
967 bws = 1e5;
968 for (vCell2D_i i = fCells2D.begin(); i != fCells2D.end(); ++i)
969 {
970 if (i->MinSize() < bws) bws = i->MinSize();
971 if (i->fSumVal > maxv) maxv = i->fSumVal;
972 }
973 bws *= 0.5f;
974 logMax = TMath::Log10(maxv + 1);
975 fValToPixel = bws / logMax;
976
977 if (rnrCtx.SecSelection())
978 {
979 // Special draw for name stack.
980
981 for (vCell2D_i i = cells2D.begin(); i != cells2D.end(); ++i)
982 {
983 glLoadName(i->fMaxSlice);
984 glPushName(i->fId);
985
986 glBegin(GL_QUADS);
987 Float_t z = fM->GetHasFixedHeightIn2DMode() ? baseOffset : i->fSumVal;
988 glVertex3f(i->fX0, i->fY0, z);
989 glVertex3f(i->fX1, i->fY0, z);
990 glVertex3f(i->fX1, i->fY1, z);
991 glVertex3f(i->fX0, i->fY1, z);
992 glEnd();
993
994 glPopName();
995 }
996 }
997 else
998 {
999 // Optimised draw without name stack.
1000
1001 if ( ! rnrCtx.Highlight())
1002 {
1003 glBegin(GL_POINTS);
1004 for (vCell2D_i i = cells2D.begin(); i != cells2D.end(); ++i)
1005 {
1007 Float_t z = fM->GetHasFixedHeightIn2DMode() ? baseOffset : i->fSumVal;
1008 glVertex3f(i->X(), i->Y() , z);
1009 }
1010 glEnd();
1011 }
1012
1013 glBegin(GL_QUADS);
1014 for (vCell2D_i i = cells2D.begin(); i != cells2D.end(); ++i)
1015 {
1017 Float_t bw = fValToPixel*TMath::Log10(i->fSumVal+1);
1018 Float_t x = i->X();
1019 Float_t y = i->Y();
1020 Float_t z = fM->GetHasFixedHeightIn2DMode() ? baseOffset : i->fSumVal;
1021 glVertex3f(x - bw, y - bw, z);
1022 glVertex3f(x + bw, y - bw, z);
1023 glVertex3f(x + bw, y + bw, z);
1024 glVertex3f(x - bw, y + bw, z);
1025 }
1026 glEnd();
1027
1029 {
1030 glPushAttrib(GL_ENABLE_BIT | GL_POLYGON_BIT);
1031 Float_t z = 0;
1032 Float_t zOff = fMaxVal*0.001 ; // avoid polygon stipling
1033 glBegin(GL_QUADS);
1034 for ( vCell2D_i i = cells2D.begin(); i != cells2D.end(); ++i) {
1035 Char_t transp = TMath::Min(100, 80 + fM->fData->GetSliceTransparency(i->fMaxSlice) / 5);
1036 TGLUtil::ColorTransparency(fM->fData->GetSliceColor(i->fMaxSlice), transp);
1037 z = fM->GetHasFixedHeightIn2DMode() ? baseOffset : i->fSumVal;
1038 z -= zOff;
1039 glVertex3f(i->fX0, i->fY0, z);
1040 glVertex3f(i->fX1, i->fY0, z);
1041 glVertex3f(i->fX1, i->fY1, z);
1042 glVertex3f(i->fX0, i->fY1, z);
1043 }
1044 glEnd();
1045
1046 glPolygonMode(GL_FRONT_AND_BACK, GL_LINE);
1047 glBegin(GL_QUADS);
1048 for ( vCell2D_i i = cells2D.begin(); i != cells2D.end(); ++i) {
1049 TGLUtil::ColorTransparency(fM->fData->GetSliceColor(i->fMaxSlice), 60);
1050 z = fM->GetHasFixedHeightIn2DMode() ? baseOffset : i->fSumVal;
1051 z += zOff;
1052 glVertex3f(i->fX0, i->fY0, z);
1053 glVertex3f(i->fX1, i->fY0, z);
1054 glVertex3f(i->fX1, i->fY1, z);
1055 glVertex3f(i->fX0, i->fY1, z);
1056 }
1057 glEnd();
1058 glPopAttrib();
1059 }
1060 }
1061 }
1062
1063 // text
1065 ! rnrCtx.Selection() && ! rnrCtx.Highlight())
1066 {
1068 TGLFont font;
1070 const char* txt;
1071 for (vCell2D_i i = cells2D.begin(); i != cells2D.end(); ++i) {
1072
1073 Float_t val = i->fSumVal;
1074 if (val > 10)
1075 txt = Form("%d", TMath::Nint(val));
1076 else if (val > 1 )
1077 txt = Form("%.1f", val);
1078 else if (val > 0.01 )
1079 txt = Form("%.2f", 0.01*TMath::Nint(val*100));
1080 else
1081 txt = Form("~1e%d", TMath::Nint(TMath::Log10(val)));
1082
1083 font.Render(txt, i->X(), i->Y(), val*1.2, TGLFont::kCenterH, TGLFont::kCenterV);
1084 }
1085 }
1086}
1087
1088////////////////////////////////////////////////////////////////////////////////
1089/// Draw highligted cells.
1090
1091void TEveCaloLegoGL::DrawHighlight(TGLRnrCtx& rnrCtx, const TGLPhysicalShape* /*pshp*/, Int_t /*lvl*/) const
1092{
1093 if (fM->fData->GetCellsSelected().empty() && fM->fData->GetCellsHighlighted().empty())
1094 {
1095 return;
1096 }
1097
1098 // modelview matrix
1099 glPushMatrix();
1100 Float_t sx, sy, sz;
1101 GetScaleForMatrix(sx, sy, sz);
1102 glScalef(sx, sy, sz);
1103 glTranslatef(-fM->GetEta(), -fM->fPhi, 0);
1104
1105 if (fCells3D)
1106 {
1107 glPushAttrib(GL_ENABLE_BIT | GL_LINE_BIT | GL_POLYGON_BIT);
1108 glDisable(GL_LIGHTING);
1109 glDisable(GL_CULL_FACE);
1110 glPolygonMode(GL_FRONT_AND_BACK, GL_LINE);
1112 }
1113
1115 if (!fM->fData->GetCellsHighlighted().empty())
1116 {
1117 glColor4ubv(rnrCtx.ColorSet().Selection(3).CArr());
1119 }
1120 if (!fM->fData->GetCellsSelected().empty())
1121 {
1122 glColor4ubv(rnrCtx.ColorSet().Selection(1).CArr());
1124 }
1126
1127 if (fCells3D)
1128 {
1129 glPopAttrib();
1130 }
1131
1132 glPopMatrix();
1133}
1134
1135
1136////////////////////////////////////////////////////////////////////////////////
1137/// Draw selected cells in highlight mode.
1138
1140{
1141 // check eta&phi range of selected cells
1142 TEveCaloData::vCellId_t cellsSelected;
1143 TEveCaloData::CellData_t cellData;
1144 for (TEveCaloData::vCellId_i i = cellsSelectedInput.begin(); i != cellsSelectedInput.end(); ++i)
1145 {
1146 fM->fData->GetCellData((*i), cellData);
1147 if (fM->CellInEtaPhiRng(cellData))
1148 cellsSelected.push_back(*i);
1149 }
1150
1151 // prepare rebin for 2D or 3D if necessary
1152 TEveCaloData::RebinData_t rebinDataSelected;
1153 if (fBinStep > 1)
1154 {
1155 fM->fData->Rebin(fEtaAxis, fPhiAxis, cellsSelected, fM->fPlotEt, rebinDataSelected);
1156 if (fM->fNormalizeRebin) {
1157 Float_t scale = 1.f / (fBinStep * fBinStep);
1158 for (std::vector<Float_t>::iterator it = rebinDataSelected.fSliceData.begin();
1159 it != rebinDataSelected.fSliceData.end(); ++it)
1160 (*it) *= scale;
1161 }
1162 }
1163
1164 if (fCells3D)
1165 {
1166 Float_t offset = 0;
1167 if (fBinStep == 1)
1168 {
1169 for (TEveCaloData::vCellId_i j = cellsSelected.begin(); j != cellsSelected.end(); ++j)
1170 {
1171 offset = 0;
1172 {
1173 Int_t orig_slice = j->fSlice;
1174 for (Int_t s = 0; s < orig_slice; ++s)
1175 {
1176 j->fSlice = s;
1177 fM->fData->GetCellData(*j, cellData);
1178 offset += cellData.Value(fM->fPlotEt);
1179 }
1180 j->fSlice = orig_slice;
1181 }
1182 fM->fData->GetCellData(*j, cellData);
1183 WrapTwoPi(cellData.fPhiMin, cellData.fPhiMax);
1184 MakeQuad(cellData.EtaMin(), cellData.PhiMin(), offset,
1185 cellData.EtaDelta(), cellData.PhiDelta(), cellData.Value(fM->fPlotEt));
1186 }
1187 }
1188 else
1189 {
1190 Float_t *vals;
1191 Float_t *valsRef;
1192 Float_t y0, y1;
1193 Int_t nSlices = fM->fData->GetNSlices();
1194 for (Int_t i = 1; i <= fEtaAxis->GetNbins(); ++i)
1195 {
1196 for (Int_t j = 1; j <= fPhiAxis->GetNbins(); ++j)
1197 {
1198 const Int_t bin = (i)+(j)*(fEtaAxis->GetNbins()+2);
1199 if (rebinDataSelected.fBinData[bin] !=-1)
1200 {
1201 offset = 0;
1202 vals = rebinDataSelected.GetSliceVals(bin);
1203 valsRef = fRebinData.GetSliceVals(bin);
1204 for (Int_t s = 0; s < nSlices; ++s)
1205 {
1206 if (vals[s] > 0)
1207 {
1208 y0 = fPhiAxis->GetBinLowEdge(j);
1209 y1 = fPhiAxis->GetBinUpEdge(j);
1210 WrapTwoPi(y0, y1);
1211 MakeQuad(fEtaAxis->GetBinLowEdge(i), y0, offset,
1212 fEtaAxis->GetBinWidth(i), y1-y0, vals[s]);
1213 }
1214 offset += valsRef[s];
1215 }
1216 }
1217 }
1218 }
1219 }
1220 }
1221 else
1222 {
1223 vCell2D_t cells2DSelected;
1224 if (fBinStep == 1)
1225 {
1226 // but is confusing since top view does not draw all slices at same time
1227 TEveCaloData::vCellId_i j = cellsSelectedInput.begin();
1228 TEveCaloData::vCellId_i jEnd = cellsSelectedInput.end();
1229 std::set<Int_t> towers;
1230 while (j != jEnd)
1231 {
1232 towers.insert(j->fTower);
1233 ++j;
1234 }
1235 for (vCell2D_i i = fCells2D.begin(); i != fCells2D.end(); ++i)
1236 {
1237 TEveCaloData::CellId_t cell = fM->fCellList[i->fId];
1238 // std::set<Int_t>::iterator ti = towers.find(cell.fTower);
1239 if (towers.find(cell.fTower) != towers.end())
1240 {
1241 cells2DSelected.push_back(*i);
1242 }
1243 }
1244 }
1245 else
1246 {
1247 PrepareCell2DDataRebin(rebinDataSelected, cells2DSelected);
1248 }
1249 DrawCells2D(rnrCtx, cells2DSelected);
1250 }
1251}
1252
1253////////////////////////////////////////////////////////////////////////////////
1254/// Draw the object.
1255
1257{
1258 if (! fM->fData || ! fM->fData->GetEtaBins() || ! fM->fData->GetPhiBins())
1259 return;
1260
1261 // projection type
1263 fCells3D = (!(rnrCtx.RefCamera().IsOrthographic() && rnrCtx.RefCamera().GetCamBase().GetBaseVec(1).Z()));
1264 else if (fM->fProjection == TEveCaloLego::k2D)
1265 fCells3D = kFALSE;
1266 else if (fM->fProjection == TEveCaloLego::k3D)
1267 fCells3D = kTRUE;
1268
1269 // rebin axsis , check limits, fix TwoPi cycling
1270 Int_t new_bin_step = GetGridStep(rnrCtx);
1271
1272 // rebin data
1273 if (fM->AssertCellIdCache() || fBinStep != new_bin_step)
1274 {
1275 fBinStep = new_bin_step;
1277 fRebinData.Clear();
1278
1281
1282 if (fBinStep > 1)
1283 {
1285
1286 fMaxVal = 0;
1287 for (UInt_t i = 0; i < fRebinData.fSliceData.size(); i += fRebinData.fNSlices)
1288 {
1289 Double_t sum = 0;
1290 for (Int_t s = 0; s < fRebinData.fNSlices; s++)
1291 {
1292 sum += fRebinData.fSliceData[i+s];
1293 }
1294 if (sum > fMaxVal) fMaxVal = sum;
1295 }
1296
1297 if (fM->fNormalizeRebin)
1298 {
1299 Float_t scale = 1.f / (fBinStep * fBinStep);
1300 for (std::vector<Float_t>::iterator it = fRebinData.fSliceData.begin(); it != fRebinData.fSliceData.end();
1301 ++it) {
1302 (*it) *= scale;
1303 }
1304 fMaxVal *= scale;
1305 }
1306 }
1307 else
1308 {
1309 fMaxVal = fM->GetMaxVal();
1310 }
1311 }
1312
1313 // modelview matrix
1314 glPushMatrix();
1315 Float_t sx, sy, sz;
1316 GetScaleForMatrix(sx, sy, sz);
1317 glScalef(sx, sy, sz);
1318 glTranslatef(-fM->GetEta(), -fM->fPhi, 0);
1319
1322 if (fGridColor < 0 || fFontColor < 0)
1323 {
1324 TColor* c1 = gROOT->GetColor(rnrCtx.ColorSet().Markup().GetColorIndex());
1325 TColor* c2 = gROOT->GetColor(rnrCtx.ColorSet().Background().GetColorIndex());
1326 Float_t f1, f2;
1327 if (fFontColor < 0) {
1328 f1 = 0.8; f2 = 0.2;
1329 fFontColor = TColor::GetColor(c1->GetRed() *f1 + c2->GetRed() *f2,
1330 c1->GetGreen()*f1 + c2->GetGreen()*f2,
1331 c1->GetBlue() *f1 + c2->GetBlue() *f2);
1332 }
1333 if (fGridColor < 0) {
1334 f1 = 0.3; f2 = 0.3;
1335 fGridColor = TColor::GetColor(c1->GetRed() *f1 + c2->GetRed() *f2,
1336 c1->GetGreen()*f1 + c2->GetGreen()*f2,
1337 c1->GetBlue() *f1 + c2->GetBlue() *f2);
1338 }
1339 }
1340
1341 glPushAttrib(GL_ENABLE_BIT | GL_LINE_BIT | GL_POLYGON_BIT);
1343 glEnable(GL_BLEND);
1344 glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
1345 if (!fM->fData->Empty())
1346 {
1347 glPushName(0);
1348 if (fCells3D)
1349 {
1350 if (fDLCacheOK == kFALSE)
1351 {
1352 if (fBinStep == 1)
1354 else
1356 fDLCacheOK = kTRUE;
1357 }
1358 glEnable(GL_NORMALIZE);
1359 glEnable(GL_POLYGON_OFFSET_FILL);
1360 glPolygonOffset(0.8, 1);
1361
1362 DrawCells3D(rnrCtx);
1363 }
1364 else
1365 {
1366 glDisable(GL_LIGHTING);
1367
1368 fCells2D.clear();
1369 if (fBinStep == 1)
1371 else
1373
1374 DrawCells2D(rnrCtx, fCells2D);
1375 }
1376 glPopName();
1377 }
1378 glPopAttrib();
1379
1380 // draw histogram base
1381 if (rnrCtx.Selection() == kFALSE && rnrCtx.IsDrawPassFilled())
1382 {
1383 glPushAttrib(GL_ENABLE_BIT | GL_LINE_BIT | GL_POLYGON_BIT);
1384 glDisable(GL_LIGHTING);
1385 DrawHistBase(rnrCtx);
1386 if (fM->fDrawHPlane) {
1387 glBlendFunc(GL_ONE, GL_ONE_MINUS_SRC_ALPHA);
1388 glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);
1389 glDisable(GL_CULL_FACE);
1391 Float_t zhp = fM->fHPlaneVal * fMaxVal;
1392 glBegin(GL_POLYGON);
1393 glVertex3f(fM->fEtaMin, fM->GetPhiMin(), zhp);
1394 glVertex3f(fM->fEtaMax, fM->GetPhiMin(), zhp);
1395 glVertex3f(fM->fEtaMax, fM->GetPhiMax(), zhp);
1396 glVertex3f(fM->fEtaMin, fM->GetPhiMax(), zhp);
1397 glEnd();
1398 }
1399 glPopAttrib();
1400 }
1401
1402 glPopMatrix();
1403}
1404
1405////////////////////////////////////////////////////////////////////////////////
1406/// Processes tower selection from TGLViewer.
1407
1409{
1411 if (rec.GetN() > 2)
1412 {
1413 Int_t slice = rec.GetItem(1);
1414 Int_t cell = rec.GetItem(2);
1415
1416 if (fBinStep == 1)
1417 {
1418 Int_t tower = fM->fCellList[cell].fTower;
1419 while (cell > 0 && tower == fM->fCellList[cell].fTower)
1420 {
1421 sel.push_back(fM->fCellList[cell]);
1422 if (fCells3D) break;
1423 --cell;
1424 }
1425 }
1426 else
1427 {
1428 if (cell > 0)
1429 {
1430 Int_t nEta = fEtaAxis->GetNbins();
1431 Int_t phiBin = Int_t(cell/(nEta+2));
1432 Int_t etaBin = cell - phiBin*(nEta+2);
1435 fPhiAxis->GetBinCenter(phiBin), fPhiAxis->GetBinWidth(phiBin),
1436 sl);
1437
1438 for (TEveCaloData::vCellId_i it = sl.begin(); it != sl.end(); ++it)
1439 {
1440 if (fCells3D) {
1441 if ((*it).fSlice == slice ) sel.push_back(*it);
1442 } else {
1443 if ((*it).fSlice <= slice ) sel.push_back(*it);
1444 }
1445 }
1446 }
1447 }
1448 }
1449 fM->fData->ProcessSelection(sel, rec);
1450}
ROOT::R::TRInterface & r
Definition: Object.C:4
#define b(i)
Definition: RSha256.hxx:100
#define h(i)
Definition: RSha256.hxx:106
#define e(i)
Definition: RSha256.hxx:103
static const double x2[5]
static const double x1[5]
int Int_t
Definition: RtypesCore.h:41
unsigned char UChar_t
Definition: RtypesCore.h:34
char Char_t
Definition: RtypesCore.h:29
unsigned int UInt_t
Definition: RtypesCore.h:42
const Bool_t kFALSE
Definition: RtypesCore.h:88
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
float Float_t
Definition: RtypesCore.h:53
const Bool_t kTRUE
Definition: RtypesCore.h:87
const char Option_t
Definition: RtypesCore.h:62
#define ClassImp(name)
Definition: Rtypes.h:365
#define gROOT
Definition: TROOT.h:414
char * Form(const char *fmt,...)
virtual Int_t GetNdivisions() const
Definition: TAttAxis.h:36
virtual Style_t GetTitleFont() const
Definition: TAttAxis.h:46
virtual Float_t GetLabelOffset() const
Definition: TAttAxis.h:40
virtual void SetAxisColor(Color_t color=1, Float_t alpha=1.)
Set color of the line axis and tick marks.
Definition: TAttAxis.cxx:163
virtual void SetTitleFont(Style_t font=62)
Set the title font.
Definition: TAttAxis.cxx:322
virtual void SetLabelOffset(Float_t offset=0.005)
Set distance between the axis and the labels The distance is expressed in per cent of the pad width.
Definition: TAttAxis.cxx:193
virtual void SetTitleSize(Float_t size=0.04)
Set size of axis title The size is expressed in per cent of the pad width.
Definition: TAttAxis.cxx:304
virtual void SetTitleColor(Color_t color=1)
Set color of axis title.
Definition: TAttAxis.cxx:313
virtual Float_t GetTitleSize() const
Definition: TAttAxis.h:43
virtual Float_t GetLabelSize() const
Definition: TAttAxis.h:41
virtual Float_t GetTickLength() const
Definition: TAttAxis.h:44
virtual void SetTickLength(Float_t length=0.03)
Set tick mark length The length is expressed in per cent of the pad width.
Definition: TAttAxis.cxx:280
virtual void SetNdivisions(Int_t n=510, Bool_t optim=kTRUE)
Set the number of divisions for this axis.
Definition: TAttAxis.cxx:229
virtual void SetLabelColor(Color_t color=1, Float_t alpha=1.)
Set color of labels.
Definition: TAttAxis.cxx:173
Class to manage histogram axis.
Definition: TAxis.h:30
virtual Double_t GetBinCenter(Int_t bin) const
Return center of bin.
Definition: TAxis.cxx:464
Double_t GetXmax() const
Definition: TAxis.h:134
virtual Int_t FindBin(Double_t x)
Find bin number corresponding to abscissa x.
Definition: TAxis.cxx:279
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition: TAxis.cxx:504
virtual void Set(Int_t nbins, Double_t xmin, Double_t xmax)
Initialize axis with fix bins.
Definition: TAxis.cxx:717
virtual void SetLimits(Double_t xmin, Double_t xmax)
Definition: TAxis.h:154
Double_t GetXmin() const
Definition: TAxis.h:133
Int_t GetNbins() const
Definition: TAxis.h:121
const char * GetTitle() const
Returns title of object.
Definition: TAxis.h:129
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width.
Definition: TAxis.cxx:526
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
Definition: TAxis.cxx:514
The color creation and management class.
Definition: TColor.h:19
static Int_t GetColor(const char *hexcolor)
Static method returning color number for color specified by hex color string of form: "#rrggbb",...
Definition: TColor.cxx:1764
virtual void GetEtaLimits(Double_t &min, Double_t &max) const =0
std::vector< CellId_t > vCellId_t
Definition: TEveCaloData.h:146
Bool_t Empty() const
Definition: TEveCaloData.h:216
virtual void GetCellList(Float_t etaMin, Float_t etaMax, Float_t phi, Float_t phiRng, vCellId_t &out) const =0
void ProcessSelection(vCellId_t &sel_cells, TGLSelectRecord &rec)
Process newly selected cells with given select-record.
virtual void GetCellData(const CellId_t &id, CellData_t &data) const =0
Char_t GetSliceTransparency(Int_t slice) const
Get transparency for given slice.
virtual TAxis * GetEtaBins() const
Definition: TEveCaloData.h:218
Color_t GetSliceColor(Int_t slice) const
Get color for given slice.
Int_t GetNSlices() const
Definition: TEveCaloData.h:202
vCellId_t & GetCellsHighlighted()
Definition: TEveCaloData.h:189
virtual void GetPhiLimits(Double_t &min, Double_t &max) const =0
virtual TAxis * GetPhiBins() const
Definition: TEveCaloData.h:221
vCellId_t & GetCellsSelected()
Definition: TEveCaloData.h:188
virtual void Rebin(TAxis *ax, TAxis *ay, vCellId_t &in, Bool_t et, RebinData_t &out) const =0
std::vector< CellId_t >::iterator vCellId_i
Definition: TEveCaloData.h:147
OpenGL renderer class for TEveCaloLego.
std::vector< Cell2D_t > vCell2D_t
void MakeQuad(Float_t x, Float_t y, Float_t z, Float_t xw, Float_t yw, Float_t zh) const
Draw an axis-aligned box using quads.
TEveVector fBackPlaneYConst[2]
virtual void DirectDraw(TGLRnrCtx &rnrCtx) const
Draw the object.
Color_t fFontColor
TEveVector fZAxisTitlePos
void WrapTwoPi(Float_t &min, Float_t &max) const
void PrepareCell2DData(TEveCaloData::vCellId_t &cellList, vCell2D_t &cells2D) const
Prepare cells 2D data non-rebinned for drawing.
Int_t fCurrentPixelsPerBin
void DrawSelectedCells(TGLRnrCtx &rnrCtx, TEveCaloData::vCellId_t cells) const
Draw selected cells in highlight mode.
void DrawCells3D(TGLRnrCtx &rnrCtx) const
Render the calo lego-plot with OpenGL.
void DrawHistBase(TGLRnrCtx &rnrCtx) const
Draw basic histogram components: x-y grid.
TEveVector fXAxisTitlePos
void DrawCells2D(TGLRnrCtx &rnrCtx, vCell2D_t &cells2D) const
Draw cells in top view.
void DrawAxis3D(TGLRnrCtx &rnrCtx) const
Draw z-axis and z-box at the appropriate grid corner-point including tick-marks and labels.
virtual void DrawHighlight(TGLRnrCtx &rnrCtx, const TGLPhysicalShape *ps, Int_t lvl=-1) const
Draw highligted cells.
void RebinAxis(TAxis *orig, TAxis *curr) const
Rebin eta, phi axis.
Int_t GetGridStep(TGLRnrCtx &rnrCtx) const
Calculate view-dependent grid density.
std::vector< Cell2D_t >::iterator vCell2D_i
vCell2D_t fCells2D
void GetScaleForMatrix(Float_t &sx, Float_t &sy, Float_t &sz) const
Get scale for matrix.
virtual void ProcessSelection(TGLRnrCtx &rnrCtx, TGLSelectRecord &rec)
Processes tower selection from TGLViewer.
TEveVector fBackPlaneXConst[2]
void Make3DDisplayListRebin(TEveCaloData::RebinData_t &rebinData, SliceDLMap_t &map, Bool_t select) const
Create display-list that draws histogram bars for rebinned data.
virtual void SetBBox()
Set bounding box.
virtual Bool_t SetModel(TObject *obj, const Option_t *opt=0)
Set model object.
Float_t fValToPixel
void SetAxis3DTitlePos(TGLRnrCtx &rnrCtx, Float_t x0, Float_t x1, Float_t y0, Float_t y1) const
Set the axis 3D title position.
TEveCaloLegoGL()
Constructor.
TEveCaloLego * fM
virtual void DLCacheDrop()
Drop all display-list definitions.
void DrawAxis2D(TGLRnrCtx &rnrCtx) const
Draw XY axis.
virtual void DLCachePurge()
Unregister all display-lists.
void Make3DDisplayList(TEveCaloData::vCellId_t &cellList, SliceDLMap_t &map, Bool_t select) const
Create display-list that draws histogram bars for non-rebinned data.
TGLAxisPainter fAxisPainter
std::map< Int_t, UInt_t > SliceDLMap_t
TEveVector fYAxisTitlePos
std::map< Int_t, UInt_t >::iterator SliceDLMap_i
virtual ~TEveCaloLegoGL()
Destructor.
void PrepareCell2DDataRebin(TEveCaloData::RebinData_t &rebinData, vCell2D_t &cells2D) const
Prepare cells 2D rebinned data for drawing.
SliceDLMap_t fDLMap
Color_t fGridColor
TEveCaloData::RebinData_t fRebinData
Visualization of calorimeter data as eta/phi histogram.
Definition: TEveCalo.h:250
EProjection_e fProjection
Definition: TEveCalo.h:278
Bool_t fNormalizeRebin
Definition: TEveCalo.h:276
Char_t fPlaneTransparency
Definition: TEveCalo.h:269
Color_t fPlaneColor
Definition: TEveCalo.h:268
@ kValSizeOutline
Definition: TEveCalo.h:256
EBoxMode_e fBoxMode
Definition: TEveCalo.h:280
Int_t fDrawNumberCellPixels
Definition: TEveCalo.h:288
E2DMode_e f2DMode
Definition: TEveCalo.h:279
Float_t fHPlaneVal
Definition: TEveCalo.h:283
float GetFixedHeightValIn2DMode() const
Definition: TEveCalo.h:333
TEveCaloData::vCellId_t fCellList
Definition: TEveCalo.h:264
Color_t fFontColor
Definition: TEveCalo.h:266
Bool_t fDrawHPlane
Definition: TEveCalo.h:282
Bool_t fAutoRebin
Definition: TEveCalo.h:274
bool GetHasFixedHeightIn2DMode() const
Definition: TEveCalo.h:330
Int_t fCellPixelFontSize
Definition: TEveCalo.h:289
Color_t fGridColor
Definition: TEveCalo.h:267
Int_t fNZSteps
Definition: TEveCalo.h:271
Int_t fPixelsPerBin
Definition: TEveCalo.h:275
Color_t GetDataSliceColor(Int_t slice) const
Get slice color from data.
Definition: TEveCalo.cxx:121
Float_t GetEtaMin() const
Definition: TEveCalo.h:136
Bool_t AssertCellIdCache() const
Assert cell id cache is ok.
Definition: TEveCalo.cxx:293
Float_t fPlotEt
Definition: TEveCalo.h:54
Float_t GetDataSliceThreshold(Int_t slice) const
Get threshold for given slice.
Definition: TEveCalo.cxx:86
Float_t GetMaxTowerH() const
Definition: TEveCalo.h:112
Float_t GetPhiRng() const
Definition: TEveCalo.h:146
TEveRGBAPalette * fPalette
Definition: TEveCalo.h:61
Double_t fEtaMax
Definition: TEveCalo.h:43
Float_t GetMaxVal() const
Definition: TEveCalo.cxx:159
TEveCaloData * GetData() const
Definition: TEveCalo.h:86
Float_t GetEta() const
Definition: TEveCalo.h:135
TEveRGBAPalette * AssertPalette()
Make sure the TEveRGBAPalette pointer is not null.
Definition: TEveCalo.cxx:378
Bool_t GetPlotEt() const
Definition: TEveCalo.h:108
Float_t GetEtaMax() const
Definition: TEveCalo.h:137
Bool_t fScaleAbs
Definition: TEveCalo.h:57
Double_t fPhi
Definition: TEveCalo.h:45
TEveCaloData * fData
Definition: TEveCalo.h:39
Float_t GetPhiMax() const
Definition: TEveCalo.h:145
Double_t fEtaMin
Definition: TEveCalo.h:42
Float_t GetEtaRng() const
Definition: TEveCalo.h:138
Float_t fMaxValAbs
Definition: TEveCalo.h:58
Bool_t CellInEtaPhiRng(TEveCaloData::CellData_t &) const
Returns true if given cell is in the ceta phi range.
Definition: TEveCalo.cxx:307
Float_t GetPhiMin() const
Definition: TEveCalo.h:144
const UChar_t * ColorFromValue(Int_t val) const
void Set(const Float_t *v)
Definition: TEveVector.h:81
TGLVector3 & RefDir()
void SetLabelAlign(TGLFont::ETextAlignH_e, TGLFont::ETextAlignV_e)
Set label align.
TGLVector3 & RefTitlePos()
void SetAttAxis(TAttAxis *a)
void SetLabelPixelFontSize(Int_t fs)
void SetTMNDim(Int_t x)
TGLVector3 & RefTMOff(Int_t i)
void SetFontMode(TGLFont::EMode m)
Int_t GetLabelPixelFontSize() const
void PaintAxis(TGLRnrCtx &ctx, TAxis *ax)
GL render TAxis.
void SetTitlePixelFontSize(Int_t fs)
Abstract base camera class - concrete classes for orthographic and perspective cameras derive from it...
Definition: TGLCamera.h:44
TGLMatrix & RefLastNoPickProjM() const
Definition: TGLCamera.h:174
const TGLMatrix & GetCamBase() const
Definition: TGLCamera.h:166
virtual Bool_t IsOrthographic() const
Definition: TGLCamera.h:118
const TGLPlane & FrustumPlane(EFrustumPlane plane) const
Definition: TGLCamera.h:219
TGLRect & RefViewport()
Definition: TGLCamera.h:128
TGLVector3 ViewportDeltaToWorld(const TGLVertex3 &worldRef, Double_t viewportXDelta, Double_t viewportYDelta, TGLMatrix *modviewMat=0) const
Apply a 2D viewport delta (shift) to the projection of worldRef onto viewport, returning the resultan...
Definition: TGLCamera.cxx:546
TGLColor & Markup()
Definition: TGLUtil.h:853
TGLColor & Selection(Int_t i)
Definition: TGLUtil.h:854
TGLColor & Background()
Definition: TGLUtil.h:850
const UChar_t * CArr() const
Definition: TGLUtil.h:799
Color_t GetColorIndex() const
Returns color-index representing the color.
Definition: TGLUtil.cxx:1218
A wrapper class for FTFont.
void Render(const char *txt, Double_t x, Double_t y, Double_t angle, Double_t mgn) const
virtual void DLCachePurge()
Purge all entries for all LODs for this drawable from the display list cache, returning the reserved ...
TObject * fExternalObj
first replica
virtual void DLCacheDrop()
Drop all entries for all LODs for this drawable from the display list cache, WITHOUT returning the re...
void PurgeDLRange(UInt_t base, Int_t size) const
External object is a fake.
Bool_t fDLCache
display-list validity bit-field
16 component (4x4) transform matrix - column MAJOR as per GL.
Definition: TGLUtil.h:597
TGLVector3 GetBaseVec(Int_t b) const
Definition: TGLUtil.h:753
const Double_t * CArr() const
Definition: TGLUtil.h:663
Base-class for direct OpenGL renderers.
Definition: TGLObject.h:22
void SetAxisAlignedBBox(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax)
Set axis-aligned bounding-box.
Definition: TGLObject.cxx:86
Concrete physical shape - a GL drawable.
Double_t D() const
Definition: TGLUtil.h:555
const Int_t * CArr() const
Definition: TGLUtil.h:442
Int_t Diagonal() const
Return the diagonal of the rectangle.
Definition: TGLUtil.cxx:287
The TGLRnrCtx class aggregates data for a given redering context as needed by various parts of the RO...
Definition: TGLRnrCtx.h:41
Short_t SceneStyle() const
Definition: TGLRnrCtx.h:184
Bool_t SecSelection() const
Definition: TGLRnrCtx.h:224
void RegisterFontNoScale(Int_t size, Int_t file, Int_t mode, TGLFont &out)
Get font in the GL rendering context.
Definition: TGLRnrCtx.cxx:367
Bool_t IsDrawPassFilled() const
Returns true if current render-pass uses filled polygon style.
Definition: TGLRnrCtx.cxx:153
TGLColorSet & ColorSet()
Return reference to current color-set (top of the stack).
Definition: TGLRnrCtx.cxx:278
TGLCamera & RefCamera()
Definition: TGLRnrCtx.h:157
Bool_t Highlight() const
Definition: TGLRnrCtx.h:218
Bool_t Selection() const
Definition: TGLRnrCtx.h:222
UInt_t GetItem(Int_t i) const
Int_t GetN() const
Standard selection record including information about containing scene and details ob out selected ob...
static void Color4ubv(const UChar_t *rgba)
Wrapper for glColor4ubv.
Definition: TGLUtil.cxx:1774
static UInt_t LockColor()
Prevent further color changes.
Definition: TGLUtil.cxx:1664
static void ColorTransparency(Color_t color_index, Char_t transparency=0)
Set color from color_index and ROOT-style transparency (default 0).
Definition: TGLUtil.cxx:1736
static UInt_t UnlockColor()
Allow color changes.
Definition: TGLUtil.cxx:1672
static void Color(const TGLColor &color)
Set color from TGLColor.
Definition: TGLUtil.cxx:1692
static Float_t LineWidth()
Get the line-width, taking the global scaling into account.
Definition: TGLUtil.cxx:1938
3 component (x/y/z) vector class.
Definition: TGLUtil.h:247
3 component (x/y/z) vertex class.
Definition: TGLUtil.h:83
Double_t X() const
Definition: TGLUtil.h:118
Double_t Z() const
Definition: TGLUtil.h:122
void Set(Double_t x, Double_t y, Double_t z)
Definition: TGLUtil.h:209
Double_t Y() const
Definition: TGLUtil.h:120
static void Optimize(Double_t A1, Double_t A2, Int_t nold, Double_t &BinLow, Double_t &BinHigh, Int_t &nbins, Double_t &BWID, Option_t *option="")
Static function to compute reasonable axis limits.
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
Mother of all ROOT objects.
Definition: TObject.h:37
return c1
Definition: legend1.C:41
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
TF1 * f1
Definition: legend1.C:11
return c2
Definition: legend2.C:14
static constexpr double s
static constexpr double mm
Int_t Nint(T x)
Round to nearest integer. Rounds half integers to the nearest even integer.
Definition: TMath.h:701
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
constexpr Double_t Sqrt2()
Definition: TMath.h:89
Int_t FloorNint(Double_t x)
Definition: TMath.h:695
Double_t Sqrt(Double_t x)
Definition: TMath.h:679
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
Double_t Hypot(Double_t x, Double_t y)
Definition: TMath.cxx:57
Double_t Log10(Double_t x)
Definition: TMath.h:752
Cell data inner structure.
Definition: TEveCaloData.h:115
Float_t Value(Bool_t) const
Return energy value associated with the cell, usually Et.
Float_t PhiDelta() const
Definition: TEveCaloData.h:98
Float_t PhiMin() const
Definition: TEveCaloData.h:95
Float_t EtaDelta() const
Definition: TEveCaloData.h:93
Float_t EtaMin() const
Definition: TEveCaloData.h:90
std::vector< Float_t > fSliceData
Definition: TEveCaloData.h:132
std::vector< Int_t > fBinData
Definition: TEveCaloData.h:133
Float_t * GetSliceVals(Int_t bin)
auto * l
Definition: textangle.C:4
static long int sum(long int i)
Definition: Factory.cxx:2258