ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
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 "TROOT.h"
13 #include "TGraph2DPainter.h"
14 #include "TMath.h"
15 #include "TList.h"
16 #include "TGraph.h"
17 #include "TGraph2D.h"
18 #include "TGraphDelaunay.h"
19 #include "TGraphDelaunay2D.h"
20 #include "TPolyLine.h"
21 #include "TPolyMarker.h"
22 #include "TVirtualPad.h"
23 #include "TView.h"
24 #include "THLimitsFinder.h"
25 #include "TStyle.h"
26 #include "Hoption.h"
27 #include "TH1.h"
28 
31 
33 
34 
35 ////////////////////////////////////////////////////////////////////////////////
36 /*! \class TGraph2DPainter
37  \ingroup Histpainter
38  \brief The TGraphDelaunay painting class.
39 
40 TGraph2DPainter paints a TGraphDelaunay using triangles or clouds of points.
41 
42 */
43 
44 
45 ////////////////////////////////////////////////////////////////////////////////
46 /// TGraph2DPainter default constructor
47 
49 {
50  fX = 0;
51  fY = 0;
52  fZ = 0;
53  fEX = 0;
54  fEY = 0;
55  fEZ = 0;
56  fXN = 0;
57  fYN = 0;
58  fPTried = 0;
59  fNTried = 0;
60  fMTried = 0;
61  fGraph2D = 0;
62  fDelaunay = 0;
63  fDelaunay2D = 0;
64  fXmin = 0.;
65  fXmax = 0.;
66  fYmin = 0.;
67  fYmax = 0.;
68  fZmin = 0.;
69  fZmax = 0.;
70  fXNmin = 0;
71  fXNmax = 0;
72  fYNmin = 0;
73  fYNmax = 0;
74  fNdt = 0;
75  fNpoints = 0;
76 }
77 
78 
79 ////////////////////////////////////////////////////////////////////////////////
80 /// TGraph2DPainter constructor using the old Delaunay algorithm
81 
83 {
84  fDelaunay = gd;
85  fDelaunay2D = 0;
87  fNpoints = fGraph2D->GetN();
88  fX = fGraph2D->GetX();
89  fY = fGraph2D->GetY();
90  fZ = fGraph2D->GetZ();
91  fEX = fGraph2D->GetEX();
92  fEY = fGraph2D->GetEY();
93  fEZ = fGraph2D->GetEZ();
94  fNdt = 0;
95  fXN = 0;
96  fYN = 0;
97  fXNmin = 0;
98  fXNmax = 0;
99  fYNmin = 0;
100  fYNmax = 0;
101  fPTried = 0;
102  fNTried = 0;
103  fMTried = 0;
104  fXmin = 0.;
105  fXmax = 0.;
106  fYmin = 0.;
107  fYmax = 0.;
108  fZmin = 0.;
109  fZmax = 0.;
110 }
111 
112 ////////////////////////////////////////////////////////////////////////////////
113 /// TGraph2DPainter constructor using the new Delaunay algorithm
114 
116 {
117  fDelaunay = 0;
118  fDelaunay2D = gd;
120  fNpoints = fGraph2D->GetN();
121  fX = fGraph2D->GetX();
122  fY = fGraph2D->GetY();
123  fZ = fGraph2D->GetZ();
124  fEX = fGraph2D->GetEX();
125  fEY = fGraph2D->GetEY();
126  fEZ = fGraph2D->GetEZ();
127  fNdt = 0;
128  fXN = 0;
129  fYN = 0;
130  fXNmin = 0;
131  fXNmax = 0;
132  fYNmin = 0;
133  fYNmax = 0;
134  fPTried = 0;
135  fNTried = 0;
136  fMTried = 0;
137  fXmin = 0.;
138  fXmax = 0.;
139  fYmin = 0.;
140  fYmax = 0.;
141  fZmin = 0.;
142  fZmax = 0.;
143 }
144 
145 
146 ////////////////////////////////////////////////////////////////////////////////
147 /// TGraph2DPainter destructor.
148 
150 {
151 }
152 
153 
154 ////////////////////////////////////////////////////////////////////////////////
155 /// Find triangles in fDelaunay and initialise the TGraph2DPainter values
156 /// needed to paint triangles or find contours.
157 
159 {
160  if (fDelaunay) {
162  fNdt = fDelaunay->GetNdt();
163  fXN = fDelaunay->GetXN();
164  fYN = fDelaunay->GetYN();
172  }
173  else if (fDelaunay2D) {
175  fNdt = fDelaunay2D->GetNdt();
180  }
181 }
182 
183 
184 ////////////////////////////////////////////////////////////////////////////////
185 /// Returns the X and Y graphs building a contour. A contour level may
186 /// consist in several parts not connected to each other. This function
187 /// finds them and returns them in a graphs' list.
188 
190 {
191  // Exit if the contour is outisde the Z range.
192  Double_t zmin = gCurrentHist->GetMinimum();
193  Double_t zmax = gCurrentHist->GetMaximum();
194  if (Hoption.Logz) {
195  if (zmin > 0) {
196  zmin = TMath::Log10(zmin);
197  zmax = TMath::Log10(zmax);
198  } else {
199  return 0;
200  }
201  }
202  if(contour<zmin || contour>zmax) {
203  Error("GetContourList", "Contour level (%g) outside the Z scope [%g,%g]",
204  contour,zmin,zmax);
205  return 0;
206  }
207 
208  if (!fNdt) FindTriangles();
209 
210  TGraph *graph = 0; // current graph
211  Int_t npg = 0; // number of points in the current graph
212  TList *list = new TList(); // list holding all the graphs
213 
214  // Find all the segments making the contour
215 
216  Double_t r21, r20, r10;
217  Int_t p0, p1, p2;
218  Double_t x0, y0, z0;
219  Double_t x1, y1, z1;
220  Double_t x2, y2, z2;
221  Int_t t[3],i,it,i0,i1,i2;
222 
223  // Allocate space to store the segments. They cannot be more than the
224  // number of triangles.
225  Double_t xs0c, ys0c, xs1c, ys1c;
226  Double_t *xs0 = new Double_t[fNdt];
227  Double_t *ys0 = new Double_t[fNdt];
228  Double_t *xs1 = new Double_t[fNdt];
229  Double_t *ys1 = new Double_t[fNdt];
230  for (i=0;i<fNdt;i++) {
231  xs0[i] = 0.;
232  ys0[i] = 0.;
233  xs1[i] = 0.;
234  ys1[i] = 0.;
235  }
236  Int_t nbSeg = 0;
237 
238  // Loop over all the triangles in order to find all the line segments
239  // making the contour.
240 
241  // old implementation
242  if (fDelaunay) {
243  for(it=0; it<fNdt; it++) {
244  t[0] = fPTried[it];
245  t[1] = fNTried[it];
246  t[2] = fMTried[it];
247  p0 = t[0]-1;
248  p1 = t[1]-1;
249  p2 = t[2]-1;
250  x0 = fX[p0]; x2 = fX[p0];
251  y0 = fY[p0]; y2 = fY[p0];
252  z0 = fZ[p0]; z2 = fZ[p0];
253 
254  // Order along Z axis the points (xi,yi,zi) where "i" belongs to {0,1,2}
255  // After this z0 < z1 < z2
256  i0=0, i1=0, i2=0;
257  if (fZ[p1]<=z0) {z0=fZ[p1]; x0=fX[p1]; y0=fY[p1]; i0=1;}
258  if (fZ[p1]>z2) {z2=fZ[p1]; x2=fX[p1]; y2=fY[p1]; i2=1;}
259  if (fZ[p2]<=z0) {z0=fZ[p2]; x0=fX[p2]; y0=fY[p2]; i0=2;}
260  if (fZ[p2]>z2) {z2=fZ[p2]; x2=fX[p2]; y2=fY[p2]; i2=2;}
261  if (i0==0 && i2==0) {
262  Error("GetContourList", "wrong vertices ordering");
263  delete [] xs0;
264  delete [] ys0;
265  delete [] xs1;
266  delete [] ys1;
267  return 0;
268  } else {
269  i1 = 3-i2-i0;
270  }
271  x1 = fX[t[i1]-1];
272  y1 = fY[t[i1]-1];
273  z1 = fZ[t[i1]-1];
274 
275  if (Hoption.Logz) {
276  z0 = TMath::Log10(z0);
277  z1 = TMath::Log10(z1);
278  z2 = TMath::Log10(z2);
279  }
280 
281  if(contour >= z0 && contour <=z2) {
282  r20 = (contour-z0)/(z2-z0);
283  xs0c = r20*(x2-x0)+x0;
284  ys0c = r20*(y2-y0)+y0;
285  if(contour >= z1 && contour <=z2) {
286  r21 = (contour-z1)/(z2-z1);
287  xs1c = r21*(x2-x1)+x1;
288  ys1c = r21*(y2-y1)+y1;
289  } else {
290  r10 = (contour-z0)/(z1-z0);
291  xs1c = r10*(x1-x0)+x0;
292  ys1c = r10*(y1-y0)+y0;
293  }
294  // do not take the segments equal to a point
295  if(xs0c != xs1c || ys0c != ys1c) {
296  nbSeg++;
297  xs0[nbSeg-1] = xs0c;
298  ys0[nbSeg-1] = ys0c;
299  xs1[nbSeg-1] = xs1c;
300  ys1[nbSeg-1] = ys1c;
301  }
302  }
303  }
304  }
305  else if (fDelaunay2D) {
306 
307  Int_t p[3];
308  for(const auto & face : *fDelaunay2D) {
309  p[0] = face.idx[0];
310  p[1] = face.idx[1];
311  p[2] = face.idx[2];
312  x0 = fX[p[0]]; x2 = fX[p[0]];
313  y0 = fY[p[0]]; y2 = fY[p[0]];
314  z0 = fZ[p[0]]; z2 = fZ[p[0]];
315 
316  // Order along Z axis the points (xi,yi,zi) where "i" belongs to {0,1,2}
317  // After this z0 < z1 < z2
318  i0=0, i1=0, i2=0;
319  if (fZ[p[1]]<=z0) {z0=fZ[p[1]]; x0=fX[p[1]]; y0=fY[p[1]]; i0=1;}
320  if (fZ[p[1]]>z2) {z2=fZ[p[1]]; x2=fX[p[1]]; y2=fY[p[1]]; i2=1;}
321  if (fZ[p[2]]<=z0) {z0=fZ[p[2]]; x0=fX[p[2]]; y0=fY[p[2]]; i0=2;}
322  if (fZ[p[2]]>z2) {z2=fZ[p[2]]; x2=fX[p[2]]; y2=fY[p[2]]; i2=2;}
323  if (i0==0 && i2==0) {
324  Error("GetContourList", "wrong vertices ordering");
325  delete [] xs0;
326  delete [] ys0;
327  delete [] xs1;
328  delete [] ys1;
329  return 0;
330  } else {
331  i1 = 3-i2-i0;
332  }
333  x1 = fX[p[i1]];
334  y1 = fY[p[i1]];
335  z1 = fZ[p[i1]];
336 
337  if (Hoption.Logz) {
338  z0 = TMath::Log10(z0);
339  z1 = TMath::Log10(z1);
340  z2 = TMath::Log10(z2);
341  }
342 
343  if(contour >= z0 && contour <=z2) {
344  r20 = (contour-z0)/(z2-z0);
345  xs0c = r20*(x2-x0)+x0;
346  ys0c = r20*(y2-y0)+y0;
347  if(contour >= z1 && contour <=z2) {
348  r21 = (contour-z1)/(z2-z1);
349  xs1c = r21*(x2-x1)+x1;
350  ys1c = r21*(y2-y1)+y1;
351  } else {
352  r10 = (contour-z0)/(z1-z0);
353  xs1c = r10*(x1-x0)+x0;
354  ys1c = r10*(y1-y0)+y0;
355  }
356  // do not take the segments equal to a point
357  if(xs0c != xs1c || ys0c != ys1c) {
358  nbSeg++;
359  xs0[nbSeg-1] = xs0c;
360  ys0[nbSeg-1] = ys0c;
361  xs1[nbSeg-1] = xs1c;
362  ys1[nbSeg-1] = ys1c;
363  }
364  }
365  }
366  }
367 
368  Bool_t *segUsed = new Bool_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 (segUsed[js] && js<nbSeg) {
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 (segUsed[js] && js<nbSeg) {
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->GetZmax();
534  fZmin = fGraph2D->GetZmin();
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;
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  }
584 }
585 
586 
587 ////////////////////////////////////////////////////////////////////////////////
588 /// Paints the 2D graph as error bars
589 
591 {
592  Double_t temp1[3],temp2[3];
593 
594  TView *view = gPad->GetView();
595  if (!view) {
596  Error("PaintErrors", "No TView in current pad");
597  return;
598  }
599 
600  Int_t it;
601 
602  Double_t *xm = new Double_t[2];
603  Double_t *ym = new Double_t[2];
604 
608  fGraph2D->TAttLine::Modify();
609 
610  for (it=0; it<fNpoints; it++) {
611  if(fX[it] < fXmin || fX[it] > fXmax) continue;
612  if(fY[it] < fYmin || fY[it] > fYmax) continue;
613  if (fEX) {
614  temp1[0] = fX[it]-fEX[it];
615  temp1[1] = fY[it];
616  temp1[2] = fZ[it];
617  temp1[0] = TMath::Max(temp1[0],fXmin);
618  temp1[1] = TMath::Max(temp1[1],fYmin);
619  temp1[2] = TMath::Max(temp1[2],fZmin);
620  temp1[2] = TMath::Min(temp1[2],fZmax);
621  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
622  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
623  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
624  view->WCtoNDC(temp1, &temp2[0]);
625  xm[0] = temp2[0];
626  ym[0] = temp2[1];
627 
628  temp1[0] = fX[it]+fEX[it];
629  temp1[0] = TMath::Max(temp1[0],fXmin);
630  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
631  view->WCtoNDC(temp1, &temp2[0]);
632  xm[1] = temp2[0];
633  ym[1] = temp2[1];
634  gPad->PaintPolyLine(2,xm,ym);
635  }
636  if (fEY) {
637  temp1[0] = fX[it];
638  temp1[1] = fY[it]-fEY[it];
639  temp1[2] = fZ[it];
640  temp1[0] = TMath::Max(temp1[0],fXmin);
641  temp1[1] = TMath::Max(temp1[1],fYmin);
642  temp1[2] = TMath::Max(temp1[2],fZmin);
643  temp1[2] = TMath::Min(temp1[2],fZmax);
644  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
645  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
646  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
647  view->WCtoNDC(temp1, &temp2[0]);
648  xm[0] = temp2[0];
649  ym[0] = temp2[1];
650 
651  temp1[1] = fY[it]+fEY[it];
652  temp1[1] = TMath::Max(temp1[1],fYmin);
653  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
654  view->WCtoNDC(temp1, &temp2[0]);
655  xm[1] = temp2[0];
656  ym[1] = temp2[1];
657  gPad->PaintPolyLine(2,xm,ym);
658  }
659  if (fEZ) {
660  temp1[0] = fX[it];
661  temp1[1] = fY[it];
662  temp1[2] = fZ[it]-fEZ[it];
663  temp1[0] = TMath::Max(temp1[0],fXmin);
664  temp1[1] = TMath::Max(temp1[1],fYmin);
665  temp1[2] = TMath::Max(temp1[2],fZmin);
666  temp1[2] = TMath::Min(temp1[2],fZmax);
667  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
668  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
669  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
670  view->WCtoNDC(temp1, &temp2[0]);
671  xm[0] = temp2[0];
672  ym[0] = temp2[1];
673 
674  temp1[2] = fZ[it]+fEZ[it];
675  temp1[2] = TMath::Max(temp1[2],fZmin);
676  temp1[2] = TMath::Min(temp1[2],fZmax);
677  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
678  view->WCtoNDC(temp1, &temp2[0]);
679  xm[1] = temp2[0];
680  ym[1] = temp2[1];
681  gPad->PaintPolyLine(2,xm,ym);
682  }
683  }
684  delete [] xm;
685  delete [] ym;
686 }
687 
688 
689 ////////////////////////////////////////////////////////////////////////////////
690 /// Paints one triangle.
691 ///
692 /// - nblev = 0 : paint the color levels
693 /// - nblev != 0 : paint the grid
694 
696  Int_t nblev, Double_t *glev)
697 {
698  Int_t i, fillColor, ncolors, theColor0, theColor2;
699 
700  Int_t p[3];
701  if (fDelaunay) {
702  p[0]=t[0]-1;
703  p[1]=t[1]-1;
704  p[2]=t[2]-1;
705  }
706  else {
707  p[0]=t[0];
708  p[1]=t[1];
709  p[2]=t[2];
710  }
711  Double_t xl[2],yl[2];
712  Double_t zl, r21, r20, r10;
713  Double_t x0 = x[0] , x2 = x[0];
714  Double_t y0 = y[0] , y2 = y[0];
715  Double_t z0 = fZ[p[0]], z2 = fZ[p[0]];
716  Double_t zmin = fGraph2D->GetMinimum();
717  Double_t zmax = fGraph2D->GetMaximum();
718  if (zmin==-1111 && zmax==-1111) {
719  zmin = fZmin;
720  zmax = fZmax;
721  }
722 
723  // Order along Z axis the points (xi,yi,zi) where "i" belongs to {0,1,2}
724  // After this z0 < z1 < z2
725  Int_t i0=0, i1=0, i2=0;
726  if (fZ[p[1]]<=z0) {z0=fZ[p[1]]; x0=x[1]; y0=y[1]; i0=1;}
727  if (fZ[p[1]]>z2) {z2=fZ[p[1]]; x2=x[1]; y2=y[1]; i2=1;}
728  if (fZ[p[2]]<=z0) {z0=fZ[p[2]]; x0=x[2]; y0=y[2]; i0=2;}
729  if (fZ[p[2]]>z2) {z2=fZ[p[2]]; x2=x[2]; y2=y[2]; i2=2;}
730  i1 = 3-i2-i0;
731  Double_t x1 = x[i1];
732  Double_t y1 = y[i1];
733  Double_t z1 = fZ[p[i1]];
734 
735  if (z0>zmax) z0 = zmax;
736  if (z2>zmax) z2 = zmax;
737  if (z0<zmin) z0 = zmin;
738  if (z2<zmin) z2 = zmin;
739  if (z1>zmax) z1 = zmax;
740  if (z1<zmin) z1 = zmin;
741 
742  if (Hoption.Logz) {
743  z0 = TMath::Log10(z0);
744  z1 = TMath::Log10(z1);
745  z2 = TMath::Log10(z2);
746  zmin = TMath::Log10(zmin);
747  zmax = TMath::Log10(zmax);
748  }
749 
750  // zi = Z values of the stripe number i
751  // zip = Previous zi
752  Double_t zi=0, zip=0;
753 
754  if (nblev <= 0) {
755  // Paint the colors levels
756 
757  // Compute the color associated to z0 (theColor0) and z2 (theColor2)
758  ncolors = gStyle->GetNumberOfColors();
759  theColor0 = (Int_t)( ((z0-zmin)/(zmax-zmin))*(ncolors-1) );
760  theColor2 = (Int_t)( ((z2-zmin)/(zmax-zmin))*(ncolors-1) );
761 
762  // The stripes drawn to fill the triangles may have up to 5 points
763  Double_t xp[5], yp[5];
764 
765  // rl = Ratio between z0 and z2 (long)
766  // rs = Ratio between z0 and z1 or z1 and z2 (short)
767  Double_t rl,rs;
768 
769  // ci = Color of the stripe number i
770  // npf = number of point needed to draw the current stripe
771  Int_t ci,npf;
772 
773  fillColor = fGraph2D->GetFillColor();
774 
775  // If the z0's color and z2's colors are the same, the whole triangle
776  // can be painted in one go.
777  if(theColor0 == theColor2) {
779  fGraph2D->TAttFill::Modify();
780  gPad->PaintFillArea(3,x,y);
781 
782  // The triangle must be painted with several colors
783  } else {
784  for(ci=theColor0; ci<=theColor2; ci++) {
786  fGraph2D->TAttFill::Modify();
787  if (ci==theColor0) {
788  zi = (((ci+1)*(zmax-zmin))/(ncolors-1))+zmin;
789  xp[0] = x0;
790  yp[0] = y0;
791  rl = (zi-z0)/(z2-z0);
792  xp[1] = rl*(x2-x0)+x0;
793  yp[1] = rl*(y2-y0)+y0;
794  if (zi>=z1 || z0==z1) {
795  rs = (zi-z1)/(z2-z1);
796  xp[2] = rs*(x2-x1)+x1;
797  yp[2] = rs*(y2-y1)+y1;
798  xp[3] = x1;
799  yp[3] = y1;
800  npf = 4;
801  } else {
802  rs = (zi-z0)/(z1-z0);
803  xp[2] = rs*(x1-x0)+x0;
804  yp[2] = rs*(y1-y0)+y0;
805  npf = 3;
806  }
807  } else if (ci==theColor2) {
808  xp[0] = xp[1];
809  yp[0] = yp[1];
810  xp[1] = x2;
811  yp[1] = y2;
812  if (zi<z1 || z2==z1) {
813  xp[3] = xp[2];
814  yp[3] = yp[2];
815  xp[2] = x1;
816  yp[2] = y1;
817  npf = 4;
818  } else {
819  npf = 3;
820  }
821  } else {
822  zi = (((ci+1)*(zmax-zmin))/(ncolors-1))+zmin;
823  xp[0] = xp[1];
824  yp[0] = yp[1];
825  rl = (zi-z0)/(z2-z0);
826  xp[1] = rl*(x2-x0)+x0;
827  yp[1] = rl*(y2-y0)+y0;
828  if ( zi>=z1 && zip<=z1) {
829  xp[3] = x1;
830  yp[3] = y1;
831  xp[4] = xp[2];
832  yp[4] = yp[2];
833  npf = 5;
834  } else {
835  xp[3] = xp[2];
836  yp[3] = yp[2];
837  npf = 4;
838  }
839  if (zi<z1) {
840  rs = (zi-z0)/(z1-z0);
841  xp[2] = rs*(x1-x0)+x0;
842  yp[2] = rs*(y1-y0)+y0;
843  } else {
844  rs = (zi-z1)/(z2-z1);
845  xp[2] = rs*(x2-x1)+x1;
846  yp[2] = rs*(y2-y1)+y1;
847  }
848  }
849  zip = zi;
850  // Paint a stripe
851  gPad->PaintFillArea(npf,xp,yp);
852  }
853  }
854  fGraph2D->SetFillColor(fillColor);
855  fGraph2D->TAttFill::Modify();
856 
857  } else {
858  // Paint the grid levels
860  fGraph2D->TAttLine::Modify();
861  for(i=0; i<nblev; i++){
862  zl=glev[i];
863  if(zl >= z0 && zl <=z2) {
864  r21=(zl-z1)/(z2-z1);
865  r20=(zl-z0)/(z2-z0);
866  r10=(zl-z0)/(z1-z0);
867  xl[0]=r20*(x2-x0)+x0;
868  yl[0]=r20*(y2-y0)+y0;
869  if(zl >= z1 && zl <=z2) {
870  xl[1]=r21*(x2-x1)+x1;
871  yl[1]=r21*(y2-y1)+y1;
872  } else {
873  xl[1]=r10*(x1-x0)+x0;
874  yl[1]=r10*(y1-y0)+y0;
875  }
876  gPad->PaintPolyLine(2,xl,yl);
877  }
878  }
880  fGraph2D->TAttLine::Modify();
881  }
882 }
883 
884 
885 ////////////////////////////////////////////////////////////////////////////////
886 /// Paints the 2D graph as PaintPolyMarker
887 
889 {
890  Double_t temp1[3],temp2[3];
891 
892  TView *view = gPad->GetView();
893  if (!view) {
894  Error("PaintPolyMarker", "No TView in current pad");
895  return;
896  }
897 
898  TString opt = option;
899  opt.ToLower();
900  Bool_t markers0 = opt.Contains("p0");
901  Bool_t colors = opt.Contains("pcol");
902  Int_t ncolors = gStyle->GetNumberOfColors();
903  Int_t it, theColor;
904 
905  // Initialize the levels on the Z axis
906  if (colors) {
907  Int_t ndiv = gCurrentHist->GetContour();
908  if (ndiv == 0 ) {
909  ndiv = gStyle->GetNumberContours();
910  gCurrentHist->SetContour(ndiv);
911  }
913  }
914 
915  Double_t *xm = new Double_t[fNpoints];
916  Double_t *ym = new Double_t[fNpoints];
917  Double_t *zm = new Double_t[fNpoints];
918  Double_t hzmin = gCurrentHist->GetMinimum();
919  Double_t hzmax = gCurrentHist->GetMaximum();
920  Int_t npd = 0;
921  for (it=0; it<fNpoints; it++) {
922  xm[it] = 0;
923  ym[it] = 0;
924  if(fX[it] < fXmin || fX[it] > fXmax) continue;
925  if(fY[it] < fYmin || fY[it] > fYmax) continue;
926  if(fZ[it] < hzmin || fZ[it] > hzmax) continue;
927  temp1[0] = fX[it];
928  temp1[1] = fY[it];
929  temp1[2] = fZ[it];
930  temp1[0] = TMath::Max(temp1[0],fXmin);
931  temp1[1] = TMath::Max(temp1[1],fYmin);
932  temp1[2] = TMath::Max(temp1[2],hzmin);
933  temp1[2] = TMath::Min(temp1[2],hzmax);
934  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
935  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
936  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
937  view->WCtoNDC(temp1, &temp2[0]);
938  xm[npd] = temp2[0];
939  ym[npd] = temp2[1];
940  zm[npd] = fZ[it];
941  npd++;
942  }
943  if (markers0) {
944  PaintPolyMarker0(npd,xm,ym);
945  } else if (colors) {
946  Int_t cols = fGraph2D->GetMarkerColor();
947  for (it=0; it<npd; it++) {
948  theColor = (Int_t)( ((zm[it]-hzmin)/(hzmax-hzmin))*(ncolors-1) );
950  fGraph2D->TAttMarker::Modify();
951  gPad->PaintPolyMarker(1,&xm[it],&ym[it]);
952  }
953  fGraph2D->SetMarkerColor(cols);
954  } else {
958  fGraph2D->TAttMarker::Modify();
959  gPad->PaintPolyMarker(npd,xm,ym);
960  }
961  delete [] xm;
962  delete [] ym;
963  delete [] zm;
964 }
965 
966 
967 ////////////////////////////////////////////////////////////////////////////////
968 /// Paints the 2D graph as PaintPolyLine
969 
971 {
972  Double_t temp1[3],temp2[3];
973 
974  TView *view = gPad->GetView();
975  if (!view) {
976  Error("PaintPolyLine", "No TView in current pad");
977  return;
978  }
979 
980  Int_t it;
981 
982  Double_t *xm = new Double_t[fNpoints];
983  Double_t *ym = new Double_t[fNpoints];
984  Int_t npd = 0;
985  for (it=0; it<fNpoints; it++) {
986  if(fX[it] < fXmin || fX[it] > fXmax) continue;
987  if(fY[it] < fYmin || fY[it] > fYmax) continue;
988  npd++;
989  temp1[0] = fX[it];
990  temp1[1] = fY[it];
991  temp1[2] = fZ[it];
992  temp1[0] = TMath::Max(temp1[0],fXmin);
993  temp1[1] = TMath::Max(temp1[1],fYmin);
994  temp1[2] = TMath::Max(temp1[2],fZmin);
995  temp1[2] = TMath::Min(temp1[2],fZmax);
996  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
997  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
998  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
999  view->WCtoNDC(temp1, &temp2[0]);
1000  xm[it] = temp2[0];
1001  ym[it] = temp2[1];
1002  }
1006  fGraph2D->TAttLine::Modify();
1007  gPad->PaintPolyLine(npd,xm,ym);
1008  delete [] xm;
1009  delete [] ym;
1010 }
1011 
1012 
1013 ////////////////////////////////////////////////////////////////////////////////
1014 /// Paints a circle at each vertex. Each circle background is white.
1015 
1017 {
1019  Int_t mc = fGraph2D->GetMarkerColor();
1020  Int_t ms = fGraph2D->GetMarkerStyle();
1021  for (Int_t i=0; i<n; i++) {
1022  fGraph2D->SetMarkerStyle(20);
1024  fGraph2D->TAttMarker::Modify();
1025  gPad->PaintPolyMarker(1,&x[i],&y[i]);
1026  fGraph2D->SetMarkerStyle(24);
1027  fGraph2D->SetMarkerColor(mc);
1028  fGraph2D->TAttMarker::Modify();
1029  gPad->PaintPolyMarker(1,&x[i],&y[i]);
1030  }
1031  fGraph2D->SetMarkerStyle(ms);
1032 }
1033 
1034 
1035 ////////////////////////////////////////////////////////////////////////////////
1036 /// Paints the 2D graph as triangles
1038 {
1039  if (fDelaunay)
1040  PaintTriangles_old(option);
1041  else if (fDelaunay2D)
1042  PaintTriangles_new(option);
1043 }
1044 
1045 
1047 {
1048  Double_t x[4], y[4], temp1[3],temp2[3];
1049  Int_t it,t[3];
1050  Int_t *order = 0;
1051  Double_t *dist = 0;
1052 
1053  TView *view = gPad->GetView();
1054  if (!view) {
1055  Error("PaintTriangles", "No TView in current pad");
1056  return;
1057  }
1058 
1059  TString opt = option;
1060  opt.ToLower();
1061  Bool_t tri1 = opt.Contains("tri1");
1062  Bool_t tri2 = opt.Contains("tri2");
1063  Bool_t markers = opt.Contains("p");
1064  Bool_t markers0 = opt.Contains("p0");
1065  Bool_t wire = opt.Contains("w");
1066 
1067  // Define the grid levels drawn on the triangles.
1068  // The grid levels are aligned on the Z axis' main tick marks.
1069  Int_t nblev=0;
1070  Double_t *glev=0;
1071  if (!tri1 && !tri2 && !wire) {
1072  Int_t ndivz = gCurrentHist->GetZaxis()->GetNdivisions()%100;
1073  Int_t nbins;
1074  Double_t binLow = 0, binHigh = 0, binWidth = 0;
1075 
1076  // Find the main tick marks positions.
1077  Double_t *r0 = view->GetRmin();
1078  Double_t *r1 = view->GetRmax();
1079  if (!r0 || !r1) return;
1080 
1081  if (ndivz > 0) {
1082  THLimitsFinder::Optimize(r0[2], r1[2], ndivz,
1083  binLow, binHigh, nbins, binWidth, " ");
1084  } else {
1085  nbins = TMath::Abs(ndivz);
1086  binLow = r0[2];
1087  binHigh = r1[2];
1088  binWidth = (binHigh-binLow)/nbins;
1089  }
1090  // Define the grid levels
1091  nblev = nbins+1;
1092  glev = new Double_t[nblev];
1093  for (Int_t i = 0; i < nblev; ++i) glev[i] = binLow+i*binWidth;
1094  }
1095 
1096  // Initialize the levels on the Z axis
1097  if (tri1 || tri2) {
1098  Int_t ndiv = gCurrentHist->GetContour();
1099  if (ndiv == 0 ) {
1100  ndiv = gStyle->GetNumberContours();
1101  gCurrentHist->SetContour(ndiv);
1102  }
1104  }
1105 
1106  // For each triangle, compute the distance between the triangle centre
1107  // and the back planes. Then these distances are sorted in order to draw
1108  // the triangles from back to front.
1109  if (!fNdt) FindTriangles();
1110  Double_t cp = TMath::Cos(view->GetLongitude()*TMath::Pi()/180.);
1111  Double_t sp = TMath::Sin(view->GetLongitude()*TMath::Pi()/180.);
1112  order = new Int_t[fNdt];
1113  dist = new Double_t[fNdt];
1114  Double_t xd,yd;
1115  Int_t p, n, m;
1116  Bool_t o = kFALSE;
1117  for (it=0; it<fNdt; it++) {
1118  p = fPTried[it];
1119  n = fNTried[it];
1120  m = fMTried[it];
1121  xd = (fXN[p]+fXN[n]+fXN[m])/3;
1122  yd = (fYN[p]+fYN[n]+fYN[m])/3;
1123  if ((cp >= 0) && (sp >= 0.)) {
1124  dist[it] = -(fXNmax-xd+fYNmax-yd);
1125  } else if ((cp <= 0) && (sp >= 0.)) {
1126  dist[it] = -(fXNmax-xd+yd-fYNmin);
1127  o = kTRUE;
1128  } else if ((cp <= 0) && (sp <= 0.)) {
1129  dist[it] = -(xd-fXNmin+yd-fYNmin);
1130  } else {
1131  dist[it] = -(xd-fXNmin+fYNmax-yd);
1132  o = kTRUE;
1133  }
1134  }
1135  TMath::Sort(fNdt, dist, order, o);
1136 
1137  // Draw the triangles and markers if requested
1139  Int_t fs = fGraph2D->GetFillStyle();
1140  fGraph2D->SetFillStyle(1001);
1141  fGraph2D->TAttFill::Modify();
1143  fGraph2D->TAttLine::Modify();
1144  int lst = fGraph2D->GetLineStyle();
1145  for (it=0; it<fNdt; it++) {
1146  t[0] = fPTried[order[it]];
1147  t[1] = fNTried[order[it]];
1148  t[2] = fMTried[order[it]];
1149  for (Int_t k=0; k<3; k++) {
1150  if(fX[t[k]-1] < fXmin || fX[t[k]-1] > fXmax) goto endloop;
1151  if(fY[t[k]-1] < fYmin || fY[t[k]-1] > fYmax) goto endloop;
1152  temp1[0] = fX[t[k]-1];
1153  temp1[1] = fY[t[k]-1];
1154  temp1[2] = fZ[t[k]-1];
1155  temp1[0] = TMath::Max(temp1[0],fXmin);
1156  temp1[1] = TMath::Max(temp1[1],fYmin);
1157  temp1[2] = TMath::Max(temp1[2],fZmin);
1158  temp1[2] = TMath::Min(temp1[2],fZmax);
1159  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
1160  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
1161  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
1162  view->WCtoNDC(temp1, &temp2[0]);
1163  x[k] = temp2[0];
1164  y[k] = temp2[1];
1165  }
1166  x[3] = x[0];
1167  y[3] = y[0];
1168  if (tri1 || tri2) PaintLevels(t,x,y);
1169  if (!tri1 && !tri2 && !wire) {
1170  gPad->PaintFillArea(3,x,y);
1171  PaintLevels(t,x,y,nblev,glev);
1172  }
1173  if (!tri2) gPad->PaintPolyLine(4,x,y);
1174  if (markers) {
1175  if (markers0) {
1176  PaintPolyMarker0(3,x,y);
1177  } else {
1181  fGraph2D->TAttMarker::Modify();
1182  gPad->PaintPolyMarker(3,x,y);
1183  }
1184  }
1185 endloop:
1186  continue;
1187  }
1188  fGraph2D->SetFillStyle(fs);
1189  fGraph2D->SetLineStyle(lst);
1190  fGraph2D->TAttLine::Modify();
1191  fGraph2D->TAttFill::Modify();
1192  delete [] order;
1193  delete [] dist;
1194  if (glev) delete [] glev;
1195 }
1196 
1197 
1198 /////////
1199 // Paints the 2D graph as triangles (new implementation)
1200 /////
1202 {
1203 
1204  Double_t x[4], y[4], temp1[3],temp2[3];
1205  Int_t p[3];
1206 
1207  TView *view = gPad->GetView();
1208  if (!view) {
1209  Error("PaintTriangles", "No TView in current pad");
1210  return;
1211  }
1212 
1213  TString opt = option;
1214  opt.ToLower();
1215  Bool_t tri1 = opt.Contains("tri1");
1216  Bool_t tri2 = opt.Contains("tri2");
1217  Bool_t markers = opt.Contains("p");
1218  Bool_t markers0 = opt.Contains("p0");
1219  Bool_t wire = opt.Contains("w");
1220 
1221  // Define the grid levels drawn on the triangles.
1222  // The grid levels are aligned on the Z axis' main tick marks.
1223  Int_t nblev=0;
1224  Double_t *glev=0;
1225  if (!tri1 && !tri2 && !wire) {
1226  Int_t ndivz = gCurrentHist->GetZaxis()->GetNdivisions()%100;
1227  Int_t nbins;
1228  Double_t binLow = 0, binHigh = 0, binWidth = 0;
1229 
1230  // Find the main tick marks positions.
1231  Double_t *r0 = view->GetRmin();
1232  Double_t *r1 = view->GetRmax();
1233  if (!r0 || !r1) return;
1234 
1235  if (ndivz > 0) {
1236  THLimitsFinder::Optimize(r0[2], r1[2], ndivz,
1237  binLow, binHigh, nbins, binWidth, " ");
1238  } else {
1239  nbins = TMath::Abs(ndivz);
1240  binLow = r0[2];
1241  binHigh = r1[2];
1242  binWidth = (binHigh-binLow)/nbins;
1243  }
1244  // Define the grid levels
1245  nblev = nbins+1;
1246  glev = new Double_t[nblev];
1247  for (Int_t i = 0; i < nblev; ++i) glev[i] = binLow+i*binWidth;
1248  }
1249 
1250  // Initialize the levels on the Z axis
1251  if (tri1 || tri2) {
1252  Int_t ndiv = gCurrentHist->GetContour();
1253  if (ndiv == 0 ) {
1254  ndiv = gStyle->GetNumberContours();
1255  gCurrentHist->SetContour(ndiv);
1256  }
1258  }
1259 
1260  // For each triangle, compute the distance between the triangle centre
1261  // and the back planes. Then these distances are sorted in order to draw
1262  // the triangles from back to front.
1263  if (!fNdt) FindTriangles();
1264  Double_t cp = TMath::Cos(view->GetLongitude()*TMath::Pi()/180.);
1265  Double_t sp = TMath::Sin(view->GetLongitude()*TMath::Pi()/180.);
1266 
1267  Bool_t reverse = kFALSE;
1268  std::function<Double_t(Double_t, Double_t)> fDist;
1269 
1270  if ((cp >= 0) && (sp >= 0.)) {
1271  fDist = [&](Double_t xd, Double_t yd) -> Double_t { return -(fXNmax-xd+fYNmax-yd);};
1272  } else if ((cp <= 0) && (sp >= 0.)) {
1273  fDist = [&](Double_t xd, Double_t yd) -> Double_t { return -(fXNmax-xd+yd-fYNmin);};
1274  reverse = kTRUE;
1275  } else if ((cp <= 0) && (sp <= 0.)) {
1276  fDist = [&](Double_t xd, Double_t yd) -> Double_t { return -(xd-fXNmin+yd-fYNmin);};
1277  } else {
1278  fDist = [&](Double_t xd, Double_t yd) -> Double_t { return -(xd-fXNmin+fYNmax-yd);};
1279  reverse = kTRUE;
1280  }
1281 
1282  typedef std::pair<Double_t, TGraphDelaunay2D::Triangles::const_iterator> DistEntry;
1283  std::vector<DistEntry> dist;
1284  for(auto it = fDelaunay2D->begin(); it != fDelaunay2D->end(); ++it){
1285  auto face = *it;
1286  Double_t xd = (face.x[0] + face.x[1] + face.x[2]) / 3;
1287  Double_t yd = (face.y[0] + face.y[1] + face.y[2]) / 3;
1288 
1289  dist.emplace_back(fDist(xd, yd), it);
1290  }
1291 
1292  std::sort(dist.begin(), dist.end(),
1293  [&](const DistEntry & a, const DistEntry & b){ return !reverse ? (a.first < b.first) : (b.first < a .first); });
1294 
1295  // Draw the triangles and markers if requested
1297  Int_t fs = fGraph2D->GetFillStyle();
1298  fGraph2D->SetFillStyle(1001);
1299  fGraph2D->TAttFill::Modify();
1301  fGraph2D->TAttLine::Modify();
1302  int lst = fGraph2D->GetLineStyle();
1303  for (const auto & it : dist) {
1304  p[0] = it.second->idx[0];
1305  p[1] = it.second->idx[1];
1306  p[2] = it.second->idx[2];
1307  for (Int_t k=0; k<3; k++) {
1308  if(fX[p[k]] < fXmin || fX[p[k]] > fXmax) goto endloop;
1309  if(fY[p[k]] < fYmin || fY[p[k]] > fYmax) goto endloop;
1310  temp1[0] = fX[p[k]];
1311  temp1[1] = fY[p[k]];
1312  temp1[2] = fZ[p[k]];
1313  temp1[0] = TMath::Max(temp1[0],fXmin);
1314  temp1[1] = TMath::Max(temp1[1],fYmin);
1315  temp1[2] = TMath::Max(temp1[2],fZmin);
1316  temp1[2] = TMath::Min(temp1[2],fZmax);
1317  if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
1318  if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
1319  if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
1320  view->WCtoNDC(temp1, &temp2[0]);
1321  x[k] = temp2[0];
1322  y[k] = temp2[1];
1323  }
1324  x[3] = x[0];
1325  y[3] = y[0];
1326  if (tri1 || tri2) PaintLevels(p,x,y);
1327  if (!tri1 && !tri2 && !wire) {
1328  gPad->PaintFillArea(3,x,y);
1329  PaintLevels(p,x,y,nblev,glev);
1330  }
1331  if (!tri2) gPad->PaintPolyLine(4,x,y);
1332  if (markers) {
1333  if (markers0) {
1334  PaintPolyMarker0(3,x,y);
1335  } else {
1339  fGraph2D->TAttMarker::Modify();
1340  gPad->PaintPolyMarker(3,x,y);
1341  }
1342  }
1343 endloop:
1344  continue;
1345  }
1346  fGraph2D->SetFillStyle(fs);
1347  fGraph2D->SetLineStyle(lst);
1348  fGraph2D->TAttLine::Modify();
1349  fGraph2D->TAttFill::Modify();
1350 
1351  if (glev) delete [] glev;
1352 }
Int_t GetFirst() const
Return first bin on the axis i.e.
Definition: TAxis.cxx:429
virtual Style_t GetLineStyle() const
Definition: TAttLine.h:48
virtual Style_t GetFillStyle() const
Definition: TAttFill.h:44
virtual void SetLineWidth(Width_t lwidth)
Definition: TAttLine.h:57
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:48
Int_t * fPTried
Equal to fDelaunay->fNdt.
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:487
ClassImp(TSeqCollection) Int_t TSeqCollection TIter next(this)
Return index of object in collection.
TLine * line
Histogram option structure.
Definition: Hoption.h:24
virtual void WCtoNDC(const Float_t *pw, Float_t *pn)=0
float Float_t
Definition: RtypesCore.h:53
return c
const char Option_t
Definition: RtypesCore.h:62
virtual Double_t * GetRmax()=0
virtual void SetContour(Int_t nlevels, const Double_t *levels=0)
Set the number and values of contour levels.
Definition: TH1.cxx:7863
R__EXTERN TStyle * gStyle
Definition: TStyle.h:423
Double_t * GetX() const
Definition: TGraph2D.h:128
int Logy
log scale in Y. Also set by histogram option
Definition: Hoption.h:68
Double_t GetZmin() const
Returns the Z minimum.
Definition: TGraph2D.cxx:1267
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
See TView3D.
Definition: TView.h:36
Int_t GetNumberContours() const
Definition: TStyle.h:249
void PaintPolyLine(Option_t *option)
Paints the 2D graph as PaintPolyLine.
virtual Double_t * GetEY() const
Definition: TGraph2D.h:132
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition: TAxis.cxx:489
virtual Double_t * GetEX() const
Definition: TGraph2D.h:131
Double_t fXNmax
Equal to fDelaunay->fXNmin.
Basic string class.
Definition: TString.h:137
void PaintPolyMarker0(Int_t n, Double_t *x, Double_t *y)
Paints a circle at each vertex. Each circle background is white.
Double_t GetYNmax() const
Double_t GetYNmax() const
void PaintLevels(Int_t *v, Double_t *x, Double_t *y, Int_t nblev=0, Double_t *glev=0)
Paints one triangle.
R__EXTERN Hoption_t Hoption
Double_t GetXNmin() const
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:170
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1075
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
Double_t GetMinimum() const
Definition: TGraph2D.h:123
const Bool_t kFALSE
Definition: Rtypes.h:92
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width.
Definition: TAxis.cxx:511
virtual void SetFillStyle(Style_t fstyle)
Definition: TAttFill.h:52
int nbins[3]
void Paint(Option_t *option)
Paint a TGraphDelaunay according to the value of "option":
int Logx
log scale in X. Also set by histogram option
Definition: Hoption.h:67
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
Double_t GetXNmax() const
TSocket * s1
Definition: hserv2.C:36
TGraphDelaunay2D * fDelaunay2D
virtual Int_t GetContour(Double_t *levels=0)
Return contour values into array levels if pointer levels is non zero.
Definition: TH1.cxx:7787
static const double x2[5]
void PaintTriangles_new(Option_t *option)
Double_t x[n]
Definition: legend1.C:17
virtual void Paint(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:1907
Double_t * fXN
Pointer to fGraph2D->fZ.
Int_t fNdt
Equal to fGraph2D->fNpoints.
Double_t * fY
Pointer to fGraph2D->fX.
Int_t GetNdt() const
Double_t * fYN
Pointer to fDelaunay->fXN.
void Class()
Definition: Class.C:29
Double_t Log10(Double_t x)
Definition: TMath.h:529
static double p2(double t, double a, double b, double c)
Double_t GetZmax() const
Returns the Z maximum.
Definition: TGraph2D.cxx:1256
virtual void SetMarkerColor(Color_t mcolor=1)
Definition: TAttMarker.h:51
Triangles::const_iterator end() const
Int_t GetNumberOfColors() const
Return number of colors in the color palette.
Definition: TStyle.cxx:796
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Definition: TMath.h:1002
Double_t * GetYN() const
Double_t * fZ
Pointer to fGraph2D->fY.
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:918
Int_t * GetMTried() const
Triangles::const_iterator begin() const
virtual Double_t * GetEZ() const
Definition: TGraph2D.h:133
Double_t GetXNmax() const
A doubly linked list.
Definition: TList.h:47
Double_t * GetXN() const
virtual void SetLineColor(Color_t lcolor)
Definition: TAttLine.h:54
Double_t fYmax
fGraph2D->fHistogram limits
virtual Size_t GetMarkerSize() const
Definition: TAttMarker.h:46
TThread * t[5]
Definition: threadsh1.C:13
Class to manage histogram axis.
Definition: TAxis.h:36
TGraph2D * fGraph2D
virtual void SetFillColor(Color_t fcolor)
Definition: TAttFill.h:50
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
Definition: TAxis.cxx:499
unsigned int r1[N_CITIES]
Definition: simanTSP.cxx:321
virtual Color_t GetFillColor() const
Definition: TAttFill.h:43
void PaintTriangles(Option_t *option)
Paints the 2D graph as triangles.
Double_t * GetZ() const
Definition: TGraph2D.h:130
Bool_t TestBit(UInt_t f) const
Definition: TObject.h:173
TMarker * m
Definition: textangle.C:8
Int_t * GetPTried() const
bool first
Definition: line3Dfit.C:48
Double_t fXNmin
Pointer to fGraph2D->fZE.
TLine * l
Definition: textangle.C:4
virtual void SetMarkerStyle(Style_t mstyle=1)
Definition: TAttMarker.h:53
TAxis * GetYaxis()
Definition: TH1.h:320
static double p1(double t, double a, double b)
void FindTriangles()
Find triangles in fDelaunay and initialise the TGraph2DPainter values needed to paint triangles or fi...
virtual Color_t GetLineColor() const
Definition: TAttLine.h:47
virtual void SetMarkerSize(Size_t msize=1)
Definition: TAttMarker.h:54
virtual Double_t * GetRmin()=0
virtual Double_t GetLongitude()=0
Double_t Cos(Double_t)
Definition: TMath.h:424
Int_t GetColorPalette(Int_t i) const
Return color number i in current palette.
Definition: TStyle.cxx:730
Double_t Pi()
Definition: TMath.h:44
void PaintErrors(Option_t *option)
Paints the 2D graph as error bars.
Int_t * fNTried
Pointer to fDelaunay->fPTried.
Color * colors
Definition: X3DBuffer.c:19
static const double x1[5]
#define ClassImp(name)
Definition: Rtypes.h:279
double Double_t
Definition: RtypesCore.h:55
Double_t fYNmin
Equal to fDelaunay->fXNmax.
Double_t fYNmax
Equal to fDelaunay->fYNmin.
Double_t * fEX
Pointer to fDelaunay->fYN.
void PaintContour(Option_t *option)
Paints the 2D graph as a contour plot.
TGraphDelaunay * fDelaunay
Pointer to fDelaunay->fMTried.
Double_t y[n]
Definition: legend1.C:17
virtual Int_t GetNdivisions() const
Definition: TAttAxis.h:50
TGraphDelaunay2D generates a Delaunay triangulation of a TGraph2D.
The TH1 histogram class.
Definition: TH1.h:80
tuple view
Definition: tornado.py:20
TList * GetContourList(Double_t contour)
Returns the X and Y graphs building a contour.
virtual void SetLineStyle(Style_t lstyle)
Definition: TAttLine.h:56
Int_t GetNdt() const
TGraph2D * GetGraph2D() const
Int_t GetLast() const
Return last bin on the axis i.e.
Definition: TAxis.cxx:440
TAxis * GetZaxis()
Definition: TH1.h:321
Mother of all ROOT objects.
Definition: TObject.h:58
TGraph2D * GetGraph2D() const
#define R__EXTERN
Definition: DllImport.h:27
virtual Color_t GetMarkerColor() const
Definition: TAttMarker.h:44
virtual void SetPoint(Int_t i, Double_t x, Double_t y)
Set x and y values for point number i.
Definition: TGraph.cxx:2127
Int_t * fMTried
Pointer to fDelaunay->fNTried.
virtual Int_t FindFixBin(Double_t x) const
Find bin number corresponding to abscissa x.
Definition: TAxis.cxx:390
Double_t * fEY
Pointer to fGraph2D->fXE.
virtual void Add(TObject *obj)
Definition: TList.h:81
void PaintTriangles_old(Option_t *option)
Int_t GetN() const
Definition: TGraph2D.h:127
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:202
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:567
A Graph is a graphics object made of two arrays X and Y with npoints each.
Definition: TGraph.h:53
Double_t Sin(Double_t)
Definition: TMath.h:421
void PaintPolyMarker(Option_t *option)
Paints the 2D graph as PaintPolyMarker.
Double_t GetYNmin() const
Double_t * GetY() const
Definition: TGraph2D.h:129
R__EXTERN TH1 * gCurrentHist
tuple xs1
Definition: hsum.py:46
#define gPad
Definition: TVirtualPad.h:288
Double_t * fEZ
Pointer to fGraph2D->fYE.
The TGraphDelaunay painting class.
TSocket * s0
Definition: hserv2.C:36
Double_t GetXNmin() const
void FindAllTriangles()
Attempt to find all the Delaunay triangles of the point set.
virtual Style_t GetMarkerStyle() const
Definition: TAttMarker.h:45
const Bool_t kTRUE
Definition: Rtypes.h:91
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:7921
virtual Width_t GetLineWidth() const
Definition: TAttLine.h:49
Double_t GetMaximum() const
Definition: TGraph2D.h:122
TObject * obj
const Int_t n
Definition: legend1.C:16
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:8006
virtual Double_t GetContourLevelPad(Int_t level) const
Return the value of contour number "level" in Pad coordinates ie: if the Pad is in log scale along Z ...
Definition: TH1.cxx:7818
int Logz
log scale in Z. Also set by histogram option
Definition: Hoption.h:69
TAxis * GetXaxis()
Definition: TH1.h:319
TGraphDelaunay generates a Delaunay triangulation of a TGraph2D.
Double_t GetYNmin() const
Int_t * GetNTried() const
Double_t fXmin
Equal to fDelaunay->fYNmax.
virtual ~TGraph2DPainter()
TGraph2DPainter destructor.