#include "sopnamsp.h" #include "machdefs.h" #include #include #include #include #include #include "pexceptions.h" #include "minuitadapt.h" #include "localmnpout.h" #define AVEC_CFORTRAN #ifdef AVEC_CFORTRAN #include "Cern/cfortran.h" #include "Cern/minuit.h" #else #include "minuitprot.h" #endif /*! \class SOPHYA::Minuit Cern Fitter adapter \ingroup MinuitAdapt */ MinuitAdapt::MinuitAdapt( void (*myfcn)(int_4 *,double *,double *,double *,int_4 *,double (*f)(double *)), double (*myfutils)(double *x) ) : mNPar(0), iErFlg(0), fcn(myfcn), futils(myfutils) { int_4 i1=5,i2=6,i3=0; #ifdef AVEC_CFORTRAN MNINIT(i1,i2,i3); #else mninit_(&i1,&i2,&i3); #endif } MinuitAdapt::~MinuitAdapt(void) { } void MinuitAdapt::SetTitle(char* title) { int_4 nc=strlen(title); if(nc>0) { #ifdef AVEC_CFORTRAN MNSETI(title); #else mnseti_(title,nc); #endif } } void MinuitAdapt::Clear(void) { char* str="CLEAR"; double arglis[1]={0.}; int_4 narg=0; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif mNPar = 0; } void MinuitAdapt::Return(void) { char* str="RETURN"; double arglis[1]={0.}; int_4 narg=0; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::PrintLevel(int_4 n) { char* str="SET PRINTOUT"; double arglis[1]={(double)n}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetWidthPage(int_4 n) { char* str="SET WIDTHPAGE"; double arglis[1]={(double)n}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetWarnings(bool w) { char* str="SET WARNINGS"; if(!w) str="SET NOWARNINGS"; double arglis[1]={0.}; int_4 narg=0; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetErrorDef(double err) // 1=xi2, 0.5=likelyhood { char* str="SET ERRORDEF"; double arglis[1]={err}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetEpsMachine(double eps) { char* str="SET EPSMACHINE"; double arglis[1]={eps}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetStrategy(int_4 n) { char* str="SET STRATEGY"; double arglis[1]={(double)n}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetGradient(int_4 n) // n<0 : gradient from finite differences computed by Minuit // n=0 : gradient from user AND finite diff.. then compare // n>0 : force gradient from user { char* str="SET NOGRADIENT"; double arglis[1]={1.}; int_4 narg=0; if(n==0) str="SET GRADIENT"; else if(n>0) {str="SET GRADIENT"; narg=1;} #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetRandom(int_4 seed) // Seed must be integer between 10000 and 900000000 { char* str="SET RANDOMGENERATOR"; double arglis[1]={(double)seed}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::DefineParameter(int_4 n,char* name,double init,double step ,double vmin,double vmax) // nom parametre, numero parametre, initial value, initial step size, parameter bounds // - ATTENTION numero de parametre [1,...] { char str1[16]; sprintf(str1,"P%d",n); char* str=str1; if(name!=NULL) str=name; #ifdef AVEC_CFORTRAN MNPARM(n,str,init,step,vmin,vmax,iErFlg); #else int_4 nc=strlen(str); mnparm_(&n,str,&init,&step,&vmin,&vmax,&iErFlg,nc); #endif mNPar++; } void MinuitAdapt::SetParameter(int_4 n,double val) { char* str="SET PARAMETER"; double arglis[2]={(double)n,val}; int_4 narg=2; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetLimits(int_4 n,double val1,double val2) { char* str="SET LIMITS"; double arglis[3]={(double)n,val1,val2}; int_4 narg=3; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::SetFix(int_4 n) { char* str="FIX"; double arglis[1]={(double)n}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::Release(int_4 n) { char* str="RELEASE"; double arglis[1]={(double)n}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::Restore(bool only_last_fixed) { char* str="RESTORE"; double arglis[1]={0.}; int_4 narg=0; if(only_last_fixed) {arglis[0]=1.; narg=1;} #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::DoFit(char *method,int_4 maxcall,double tol) { char* str=method; double arglis[2]={(double)maxcall,tol}; int_4 narg=2; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::Improve(int_4 maxcall) { char* str="IMPROVE"; double arglis[1]={(double)maxcall}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::Hesse(int_4 maxcall) { char* str="HESSE"; double arglis[1]={(double)maxcall}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::Minos(int_4 maxcall) { char* str="MINOS"; double arglis[1]={(double)maxcall}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::Seek(int_4 maxcall,double stdev) { char* str="SEEK"; double arglis[2]={(double)maxcall,stdev}; int_4 narg=2; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::Call(int_4 iflag) // See meaning in fcn exemple { char* str="CALL"; double arglis[1]={(double)iflag}; int_4 narg=1; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } void MinuitAdapt::DrawContour(int_4 n1,int_4 n2,int_4 npts) { char* str="MNCONTOUR"; double arglis[3]={(double)n1,(double)n2,(double)npts}; int_4 narg=3; #ifdef AVEC_CFORTRAN MNEXCM(fcn,str,arglis,narg,iErFlg,futils); #else int_4 nc=strlen(str); mnexcm_(fcn,str,arglis,&narg,&iErFlg,futils,nc); #endif } int_4 MinuitAdapt::GetContour(int_4 n1,int_4 n2,int_4 npts ,TVector& xpt,TVector& ypt) { if(npts<5) npts=5; int_4 ncontok = 0; double* xcont = new double[npts]; double* ycont = new double[npts]; #ifdef AVEC_CFORTRAN MNCONT(fcn,n1,n2,npts,xcont[0],ycont[0],ncontok,futils); #else mncont_(fcn,&n1,&n2,&npts,xcont,ycont,&ncontok,futils); #endif if(ncontok<1) return 0; xpt.ReSize(ncontok); ypt.ReSize(ncontok); for(int_4 i=0;i0) for(unsigned int i=0;i MinuitAdapt::GetErrorsMatrix(void) { TMatrix mat; int_4 npar=NPar(); if(npar<=0) return mat; double* mcov = new double[npar*npar]; #ifdef AVEC_CFORTRAN MNEMAT(mcov[0],npar); #else mnemat_(mcov,&npar); #endif mat.ReSize(npar,npar); for(int i=0;i