GSI Object Oriented Online Offline (Go4)  GO4-6.2.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
TGo4Fitter.cxx
Go to the documentation of this file.
1 // $Id: TGo4Fitter.cxx 3494 2022-01-21 15:54:00Z linev $
2 //-----------------------------------------------------------------------
3 // The GSI Online Offline Object Oriented (Go4) Project
4 // Experiment Data Processing at EE department, GSI
5 //-----------------------------------------------------------------------
6 // Copyright (C) 2000- GSI Helmholtzzentrum fuer Schwerionenforschung GmbH
7 // Planckstr. 1, 64291 Darmstadt, Germany
8 // Contact: http://go4.gsi.de
9 //-----------------------------------------------------------------------
10 // This software can be used under the license agreements as stated
11 // in Go4License.txt file which is part of the distribution.
12 //-----------------------------------------------------------------------
13 
14 #include "TGo4Fitter.h"
15 
16 #include <cstdlib>
17 #include <iostream>
18 #include <iomanip>
19 
20 #include "TClass.h"
21 #include "TMath.h"
22 #include "TH1.h"
23 #include "TGraph.h"
24 #include "TCanvas.h"
25 #include "TObjArray.h"
26 #include "TObjString.h"
27 #include "TROOT.h"
28 #include "THStack.h"
29 #include "TArrayD.h"
30 #include "TBuffer.h"
31 
32 #include "TGo4FitDataHistogram.h"
33 #include "TGo4FitDataGraph.h"
34 #include "TGo4FitModelGauss1.h"
35 #include "TGo4FitModelPolynom.h"
36 #include "TGo4FitAmplEstimation.h"
37 #include "TGo4FitSlot.h"
38 
41  fxDatas(), fxModels(), fiFitFunctionType(0), fiMemoryUsage(0),
42  fxUserFitFunction(0), fxDrawObjs(0)
43 {
44 }
45 
46 TGo4Fitter::TGo4Fitter(const char* iName, const char* iTitle) :
47  TGo4FitterAbstract(iName, iTitle),
48  fxDatas(), fxModels(), fiFitFunctionType(0), fiMemoryUsage(0), fxUserFitFunction(0), fxDrawObjs(0)
49 {
50  fxDatas.SetOwner(kTRUE);
51  fxModels.SetOwner(kTRUE);
52 }
53 
54 TGo4Fitter::TGo4Fitter(const char* iName, Int_t iFitFunctionType, Bool_t IsAddStandardActions) :
55  TGo4FitterAbstract(iName,"TGo4Fitter object"),
56  fxDatas(), fxModels(), fiFitFunctionType(0), fiMemoryUsage(100), fxUserFitFunction(0), fxDrawObjs(0)
57 {
58  fxDatas.SetOwner(kTRUE);
59  fxModels.SetOwner(kTRUE);
60  SetFitFunctionType(iFitFunctionType);
61  if (IsAddStandardActions) AddStandardActions();
62 }
63 
65 {
68 }
69 
70 void TGo4Fitter::SetMemoryUsage(Int_t iMemoryUsage)
71 {
72  if (iMemoryUsage<0) fiMemoryUsage=0; else
73  if (iMemoryUsage>100) fiMemoryUsage=100; else
74  fiMemoryUsage = iMemoryUsage;
75 }
76 
78 {
80 
81  for(Int_t ndata=0;ndata<GetNumData();ndata++)
82  GetData(ndata)->CollectParsTo(*this);
83 
84  for(Int_t nmodel=0;nmodel<GetNumModel();nmodel++) {
85  TGo4FitModel* model = GetModel(nmodel);
86  for (Int_t n=0;n<model->NumAssigments();n++)
87  if (FindData(model->AssignmentName(n))) {
88  model->CollectParsTo(*this);
89  break;
90  }
91  }
92 }
93 
94 void TGo4Fitter::Clear(Option_t* option)
95 {
97  DeleteAllData();
99 }
100 
102 {
103  return (n>=0) && (n<GetNumData()) ? dynamic_cast<TGo4FitData*> (fxDatas[n]) : 0;
104 }
105 
106 const char* TGo4Fitter::GetDataName(Int_t n)
107 {
108  return GetData(n) ? GetData(n)->GetName() : 0;
109 }
110 
111 TGo4FitData* TGo4Fitter::FindData(const char* DataName)
112 {
113  return (DataName==0) ? 0 : dynamic_cast<TGo4FitData*> (fxDatas.FindObject(DataName));
114 }
115 
117 {
118  fxDatas.Add(data);
121  return data;
122 }
123 
124 TGo4FitDataHistogram* TGo4Fitter::AddH1(const char* DataName, TH1* histo, Bool_t Owned, Double_t lrange, Double_t rrange)
125 {
126  TGo4FitDataHistogram *data = new TGo4FitDataHistogram(DataName, histo, Owned);
127  if ((lrange<rrange) || (rrange!=0.)) data->SetRange(0,lrange,rrange);
128  AddData(data);
129  return data;
130 }
131 
132 TGo4FitDataHistogram* TGo4Fitter::SetH1(const char* DataName, TH1* histo, Bool_t Owned)
133 {
134  TGo4FitDataHistogram* data = dynamic_cast<TGo4FitDataHistogram*> (FindData(DataName));
135  if (data!=0) data->SetHistogram(histo, Owned);
136  return data;
137 }
138 
139 
140 TGo4FitDataGraph* TGo4Fitter::AddGraph(const char* DataName, TGraph* gr, Bool_t Owned, Double_t lrange, Double_t rrange)
141 {
142  TGo4FitDataGraph *data = new TGo4FitDataGraph(DataName, gr, Owned);
143  if ((lrange<rrange) || (rrange!=0.)) data->SetRange(0,lrange,rrange);
144  AddData(data);
145  return data;
146 }
147 
148 TGo4FitDataGraph* TGo4Fitter::SetGraph(const char* DataName, TGraph* gr, Bool_t Owned)
149 {
150  TGo4FitDataGraph *data = dynamic_cast<TGo4FitDataGraph*> (FindData(DataName));
151  if (data!=0) data->SetGraph(gr, Owned);
152  return data;
153 }
154 
156 {
157  if (comp)
158  for (Int_t n=0;n<comp->NumSlots();n++) {
159  TGo4FitSlot* slot = comp->GetSlot(n);
160  ClearSlot(slot, kFALSE);
161  for(Int_t n2=0;n2<NumSlots();n2++) {
162  TGo4FitSlot* slot2 = GetSlot(n2);
163  if (slot2->GetConnectedSlot()==slot)
164  slot2->ClearConnectionToSlot();
165  }
166 
167  }
168 }
169 
170 TGo4FitData* TGo4Fitter::RemoveData(const char* DataName, Bool_t IsDel)
171 {
172  TGo4FitData* dat = FindData(DataName);
173  if (dat) {
176  fxDatas.Remove(dat);
177  if(IsDel) { CheckSlotsBeforeDelete(dat); delete dat; dat = 0; }
178  fxDatas.Compress();
179  }
180  return dat;
181 }
182 
184 {
185  fxDatas.Delete();
186  fxDatas.Compress();
189 }
190 
192 {
193  return (n>=0) && (n<GetNumModel()) ? dynamic_cast<TGo4FitModel*> (fxModels[n]) : 0;
194 }
195 
196 TGo4FitModel* TGo4Fitter::FindModel(const char* ModelName)
197 {
198  return (ModelName==0) ? 0 : dynamic_cast<TGo4FitModel*> (fxModels.FindObject(ModelName));
199 }
200 
202 {
203  fxModels.Add(model);
206  return model;
207 }
208 
209 TGo4FitModel* TGo4Fitter::AddModel(const char* DataName, TGo4FitModel* model)
210 {
211  model->AssignToData(DataName);
212  fxModels.Add(model);
215  return model;
216 }
217 
219 {
220  if(n<GetNumData()) return GetData(n); else return GetModel(n-GetNumData());
221 }
222 
223 void TGo4Fitter::AddPolynomX(const char* DataName, const char* NamePrefix, Int_t MaxOrder, Int_t GroupIndex, Double_t lrange, Double_t rrange)
224 {
225  if (DataName==0) return;
226 
227  Bool_t flag = kFALSE;
228  Int_t NumTry = 0;
229  Bool_t createmodel = kFALSE;
230 
231  do {
232  TString Prefix(NamePrefix);
233  if (NumTry>0) Prefix+=NumTry;
234  flag = kFALSE;
235  Bool_t findsame = kFALSE;
236 
237  for(Int_t Order=0; Order<=MaxOrder; Order++) {
238  TString Name(Prefix);
239  Name += "_";
240  Name += Order;
241  findsame = FindModel(Name.Data());
242  if (findsame && !createmodel) break;
243 
244  if (createmodel) {
245  TGo4FitModelPolynom* comp = new TGo4FitModelPolynom(Name, Order);
246  comp->SetGroupIndex(GroupIndex);
247  if ((lrange<rrange) || (rrange!=0.)) comp->SetRange(0,lrange,rrange);
248  AddModel(DataName, comp);
249  }
250  }
251 
252  flag = kTRUE;
253  if (findsame) { NumTry++; createmodel = kFALSE; } else
254  if (createmodel) flag = kFALSE;
255  else createmodel = kTRUE;
256 
257  } while (flag);
258 }
259 
260 void TGo4Fitter::AddPolynomX(const char* DataName, const char* NamePrefix, TArrayD& Coef, Int_t GroupIndex)
261 {
262  if (DataName==0) return;
263 
264  Bool_t flag = kFALSE;
265  Int_t NumTry = 0;
266  Bool_t createmodel = kFALSE;
267 
268  do {
269  Bool_t findsame = kFALSE;
270 
271  for (Int_t n=0;n<Coef.GetSize();n++) {
272  TString Name(NamePrefix);
273  if (NumTry>0) Name+=NumTry;
274  Name+="_";
275  Name+=n;
276  findsame = FindModel(Name.Data());
277  if (findsame && !createmodel) break;
278 
279  if (createmodel) {
280  TGo4FitModelPolynom* model = new TGo4FitModelPolynom(Name, n);
281  model->SetAmplValue(Coef[n]);
282  model->SetGroupIndex(GroupIndex);
283  AddModel(DataName, model);
284  }
285  }
286 
287  flag = kTRUE;
288  if (findsame) { NumTry++; createmodel = kFALSE; } else
289  if (createmodel) flag = kFALSE;
290  else createmodel = kTRUE;
291 
292  } while(flag);
293 }
294 
295 void TGo4Fitter::AddPolynoms(const char* DataName, const char* NamePrefix, Int_t MaxOrder, Int_t NumAxis, Int_t GroupIndex)
296 {
297  if (DataName==0) return;
298  TArrayD Orders(NumAxis);
299  Orders.Reset(0.);
300 
301  Bool_t flag = kFALSE;
302  Int_t NumTry = 0;
303  Bool_t createmodel = kFALSE;
304 
305  do {
306 
307  TString Prefix(NamePrefix);
308  if (NumTry>0) Prefix+=NumTry;
309  flag = kFALSE;
310  Bool_t findsame = kFALSE;
311 
312  do {
313  TString Name(Prefix);
314  for(Int_t n=0; n<NumAxis; n++) {
315  Name+="_";
316  Name+=Int_t(Orders[n]);
317  }
318  findsame = FindModel(Name.Data());
319  if (findsame && !createmodel) break;
320 
321  if (createmodel) {
322  TGo4FitModelPolynom* comp = new TGo4FitModelPolynom(Name, Orders);
323  comp->SetGroupIndex(GroupIndex);
324  AddModel(DataName, comp);
325  }
326 
327  Int_t nn = 0;
328  do {
329  Orders[nn] += 1.;
330  if (Orders[nn]<=MaxOrder) break;
331  Orders[nn]=0;
332  nn++;
333  } while (nn<NumAxis);
334  flag = (nn<NumAxis);
335  } while (flag);
336 
337  flag = kTRUE;
338  if (findsame) { NumTry++; createmodel = kFALSE; } else
339  if (createmodel) flag = kFALSE;
340  else createmodel = kTRUE;
341 
342  } while (flag);
343 }
344 
345 TGo4FitModelGauss1* TGo4Fitter::AddGauss1(const char* DataName, const char* ModelName, Double_t iPosition, Double_t iWidth, Double_t iAmpl, Int_t Axis)
346 {
347  TGo4FitModelGauss1* gauss = new TGo4FitModelGauss1(ModelName, iPosition, iWidth, Axis);
348  gauss->SetAmplValue(iAmpl);
349  AddModel(DataName, gauss);
350  return gauss;
351 }
352 
353 TGo4FitModel* TGo4Fitter::CloneModel(const char* ModelName, const char* NewName)
354 {
355  TGo4FitModel* mod = FindModel(ModelName);
356  if (mod==0) return 0;
357 
358  TString newname;
359  int cnt = 0;
360 
361  do {
362  newname = NewName ? NewName : mod->GetName();
363  if (cnt>0) {
364  newname += "_";
365  newname += cnt;
366  }
367  cnt++;
368  } while (FindModel(newname.Data()));
369 
370  TGo4FitModel* newmod = (TGo4FitModel*) mod->Clone(newname.Data());
371 
372  return AddModel(newmod);
373 }
374 
375 TGo4FitModel* TGo4Fitter::RemoveModel(const char * ModelName, Bool_t IsDel)
376 {
377  TGo4FitModel* mod = FindModel(ModelName);
378  if (mod) {
381  fxModels.Remove(mod);
382  if(IsDel) { CheckSlotsBeforeDelete(mod); delete mod; mod = 0; }
383  fxModels.Compress();
384  }
385  return mod;
386 }
387 
388 Int_t TGo4Fitter::NumModelsAssosiatedTo(const char* DataName)
389 {
390  Int_t res = 0;
391  for (Int_t n=0;n<GetNumModel();n++)
392  if (GetModel(n)->IsAssignTo(DataName)) res++;
393  return res;
394 }
395 
396 void TGo4Fitter::DeleteModelsAssosiatedTo(const char* DataName)
397 {
398  Int_t n=0;
399  while (n<GetNumModel()) {
400  TGo4FitModel* model = GetModel(n++);
401  if (model->IsAssignTo(DataName)){
402  if (model->NumAssigments()==1) {
405  fxModels.Remove(model);
406  delete model;
407  fxModels.Compress();
408  n--;
409  } else model->ClearAssignmentTo(DataName);
410  }
411  }
412 }
413 
415 {
416  fxModels.Delete();
417  fxModels.Compress();
420 }
421 
422 void TGo4Fitter::AssignModelTo(const char* ModelName, const char* DataName, Double_t RatioValue, Bool_t FixRatio)
423 {
424  TGo4FitModel* model = FindModel(ModelName);
425  if (model) {
426  if (DataName!=0)
427  model->AssignToData(DataName, RatioValue, FixRatio);
428  else {
429  model->ClearAssignments();
430  for (int n=0; n<GetNumData(); n++)
431  model->AssignToData(GetData(n)->GetName(), RatioValue, FixRatio);
432  }
434  }
435 }
436 
437 void TGo4Fitter::ClearModelAssignmentTo(const char* ModelName, const char* DataName)
438 {
439  TGo4FitModel* model = FindModel(ModelName);
440  if (model==0) return;
441  if (DataName!=0) model->ClearAssignmentTo(DataName);
442  else model->ClearAssignments();
444 }
445 
446 void TGo4Fitter::ChangeDataNameInAssignments(const char* oldname, const char* newname)
447 {
448  for(Int_t n=0;n<GetNumModel();n++)
449  GetModel(n)->ChangeDataNameInAssignments(oldname,newname);
450 }
451 
453 {
454  Int_t dbuf = -1, mbuf = -1;
455  switch (GetMemoryUsage()) {
456  case 0: dbuf = 0; mbuf = 0; break;
457  case 1: dbuf = 1; mbuf = 0; break;
458  case 2: dbuf = 1; mbuf = 1; break;
459  default: dbuf = -1; mbuf = -1;
460  }
461 
462  for(Int_t i1=0;i1<GetNumData();i1++) {
463  TGo4FitData *data = GetData(i1);
464  if (!data->Initialize(dbuf)) return kFALSE;
465  for(Int_t i2=0;i2<GetNumModel();i2++)
467  }
468 
469  for(Int_t i2=0;i2<GetNumModel();i2++) {
470  TGo4FitModel *fModel = GetModel(i2);
471  if (!fModel->Initialize(mbuf)) return kFALSE;
472  }
473 
474  if( (fiFitFunctionType == ff_user) && (fxUserFitFunction==0) ) {
475  std::cout << " User fit function not set. Switch to least squares " << std::endl;
477  }
478 
479  return kTRUE;
480 }
481 
483 {
484  for(Int_t i=0;i<GetNumData();i++) GetData(i)->Finalize();
485 
486  for(Int_t i=0;i<GetNumModel();i++) GetModel(i)->Finalize();
487 }
488 
489 Double_t TGo4Fitter::PointFitFunction(Int_t FitFunctionType, Double_t value, Double_t modelvalue, Double_t standdev)
490 {
491  switch (FitFunctionType) {
492  case ff_least_squares: {
493  Double_t zn1 = (value-modelvalue);
494  return zn1*zn1; }
495  case ff_chi_square: {
496  if (standdev<=0.) return 0.;
497  Double_t zn2 = (value-modelvalue);
498  return zn2*zn2/standdev; }
499  case ff_chi_Pearson: {
500  if (modelvalue<=0.) return 0.;
501  Double_t zn3 = (value-modelvalue);
502  return zn3*zn3/modelvalue; }
503  case ff_chi_Neyman: {
504  Double_t zn4 = (value-modelvalue);
505  return zn4*zn4/((value<1.) ? 1. : value); }
506  case ff_chi_gamma: {
507  if (value<0.) return 0.;
508  Double_t zn5 = (value+((value<1.) ? 0. : 1.)-modelvalue);
509  return zn5*zn5/(value+1.); }
510  case ff_ML_Poisson:
511  if (modelvalue<=0.) return 0.;
512  return modelvalue - value*TMath::Log(modelvalue);
513  case ff_user:
514  if (fxUserFitFunction==0) return 0;
515  return fxUserFitFunction(value, modelvalue, standdev);
516  default:
517  return (value-modelvalue)*(value-modelvalue);
518  }
519 }
520 
521 Double_t TGo4Fitter::CalculateFCN(Int_t FitFunctionType, TGo4FitData* selectdata)
522 {
523  if (GetMemoryUsage()>0) RebuildAll();
524 
525  Double_t fSum = 0.;
526 
527  for(Int_t n=0;n<GetNumData();n++) {
528  TGo4FitData* dat = GetData(n);
529  if (selectdata && (dat!=selectdata)) continue;
530  Double_t DataAmpl = dat->GetAmplValue();
531 
532  if (dat->BuffersAllocated()) {
533  Int_t size = dat->GetBinsSize();
534  Double_t* values = dat->GetBinsValues();
535  Double_t* res = dat->GetBinsResult();
536  Double_t* devs = dat->GetBinsDevs();
537 
538  for(Int_t nbin=0;nbin<size;nbin++)
539  fSum += PointFitFunction(FitFunctionType, DataAmpl*values[nbin], res[nbin], DataAmpl*DataAmpl*devs[nbin] );
540  } else {
541 
542  TGo4FitDataIter* iter = dat->MakeIter();
543  if (!iter->Reset()) { delete iter; continue; }
544 
545  TObjArray Models;
546 
547  TArrayD Ampls(GetNumModel());
548  for(Int_t nm=0;nm<GetNumModel();nm++)
549  if (GetModel(nm)->IsAssignTo(dat->GetName())) {
550  TGo4FitModel* model = GetModel(nm);
551  Models.Add(model);
552  model->BeforeEval(iter->ScalesSize());
553  Ampls[Models.GetLast()] = model->GetAmplValue() * model->GetRatioValueFor(dat->GetName());
554  }
555 
556  do {
557  Double_t value = DataAmpl * iter->Value();
558  Double_t deviat = DataAmpl * DataAmpl * iter->StandardDeviation();
559 
560  Double_t modelvalue = 0.;
561  for(Int_t nm=0;nm<=Models.GetLast();nm++) {
562  TGo4FitModel* model = dynamic_cast<TGo4FitModel*> (Models.At(nm));
563  modelvalue += Ampls[nm] * model->EvaluateAtPoint(iter);
564  }
565 
566  fSum += PointFitFunction(FitFunctionType, value, modelvalue, deviat);
567 
568  } while (iter->Next());
569 
570  for(Int_t nm=0;nm<=Models.GetLast();nm++)
571  ((TGo4FitModel*) Models.At(nm))->AfterEval();
572 
573  delete iter;
574  }
575  }
576 
577  return fSum;
578 }
579 
580 Double_t TGo4Fitter::CalculateFitFunction(Double_t* pars, Int_t FitFunctionType, const char* DataName)
581 {
582 
583  if (FitFunctionType<0) FitFunctionType = GetFitFunctionType();
584 
585  if (pars) {
586  ExecuteDependencies(pars);
587  SetParsValues(pars);
588  }
589 
590  return CalculateFCN(FitFunctionType, FindData(DataName));
591 }
592 
593 Int_t TGo4Fitter::CalculateNDF(const char* DataName)
594 {
595  TGo4FitData* selectdata = FindData(DataName);
596 
597  Int_t NDF = 0;
598 
599  if (selectdata!=0) {
600  NDF = selectdata->DefineBinsSize();
601  for(Int_t nm=0;nm<GetNumModel();nm++) {
602  TGo4FitModel* model = GetModel(nm);
603  if (model->IsAssignTo(selectdata->GetName()))
604  NDF -= model->NumFreePars();
605  }
606  } else {
607  for (Int_t n=0;n<GetNumData();n++)
608  NDF += GetData(n)->DefineBinsSize();
609  NDF -= NumFreePars();
610  }
611 
612  return NDF;
613 }
614 
616 {
618 }
619 
621 {
622  return CalculateNDF(0);
623 }
624 
625 void TGo4Fitter::RebuildAll(Bool_t ForceBuild)
626 {
627  for(Int_t i2=0;i2<GetNumModel();i2++)
628  GetModel(i2)->RebuildShape(ForceBuild);
629 
630  for(Int_t i1=0;i1<GetNumData();i1++) {
631  TGo4FitData* data = GetData(i1);
632  if (!data->BuffersAllocated()) continue;
633 
634  Int_t size = data->GetBinsSize();
635  Double_t* result = data->GetBinsResult();
636  for (Int_t nbin=0;nbin<size;nbin++) result[nbin] = 0.;
637 
638  for(Int_t i2=0;i2<GetNumModel();i2++) {
639  TGo4FitModel *model = GetModel(i2);
640  if (model->IsAssignTo(data->GetName()))
641  model->AddModelToDataResult(data);
642  }
643  }
644 }
645 
646 void TGo4Fitter::EstimateAmplitudes(Int_t NumIters)
647 {
648  TGo4FitAmplEstimation abc("this",NumIters);
649  abc.DoAction(this);
650 }
651 
652 void TGo4Fitter::AddAmplEstimation(Int_t NumIters)
653 {
654  AddAction(new TGo4FitAmplEstimation("AmplEstim",NumIters));
655 }
656 
658 {
660  AddSimpleMinuit();
661 }
662 
663 void TGo4Fitter::FillSlotList(TSeqCollection* list)
664 {
666  for(Int_t i=0;i<GetNumComp();i++)
667  GetComp(i)->FillSlotList(list);
668 }
669 
671 {
672  return model==0 ? kFALSE : model->BuffersAllocated();
673 }
674 
676 {
677  return data==0 ? kFALSE : data->BuffersAllocated();
678 }
679 
681 {
682  return (data==0) ? 0 : data->GetBinsSize();
683 }
684 
686 {
687  return (data==0) ? 0 : data->GetBinsValues();
688 }
689 
691 {
692  return (data==0) ? 0 : data->GetBinsDevs();
693 }
694 
696 {
697  return (data==0) ? 0 : data->GetBinsResult();
698 }
699 
700 Double_t* TGo4Fitter::GetModelBinsValues(TGo4FitModel* model, const char* DataName)
701 {
702  return (model==0) ? 0 : model->GetModelBins(DataName);
703 }
704 
705 void TGo4Fitter::Print(Option_t* option) const
706 {
707  TString Opt(option);
708  if (Opt=="Ampls") { PrintAmpls(); return; } else
709  if (Opt=="Pars") { PrintPars(); return; } else
710  if (Opt=="Results") { PrintResults(); return; } else
711  if (Opt=="Lines") { PrintLines(); return; }
712 
714  std::cout << "Fitiing function type: ";
715  switch (fiFitFunctionType) {
716  case ff_chi_square : std::cout << "ff_chi_square" << std::endl; break;
717  case ff_chi_Pearson : std::cout << "ff_chi_Pearson" << std::endl; break;
718  case ff_chi_Neyman : std::cout << "ff_chi_Neyman" << std::endl; break;
719  case ff_chi_gamma : std::cout << "ff_chi_gamma" << std::endl; break;
720  case ff_ML_Poisson : std::cout << "ff_ML_Poisson" << std::endl; break;
721  case ff_user : std::cout << "user defined" << std::endl; break;
722  default: std::cout << "ff_least_squares" << std::endl;
723  }
724  std::cout << std::endl << " LIST OF DATA OBJECTS" << std::endl;
725  fxDatas.Print(option);
726  std::cout << std::endl << " LIST OF MODEL OBJECTS" << std::endl;
727  fxModels.Print(option);
728 }
729 
730 Bool_t TGo4Fitter::CalculatesMomentums(const char* DataName, Bool_t UseRanges, Bool_t SubstractModels, Double_t& first, Double_t& second)
731 {
732  TGo4FitData* data = FindData(DataName);
733  if (data==0) return kFALSE;
734 
735  TGo4FitDataIter* iter = data->MakeIter();
736  if (iter==0) return kFALSE;
737 
738  Int_t size = iter->CountPoints(UseRanges);
739  if (size==0) { delete iter; return kFALSE; }
740 
741  TObjArray Models;
742  if (SubstractModels)
743  for(Int_t nm=0; nm<GetNumModel(); nm++)
744  if (GetModel(nm)->IsAssignTo(DataName))
745  Models.Add(GetModel(nm));
746 
747  TArrayD Ampls(Models.GetLast()+1);
748  for(Int_t n=0;n<=Models.GetLast();n++) {
749  TGo4FitModel* model = (TGo4FitModel*) Models[n];
750  model->BeforeEval(iter->ScalesSize());
751  Ampls[n] = model->GetAmplValue() * model->GetRatioValueFor(DataName);
752  }
753 
754  TArrayD bins(size), scales(size);
755 
756  Int_t ppnt=0;
757 
758  if (iter->Reset(UseRanges)) do {
759  Double_t value = iter->Value();
760  for(Int_t n=0;n<=Models.GetLast();n++) {
761  TGo4FitModel* model = (TGo4FitModel*) Models[n];
762  value -= Ampls[n] * model->EvaluateAtPoint(iter);
763  }
764  value = TMath::Abs(value);
765  bins[ppnt] = value;
766  scales[ppnt] = iter->x();
767  ppnt++;
768  } while (iter->Next(UseRanges));
769 
770  delete iter;
771 
772  for(Int_t n=0;n<=Models.GetLast();n++) {
773  TGo4FitModel* model = (TGo4FitModel*) Models[n];
774  model->AfterEval();
775  }
776 
777  Int_t niter=0;
778 
779  do {
780  Double_t sum00 = 0., sum11 = 0., sum22 = 0.;
781  for (Int_t pnt=0;pnt<size;pnt++)
782  if ((bins[pnt]>0.) && ((niter==0) || (TMath::Abs(scales[pnt]-first)<second*2.))) {
783  sum00 += bins[pnt];
784  sum11 += bins[pnt]*scales[pnt];
785  sum22 += bins[pnt]*scales[pnt]*scales[pnt];
786  }
787 
788  if (sum00>0.) {
789  Double_t mid = sum11/sum00;
790  Double_t dev = TMath::Sqrt(sum22/sum00-mid*mid);
791 
792  if (niter>0)
793  if ((dev/second>0.8) && (dev/second<1.2)) niter=10;
794 
795  first = mid; second = dev;
796 
797  } else niter=10;
798 
799  } while (niter++<8);
800 
801  return kTRUE;
802 }
803 
804 Double_t TGo4Fitter::CalculatesIntegral(const char* DataName, const char* ModelName, Bool_t onlycounts)
805 {
806  TGo4FitData* data = FindData(DataName);
807  if (data==0) return 0.;
808 
809  TGo4FitModel* model = ModelName!=0 ? FindModel(ModelName) : 0;
810  if ((ModelName!=0) && (model==0)) return 0.;
811 
812  TGo4FitDataIter* iter = data->MakeIter();
813  if (iter==0) return 0.;
814 
815  double sum = 0.;
816  double ampl = 1.;
817 
818  if (!iter->Reset(kTRUE)) {
819  delete iter;
820  return 0.;
821  }
822 
823  if (model!=0) {
824  model->BeforeEval(iter->ScalesSize());
825  ampl = model->GetAmplValue() * model->GetRatioValueFor(DataName);
826  }
827 
828  do {
829  double dx = onlycounts ? 1. : iter->xWidths();
830  double value = model==0 ? iter->Value() : model->EvaluateAtPoint(iter);
831  sum += ampl*value*dx;
832  } while (iter->Next(kTRUE));
833 
834  if (model!=0)
835  model->AfterEval();
836 
837  delete iter;
838 
839  return sum;
840 }
841 
842 Double_t TGo4Fitter::CalculatesModelIntegral(const char* ModelName, Bool_t OnlyCounts)
843 {
844  TGo4FitModel* model = FindModel(ModelName);
845  if (model==0) return 0.;
846  return CalculatesIntegral(model->AssignmentName(0), ModelName, OnlyCounts);
847 }
848 
849 TObject* TGo4Fitter::CreateDrawObject(const char* ResName, const char* DataName, Bool_t IsModel, const char* ModelName)
850 {
851  TGo4FitData* data = FindData(DataName);
852  if (data==0) return 0;
853 
854  TObjArray Models;
855  if (IsModel) {
856  TGo4FitModel* model = FindModel(ModelName);
857  if (model) {
858  Models.Add(model);
859  } else {
860  Int_t groupindex = -1;
861 
862  if (ModelName!=0) {
863  TString modelname(ModelName);
864  if (modelname=="Background") groupindex = 0; else
865  if (modelname.Index("Group",0,TString::kExact)==0) {
866  modelname.Remove(0,5);
867  char* err = 0;
868  groupindex = strtol(modelname.Data(),&err,10);
869  if (err && (*err!=0)) groupindex=-1;
870  }
871  }
872 
873  for(Int_t nm=0; nm<GetNumModel(); nm++) {
874  model = GetModel(nm);
875  if (model->IsAssignTo(DataName))
876  if ((groupindex<0) || (model->GetGroupIndex()==groupindex))
877  Models.Add(model);
878  }
879  }
880  if (Models.GetLast()<0) return 0;
881  }
882 
883  TGo4FitDataIter* iter = data->MakeIter();
884  if (iter==0) return 0;
885 
886  if (!iter->Reset(kFALSE)) { delete iter; return 0; }
887 
888  TH1* histo = 0;
889  Int_t ndim = 0;
890  TGraph* gr = 0;
891  Bool_t UseRanges = kTRUE;
892 
893  if (iter->HasIndexes() && (iter->IndexesSize()==iter->ScalesSize()) && iter->HasWidths()) {
894  histo = iter->CreateHistogram(ResName, kFALSE, !IsModel);
895 
896  if (histo) {
897  ndim = histo->GetDimension();
898  UseRanges = Models.GetLast() >= 0;
899  }
900  } else {
901  gr = iter->CreateGraph(ResName, kTRUE, !IsModel);
902  UseRanges = kTRUE;
903  }
904 
905  if (((histo!=0) || (gr!=0)) && IsModel) {
906  TArrayD Ampls(Models.GetLast()+1);
907 
908  for(Int_t n=0;n<=Models.GetLast();n++) {
909  TGo4FitModel* model = (TGo4FitModel*) Models[n];
910  model->BeforeEval(iter->ScalesSize());
911  Ampls[n] = model->GetAmplValue() * model->GetRatioValueFor(DataName);
912  }
913 
914  if (iter->Reset(UseRanges)) do {
915  Double_t zn = 0;
916  for(Int_t n=0;n<=Models.GetLast();n++) {
917  TGo4FitModel* model = (TGo4FitModel*) Models[n];
918  zn += Ampls[n] * model->EvaluateAtPoint(iter);
919  }
920  if (histo)
921  switch(ndim) {
922  case 1: histo->SetBinContent(iter->Indexes()[0]+1,zn); break;
923  case 2: histo->SetBinContent(iter->Indexes()[0]+1,iter->Indexes()[1]+1,zn); break;
924  case 3: histo->SetBinContent(iter->Indexes()[0]+1,iter->Indexes()[1]+1,iter->Indexes()[2]+1,zn); break;
925  }
926  if(gr) {
927  (gr->GetX())[iter->Point()] = iter->x();
928  (gr->GetY())[iter->Point()] = zn;
929  }
930 
931  } while (iter->Next(UseRanges));
932 
933  for(Int_t n=0;n<=Models.GetLast();n++) {
934  TGo4FitModel* model = (TGo4FitModel*) Models[n];
935  model->AfterEval();
936  }
937  }
938 
939  delete iter;
940 
941  if (gr && IsModel)
942  for(Int_t n1=0;n1<gr->GetN()-1;n1++)
943  for(Int_t n2=n1+1;n2<gr->GetN();n2++)
944  if ((gr->GetX())[n1]>(gr->GetX())[n2]) {
945  Double_t xx = (gr->GetX())[n1];
946  (gr->GetX())[n1] = (gr->GetX())[n2];
947  (gr->GetX())[n2] = xx;
948  Double_t yy = (gr->GetY())[n1];
949  (gr->GetY())[n1] = (gr->GetY())[n2];
950  (gr->GetY())[n2] = yy;
951  }
952 
953  TNamed* res = 0;
954  if (histo) res = histo; else res = gr;
955  if (res) {
956  TString title;
957  if (IsModel)
958  if (ModelName) { title = "Draw of model "; title+=ModelName; }
959  else { title = "Draw of full model of "; title+=DataName; }
960  else { title = "Draw of data "; title+=DataName; }
961  res->SetTitle(title.Data());
962  }
963 
964  return res;
965 }
966 
967 
968 /* valid format for draw options
969 
970  "" - draw all data
971  "*" - draw all data with models
972  "**" - draw all data with models and all components
973  "DataName" - draw select data with model
974  "DataName-" - draw select data without model
975  "DataName*" - draw select data with model and components
976  "DataName-, Background" - draw data and sum of components, belongs to background group (index=0)
977  "DataName-, Group1" - draw data and sum of components, belongs to group 1
978  "ModelName" - draw model component and correspondent data
979 */
980 
981 void TGo4Fitter::Draw(Option_t* option)
982 {
983  TString opt(option);
984 
985  TCanvas *fCanvas = 0;
986 
987  if ((opt.Length()>0) && (opt[0]=='#')) {
988  opt.Remove(0,1);
989  TString CanvasName;
990  Int_t n=0;
991  do {
992  CanvasName = "Canvas";
993  CanvasName+=n++;
994  } while (gROOT->FindObject(CanvasName.Data()));
995  fCanvas = new TCanvas(CanvasName,TString("Draw of fitter ")+GetName()+" "+opt,3);
996  fCanvas->cd();
997  }
998 
999  Bool_t drawdata = kFALSE;
1000  TGo4FitData *selectdata = 0;
1001  TObjArray selectmodels;
1002  selectmodels.SetOwner(kTRUE);
1003  Bool_t drawcomp = kFALSE;
1004  Bool_t drawmodel = kFALSE;
1005 
1006  if (opt=="*") { opt = ""; drawdata = kTRUE; }
1007 
1008  while (opt.Length()>0) {
1009  Int_t len = opt.Index(",",0,TString::kExact);
1010 
1011  if (len<0) len = opt.Length();
1012  if (len==0) break;
1013 
1014  TString optpart(opt.Data(), len);
1015  while ((optpart.Length()>0) && (optpart[0]==' ')) optpart.Remove(0,1);
1016 
1017  Bool_t find = kFALSE;
1018 
1019  for(Int_t n=0;n<GetNumData();n++)
1020  if (optpart.Index(GetDataName(n),0,TString::kExact)==0) {
1021  selectdata = GetData(n);
1022  drawdata = kTRUE;
1023  drawmodel = kTRUE;
1024  find = kTRUE;
1025  optpart.Remove(0,strlen(GetDataName(n)));
1026  if (optpart=="*") drawcomp = kTRUE; else
1027  if (optpart=="-") drawmodel = kFALSE;
1028  break;
1029  }
1030 
1031  if (!find) {
1032  selectmodels.Add(new TObjString(optpart));
1033  TGo4FitModel* model = FindModel(optpart.Data());
1034  if (model && (selectdata==0))
1035  for(Int_t n=0;n<GetNumData();n++)
1036  if (model->IsAssignTo(GetDataName(n)))
1037  selectdata = GetData(n);
1038  }
1039 
1040  opt.Remove(0,len);
1041  if (opt.Length()>0) opt.Remove(0,1);
1042  }
1043 
1044  if ((selectdata==0) && !drawdata)
1045  if (GetNumData()>0) selectdata = GetData(0);
1046 
1048  fxDrawObjs = new TObjArray();
1049 
1050  for(Int_t nn = 0; nn < GetNumData(); nn++) {
1051  TGo4FitData *data = GetData(nn);
1052  if (selectdata && (data!=selectdata)) continue;
1053 
1054  if (drawdata) {
1055  TObject* obj = CreateDrawObject(TString(data->GetName())+"_bins", data->GetName(), kFALSE);
1056 
1057  TAttLine* line = dynamic_cast<TAttLine*> (obj);
1058  if (line) {
1059  line->SetLineColor(1);
1060  line->SetLineWidth(1);
1061  }
1062  if (obj) fxDrawObjs->Add(obj);
1063  }
1064 
1065  if (drawmodel) {
1066  TObject* mobj = CreateDrawObject(TString(data->GetName())+"_fullmodel", data->GetName(), kTRUE);
1067 
1068  TAttLine* line = dynamic_cast<TAttLine*> (mobj);
1069  if (line) {
1070  line->SetLineColor(4);
1071  line->SetLineWidth(2);
1072  }
1073  if (mobj) fxDrawObjs->Add(mobj);
1074  }
1075 
1076  if (drawcomp)
1077  for(Int_t nmodel=0;nmodel<GetNumModel();nmodel++) {
1078  TGo4FitModel *model = GetModel(nmodel);
1079  if ( !model->IsAssignTo(data->GetName())) continue;
1080 
1081  TObject* cobj = CreateDrawObject(TString(data->GetName())+"_"+model->GetName(), data->GetName(), kTRUE, model->GetName());
1082 
1083  TAttLine* line = dynamic_cast<TAttLine*> (cobj);
1084  if (line) {
1085  line->SetLineColor(6);
1086  line->SetLineWidth(1);
1087  }
1088  if (cobj) fxDrawObjs->Add(cobj);
1089  }
1090 
1091  for (Int_t n=0;n<=selectmodels.GetLast();n++) {
1092  TString name = ((TObjString*) (selectmodels[n]))->String();
1093 
1094  TObject* cobj = CreateDrawObject(TString(data->GetName())+"_" +name, data->GetName(), kTRUE, name.Data());
1095 
1096  TAttLine* line = dynamic_cast<TAttLine*> (cobj);
1097  if (line) {
1098  line->SetLineColor(6);
1099  line->SetLineWidth(1);
1100  }
1101  if (cobj) fxDrawObjs->Add(cobj);
1102  }
1103  }
1104 
1105  if (fxDrawObjs->GetLast()==0) fxDrawObjs->At(0)->Draw(); else
1106  if (fxDrawObjs->GetLast()>0) {
1107  Bool_t allhisto = kTRUE;
1108  for (Int_t n=0;n<=fxDrawObjs->GetLast();n++)
1109  if (!(fxDrawObjs->At(n)->InheritsFrom(TH1::Class()))) allhisto = kFALSE;
1110  if (allhisto) {
1111  THStack* stack = new THStack(TString("Stack")+"_"+fxDrawObjs->At(0)->GetName(),fxDrawObjs->At(0)->GetName());
1112  for (Int_t n=0;n<=fxDrawObjs->GetLast();n++) stack->Add((TH1*)fxDrawObjs->At(n));
1113  fxDrawObjs->Clear();
1114  fxDrawObjs->Add(stack);
1115  stack->Draw("nostack");
1116  } else
1117  for (Int_t n=0;n<=fxDrawObjs->GetLast();n++)
1118  if (n==0) fxDrawObjs->At(n)->Draw("A*");
1119  else fxDrawObjs->At(n)->Draw("L");
1120  }
1121 
1122  if (fCanvas!=0) fCanvas->Update();
1123 }
1124 
1126 {
1127  std::cout << std::endl << "*** LIST OF AMPLITUDES VALUE ***" << std::endl;
1128  for(Int_t n=0;n<GetNumComp();n++) {
1129  TGo4FitComponent* comp = ((TGo4Fitter*) this)->GetComp(n);
1130  if (comp->GetAmplPar() != 0)
1131  std::cout << " " << comp->GetAmplFullName() << " " << comp->GetAmplValue() << " " << comp->GetAmplError() << std::endl;
1132  }
1133 }
1134 
1136 {
1137  std::cout << std::endl << " *** LIST OF LINES PARAMETERS ***" << std::endl;
1138 
1139  int MaxAxis = 0;
1140  for (Int_t n=0; n<GetNumModel();n++) {
1141  TGo4FitModel* m = ((TGo4Fitter*) this)->GetModel(n);
1142  if (m==0) continue;
1143  Double_t zn;
1144  for (int naxis=0;naxis<3;naxis++)
1145  if (m->GetPosition(naxis,zn) || m->GetWidth(naxis,zn)) MaxAxis = naxis;
1146  }
1147 
1148  std::cout << std::setw(10) << "Name" << std::setw(12) << "Ampl";
1149  for(Int_t naxis=0;naxis<=MaxAxis;naxis++)
1150  std::cout << std::setw(11) << "Pos" << naxis << std::setw(11) << "Width" << naxis;
1151  std::cout << std::endl;
1152 
1153  for (Int_t n=0; n<GetNumModel();n++) {
1154  TGo4FitModel* m = ((TGo4Fitter*) this)->GetModel(n);
1155  if (m==0) continue;
1156  std::cout << std::setw(10) << m->GetName() << std::setw(12) << m->GetAmplValue();
1157 
1158  for (int naxis=0;naxis<=MaxAxis;naxis++) {
1159  Double_t pos, width;
1160  std::cout << std::setw(12);
1161  if (m->GetPosition(naxis,pos)) std::cout << pos;
1162  else std::cout << "---";
1163  std::cout << std::setw(12);
1164  if (m->GetWidth(naxis,width)) std::cout << width;
1165  else std::cout << "---";
1166  }
1167  std::cout << std::endl;
1168  }
1169 }
1170 
1171 
1173 {
1174  if (fxDrawObjs) {
1175  lst.AddAll(fxDrawObjs);
1176  delete fxDrawObjs;
1177  fxDrawObjs = 0;
1178  }
1179 }
1180 
1182 {
1183  if (fxDrawObjs) {
1184  for(Int_t n=0;n<=fxDrawObjs->GetLast();n++)
1185  gROOT->Add(fxDrawObjs->At(n));
1186  delete fxDrawObjs;
1187  fxDrawObjs = 0;
1188  }
1189 }
1190 
1191 TString TGo4Fitter::FindNextName(const char* Head, Int_t start, Bool_t isModel)
1192 {
1193  TString name;
1194  Int_t n = start;
1195  do {
1196  name.Form("%s%d", Head, n++);
1197  } while (isModel ? FindModel(name.Data())!=0 : FindData(name.Data())!=0 );
1198  return name;
1199 }
1200 
1201 void TGo4Fitter::Streamer(TBuffer& b)
1202 {
1203  if (b.IsReading()) {
1204  TGo4Fitter::Class()->ReadBuffer(b, this);
1208  } else {
1210  TGo4Fitter::Class()->WriteBuffer(b, this);
1211  }
1212 }
virtual void DoAction(TGo4FitterAbstract *Fitter)
Double_t xWidths() const
TH1 * CreateHistogram(const char *HistoName, Bool_t UseRanges=kFALSE, Bool_t SetBins=kFALSE)
void SetUpdateSlotList()
virtual void Print(Option_t *option) const
Definition: TGo4Fitter.cxx:705
friend class TGo4FitAmplEstimation
Definition: TGo4Fitter.h:385
TGo4FitDataGraph * AddGraph(const char *DataName, TGraph *gr, Bool_t Owned=kFALSE, Double_t lrange=0., Double_t rrange=0.)
Definition: TGo4Fitter.cxx:140
void PrintLines() const
void RebuildAll(Bool_t ForceBuild=kFALSE)
Definition: TGo4Fitter.cxx:625
virtual Double_t DoCalculation()
Definition: TGo4Fitter.cxx:615
Bool_t HasWidths() const
Definition: TGo4FitData.h:499
virtual void Clear(Option_t *option="")
void SetMemoryUsage(Int_t iMemoryUsage)
Definition: TGo4Fitter.cxx:70
TObjArray * fxDrawObjs
Definition: TGo4Fitter.h:446
void SetHistogram(TH1 *iHistogram, Bool_t iHistogramOwned=kFALSE)
virtual void Finalize()
virtual Int_t DoNDFCalculation()
Definition: TGo4Fitter.cxx:620
TGo4FitSlot * GetSlot(Int_t nslot)
Double_t x() const
Definition: TGo4FitData.h:484
void ClearAssignments()
TGo4FitComponent * GetComp(Int_t n)
Definition: TGo4Fitter.cxx:218
Double_t * GetBinsValues()
Definition: TGo4FitData.h:255
void MoveDrawObjectsToROOT()
void SetFitFunctionType(Int_t iFitFunctionType)
Definition: TGo4Fitter.h:80
virtual void Draw(Option_t *option)
Definition: TGo4Fitter.cxx:981
void AddAmplEstimation(Int_t NumIters=1)
Definition: TGo4Fitter.cxx:652
void ClearModelAssignmentTo(const char *ModelName, const char *DataName=0)
Definition: TGo4Fitter.cxx:437
virtual void CollectAllPars()
Definition: TGo4Fitter.cxx:77
const char * GetAmplFullName()
virtual TGo4FitDataIter * MakeIter()
Definition: TGo4FitData.h:165
virtual void Finalize()
virtual void FinalizeFitterData()
Definition: TGo4Fitter.cxx:482
Int_t GetDataBinsSize(TGo4FitData *data)
Definition: TGo4Fitter.cxx:680
void EstimateAmplitudes(Int_t NumIters=1)
Definition: TGo4Fitter.cxx:646
TGo4FitModel * FindModel(const char *ModelName)
Definition: TGo4Fitter.cxx:196
Int_t GetGroupIndex() const
Definition: TGo4FitModel.h:119
void AddPolynoms(const char *DataName, const char *NamePrefix, Int_t MaxOrder=1, Int_t NumAxis=1, Int_t GroupIndex=0)
Definition: TGo4Fitter.cxx:295
void ClearConnectionToSlot()
Definition: TGo4FitSlot.cxx:93
Double_t * GetDataBinsResult(TGo4FitData *data)
Definition: TGo4Fitter.cxx:695
void AssignToData(const char *DataName, Double_t RatioValue=1., Bool_t FixRatio=kFALSE)
Double_t * GetModelBins(const char *DataName) const
Int_t NumModelsAssosiatedTo(const char *DataName)
Definition: TGo4Fitter.cxx:388
void AddAction(TGo4FitterAction *Action)
void CheckSlotsBeforeDelete(TGo4FitComponent *comp)
Definition: TGo4Fitter.cxx:155
TGraph * CreateGraph(const char *GraphName, Bool_t UseRanges=kFALSE, Bool_t SetBins=kFALSE)
void AddPolynomX(const char *DataName, const char *NamePrefix, Int_t MaxOrder=1, Int_t GroupIndex=0, Double_t lrange=0., Double_t rrange=0.)
Definition: TGo4Fitter.cxx:223
virtual Bool_t Initialize(Int_t UseBuffers=-1)
void ClearAssignmentTo(const char *DataName)
Double_t PointFitFunction(Int_t FitFunctionType, Double_t value, Double_t modelvalue, Double_t standdev)
Definition: TGo4Fitter.cxx:489
void AssignModelTo(const char *ModelName, const char *DataName, Double_t RatioValue=1., Bool_t FixRatio=kFALSE)
Definition: TGo4Fitter.cxx:422
void ClearSlot(TGo4FitSlot *slot, Bool_t NonOwned)
void SetAmplValue(Double_t iAmpl)
const char * AssignmentName(Int_t n)
Definition: TGo4FitModel.h:129
void SetOwner(TNamed *iOwner)
Definition: TGo4FitNamed.h:58
TObjArray fxModels
Definition: TGo4Fitter.h:424
TGo4FitData * RemoveData(const char *DataName, Bool_t IsDel=kFALSE)
Definition: TGo4Fitter.cxx:170
void DeleteAllData()
Definition: TGo4Fitter.cxx:183
TGo4FitModel * RemoveModel(const char *ModelName, Bool_t IsDel=kFALSE)
Definition: TGo4Fitter.cxx:375
Int_t fiMemoryUsage
Definition: TGo4Fitter.h:434
TUserFitFunction fxUserFitFunction
Definition: TGo4Fitter.h:444
Double_t * GetDataBinsDevs(TGo4FitData *data)
Definition: TGo4Fitter.cxx:690
virtual Bool_t Reset(Bool_t UseRanges=kTRUE)
Int_t DefineBinsSize()
virtual Double_t EvaluateAtPoint(TGo4FitData *data, Int_t nbin, Bool_t UseRanges=kTRUE)
Int_t fiFitFunctionType
Definition: TGo4Fitter.h:429
void DeleteModelsAssosiatedTo(const char *DataName)
Definition: TGo4Fitter.cxx:396
virtual void FillSlotList(TSeqCollection *list)
Definition: TGo4Fitter.cxx:663
Double_t CalculatesIntegral(const char *DataName, const char *ModelName=0, Bool_t OnlyCounts=kFALSE)
Definition: TGo4Fitter.cxx:804
Int_t Point() const
Definition: TGo4FitData.h:524
Int_t GetMemoryUsage()
Definition: TGo4Fitter.h:113
void PrepareSlotsForWriting()
TGo4FitData * AddData(TGo4FitData *d)
Definition: TGo4Fitter.cxx:116
virtual Bool_t GetPosition(Int_t naxis, Double_t &pos)
virtual void Print(Option_t *option) const
Double_t * GetBinsResult()
Definition: TGo4FitData.h:267
void ChangeDataNameInAssignments(const char *oldname, const char *newname)
void SetGraph(TGraph *iGraph, Bool_t iGraphOwned=kFALSE)
Bool_t BuffersAllocated() const
Definition: TGo4FitData.h:238
TGo4FitSlot * GetConnectedSlot() const
Definition: TGo4FitSlot.h:124
void CheckDuplicatesOnSlot()
TGo4FitModelGauss1 * AddGauss1(const char *DataName, const char *ModelName, Double_t iPosition, Double_t iWidth, Double_t iAmpl=1., Int_t Axis=0)
Definition: TGo4Fitter.cxx:345
void ChangeDataNameInAssignments(const char *oldname, const char *newname)
Definition: TGo4Fitter.cxx:446
void SetRange(Int_t naxis, Double_t min, Double_t max)
virtual void CollectAllPars()
Int_t CalculateNDF(const char *DataName=0)
Definition: TGo4Fitter.cxx:593
TGo4FitDataGraph * SetGraph(const char *DataName, TGraph *gr, Bool_t Owned=kFALSE)
Definition: TGo4Fitter.cxx:148
TGo4FitModel * GetModel(Int_t n)
Definition: TGo4Fitter.cxx:191
Double_t CalculatesModelIntegral(const char *ModelName, Bool_t OnlyCounts=kFALSE)
Definition: TGo4Fitter.cxx:842
Double_t * GetBinsDevs()
Definition: TGo4FitData.h:261
TObjArray fxDatas
Definition: TGo4Fitter.h:419
Int_t GetNumComp() const
Definition: TGo4Fitter.h:281
Int_t CountPoints(Bool_t UseRanges=kTRUE)
virtual Bool_t Initialize(Int_t UseBuffers=-1)
Int_t GetBinsSize() const
Definition: TGo4FitData.h:243
virtual Bool_t Next(Bool_t UseRanges=kTRUE)
virtual Bool_t BeforeEval(Int_t ndim)
void ConnectToDataIfAssigned(TGo4FitData *data)
TString FindNextName(const char *Head, Int_t start, Bool_t isModel=kTRUE)
virtual void FillSlotList(TSeqCollection *lst)
void PrintAmpls() const
virtual void AfterEval()
Definition: TGo4FitModel.h:253
TGo4FitData * GetData(Int_t n)
Definition: TGo4Fitter.cxx:101
void SetGroupIndex(Int_t index=-1)
Definition: TGo4FitModel.h:109
virtual Bool_t GetWidth(Int_t naxis, Double_t &width)
virtual ~TGo4Fitter()
Definition: TGo4Fitter.cxx:64
Int_t NumAssigments() const
Definition: TGo4FitModel.h:124
Int_t GetFitFunctionType()
Definition: TGo4Fitter.h:85
Double_t * GetModelBinsValues(TGo4FitModel *model, const char *DataName)
Definition: TGo4Fitter.cxx:700
Double_t CalculateFitFunction(Double_t *pars=0, Int_t FitFunctionType=-1, const char *DataName=0)
Definition: TGo4Fitter.cxx:580
virtual void CollectParsTo(TGo4FitParsList &list)
Bool_t DataBuffersAllocated(TGo4FitData *data)
Definition: TGo4Fitter.cxx:675
Bool_t HasIndexes() const
Definition: TGo4FitData.h:454
TGo4FitModel * CloneModel(const char *ModelName, const char *NewName=0)
Definition: TGo4Fitter.cxx:353
const char * GetDataName(Int_t n)
Definition: TGo4Fitter.cxx:106
Bool_t BuffersAllocated() const
Double_t Value() const
Definition: TGo4FitData.h:514
Int_t GetNumModel() const
Definition: TGo4Fitter.h:178
TGo4FitModel * AddModel(TGo4FitModel *m)
Definition: TGo4Fitter.cxx:201
const Int_t * Indexes() const
Definition: TGo4FitData.h:464
void PrintPars() const
TObject * CreateDrawObject(const char *ObjName, const char *DataName, Bool_t IsModel=kFALSE, const char *ModelName=0)
Definition: TGo4Fitter.cxx:849
virtual Bool_t InitFitterData()
Definition: TGo4Fitter.cxx:452
void RebuildShape(Bool_t ForceBuild=kFALSE)
Bool_t IsAssignTo(const char *DataName) const
Definition: TGo4FitModel.h:141
TGo4FitDataHistogram * SetH1(const char *DataName, TH1 *histo, Bool_t Owned=kFALSE)
Definition: TGo4Fitter.cxx:132
Bool_t CalculatesMomentums(const char *DataName, Bool_t UseRanges, Bool_t SubstractModels, Double_t &first, Double_t &second)
Definition: TGo4Fitter.cxx:730
Double_t * GetDataBinsValues(TGo4FitData *data)
Definition: TGo4Fitter.cxx:685
virtual void Clear(Option_t *option=0)
Definition: TGo4Fitter.cxx:94
TGo4FitData * FindData(const char *DataName)
Definition: TGo4Fitter.cxx:111
Double_t CalculateFCN(Int_t FitFunctionType, TGo4FitData *selectdata=0)
Definition: TGo4Fitter.cxx:521
Bool_t AddModelToDataResult(TGo4FitData *data)
Int_t IndexesSize() const
Definition: TGo4FitData.h:459
Double_t StandardDeviation() const
Definition: TGo4FitData.h:519
Double_t GetRatioValueFor(const char *DataName)
void ProvideLastDrawObjects(TObjArray &lst)
Bool_t ModelBuffersAllocated(TGo4FitModel *model)
Definition: TGo4Fitter.cxx:670
void AddStandardActions()
Definition: TGo4Fitter.cxx:657
Int_t ScalesSize() const
Definition: TGo4FitData.h:474
TGo4FitParameter * GetAmplPar()
void ExecuteDependencies(Double_t *pars)
void DeleteAllModels()
Definition: TGo4Fitter.cxx:414
Int_t GetNumData() const
Definition: TGo4Fitter.h:123
void SetParsValues(Double_t *pars)
TGo4FitDataHistogram * AddH1(const char *DataName, TH1 *histo, Bool_t Owned=kFALSE, Double_t lrange=0., Double_t rrange=0.)
Definition: TGo4Fitter.cxx:124