Logo ROOT   6.16/01
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 "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
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 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;
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;
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) {
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 outside the Z range.
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 Int_t *segUsed = new Int_t[fNdt];
369 for(i=0; i<fNdt; i++) segUsed[i]=kFALSE;
370
371 // Find all the graphs making the contour. There is two kind of graphs,
372 // either they are "opened" or they are "closed"
373
374 // Find the opened graphs
375 Double_t xc=0, yc=0, xnc=0, ync=0;
376 Bool_t findNew;
377 Bool_t s0, s1;
378 Int_t is, js;
379 for (is=0; is<nbSeg; is++) {
380 if (segUsed[is]) continue;
381 s0 = s1 = kFALSE;
382
383 // Find to which segment is is connected. It can be connected
384 // via 0, 1 or 2 vertices.
385 for (js=0; js<nbSeg; js++) {
386 if (is==js) continue;
387 if (xs0[is]==xs0[js] && ys0[is]==ys0[js]) s0 = kTRUE;
388 if (xs0[is]==xs1[js] && ys0[is]==ys1[js]) s0 = kTRUE;
389 if (xs1[is]==xs0[js] && ys1[is]==ys0[js]) s1 = kTRUE;
390 if (xs1[is]==xs1[js] && ys1[is]==ys1[js]) s1 = kTRUE;
391 }
392
393 // Segment is is alone, not connected. It is stored in the
394 // list and the next segment is examined.
395 if (!s0 && !s1) {
396 graph = new TGraph();
397 graph->SetPoint(npg,xs0[is],ys0[is]); npg++;
398 graph->SetPoint(npg,xs1[is],ys1[is]); npg++;
399 segUsed[is] = kTRUE;
400 list->Add(graph); npg = 0;
401 continue;
402 }
403
404 // Segment is is connected via 1 vertex only and can be considered
405 // as the starting point of an opened contour.
406 if (!s0 || !s1) {
407 // Find all the segments connected to segment is
408 graph = new TGraph();
409 if (s0) {xc = xs0[is]; yc = ys0[is]; xnc = xs1[is]; ync = ys1[is];}
410 if (s1) {xc = xs1[is]; yc = ys1[is]; xnc = xs0[is]; ync = ys0[is];}
411 graph->SetPoint(npg,xnc,ync); npg++;
412 segUsed[is] = kTRUE;
413 js = 0;
414L01:
415 findNew = kFALSE;
416 if (js < nbSeg && segUsed[js]) {
417 js++;
418 goto L01;
419 } else if (xc==xs0[js] && yc==ys0[js]) {
420 xc = xs1[js];
421 yc = ys1[js];
422 findNew = kTRUE;
423 } else if (xc==xs1[js] && yc==ys1[js]) {
424 xc = xs0[js];
425 yc = ys0[js];
426 findNew = kTRUE;
427 }
428 if (findNew) {
429 segUsed[js] = kTRUE;
430 graph->SetPoint(npg,xc,yc); npg++;
431 js = 0;
432 goto L01;
433 }
434 js++;
435 if (js<nbSeg) goto L01;
436 list->Add(graph); npg = 0;
437 }
438 }
439
440
441 // Find the closed graphs. At this point all the remaining graphs
442 // are closed. Any segment can be used to start the search.
443 for (is=0; is<nbSeg; is++) {
444 if (segUsed[is]) continue;
445
446 // Find all the segments connected to segment is
447 graph = new TGraph();
448 segUsed[is] = kTRUE;
449 xc = xs0[is];
450 yc = ys0[is];
451 js = 0;
452 graph->SetPoint(npg,xc,yc); npg++;
453L02:
454 findNew = kFALSE;
455 if (js < nbSeg && segUsed[js]) {
456 js++;
457 goto L02;
458 } else if (xc==xs0[js] && yc==ys0[js]) {
459 xc = xs1[js];
460 yc = ys1[js];
461 findNew = kTRUE;
462 } else if (xc==xs1[js] && yc==ys1[js]) {
463 xc = xs0[js];
464 yc = ys0[js];
465 findNew = kTRUE;
466 }
467 if (findNew) {
468 segUsed[js] = kTRUE;
469 graph->SetPoint(npg,xc,yc); npg++;
470 js = 0;
471 goto L02;
472 }
473 js++;
474 if (js<nbSeg) goto L02;
475 // Close the contour
476 graph->SetPoint(npg,xs0[is],ys0[is]); npg++;
477 list->Add(graph); npg = 0;
478 }
479
480 delete [] xs0;
481 delete [] ys0;
482 delete [] xs1;
483 delete [] ys1;
484 delete [] segUsed;
485 return list;
486}
487
488
489////////////////////////////////////////////////////////////////////////////////
490/// Paint a TGraphDelaunay according to the value of "option":
491///
492/// | Option | Description |
493/// |----------|---------------------------------------------------------------|
494/// | "TRI" | The Delaunay triangles are drawn using filled area. An hidden surface drawing technique is used. The surface is painted with the current fill area color. The edges of each triangles are painted with the current line color. |
495/// | "TRIW" | The Delaunay triangles are drawn as wire frame. |
496/// | "TRI1" | The Delaunay triangles are painted with color levels. The edges of each triangles are painted with the current line color. |
497/// | "TRI2" | the Delaunay triangles are painted with color levels. |
498/// | "P" | Draw a marker at each vertex. |
499/// | "P0" | Draw a circle at each vertex. Each circle background is white. |
500/// | "PCOL" | Draw a marker at each vertex. The color of each marker is defined according to its Z position. |
501/// | "CONT" | Draw contours. |
502/// | "LINE" | Draw a 3D polyline. |
503
505{
506 TString opt = option;
507 opt.ToLower();
508 Bool_t triangles = opt.Contains("tri") ||
509 opt.Contains("tri1") ||
510 opt.Contains("tri2");
511 if (opt.Contains("tri0")) triangles = kFALSE;
512
513 Bool_t markers = opt.Contains("p") && !triangles;
514 Bool_t contour = opt.Contains("cont");
515 Bool_t line = opt.Contains("line");
516 Bool_t err = opt.Contains("err");
517
518 fGraph2D->TAttLine::Modify();
519 fGraph2D->TAttFill::Modify();
520 fGraph2D->TAttMarker::Modify();
521
522 // Compute minimums and maximums
523 TAxis *xaxis = gCurrentHist->GetXaxis();
524 Int_t first = xaxis->GetFirst();
525 fXmin = xaxis->GetBinLowEdge(first);
526 if (Hoption.Logx && fXmin <= 0) fXmin = xaxis->GetBinUpEdge(xaxis->FindFixBin(0.01*xaxis->GetBinWidth(first)));
527 fXmax = xaxis->GetBinUpEdge(xaxis->GetLast());
528 TAxis *yaxis = gCurrentHist->GetYaxis();
529 first = yaxis->GetFirst();
530 fYmin = yaxis->GetBinLowEdge(first);
531 if (Hoption.Logy && fYmin <= 0) fYmin = yaxis->GetBinUpEdge(yaxis->FindFixBin(0.01*yaxis->GetBinWidth(first)));
532 fYmax = yaxis->GetBinUpEdge(yaxis->GetLast());
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();
557 }
558 Int_t ndivz = TMath::Abs(ndiv);
560
561 Int_t theColor;
562 TList *l;
563 TGraph *g;
564 TObject *obj;
565 Double_t c;
566
567 if (!fNdt) FindTriangles();
568
569 for (Int_t k=0; k<ndiv; k++) {
571 l = GetContourList(c);
572 TIter next(l);
573 while ((obj = next())) {
574 if(obj->InheritsFrom(TGraph::Class()) ) {
575 g=(TGraph*)obj;
576 g->SetLineWidth(fGraph2D->GetLineWidth());
577 g->SetLineStyle(fGraph2D->GetLineStyle());
578 theColor = Int_t((k+0.99)*Float_t(ncolors)/Float_t(ndivz));
579 g->SetLineColor(gStyle->GetColorPalette(theColor));
580 g->Paint("l");
581 }
582 }
583 }
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 = TMath::Min(fZmin, 0.);
720 if (Hoption.Logz && zmin <= 0) zmin = TMath::Min((Double_t)1, (Double_t)0.001*fGraph2D->GetZmax());
721 zmax = fZmax;
722 }
723
724 // Order along Z axis the points (xi,yi,zi) where "i" belongs to {0,1,2}
725 // After this z0 < z1 < z2
726 Int_t i0=0, i1=0, i2=0;
727 if (fZ[p[1]]<=z0) {z0=fZ[p[1]]; x0=x[1]; y0=y[1]; i0=1;}
728 if (fZ[p[1]]>z2) {z2=fZ[p[1]]; x2=x[1]; y2=y[1]; i2=1;}
729 if (fZ[p[2]]<=z0) {z0=fZ[p[2]]; x0=x[2]; y0=y[2]; i0=2;}
730 if (fZ[p[2]]>z2) {z2=fZ[p[2]]; x2=x[2]; y2=y[2]; i2=2;}
731 i1 = 3-i2-i0;
732 Double_t x1 = x[i1];
733 Double_t y1 = y[i1];
734 Double_t z1 = fZ[p[i1]];
735
736 if (z0>zmax) z0 = zmax;
737 if (z2>zmax) z2 = zmax;
738 if (z0<zmin) z0 = zmin;
739 if (z2<zmin) z2 = zmin;
740 if (z1>zmax) z1 = zmax;
741 if (z1<zmin) z1 = zmin;
742
743 if (Hoption.Logz) {
744 z0 = TMath::Log10(z0);
745 z1 = TMath::Log10(z1);
746 z2 = TMath::Log10(z2);
747 zmin = TMath::Log10(zmin);
748 zmax = TMath::Log10(zmax);
749 }
750
751 // zi = Z values of the stripe number i
752 // zip = Previous zi
753 Double_t zi=0, zip=0;
754
755 if (nblev <= 0) {
756 // Paint the colors levels
757
758 // Compute the color associated to z0 (theColor0) and z2 (theColor2)
759 ncolors = gStyle->GetNumberOfColors();
760 theColor0 = (Int_t)( ((z0-zmin)/(zmax-zmin))*(ncolors-1) );
761 theColor2 = (Int_t)( ((z2-zmin)/(zmax-zmin))*(ncolors-1) );
762
763 // The stripes drawn to fill the triangles may have up to 5 points
764 Double_t xp[5], yp[5];
765
766 // rl = Ratio between z0 and z2 (long)
767 // rs = Ratio between z0 and z1 or z1 and z2 (short)
768 Double_t rl,rs;
769
770 // ci = Color of the stripe number i
771 // npf = number of point needed to draw the current stripe
772 Int_t ci,npf;
773
774 fillColor = fGraph2D->GetFillColor();
775
776 // If the z0's color and z2's colors are the same, the whole triangle
777 // can be painted in one go.
778 if(theColor0 == theColor2) {
780 fGraph2D->TAttFill::Modify();
781 gPad->PaintFillArea(3,x,y);
782
783 // The triangle must be painted with several colors
784 } else {
785 for(ci=theColor0; ci<=theColor2; ci++) {
787 fGraph2D->TAttFill::Modify();
788 if (ci==theColor0) {
789 zi = (((ci+1)*(zmax-zmin))/(ncolors-1))+zmin;
790 xp[0] = x0;
791 yp[0] = y0;
792 rl = (zi-z0)/(z2-z0);
793 xp[1] = rl*(x2-x0)+x0;
794 yp[1] = rl*(y2-y0)+y0;
795 if (zi>=z1 || z0==z1) {
796 rs = (zi-z1)/(z2-z1);
797 xp[2] = rs*(x2-x1)+x1;
798 yp[2] = rs*(y2-y1)+y1;
799 xp[3] = x1;
800 yp[3] = y1;
801 npf = 4;
802 } else {
803 rs = (zi-z0)/(z1-z0);
804 xp[2] = rs*(x1-x0)+x0;
805 yp[2] = rs*(y1-y0)+y0;
806 npf = 3;
807 }
808 } else if (ci==theColor2) {
809 xp[0] = xp[1];
810 yp[0] = yp[1];
811 xp[1] = x2;
812 yp[1] = y2;
813 if (zi<z1 || z2==z1) {
814 xp[3] = xp[2];
815 yp[3] = yp[2];
816 xp[2] = x1;
817 yp[2] = y1;
818 npf = 4;
819 } else {
820 npf = 3;
821 }
822 } else {
823 zi = (((ci+1)*(zmax-zmin))/(ncolors-1))+zmin;
824 xp[0] = xp[1];
825 yp[0] = yp[1];
826 rl = (zi-z0)/(z2-z0);
827 xp[1] = rl*(x2-x0)+x0;
828 yp[1] = rl*(y2-y0)+y0;
829 if ( zi>=z1 && zip<=z1) {
830 xp[3] = x1;
831 yp[3] = y1;
832 xp[4] = xp[2];
833 yp[4] = yp[2];
834 npf = 5;
835 } else {
836 xp[3] = xp[2];
837 yp[3] = yp[2];
838 npf = 4;
839 }
840 if (zi<z1) {
841 rs = (zi-z0)/(z1-z0);
842 xp[2] = rs*(x1-x0)+x0;
843 yp[2] = rs*(y1-y0)+y0;
844 } else {
845 rs = (zi-z1)/(z2-z1);
846 xp[2] = rs*(x2-x1)+x1;
847 yp[2] = rs*(y2-y1)+y1;
848 }
849 }
850 zip = zi;
851 // Paint a stripe
852 gPad->PaintFillArea(npf,xp,yp);
853 }
854 }
855 fGraph2D->SetFillColor(fillColor);
856 fGraph2D->TAttFill::Modify();
857
858 } else {
859 // Paint the grid levels
861 fGraph2D->TAttLine::Modify();
862 for(i=0; i<nblev; i++){
863 zl=glev[i];
864 if(zl >= z0 && zl <=z2) {
865 r21=(zl-z1)/(z2-z1);
866 r20=(zl-z0)/(z2-z0);
867 r10=(zl-z0)/(z1-z0);
868 xl[0]=r20*(x2-x0)+x0;
869 yl[0]=r20*(y2-y0)+y0;
870 if(zl >= z1 && zl <=z2) {
871 xl[1]=r21*(x2-x1)+x1;
872 yl[1]=r21*(y2-y1)+y1;
873 } else {
874 xl[1]=r10*(x1-x0)+x0;
875 yl[1]=r10*(y1-y0)+y0;
876 }
877 gPad->PaintPolyLine(2,xl,yl);
878 }
879 }
881 fGraph2D->TAttLine::Modify();
882 }
883}
884
885
886////////////////////////////////////////////////////////////////////////////////
887/// Paints the 2D graph as PaintPolyMarker
888
890{
891 Double_t temp1[3],temp2[3];
892
893 TView *view = gPad->GetView();
894 if (!view) {
895 Error("PaintPolyMarker", "No TView in current pad");
896 return;
897 }
898
899 TString opt = option;
900 opt.ToLower();
901 Bool_t markers0 = opt.Contains("p0");
902 Bool_t colors = opt.Contains("pcol");
903 Int_t ncolors = gStyle->GetNumberOfColors();
904 Int_t it, theColor;
905
906 // Initialize the levels on the Z axis
907 if (colors) {
908 Int_t ndiv = gCurrentHist->GetContour();
909 if (ndiv == 0 ) {
910 ndiv = gStyle->GetNumberContours();
912 }
914 }
915
916 Double_t *xm = new Double_t[fNpoints];
917 Double_t *ym = new Double_t[fNpoints];
918 Double_t *zm = new Double_t[fNpoints];
921
922 // min and max for colors
923 Double_t hzmincol = hzmin;
924 Double_t hzmaxcol = hzmax;
925 if (hzmincol==-1111 && hzmaxcol==-1111) {
926 hzmincol = TMath::Min(hzmincol, 0.);
927 if (Hoption.Logz && hzmincol <= 0) hzmincol = TMath::Min((Double_t)1, (Double_t)0.001*fGraph2D->GetZmax());
928 hzmaxcol = fZmax;
929 }
930 if (Hoption.Logz) {
931 hzmincol = TMath::Log10(hzmincol);
932 hzmaxcol = TMath::Log10(hzmaxcol);
933 }
934
935 Double_t Xeps = (fXmax-fXmin)*0.0001;
936 Double_t Yeps = (fYmax-fYmin)*0.0001;
937 Double_t Zeps = (hzmax-hzmin)*0.0001;
938
939 Int_t npd = 0;
940 for (it=0; it<fNpoints; it++) {
941 xm[it] = 0;
942 ym[it] = 0;
943 if(fXmin - fX[it] > Xeps || fX[it] - fXmax > Xeps) continue;
944 if(fYmin - fY[it] > Yeps || fY[it] - fYmax > Yeps) continue;
945 if(hzmin - fZ[it] > Zeps || fZ[it] - hzmax > Zeps) continue;
946 temp1[0] = fX[it];
947 temp1[1] = fY[it];
948 temp1[2] = fZ[it];
949 temp1[0] = TMath::Max(temp1[0],fXmin);
950 temp1[1] = TMath::Max(temp1[1],fYmin);
951 temp1[2] = TMath::Max(temp1[2],hzmin);
952 temp1[2] = TMath::Min(temp1[2],hzmax);
953 if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
954 if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
955 if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
956 view->WCtoNDC(temp1, &temp2[0]);
957 xm[npd] = temp2[0];
958 ym[npd] = temp2[1];
959 zm[npd] = temp1[2];
960 npd++;
961 }
962 if (markers0) {
963 PaintPolyMarker0(npd,xm,ym);
964 } else if (colors) {
965 Int_t cols = fGraph2D->GetMarkerColor();
966 for (it=0; it<npd; it++) {
967 theColor = (Int_t)( ((zm[it]-hzmincol)/(hzmaxcol-hzmincol))*(ncolors-1) );
969 fGraph2D->TAttMarker::Modify();
970 gPad->PaintPolyMarker(1,&xm[it],&ym[it]);
971 }
973 } else {
977 fGraph2D->TAttMarker::Modify();
978 gPad->PaintPolyMarker(npd,xm,ym);
979 }
980 delete [] xm;
981 delete [] ym;
982 delete [] zm;
983}
984
985
986////////////////////////////////////////////////////////////////////////////////
987/// Paints the 2D graph as PaintPolyLine
988
990{
991 Double_t temp1[3],temp2[3];
992
993 TView *view = gPad->GetView();
994 if (!view) {
995 Error("PaintPolyLine", "No TView in current pad");
996 return;
997 }
998
999 Int_t it;
1000 Double_t Xeps = (fXmax - fXmin) * 0.0001;
1001 Double_t Yeps = (fYmax - fYmin) * 0.0001;
1002
1003 Double_t *xm = new Double_t[fNpoints];
1004 Double_t *ym = new Double_t[fNpoints];
1005 Int_t npd = 0;
1006
1007 for (it=0; it<fNpoints; it++) {
1008 if (fXmin - fX[it] > Xeps || fX[it] - fXmax > Xeps)
1009 continue;
1010 if (fYmin - fY[it] > Yeps || fY[it] - fYmax > Yeps)
1011 continue;
1012 npd++;
1013 temp1[0] = fX[it];
1014 temp1[1] = fY[it];
1015 temp1[2] = fZ[it];
1016 temp1[0] = TMath::Max(temp1[0],fXmin);
1017 temp1[1] = TMath::Max(temp1[1],fYmin);
1018 temp1[2] = TMath::Max(temp1[2],fZmin);
1019 temp1[2] = TMath::Min(temp1[2],fZmax);
1020 if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
1021 if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
1022 if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
1023 view->WCtoNDC(temp1, &temp2[0]);
1024 xm[npd - 1] = temp2[0];
1025 ym[npd - 1] = temp2[1];
1026 }
1030 fGraph2D->TAttLine::Modify();
1031 gPad->PaintPolyLine(npd,xm,ym);
1032 delete [] xm;
1033 delete [] ym;
1034}
1035
1036
1037////////////////////////////////////////////////////////////////////////////////
1038/// Paints a circle at each vertex. Each circle background is white.
1039
1041{
1045 for (Int_t i=0; i<n; i++) {
1048 fGraph2D->TAttMarker::Modify();
1049 gPad->PaintPolyMarker(1,&x[i],&y[i]);
1052 fGraph2D->TAttMarker::Modify();
1053 gPad->PaintPolyMarker(1,&x[i],&y[i]);
1054 }
1056}
1057
1058
1059////////////////////////////////////////////////////////////////////////////////
1060/// Paints the 2D graph as triangles
1061
1063{
1064 if (fDelaunay)
1065 PaintTriangles_old(option);
1066 else if (fDelaunay2D)
1067 PaintTriangles_new(option);
1068}
1069
1070
1071////////////////////////////////////////////////////////////////////////////////
1072/// Paints the 2D graph as triangles (old implementation)
1073
1075{
1076 Double_t x[4], y[4], temp1[3],temp2[3];
1077 Int_t it,t[3];
1078 Int_t *order = 0;
1079 Double_t *dist = 0;
1080
1081 TView *view = gPad->GetView();
1082 if (!view) {
1083 Error("PaintTriangles", "No TView in current pad");
1084 return;
1085 }
1086
1087 TString opt = option;
1088 opt.ToLower();
1089 Bool_t tri1 = opt.Contains("tri1");
1090 Bool_t tri2 = opt.Contains("tri2");
1091 Bool_t markers = opt.Contains("p");
1092 Bool_t markers0 = opt.Contains("p0");
1093 Bool_t wire = opt.Contains("w");
1094
1095 // Define the grid levels drawn on the triangles.
1096 // The grid levels are aligned on the Z axis' main tick marks.
1097 Int_t nblev=0;
1098 Double_t *glev=0;
1099 if (!tri1 && !tri2 && !wire) {
1100 Int_t ndivz = gCurrentHist->GetZaxis()->GetNdivisions()%100;
1101 Int_t nbins;
1102 Double_t binLow = 0, binHigh = 0, binWidth = 0;
1103
1104 // Find the main tick marks positions.
1105 Double_t *r0 = view->GetRmin();
1106 Double_t *r1 = view->GetRmax();
1107 if (!r0 || !r1) return;
1108
1109 if (ndivz > 0) {
1110 THLimitsFinder::Optimize(r0[2], r1[2], ndivz,
1111 binLow, binHigh, nbins, binWidth, " ");
1112 } else {
1113 nbins = TMath::Abs(ndivz);
1114 binLow = r0[2];
1115 binHigh = r1[2];
1116 binWidth = (binHigh-binLow)/nbins;
1117 }
1118 // Define the grid levels
1119 nblev = nbins+1;
1120 glev = new Double_t[nblev];
1121 for (Int_t i = 0; i < nblev; ++i) glev[i] = binLow+i*binWidth;
1122 }
1123
1124 // Initialize the levels on the Z axis
1125 if (tri1 || tri2) {
1126 Int_t ndiv = gCurrentHist->GetContour();
1127 if (ndiv == 0 ) {
1128 ndiv = gStyle->GetNumberContours();
1129 gCurrentHist->SetContour(ndiv);
1130 }
1132 }
1133
1134 // For each triangle, compute the distance between the triangle centre
1135 // and the back planes. Then these distances are sorted in order to draw
1136 // the triangles from back to front.
1137 if (!fNdt) FindTriangles();
1138 Double_t cp = TMath::Cos(view->GetLongitude()*TMath::Pi()/180.);
1139 Double_t sp = TMath::Sin(view->GetLongitude()*TMath::Pi()/180.);
1140 order = new Int_t[fNdt];
1141 dist = new Double_t[fNdt];
1142 Double_t xd,yd;
1143 Int_t p, n, m;
1144 Bool_t o = kFALSE;
1145 for (it=0; it<fNdt; it++) {
1146 p = fPTried[it];
1147 n = fNTried[it];
1148 m = fMTried[it];
1149 xd = (fXN[p]+fXN[n]+fXN[m])/3;
1150 yd = (fYN[p]+fYN[n]+fYN[m])/3;
1151 if ((cp >= 0) && (sp >= 0.)) {
1152 dist[it] = -(fXNmax-xd+fYNmax-yd);
1153 } else if ((cp <= 0) && (sp >= 0.)) {
1154 dist[it] = -(fXNmax-xd+yd-fYNmin);
1155 o = kTRUE;
1156 } else if ((cp <= 0) && (sp <= 0.)) {
1157 dist[it] = -(xd-fXNmin+yd-fYNmin);
1158 } else {
1159 dist[it] = -(xd-fXNmin+fYNmax-yd);
1160 o = kTRUE;
1161 }
1162 }
1163 TMath::Sort(fNdt, dist, order, o);
1164
1165 // Draw the triangles and markers if requested
1167 Int_t fs = fGraph2D->GetFillStyle();
1168 fGraph2D->SetFillStyle(1001);
1169 fGraph2D->TAttFill::Modify();
1171 fGraph2D->TAttLine::Modify();
1172 int lst = fGraph2D->GetLineStyle();
1173 for (it=0; it<fNdt; it++) {
1174 t[0] = fPTried[order[it]];
1175 t[1] = fNTried[order[it]];
1176 t[2] = fMTried[order[it]];
1177 for (Int_t k=0; k<3; k++) {
1178 if(fX[t[k]-1] < fXmin || fX[t[k]-1] > fXmax) goto endloop;
1179 if(fY[t[k]-1] < fYmin || fY[t[k]-1] > fYmax) goto endloop;
1180 temp1[0] = fX[t[k]-1];
1181 temp1[1] = fY[t[k]-1];
1182 temp1[2] = fZ[t[k]-1];
1183 temp1[0] = TMath::Max(temp1[0],fXmin);
1184 temp1[1] = TMath::Max(temp1[1],fYmin);
1185 temp1[2] = TMath::Max(temp1[2],fZmin);
1186 temp1[2] = TMath::Min(temp1[2],fZmax);
1187 if (Hoption.Logx) temp1[0] = TMath::Log10(temp1[0]);
1188 if (Hoption.Logy) temp1[1] = TMath::Log10(temp1[1]);
1189 if (Hoption.Logz) temp1[2] = TMath::Log10(temp1[2]);
1190 view->WCtoNDC(temp1, &temp2[0]);
1191 x[k] = temp2[0];
1192 y[k] = temp2[1];
1193 }
1194 x[3] = x[0];
1195 y[3] = y[0];
1196 if (tri1 || tri2) PaintLevels(t,x,y);
1197 if (!tri1 && !tri2 && !wire) {
1198 gPad->PaintFillArea(3,x,y);
1199 PaintLevels(t,x,y,nblev,glev);
1200 }
1201 if (!tri2) gPad->PaintPolyLine(4,x,y);
1202 if (markers) {
1203 if (markers0) {
1204 PaintPolyMarker0(3,x,y);
1205 } else {
1209 fGraph2D->TAttMarker::Modify();
1210 gPad->PaintPolyMarker(3,x,y);
1211 }
1212 }
1213endloop:
1214 continue;
1215 }
1217 fGraph2D->SetLineStyle(lst);
1218 fGraph2D->TAttLine::Modify();
1219 fGraph2D->TAttFill::Modify();
1220 delete [] order;
1221 delete [] dist;
1222 if (glev) delete [] glev;
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->GetView();
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 Double_t *glev=0;
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 = new Double_t[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
1325 Int_t fs = fGraph2D->GetFillStyle();
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);
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
1379 if (glev) delete [] glev;
1380}
void Class()
Definition: Class.C:29
#define R__EXTERN
Definition: DllImport.h:27
#define b(i)
Definition: RSha256.hxx:100
#define c(i)
Definition: RSha256.hxx:101
#define g(i)
Definition: RSha256.hxx:105
#define s0(x)
Definition: RSha256.hxx:90
#define s1(x)
Definition: RSha256.hxx:91
static double p1(double t, double a, double b)
static double p2(double t, double a, double b, double c)
static const double x2[5]
static const double x1[5]
int Int_t
Definition: RtypesCore.h:41
const Bool_t kFALSE
Definition: RtypesCore.h:88
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
float Float_t
Definition: RtypesCore.h:53
const Bool_t kTRUE
Definition: RtypesCore.h:87
const char Option_t
Definition: RtypesCore.h:62
#define ClassImp(name)
Definition: Rtypes.h:363
R__EXTERN TH1 * gCurrentHist
R__EXTERN Hoption_t Hoption
R__EXTERN TStyle * gStyle
Definition: TStyle.h:406
#define gPad
Definition: TVirtualPad.h:286
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:41
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:504
virtual Int_t FindFixBin(Double_t x) const
Find bin number corresponding to abscissa x.
Definition: TAxis.cxx:405
Int_t GetLast() const
Return last bin on the axis i.e.
Definition: TAxis.cxx:455
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width.
Definition: TAxis.cxx:526
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
Definition: TAxis.cxx:514
Int_t GetFirst() const
Return first bin on the axis i.e.
Definition: TAxis.cxx:444
The TGraphDelaunay painting class.
void PaintContour(Option_t *option)
Paints the 2D graph as a contour plot.
void Paint(Option_t *option)
Paint a TGraphDelaunay according to the value of "option":
Int_t * fNTried
Pointer to fDelaunay->fPTried.
Double_t * fYN
Pointer to fDelaunay->fXN.
TGraphDelaunay * fDelaunay
Pointer to fDelaunay->fMTried.
Double_t * fZ
Pointer to fGraph2D->fY.
TGraphDelaunay2D * fDelaunay2D
virtual ~TGraph2DPainter()
TGraph2DPainter destructor.
TGraph2D * fGraph2D
void PaintErrors(Option_t *option)
Paints the 2D graph as error bars.
Double_t * fXN
Pointer to fGraph2D->fZ.
Int_t * fMTried
Pointer to fDelaunay->fNTried.
void FindTriangles()
Find triangles in fDelaunay and initialise the TGraph2DPainter values needed to paint triangles or fi...
Int_t * fPTried
Equal to fDelaunay->fNdt.
Double_t fXNmin
Pointer to fGraph2D->fZE.
void PaintLevels(Int_t *v, Double_t *x, Double_t *y, Int_t nblev=0, Double_t *glev=0)
Paints one triangle.
Int_t fNdt
Equal to fGraph2D->fNpoints.
Double_t * fEZ
Pointer to fGraph2D->fYE.
TList * GetContourList(Double_t contour)
Returns the X and Y graphs building a contour.
Double_t fYmax
fGraph2D->fHistogram limits
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->fX.
Double_t fYNmin
Equal to fDelaunay->fXNmax.
Double_t * fEY
Pointer to fGraph2D->fXE.
void PaintPolyLine(Option_t *option)
Paints the 2D graph as PaintPolyLine.
TGraph2DPainter()
TGraph2DPainter default constructor.
void PaintTriangles_old(Option_t *option)
Paints the 2D graph as triangles (old implementation)
Double_t fYNmax
Equal to fDelaunay->fYNmin.
Double_t fXmin
Equal to fDelaunay->fYNmax.
Double_t * fEX
Pointer to fDelaunay->fYN.
Double_t fXNmax
Equal to fDelaunay->fXNmin.
void PaintTriangles_new(Option_t *option)
Paints the 2D graph as triangles (new implementation)
Double_t GetMaximum() const
Definition: TGraph2D.h:112
Double_t GetMinimum() const
Definition: TGraph2D.h:113
Double_t * GetY() const
Definition: TGraph2D.h:119
Double_t GetZmin() const
Returns the Z minimum.
Definition: TGraph2D.cxx:1257
Double_t * GetX() const
Definition: TGraph2D.h:118
virtual Double_t * GetEZ() const
Definition: TGraph2D.h:123
Double_t GetZmax() const
Returns the Z maximum.
Definition: TGraph2D.cxx:1246
virtual Double_t * GetEY() const
Definition: TGraph2D.h:122
Int_t GetN() const
Definition: TGraph2D.h:117
virtual Double_t * GetEX() const
Definition: TGraph2D.h:121
Double_t * GetZ() const
Definition: TGraph2D.h:120
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 Graph is a graphics object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
The TH1 histogram class.
Definition: TH1.h:56
TAxis * GetZaxis()
Definition: TH1.h:318
virtual Double_t GetContourLevelPad(Int_t level) const
Return the value of contour number "level" in Pad coordinates.
Definition: TH1.cxx:7770
@ kUserContour
user specified contour levels
Definition: TH1.h:161
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:316
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:7872
TAxis * GetYaxis()
Definition: TH1.h:317
virtual void SetContour(Int_t nlevels, const Double_t *levels=0)
Set the number and values of contour levels.
Definition: TH1.cxx:7813
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:7957
virtual Int_t GetContour(Double_t *levels=0)
Return contour values into array levels if pointer levels is non zero.
Definition: TH1.cxx:7741
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:44
virtual void Add(TObject *obj)
Definition: TList.h:87
Mother of all ROOT objects.
Definition: TObject.h:37
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:443
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:880
Basic string class.
Definition: TString.h:131
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1100
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:619
Int_t GetColorPalette(Int_t i) const
Return color number i in current palette.
Definition: TStyle.cxx:913
Int_t GetNumberOfColors() const
Return number of colors in the color palette.
Definition: TStyle.cxx:979
Int_t GetNumberContours() const
Definition: TStyle.h:228
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:151
static constexpr double ms
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
Double_t Cos(Double_t)
Definition: TMath.h:629
constexpr Double_t Pi()
Definition: TMath.h:38
Double_t Sin(Double_t)
Definition: TMath.h:625
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Definition: TMathBase.h:362
Double_t Log10(Double_t x)
Definition: TMath.h:752
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
Definition: first.py:1
Definition: graph.py:1
Histogram option structure.
Definition: Hoption.h:24
int Logx
log scale in X. Also set by histogram option
Definition: Hoption.h:66
int Logz
log scale in Z. Also set by histogram option
Definition: Hoption.h:68
int Logy
log scale in Y. Also set by histogram option
Definition: Hoption.h:67
auto * m
Definition: textangle.C:8
auto * l
Definition: textangle.C:4
auto * a
Definition: textangle.C:12