Logo ROOT  
Reference Guide
TGraph2DPainter.cxx
Go to the documentation of this file.
1 // @(#)root/histpainter:$Id: TGraph2DPainter.cxx,v 1.00
2 // Author: Olivier Couet
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2000, 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 "TGraph2DPainter.h"
13 #include "TMath.h"
14 #include "TList.h"
15 #include "TGraph.h"
16 #include "TGraph2D.h"
17 #include "TGraphDelaunay.h"
18 #include "TGraphDelaunay2D.h"
19 #include "TPolyLine.h"
20 #include "TPolyMarker.h"
21 #include "TVirtualPad.h"
22 #include "TView.h"
23 #include "THLimitsFinder.h"
24 #include "TStyle.h"
25 #include "Hoption.h"
26 #include "TH1.h"
27 
30 
32 
33 
34 ////////////////////////////////////////////////////////////////////////////////
35 /*! \class TGraph2DPainter
36  \ingroup Histpainter
37  \brief The TGraphDelaunay painting class.
38 
39 TGraph2DPainter paints a TGraphDelaunay using triangles or clouds of points.
40 
41 */
42 
43 
44 ////////////////////////////////////////////////////////////////////////////////
45 /// TGraph2DPainter default constructor
46 
48 {
49  fX = 0;
50  fY = 0;
51  fZ = 0;
52  fEX = 0;
53  fEY = 0;
54  fEZ = 0;
55  fXN = 0;
56  fYN = 0;
57  fPTried = 0;
58  fNTried = 0;
59  fMTried = 0;
60  fGraph2D = 0;
61  fDelaunay = 0;
62  fDelaunay2D = 0;
63  fXmin = 0.;
64  fXmax = 0.;
65  fYmin = 0.;
66  fYmax = 0.;
67  fZmin = 0.;
68  fZmax = 0.;
69  fXNmin = 0;
70  fXNmax = 0;
71  fYNmin = 0;
72  fYNmax = 0;
73  fNdt = 0;
74  fNpoints = 0;
75 }
76 
77 
78 ////////////////////////////////////////////////////////////////////////////////
79 /// TGraph2DPainter constructor using the old Delaunay algorithm
80 
82 {
83  fDelaunay = gd;
84  fDelaunay2D = 0;
86  fNpoints = fGraph2D->GetN();
87  fX = fGraph2D->GetX();
88  fY = fGraph2D->GetY();
89  fZ = fGraph2D->GetZ();
90  fEX = fGraph2D->GetEX();
91  fEY = fGraph2D->GetEY();
92  fEZ = fGraph2D->GetEZ();
93  fNdt = 0;
94  fXN = 0;
95  fYN = 0;
96  fXNmin = 0;
97  fXNmax = 0;
98  fYNmin = 0;
99  fYNmax = 0;
100  fPTried = 0;
101  fNTried = 0;
102  fMTried = 0;
103  fXmin = 0.;
104  fXmax = 0.;
105  fYmin = 0.;
106  fYmax = 0.;
107  fZmin = 0.;
108  fZmax = 0.;
109 }
110 
111 ////////////////////////////////////////////////////////////////////////////////
112 /// TGraph2DPainter constructor using the new Delaunay algorithm
113 
115 {
116  fDelaunay = 0;
117  fDelaunay2D = gd;
119  fNpoints = fGraph2D->GetN();
120  fX = fGraph2D->GetX();
121  fY = fGraph2D->GetY();
122  fZ = fGraph2D->GetZ();
123  fEX = fGraph2D->GetEX();
124  fEY = fGraph2D->GetEY();
125  fEZ = fGraph2D->GetEZ();
126  fNdt = 0;
127  fXN = 0;
128  fYN = 0;
129  fXNmin = 0;
130  fXNmax = 0;
131  fYNmin = 0;
132  fYNmax = 0;
133  fPTried = 0;
134  fNTried = 0;
135  fMTried = 0;
136  fXmin = 0.;
137  fXmax = 0.;
138  fYmin = 0.;
139  fYmax = 0.;
140  fZmin = 0.;
141  fZmax = 0.;
142 }
143 
144 
145 ////////////////////////////////////////////////////////////////////////////////
146 /// TGraph2DPainter destructor.
147 
149 {
150 }
151 
152 
153 ////////////////////////////////////////////////////////////////////////////////
154 /// Find triangles in fDelaunay and initialise the TGraph2DPainter values
155 /// needed to paint triangles or find contours.
156 
158 {
159  if (fDelaunay) {
161  fNdt = fDelaunay->GetNdt();
162  fXN = fDelaunay->GetXN();
163  fYN = fDelaunay->GetYN();
171  }
172  else if (fDelaunay2D) {
174  fNdt = fDelaunay2D->GetNdt();
179  }
180 }
181 
182 
183 ////////////////////////////////////////////////////////////////////////////////
184 /// Returns the X and Y graphs building a contour. A contour level may
185 /// consist in several parts not connected to each other. This function
186 /// finds them and returns them in a graphs' list.
187 
189 {
190  // Exit if the contour is outside the Z range.
191  Double_t zmin = gCurrentHist->GetMinimum();
192  Double_t zmax = gCurrentHist->GetMaximum();
193  if (Hoption.Logz) {
194  if (zmin > 0) {
195  zmin = TMath::Log10(zmin);
196  zmax = TMath::Log10(zmax);
197  } else {
198  return 0;
199  }
200  }
201  if(contour<zmin || contour>zmax) {
202  Error("GetContourList", "Contour level (%g) outside the Z scope [%g,%g]",
203  contour,zmin,zmax);
204  return 0;
205  }
206 
207  if (!fNdt) FindTriangles();
208 
209  TGraph *graph = nullptr; // current graph
210  Int_t npg = 0; // number of points in the current graph
211 
212  // Find all the segments making the contour
213 
214  Double_t r21, r20, r10;
215  Int_t p0, p1, p2;
216  Double_t x0, y0, z0;
217  Double_t x1, y1, z1;
218  Double_t x2, y2, z2;
219  Int_t t[3],i,it,i0,i1,i2;
220 
221  // Allocate space to store the segments. They cannot be more than the
222  // number of triangles.
223  Double_t xs0c, ys0c, xs1c, ys1c;
224  Double_t *xs0 = new Double_t[fNdt];
225  Double_t *ys0 = new Double_t[fNdt];
226  Double_t *xs1 = new Double_t[fNdt];
227  Double_t *ys1 = new Double_t[fNdt];
228  for (i=0;i<fNdt;i++) {
229  xs0[i] = 0.;
230  ys0[i] = 0.;
231  xs1[i] = 0.;
232  ys1[i] = 0.;
233  }
234  Int_t nbSeg = 0;
235 
236  // Loop over all the triangles in order to find all the line segments
237  // making the contour.
238 
239  // old implementation
240  if (fDelaunay) {
241  for(it=0; it<fNdt; it++) {
242  t[0] = fPTried[it];
243  t[1] = fNTried[it];
244  t[2] = fMTried[it];
245  p0 = t[0]-1;
246  p1 = t[1]-1;
247  p2 = t[2]-1;
248  x0 = fX[p0]; x2 = fX[p0];
249  y0 = fY[p0]; y2 = fY[p0];
250  z0 = fZ[p0]; z2 = fZ[p0];
251 
252  // Order along Z axis the points (xi,yi,zi) where "i" belongs to {0,1,2}
253  // After this z0 < z1 < z2
254  i0 = i1 = i2 = 0;
255  if (fZ[p1]<=z0) {z0=fZ[p1]; x0=fX[p1]; y0=fY[p1]; i0=1;}
256  if (fZ[p1]>z2) {z2=fZ[p1]; x2=fX[p1]; y2=fY[p1]; i2=1;}
257  if (fZ[p2]<=z0) {z0=fZ[p2]; x0=fX[p2]; y0=fY[p2]; i0=2;}
258  if (fZ[p2]>z2) {z2=fZ[p2]; x2=fX[p2]; y2=fY[p2]; i2=2;}
259  if (i0==0 && i2==0) {
260  Error("GetContourList", "wrong vertices ordering");
261  delete [] xs0;
262  delete [] ys0;
263  delete [] xs1;
264  delete [] ys1;
265  return nullptr;
266  } else {
267  i1 = 3-i2-i0;
268  }
269  x1 = fX[t[i1]-1];
270  y1 = fY[t[i1]-1];
271  z1 = fZ[t[i1]-1];
272 
273  if (Hoption.Logz) {
274  z0 = TMath::Log10(z0);
275  z1 = TMath::Log10(z1);
276  z2 = TMath::Log10(z2);
277  }
278 
279  if(contour >= z0 && contour <=z2) {
280  r20 = (contour-z0)/(z2-z0);
281  xs0c = r20*(x2-x0)+x0;
282  ys0c = r20*(y2-y0)+y0;
283  if(contour >= z1 && contour <=z2) {
284  r21 = (contour-z1)/(z2-z1);
285  xs1c = r21*(x2-x1)+x1;
286  ys1c = r21*(y2-y1)+y1;
287  } else {
288  r10 = (contour-z0)/(z1-z0);
289  xs1c = r10*(x1-x0)+x0;
290  ys1c = r10*(y1-y0)+y0;
291  }
292  // do not take the segments equal to a point
293  if(xs0c != xs1c || ys0c != ys1c) {
294  nbSeg++;
295  xs0[nbSeg-1] = xs0c;
296  ys0[nbSeg-1] = ys0c;
297  xs1[nbSeg-1] = xs1c;
298  ys1[nbSeg-1] = ys1c;
299  }
300  }
301  }
302  }
303  else if (fDelaunay2D) {
304 
305  Int_t p[3];
306  for(const auto & face : *fDelaunay2D) {
307  p[0] = face.idx[0];
308  p[1] = face.idx[1];
309  p[2] = face.idx[2];
310  x0 = fX[p[0]]; x2 = fX[p[0]];
311  y0 = fY[p[0]]; y2 = fY[p[0]];
312  z0 = fZ[p[0]]; z2 = fZ[p[0]];
313 
314  // Order along Z axis the points (xi,yi,zi) where "i" belongs to {0,1,2}
315  // After this z0 < z1 < z2
316  i0 = i1 = i2 = 0;
317  if (fZ[p[1]]<=z0) {z0=fZ[p[1]]; x0=fX[p[1]]; y0=fY[p[1]]; i0=1;}
318  if (fZ[p[1]]>z2) {z2=fZ[p[1]]; x2=fX[p[1]]; y2=fY[p[1]]; i2=1;}
319  if (fZ[p[2]]<=z0) {z0=fZ[p[2]]; x0=fX[p[2]]; y0=fY[p[2]]; i0=2;}
320  if (fZ[p[2]]>z2) {z2=fZ[p[2]]; x2=fX[p[2]]; y2=fY[p[2]]; i2=2;}
321  if (i0==0 && i2==0) {
322  Error("GetContourList", "wrong vertices ordering");
323  delete [] xs0;
324  delete [] ys0;
325  delete [] xs1;
326  delete [] ys1;
327  return nullptr;
328  } else {
329  i1 = 3-i2-i0;
330  }
331  x1 = fX[p[i1]];
332  y1 = fY[p[i1]];
333  z1 = fZ[p[i1]];
334 
335  if (Hoption.Logz) {
336  z0 = TMath::Log10(z0);
337  z1 = TMath::Log10(z1);
338  z2 = TMath::Log10(z2);
339  }
340 
341  if(contour >= z0 && contour <=z2) {
342  r20 = (contour-z0)/(z2-z0);
343  xs0c = r20*(x2-x0)+x0;
344  ys0c = r20*(y2-y0)+y0;
345  if(contour >= z1 && contour <=z2) {
346  r21 = (contour-z1)/(z2-z1);
347  xs1c = r21*(x2-x1)+x1;
348  ys1c = r21*(y2-y1)+y1;
349  } else {
350  r10 = (contour-z0)/(z1-z0);
351  xs1c = r10*(x1-x0)+x0;
352  ys1c = r10*(y1-y0)+y0;
353  }
354  // do not take the segments equal to a point
355  if(xs0c != xs1c || ys0c != ys1c) {
356  nbSeg++;
357  xs0[nbSeg-1] = xs0c;
358  ys0[nbSeg-1] = ys0c;
359  xs1[nbSeg-1] = xs1c;
360  ys1[nbSeg-1] = ys1c;
361  }
362  }
363  }
364  }
365 
366  TList *list = new TList(); // list holding all the graphs
367 
368  Int_t *segUsed = new Int_t[fNdt];
369  for(i=0; i<fNdt; i++) segUsed[i]=kFALSE;
370 
371  // Find all the graphs making the contour. There is two kind of graphs,
372  // either they are "opened" or they are "closed"
373 
374  // Find the opened graphs
375  Double_t xc=0, yc=0, xnc=0, ync=0;
376  Bool_t findNew;
377  Bool_t s0, s1;
378  Int_t is, js;
379  for (is=0; is<nbSeg; is++) {
380  if (segUsed[is]) continue;
381  s0 = s1 = kFALSE;
382 
383  // Find to which segment is is connected. It can be connected
384  // via 0, 1 or 2 vertices.
385  for (js=0; js<nbSeg; js++) {
386  if (is==js) continue;
387  if (xs0[is]==xs0[js] && ys0[is]==ys0[js]) s0 = kTRUE;
388  if (xs0[is]==xs1[js] && ys0[is]==ys1[js]) s0 = kTRUE;
389  if (xs1[is]==xs0[js] && ys1[is]==ys0[js]) s1 = kTRUE;
390  if (xs1[is]==xs1[js] && ys1[is]==ys1[js]) s1 = kTRUE;
391  }
392 
393  // Segment is is alone, not connected. It is stored in the
394  // list and the next segment is examined.
395  if (!s0 && !s1) {
396  graph = new TGraph();
397  graph->SetPoint(npg,xs0[is],ys0[is]); npg++;
398  graph->SetPoint(npg,xs1[is],ys1[is]); npg++;
399  segUsed[is] = kTRUE;
400  list->Add(graph); npg = 0;
401  continue;
402  }
403 
404  // Segment is is connected via 1 vertex only and can be considered
405  // as the starting point of an opened contour.
406  if (!s0 || !s1) {
407  // Find all the segments connected to segment is
408  graph = new TGraph();
409  if (s0) {xc = xs0[is]; yc = ys0[is]; xnc = xs1[is]; ync = ys1[is];}
410  if (s1) {xc = xs1[is]; yc = ys1[is]; xnc = xs0[is]; ync = ys0[is];}
411  graph->SetPoint(npg,xnc,ync); npg++;
412  segUsed[is] = kTRUE;
413  js = 0;
414 L01:
415  findNew = kFALSE;
416  if (js < nbSeg && segUsed[js]) {
417  js++;
418  goto L01;
419  } else if (xc==xs0[js] && yc==ys0[js]) {
420  xc = xs1[js];
421  yc = ys1[js];
422  findNew = kTRUE;
423  } else if (xc==xs1[js] && yc==ys1[js]) {
424  xc = xs0[js];
425  yc = ys0[js];
426  findNew = kTRUE;
427  }
428  if (findNew) {
429  segUsed[js] = kTRUE;
430  graph->SetPoint(npg,xc,yc); npg++;
431  js = 0;
432  goto L01;
433  }
434  js++;
435  if (js<nbSeg) goto L01;
436  list->Add(graph); npg = 0;
437  }
438  }
439 
440 
441  // Find the closed graphs. At this point all the remaining graphs
442  // are closed. Any segment can be used to start the search.
443  for (is=0; is<nbSeg; is++) {
444  if (segUsed[is]) continue;
445 
446  // Find all the segments connected to segment is
447  graph = new TGraph();
448  segUsed[is] = kTRUE;
449  xc = xs0[is];
450  yc = ys0[is];
451  js = 0;
452  graph->SetPoint(npg,xc,yc); npg++;
453 L02:
454  findNew = kFALSE;
455  if (js < nbSeg && segUsed[js]) {
456  js++;
457  goto L02;
458  } else if (xc==xs0[js] && yc==ys0[js]) {
459  xc = xs1[js];
460  yc = ys1[js];
461  findNew = kTRUE;
462  } else if (xc==xs1[js] && yc==ys1[js]) {
463  xc = xs0[js];
464  yc = ys0[js];
465  findNew = kTRUE;
466  }
467  if (findNew) {
468  segUsed[js] = kTRUE;
469  graph->SetPoint(npg,xc,yc); npg++;
470  js = 0;
471  goto L02;
472  }
473  js++;
474  if (js<nbSeg) goto L02;
475  // Close the contour
476  graph->SetPoint(npg,xs0[is],ys0[is]); npg++;
477  list->Add(graph); npg = 0;
478  }
479 
480  delete [] xs0;
481  delete [] ys0;
482  delete [] xs1;
483  delete [] ys1;
484  delete [] segUsed;
485  return list;
486 }
487 
488 
489 ////////////////////////////////////////////////////////////////////////////////
490 /// Paint a TGraphDelaunay according to the value of "option":
491 ///
492 /// | Option | Description |
493 /// |----------|---------------------------------------------------------------|
494 /// | "TRI" | The Delaunay triangles are drawn using filled area. An hidden surface drawing technique is used. The surface is painted with the current fill area color. The edges of each triangles are painted with the current line color. |
495 /// | "TRIW" | The Delaunay triangles are drawn as wire frame. |
496 /// | "TRI1" | The Delaunay triangles are painted with color levels. The edges of each triangles are painted with the current line color. |
497 /// | "TRI2" | the Delaunay triangles are painted with color levels. |
498 /// | "P" | Draw a marker at each vertex. |
499 /// | "P0" | Draw a circle at each vertex. Each circle background is white. |
500 /// | "PCOL" | Draw a marker at each vertex. The color of each marker is defined according to its Z position. |
501 /// | "CONT" | Draw contours. |
502 /// | "LINE" | Draw a 3D polyline. |
503 
505 {
506  TString opt = option;
507  opt.ToLower();
508  Bool_t triangles = opt.Contains("tri") ||
509  opt.Contains("tri1") ||
510  opt.Contains("tri2");
511  if (opt.Contains("tri0")) triangles = kFALSE;
512 
513  Bool_t markers = opt.Contains("p") && !triangles;
514  Bool_t contour = opt.Contains("cont");
515  Bool_t line = opt.Contains("line");
516  Bool_t err = opt.Contains("err");
517 
518  fGraph2D->TAttLine::Modify();
519  fGraph2D->TAttFill::Modify();
520  fGraph2D->TAttMarker::Modify();
521 
522  // Compute minimums and maximums
523  TAxis *xaxis = gCurrentHist->GetXaxis();
524  Int_t first = xaxis->GetFirst();
525  fXmin = xaxis->GetBinLowEdge(first);
526  if (Hoption.Logx && fXmin <= 0) fXmin = xaxis->GetBinUpEdge(xaxis->FindFixBin(0.01*xaxis->GetBinWidth(first)));
527  fXmax = xaxis->GetBinUpEdge(xaxis->GetLast());
528  TAxis *yaxis = gCurrentHist->GetYaxis();
529  first = yaxis->GetFirst();
530  fYmin = yaxis->GetBinLowEdge(first);
531  if (Hoption.Logy && fYmin <= 0) fYmin = yaxis->GetBinUpEdge(yaxis->FindFixBin(0.01*yaxis->GetBinWidth(first)));
532  fYmax = yaxis->GetBinUpEdge(yaxis->GetLast());
533  fZmax = fGraph2D->GetZmaxE();
534  fZmin = fGraph2D->GetZminE();
535  if (Hoption.Logz && fZmin <= 0) fZmin = TMath::Min((Double_t)1, (Double_t)0.001*fGraph2D->GetZmax());
536 
537  if (triangles) PaintTriangles(option);
538  if (contour) PaintContour(option);
539  if (line) PaintPolyLine(option);
540  if (err) PaintErrors(option);
541  if (markers) PaintPolyMarker(option);
542 }
543 
544 
545 ////////////////////////////////////////////////////////////////////////////////
546 /// Paints the 2D graph as a contour plot. Delaunay triangles are used
547 /// to compute the contours.
548 
550 {
551  // Initialize the levels on the Z axis
552  Int_t ncolors = gStyle->GetNumberOfColors();
553  Int_t ndiv = gCurrentHist->GetContour();
554  if (ndiv == 0 ) {
555  ndiv = gStyle->GetNumberContours();
556  gCurrentHist->SetContour(ndiv);
557  }
558  Int_t ndivz = TMath::Abs(ndiv);
560 
561  Int_t theColor;
562  TList *l;
563  TGraph *g;
564  TObject *obj;
565  Double_t c;
566 
567  if (!fNdt) FindTriangles();
568 
569  for (Int_t k=0; k<ndiv; k++) {
571  l = GetContourList(c);
572  TIter next(l);
573  while ((obj = next())) {
574  if(obj->InheritsFrom(TGraph::Class()) ) {
575  g=(TGraph*)obj;
576  g->SetLineWidth(fGraph2D->GetLineWidth());
577  g->SetLineStyle(fGraph2D->GetLineStyle());
578  theColor = Int_t((k+0.99)*Float_t(ncolors)/Float_t(ndivz));
579  g->SetLineColor(gStyle->GetColorPalette(theColor));
580  g->Paint("l");
581  }
582  }
583  if (l) { l->Delete(); delete l; }
584  }
585 }
586 
587 
588 ////////////////////////////////////////////////////////////////////////////////
589 /// Paints the 2D graph as error bars
590 
592 {
593  Double_t temp1[3],temp2[3];
594 
595  TView *view = gPad->GetView();
596  if (!view) {
597  Error("PaintErrors", "No TView in current pad");
598  return;
599  }
600 
601  Int_t it;
602 
603  Double_t *xm = new Double_t[2];
604  Double_t *ym = new Double_t[2];
605 
606  fGraph2D->SetLineStyle(fGraph2D->GetLineStyle());
607  fGraph2D->SetLineWidth(fGraph2D->GetLineWidth());
608  fGraph2D->SetLineColor(fGraph2D->GetLineColor());
609  fGraph2D->TAttLine::Modify();
610 
611  for (it=0; it<fNpoints; it++) {
612  if(fX[it] < fXmin || fX[it] > fXmax) continue;
613  if(fY[it] < fYmin || fY[it] > fYmax) continue;
614  if (fEX) {
615  temp1[0] = fX[it]-fEX[it];
616  temp1[1] = fY[it];
617  temp1[2] = fZ[it];
618  temp1[0] = TMath::Max(temp1[0],fXmin);
619  temp1[1] = TMath::Max(temp1[1],fYmin);
620  temp1[2] = TMath::Max(temp1[2],fZmin);
621  temp1[2] = TMath::Min(temp1[2],fZmax);
622  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
623  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
624  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
625  view->WCtoNDC(temp1, &temp2[0]);
626  xm[0] = temp2[0];
627  ym[0] = temp2[1];
628 
629  temp1[0] = fX[it]+fEX[it];
630  temp1[0] = TMath::Max(temp1[0],fXmin);
631  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
632  view->WCtoNDC(temp1, &temp2[0]);
633  xm[1] = temp2[0];
634  ym[1] = temp2[1];
635  gPad->PaintPolyLine(2,xm,ym);
636  }
637  if (fEY) {
638  temp1[0] = fX[it];
639  temp1[1] = fY[it]-fEY[it];
640  temp1[2] = fZ[it];
641  temp1[0] = TMath::Max(temp1[0],fXmin);
642  temp1[1] = TMath::Max(temp1[1],fYmin);
643  temp1[2] = TMath::Max(temp1[2],fZmin);
644  temp1[2] = TMath::Min(temp1[2],fZmax);
645  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
646  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
647  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
648  view->WCtoNDC(temp1, &temp2[0]);
649  xm[0] = temp2[0];
650  ym[0] = temp2[1];
651 
652  temp1[1] = fY[it]+fEY[it];
653  temp1[1] = TMath::Max(temp1[1],fYmin);
654  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
655  view->WCtoNDC(temp1, &temp2[0]);
656  xm[1] = temp2[0];
657  ym[1] = temp2[1];
658  gPad->PaintPolyLine(2,xm,ym);
659  }
660  if (fEZ) {
661  temp1[0] = fX[it];
662  temp1[1] = fY[it];
663  temp1[2] = fZ[it]-fEZ[it];
664  temp1[0] = TMath::Max(temp1[0],fXmin);
665  temp1[1] = TMath::Max(temp1[1],fYmin);
666  temp1[2] = TMath::Max(temp1[2],fZmin);
667  temp1[2] = TMath::Min(temp1[2],fZmax);
668  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
669  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
670  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
671  view->WCtoNDC(temp1, &temp2[0]);
672  xm[0] = temp2[0];
673  ym[0] = temp2[1];
674 
675  temp1[2] = fZ[it]+fEZ[it];
676  temp1[2] = TMath::Max(temp1[2],fZmin);
677  temp1[2] = TMath::Min(temp1[2],fZmax);
678  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
679  view->WCtoNDC(temp1, &temp2[0]);
680  xm[1] = temp2[0];
681  ym[1] = temp2[1];
682  gPad->PaintPolyLine(2,xm,ym);
683  }
684  }
685  delete [] xm;
686  delete [] ym;
687 }
688 
689 
690 ////////////////////////////////////////////////////////////////////////////////
691 /// Paints one triangle.
692 ///
693 /// - nblev = 0 : paint the color levels
694 /// - nblev != 0 : paint the grid
695 
697  Int_t nblev, Double_t *glev)
698 {
699  Int_t i, fillColor, ncolors, theColor0, theColor2;
700 
701  Int_t p[3];
702  if (fDelaunay) {
703  p[0]=t[0]-1;
704  p[1]=t[1]-1;
705  p[2]=t[2]-1;
706  }
707  else {
708  p[0]=t[0];
709  p[1]=t[1];
710  p[2]=t[2];
711  }
712  Double_t xl[2],yl[2];
713  Double_t zl, r21, r20, r10;
714  Double_t x0 = x[0] , x2 = x[0];
715  Double_t y0 = y[0] , y2 = y[0];
716  Double_t z0 = fZ[p[0]], z2 = fZ[p[0]];
717  Double_t zmin = fGraph2D->GetMinimum();
718  Double_t zmax = fGraph2D->GetMaximum();
719  if (zmin==-1111 && zmax==-1111) {
720  zmin = TMath::Min(fZmin, 0.);
721  if (Hoption.Logz && zmin <= 0) zmin = TMath::Min((Double_t)1, (Double_t)0.001*fGraph2D->GetZmax());
722  zmax = fZmax;
723  }
724 
725  // Order along Z axis the points (xi,yi,zi) where "i" belongs to {0,1,2}
726  // After this z0 < z1 < z2
727  Int_t i0=0, i1=0, i2=0;
728  if (fZ[p[1]]<=z0) {z0=fZ[p[1]]; x0=x[1]; y0=y[1]; i0=1;}
729  if (fZ[p[1]]>z2) {z2=fZ[p[1]]; x2=x[1]; y2=y[1]; i2=1;}
730  if (fZ[p[2]]<=z0) {z0=fZ[p[2]]; x0=x[2]; y0=y[2]; i0=2;}
731  if (fZ[p[2]]>z2) {z2=fZ[p[2]]; x2=x[2]; y2=y[2]; i2=2;}
732  i1 = 3-i2-i0;
733  Double_t x1 = x[i1];
734  Double_t y1 = y[i1];
735  Double_t z1 = fZ[p[i1]];
736 
737  if (z0>zmax) z0 = zmax;
738  if (z2>zmax) z2 = zmax;
739  if (z0<zmin) z0 = zmin;
740  if (z2<zmin) z2 = zmin;
741  if (z1>zmax) z1 = zmax;
742  if (z1<zmin) z1 = zmin;
743 
744  if (Hoption.Logz) {
745  z0 = TMath::Log10(z0);
746  z1 = TMath::Log10(z1);
747  z2 = TMath::Log10(z2);
748  zmin = TMath::Log10(zmin);
749  zmax = TMath::Log10(zmax);
750  }
751 
752  // zi = Z values of the stripe number i
753  // zip = Previous zi
754  Double_t zi=0, zip=0;
755 
756  if (nblev <= 0) {
757  // Paint the colors levels
758 
759  // Compute the color associated to z0 (theColor0) and z2 (theColor2)
760  ncolors = gStyle->GetNumberOfColors();
761  theColor0 = (Int_t)( ((z0-zmin)/(zmax-zmin))*(ncolors-1) );
762  theColor2 = (Int_t)( ((z2-zmin)/(zmax-zmin))*(ncolors-1) );
763 
764  // The stripes drawn to fill the triangles may have up to 5 points
765  Double_t xp[5], yp[5];
766 
767  // rl = Ratio between z0 and z2 (long)
768  // rs = Ratio between z0 and z1 or z1 and z2 (short)
769  Double_t rl,rs;
770 
771  // ci = Color of the stripe number i
772  // npf = number of point needed to draw the current stripe
773  Int_t ci,npf;
774 
775  fillColor = fGraph2D->GetFillColor();
776 
777  // If the z0's color and z2's colors are the same, the whole triangle
778  // can be painted in one go.
779  if(theColor0 == theColor2) {
780  fGraph2D->SetFillColor(gStyle->GetColorPalette(theColor0));
781  fGraph2D->TAttFill::Modify();
782  gPad->PaintFillArea(3,x,y);
783 
784  // The triangle must be painted with several colors
785  } else {
786  for(ci=theColor0; ci<=theColor2; ci++) {
787  fGraph2D->SetFillColor(gStyle->GetColorPalette(ci));
788  fGraph2D->TAttFill::Modify();
789  if (ci==theColor0) {
790  zi = (((ci+1)*(zmax-zmin))/(ncolors-1))+zmin;
791  xp[0] = x0;
792  yp[0] = y0;
793  rl = (zi-z0)/(z2-z0);
794  xp[1] = rl*(x2-x0)+x0;
795  yp[1] = rl*(y2-y0)+y0;
796  if (zi>=z1 || z0==z1) {
797  rs = (zi-z1)/(z2-z1);
798  xp[2] = rs*(x2-x1)+x1;
799  yp[2] = rs*(y2-y1)+y1;
800  xp[3] = x1;
801  yp[3] = y1;
802  npf = 4;
803  } else {
804  rs = (zi-z0)/(z1-z0);
805  xp[2] = rs*(x1-x0)+x0;
806  yp[2] = rs*(y1-y0)+y0;
807  npf = 3;
808  }
809  } else if (ci==theColor2) {
810  xp[0] = xp[1];
811  yp[0] = yp[1];
812  xp[1] = x2;
813  yp[1] = y2;
814  if (zi<z1 || z2==z1) {
815  xp[3] = xp[2];
816  yp[3] = yp[2];
817  xp[2] = x1;
818  yp[2] = y1;
819  npf = 4;
820  } else {
821  npf = 3;
822  }
823  } else {
824  zi = (((ci+1)*(zmax-zmin))/(ncolors-1))+zmin;
825  xp[0] = xp[1];
826  yp[0] = yp[1];
827  rl = (zi-z0)/(z2-z0);
828  xp[1] = rl*(x2-x0)+x0;
829  yp[1] = rl*(y2-y0)+y0;
830  if ( zi>=z1 && zip<=z1) {
831  xp[3] = x1;
832  yp[3] = y1;
833  xp[4] = xp[2];
834  yp[4] = yp[2];
835  npf = 5;
836  } else {
837  xp[3] = xp[2];
838  yp[3] = yp[2];
839  npf = 4;
840  }
841  if (zi<z1) {
842  rs = (zi-z0)/(z1-z0);
843  xp[2] = rs*(x1-x0)+x0;
844  yp[2] = rs*(y1-y0)+y0;
845  } else {
846  rs = (zi-z1)/(z2-z1);
847  xp[2] = rs*(x2-x1)+x1;
848  yp[2] = rs*(y2-y1)+y1;
849  }
850  }
851  zip = zi;
852  // Paint a stripe
853  gPad->PaintFillArea(npf,xp,yp);
854  }
855  }
856  fGraph2D->SetFillColor(fillColor);
857  fGraph2D->TAttFill::Modify();
858 
859  } else {
860  // Paint the grid levels
861  fGraph2D->SetLineStyle(3);
862  fGraph2D->TAttLine::Modify();
863  for(i=0; i<nblev; i++){
864  zl=glev[i];
865  if(zl >= z0 && zl <=z2) {
866  r21=(zl-z1)/(z2-z1);
867  r20=(zl-z0)/(z2-z0);
868  r10=(zl-z0)/(z1-z0);
869  xl[0]=r20*(x2-x0)+x0;
870  yl[0]=r20*(y2-y0)+y0;
871  if(zl >= z1 && zl <=z2) {
872  xl[1]=r21*(x2-x1)+x1;
873  yl[1]=r21*(y2-y1)+y1;
874  } else {
875  xl[1]=r10*(x1-x0)+x0;
876  yl[1]=r10*(y1-y0)+y0;
877  }
878  gPad->PaintPolyLine(2,xl,yl);
879  }
880  }
881  fGraph2D->SetLineStyle(1);
882  fGraph2D->TAttLine::Modify();
883  }
884 }
885 
886 
887 ////////////////////////////////////////////////////////////////////////////////
888 /// Paints the 2D graph as PaintPolyMarker
889 
891 {
892  Double_t temp1[3],temp2[3];
893 
894  TView *view = gPad->GetView();
895  if (!view) {
896  Error("PaintPolyMarker", "No TView in current pad");
897  return;
898  }
899 
900  TString opt = option;
901  opt.ToLower();
902  Bool_t markers0 = opt.Contains("p0");
903  Bool_t colors = opt.Contains("pcol");
904  Int_t ncolors = gStyle->GetNumberOfColors();
905  Int_t it, theColor;
906 
907  // Initialize the levels on the Z axis
908  if (colors) {
909  Int_t ndiv = gCurrentHist->GetContour();
910  if (ndiv == 0 ) {
911  ndiv = gStyle->GetNumberContours();
912  gCurrentHist->SetContour(ndiv);
913  }
915  }
916 
917  Double_t *xm = new Double_t[fNpoints];
918  Double_t *ym = new Double_t[fNpoints];
919  Double_t *zm = new Double_t[fNpoints];
920  Double_t hzmin = gCurrentHist->GetMinimum();
921  Double_t hzmax = gCurrentHist->GetMaximum();
922 
923  // min and max for colors
924  Double_t hzmincol = hzmin;
925  Double_t hzmaxcol = hzmax;
926  if (hzmincol==-1111 && hzmaxcol==-1111) {
927  hzmincol = TMath::Min(hzmincol, 0.);
928  if (Hoption.Logz && hzmincol <= 0) hzmincol = TMath::Min((Double_t)1, (Double_t)0.001*fGraph2D->GetZmax());
929  hzmaxcol = fZmax;
930  }
931  if (Hoption.Logz) {
932  hzmincol = TMath::Log10(hzmincol);
933  hzmaxcol = TMath::Log10(hzmaxcol);
934  }
935 
936  Double_t Xeps = (fXmax-fXmin)*0.0001;
937  Double_t Yeps = (fYmax-fYmin)*0.0001;
938  Double_t Zeps = (hzmax-hzmin)*0.0001;
939 
940  Int_t npd = 0;
941  for (it=0; it<fNpoints; it++) {
942  xm[it] = 0;
943  ym[it] = 0;
944  if(fXmin - fX[it] > Xeps || fX[it] - fXmax > Xeps) continue;
945  if(fYmin - fY[it] > Yeps || fY[it] - fYmax > Yeps) continue;
946  if(hzmin - fZ[it] > Zeps || fZ[it] - hzmax > Zeps) continue;
947  temp1[0] = fX[it];
948  temp1[1] = fY[it];
949  temp1[2] = fZ[it];
950  temp1[0] = TMath::Max(temp1[0],fXmin);
951  temp1[1] = TMath::Max(temp1[1],fYmin);
952  temp1[2] = TMath::Max(temp1[2],hzmin);
953  temp1[2] = TMath::Min(temp1[2],hzmax);
954  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
955  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
956  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
957  view->WCtoNDC(temp1, &temp2[0]);
958  xm[npd] = temp2[0];
959  ym[npd] = temp2[1];
960  zm[npd] = temp1[2];
961  npd++;
962  }
963  if (markers0) {
964  PaintPolyMarker0(npd,xm,ym);
965  } else if (colors) {
966  Int_t cols = fGraph2D->GetMarkerColor();
967  for (it=0; it<npd; it++) {
968  theColor = (Int_t)( ((zm[it]-hzmincol)/(hzmaxcol-hzmincol))*(ncolors-1) );
969  fGraph2D->SetMarkerColor(gStyle->GetColorPalette(theColor));
970  fGraph2D->TAttMarker::Modify();
971  gPad->PaintPolyMarker(1,&xm[it],&ym[it]);
972  }
973  fGraph2D->SetMarkerColor(cols);
974  } else {
975  fGraph2D->SetMarkerStyle(fGraph2D->GetMarkerStyle());
976  fGraph2D->SetMarkerSize(fGraph2D->GetMarkerSize());
977  fGraph2D->SetMarkerColor(fGraph2D->GetMarkerColor());
978  fGraph2D->TAttMarker::Modify();
979  gPad->PaintPolyMarker(npd,xm,ym);
980  }
981  delete [] xm;
982  delete [] ym;
983  delete [] zm;
984 }
985 
986 
987 ////////////////////////////////////////////////////////////////////////////////
988 /// Paints the 2D graph as PaintPolyLine
989 
991 {
992  Double_t temp1[3],temp2[3];
993 
994  TView *view = gPad->GetView();
995  if (!view) {
996  Error("PaintPolyLine", "No TView in current pad");
997  return;
998  }
999 
1000  Int_t it;
1001  Double_t Xeps = (fXmax - fXmin) * 0.0001;
1002  Double_t Yeps = (fYmax - fYmin) * 0.0001;
1003 
1004  Double_t *xm = new Double_t[fNpoints];
1005  Double_t *ym = new Double_t[fNpoints];
1006  Int_t npd = 0;
1007 
1008  for (it=0; it<fNpoints; it++) {
1009  if (fXmin - fX[it] > Xeps || fX[it] - fXmax > Xeps)
1010  continue;
1011  if (fYmin - fY[it] > Yeps || fY[it] - fYmax > Yeps)
1012  continue;
1013  npd++;
1014  temp1[0] = fX[it];
1015  temp1[1] = fY[it];
1016  temp1[2] = fZ[it];
1017  temp1[0] = TMath::Max(temp1[0],fXmin);
1018  temp1[1] = TMath::Max(temp1[1],fYmin);
1019  temp1[2] = TMath::Max(temp1[2],fZmin);
1020  temp1[2] = TMath::Min(temp1[2],fZmax);
1021  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
1022  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
1023  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
1024  view->WCtoNDC(temp1, &temp2[0]);
1025  xm[npd - 1] = temp2[0];
1026  ym[npd - 1] = temp2[1];
1027  }
1028  fGraph2D->SetLineStyle(fGraph2D->GetLineStyle());
1029  fGraph2D->SetLineWidth(fGraph2D->GetLineWidth());
1030  fGraph2D->SetLineColor(fGraph2D->GetLineColor());
1031  fGraph2D->TAttLine::Modify();
1032  gPad->PaintPolyLine(npd,xm,ym);
1033  delete [] xm;
1034  delete [] ym;
1035 }
1036 
1037 
1038 ////////////////////////////////////////////////////////////////////////////////
1039 /// Paints a circle at each vertex. Each circle background is white.
1040 
1042 {
1043  fGraph2D->SetMarkerSize(fGraph2D->GetMarkerSize());
1044  Int_t mc = fGraph2D->GetMarkerColor();
1045  Int_t ms = fGraph2D->GetMarkerStyle();
1046  for (Int_t i=0; i<n; i++) {
1047  fGraph2D->SetMarkerStyle(20);
1048  fGraph2D->SetMarkerColor(0);
1049  fGraph2D->TAttMarker::Modify();
1050  gPad->PaintPolyMarker(1,&x[i],&y[i]);
1051  fGraph2D->SetMarkerStyle(24);
1052  fGraph2D->SetMarkerColor(mc);
1053  fGraph2D->TAttMarker::Modify();
1054  gPad->PaintPolyMarker(1,&x[i],&y[i]);
1055  }
1056  fGraph2D->SetMarkerStyle(ms);
1057 }
1058 
1059 
1060 ////////////////////////////////////////////////////////////////////////////////
1061 /// Paints the 2D graph as triangles
1062 
1064 {
1065  if (fDelaunay)
1066  PaintTriangles_old(option);
1067  else if (fDelaunay2D)
1068  PaintTriangles_new(option);
1069 }
1070 
1071 
1072 ////////////////////////////////////////////////////////////////////////////////
1073 /// Paints the 2D graph as triangles (old implementation)
1074 
1076 {
1077  Double_t x[4], y[4], temp1[3],temp2[3];
1078  Int_t it,t[3];
1079  Int_t *order = 0;
1080  Double_t *dist = 0;
1081 
1082  TView *view = gPad->GetView();
1083  if (!view) {
1084  Error("PaintTriangles", "No TView in current pad");
1085  return;
1086  }
1087 
1088  TString opt = option;
1089  opt.ToLower();
1090  Bool_t tri1 = opt.Contains("tri1");
1091  Bool_t tri2 = opt.Contains("tri2");
1092  Bool_t markers = opt.Contains("p");
1093  Bool_t markers0 = opt.Contains("p0");
1094  Bool_t wire = opt.Contains("w");
1095 
1096  // Define the grid levels drawn on the triangles.
1097  // The grid levels are aligned on the Z axis' main tick marks.
1098  Int_t nblev=0;
1099  Double_t *glev=0;
1100  if (!tri1 && !tri2 && !wire) {
1101  Int_t ndivz = gCurrentHist->GetZaxis()->GetNdivisions()%100;
1102  Int_t nbins;
1103  Double_t binLow = 0, binHigh = 0, binWidth = 0;
1104 
1105  // Find the main tick marks positions.
1106  Double_t *r0 = view->GetRmin();
1107  Double_t *r1 = view->GetRmax();
1108  if (!r0 || !r1) return;
1109 
1110  if (ndivz > 0) {
1111  THLimitsFinder::Optimize(r0[2], r1[2], ndivz,
1112  binLow, binHigh, nbins, binWidth, " ");
1113  } else {
1114  nbins = TMath::Abs(ndivz);
1115  binLow = r0[2];
1116  binHigh = r1[2];
1117  binWidth = (binHigh-binLow)/nbins;
1118  }
1119  // Define the grid levels
1120  nblev = nbins+1;
1121  glev = new Double_t[nblev];
1122  for (Int_t i = 0; i < nblev; ++i) glev[i] = binLow+i*binWidth;
1123  }
1124 
1125  // Initialize the levels on the Z axis
1126  if (tri1 || tri2) {
1127  Int_t ndiv = gCurrentHist->GetContour();
1128  if (ndiv == 0 ) {
1129  ndiv = gStyle->GetNumberContours();
1130  gCurrentHist->SetContour(ndiv);
1131  }
1133  }
1134 
1135  // For each triangle, compute the distance between the triangle centre
1136  // and the back planes. Then these distances are sorted in order to draw
1137  // the triangles from back to front.
1138  if (!fNdt) FindTriangles();
1139  Double_t cp = TMath::Cos(view->GetLongitude()*TMath::Pi()/180.);
1140  Double_t sp = TMath::Sin(view->GetLongitude()*TMath::Pi()/180.);
1141  order = new Int_t[fNdt];
1142  dist = new Double_t[fNdt];
1143  Double_t xd,yd;
1144  Int_t p, n, m;
1145  Bool_t o = kFALSE;
1146  for (it=0; it<fNdt; it++) {
1147  p = fPTried[it];
1148  n = fNTried[it];
1149  m = fMTried[it];
1150  xd = (fXN[p]+fXN[n]+fXN[m])/3;
1151  yd = (fYN[p]+fYN[n]+fYN[m])/3;
1152  if ((cp >= 0) && (sp >= 0.)) {
1153  dist[it] = -(fXNmax-xd+fYNmax-yd);
1154  } else if ((cp <= 0) && (sp >= 0.)) {
1155  dist[it] = -(fXNmax-xd+yd-fYNmin);
1156  o = kTRUE;
1157  } else if ((cp <= 0) && (sp <= 0.)) {
1158  dist[it] = -(xd-fXNmin+yd-fYNmin);
1159  } else {
1160  dist[it] = -(xd-fXNmin+fYNmax-yd);
1161  o = kTRUE;
1162  }
1163  }
1164  TMath::Sort(fNdt, dist, order, o);
1165 
1166  // Draw the triangles and markers if requested
1167  fGraph2D->SetFillColor(fGraph2D->GetFillColor());
1168  Int_t fs = fGraph2D->GetFillStyle();
1169  fGraph2D->SetFillStyle(1001);
1170  fGraph2D->TAttFill::Modify();
1171  fGraph2D->SetLineColor(fGraph2D->GetLineColor());
1172  fGraph2D->TAttLine::Modify();
1173  int lst = fGraph2D->GetLineStyle();
1174  for (it=0; it<fNdt; it++) {
1175  t[0] = fPTried[order[it]];
1176  t[1] = fNTried[order[it]];
1177  t[2] = fMTried[order[it]];
1178  for (Int_t k=0; k<3; k++) {
1179  if(fX[t[k]-1] < fXmin || fX[t[k]-1] > fXmax) goto endloop;
1180  if(fY[t[k]-1] < fYmin || fY[t[k]-1] > fYmax) goto endloop;
1181  temp1[0] = fX[t[k]-1];
1182  temp1[1] = fY[t[k]-1];
1183  temp1[2] = fZ[t[k]-1];
1184  temp1[0] = TMath::Max(temp1[0],fXmin);
1185  temp1[1] = TMath::Max(temp1[1],fYmin);
1186  temp1[2] = TMath::Max(temp1[2],fZmin);
1187  temp1[2] = TMath::Min(temp1[2],fZmax);
1188  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
1189  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
1190  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
1191  view->WCtoNDC(temp1, &temp2[0]);
1192  x[k] = temp2[0];
1193  y[k] = temp2[1];
1194  }
1195  x[3] = x[0];
1196  y[3] = y[0];
1197  if (tri1 || tri2) PaintLevels(t,x,y);
1198  if (!tri1 && !tri2 && !wire) {
1199  gPad->PaintFillArea(3,x,y);
1200  PaintLevels(t,x,y,nblev,glev);
1201  }
1202  if (!tri2) gPad->PaintPolyLine(4,x,y);
1203  if (markers) {
1204  if (markers0) {
1205  PaintPolyMarker0(3,x,y);
1206  } else {
1207  fGraph2D->SetMarkerStyle(fGraph2D->GetMarkerStyle());
1208  fGraph2D->SetMarkerSize(fGraph2D->GetMarkerSize());
1209  fGraph2D->SetMarkerColor(fGraph2D->GetMarkerColor());
1210  fGraph2D->TAttMarker::Modify();
1211  gPad->PaintPolyMarker(3,x,y);
1212  }
1213  }
1214 endloop:
1215  continue;
1216  }
1217  fGraph2D->SetFillStyle(fs);
1218  fGraph2D->SetLineStyle(lst);
1219  fGraph2D->TAttLine::Modify();
1220  fGraph2D->TAttFill::Modify();
1221  delete [] order;
1222  delete [] dist;
1223  if (glev) delete [] glev;
1224 }
1225 
1226 
1227 ////////////////////////////////////////////////////////////////////////////////
1228 /// Paints the 2D graph as triangles (new implementation)
1229 
1231 {
1232 
1233  Double_t x[4], y[4], temp1[3],temp2[3];
1234  Int_t p[3];
1235 
1236  TView *view = gPad->GetView();
1237  if (!view) {
1238  Error("PaintTriangles", "No TView in current pad");
1239  return;
1240  }
1241 
1242  TString opt = option;
1243  opt.ToLower();
1244  Bool_t tri1 = opt.Contains("tri1");
1245  Bool_t tri2 = opt.Contains("tri2");
1246  Bool_t markers = opt.Contains("p");
1247  Bool_t markers0 = opt.Contains("p0");
1248  Bool_t wire = opt.Contains("w");
1249 
1250  // Define the grid levels drawn on the triangles.
1251  // The grid levels are aligned on the Z axis' main tick marks.
1252  Int_t nblev=0;
1253  Double_t *glev=0;
1254  if (!tri1 && !tri2 && !wire) {
1255  Int_t ndivz = gCurrentHist->GetZaxis()->GetNdivisions()%100;
1256  Int_t nbins;
1257  Double_t binLow = 0, binHigh = 0, binWidth = 0;
1258 
1259  // Find the main tick marks positions.
1260  Double_t *r0 = view->GetRmin();
1261  Double_t *r1 = view->GetRmax();
1262  if (!r0 || !r1) return;
1263 
1264  if (ndivz > 0) {
1265  THLimitsFinder::Optimize(r0[2], r1[2], ndivz,
1266  binLow, binHigh, nbins, binWidth, " ");
1267  } else {
1268  nbins = TMath::Abs(ndivz);
1269  binLow = r0[2];
1270  binHigh = r1[2];
1271  binWidth = (binHigh-binLow)/nbins;
1272  }
1273  // Define the grid levels
1274  nblev = nbins+1;
1275  glev = new Double_t[nblev];
1276  for (Int_t i = 0; i < nblev; ++i) glev[i] = binLow+i*binWidth;
1277  }
1278 
1279  // Initialize the levels on the Z axis
1280  if (tri1 || tri2) {
1281  Int_t ndiv = gCurrentHist->GetContour();
1282  if (ndiv == 0 ) {
1283  ndiv = gStyle->GetNumberContours();
1284  gCurrentHist->SetContour(ndiv);
1285  }
1287  }
1288 
1289  // For each triangle, compute the distance between the triangle centre
1290  // and the back planes. Then these distances are sorted in order to draw
1291  // the triangles from back to front.
1292  if (!fNdt) FindTriangles();
1293  Double_t cp = TMath::Cos(view->GetLongitude()*TMath::Pi()/180.);
1294  Double_t sp = TMath::Sin(view->GetLongitude()*TMath::Pi()/180.);
1295 
1296  Bool_t reverse = kFALSE;
1298 
1299  if ((cp >= 0) && (sp >= 0.)) {
1300  fDist = [&](Double_t xd, Double_t yd) -> Double_t { return -(fXNmax-xd+fYNmax-yd);};
1301  } else if ((cp <= 0) && (sp >= 0.)) {
1302  fDist = [&](Double_t xd, Double_t yd) -> Double_t { return -(fXNmax-xd+yd-fYNmin);};
1303  reverse = kTRUE;
1304  } else if ((cp <= 0) && (sp <= 0.)) {
1305  fDist = [&](Double_t xd, Double_t yd) -> Double_t { return -(xd-fXNmin+yd-fYNmin);};
1306  } else {
1307  fDist = [&](Double_t xd, Double_t yd) -> Double_t { return -(xd-fXNmin+fYNmax-yd);};
1308  reverse = kTRUE;
1309  }
1310 
1311  typedef std::pair<Double_t, TGraphDelaunay2D::Triangles::const_iterator> DistEntry;
1312  std::vector<DistEntry> dist;
1313  for(auto it = fDelaunay2D->begin(); it != fDelaunay2D->end(); ++it){
1314  auto face = *it;
1315  Double_t xd = (face.x[0] + face.x[1] + face.x[2]) / 3;
1316  Double_t yd = (face.y[0] + face.y[1] + face.y[2]) / 3;
1317 
1318  dist.emplace_back(fDist(xd, yd), it);
1319  }
1320 
1321  std::sort(dist.begin(), dist.end(),
1322  [&](const DistEntry & a, const DistEntry & b){ return !reverse ? (a.first < b.first) : (b.first < a .first); });
1323 
1324  // Draw the triangles and markers if requested
1325  fGraph2D->SetFillColor(fGraph2D->GetFillColor());
1326  Int_t fs = fGraph2D->GetFillStyle();
1327  fGraph2D->SetFillStyle(1001);
1328  fGraph2D->TAttFill::Modify();
1329  fGraph2D->SetLineColor(fGraph2D->GetLineColor());
1330  fGraph2D->TAttLine::Modify();
1331  int lst = fGraph2D->GetLineStyle();
1332  for (const auto & it : dist) {
1333  p[0] = it.second->idx[0];
1334  p[1] = it.second->idx[1];
1335  p[2] = it.second->idx[2];
1336  for (Int_t k=0; k<3; k++) {
1337  if(fX[p[k]] < fXmin || fX[p[k]] > fXmax) goto endloop;
1338  if(fY[p[k]] < fYmin || fY[p[k]] > fYmax) goto endloop;
1339  temp1[0] = fX[p[k]];
1340  temp1[1] = fY[p[k]];
1341  temp1[2] = fZ[p[k]];
1342  temp1[0] = TMath::Max(temp1[0],fXmin);
1343  temp1[1] = TMath::Max(temp1[1],fYmin);
1344  temp1[2] = TMath::Max(temp1[2],fZmin);
1345  temp1[2] = TMath::Min(temp1[2],fZmax);
1346  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
1347  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
1348  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
1349  view->WCtoNDC(temp1, &temp2[0]);
1350  x[k] = temp2[0];
1351  y[k] = temp2[1];
1352  }
1353  x[3] = x[0];
1354  y[3] = y[0];
1355  if (tri1 || tri2) PaintLevels(p,x,y);
1356  if (!tri1 && !tri2 && !wire) {
1357  gPad->PaintFillArea(3,x,y);
1358  PaintLevels(p,x,y,nblev,glev);
1359  }
1360  if (!tri2) gPad->PaintPolyLine(4,x,y);
1361  if (markers) {
1362  if (markers0) {
1363  PaintPolyMarker0(3,x,y);
1364  } else {
1365  fGraph2D->SetMarkerStyle(fGraph2D->GetMarkerStyle());
1366  fGraph2D->SetMarkerSize(fGraph2D->GetMarkerSize());
1367  fGraph2D->SetMarkerColor(fGraph2D->GetMarkerColor());
1368  fGraph2D->TAttMarker::Modify();
1369  gPad->PaintPolyMarker(3,x,y);
1370  }
1371  }
1372 endloop:
1373  continue;
1374  }
1375  fGraph2D->SetFillStyle(fs);
1376  fGraph2D->SetLineStyle(lst);
1377  fGraph2D->TAttLine::Modify();
1378  fGraph2D->TAttFill::Modify();
1379 
1380  if (glev) delete [] glev;
1381 }
c
#define c(i)
Definition: RSha256.hxx:119
l
auto * l
Definition: textangle.C:4
Hoption_t::Logz
int Logz
log scale in Z. Also set by histogram option
Definition: Hoption.h:69
m
auto * m
Definition: textangle.C:8
n
const Int_t n
Definition: legend1.C:16
TAxis
Class to manage histogram axis.
Definition: TAxis.h:30
TGraph2D::GetX
Double_t * GetX() const
Definition: TGraph2D.h:120
first
Definition: first.py:1
Hoption_t
Histogram option structure.
Definition: Hoption.h:24
TGraph2DPainter::fXNmin
Double_t fXNmin
Pointer to fGraph2D->fZE.
Definition: TGraph2DPainter.h:43
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
TObject::TestBit
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
TGraphDelaunay2D::GetGraph2D
TGraph2D * GetGraph2D() const
Definition: TGraphDelaunay2D.h:55
TAxis::GetBinLowEdge
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition: TAxis.cxx:518
TH1::SetContour
virtual void SetContour(Int_t nlevels, const Double_t *levels=0)
Set the number and values of contour levels.
Definition: TH1.cxx:7946
TGraphDelaunay::GetNdt
Int_t GetNdt() const
Definition: TGraphDelaunay.h:86
TMath::Max
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
TGraph2DPainter::fYN
Double_t * fYN
Pointer to fDelaunay->fXN.
Definition: TGraph2DPainter.h:39
TGraph2DPainter::GetContourList
TList * GetContourList(Double_t contour)
Returns the X and Y graphs building a contour.
Definition: TGraph2DPainter.cxx:188
TMath::Cos
Double_t Cos(Double_t)
Definition: TMath.h:630
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
TPolyLine.h
TGraph2D.h
TGraphDelaunay::FindAllTriangles
void FindAllTriangles()
Attempt to find all the Delaunay triangles of the point set.
Definition: TGraphDelaunay.cxx:290
TGraph.h
TGraph2DPainter::PaintContour
void PaintContour(Option_t *option)
Paints the 2D graph as a contour plot.
Definition: TGraph2DPainter.cxx:549
TGraphDelaunay::GetGraph2D
TGraph2D * GetGraph2D() const
Definition: TGraphDelaunay.h:84
TView.h
TObject::Error
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:890
colors
Color * colors
Definition: X3DBuffer.c:21
TGraph2D::GetEZ
virtual Double_t * GetEZ() const
Definition: TGraph2D.h:125
Hoption_t::Logy
int Logy
log scale in Y. Also set by histogram option
Definition: Hoption.h:68
TGraphDelaunay2D::GetNdt
Int_t GetNdt() const
Definition: TGraphDelaunay2D.h:57
TGraph2D::GetN
Int_t GetN() const
Definition: TGraph2D.h:119
TStyle::GetNumberOfColors
Int_t GetNumberOfColors() const
Return number of colors in the color palette.
Definition: TStyle.cxx:1124
TGraphDelaunay2D
Definition: TGraphDelaunay2D.h:32
TAxis::GetFirst
Int_t GetFirst() const
Return first bin on the axis i.e.
Definition: TAxis.cxx:458
TGraph2DPainter::PaintTriangles_new
void PaintTriangles_new(Option_t *option)
Paints the 2D graph as triangles (new implementation)
Definition: TGraph2DPainter.cxx:1230
Float_t
float Float_t
Definition: RtypesCore.h:57
TStyle.h
Int_t
int Int_t
Definition: RtypesCore.h:45
TAxis::GetBinUpEdge
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
Definition: TAxis.cxx:528
TString::Contains
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
x
Double_t x[n]
Definition: legend1.C:17
TGraph2DPainter::PaintPolyMarker0
void PaintPolyMarker0(Int_t n, Double_t *x, Double_t *y)
Paints a circle at each vertex. Each circle background is white.
Definition: TGraph2DPainter.cxx:1041
TList.h
TGraph2DPainter::fEY
Double_t * fEY
Pointer to fGraph2D->fXE.
Definition: TGraph2DPainter.h:41
TGraph2DPainter::fGraph2D
TGraph2D * fGraph2D
Pointer to the TGraphDelaunay2D to be painted.
Definition: TGraph2DPainter.h:62
TGraph2DPainter::fNdt
Int_t fNdt
Equal to fGraph2D->fNpoints.
Definition: TGraph2DPainter.h:54
TMath::Abs
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
TGraphDelaunay::GetXN
Double_t * GetXN() const
Definition: TGraphDelaunay.h:90
TH1::GetContourLevelPad
virtual Double_t GetContourLevelPad(Int_t level) const
Return the value of contour number "level" in Pad coordinates.
Definition: TH1.cxx:7903
TGraph2DPainter::~TGraph2DPainter
virtual ~TGraph2DPainter()
TGraph2DPainter destructor.
Definition: TGraph2DPainter.cxx:148
Hoption_t::Logx
int Logx
log scale in X. Also set by histogram option
Definition: Hoption.h:67
TMath::Sort
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Definition: TMathBase.h:362
TGraph2D::GetMinimum
Double_t GetMinimum() const
Definition: TGraph2D.h:115
TGraphDelaunay::GetYNmax
Double_t GetYNmax() const
Definition: TGraphDelaunay.h:95
TGraph2D::GetZmax
Double_t GetZmax() const
Returns the Z maximum.
Definition: TGraph2D.cxx:1269
TGraph2DPainter::PaintErrors
void PaintErrors(Option_t *option)
Paints the 2D graph as error bars.
Definition: TGraph2DPainter.cxx:591
TGraph2D::GetZminE
virtual Double_t GetZminE() const
Definition: TGraph2D.h:137
TString
Definition: TString.h:136
THLimitsFinder::Optimize
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.
Definition: THLimitsFinder.cxx:184
TGraph2DPainter::PaintTriangles_old
void PaintTriangles_old(Option_t *option)
Paints the 2D graph as triangles (old implementation)
Definition: TGraph2DPainter.cxx:1075
TGraph2D::GetEY
virtual Double_t * GetEY() const
Definition: TGraph2D.h:124
TH1::GetZaxis
TAxis * GetZaxis()
Definition: TH1.h:319
TObject::InheritsFrom
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
TGraphDelaunay2D::GetXNmax
Double_t GetXNmax() const
Definition: TGraphDelaunay2D.h:59
TGraph2DPainter::fZ
Double_t * fZ
Pointer to fGraph2D->fY.
Definition: TGraph2DPainter.h:37
b
#define b(i)
Definition: RSha256.hxx:118
TGraphDelaunay::GetMTried
Int_t * GetMTried() const
Definition: TGraphDelaunay.h:89
TGraph2DPainter::PaintPolyMarker
void PaintPolyMarker(Option_t *option)
Paints the 2D graph as PaintPolyMarker.
Definition: TGraph2DPainter.cxx:890
TGraph2DPainter::fYNmax
Double_t fYNmax
Equal to fDelaunay->fYNmin.
Definition: TGraph2DPainter.h:46
bool
TGraphDelaunay2D::GetYNmax
Double_t GetYNmax() const
Definition: TGraphDelaunay2D.h:61
TGraph2DPainter::fNTried
Int_t * fNTried
Pointer to fDelaunay->fPTried.
Definition: TGraph2DPainter.h:56
THLimitsFinder.h
TGraph2DPainter::fDelaunay
TGraphDelaunay * fDelaunay
Pointer to fDelaunay->fMTried.
Definition: TGraph2DPainter.h:60
TGraphDelaunay2D::GetXNmin
Double_t GetXNmin() const
Definition: TGraphDelaunay2D.h:58
x1
static const double x1[5]
Definition: RooGaussKronrodIntegrator1D.cxx:346
Hoption.h
TH1::GetContour
virtual Int_t GetContour(Double_t *levels=0)
Return contour values into array levels if pointer levels is non zero.
Definition: TH1.cxx:7874
TGraph2DPainter::fEX
Double_t * fEX
Pointer to fDelaunay->fYN.
Definition: TGraph2DPainter.h:40
TGraph2DPainter::fPTried
Int_t * fPTried
Equal to fDelaunay->fNdt.
Definition: TGraph2DPainter.h:55
ROOT::Math::gv_detail::dist
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:63
TGraphDelaunay::GetXNmin
Double_t GetXNmin() const
Definition: TGraphDelaunay.h:92
TGraph2DPainter::fDelaunay2D
TGraphDelaunay2D * fDelaunay2D
Pointer to the TGraphDelaunay2D to be painted.
Definition: TGraph2DPainter.h:61
s0
#define s0(x)
Definition: RSha256.hxx:108
TGraphDelaunay.h
TGraphDelaunay::GetPTried
Int_t * GetPTried() const
Definition: TGraphDelaunay.h:87
TMath::Pi
constexpr Double_t Pi()
Definition: TMath.h:43
Option_t
const typedef char Option_t
Definition: RtypesCore.h:66
TGraphDelaunay::GetYNmin
Double_t GetYNmin() const
Definition: TGraphDelaunay.h:94
TMath::Log10
Double_t Log10(Double_t x)
Definition: TMath.h:753
TView
Definition: TView.h:25
gStyle
R__EXTERN TStyle * gStyle
Definition: TStyle.h:412
TH1::GetYaxis
TAxis * GetYaxis()
Definition: TH1.h:318
TGraph2DPainter::fXmin
Double_t fXmin
Equal to fDelaunay->fYNmax.
Definition: TGraph2DPainter.h:47
TGraph2D::GetZmaxE
virtual Double_t GetZmaxE() const
Definition: TGraph2D.h:136
TGraph2DPainter::fXNmax
Double_t fXNmax
Equal to fDelaunay->fXNmin.
Definition: TGraph2DPainter.h:44
TGraph2DPainter::PaintTriangles
void PaintTriangles(Option_t *option)
Paints the 2D graph as triangles.
Definition: TGraph2DPainter.cxx:1063
TGraph2D::GetMaximum
Double_t GetMaximum() const
Definition: TGraph2D.h:114
TGraph2DPainter::fYmin
Double_t fYmin
Definition: TGraph2DPainter.h:49
TGraphDelaunay
Definition: TGraphDelaunay.h:30
TGraph2DPainter::fX
Double_t * fX
Definition: TGraph2DPainter.h:35
a
auto * a
Definition: textangle.C:12
TGraph2DPainter::TGraph2DPainter
TGraph2DPainter()
TGraph2DPainter default constructor.
Definition: TGraph2DPainter.cxx:47
TGraph2DPainter::fZmax
Double_t fZmax
Definition: TGraph2DPainter.h:52
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
s1
#define s1(x)
Definition: RSha256.hxx:109
ROOT::R::function
void function(const Char_t *name_, T fun, const Char_t *docstring=0)
Definition: RExports.h:151
TGraph2DPainter::fYmax
Double_t fYmax
fGraph2D->fHistogram limits
Definition: TGraph2DPainter.h:50
TGraph2DPainter::FindTriangles
void FindTriangles()
Pointer to the TGraph2D in fDelaunay.
Definition: TGraph2DPainter.cxx:157
TAxis::GetLast
Int_t GetLast() const
Return last bin on the axis i.e.
Definition: TAxis.cxx:469
TH1::GetMinimum
virtual Double_t GetMinimum(Double_t minval=-FLT_MAX) const
Return minimum value larger than minval of bins in the range, unless the value has been overridden by...
Definition: TH1.cxx:8090
TGraphDelaunay2D::end
Triangles::const_iterator end() const
Definition: TGraphDelaunay2D.h:66
TStyle::GetColorPalette
Int_t GetColorPalette(Int_t i) const
Return color number i in current palette.
Definition: TStyle.cxx:1058
TMath::Sin
Double_t Sin(Double_t)
Definition: TMath.h:626
TView::GetRmin
virtual Double_t * GetRmin()=0
TView::GetLongitude
virtual Double_t GetLongitude()=0
gCurrentHist
R__EXTERN TH1 * gCurrentHist
Definition: TGraph2DPainter.cxx:28
TVirtualPad.h
TGraph2DPainter::fXmax
Double_t fXmax
Definition: TGraph2DPainter.h:48
TStyle::GetNumberContours
Int_t GetNumberContours() const
Definition: TStyle.h:232
y
Double_t y[n]
Definition: legend1.C:17
TGraph2DPainter::fY
Double_t * fY
Pointer to fGraph2D->fX.
Definition: TGraph2DPainter.h:36
TGraph2DPainter::PaintPolyLine
void PaintPolyLine(Option_t *option)
Paints the 2D graph as PaintPolyLine.
Definition: TGraph2DPainter.cxx:990
line
TLine * line
Definition: entrylistblock_figure1.C:235
TGraph2D::GetEX
virtual Double_t * GetEX() const
Definition: TGraph2D.h:123
TGraph2DPainter::fXN
Double_t * fXN
Pointer to fGraph2D->fZ.
Definition: TGraph2DPainter.h:38
TMath::Min
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
TGraphDelaunay2D::begin
Triangles::const_iterator begin() const
Definition: TGraphDelaunay2D.h:65
TGraph2DPainter.h
TGraphDelaunay::GetNTried
Int_t * GetNTried() const
Definition: TGraphDelaunay.h:88
TGraph2DPainter::fZmin
Double_t fZmin
Definition: TGraph2DPainter.h:51
TGraph2D::GetY
Double_t * GetY() const
Definition: TGraph2D.h:121
Double_t
double Double_t
Definition: RtypesCore.h:59
TGraph
Definition: TGraph.h:41
TGraph2DPainter::Paint
void Paint(Option_t *option)
Paint a TGraphDelaunay according to the value of "option":
Definition: TGraph2DPainter.cxx:504
TGraph2DPainter::fMTried
Int_t * fMTried
Pointer to fDelaunay->fNTried.
Definition: TGraph2DPainter.h:57
TH1::kUserContour
@ kUserContour
user specified contour levels
Definition: TH1.h:162
TGraph2D::GetZ
Double_t * GetZ() const
Definition: TGraph2D.h:122
TGraphDelaunay2D.h
graph
Definition: graph.py:1
TGraph2DPainter::fYNmin
Double_t fYNmin
Equal to fDelaunay->fXNmax.
Definition: TGraph2DPainter.h:45
TList::Add
virtual void Add(TObject *obj)
Definition: TList.h:87
TObject
Definition: TObject.h:37
TH1
Definition: TH1.h:57
TGraph2DPainter::fEZ
Double_t * fEZ
Pointer to fGraph2D->fYE.
Definition: TGraph2DPainter.h:42
x2
static const double x2[5]
Definition: RooGaussKronrodIntegrator1D.cxx:364
TPolyMarker.h
gPad
#define gPad
Definition: TVirtualPad.h:287
TGraph2DPainter
The TGraphDelaunay painting class.
Definition: TGraph2DPainter.h:31
TIter
Definition: TCollection.h:233
Hoption
R__EXTERN Hoption_t Hoption
Definition: TGraph2DPainter.cxx:29
TAxis::FindFixBin
virtual Int_t FindFixBin(Double_t x) const
Find bin number corresponding to abscissa x.
Definition: TAxis.cxx:419
TGraphDelaunay2D::FindAllTriangles
void FindAllTriangles()
Definition: TGraphDelaunay2D.h:53
TGraph2DPainter::PaintLevels
void PaintLevels(Int_t *v, Double_t *x, Double_t *y, Int_t nblev=0, Double_t *glev=0)
Paints one triangle.
Definition: TGraph2DPainter.cxx:696
R__EXTERN
#define R__EXTERN
Definition: DllImport.h:27
TGraphDelaunay::GetXNmax
Double_t GetXNmax() const
Definition: TGraphDelaunay.h:93
TGeant4Unit::ms
static constexpr double ms
Definition: TGeant4SystemOfUnits.h:169
TGraph2DPainter::fNpoints
Int_t fNpoints
Definition: TGraph2DPainter.h:53
Class
void Class()
Definition: Class.C:29
TString::ToLower
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1125
TH1::GetXaxis
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:317
TH1::GetMaximum
virtual Double_t GetMaximum(Double_t maxval=FLT_MAX) const
Return maximum value smaller than maxval of bins in the range, unless the value has been overridden b...
Definition: TH1.cxx:8005
TGraphDelaunay2D::GetYNmin
Double_t GetYNmin() const
Definition: TGraphDelaunay2D.h:60
TH1.h
TView::WCtoNDC
virtual void WCtoNDC(const Float_t *pw, Float_t *pn)=0
TAttAxis::GetNdivisions
virtual Int_t GetNdivisions() const
Definition: TAttAxis.h:42
TGraphDelaunay::GetYN
Double_t * GetYN() const
Definition: TGraphDelaunay.h:91
TView::GetRmax
virtual Double_t * GetRmax()=0
TList
Definition: TList.h:44
TAxis::GetBinWidth
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width.
Definition: TAxis.cxx:540
TMath.h
int
g
#define g(i)
Definition: RSha256.hxx:123