Hi Topher,
I assume that Wim Lavrijsen will answer your points concerning directly pyroot.
Concerning access to the covariance matrix following a fit,
I have added a new method in TVirtualFitter, TFitter, TFumili, etc
TVirtualFitter::GetCovarianceMatrixElement(int i, int j);
you can access it as indicated in the comments of TH1::Fit
// Access to the fit covariance matrix
// ===================================
// Example1:
// TH1F h("h","test",100,-2,2);
// h.FillRandom("gaus",1000);
// h.Fit("gaus");
// Double_t matrix[3][3];
// gMinuit->mnemat(&matrix[0][0],3);
// Example2:
// TH1F h("h","test",100,-2,2);
// h.FillRandom("gaus",1000);
// h.Fit("gaus");
// TMatrixD matrix(npar,npar,gMinuit->GetCovarianceMatrix());
// matrix.Print();
// matrix.Draw("text");
// Example3: via the TVirtualFitter interface (works with TMinuit, TFumili,
etc)
// TH1F h("h","test",100,-2,2);
// h.FillRandom("gaus",1000);
// h.Fit("gaus");
// TVirtualfitter *fitter = TVirtualFitter::GetFitter();
// TMatrixD matrix(npar,npar,fitter->GetCovarianceMatrix());
// Double_t errorFirstPar = fitter->GetCovarianceMatrixElement(0,0);
Also note that in the class TmatrixD, you can use
TMatrixD::GetMatrixArray
to get a pointer to the matrix elements directly.
There is no point in adding a function to TF1 to return the covariance matrix
or one of its elements. TF1 does not store a copy of the covariance matrix.
Accesd to the matrix can only be done (via TMinuit or TVirtualFitter)
immediatly after the fit as shown in the comments above.
Rene Brun
Topher Cawlfield wrote:
>
> Hi,
>
> I haven't seen much discussion on roottalk about pyroot, so I'm not sure
> if this is the best place to ask my questions. Please let me know if I
> should direct my queries elsewhere.
>
> I'm having a hard time getting at Minuit's error matrix from a python
> script. Two methods to access the error matrix (from C++) are given in
> the TH1 documentation:
>
> > Access to the fit covariance matrix
> > ===================================
> > Example1:
> > TH1F <http://root.cern.ch/root/html/TH1F.html> h("h","test",100,-2,2);
> > h.FillRandom <http://root.cern.ch/root/html/TH1.html#TH1:FillRandom>("gaus",1000);
> > h.Fit <http://root.cern.ch/root/html/TH1.html#TH1:Fit>("gaus");
> > Double_t <http://root.cern.ch/root/html/ListOfTypes.html#Double_t> matrix[3][3];
> > gMinuit <http://root.cern.ch/root/html/TMinuit.html>->mnemat <http://root.cern.ch/root/html/TMinuit.html#TMinuit:mnemat>(&matrix[0][0],3);
> > Example2:
> > TH1F <http://root.cern.ch/root/html/TH1F.html> h("h","test",100,-2,2);
> > h.FillRandom <http://root.cern.ch/root/html/TH1.html#TH1:FillRandom>("gaus",1000);
> > h.Fit <http://root.cern.ch/root/html/TH1.html#TH1:Fit>("gaus");
> > TMatrixD <http://root.cern.ch/root/html/TMatrixD.html> matrix(npar,npar);
> > gMinuit <http://root.cern.ch/root/html/TMinuit.html>->mnemat <http://root.cern.ch/root/html/TMinuit.html#TMinuit:mnemat>(matrix.GetElements(),npar);
> > matrix.Print <http://root.cern.ch/root/html/TH1.html#TH1:Print>();
> > matrix.Draw <http://root.cern.ch/root/html/TH1.html#TH1:Draw>("text");
>
> I'm having two problems with this approach. First, I can't find gMinuit
> at all. I can get gROOT, gGeometry, and gDirectory like this:
> >>> not not ROOT.gROOT.GetGlobal( 'gDirectory' )
> True
>
> <>But not gMinuit:
> >>> not not ROOT.gROOT.GetGlobal( 'gMinuit' )
> False
>
> I can get fits to work, so the object does exist somewhere.
>
> My second problem is that I don't see how to implement either of these
> methods in python. I can't make a 2-D matrix with "Double_t
> matrix[3][3]" and pass it to mnemat() through pyroot. I CAN create a
> TMatrixD, but the only methods to access the elements are overloaded
> array-lookup operators (irritating!), and I can't find a way to use
> these through pyroot.
>
> It would be very nice (feature request) if a GetErrorMatrixElem(ipar,
> jpar) method would be added to TF1 which does all the work of mapping
> parameters to Minuit's reduced free-parameter list, and returns the
> corresponding element in the error matrix. It doesn't sound like that
> should be too hard -- should I attempt it myself?
>
> Thanks in advance for any help on these problems!
> Topher Cawlfield
>
> p.s. I really appreciate pyROOT! Thanks! This was the final
> ingredient that I was hoping for to create a "timeline" system for
> CLEO. What I needed was a way to access histogram files, initially in
> hbook (zebra) format, through a high-level language like Perl, Python,
> Ruby, etc. The combination of h2root and pyroot was just perfect, and
> Python is proving itself to be the ideal language for the project in
> other ways as well. I'm implementing a general-purpose fitting
> capability, and ran into trouble computing the error in the yield of a
> fit of gaus + background. Root (inheriting bad genes from PAW) uses a
> most unfortunate parameterization for Gaussians, picking amplitude
> instead of area for the first parameter. So to get the error in the
> yield one needs the error in the amplitude, error in the width, and the
> correlation coefficient between the two. Really tragic since we all
> like to fit Gaussians for yields twice each morning before breakfast.
This archive was generated by hypermail 2b29 : Sun Jan 02 2005 - 05:50:09 MET