Logo ROOT  
Reference Guide
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 "strlcpy.h"
20 
25 
27 
28 /** \class TCandle
29 \ingroup BasicGraphics
30 
31 The candle plot painter class.
32 
33 Instances of this class are generated by the histograms painting
34 classes (THistPainter and THStack) when an candle plot (box plot) is drawn.
35 TCandle is the "painter class" of the box plots. Therefore it is never used
36 directly to draw a candle.
37 */
38 
39 ////////////////////////////////////////////////////////////////////////////////
40 /// TCandle default constructor.
41 
43 {
44  fIsCalculated = 0;
45  fIsRaw = 0;
46  fPosCandleAxis = 0.;
47  fCandleWidth = 1.0;
48  fHistoWidth = 1.0;
49  fMean = 0.;
50  fMedian = 0.;
51  fMedianErr = 0;
52  fBoxUp = 0.;
53  fBoxDown = 0.;
54  fWhiskerUp = 0.;
55  fWhiskerDown = 0.;
56  fNDatapoints = 0;
57  fDismiss = 0;
58  fLogX = 0;
59  fLogY = 0;
60  fLogZ = 0;
61  fNDrawPoints = 0;
62  fNHistoPoints = 0;
63  fAxisMin = 0.;
64  fAxisMax = 0.;
66  fProj = NULL;
67  fDatapoints = 0;
68 
69 }
70 
71 ////////////////////////////////////////////////////////////////////////////////
72 /// TCandle constructor passing a draw-option.
73 
74 TCandle::TCandle(const char *opt)
75 {
76  fIsCalculated = 0;
77  fIsRaw = 0;
78  fPosCandleAxis = 0.;
79  fCandleWidth = 1.0;
80  fHistoWidth = 1.0;
81  fMean = 0.;
82  fMedian = 0.;
83  fMedianErr = 0;
84  fBoxUp = 0.;
85  fBoxDown = 0.;
86  fWhiskerUp = 0.;
87  fWhiskerDown = 0.;
88  fNDatapoints = 0;
89  fDismiss = 0;
90  fLogX = 0;
91  fLogY = 0;
92  fLogZ = 0;
93  fNDrawPoints = 0;
94  fNHistoPoints = 0;
95  fAxisMin = 0.;
96  fAxisMax = 0.;
98  fProj = NULL;
99  fDatapoints = 0;
100 
101 
102  // Conversion necessary in order to cast from const char* to char*
103  char myopt[128];
104  strlcpy(myopt,opt,128);
105 
106 
107  ParseOption(myopt);
108 }
109 
110 
111 ////////////////////////////////////////////////////////////////////////////////
112 /// TCandle constructor for raw-data candles.
113 
114 TCandle::TCandle(const Double_t candlePos, const Double_t candleWidth, Long64_t n, Double_t * points)
115  : TAttLine(), TAttFill(), TAttMarker()
116 {
117  //Preliminary values only, need to be calculated before paint
118  fMean = 0;
119  fMedian = 0;
120  fMedianErr = 0;
121  fBoxUp = 0;
122  fBoxDown = 0;
123  fWhiskerUp = 0;
124  fWhiskerDown = 0;
125  fNDatapoints = n;
126  fIsCalculated = 0;
127  fIsRaw = true;
128  fPosCandleAxis = candlePos;
129  fCandleWidth = candleWidth;
130  fHistoWidth = candleWidth;
132  fProj = NULL;
133  fDismiss = 0;
134  fOption = kNoOption;
135  fLogX = 0;
136  fLogY = 0;
137  fLogZ = 0;
138  fNDrawPoints = 0;
139  fNHistoPoints = 0;
140  fAxisMin = 0.;
141  fAxisMax = 0.;
142  sprintf(fOptionStr," ");
143 }
144 
145 ////////////////////////////////////////////////////////////////////////////////
146 /// TCandle TH1 data constructor.
147 
148 TCandle::TCandle(const Double_t candlePos, const Double_t candleWidth, TH1D *proj)
149  : TAttLine(), TAttFill(), TAttMarker()
150 {
151  //Preliminary values only, need to be calculated before paint
152  fMean = 0;
153  fMedian = 0;
154  fMedianErr = 0;
155  fBoxUp = 0;
156  fBoxDown = 0;
157  fWhiskerUp = 0;
158  fWhiskerDown = 0;
159  fNDatapoints = 0;
160  fIsCalculated = 0;
161  fIsRaw = 0;
162  fPosCandleAxis = candlePos;
163  fCandleWidth = candleWidth;
164  fHistoWidth = candleWidth;
165  fDatapoints = 0;
166  fProj = proj;
167  fDismiss = 0;
168  fOption = kNoOption;
169  fLogX = 0;
170  fLogY = 0;
171  fLogZ = 0;
172  fNDrawPoints = 0;
173  fNHistoPoints = 0;
174  fAxisMin = 0.;
175  fAxisMax = 0.;
176  sprintf(fOptionStr," ");
177 }
178 
179 ////////////////////////////////////////////////////////////////////////////////
180 /// TCandle default destructor.
181 
183  if (fIsRaw && fProj) delete fProj;
184 }
185 
187 {
188  return fScaledCandle;
189 }
190 
192 {
193  return fScaledViolin;
194 }
195 
196 ////////////////////////////////////////////////////////////////////////////////
197 /// Static function to set fWhiskerRange, by setting whisker-range, one can force
198 /// the whiskers to cover the fraction of the distribution.
199 /// Set wRange between 0 and 1. Default is 1
200 /// TCandle::SetWhiskerRange(0.95) will set all candle-charts to cover 95% of
201 /// the distribution with the whiskers.
202 /// Can only be used with the standard-whisker definition
203 
204 void TCandle::SetWhiskerRange(const Double_t wRange) {
205  if (wRange < 0) fWhiskerRange = 0;
206  else if (wRange > 1) fWhiskerRange = 1;
207  else fWhiskerRange = wRange;
208 
209 }
210 
211 ////////////////////////////////////////////////////////////////////////////////
212 /// Static function to set fBoxRange, by setting whisker-range, one can force the
213 /// box of the candle-chart to cover that given fraction of the distribution.
214 /// Set bRange between 0 and 1. Default is 0.5
215 /// TCandle::SetBoxRange(0.68) will set all candle-charts to cover 68% of the
216 /// distribution by the box
217 
218 void TCandle::SetBoxRange(const Double_t bRange) {
219  if (bRange < 0) fBoxRange = 0;
220  else if (bRange > 1) fBoxRange = 1;
221  else fBoxRange = bRange;
222 }
223 
224 ////////////////////////////////////////////////////////////////////////////////
225 /// Static function to set scaling between candles-withs. A candle containing
226 /// 100 entries with be two times wider than a candle containing 50 entries
227 
228 void TCandle::SetScaledCandle(const Bool_t cScale) {
229  fScaledCandle = cScale;
230 }
231 
232 ////////////////////////////////////////////////////////////////////////////////
233 /// Static function to set scaling between violin-withs. A violin or histo chart
234 /// with a maximum bin content to 100 will be two times as high as a violin with
235 /// a maximum bin content of 50
236 
237 void TCandle::SetScaledViolin(const Bool_t vScale) {
238  fScaledViolin = vScale;
239 }
240 
241 ////////////////////////////////////////////////////////////////////////////////
242 /// Parsing of the option-string.
243 /// The option-string will be empty at the end (by-reference).
244 
245 int TCandle::ParseOption(char * opt) {
246  fOption = kNoOption;
247  char *l = strstr(opt,"CANDLE");
248 
249  if (l) {
250  const CandleOption fallbackCandle = (CandleOption)(kBox + kMedianLine + kMeanCircle + kWhiskerAll + kAnchor);
251 
252  char direction = ' ';
253  char preset = ' ';
254 
255  if (l[6] >= 'A' && l[6] <= 'Z') direction = l[6];
256  if (l[6] >= '1' && l[6] <= '9') preset = l[6];
257  if (l[7] >= 'A' && l[7] <= 'Z' && preset != ' ') direction = l[7];
258  if (l[7] >= '1' && l[7] <= '9' && direction != ' ') preset = l[7];
259 
260  if (direction == 'X' || direction == 'V') { /* nothing */ }
261  if (direction == 'Y' || direction == 'H') { fOption = (CandleOption)(fOption + kHorizontal); }
262  if (preset == '1') //Standard candle using old candle-definition
263  fOption = (CandleOption)(fOption + fallbackCandle);
264  else if (preset == '2') //New standard candle with better whisker definition + outlier
266  else if (preset == '3') //Like candle2 but with a fMean as a circle
268  else if (preset == '4') //Like candle3 but showing the uncertainty of the fMedian as well
270  else if (preset == '5') //Like candle2 but showing all datapoints
272  else if (preset == '6') //Like candle2 but showing all datapoints scattered
274  else if (preset != ' ') //For all other presets not implemented yet used fallback candle
275  fOption = (CandleOption)(fOption + fallbackCandle);
276 
277  if (preset != ' ' && direction != ' ')
278  memcpy(l," ",8);
279  else if (preset != ' ' || direction != ' ')
280  memcpy(l," ",7);
281  else
282  memcpy(l," ",6);
283 
284  Bool_t useIndivOption = false;
285 
286  if (direction == ' ') direction = 'X';
287  if (preset == ' ') { // Check if the user wants to set the properties individually
288  char *brOpen = strstr(opt,"(");
289  char *brClose = strstr(opt,")");
290  char indivOption[32];
291  if (brOpen && brClose) {
292  useIndivOption = true;
293  bool isHorizontal = IsHorizontal();
294  strlcpy(indivOption, brOpen, brClose-brOpen+2); //Now the string "(....)" including brackets is in this array
295  sscanf(indivOption,"(%d)", (int*) &fOption);
296  if (isHorizontal) {fOption = (CandleOption)(fOption + kHorizontal);}
297  memcpy(brOpen," ",brClose-brOpen+1); //Cleanup
298 
299  snprintf(fOptionStr, sizeof(fOptionStr), "CANDLE%c(%ld)",direction,(long)fOption);
300  } else {
301  fOption = (CandleOption)(fOption + fallbackCandle);
302  }
303  } else {
304  snprintf(fOptionStr, sizeof(fOptionStr), "CANDLE%c%c",direction,preset);
305  }
306  //Handle option "CANDLE" ,"CANDLEX" or "CANDLEY" to behave like "CANDLEX1" or "CANDLEY1"
307  if (!useIndivOption && !fOption ) {
308  fOption = fallbackCandle;
309  snprintf(fOptionStr, sizeof(fOptionStr), "CANDLE%c2",direction);
310  }
311  }
312 
313  l = strstr(opt,"VIOLIN");
314  if (l) {
316 
317  char direction = ' ';
318  char preset = ' ';
319 
320  if (l[6] >= 'A' && l[6] <= 'Z') direction = l[6];
321  if (l[6] >= '1' && l[6] <= '9') preset = l[6];
322  if (l[7] >= 'A' && l[7] <= 'Z' && preset != ' ') direction = l[7];
323  if (l[7] >= '1' && l[7] <= '9' && direction != ' ') preset = l[7];
324 
325  if (direction == 'X' || direction == 'V') { /* nothing */ }
326  if (direction == 'Y' || direction == 'H') { fOption = (CandleOption)(fOption + kHorizontal); }
327  if (preset == '1') //Standard candle using old candle-definition
328  fOption = (CandleOption)(fOption + fallbackCandle);
329  else if (preset == '2') //New standard candle with better whisker definition + outlier
331  else if (preset != ' ') //For all other presets not implemented yet used fallback candle
332  fOption = (CandleOption)(fOption + fallbackCandle);
333 
334  if (preset != ' ' && direction != ' ')
335  memcpy(l," ",8);
336  else if (preset != ' ' || direction != ' ')
337  memcpy(l," ",7);
338  else
339  memcpy(l," ",6);
340 
341  Bool_t useIndivOption = false;
342 
343  if (direction == ' ') direction = 'X';
344  if (preset == ' ') { // Check if the user wants to set the properties individually
345  char *brOpen = strstr(opt,"(");
346  char *brClose = strstr(opt,")");
347  char indivOption[32];
348  if (brOpen && brClose) {
349  useIndivOption = true;
350  bool isHorizontal = IsHorizontal();
351  strlcpy(indivOption, brOpen, brClose-brOpen +2); //Now the string "(....)" including brackets is in this array
352  sscanf(indivOption,"(%d)", (int*) &fOption);
353  if (isHorizontal) {fOption = (CandleOption)(fOption + kHorizontal);}
354  memcpy(brOpen," ",brClose-brOpen+1); //Cleanup
355 
356  snprintf(fOptionStr, sizeof(fOptionStr), "VIOLIN%c(%ld)",direction,(long)fOption);
357  } else {
358  fOption = (CandleOption)(fOption + fallbackCandle);
359  }
360  } else {
361  snprintf(fOptionStr, sizeof(fOptionStr), "VIOLIN%c%c",direction,preset);
362  }
363  //Handle option "VIOLIN" ,"VIOLINX" or "VIOLINY" to behave like "VIOLINX1" or "VIOLINY1"
364  if (!useIndivOption && !fOption ) {
365  fOption = fallbackCandle;
366  snprintf(fOptionStr, sizeof(fOptionStr), "VIOLIN%c1",direction);
367  }
368  }
369 
370  fIsCalculated = false;
371 
372  return fOption;
373 
374 }
375 
376 ////////////////////////////////////////////////////////////////////////////////
377 /// Calculates all values needed by the candle definition depending on the
378 /// candle options.
379 
381  //Reset everything
382  fNDrawPoints = 0;
383  fNHistoPoints = 0;
384 
385  Bool_t swapXY = IsOption(kHorizontal);
386  Bool_t doLogY = (!(swapXY) && fLogY) || (swapXY && fLogX);
387  Bool_t doLogX = (!(swapXY) && fLogX) || (swapXY && fLogY);
388  Bool_t doLogZ = fLogZ;
389 
390  //Will be min and max values of raw-data
391  Double_t min = 1e15;
392  Double_t max = -1e15;
393 
394  // Determining the quantiles
395  Double_t *prob = new Double_t[5];
396 
397  if (fWhiskerRange >= 1) {
398  prob[0] = 1e-15;
399  prob[4] = 1-1e-15;
400  } else {
401  prob[0] = 0.5 - fWhiskerRange/2.;
402  prob[4] = 0.5 + fWhiskerRange/2.;
403  }
404 
405 
406  if (fBoxRange >= 1) {
407  prob[1] = 1E-14;
408  prob[3] = 1-1E-14;
409  } else {
410  prob[1] = 0.5 - fBoxRange/2.;
411  prob[3] = 0.5 + fBoxRange/2.;
412  }
413 
414  prob[2]=0.5;
415  Double_t *quantiles = new Double_t[5];
416  quantiles[0]=0.; quantiles[1]=0.; quantiles[2] = 0.; quantiles[3] = 0.; quantiles[4] = 0.;
417  if (!fIsRaw && fProj) { //Need a calculation for a projected histo
418  if (((IsOption(kHistoLeft)) || (IsOption(kHistoRight)) || (IsOption(kHistoViolin))) && fProj->GetNbinsX() > 500) {
419  // When using the histooption the number of bins of the projection is
420  // limited because of the array space defined by kNMAXPOINTS.
421  // So the histo is rebinned, that it can be displayed at any time.
422  // Finer granularity is not useful anyhow
423  int divideBy = ((fProj->GetNbinsX() - 1)/((kNMAXPOINTS-10)/4))+1;
424  fProj->RebinX(divideBy);
425  }
426  fProj->GetQuantiles(5, quantiles, prob);
427  } else { //Need a calculation for a raw-data candle
429  }
430 
431  // Check if the quantiles are valid, seems the under- and overflow is taken
432  // into account as well, we need to ignore this!
433  if (quantiles[0] >= quantiles[4] ||
434  quantiles[1] >= quantiles[3]) {
435  delete [] prob;
436  delete [] quantiles;
437  return;
438  }
439 
440  // Definition of the candle in the standard case
441  fBoxUp = quantiles[3];
442  fBoxDown = quantiles[1];
443  fWhiskerUp = quantiles[4]; //Standard case
444  fWhiskerDown = quantiles[0]; //Standard case
445  fMedian = quantiles[2];
446  Double_t iqr = fBoxUp-fBoxDown;
447  Int_t nOutliers = 0;
448 
449  if (IsOption(kWhisker15)) { // Improved whisker definition, with 1.5*iqr
450  if (!fIsRaw && fProj) { //Need a calculation for a projected histo
451  int bin = fProj->FindBin(fBoxDown-1.5*iqr);
452  // extending only to the lowest data value within this range
453  while (fProj->GetBinContent(bin) == 0 && bin <= fProj->GetNbinsX()) bin++;
455 
456  bin = fProj->FindBin(fBoxUp+1.5*iqr);
457  while (fProj->GetBinContent(bin) == 0 && bin >= 1) bin--;
458  fWhiskerUp = fProj->GetBinCenter(bin);
459  } else { //Need a calculation for a raw-data candle
462 
463  //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
464  for (Long64_t i = 0; i < fNDatapoints; ++i) {
465  Double_t myData = fDatapoints[i];
466  if (myData > fWhiskerUp && myData <= fBoxUp + 1.5*iqr) fWhiskerUp = myData;
467  if (myData < fWhiskerDown && myData >= fBoxDown - 1.5*iqr) fWhiskerDown = myData;
468  }
469  }
470  }
471 
472  if (!fIsRaw && fProj) { //Need a calculation for a projected histo
473  fMean = fProj->GetMean();
474  fMedianErr = 1.57*iqr/sqrt(fProj->GetEntries());
475  fAxisMin = fProj->GetXaxis()->GetXmin();
476  fAxisMax = fProj->GetXaxis()->GetXmax();
477  } else { //Need a calculation for a raw-data candle
478  //Calculate the Mean
479  fMean = 0;
480  for (Long64_t i = 0; i < fNDatapoints; ++i) {
481  fMean += fDatapoints[i];
482  if (fDatapoints[i] < min) min = fDatapoints[i];
483  if (fDatapoints[i] > max) max = fDatapoints[i];
484  if (fDatapoints[i] < fWhiskerDown || fDatapoints[i] > fWhiskerUp) nOutliers++;
485  }
486  fMean /= fNDatapoints;
487  fMedianErr = 1.57*iqr/sqrt(fNDatapoints);
488  }
489 
490  delete [] prob;
491  delete [] quantiles;
492 
493  //Doing the outliers and other single points to show
494  if (GetCandleOption(5) > 0) { //Draw outliers
495  TRandom2 random;
496  const int maxOutliers = kNMAXPOINTS;
497  Double_t myScale = 1.;
498  if (!fIsRaw && fProj) { //Need a calculation for a projected histo
499  if (fProj->GetEntries() > maxOutliers/2) myScale = fProj->GetEntries()/(maxOutliers/2.);
500  fNDrawPoints = 0;
501  for (int bin = 0; bin < fProj->GetNbinsX(); bin++) {
502  // Either show them only outside the whiskers, or all of them
503  if (fProj->GetBinContent(bin) > 0 && (fProj->GetBinCenter(bin) < fWhiskerDown || fProj->GetBinCenter(bin) > fWhiskerUp || (GetCandleOption(5) > 1)) ) {
504  Double_t scaledBinContent = fProj->GetBinContent(bin)/myScale;
505  if (scaledBinContent >0 && scaledBinContent < 1) scaledBinContent = 1; //Outliers have a typical bin content between 0 and 1, when scaling they would disappear
506  for (int j=0; j < (int)scaledBinContent; j++) {
507  if (fNDrawPoints > maxOutliers) break;
508  if (IsOption(kPointsAllScat)) { //Draw outliers and "all" values scattered
511  } else { //Draw them in the "candle line"
513  if ((int)scaledBinContent == 1) //If there is only one datapoint available put it in the middle of the bin
515  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
517  }
518  if (swapXY) {
519  //Swap X and Y
520  Double_t keepCurrently;
521  keepCurrently = fDrawPointsX[fNDrawPoints];
523  fDrawPointsY[fNDrawPoints] = keepCurrently;
524  }
525  // Continue fMeans, that fNDrawPoints is not increased, so that value will not be shown
526  if (doLogX) {
528  }
529  if (doLogY) {
531  }
532  fNDrawPoints++;
533  }
534  }
535  if (fNDrawPoints > maxOutliers) { //Should never happen, due to myScale!!!
536  Error ("PaintCandlePlot","Not possible to draw all outliers.");
537  break;
538  }
539  }
540  } else { //Raw data candle
541  //If only outliers are shown, calculate myScale only based on nOutliers, use fNDatapoints (all) instead
542  if (IsOption(kPointsOutliers) && nOutliers > maxOutliers/2) {
543  myScale = nOutliers/(maxOutliers/2.);
544  } else {
545  if (fNDatapoints > maxOutliers/2) myScale = fNDatapoints/(maxOutliers/2.);
546  }
547  fNDrawPoints = 0;
548  for (int i = 0; i < fNDatapoints; i++ ) {
549  Double_t myData = fDatapoints[i];
550  Double_t maxScatter = (fWhiskerUp-fWhiskerDown)/100;
551  if (!(i % (int) myScale == 0 )) continue; //If the amount of data is too large take only every 2nd or 3rd to reduce the amount
552  // Either show them only outside the whiskers, or all of them
553  if (myData < fWhiskerDown || myData > fWhiskerUp || (GetCandleOption(5) > 1)) {
554  if (IsOption(kPointsAllScat)) { //Draw outliers and "all" values scattered
556  fDrawPointsY[fNDrawPoints] = myData + (random.Rndm() - 0.5)*maxScatter; //random +- 0.5 of candle-height
557  } else { //Draw them in the "candle line"
559  fDrawPointsY[fNDrawPoints] = myData + (random.Rndm() - 0.5)*maxScatter; //random +- 0.5 of candle-height
560  }
561  if (swapXY) {
562  //Swap X and Y
563  Double_t keepCurrently;
564  keepCurrently = fDrawPointsX[fNDrawPoints];
566  fDrawPointsY[fNDrawPoints] = keepCurrently;
567  }
568  // Continue fMeans, that fNDrawPoints is not increased, so that value will not be shown
569  if (doLogX) {
571  else continue;
572  }
573  if (doLogY) {
575  else continue;
576  }
577  fNDrawPoints++;
578  if (fNDrawPoints > maxOutliers) { //Should never happen, due to myScale!!!
579  Error ("PaintCandlePlotRaw","Not possible to draw all outliers.");
580  break;
581  }
582  }
583  }
584  }
585  }
587  //We are starting with kHistoRight, left will be modified from right later
588  if (fIsRaw) { //This is a raw-data candle
589  if (!fProj) {
590  fProj = new TH1D("hpa","hpa",100,min,max+0.0001*(max-min));
591  for (Long64_t i = 0; i < fNDatapoints; ++i) {
592  fProj->Fill(fDatapoints[i]);
593  }
594  }
595  }
596 
597  fNHistoPoints = 0;
598  Double_t maxContent = fProj->GetMaximum();
599  Double_t maxHistoHeight = fHistoWidth;
600  if (IsOption(kHistoViolin)) maxHistoHeight *= 0.5;
601 
602  bool isFirst = true;
603  int lastNonZero = 0;
604  for (int bin = 1; bin <= fProj->GetNbinsX(); bin++) {
605  if (isFirst) {
606  if (fProj->GetBinContent(bin) > 0) {
609  if (doLogX) {
611  }
612  if (doLogY) {
614  }
615  fNHistoPoints++;
616  isFirst = false;
617  } else {
618  continue;
619  }
620  }
621 
622  Double_t myBinValue = fProj->GetBinContent(bin);
623  if (doLogZ) {
624  if (myBinValue > 0) myBinValue = TMath::Log10(myBinValue); else myBinValue = 0;
625  }
626  fHistoPointsX[fNHistoPoints] = fPosCandleAxis + myBinValue/maxContent*maxHistoHeight;
628  fNHistoPoints++;
629  fHistoPointsX[fNHistoPoints] = fPosCandleAxis + myBinValue/maxContent*maxHistoHeight;
631  if (doLogX) {
634  }
635  if (doLogY) {
638  }
639 
640  fNHistoPoints++;
641  if (fProj->GetBinContent(bin) > 0) lastNonZero = fNHistoPoints;
642  }
643 
646  fNHistoPoints = lastNonZero+1; //+1 so that the line down to 0 is added as well
647 
648  if (IsOption(kHistoLeft)) {
649  for (int i = 0; i < fNHistoPoints; i++) {
651  }
652  }
653  if (IsOption(kHistoViolin)) {
654  for (int i = 0; i < fNHistoPoints; i++) {
657  }
658  fNHistoPoints *= 2;
659  }
660  }
661 
662 
663  fIsCalculated = true;
664 }
665 
666 ////////////////////////////////////////////////////////////////////////////////
667 /// Paint one candle with its current attributes.
668 
670 {
671  //If something was changed before, we need to recalculate some values
672  if (!fIsCalculated) Calculate();
673 
674  // Save the attributes as they were set originally
675  Style_t saveLine = GetLineStyle();
676  Style_t saveMarker = GetMarkerStyle();
677  Style_t saveFillStyle = GetFillStyle();
678  Style_t saveFillColor = GetFillColor();
679  Style_t saveLineColor = GetLineColor();
680 
681  Double_t dimLeft = fPosCandleAxis-0.5*fCandleWidth;
682  Double_t dimRight = fPosCandleAxis+0.5*fCandleWidth;
683 
687 
688  Bool_t swapXY = IsOption(kHorizontal);
689  Bool_t doLogY = (!(swapXY) && fLogY) || (swapXY && fLogX);
690  Bool_t doLogX = (!(swapXY) && fLogX) || (swapXY && fLogY);
691 
692  // From now on this is real painting only, no calculations anymore
693 
695  SetLineColor(saveFillColor);
698  SetLineColor(saveLineColor);
700  }
701 
702 
704  if (IsOption(kHistoZeroIndicator) && (saveFillStyle != 0)) {
705  SetLineColor(saveFillColor);
707  }
708  if (!swapXY) {
709  gPad->PaintFillArea(fNHistoPoints, fHistoPointsX, fHistoPointsY);
710  gPad->PaintPolyLine(fNHistoPoints, fHistoPointsX, fHistoPointsY);
711  } else {
712  gPad->PaintFillArea(fNHistoPoints, fHistoPointsY, fHistoPointsX);
713  gPad->PaintPolyLine(fNHistoPoints, fHistoPointsY, fHistoPointsX);
714  }
715  if (IsOption(kHistoZeroIndicator) && (saveFillStyle != 0)) {
716  SetLineColor(saveLineColor);
718  }
719  }
720 
721  if (IsOption(kBox)) { // Draw a simple box
722  if (IsOption(kMedianNotched)) { // Check if we have to draw a box with notches
723  Double_t x[] = {dimLeft, dimLeft, dimLeft+fCandleWidth/3., dimLeft, dimLeft, dimRight,
724  dimRight, dimRight-fCandleWidth/3., dimRight, dimRight, dimLeft};
727  PaintBox(11, x, y, swapXY);
728  } else { // draw a simple box
729  Double_t x[] = {dimLeft, dimLeft, dimRight, dimRight, dimLeft};
731  PaintBox(5, x, y, swapXY);
732  }
733  }
734 
735  if (IsOption(kAnchor)) { // Draw the anchor line
736  PaintLine(dimLeft, fWhiskerUp, dimRight, fWhiskerUp, swapXY);
737  PaintLine(dimLeft, fWhiskerDown, dimRight, fWhiskerDown, swapXY);
738  }
739 
740  if (IsOption(kWhiskerAll) && !IsOption(kHistoZeroIndicator)) { // Whiskers are dashed
741  SetLineStyle(2);
745  SetLineStyle(saveLine);
747  } else if ((IsOption(kWhiskerAll) && IsOption(kHistoZeroIndicator)) || IsOption(kWhisker15) ) { // Whiskers without dashing, better whisker definition, or forced when using zero line
750  }
751 
752  if (IsOption(kMedianLine)) { // Paint fMedian as a line
753  PaintLine(dimLeft, fMedian, dimRight, fMedian, swapXY);
754  } else if (IsOption(kMedianNotched)) { // Paint fMedian as a line (using notches, fMedian line is shorter)
755  PaintLine(dimLeft+fCandleWidth/3, fMedian, dimRight-fCandleWidth/3., fMedian, swapXY);
756  } else if (IsOption(kMedianCircle)) { // Paint fMedian circle
757  Double_t myMedianX[1], myMedianY[1];
758  if (!swapXY) {
759  myMedianX[0] = fPosCandleAxis;
760  myMedianY[0] = fMedian;
761  } else {
762  myMedianX[0] = fMedian;
763  myMedianY[0] = fPosCandleAxis;
764  }
765 
766  Bool_t isValid = true;
767  if (doLogX) {
768  if (myMedianX[0] > 0) myMedianX[0] = TMath::Log10(myMedianX[0]); else isValid = false;
769  }
770  if (doLogY) {
771  if (myMedianY[0] > 0) myMedianY[0] = TMath::Log10(myMedianY[0]); else isValid = false;
772  }
773 
774  SetMarkerStyle(24);
776 
777  if (isValid) gPad->PaintPolyMarker(1,myMedianX,myMedianY); // A circle for the fMedian
778 
779  SetMarkerStyle(saveMarker);
781 
782  }
783 
784  if (IsOption(kMeanCircle)) { // Paint fMean as a circle
785  Double_t myMeanX[1], myMeanY[1];
786  if (!swapXY) {
787  myMeanX[0] = fPosCandleAxis;
788  myMeanY[0] = fMean;
789  } else {
790  myMeanX[0] = fMean;
791  myMeanY[0] = fPosCandleAxis;
792  }
793 
794  Bool_t isValid = true;
795  if (doLogX) {
796  if (myMeanX[0] > 0) myMeanX[0] = TMath::Log10(myMeanX[0]); else isValid = false;
797  }
798  if (doLogY) {
799  if (myMeanY[0] > 0) myMeanY[0] = TMath::Log10(myMeanY[0]); else isValid = false;
800  }
801 
802  SetMarkerStyle(24);
804 
805  if (isValid) gPad->PaintPolyMarker(1,myMeanX,myMeanY); // A circle for the fMean
806 
807  SetMarkerStyle(saveMarker);
809 
810  } else if (IsOption(kMeanLine)) { // Paint fMean as a dashed line
811  SetLineStyle(2);
813 
814  PaintLine(dimLeft, fMean, dimRight, fMean, swapXY);
815  SetLineStyle(saveLine);
817 
818  }
819 
820  if (IsOption(kAnchor)) { //Draw standard anchor
821  PaintLine(dimLeft, fWhiskerDown, dimRight, fWhiskerDown, swapXY); // the lower anchor line
822  PaintLine(dimLeft, fWhiskerUp, dimRight, fWhiskerUp, swapXY); // the upper anchor line
823  }
824 
825  // This is a bit complex. All values here are handled as outliers. Usually
826  // only the datapoints outside the whiskers are shown.
827  // One can show them in one row as crosses, or scattered randomly. If activated
828  // all datapoint are shown in the same way
829 
830  if (GetCandleOption(5) > 0) { //Draw outliers
831  if (IsOption(kPointsAllScat)) { //Draw outliers and "all" values scattered
832  SetMarkerStyle(0);
833  } else {
834  SetMarkerStyle(5);
835  }
837  gPad->PaintPolyMarker(fNDrawPoints,fDrawPointsX, fDrawPointsY);
838  }
839 }
840 
841 ////////////////////////////////////////////////////////////////////////////////
842 /// Return true is this option is activated in fOption
843 
845  long myOpt = 9;
846  int pos = 0;
847  for (pos = 0; pos < 16; pos++) {
848  if (myOpt > opt) break;
849  else myOpt *=10;
850  }
851  myOpt /= 9;
852  int thisOpt = GetCandleOption(pos);
853 
854  return ((thisOpt * myOpt) == opt);
855 }
856 
857 ////////////////////////////////////////////////////////////////////////////////
858 /// Paint a box for candle.
859 
860 void TCandle::PaintBox(Int_t nPoints, Double_t *x, Double_t *y, Bool_t swapXY)
861 {
862  Bool_t doLogY = (!(swapXY) && fLogY) || (swapXY && fLogX);
863  Bool_t doLogX = (!(swapXY) && fLogX) || (swapXY && fLogY);
864  if (doLogY) {
865  for (int i=0; i<nPoints; i++) {
866  if (y[i] > 0) y[i] = TMath::Log10(y[i]);
867  else return;
868  }
869  }
870  if (doLogX) {
871  for (int i=0; i<nPoints; i++) {
872  if (x[i] > 0) x[i] = TMath::Log10(x[i]);
873  else return;
874  }
875  }
876  if (!swapXY) {
877  gPad->PaintFillArea(nPoints, x, y);
878  gPad->PaintPolyLine(nPoints, x, y);
879  } else {
880  gPad->PaintFillArea(nPoints, y, x);
881  gPad->PaintPolyLine(nPoints, y, x);
882  }
883 }
884 
885 ////////////////////////////////////////////////////////////////////////////////
886 /// Paint a line for candle.
887 
889 {
890  Bool_t doLogY = (!(swapXY) && fLogY) || (swapXY && fLogX);
891  Bool_t doLogX = (!(swapXY) && fLogX) || (swapXY && fLogY);
892  if (doLogY) {
893  if (y1 > 0) y1 = TMath::Log10(y1); else return;
894  if (y2 > 0) y2 = TMath::Log10(y2); else return;
895  }
896  if (doLogX) {
897  if (x1 > 0) x1 = TMath::Log10(x1); else return;
898  if (x2 > 0) x2 = TMath::Log10(x2); else return;
899  }
900  if (!swapXY) {
901  gPad->PaintLine(x1, y1, x2, y2);
902  } else {
903  gPad->PaintLine(y1, x1, y2, x2);
904  }
905 }
906 
907 ////////////////////////////////////////////////////////////////////////////////
908 /// Stream an object of class TCandle.
909 
910 void TCandle::Streamer(TBuffer &R__b)
911 {
912  if (R__b.IsReading()) {
913  UInt_t R__s, R__c;
914  Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
915  if (R__v > 3) {
916  R__b.ReadClassBuffer(TCandle::Class(), this, R__v, R__s, R__c);
917  return;
918  }
919  } else {
920  R__b.WriteClassBuffer(TCandle::Class(),this);
921  }
922 }
923 
924 ////////////////////////////////////////////////////////////////////////////////
925 /// The coordinates in the TParallelCoordVar-class are in Pad-Coordinates, so we need to convert them
926 
927 void TCandle::ConvertToPadCoords(Double_t minAxis, Double_t maxAxis, Double_t axisMinCoord, Double_t axisMaxCoord)
928 {
929  if (!fIsCalculated) Calculate();
930  Double_t a,b;
931  if (fLogY) {
932  a = TMath::Log10(minAxis);
933  b = TMath::Log10(maxAxis/minAxis);
934  } else {
935  a = minAxis;
936  b = maxAxis-minAxis;
937  }
938 
939  fMean = axisMinCoord + ((fMean-a)/b)*(axisMaxCoord-axisMinCoord);
940  fMedian = axisMinCoord + ((fMedian-a)/b)*(axisMaxCoord-axisMinCoord);
941  fMedianErr = axisMinCoord + ((fMedianErr-a)/b)*(axisMaxCoord-axisMinCoord);
942  fBoxUp = axisMinCoord + ((fBoxUp-a)/b)*(axisMaxCoord-axisMinCoord);
943  fBoxDown = axisMinCoord + ((fBoxDown-a)/b)*(axisMaxCoord-axisMinCoord);
944  fWhiskerUp = axisMinCoord + ((fWhiskerUp-a)/b)*(axisMaxCoord-axisMinCoord);
945  fWhiskerDown = axisMinCoord + ((fWhiskerDown-a)/b)*(axisMaxCoord-axisMinCoord);
946 
947  for (int i = 0; i < fNDrawPoints; i++) {
948  fDrawPointsY[i] = axisMinCoord + ((fDrawPointsY[i]-a)/b)*(axisMaxCoord-axisMinCoord);
949  }
950  for (int i = 0; i < fNHistoPoints; i++) {
951  fHistoPointsY[i] = axisMinCoord + ((fHistoPointsY[i]-a)/b)*(axisMaxCoord-axisMinCoord);
952  }
953 }
TCandle::IsViolinScaled
Bool_t IsViolinScaled()
Definition: TCandle.cxx:191
l
auto * l
Definition: textangle.C:4
TCandle::fHistoWidth
Double_t fHistoWidth
The histo width (the height of the max bin)
Definition: TCandle.h:60
TCandle::kMedianLine
@ kMedianLine
Definition: TCandle.h:32
n
const Int_t n
Definition: legend1.C:16
TCandle::fCandleWidth
Double_t fCandleWidth
The candle width.
Definition: TCandle.h:59
TCandle::fHistoPointsY
Double_t fHistoPointsY[kNMAXPOINTS]
y-coord for the polyline of the histo
Definition: TCandle.h:78
e
#define e(i)
Definition: RSha256.hxx:103
Style_t
short Style_t
Definition: RtypesCore.h:84
Version_t
short Version_t
Definition: RtypesCore.h:65
snprintf
#define snprintf
Definition: civetweb.c:1540
Option_t
const char Option_t
Definition: RtypesCore.h:66
TCandle::fAxisMin
Double_t fAxisMin
The Minimum which is visible by the axis (used by zero indicator)
Definition: TCandle.h:87
TCandle::fDrawPointsX
Double_t fDrawPointsX[kNMAXPOINTS]
x-coord for every outlier, ..
Definition: TCandle.h:73
TAttFill::Modify
virtual void Modify()
Change current fill area attributes if necessary.
Definition: TAttFill.cxx:211
TCandle::fWhiskerDown
Double_t fWhiskerDown
Position of the lower whisker end.
Definition: TCandle.h:68
TCandle::Paint
virtual void Paint(Option_t *option="")
Paint one candle with its current attributes.
Definition: TCandle.cxx:669
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
TRandom2.h
TCandle::ConvertToPadCoords
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:927
TCandle::fHistoPointsX
Double_t fHistoPointsX[kNMAXPOINTS]
x-coord for the polyline of the histo
Definition: TCandle.h:77
TCandle::fDismiss
bool fDismiss
True if the candle cannot be painted.
Definition: TCandle.h:56
TBuffer::ReadClassBuffer
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
Long64_t
long long Long64_t
Definition: RtypesCore.h:75
TCandle::fBoxUp
Double_t fBoxUp
Position of the upper box end.
Definition: TCandle.h:65
TCandle::kPointsOutliers
@ kPointsOutliers
Definition: TCandle.h:40
TCandle::fProj
TH1D * fProj
Definition: TCandle.h:55
TH1D
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:618
TCandle::fLogZ
int fLogZ
make the candle appear logz-like
Definition: TCandle.h:85
TCandle::fLogX
int fLogX
make the candle appear logx-like
Definition: TCandle.h:83
TCandle::IsOption
bool IsOption(CandleOption opt)
Return true is this option is activated in fOption.
Definition: TCandle.cxx:844
TCandle::fNHistoPoints
int fNHistoPoints
Definition: TCandle.h:79
TCandle::PaintBox
void PaintBox(Int_t nPoints, Double_t *x, Double_t *y, Bool_t swapXY)
Paint a box for candle.
Definition: TCandle.cxx:860
TH1::GetBinWidth
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width for 1D histogram.
Definition: TH1.cxx:9003
TCandle::fLogY
int fLogY
make the candle appear logy-like
Definition: TCandle.h:84
TCandle::kHistoViolin
@ kHistoViolin
Definition: TCandle.h:45
TH1::GetEntries
virtual Double_t GetEntries() const
Return the current number of entries.
Definition: TH1.cxx:4386
x
Double_t x[n]
Definition: legend1.C:17
TAttMarker::GetMarkerStyle
virtual Style_t GetMarkerStyle() const
Return the marker style.
Definition: TAttMarker.h:32
TCandle::fAxisMax
Double_t fAxisMax
The Maximum which is visible by the axis (used by zero indicator)
Definition: TCandle.h:88
TCandle.h
TBuffer
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
TAttLine::SetLineColor
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
TCandle::kHistoLeft
@ kHistoLeft
Definition: TCandle.h:43
TCandle::kHistoRight
@ kHistoRight
Definition: TCandle.h:44
Bool_t
bool Bool_t
Definition: RtypesCore.h:63
TCandle::SetScaledViolin
static void SetScaledViolin(const Bool_t vScale=true)
Static function to set scaling between violin-withs.
Definition: TCandle.cxx:237
TCandle::kMeanCircle
@ kMeanCircle
Definition: TCandle.h:36
TRandom2::Rndm
virtual Double_t Rndm()
TausWorth generator from L'Ecuyer, uses as seed 3x32bits integers Use a mask of 0xffffffffUL to make ...
Definition: TRandom2.cxx:56
TH1::RebinX
virtual TH1 * RebinX(Int_t ngroup=2, const char *newname="")
Definition: TH1.h:350
b
#define b(i)
Definition: RSha256.hxx:100
TCandle::fWhiskerRange
static Double_t fWhiskerRange
The fraction which is covered by the whiskers (0 < x < 1), default 1.
Definition: TCandle.h:90
bool
TCandle::SetScaledCandle
static void SetScaledCandle(const Bool_t cScale=true)
Static function to set scaling between candles-withs.
Definition: TCandle.cxx:228
TCandle::~TCandle
virtual ~TCandle()
TCandle default destructor.
Definition: TCandle.cxx:182
x1
static const double x1[5]
Definition: RooGaussKronrodIntegrator1D.cxx:346
TAttLine::Modify
virtual void Modify()
Change current line attributes if necessary.
Definition: TAttLine.cxx:242
TCandle::fBoxRange
static Double_t fBoxRange
The fraction which is covered by the box (0 < x < 1), default 0.5.
Definition: TCandle.h:91
TCandle::fOption
CandleOption fOption
Setting the style of the candle.
Definition: TCandle.h:81
TCandle
The candle plot painter class.
Definition: TCandle.h:26
TH1::GetMean
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:7428
TCandle::CandleOption
CandleOption
Definition: TCandle.h:29
TCandle::Calculate
void Calculate()
Calculates all values needed by the candle definition depending on the candle options.
Definition: TCandle.cxx:380
TH1::GetBinContent
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH1.cxx:4993
TBuffer.h
TAttLine
Line Attributes class.
Definition: TAttLine.h:18
TMath::Log10
Double_t Log10(Double_t x)
Definition: TMath.h:764
TAxis::GetXmin
Double_t GetXmin() const
Definition: TAxis.h:133
TH2D.h
TCandle::kWhisker15
@ kWhisker15
Definition: TCandle.h:38
TCandle::SetWhiskerRange
static void SetWhiskerRange(const Double_t wRange)
Static function to set fWhiskerRange, by setting whisker-range, one can force the whiskers to cover t...
Definition: TCandle.cxx:204
TCandle::ParseOption
int ParseOption(char *optin)
Parsing of the option-string.
Definition: TCandle.cxx:245
TAttLine::GetLineColor
virtual Color_t GetLineColor() const
Return the line color.
Definition: TAttLine.h:33
TCandle::kMeanLine
@ kMeanLine
Definition: TCandle.h:35
TCandle::TCandle
TCandle()
TCandle default constructor.
Definition: TCandle.cxx:42
a
auto * a
Definition: textangle.C:12
TCandle::SetBoxRange
static void SetBoxRange(const Double_t bRange)
Static function to set fBoxRange, by setting whisker-range, one can force the box of the candle-chart...
Definition: TCandle.cxx:218
TAttLine::GetLineStyle
virtual Style_t GetLineStyle() const
Return the line style.
Definition: TAttLine.h:34
TBuffer::WriteClassBuffer
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
TCandle::fScaledViolin
static Bool_t fScaledViolin
shall the violin or histos be scaled to each other by the maximum height?
Definition: TCandle.h:94
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:96
TCandle::fMedian
Double_t fMedian
Position of the median.
Definition: TCandle.h:63
TH1::GetBinCenter
virtual Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
Definition: TH1.cxx:8981
TH1::Fill
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3350
TAttMarker
Marker Attributes class.
Definition: TAttMarker.h:19
TCandle::fIsRaw
bool fIsRaw
0: for TH1 projection, 1: using raw data
Definition: TCandle.h:53
TRandom2
Random number generator class based on the maximally quidistributed combined Tausworthe generator by ...
Definition: TRandom2.h:27
TCandle::fWhiskerUp
Double_t fWhiskerUp
Position of the upper whisker end.
Definition: TCandle.h:67
TVirtualPad.h
UInt_t
unsigned int UInt_t
Definition: RtypesCore.h:46
TAttFill::GetFillStyle
virtual Style_t GetFillStyle() const
Return the fill area style.
Definition: TAttFill.h:31
y
Double_t y[n]
Definition: legend1.C:17
sqrt
double sqrt(double)
TAttMarker::Modify
virtual void Modify()
Change current marker attributes if necessary.
Definition: TAttMarker.cxx:314
TBuffer::ReadVersion
virtual Version_t ReadVersion(UInt_t *start=0, UInt_t *bcnt=0, const TClass *cl=0)=0
TCandle::kHorizontal
@ kHorizontal
If this bit is not set it is vertical.
Definition: TCandle.h:47
TCandle::fIsCalculated
bool fIsCalculated
Definition: TCandle.h:54
TH1::GetQuantiles
virtual Int_t GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum=0)
Compute Quantiles for this histogram Quantile x_q of a probability distribution Function F is defined...
Definition: TH1.cxx:4543
TCandle::fNDatapoints
Long64_t fNDatapoints
Number of Datapoints within this candle.
Definition: TCandle.h:71
TH1::GetBinLowEdge
virtual Double_t GetBinLowEdge(Int_t bin) const
Return bin lower edge for 1D histogram.
Definition: TH1.cxx:8992
kNMAXPOINTS
const Int_t kNMAXPOINTS
Definition: TCandle.h:24
TBuffer::IsReading
Bool_t IsReading() const
Definition: TBuffer.h:86
TCandle::fPosCandleAxis
Double_t fPosCandleAxis
x-pos for a vertical candle
Definition: TCandle.h:58
Double_t
double Double_t
Definition: RtypesCore.h:59
TCandle::fScaledCandle
static Bool_t fScaledCandle
shall the box-width be scaled to each other by the integral of a box?
Definition: TCandle.h:93
TH1::FindBin
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:3680
TCandle::IsHorizontal
Bool_t IsHorizontal()
Definition: TCandle.h:117
TCandle::kAnchor
@ kAnchor
Definition: TCandle.h:39
points
point * points
Definition: X3DBuffer.c:22
TCandle::PaintLine
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:888
TCandle::kPointsAll
@ kPointsAll
Definition: TCandle.h:41
TCandle::kWhiskerAll
@ kWhiskerAll
Definition: TCandle.h:37
TCandle::kPointsAllScat
@ kPointsAllScat
Definition: TCandle.h:42
x2
static const double x2[5]
Definition: RooGaussKronrodIntegrator1D.cxx:364
gPad
#define gPad
Definition: TVirtualPad.h:287
TAttFill::GetFillColor
virtual Color_t GetFillColor() const
Return the fill area color.
Definition: TAttFill.h:30
TCandle::kMedianCircle
@ kMedianCircle
Definition: TCandle.h:34
TCandle::fOptionStr
char fOptionStr[128]
String to draw the candle.
Definition: TCandle.h:82
TCandle::kBox
@ kBox
Definition: TCandle.h:31
TCandle::fMedianErr
Double_t fMedianErr
The size of the notch.
Definition: TCandle.h:64
TAttFill
Fill Area Attributes class.
Definition: TAttFill.h:19
TAttMarker::SetMarkerStyle
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
TCandle::GetCandleOption
int GetCandleOption(const int pos)
Definition: TCandle.h:98
TAxis::GetXmax
Double_t GetXmax() const
Definition: TAxis.h:134
TCandle::kMedianNotched
@ kMedianNotched
Definition: TCandle.h:33
TCandle::IsCandleScaled
Bool_t IsCandleScaled()
Definition: TCandle.cxx:186
Class
void Class()
Definition: Class.C:29
TCandle::fNDrawPoints
Long64_t fNDrawPoints
max number of outliers or other point to be shown
Definition: TCandle.h:75
TH1::GetXaxis
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:320
TAttLine::SetLineStyle
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition: TAttLine.h:42
TH1::GetMaximum
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:8390
TMath::Quantiles
void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob, Bool_t isSorted=kTRUE, Int_t *index=0, Int_t type=7)
Computes sample quantiles, corresponding to the given probabilities.
Definition: TMath.cxx:1183
TCandle::fDrawPointsY
Double_t fDrawPointsY[kNMAXPOINTS]
y-coord for every outlier, ..
Definition: TCandle.h:74
TCandle::kNoOption
@ kNoOption
Definition: TCandle.h:30
TCandle::fDatapoints
Double_t * fDatapoints
position of all Datapoints within this candle
Definition: TCandle.h:70
TMath::E
constexpr Double_t E()
Base of natural log:
Definition: TMath.h:96
TCandle::fMean
Double_t fMean
Position of the mean.
Definition: TCandle.h:62
TH1::GetNbinsX
virtual Int_t GetNbinsX() const
Definition: TH1.h:296
TCandle::kHistoZeroIndicator
@ kHistoZeroIndicator
Definition: TCandle.h:46
TCandle::fBoxDown
Double_t fBoxDown
Position of the lower box end.
Definition: TCandle.h:66
int
Error
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Definition: TError.cxx:187