1
0
mirror of https://github.com/aportelli/LatAnalyze.git synced 2024-11-10 00:45:36 +00:00

Fit: various fixes in indexing for error-on-variable chi^2

This commit is contained in:
Antonin Portelli 2015-06-26 20:34:48 +01:00
parent 1adf182caf
commit 57d4bfb4f8
2 changed files with 23 additions and 18 deletions

View File

@ -253,25 +253,25 @@ double Chi2Function::operator()(const double *arg) const
Index is;
ConstMap<DVec> xi(arg + nPar_, getNArg() - nPar_, 1);
double res;
// initialize buffers if necessary
if (!buffer_->isInit)
{
initBuffer();
}
// set the upper part of v: f_j(xi, p) - y_j
is = 0;
for (Index j = 0; j < yDim; ++j)
FOR_VEC(buffer_->dInd, k)
{
const DoubleModel *f = model_[static_cast<size_type>(j)];
double f_jk, y_jk = data_.y(j, buffer_->dInd(k));
if (!f)
{
LATAN_ERROR(Memory, "null model");
}
is = 0;
for (Index i = 0; i < xDim; ++i)
{
if (data_.isXExact(i))
@ -288,7 +288,7 @@ double Chi2Function::operator()(const double *arg) const
f_jk = (*f)(buffer_->x.data(), arg);
buffer_->v(j*nPoint + k) = f_jk - y_jk;
}
// set the down part of v: xi_ik - x_ik
FOR_VEC(buffer_->xInd, i)
FOR_VEC(buffer_->dInd, k)
@ -296,9 +296,9 @@ double Chi2Function::operator()(const double *arg) const
double x_ik = data_.x(buffer_->xInd(i), buffer_->dInd(k));
double xi_ik = xi(i*nPoint + k);
buffer_->v(i*nPoint + k) = xi_ik - x_ik;
buffer_->v(yDim*nPoint + i*nPoint + k) = xi_ik - x_ik;
}
// compute result
res = buffer_->v.dot(buffer_->invVar*buffer_->v);
@ -314,11 +314,11 @@ DoubleFunction Chi2Function::makeFunction(const bool makeHardCopy) const
{
Chi2Function copy(*this);
res.setFunction([copy](const double *p){return copy(p);}, getNPar());
res.setFunction([copy](const double *p){return copy(p);}, getNArg());
}
else
{
res.setFunction([this](const double *p){return (*this)(p);}, getNPar());
res.setFunction([this](const double *p){return (*this)(p);}, getNArg());
}
return res;

View File

@ -233,25 +233,30 @@ FitResult XYStatData::fit(Minimizer &minimizer, const DVec &init,
const Index nPoint = getNFitPoint();
DVec fullInit = init;
Index is = 0, kf = 0;
fullInit.conservativeResize(chi2_.getNArg());
for (Index i = 0; i < getXDim(); ++i)
if (!isXExact(i))
{
for (Index k = 0; k < getNData(); ++k)
if (isFitPoint(k))
if (!isXExact(i))
{
fullInit(chi2_.getNPar() + nPoint*is + kf) = x(i, k);
kf++;
kf = 0;
for (Index k = 0; k < getNData(); ++k)
{
if (isFitPoint(k))
{
fullInit(chi2_.getNPar() + nPoint*is + kf) = x(i, k);
kf++;
}
}
is++;
}
is++;
}
minimizer.setInit(fullInit);
// fit
DoubleFunction chi2 = chi2_.makeFunction(false);
FitResult result;
result = minimizer(chi2);
result.chi2_ = chi2(result);
result.nDof_ = chi2_.getNDof();