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