00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016 #include "TGo4FitModel.h"
00017
00018 #include <iostream.h>
00019
00020 #include "TArrayD.h"
00021 #include "TArrayI.h"
00022
00023 TGo4FitAssignment::TGo4FitAssignment() :
00024 TNamed(), fxRatio(0), fxData(0), fxModelMask(0), fxModelBins(0) {
00025 }
00026
00027 TGo4FitAssignment::TGo4FitAssignment(const char* DataName) :
00028 TNamed(DataName,""),
00029 fxRatio(0), fxData(0), fxModelMask(0), fxModelBins(0) {
00030 }
00031
00032 TGo4FitAssignment::~TGo4FitAssignment() {
00033 if (fxRatio) delete fxRatio;
00034 if (fxModelMask) delete[] fxModelMask;
00035 if (fxModelBins) delete[] fxModelBins;
00036 }
00037
00038 Double_t TGo4FitAssignment::RatioValue() {
00039 return (fxRatio==0) ? 1. : fxRatio->GetValue();
00040 }
00041
00042 void TGo4FitAssignment::Print(Option_t* option) const {
00043 cout << " " << GetName();
00044 }
00045
00046 ClassImp(TGo4FitAssignment)
00047
00048
00049
00050
00051
00052 TGo4FitModel::TGo4FitModel() : TGo4FitComponent(),
00053 fiMinIntegrDepth(0), fiMaxIntegrDepth(0), fdIntegrEps(0.), fbAbsoluteEps(kFALSE), fbIntegrScaling(kFALSE),
00054 fxAssigments(), fiGroupIndex(-1),
00055 fxCurrentPars(), fxCurrentParsArray(0),
00056 fbNeedToRebuild(kFALSE), fxAllPars(0), fxAllParsValues(0) {
00057 }
00058
00059 TGo4FitModel::TGo4FitModel(const char* iName, const char* iTitle, Bool_t MakeAmplitude) :
00060 TGo4FitComponent(iName,iTitle),
00061 fiMinIntegrDepth(0), fiMaxIntegrDepth(0), fdIntegrEps(0.), fbAbsoluteEps(kFALSE), fbIntegrScaling(kFALSE),
00062 fxAssigments(), fiGroupIndex(-1),
00063 fxCurrentPars(), fxCurrentParsArray(0),
00064 fbNeedToRebuild(kFALSE), fxAllPars(0), fxAllParsValues(0) {
00065 fxAssigments.SetOwner(kTRUE);
00066 if (MakeAmplitude) NewAmplitude();
00067 }
00068
00069 TGo4FitModel::~TGo4FitModel() {
00070 RemoveAllPars();
00071 }
00072
00073 Int_t TGo4FitModel::NumPars() {
00074 Int_t res = TGo4FitParsList::NumPars();
00075 if (NumAssigments()>0) res += NumAssigments()-1;
00076 return res;
00077 }
00078
00079 TGo4FitParameter* TGo4FitModel::Get(Int_t n) {
00080 Int_t numpars = TGo4FitParsList::NumPars();
00081 if (n<numpars) return TGo4FitParsList::Get(n);
00082 n-=numpars-1;
00083 return (n>=0) && (n<NumAssigments()) ? GetAssigment(n)->fxRatio : 0;
00084 }
00085
00086 void TGo4FitModel::AssignToData(const char* DataName, Double_t RatioValue, Bool_t FixRatio) {
00087 if (FindAssigment(DataName)==0) {
00088 TGo4FitAssignment* ass = new TGo4FitAssignment(DataName);
00089 if (NumAssigments()>0) {
00090 ass->fxRatio = new TGo4FitParameter(GetRatioName(NumAssigments()),
00091 "Amplitude ratio to first data",
00092 RatioValue);
00093 ass->fxRatio->SetOwner(this);
00094 if (FixRatio) ass->fxRatio->SetFixed(kTRUE);
00095 }
00096 fxAssigments.Add(ass);
00097 }
00098 }
00099
00100 void TGo4FitModel::ChangeDataNameInAssignments(const char* oldname, const char* newname) {
00101 TGo4FitAssignment* ass = FindAssigment(oldname);
00102 if (ass) ass->SetName(newname);
00103 }
00104
00105 void TGo4FitModel::ClearAssignmentTo(const char* DataName) {
00106 TGo4FitAssignment* ass = FindAssigment(DataName);
00107 if (ass==0) return;
00108
00109 fxAssigments.Remove(ass);
00110 delete ass;
00111 fxAssigments.Compress();
00112
00113 for (Int_t n=0;n<NumAssigments();n++) {
00114 ass = GetAssigment(n);
00115 if (ass->fxRatio==0) continue;
00116 if (n==0) {
00117 delete ass->fxRatio;
00118 ass->fxRatio = 0;
00119 } else
00120 ass->fxRatio->SetName(GetRatioName(n+1));
00121 }
00122 }
00123
00124 void TGo4FitModel::ClearAssignments() {
00125 fxAssigments.Clear();
00126 fxAssigments.Compress();
00127 }
00128
00129 void TGo4FitModel::ConnectToDataIfAssigned(TGo4FitData* data) {
00130 if (data==0) return;
00131 TGo4FitAssignment* ass = FindAssigment(data->GetName());
00132 if ( ass != 0 ) ass->fxData = data;
00133 }
00134
00135 TGo4FitData* TGo4FitModel::GetAssignedConnection(Int_t n) {
00136 TGo4FitAssignment* ass = GetAssigment(n);
00137 if(ass) return ass->fxData;
00138 else return 0;
00139 }
00140
00141 Bool_t TGo4FitModel::GetPosition(Int_t naxis, Double_t& pos) {
00142 if(GetPosPar(naxis)) {
00143 pos = GetPosPar(naxis)->GetValue();
00144 return kTRUE;
00145 } else return kFALSE;
00146 }
00147
00148 Bool_t TGo4FitModel::SetPosition(Int_t naxis, Double_t pos) {
00149 if(GetPosPar(naxis)) {
00150 GetPosPar(naxis)->SetValue(pos);
00151 return kTRUE;
00152 } else return kFALSE;
00153 }
00154
00155 Bool_t TGo4FitModel::GetWidth(Int_t naxis, Double_t& width) {
00156 if(GetWidthPar(naxis)) {
00157 width = GetWidthPar(naxis)->GetValue();
00158 return kTRUE;
00159 } else return kFALSE;
00160 }
00161
00162 Bool_t TGo4FitModel::SetWidth(Int_t naxis, Double_t width) {
00163 if(GetWidthPar(naxis)) {
00164 GetWidthPar(naxis)->SetValue(width);
00165 return kTRUE;
00166 } else return kFALSE;
00167 }
00168
00169 void TGo4FitModel::SetIntegrationsProperty(Int_t iMinIntegrDepth, Int_t iMaxIntegrDepth, Double_t iIntegrEps, Bool_t iAbsoluteEps, Bool_t iIntegrScaling) {
00170 fiMinIntegrDepth = iMinIntegrDepth;
00171 fiMaxIntegrDepth = iMaxIntegrDepth;
00172 fdIntegrEps = iIntegrEps;
00173 fbAbsoluteEps = iAbsoluteEps;
00174 fbIntegrScaling = iIntegrScaling;
00175 if (fdIntegrEps<=0.)
00176 fiMaxIntegrDepth = fiMinIntegrDepth;
00177 }
00178
00179 void TGo4FitModel::RemoveAllPars() {
00180 if (fxAllPars) { delete fxAllPars; fxAllPars=0; }
00181 if (fxAllParsValues) { delete fxAllParsValues; fxAllParsValues=0; }
00182 }
00183
00184 Bool_t TGo4FitModel::Initialize(Int_t UseBuffers) {
00185 Bool_t use = ((UseBuffers<0) && GetUseBuffers()) || (UseBuffers>0);
00186 for(Int_t n=0;n<NumAssigments();n++) {
00187 TGo4FitAssignment* ass = GetAssigment(n);
00188 if (ass->fxData==0) {
00189 cout << "Data " << ass->GetName() << " not assigned to model " << GetName() << endl;
00190 continue;
00191 }
00192
00193 if (use && ass->fxData->BuffersAllocated()) {
00194 if (IsAnyRangeLimits()) {
00195 ass->fxModelMask = new Char_t[ass->fxData->GetBinsSize()];
00196 ass->fxData->ApplyRangesForModelMask(this,ass->fxModelMask);
00197 }
00198
00199 ass->fxModelBins = new Double_t[ass->fxData->GetBinsSize()];
00200 }
00201 }
00202
00203 SetNeedToRebuild();
00204
00205 if (use) {
00206 RemoveAllPars();
00207
00208 fxAllPars = new TGo4FitParsList(kFALSE);
00209 CollectParsTo(*fxAllPars);
00210
00211 fxAllParsValues = new TArrayD(fxAllPars->NumPars());
00212 fxAllPars->GetParsValues(fxAllParsValues->GetArray());
00213 }
00214
00215 return kTRUE;
00216 }
00217
00218 void TGo4FitModel::Finalize() {
00219 RemoveAllPars();
00220 fxCurrentPars.Set(0);
00221 for(Int_t n=0;n<NumAssigments();n++) {
00222 TGo4FitAssignment* ass = GetAssigment(n);
00223 if (ass->fxModelMask) { delete[] ass->fxModelMask; ass->fxModelMask = 0; }
00224 if (ass->fxModelBins) { delete[] ass->fxModelBins; ass->fxModelBins = 0; }
00225 }
00226 }
00227
00228 Bool_t TGo4FitModel::BuffersAllocated() const {
00229 Bool_t res = (NumAssigments()>0);
00230 for(Int_t n=0;n<NumAssigments();n++)
00231 res = res && (GetAssigment(n)->fxModelBins!=0);
00232 return res;
00233 }
00234
00235 Double_t* TGo4FitModel::GetModelBins(const char* DataName) const {
00236 TGo4FitAssignment* ass = FindAssigment(DataName);
00237 return ass != 0 ? ass->fxModelBins : 0;
00238 }
00239
00240 Double_t TGo4FitModel::EvaluateAndIntegrate(Int_t NumScales, const Double_t* Scales, const Double_t* Widths) {
00241 if ((NumScales<1) || (Scales==0)) return 0.;
00242
00243 if ( (fiMinIntegrDepth>0) && (fiMaxIntegrDepth>0) && (Widths!=0)) {
00244 TArrayI IntegrIndexes(NumScales);
00245 TArrayD ScaleValues(NumScales, Scales);
00246 TArrayD WidthValues(NumScales, Widths);
00247 TArrayD dScaleValues(NumScales);
00248
00249 Int_t* dindx = IntegrIndexes.GetArray();
00250 Double_t* vector = ScaleValues.GetArray();
00251 Double_t* width = WidthValues.GetArray();
00252 Double_t* dvector = dScaleValues.GetArray();
00253
00254 Int_t n=0;
00255
00256
00257 for(n=0;n<NumScales;n++)
00258 vector[n] -= 0.5*width[n];
00259
00260
00261 Double_t TotalSum = EvalN(vector);
00262 Int_t TotalNumPnt = 1;
00263
00264 Int_t IntegrDepth = 0;
00265 Int_t NumStep = 1;
00266 Double_t Step = 1.0;
00267 Bool_t stopcondition = kFALSE;
00268
00269
00270 do {
00271
00272 IntegrIndexes.Reset(0);
00273 Double_t Sum = 0.;
00274 Int_t NumPnt = 0;
00275
00276
00277 do {
00278 for(n=0;n<NumScales;n++)
00279 dvector[n] = vector[n]+Step*width[n]*(dindx[n]+0.5);
00280
00281 Sum+=EvalN(dvector);
00282 NumPnt++;
00283
00284 n=0;
00285 do {
00286 dindx[n]++;
00287 if (dindx[n]<NumStep) break;
00288 dindx[n]=0; n++;
00289 } while (n<NumScales);
00290 } while (n<NumScales);
00291
00292
00293 stopcondition = kFALSE;
00294 if (IntegrDepth>=fiMinIntegrDepth)
00295 if (IntegrDepth>=fiMaxIntegrDepth) stopcondition = kTRUE; else {
00296 Double_t v1 = TotalSum/TotalNumPnt;
00297 Double_t v2 = Sum/NumPnt;
00298 Double_t v = TMath::Abs(v1-v2);
00299 if (!fbAbsoluteEps) {
00300 if ((v1!=0) || (v2!=0)) v=v/(TMath::Abs(v1)+TMath::Abs(v2));
00301 else v=0.;
00302 }
00303 if (v<=fdIntegrEps) stopcondition = kTRUE;
00304 }
00305
00306 TotalSum+=Sum;
00307 TotalNumPnt+=NumPnt;
00308 IntegrDepth++;
00309 NumStep*=2;
00310 Step = Step/2.;
00311 } while (!stopcondition);
00312
00313 Double_t value = 0;
00314 if (TotalNumPnt>0) value = TotalSum / TotalNumPnt;
00315 if (fbIntegrScaling)
00316 for(n=0;n<NumScales;n++) value*=width[n];
00317 return value;
00318 } else return EvalN(Scales);
00319 }
00320
00321 Double_t TGo4FitModel::EvaluateAtPoint(TGo4FitData* data, Int_t nbin, Bool_t UseRanges) {
00322 if (data==0) return 0.;
00323 const Double_t* scales = data->GetScaleValues(nbin);
00324 Int_t scalessize = data->GetScalesSize();
00325 if (UseRanges && !CheckRangeConditions(scales, scalessize)) return 0.;
00326 return EvaluateAndIntegrate(scalessize, scales, data->GetWidthValues(nbin));
00327 }
00328
00329 Double_t TGo4FitModel::EvaluateAtPoint(TGo4FitDataIter* iter, Bool_t UseRanges) {
00330 if (iter==0) return 0;
00331 const Double_t* scales = iter->Scales();
00332 Int_t scalessize = iter->ScalesSize();
00333 if (UseRanges && !CheckRangeConditions(scales, scalessize)) return 0.;
00334 return EvaluateAndIntegrate(scalessize, scales, iter->Widths());
00335 }
00336
00337 void TGo4FitModel::RebuildShape(Bool_t ForceBuild) {
00338 if ((fxAllPars==0) || (fxAllParsValues==0)) return;
00339
00340 Bool_t fill = ForceBuild || fbNeedToRebuild;
00341 if (!fill)
00342 for(Int_t n=0;n<fxAllPars->NumPars();n++) {
00343 if (n==GetAmplIndex()) continue;
00344 if ((*fxAllParsValues)[n] != fxAllPars->GetPar(n)->GetValue()) {
00345 fill = kTRUE;
00346 break;
00347 }
00348 }
00349
00350 if(fill)
00351 for(Int_t n=0;n<NumAssigments();n++) {
00352 TGo4FitAssignment* ass = GetAssigment(n);
00353
00354 if ((ass->fxData==0) || (ass->fxModelBins==0)) continue;
00355
00356 TGo4FitData* data = ass->fxData;
00357
00358 if (!BeforeEval(data->GetScalesSize())) continue;
00359
00360 Int_t size = data->GetBinsSize();
00361 Double_t* bins = ass->fxModelBins;
00362 Char_t* mask = ass->fxModelMask;
00363 Double_t ratio = ass->RatioValue();
00364
00365 for(Int_t nbin=0;nbin<size;nbin++) {
00366 if ((mask!=0) && (mask[nbin]==0)) continue;
00367 bins[nbin] = ratio * EvaluateAtPoint(data, nbin, kFALSE);
00368 }
00369
00370 AfterEval();
00371 }
00372
00373
00374 for(Int_t n=0;n<fxAllPars->NumPars();n++)
00375 (*fxAllParsValues)[n] = fxAllPars->GetPar(n)->GetValue();
00376
00377 fbNeedToRebuild = kFALSE;
00378 }
00379
00380 Bool_t TGo4FitModel::AddModelToDataResult(TGo4FitData* data) {
00381 if ((data==0) || (!data->BuffersAllocated())) return kFALSE;
00382
00383 Double_t* result = data->GetBinsResult();
00384 if (result==0) return kFALSE;
00385 Int_t size = data->GetBinsSize();
00386
00387 Double_t* modelbins = GetModelBins(data->GetName());
00388 if (modelbins) {
00389 Double_t ampl = GetAmplValue();
00390 for (Int_t nbin=0;nbin<size;nbin++)
00391 result[nbin] += ampl * modelbins[nbin];
00392 return kTRUE;
00393 }
00394
00395 if (!BeforeEval(data->GetScalesSize())) return kFALSE;
00396
00397 Double_t ampl = GetAmplValue() * GetRatioValueFor(data->GetName());
00398
00399 for(Int_t nbin=0;nbin<size;nbin++)
00400 result[nbin] += ampl * EvaluateAtPoint(data, nbin);
00401
00402 AfterEval();
00403
00404 return kTRUE;
00405 }
00406
00407 Bool_t TGo4FitModel::BeforeEval(Int_t ndim) {
00408 fxCurrentPars.Set(NumPars());
00409 GetParsValues(fxCurrentPars.GetArray());
00410 fxCurrentParsArray = fxCurrentPars.GetArray();
00411 if (GetAmplPar()!=0) fxCurrentParsArray++;
00412 return kTRUE;
00413 }
00414
00415 Double_t TGo4FitModel::EvalN(const Double_t* v) {
00416 return UserFunction((Double_t*)v, fxCurrentParsArray);
00417 }
00418
00419 Double_t TGo4FitModel::Evaluate(Double_t x) {
00420 Double_t zn = 0.;
00421 if (BeforeEval(1)) {
00422 zn = GetAmplValue() * EvalN(&x);
00423 AfterEval();
00424 }
00425 return zn;
00426 }
00427
00428 Double_t TGo4FitModel::Evaluate(Double_t x, Double_t y) {
00429 Double_t zn = 0.;
00430 if (BeforeEval(2)) {
00431 Double_t vector[2] = {x,y};
00432 zn = GetAmplValue() * EvalN(vector);
00433 AfterEval();
00434 }
00435 return zn;
00436 }
00437
00438 Double_t TGo4FitModel::Evaluate(Double_t x, Double_t y, Double_t z) {
00439 Double_t zn = 0.;
00440 if (BeforeEval(3)) {
00441 Double_t vector[3] = {x,y,z};
00442 zn = GetAmplValue() * EvalN(vector);
00443 AfterEval();
00444 }
00445 return zn;}
00446
00447 Double_t TGo4FitModel::Evaluate(Double_t* v, Int_t ndim) {
00448 Double_t zn = 0.;
00449 if (BeforeEval(ndim)) {
00450 zn = GetAmplValue() * EvalN(v);
00451 AfterEval();
00452 }
00453 return zn;
00454 }
00455
00456 Double_t TGo4FitModel::Integral() {
00457 return 0;
00458 }
00459
00460
00461 const Int_t* TGo4FitModel::GetDataFullIndex(TGo4FitData* data, Int_t nbin) {
00462 return data==0 ? 0 : data->GetFullIndex(nbin);
00463 }
00464
00465 Int_t TGo4FitModel::GetDataIndexesSize(TGo4FitData* data) {
00466 return data==0 ? 0 : data->GetIndexesSize();
00467 }
00468
00469 TGo4FitAssignment* TGo4FitModel::FindAssigment(const char* DataName) const {
00470 for (Int_t n=0;n<NumAssigments();n++)
00471 if (strcmp(DataName, GetAssigment(n)->GetName()) == 0)
00472 return GetAssigment(n);
00473 return 0;
00474 }
00475
00476 TString TGo4FitModel::GetRatioName(Int_t n) {
00477 char str[20];
00478 snprintf(str,19,"Ratio%d",n);
00479 return TString(str);
00480 }
00481
00482 Double_t TGo4FitModel::GetRatioValueFor(const char* DataName) {
00483 TGo4FitAssignment* ass = FindAssigment(DataName);
00484 return (ass==0) ? 1. : ass->RatioValue();
00485 }
00486
00487 void TGo4FitModel::Print(Option_t* option) const {
00488 TGo4FitComponent::Print(option);
00489 cout << " Assigned to: ";
00490 fxAssigments.Print(option);
00491 cout << endl;
00492 if ( (fiMinIntegrDepth>0) && (fiMaxIntegrDepth>0) ) {
00493 cout << " Integration property: depths from " << fiMinIntegrDepth << " to " << fiMaxIntegrDepth;
00494 if (fbAbsoluteEps) cout << " absolute"; else cout << " relative";
00495 cout << " error " << fdIntegrEps << endl;
00496 }
00497 }
00498
00499 void TGo4FitModel::Streamer(TBuffer& b) {
00500 if (b.IsReading()) {
00501
00502 TGo4FitModel::Class()->ReadBuffer(b, this);
00503
00504 for (Int_t n=0;n<NumAssigments();n++)
00505 if (GetAssigment(n)->fxRatio)
00506 GetAssigment(n)->fxRatio->SetOwner(this);
00507
00508 } else {
00509 TGo4FitModel::Class()->WriteBuffer(b, this);
00510 }
00511 }
00512
00513 ClassImp(TGo4FitModel)
00514
00515
