diff --git a/ChangeLog b/ChangeLog index 91d71410..e7e02cc6 100644 --- a/ChangeLog +++ b/ChangeLog @@ -53,6 +53,9 @@ FIXED 2012-09-23 fixed wrong chisq output in musrview if expected chisq is present. FIXED 2012-05-30 fixed RRF bug in single histo plotting. FIXED 2012-05-18 fixed wrong forward/backward tag for ROOT-PPC (MUSR-215) +CHANGED 2013-11-12 changed normalization of log max likelihood according to S. +Backer and R.D. Cousins NIM 221, 437 (1984), in order to have a +"goodness-of-fit" criteria. CHANGED 2012-12-11 if multiple SAVE are present in the COMMAND block, append MINUIT2.OUTPUT file. Added docu for PRINT_LEVEL (MUSR-244). CHANGED 2012-11-19 replaced hard coded gyromagnetic ratio of the muon in the diff --git a/src/classes/PRunMuMinus.cpp b/src/classes/PRunMuMinus.cpp index 4e4f95fa..e61d5d55 100644 --- a/src/classes/PRunMuMinus.cpp +++ b/src/classes/PRunMuMinus.cpp @@ -261,6 +261,11 @@ Double_t PRunMuMinus::CalcMaxLikelihood(const std::vector& par) data = fData.GetValue()->at(i); + if (theo <= 0.0) { + cerr << ">> PRunMuMinus::CalcMaxLikelihood: **WARNING** NEGATIVE theory!!" << endl; + continue; + } + if (data > 1.0e-9) { mllh += (theo-data) + data*log(data/theo); } else { @@ -268,7 +273,7 @@ Double_t PRunMuMinus::CalcMaxLikelihood(const std::vector& par) } } - return mllh; + return 2.0*mllh; } //-------------------------------------------------------------------------- diff --git a/src/classes/PRunSingleHisto.cpp b/src/classes/PRunSingleHisto.cpp index 7dc81d05..b18092f3 100644 --- a/src/classes/PRunSingleHisto.cpp +++ b/src/classes/PRunSingleHisto.cpp @@ -376,14 +376,19 @@ Double_t PRunSingleHisto::CalcMaxLikelihood(const std::vector& par) chunk = 10; #pragma omp parallel for default(shared) private(i,time,theo,data) schedule(dynamic,chunk) reduction(-:mllh) #endif - for (i=startTimeBin; i < endTimeBin; ++i) { + for (i=startTimeBin; iFunc(time, par, fFuncValues))+bkg; + theo = N0*TMath::Exp(-time/tau)*(1.0+fTheory->Func(time, par, fFuncValues))+bkg; theo *= normalizer; data = normalizer*fData.GetValue()->at(i); + if (theo <= 0.0) { + cerr << ">> PRunSingleHisto::CalcMaxLikelihood: **WARNING** NEGATIVE theory!!" << endl; + continue; + } + if (data > 1.0e-9) { mllh += (theo-data) + data*log(data/theo); } else { @@ -391,7 +396,7 @@ Double_t PRunSingleHisto::CalcMaxLikelihood(const std::vector& par) } } - return mllh; + return 2.0*mllh; } //--------------------------------------------------------------------------