Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TSelectorDraw.cxx
Go to the documentation of this file.
1// @(#)root/treeplayer:$Id$
2// Author: Rene Brun 08/01/2003
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 TSelectorDraw
13A specialized TSelector for TTree::Draw.
14*/
15
16#include "TSelectorDraw.h"
17#include "TROOT.h"
18#include "TH2.h"
19#include "TH3.h"
20#include "TView.h"
21#include "TGraph.h"
22#include "TPolyMarker3D.h"
23#include "TDirectory.h"
24#include "TVirtualPad.h"
25#include "TProfile.h"
26#include "TProfile2D.h"
27#include "TTreeFormulaManager.h"
28#include "TEnv.h"
29#include "TTree.h"
30#include "TCut.h"
31#include "TEntryList.h"
32#include "TEventList.h"
33#include "TEntryListArray.h"
34#include "THLimitsFinder.h"
35#include "TStyle.h"
36#include "TClass.h"
37#include "TColor.h"
38#include "strlcpy.h"
39
41
43
44////////////////////////////////////////////////////////////////////////////////
45/// Default selector constructor.
46
48{
49 fTree = nullptr;
50 fW = nullptr;
51 fValSize = 4;
52 fVal = new Double_t*[fValSize];
53 fVmin = new Double_t[fValSize];
54 fVmax = new Double_t[fValSize];
55 fNbins = new Int_t[fValSize];
56 fVarMultiple = new bool[fValSize];
58 for (Int_t i = 0; i < fValSize; ++i) {
59 fVal[i] = nullptr;
60 fVar[i] = nullptr;
61 }
62 fManager = nullptr;
63 fMultiplicity = 0;
64 fSelect = nullptr;
65 fSelectedRows = 0;
66 fDraw = 0;
67 fObject = nullptr;
68 fOldHistogram = nullptr;
69 fObjEval = false;
70 fSelectMultiple = false;
71 fCleanElist = false;
72 fTreeElist = nullptr;
73 fAction = 0;
74 fNfill = 0;
75 fDimension = 0;
76 fOldEstimate = 0;
77 fForceRead = 0;
78 fWeight = 1;
80 fTreeElistArray = nullptr;
81}
82
83////////////////////////////////////////////////////////////////////////////////
84/// Selector destructor.
85
87{
89 delete [] fVar;
90 if (fVal) {
91 for (Int_t i = 0; i < fValSize; ++i)
92 delete [] fVal[i];
93 delete [] fVal;
94 }
95 if (fVmin) delete [] fVmin;
96 if (fVmax) delete [] fVmax;
97 if (fNbins) delete [] fNbins;
98 if (fVarMultiple) delete [] fVarMultiple;
99 if (fW) delete [] fW;
100}
101
102////////////////////////////////////////////////////////////////////////////////
103/// Called every time a loop on the tree(s) starts.
104
106{
107 SetStatus(0);
108 ResetAbort();
110 fSelectedRows = 0;
111 fTree = tree;
112 fDimension = 0;
113 fAction = 0;
114
115 TObject *obj = fInput->FindObject("varexp");
116 const char *varexp0 = obj ? obj->GetTitle() : "";
117 obj = fInput->FindObject("selection");
118 const char *selection = obj ? obj->GetTitle() : "";
119 const char *option = GetOption();
120
121 TString opt, abrt;
122 char *hdefault = (char *)"htemp";
123 char *varexp = nullptr;
124 Int_t i, j, hkeep;
125 opt = option;
126 opt.ToLower();
127 fOldHistogram = nullptr;
128 TEntryList *enlist = nullptr;
129 TEventList *evlist = nullptr;
131 bool profile = false;
132 bool optSame = false;
133 bool optEnlist = false;
134 bool optEnlistArray = false;
135 bool optpara = false;
136 bool optcandle = false;
137 bool opt5d = false;
138 if (opt.Contains("same")) {
139 optSame = true;
140 opt.ReplaceAll("same", "");
141 }
142 if (opt.Contains("entrylist")) {
143 optEnlist = true;
144 if (opt.Contains("entrylistarray")) {
145 optEnlistArray = true;
146 opt.ReplaceAll("entrylistarray", "");
147 } else {
148 opt.ReplaceAll("entrylist", "");
149 }
150 }
151 if (opt.Contains("para")) {
152 optpara = true;
153 opt.ReplaceAll("para", "");
154 }
155 if (opt.Contains("candle")) {
156 optcandle = true;
157 opt.ReplaceAll("candle", "");
158 }
159 if (opt.Contains("gl5d")) {
160 opt5d = true;
161 opt.ReplaceAll("gl5d", "");
162 }
164 //input list - only TEntryList
167 if (evlist && inElist) {
168 //this is needed because the input entry list was created
169 //by the fTree from the input TEventList and is owned by the fTree.
170 //Calling GetEntryList function changes ownership and here
171 //we want fTree to still delete this entry list
172
173 inElist->SetBit(kCanDelete, true);
174 }
175 fCleanElist = false;
177
178 fTreeElistArray = inElist ? dynamic_cast<TEntryListArray*>(fTreeElist) : nullptr;
179
180
181 if (inElist && inElist->GetReapplyCut()) {
182 realSelection *= inElist->GetTitle();
183 }
184
185 // what each variable should contain:
186 // varexp0 - original expression eg "a:b>>htest"
187 // hname - name of new or old histogram
188 // hkeep - flag if to keep new histogram
189 // hnameplus - flag if to add to current histo
190 // i - length of variable expression stripped of everything including/after ">>"
191 // varexp - variable expression stripped of everything including/after ">>"
192 // fOldHistogram - pointer to hist hname
193 // elist - pointer to selection list of hname
194
195 bool canExtend = true;
196 if (optSame) canExtend = false;
197
198 Int_t nbinsx = 0, nbinsy = 0, nbinsz = 0;
199 Double_t xmin = 0, xmax = 0, ymin = 0, ymax = 0, zmin = 0, zmax = 0;
200
201 fObject = nullptr;
202 char *hname = nullptr;
204 i = 0;
205 if (varexp0 && strlen(varexp0)) {
206 for (UInt_t k = strlen(varexp0) - 1; k > 0; k--) {
207 if (varexp0[k] == '>' && varexp0[k-1] == '>') {
208 i = (int)(&(varexp0[k-1]) - varexp0); // length of varexp0 before ">>"
209 hnamealloc = &(varexp0[k+1]); // use TString to copy value of the
210 hname = (char *) hnamealloc.Data();
211 break;
212 }
213 }
214 }
215 // char *hname = (char*)strstr(varexp0,">>");
216 if (hname) {
217 hkeep = 1;
218 varexp = new char[i+1];
219 varexp[0] = 0; //necessary if i=0
220 bool hnameplus = false;
221 while (*hname == ' ') hname++;
222 if (*hname == '+') {
223 hnameplus = true;
224 hname++;
225 while (*hname == ' ') hname++; //skip ' '
226 }
227 j = strlen(hname) - 1; // skip ' ' at the end
228 while (j) {
229 if (hname[j] != ' ') break;
230 hname[j] = 0;
231 j--;
232 }
233
234 if (i) {
236
237 Int_t mustdelete = 0;
239
240 // parse things that follow the name of the histo between '(' and ')'.
241 // At this point hname contains the name of the specified histogram.
242 // Now the syntax is extended to handle an hname of the following format
243 // hname(nBIN [[,[xlow]][,xhigh]],...)
244 // so enclosed in brackets is the binning information, xlow, xhigh, and
245 // the same for the other dimensions
246
247 char *pstart; // pointer to '('
248 char *pend; // pointer to ')'
249 char *cdummy; // dummy pointer
250 int ncomma; // number of commas between '(' and ')', later number of arguments
251 int ncols; // number of columns in varexpr
252 Double_t value; // parsed value (by sscanf)
253
254 const Int_t maxvalues = 9;
255
256 pstart = strchr(hname, '(');
257 pend = strchr(hname, ')');
258 if (pstart != nullptr) { // found the bracket
259
260 mustdelete = 1;
261
262 // check that there is only one open and close bracket
263 if (pstart == strrchr(hname, '(') && pend == strrchr(hname, ')')) {
264
265 // count number of ',' between '(' and ')'
266 ncomma = 0;
267 cdummy = pstart;
268 cdummy = strchr(&cdummy[1], ',');
269 while (cdummy != nullptr) {
270 cdummy = strchr(&cdummy[1], ',');
271 ncomma++;
272 }
273
274 if (ncomma + 1 > maxvalues) {
275 Error("DrawSelect", "ncomma+1>maxvalues, ncomma=%d, maxvalues=%d", ncomma, maxvalues);
276 ncomma = maxvalues - 1;
277 }
278
279 ncomma++; // number of arguments
280 cdummy = pstart;
281
282 // number of columns
283 ncols = 1;
284 for (j = 0; j < i; j++) {
285 if (varexp[j] == ':'
286 && !((j > 0 && varexp[j-1] == ':') || varexp[j+1] == ':')
287 ) {
288 ncols++;
289 } else if (varexp[j] == '?') {
290 ncols--;
291 // We substract to compensate the later `:` of the ternary operator
292 }
293 }
294 if (ncols > 3) { // max 3 columns
295 Error("DrawSelect", "ncols > 3, ncols=%d", ncols);
296 ncols = 0;
297 }
298
299 // check dimensions before and after ">>"
300 if (ncols * 3 < ncomma) {
301 Error("DrawSelect", "ncols*3 < ncomma ncols=%d, ncomma=%d", ncols, ncomma);
302 ncomma = ncols * 3;
303 }
304
305 // scan the values one after the other
306 for (j = 0; j < ncomma; j++) {
307 cdummy++; // skip '(' or ','
308 if (sscanf(cdummy, " %lf ", &value) == 1) {
309 cdummy = strchr(&cdummy[1], ',');
310
311 switch (j) { // do certain settings depending on position of argument
312 case 0: // binning x-axis
313 nbinsx = (Int_t)value;
314 if (ncols < 2) {
315 gEnv->SetValue("Hist.Binning.1D.x", nbinsx);
316 } else if (ncols < 3) {
317 gEnv->SetValue("Hist.Binning.2D.x", nbinsx);
318 gEnv->SetValue("Hist.Binning.2D.Prof", nbinsx);
319 } else {
320 gEnv->SetValue("Hist.Binning.3D.x", nbinsx);
321 gEnv->SetValue("Hist.Binning.3D.Profx", nbinsx);
322 }
323
324 break;
325 case 1: // lower limit x-axis
326 xmin = value;
327 break;
328 case 2: // upper limit x-axis
329 xmax = value;
330 break;
331 case 3: // binning y-axis
332 nbinsy = (Int_t)value;
333 if (ncols < 3) gEnv->SetValue("Hist.Binning.2D.y", nbinsy);
334 else {
335 gEnv->SetValue("Hist.Binning.3D.y", nbinsy);
336 gEnv->SetValue("Hist.Binning.3D.Profy", nbinsy);
337 }
338 break;
339 case 4: // lower limit y-axis
340 ymin = value;
341 break;
342 case 5: // upper limit y-axis
343 ymax = value;
344 break;
345 case 6: // binning z-axis
346 nbinsz = (Int_t)value;
347 gEnv->SetValue("Hist.Binning.3D.z", nbinsz);
348 break;
349 case 7: // lower limit z-axis
350 zmin = value;
351 break;
352 case 8: // upper limit z-axis
353 zmax = value;
354 break;
355 default:
356 Error("DrawSelect", "j>8");
357 break;
358 }
359 } // if sscanf == 1
360 } // for j=0;j<ncomma;j++
361 } else {
362 Error("Begin", "Two open or close brackets found, hname=%s", hname);
363 }
364
365 // fix up hname
366 pstart[0] = '\0'; // removes things after (and including) '('
367 } // if '(' is found
368
369 j = strlen(hname) - 1; // skip ' ' at the end
370 while (j > 0) {
371 if (hname[j] != ' ') break; // skip ' ' at the end
372 hname[j] = 0;
373 j--;
374 }
375
376 TObject *oldObject = gDirectory->Get(hname); // if hname contains '(...)' the return values is NULL, which is what we want
377 fOldHistogram = oldObject ? dynamic_cast<TH1*>(oldObject) : nullptr;
378
379 if (!fOldHistogram && oldObject && !oldObject->InheritsFrom(TH1::Class())) {
380 abrt.Form("An object of type '%s' has the same name as the requested histo (%s)", oldObject->IsA()->GetName(), hname);
381 Abort(abrt);
382 delete[] varexp;
383 return;
384 }
385 if (fOldHistogram && !hnameplus) fOldHistogram->Reset(); // reset unless adding is wanted
386
387 if (mustdelete) {
388 if (gDebug) {
389 Warning("Begin", "Deleting old histogram, since (possibly new) limits and binnings have been given");
390 }
391 delete fOldHistogram; fOldHistogram=nullptr;
392 }
393
394 } else {
395 // make selection list (i.e. varexp0 starts with ">>")
397 if (optEnlist) {
398 //write into a TEntryList
399 enlist = oldObject ? dynamic_cast<TEntryList*>(oldObject) : nullptr;
400
401 if (!enlist && oldObject) {
402 abrt.Form("An object of type '%s' has the same name as the requested event list (%s)",
403 oldObject->IsA()->GetName(), hname);
404 Abort(abrt);
405 delete[] varexp;
406 return;
407 }
408 if (!enlist) {
409 if (optEnlistArray) {
410 enlist = new TEntryListArray(hname, realSelection.GetTitle());
411 } else {
412 enlist = new TEntryList(hname, realSelection.GetTitle());
413 }
414 }
415 if (enlist) {
416 if (!hnameplus) {
417 if (enlist == inElist) {
418 // We have been asked to reset the input list!!
419 // Let's set it aside for now ...
420 if (optEnlistArray) {
422 } else {
423 inElist = new TEntryList(*enlist);
424 }
425 fCleanElist = true;
427 }
428 enlist->Reset();
429 enlist->SetTitle(realSelection.GetTitle());
430 } else {
431 TCut old = enlist->GetTitle();
432 TCut upd = old || realSelection.GetTitle();
433 enlist->SetTitle(upd.GetTitle());
434 }
435 }
436 } else {
437 //write into a TEventList
438 evlist = oldObject ? dynamic_cast<TEventList*>(oldObject) : nullptr;
439
440 if (!evlist && oldObject) {
441 abrt.Form("An object of type '%s' has the same name as the requested event list (%s)",
442 oldObject->IsA()->GetName(), hname);
443 Abort(abrt);
444 delete[] varexp;
445 return;
446 }
447 if (!evlist) {
448 evlist = new TEventList(hname, realSelection.GetTitle(), 1000, 0);
449 }
450 if (evlist) {
451 if (!hnameplus) {
452 if (evlist == fTree->GetEventList()) {
453 // We have been asked to reset the input list!!
454 // Let's set it aside for now ...
455 Abort("Input and output lists are the same!");
456 delete[] varexp;
457 return;
458 }
459 evlist->Reset();
460 evlist->SetTitle(realSelection.GetTitle());
461 } else {
462 TCut old = evlist->GetTitle();
463 TCut upd = old || realSelection.GetTitle();
464 evlist->SetTitle(upd.GetTitle());
465 }
466 }
467 }
468
469 } // if (i)
470 } else { // if (hname)
471 hname = hdefault;
472 hkeep = 0;
473 const size_t varexpLen = strlen(varexp0) + 1;
474 varexp = new char[varexpLen];
476 if (gDirectory) {
478 if (fOldHistogram) { fOldHistogram->Delete(); fOldHistogram = nullptr;}
479 }
480 }
481
482 // Decode varexp and selection
483 if (!CompileVariables(varexp, realSelection.GetTitle())) {
484 abrt.Form("Variable compilation failed: {%s,%s}", varexp, realSelection.GetTitle());
485 Abort(abrt);
486 delete[] varexp;
487 return;
488 }
489 if (fDimension > 4 && !(optpara || optcandle || opt5d || opt.Contains("goff"))) {
490 Abort("Too many variables. Use the option \"para\", \"gl5d\" or \"candle\" to display more than 4 variables.");
491 delete[] varexp;
492 return;
493 }
494 if (fDimension < 2 && (optpara || optcandle)) {
495 Abort("The options \"para\" and \"candle\" require at least 2 variables.");
496 delete[] varexp;
497 return;
498 }
499
500 // In case fOldHistogram exists, check dimensionality
502 if (nsel > 1) {
503 htitle.Form("%s {%s}", varexp, selection);
504 } else {
505 htitle = varexp;
506 }
507 if (fOldHistogram) {
509 Int_t mustdelete = 0;
511 profile = true;
512 olddim = 2;
513 }
515 profile = true;
516 olddim = 3;
517 }
518 if (opt.Contains("prof") && fDimension > 1) {
519 // ignore "prof" for 1D.
520 if (!profile || olddim != fDimension) mustdelete = 1;
521 } else if (opt.Contains("col") && fDimension>2) {
522 if (olddim+1 != fDimension) mustdelete = 1;
523 } else {
524 if (olddim != fDimension) mustdelete = 1;
525 }
526 if (mustdelete) {
527 Warning("Begin", "Deleting old histogram with different dimensions");
528 delete fOldHistogram;
529 fOldHistogram = nullptr;
530 }
531 }
532
533 // Create a default canvas if none exists
534 fDraw = 0;
535 if (!gPad && !opt.Contains("goff") && fDimension > 0) {
536 gROOT->MakeDefCanvas();
537 if (!gPad) {
538 Abort("Creation of default canvas failed");
539 delete[] varexp;
540 return;
541 }
542 }
543
544 // 1-D distribution
545 TH1 *hist;
546 if (fDimension == 1) {
547 fAction = 1;
548 if (!fOldHistogram) {
549 fNbins[0] = gEnv->GetValue("Hist.Binning.1D.x", 100);
550 if (gPad && optSame) {
551 TListIter np(gPad->GetListOfPrimitives());
552 TObject *op;
553 TH1 *oldhtemp = nullptr;
554 while ((op = np()) && !oldhtemp) {
555 if (op->InheritsFrom(TH1::Class())) oldhtemp = (TH1 *)op;
556 }
557 if (oldhtemp) {
558 fNbins[0] = oldhtemp->GetXaxis()->GetNbins();
559 fVmin[0] = oldhtemp->GetXaxis()->GetXmin();
560 fVmax[0] = oldhtemp->GetXaxis()->GetXmax();
561 } else {
562 fVmin[0] = gPad->GetUxmin();
563 fVmax[0] = gPad->GetUxmax();
564 }
565 } else {
566 fAction = -1;
567 fVmin[0] = xmin;
568 fVmax[0] = xmax;
569 if (xmin < xmax) canExtend = false;
570 }
571 }
572 if (fOldHistogram) {
573 hist = fOldHistogram;
574 fNbins[0] = hist->GetXaxis()->GetNbins();
575 } else {
576 TString precision = gEnv->GetValue("Hist.Precision.1D", "float");
577 if (precision.Contains("float")) {
578 hist = new TH1F(hname, htitle.Data(), fNbins[0], fVmin[0], fVmax[0]);
579 } else {
580 hist = new TH1D(hname, htitle.Data(), fNbins[0], fVmin[0], fVmax[0]);
581 }
591 if (!hkeep) {
592 hist->GetXaxis()->SetTitle(fVar[0]->GetTitle());
593 hist->SetBit(kCanDelete);
594 if (!opt.Contains("goff")) hist->SetDirectory(nullptr);
595 }
596 if (opt.Length() && opt.Contains("e")) hist->Sumw2();
597 }
598 fVar[0]->SetAxis(hist->GetXaxis());
599 fObject = hist;
600
601 // 2-D distribution
602 } else if (fDimension == 2 && !(optpara || optcandle)) {
603 fAction = 2;
604 if (!fOldHistogram || !optSame) {
605 fNbins[0] = gEnv->GetValue("Hist.Binning.2D.y", 40);
606 fNbins[1] = gEnv->GetValue("Hist.Binning.2D.x", 40);
607 if (opt.Contains("prof")) fNbins[1] = gEnv->GetValue("Hist.Binning.2D.Prof", 100);
608 if (optSame) {
609 TH1 *oldhtemp = (TH1*)gPad->FindObject(hdefault);
610 if (oldhtemp) {
611 fNbins[1] = oldhtemp->GetXaxis()->GetNbins();
612 fVmin[1] = oldhtemp->GetXaxis()->GetXmin();
613 fVmax[1] = oldhtemp->GetXaxis()->GetXmax();
614 fNbins[0] = oldhtemp->GetYaxis()->GetNbins();
615 fVmin[0] = oldhtemp->GetYaxis()->GetXmin();
616 fVmax[0] = oldhtemp->GetYaxis()->GetXmax();
617 } else {
618 fNbins[1] = gEnv->GetValue("Hist.Binning.2D.x", 40);
619 fVmin[1] = gPad->GetUxmin();
620 fVmax[1] = gPad->GetUxmax();
621 fNbins[0] = gEnv->GetValue("Hist.Binning.2D.y", 40);
622 fVmin[0] = gPad->GetUymin();
623 fVmax[0] = gPad->GetUymax();
624 }
625 } else {
626 if (!fOldHistogram) fAction = -2;
627 fVmin[1] = xmin;
628 fVmax[1] = xmax;
629 fVmin[0] = ymin;
630 fVmax[0] = ymax;
631 if (xmin < xmax && ymin < ymax) canExtend = false;
632 }
633 }
634 if (profile || opt.Contains("prof")) {
635 TProfile *hp;
636 if (fOldHistogram) {
637 fAction = 4;
639 } else {
640 if (fAction < 0) {
641 fAction = -4;
642 fVmin[1] = xmin;
643 fVmax[1] = xmax;
644 if (xmin < xmax) canExtend = false;
645 }
646 if (fAction == 2) {
647 //we come here when option = "same prof"
648 fAction = -4;
649 TH1 *oldhtemp = (TH1*)gPad->FindObject(hdefault);
650 if (oldhtemp) {
651 fNbins[1] = oldhtemp->GetXaxis()->GetNbins();
652 fVmin[1] = oldhtemp->GetXaxis()->GetXmin();
653 fVmax[1] = oldhtemp->GetXaxis()->GetXmax();
654 }
655 }
656 if (opt.Contains("profs")) {
657 hp = new TProfile(hname, htitle.Data(), fNbins[1], fVmin[1], fVmax[1], "s");
658 } else if (opt.Contains("profi")) {
659 hp = new TProfile(hname, htitle.Data(), fNbins[1], fVmin[1], fVmax[1], "i");
660 } else if (opt.Contains("profg")) {
661 hp = new TProfile(hname, htitle.Data(), fNbins[1], fVmin[1], fVmax[1], "g");
662 } else {
663 hp = new TProfile(hname, htitle.Data(), fNbins[1], fVmin[1], fVmax[1], "");
664 }
665 if (!hkeep) {
666 hp->SetBit(kCanDelete);
667 if (!opt.Contains("goff")) hp->SetDirectory(nullptr);
668 }
669 hp->SetLineColor(fTree->GetLineColor());
670 hp->SetLineWidth(fTree->GetLineWidth());
671 hp->SetLineStyle(fTree->GetLineStyle());
672 hp->SetFillColor(fTree->GetFillColor());
673 hp->SetFillStyle(fTree->GetFillStyle());
674 hp->SetMarkerStyle(fTree->GetMarkerStyle());
675 hp->SetMarkerColor(fTree->GetMarkerColor());
676 hp->SetMarkerSize(fTree->GetMarkerSize());
677 if (canExtend) hp->SetCanExtend(TH1::kAllAxes);
678 }
679 fVar[1]->SetAxis(hp->GetXaxis());
680 fObject = hp;
681
682 } else {
683 TH2 *h2;
684 if (fOldHistogram) {
685 h2 = (TH2F*)fOldHistogram;
686 } else {
687 TString precision = gEnv->GetValue("Hist.Precision.2D", "float");
688 if (precision.Contains("float")) {
689 h2 = new TH2F(hname, htitle.Data(), fNbins[1], fVmin[1], fVmax[1], fNbins[0], fVmin[0], fVmax[0]);
690 } else {
691 h2 = new TH2D(hname, htitle.Data(), fNbins[1], fVmin[1], fVmax[1], fNbins[0], fVmin[0], fVmax[0]);
692 }
693 h2->SetLineColor(fTree->GetLineColor());
694 h2->SetLineWidth(fTree->GetLineWidth());
695 h2->SetLineStyle(fTree->GetLineStyle());
696 h2->SetFillColor(fTree->GetFillColor());
697 h2->SetFillStyle(fTree->GetFillStyle());
698 h2->SetMarkerStyle(fTree->GetMarkerStyle());
699 h2->SetMarkerColor(fTree->GetMarkerColor());
700 h2->SetMarkerSize(fTree->GetMarkerSize());
701 if (canExtend) h2->SetCanExtend(TH1::kAllAxes);
702 if (!hkeep) {
703 h2->GetXaxis()->SetTitle(fVar[1]->GetTitle());
704 h2->GetYaxis()->SetTitle(fVar[0]->GetTitle());
705 h2->SetBit(TH1::kNoStats);
706 h2->SetBit(kCanDelete);
707 if (!opt.Contains("goff")) h2->SetDirectory(nullptr);
708 }
709 }
710 fVar[0]->SetAxis(h2->GetYaxis());
711 fVar[1]->SetAxis(h2->GetXaxis());
712 bool graph = false;
713 Int_t l = opt.Length();
714 if (l == 0 || optSame) graph = true;
715 if (opt.Contains("p") || opt.Contains("*") || opt.Contains("l")) graph = true;
716 if (opt.Contains("surf") || opt.Contains("lego") || opt.Contains("cont")) graph = false;
717 if (opt.Contains("col") || opt.Contains("hist") || opt.Contains("scat")) graph = false;
718 if (opt.Contains("box")) graph = false;
719 fObject = h2;
720 if (graph) {
721 fAction = 12;
722 if (!fOldHistogram && !optSame) fAction = -12;
723 }
724 }
725
726 // 3-D distribution
727 } else if ((fDimension == 3 || fDimension == 4) && !(optpara || optcandle)) {
728 fAction = 3;
729 if (fDimension == 4) fAction = 40;
730 if (!fOldHistogram || !optSame) {
731 fNbins[0] = gEnv->GetValue("Hist.Binning.3D.z", 20);
732 fNbins[1] = gEnv->GetValue("Hist.Binning.3D.y", 20);
733 fNbins[2] = gEnv->GetValue("Hist.Binning.3D.x", 20);
734 if (fDimension == 3 && opt.Contains("prof")) {
735 fNbins[1] = gEnv->GetValue("Hist.Binning.3D.Profy", 20);
736 fNbins[2] = gEnv->GetValue("Hist.Binning.3D.Profx", 20);
737 } else if (fDimension == 3 && opt.Contains("col")) {
738 fNbins[0] = gEnv->GetValue("Hist.Binning.2D.y", 40);
739 fNbins[1] = gEnv->GetValue("Hist.Binning.2D.x", 40);
740 }
741 if (optSame) {
742 TH1 *oldhtemp = (TH1*)gPad->FindObject(hdefault);
743 if (oldhtemp) {
744 fNbins[2] = oldhtemp->GetXaxis()->GetNbins();
745 fVmin[2] = oldhtemp->GetXaxis()->GetXmin();
746 fVmax[2] = oldhtemp->GetXaxis()->GetXmax();
747 fNbins[1] = oldhtemp->GetYaxis()->GetNbins();
748 fVmin[1] = oldhtemp->GetYaxis()->GetXmin();
749 fVmax[1] = oldhtemp->GetYaxis()->GetXmax();
750 fNbins[0] = oldhtemp->GetZaxis()->GetNbins();
751 fVmin[0] = oldhtemp->GetZaxis()->GetXmin();
752 fVmax[0] = oldhtemp->GetZaxis()->GetXmax();
753 } else {
754 TView *view = gPad->GetView();
755 if (!view) {
756 Error("Begin", "You cannot use option same when no 3D view exists");
757 fVmin[0] = fVmin[1] = fVmin[2] = -1;
758 fVmax[0] = fVmax[1] = fVmax[2] = 1;
759 view = TView::CreateView(1, fVmin, fVmax);
760 }
761 Double_t *rmin = view->GetRmin();
762 Double_t *rmax = view->GetRmax();
763 fNbins[2] = gEnv->GetValue("Hist.Binning.3D.z", 20);
764 fVmin[2] = rmin[0];
765 fVmax[2] = rmax[0];
766 fNbins[1] = gEnv->GetValue("Hist.Binning.3D.y", 20);
767 fVmin[1] = rmin[1];
768 fVmax[1] = rmax[1];
769 fNbins[0] = gEnv->GetValue("Hist.Binning.3D.x", 20);
770 fVmin[0] = rmin[2];
771 fVmax[0] = rmax[2];
772 }
773 } else {
774 if (!fOldHistogram && fDimension == 3) fAction = -3;
775 fVmin[2] = xmin;
776 fVmax[2] = xmax;
777 fVmin[1] = ymin;
778 fVmax[1] = ymax;
779 fVmin[0] = zmin;
780 fVmax[0] = zmax;
781 if (xmin < xmax && ymin < ymax && zmin < zmax) canExtend = false;
782 }
783 }
784 if ((fDimension == 3) && (profile || opt.Contains("prof"))) {
785 TProfile2D *hp;
786 if (fOldHistogram) {
787 fAction = 23;
789 } else {
790 if (fAction < 0) {
791 fAction = -23;
792 fVmin[2] = xmin;
793 fVmax[2] = xmax;
794 fVmin[1] = ymin;
795 fVmax[1] = ymax;
796 if (xmin < xmax && ymin < ymax) canExtend = false;
797 }
798 if (opt.Contains("profs")) {
799 hp = new TProfile2D(hname, htitle.Data(), fNbins[2], fVmin[2], fVmax[2], fNbins[1], fVmin[1], fVmax[1], "s");
800 } else if (opt.Contains("profi")) {
801 hp = new TProfile2D(hname, htitle.Data(), fNbins[2], fVmin[2], fVmax[2], fNbins[1], fVmin[1], fVmax[1], "i");
802 } else if (opt.Contains("profg")) {
803 hp = new TProfile2D(hname, htitle.Data(), fNbins[2], fVmin[2], fVmax[2], fNbins[1], fVmin[1], fVmax[1], "g");
804 } else {
805 hp = new TProfile2D(hname, htitle.Data(), fNbins[2], fVmin[2], fVmax[2], fNbins[1], fVmin[1], fVmax[1], "");
806 }
807 if (!hkeep) {
808 hp->SetBit(kCanDelete);
809 if (!opt.Contains("goff")) hp->SetDirectory(nullptr);
810 }
811 hp->SetLineColor(fTree->GetLineColor());
812 hp->SetLineWidth(fTree->GetLineWidth());
813 hp->SetLineStyle(fTree->GetLineStyle());
814 hp->SetFillColor(fTree->GetFillColor());
815 hp->SetFillStyle(fTree->GetFillStyle());
816 hp->SetMarkerStyle(fTree->GetMarkerStyle());
817 hp->SetMarkerColor(fTree->GetMarkerColor());
818 hp->SetMarkerSize(fTree->GetMarkerSize());
819 if (canExtend) hp->SetCanExtend(TH1::kAllAxes);
820 }
821 fVar[1]->SetAxis(hp->GetYaxis());
822 fVar[2]->SetAxis(hp->GetXaxis());
823 fObject = hp;
824 } else if (fDimension == 3 && opt.Contains("col")) {
825 TH2F *h2;
826 if (fOldHistogram) {
827 h2 = (TH2F*)fOldHistogram;
828 } else {
829 h2 = new TH2F(hname, htitle.Data(), fNbins[1], fVmin[1], fVmax[1], fNbins[0], fVmin[0], fVmax[0]);
830 h2->SetLineColor(fTree->GetLineColor());
831 h2->SetLineWidth(fTree->GetLineWidth());
832 h2->SetLineStyle(fTree->GetLineStyle());
833 h2->SetFillColor(fTree->GetFillColor());
834 h2->SetFillStyle(fTree->GetFillStyle());
835 h2->SetMarkerStyle(fTree->GetMarkerStyle());
836 h2->SetMarkerColor(fTree->GetMarkerColor());
837 h2->SetMarkerSize(fTree->GetMarkerSize());
838 if (canExtend) h2->SetCanExtend(TH1::kAllAxes);
839 if (!hkeep) {
840 h2->GetXaxis()->SetTitle(fVar[1]->GetTitle());
841 h2->GetYaxis()->SetTitle(fVar[0]->GetTitle());
842 h2->GetZaxis()->SetTitle(fVar[2]->GetTitle());
843 h2->SetBit(TH1::kNoStats);
844 h2->SetBit(kCanDelete);
845 if (!opt.Contains("goff")) h2->SetDirectory(nullptr);
846 }
847 }
848 fVar[0]->SetAxis(h2->GetYaxis());
849 fVar[1]->SetAxis(h2->GetXaxis());
850 fObject = h2;
851 fAction = 33;
852 } else {
853 TH3 *h3;
854 if (fOldHistogram) {
855 h3 = (TH3F*)fOldHistogram;
856 } else {
857 TString precision = gEnv->GetValue("Hist.Precision.3D", "float");
858 if (precision.Contains("float")) {
859 h3 = new TH3F(hname, htitle.Data(), fNbins[2], fVmin[2], fVmax[2], fNbins[1], fVmin[1], fVmax[1], fNbins[0], fVmin[0], fVmax[0]);
860 } else {
861 h3 = new TH3D(hname, htitle.Data(), fNbins[2], fVmin[2], fVmax[2], fNbins[1], fVmin[1], fVmax[1], fNbins[0], fVmin[0], fVmax[0]);
862 }
863 h3->SetLineColor(fTree->GetLineColor());
864 h3->SetLineWidth(fTree->GetLineWidth());
865 h3->SetLineStyle(fTree->GetLineStyle());
866 h3->SetFillColor(fTree->GetFillColor());
867 h3->SetFillStyle(fTree->GetFillStyle());
868 h3->SetMarkerStyle(fTree->GetMarkerStyle());
869 h3->SetMarkerColor(fTree->GetMarkerColor());
870 h3->SetMarkerSize(fTree->GetMarkerSize());
871 if (canExtend) h3->SetCanExtend(TH1::kAllAxes);
872 if (!hkeep) {
873 //small correction for the title offsets in x,y to take into account the angles
874 Double_t xoffset = h3->GetXaxis()->GetTitleOffset();
875 Double_t yoffset = h3->GetYaxis()->GetTitleOffset();
876 h3->GetXaxis()->SetTitleOffset(1.2 * xoffset);
877 h3->GetYaxis()->SetTitleOffset(1.2 * yoffset);
878 h3->GetXaxis()->SetTitle(fVar[2]->GetTitle());
879 h3->GetYaxis()->SetTitle(fVar[1]->GetTitle());
880 h3->GetZaxis()->SetTitle(fVar[0]->GetTitle());
881 h3->SetBit(kCanDelete);
882 h3->SetBit(TH1::kNoStats);
883 if (!opt.Contains("goff")) h3->SetDirectory(nullptr);
884 }
885 }
886 fVar[0]->SetAxis(h3->GetZaxis());
887 fVar[1]->SetAxis(h3->GetYaxis());
888 fVar[2]->SetAxis(h3->GetXaxis());
889 fObject = h3;
891 if (optSame) noscat -= 4;
892 if (!noscat && fDimension == 3) {
893 fAction = 13;
894 if (!fOldHistogram && !optSame) fAction = -13;
895 }
896 }
897 // An Event List
898 } else if (enlist) {
899 fAction = 5;
901 fTree->SetEstimate(1);
902 fObject = enlist;
903 } else if (evlist) {
904 fAction = 5;
906 fTree->SetEstimate(1);
907 fObject = evlist;
908 } else if (optcandle || optpara || opt5d) {
909 if (optcandle) fAction = 7;
910 else if (opt5d) fAction = 8;
911 else fAction = 6;
912 }
913 if (varexp) delete[] varexp;
914 for (i = 0; i < fValSize; ++i)
915 fVarMultiple[i] = false;
916 fSelectMultiple = false;
917 for (i = 0; i < fDimension; ++i) {
918 if (fVar[i] && fVar[i]->GetMultiplicity()) fVarMultiple[i] = true;
919 }
920
922
925 fNfill = 0;
926
927 for (i = 0; i < fDimension; ++i) {
928 if (!fVal[i] && fVar[i]) {
929 fVal[i] = new Double_t[(Int_t)fTree->GetEstimate()];
930 }
931 }
932
933 if (!fW) fW = new Double_t[(Int_t)fTree->GetEstimate()];
934
935 for (i = 0; i < fValSize; ++i) {
936 fVmin[i] = DBL_MAX;
937 fVmax[i] = -DBL_MAX;
938 }
939}
940
941////////////////////////////////////////////////////////////////////////////////
942/// Delete internal buffers.
943
945{
947 for (Int_t i = 0; i < fValSize; ++i) {
948 delete fVar[i];
949 fVar[i] = nullptr;
950 }
951 delete fSelect; fSelect = nullptr;
952 fManager = nullptr;
953 fMultiplicity = 0;
954}
955
956////////////////////////////////////////////////////////////////////////////////
957/// Compile input variables and selection expression.
958///
959/// varexp is an expression of the general form e1:e2:e3
960/// where e1,etc is a formula referencing a combination of the columns
961///
962/// Example:
963///
964/// varexp = x simplest case: draw a 1-Dim distribution of column named x
965/// = sqrt(x) : draw distribution of sqrt(x)
966/// = x*y/z
967/// = y:sqrt(x) 2-Dim distribution of y versus sqrt(x)
968///
969/// selection is an expression with a combination of the columns
970///
971/// Example:
972///
973/// selection = "x<y && sqrt(z)>3.2"
974///
975/// in a selection all the C++ operators are authorized
976///
977/// Return false if any of the variable is not compilable.
978
980{
981 Int_t i, nch, ncols;
982
983 // Compile selection expression if there is one
984 fDimension = 0;
985 ClearFormula();
986 fMultiplicity = 0;
987 fObjEval = false;
988
989 if (strlen(selection)) {
990 fSelect = new TTreeFormula("Selection", selection, fTree);
991 fSelect->SetQuickLoad(true);
992 if (!fSelect->GetNdim()) {
993 delete fSelect;
994 fSelect = nullptr;
995 return false;
996 }
997 }
998
999 // if varexp is empty, take first column by default
1000 nch = strlen(varexp);
1001 if (nch == 0) {
1002 fDimension = 0;
1003 if (fSelect) {
1005 }
1007
1008 if (fManager) {
1009 fManager->Sync();
1010
1013 }
1014
1015 return true;
1016 }
1017
1018 // otherwise select only the specified columns
1019 std::vector<TString> varnames;
1021
1023
1025 if (fSelect) fManager->Add(fSelect);
1027 for (i = 0; i < ncols; ++i) {
1028 fVar[i] = new TTreeFormula(TString::Format("Var%i", i + 1), varnames[i].Data(), fTree);
1029 fVar[i]->SetQuickLoad(true);
1030 if(!fVar[i]->GetNdim()) { ClearFormula(); return false; }
1031 fManager->Add(fVar[i]);
1032 }
1033 fManager->Sync();
1034
1037
1038 fDimension = ncols;
1039
1040 if (ncols == 1) {
1041 TClass *cl = fVar[0]->EvalClass();
1042 if (cl) {
1043 fObjEval = true;
1044 }
1045 }
1046 return true;
1047}
1048
1049////////////////////////////////////////////////////////////////////////////////
1050/// Return the last values corresponding to the i-th component
1051/// of the formula being processed (where the component are ':' separated).
1052/// The actual number of entries is:
1053///
1054/// GetSelectedRows() % tree->GetEstimate()
1055///
1056/// Note GetSelectedRows currently returns the actual number of values plotted
1057/// and thus if the formula contains arrays, this number might be greater than
1058/// the number of entries in the trees.
1059///
1060/// By default TTree::Draw creates the arrays obtained
1061/// with all GetVal and GetW with a length corresponding to the
1062/// parameter fEstimate. By default fEstimate=10000 and can be modified
1063/// via TTree::SetEstimate. A possible recipe is to do
1064///
1065/// tree->SetEstimate(tree->GetEntries());
1066///
1067/// You must call SetEstimate if the expected number of selected rows
1068/// is greater than 10000.
1069///
1070/// See TTree::Draw for additional details.
1071
1073{
1075 return nullptr;
1076 else
1077 return fVal[i];
1078}
1079
1080////////////////////////////////////////////////////////////////////////////////
1081/// Return the TTreeFormula corresponding to the i-th component
1082/// of the request formula (where the component are ':' separated).
1083
1085{
1087 return nullptr;
1088 else
1089 return fVar[i];
1090}
1091
1092////////////////////////////////////////////////////////////////////////////////
1093/// Initialization of the primitive type arrays if the new size is bigger than the available space.
1094
1096{
1097 if (newsize > fValSize) {
1099 while (fValSize < newsize)
1100 fValSize *= 2; // Double the available space until it matches the new size.
1101 if (fNbins) delete [] fNbins;
1102 if (fVmin) delete [] fVmin;
1103 if (fVmax) delete [] fVmax;
1104 if (fVarMultiple) delete [] fVarMultiple;
1105
1106 fNbins = new Int_t[fValSize];
1107 fVmin = new Double_t[fValSize];
1108 fVmax = new Double_t[fValSize];
1109 fVarMultiple = new bool[fValSize];
1110
1111 for (Int_t i = 0; i < oldsize; ++i)
1112 delete [] fVal[i];
1113 delete [] fVal;
1114 delete [] fVar;
1115 fVal = new Double_t*[fValSize];
1116 fVar = new TTreeFormula*[fValSize];
1117 for (Int_t i = 0; i < fValSize; ++i) {
1118 fVal[i] = nullptr;
1119 fVar[i] = nullptr;
1120 }
1121 }
1122}
1123
1124////////////////////////////////////////////////////////////////////////////////
1125/// Build Index array for names in varexp.
1126/// This will allocated a C style array of TString and Ints
1127
1128UInt_t TSelectorDraw::SplitNames(const TString &varexp, std::vector<TString> &names)
1129{
1130 names.clear();
1131
1132 bool ternary = false;
1133 Int_t prev = 0;
1134 for (Int_t i = 0; i < varexp.Length(); i++) {
1135 if (varexp[i] == ':'
1136 && !((i > 0 && varexp[i-1] == ':') || varexp[i+1] == ':')
1137 ) {
1138 if (ternary) {
1139 ternary = false;
1140 } else {
1141 names.push_back(varexp(prev, i - prev));
1142 prev = i + 1;
1143 }
1144 }
1145 if (varexp[i] == '?') {
1146 ternary = true;
1147 }
1148 }
1149 names.push_back(varexp(prev, varexp.Length() - prev));
1150 return names.size();
1151}
1152
1153
1154////////////////////////////////////////////////////////////////////////////////
1155/// This function is called at the first entry of a new tree in a chain.
1156
1158{
1159 if (fTree) fWeight = fTree->GetWeight();
1160 if (fVar) {
1161 for (Int_t i = 0; i < fDimension; ++i) {
1162 if (fVar[i]) fVar[i]->UpdateFormulaLeaves();
1163 }
1164 }
1166 return true;
1167}
1168
1169////////////////////////////////////////////////////////////////////////////////
1170/// Called in the entry loop for all entries accepted by Select.
1171
1173{
1174 if (fObjEval) {
1176 return;
1177 }
1178
1179 if (fMultiplicity) {
1181 return;
1182 }
1183
1184 if (fNfill >= fTree->GetEstimate())
1185 fNfill = 0;
1186
1187 // simple case with no multiplicity
1188 if (fForceRead && fManager->GetNdata() <= 0) return;
1189
1190 if (fSelect) {
1192 if (!fW[fNfill]) return;
1193 } else fW[fNfill] = fWeight;
1194 if (fVal) {
1195 for (Int_t i = 0; i < fDimension; ++i) {
1196 if (fVar[i]) fVal[i][fNfill] = fVar[i]->EvalInstance(0);
1197 }
1198 }
1199 fNfill++;
1200 if (fNfill >= fTree->GetEstimate()) {
1201 TakeAction();
1202 }
1203}
1204
1205////////////////////////////////////////////////////////////////////////////////
1206/// Called in the entry loop for all entries accepted by Select.
1207/// Complex case with multiplicity.
1208
1210{
1211 if (fNfill >= fTree->GetEstimate())
1212 fNfill = 0;
1213
1214 // Grab the array size of the formulas for this entry
1216
1217 // No data at all, let's move on to the next entry.
1218 if (!ndata) return;
1219
1220 // If the entry list is a TEntryListArray, get the selected subentries for this entry
1221 TEntryList *subList = nullptr;
1222 if (fTreeElistArray) {
1224 }
1225
1227
1228 // Calculate the first values
1229 if (fSelect) {
1230 // coverity[var_deref_model] fSelectMultiple==true => fSelect != 0
1232 if (!fW[fNfill] && !fSelectMultiple) return;
1233 } else fW[fNfill] = fWeight;
1234
1235 // Always call EvalInstance(0) to insure the loading
1236 // of the branches.
1237 if (fW[fNfill] && (!subList || subList->Contains(0))) {
1238 if (fDimension == 0 && fSelectMultiple) fCurrentSubEntry = (Long64_t) 0; // to fill TEntryListArray
1239 for (Int_t i = 0; i < fDimension; ++i) {
1240 if (fVar[i]) fVal[i][fNfill] = fVar[i]->EvalInstance(0);
1241 }
1242 fNfill++;
1243 if (fNfill >= fTree->GetEstimate()) {
1244 TakeAction();
1245 }
1246 } else {
1247 for (Int_t i = 0; i < fDimension; ++i) {
1248 if (fVar[i]) fVar[i]->ResetLoading();
1249 }
1250 }
1251 Double_t ww = fW[nfill0];
1252
1253 for (Int_t i = 1; i < ndata; i++) {
1254 if (fNfill >= fTree->GetEstimate())
1255 fNfill = 0;
1256 if (subList && !subList->Contains(i))
1257 continue;
1258 if (fSelectMultiple) {
1259 // coverity[var_deref_model] fSelectMultiple==true => fSelect != 0
1260 ww = fWeight * fSelect->EvalInstance(i);
1261 if (ww == 0) continue;
1262 if (fNfill == nfill0) {
1263 for (Int_t k = 0; k < fDimension; ++k) {
1264 if (!fVarMultiple[k]) fVal[k][fNfill] = fVar[k]->EvalInstance(0);
1265 }
1266 }
1267 if (fDimension == 0) fCurrentSubEntry = (Long64_t) i; // to fill TEntryListArray
1268 }
1269 for (Int_t k = 0; k < fDimension; ++k) {
1270 if (fVarMultiple[k]) fVal[k][fNfill] = fVar[k]->EvalInstance(i);
1271 else fVal[k][fNfill] = fVal[k][nfill0];
1272 }
1273 fW[fNfill] = ww;
1274
1275 fNfill++;
1276 if (fNfill >= fTree->GetEstimate()) {
1277 TakeAction();
1278 }
1279 }
1280}
1281
1282////////////////////////////////////////////////////////////////////////////////
1283/// Called in the entry loop for all entries accepted by Select.
1284/// Case where the only variable returns an object (or pointer to).
1285
1287{
1288 // Complex case with multiplicity.
1289
1290 if (fNfill >= fTree->GetEstimate())
1291 fNfill = 0;
1292
1293 // Grab the array size of the formulas for this entry
1295
1296 // No data at all, let's move on to the next entry.
1297 if (!ndata) return;
1298
1300 Double_t ww = 0;
1301
1302 for (Int_t i = 0; i < ndata; i++) {
1303 if (i == 0) {
1304 if (fSelect) {
1306 if (!fW[fNfill] && !fSelectMultiple) return;
1307 } else fW[fNfill] = fWeight;
1308 ww = fW[nfill0];
1309 } else if (fSelectMultiple) {
1310 ww = fWeight * fSelect->EvalInstance(i);
1311 if (ww == 0) continue;
1312 }
1313 if (fDimension >= 1 && fVar[0]) {
1314 TClass *cl = fVar[0]->EvalClass();
1315 if (cl == TBits::Class()) {
1316
1317 void *obj = fVar[0]->EvalObject(i);
1318
1319 if (obj) {
1320 TBits *bits = (TBits*)obj;
1321 Int_t nbits = bits->GetNbits();
1322
1323 Int_t nextbit = -1;
1324 while (true) {
1325 nextbit = bits->FirstSetBit(nextbit + 1);
1326 if (nextbit >= nbits) break;
1327 fVal[0][fNfill] = nextbit;
1328 fW[fNfill] = ww;
1329 fNfill++;
1330 }
1331 }
1332
1333 } else {
1334
1335 if (!TestBit(kWarn)) {
1336 Warning("ProcessFillObject",
1337 "Not implemented for %s",
1338 cl ? cl->GetName() : "unknown class");
1339 SetBit(kWarn);
1340 }
1341
1342 }
1343 }
1344 if (fNfill >= fTree->GetEstimate()) {
1345 TakeAction();
1346 }
1347 }
1348
1349}
1350
1351////////////////////////////////////////////////////////////////////////////////
1352/// Set number of entries to estimate variable limits.
1353
1355{
1356 if (fVal) {
1357 for (Int_t i = 0; i < fValSize; ++i) {
1358 delete [] fVal[i];
1359 fVal[i] = nullptr;
1360 }
1361 }
1362 delete [] fW; fW = nullptr;
1363}
1364
1365////////////////////////////////////////////////////////////////////////////////
1366/// Execute action for object obj fNfill times.
1367
1369{
1370 Int_t i;
1371 //__________________________1-D histogram_______________________
1372 if (fAction == 1)((TH1*)fObject)->FillN(fNfill, fVal[0], fW);
1373 //__________________________2-D histogram_______________________
1374 else if (fAction == 2) {
1375 TH2 *h2 = (TH2*)fObject;
1376 for (i = 0; i < fNfill; i++) h2->Fill(fVal[1][i], fVal[0][i], fW[i]);
1377 }
1378 //__________________________Profile histogram_______________________
1379 else if (fAction == 4)((TProfile*)fObject)->FillN(fNfill, fVal[1], fVal[0], fW);
1380 //__________________________Event List______________________________
1381 else if (fAction == 5) {
1384 Long64_t enumb = fTree->GetTree()->GetReadEntry();
1385 enlistarray->Enter(enumb, nullptr, fCurrentSubEntry);
1386 } else if (fObject->InheritsFrom(TEntryList::Class())) {
1388 Long64_t enumb = fTree->GetTree()->GetReadEntry();
1389 enlist->Enter(enumb);
1390 } else {
1392 Long64_t enumb = fTree->GetChainOffset() + fTree->GetTree()->GetReadEntry();
1393 if (evlist->GetIndex(enumb) < 0) evlist->Enter(enumb);
1394 }
1395 }
1396 //__________________________2D scatter plot_______________________
1397 else if (fAction == 12) {
1398 TH2 *h2 = (TH2*)fObject;
1399 if (h2->CanExtendAllAxes() && h2->TestBit(kCanDelete)) {
1400 for (i = 0; i < fNfill; i++) {
1401 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1402 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1403 if (fVmin[1] > fVal[1][i]) fVmin[1] = fVal[1][i];
1404 if (fVmax[1] < fVal[1][i]) fVmax[1] = fVal[1][i];
1405 }
1406 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(h2, fVmin[1], fVmax[1], fVmin[0], fVmax[0]);
1407 }
1408 TGraph *pm = new TGraph(fNfill, fVal[1], fVal[0]);
1409 pm->SetEditable(false);
1410 pm->SetBit(kCanDelete);
1411 pm->SetMarkerStyle(fTree->GetMarkerStyle());
1412 pm->SetMarkerColor(fTree->GetMarkerColor());
1413 pm->SetMarkerSize(fTree->GetMarkerSize());
1414 pm->SetLineColor(fTree->GetLineColor());
1415 pm->SetLineWidth(fTree->GetLineWidth());
1416 pm->SetLineStyle(fTree->GetLineStyle());
1417 pm->SetFillColor(fTree->GetFillColor());
1418 pm->SetFillStyle(fTree->GetFillStyle());
1419
1420 if (!fDraw && !strstr(fOption.Data(), "goff")) {
1421 if (fOption.Length() == 0 || strcasecmp(fOption.Data(), "same") == 0) pm->Draw("p");
1422 else pm->Draw(fOption.Data());
1423 }
1424 if (!h2->TestBit(kCanDelete)) {
1425 for (i = 0; i < fNfill; i++) h2->Fill(fVal[1][i], fVal[0][i], fW[i]);
1426 }
1427 }
1428 //__________________________3D scatter plot_______________________
1429 else if (fAction == 3) {
1430 TH3 *h3 = (TH3*)fObject;
1431 if (!h3->TestBit(kCanDelete)) {
1432 for (i = 0; i < fNfill; i++) h3->Fill(fVal[2][i], fVal[1][i], fVal[0][i], fW[i]);
1433 }
1434 } else if (fAction == 13) {
1436 pm3d->SetMarkerStyle(fTree->GetMarkerStyle());
1437 pm3d->SetMarkerColor(fTree->GetMarkerColor());
1438 pm3d->SetMarkerSize(fTree->GetMarkerSize());
1439 for (i = 0; i < fNfill; i++) {
1440 pm3d->SetPoint(i, fVal[2][i], fVal[1][i], fVal[0][i]);
1441 }
1442 pm3d->Draw();
1443 TH3 *h3 = (TH3*)fObject;
1444 if (!h3->TestBit(kCanDelete)) {
1445 for (i = 0; i < fNfill; i++) h3->Fill(fVal[2][i], fVal[1][i], fVal[0][i], fW[i]);
1446 }
1447 }
1448 //__________________________3D scatter plot (3rd variable = col)__
1449 else if (fAction == 33) {
1450 TH2 *h2 = (TH2*)fObject;
1451 TakeEstimate();
1452 Int_t ncolors = gStyle->GetNumberOfColors();
1453 TObjArray *grs = (TObjArray*)h2->GetListOfFunctions()->FindObject("graphs");
1454 Int_t col;
1455 TGraph *gr;
1456 if (!grs) {
1457 grs = new TObjArray(ncolors);
1458 grs->SetOwner();
1459 grs->SetName("graphs");
1460 h2->GetListOfFunctions()->Add(grs, "P");
1461 for (col = 0; col < ncolors; col++) {
1462 gr = new TGraph();
1466 grs->AddAt(gr, col);
1467 }
1468 }
1469 h2->SetEntries(fNfill);
1470 h2->SetMinimum(fVmin[2]);
1471 h2->SetMaximum(fVmax[2]);
1472 // Fill the graphs according to the color
1473 for (i = 0; i < fNfill; i++) {
1474 col = Int_t(ncolors * ((fVal[2][i] - fVmin[2]) / (fVmax[2] - fVmin[2])));
1475 if (col < 0) col = 0;
1476 if (col > ncolors - 1) col = ncolors - 1;
1477 gr = (TGraph*)grs->UncheckedAt(col);
1478 if (gr) gr->SetPoint(gr->GetN(), fVal[1][i], fVal[0][i]);
1479 }
1480 // Remove potential empty graphs
1481 for (col = 0; col < ncolors; col++) {
1482 gr = (TGraph*)grs->At(col);
1483 if (gr && gr->GetN() <= 0) grs->Remove(gr);
1484 }
1485 }
1486 //__________________________2D Profile Histogram__________________
1487 else if (fAction == 23) {
1489 for (i = 0; i < fNfill; i++) hp2->Fill(fVal[2][i], fVal[1][i], fVal[0][i], fW[i]);
1490 }
1491 //__________________________4D scatter plot_______________________
1492 else if (fAction == 40) {
1493 TakeEstimate();
1494 TH3 *h3 = (TH3*)fObject;
1495 Int_t ncolors = gStyle->GetNumberOfColors();
1496 if (ncolors == 0) {
1498 ncolors = gStyle->GetNumberOfColors();
1499 }
1500 TObjArray *pms = (TObjArray*)h3->GetListOfFunctions()->FindObject("polymarkers");
1501 Int_t col;
1503 if (!pms) {
1504 pms = new TObjArray(ncolors);
1505 pms->SetOwner();
1506 pms->SetName("polymarkers");
1507 h3->GetListOfFunctions()->Add(pms);
1508 for (col = 0; col < ncolors; col++) {
1509 pm3d = new TPolyMarker3D();
1510 pm3d->SetMarkerColor(gStyle->GetColorPalette(col));
1511 pm3d->SetMarkerStyle(fTree->GetMarkerStyle());
1512 pm3d->SetMarkerSize(fTree->GetMarkerSize());
1513 pms->AddAt(pm3d, col);
1514 }
1515 }
1516 h3->SetEntries(fNfill);
1517 h3->SetMinimum(fVmin[3]);
1518 h3->SetMaximum(fVmax[3]);
1519 for (i = 0; i < fNfill; i++) {
1520 col = Int_t(ncolors * ((fVal[3][i] - fVmin[3]) / (fVmax[3] - fVmin[3])));
1521 if (col > ncolors-1) col = ncolors-1;
1522 if (col < 0) col = 0;
1523 pm3d = (TPolyMarker3D*)pms->UncheckedAt(col);
1524 pm3d->SetPoint(pm3d->GetLastPoint() + 1, fVal[2][i], fVal[1][i], fVal[0][i]);
1525 }
1526 }
1527 //__________________________Parallel coordinates / candle chart_______________________
1528 else if (fAction == 6 || fAction == 7) {
1529 TakeEstimate();
1530 bool candle = (fAction == 7);
1531 // Using CINT to avoid a dependency in TParallelCoord
1532 if (!fOption.Contains("goff"))
1533 gROOT->ProcessLine(TString::Format("TParallelCoord::BuildParallelCoord((TSelectorDraw*)0x%zx,0x%zx)",
1534 (size_t)this, (size_t)candle));
1535 } else if (fAction == 8) {
1536 //gROOT->ProcessLineFast(TString::Format("(new TGL5DDataSet((TTree *)0x%1x))->Draw(\"%s\");", fTree, fOption.Data()));
1537 }
1538 //__________________________something else_______________________
1539 else if (fAction < 0) {
1540 fAction = -fAction;
1541 TakeEstimate();
1542 }
1543
1544 // Do we need to update screen?
1546 if (!fTree->GetUpdate()) return;
1547 if (fSelectedRows > fDraw + fTree->GetUpdate()) {
1548 if (fDraw) gPad->Modified();
1549 else fObject->Draw(fOption.Data());
1550 gPad->Update();
1552 }
1553}
1554
1555////////////////////////////////////////////////////////////////////////////////
1556/// Estimate limits for 1-D, 2-D or 3-D objects.
1557
1559{
1560 Int_t i;
1561 Double_t rmin[3], rmax[3];
1562 Double_t vminOld[4], vmaxOld[4];
1563 for (i = 0; i < fValSize && i < 4; i++) {
1564 vminOld[i] = fVmin[i];
1565 vmaxOld[i] = fVmax[i];
1566 }
1567 for (i = 0; i < fValSize; ++i) {
1568 fVmin[i] = DBL_MAX;
1569 fVmax[i] = - DBL_MAX;
1570 }
1571 //__________________________1-D histogram_______________________
1572 if (fAction == 1) {
1573 TH1 *h1 = (TH1*)fObject;
1574 if (h1->CanExtendAllAxes()) {
1575 for (i = 0; i < fNfill; i++) {
1576 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1577 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1578 }
1579 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(h1, fVmin[0], fVmax[0]);
1580 }
1581 h1->FillN(fNfill, fVal[0], fW);
1582 //__________________________2-D histogram_______________________
1583 } else if (fAction == 2) {
1584 TH2 *h2 = (TH2*)fObject;
1585 if (h2->CanExtendAllAxes()) {
1586 for (i = 0; i < fNfill; i++) {
1587 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1588 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1589 if (fVmin[1] > fVal[1][i]) fVmin[1] = fVal[1][i];
1590 if (fVmax[1] < fVal[1][i]) fVmax[1] = fVal[1][i];
1591 }
1592 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(h2, fVmin[1], fVmax[1], fVmin[0], fVmax[0]);
1593 }
1594 for (i = 0; i < fNfill; i++) h2->Fill(fVal[1][i], fVal[0][i], fW[i]);
1595 //__________________________Profile histogram_______________________
1596 } else if (fAction == 4) {
1598 if (hp->CanExtendAllAxes()) {
1599 for (i = 0; i < fNfill; i++) {
1600 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1601 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1602 if (fVmin[1] > fVal[1][i]) fVmin[1] = fVal[1][i];
1603 if (fVmax[1] < fVal[1][i]) fVmax[1] = fVal[1][i];
1604 }
1605 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(hp, fVmin[1], fVmax[1]);
1606 }
1607 hp->FillN(fNfill, fVal[1], fVal[0], fW);
1608 //__________________________2D scatter plot_______________________
1609 } else if (fAction == 12) {
1610 TH2 *h2 = (TH2*)fObject;
1611 if (h2->CanExtendAllAxes()) {
1612 for (i = 0; i < fNfill; i++) {
1613 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1614 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1615 if (fVmin[1] > fVal[1][i]) fVmin[1] = fVal[1][i];
1616 if (fVmax[1] < fVal[1][i]) fVmax[1] = fVal[1][i];
1617 }
1618 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(h2, fVmin[1], fVmax[1], fVmin[0], fVmax[0]);
1619 // In case the new lower limits of h2 axis are 0, it is better to set them to the minimum of
1620 // the data set (which should be >0) to avoid data cut when plotting in log scale.
1621 TAxis *aX = h2->GetXaxis();
1622 TAxis *aY = h2->GetYaxis();
1623 Double_t xmin = aX->GetXmin();
1624 Double_t ymin = aY->GetXmin();
1625 if (xmin == 0 || ymin == 0) {
1626 if (aX->GetBinUpEdge(aX->FindFixBin(0.01*aX->GetBinWidth(aX->GetFirst()))) > fVmin[1]) xmin = fVmin[1];
1627 if (aY->GetBinUpEdge(aY->FindFixBin(0.01*aY->GetBinWidth(aY->GetFirst()))) > fVmin[0]) ymin = fVmin[0];
1628 h2->SetBins(aX->GetNbins(), xmin, aX->GetXmax(), aY->GetNbins(), ymin, aY->GetXmax());
1629 }
1630 }
1631
1632 if (!strstr(fOption.Data(), "same") && !strstr(fOption.Data(), "goff")) {
1633 if (!h2->TestBit(kCanDelete)) {
1634 // case like: T.Draw("y:x>>myhist")
1635 // we must draw a copy before filling the histogram h2=myhist
1636 // because h2 will be filled below and we do not want to show
1637 // the binned scatter-plot, the TGraph being better.
1638 TH1 *h2c = h2->DrawCopy(fOption.Data(),"");
1639 if (h2c) h2c->SetStats(false);
1640 } else {
1641 // case like: T.Draw("y:x")
1642 // h2 is a temporary histogram (htemp). This histogram
1643 // will be automatically deleted by TPad::Clear
1644 h2->Draw();
1645 }
1646 gPad->Update();
1647 }
1648 TGraph *pm = new TGraph(fNfill, fVal[1], fVal[0]);
1649 pm->SetEditable(false);
1650 pm->SetBit(kCanDelete);
1651 pm->SetMarkerStyle(fTree->GetMarkerStyle());
1652 pm->SetMarkerColor(fTree->GetMarkerColor());
1653 pm->SetMarkerSize(fTree->GetMarkerSize());
1654 pm->SetLineColor(fTree->GetLineColor());
1655 pm->SetLineWidth(fTree->GetLineWidth());
1656 pm->SetLineStyle(fTree->GetLineStyle());
1657 pm->SetFillColor(fTree->GetFillColor());
1658 pm->SetFillStyle(fTree->GetFillStyle());
1659 if (!fDraw && !strstr(fOption.Data(),"goff")) {
1660 if (fOption.Length() == 0 || strcasecmp(fOption.Data(),"same")==0) {
1661 pm->Draw("p");
1662 }
1663 else {
1664 TString opt = fOption;
1665 opt.ToLower();
1666 if (opt.Contains("a")) {
1667 TString temp(opt);
1668 temp.ReplaceAll("same","");
1669 if (temp.Contains("a")) {
1670 if (h2->TestBit(kCanDelete)) {
1671 // h2 will be deleted, the axis setting is delegated to only
1672 // the TGraph.
1673 h2 = nullptr;
1674 fObject = pm->GetHistogram();
1675 }
1676 }
1677 }
1678 pm->Draw(fOption.Data());
1679 }
1680 }
1681 if (h2 && !h2->TestBit(kCanDelete)) {
1682 for (i = 0; i < fNfill; i++) h2->Fill(fVal[1][i], fVal[0][i], fW[i]);
1683 }
1684 //__________________________3D scatter plot with option col_______________________
1685 } else if (fAction == 33) {
1686 TH2 *h2 = (TH2*)fObject;
1687 bool process2 = false;
1688 if (h2->CanExtendAllAxes()) {
1689 if (vminOld[2] == DBL_MAX)
1690 process2 = true;
1691 for (i = 0; i < fValSize && i < 4; i++) {
1692 fVmin[i] = vminOld[i];
1693 fVmax[i] = vmaxOld[i];
1694 }
1695 for (i = 0; i < fNfill; i++) {
1696 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1697 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1698 if (fVmin[1] > fVal[1][i]) fVmin[1] = fVal[1][i];
1699 if (fVmax[1] < fVal[1][i]) fVmax[1] = fVal[1][i];
1700 if (process2) {
1701 if (fVmin[2] > fVal[2][i]) fVmin[2] = fVal[2][i];
1702 if (fVmax[2] < fVal[2][i]) fVmax[2] = fVal[2][i];
1703 }
1704 }
1705 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(h2, fVmin[1], fVmax[1], fVmin[0], fVmax[0]);
1706 // In case the new lower limits of h2 axis are 0, it is better to set them to the minimum of
1707 // the data set (which should be >0) to avoid data cut when plotting in log scale.
1708 TAxis *aX = h2->GetXaxis();
1709 TAxis *aY = h2->GetYaxis();
1710 Double_t xmin = aX->GetXmin();
1711 Double_t ymin = aY->GetXmin();
1712 if (xmin == 0 || ymin == 0) {
1713 if (aX->GetBinUpEdge(aX->FindFixBin(0.01*aX->GetBinWidth(aX->GetFirst()))) > fVmin[1]) xmin = fVmin[1];
1714 if (aY->GetBinUpEdge(aY->FindFixBin(0.01*aY->GetBinWidth(aY->GetFirst()))) > fVmin[0]) ymin = fVmin[0];
1715 h2->SetBins(aX->GetNbins(), xmin, aX->GetXmax(), aY->GetNbins(), ymin, aY->GetXmax());
1716 }
1717 } else {
1718 for (i = 0; i < fNfill; i++) {
1719 if (fVmin[2] > fVal[2][i]) fVmin[2] = fVal[2][i];
1720 if (fVmax[2] < fVal[2][i]) fVmax[2] = fVal[2][i];
1721 }
1722 }
1723 //__________________________3D scatter plot_______________________
1724 } else if (fAction == 3 || fAction == 13) {
1725 TH3 *h3 = (TH3*)fObject;
1726 if (h3->CanExtendAllAxes()) {
1727 for (i = 0; i < fNfill; i++) {
1728 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1729 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1730 if (fVmin[1] > fVal[1][i]) fVmin[1] = fVal[1][i];
1731 if (fVmax[1] < fVal[1][i]) fVmax[1] = fVal[1][i];
1732 if (fVmin[2] > fVal[2][i]) fVmin[2] = fVal[2][i];
1733 if (fVmax[2] < fVal[2][i]) fVmax[2] = fVal[2][i];
1734 }
1735 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(h3, fVmin[2], fVmax[2], fVmin[1], fVmax[1], fVmin[0], fVmax[0]);
1736 }
1737 if (fAction == 3) {
1738 for (i = 0; i < fNfill; i++) h3->Fill(fVal[2][i], fVal[1][i], fVal[0][i], fW[i]);
1739 return;
1740 }
1741 if (!strstr(fOption.Data(), "same") && !strstr(fOption.Data(), "goff")) {
1742 if (!h3->TestBit(kCanDelete)) {
1743 // case like: T.Draw("y:x>>myhist")
1744 // we must draw a copy before filling the histogram h3=myhist
1745 // because h3 will be filled below and we do not want to show
1746 // the binned scatter-plot, the TGraph being better.
1747 TH1 *h3c = h3->DrawCopy(fOption.Data(),"");
1748 if (h3c) h3c->SetStats(false);
1749 } else {
1750 // case like: T.Draw("y:x")
1751 // h3 is a temporary histogram (htemp). This histogram
1752 // will be automatically deleted by TPad::Clear
1753 h3->Draw(fOption.Data());
1754 }
1755 gPad->Update();
1756 } else {
1757 rmin[0] = fVmin[2]; rmin[1] = fVmin[1]; rmin[2] = fVmin[0];
1758 rmax[0] = fVmax[2]; rmax[1] = fVmax[1]; rmax[2] = fVmax[0];
1759 gPad->Clear();
1760 gPad->Range(-1, -1, 1, 1);
1762 }
1764 pm3d->SetMarkerStyle(fTree->GetMarkerStyle());
1765 pm3d->SetMarkerColor(fTree->GetMarkerColor());
1766 pm3d->SetMarkerSize(fTree->GetMarkerSize());
1767 for (i = 0; i < fNfill; i++) {
1768 pm3d->SetPoint(i, fVal[2][i], fVal[1][i], fVal[0][i]);
1769 }
1770 if (!fDraw && !strstr(fOption.Data(), "goff")) pm3d->Draw();
1771 if (!h3->TestBit(kCanDelete)) {
1772 for (i = 0; i < fNfill; i++) h3->Fill(fVal[2][i], fVal[1][i], fVal[0][i], fW[i]);
1773 }
1774
1775 //__________________________2D Profile Histogram__________________
1776 } else if (fAction == 23) {
1778 if (hp->CanExtendAllAxes()) {
1779 for (i = 0; i < fNfill; i++) {
1780 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1781 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1782 if (fVmin[1] > fVal[1][i]) fVmin[1] = fVal[1][i];
1783 if (fVmax[1] < fVal[1][i]) fVmax[1] = fVal[1][i];
1784 if (fVmin[2] > fVal[2][i]) fVmin[2] = fVal[2][i];
1785 if (fVmax[2] < fVal[2][i]) fVmax[2] = fVal[2][i];
1786 }
1787 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(hp, fVmin[2], fVmax[2], fVmin[1], fVmax[1]);
1788 }
1789 for (i = 0; i < fNfill; i++) hp->Fill(fVal[2][i], fVal[1][i], fVal[0][i], fW[i]);
1790 //__________________________4D scatter plot_______________________
1791 } else if (fAction == 40) {
1792 TH3 *h3 = (TH3*)fObject;
1793 if (h3->CanExtendAllAxes()) {
1794 for (i = 0; i < fValSize && i < 4; i++) {
1795 fVmin[i] = vminOld[i];
1796 fVmax[i] = vmaxOld[i];
1797 }
1798 for (i = 0; i < fNfill; i++) {
1799 if (fVmin[0] > fVal[0][i]) fVmin[0] = fVal[0][i];
1800 if (fVmax[0] < fVal[0][i]) fVmax[0] = fVal[0][i];
1801 if (fVmin[1] > fVal[1][i]) fVmin[1] = fVal[1][i];
1802 if (fVmax[1] < fVal[1][i]) fVmax[1] = fVal[1][i];
1803 if (fVmin[2] > fVal[2][i]) fVmin[2] = fVal[2][i];
1804 if (fVmax[2] < fVal[2][i]) fVmax[2] = fVal[2][i];
1805 if (fVmin[3] > fVal[3][i]) fVmin[3] = fVal[3][i];
1806 if (fVmax[3] < fVal[3][i]) fVmax[3] = fVal[3][i];
1807 }
1808 THLimitsFinder::GetLimitsFinder()->FindGoodLimits(h3, fVmin[2], fVmax[2], fVmin[1], fVmax[1], fVmin[0], fVmax[0]);
1809 } else {
1810 for (i = 0; i < fNfill; i++) {
1811 if (fVmin[3] > fVal[3][i]) fVmin[3] = fVal[3][i];
1812 if (fVmax[3] < fVal[3][i]) fVmax[3] = fVal[3][i];
1813 }
1814 }
1815 }
1816 //__________________________Parallel coordinates plot / candle chart_______________________
1817 else if (fAction == 6 || fAction == 7) {
1818 for (i = 0; i < fDimension; ++i) {
1819 for (Long64_t entry = 0; entry < fNfill; entry++) {
1820 if (fVmin[i] > fVal[i][entry]) fVmin[i] = fVal[i][entry];
1821 if (fVmax[i] < fVal[i][entry]) fVmax[i] = fVal[i][entry];
1822 }
1823 }
1824 }
1825}
1826
1827////////////////////////////////////////////////////////////////////////////////
1828/// Called at the end of a loop on a TTree.
1829
1831{
1832 // We take action (in Process) when we reach GetEstimate but
1833 // we reset at the beginning of Process, so
1834 // if fNfill == GetEstimate(), we just took action.
1835 if (fNfill && fNfill < fTree->GetEstimate())
1836 TakeAction();
1837
1838 if ((fSelectedRows == 0) && (TestBit(kCustomHistogram) == 0)) fDraw = 1; // do not draw
1839
1841}
int Int_t
Definition RtypesCore.h:45
long long Long64_t
Definition RtypesCore.h:69
#define BIT(n)
Definition Rtypes.h:90
#define ClassImp(name)
Definition Rtypes.h:374
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
#define gDirectory
Definition TDirectory.h:384
R__EXTERN TEnv * gEnv
Definition TEnv.h:170
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t np
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
float xmin
float ymin
float xmax
float ymax
Int_t gDebug
Definition TROOT.cxx:622
#define gROOT
Definition TROOT.h:414
const Int_t kCustomHistogram
R__EXTERN TStyle * gStyle
Definition TStyle.h:442
#define gPad
virtual Int_t GetNdim() const
Definition TFormula.h:237
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 SetFillColor(Color_t fcolor)
Set the fill area color.
Definition TAttFill.h:38
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition TAttFill.h:40
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 Width_t GetLineWidth() const
Return the line width.
Definition TAttLine.h:37
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition TAttLine.h:45
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:42
virtual Style_t GetLineStyle() const
Return the line style.
Definition TAttLine.h:36
virtual Style_t GetMarkerStyle() const
Return the marker style.
Definition TAttMarker.h:33
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition TAttMarker.h:39
virtual Color_t GetMarkerColor() const
Return the marker color.
Definition TAttMarker.h:32
virtual Size_t GetMarkerSize() const
Return the marker size.
Definition TAttMarker.h:34
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition TAttMarker.h:41
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition TAttMarker.h:46
Class to manage histogram axis.
Definition TAxis.h:32
Int_t GetNbins() const
Definition TAxis.h:127
Container of bits.
Definition TBits.h:26
UInt_t GetNbits() const
Definition TBits.h:134
UInt_t FirstSetBit(UInt_t startBit=0) const
Return position of first non null bit (starting from position 0 and up)
Definition TBits.cxx:359
static TClass * Class()
TClass instances represent classes, structs and namespaces in the ROOT type system.
Definition TClass.h:84
static void InitializeColors()
Initialize colors used by the TCanvas based graphics (via TColor objects).
Definition TColor.cxx:1173
A specialized string object used for TTree selections.
Definition TCut.h:25
A list of entries and subentries in a TTree or TChain.
virtual TEntryListArray * GetSubListForEntry(Long64_t entry, TTree *tree=nullptr)
Return the list holding the subentries for the given entry or 0.
static TClass * Class()
A List of entry numbers in a TTree or TChain.
Definition TEntryList.h:26
static TClass * Class()
virtual Int_t GetValue(const char *name, Int_t dflt) const
Returns the integer value for a resource.
Definition TEnv.cxx:491
virtual void SetValue(const char *name, const char *value, EEnvLevel level=kEnvChange, const char *type=nullptr)
Set the value of a resource or create a new resource.
Definition TEnv.cxx:736
<div class="legacybox"><h2>Legacy Code</h2> TEventList is a legacy interface: there will be no bug fi...
Definition TEventList.h:31
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
virtual void SetPoint(Int_t i, Double_t x, Double_t y)
Set x and y values for point number i.
Definition TGraph.cxx:2291
Int_t GetN() const
Definition TGraph.h:131
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:698
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:650
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:59
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:8933
@ kAllAxes
Definition TH1.h:76
static TClass * Class()
virtual Int_t GetDimension() const
Definition TH1.h:300
@ kNoStats
Don't draw stats box.
Definition TH1.h:176
virtual Bool_t CanExtendAllAxes() const
Returns true if all axes are extendable.
Definition TH1.cxx:6639
virtual void Reset(Option_t *option="")
Reset this histogram: contents, errors, etc.
Definition TH1.cxx:7106
TAxis * GetXaxis()
Definition TH1.h:344
virtual UInt_t SetCanExtend(UInt_t extendBitMask)
Make the histogram axes extendable / not extendable according to the bit mask returns the previous bi...
Definition TH1.cxx:6652
virtual void FillN(Int_t ntimes, const Double_t *x, const Double_t *w, Int_t stride=1)
Fill this histogram with an array x and weights w.
Definition TH1.cxx:3419
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition TH1.cxx:9016
2-D histogram with a double per channel (see TH1 documentation)
Definition TH2.h:356
2-D histogram with a float per channel (see TH1 documentation)
Definition TH2.h:307
Service class for 2-D histogram classes.
Definition TH2.h:30
3-D histogram with a double per channel (see TH1 documentation)
Definition TH3.h:371
3-D histogram with a float per channel (see TH1 documentation)
Definition TH3.h:326
The 3-D histogram classes derived from the 1-D histogram classes.
Definition TH3.h:40
static THLimitsFinder * GetLimitsFinder()
Return pointer to the current finder.
Iterator of linked list.
Definition TList.h:191
TObject * FindObject(const char *name) const override
Find an object in this list using its name.
Definition TList.cxx:576
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:174
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
An array of TObjects.
Definition TObjArray.h:31
Mother of all ROOT objects.
Definition TObject.h:41
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:205
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
virtual void Delete(Option_t *option="")
Delete this object.
Definition TObject.cxx:267
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:864
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition TObject.cxx:543
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
virtual const char * GetTitle() const
Returns title of object.
Definition TObject.cxx:501
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition TObject.cxx:293
void ResetBit(UInt_t f)
Definition TObject.h:204
@ kCanDelete
if object in a list can be deleted
Definition TObject.h:68
A 3D polymarker.
Profile2D histograms are used to display the mean value of Z and its error for each cell in X,...
Definition TProfile2D.h:27
static TClass * Class()
Profile Histogram.
Definition TProfile.h:32
static TClass * Class()
A specialized TSelector for TTree::Draw.
void ProcessFill(Long64_t entry) override
Called in the entry loop for all entries accepted by Select.
TEntryListArray * fTreeElistArray
! Pointer to Tree Event list array
virtual void SetEstimate(Long64_t n)
Set number of entries to estimate variable limits.
virtual void InitArrays(Int_t newsize)
Initialization of the primitive type arrays if the new size is bigger than the available space.
void Terminate() override
Called at the end of a loop on a TTree.
Int_t fAction
! Action type
Int_t GetMultiplicity() const
TTreeFormulaManager * fManager
Pointer to the formula manager.
TTreeFormula * GetVar(Int_t i) const
Return the TTreeFormula corresponding to the i-th component of the request formula (where the compone...
bool * fVarMultiple
![fDimension] True if fVar[i] has a variable index
~TSelectorDraw() override
Selector destructor.
TSelectorDraw()
Default selector constructor.
Double_t * fW
![fSelectedRows]Local buffer for weights
virtual UInt_t SplitNames(const TString &varexp, std::vector< TString > &names)
Build Index array for names in varexp.
Long64_t fCurrentSubEntry
Current subentry when fSelectMultiple is true. Used to fill TEntryListArray.
TTreeFormula * fSelect
Pointer to selection formula.
Long64_t fSelectedRows
Number of selected entries.
Int_t fForceRead
Force Read flag.
virtual void ClearFormula()
Delete internal buffers.
virtual void TakeAction()
Execute action for object obj fNfill times.
Long64_t fOldEstimate
Value of Tree fEstimate when selector is called.
Double_t fWeight
Tree weight (see TTree::SetWeight)
Int_t fNfill
! Total number of histogram fills
TH1 * fOldHistogram
! Pointer to previously used histogram
virtual bool CompileVariables(const char *varexp="", const char *selection="")
Compile input variables and selection expression.
Double_t * fVmax
![fDimension] Maxima of varexp columns
Int_t fMultiplicity
Indicator of the variability of the size of entries.
TTreeFormula ** fVar
![fDimension] Array of pointers to variables formula
virtual void ProcessFillMultiple(Long64_t entry)
Called in the entry loop for all entries accepted by Select.
bool fSelectMultiple
True if selection has a variable index.
TTree * fTree
Pointer to current Tree.
Int_t fDimension
Dimension of the current expression.
TObject * fTreeElist
Pointer to Tree Event list.
virtual Double_t * GetVal(Int_t i) const
Return the last values corresponding to the i-th component of the formula being processed (where the ...
Int_t * fNbins
![fDimension] Number of bins per dimension
Double_t ** fVal
![fSelectedRows][fDimension] Local buffer for the variables
Long64_t fDraw
! Last entry loop number when object was drawn
virtual void TakeEstimate()
Estimate limits for 1-D, 2-D or 3-D objects.
bool Notify() override
This function is called at the first entry of a new tree in a chain.
bool fCleanElist
True if original Tree elist must be saved.
bool fObjEval
True if fVar1 returns an object (or pointer to).
void Begin(TTree *tree) override
Called every time a loop on the tree(s) starts.
virtual void ProcessFillObject(Long64_t entry)
Called in the entry loop for all entries accepted by Select.
Double_t * fVmin
![fDimension] Minima of varexp columns
virtual void SetStatus(Long64_t status)
Definition TSelector.h:67
TList * fInput
List of objects available during processing.
Definition TSelector.h:41
TString fOption
Option given to TTree::Process.
Definition TSelector.h:39
virtual void Abort(const char *why, EAbort what=kAbortProcess)
Abort processing.
const char * GetOption() const override
Definition TSelector.h:57
TObject * fObject
! Current object if processing object (vs. TTree)
Definition TSelector.h:40
virtual void ResetAbort()
Definition TSelector.h:74
Basic string class.
Definition TString.h:139
Ssiz_t Length() const
Definition TString.h:417
void ToLower()
Change string to lower-case.
Definition TString.cxx:1182
const char * Data() const
Definition TString.h:376
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition TString.h:704
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
Definition TString.cxx:2378
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:632
Int_t GetColorPalette(Int_t i) const
Return color number i in current palette.
Definition TStyle.cxx:1103
Int_t GetNumberOfColors() const
Return number of colors in the color palette.
Definition TStyle.cxx:1177
Used to coordinate one or more TTreeFormula objects.
virtual Int_t GetNdata(bool forceLoadDim=false)
Return number of available instances in the formulas.
virtual void Add(TTreeFormula *)
Add a new formula to the list of formulas managed The manager of the formula will be changed and the ...
virtual bool Sync()
Synchronize all the formulae.
virtual Int_t GetMultiplicity() const
Used to pass a selection expression to the Tree drawing routine.
virtual void ResetLoading()
Tell the formula that we are going to request a new entry.
TTreeFormulaManager * GetManager() const
virtual Int_t GetMultiplicity() const
T EvalInstance(Int_t i=0, const char *stringStack[]=nullptr)
Evaluate this treeformula.
void SetQuickLoad(bool quick)
virtual void UpdateFormulaLeaves()
This function is called TTreePlayer::UpdateFormulaLeaves, itself called by TChain::LoadTree when a ne...
virtual void * EvalObject(Int_t i=0)
Evaluate this treeformula.
virtual void SetAxis(TAxis *axis=nullptr)
Set the axis (in particular get the type).
virtual TClass * EvalClass(Int_t oper) const
Evaluate the class of the operation oper.
A TTree represents a columnar dataset.
Definition TTree.h:84
virtual Long64_t GetEstimate() const
Definition TTree.h:519
virtual Double_t GetWeight() const
Definition TTree.h:596
virtual TEntryList * GetEntryList()
Returns the entry list assigned to this tree.
Definition TTree.cxx:5874
virtual void SetEstimate(Long64_t nentries=1000000)
Set number of entries to estimate variable limits.
Definition TTree.cxx:9221
virtual TTree * GetTree() const
Definition TTree.h:569
virtual void SetEntryList(TEntryList *list, Option_t *opt="")
Set an EntryList.
Definition TTree.cxx:9157
TEventList * GetEventList() const
Definition TTree.h:525
@ kForceRead
Definition TTree.h:263
virtual Int_t GetUpdate() const
Definition TTree.h:573
virtual Long64_t GetChainOffset() const
Definition TTree.h:468
See TView3D.
Definition TView.h:25
virtual Double_t * GetRmax()=0
virtual Double_t * GetRmin()=0
static TView * CreateView(Int_t system=1, const Double_t *rmin=nullptr, const Double_t *rmax=nullptr)
Create a concrete default 3-d view via the plug-in manager.
Definition TView.cxx:27
TGraphErrors * gr
Definition legend1.C:25
TH1F * h1
Definition legend1.C:5
TLine l
Definition textangle.C:4