loon: NC::Fitter::MarginalizeHybrid Class Reference

00342 {}//{pminu = 0;}

01532     :fFunc(f), fParams(p), fTestMinimum(testmin), fXpos(xcoord), fYpos(ycoord)
01533   {
01534     assert(fXpos < int(fParams.size()) && fYpos < int(fParams.size()));
01535     assert(fYpos > fXpos || fYpos == -1);
01536     // Recycling the old Minuit instance each time is a possible cause of
01537     // getting stuck in the wrong minimum.
01538     /*
01539     pminu = new MinuitForCallable(fFunc, fParams.size());
01540 
01541     pminu->SetPrintLevel(-1);
01542     //  minu.Command("set strategy 2"); // Be careful, instead of fast.
01543 
01544     for(int n = 0; n < int(fParams.size()); ++n){
01545       pminu->DefineParameter(n, fParams[n].ShortName().c_str(),
01546                              (fParams[n].min+fParams[n].Max())/2,
01547                              (fParams[n].Max()-fParams[n].Min())/4, 0, 0);
01548     }
01549     */
01550     fLastMinPos.resize(fParams.size());
01551   }

00349     {
00350       //      delete pminu;
00351     }

01630   {
01631     CoordNDim p = fLastMinPos;   // Last place we were
01632     // Move to where we are now on one axis of the contour
01633     p[fXpos] = r.x;
01634     if(fYpos > 0) p[fYpos] = r.y; // If we have a second axis move there too.
01635 
01636     CoordNDim junk;
01637     valAtMin = fFunc->EvalAtEx(p, &junk);
01638 
01639     // What direction to take a step in.
01640     for(int n = 0; n < int(fParams.size()); ++n){
01641       // Not in the contour axis directions.
01642       if(n == fXpos || n == fYpos) continue;
01643       CoordNDim q = p;
01644       q[n] += fParams[n].Precision();   // Move one step positive in this axis
01645       if(fFunc->EvalAtEx(q, &junk) < valAtMin) return false;
01646       q[n] -= 2*fParams[n].Precision(); // Also try one step negative.
01647       if(fFunc->EvalAtEx(q, &junk) < valAtMin) return false;
01648     }
01649     MSG("NCFitter", Msg::kVerbose) << "MarginalizeHybrid: "
01650                                    << "Already at minimum" << endl;
01651     return true;
01652   }

01556   {
01557     //assert(pminu);
01558 
01559     double valAtMin;
01560     //    if(fTestMinimum && AlreadyAtMinimum(r, valAtMin)) return valAtMin;
01561 
01562     MinuitForCallable* pminu = new MinuitForCallable(fFunc, fParams.size());
01563 
01564     pminu->SetPrintLevel(-1);
01565     //    pminu->Command("set strategy 2"); // Be careful, instead of fast.
01566 
01567     for(int n = 0; n < int(fParams.size()); ++n){
01568       pminu->DefineParameter(n, fParams[n].ShortName().c_str(),
01569                              (fParams[n].Min()+fParams[n].Max())/2,
01570                              (fParams[n].Max()-fParams[n].Min())/4, 0, 0);
01571     }
01572 
01573     // Minuit tries to return parameters we don't want to keep.
01574     double junkd; int junki;
01575 
01576     // Add one to the index because of ridiculous fortran numbering scheme.
01577     double args1[] = {fXpos+1, r.x};
01578     pminu->mnexcm("set param", args1, 2, junki);
01579     if(fYpos > 0){
01580       double args2[] = {fYpos+1, r.y};
01581       pminu->mnexcm("set param", args2, 2, junki);
01582     }
01583 
01584     pminu->mnfixp(fXpos, junki);
01585     if(fYpos > 0){
01586       // The parameters all got renumbered when we fixed x.
01587       pminu->mnfixp(fYpos-1, junki);
01588     }
01589 
01590     // Simplex seems to be much faster than Migrad.
01591     // I think Migrad can fail when put in a bad place, and spends
01592     // a long time doing it.
01593     // However - I have had bad results using Simplex, so stick with Migrad.
01594     //  minu.mnsimp();
01595     if(pminu->Migrad()){
01596       //return FLT_MAX; // this messes up the contours...
01597       MSG("NCFitter", Msg::kVerbose) << "MinWRT: Migrad didn't converge"
01598                                      << endl;
01599     }
01600 
01601     pminu->mnstat(valAtMin, junkd, junkd, junki, junki, junki);
01602     for(unsigned int n = 0; n < fParams.size(); ++n){
01603       pminu->GetParameter(n, fLastMinPos[n], junkd);
01604       if(int(n) != fXpos && int(n) != fYpos)
01605             ret.push_back(fLastMinPos[n]);
01606     }
01607 
01608     TString txt = "MarginalizeHybrid: min at ";
01609     for(unsigned int n = 0; n < fParams.size(); ++n){
01610       txt += fParams[n].ShortName() + "=";
01611       txt += fLastMinPos[n];
01612       txt += +"\t";
01613     }
01614     MSG("NCFitter", Msg::kVerbose) << txt << " chisq = "
01615                                    << valAtMin << endl;
01616 
01617     // Free all the parameters, we will refix them next time round.
01618     pminu->mnfree(0);
01619 
01620     delete pminu;
01621 
01622     return valAtMin;
01623   }


Public Member Functions
	MarginalizeHybrid ()
	MarginalizeHybrid (const ICallableND *f, const std::vector< NCParameter > params, const bool testmin, const int xparam, const int yparam=-1)
	Minimize callable f with respect to some of its arguments.
virtual	~MarginalizeHybrid ()
virtual double	EvalAtEx (const Coord &r, CoordNDim &ret) const
Private Member Functions
MarginalizeHybrid &	operator= (const MarginalizeHybrid &)
bool	AlreadyAtMinimum (const Coord &r, double &valAtMin) const
Private Attributes
const ICallableND *	fFunc
std::vector< NCParameter >	fParams
bool	fTestMinimum
int	fXpos
int	fYpos
CoordNDim	fLastMinPos

NC::Fitter::MarginalizeHybrid Class Reference

Public Member Functions

Private Member Functions

Private Attributes

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation