Mandeep,
On Tue, 12 Feb 2002 14:34:34 -0800 (PST)
"Mandeep S. Gill" <msgill@SLAC.stanford.edu> wrote
concerning "[ROOT] Chisq for 4D histos..":
>
> Dear Rootfolk: i hope this is the right place to send this to, it's my
> first roottalk posting so: has anyone implemented anything higher than a
> 3D histo in Root?
It's hard to visualise more than 3D histos (heck even a 3D histo is
difficult), so I don't think that's something ROOT really needs.
> (The reason i ask is because i wanted to calculate a chisq/dof value
> for some 2D fits i'm doing (in RooFitCore, based on Root), and since
> chisq is only pre-implemented for 1D in RFC, i wrote the code for
> doing this myself based on the GetBinContents function for a TH2F
> object, but there are only histos predefined up to a TH3F on the
> Root class doc page so i'm wondering if i'll have to either use a
> different method for getting the chisq, or extend the functionality
> myself to higher dim'l TH objects..)
I think a better way to fit a R^4->R mapping is to use the TFitter
class (which uses TMinuit) to do your fit. Suppose you have M tuples
of (v, x, y, z, d, e) where (v, x, y, z) is the independent variables,
and (d, e) is the dependent variable, and it's error, then you can fit
the function f: (v,x,y,z)->d to the M data points (which we'll assume
we can get from some object data, using the method GetEntry(i)) to the
function f like:
// The function we need to minimise, here a Chi^2 function
void chiSquare(Int_t& nPar, Double_t* in, Double_t& f,
Double_t* p, Int_t flag)
{
f = Data::Instance()->ChiSquare(p);
}
class Data {
private:
Int_t fN; // Num entries
Double_t* fX; // [4 * fN] independent points
Double_t* fD; // [fN] dependent points
Double_t* fE; // [fN] error dependent points
static Data* fgInstance;
public:
Data() : fN(0), fX(0), fD(0), fE(0) {};
static Data* Instance() {
if (!fgInstance) fgInstance = new Data();
return fgInstance;
}
Double_t Function(Int_t i, Double_t* p) {
// The parameters are in the array pp.
// Here we need to evaluate the function at the point (v,x,y,z). As
// an example we'll just add the numbers
return p[0]
+ p[1] * fX[4 * i]
+ p[2] * fX[4 * i + 1]
+ p[3] * fX[4 * i + 2]
+ p[4] * fX[4 * i + 3];
}
Double_t ChiSquare(Double_t* p)
Double_t chi2 = 0;
for (Int_t i = 0; i < fN; i++) {
chi2 += TMath::Power((fD[i] - f(i, p))/ fE[i], 2);
return chi2;
}
void Fit() {
TFitter* fitter = new TFitter(5); // Arg is number of params
fitter->SetFCN(chiSquare);
Double_t arguments[] = { 500, 1 }
fitter->ExecuteCommand("MIGRAD", arguments, 2);
Double_t amin, edm, errdef,
Int_t nvpar, nparx;
fitter->GetStats(amin, edm, errdef, nvpar, nparx);
fitter->PrintResults(3, amin);
}
}
If you want to bin that data, you need to implement that.
Please also refer to the tutorials, HOWTOs, and the Users guide. A
search on roottalk will probably also give a good number of hits.
Rene and Fons, couldn't the TFitter have a more intuitive interface,
like for example
TFitter::Fit(Option_t* option);
TFitter::SetFunction((*)(Int_t npar, Double_t* par, Double_t& f));
TFitter::Print(Option_t*)
And maybe a special base class, like
class TFitObject {
protected:
Int_t fNParameters;
static TFitObject* fgInstance;
public:
TFitObject(Int_t nPar) : fNParameters(nPar) { fgInstance = this; }
static TFitObject* Instance() { return fgInstance; }
virtual ~TFitObject() {}
virtual Double_t MinimisationFunction(Double_t* par) = 0;
};
void FitClassFunction(Int_t&, Double_t*, Double_t&f, Double_t*, Int_t)
{
TFitObject* fit = TFitObject::Instance();
if (!fit) {
Error("FitClassFunction", "no fit object set");
return;
}
f = fit->MinimisationFunction(p);
}
so that one could sub-class the TFitObject class to do various kinds
of fit:
class TChiSquare : public TFitObject {
protected:
Int_t fNSample;
Int_t fNVariables;
public:
virtual Double_t Function(Double_t* vars, Double_t* pars) = 0;
virtual Double_t* GetX(Int_t i) = 0;
virtual Double_t* GetD(Int_t i) = 0;
virtual Double_t* GetE(Int_t i) = 0;
Double_t MinimisationFunction(Double_t* p) {
Double_t chi2 = 0;
for (Int_t i = 0; i < GetNDataPoints(); i++)
chi2 +=
TMath::Power((GetD() - Function(GetX(i), p)) / GetE(), 2);
return chi2;
}
};
and then individual users can subclass these sub-classes, implementing
the abstract member functions. Just my two-pennies worth.
Yours,
Christian Holm Christensen -------------------------------------------
Address: Sankt Hansgade 23, 1. th. Phone: (+45) 35 35 96 91
DK-2200 Copenhagen N Cell: (+45) 28 82 16 23
Denmark Office: (+45) 353 25 305
Email: cholm@nbi.dk Web: www.nbi.dk/~cholm
This archive was generated by hypermail 2b29 : Sat Jan 04 2003 - 23:50:41 MET