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