mirror of
https://github.com/aportelli/LatAnalyze.git
synced 2025-04-11 03:20:46 +01:00
first commit for the new fit interface
This commit is contained in:
parent
352caf30b6
commit
917031747f
@ -26,7 +26,7 @@ noinst_PROGRAMS = \
|
|||||||
exRootFinder
|
exRootFinder
|
||||||
|
|
||||||
if HAVE_MINUIT
|
if HAVE_MINUIT
|
||||||
noinst_PROGRAMS += exMin
|
noinst_PROGRAMS += exFit exMin
|
||||||
endif
|
endif
|
||||||
|
|
||||||
exCompiledDoubleFunction_SOURCES = exCompiledDoubleFunction.cpp
|
exCompiledDoubleFunction_SOURCES = exCompiledDoubleFunction.cpp
|
||||||
@ -37,6 +37,12 @@ exDerivative_SOURCES = exDerivative.cpp
|
|||||||
exDerivative_CFLAGS = -g -O2
|
exDerivative_CFLAGS = -g -O2
|
||||||
exDerivative_LDFLAGS = -L../lib/.libs -lLatAnalyze
|
exDerivative_LDFLAGS = -L../lib/.libs -lLatAnalyze
|
||||||
|
|
||||||
|
if HAVE_MINUIT
|
||||||
|
exFit_SOURCES = exFit.cpp
|
||||||
|
exFit_CFLAGS = -g -O2
|
||||||
|
exFit_LDFLAGS = -L../lib/.libs -lLatAnalyze
|
||||||
|
endif
|
||||||
|
|
||||||
exInterp_SOURCES = exInterp.cpp
|
exInterp_SOURCES = exInterp.cpp
|
||||||
exInterp_CFLAGS = -g -O2
|
exInterp_CFLAGS = -g -O2
|
||||||
exInterp_LDFLAGS = -L../lib/.libs -lLatAnalyze
|
exInterp_LDFLAGS = -L../lib/.libs -lLatAnalyze
|
||||||
|
@ -1,54 +1,17 @@
|
|||||||
#include <iostream>
|
#include <iostream>
|
||||||
#include <cmath>
|
#include <LatAnalyze/FitInterface.hpp>
|
||||||
#include <LatAnalyze/CompiledModel.hpp>
|
|
||||||
#include <LatAnalyze/MinuitMinimizer.hpp>
|
|
||||||
#include <LatAnalyze/Plot.hpp>
|
|
||||||
#include <LatAnalyze/RandGen.hpp>
|
|
||||||
#include <LatAnalyze/XYStatData.hpp>
|
|
||||||
|
|
||||||
using namespace std;
|
using namespace std;
|
||||||
using namespace Latan;
|
using namespace Latan;
|
||||||
|
|
||||||
const Index nPoint = 20;
|
|
||||||
const double exactPar[2] = {0.5,5.0}, dx = 10.0/static_cast<double>(nPoint);
|
|
||||||
|
|
||||||
int main(void)
|
int main(void)
|
||||||
{
|
{
|
||||||
// generate fake data
|
FitInterface f;
|
||||||
XYStatData data(nPoint, 1, 1);
|
|
||||||
RandGen rg;
|
|
||||||
double x_k, y_k;
|
|
||||||
DoubleModel f = compile("return p_1*exp(-x_0*p_0);", 1, 2);
|
|
||||||
|
|
||||||
for (Index k = 0; k < nPoint; ++k)
|
|
||||||
{
|
|
||||||
x_k = k*dx;
|
|
||||||
y_k = f(&x_k, exactPar) + rg.gaussian(0.0, 0.1);
|
|
||||||
cout << x_k << " " << y_k << " " << 0.1 << endl;
|
|
||||||
data.x(0, k) = x_k;
|
|
||||||
data.y(0, k) = y_k;
|
|
||||||
}
|
|
||||||
data.yyVar(0, 0).diagonal() = DMat::Constant(nPoint, 1, 0.1*0.1);
|
|
||||||
data.assumeXExact(0);
|
|
||||||
|
|
||||||
// fit
|
|
||||||
DVec init = DVec::Constant(2, 0.5);
|
|
||||||
FitResult p;
|
|
||||||
MinuitMinimizer minimizer;
|
|
||||||
|
|
||||||
data.fitAllPoints();
|
f.addXDim("a", 3);
|
||||||
p = data.fit(minimizer, init, f);
|
f.addXDim("b", 4);
|
||||||
|
f.addXDim("c", 3);
|
||||||
cout << "a= " << p(0) << " b= " << p(1);
|
f.addYDim("y");
|
||||||
cout << " chi^2/ndof= " << p.getChi2PerDof();
|
|
||||||
cout << " p-value= " << p.getPValue() <<endl;
|
|
||||||
|
|
||||||
// plot result
|
|
||||||
Plot plot;
|
|
||||||
|
|
||||||
plot << LogScale(Axis::y) << PlotData(data);
|
|
||||||
plot << Color("rgb 'blue'") << PlotFunction(p.getModel(), 0.0, 10.0);
|
|
||||||
plot.display();
|
|
||||||
|
|
||||||
return EXIT_SUCCESS;
|
return EXIT_SUCCESS;
|
||||||
}
|
}
|
||||||
|
54
examples/exFit.old.cpp
Normal file
54
examples/exFit.old.cpp
Normal file
@ -0,0 +1,54 @@
|
|||||||
|
#include <iostream>
|
||||||
|
#include <cmath>
|
||||||
|
#include <LatAnalyze/CompiledModel.hpp>
|
||||||
|
#include <LatAnalyze/MinuitMinimizer.hpp>
|
||||||
|
#include <LatAnalyze/Plot.hpp>
|
||||||
|
#include <LatAnalyze/RandGen.hpp>
|
||||||
|
#include <LatAnalyze/XYStatData.hpp>
|
||||||
|
|
||||||
|
using namespace std;
|
||||||
|
using namespace Latan;
|
||||||
|
|
||||||
|
const Index nPoint = 20;
|
||||||
|
const double exactPar[2] = {0.5,5.0}, dx = 10.0/static_cast<double>(nPoint);
|
||||||
|
|
||||||
|
int main(void)
|
||||||
|
{
|
||||||
|
// generate fake data
|
||||||
|
XYStatData data(nPoint, 1, 1);
|
||||||
|
RandGen rg;
|
||||||
|
double x_k, y_k;
|
||||||
|
DoubleModel f = compile("return p_1*exp(-x_0*p_0);", 1, 2);
|
||||||
|
|
||||||
|
for (Index k = 0; k < nPoint; ++k)
|
||||||
|
{
|
||||||
|
x_k = k*dx;
|
||||||
|
y_k = f(&x_k, exactPar) + rg.gaussian(0.0, 0.1);
|
||||||
|
cout << x_k << " " << y_k << " " << 0.1 << endl;
|
||||||
|
data.x(0, k) = x_k;
|
||||||
|
data.y(0, k) = y_k;
|
||||||
|
}
|
||||||
|
data.yyVar(0, 0).diagonal() = DMat::Constant(nPoint, 1, 0.1*0.1);
|
||||||
|
data.assumeXExact(0);
|
||||||
|
|
||||||
|
// fit
|
||||||
|
DVec init = DVec::Constant(2, 0.5);
|
||||||
|
FitResult p;
|
||||||
|
MinuitMinimizer minimizer;
|
||||||
|
|
||||||
|
data.fitAllPoints();
|
||||||
|
p = data.fit(minimizer, init, f);
|
||||||
|
|
||||||
|
cout << "a= " << p(0) << " b= " << p(1);
|
||||||
|
cout << " chi^2/ndof= " << p.getChi2PerDof();
|
||||||
|
cout << " p-value= " << p.getPValue() <<endl;
|
||||||
|
|
||||||
|
// plot result
|
||||||
|
Plot plot;
|
||||||
|
|
||||||
|
plot << LogScale(Axis::y) << PlotData(data);
|
||||||
|
plot << Color("rgb 'blue'") << PlotFunction(p.getModel(), 0.0, 10.0);
|
||||||
|
plot.display();
|
||||||
|
|
||||||
|
return EXIT_SUCCESS;
|
||||||
|
}
|
74
lib/FitInterface.cpp
Normal file
74
lib/FitInterface.cpp
Normal file
@ -0,0 +1,74 @@
|
|||||||
|
/*
|
||||||
|
* FitInterface.cpp, part of LatAnalyze 3
|
||||||
|
*
|
||||||
|
* Copyright (C) 2013 - 2016 Antonin Portelli
|
||||||
|
*
|
||||||
|
* LatAnalyze 3 is free software: you can redistribute it and/or modify
|
||||||
|
* it under the terms of the GNU General Public License as published by
|
||||||
|
* the Free Software Foundation, either version 3 of the License, or
|
||||||
|
* (at your option) any later version.
|
||||||
|
*
|
||||||
|
* LatAnalyze 3 is distributed in the hope that it will be useful,
|
||||||
|
* but WITHOUT ANY WARRANTY; without even the implied warranty of
|
||||||
|
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
||||||
|
* GNU General Public License for more details.
|
||||||
|
*
|
||||||
|
* You should have received a copy of the GNU General Public License
|
||||||
|
* along with LatAnalyze 3. If not, see <http://www.gnu.org/licenses/>.
|
||||||
|
*/
|
||||||
|
|
||||||
|
#include <LatAnalyze/FitInterface.hpp>
|
||||||
|
#include <LatAnalyze/includes.hpp>
|
||||||
|
|
||||||
|
using namespace std;
|
||||||
|
using namespace Latan;
|
||||||
|
|
||||||
|
/******************************************************************************
|
||||||
|
* FitInterface implementation *
|
||||||
|
******************************************************************************/
|
||||||
|
// constructor /////////////////////////////////////////////////////////////////
|
||||||
|
FitInterface::FitInterface(void)
|
||||||
|
{}
|
||||||
|
|
||||||
|
// add dimensions //////////////////////////////////////////////////////////////
|
||||||
|
void FitInterface::addXDim(const string name, const Index nData)
|
||||||
|
{
|
||||||
|
xDimName_.push_back(name);
|
||||||
|
xSize_.push_back(nData);
|
||||||
|
xIndex_[name] = xDimName_.size();
|
||||||
|
updateYSize();
|
||||||
|
}
|
||||||
|
|
||||||
|
void FitInterface::addYDim(const string name)
|
||||||
|
{
|
||||||
|
yDimName_.push_back(name);
|
||||||
|
yIndex_[name] = yDimName_.size();
|
||||||
|
}
|
||||||
|
|
||||||
|
unsigned int FitInterface::getNXDim(void)
|
||||||
|
{
|
||||||
|
return xDimName_.size();
|
||||||
|
}
|
||||||
|
|
||||||
|
unsigned int FitInterface::getNYDim(void)
|
||||||
|
{
|
||||||
|
return yDimName_.size();
|
||||||
|
}
|
||||||
|
|
||||||
|
void FitInterface::updateYSize(void)
|
||||||
|
{
|
||||||
|
ySize_ = 1;
|
||||||
|
for (Index size: xSize_)
|
||||||
|
{
|
||||||
|
ySize_ *= size;
|
||||||
|
}
|
||||||
|
isFitPoint_.resize(ySize_);
|
||||||
|
isFitPoint_.fill(1);
|
||||||
|
}
|
||||||
|
|
||||||
|
// enable fit points ///////////////////////////////////////////////////////////
|
||||||
|
void FitInterface::fitPoint(const bool isFitPoint, const Index j)
|
||||||
|
{
|
||||||
|
DEBUG_VAR(j);
|
||||||
|
isFitPoint_(j) = isFitPoint ? 1 : 0;
|
||||||
|
}
|
100
lib/FitInterface.hpp
Normal file
100
lib/FitInterface.hpp
Normal file
@ -0,0 +1,100 @@
|
|||||||
|
/*
|
||||||
|
* FitInterface.hpp, part of LatAnalyze 3
|
||||||
|
*
|
||||||
|
* Copyright (C) 2013 - 2016 Antonin Portelli
|
||||||
|
*
|
||||||
|
* LatAnalyze 3 is free software: you can redistribute it and/or modify
|
||||||
|
* it under the terms of the GNU General Public License as published by
|
||||||
|
* the Free Software Foundation, either version 3 of the License, or
|
||||||
|
* (at your option) any later version.
|
||||||
|
*
|
||||||
|
* LatAnalyze 3 is distributed in the hope that it will be useful,
|
||||||
|
* but WITHOUT ANY WARRANTY; without even the implied warranty of
|
||||||
|
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
||||||
|
* GNU General Public License for more details.
|
||||||
|
*
|
||||||
|
* You should have received a copy of the GNU General Public License
|
||||||
|
* along with LatAnalyze 3. If not, see <http://www.gnu.org/licenses/>.
|
||||||
|
*/
|
||||||
|
|
||||||
|
#ifndef Latan_FitInterface_hpp_
|
||||||
|
#define Latan_FitInterface_hpp_
|
||||||
|
|
||||||
|
#include <LatAnalyze/Global.hpp>
|
||||||
|
|
||||||
|
BEGIN_LATAN_NAMESPACE
|
||||||
|
|
||||||
|
/******************************************************************************
|
||||||
|
* FitInterface *
|
||||||
|
******************************************************************************/
|
||||||
|
|
||||||
|
class FitInterface
|
||||||
|
{
|
||||||
|
public:
|
||||||
|
// constructor
|
||||||
|
FitInterface(void);
|
||||||
|
// destructor
|
||||||
|
virtual ~FitInterface(void) = default;
|
||||||
|
// add dimensions
|
||||||
|
void addXDim(const std::string name, const Index nData);
|
||||||
|
void addYDim(const std::string name);
|
||||||
|
unsigned int getNXDim(void);
|
||||||
|
unsigned int getNYDim(void);
|
||||||
|
// Y dimension index helper
|
||||||
|
template <typename... Ts>
|
||||||
|
Index yIndex(const Ts... is);
|
||||||
|
// enable fit points
|
||||||
|
void fitPoint(const bool isFitPoint, const Index j);
|
||||||
|
private:
|
||||||
|
void updateYSize(void);
|
||||||
|
private:
|
||||||
|
std::vector<std::string> xDimName_, yDimName_;
|
||||||
|
std::map<std::string, Index> xIndex_, yIndex_;
|
||||||
|
std::vector<Index> xSize_;
|
||||||
|
Index ySize_;
|
||||||
|
public:
|
||||||
|
IVec isFitPoint_;
|
||||||
|
};
|
||||||
|
|
||||||
|
/******************************************************************************
|
||||||
|
* FitInterface template implementation *
|
||||||
|
******************************************************************************/
|
||||||
|
// Y dimension index helper ////////////////////////////////////////////////////
|
||||||
|
template <typename... Ts>
|
||||||
|
Index FitInterface::yIndex(const Ts... is)
|
||||||
|
{
|
||||||
|
static_assert(static_or<std::is_convertible<Index, Ts>::value...>::value,
|
||||||
|
"fitPoint arguments are not compatible with Index");
|
||||||
|
|
||||||
|
constexpr size_t n = sizeof...(is);
|
||||||
|
|
||||||
|
if (n != getNXDim())
|
||||||
|
{
|
||||||
|
LATAN_ERROR(Size, "number of arguments and number of X dimensions "
|
||||||
|
"mismatch");
|
||||||
|
}
|
||||||
|
|
||||||
|
constexpr std::array<Index, sizeof...(is)> i = {is...};
|
||||||
|
Index j;
|
||||||
|
|
||||||
|
for (unsigned int d = 0; d < n; ++d)
|
||||||
|
{
|
||||||
|
if (i[d] >= xSize_[d])
|
||||||
|
{
|
||||||
|
LATAN_ERROR(Range, "index in X dimension " + strFrom(d)
|
||||||
|
+ " out of range");
|
||||||
|
}
|
||||||
|
}
|
||||||
|
j = xSize_[1]*i[0];
|
||||||
|
for (unsigned int d = 1; d < n-1; ++d)
|
||||||
|
{
|
||||||
|
j = xSize_[d+1]*(i[d] + j);
|
||||||
|
}
|
||||||
|
j += i[n-1];
|
||||||
|
|
||||||
|
return j;
|
||||||
|
}
|
||||||
|
|
||||||
|
END_LATAN_NAMESPACE
|
||||||
|
|
||||||
|
#endif // Latan_FitInterface_hpp_
|
@ -23,6 +23,7 @@ libLatAnalyze_la_SOURCES = \
|
|||||||
Derivative.cpp \
|
Derivative.cpp \
|
||||||
Exceptions.cpp \
|
Exceptions.cpp \
|
||||||
File.cpp \
|
File.cpp \
|
||||||
|
FitInterface.cpp \
|
||||||
Function.cpp \
|
Function.cpp \
|
||||||
Global.cpp \
|
Global.cpp \
|
||||||
GslHybridRootFinder.cpp\
|
GslHybridRootFinder.cpp\
|
||||||
@ -55,6 +56,7 @@ libLatAnalyze_la_HEADERS = \
|
|||||||
Exceptions.hpp \
|
Exceptions.hpp \
|
||||||
Function.hpp \
|
Function.hpp \
|
||||||
File.hpp \
|
File.hpp \
|
||||||
|
FitInterface.hpp \
|
||||||
Global.hpp \
|
Global.hpp \
|
||||||
GslHybridRootFinder.hpp\
|
GslHybridRootFinder.hpp\
|
||||||
GslQagsIntegrator.hpp \
|
GslQagsIntegrator.hpp \
|
||||||
|
Loading…
x
Reference in New Issue
Block a user