Logo ROOT  
Reference Guide
RFitPanel.cxx
Go to the documentation of this file.
1// Authors: Sergey Linev <S.Linev@gsi.de> Iliana Betsou <Iliana.Betsou@cern.ch>
2// Date: 2019-04-11
3// Warning: This is part of the ROOT 7 prototype! It will change without notice. It might trigger earthquakes. Feedback is welcome!
4
5/*************************************************************************
6 * Copyright (C) 1995-2019, Rene Brun and Fons Rademakers. *
7 * All rights reserved. *
8 * *
9 * For the licensing terms see $ROOTSYS/LICENSE. *
10 * For the list of contributors see $ROOTSYS/README/CREDITS. *
11 *************************************************************************/
12
13#include <ROOT/RFitPanel.hxx>
14
15#include <ROOT/RLogger.hxx>
16
17#include "Fit/BinData.h"
18#include "Fit/Fitter.h"
19
20#include "TString.h"
21#include "TGraph.h"
22#include "TGraphErrors.h"
23#include "TGraph2D.h"
24#include "TGraph2DErrors.h"
25#include "TMultiGraph.h"
26#include "THStack.h"
27#include "TROOT.h"
28#include "TH1.h"
29#include "TH2.h"
30#include "TF1.h"
31#include "TF2.h"
32#include "TF3.h"
33#include "TFitResult.h"
34#include "TList.h"
35#include "TVirtualPad.h"
36#include "TCanvas.h"
37#include "TDirectory.h"
38#include "TBufferJSON.h"
39#include "Math/Minimizer.h"
40#include "HFitInterface.h"
41#include "TColor.h"
42
43#include <iostream>
44#include <iomanip>
45#include <memory>
46#include <sstream>
47
48using namespace std::string_literals;
49
50using namespace ROOT::Experimental;
51
52/** \class RFitPanel
53\ingroup webwidgets
54
55web-based FitPanel prototype.
56*/
57
58////////////////////////////////////////////////////////////////////////////////
59/// Constructor
60
61RFitPanel::FitRes::FitRes(const std::string &_objid, std::unique_ptr<TF1> &_func, TFitResultPtr &_res)
62 : objid(_objid), res(_res)
63{
64 std::swap(func, _func);
65}
66
67////////////////////////////////////////////////////////////////////////////////
68/// Destructor
69
71{
72 // to avoid dependency from TF1
73
74 // if TF1 object deleted before - prevent second delete
75 if (func && func->IsZombie())
76 func.release();
77
78}
79
80////////////////////////////////////////////////////////////////////////////////
81/// Constructor
82
83RFitPanel::RFitPanel(const std::string &title)
84{
85 model().fTitle = title;
86
88}
89
90////////////////////////////////////////////////////////////////////////////////
91/// Destructor
92
94{
95 // to avoid dependency from TF1
96}
97
98////////////////////////////////////////////////////////////////////////////////
99/// Returns RWebWindow instance, used to display FitPanel
100
101std::shared_ptr<RWebWindow> RFitPanel::GetWindow()
102{
103 if (!fWindow) {
105
106 fWindow->SetPanelName("rootui5.fitpanel.view.FitPanel");
107
108 fWindow->SetCallBacks(
109 [this](unsigned connid)
110 {
111 fConnId = connid;
112 fWindow->Send(fConnId, "INITDONE");
113 if (!model().fInitialized)
114 SelectObject("$$$");
115 SendModel();
116 },
117 [this](unsigned connid, const std::string &arg) { ProcessData(connid, arg); },
118 [this](unsigned) { fConnId = 0; });
119
120 fWindow->SetGeometry(400, 650); // configure predefined geometry
121 }
122
123 return fWindow;
124}
125
126////////////////////////////////////////////////////////////////////////////////
127/// Update list of available data
128
130{
131 auto &m = model();
132
133 m.fDataSet.clear();
134
135 for (auto &obj : fObjects)
136 m.fDataSet.emplace_back("Panel", "panel::"s + obj->GetName(), Form("%s::%s", obj->ClassName(), obj->GetName()));
137
138 if (gDirectory) {
139 TIter iter(gDirectory->GetList());
140 TObject *obj = nullptr;
141 while ((obj = iter()) != nullptr) {
143 m.fDataSet.emplace_back("gDirectory", "gdir::"s + obj->GetName(), Form("%s::%s", obj->ClassName(), obj->GetName()));
144 }
145 }
146}
147
148////////////////////////////////////////////////////////////////////////////////
149/// Select object for fitting
150
151void RFitPanel::SelectObject(const std::string &objid)
152{
154
155 auto &m = model();
156
157 std::string id = objid;
158 if (id.compare("$$$") == 0) {
159 if (m.fDataSet.size() > 0)
160 id = m.fDataSet[0].id;
161 else
162 id.clear();
163 }
164
165 TObject *obj = GetSelectedObject(id);
166 auto kind = GetFitObjectType(obj);
167
168 m.SetObjectKind(kind);
169
170 TH1 *hist = nullptr;
171 switch (kind) {
173 hist = (TH1*)obj;
174 break;
175
177 hist = ((TGraph*)obj)->GetHistogram();
178 break;
179
181 hist = ((TMultiGraph*)obj)->GetHistogram();
182 break;
183
185 hist = ((TGraph2D*)obj)->GetHistogram("empty");
186 break;
187
189 hist = (TH1 *)((THStack *)obj)->GetHists()->First();
190 break;
191
192 //case RFitPanelModel::kObjectTree:
193 // m.fFitMethods = {{ kFP_MUBIN, "Unbinned Likelihood" }};
194 // m.fFitMethod = kFP_MUBIN;
195 // break;
196
197 default:
198 break;
199 }
200
201 if (!obj)
202 m.fSelectedData = "";
203 else
204 m.fSelectedData = id;
205
206 m.fInitialized = true;
207
208 // update list of data
209
210 m.UpdateRange(hist);
211
213
214 std::string selfunc = m.fSelectedFunc;
215
216 if (!m.HasFunction(selfunc)) {
217 if (m.fFuncList.size() > 0)
218 selfunc = m.fFuncList[0].id;
219 else
220 selfunc.clear();
221 }
222
223 SelectFunction(selfunc);
224}
225
226////////////////////////////////////////////////////////////////////////////////
227/// Returns object based on it string id
228/// Searches either in gDirectory or in internal panel list
229
230TObject *RFitPanel::GetSelectedObject(const std::string &objid)
231{
232 if (objid.compare(0,6,"gdir::") == 0) {
233 std::string name = objid.substr(6);
234 if (gDirectory)
235 return gDirectory->GetList()->FindObject(name.c_str());
236 } else if (objid.compare(0,7,"panel::") == 0) {
237 std::string name = objid.substr(7);
238 for (auto &item : fObjects)
239 if (name.compare(item->GetName()) == 0)
240 return item;
241 }
242
243 return nullptr;
244}
245
246////////////////////////////////////////////////////////////////////////////////
247/// Returns kind of object
248
250{
251 if (!obj)
253
254 if (obj->InheritsFrom(TH1::Class()))
256
257 if (obj->InheritsFrom(TGraph::Class()))
259
260 if (obj->InheritsFrom(TGraph2D::Class()))
262
263 if (obj->InheritsFrom(THStack::Class()))
265
268
270}
271
272////////////////////////////////////////////////////////////////////////////////
273/// Update list of available functions
274
276{
277 auto &m = model();
278
279 m.fFuncList.clear();
280
281 if (m.fDim == 1) {
282 m.fFuncList = { {"gaus"}, {"gausn"}, {"expo"}, {"landau"}, {"landaun"},
283 {"pol0"},{"pol1"},{"pol2"},{"pol3"},{"pol4"},{"pol5"},{"pol6"},{"pol7"},{"pol8"},{"pol9"},
284 {"cheb0"}, {"cheb1"}, {"cheb2"}, {"cheb3"}, {"cheb4"}, {"cheb5"}, {"cheb6"}, {"cheb7"}, {"cheb8"}, {"cheb9"} };
285 } else if (m.fDim == 2) {
286 m.fFuncList = { {"xygaus"}, {"bigaus"}, {"xyexpo"}, {"xylandau"}, {"xylandaun"} };
287 }
288
289 for (auto &func : fSystemFuncs) {
290 m.fFuncList.emplace_back("System", "system::"s + func->GetName(), func->GetName());
291 }
292
293 for (auto &entry : fPrevRes) {
294 if (entry.objid == m.fSelectedData)
295 m.fFuncList.emplace_back("Previous", "previous::"s + entry.func->GetName(), entry.func->GetName());
296 }
297}
298
299////////////////////////////////////////////////////////////////////////////////
300/// Select fit function
301
302void RFitPanel::SelectFunction(const std::string &funcid)
303{
304 model().SelectedFunc(funcid, FindFunction(funcid));
305
307}
308
309////////////////////////////////////////////////////////////////////////////////
310/// Assign histogram to use with fit panel - without ownership
311
313{
314 fObjects.emplace_back(hist);
315 SelectObject("panel::"s + hist->GetName());
316 SendModel();
317}
318
319/// Assign histogram name to use with fit panel - it should be available in gDirectory
320
321void RFitPanel::AssignHistogram(const std::string &hname)
322{
323 SelectObject("gdir::" + hname);
324 SendModel();
325}
326
327////////////////////////////////////////////////////////////////////////////////
328/// assign canvas to use for drawing results of fitting or showing fitpanel itself
329
330void RFitPanel::AssignCanvas(std::shared_ptr<RCanvas> &canv)
331{
332 if (!fCanvas) {
333 fCanvas = canv;
334 } else {
335 R__LOG_ERROR(FitPanelLog()) << "FitPanel already bound to the canvas - change is not yet supported";
336 }
337}
338
339////////////////////////////////////////////////////////////////////////////////
340/// assign histogram for fitting
341
342void RFitPanel::AssignHistogram(std::shared_ptr<RH1D> &hist)
343{
344 fFitHist = hist;
345}
346
347////////////////////////////////////////////////////////////////////////////////
348/// Show FitPanel
349
350void RFitPanel::Show(const std::string &where)
351{
352 GetWindow()->Show(where);
353}
354
355////////////////////////////////////////////////////////////////////////////////
356/// Hide FitPanel
357
359{
360 if (fWindow)
361 fWindow->CloseConnections();
362}
363
364////////////////////////////////////////////////////////////////////////////////
365/// Return reference on model object
366/// Model created if was not exists before
367
368
370{
371 if (!fModel) {
372 fModel = std::make_unique<RFitPanelModel>();
373 fModel->Initialize();
374 }
375
376 return *fModel.get();
377}
378
379////////////////////////////////////////////////////////////////////////////////
380/// Send model object to the client
381
383{
384 if (fWindow && (fConnId > 0)) {
386 fWindow->Send(fConnId, "MODEL:"s + json.Data());
387 }
388}
389
390////////////////////////////////////////////////////////////////////////////////
391/// Process data from FitPanel
392/// OpenUI5-based FitPanel sends commands or status changes
393
394void RFitPanel::ProcessData(unsigned, const std::string &arg)
395{
396 if (arg == "RELOAD") {
397
399
402
403 SendModel();
404
405 } else if (arg.compare(0, 7, "UPDATE:") == 0) {
406
407 if (UpdateModel(arg.substr(7)) > 0)
408 SendModel();
409
410 } else if (arg.compare(0, 6, "DOFIT:") == 0) {
411
412 if (UpdateModel(arg.substr(6)) >= 0)
413 if (DoFit())
414 SendModel();
415
416 } else if (arg.compare(0, 7, "DODRAW:") == 0) {
417
418 if (UpdateModel(arg.substr(7)) >= 0)
419 if (DoDraw())
420 SendModel();
421
422 } else if (arg.compare(0, 8, "SETPARS:") == 0) {
423
424 auto info = TBufferJSON::FromJSON<RFitPanelModel::RFuncParsList>(arg.substr(8));
425
426 if (info) {
427 TF1 *func = FindFunction(info->id);
428
429 // copy all parameters back to the function
430 if (func)
431 info->SetParameters(func);
432 }
433
434 }
435}
436
437////////////////////////////////////////////////////////////////////////////////
438/// Search for existing functions, ownership still belongs to FitPanel or global lists
439
440TF1 *RFitPanel::FindFunction(const std::string &id)
441{
442 if (id.compare(0,8,"system::") == 0) {
443 std::string name = id.substr(8);
444
445 for (auto &item : fSystemFuncs)
446 if (name.compare(item->GetName()) == 0)
447 return item.get();
448 }
449
450 if (id.compare(0,10,"previous::") == 0) {
451 std::string name = id.substr(10);
452
453 for (auto &entry : fPrevRes)
454 if (name.compare(entry.func->GetName()) == 0)
455 return entry.func.get();
456 }
457
458 return nullptr;
459}
460
461////////////////////////////////////////////////////////////////////////////////
462/// Creates new instance to make fitting
463
465{
466 if (id.compare(0,10,"previous::") == 0) {
467 std::string name = id.substr(10);
468
469 for (auto &entry : fPrevRes)
470 if (name.compare(entry.func->GetName()) == 0)
471 return entry.res.Get();
472 }
473
474 return nullptr;
475}
476
477////////////////////////////////////////////////////////////////////////////////
478/// Creates new instance to make fitting
479
480std::unique_ptr<TF1> RFitPanel::GetFitFunction(const std::string &funcname)
481{
482 std::unique_ptr<TF1> res;
483
484 TF1 *func = FindFunction(funcname);
485
486 if (func) {
487 // Now we make a copy.
488 res.reset((TF1*) func->IsA()->New());
489 func->Copy(*res);
490 } else if (funcname.compare(0,6,"dflt::") == 0) {
491
492 std::string formula = funcname.substr(6);
493
495
496 double xmin, xmax, ymin, ymax, zmin, zmax;
497 drange.GetRange(xmin, xmax, ymin, ymax, zmin, zmax);
498
499 if ( model().fDim == 1 || model().fDim == 0 ) {
500 res.reset(new TF1(formula.c_str(), formula.c_str(), xmin, xmax));
501 } else if ( model().fDim == 2 ) {
502 res.reset(new TF2(formula.c_str(), formula.c_str(), xmin, xmax, ymin, ymax));
503 } else if ( model().fDim == 3 ) {
504 res.reset(new TF3(formula.c_str(), formula.c_str(), xmin, xmax, ymin, ymax, zmin, zmax));
505 }
506 }
507
508 return res;
509}
510
511////////////////////////////////////////////////////////////////////////////////
512/// Update fit model
513/// returns -1 if JSON fails
514/// return 0 if nothing large changed
515/// return 1 if important selection are changed and client need to be updated
516
517int RFitPanel::UpdateModel(const std::string &json)
518{
519 auto m = TBufferJSON::FromJSON<RFitPanelModel>(json);
520
521 if (!m) {
522 R__LOG_ERROR(FitPanelLog()) << "Fail to parse JSON for RFitPanelModel";
523 return -1;
524 }
525
526 m->fInitialized = true;
527
528 int res = 0; // nothing changed
529
530 if (model().fSelectedData != m->fSelectedData) {
531 res |= 1;
532 }
533
534 if (model().fSelectedFunc != m->fSelectedFunc) {
535 res |= 2;
536 }
537
538 std::swap(fModel, m); // finally replace model
539
540 if (res & 1)
541 SelectObject(model().fSelectedData);
542
543 if (res != 0)
544 SelectFunction(model().fSelectedFunc);
545
546 return res;
547}
548
549////////////////////////////////////////////////////////////////////////////////
550/// Copies f into a new TF1 to be stored in the fitpanel with it's
551/// own ownership. This is taken from Fit::StoreAndDrawFitFunction in
552/// HFitImpl.cxx
553
555{
556 double xmin = 0, xmax = 0, ymin = 0, ymax = 0, zmin = 0, zmax = 0;
557
558 // no need to use kNotGlobal bit. TF1::Copy does not add in the list by default
559 if ( dynamic_cast<TF3*>(f)) {
560 TF3* fnew = (TF3*)f->IsA()->New();
561 f->Copy(*fnew);
562 f->GetRange(xmin,ymin,zmin,xmax,ymax,zmax);
563 fnew->SetRange(xmin,ymin,zmin,xmax,ymax,zmax);
564 fnew->SetParent( nullptr );
565 fnew->AddToGlobalList(false);
566 return fnew;
567 } else if ( dynamic_cast<TF2*>(f) != 0 ) {
568 TF2* fnew = (TF2*)f->IsA()->New();
569 f->Copy(*fnew);
570 f->GetRange(xmin,ymin,xmax,ymax);
571 fnew->SetRange(xmin,ymin,xmax,ymax);
572 fnew->Save(xmin,xmax,ymin,ymax,0,0);
573 fnew->SetParent( nullptr );
574 fnew->AddToGlobalList(false);
575 return fnew;
576 }
577
578 TF1* fnew = (TF1*)f->IsA()->New();
579 f->Copy(*fnew);
580 f->GetRange(xmin,xmax);
581 fnew->SetRange(xmin,xmax);
582 // This next line is added, as fnew-Save fails with gausND! As
583 // the number of dimensions is unknown...
584 if ( '\0' != fnew->GetExpFormula()[0] )
585 fnew->Save(xmin,xmax,0,0,0,0);
586 fnew->SetParent( nullptr );
587 fnew->AddToGlobalList(false);
588 return fnew;
589}
590
591////////////////////////////////////////////////////////////////////////////////
592/// Looks for all the functions registered in the current ROOT
593/// session.
594
596{
597
598 fSystemFuncs.clear();
599
600 // Be carefull not to store functions that will be in the
601 // predefined section
602 std::vector<std::string> fnames = { "gaus" , "gausn", "expo", "landau",
603 "landaun", "pol0", "pol1", "pol2",
604 "pol3", "pol4", "pol5", "pol6",
605 "pol7", "pol8", "pol9", "user" };
606
607 // No go through all the objects registered in gROOT
608 TIter functionsIter(gROOT->GetListOfFunctions());
609 TObject* obj;
610 while( (obj = functionsIter()) != nullptr ) {
611 // And if they are TF1s
612 if ( TF1* func = dynamic_cast<TF1*>(obj) ) {
613 bool addFunction = true;
614 // And they are not already registered in fSystemFunc
615 for ( auto &name : fnames) {
616 if ( name.compare(func->GetName()) == 0 ) {
617 addFunction = false;
618 break;
619 }
620 }
621 // Add them.
622 if ( addFunction )
623 fSystemFuncs.emplace_back( copyTF1(func) );
624 }
625 }
626}
627
628////////////////////////////////////////////////////////////////////////////////
629/// Returns pad where histogram is drawn
630/// If canvas not exists, create new one
631
633{
634 if (!obj || (!force && (model().fNoDrawing || model().fNoStoreDraw)))
635 return nullptr;
636
637 std::function<TPad*(TPad*)> check = [&](TPad *pad) {
638 TPad *res = nullptr;
639 if (!pad) return res;
640 if (!fPadName.empty() && (fPadName.compare(pad->GetName()) == 0)) return pad;
641 TIter next(pad->GetListOfPrimitives());
642 TObject *prim = nullptr;
643 while (!res && (prim = next())) {
644 res = (prim == obj) ? pad : check(dynamic_cast<TPad *>(prim));
645 }
646 return res;
647 };
648
649 if (!fCanvName.empty()) {
650 auto drawcanv = dynamic_cast<TCanvas *> (gROOT->GetListOfCanvases()->FindObject(fCanvName.c_str()));
651 auto drawpad = check(drawcanv);
652 if (drawpad) {
653 drawpad->cd();
654 return drawpad;
655 }
656 if (drawcanv) {
657 drawcanv->Clear();
658 drawcanv->cd();
659 obj->Draw();
660 return drawcanv;
661 }
662 fCanvName.clear();
663 fPadName.clear();
664 }
665
666 TObject *c = nullptr;
667 TIter nextc(gROOT->GetListOfCanvases());
668 while ((c = nextc())) {
669 auto drawpad = check(dynamic_cast<TCanvas* >(c));
670 if (drawpad) {
671 drawpad->cd();
672 fCanvName = c->GetName();
673 fPadName = drawpad->GetName();
674 return drawpad;
675 }
676 }
677
678 auto canv = gROOT->MakeDefCanvas();
679 canv->SetName("fpc");
680 canv->SetTitle("Fit panel drawings");
681 fPadName = fCanvName = canv->GetName();
682
683 canv->cd();
684 obj->Draw();
685
686 return canv;
687}
688
689////////////////////////////////////////////////////////////////////////////////
690/// Perform fitting using current model settings
691/// Returns true if any action was done
692
694{
695 auto &m = model();
696
697 TObject *obj = GetSelectedObject(m.fSelectedData);
698 if (!obj) return false;
699 auto kind = GetFitObjectType(obj);
700
701 auto f1 = GetFitFunction(m.fSelectedFunc);
702 if (!f1) return false;
703
704 auto drange = m.GetRanges();
705 auto minOption = m.GetMinimizerOptions();
706 auto fitOpts = m.GetFitOptions();
707 auto drawOpts = m.GetDrawOption();
708
709 fitOpts.StoreResult = 1;
710
712
713 auto pad = GetDrawPad(obj);
714
715 TFitResultPtr res;
716
717 switch (kind) {
719
720 TH1 *hist = dynamic_cast<TH1*>(obj);
721 if (hist)
722 res = ROOT::Fit::FitObject(hist, f1.get(), fitOpts, minOption, drawOpts, drange);
723
724 break;
725 }
727
728 TGraph *gr = dynamic_cast<TGraph*>(obj);
729 if (gr)
730 res = ROOT::Fit::FitObject(gr, f1.get(), fitOpts, minOption, drawOpts, drange);
731 break;
732 }
734
735 TMultiGraph *mg = dynamic_cast<TMultiGraph*>(obj);
736 if (mg)
737 res = ROOT::Fit::FitObject(mg, f1.get(), fitOpts, minOption, drawOpts, drange);
738
739 break;
740 }
742
743 TGraph2D *g2d = dynamic_cast<TGraph2D*>(obj);
744 if (g2d)
745 res = ROOT::Fit::FitObject(g2d, f1.get(), fitOpts, minOption, drawOpts, drange);
746
747 break;
748 }
750 // N/A
751 break;
752 }
753
754 default: {
755 // N/A
756 break;
757 }
758 }
759
760 // After fitting function appears in global list
761 if (f1 && gROOT->GetListOfFunctions()->FindObject(f1.get()))
762 gROOT->GetListOfFunctions()->Remove(f1.get());
763
764 if (m.fSame && f1 && pad) {
765 TF1 *copy = copyTF1(f1.get());
766 copy->SetBit(kCanDelete);
767 copy->Draw("same");
768 }
769
770 if (pad) {
771 pad->Modified();
772 pad->Update();
773 }
774
775 std::string funcname = f1->GetName();
776 if ((funcname.compare(0,4,"prev") == 0) && (funcname.find("-") > 4))
777 funcname.erase(0, funcname.find("-") + 1);
778 funcname = "prev"s + std::to_string(fPrevRes.size() + 1) + "-"s + funcname;
779 f1->SetName(funcname.c_str());
780
781 fPrevRes.emplace_back(m.fSelectedData, f1, res);
782
784
785 SelectFunction("previous::"s + funcname);
786
787 return true; // provide client with latest changes
788}
789
790////////////////////////////////////////////////////////////////////////////////
791/// Extract color from string
792/// Should be coded as #%ff00ff string
793
794Color_t RFitPanel::GetColor(const std::string &colorid)
795{
796 if ((colorid.length() != 7) || (colorid.compare(0,1,"#") != 0)) return 0;
797
798 return TColor::GetColor(colorid.c_str());
799}
800
801////////////////////////////////////////////////////////////////////////////////
802/// Create confidence levels drawing
803/// tab. Then it call Virtual Fitter to perform it.
804
806{
807 if (!result)
808 return nullptr;
809
810 // try to use provided method
811 // auto conf = result->GetConfidenceIntervals();
812 // printf("GET INTERVALS %d\n", (int) conf.size());
813
814 const auto *function = result->FittedFunction();
815 if (!function) {
816 R__LOG_ERROR(FitPanelLog()) << "Fit Function does not exist!";
817 return nullptr;
818 }
819
820 const auto *data = result->FittedBinData();
821 if (!data) {
822 R__LOG_ERROR(FitPanelLog()) << "Unbinned data set cannot draw confidence levels.";
823 return nullptr;
824 }
825
826 std::vector<Double_t> ci(data->Size());
827 result->GetConfidenceIntervals(*data, &ci[0], model().fConfidenceLevel);
828
829 if (model().fDim == 1) {
830 TGraphErrors *g = new TGraphErrors(ci.size());
831 for (unsigned int i = 0; i < ci.size(); ++i) {
832 const Double_t *x = data->Coords(i);
833 const Double_t y = (*function)(x);
834 g->SetPoint(i, *x, y);
835 g->SetPointError(i, 0, ci[i]);
836 }
837 // std::ostringstream os;
838 // os << "Confidence Intervals with " << conflvl << " conf. band.";
839 // g->SetTitle(os.str().c_str());
840 g->SetTitle("Confidence Intervals with");
841
842 auto icol = GetColor(model().fConfidenceColor);
843 g->SetLineColor(icol);
844 g->SetFillColor(icol);
845 g->SetFillStyle(3001);
846 return g;
847 } else if (model().fDim == 2) {
848 TGraph2DErrors *g = new TGraph2DErrors(ci.size());
849 for (unsigned int i = 0; i < ci.size(); ++i) {
850 const Double_t *x = data->Coords(i);
851 const Double_t y = (*function)(x);
852 g->SetPoint(i, x[0], x[1], y);
853 g->SetPointError(i, 0, 0, ci[i]);
854 }
855 // std::ostringstream os;
856 // os << "Confidence Intervals with " << fConfLevel->GetNumber() << " conf. band.";
857 // g->SetTitle(os.str().c_str());
858
859 g->SetTitle("Confidence Intervals with");
860
861 auto icol = GetColor(model().fConfidenceColor);
862 g->SetLineColor(icol);
863 g->SetFillColor(icol);
864 g->SetFillStyle(3001);
865 return g;
866 }
867
868 return nullptr;
869}
870
871////////////////////////////////////////////////////////////////////////////////
872/// Perform drawing using current model settings
873/// Returns true if any action was done
874
876{
877 auto &m = model();
878
879 TObject *obj = GetSelectedObject(m.fSelectedData);
880 if (!obj) return false;
881
882 TObject *drawobj = nullptr;
883 std::string drawopt;
884 bool superimpose = true, objowner = true;
885
886 if (m.fHasAdvanced && (m.fSelectedTab == "Advanced")) {
887
888 TFitResult *res = FindFitResult(m.fSelectedFunc);
889 if (!res) return false;
890
891 if (m.fAdvancedTab == "Contour") {
892
893 superimpose = m.fContourSuperImpose;
894 int par1 = std::stoi(m.fContourPar1Id);
895 int par2 = std::stoi(m.fContourPar2Id);
896
897 TGraph *graph = new TGraph(m.fContourPoints);
898
899 if (!res->Contour(par1, par2, graph, m.fConfidenceLevel)) {
900 delete graph;
901 return false;
902 }
903
904 auto fillcolor = GetColor(m.fContourColor);
905 graph->SetFillColor(fillcolor);
906 graph->GetXaxis()->SetTitle( res->ParName(par1).c_str() );
907 graph->GetYaxis()->SetTitle( res->ParName(par2).c_str() );
908
909 drawobj = graph;
910 drawopt = superimpose ? "LF" : "ALF";
911
912 } else if (m.fAdvancedTab == "Scan") {
913
914 int par = std::stoi(m.fScanId);
915 TGraph *graph = new TGraph( m.fScanPoints);
916 if (!res->Scan( par, graph, m.fScanMin, m.fScanMax)) {
917 delete graph;
918 return false;
919 }
920 auto linecolor = GetColor(m.fScanColor);
921 if (!linecolor) linecolor = kBlue;
922 graph->SetLineColor(linecolor);
923 graph->SetLineWidth(2);
924 graph->GetXaxis()->SetTitle(res->ParName(par).c_str());
925 graph->GetYaxis()->SetTitle("FCN" );
926
927 superimpose = false;
928 drawobj = graph;
929 drawopt = "ALF";
930
931 } else if (m.fAdvancedTab == "Confidence") {
932
933 drawobj = MakeConfidenceLevels(res);
934 drawopt = "C3same";
935
936 } else {
937 return false;
938 }
939 } else {
940
941 // find already existing functions, not try to create something new
942 TF1 *func = FindFunction(m.fSelectedFunc);
943
944 // when "Pars" tab is selected, automatically update function parameters
945 if (func && (m.fSelectedTab.compare("Pars") == 0) && (m.fSelectedFunc == m.fFuncPars.id))
946 m.fFuncPars.SetParameters(func);
947
948 drawobj = func;
949 drawopt = "same";
950 objowner = true;
951 }
952
953 if (!drawobj)
954 return false;
955
956 auto pad = GetDrawPad(obj, true);
957 if (!pad) {
958 if (objowner)
959 delete drawobj;
960 return false;
961 }
962
963 if (!superimpose)
964 pad->Clear();
965
966 if (objowner)
967 drawobj->SetBit(kCanDelete);
968
969 drawobj->Draw(drawopt.c_str());
970
971 pad->Modified();
972 pad->Update();
973
974 return true;
975}
nlohmann::json json
#define R__LOG_ERROR(...)
Definition: RLogger.hxx:362
#define f(i)
Definition: RSha256.hxx:104
#define c(i)
Definition: RSha256.hxx:101
short Color_t
Definition: RtypesCore.h:92
const Bool_t kFALSE
Definition: RtypesCore.h:101
@ kBlue
Definition: Rtypes.h:66
#define gDirectory
Definition: TDirectory.h:385
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
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 result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize id
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 g
char name[80]
Definition: TGX11.cxx:110
float xmin
Definition: THbookFile.cxx:95
float ymin
Definition: THbookFile.cxx:95
float xmax
Definition: THbookFile.cxx:95
float ymax
Definition: THbookFile.cxx:95
@ kCanDelete
Definition: TObject.h:369
#define gROOT
Definition: TROOT.h:405
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2468
RFitPanel(const std::string &title="Fit panel")
Constructor.
Definition: RFitPanel.cxx:83
TPad * GetDrawPad(TObject *obj, bool force=false)
Returns pad where histogram is drawn If canvas not exists, create new one.
Definition: RFitPanel.cxx:632
void UpdateDataSet()
Update list of available data.
Definition: RFitPanel.cxx:129
void UpdateFunctionsList()
Update list of available functions.
Definition: RFitPanel.cxx:275
void SelectFunction(const std::string &funcid)
Select fit function.
Definition: RFitPanel.cxx:302
void GetFunctionsFromSystem()
Looks for all the functions registered in the current ROOT session.
Definition: RFitPanel.cxx:595
RFitPanelModel::EFitObjectType GetFitObjectType(TObject *obj)
Returns kind of object.
Definition: RFitPanel.cxx:249
std::shared_ptr< RCanvas > fCanvas
! v7 canvas used to display results
Definition: RFitPanel.hxx:45
TF1 * FindFunction(const std::string &funcid)
Search for existing functions, ownership still belongs to FitPanel or global lists.
Definition: RFitPanel.cxx:440
std::shared_ptr< RWebWindow > fWindow
! configured display
Definition: RFitPanel.hxx:48
std::string fCanvName
! v6 canvas name used to display fit, will be created if not exists
Definition: RFitPanel.hxx:42
std::unique_ptr< RFitPanelModel > fModel
Definition: RFitPanel.hxx:39
void Show(const std::string &where="")
show FitPanel in specified place
Definition: RFitPanel.cxx:350
std::shared_ptr< RWebWindow > GetWindow()
Returns RWebWindow instance, used to display FitPanel.
Definition: RFitPanel.cxx:101
Color_t GetColor(const std::string &colorid)
Extract color from string Should be coded as #ff00ff string.
Definition: RFitPanel.cxx:794
unsigned fConnId
! client connection id
Definition: RFitPanel.hxx:49
std::string fPadName
! v6 pad name in the canvas, where object is (was) drawn
Definition: RFitPanel.hxx:43
bool DoFit()
Perform fitting using current model settings Returns true if any action was done.
Definition: RFitPanel.cxx:693
TObject * MakeConfidenceLevels(TFitResult *res)
Create confidence levels drawing tab.
Definition: RFitPanel.cxx:805
std::vector< TObject * > fObjects
! objects provided directly to panel for fitting
Definition: RFitPanel.hxx:41
std::unique_ptr< TF1 > GetFitFunction(const std::string &funcid)
Creates new instance to make fitting.
Definition: RFitPanel.cxx:480
void AssignHistogram(TH1 *hist)
Assign histogram to use with fit panel - without ownership.
Definition: RFitPanel.cxx:312
TFitResult * FindFitResult(const std::string &funcid)
Creates new instance to make fitting.
Definition: RFitPanel.cxx:464
TObject * GetSelectedObject(const std::string &objid)
Returns object based on it string id Searches either in gDirectory or in internal panel list.
Definition: RFitPanel.cxx:230
void SendModel()
Send model object to the client.
Definition: RFitPanel.cxx:382
void Hide()
hide FitPanel
Definition: RFitPanel.cxx:358
int UpdateModel(const std::string &json)
Update fit model returns -1 if JSON fails return 0 if nothing large changed return 1 if important sel...
Definition: RFitPanel.cxx:517
std::shared_ptr< RH1D > fFitHist
! v7 histogram for fitting
Definition: RFitPanel.hxx:46
RFitPanelModel & model()
Return reference on model object Model created if was not exists before.
Definition: RFitPanel.cxx:369
TF1 * copyTF1(TF1 *f)
Copies f into a new TF1 to be stored in the fitpanel with it's own ownership.
Definition: RFitPanel.cxx:554
std::vector< std::unique_ptr< TF1 > > fSystemFuncs
! local copy of all internal system funcs
Definition: RFitPanel.hxx:51
std::list< FitRes > fPrevRes
! all previous functions used for fitting
Definition: RFitPanel.hxx:62
void ProcessData(unsigned connid, const std::string &arg)
Process data from FitPanel OpenUI5-based FitPanel sends commands or status changes.
Definition: RFitPanel.cxx:394
void AssignCanvas(const std::string &cname)
Definition: RFitPanel.hxx:107
bool DoDraw()
Perform drawing using current model settings Returns true if any action was done.
Definition: RFitPanel.cxx:875
void SelectObject(const std::string &objid)
Select object for fitting.
Definition: RFitPanel.cxx:151
static std::shared_ptr< RWebWindow > Create()
Create new RWebWindow Using default RWebWindowsManager.
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition: DataRange.h:35
void GetRange(unsigned int irange, unsigned int icoord, double &xmin, double &xmax) const
get the i-th range for given coordinate.
Definition: DataRange.h:104
std::string ParName(unsigned int i) const
name of the parameter
Definition: FitResult.cxx:371
static TString ToJSON(const T *obj, Int_t compact=0, const char *member_name=nullptr)
Definition: TBufferJSON.h:75
The Canvas class.
Definition: TCanvas.h:23
void * New(ENewType defConstructor=kClassNew, Bool_t quiet=kFALSE) const
Return a pointer to a newly allocated object of this class.
Definition: TClass.cxx:4967
static Int_t GetColor(const char *hexcolor)
Static method returning color number for color specified by hex color string of form: "#rrggbb",...
Definition: TColor.cxx:1823
1-Dim function class
Definition: TF1.h:213
virtual void SetRange(Double_t xmin, Double_t xmax)
Initialize the upper and lower bounds to draw the function.
Definition: TF1.cxx:3520
virtual TString GetExpFormula(Option_t *option="") const
Definition: TF1.h:466
virtual void SetParent(TObject *p=nullptr)
Definition: TF1.h:670
void Copy(TObject &f1) const override
Copy this F1 to a new F1.
Definition: TF1.cxx:1006
void Draw(Option_t *option="") override
Draw this function with its current attributes.
Definition: TF1.cxx:1334
virtual void Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax)
Save values of function in array fSave.
Definition: TF1.cxx:3163
TClass * IsA() const override
Definition: TF1.h:719
virtual Bool_t AddToGlobalList(Bool_t on=kTRUE)
Add to global list of functions (gROOT->GetListOfFunctions() ) return previous status (true if the fu...
Definition: TF1.cxx:848
A 2-Dim function with parameters.
Definition: TF2.h:29
void Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax) override
Save values of function in array fSave.
Definition: TF2.cxx:798
void SetRange(Double_t xmin, Double_t xmax) override
Initialize the upper and lower bounds to draw the function.
Definition: TF2.h:148
A 3-Dim function with parameters.
Definition: TF3.h:28
void SetRange(Double_t xmin, Double_t xmax) override
Initialize the upper and lower bounds to draw the function.
Definition: TF3.h:145
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
Definition: TFitResultPtr.h:32
Extends the ROOT::Fit::Result class with a TNamed inheritance providing easy possibility for I/O.
Definition: TFitResult.h:34
bool Contour(unsigned int ipar, unsigned int jpar, TGraph *gr, double confLevel=0.683)
Create a 2D contour around the minimum for the parameter ipar and jpar if a minimum does not exist or...
Definition: TFitResult.cxx:119
bool Scan(unsigned int ipar, TGraph *gr, double xmin=0, double xmax=0)
Scan parameter ipar between value of xmin and xmax A graph must be given which will be on return fill...
Definition: TFitResult.cxx:93
Graph 2D class with errors.
Graphics object made of three arrays X, Y and Z with the same number of points each.
Definition: TGraph2D.h:41
static TClass * Class()
A TGraphErrors is a TGraph with error bars.
Definition: TGraphErrors.h:26
A TGraph is an object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
static TClass * Class()
TH1 is the base class of all histogram classes in ROOT.
Definition: TH1.h:58
static TClass * Class()
The Histogram stack class.
Definition: THStack.h:38
static TClass * Class()
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition: TMultiGraph.h:34
static TClass * Class()
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:140
Mother of all ROOT objects.
Definition: TObject.h:41
virtual void Clear(Option_t *="")
Definition: TObject.h:119
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:439
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:207
virtual Option_t * GetDrawOption() const
Get option used by the graphics system to draw this object.
Definition: TObject.cxx:423
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:774
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:525
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition: TObject.cxx:274
The most important graphics class in the ROOT system.
Definition: TPad.h:28
TList * GetListOfPrimitives() const override
Definition: TPad.h:243
const char * GetName() const override
Returns name of object.
Definition: TPad.h:258
Basic string class.
Definition: TString.h:136
small helper class to store/restore gPad context in TPad methods
Definition: TVirtualPad.h:61
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
TGraphErrors * gr
Definition: legend1.C:25
TF1 * f1
Definition: legend1.C:11
void swap(RDirectoryEntry &e1, RDirectoryEntry &e2) noexcept
RLogChannel & FitPanelLog()
Log channel for FitPanel diagnostics.
TFitResultPtr FitObject(TH1 *h1, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption, const char *goption, ROOT::Fit::DataRange &range)
fitting function for a TH1 (called from TH1::Fit)
Definition: HFitImpl.cxx:965
void function(const Char_t *name_, T fun, const Char_t *docstring=0)
Definition: RExports.h:167
static constexpr double s
static constexpr double mg
Definition: graph.py:1
Data structure for the fit panel.
void UpdateAdvanced(TFitResult *res)
Update advanced parameters associated with fit function.
void SelectedFunc(const std::string &name, TF1 *func)
Select function.
std::string fTitle
title of the fit panel
TMarker m
Definition: textangle.C:8