Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TCandle.cxx
Go to the documentation of this file.
1// @(#)root/graf:$Id$
2// Author: Georg Troska 19/05/16
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 <cstdlib>
13
14#include "TBuffer.h"
15#include "TCandle.h"
16#include "TVirtualPad.h"
17#include "TH2D.h"
18#include "TRandom2.h"
19#include "TStyle.h"
20#include "strlcpy.h"
21
22
23/** \class TCandle
24\ingroup BasicGraphics
25
26The candle plot painter class.
27
28Instances of this class are generated by the histograms painting
29classes (THistPainter and THStack) when an candle plot (box plot) is drawn.
30TCandle is the "painter class" of the box plots. Therefore it is never used
31directly to draw a candle.
32*/
33
34////////////////////////////////////////////////////////////////////////////////
35/// TCandle default constructor.
36
38{
39 fIsCalculated = false;
40 fIsRaw = false;
41 fPosCandleAxis = 0.;
42 fCandleWidth = 1.0;
43 fHistoWidth = 1.0;
44 fMean = 0.;
45 fMedian = 0.;
46 fMedianErr = 0;
47 fBoxUp = 0.;
48 fBoxDown = 0.;
49 fWhiskerUp = 0.;
50 fWhiskerDown = 0.;
51 fNDatapoints = 0;
52 fDismiss = false;
53 fLogX = 0;
54 fLogY = 0;
55 fLogZ = 0;
56 fNDrawPoints = 0;
57 fNHistoPoints = 0;
58 fAxisMin = 0.;
59 fAxisMax = 0.;
61 fProj = nullptr;
62 fDatapoints = nullptr;
63
64}
65
66////////////////////////////////////////////////////////////////////////////////
67/// TCandle constructor passing a draw-option.
68
69TCandle::TCandle(const char *opt)
70{
71 fIsCalculated = false;
72 fIsRaw = false;
73 fPosCandleAxis = 0.;
74 fCandleWidth = 1.0;
75 fHistoWidth = 1.0;
76 fMean = 0.;
77 fMedian = 0.;
78 fMedianErr = 0;
79 fBoxUp = 0.;
80 fBoxDown = 0.;
81 fWhiskerUp = 0.;
82 fWhiskerDown = 0.;
83 fNDatapoints = 0;
84 fDismiss = false;
85 fLogX = 0;
86 fLogY = 0;
87 fLogZ = 0;
88 fNDrawPoints = 0;
89 fNHistoPoints = 0;
90 fAxisMin = 0.;
91 fAxisMax = 0.;
93 fProj = nullptr;
94 fDatapoints = nullptr;
95
96 // Conversion necessary in order to cast from const char* to char*
97 char myopt[128];
98 strlcpy(myopt,opt,128);
99
100
102}
103
104
105////////////////////////////////////////////////////////////////////////////////
106/// TCandle constructor for raw-data candles.
107
110{
111 //Preliminary values only, need to be calculated before paint
112 fMean = 0;
113 fMedian = 0;
114 fMedianErr = 0;
115 fBoxUp = 0;
116 fBoxDown = 0;
117 fWhiskerUp = 0;
118 fWhiskerDown = 0;
119 fNDatapoints = n;
120 fIsCalculated = false;
121 fIsRaw = true;
126 fProj = nullptr;
127 fDismiss = false;
129 fLogX = 0;
130 fLogY = 0;
131 fLogZ = 0;
132 fNDrawPoints = 0;
133 fNHistoPoints = 0;
134 fAxisMin = 0.;
135 fAxisMax = 0.;
136}
137
138////////////////////////////////////////////////////////////////////////////////
139/// TCandle TH1 data constructor.
140
143{
144 //Preliminary values only, need to be calculated before paint
145 fMean = 0;
146 fMedian = 0;
147 fMedianErr = 0;
148 fBoxUp = 0;
149 fBoxDown = 0;
150 fWhiskerUp = 0;
151 fWhiskerDown = 0;
152 fNDatapoints = 0;
153 fIsCalculated = false;
154 fIsRaw = false;
158 fDatapoints = nullptr;
159 fProj = proj;
160 fDismiss = false;
162 fLogX = 0;
163 fLogY = 0;
164 fLogZ = 0;
165 fNDrawPoints = 0;
166 fNHistoPoints = 0;
167 fAxisMin = 0.;
168 fAxisMax = 0.;
169}
170
171////////////////////////////////////////////////////////////////////////////////
172/// TCandle default destructor.
173
175{
176 if (fIsRaw && fProj) delete fProj;
177}
178
179////////////////////////////////////////////////////////////////////////////////
180/// Returns true if candle plot should be scaled
181
183{
184 return gStyle->GetCandleScaled();
185}
186
187////////////////////////////////////////////////////////////////////////////////
188/// Returns true if violin plot should be scaled
189
191{
192 return gStyle->GetViolinScaled();
193}
194
195////////////////////////////////////////////////////////////////////////////////
196/// Static function to set WhiskerRange, by setting whisker-range, one can force
197/// the whiskers to cover the fraction of the distribution.
198/// Set wRange between 0 and 1. Default is 1
199/// TCandle::SetWhiskerRange(0.95) will set all candle-charts to cover 95% of
200/// the distribution with the whiskers.
201/// Can only be used with the standard-whisker definition
202
207
208////////////////////////////////////////////////////////////////////////////////
209/// Static function to set BoxRange, by setting box-range, one can force the
210/// box of the candle-chart to cover that given fraction of the distribution.
211/// Set bRange between 0 and 1. Default is 0.5
212/// TCandle::SetBoxRange(0.68) will set all candle-charts to cover 68% of the
213/// distribution by the box
214
219
220////////////////////////////////////////////////////////////////////////////////
221/// Static function to set scaling between candles-withs. A candle containing
222/// 100 entries will be two times wider than a candle containing 50 entries
223
228
229////////////////////////////////////////////////////////////////////////////////
230/// Static function to set scaling between violin-withs. A violin or histo chart
231/// with a maximum bin content to 100 will be two times as high as a violin with
232/// a maximum bin content of 50
233
238
239////////////////////////////////////////////////////////////////////////////////
240/// Parsing of the option-string.
241/// The option-string will be empty at the end (by-reference).
242
243int TCandle::ParseOption(char * opt) {
245 char *l = strstr(opt,"CANDLE");
246
247 if (l) {
249
250 char direction = ' ';
251 char preset = ' ';
252
253 if (l[6] >= 'A' && l[6] <= 'Z') direction = l[6];
254 if (l[6] >= '1' && l[6] <= '9') preset = l[6];
255 if (l[7] >= 'A' && l[7] <= 'Z' && preset != ' ') direction = l[7];
256 if (l[7] >= '1' && l[7] <= '9' && direction != ' ') preset = l[7];
257
258 if (direction == 'X' || direction == 'V') { /* nothing */ }
259 if (direction == 'Y' || direction == 'H') { fOption = (CandleOption)(fOption + kHorizontal); }
260 if (preset == '1') //Standard candle using old candle-definition
262 else if (preset == '2') //New standard candle with better whisker definition + outlier
264 else if (preset == '3') //Like candle2 but with a fMean as a circle
266 else if (preset == '4') //Like candle3 but showing the uncertainty of the fMedian as well
268 else if (preset == '5') //Like candle2 but showing all datapoints
270 else if (preset == '6') //Like candle2 but showing all datapoints scattered
272 else if (preset != ' ') //For all other presets not implemented yet used fallback candle
274
275 if (preset != ' ' && direction != ' ')
276 memcpy(l," ",8);
277 else if (preset != ' ' || direction != ' ')
278 memcpy(l," ",7);
279 else
280 memcpy(l," ",6);
281
282 Bool_t useIndivOption = false;
283
284 if (direction == ' ') direction = 'X';
285 if (preset == ' ') { // Check if the user wants to set the properties individually
286 char *brOpen = strstr(opt,"(");
287 char *brClose = strstr(opt,")");
288 char indivOption[32];
289 if (brOpen && brClose) {
290 useIndivOption = true;
292 strlcpy(indivOption, brOpen, brClose-brOpen+2); //Now the string "(....)" including brackets is in this array
293 sscanf(indivOption,"(%d)", (int*) &fOption);
295 memcpy(brOpen," ",brClose-brOpen+1); //Cleanup
296
297 fOptionStr.Form("CANDLE%c(%ld)", direction, (long)fOption);
298 } else {
300 }
301 } else {
302 fOptionStr.Form("CANDLE%c%c",direction,preset);
303 }
304 //Handle option "CANDLE" ,"CANDLEX" or "CANDLEY" to behave like "CANDLEX1" or "CANDLEY1"
305 if (!useIndivOption && !fOption ) {
307 fOptionStr.Form("CANDLE%c2",direction);
308 }
309 }
310
311 l = strstr(opt,"VIOLIN");
312 if (l) {
314
315 char direction = ' ';
316 char preset = ' ';
317
318 if (l[6] >= 'A' && l[6] <= 'Z') direction = l[6];
319 if (l[6] >= '1' && l[6] <= '9') preset = l[6];
320 if (l[7] >= 'A' && l[7] <= 'Z' && preset != ' ') direction = l[7];
321 if (l[7] >= '1' && l[7] <= '9' && direction != ' ') preset = l[7];
322
323 if (direction == 'X' || direction == 'V') { /* nothing */ }
324 if (direction == 'Y' || direction == 'H') { fOption = (CandleOption)(fOption + kHorizontal); }
325 if (preset == '1') //Standard candle using old candle-definition
327 else if (preset == '2') //New standard candle with better whisker definition + outlier
329 else if (preset != ' ') //For all other presets not implemented yet used fallback candle
331
332 if (preset != ' ' && direction != ' ')
333 memcpy(l," ",8);
334 else if (preset != ' ' || direction != ' ')
335 memcpy(l," ",7);
336 else
337 memcpy(l," ",6);
338
339 Bool_t useIndivOption = false;
340
341 if (direction == ' ') direction = 'X';
342 if (preset == ' ') { // Check if the user wants to set the properties individually
343 char *brOpen = strstr(opt,"(");
344 char *brClose = strstr(opt,")");
345 char indivOption[32];
346 if (brOpen && brClose) {
347 useIndivOption = true;
349 strlcpy(indivOption, brOpen, brClose-brOpen +2); //Now the string "(....)" including brackets is in this array
350 sscanf(indivOption,"(%d)", (int*) &fOption);
352 memcpy(brOpen," ",brClose-brOpen+1); //Cleanup
353
354 fOptionStr.Form("VIOLIN%c(%ld)", direction, (long)fOption);
355 } else {
357 }
358 } else {
359 fOptionStr.Form("VIOLIN%c%c", direction, preset);
360 }
361 //Handle option "VIOLIN" ,"VIOLINX" or "VIOLINY" to behave like "VIOLINX1" or "VIOLINY1"
362 if (!useIndivOption && !fOption ) {
364 fOptionStr.Form("VIOLIN%c1", direction);
365 }
366 }
367
368 fIsCalculated = false;
369
370 return fOption;
371
372}
373
374////////////////////////////////////////////////////////////////////////////////
375/// Calculates all values needed by the candle definition depending on the
376/// candle options.
377
379 //Reset everything
380 fNDrawPoints = 0;
381 fNHistoPoints = 0;
382
384 Bool_t doLogY = (!(swapXY) && fLogY) || (swapXY && fLogX);
385 Bool_t doLogX = (!(swapXY) && fLogX) || (swapXY && fLogY);
387
388 //Will be min and max values of raw-data
389 Double_t min = 1e15;
390 Double_t max = -1e15;
391
392 // Determining the quantiles
393 Double_t prob[5];
394
397
398
399 if (wRange >= 1) {
400 prob[0] = 1e-15;
401 prob[4] = 1-1e-15;
402 } else {
403 prob[0] = 0.5 - wRange/2.;
404 prob[4] = 0.5 + wRange/2.;
405 }
406
407 if (bRange >= 1) {
408 prob[1] = 1E-14;
409 prob[3] = 1-1E-14;
410 } else {
411 prob[1] = 0.5 - bRange/2.;
412 prob[3] = 0.5 + bRange/2.;
413 }
414
415 prob[2] = 0.5;
416
417 Double_t quantiles[5] = { 0., 0., 0., 0., 0. };
418 if (!fIsRaw && fProj) { //Need a calculation for a projected histo
419 if (((IsOption(kHistoLeft)) || (IsOption(kHistoRight)) || (IsOption(kHistoViolin))) && fProj->GetNbinsX() > 500) {
420 // When using the histooption the number of bins of the projection is
421 // limited because of the array space defined by kNMAXPOINTS.
422 // So the histo is rebinned, that it can be displayed at any time.
423 // Finer granularity is not useful anyhow
424 int divideBy = ((fProj->GetNbinsX() - 1)/((kNMAXPOINTS-10)/4))+1;
426 }
428 } else { //Need a calculation for a raw-data candle
430 }
431
432 // Check if the quantiles are valid, seems the under- and overflow is taken
433 // into account as well, we need to ignore this!
434 if (quantiles[0] >= quantiles[4] ||
435 quantiles[1] >= quantiles[3]) {
436 return;
437 }
438
439 // Definition of the candle in the standard case
440 fBoxUp = quantiles[3];
441 fBoxDown = quantiles[1];
442 fWhiskerUp = quantiles[4]; //Standard case
443 fWhiskerDown = quantiles[0]; //Standard case
444 fMedian = quantiles[2];
446 Int_t nOutliers = 0;
447
448 if (IsOption(kWhisker15)) { // Improved whisker definition, with 1.5*iqr
449 if (!fIsRaw && fProj) { //Need a calculation for a projected histo
450 int bin = fProj->FindBin(fBoxDown-1.5*iqr);
451 // extending only to the lowest data value within this range
452 while (fProj->GetBinContent(bin) == 0 && bin <= fProj->GetNbinsX()) bin++;
454
455 bin = fProj->FindBin(fBoxUp+1.5*iqr);
456 while (fProj->GetBinContent(bin) == 0 && bin >= 1) bin--;
458 } else { //Need a calculation for a raw-data candle
461
462 //Need to find highest value up to 1.5*iqr from the BoxUp-pos, and the lowest value up to -1.5*iqr from the boxLow-pos
463 for (Long64_t i = 0; i < fNDatapoints; ++i) {
465 if (myData > fWhiskerUp && myData <= fBoxUp + 1.5*iqr) fWhiskerUp = myData;
467 }
468 }
469 }
470
471 if (!fIsRaw && fProj) { //Need a calculation for a projected histo
472 fMean = fProj->GetMean();
473 fMedianErr = 1.57*iqr/sqrt(fProj->GetEntries());
476 } else { //Need a calculation for a raw-data candle
477 //Calculate the Mean
478 fMean = 0;
479 for (Long64_t i = 0; i < fNDatapoints; ++i) {
480 fMean += fDatapoints[i];
481 if (fDatapoints[i] < min) min = fDatapoints[i];
482 if (fDatapoints[i] > max) max = fDatapoints[i];
484 }
486 fMedianErr = 1.57*iqr/sqrt(fNDatapoints);
487 }
488
489 //Doing the outliers and other single points to show
490 if (GetCandleOption(5) > 0) { //Draw outliers
492 const int maxOutliers = kNMAXPOINTS;
493 Double_t myScale = 1.;
494 if (!fIsRaw && fProj) { //Need a calculation for a projected histo
496 fNDrawPoints = 0;
497 for (int bin = 0; bin < fProj->GetNbinsX(); bin++) {
498 // Either show them only outside the whiskers, or all of them
499 if (fProj->GetBinContent(bin) > 0 && (fProj->GetBinCenter(bin) < fWhiskerDown || fProj->GetBinCenter(bin) > fWhiskerUp || (GetCandleOption(5) > 1)) ) {
501 if (scaledBinContent >0 && scaledBinContent < 1) scaledBinContent = 1; //Outliers have a typical bin content between 0 and 1, when scaling they would disappear
502 for (int j=0; j < (int)scaledBinContent; j++) {
503 if (fNDrawPoints > maxOutliers) break;
504 if (IsOption(kPointsAllScat)) { //Draw outliers and "all" values scattered
507 } else { //Draw them in the "candle line"
509 if ((int)scaledBinContent == 1) //If there is only one datapoint available put it in the middle of the bin
511 else //If there is more than one datapoint scatter it along the bin, otherwise all marker would be (invisibly) stacked on top of each other
513 }
514 if (swapXY) {
515 //Swap X and Y
520 }
521 // Continue fMeans, that fNDrawPoints is not increased, so that value will not be shown
522 if (doLogX) {
524 }
525 if (doLogY) {
527 }
528 fNDrawPoints++;
529 }
530 }
531 if (fNDrawPoints > maxOutliers) { //Should never happen, due to myScale!!!
532 Error ("PaintCandlePlot","Not possible to draw all outliers.");
533 break;
534 }
535 }
536 } else { //Raw data candle
537 //If only outliers are shown, calculate myScale only based on nOutliers, use fNDatapoints (all) instead
540 } else {
542 }
543 fNDrawPoints = 0;
544 for (int i = 0; i < fNDatapoints; i++ ) {
547 if (!(i % (int) myScale == 0 )) continue; //If the amount of data is too large take only every 2nd or 3rd to reduce the amount
548 // Either show them only outside the whiskers, or all of them
550 if (IsOption(kPointsAllScat)) { //Draw outliers and "all" values scattered
552 fDrawPointsY[fNDrawPoints] = myData + (random.Rndm() - 0.5)*maxScatter; //random +- 0.5 of candle-height
553 } else { //Draw them in the "candle line"
555 fDrawPointsY[fNDrawPoints] = myData + (random.Rndm() - 0.5)*maxScatter; //random +- 0.5 of candle-height
556 }
557 if (swapXY) {
558 //Swap X and Y
563 }
564 // Continue fMeans, that fNDrawPoints is not increased, so that value will not be shown
565 if (doLogX) {
567 else continue;
568 }
569 if (doLogY) {
571 else continue;
572 }
573 fNDrawPoints++;
574 if (fNDrawPoints > maxOutliers) { //Should never happen, due to myScale!!!
575 Error ("PaintCandlePlotRaw","Not possible to draw all outliers.");
576 break;
577 }
578 }
579 }
580 }
581 }
583 //We are starting with kHistoRight, left will be modified from right later
584 if (fIsRaw) { //This is a raw-data candle
585 if (!fProj) {
586 fProj = new TH1D("hpa","hpa",100,min,max+0.0001*(max-min));
587 fProj->SetDirectory(nullptr);
588 for (Long64_t i = 0; i < fNDatapoints; ++i) {
590 }
591 }
592 }
593
594 fNHistoPoints = 0;
598
599 bool isFirst = true;
600 int lastNonZero = 0;
601 for (int bin = 1; bin <= fProj->GetNbinsX(); bin++) {
602 if (isFirst) {
603 if (fProj->GetBinContent(bin) > 0) {
606 if (doLogX) {
608 }
609 if (doLogY) {
611 }
613 isFirst = false;
614 } else {
615 continue;
616 }
617 }
618
620 if (doLogZ) {
622 }
628 if (doLogX) {
631 }
632 if (doLogY) {
635 }
636
639 }
640
643 fNHistoPoints = lastNonZero+1; //+1 so that the line down to 0 is added as well
644
645 if (IsOption(kHistoLeft)) {
646 for (int i = 0; i < fNHistoPoints; i++) {
648 }
649 }
650 if (IsOption(kHistoViolin)) {
651 for (int i = 0; i < fNHistoPoints; i++) {
654 }
655 fNHistoPoints *= 2;
656 }
657 }
658
659 fIsCalculated = true;
660}
661
662////////////////////////////////////////////////////////////////////////////////
663/// Paint one candle with its current attributes.
664
666{
667 if (!gPad) return;
668
669 //If something was changed before, we need to recalculate some values
670 if (!fIsCalculated) Calculate();
671
672 // Save the attributes as they were set originally
678
681
685
687 Bool_t doLogY = (!(swapXY) && fLogY) || (swapXY && fLogX);
688 Bool_t doLogX = (!(swapXY) && fLogX) || (swapXY && fLogY);
689
690 // From now on this is real painting only, no calculations anymore
691
698 }
699
700
705 }
706 if (!swapXY) {
709 } else {
712 }
716 }
717 }
718
719 if (IsOption(kBox)) { // Draw a simple box
720 if (IsOption(kMedianNotched)) { // Check if we have to draw a box with notches
725 PaintBox(11, x, y, swapXY);
726 } else { // draw a simple box
729 PaintBox(5, x, y, swapXY);
730 }
731 }
732
733 if (IsOption(kAnchor)) { // Draw the anchor line
736 }
737
738 if (IsOption(kWhiskerAll) && !IsOption(kHistoZeroIndicator)) { // Whiskers are dashed
739 SetLineStyle(2);
745 } else if ((IsOption(kWhiskerAll) && IsOption(kHistoZeroIndicator)) || IsOption(kWhisker15) ) { // Whiskers without dashing, better whisker definition, or forced when using zero line
748 }
749
750 if (IsOption(kMedianLine)) { // Paint fMedian as a line
752 } else if (IsOption(kMedianNotched)) { // Paint fMedian as a line (using notches, fMedian line is shorter)
754 } else if (IsOption(kMedianCircle)) { // Paint fMedian circle
756 if (!swapXY) {
758 myMedianY[0] = fMedian;
759 } else {
760 myMedianX[0] = fMedian;
762 }
763
764 Bool_t isValid = true;
765 if (doLogX) {
766 if (myMedianX[0] > 0) myMedianX[0] = TMath::Log10(myMedianX[0]); else isValid = false;
767 }
768 if (doLogY) {
769 if (myMedianY[0] > 0) myMedianY[0] = TMath::Log10(myMedianY[0]); else isValid = false;
770 }
772 if (mw == 1) SetMarkerStyle(24);
773 else SetMarkerStyle(18*mw+17);
775
776 if (isValid) gPad->PaintPolyMarker(1,myMedianX,myMedianY); // A circle for the fMedian
777
780
781 }
782
783 if (IsOption(kMeanCircle)) { // Paint fMean as a circle
784 Double_t myMeanX[1], myMeanY[1];
785 if (!swapXY) {
787 myMeanY[0] = fMean;
788 } else {
789 myMeanX[0] = fMean;
791 }
792
793 Bool_t isValid = true;
794 if (doLogX) {
795 if (myMeanX[0] > 0) myMeanX[0] = TMath::Log10(myMeanX[0]); else isValid = false;
796 }
797 if (doLogY) {
798 if (myMeanY[0] > 0) myMeanY[0] = TMath::Log10(myMeanY[0]); else isValid = false;
799 }
800
802 if (mw == 1) SetMarkerStyle(24);
803 else SetMarkerStyle(18*mw+17);
805
806 if (isValid) gPad->PaintPolyMarker(1,myMeanX,myMeanY); // A circle for the fMean
807
810
811 } else if (IsOption(kMeanLine)) { // Paint fMean as a dashed line
812 SetLineStyle(2);
814
818
819 }
820
821 if (IsOption(kAnchor)) { //Draw standard anchor
822 PaintLine(dimLeft, fWhiskerDown, dimRight, fWhiskerDown, swapXY); // the lower anchor line
823 PaintLine(dimLeft, fWhiskerUp, dimRight, fWhiskerUp, swapXY); // the upper anchor line
824 }
825
826 // This is a bit complex. All values here are handled as outliers. Usually
827 // only the datapoints outside the whiskers are shown.
828 // One can show them in one row as crosses, or scattered randomly. If activated
829 // all datapoint are shown in the same way
830
831 if (GetCandleOption(5) > 0) { //Draw outliers
832 if (IsOption(kPointsAllScat)) { //Draw outliers and "all" values scattered
834 } else {
836 if (mw == 1) SetMarkerStyle(5);
837 else SetMarkerStyle(18*mw+16);
838 }
840 gPad->PaintPolyMarker(fNDrawPoints,fDrawPointsX, fDrawPointsY);
841 }
842}
843
844////////////////////////////////////////////////////////////////////////////////
845/// Return true is this option is activated in fOption
846
848 long myOpt = 9;
849 int pos = 0;
850 for (pos = 0; pos < 16; pos++) {
851 if (myOpt > opt) break;
852 else myOpt *=10;
853 }
854 myOpt /= 9;
855 int thisOpt = GetCandleOption(pos);
856
857 return (thisOpt * myOpt) == opt;
858}
859
860////////////////////////////////////////////////////////////////////////////////
861/// Paint a box for candle.
862
864{
865 if (!gPad) return;
866
867 Bool_t doLogY = (!(swapXY) && fLogY) || (swapXY && fLogX);
868 Bool_t doLogX = (!(swapXY) && fLogX) || (swapXY && fLogY);
869 if (doLogY) {
870 for (int i=0; i<nPoints; i++) {
871 if (y[i] > 0) y[i] = TMath::Log10(y[i]);
872 else return;
873 }
874 }
875 if (doLogX) {
876 for (int i=0; i<nPoints; i++) {
877 if (x[i] > 0) x[i] = TMath::Log10(x[i]);
878 else return;
879 }
880 }
881 if (!swapXY) {
882 gPad->PaintFillArea(nPoints, x, y);
883 gPad->PaintPolyLine(nPoints, x, y);
884 } else {
885 gPad->PaintFillArea(nPoints, y, x);
886 gPad->PaintPolyLine(nPoints, y, x);
887 }
888}
889
890////////////////////////////////////////////////////////////////////////////////
891/// Paint a line for candle.
892
894{
895 if (!gPad) return;
896
897 Bool_t doLogY = (!(swapXY) && fLogY) || (swapXY && fLogX);
898 Bool_t doLogX = (!(swapXY) && fLogX) || (swapXY && fLogY);
899 if (doLogY) {
900 if (y1 > 0) y1 = TMath::Log10(y1); else return;
901 if (y2 > 0) y2 = TMath::Log10(y2); else return;
902 }
903 if (doLogX) {
904 if (x1 > 0) x1 = TMath::Log10(x1); else return;
905 if (x2 > 0) x2 = TMath::Log10(x2); else return;
906 }
907 if (!swapXY) {
908 gPad->PaintLine(x1, y1, x2, y2);
909 } else {
910 gPad->PaintLine(y1, x1, y2, x2);
911 }
912}
913
914////////////////////////////////////////////////////////////////////////////////
915/// Stream an object of class TCandle.
916
918{
919 if (R__b.IsReading()) {
921 Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
922 if (R__v > 3) {
923 R__b.ReadClassBuffer(TCandle::Class(), this, R__v, R__s, R__c);
924 return;
925 }
926 } else {
927 R__b.WriteClassBuffer(TCandle::Class(),this);
928 }
929}
930
931////////////////////////////////////////////////////////////////////////////////
932/// The coordinates in the TParallelCoordVar-class are in Pad-Coordinates, so we need to convert them
933
#define b(i)
Definition RSha256.hxx:100
#define a(i)
Definition RSha256.hxx:99
#define e(i)
Definition RSha256.hxx:103
short Style_t
Style number (short)
Definition RtypesCore.h:96
short Version_t
Class version identifier (short)
Definition RtypesCore.h:79
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
long long Long64_t
Portable signed long integer 8 bytes.
Definition RtypesCore.h:83
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
const Int_t kNMAXPOINTS
Definition TCandle.h:25
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Definition TError.cxx:208
Option_t Option_t SetLineColor
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
Option_t Option_t TPoint TPoint const char y2
Option_t Option_t SetMarkerStyle
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t points
Option_t Option_t TPoint TPoint const char y1
R__EXTERN TStyle * gStyle
Definition TStyle.h:442
#define gPad
Fill Area Attributes class.
Definition TAttFill.h:20
virtual Color_t GetFillColor() const
Return the fill area color.
Definition TAttFill.h:31
virtual Style_t GetFillStyle() const
Return the fill area style.
Definition TAttFill.h:32
virtual void Modify()
Change current fill area attributes if necessary.
Definition TAttFill.cxx:215
Line Attributes class.
Definition TAttLine.h:20
virtual Color_t GetLineColor() const
Return the line color.
Definition TAttLine.h:35
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition TAttLine.h:44
virtual Style_t GetLineStyle() const
Return the line style.
Definition TAttLine.h:36
virtual void Modify()
Change current line attributes if necessary.
Definition TAttLine.cxx:246
Marker Attributes class.
Definition TAttMarker.h:20
virtual void Modify()
Change current marker attributes if necessary.
virtual Style_t GetMarkerStyle() const
Return the marker style.
Definition TAttMarker.h:33
Double_t GetXmax() const
Definition TAxis.h:142
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition TAxis.cxx:521
Double_t GetXmin() const
Definition TAxis.h:141
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
Definition TAxis.cxx:531
Buffer base class used for serializing objects.
Definition TBuffer.h:43
Double_t fMedian
Position of the median.
Definition TCandle.h:64
CandleOption
Definition TCandle.h:30
@ kHorizontal
If this bit is not set it is vertical.
Definition TCandle.h:48
@ kAnchor
Definition TCandle.h:40
@ kWhiskerAll
Definition TCandle.h:38
@ kMeanCircle
Definition TCandle.h:37
@ kHistoZeroIndicator
Definition TCandle.h:47
@ kPointsAllScat
Definition TCandle.h:43
@ kNoOption
Definition TCandle.h:31
@ kHistoRight
Definition TCandle.h:45
@ kMeanLine
Definition TCandle.h:36
@ kBox
Definition TCandle.h:32
@ kWhisker15
Definition TCandle.h:39
@ kMedianLine
Definition TCandle.h:33
@ kPointsAll
Definition TCandle.h:42
@ kPointsOutliers
Definition TCandle.h:41
@ kMedianNotched
Definition TCandle.h:34
@ kHistoLeft
Definition TCandle.h:44
@ kMedianCircle
Definition TCandle.h:35
@ kHistoViolin
Definition TCandle.h:46
TH1D * fProj
Definition TCandle.h:56
static TClass * Class()
void Streamer(TBuffer &) override
Stream an object of class TCandle.
Definition TCandle.cxx:917
static void SetScaledViolin(const Bool_t vScale=true)
Static function to set scaling between violin-withs.
Definition TCandle.cxx:234
bool fIsRaw
0: for TH1 projection, 1: using raw data
Definition TCandle.h:54
CandleOption fOption
Setting the style of the candle.
Definition TCandle.h:82
TString fOptionStr
String to draw the candle.
Definition TCandle.h:83
Double_t fAxisMax
The Maximum which is visible by the axis (used by zero indicator)
Definition TCandle.h:89
~TCandle() override
TCandle default destructor.
Definition TCandle.cxx:174
int fLogZ
make the candle appear logz-like
Definition TCandle.h:86
Double_t fPosCandleAxis
x-pos for a vertical candle
Definition TCandle.h:59
int GetCandleOption(const int pos) const
Definition TCandle.h:93
int fNHistoPoints
Definition TCandle.h:80
Double_t fMean
Position of the mean.
Definition TCandle.h:63
Double_t fWhiskerDown
Position of the lower whisker end.
Definition TCandle.h:69
Long64_t fNDrawPoints
max number of outliers or other point to be shown
Definition TCandle.h:76
virtual void Paint(Option_t *option="")
Paint one candle with its current attributes.
Definition TCandle.cxx:665
Double_t fAxisMin
The Minimum which is visible by the axis (used by zero indicator)
Definition TCandle.h:88
void ConvertToPadCoords(Double_t minAxis, Double_t maxAxis, Double_t axisMinCoord, Double_t axisMaxCoord)
The coordinates in the TParallelCoordVar-class are in Pad-Coordinates, so we need to convert them.
Definition TCandle.cxx:934
Double_t fDrawPointsY[kNMAXPOINTS]
y-coord for every outlier, ..
Definition TCandle.h:75
Long64_t fNDatapoints
Number of Datapoints within this candle.
Definition TCandle.h:72
Bool_t IsViolinScaled() const
Returns true if violin plot should be scaled.
Definition TCandle.cxx:190
Double_t fHistoPointsY[kNMAXPOINTS]
y-coord for the polyline of the histo
Definition TCandle.h:79
Bool_t IsHorizontal() const
Definition TCandle.h:112
Double_t fHistoPointsX[kNMAXPOINTS]
x-coord for the polyline of the histo
Definition TCandle.h:78
static void SetScaledCandle(const Bool_t cScale=true)
Static function to set scaling between candles-withs.
Definition TCandle.cxx:224
Double_t * fDatapoints
position of all Datapoints within this candle
Definition TCandle.h:71
int ParseOption(char *optin)
Parsing of the option-string.
Definition TCandle.cxx:243
bool fDismiss
True if the candle cannot be painted.
Definition TCandle.h:57
Double_t fWhiskerUp
Position of the upper whisker end.
Definition TCandle.h:68
static void SetBoxRange(const Double_t bRange)
Static function to set BoxRange, by setting box-range, one can force the box of the candle-chart to c...
Definition TCandle.cxx:215
bool IsOption(CandleOption opt) const
Return true is this option is activated in fOption.
Definition TCandle.cxx:847
int fLogY
make the candle appear logy-like
Definition TCandle.h:85
Double_t fMedianErr
The size of the notch.
Definition TCandle.h:65
void PaintBox(Int_t nPoints, Double_t *x, Double_t *y, Bool_t swapXY)
Paint a box for candle.
Definition TCandle.cxx:863
Double_t fBoxUp
Position of the upper box end.
Definition TCandle.h:66
int fLogX
make the candle appear logx-like
Definition TCandle.h:84
Double_t fDrawPointsX[kNMAXPOINTS]
x-coord for every outlier, ..
Definition TCandle.h:74
TCandle()
TCandle default constructor.
Definition TCandle.cxx:37
Double_t fBoxDown
Position of the lower box end.
Definition TCandle.h:67
void Calculate()
Calculates all values needed by the candle definition depending on the candle options.
Definition TCandle.cxx:378
void PaintLine(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Bool_t swapXY)
Paint a line for candle.
Definition TCandle.cxx:893
Double_t fHistoWidth
The histo width (the height of the max bin)
Definition TCandle.h:61
Bool_t IsCandleScaled() const
Returns true if candle plot should be scaled.
Definition TCandle.cxx:182
bool fIsCalculated
Definition TCandle.h:55
static void SetWhiskerRange(const Double_t wRange)
Static function to set WhiskerRange, by setting whisker-range, one can force the whiskers to cover th...
Definition TCandle.cxx:203
Double_t fCandleWidth
The candle width.
Definition TCandle.h:60
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:927
virtual void SetDirectory(TDirectory *dir)
By default, when a histogram is created, it is added to the list of histogram objects in the current ...
Definition TH1.cxx:8978
virtual Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
Definition TH1.cxx:9179
virtual Double_t GetMean(Int_t axis=1) const
For axis = 1,2 or 3 returns the mean value of the histogram along X,Y or Z axis.
Definition TH1.cxx:7571
TAxis * GetXaxis()
Definition TH1.h:572
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:8586
virtual Int_t GetNbinsX() const
Definition TH1.h:542
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3326
virtual TH1 * RebinX(Int_t ngroup=2, const char *newname="")
Definition TH1.h:603
virtual Double_t GetBinLowEdge(Int_t bin) const
Return bin lower edge for 1D histogram.
Definition TH1.cxx:9190
virtual Double_t GetEntries() const
Return the current number of entries.
Definition TH1.cxx:4411
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition TH1.cxx:5076
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width for 1D histogram.
Definition TH1.cxx:9201
virtual Int_t FindBin(Double_t x, Double_t y=0, Double_t z=0)
Return Global bin number corresponding to x,y,z.
Definition TH1.cxx:3660
virtual Int_t GetQuantiles(Int_t n, Double_t *xp, const Double_t *p=nullptr)
Compute Quantiles for this histogram.
Definition TH1.cxx:4603
Random number generator class based on the maximally quidistributed combined Tausworthe generator by ...
Definition TRandom2.h:27
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition TString.cxx:2362
Bool_t GetViolinScaled() const
Definition TStyle.h:295
Bool_t GetCandleScaled() const
Definition TStyle.h:294
Double_t GetCandleWhiskerRange() const
Definition TStyle.h:292
void SetCandleScaled(Bool_t on=kFALSE)
Definition TStyle.h:429
Double_t GetCandleBoxRange() const
Definition TStyle.h:293
Int_t GetCandleCircleLineWidth() const
Definition TStyle.h:296
Int_t GetCandleCrossLineWidth() const
Definition TStyle.h:297
void SetViolinScaled(Bool_t on=kTRUE)
Definition TStyle.h:430
void SetCandleWhiskerRange(Double_t wRange=1.0)
By setting whisker-range for candle plot, one can force the whiskers to cover the fraction of the dis...
Definition TStyle.cxx:1940
void SetCandleBoxRange(Double_t bRange=0.5)
By setting box-range for candle plot, one can force the box of the candle-chart to cover that given f...
Definition TStyle.cxx:1957
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob, Bool_t isSorted=kTRUE, Int_t *index=nullptr, Int_t type=7)
Computes sample quantiles, corresponding to the given probabilities.
Definition TMath.cxx:1207
Double_t Log10(Double_t x)
Returns the common (base-10) logarithm of x.
Definition TMath.h:773
TLine l
Definition textangle.C:4