GSI Object Oriented Online Offline (Go4) GO4-6.4.0
Loading...
Searching...
No Matches
TGo4FitData.cxx
Go to the documentation of this file.
1// $Id$
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 "TGo4FitData.h"
15
16#include <iostream>
17
18#include "TBuffer.h"
19#include "TMath.h"
20#include "TH1.h"
21#include "TH2.h"
22#include "TH3.h"
23#include "TGraph.h"
24#include "TClass.h"
25
26#include "TGo4FitAxisTrans.h"
27
34
35TGo4FitData::TGo4FitData(const char *iName, const char *iTitle, Int_t iDataType, Bool_t AddAmpl)
36 : TGo4FitComponent(iName, iTitle), fiDataType(iDataType), fbUseBinScale(kFALSE), fiTakeSigmasFrom(1),
38{
39
41
42 if (AddAmpl)
43 NewAmplitude("Ampl", 1.0, kTRUE);
44
45 fxAxisTrans.SetOwner(kTRUE);
46}
47
52
54{
55 Int_t oldnum = GetNumberOfTransSlots();
56 if ((nslots < 0) || (nslots == oldnum))
57 return kFALSE;
58
59 if (oldnum < nslots)
60 for (Int_t n = oldnum; n < nslots; n++) {
61 TString name("Trans");
62 name += n;
63 fxAxisTrans.Add(new TGo4FitSlot(name.Data(), "Axis transformation", this, TGo4FitAxisTrans::Class(), kFALSE));
64 }
65 else
66 for (Int_t n = oldnum; n > nslots; n--) {
67 TObject *slot = fxAxisTrans.Last();
68 fxAxisTrans.Remove(slot);
69 fxAxisTrans.Compress();
70 delete slot;
71 }
73 return kTRUE;
74}
75
77{
78 return (nslot >= 0) && (nslot <= fxAxisTrans.GetLast()) ? dynamic_cast<TGo4FitSlot *>(fxAxisTrans[nslot]) : nullptr;
79}
80
82{
83 TGo4FitSlot *slot = GetAxisTransSlot(nslot);
84 return !slot ? nullptr : dynamic_cast<TGo4FitAxisTrans *>(slot->GetObject());
85}
86
87void TGo4FitData::SetAxisTrans(Int_t nslot, TGo4FitAxisTrans *Trans, Bool_t TransOwned)
88{
89 if (nslot < 0)
90 return;
91 if (nslot >= GetNumberOfTransSlots())
92 SetNumberOfTransSlots(nslot + 1);
93 ((TGo4FitSlot *)(fxAxisTrans[nslot]))->SetObject(Trans, TransOwned);
94}
95
96void TGo4FitData::AddAxisTrans(TGo4FitAxisTrans *Trans, Bool_t TransOwned)
97{
98 Int_t nslot = GetNumberOfTransSlots();
99 SetNumberOfTransSlots(nslot + 1);
100 ((TGo4FitSlot *)(fxAxisTrans[nslot]))->SetObject(Trans, TransOwned);
101}
102
103void TGo4FitData::SetAxisTransNeeded(Int_t nslot, Bool_t iNeeded)
104{
105 if (iNeeded && (nslot >= GetNumberOfTransSlots()))
106 SetNumberOfTransSlots(nslot + 1);
107 if ((nslot >= 0) && (nslot < GetNumberOfTransSlots()))
108 ((TGo4FitSlot *)(fxAxisTrans[nslot]))->SetNeeded(iNeeded);
109}
110
112{
113 if (GetUseBinScale() || (GetExcludeLessThen() > 0))
114 return kTRUE;
115 for (Int_t n = 0; n < GetNumberOfTransSlots(); n++)
116 if (GetAxisTrans(n))
117 return kTRUE;
118 return kFALSE;
119}
120
121TObject *TGo4FitData::CreateDrawObject(const char *ObjName)
122{
123 auto iter = MakeIter();
124 if (!iter)
125 return nullptr;
126 TObject *obj = iter->CreateDrawObject(ObjName);
127 return obj;
128}
129
130Bool_t TGo4FitData::Initialize(Int_t UseBuffers)
131{
132 auto iter = MakeIter();
133 if (!iter)
134 return kFALSE;
135
136 fiBinsSize = iter->CountPoints(kTRUE);
137
138 fiIndexesSize = iter->IndexesSize();
139 fiScalesSize = iter->ScalesSize();
140
141 Bool_t use = ((UseBuffers < 0) && GetUseBuffers()) || (UseBuffers > 0);
142
143 if (use)
144 for (Int_t n = 0; n < GetNumberOfTransSlots(); n++) {
145 TGo4FitAxisTrans *trans = GetAxisTrans(n);
146 if (trans && !trans->IsAllParsFixed()) {
147 use = kFALSE;
148 break;
149 }
150 }
151
152 if (use) {
153
154 fxValues = new Double_t[fiBinsSize];
155 fxStandDev = new Double_t[fiBinsSize];
156 fxBinsResult = new Double_t[fiBinsSize];
157
158 if (iter->HasIndexes())
159 fxFullIndex = new Int_t[fiBinsSize * fiIndexesSize];
160 fxFullScale = new Double_t[fiBinsSize * fiScalesSize];
161 if (iter->HasWidths())
162 fxFullWidth = new Double_t[fiBinsSize * fiScalesSize];
163
164 Int_t nbin = 0;
165 if (iter->Reset())
166 do {
167
168 fxValues[nbin] = iter->Value();
169 fxStandDev[nbin] = iter->StandardDeviation();
170
171 if (fxFullIndex)
172 for (Int_t n = 0; n < fiIndexesSize; n++)
173 fxFullIndex[nbin * fiIndexesSize + n] = iter->Indexes()[n];
174
175 if (fxFullScale)
176 for (Int_t naxis = 0; naxis < fiScalesSize; naxis++)
177 fxFullScale[nbin * fiScalesSize + naxis] = iter->Scales()[naxis];
178
179 if (fxFullWidth && iter->HasWidths())
180 for (Int_t naxis = 0; naxis < fiScalesSize; naxis++)
181 fxFullWidth[nbin * fiScalesSize + naxis] = iter->Widths()[naxis];
182
183 nbin++;
184 } while (iter->Next());
185 }
186
187 return kTRUE;
188}
189
194
196{
197 fiBinsSize = 0;
198 fiIndexesSize = 0;
199 fiScalesSize = 0;
200
201 fxValues = nullptr;
202 fxStandDev = nullptr;
203 fxBinsResult = nullptr;
204
205 fxFullScale = nullptr;
206 fxFullWidth = nullptr;
207
208 fxFullIndex = nullptr;
209}
210
212{
213 if (fxValues)
214 delete[] fxValues;
215
216 if (fxStandDev)
217 delete[] fxStandDev;
218
219 if (fxBinsResult)
220 delete[] fxBinsResult;
221
222 if (fxFullIndex)
223 delete[] fxFullIndex;
224
225 if (fxFullScale)
226 delete[] fxFullScale;
227
228 if (fxFullWidth)
229 delete[] fxFullWidth;
230
232}
233
234Bool_t TGo4FitData::DefineScaleMinMax(Int_t naxis, Double_t &min, Double_t &max)
235{
236 auto iter = MakeIter();
237 if (!iter)
238 return kFALSE;
239 Bool_t res = kFALSE;
240 if (iter->Reset(kFALSE) && (iter->ScalesSize() <= naxis)) {
241 min = iter->Scales()[naxis];
242 max = min;
243 do {
244 Double_t value = iter->Scales()[naxis];
245 if (value < min)
246 min = value;
247 else if (value > max)
248 max = value;
249 } while (iter->Next(kFALSE));
250 res = kTRUE;
251 }
252
253 return res;
254}
255
257{
258 auto iter = MakeIter();
259 if (!iter)
260 return 0;
261 Int_t res = 0;
262 if (iter->Reset(kFALSE))
263 res = iter->IndexesSize();
264 return res;
265}
266
268{
269 auto iter = MakeIter();
270
271 return iter ? iter->CountPoints(kTRUE) : 0;
272}
273
274const Double_t *TGo4FitData::GetScaleValues(Int_t nbin) const
275{
276 if (fxFullScale)
277 return &(fxFullScale[nbin * GetScalesSize()]);
278
279 return nullptr;
280}
281
282const Double_t *TGo4FitData::GetWidthValues(Int_t nbin) const
283{
284 if (fxFullWidth)
285 return &(fxFullWidth[nbin * GetScalesSize()]);
286 return nullptr;
287}
288
289const Int_t *TGo4FitData::GetFullIndex(Int_t nbin) const
290{
291 if (fxFullIndex)
292 return &(fxFullIndex[nbin * GetIndexesSize()]);
293 return nullptr;
294}
295
297{
298 if (!data)
299 return kFALSE;
300 auto iter = data->MakeIter();
301 if (!iter)
302 return kFALSE;
303
304 Bool_t res = kFALSE;
305 if (iter->Reset(kFALSE))
306 res = (iter->IndexesSize() == GetIndexesSize()) && (GetIndexesSize() > 0);
307
308 return res;
309}
310
312{
313 if (!ModelMask)
314 return;
315
316 if (BuffersAllocated()) {
317 for (Int_t nbin = 0; nbin < GetBinsSize(); nbin++) {
318 const Double_t *values = GetScaleValues(nbin);
319
320 Bool_t res = model->CheckRangeConditions(values, GetScalesSize());
321
322 ModelMask[nbin] = res ? 1 : 0;
323 }
324 } else {
325 auto iter = MakeIter();
326 Int_t nbin = 0;
327 if (iter->Reset())
328 do {
329 Bool_t res = model->CheckRangeConditions(iter->Scales(), iter->ScalesSize());
330 ModelMask[nbin] = res ? 1 : 0;
331 nbin++;
332 } while (iter->Next());
333 }
334}
335
336void TGo4FitData::FillSlotList(TSeqCollection *list)
337{
339 for (Int_t n = 0; n <= fxAxisTrans.GetLast(); n++)
340 list->Add(fxAxisTrans[n]);
341}
342
343void TGo4FitData::Print(Option_t *option) const
344{
346 std::cout << " Data type: ";
347 switch (fiDataType) {
348 case 1: std::cout << "histogram" << std::endl; break;
349 case 2: std::cout << "graph" << std::endl; break;
350 default: std::cout << fiDataType << std::endl;
351 }
352 std::cout << " Use bin scale: " << fbUseBinScale << std::endl;
353 std::cout << " Take sigmas from: ";
354 switch (GetSigmaSource()) {
355 case 0: std::cout << "none" << std::endl; break;
356 case 1: std::cout << "data" << std::endl; break;
357 case 2: std::cout << "const value " << GetSigmaValue() << std::endl; break;
358 }
359 std::cout << " Exclude bins less then: " << GetExcludeLessThen() << std::endl;
360 std::cout << " Axis transformation data: " << std::endl;
361 fxAxisTrans.Print(option);
362}
363
364void TGo4FitData::Streamer(TBuffer &b)
365{
366 if (b.IsReading()) {
367
368 TGo4FitData::Class()->ReadBuffer(b, this);
369
370 for (Int_t n = 0; n <= fxAxisTrans.GetLast(); n++) {
372 dc->SetDefaults(this, TGo4FitAxisTrans::Class());
373 }
374
375 } else {
376 TGo4FitData::Class()->WriteBuffer(b, this);
377 }
378}
379
380// *******************************************************************************************
381
387
389
390Bool_t TGo4FitDataIter::ReserveArrays(Int_t NumDimen, Int_t NumOwnAxis, Bool_t HasWidth)
391{
392 TGo4FitData *data = GetData();
393 if (!data)
394 return kFALSE;
395
396 fxIndexes.Set(NumDimen);
397 fxIndexes.Reset(0);
398
399 Int_t size = 0;
400 if (data->GetUseBinScale())
401 size = NumDimen;
402 else
403 size = NumOwnAxis;
404
405 if (size <= 0)
406 return kFALSE;
407
408 fxScales.Set(size);
409 fxScales.Reset(0.);
410 if (HasWidth) {
411 fxWidths.Set(size);
412 fxWidths.Reset(1.);
413 } else
414 fxWidths.Set(0);
415
416 return kTRUE;
417}
418
420{
421 TGo4FitData *data = GetData();
422 for (Int_t nslot = 0; nslot < data->GetNumberOfTransSlots(); nslot++) {
423 TGo4FitAxisTrans *trans = data->GetAxisTrans(nslot);
424 if (trans)
425 trans->Transformation(scales, ScalesSize());
426 }
427}
428
429Bool_t TGo4FitDataIter::ProduceScales(const Int_t *index, const Double_t *ownscales, const Double_t *ownwidths)
430{
431 TGo4FitData *data = GetData();
432 if (!data)
433 return kFALSE;
434
435 if (data->GetUseBinScale() || !ownscales) {
436 if (!index)
437 return kFALSE;
438 Double_t add = (data->GetDataType() == TGo4FitData::dtHistogram) ? .5 : 0.;
439 for (Int_t n = 0; n < fxScales.GetSize(); n++)
440 fxScales[n] = index[n] + add;
441 fxWidths.Reset(1.);
442 } else {
443 for (Int_t n = 0; n < fxScales.GetSize(); n++)
444 fxScales[n] = ownscales[n];
445 if (ownwidths)
446 for (Int_t n = 0; n < fxWidths.GetSize(); n++)
447 fxWidths[n] = ownwidths[n];
448 }
449
450 if (data->GetNumberOfTransSlots() > 0) {
451 if (fxWidths.GetSize() == ScalesSize()) {
452 TArrayD arr1(ScalesSize()), arr2(ScalesSize());
453 for (Int_t n = 0; n < ScalesSize(); n++) {
454 arr1[n] = fxScales[n] - fxWidths[n] / 2.;
455 arr2[n] = fxScales[n] + fxWidths[n] / 2.;
456 }
457 TransformScales(arr1.GetArray());
458 TransformScales(arr2.GetArray());
459 for (Int_t n = 0; n < ScalesSize(); n++)
460 fxWidths[n] = TMath::Abs(arr2[n] - arr1[n]);
461 }
462
463 TransformScales(fxScales.GetArray());
464 }
465
466 return kTRUE;
467}
468
469Bool_t TGo4FitDataIter::NextIndex(TArrayI &Index, TArrayI &Limits)
470{
471 Int_t n = 0;
472 while (n < Index.GetSize()) {
473 Index[n]++;
474 if (Index[n] < Limits[n])
475 return kTRUE;
476 Index[n] = 0;
477 n++;
478 }
479 return kFALSE;
480}
481
483{
484 TGo4FitData *data = GetData();
485 if (!data)
486 return kFALSE;
487 if (data->GetSigmaSource() == 2) {
489 return kTRUE;
490 } else
491 return kFALSE;
492}
493
495{
496 TGo4FitData *data = GetData();
497 if (!data)
498 return kFALSE;
499 if (Value() < data->GetExcludeLessThen())
500 return kFALSE;
501 return data->CheckRangeConditions(Scales(), ScalesSize());
502}
503
504Bool_t TGo4FitDataIter::Reset(Bool_t UseRanges)
505{
506 fbReachEnd = kTRUE;
507
508 if (!StartReset())
509 return kFALSE;
510
511 fiNumPoint = 0;
512
513 if (!ReadCurrentPoint())
514 return kFALSE;
515 if (!UseRanges) {
516 fbReachEnd = kFALSE;
517 return kTRUE;
518 }
519
520 while (!CheckPointForRange()) {
521 if (!ShiftToNextPoint())
522 return kFALSE;
523 if (!ReadCurrentPoint())
524 return kFALSE;
525 }
526
527 fbReachEnd = kFALSE;
528 return kTRUE;
529}
530
531Bool_t TGo4FitDataIter::Next(Bool_t UseRanges)
532{
533 fiNumPoint++;
534
535 if (fbReachEnd || !GetData()) {
536 fbReachEnd = kTRUE;
537 return kFALSE;
538 }
539
540 do {
541 if (!ShiftToNextPoint()) {
542 fbReachEnd = kTRUE;
543 return kFALSE;
544 }
545
546 if (!ReadCurrentPoint()) {
547 fbReachEnd = kTRUE;
548 return kFALSE;
549 }
550
551 if (!UseRanges)
552 return kTRUE;
553
554 } while (!CheckPointForRange());
555
556 return kTRUE;
557}
558
560{
561 double res = 1.;
562 if (HasWidths())
563 for (int n = 0; n < fxWidths.GetSize(); n++)
564 res = res * fxWidths[n];
565 return res;
566}
567
568Int_t TGo4FitDataIter::CountPoints(Bool_t UseRanges)
569{
570 if (!Reset(UseRanges))
571 return 0;
572 Int_t cnt = 0;
573 do {
574 cnt += 1;
575 } while (Next(UseRanges));
576 return cnt;
577}
578
580{
581 if (!Reset(kFALSE))
582 return kFALSE;
583 if (IndexesSize() <= 0)
584 return kFALSE;
585 Limits.Set(IndexesSize());
586 Limits.Reset(0);
587 do {
588 for (Int_t n = 0; n < IndexesSize(); n++)
589 if (Indexes()[n] > Limits[n])
590 Limits[n] = Indexes()[n];
591 } while (Next(kFALSE));
592 return kTRUE;
593}
594
595TH1 *TGo4FitDataIter::CreateHistogram(const char *HistoName, Bool_t UseRanges, Bool_t SetBins)
596{
597 TArrayI Limits;
598 if (!DefineIndexesLimits(Limits))
599 return nullptr;
600 if (!HasIndexes() || (IndexesSize() != ScalesSize()) || !HasWidths())
601 return nullptr;
602
603 Int_t NumDim = IndexesSize();
604 if (NumDim > 3)
605 NumDim = 3;
606
607 Double_t *dummy = nullptr;
608 TH1 *histo = nullptr;
609 switch (NumDim) {
610 case 1: histo = new TH1D(HistoName, "result", Limits[0] + 1, dummy); break;
611 case 2: histo = new TH2D(HistoName, "result", Limits[0] + 1, dummy, Limits[1] + 1, dummy); break;
612 case 3:
613 histo = new TH3D(HistoName, "result", Limits[0] + 1, dummy, Limits[1] + 1, dummy, Limits[2] + 1, dummy);
614 break;
615 default: return nullptr;
616 }
617
618 histo->SetDirectory(nullptr);
619
620 Double_t *Axises[3];
621 for (Int_t n = 0; n < NumDim; n++)
622 Axises[n] = new Double_t[Limits[n] + 2];
623
624 Double_t ampl = GetData()->GetAmplValue();
625
626 if (Reset(UseRanges))
627 do {
628 if (SetBins)
629 switch (NumDim) {
630 case 1: histo->SetBinContent(Indexes()[0] + 1, ampl * Value()); break;
631 case 2: histo->SetBinContent(Indexes()[0] + 1, Indexes()[1] + 1, ampl * Value()); break;
632 case 3: histo->SetBinContent(Indexes()[0] + 1, Indexes()[1] + 1, Indexes()[2] + 1, ampl * Value()); break;
633 }
634 for (Int_t n = 0; n < NumDim; n++) {
635 Int_t indx = Indexes()[n];
636 Axises[n][indx] = Scales()[n] - Widths()[n] / 2.;
637 Axises[n][indx + 1] = Scales()[n] + Widths()[n] / 2.;
638 }
639 } while (Next(UseRanges));
640
641 histo->GetXaxis()->Set(Limits[0] + 1, Axises[0]);
642 if (NumDim > 1)
643 histo->GetYaxis()->Set(Limits[1] + 1, Axises[1]);
644 if (NumDim > 2)
645 histo->GetZaxis()->Set(Limits[2] + 1, Axises[2]);
646
647 for (Int_t n = 0; n < NumDim; n++)
648 delete[] Axises[n];
649
650 return histo;
651}
652
653TGraph *TGo4FitDataIter::CreateGraph(const char *GraphName, Bool_t UseRanges, Bool_t SetBins)
654{
655 Int_t NumPoints = CountPoints(UseRanges);
656 if ((NumPoints <= 0) || (ScalesSize() < 1))
657 return nullptr;
658
659 TGraph *gr = new TGraph(NumPoints);
660 gr->SetName(GraphName);
661 if (Reset(UseRanges))
662 do {
663 (gr->GetX())[Point()] = x();
664 if (SetBins)
665 (gr->GetY())[Point()] = GetData()->GetAmplValue() * Value();
666 } while (Next(UseRanges));
667
668 return gr;
669}
670
671TObject *TGo4FitDataIter::CreateDrawObject(const char *ObjName)
672{
673 if (!Reset(kFALSE))
674 return nullptr;
675
676 if (HasIndexes() && (IndexesSize() == ScalesSize()) && HasWidths())
677 return CreateHistogram(ObjName, kFALSE, kTRUE);
678
679 return CreateGraph(ObjName, kFALSE, kTRUE);
680}
Base class for axis transformation objects.
virtual void Transformation(Double_t *scales, Int_t naxis)=0
Bool_t CheckRangeConditions(const Double_t *values, Int_t numaxis)
Check all range conditions for specified point.
Bool_t GetUseBuffers() const
Returns flag of usage of additional buffers.
TGo4FitParameter * NewAmplitude(const char *Name=nullptr, Double_t iValue=0., Bool_t IsFixed=kFALSE, Int_t AtIndx=0)
Create amplitude parameter with specified properties.
TGo4FitComponent()
Default constructor.
Double_t GetAmplValue()
Return value of amplitude parameter.
void Print(Option_t *option="") const override
Print info about object on standard output.
TGraph * CreateGraph(const char *GraphName, Bool_t UseRanges=kFALSE, Bool_t SetBins=kFALSE)
Create TGraph object with appropriate to data object size.
Double_t xWidths() const
Return production of all width parameters (1 if no widths)
const Double_t * Widths() const
Return scales widths values.
Double_t fdStandardDeviation
virtual ~TGo4FitDataIter()
Destroys TGo4FitDataIter object.
virtual Bool_t Next(Bool_t UseRanges=kTRUE)
Shift pointer to next data point.
virtual Bool_t StartReset()=0
Reset pointer and other specific values to the beginning of data.
Bool_t ProduceScales(const Int_t *index, const Double_t *ownscales, const Double_t *ownwidths)
Converts scale values.
Bool_t GetDeviation()
Calculates standard deviation from GetSigmaValue() of data object.
Bool_t DefineIndexesLimits(TArrayI &Limits)
Iterate over all data points and returns maximum value for indexes.
Bool_t NextIndex(TArrayI &Index, TArrayI &Limits)
Producing next indexes set according limits values.
virtual Bool_t Reset(Bool_t UseRanges=kTRUE)
Initialize iterator and positioning pointer on first point.
Bool_t HasWidths() const
Return kTRUE, if scale widths values exists.
Bool_t HasIndexes() const
Return kTRUE, if data object has indexes.
const Int_t * Indexes() const
Return indexes for current data point.
Int_t CountPoints(Bool_t UseRanges=kTRUE)
Counts total number of points in data object.
Double_t Value() const
Return bin content (Value) for current point.
virtual Bool_t ReadCurrentPoint()=0
Perform specific actions to read all values from data object.
Bool_t CheckPointForRange()
Check range conditions and amplitude threshold for current point.
void TransformScales(Double_t *scales)
Transform scales values, using transformation objects in data.
Double_t x() const
Return current x coordinate if exists, otherwise 0.
TGo4FitDataIter()
Default constructor.
const Double_t * Scales() const
Return scale values for current data points.
TH1 * CreateHistogram(const char *HistoName, Bool_t UseRanges=kFALSE, Bool_t SetBins=kFALSE)
Create histogram (if possible) with appropriate to data object dimensions number and size.
Bool_t ReserveArrays(Int_t NumDimen, Int_t NumOwnAxis, Bool_t HasWidth)
Reserve buffers for indexes, scales and width values.
virtual Bool_t ShiftToNextPoint()=0
Move pointer to following data point.
virtual TGo4FitData * GetData() const =0
Return pointer on correspondent TGo4FitData object, which create iterator.
Int_t ScalesSize() const
Return size (number) of scale values for each data point.
Int_t Point() const
Return number of current point, starting from 0.
Int_t IndexesSize() const
Return size (number of dimensions) of data indexes.
TObject * CreateDrawObject(const char *ObjName)
Create either histogram or graph object.
Basic abstract class for representing data, which should be fitted.
Definition TGo4FitData.h:39
Double_t * fxStandDev
Buffer for standard deviations of bins values.
const Double_t * GetWidthValues(Int_t nbin) const
Return scales width values for specified index from buffer.
const Double_t * GetScaleValues(Int_t nbin) const
Return scale values for specified index from buffer.
Bool_t fbUseBinScale
Use binary numbers as scale values.
Int_t GetNumberOfTransSlots() const
Return number of slots for scale transformation objects.
TGo4FitAxisTrans * GetAxisTrans(Int_t nslot) const
Return transformation object for given slot.
Int_t DefineDimensions()
Define dimension number of data Create iterator and checks number of dimension.
Double_t * fxBinsResult
Buffer for complete model of bins values.
TObject * CreateDrawObject(const char *ObjName)
Creates object, which can be drawn on canvas by ROOT.
Int_t * fxFullIndex
Store combination of indexes for each data bins.
Double_t fdSigmaValue
Value of sigma when fiTakeSigmasFrom = 2.
Bool_t BuffersAllocated() const
Checks, if buffers allocated for data.
Int_t GetBinsSize() const
Return number of data bins in buffers.
Int_t fiDataType
Specified type of data: 0 - Histogram, 1 - graphics, 2 and so on - user defined.
void SetAxisTrans(Int_t nslot, TGo4FitAxisTrans *Trans, Bool_t TransOwned=kFALSE)
Sets transformation object for given slot.
void SetAxisTransNeeded(Int_t nslot, Bool_t iNeeded=kFALSE)
Specified, when iNeeded = kTRUE, that transformation object should always be provided to data.
Int_t fiTakeSigmasFrom
Specify sigma source.
Double_t * fxFullScale
Array of axis values for each bins.
void ResetAllPoinetrs()
Clears (sets to 0) all pointers, used for buffer allocations.
Int_t GetSigmaSource() const
Return source of sigma values.
Definition TGo4FitData.h:90
Int_t GetIndexesSize() const
Returns dimension of indexes arrays.
Double_t * fxValues
Buffer for bins values.
const Int_t * GetFullIndex(Int_t nbin) const
Return indexes values for specified data bin from buffer.
virtual ~TGo4FitData()
Destroys TGo4FitData object.
Int_t fiScalesSize
Number of scales values for each point.
Int_t GetDataType() const
Returns type of data source.
Definition TGo4FitData.h:64
Int_t fiIndexesSize
Number of indexes for each point.
Bool_t SetNumberOfTransSlots(Int_t nslots)
Sets number of slots for scale transformation objects.
Bool_t DefineScaleMinMax(Int_t naxis, Double_t &min, Double_t &max)
Return scales minimum and maximum for specified axis.
Int_t fiBinsSize
Number of entries in buffers.
void FillSlotList(TSeqCollection *list) override
Copy pointers on all slots of data object to list.
TGo4FitData()
Default constructor.
virtual void Finalize()
Removes all buffers, created in initialize() routine.
virtual Bool_t Initialize(Int_t UseBuffers=-1)
Initialize data object.
Bool_t GetUseBinScale() const
Returns kTRUE, if binary numbers used as scale values.
Definition TGo4FitData.h:69
void AddAxisTrans(TGo4FitAxisTrans *Trans, Bool_t TransOwned=kFALSE)
Add transformation object to data.
Double_t GetSigmaValue() const
Return constant sigma value.
Definition TGo4FitData.h:96
Double_t GetExcludeLessThen() const
Returns limit, which uses to exclude bins, less then this limit.
virtual Bool_t IsAnyDataTransform() const
Return kTRUE, if either initial data axis or data bins are transformed by TGo4FitData object.
Double_t fdExcludeLessThen
Sets limit for exclude bins, which less then this limit.
void ReleaseAllPointers()
Release all memory, allocated for buffers.
virtual std::unique_ptr< TGo4FitDataIter > MakeIter()
Creates iterator for data object.
Double_t * fxFullWidth
Array of width values for each bin.
void Print(Option_t *option="") const override
Display information about data object on standard output.
Int_t DefineBinsSize()
Defines number of selected bins in data Creates iterator and sequentially checks all points.
TGo4FitSlot * GetAxisTransSlot(Int_t nslot) const
Return slot for transformation object.
TObjArray fxAxisTrans
Array of slots for scale transformation objects.
Int_t GetScalesSize() const
Returns number of axis values for each point.
void ApplyRangesForModelMask(TGo4FitComponent *model, Char_t *ModelMask)
Exclude points from model according model range conditions.
Bool_t IsCompatibleData(TGo4FitData *data)
Checks, if data has same dimensions number and size of each dimensions.
Bool_t IsAllParsFixed()
Returns true, if all parameters in list fixed;.
virtual void FillSlotList(TSeqCollection *lst)
Fill list of TGo4FitSlot objects to provided collection.
void SetUpdateSlotList()
Update internal list of slots (if exists).
Managing pointers on specific objects.
Definition TGo4FitSlot.h:28
void SetDefaults(TNamed *iOwner, TClass *iClass)
Set basic fields of slot.
TObject * GetObject() const
Return pointer on assigned object.