REST-for-Physics  v2.3
Rare Event Searches ToolKit for Physics
TRestDataSetPlot.cxx
1/*************************************************************************
2 * This file is part of the REST software framework. *
3 * *
4 * Copyright (C) 2016 GIFNA/TREX (University of Zaragoza) *
5 * For more information see https://gifna.unizar.es/trex *
6 * *
7 * REST is free software: you can redistribute it and/or modify *
8 * it under the terms of the GNU General Public License as published by *
9 * the Free Software Foundation, either version 3 of the License, or *
10 * (at your option) any later version. *
11 * *
12 * REST is distributed in the hope that it will be useful, *
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15 * GNU General Public License for more details. *
16 * *
17 * You should have a copy of the GNU General Public License along with *
18 * REST in $REST_PATH/LICENSE. *
19 * If not, see https://www.gnu.org/licenses/. *
20 * For the list of contributors see $REST_PATH/CREDITS. *
21 *************************************************************************/
22
272
273#include "TRestDataSetPlot.h"
274
275#include "TCanvas.h"
276#include "TDirectory.h"
277#include "TStyle.h"
278
279ClassImp(TRestDataSetPlot);
280
285
289TRestDataSetPlot::TRestDataSetPlot(const char* configFilename, std::string name)
290 : TRestMetadata(configFilename) {
291 Initialize();
294}
295
300
305void TRestDataSetPlot::Initialize() { SetSectionName(this->ClassName()); }
306
312
313 if (fDataSetName == "") fDataSetName = GetParameter("inputFileName", "");
314 if (fOutputFileName == "") fOutputFileName = GetParameter("outputFileName", "");
315
316 fCut = ReadCut(fCut);
317
318 ReadPlotInfo();
320}
321
328TRestCut* TRestDataSetPlot::ReadCut(TRestCut* cut, TiXmlElement* ele) {
329 TiXmlElement* cutele = GetElement("addCut", ele);
330 while (cutele != nullptr) {
331 std::string cutName = GetParameter("name", cutele, "");
332 if (!cutName.empty()) {
333 if (cut == nullptr) {
334 cut = (TRestCut*)InstantiateChildMetadata("TRestCut", cutName);
335 } else {
336 cut->AddCut((TRestCut*)InstantiateChildMetadata("TRestCut", cutName));
337 }
338 }
339 cutele = GetNextElement(cutele);
340 }
341
342 return cut;
343}
344
350 if (!fPanels.empty()) {
351 RESTWarning << "Plot metadata already initialized" << RESTendl;
352 }
353
354 TiXmlElement* panelele = GetElement("panel");
355 while (panelele != nullptr) {
356 std::string active = GetParameter("value", panelele, "ON");
357 if (ToUpper(active) != "ON") continue;
358
359 PanelInfo panel;
360 panel.font_size = StringToDouble(GetParameter("font_size", panelele, "0.1"));
361 panel.precision = StringToInteger(GetParameter("precision", panelele, "2"));
362
363 panel.panelCut = ReadCut(panel.panelCut, panelele);
364
365 TiXmlElement* labelele = GetElement("variable", panelele);
366 while (labelele != nullptr) {
367 std::array<std::string, 3> label;
368 label[0] = GetParameter("value", labelele, "");
369 label[1] = GetParameter("label", labelele, "");
370 label[2] = GetParameter("units", labelele, "");
371 double posX = StringToDouble(GetParameter("x", labelele, "0.1"));
372 double posY = StringToDouble(GetParameter("y", labelele, "0.1"));
373
374 panel.variablePos.push_back(std::make_pair(label, TVector2(posX, posY)));
375
376 labelele = GetNextElement(labelele);
377 }
378 TiXmlElement* metadata = GetElement("metadata", panelele);
379 while (metadata != nullptr) {
380 std::array<std::string, 3> label;
381 label[0] = GetParameter("value", metadata, "");
382 label[1] = GetParameter("label", metadata, "");
383 label[2] = GetParameter("units", metadata, "");
384 double posX = StringToDouble(GetParameter("x", metadata, "0.1"));
385 double posY = StringToDouble(GetParameter("y", metadata, "0.1"));
386
387 panel.metadataPos.push_back(std::make_pair(label, TVector2(posX, posY)));
388
389 metadata = GetNextElement(metadata);
390 }
391 TiXmlElement* observable = GetElement("observable", panelele);
392 while (observable != nullptr) {
393 std::array<std::string, 3> label;
394 label[0] = GetParameter("value", observable, "");
395 label[1] = GetParameter("label", observable, "");
396 label[2] = GetParameter("units", observable, "");
397 double posX = StringToDouble(GetParameter("x", observable, "0.1"));
398 double posY = StringToDouble(GetParameter("y", observable, "0.1"));
399
400 panel.obsPos.push_back(std::make_pair(label, TVector2(posX, posY)));
401
402 observable = GetNextElement(observable);
403 }
404 TiXmlElement* expression = GetElement("expression", panelele);
405 while (expression != nullptr) {
406 std::array<std::string, 3> label;
407 label[0] = GetParameter("value", expression, "");
408 label[1] = GetParameter("label", expression, "");
409 label[2] = GetParameter("units", expression, "");
410 double posX = StringToDouble(GetParameter("x", expression, "0.1"));
411 double posY = StringToDouble(GetParameter("y", expression, "0.1"));
412
413 panel.expPos.push_back(std::make_pair(label, TVector2(posX, posY)));
414
415 expression = GetNextElement(expression);
416 }
417
418 fPanels.push_back(panel);
419 panelele = GetNextElement(panelele);
420 }
421}
422
428 if (!fPlots.empty()) {
429 RESTWarning << "Plot metadata already initialized" << RESTendl;
430 return;
431 }
432
433 TiXmlElement* plotele = GetElement("plot");
434 while (plotele != nullptr) {
435 std::string active = GetParameter("value", plotele, "ON");
436 if (ToUpper(active) != "ON") continue;
437 int N = fPlots.size();
438 PlotInfo plot;
439 plot.name = RemoveWhiteSpaces(GetParameter("name", plotele, "plot_" + ToString(N)));
440 plot.title = GetParameter("title", plotele, plot.name);
441 plot.logX = StringToBool(GetParameter("logX", plotele, "false"));
442 plot.logY = StringToBool(GetParameter("logY", plotele, "false"));
443 plot.logZ = StringToBool(GetParameter("logZ", plotele, "false"));
444 plot.gridY = StringToBool(GetParameter("gridY", plotele, "false"));
445 plot.gridX = StringToBool(GetParameter("gridX", plotele, "false"));
446 plot.normalize = StringToDouble(GetParameter("norm", plotele, ""));
447 plot.scale = GetParameter("scale", plotele, "");
448 plot.labelX = GetParameter("xlabel", plotele, "");
449 plot.labelY = GetParameter("ylabel", plotele, "");
450 plot.marginBottom = StringToDouble(GetParameter("marginBottom", plotele, "0.15"));
451 plot.marginTop = StringToDouble(GetParameter("marginTop", plotele, "0.07"));
452 plot.marginLeft = StringToDouble(GetParameter("marginLeft", plotele, "0.25"));
453 plot.marginRight = StringToDouble(GetParameter("marginRight", plotele, "0.1"));
454 plot.legendOn = StringToBool(GetParameter("legend", plotele, "OFF"));
455 plot.stackDrawOption = GetParameter("stackOption", plotele, "nostack");
456 // plot.annotationOn = StringToBool(GetParameter("annotation", plotele, "OFF"));
457 plot.xOffset = StringToDouble(GetParameter("xOffset", plotele, "0"));
458 plot.yOffset = StringToDouble(GetParameter("yOffset", plotele, "0"));
459 plot.timeDisplay = StringToBool(GetParameter("timeDisplay", plotele, "OFF"));
460 plot.save = RemoveWhiteSpaces(GetParameter("save", plotele, ""));
461
462 TiXmlElement* histele = GetElement("histo", plotele);
463 if (histele == nullptr) {
464 histele = plotele;
465 }
466 while (histele != nullptr) {
467 HistoInfo hist;
468 hist.name = RemoveWhiteSpaces(GetParameter("name", histele, plot.name));
469 hist.drawOption = GetParameter("option", histele, "colz");
470 TiXmlElement* varele = GetElement("variable", histele);
471 while (varele != nullptr) {
472 hist.variable.push_back(GetParameter("name", varele));
473 std::string rangeStr = GetParameter("range", varele);
474 hist.range.push_back(StringTo2DVector(rangeStr));
475 hist.nBins.push_back(StringToInteger(GetParameter("nbins", varele)));
476 varele = GetNextElement(varele);
477 }
478
479 hist.lineColor = GetIDFromMapString(ColorIdMap, GetParameter("lineColor", histele, "602"));
480 hist.lineWidth = StringToInteger(GetParameter("lineWidth", histele, "1"));
481 hist.lineStyle = GetIDFromMapString(LineStyleMap, GetParameter("lineStyle", histele, "1"));
482 hist.fillStyle = GetIDFromMapString(FillStyleMap, GetParameter("fillStyle", histele, "1001"));
483 hist.fillColor = GetIDFromMapString(ColorIdMap, GetParameter("fillColor", histele, "0"));
484 hist.statistics = StringToBool(GetParameter("stats", histele, "OFF"));
485 // hist.weight = GetParameter("weight", histele, "");
486 hist.histoCut = ReadCut(hist.histoCut, histele);
487 plot.histos.push_back(hist);
488
489 if (histele == plotele) {
490 break;
491 }
492 histele = GetNextElement(histele);
493 }
494
495 fPlots.push_back(plot);
496 plotele = GetNextElement(plotele);
497 }
498}
499
506 std::vector<std::string> obsList;
507
508 // Add obserbables from global cut
509 if (fCut != nullptr) {
510 const auto paramCut = fCut->GetParamCut();
511 for (const auto& [param, condition] : paramCut) {
512 obsList.push_back(param);
513 }
514 }
515
516 // Add obserbables from plot info, both variables and cuts
517 for (const auto& plots : fPlots) {
518 for (const auto& hist : plots.histos) {
519 for (const auto& var : hist.variable) {
520 obsList.push_back(var);
521 }
522 if (hist.histoCut == nullptr) continue;
523 const auto paramCut = hist.histoCut->GetParamCut();
524 for (const auto& [param, condition] : paramCut) {
525 obsList.push_back(param);
526 }
527 }
528 }
529
530 std::map<std::string, TRestDataSet::RelevantQuantity> quantity;
531
532 for (auto& panel : fPanels) {
533 // Add metadata and observables from expression key from panel info
534 for (auto& [key, posLabel] : panel.expPos) {
535 // look for metadata which are surrounded by [] but not [[]] (variables)
536 auto&& [exp, label, units] = key;
537 std::string text = exp;
538 while (text.find_last_of('[') != std::string::npos) {
539 int squareBracketCorrector = 0;
540 size_t posOpen = text.find_last_of('[');
541 size_t posClose = text.find_first_of(']', posOpen);
542 if (posOpen != 0) {
543 if (text[posOpen - 1] == '[') {
544 squareBracketCorrector = 1;
545 }
546 }
547 std::string varOrMeta = text.substr(posOpen - squareBracketCorrector,
548 posClose + 1 - posOpen + 2 * squareBracketCorrector);
549 if (squareBracketCorrector == 0) {
550 // Here varOrMeta is a metadata
551 // Add metadata to relevant quantity from the panel info
553 quant.metadata = varOrMeta;
554 quant.strategy = "unique";
555 quantity[label] = quant;
556 }
557 text = Replace(text, varOrMeta, "1");
558 }
559
560 // look for observables (characterized by having a _ in the name)
561 while (text.find("_") != std::string::npos) {
562 size_t pos_ = text.find("_");
563 size_t beginning = text.find_last_of(" -+*/)(^%", pos_) + 1;
564 size_t end = text.find_first_of(" -+*/)(^%", pos_);
565 std::string obs = text.substr(beginning, end - beginning);
566 text = Replace(text, obs, "1");
567 obsList.push_back(obs);
568 }
569
570 if (!(isAExpression(text) || isANumber(text)))
571 RESTWarning << "The expression " << exp
572 << " has not been correctly parsed into variables, metadata and observables"
573 << RESTendl;
574 }
575
576 // Add observables from observable key of panel info
577 for (auto& [key, posLabel] : panel.obsPos) {
578 auto&& [obs, label, units] = key;
579 obsList.push_back(obs);
580 }
581 // Add relevant quantity to metadata from the panel info
582 for (auto& [key, posLabel] : panel.metadataPos) {
583 auto&& [metadata, label, units] = key;
585 quant.metadata = metadata;
586 quant.strategy = "unique";
587 quantity[label] = quant;
588 }
589 }
590
591 // Remove duplicated observables if any
592 std::sort(obsList.begin(), obsList.end());
593 obsList.erase(std::unique(obsList.begin(), obsList.end()), obsList.end());
594 dataSet.SetFilePattern(fDataSetName);
595 dataSet.SetObservablesList(obsList);
596 dataSet.SetQuantity(quantity);
597 dataSet.GenerateDataSet();
598}
599
605 TRestDataSet dataSet;
606 dataSet.EnableMultiThreading(true);
607
608 // Import dataSet
609 dataSet.Import(fDataSetName);
610
611 // If dataSet is not valid, try to generate it, but note that this is deprecated
612 if (dataSet.GetTree() == nullptr) {
613 RESTWarning << "Cannot import dataSet, trying to generate it with pattern " << fDataSetName
614 << RESTendl;
615 RESTWarning << "Note that the generation of a dataSet inside TRestDataSetPlot is deplecated. Check "
616 "TRestDataSet documentation to generate a dataSet"
617 << RESTendl;
619 if (dataSet.GetTree() == nullptr) {
620 RESTError << "Cannot generate dataSet " << RESTendl;
621 exit(1);
622 }
623 }
624
625 // Perform the global cut over the dataSet
626 dataSet.SetDataFrame(dataSet.MakeCut(fCut));
627
628 TCanvas combinedCanvas(this->GetName(), this->GetName(), 0, 0, fCanvasSize.X(), fCanvasSize.Y());
629 combinedCanvas.Divide((Int_t)fCanvasDivisions.X(), (Int_t)fCanvasDivisions.Y(),
630 fCanvasDivisionMargins.X(), fCanvasDivisionMargins.Y());
631
632 gStyle->SetPalette(fPaletteStyle);
633
634 const double duration = dataSet.GetTotalTimeInSeconds();
635 const double startTime = dataSet.GetStartTime();
636 const double endTime = dataSet.GetEndTime();
637
638 // Fill general variables for the panel
639 std::map<std::string, std::string> paramMap;
640 paramMap["[[startTime]]"] = ToDateTimeString(startTime);
641 paramMap["[[endTime]]"] = ToDateTimeString(endTime);
642
643 // DataSet quantity is used to replace metadata parameters
644 const auto quantity = dataSet.GetQuantity();
645
646 int canvasIndex = 1;
647
648 for (auto& panel : fPanels) {
649 combinedCanvas.cd(canvasIndex);
650 // Gets a dataFrame with the panel cut
651 auto dataFrame = dataSet.MakeCut(panel.panelCut);
652 const int entries = *dataFrame.Count();
653 const double meanRate = entries / duration;
654 const double runLength = duration / 3600.;
655 paramMap["[[runLength]]"] = StringWithPrecision(runLength, panel.precision);
656 paramMap["[[entries]]"] = StringWithPrecision(entries, panel.precision);
657 paramMap["[[meanRate]]"] = StringWithPrecision(meanRate, panel.precision);
658
659 paramMap["[[cutNames]]"] = "";
660 paramMap["[[cuts]]"] = "";
661 if (fCut) {
662 for (const auto& cut : fCut->GetCuts()) {
663 if (paramMap["[[cutNames]]"].empty())
664 paramMap["[[cutNames]]"] += cut.GetName();
665 else
666 paramMap["[[cutNames]]"] += "," + (std::string)cut.GetName();
667 if (paramMap["[[cuts]]"].empty())
668 paramMap["[[cuts]]"] += cut.GetTitle();
669 else
670 paramMap["[[cuts]]"] += " && " + (std::string)cut.GetTitle();
671 }
672 }
673
674 paramMap["[[panelCutNames]]"] = "";
675 paramMap["[[panelCuts]]"] = "";
676 if (panel.panelCut) {
677 for (const auto& cut : panel.panelCut->GetCuts()) {
678 if (paramMap["[[panelCutNames]]"].empty())
679 paramMap["[[panelCutNames]]"] += cut.GetName();
680 else
681 paramMap["[[panelCutNames]]"] += "," + (std::string)cut.GetName();
682 if (paramMap["[[panelCuts]]"].empty())
683 paramMap["[[panelCuts]]"] += cut.GetTitle();
684 else
685 paramMap["[[panelCuts]]"] += " && " + (std::string)cut.GetTitle();
686 }
687 }
688
689 RESTInfo << "Global cuts: " << paramMap["[[cuts]]"] << RESTendl;
690 if (!paramMap["[[panelCuts]]"].empty())
691 RESTInfo << "Additional panel cuts: " << paramMap["[[panelCuts]]"] << RESTendl;
692
693 // Replace panel variables and generate a TLatex label
694 for (const auto& [key, posLabel] : panel.variablePos) {
695 auto&& [variable, label, units] = key;
696 bool found = false;
697 std::string var = variable;
698 if (!var.empty()) {
699 for (const auto& [param, val] : paramMap) {
700 if (var == param) {
701 size_t pos = 0;
702 var = Replace(var, param, val, pos);
703 found = true;
704 break;
705 }
706 }
707 if (!found) RESTWarning << "Variable " << variable << " not found" << RESTendl;
708 }
709 std::string lab = label + ": " + StringWithPrecision(var, panel.precision) + " " + units;
710 panel.text.emplace_back(new TLatex(posLabel.X(), posLabel.Y(), lab.c_str()));
711 }
712
713 // Replace metadata variables and generate a TLatex label
714 for (const auto& [key, posLabel] : panel.metadataPos) {
715 auto&& [metadata, label, units] = key;
716 std::string value = "";
717
718 for (const auto& [name, quant] : quantity) {
719 if (quant.metadata == metadata) value = quant.value;
720 }
721
722 if (value.empty()) {
723 RESTWarning << "Metadata quantity " << metadata << " not found in dataSet" << RESTendl;
724 continue;
725 }
726
727 std::string lab = label + ": " + StringWithPrecision(value, panel.precision) + " " + units;
728 panel.text.emplace_back(new TLatex(posLabel.X(), posLabel.Y(), lab.c_str()));
729 }
730
731 // Replace observable variables and generate a TLatex label
732 for (const auto& [key, posLabel] : panel.obsPos) {
733 auto&& [obs, label, units] = key;
734 auto value = *dataFrame.Mean(obs);
735
736 std::string lab = label + ": " + StringWithPrecision(value, panel.precision) + " " + units;
737 panel.text.emplace_back(new TLatex(posLabel.X(), posLabel.Y(), lab.c_str()));
738 }
739
740 // Replace any expression and generate TLatex label
741 for (const auto& [key, posLabel] : panel.expPos) {
742 auto&& [text, label, units] = key;
743 std::string var = text;
744
745 // replace variables
746 for (const auto& [param, val] : paramMap) {
747 var = Replace(var, param, val);
748 }
749 // replace metadata
750 for (const auto& [name, quant] : quantity) {
751 var = Replace(var, name, quant.value);
752 }
753 // replace observables
754 for (const auto& obs : dataFrame.GetColumnNames()) {
755 if (var.find(obs) == std::string::npos) continue;
756 // here there should be a checking that the mean(obs) can be calculated
757 // (checking obs data type?)
758 double value = *dataFrame.Mean(obs);
759 var = Replace(var, obs, DoubleToString(value));
760 }
761 var = Replace(var, "[", "(");
762 var = Replace(var, "]", ")");
763 var = EvaluateExpression(var);
764
765 std::string lab = label + ": " + var + " " + units;
766 panel.text.emplace_back(new TLatex(posLabel.X(), posLabel.Y(), lab.c_str()));
767 }
768
769 // Draw the labels inside the pad
770 for (const auto& text : panel.text) {
771 text->SetTextColor(1);
772 text->SetTextSize(panel.font_size);
773 text->Draw("same");
774 }
775 canvasIndex++;
776 }
777
778 for (auto& plots : fPlots) {
779 // Histograms are added to a THStack and will be ploted later on
780 combinedCanvas.cd(canvasIndex);
781 plots.hs = new THStack(plots.name.c_str(), plots.title.c_str());
782 if (plots.legendOn) plots.legend = new TLegend(fLegendX1, fLegendY1, fLegendX2, fLegendY2);
784 for (auto& hist : plots.histos) {
785 auto dataFrame = dataSet.MakeCut(hist.histoCut);
786 if (hist.variable.front() == "timeStamp") {
787 hist.range.front().SetX(startTime);
788 hist.range.front().SetY(endTime);
789 }
790 // 1-D Histograms
791 if (hist.variable.size() == 1) {
792 auto histo = dataFrame.Histo1D({hist.name.c_str(), hist.name.c_str(), hist.nBins.front(),
793 hist.range.front().X(), hist.range.front().Y()},
794 hist.variable.front());
795 hist.histo = static_cast<TH1*>(histo->DrawClone());
796 // 2-D Histograms
797 } else if (hist.variable.size() == 2) {
798 auto histo = dataFrame.Histo2D(
799 {hist.name.c_str(), hist.name.c_str(), hist.nBins.front(), hist.range.front().X(),
800 hist.range.front().Y(), hist.nBins.back(), hist.range.back().X(), hist.range.back().Y()},
801 hist.variable.front(), hist.variable.back());
802 hist.histo = static_cast<TH1*>(histo->DrawClone());
803 } else {
804 RESTError << "Only 1D or 2D histograms are supported " << RESTendl;
805 continue;
806 }
807 hist.histo->SetLineColor(hist.lineColor);
808 hist.histo->SetLineWidth(hist.lineWidth);
809 hist.histo->SetLineStyle(hist.lineStyle);
810 hist.histo->SetFillColor(hist.fillColor);
811 hist.histo->SetFillStyle(hist.fillStyle);
812 // If stats are on histos must we drawn
813 if (hist.statistics) {
814 hist.histo->SetStats(true);
815 hist.histo->Draw();
816 combinedCanvas.Update();
817 } else {
818 hist.histo->SetStats(false);
819 }
820 // Normalize histos
821 if (plots.normalize > 0) {
822 const double integral = hist.histo->Integral();
823 if (integral > 0) hist.histo->Scale(plots.normalize / integral);
824 }
825 // Scale histos
826 if (plots.scale != "") {
827 std::string inputScale = plots.scale;
828 double binSize = hist.histo->GetXaxis()->GetBinWidth(1);
829 double entries = hist.histo->GetEntries();
830 double runLength = dataSet.GetTotalTimeInSeconds() / 3600.; // in hours
831 double integral = hist.histo->Integral("width");
832
833 inputScale = Replace(inputScale, "binSize", DoubleToString(binSize));
834 inputScale = Replace(inputScale, "entries", DoubleToString(entries));
835 inputScale = Replace(inputScale, "runLength", DoubleToString(runLength));
836 inputScale = Replace(inputScale, "integral", DoubleToString(integral));
837
838 std::string scale = "1./(" + inputScale + ")";
839 hist.histo->Scale(StringToDouble(EvaluateExpression(scale))); // -1 if 'scale' isn't valid
840 }
841
842 // Add histos to the THStack
843 plots.hs->Add(hist.histo, hist.drawOption.c_str());
844 // Add histos to the legend
845 if (plots.legend != nullptr) plots.legend->AddEntry(hist.histo, hist.histo->GetName(), "lf");
846 }
847 }
848
849 // This function do the actual drawing of the THStack with the different options
850 for (auto& plots : fPlots) {
851 if (plots.hs == nullptr) continue;
852 // TPad parameters
853 TPad* targetPad = (TPad*)combinedCanvas.cd(canvasIndex);
854 targetPad->SetLogx(plots.logX);
855 targetPad->SetLogy(plots.logY);
856 targetPad->SetLogz(plots.logZ);
857 targetPad->SetGridx(plots.gridX);
858 targetPad->SetGridy(plots.gridY);
859 targetPad->SetLeftMargin(plots.marginLeft);
860 targetPad->SetRightMargin(plots.marginRight);
861 targetPad->SetBottomMargin(plots.marginBottom);
862 targetPad->SetTopMargin(plots.marginTop);
863
864 // HStack draw parameters
865 plots.hs->Draw(plots.stackDrawOption.c_str());
866 plots.hs->GetXaxis()->SetTitle(plots.labelX.c_str());
867 plots.hs->GetYaxis()->SetTitle(plots.labelY.c_str());
868 plots.hs->GetXaxis()->SetLabelSize(1.1 * plots.hs->GetXaxis()->GetLabelSize());
869 plots.hs->GetYaxis()->SetLabelSize(1.1 * plots.hs->GetYaxis()->GetLabelSize());
870 plots.hs->GetXaxis()->SetTitleSize(1.1 * plots.hs->GetXaxis()->GetTitleSize());
871 plots.hs->GetYaxis()->SetTitleSize(1.1 * plots.hs->GetYaxis()->GetTitleSize());
872
873 if (plots.timeDisplay) plots.hs->GetXaxis()->SetTimeDisplay(1);
874 if (plots.legend != nullptr) plots.legend->Draw();
875
876 targetPad->Update();
877 combinedCanvas.Update();
878 canvasIndex++;
879 }
880
881 // Preview plot. User can make some changed before saving
882 if (!REST_Display_CompatibilityMode && fPreviewPlot) {
883 combinedCanvas.Resize();
884 GetChar();
885 }
886
887 // Save single pads if save is marked
888 for (auto& plots : fPlots) {
889 if (plots.save.empty()) continue;
890 std::unique_ptr<TCanvas> canvas(new TCanvas());
891 canvas->SetLogx(plots.logX);
892 canvas->SetLogy(plots.logY);
893 canvas->SetLogz(plots.logZ);
894 canvas->SetGridx(plots.gridX);
895 canvas->SetGridy(plots.gridY);
896 canvas->SetLeftMargin(plots.marginLeft);
897 canvas->SetRightMargin(plots.marginRight);
898 canvas->SetBottomMargin(plots.marginBottom);
899 canvas->SetTopMargin(plots.marginTop);
900 plots.hs->Draw(plots.stackDrawOption.c_str());
901 canvas->Print(plots.save.c_str());
902 }
903
904 // Save combined canvas
905 if (!fOutputFileName.empty()) {
906 for (const auto& [name, quant] : quantity) {
907 size_t pos = 0;
908 fOutputFileName = Replace(fOutputFileName, quant.metadata, quant.value, pos);
909 }
910 combinedCanvas.Print(fOutputFileName.c_str());
911 // In case of root file save also the histograms
913 std::unique_ptr<TFile> f(TFile::Open(fOutputFileName.c_str(), "UPDATE"));
914 for (auto& plots : fPlots) {
915 for (auto& hist : plots.histos) {
916 hist.histo->Write();
917 }
918 }
919 this->Write();
920 }
921 }
922
923 CleanUp();
924}
925
931 for (auto& plots : fPlots) {
932 for (auto& hist : plots.histos) {
933 delete hist.histo;
934 }
935 delete plots.hs;
936 delete plots.legend;
937 }
938
939 for (auto& panel : fPanels) {
940 for (auto& text : panel.text) {
941 delete text;
942 }
943 panel.text.clear();
944 }
945}
946
952Int_t TRestDataSetPlot::GetIDFromMapString(const std::map<std::string, int>& mapStr, const std::string& in) {
953 if (in.find_first_not_of("0123456789") == std::string::npos) {
954 return StringToInteger(in);
955 }
956 auto it = mapStr.find(in);
957 if (it != mapStr.end()) {
958 return it->second;
959 } else {
960 RESTWarning << "cannot find ID with name \"" << in << "\"" << RESTendl;
961 }
962 return -1;
963}
964
970
971 RESTMetadata << "DataSet name: " << fDataSetName << RESTendl;
972 RESTMetadata << "PaletteStyle: " << fPaletteStyle << RESTendl;
973 if (fPreviewPlot) RESTMetadata << "Preview plot is ACTIVE" << RESTendl;
974 RESTMetadata << "Canvas size: (" << fCanvasSize.X() << " ," << fCanvasSize.Y() << ")" << RESTendl;
975 RESTMetadata << "Canvas divisions: (" << fCanvasDivisions.X() << " ," << fCanvasDivisions.Y() << ")"
976 << RESTendl;
977 RESTMetadata << "-------------------" << RESTendl;
978 for (const auto& plot : fPlots) {
979 RESTMetadata << "-------------------" << RESTendl;
980 RESTMetadata << "Plot name/title: " << plot.name << " " << plot.title << RESTendl;
981 RESTMetadata << "Save string: " << plot.save << RESTendl;
982 RESTMetadata << "Set log X,Y,Z: " << plot.logX << ", " << plot.logY << ", " << plot.logZ << RESTendl;
983 RESTMetadata << "Stack draw Option: " << plot.stackDrawOption << RESTendl;
984 if (plot.legend) RESTMetadata << "Legend is ON" << RESTendl;
985 if (plot.timeDisplay) RESTMetadata << "Time display is ON" << RESTendl;
986 RESTMetadata << "Labels X,Y: " << plot.labelX << ", " << plot.labelY << RESTendl;
987 for (const auto& hist : plot.histos) {
988 RESTMetadata << "****************" << RESTendl;
989 RESTMetadata << "Histo name: " << hist.name << RESTendl;
990 RESTMetadata << "Draw Option: " << hist.drawOption << RESTendl;
991 RESTMetadata << "Histogram size: " << hist.variable.size() << " with parameters:" << RESTendl;
992 for (size_t i = 0; i < hist.variable.size(); i++) {
993 RESTMetadata << "\t" << i << " " << hist.variable[i] << ", " << hist.nBins[i] << ", "
994 << hist.range[i].X() << ", " << hist.range[i].Y() << RESTendl;
995 }
996 RESTMetadata << "****************" << RESTendl;
997 }
998 }
999 RESTMetadata << "-------------------" << RESTendl;
1000 for (auto& panel : fPanels) {
1001 RESTMetadata << "-------------------" << RESTendl;
1002 RESTMetadata << "Panel font size/precision " << panel.font_size << ", " << panel.precision
1003 << RESTendl;
1004 RESTMetadata << "****************" << RESTendl;
1005 for (auto& [key, posLabel] : panel.variablePos) {
1006 auto&& [obs, label, units] = key;
1007 RESTMetadata << "Label variable " << obs << ", label " << label << ", units " << units << " Pos ("
1008 << posLabel.X() << ", " << posLabel.Y() << ")" << RESTendl;
1009 }
1010 RESTMetadata << "****************" << RESTendl;
1011 for (auto& [key, posLabel] : panel.metadataPos) {
1012 auto&& [obs, label, units] = key;
1013 RESTMetadata << "Label metadata " << obs << ", label " << label << ", units " << units << " Pos ("
1014 << posLabel.X() << ", " << posLabel.Y() << ")" << RESTendl;
1015 }
1016 RESTMetadata << "****************" << RESTendl;
1017 for (auto& [key, posLabel] : panel.obsPos) {
1018 auto&& [obs, label, units] = key;
1019 RESTMetadata << "Label Observable " << obs << ", label " << label << ", units " << units
1020 << " Pos (" << posLabel.X() << ", " << posLabel.Y() << ")" << RESTendl;
1021 }
1022 RESTMetadata << "****************" << RESTendl;
1023 for (auto& [key, posLabel] : panel.expPos) {
1024 auto&& [obs, label, units] = key;
1025 RESTMetadata << "Label Expression " << obs << ", label " << label << ", units " << units
1026 << " Pos (" << posLabel.X() << ", " << posLabel.Y() << ")" << RESTendl;
1027 }
1028 RESTMetadata << "****************" << RESTendl;
1029 }
1030 RESTMetadata << "-------------------" << RESTendl;
1031
1032 RESTMetadata << RESTendl;
1033}
A class to help on cuts definitions. To be used with TRestAnalysisTree.
Definition: TRestCut.h:31
Perform the plot over datasets.
Int_t GetIDFromMapString(const std::map< std::string, int > &mapStr, const std::string &in)
This functions gets the ID from a map string that is passed by reference. It is used to translate col...
std::vector< PanelInfo > fPanels
Vector with panels/label options.
TRestCut * ReadCut(TRestCut *cut, TiXmlElement *ele=nullptr)
this function is used to add the different cuts provided in different metadata sections,...
void Initialize() override
Function to initialize input/output event members and define the section name.
void PlotCombinedCanvas()
This functions performs the plot of the combined canvas with the different panels and plots.
const std::map< std::string, int > LineStyleMap
LineStyleMap as enum "ELineStyle" defined in TAttLine.h.
std::string fDataSetName
Name of the dataset to be imported.
void ReadPlotInfo()
This function reads the config file plot info and stores it in a vector of PlotInfo.
void ReadPanelInfo()
This function reads the config file panel info and stores it in a vector of PanelInfo.
void GenerateDataSetFromFilePattern(TRestDataSet &dataSet)
This functions generates a dataSet based on the information of the rml file. A TRestDataSet is pased ...
Int_t fPaletteStyle
Palette style.
const std::map< std::string, int > FillStyleMap
FillStyleMap as enum "EFillStyle" defined in TAttFill.h.
TRestDataSetPlot()
Default constructor.
TRestCut * fCut
Global cut for the entire dataSet.
Bool_t fPreviewPlot
Preview plot.
const std::map< std::string, int > ColorIdMap
Maps for internal use only.
void PrintMetadata() override
Prints on screen the information about the metadata members of TRestDataSetPlot.
TVector2 fCanvasSize
Canvas options, size, divisions and margins.
std::string fOutputFileName
OutputFileName.
~TRestDataSetPlot()
Default destructor.
Double_t fLegendX1
Legend position and size.
std::vector< PlotInfo > fPlots
Vector with plots/pads options.
void InitFromConfigFile() override
Initialization of specific TRestDataSetPlot members through an RML file.
void CleanUp()
Clean up histos and text but note that the metadata is unchanged.
It allows to group a number of runs that satisfy given metadata conditions.
Definition: TRestDataSet.h:34
void Import(const std::string &fileName)
This function imports metadata from a root file it import metadata info from the previous dataSet whi...
Double_t GetTotalTimeInSeconds() const
It returns the accumulated run time in seconds.
Definition: TRestDataSet.h:166
ROOT::RDF::RNode MakeCut(const TRestCut *cut)
This function applies a TRestCut to the dataframe and returns a dataframe with the applied cuts....
void GenerateDataSet()
This function generates the data frame with the filelist and column names (or observables) that have ...
TTree * GetTree() const
Gives access to the tree.
Definition: TRestDataSet.h:138
A base class for any REST metadata class.
Definition: TRestMetadata.h:74
virtual void PrintMetadata()
Implemented it in the derived metadata class to print out specific metadata information.
endl_t RESTendl
Termination flag object for TRestStringOutput.
TiXmlElement * GetElement(std::string eleDeclare, TiXmlElement *e=nullptr)
Get an xml element from a given parent element, according to its declaration.
Int_t LoadConfigFromFile(const std::string &configFilename, const std::string &sectionName="")
Give the file name, find out the corresponding section. Then call the main starter.
TRestMetadata * InstantiateChildMetadata(int index, std::string pattern="")
This method will retrieve a new TRestMetadata instance of a child element of the present TRestMetadat...
virtual void InitFromConfigFile()
To make settings from rml file. This method must be implemented in the derived class.
TRestStringOutput::REST_Verbose_Level GetVerboseLevel()
returns the verboselevel in type of REST_Verbose_Level enumerator
void SetSectionName(std::string sName)
set the section name, clear the section content
std::string fConfigFileName
Full name of the rml file.
virtual Int_t Write(const char *name=nullptr, Int_t option=0, Int_t bufsize=0)
overwriting the write() method with fStore considered
TiXmlElement * GetNextElement(TiXmlElement *e)
Get the next sibling xml element of this element, with same eleDeclare.
std::string GetParameter(std::string parName, TiXmlElement *e, TString defaultValue=PARAMETER_NOT_FOUND_STR)
Returns the value for the parameter named parName in the given section.
@ REST_Info
+show most of the information for each steps
static std::string GetFileNameExtension(const std::string &fullname)
Gets the file extension as the substring found after the latest ".".
Definition: TRestTools.cxx:823
Int_t GetChar(std::string hint="Press a KEY to continue ...")
Helps to pause the program, printing a message before pausing.
Double_t StringToDouble(std::string in)
Gets a double from a string.
std::string ToUpper(std::string in)
Convert string to its upper case. Alternative of TString::ToUpper.
Int_t StringToInteger(std::string in)
Gets an integer from a string.
Int_t isAExpression(const std::string &in)
Returns 1 only if valid mathematical expression keywords (or numbers) are found in the string in....
std::string DoubleToString(Double_t d, std::string format="%8.6e")
Gets a string from a double.
TVector2 StringTo2DVector(std::string in)
Gets a 2D-vector from a string.
std::string EvaluateExpression(std::string exp)
Evaluates a complex numerical expression and returns the resulting value using TFormula.
Int_t isANumber(std::string in)
Returns 1 only if a valid number is found in the string in. If not it returns 0.
std::string RemoveWhiteSpaces(std::string in)
Returns the input string removing all white spaces.
std::string ToDateTimeString(time_t time)
Format time_t into string.
std::string Replace(std::string in, std::string thisString, std::string byThisString, size_t fromPosition=0, Int_t N=0)
Replace any occurences of thisSring by byThisString inside string in.
Nested classes for internal use only.
Auxiliary class for panels/labels.
Auxiliary struct for plots/pads.
std::string metadata
The associated metadata member used to register the relevant quantity.
Definition: TRestDataSet.h:38
std::string strategy
It determines how to produce the relevant quantity (accumulate/unique/last/max/min)
Definition: TRestDataSet.h:41