Logo ROOT  
Reference Guide
TGraphPolargram.cxx
Go to the documentation of this file.
1// @(#)root/graf:$Id$
2// Author: Sebastian Boser, Mathieu Demaret 02/02/06
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/** \class TGraphPolargram
13\ingroup BasicGraphics
14
15To draw polar axis
16
17TGraphPolargram draw the polar axis of the TGraphPolar.
18
19Example:
20
21Begin_Macro(source)
22{
23 TCanvas * CPol = new TCanvas("CPol","TGraphPolar Examples",500,500);
24
25 Double_t rmin=0;
26 Double_t rmax=TMath::Pi()*2;
27 Double_t r[1000];
28 Double_t theta[1000];
29
30 TF1 * fp1 = new TF1("fplot","cos(x)",rmin,rmax);
31 for (Int_t ipt = 0; ipt < 1000; ipt++) {
32 r[ipt] = ipt*(rmax-rmin)/1000+rmin;
33 theta[ipt] = fp1->Eval(r[ipt]);
34 }
35 TGraphPolar * grP1 = new TGraphPolar(1000,r,theta);
36 grP1->SetTitle("");
37 grP1->SetLineColor(2);
38 grP1->Draw("AOL");
39}
40End_Macro
41*/
42
43#include "TGraphPolar.h"
44#include "TGraphPolargram.h"
45#include "TGaxis.h"
46#include "THLimitsFinder.h"
47#include "TVirtualPad.h"
48#include "TROOT.h"
49#include "TLatex.h"
50#include "TEllipse.h"
51#include "TMath.h"
52
54
55////////////////////////////////////////////////////////////////////////////////
56/// TGraphPolargram Constructor.
57
59 Double_t tmin, Double_t tmax):
60 TNamed(name,"Polargram")
61{
62 Init();
63 fNdivRad = 508;
64 fNdivPol = 508;
65 fPolarLabels = NULL;
66 fRwrmax = rmax;
67 fRwrmin = rmin;
68 fRwtmin = tmin;
69 fRwtmax = tmax;
70}
71
72////////////////////////////////////////////////////////////////////////////////
73/// Short constructor used in the case of a spider plot.
74
76 TNamed(name,"Polargram")
77{
78 Init();
79 fNdivRad = 0;
80 fNdivPol = 0;
81 fPolarLabels = NULL;
82 fRwrmax = 1;
83 fRwrmin = 0;
84 fRwtmax = 0;
85 fRwtmin = 0;
86}
87
88////////////////////////////////////////////////////////////////////////////////
89/// TGraphPolargram destructor.
90
92{
93 if (fPolarLabels != NULL) delete [] fPolarLabels;
94}
95
96////////////////////////////////////////////////////////////////////////////////
97/// Set the Polar range.
98/// \param[in] tmin the start number.
99/// \param[in] tmax the end number.
100
102{
103 if (tmin < tmax) {
104 fRwtmin = tmin;
105 fRwtmax = tmax;
106 }
107 if (gPad) gPad->Modified();
108}
109
110////////////////////////////////////////////////////////////////////////////////
111/// Everything within the circle belongs to the TGraphPolargram.
112
114{
115 Int_t i;
116 Double_t x = gPad->AbsPixeltoX(px);
117 Double_t y = gPad->AbsPixeltoY(py);
118
119 // Check if close to a (major) radial line.
121 Int_t div = (Int_t)rad*(fNdivRad%100);
122 Double_t dr = TMath::Min(TMath::Abs(rad-div*1./(fNdivRad%100)),
123 TMath::Abs(rad-(div+1)*1./(fNdivRad%100)));
124 Int_t drad = gPad->XtoPixel(dr)-gPad->XtoPixel(0);
125
126 // Check if close to a (major) Polar line.
127 // This is not a proper calculation, but rather fast.
128 Int_t dt = kMaxPixel;
129 for (i=0; i<(fNdivPol%100); i++) {
130 Double_t theta = i*2*TMath::Pi()/(fNdivPol%100);
131
132 // Attention: px,py in pixel units, line given in user coordinates.
133 Int_t dthis = DistancetoLine(px,py,0.,0.,TMath::Cos(theta),
134 TMath::Sin(theta));
135
136 // Fails if we are outside box described by the line.
137 // (i.e for all hor/vert lines)
138 if (dthis==9999) {
139
140 // Outside -> Get distance to endpoint of line.
141 if (rad>1) {
142 dthis = (Int_t)TMath::Sqrt(
143 TMath::Power(px-gPad->XtoPixel(TMath::Cos(theta)),2)+
144 TMath::Power(py-gPad->YtoPixel(TMath::Sin(theta)),2));
145 } else {
146
147 // Check for horizontal line.
148 if (((TMath::Abs(theta-TMath::Pi())<0.1) &&
149 ((px-gPad->XtoPixel(0))<0)) ||
150 ((TMath::Abs(theta)<0.1) &&
151 ((px-gPad->XtoPixel(0))>0))) {
152 dthis = TMath::Abs(py-gPad->YtoPixel(0.));
153 }
154
155 //Check for vertical line.
156 if (((TMath::Abs(theta-TMath::PiOver2())<0.1) &&
157 ((py-gPad->YtoPixel(0))>0)) ||
158 ((TMath::Abs(theta-3*TMath::PiOver2())<0.1) &&
159 (py-gPad->YtoPixel(0))<0)) {
160 dthis = TMath::Abs(px-gPad->XtoPixel(0.));
161 }
162 if (dthis==9999) {
163
164 // Inside, but out of box for nonorthogonal line ->
165 // get distance to start point.
166 dthis = (Int_t)TMath::Sqrt(
167 TMath::Power(px-gPad->XtoPixel(0.),2)+
168 TMath::Power(py-gPad->YtoPixel(0.),2));
169 }
170 }
171 }
172
173 // Take distance to closes line.
174 dt = TMath::Min(dthis,dt);
175 }
176 return TMath::Min(drad, dt);
177}
178
179////////////////////////////////////////////////////////////////////////////////
180/// Draw Polargram.
181
183{
184 Paint(options);
185 AppendPad(options);
186}
187
188////////////////////////////////////////////////////////////////////////////////
189/// Indicate that there is something to click here.
190
192{
193 if (!gPad) return;
194
195 Int_t kMaxDiff = 20;
196 static Int_t d1, d2, d3, px1, py1, px3, py3;
197 static Bool_t p1, p2, p3, p4, p5, p6, p7, p8;
198 Double_t px2, py2;
199 p2 = p3 = p4 = p5 = p6 = p7 = p8 = kFALSE;
200 if (!gPad->IsEditable()) return;
201 switch (event) {
202 case kMouseMotion:
203 px1 = gPad->XtoAbsPixel(TMath::Cos(GetAngle()));
204 py1 = gPad->YtoAbsPixel(TMath::Sin(GetAngle()));
205 d1 = TMath::Abs(px1 - px) + TMath::Abs(py1-py); //simply take sum of pixels differences
206 p1 = kFALSE;
207 px2 = gPad->XtoAbsPixel(-1);
208 py2 = gPad->YtoAbsPixel(1);
209 d2 = (Int_t)(TMath::Abs(px2 - px) + TMath::Abs(py2 - py)) ;
210 px3 = gPad->XtoAbsPixel(-1);
211 py3 = gPad->YtoAbsPixel(-1);
212 d3 = TMath::Abs(px3 - px) + TMath::Abs(py3 - py) ; //simply take sum of pixels differences
213 // check if point is close to the radial axis
214 if (d1 < kMaxDiff) {
215 gPad->SetCursor(kMove);
216 p1 = kTRUE;
217 }
218 // check if point is close to the left high axis
219 if ( d2 < kMaxDiff) {
220 gPad->SetCursor(kHand);
221 p7 = kTRUE;
222 }
223 // check if point is close to the left down axis
224 if ( d3 < kMaxDiff) {
225 gPad->SetCursor(kHand);
226 p8 = kTRUE;
227 }
228 // check if point is close to a main circle
229 if (!p1 && !p7 ) {
230 p6 = kTRUE;
231 gPad->SetCursor(kHand);
232 }
233 break;
234
235 case kButton1Down:
236 // Record initial coordinates
237 //px4 = px;
238 //py4 = py;
239
240 case kButton1Motion:
241 if (p1) {
242 px2 = gPad->AbsPixeltoX(px);
243 py2 = gPad->AbsPixeltoY(py);
244 if ( px2 < 0 && py2 < 0) {p2 = kTRUE;};
245 if ( px2 < 0 && py2 > 0 ) {p3 = kTRUE;};
246 if ( px2 > 0 && py2 > 0 ) {p4 = kTRUE;};
247 if ( px2 > 0 && py2 < 0 ) {p5 = kTRUE;};
248 px2 = TMath::ACos(TMath::Abs(px2));
249 py2 = TMath::ASin(TMath::Abs(py2));
250 if (p2) {
251 fAxisAngle = TMath::Pi()+(px2+py2)/2;
252 p2 = kFALSE;
253 };
254 if (p3) {
255 fAxisAngle = TMath::Pi()-(px2+py2)/2;
256 p3 = kFALSE;
257 };
258 if (p4) {
259 fAxisAngle = (px2+py2)/2;
260 p4 = kFALSE;
261 };
262 if (p5) {
263 fAxisAngle = -(px2+py2)/2;
264 p5 = kFALSE;
265 };
266 }
267 break;
268
269 case kButton1Up:
270 Paint();
271 }
272}
273
274////////////////////////////////////////////////////////////////////////////////
275/// Find the alignement rule to apply for TText::SetTextAlign(Short_t).
276
278{
280
281 while(angle < 0 || angle > 2*pi){
282 if(angle < 0) angle+=2*pi;
283 if(angle > 2*pi) angle-=2*pi;
284 }
286 if(angle > 0 && angle < pi/2) return 11;
287 else if(angle > pi/2 && angle < pi) return 31;
288 else if(angle > pi && angle < 3*pi/2) return 33;
289 else if(angle > 3*pi/2 && angle < 2*pi) return 13;
290 else if(angle == 0 || angle == 2*pi) return 12;
291 else if(angle == pi/2) return 21;
292 else if(angle == pi) return 32;
293 else if(angle == 3*pi/2) return 23;
294 else return 0;
295 }
296 else{
297 if(angle >= 0 && angle <= pi/2) return 12;
298 else if((angle > pi/2 && angle <= pi) || (angle > pi && angle <= 3*pi/2)) return 32;
299 else if(angle > 3*pi/2 && angle <= 2*pi) return 12;
300 else return 0;
301 }
302}
303
304////////////////////////////////////////////////////////////////////////////////
305/// Determine the orientation of the polar labels according to their angle.
306
308{
310 Double_t convraddeg = 180.0/pi;
311
312 while(angle < 0 || angle > 2*pi){
313 if(angle < 0) angle+=2*pi;
314 if(angle > 2*pi) angle-=2*pi;
315 }
316
317 if(angle >= 0 && angle <= pi/2) return angle*convraddeg;
318 else if(angle > pi/2 && angle <= pi) return (angle + pi)*convraddeg;
319 else if(angle > pi && angle <= 3*pi/2) return (angle - pi)*convraddeg;
320 else if(angle > 3*pi/2 && angle <= 2*pi) return angle*convraddeg;
321 else return 0;
322}
323
324////////////////////////////////////////////////////////////////////////////////
325/// Initialize some of the fields of TGraphPolargram.
326
328{
329 fAxisAngle = 0;
330 fCutRadial = 0;
331 fDegree = kFALSE;
332 fGrad = kFALSE;
333 fLineStyle = 3;
335 fPolarLabelFont = 62;
336 fPolarOffset = 0.04;
337 fPolarTextSize = 0.04;
338 fRadialOffset = 0.025;
339 fRadian = kTRUE;
341 fRadialLabelFont = 62;
342 fRadialTextSize = 0.035;
343 fTickpolarSize = 0.02;
344}
345
346////////////////////////////////////////////////////////////////////////////////
347/// Paint TGraphPolargram.
348
350{
351 Int_t optionpoldiv, optionraddiv;
352 Bool_t optionLabels = kTRUE;
353
354 TString opt = chopt;
355 opt.ToUpper();
356
357 if(opt.Contains('P')) optionpoldiv=1; else optionpoldiv=0;
358 if(opt.Contains('R')) optionraddiv=1; else optionraddiv=0;
361 if(!opt.Contains('P') && !opt.Contains('R')) optionpoldiv=optionraddiv=1;
362 if(opt.Contains('N')) optionLabels = kFALSE;
363
364 if(optionraddiv) PaintRadialDivisions(kTRUE);
366 if(optionpoldiv) PaintPolarDivisions(optionLabels);
367}
368
369////////////////////////////////////////////////////////////////////////////////
370/// This is simplified from TEllipse::PaintEllipse.
371/// Draw this ellipse with new coordinates.
372
374 Double_t phimin, Double_t phimax, Double_t theta)
375{
376 Int_t i;
377 const Int_t np = 200; // Number of point to draw circle
378 static Double_t x[np+3], y[np+3];
379
380 // Set number of points approximatively proportional to the ellipse
381 // circumference.
382
383 Double_t circ = TMath::Pi()*2*r*(phimax-phimin)/36;
384 Int_t n = (Int_t)(np*circ/((gPad->GetX2()-gPad->GetX1())+
385 (gPad->GetY2()-gPad->GetY1())));
386 if (n < 8) n = 8;
387 if (n > np) n = np;
388 Double_t angle,dx,dy;
389 Double_t dphi = (phimax-phimin)*TMath::Pi()/(180*n);
390 Double_t ct = TMath::Cos(TMath::Pi()*theta/180);
391 Double_t st = TMath::Sin(TMath::Pi()*theta/180);
392 for (i=0; i<=n; i++) {
393 angle = phimin*TMath::Pi()/180 + Double_t(i)*dphi;
394 dx = r*TMath::Cos(angle);
395 dy = r*TMath::Sin(angle);
396 x[i] = x1 + dx*ct - dy*st;
397 y[i] = y1 + dx*st + dy*ct;
398 }
399 gPad->PaintPolyLine(n+1,x,y);
400}
401
402////////////////////////////////////////////////////////////////////////////////
403/// Draw Polar divisions.
404/// Check for editable pad or create default.
405
407{
408 Int_t i, j, rnum, rden, first, last;
409 if (!gPad) return ;
410
411 gPad->RangeAxis(-1,-1,1,1);
412 gPad->Range(-1.25,-1.25,1.25,1.25);
413 Int_t ndivMajor = fNdivPol%100;
414 Int_t ndivMinor = fNdivPol/100;
415
416 if (!gPad->GetLogy()) {
417 for (i=0; i<ndivMajor; i++) {
418 Double_t txtval = fRwtmin + i*(fRwtmax-fRwtmin)/ndivMajor;
419 Double_t theta = i*2*TMath::Pi()/ndivMajor;
420 Double_t costheta = TMath::Cos(theta);
421 Double_t sintheta = TMath::Sin(theta);
422 Double_t tantheta = TMath::Tan(theta);
423 Double_t costhetas = (1+fPolarOffset)*costheta;
424 Double_t sinthetas = (1+fPolarOffset)*sintheta;
425 Double_t corr = 0.01;
426
427 TLatex *textangular = new TLatex();
428 textangular->SetTextColor(GetPolarColorLabel());
429 textangular->SetTextFont(GetPolarLabelFont());
430
431 const char* form = (char *)" ";
432 TGaxis axis;
434 // Polar numbers are aligned with their axis.
435 if(fPolarLabels == NULL && optionLabels){;
436 if (fRadian) {
437 // Radian case.
438 ReduceFraction(2*i, ndivMajor, rnum, rden); // Reduces the fraction.
439 if (rnum == 0) form = Form("%d",rnum);
440 if (rnum == 1 && rden == 1) form = Form("#pi");
441 if (rnum == 1 && rden != 1) form = Form("#frac{#pi}{%d}",rden);
442 if (rnum != 1 && rden == 1 && i !=0) form= Form("%d#pi",rnum);
443 if (rnum != 1 && rden != 1) form = Form("#frac{%d#pi}{%d}",rnum,rden);
444 textangular->SetTextAlign(FindAlign(theta));
445 textangular->PaintLatex(costhetas,
446 sinthetas, FindTextAngle(theta),
447 GetPolarLabelSize(), form);
448 } else {
449 // Any other cases: numbers are aligned with their axis.
450 form = Form("%5.3g",txtval);
451 axis.LabelsLimits(form,first,last);
452 TString s = Form("%s",form);
453 if (first != 0) s.Remove(0, first);
454 textangular->SetTextAlign(FindAlign(theta));
455 textangular->PaintLatex(costhetas,
456 sinthetas, FindTextAngle(theta),
458 }
459 } else if (fPolarLabels){
460 // print the specified polar labels
461 textangular->SetTextAlign(FindAlign(theta));
462 textangular->PaintLatex(costhetas,sinthetas,FindTextAngle(theta),
464 }
465 } else {
466 // Polar numbers are shown horizontally.
467 if(fPolarLabels == NULL && optionLabels){
468 if (fRadian) {
469 // Radian case
470 ReduceFraction(2*i, ndivMajor, rnum, rden);
471 if (rnum == 0) form = Form("%d",rnum);
472 if (rnum == 1 && rden == 1) form = Form("#pi");
473 if (rnum == 1 && rden != 1) form = Form("#frac{#pi}{%d}",rden);
474 if (rnum != 1 && rden == 1 && i !=0) form = Form("%d#pi",rnum);
475 if (rnum != 1 && rden != 1) form = Form("#frac{%d#pi}{%d}",rnum,rden);
476 if(theta >= 3*TMath::Pi()/12.0 && theta < 2*TMath::Pi()/3.0) corr=0.04;
477 textangular->SetTextAlign(FindAlign(theta));
478 textangular->PaintLatex(costhetas,corr+sinthetas,0,
479 GetPolarLabelSize(),form);
480 } else {
481 // Any other cases where numbers are shown horizontally.
482 form = Form("%5.3g",txtval);
483 axis.LabelsLimits(form,first,last);
484 TString s = Form("%s",form);
485 if (first != 0) s.Remove(0, first);
486 if(theta >= 3*TMath::Pi()/12.0 && theta < 2*TMath::Pi()/3.0) corr=0.04;
487 textangular->SetTextAlign(FindAlign(theta));
488 textangular->PaintLatex(costhetas, //j'ai efface des offset la
489 corr+sinthetas,0,GetPolarLabelSize(),s);
490 }
491 } else if (fPolarLabels) {
492 // print the specified polar labels
493 textangular->SetTextAlign(FindAlign(theta));
494 textangular->PaintText(costhetas,sinthetas,fPolarLabels[i]);
495 }
496 }
498 //Check if SetTickPolar is activated, and draw tick marks
499 Bool_t issettickpolar = gPad->GetTicky();
500
501 if (issettickpolar) {
502 if (theta != 0 && theta !=TMath::Pi()) {
503 gPad->PaintLine((sintheta-GetTickpolarSize())/tantheta,sintheta-GetTickpolarSize(),
504 (sintheta+GetTickpolarSize())/tantheta,sintheta+GetTickpolarSize());
505 }
506 if (theta == 0 || theta ==TMath::Pi()) {
507 gPad->PaintLine(1-GetTickpolarSize(),0,1+GetTickpolarSize(),0);
508 gPad->PaintLine(-1+GetTickpolarSize(),0,-1-GetTickpolarSize(),0);
509 }
510 }
513 gPad->PaintLine(0.,0.,costheta,sintheta);
514 delete textangular;
515 // Add minor lines w/o text.
516 Int_t oldLineStyle = GetLineStyle();
517 TAttLine::SetLineStyle(2); //Minor lines always in this style.
518 TAttLine::Modify(); //Changes line attributes apart from style.
519 for (j=1; j<ndivMinor; j++) {
520 Double_t thetamin = theta+j*2*TMath::Pi()/(ndivMajor*ndivMinor);
521 gPad->PaintLine(0.,0.,TMath::Cos(thetamin),TMath::Sin(thetamin));
522 }
523 TAttLine::SetLineStyle(oldLineStyle);
525 }
526 } else {
527 Int_t big = (Int_t)fRwtmax;
528 Int_t test= 1;
529 while (big >= 10) {
530 big = big/10;
531 test++;
532 }
533 for (i=1; i<=test; i++) {
534 Double_t txtval = pow((double)10,(double)(i-1));
535 Double_t theta = (i-1)*2*TMath::Pi()/(double)(test);
536 Double_t costheta = TMath::Cos(theta);
537 Double_t sintheta = TMath::Sin(theta);
538 Double_t tantheta = TMath::Tan(theta);
539 Double_t costhetas = (1+fPolarOffset)*costheta;
540 Double_t sinthetas = (1+fPolarOffset)*sintheta;
541 Double_t corr = 0.01;
542
543 TLatex *textangular = new TLatex();
544 textangular->SetTextColor(GetPolarColorLabel());
545 textangular->SetTextFont(GetPolarLabelFont());
546
547 const char* form = (char *)" ";
548 TGaxis axis;
549
551 if(fPolarLabels==NULL && optionLabels){
552 // Polar numbers are aligned with their axis.
553 form = Form("%5.3g",txtval);
554 axis.LabelsLimits(form,first,last);
555 TString s = Form("%s",form);
556 if (first != 0) s.Remove(0, first);
557 textangular->SetTextAlign(FindAlign(theta));
558 textangular->PaintLatex(costhetas,
559 sinthetas, FindTextAngle(theta), GetPolarLabelSize(), s);
560 }
561 else if (fPolarLabels){
562 // print the specified polar labels
563 textangular->SetTextAlign(FindAlign(theta));
564 textangular->PaintText(costhetas,sinthetas,fPolarLabels[i]);
565 }
566
567 } else {
568 if(fPolarLabels==NULL && optionLabels){
569 // Polar numbers are shown horizontally.
570 form = Form("%5.3g",txtval);
571 axis.LabelsLimits(form,first,last);
572 TString s = Form("%s",form);
573 if (first != 0) s.Remove(0, first);
574 if(theta >= 3*TMath::Pi()/12.0 && theta < 2*TMath::Pi()/3.0) corr=0.04;
575 textangular->SetTextAlign(FindAlign(theta));
576 textangular->PaintLatex(costhetas,
577 corr+sinthetas,0,GetPolarLabelSize(),s);
578 } else if (fPolarLabels){
579 // print the specified polar labels
580 textangular->SetTextAlign(FindAlign(theta));
581 textangular->PaintText(costhetas,sinthetas,fPolarLabels[i]);
582 }
583 }
584
586 //Check if SetTickPolar is activated, and draw tick marks
587 Bool_t issettickpolar = gPad->GetTicky();
588 if (issettickpolar) {
589 if (theta != 0 && theta !=TMath::Pi()) {
590 gPad->PaintLine((sintheta-GetTickpolarSize())/tantheta,sintheta-GetTickpolarSize(),
591 (sintheta+GetTickpolarSize())/tantheta,sintheta+GetTickpolarSize());
592 }
593 if (theta == 0 || theta ==TMath::Pi()) {
594 gPad->PaintLine(1-GetTickpolarSize(),0,1+GetTickpolarSize(),0);
595 gPad->PaintLine(-1+GetTickpolarSize(),0,-1-GetTickpolarSize(),0);
596 }
597 }
600 gPad->PaintLine(0.,0.,costheta,sintheta);
601 delete textangular;
602 // Add minor lines w/o text.
603 Int_t oldLineStyle = GetLineStyle();
604 TAttLine::SetLineStyle(2); //Minor lines always in this style.
605 TAttLine::Modify(); //Changes line attributes apart from style.
606 Double_t a=0;
607 Double_t b,c,d;
608 b = TMath::Log(10)*test;
609 d= 2*TMath::Pi()/(double)test;
610 for (j=1; j<9; j++) {
611 a=TMath::Log(j+1)-TMath::Log(j)+a;
612 c=a/b*6.28+d*(i-1);
613 gPad->PaintLine(0.,0.,TMath::Cos(c),TMath::Sin(c));
614 }
615 TAttLine::SetLineStyle(oldLineStyle);
617 }
618 }
619}
620
621////////////////////////////////////////////////////////////////////////////////
622/// Paint radial divisions.
623/// Check for editable pad or create default.
624
626{
627 static char chopt[8] = "";
628 Int_t i,j;
629 Int_t ndiv = TMath::Abs(fNdivRad);
630 Int_t ndivMajor = ndiv%100;
631 Int_t ndivMinor = ndiv/100;
632 Int_t ndivmajor = 0;
633 Double_t frwrmin = 0., frwrmax = 0., binWidth = 0;
634
635 THLimitsFinder::Optimize(fRwrmin,fRwrmax,ndivMajor,frwrmin,
636 frwrmax, ndivmajor,binWidth,"");
637
638 if (!gPad) return ;
639 if (!gPad->GetLogx()) {
640 gPad->RangeAxis(-1,-1,1,1);
641 gPad->Range(-1.25,-1.25,1.25,1.25);
642 Double_t umin = fRwrmin;
643 Double_t umax = fRwrmax;
644 Double_t rmajmin = (frwrmin-fRwrmin)/(fRwrmax-fRwrmin);
645 Double_t rmajmax = (frwrmax-fRwrmin)/(fRwrmax-fRwrmin);
646 Double_t dist = (rmajmax-rmajmin)/ndivmajor;
647 Int_t ndivminor = 0;
648
649 chopt[0] = 0;
650 strncat(chopt, "SDH", 4);
651 if (fNdivRad < 0) strncat(chopt, "N",2);
652 if(drawaxis){
653 // Paint axis.
654 TGaxis axis;
660 umin, umax, ndiv, chopt, 0., kFALSE);
661 }
662
663 // Paint Circles.
664 // First paint main circle.
665 PaintCircle(0.,0.,1,0.,360,0);
666 // Optimised case.
667 if (fNdivRad>0 ) {
668 Double_t frwrmini = 0., frwrmaxi = 0., binWidth2 =0;
669 THLimitsFinder::Optimize(frwrmin,frwrmin+binWidth,ndivMinor,frwrmini,
670 frwrmaxi, ndivminor,binWidth2,"");
671 Double_t dist2 = dist/(ndivminor);
672 // Paint major circles.
673 for (i=1; i<=ndivmajor+2; i++) {
676 PaintCircle(0.,0.,rmajmin,0.,360,0);
677
678 //Paint minor circles.
681 for (j=1; j<ndivminor+1; j++) {
682 if (rmajmin+j*dist2<=1) PaintCircle(0.,0.,rmajmin+j*dist2,0.,360,0);
683 }
684 rmajmin = (frwrmin-fRwrmin)/(fRwrmax-fRwrmin)+(i-1)*dist;
685 }
686 // Non-optimized case.
687 } else {
688
689 // Paint major circles.
690 for (i=1; i<=ndivMajor; i++) {
693 Double_t rmaj = i*1./ndivMajor;
694 PaintCircle(0.,0.,rmaj,0.,360,0);
695
696 // Paint minor circles.
697 for (j=1; j<ndivMinor; j++) {
700 PaintCircle(0.,0.,rmaj- j*1./(ndivMajor*ndivMinor),0.,360,0);
701 }
702 }
703 }
704 } else {
705 // Draw Log scale on radial axis if option activated.
706 Int_t big = (Int_t)fRwrmax;
707 Int_t test= 0;
708 while (big >= 10) {
709 big = big/10;
710 test++;
711 }
712 for (i=1; i<=test; i++) {
715 Double_t ecart;
716 ecart = ((double) i)/ ((double) test);
717 PaintCircle(0.,0.,ecart,0,360,0);
720 Double_t a=0;
721 Double_t b,c,d;
722 b = TMath::Log(10)*test;
723 d = 1/(double)test;
724 for (j=1; j<9; j++) {
725 a = TMath::Log(j+1)-TMath::Log(j)+a;
726 c = a/b+d*(i-1);
727 PaintCircle(0,0.,c,0.,360,0);
728 }
729 }
730 }
733}
734
735////////////////////////////////////////////////////////////////////////////////
736/// Reduce fractions.
737
739{
740 Int_t a = 0;
741 Int_t b = 0;
742 Int_t i = 0;
743 Int_t j = 0;
744 a = den;
745 b = num;
746 if (b > a) {
747 j = b;
748 } else {
749 j = a;
750 }
751 for (i=j; i > 1; i--) {
752 if ((a % i == 0) && (b % i == 0)) {
753 a = a/i;
754 b = b/i;
755 }
756 }
757 rden = a;
758 rnum = b;
759}
760
761////////////////////////////////////////////////////////////////////////////////
762/// Set axis angle.
763
765{
766 fAxisAngle = angle/180*TMath::Pi();
767}
768
769////////////////////////////////////////////////////////////////////////////////
770/// Set the number of Polar divisions: enter a number ij with 0<i<99 and 0<j<99
771/// - i sets the major Polar divisions.
772/// - j sets the minor Polar divisions.
773
775{
776 if (ndiv > 0)
777 fNdivPol = ndiv;
778 if (gPad) gPad->Modified();
779}
780
781////////////////////////////////////////////////////////////////////////////////
782/// Set the number of radial divisions: enter a number ij with 0<i<99 and 0<j<99
783/// - i sets the major radial divisions.
784/// - j sets the minor radial divisions.
785
787{
788 fNdivRad = ndiv;
789 if (gPad) gPad->Modified();
790}
791
792////////////////////////////////////////////////////////////////////////////////
793/// Set some specified polar labels, used in the case of a spider plot.
794
796{
797 if(fPolarLabels == NULL)
799 fPolarLabels[div]=label;
800 if (gPad) gPad->Modified();
801}
802
803////////////////////////////////////////////////////////////////////////////////
804/// Set Polar labels color.
805
807{
808 fPolarLabelColor = tcolorangular;
809}
810
811////////////////////////////////////////////////////////////////////////////////
812
814{
815 // Set Polar label font.
816
817 fPolarLabelFont = tfontangular;
818}
819
820////////////////////////////////////////////////////////////////////////////////
821/// Set angular labels size.
822
824{
825 fPolarTextSize = angularsize;
826}
827
828////////////////////////////////////////////////////////////////////////////////
829/// Set the labels offset.
830
832{
833 fPolarOffset = angularOffset;
834 if (gPad) gPad->Modified();
835}
836
837////////////////////////////////////////////////////////////////////////////////
838/// Set radial labels color.
839
841{
842 fRadialLabelColor = tcolorradial;
843}
844
845////////////////////////////////////////////////////////////////////////////////
846/// Set radial label font.
847
849{
850 fRadialLabelFont = tfontradial;
851}
852
853////////////////////////////////////////////////////////////////////////////////
854/// Set radial labels size.
855
857{
858 fRadialTextSize = radialsize;
859}
860
861////////////////////////////////////////////////////////////////////////////////
862/// Set the labels offset.
863
865{
866 fRadialOffset = radialOffset;
867 if (gPad) gPad->Modified();
868}
869
870////////////////////////////////////////////////////////////////////////////////
871/// Allows to change range Polar.
872/// \param[in] tmin the start number.
873/// \param[in] tmax the end number.
874
876{
877 fDegree = kFALSE;
878 fGrad = kFALSE;
879 fRadian = kFALSE;
880
881 if (tmin < tmax) {
882 fRwtmin = tmin;
883 fRwtmax = tmax;
884 }
885 if (gPad) gPad->Modified();
886}
887
888////////////////////////////////////////////////////////////////////////////////
889/// Set the radial range.
890/// \param[in] rmin radius at center of the circle.
891/// \param[in] rmax radius at the intersection of the right X axis part and the circle.
892
894{
895 if (rmin < rmax) {
896 fRwrmin = rmin;
897 fRwrmax = rmax;
898 }
899 if (gPad) gPad->Modified();
900}
901
902////////////////////////////////////////////////////////////////////////////////
903/// Set polar ticks size.
904
906{
907 fTickpolarSize = tickpolarsize;
908}
909
910////////////////////////////////////////////////////////////////////////////////
911/// The Polar circle is labelled using degree.
912
914{
915 fDegree = kTRUE;
916 fGrad = kFALSE;
917 fRadian = kFALSE;
918 ChangeRangePolar(0,360);
919}
920
921////////////////////////////////////////////////////////////////////////////////
922/// The Polar circle is labelled using gradian.
923
925{
926 fGrad = kTRUE;
927 fRadian = kFALSE;
928 fDegree = kFALSE;
929 ChangeRangePolar(0,200);
930}
931
932////////////////////////////////////////////////////////////////////////////////
933/// The Polar circle is labelled using radian.
934
936{
937 fRadian = kTRUE;
938 fGrad = kFALSE;
939 fDegree = kFALSE;
941}
942
943////////////////////////////////////////////////////////////////////////////////
944/// Set range from 0 to 2*pi
945
947{
949}
@ kMouseMotion
Definition: Buttons.h:23
@ kButton1Motion
Definition: Buttons.h:20
@ kButton1Up
Definition: Buttons.h:19
@ kButton1Down
Definition: Buttons.h:17
const Int_t kMaxPixel
Definition: GuiTypes.h:368
ROOT::R::TRInterface & r
Definition: Object.C:4
#define d(i)
Definition: RSha256.hxx:102
#define b(i)
Definition: RSha256.hxx:100
#define c(i)
Definition: RSha256.hxx:101
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
short Font_t
Definition: RtypesCore.h:75
double Double_t
Definition: RtypesCore.h:55
short Color_t
Definition: RtypesCore.h:79
const Bool_t kTRUE
Definition: RtypesCore.h:87
const char Option_t
Definition: RtypesCore.h:62
#define ClassImp(name)
Definition: Rtypes.h:365
char name[80]
Definition: TGX11.cxx:109
double pow(double, double)
char * Form(const char *fmt,...)
#define gPad
Definition: TVirtualPad.h:286
@ kMove
Definition: TVirtualX.h:46
@ kHand
Definition: TVirtualX.h:46
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition: TAttLine.h:42
virtual Style_t GetLineStyle() const
Return the line style.
Definition: TAttLine.h:34
virtual void Modify()
Change current line attributes if necessary.
Definition: TAttLine.cxx:242
Style_t fLineStyle
Line style.
Definition: TAttLine.h:22
Int_t DistancetoLine(Int_t px, Int_t py, Double_t xp1, Double_t yp1, Double_t xp2, Double_t yp2)
Compute distance from point px,py to a line.
Definition: TAttLine.cxx:206
virtual void SetTextAlign(Short_t align=11)
Set the text alignment.
Definition: TAttText.h:41
virtual void SetTextColor(Color_t tcolor=1)
Set the text color.
Definition: TAttText.h:43
virtual void SetTextFont(Font_t tfont=62)
Set the text font.
Definition: TAttText.h:45
The axis painter class.
Definition: TGaxis.h:24
virtual void PaintAxis(Double_t xmin, Double_t ymin, Double_t xmax, Double_t ymax, Double_t &wmin, Double_t &wmax, Int_t &ndiv, Option_t *chopt="", Double_t gridlength=0, Bool_t drawGridOnly=kFALSE)
Control function to draw an axis.
Definition: TGaxis.cxx:956
void SetLabelFont(Int_t labelfont)
Definition: TGaxis.h:106
void SetLabelOffset(Float_t labeloffset)
Definition: TGaxis.h:107
void SetLabelColor(Int_t labelcolor)
Definition: TGaxis.h:105
void SetLabelSize(Float_t labelsize)
Definition: TGaxis.h:108
void LabelsLimits(const char *label, Int_t &first, Int_t &last)
Internal method to find first and last character of a label.
Definition: TGaxis.cxx:2411
To draw polar axis.
Double_t fRwrmin
Minimal radial value (real world)
void ReduceFraction(Int_t Num, Int_t Denom, Int_t &rnum, Int_t &rden)
Reduce fractions.
Double_t GetPolarLabelSize()
Color_t fRadialLabelColor
Set color of the radial labels.
void SetPolarLabelColor(Color_t tcolorangular=1)
Set Polar labels color.
void PaintPolarDivisions(Bool_t noLabels)
Draw Polar divisions.
void SetRangePolar(Double_t tmin, Double_t tmax)
Allows to change range Polar.
Double_t fPolarTextSize
Set Polar text size.
Double_t fRwtmin
Minimal angular value (real world)
Double_t GetRadialOffset()
void Init()
Initialize some of the fields of TGraphPolargram.
Double_t fRadialTextSize
Double_t GetTickpolarSize()
void SetPolarOffset(Double_t PolarOffset=0.04)
Set the labels offset.
void SetTwoPi()
Set range from 0 to 2*pi.
Double_t FindTextAngle(Double_t theta)
Determine the orientation of the polar labels according to their angle.
void SetRadialLabelColor(Color_t tcolorradial=1)
Set radial labels color.
void PaintCircle(Double_t x, Double_t y, Double_t r, Double_t phimin, Double_t phimax, Double_t theta)
This is simplified from TEllipse::PaintEllipse.
virtual ~TGraphPolargram()
TGraphPolargram destructor.
void SetRangeRadial(Double_t rmin, Double_t rmax)
Set the radial range.
void SetTickpolarSize(Double_t tickpolarsize=0.02)
Set polar ticks size.
void SetAxisAngle(Double_t angle=0)
Set axis angle.
void SetToDegree()
The Polar circle is labelled using degree.
void SetNdivPolar(Int_t Ndiv=508)
Set the number of Polar divisions: enter a number ij with 0<i<99 and 0<j<99.
Int_t fNdivRad
Number of radial divisions.
Font_t GetPolarLabelFont()
void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Indicate that there is something to click here.
Double_t fTickpolarSize
Set size of Tickmarks.
void SetRadialLabelSize(Double_t radialsize=0.035)
Set radial labels size.
void Paint(Option_t *options="")
Paint TGraphPolargram.
Font_t fPolarLabelFont
Set font of angular labels.
TString * fPolarLabels
! [fNdivPol] Specified polar labels
void PaintRadialDivisions(Bool_t drawaxis)
Paint radial divisions.
void SetToGrad()
The Polar circle is labelled using gradian.
Double_t GetAngle()
void ChangeRangePolar(Double_t tmin, Double_t tmax)
Set the Polar range.
void SetPolarLabelSize(Double_t angularsize=0.04)
Set angular labels size.
void SetRadialLabelFont(Font_t tfontradial=62)
Set radial label font.
Double_t fAxisAngle
Set angle of the radial axis.
void SetToRadian()
The Polar circle is labelled using radian.
Double_t fRwtmax
Minimal angular value (real world)
Double_t fRwrmax
Maximal radial value (real world)
void SetPolarLabel(Int_t div, const TString &label)
Set some specified polar labels, used in the case of a spider plot.
TGraphPolargram(const char *name, Double_t rmin, Double_t rmax, Double_t tmin, Double_t tmax)
TGraphPolargram Constructor.
Font_t fRadialLabelFont
Set font of radial labels.
Int_t fNdivPol
Number of polar divisions.
Font_t GetRadialLabelFont()
Int_t fCutRadial
if fCutRadial = 0, circles are cut by radial axis if fCutRadial = 1, circles are not cut
Color_t fPolarLabelColor
Set color of the angular labels.
Double_t fPolarOffset
Offset for Polar labels.
Int_t DistancetoPrimitive(Int_t px, Int_t py)
Everything within the circle belongs to the TGraphPolargram.
void Draw(Option_t *options="")
Draw Polargram.
void SetPolarLabelFont(Font_t tfontangular=62)
void SetRadialOffset(Double_t RadialOffset=0.025)
Set the labels offset.
Color_t GetRadialColorLabel()
Int_t FindAlign(Double_t angle)
Find the alignement rule to apply for TText::SetTextAlign(Short_t).
void SetNdivRadial(Int_t Ndiv=508)
Set the number of radial divisions: enter a number ij with 0<i<99 and 0<j<99.
Double_t GetRadialLabelSize()
Color_t GetPolarColorLabel()
Double_t fRadialOffset
Offset for radial labels.
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.
To draw Mathematical Formula.
Definition: TLatex.h:18
virtual void PaintLatex(Double_t x, Double_t y, Double_t angle, Double_t size, const char *text)
Main drawing function.
Definition: TLatex.cxx:2038
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:105
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:694
void ResetBit(UInt_t f)
Definition: TObject.h:171
Basic string class.
Definition: TString.h:131
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1138
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:619
virtual void PaintText(Double_t x, Double_t y, const char *text)
Draw this text with new coordinates.
Definition: TText.cxx:743
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
static constexpr double rad
static constexpr double s
static constexpr double pi
Double_t ACos(Double_t)
Definition: TMath.h:658
Double_t ASin(Double_t)
Definition: TMath.h:651
constexpr Double_t PiOver2()
Definition: TMath.h:52
Double_t Log(Double_t x)
Definition: TMath.h:750
Double_t Sqrt(Double_t x)
Definition: TMath.h:681
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:725
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
Double_t Tan(Double_t)
Definition: TMath.h:635
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
Definition: first.py:1
Definition: test.py:1
auto * a
Definition: textangle.C:12