ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
TEveProjectionAxesGL.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 "TEveProjectionAxesGL.h"
13 #include "TEveProjectionAxes.h"
14 #include "TEveProjectionManager.h"
15 #include "THLimitsFinder.h"
16 
17 #include "TGLIncludes.h"
18 #include "TGLRnrCtx.h"
19 #include "TGLFontManager.h"
20 #include "TGLCamera.h"
21 
22 #include "TMath.h"
23 
24 /** \class TEveProjectionAxesGL
25 \ingroup TEve
26 OpenGL renderer class for TEveProjectionAxes.
27 */
28 
30 
31 ////////////////////////////////////////////////////////////////////////////////
32 /// Constructor.
33 
35  TGLObject(),
36  fM(0),
37  fProjection(0)
38 {
39  fDLCache = kFALSE; // Disable display list.
40 }
41 
42 ////////////////////////////////////////////////////////////////////////////////
43 /// Set model object.
44 /// Virtual from TGLObject.
45 
47 {
48  fM = SetModelDynCast<TEveProjectionAxes>(obj);
50  return fM->GetManager() ? kTRUE : kFALSE;
51 }
52 
53 ////////////////////////////////////////////////////////////////////////////////
54 /// Fill the bounding-box data of the logical-shape.
55 /// Virtual from TGLObject.
56 
58 {
60 }
61 
62 ////////////////////////////////////////////////////////////////////////////////
63 /// Filter overlapping labels.
64 
66 {
68  if (orig.size() == 0) return;
69 
70  Float_t center = fM->GetManager()->GetProjection()->GetProjectedCenter()[idx];
71 
72  // Get index of label closest to the distortion center.
73  // Needed to keep symmetry around center.
74  Int_t minIdx = 0;
75  Int_t cnt = 0;
76  Float_t currD = 0;
77  Float_t minD = TMath::Abs(orig[0].first -center);
78  for (TGLAxisPainter::LabVec_t::iterator it = orig.begin(); it != orig.end(); ++it)
79  {
80  currD = TMath::Abs((*it).first - center);
81  if (minD > currD)
82  {
83  minD = currD;
84  minIdx = cnt;
85  }
86  cnt++;
87  }
88 
89  // Minimum allowed distance 4* font size.
90  TGLAxisPainter::LabVec_t filtered;
91  filtered.push_back(orig[minIdx]);
92  Int_t size = orig.size();
93  Float_t minDist = 4*fM->GetLabelSize()*ref;
94  Float_t pos = 0;
95 
96  // Go from center to minimum.
97  if (minIdx > 0)
98  {
99  pos = orig[minIdx].first;
100  for (Int_t i=minIdx-1; i>=0; --i)
101  {
102  if (TMath::Abs(pos - orig[i].first) > minDist)
103  {
104  filtered.push_back(orig[i]);
105  pos = orig[i].first;
106  }
107  }
108  }
109 
110  // Go from center to maximum.
111  if (minIdx < (size -1))
112  {
113  pos = orig[minIdx].first;
114  for (Int_t i=minIdx+1; i<size; ++i)
115  {
116  if (TMath::Abs(orig[i].first - pos) > minDist)
117  {
118  filtered.push_back(orig[i]);
119  pos = orig[i].first;
120  }
121  }
122  }
123 
124  // Set labels list and text format.
125  if (filtered.size() >= 2)
126  {
127  if ( minIdx > 0 )
128  fAxisPainter.SetTextFormat(orig.front().second, orig.back().second, orig[minIdx].second - orig[minIdx-1].second);
129  else
130  fAxisPainter.SetTextFormat(orig.front().second, orig.back().second, orig[minIdx+1].second - orig[minIdx].second);
131 
132  fAxisPainter.RefLabVec().swap(filtered);
133  }
134  else
135  {
136  fAxisPainter.SetTextFormat(orig.front().second, orig.back().second, orig.back().second - orig.front().second);
137  }
138 }
139 
140 ////////////////////////////////////////////////////////////////////////////////
141 /// Build an array of tick-mark position-value pairs.
142 
144 {
145  fAxisPainter.RefLabVec().clear();
146  fAxisPainter.RefTMVec().clear();
147 
148  // Get list of label position-value pairs.
149 
150 
151  // Minimum/maximum are defined at the front/back element of list.
152  fAxisPainter.RefTMVec().push_back(TGLAxisPainter::TM_t(p1, -1));
153 
155  {
156  SplitIntervalByVal(p1, p2, ax);
157  }
159  {
160  SplitIntervalByPos(p1, p2, ax);
161  }
162 }
163 
164 ////////////////////////////////////////////////////////////////////////////////
165 /// Add tick-marks at equidistant position.
166 
168 {
169  // Limits.
170  Int_t n1a = TMath::FloorNint(fM->GetNdivisions() / 100);
171  Int_t n2a = fM->GetNdivisions() - n1a * 100;
172  Int_t bn1, bn2;
173  Double_t bw1, bw2; // bin with first second order
174  Double_t bl1=0, bh1=0, bl2=0, bh2=0; // bin low, high first second order
175  THLimitsFinder::Optimize(p1, p2, n1a, bl1, bh1, bn1, bw1);
176  THLimitsFinder::Optimize(bl1, bl1+bw1, n2a, bl2, bh2, bn2, bw2);
177 
178  Int_t n1=TMath::CeilNint(p1/bw1);
179  Int_t n2=TMath::FloorNint(p2/bw1);
180 
183 
184  Float_t p = n1*bw1;
185  Float_t pMinor = p;
186  for (Int_t l=n1; l<=n2; l++)
187  {
188  // Labels.
189  labVec.push_back( TGLAxisPainter::Lab_t(p , fProjection->GetValForScreenPos(ax, p)));
190 
191  // Tick-marks.
192  tmVec.push_back(TGLAxisPainter::TM_t(p, 0));
193  pMinor = p+bw2;
194  for (Int_t i=1; i<bn2; i++)
195  {
196  if (pMinor > p2) break;
197  tmVec.push_back( TGLAxisPainter::TM_t(pMinor, 1));
198  pMinor += bw2;
199  }
200  p += bw1;
201  }
202 
203  // Complete second order tick-marks.
204  pMinor = n1*bw1 -bw2;
205  while ( pMinor > p1)
206  {
207  tmVec.push_back(TGLAxisPainter::TM_t(pMinor, 1));
208  pMinor -=bw2;
209  }
210 }
211 
212 ////////////////////////////////////////////////////////////////////////////////
213 /// Add tick-marks on fixed value step.
214 
216 {
217 
220 
221  // Limits
222  Int_t n1a = TMath::FloorNint(fM->GetNdivisions() / 100);
223  Int_t n2a = fM->GetNdivisions() - n1a * 100;
224  Int_t bn1, bn2;
225  Double_t bw1, bw2; // bin width first / second order
226  Double_t bl1=0, bh1=0, bl2=0, bh2=0; // bin low, high first / second order
228  Float_t v2 = fProjection->GetValForScreenPos(ax, p2);
229  THLimitsFinder::Optimize(v1, v2, n1a, bl1, bh1, bn1, bw1);
230  THLimitsFinder::Optimize(bl1, bl1+bw1, n2a, bl2, bh2, bn2, bw2);
231 
232  Float_t pFirst, pSecond; // position of first, second order of tick-marks
233  Float_t v = bl1;
234 
235  // cache values here
236  TEveVector dirVec;
237  fProjection->SetDirectionalVector(ax, dirVec);
238  TEveVector oCenter;
239  fProjection->GetOrthogonalCenter(ax, oCenter);
240 
241  // step
242  for (Int_t l=0; l<=bn1; l++)
243  {
244  // Labels.
245  pFirst = fProjection->GetScreenVal(ax, v);
246  labVec.push_back(TGLAxisPainter::Lab_t(pFirst , v));
247  tmVec.push_back(TGLAxisPainter::TM_t(pFirst, 0));
248 
249  // Tick-marks.
250  for (Int_t k=1; k<bn2; k++)
251  {
252  pSecond = fProjection->GetScreenVal(ax, v+k*bw2, dirVec, oCenter);
253  if (pSecond > p2) break;
254  tmVec.push_back(TGLAxisPainter::TM_t(pSecond, 1));
255  }
256  v += bw1;
257  }
258 
259  // Complete second order tick-marks.
260  v = bl1 -bw2;
261  while ( v > v1)
262  {
263  pSecond = fProjection->GetScreenVal(ax, v, dirVec, oCenter);
264  if (pSecond < p1) break;
265  tmVec.push_back(TGLAxisPainter::TM_t(pSecond, 1));
266  v -= bw2;
267  }
268 }
269 
270 ////////////////////////////////////////////////////////////////////////////////
271 /// Get range from bounding box of projection manager and frustum size.
272 
274 {
275  Float_t* bb = fM->fManager->GetBBox();
276  // enlarge bbox times 2
277  Float_t bbMin = bb[ax*2];
278  Float_t bbMax = bb[ax*2 + 1];
279  Float_t off = ( bb[ax*2 + 1] - bb[ax*2]) * 0.5;
280  bbMin -= off;
281  bbMax += off;
282 
283 
284  // minimum
285  if (frustMin > bbMin) {
286  min = frustMin;
287  min += (frustMax - frustMin) * 0.1;
288  }
289  else {
290  min = bbMin;
291  }
292 
293  // maximum
294  if (frustMax < bbMax) {
295  max = frustMax;
296  max -= (frustMax - frustMin) * 0.1;
297  }
298  else {
299  max = bbMax;
300  }
301 }
302 
303 ////////////////////////////////////////////////////////////////////////////////
304 /// Draw function for TEveProjectionAxesGL. Skips line-pass of outline mode.
305 
307 {
308  if (rnrCtx.IsDrawPassOutlineLine())
309  return;
310 
311  TGLObject::Draw(rnrCtx);
312 }
313 
314 ////////////////////////////////////////////////////////////////////////////////
315 /// Actual rendering code.
316 /// Virtual from TGLLogicalShape.
317 
319 {
320  if (rnrCtx.Selection() || rnrCtx.Highlight() || fM->fManager->GetBBox() == 0) return;
321 
322  glPushAttrib(GL_ENABLE_BIT | GL_POLYGON_BIT);
323 
324  glDisable(GL_LIGHTING);
325  glColorMaterial(GL_FRONT_AND_BACK, GL_DIFFUSE);
326  glEnable(GL_COLOR_MATERIAL);
327  glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);
328  glDisable(GL_CULL_FACE);
329 
330  // Draw on front-clipping plane.
331  Float_t old_depth_range[2];
332  glGetFloatv(GL_DEPTH_RANGE, old_depth_range);
333  glDepthRange(0, 0.001);
334 
335  // Frustum size.
336  TGLCamera &camera = rnrCtx.RefCamera();
337  Float_t l = -camera.FrustumPlane(TGLCamera::kLeft).D();
339  Float_t t = camera.FrustumPlane(TGLCamera::kTop).D();
340  Float_t b = -camera.FrustumPlane(TGLCamera::kBottom).D();
341 
342  if (fM->fUseColorSet)
343  {
344  TGLUtil::Color(rnrCtx.ColorSet().Markup());
346  }
347 
349  glDisable(GL_LIGHTING);
350  // Projection center and origin marker.
351  {
352  Float_t d = ((r-l) > (b-t)) ? (b-t) : (r-l);
353  d *= 0.02f;
354  if (fM->GetDrawCenter())
355  {
358  glBegin(GL_LINES);
359  glVertex3f(c[0] + d, c[1], c[2]); glVertex3f(c[0] - d, c[1], c[2]);
360  glVertex3f(c[0], c[1] + d, c[2]); glVertex3f(c[0], c[1] - d, c[2]);
361  glVertex3f(c[0], c[1], c[2] + d); glVertex3f(c[0], c[1], c[2] - d);
362  glEnd();
363  }
364  if (fM->GetDrawOrigin())
365  {
366  TEveVector zero;
367  fProjection->ProjectVector(zero, 0);
369  glBegin(GL_LINES);
370  glVertex3f(zero[0] + d, zero[1], zero[2]); glVertex3f(zero[0] - d, zero[1], zero[2]);
371  glVertex3f(zero[0], zero[1] + d, zero[2]); glVertex3f(zero[0], zero[1] - d, zero[2]);
372  glVertex3f(zero[0], zero[1], zero[2] + d); glVertex3f(zero[0], zero[1], zero[2] - d);
373  glEnd();
374  }
375  }
376 
377  //
378  // Axes.
379  try {
380  using namespace TMath;
381  GLint vp[4];
382  glGetIntegerv(GL_VIEWPORT, vp);
383  Float_t refLength = TMath::Sqrt((TMath::Power(vp[2]-vp[0], 2) + TMath::Power(vp[3]-vp[1], 2)));
384  Float_t tickLength = TMath::Sqrt((TMath::Power(r-l, 2) + TMath::Power(t-b, 2)));
387 
388  Float_t min, max;
389  // X-axis.
392  {
393  GetRange(0, l, r, min, max);
394  SplitInterval(min, max, 0);
395 
396  FilterOverlappingLabels(0, r-l);
397  fAxisPainter.RefTMVec().push_back(TGLAxisPainter::TM_t(max, -1));
398 
399  fAxisPainter.RefDir().Set(1, 0, 0);
400  fAxisPainter.RefTMOff(0).Set(0, tickLength, 0);
401 
402  // Bottom.
403  glPushMatrix();
404  glTranslatef( 0, b, 0);
408  glPopMatrix();
409 
410  // Top.
411  glPushMatrix();
412  glTranslatef( 0, t, 0);
417  glPopMatrix();
418  }
419 
420  // Y-axis.
423  {
424  GetRange(1, b, t, min, max);
425  SplitInterval(min, max, 1);
426 
427  FilterOverlappingLabels(1, t-b);
428  fAxisPainter.RefTMVec().push_back(TGLAxisPainter::TM_t(max, -1));
429 
430  fAxisPainter.RefDir().Set(0, 1, 0);
431  fAxisPainter.RefTMOff(0).Set(tickLength, 0 , 0);
432 
433  // Left.
434  glPushMatrix();
435  glTranslatef(l, 0, 0);
439  glPopMatrix();
440 
441  // Right.
442  glPushMatrix();
443  glTranslatef(r, 0, 0);
448  glPopMatrix();
449  }
450  }
451  catch (TEveException& exc)
452  {
453  Warning("TEveProjectionAxesGL::DirectDraw", "caught exception: '%s'.", exc.Data());
454  }
455 
456  glDepthRange(old_depth_range[0], old_depth_range[1]);
457 
458  glPopAttrib();
459 }
The TGLRnrCtx class aggregates data for a given redering context as needed by various parts of the RO...
Definition: TGLRnrCtx.h:40
void GetRange(Int_t ax, Float_t frustMin, Float_t frustMax, Float_t &start, Float_t &en) const
Get range from bounding box of projection manager and frustum size.
static Vc_ALWAYS_INLINE int_v min(const int_v &x, const int_v &y)
Definition: vector.h:433
Abstract base camera class - concrete classes for orthographic and perspective cameras derive from it...
Definition: TGLCamera.h:43
virtual void SetDirectionalVector(Int_t screenAxis, TEveVector &vec)
Get vector for axis in a projected space.
void SetUseAxisColors(Bool_t x)
static void Color(const TGLColor &color)
Set color from TGLColor.
Definition: TGLUtil.cxx:1658
const Double_t * v1
Definition: TArcBall.cxx:33
TEveProjectionManager * fManager
static const char * GetFontNameFromId(Int_t)
Get font name from TAttAxis font id.
float Float_t
Definition: RtypesCore.h:53
Double_t D() const
Definition: TGLUtil.h:559
return c
const char Option_t
Definition: RtypesCore.h:62
virtual void SetBBox()
Fill the bounding-box data of the logical-shape.
std::vector< TM_t > TMVec_t
virtual void Draw(TGLRnrCtx &rnrCtx) const
Draw function for TEveProjectionAxesGL. Skips line-pass of outline mode.
TGLCamera & RefCamera()
Definition: TGLRnrCtx.h:157
Bool_t Selection() const
Definition: TGLRnrCtx.h:222
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
TEveProjectionAxes * fM
void SetFontMode(TGLFont::EMode m)
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
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
const Bool_t kFALSE
Definition: Rtypes.h:92
TGLVector3 & RefDir()
virtual Float_t GetScreenVal(Int_t ax, Float_t value)
Project point on given axis and return projected value.
std::pair< Float_t, Int_t > TM_t
Int_t FloorNint(Double_t x)
Definition: TMath.h:476
TGLColor & Markup()
Definition: TGLUtil.h:856
Axes for non-linear projections.
void Set(Double_t x, Double_t y, Double_t z)
Definition: TGLUtil.h:213
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:501
TGLViewer::ECameraType camera
TMVec_t & RefTMVec()
TGLVector3 & RefTMOff(Int_t i)
const char * Data() const
Definition: TString.h:349
std::pair< Float_t, Float_t > Lab_t
std::vector< Lab_t > LabVec_t
TEveProjection * fProjection
Bool_t Highlight() const
Definition: TGLRnrCtx.h:218
void SplitInterval(Float_t x1, Float_t x2, Int_t axis) const
Build an array of tick-mark position-value pairs.
int d
Definition: tornado.py:11
static double p2(double t, double a, double b, double c)
Base-class for direct OpenGL renderers.
Definition: TGLObject.h:21
void ProjectVector(TEveVector &v, Float_t d)
Project TEveVector.
void FilterOverlappingLabels(Int_t idx, Float_t ref) const
Filter overlapping labels.
virtual void DirectDraw(TGLRnrCtx &rnrCtx) const
Actual rendering code.
TObject * fExternalObj
first replica
TThread * t[5]
Definition: threadsh1.C:13
virtual Bool_t SetModel(TObject *obj, const Option_t *opt=0)
Set model object.
void SplitIntervalByPos(Float_t min, Float_t max, Int_t axis) const
Add tick-marks at equidistant position.
Bool_t GetDrawCenter() const
void RnrLabels() const
Render label reading prepared list ov value-pos pairs.
ROOT::R::TRInterface & r
Definition: Object.C:4
SVector< double, 2 > v
Definition: Dict.h:5
void RnrLines() const
Render axis main line and tick-marks.
virtual Float_t GetValForScreenPos(Int_t ax, Float_t value)
Inverse projection.
ClassImp(TEveProjectionAxesGL)
virtual void Draw(TGLRnrCtx &rnrCtx) const
Draw the GL drawable, using draw flags.
bool first
Definition: line3Dfit.C:48
TLine * l
Definition: textangle.C:4
Bool_t IsDrawPassOutlineLine() const
Definition: TGLRnrCtx.h:207
static double p1(double t, double a, double b)
void Warning(const char *location, const char *msgfmt,...)
void SetLabelAlign(TGLFont::ETextAlignH_e, TGLFont::ETextAlignV_e)
Set label align.
TEveProjection * GetProjection()
virtual Float_t * GetProjectedCenter()
Get projected center.
Float_t * GetBBox()
Definition: TAttBBox.h:46
double Double_t
Definition: RtypesCore.h:55
Bool_t GetDrawOrigin() const
virtual Style_t GetLabelFont() const
Definition: TAttAxis.h:53
virtual Int_t GetNdivisions() const
Definition: TAttAxis.h:50
void Negate()
Definition: TGLUtil.h:144
static Vc_ALWAYS_INLINE int_v max(const int_v &x, const int_v &y)
Definition: vector.h:440
TEveVector GetOrthogonalCenter(int idx, TEveVector &out)
Get center ortogonal to given axis index.
LabVec_t & RefLabVec()
Mother of all ROOT objects.
Definition: TObject.h:58
virtual Float_t GetLabelSize() const
Definition: TAttAxis.h:55
TEveProjectionManager * GetManager()
TEveProjectionAxesGL()
Constructor.
Exception class thrown by TEve classes and macros.
Definition: TEveUtil.h:102
const TGLPlane & FrustumPlane(EFrustumPlane plane) const
Definition: TGLCamera.h:219
void SetLabelFont(TGLRnrCtx &rnrCtx, const char *fontName, Int_t pixelSize=64, Double_t font3DSize=-1)
Set label font derived from TAttAxis.
static Float_t LineWidth()
Get the line-width, taking the global scaling into account.
Definition: TGLUtil.cxx:1904
Bool_t fDLCache
display-list validity bit-field
Double_t Sqrt(Double_t x)
Definition: TMath.h:464
void SetTextFormat(Double_t min, Double_t max, Double_t binWidth)
Construct print format from given primary bin width.
TGLColorSet & ColorSet()
Return reference to current color-set (top of the stack).
Definition: TGLRnrCtx.cxx:278
void SplitIntervalByVal(Float_t min, Float_t max, Int_t axis) const
Add tick-marks on fixed value step.
const Bool_t kTRUE
Definition: Rtypes.h:91
ELabMode GetLabMode() const
TObject * obj
OpenGL renderer class for TEveProjectionAxes.
Int_t CeilNint(Double_t x)
Definition: TMath.h:470
const char * cnt
Definition: TXMLSetup.cxx:75
void SetAttAxis(TAttAxis *a)