219 double ampMaxAlphaBeta = 0;
220 double tMaxAlphaBeta = 5;
221 double tMaxErrorAlphaBeta = 999;
223 double tMaxErrorRatio = 999;
244 NullChi2 = (sumAA - sumA * sumA / sum1) / sum0;
256 double alphabeta = amplitudeFitParameters[0] * amplitudeFitParameters[1];
257 double invalphabeta = 1. / alphabeta;
258 double alpha = amplitudeFitParameters[0];
259 double beta = amplitudeFitParameters[1];
261 Ratio ratios_[C::MAXSAMPLES * (C::MAXSAMPLES - 1) / 2];
262 unsigned int ratios_size = 0;
289 double err1 = Rtmp * Rtmp * (relErr2[
i] + relErr2[
j]);
303 double err3 = (0.289 * 0.289) * (invampl[
j] * invampl[
j]);
305 double totalError =
std::sqrt(err1 + err2 + err3);
308 if (totalError < 1.0 && Rtmp > 0.001 && Rtmp < Rlim[
j -
i]) {
309 Ratio currentRatio = {
i, (
j -
i), Rtmp, totalError};
310 ratios_[ratios_size++] = currentRatio;
317 if (0 == ratios_size)
321 Tmax timesAB_[C::MAXSAMPLES * (C::MAXSAMPLES - 1) / 2];
322 unsigned int timesAB_size = 0;
327 for (
unsigned int i = 0;
i < ratios_size;
i++) {
328 double stepOverBeta = double(ratios_[
i].
step) /
beta;
331 double Rmin = ratios_[
i].value - ratios_[
i].error;
335 double Rmax = ratios_[
i].value + ratios_[
i].error;
336 double RLimit = Rlim[ratios_[
i].step];
346 double tmax = 0.5 * (time1 + time2);
347 double tmaxerr = 0.5 *
std::sqrt((time1 - time2) * (time1 - time2));
352 int itmin =
std::max(-1,
int(std::floor(
tmax - alphabeta)));
353 double loffset = (double(itmin) -
tmax) * invalphabeta;
354 for (
unsigned int it = itmin + 1; it <
amplitudes_.size(); it++) {
355 loffset += invalphabeta;
360 double term1 = 1.0 + loffset;
364 sumff +=
f * (
f * inverr2);
370 chi2 = sumAA - sumAf * (sumAf / sumff);
371 amp = (sumAf / sumff);
377 if (
chi2 > 0 && tmaxerr > 0 &&
tmax > 0) {
378 Tmax currentTmax = {ratios_[
i].index, ratios_[
i].step,
tmax, tmaxerr, amp,
chi2};
379 timesAB_[timesAB_size++] = currentTmax;
384 if (0 == timesAB_size)
388 double chi2min = 1.0e+9;
391 for (
unsigned int i = 0;
i < timesAB_size;
i++) {
392 if (timesAB_[
i].
chi2 < chi2min) {
393 chi2min = timesAB_[
i].chi2;
401 double chi2Limit = chi2min + 1.0;
404 for (
unsigned int i = 0;
i < timesAB_size;
i++) {
405 if (timesAB_[
i].
chi2 < chi2Limit) {
406 double inverseSigmaSquared = 1.0 / (timesAB_[
i].error * timesAB_[
i].error);
407 time_wgt += inverseSigmaSquared;
408 time_max += timesAB_[
i].value * inverseSigmaSquared;
412 tMaxAlphaBeta = time_max / time_wgt;
413 tMaxErrorAlphaBeta = 1.0 /
sqrt(time_wgt);
422 double offset = (double(
i) - tMaxAlphaBeta) * invalphabeta;
423 double term1 = 1.0 +
offset;
427 sumff +=
f * (
f * inverr2);
432 ampMaxAlphaBeta = sumAf / sumff;
433 double chi2AlphaBeta = (sumAA - sumAf * sumAf / sumff) / sum0;
434 if (chi2AlphaBeta > NullChi2) {
461 for (
unsigned int i = 0;
i < ratios_size;
i++) {
462 if (ratios_[
i].
step == 1 && ratios_[
i].
value >= timeFitLimits.first && ratios_[
i].value <= timeFitLimits.second) {
463 double time_max_i = ratios_[
i].index;
467 double u = timeFitParameters[timeFitParameters.size() - 1];
468 for (
int k = timeFitParameters.size() - 2;
k >= 0;
k--) {
469 u = u * ratios_[
i].value + timeFitParameters[
k];
473 double du = (timeFitParameters.size() - 1) * timeFitParameters[timeFitParameters.size() - 1];
474 for (
int k = timeFitParameters.size() - 2;
k >= 1;
k--) {
475 du = du * ratios_[
i].value +
k * timeFitParameters[
k];
479 double errorsquared = ratios_[
i].error * ratios_[
i].error * du * du;
480 if (errorsquared > 0) {
481 time_max += (time_max_i - u) / errorsquared;
482 time_wgt += 1.0 / errorsquared;
494 tMaxErrorRatio = 1.0 /
sqrt(time_wgt);
509 tMaxErrorRatio * (ampMaxAlphaBeta /
ampMaxError_ - 5.0)) /