ROOT logo
// @(#)root/fumili:$Id: TFumili.cxx 23472 2008-04-23 14:13:30Z brun $
// Author: Stanislav Nesterov  07/05/2003

//______________________________________________________________________________
//         FUMILI  
//  Based on ideas, proposed by I.N. Silin
//    [See NIM A440, 2000 (p431)]
// converted from FORTRAN to C  by
//     Sergey Yaschenko <s.yaschenko@fz-juelich.de>
//
//______________________________________________________________________________
//BEGIN_HTML <!--
/* -->
<H2>FUMILI minimization package</H2>
<p>FUMILI is used to minimize Chi-square function or to search maximum of
likelihood function.

<p>Experimentally measured values $F_i$ are fitted with theoretical
functions $f_i({\vec x}_i,\vec\theta\,\,)$, where ${\vec x}_i$ are
coordinates, and $\vec\theta$ -- vector of parameters.

<p>For better convergence Chi-square function has to be the following form

<p>$$
{\chi^2\over2}={1\over2}\sum^n_{i=1}\left(f_i(\vec
x_i,\vec\theta\,\,)-F_i\over\sigma_i\right)^2 \eqno(1)
$$
<p>where $\sigma_i$ are errors of measured function.

<p>The minimum condition is
<p>$$
{\partial\chi^2\over\partial\theta_i}=\sum^n_{j=1}{1\over\sigma^2_j}\cdot
{\partial f_j\over\partial\theta_i}\left[f_j(\vec
x_j,\vec\theta\,\,)-F_j\right]=0,\qquad i=1\ldots m\eqno(2)
$$
<p>where m is the quantity of parameters.

<p>Expanding left part of (2) over parameter increments and
retaining only linear terms one gets
<p>$$
\left(\partial\chi^2\over\theta_i\right)_{\vec\theta={\vec\theta}^0}
+\sum_k\left(\partial^2\chi^2\over\partial\theta_i\partial\theta_k\right)_{
\vec\theta={\vec\theta}^0}\cdot(\theta_k-\theta_k^0)
= 0\eqno(3)
$$

 <p>Here ${\vec\theta}_0$ is some initial value of parameters. In general
case:
<p>$$
{\partial^2\chi^2\over\partial\theta_i\partial\theta_k}=
\sum^n_{j=1}{1\over\sigma^2_j}{\partial f_j\over\theta_i}
{\partial f_j\over\theta_k} + 
\sum^n_{j=1}{(f_j - F_j)\over\sigma^2_j}\cdot 
{\partial^2f_j\over\partial\theta_i\partial\theta_k}\eqno(4)
$$

<p>In FUMILI algorithm for second derivatives of Chi-square approximate
expression is used when last term in (4) is discarded. It is often
done, not always wittingly, and sometimes causes troubles, for example,
if user wants to limit parameters with positive values by writing down
$\theta_i^2$ instead of $\theta_i$. FUMILI will fail if one tries
minimize $\chi^2 = g^2(\vec\theta)$ where g is arbitrary function.

<p>Approximate value is:
<p>$${\partial^2\chi^2\over\partial\theta_i\partial\theta_k}\approx
Z_{ik}=
\sum^n_{j=1}{1\over\sigma^2_j}{\partial f_j\over\theta_i}
{\partial f_j\over\theta_k}\eqno(5)
$$

<p>Then the equations for parameter increments are
<p>$$\left(\partial\chi^2\over\partial\theta_i\right)_{\vec\theta={\vec\theta}^0}
+\sum_k Z_{ik}\cdot(\theta_k-\theta^0_k) = 0, 
\qquad i=1\ldots m\eqno(6)
$$

<p>Remarkable feature of algorithm is the technique for step
restriction. For an initial value of parameter ${\vec\theta}^0$ a
parallelepiped $P_0$ is built with the center at ${\vec\theta}^0$ and
axes parallel to coordinate axes $\theta_i$. The lengths of
parallelepiped sides along i-th axis is $2b_i$, where $b_i$ is such a
value that the functions $f_j(\vec\theta)$ are quasi-linear all over
the parallelepiped. 

<p>FUMILI takes into account simple linear inequalities in the form:
$$
\theta_i^{\rm min}\le\theta_i\le\theta^{\rm max}_i\eqno(7)
$$

<p>They form parallelepiped $P$ ($P_0$ may be deformed by $P$). 
Very similar step formulae are used in FUMILI for negative logarithm
of the likelihood function with the same idea - linearization of
function argument.

<!--*/
// -->END_HTML
//______________________________________________________________________________



#include "TROOT.h"
#include "TFumili.h"
#include "TMath.h"
#include "TF1.h"
#include "TH1.h"
#include "TGraphAsymmErrors.h"

#include "Riostream.h"


extern void H1FitChisquareFumili(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
extern void H1FitLikelihoodFumili(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
extern void GraphFitChisquareFumili(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);


ClassImp(TFumili);

TFumili *gFumili=0;
// Machine dependent values  fiXME!!
// But don't set min=max=0 if param is unlimited
static const Double_t gMAXDOUBLE=1e300;
static const Double_t gMINDOUBLE=-1e300;


//______________________________________________________________________________
TFumili::TFumili(Int_t maxpar) 
{//----------- FUMILI constructor ---------
   // maxpar is the maximum number of parameters used with TFumili object
   //
   fMaxParam = TMath::Max(maxpar,25);
   if (fMaxParam>200) fMaxParam=200;
   fMaxParam2 *= fMaxParam;
   BuildArrays();
  
   fNumericDerivatives = true;
   fLogLike = false;
   fNpar    = fMaxParam;
   fGRAD    = false;
   fWARN    = true;
   fDEBUG   = false;
   fNlog    = 0;
   fSumLog  = 0;
   fNED1    = 0;
   fNED2    = 0;
   fNED12   = fNED1+fNED2;
   fEXDA    = 0;
   fFCN     = 0;
   fNfcn    = 0;
   fRP      = 1.e-15; //precision
   fS       = 1e10;
   fEPS     =0.01;
   fENDFLG  = 0;
   fNlimMul = 2;
   fNmaxIter= 150;
   fNstepDec= 3;
   fLastFixed = -1;
  
   SetName("Fumili");
   gFumili = this;
   gROOT->GetListOfSpecials()->Add(gFumili);
}

//______________________________________________________________________________
void TFumili::BuildArrays(){
   //
   //   Allocates memory for internal arrays. Called by TFumili::TFumili
   //
   fCmPar      = new Double_t[fMaxParam];
   fA          = new Double_t[fMaxParam];
   fPL0        = new Double_t[fMaxParam];
   fPL         = new Double_t[fMaxParam];
   fParamError = new Double_t[fMaxParam];
   fDA         = new Double_t[fMaxParam];
   fAMX        = new Double_t[fMaxParam];
   fAMN        = new Double_t[fMaxParam];
   fR          = new Double_t[fMaxParam];
   fDF         = new Double_t[fMaxParam];
   fGr         = new Double_t[fMaxParam]; 
   fANames     = new TString[fMaxParam];
  
   //   fX = new Double_t[10];

   Int_t zSize = fMaxParam*(fMaxParam+1)/2;
   fZ0 = new Double_t[zSize];
   fZ  = new Double_t[zSize];

   for (Int_t i=0;i<fMaxParam;i++) {
      fA[i] =0.;
      fDF[i]=0.;
      fAMN[i]=gMINDOUBLE;
      fAMX[i]=gMAXDOUBLE;
      fPL0[i]=.1;
      fPL[i] =.1;
      fParamError[i]=0.;
      fANames[i]=Form("%d",i);
   }
}


//______________________________________________________________________________
TFumili::~TFumili() {
   // 
   // TFumili destructor
   //
   DeleteArrays();
   gROOT->GetListOfSpecials()->Remove(this);
   if (gFumili == this) gFumili = 0; 
}

//______________________________________________________________________________
Double_t TFumili::Chisquare(Int_t npar, Double_t *params) const
{
   // return a chisquare equivalent
   
   Double_t amin = 0;
   H1FitChisquareFumili(npar,params,amin,params,1);
   return 2*amin;
}


//______________________________________________________________________________
void TFumili::Clear(Option_t *)
{
   //
   // Resets all parameter names, values and errors to zero
   // 
   // Argument opt is ignored
   //
   // NB: this procedure doesn't reset parameter limits 
   //
   fNpar = fMaxParam;
   fNfcn = 0;
   for (Int_t i=0;i<fNpar;i++) {
      fA[i]   =0.;
      fDF[i]  =0.;
      fPL0[i] =.1;
      fPL[i]  =.1;
      fAMN[i] = gMINDOUBLE;
      fAMX[i] = gMAXDOUBLE;
      fParamError[i]=0.;
      fANames[i]=Form("%d",i);
   }
}


//______________________________________________________________________________
void TFumili::DeleteArrays(){
   //
   // Deallocates memory. Called from destructor TFumili::~TFumili
   //
   delete[] fCmPar;
   delete[] fANames;
   delete[] fDF;
   // delete[] fX;
   delete[] fZ0;
   delete[] fZ;
   delete[] fGr;
   delete[] fA;
   delete[] fPL0;
   delete[] fPL;
   delete[] fDA;
   delete[] fAMN;
   delete[] fAMX;
   delete[] fParamError;
   delete[] fR;
}


//______________________________________________________________________________
void TFumili::Derivatives(Double_t *df,Double_t *fX){
   //  
   // Calculates partial derivatives of theoretical function
   //
   // Input:
   //    fX  - vector of data point
   // Output:
   //    DF - array of derivatives
   //
   // ARITHM.F 
   // Converted from CERNLIB
   //
   Double_t ff,ai,hi,y,pi;
   y = EvalTFN(df,fX);
   for (Int_t i=0;i<fNpar;i++) {
      df[i]=0;
      if(fPL0[i]>0.) {
         ai = fA[i]; // save current parameter value
         hi = 0.01*fPL0[i]; // diff step 
         pi = fRP*TMath::Abs(ai);
         if (hi<pi) hi = pi; // if diff step is less than precision
         fA[i] = ai+hi;
   
         if (fA[i]>fAMX[i]) { // if param is out of limits
            fA[i] = ai-hi;
            hi = -hi;
            if (fA[i]<fAMN[i]) { // again out of bounds
               fA[i] = fAMX[i];   // set param to high limit
               hi = fAMX[i]-ai;
               if (fAMN[i]-ai+hi<0) { // if hi < (ai-fAMN)
                  fA[i]=fAMN[i];
                  hi=fAMN[i]-ai;
               }
            }
         }
         ff = EvalTFN(df,fX);
         df[i] = (ff-y)/hi;
         fA[i] = ai;
      }
   }
}



//______________________________________________________________________________
Int_t TFumili::Eval(Int_t& npar, Double_t *grad, Double_t &fval, Double_t *par, Int_t flag)
{  
   // Evaluate the minimisation function
   //  Input parameters:
   //    npar:    number of currently variable parameters
   //    par:     array of (constant and variable) parameters
   //    flag:    Indicates what is to be calculated
   //    grad:    array of gradients
   //  Output parameters:
   //    fval:    The calculated function value. 
   //    grad:    The vector of first derivatives.
   // 
   // The meaning of the parameters par is of course defined by the user, 
   // who uses the values of those parameters to calculate his function value. 
   // The starting values must be specified by the user.
   // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
   // Inside FCN user has to define Z-matrix by means TFumili::GetZ 
   //  and TFumili::Derivatives,
   // set theoretical function by means of TFumili::SetUserFunc, 
   // but first - pass number of parameters by TFumili::SetParNumber
   // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
   // Later values are determined by Fumili as it searches for the minimum 
   // or performs whatever analysis is requested by the user.
   //
   // The default function calls the function specified in SetFCN
   //
   
   if (fFCN) (*fFCN)(npar,grad,fval,par,flag);
   return npar;
}


//______________________________________________________________________________
Double_t TFumili::EvalTFN(Double_t * /*df*/, Double_t *X)
{
   // Evaluate theoretical function
   // df: array of partial derivatives
   // X:  vector of theoretical function argument
   
   // for the time being disable possibility to compute derivatives
   //if(fTFN) 
   //  return (*fTFN)(df,X,fA);
   //else if(fTFNF1) {

   TF1 *f1 = (TF1*)fUserFunc;
   return f1->EvalPar(X,fA);
   //}
   //return 0.;
}

//______________________________________________________________________________
Int_t TFumili::ExecuteCommand(const char *command, Double_t *args, Int_t nargs){
   // 
   //  Execute MINUIT commands. MINImize, SIMplex, MIGrad and FUMili all
   //  will call TFumili::Minimize method.
   // 
   //  For full command list see 
   //  MINUIT. Reference Manual. CERN Program Library Long Writeup D506.
   //
   //  Improvement and errors calculation are not yet implemented as well
   //  as Monte-Carlo seeking and minimization. 
   //  Contour commands are also unsupported.
   //
   //  command   : command string
   //  args      : array of arguments
   //  nargs     : number of arguments
   //
   TString comand = command;
   static TString clower = "abcdefghijklmnopqrstuvwxyz";
   static TString cupper = "ABCDEFGHIJKLMNOPQRSTUVWXYZ";
   const Int_t nntot = 40;
   const char *cname[nntot] = {
      "MINImize  ",    //  0    checked
      "SEEk      ",    //  1    none 
      "SIMplex   ",    //  2    checked same as 0
      "MIGrad    ",    //  3    checked  same as 0
      "MINOs     ",    //  4    none 
      "SET xxx   ",    //  5 lot of stuff
      "SHOw xxx  ",    //  6 -----------
      "TOP of pag",    //  7 .
      "fiX       ",   //   8 . 
      "REStore   ",   //   9 .
      "RELease   ",   //   10 .
      "SCAn      ",   //   11  not yet implemented
      "CONtour   ",   //   12  not yet implemented
      "HESse     ",   //   13  not yet implemented
      "SAVe      ",   //   14  obsolete
      "IMProve   ",   //   15  not yet implemented
      "CALl fcn  ",   //   16 .  
      "STAndard  ",   //   17 .
      "END       ",   //   18 .
      "EXIt      ",   //   19 .
      "RETurn    ",   //   20 .
      "CLEar     ",   //   21 .
      "HELP      ",   //   22 not yet implemented
      "MNContour ",   //   23 not yet implemented
      "STOp      ",   //   24 .
      "JUMp      ",   //   25 not yet implemented
      "          ",   //   
      "          ",   // 
      "FUMili    ",   //    28 checked same as 0
      "          ",   //
      "          ",  //
      "          ",  //
      "          ",  //
      "COVARIANCE",  // 33
      "PRINTOUT  ",  // 34
      "GRADIENT  ",  // 35
      "MATOUT    ",  // 36
      "ERROR DEF ",  // 37
      "LIMITS    ",  // 38
      "PUNCH     "}; // 39

  
   fCword = comand;
   fCword.ToUpper();
   if (nargs<=0) fCmPar[0] = 0;
   Int_t i;
   for(i=0;i<fMaxParam;i++) {
      if(i<=nargs) fCmPar[i] = args[i];
   }
   /*
   fNmaxIter = int(fCmPar[0]);
   if (fNmaxIter <= 0) {
       fNmaxIter = fNpar*10 + 20 + fNpar*M*5;
   }
   fEPS = fCmPar[1];
   */
   //*-*-               look for command in list CNAME . . . . . . . . . .
   TString ctemp = fCword(0,3);
   Int_t ind;
   for (ind = 0; ind < nntot; ++ind) {
      if (strncmp(ctemp.Data(),cname[ind],3) == 0) break;
   }
   if (ind==nntot) return -3; // Unknow command - input ignored
   if (fCword(0,4) == "MINO") ind=3;
   switch (ind) {
      case 0:  case 3: case 2: case 28:
         // MINImize [maxcalls] [tolerance]
         // also SIMplex, MIGrad  and  FUMili 
         if(nargs>=1)
         fNmaxIter=TMath::Max(Int_t(fCmPar[0]),fNmaxIter); // fiXME!!
         if(nargs==2) 
            fEPS=fCmPar[1];
         return Minimize();
      case 1:
         // SEEk not implemented in this package
         return -10;

      case 4: // MINos errors analysis not implemented
         return -10;

      case 5: case 6: // SET xxx & SHOW xxx
         return ExecuteSetCommand(nargs);

      case 7: // Obsolete command
         Printf("1");
         return 0;
      case 8: // fiX <parno> ....
         if (nargs<1) return -1; // No parameters specified
         for (i=0;i<nargs;i++) {
            Int_t parnum = Int_t(fCmPar[i])-1;
            FixParameter(parnum);
         }
         return 0;
      case 9: // REStore <code>
         if (nargs<1) return 0;
         if(fCmPar[0]==0.) 
            for (i=0;i<fNpar;i++)
               ReleaseParameter(i);
         else
            if(fCmPar[0]==1.) {
               ReleaseParameter(fLastFixed);
               cout <<fLastFixed<<endl;
            }
         return 0;
      case 10: // RELease <parno> ...
         if (nargs<1) return -1; // No parameters specified
         for (i=0;i<nargs;i++) {
            Int_t parnum = Int_t(fCmPar[i])-1;
            ReleaseParameter(parnum);
         }
         return 0;
      case 11: // SCAn not implemented
         return -10;
      case 12: // CONt not implemented 
         return -10;
  
      case 13: // HESSe not implemented
         return -10;
      case 14: // SAVe
         Printf("SAVe command is obsolete");
         return -10;
      case 15: // IMProve not implemented
         return -10;
      case 16: // CALl fcn <iflag>
         {if(nargs<1) return -1;
         Int_t flag = Int_t(fCmPar[0]);
         Double_t fval;
         Eval(fNpar,fGr,fval,fA,flag);
         return 0;}
      case 17: // STAndard must call function STAND
         return 0;
      case 18:   case 19:
      case 20:  case 24: {
         Double_t fval;
         Int_t flag = 3;
         Eval(fNpar,fGr,fval,fA,flag);
         return 0;
      }
      case 21:
         Clear();
         return 0;
      case 22: //HELp not implemented
      case 23: //MNContour not implemented
      case 25: // JUMp not implemented
         return -10;
      case 26:   case 27:  case 29:  case 30:  case 31:  case 32: 
         return 0; // blank commands
      case 33:   case 34:   case 35:  case 36:   case 37:  case 38: 
      case 39:
         Printf("Obsolete command. Use corresponding SET command instead");
         return -10;
      default:
         break;
   }
   return 0;
}



//______________________________________________________________________________
Int_t TFumili::ExecuteSetCommand(Int_t nargs){
   //
   // Called from TFumili::ExecuteCommand in case 
   // of "SET xxx" and "SHOW xxx".
   //
   static Int_t nntot = 30;
   static const char *cname[30] = {
      "FCN value ", // 0 .
      "PARameters", // 1 .
      "LIMits    ", // 2 .
      "COVariance", // 3 .
      "CORrelatio", // 4 .
      "PRInt levl", // 5 not implemented yet
      "NOGradient", // 6 .
      "GRAdient  ", // 7 .
      "ERRor def ", // 8 not sure how to implement - by time being ignored
      "INPut file", // 9 not implemented 
      "WIDth page", // 10 not implemented yet
      "LINes page", // 11 not implemented yet
      "NOWarnings", // 12 .
      "WARnings  ", // 13 .
      "RANdom gen", // 14 not implemented
      "TITle     ", // 15 ignored
      "STRategy  ", // 16 ignored
      "EIGenvalue", // 17 not implemented yet 
      "PAGe throw", // 18 ignored
      "MINos errs", // 19 not implemented yet
      "EPSmachine", // 20 .
      "OUTputfile", // 21 not implemented
      "BATch     ", // 22 ignored
      "INTeractiv", // 23 ignored
      "VERsion   ", // 24 .
      "reserve   ", // 25 .
      "NODebug   ", // 26 .
      "DEBug     ", // 27 .
      "SHOw      ", // 28 err
      "SET       "};// 29 err

   TString  cfname, cmode, ckind,  cwarn, copt, ctemp, ctemp2;
   Int_t i, ind;
   Bool_t setCommand=kFALSE;
   for (ind = 0; ind < nntot; ++ind) {
      ctemp  = cname[ind];
      ckind  = ctemp(0,3);
      ctemp2 = fCword(4,6);
      if (strstr(ctemp2.Data(),ckind.Data())) break;
   }
   ctemp2 = fCword(0,3);
   if(ctemp2.Contains("SET")) setCommand=true;
   if(ctemp2.Contains("HEL") || ctemp2.Contains("SHO")) setCommand=false;
  
   if (ind>=nntot) return -3;

   switch(ind) {
      case 0: // SET FCN value illegial // SHOw only
         if(!setCommand) Printf("FCN=%f",fS);
         return 0;
      case 1: // PARameter <parno> <value> 
         {  
         if (nargs<2 && setCommand) return -1;
         Int_t parnum;
         Double_t val;
         if(setCommand) {
            parnum = Int_t(fCmPar[0])-1;
            val= fCmPar[1];
            if(parnum<0 || parnum>=fNpar) return -2; //no such parameter
            fA[parnum] = val;
         } else {
            if (nargs>0) {
               parnum = Int_t(fCmPar[0])-1;
               if(parnum<0 || parnum>=fNpar) return -2; //no such parameter
               Printf("Parameter %s = %E",fANames[parnum].Data(),fA[parnum]);
            } else
               for (i=0;i<fNpar;i++)
               Printf("Parameter %s = %E",fANames[i].Data(),fA[i]);
        
         }
         return 0;
         }
      case 2: // LIMits [parno] [ <lolim> <uplim> ]
         {
         Int_t parnum;
         Double_t lolim,uplim;
         if (nargs<1) {
            for(i=0;i<fNpar;i++) 
               if(setCommand) {
                  fAMN[i] = gMINDOUBLE;
                  fAMX[i] = gMAXDOUBLE;
               } else 
                  Printf("Limits for param %s: Low=%E, High=%E",
                         fANames[i].Data(),fAMN[i],fAMX[i]);
         } else {
            parnum = Int_t(fCmPar[0])-1;
            if(parnum<0 || parnum>=fNpar)return -1;
            if(setCommand) {
               if(nargs>2) {
                  lolim = fCmPar[1];
                  uplim = fCmPar[2];
                  if(uplim==lolim) return -1;
                  if(lolim>uplim) {
                     Double_t tmp = lolim;
                     lolim = uplim;
                     uplim = tmp;
                  }
               } else {
                  lolim = gMINDOUBLE;
                  uplim = gMAXDOUBLE;
               }
               fAMN[parnum] = lolim;
               fAMX[parnum] = uplim;
            } else 
               Printf("Limits for param %s Low=%E, High=%E",
                    fANames[parnum].Data(),fAMN[parnum],fAMX[parnum]);
            }
            return 0;   
         }
      case 3:
         {
         if(setCommand) return 0;
         Printf("\nCovariant matrix ");
         Int_t l = 0,nn=0,nnn=0;
         for (i=0;i<fNpar;i++) if(fPL0[i]>0.) nn++;
         for (i=0;i<nn;i++) {
            for(;fPL0[nnn]<=0.;nnn++) { }
            printf("%5s: ",fANames[nnn++].Data());
            for (Int_t j=0;j<=i;j++) 
               printf("%11.2E",fZ[l++]);
            cout<<endl;
         }
         cout<<endl;
         return 0;
         }
      case 4:
         if(setCommand) return 0;
         Printf("\nGlobal correlation factors (maximum correlation of the parameter\n  with arbitrary linear combination of other parameters)");
         for(i=0;i<fNpar;i++) {
            printf("%5s: ",fANames[i].Data());
            printf("%11.3E\n",TMath::Sqrt(1-1/((fR[i]!=0.)?fR[i]:1.)) );
         }
         cout<<endl;
         return 0;
      case 5:   // PRIntout not implemented
         return -10;
      case 6: // NOGradient
         if(!setCommand) return 0;
         fGRAD = false;
         return 0;
      case 7: // GRAdient
         if(!setCommand) return 0;
         fGRAD = true;
         return 0;
      case 8: // ERRordef - now ignored
         return 0;
      case 9: // INPut - not implemented
         return -10;
      case 10: // WIDthpage - not implemented
         return -10;
      case 11: // LINesperpage - not implemented
         return -10;
      case 12: //NOWarnings
         if(!setCommand) return 0;
         fWARN = false;
         return 0;
      case 13: // WARnings 
         if(!setCommand) return 0;
         fWARN = true;
         return 0;
      case 14: // RANdomgenerator - not implemented
         return -10;
      case 15: // TITle - ignored
         return 0; 
      case 16: // STRategy - ignored
         return 0;
      case 17: // EIGenvalues - not implemented
         return -10;
      case 18: // PAGethrow - ignored
         return 0;
      case 19: // MINos errors - not implemented
         return -10;
      case 20: //EPSmachine
         if(!setCommand) {
            Printf("Relative floating point presicion RP=%E",fRP);
         } else 
            if (nargs>0) {
               Double_t pres=fCmPar[0];
               if (pres<1e-5 && pres>1e-34) fRP=pres;
            }
         return 0;
      case 21: // OUTputfile - not implemented
         return -10;
      case 22: // BATch - ignored
         return 0;
      case 23: // INTerative - ignored
         return 0;
      case 24: // VERsion
         if(setCommand) return 0;
         Printf("FUMILI-ROOT version 0.1");
         return 0;
      case 25: // reserved
         return 0;
      case 26: // NODebug
         if(!setCommand) return 0;
         fDEBUG = false;
         return 0;
      case 27: // DEBug
         if(!setCommand) return 0;
         fDEBUG = true;
         return 0;
      case 28:
      case 29:
         return -3;
      default:
         break;
   }
   return -3;
}

//______________________________________________________________________________
void TFumili::FixParameter(Int_t ipar) { 
   // Fixes parameter number ipar

   if(ipar>=0 && ipar<fNpar && fPL0[ipar]>0.) {
      fPL0[ipar] = -fPL0[ipar]; 
      fLastFixed = ipar;
   }
}

//______________________________________________________________________________
Double_t *TFumili::GetCovarianceMatrix() const
{
   // return a pointer to the covariance matrix

   return fZ;

}

//______________________________________________________________________________
Double_t TFumili::GetCovarianceMatrixElement(Int_t i, Int_t j) const
{
   // return element i,j from the covariance matrix

   if (!fZ) return 0;
   if (i < 0 || i >= fNpar || j < 0 || j >= fNpar) {
      Error("GetCovarianceMatrixElement","Illegal arguments i=%d, j=%d",i,j);
      return 0;
   }
   return fZ[j+fNpar*i];
}


//______________________________________________________________________________
Int_t TFumili::GetNumberTotalParameters() const
{
   // return the total number of parameters (free + fixed)

   return fNpar;
}

//______________________________________________________________________________
Int_t TFumili::GetNumberFreeParameters() const
{
   // return the number of free parameters

   Int_t nfree = fNpar;
   for (Int_t i=0;i<fNpar;i++) {
      if (IsFixed(i)) nfree--;
   }
   return nfree;
}

//______________________________________________________________________________
Double_t TFumili::GetParError(Int_t ipar) const
{
   // return error of parameter ipar

   if (ipar<0 || ipar>=fNpar) return 0;
   else return fParamError[ipar];
}

//______________________________________________________________________________
Double_t TFumili::GetParameter(Int_t ipar) const
{
   // return current value of parameter ipar

   if (ipar<0 || ipar>=fNpar) return 0;
   else return fA[ipar];
}


//______________________________________________________________________________
Int_t TFumili::GetParameter(Int_t ipar,char *cname,Double_t &value,Double_t &verr,Double_t &vlow, Double_t &vhigh) const
{
   // Get various ipar parameter attributs:
   // 
   // cname:    parameter name
   // value:    parameter value
   // verr:     parameter error
   // vlow:     lower limit
   // vhigh:    upper limit
   // WARNING! parname must be suitably dimensionned in the calling function.

   if (ipar<0 || ipar>=fNpar) {
      value = 0;
      verr  = 0;
      vlow  = 0;
      vhigh = 0;
      return -1;
   }
   strcpy(cname,fANames[ipar].Data());
   value = fA[ipar];
   verr  = fParamError[ipar];
   vlow  = fAMN[ipar];
   vhigh = fAMX[ipar];
   return 0;
}

//______________________________________________________________________________
const char *TFumili::GetParName(Int_t ipar) const
{
   // return name of parameter ipar
   
   if (ipar < 0 || ipar > fNpar) return "";
   return fANames[ipar].Data();
}

//______________________________________________________________________________
Int_t TFumili::GetErrors(Int_t ipar,Double_t &eplus, Double_t &eminus, Double_t &eparab, Double_t &globcc) const
{
   // Return errors after MINOs
   // not implemented
   eparab = 0;
   globcc = 0;  
   if (ipar<0 || ipar>=fNpar) {
      eplus  = 0;
      eminus = 0;
      return -1;
   }
   eplus=fParamError[ipar];
   eminus=-eplus;
   return 0;
}

//______________________________________________________________________________
Int_t TFumili::GetStats(Double_t &amin, Double_t &edm, Double_t &errdef, Int_t &nvpar, Int_t &nparx) const
{
   // return global fit parameters
   //   amin     : chisquare
   //   edm      : estimated distance to minimum
   //   errdef
   //   nvpar    : number of variable parameters
   //   nparx    : total number of parameters
   amin   = 2*fS;
   edm    = fGT; // 
   errdef = 0; // ??
   nparx  = fNpar;
   nvpar  = 0;
   for(Int_t ii=0; ii<fNpar; ii++) {
      if(fPL0[ii]>0.) nvpar++;
   }  
   return 0;
}



//______________________________________________________________________________
Double_t TFumili::GetSumLog(Int_t n)
{
   // return Sum(log(i) i=0,n
   // used by log likelihood fits

   if (n < 0) return 0;
   if (n > fNlog) {
      if (fSumLog) delete [] fSumLog;
      fNlog = 2*n+1000;
      fSumLog = new Double_t[fNlog+1];
      Double_t fobs = 0;
      for (Int_t j=0;j<=fNlog;j++) {
         if (j > 1) fobs += TMath::Log(j);
         fSumLog[j] = fobs;
      }
   }
   if (fSumLog) return fSumLog[n];
   return 0;
}



//______________________________________________________________________________
void TFumili::InvertZ(Int_t n)
{
   // Inverts packed diagonal matrix Z by square-root method.
   //  Matrix elements corresponding to 
   // fix parameters are removed.
   //
   // n: number of variable parameters
   //
   static Double_t am = 3.4e138;
   static Double_t rp = 5.0e-14;
   Double_t  ap, aps, c, d;
   Double_t *r_1=fR;
   Double_t *pl_1=fPL;
   Double_t *z_1=fZ;
   Int_t i, k, l, ii, ki, li, kk, ni, ll, nk, nl, ir, lk;
   if (n < 1) {
      return;
   }
   --pl_1;
   --r_1;
   --z_1;
   aps = am / n;
   aps = sqrt(aps);
   ap = 1.0e0 / (aps * aps);
   ir = 0;
   for (i = 1; i <= n; ++i) {
      L1:
         ++ir;
         if (pl_1[ir] <= 0.0e0) goto L1;
         else                   goto L2;
      L2:
         ni = i * (i - 1) / 2;
         ii = ni + i;
         k = n + 1;
         if (z_1[ii] <= rp * TMath::Abs(r_1[ir]) || z_1[ii] <= ap) {
            goto L19;
         }
         z_1[ii] = 1.0e0 / sqrt(z_1[ii]);
         nl = ii - 1;
      L3:
         if (nl - ni <= 0) goto L5;
         else              goto L4;
         L4:
         z_1[nl] *= z_1[ii];
         if (TMath::Abs(z_1[nl]) >= aps) {
            goto L16;
         }
         --nl;
         goto L3;
      L5:
         if (i - n >= 0) goto L12;
         else            goto L6;
      L6:
         --k;
         nk = k * (k - 1) / 2;
         nl = nk;
         kk = nk + i;
         d = z_1[kk] * z_1[ii];
         c = d * z_1[ii];
         l = k;
      L7:
         ll = nk + l;
         li = nl + i;
         z_1[ll] -= z_1[li] * c;
         --l;
         nl -= l;
         if (l - i <= 0) goto L9;
         else            goto L7;
      L8:
         ll = nk + l;
         li = ni + l;
         z_1[ll] -= z_1[li] * d;
      L9:
         --l;
         if (l <= 0) goto L10;
         else        goto L8;
      L10:
         z_1[kk] = -c;
         if (k - i - 1 <= 0) goto L11;
         else                goto L6;
      L11:
         ;
   }
   L12:
      for (i = 1; i <= n; ++i) {
         for (k = i; k <= n; ++k) {
            nl = k * (k - 1) / 2;
            ki = nl + i;
            d = 0.0e0;
            for (l = k; l <= n; ++l) {
               li = nl + i;
               lk = nl + k;
               d += z_1[li] * z_1[lk];
               nl += l;
            }
            ki = k * (k - 1) / 2 + i;
            z_1[ki] = d;
         }
      }
   L15:
      return;
   L16:
      k = i + nl - ii;
      ir = 0;
      for (i = 1; i <= k; ++i) {
         L17:
            ++ir;
            if (pl_1[ir] <= 0.0e0) {
               goto L17;
            }
      }
   L19:
      pl_1[ir] = -2.0e0;
      r_1[ir] = 0.0e0;
      fINDFLG[0] = ir - 1;
      goto L15;
}




//______________________________________________________________________________
Bool_t TFumili::IsFixed(Int_t ipar) const
{
   //return kTRUE if parameter ipar is fixed, kFALSE othersise)
   
   if(ipar < 0 || ipar >= fNpar) {
      Warning("IsFixed","Illegal parameter number :%d",ipar);
      return kFALSE;
   }
   if (fPL0[ipar] < 0) return kTRUE;
   else                return kFALSE;
}


//______________________________________________________________________________
Int_t TFumili::Minimize()
{// Main minimization procedure
   //*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*//
   //         FUMILI  
   //  Based on ideas, proposed by I.N. Silin
   //    [See NIM A440, 2000 (p431)]
   // converted from FORTRAN to C  by
   //     Sergey Yaschenko <s.yaschenko@fz-juelich.de>
   //
   //*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*//
   //
   // This function is called after setting theoretical function 
   // by means of TFumili::SetUserFunc and initializing parameters.
   // Optionally one can set FCN function (see TFumili::SetFCN and TFumili::Eval)
   // If FCN is undefined then user has to provide data arrays by calling
   //  TFumili::SetData procedure.
   //
   // TFumili::Minimize return following values:
   //    0  - fit is converged
   //   -2  - function is not decreasing (or bad derivatives)
   //   -3  - error estimations are infinite
   //   -4  - maximum number of iterations is exceeded
   //
   Int_t i;
   // Flag3 - is fit is chi2 or likelihood? 0 - chi2, 1 - likelihood
   fINDFLG[2]=0;
   //
   // Are the parameters outside of the boundaries ?
   //
   Int_t parn;

   if(fFCN) {
      Eval(parn,fGr,fS,fA,9); fNfcn++;
   }
   for( i = 0; i < fNpar; i++) {
      if(fA[i] > fAMX[i]) fA[i] = fAMX[i];
      if(fA[i] < fAMN[i]) fA[i] = fAMN[i];
   }

   Int_t nn2, n, fixFLG,  ifix1, fi, nn3, nn1, n0;
   Double_t t1;
   Double_t sp, t, olds=0;
   Double_t bi, aiMAX=0, amb;
   Double_t afix, sigi, akap;
   Double_t alambd, al, bm, abi, abm;
   Int_t l1, k, ifix;
  
   nn2=0;

   // Number of parameters;
   n=fNpar;
   fixFLG=0;

   // Exit flag
   fENDFLG=0;

   // Flag2
   fINDFLG[1] = 0;
   ifix1=-1;
   fi=0;
   nn3=0;
  
   // Initialize param.step limits
   for( i=0; i < n; i++) {
      fR[i]=0.;
      if ( fEPS > 0.) fParamError[i] = 0.;
      fPL[i] = fPL0[i];
   }

L3: // Start Iteration
 
   nn1 = 1;
   t1 = 1.;
 
L4: // New iteration
 
   // fS - objective function value - zero first
   fS = 0.;
   // n0 - number of variable parameters in fit
   n0 = 0;
   for( i = 0; i < n; i++) {
      fGr[i]=0.; // zero gradients
      if (fPL0[i] > .0) {
         n0=n0+1; 
         // new iteration - new parallelepiped
         if (fPL[i] > .0) fPL0[i]=fPL[i];
      }
   }
   Int_t nn0;
   // Calculate number of fZ-matrix elements as nn0=1+2+..+n0 
   nn0 = n0*(n0+1)/2;
   // if (nn0 >= 1) ????
   // fZ-matrix is initialized
   for( i=0; i < nn0; i++) fZ[i]=0.;

   // Flag1
   fINDFLG[0] = 0;
   Int_t ijkl=1;

   // Calculate fS - objective function, fGr - gradients, fZ - fZ-matrix
   if(fFCN) {
      Eval(parn,fGr,fS,fA,2);
      fNfcn++;
   } else
      ijkl = SGZ();
   if(!ijkl) return 10; 
   if (ijkl == -1) fINDFLG[0]=1;

   // sp - scaled on fS machine precision
   sp=fRP*TMath::Abs(fS);

   // save fZ-matrix
   for( i=0; i < nn0; i++) fZ0[i] = fZ[i];
   if (nn3 > 0) {
      if (nn1 <= fNstepDec) {
         t=2.*(fS-olds-fGT);
         if (fINDFLG[0] == 0) {
            if (TMath::Abs(fS-olds) <= sp && -fGT <= sp) goto L19;
            if(        0.59*t < -fGT) goto L19;
            t = -fGT/t;
            if (t < 0.25 ) t = 0.25;
         }
         else   t = 0.25;
         fGT = fGT*t;
         t1 = t1*t;
         nn2=0;
         for( i = 0; i < n; i++) {
            if (fPL[i] > 0.) {
               fA[i]=fA[i]-fDA[i];
               fPL[i]=fPL[i]*t;
               fDA[i]=fDA[i]*t;
               fA[i]=fA[i]+fDA[i];
            }
         }
         nn1=nn1+1;
         goto L4;
      }
   }
 
L19:
 
   if(fINDFLG[0] != 0) {
      fENDFLG=-4;
      printf("trying to execute an illegal junp at L85\n");
      //goto L85;
   } 

 
   Int_t k1, k2, i1, j, l;
   k1 = 1;
   k2 = 1;
   i1 = 1;
   // In this cycle we removed from fZ contributions from fixed parameters
   // We'll get fixed parameters after boudary check
   for( i = 0; i < n; i++) {
      if (fPL0[i] > .0) { 
         // if parameter was fixed - release it
         if (fPL[i] == 0.) fPL[i]=fPL0[i];
         if (fPL[i] > .0) { // ??? it is already non-zero
            // if derivative is negative and we above maximum
            // or vice versa then fix parameter again and increment k1 by i1
            if ((fA[i] >= fAMX[i] && fGr[i] < 0.) ||
                   (fA[i] <= fAMN[i] && fGr[i] > 0.)) {
               fPL[i] = 0.;
               k1 = k1 + i1; // i1 stands for fZ-matrix row-number multiplier
               ///  - skip this row
               //  in case we are fixing parameter number i
            } else {
               for( j=0; j <= i; j++) {// cycle on columns of fZ-matrix
                  if (fPL0[j] > .0) {
                     // if parameter is not fixed then fZ = fZ0 
                     // Now matrix fZ of other dimension
                     if (fPL[j] > .0) {
                        fZ[k2 -1] = fZ0[k1 -1];
                        k2=k2+1;
                     }  
                     k1=k1+1;
                  }
               }
            }  
         }  
         else k1 = k1 + i1; // In case of negative fPL[i] - after mconvd
         i1=i1+1;  // Next row of fZ0
      }
   }

   // INVERT fZ-matrix (mconvd() procedure)
   i1 = 1;
   l  = 1;
   for( i = 0; i < n; i++) {// extract diagonal elements to fR-vector
      if (fPL[i] > .0) {
         fR[i] = fZ[l - 1];
         i1 = i1+1;
         l = l + i1;
      }
   }

   n0 = i1 - 1;
   InvertZ(n0);

   // fZ matrix now is inversed
   if (fINDFLG[0] != 0) { // problems
      // some PLs now have negative values, try to reduce fZ-matrix again
      fINDFLG[0] = 0;
      fINDFLG[1] = 1; // errors can be infinite
      fixFLG = fixFLG + 1;
      fi = 0;
      goto L19;
   }

   // ... CALCULATE THEORETICAL STEP TO MINIMUM
   i1 = 1;
   for( i = 0; i < n; i++) {
      fDA[i]=0.; // initial step is zero
      if (fPL[i] > .0) {   // for non-fixed parameters
         l1=1;
         for( l = 0; l < n; l++) {
            if (fPL[l] > .0) {
               // Caluclate offset of Z^-1(i1,l1) element in packed matrix
               // because we skip fixed param numbers we need also i,l
               if (i1 <= l1 ) k=l1*(l1-1)/2+i1;
               else k=i1*(i1-1)/2+l1;
               // dA_i = \sum (-Z^{-1}_{il}*grad(fS)_l)
               fDA[i]=fDA[i]-fGr[l]*fZ[k - 1];
               l1=l1+1;
            }
         }
         i1=i1+1;
      }
   }
   //          ... CHECK FOR PARAMETERS ON BOUNDARY

   afix=0.;
   ifix = -1;
   i1 = 1;
   l = i1;
   for( i = 0; i < n; i++)
      if (fPL[i] > .0) {
         sigi = TMath::Sqrt(TMath::Abs(fZ[l - 1])); // calculate \sqrt{Z^{-1}_{ii}} 
         fR[i] = fR[i]*fZ[l - 1];      // Z_ii * Z^-1_ii
         if (fEPS > .0) fParamError[i]=sigi;
         if ((fA[i] >= fAMX[i] && fDA[i] > 0.) || (fA[i] <= fAMN[i]
                                               && fDA[i] < .0)) {
            // if parameter out of bounds and if step is making things worse
      
            akap = TMath::Abs(fDA[i]/sigi);
            // let's found maximum of dA/sigi - the worst of parameter steps
            if (akap > afix) {
               afix=akap;
               ifix=i;
               ifix1=i;
            }
         }
         i1=i1+1;
         l=l+i1;
      }
   if (ifix != -1) {
      // so the worst parameter is found - fix it and exclude,
      //  reduce fZ-matrix again
      fPL[ifix] = -1.;
      fixFLG = fixFLG + 1;
      fi = 0;
      //.. REPEAT CALCULATION OF THEORETICAL STEP AFTER fiXING EACH PARAMETER
      goto L19;
   }

   //... CALCULATE STEP CORRECTION FACTOR

   alambd = 1.;
   fAKAPPA = 0.;
   Int_t imax;
   imax = -1;


   for( i = 0; i < n; i++) {
      if (fPL[i] > .0) {
         bm = fAMX[i] - fA[i];  
         abi = fA[i] + fPL[i]; // upper  parameter limit
         abm = fAMX[i];
         if (fDA[i] <= .0) {
            bm = fA[i] - fAMN[i];
            abi = fA[i] - fPL[i]; // lower parameter limit
            abm = fAMN[i];
         }
         bi = fPL[i];
         // if parallelepiped boundary is crossing limits
         // then reduce it (deforming)
         if ( bi > bm) {
            bi = bm;
            abi = abm;
         }
         // if calculated step is out of bounds
         if ( TMath::Abs(fDA[i]) > bi) {
            // derease step splitter alambdA if needed
            al = TMath::Abs(bi/fDA[i]);
            if (alambd > al) {
               imax=i;
               aiMAX=abi;
               alambd=al;
            }
         }
         // fAKAPPA - parameter will be <fEPS if fit is converged
         akap = TMath::Abs(fDA[i]/fParamError[i]); 
         if (akap > fAKAPPA) fAKAPPA=akap;
      }
   }
   //... CALCULATE NEW CORRECTED STEP
   fGT = 0.;
   amb = 1.e18;
   // alambd - multiplier to split teoretical step dA
   if (alambd > .0) amb = 0.25/alambd;
   for( i = 0; i < n; i++) {
      if (fPL[i] > .0) {
         if (nn2 > fNlimMul ) {
            if (TMath::Abs(fDA[i]/fPL[i]) > amb ) {
               fPL[i] = 4.*fPL[i]; // increase parallelepiped
               t1=4.; // flag - that fPL was increased
            }
         }
         // cut step
         fDA[i] = fDA[i]*alambd;
         // expected function value change in next iteration
         fGT = fGT + fDA[i]*fGr[i];
      }
   }

   //.. CHECK IF MINIMUM ATTAINED AND SET EXIT MODE
   // if expected fGT smaller than precision
   // and other stuff
   if (-fGT <= sp && t1 < 1. && alambd < 1.)fENDFLG = -1; // function is not decreasing 

   if (fENDFLG >= 0) {
      if (fAKAPPA < TMath::Abs(fEPS)) { // fit is converging
         if (fixFLG == 0) 
            fENDFLG=1; // successful fit
         else {// we have fixed parameters
            if (fENDFLG == 0) {
               //... CHECK IF fiXING ON BOUND IS CORRECT
               fENDFLG = 1;
               fixFLG = 0;
               ifix1=-1;
               // release fixed parameters
               for( i = 0; i < fNpar; i++) fPL[i] = fPL0[i];
               fINDFLG[1] = 0;
               // and repeat iteration
               goto L19;
            } else {
               if( ifix1 >= 0) {
                  fi = fi + 1;
                  fENDFLG = 0;
               }
            }
         }
      } else { // fit does not converge
         if( fixFLG != 0) {
            if( fi > fixFLG ) {
               //... CHECK IF fiXING ON BOUND IS CORRECT
               fENDFLG = 1;
               fixFLG = 0;
               ifix1=-1;
               for( i = 0; i < fNpar; i++) fPL[i] = fPL0[i];
               fINDFLG[1] = 0;
               goto L19;
            } else {
               fi = fi + 1;
               fENDFLG = 0;
            }
         } else {
            fi = fi + 1;
            fENDFLG = 0;
         }
      }
   }

// L85:
   // iteration number limit is exceeded
   if(fENDFLG == 0 && nn3 >= fNmaxIter) fENDFLG=-3;
  
   // fit errors are infinite;
   if(fENDFLG > 0 && fINDFLG[1] > 0) fENDFLG=-2;
  
   //MONITO (fS,fNpar,nn3,IT,fEPS,fGT,fAKAPPA,alambd);
   if (fENDFLG == 0) { // make step
      for ( i = 0; i < n; i++) fA[i] = fA[i] + fDA[i];
      if (imax >= 0) fA[imax] = aiMAX;
      olds=fS;
      nn2=nn2+1;
      nn3=nn3+1;
   } else { 
      // fill covariant matrix VL
      // fill parameter error matrix up
      Int_t il;
      il = 0;
      for( Int_t ip = 0; ip < fNpar; ip++) { 
         if( fPL0[ip] > .0) {
            for( Int_t jp = 0; jp <= ip; jp++) {
               if(fPL0[jp] > .0) {
                  //         VL[ind(ip,jp)] = fZ[il];
                  il = il + 1;
               }
            }
         }
      }
      return fENDFLG - 1;
   }
   goto L3;
}

//______________________________________________________________________________
void TFumili::PrintResults(Int_t ikode,Double_t p) const
{
   // Prints fit results. 
   //  
   // ikode is the type of printing parameters
   // p is function value
   //
   //  ikode = 1   - print values, errors and limits
   //  ikode = 2   - print values, errors and steps
   //  ikode = 3   - print values, errors, steps and derivatives
   //  ikode = 4   - print only values and errors
   //
   TString exitStatus="";
   TString xsexpl="";
   TString colhdu[3],colhdl[3],cx2,cx3;
   switch (fENDFLG) {
   case 1:
      exitStatus="CONVERGED";
      break;
   case -1:
      exitStatus="CONST FCN";
      xsexpl="****\n* FUNCTION IS NOT DECREASING OR BAD DERIVATIVES\n****";
      break;
   case -2:
      exitStatus="ERRORS INF";
      xsexpl="****\n* ESTIMATED ERRORS ARE INfiNITE\n****";
      break;
   case -3:
      exitStatus="MAX ITER.";
      xsexpl="****\n* MAXIMUM NUMBER OF ITERATIONS IS EXCEEDED\n****";
      break;
   case -4:
      exitStatus="ZERO PROBAB";
      xsexpl="****\n* PROBABILITY OF LIKLIHOOD FUNCTION IS NEGATIVE OR ZERO\n****";
      break;
   default:
      exitStatus="UNDEfiNED";
      xsexpl="****\n* fiT IS IN PROGRESS\n****";
      break;
   }
   if (ikode == 1) {
      colhdu[0] = "              ";
      colhdl[0] = "      ERROR   ";
      colhdu[1] = "      PHYSICAL";
      colhdu[2] = " LIMITS       ";
      colhdl[1] = "    NEGATIVE  ";
      colhdl[2] = "    POSITIVE  ";
   }
   if (ikode == 2) {
      colhdu[0] = "              ";
      colhdl[0] = "      ERROR   ";
      colhdu[1] = "    INTERNAL  ";
      colhdl[1] = "    STEP SIZE ";
      colhdu[2] = "    INTERNAL  ";
      colhdl[2] = "      VALUE   ";
   }
   if (ikode == 3) {
      colhdu[0] = "              ";
      colhdl[0] = "      ERROR   ";
      colhdu[1] = "       STEP   ";
      colhdl[1] = "       SIZE   ";
      colhdu[2] = "       fiRST  ";
      colhdl[2] = "    DERIVATIVE";
   }
   if (ikode == 4) {
      colhdu[0] = "    PARABOLIC ";
      colhdl[0] = "      ERROR   ";
      colhdu[1] = "        MINOS ";
      colhdu[2] = "ERRORS        ";
      colhdl[1] = "   NEGATIVE   ";
      colhdl[2] = "   POSITIVE   ";
   }
   if(fENDFLG<1)Printf((const char*)xsexpl.Data());
   Printf(" FCN=%g FROM FUMILI  STATUS=%-10s %9d CALLS OF FCN",
         p,exitStatus.Data(),fNfcn);
   Printf(" EDM=%g ",-fGT);
   Printf("  EXT PARAMETER              %-14s%-14s%-14s",
         (const char*)colhdu[0].Data()
         ,(const char*)colhdu[1].Data()
         ,(const char*)colhdu[2].Data());
   Printf("  NO.   NAME          VALUE  %-14s%-14s%-14s",
         (const char*)colhdl[0].Data()
         ,(const char*)colhdl[1].Data()
         ,(const char*)colhdl[2].Data());

   for (Int_t i=0;i<fNpar;i++) { 
      if (ikode==3) { 
         cx2 = Form("%14.5e",fDA[i]);
         cx3 = Form("%14.5e",fGr[i]);

      }
      if (ikode==1) {
         cx2 = Form("%14.5e",fAMN[i]);
         cx3 = Form("%14.5e",fAMX[i]);
      }
      if (ikode==2) {
         cx2 = Form("%14.5e",fDA[i]);
         cx3 = Form("%14.5e",fA[i]);
      }
      if(ikode==4){
         cx2 = " *undefined*  ";
         cx3 = " *undefined*  ";
      }
      if(fPL0[i]<=0.) { cx2="    *fixed*   ";cx3=""; }
      Printf("%4d %-11s%14.5e%14.5e%-14s%-14s",i+1
           ,fANames[i].Data(),fA[i],fParamError[i]
           ,cx2.Data(),cx3.Data());
   }
}


//______________________________________________________________________________
void TFumili::ReleaseParameter(Int_t ipar) {
   // Releases parameter number ipar

   if(ipar>=0 && ipar<fNpar && fPL0[ipar]<=0.) {
      fPL0[ipar] = -fPL0[ipar]; 
      if (fPL0[ipar] == 0. || fPL0[ipar]>=1.) fPL0[ipar]=.1;
   }
}


//______________________________________________________________________________
void TFumili::SetData(Double_t *exdata,Int_t numpoints,Int_t vecsize){
   // Sets pointer to data array provided by user. 
   // Necessary if SetFCN is not called.
   // 
   // numpoints:    number of experimental points
   // vecsize:      size of data point vector + 2 
   //               (for N-dimensional fit vecsize=N+2)
   // exdata:       data array with following format
   //
   //   exdata[0] = ExpValue_0     - experimental data value number 0
   //   exdata[1] = ExpSigma_0     - error of value number 0
   //   exdata[2] = X_0[0]        
   //   exdata[3] = X_0[1]
   //       .........
   //   exdata[vecsize-1] = X_0[vecsize-3]
   //   exdata[vecsize]   = ExpValue_1
   //   exdata[vecsize+1] = ExpSigma_1
   //   exdata[vecsize+2] = X_1[0]
   //       .........
   //   exdata[vecsize*(numpoints-1)] = ExpValue_(numpoints-1)
   //       .........
   //   exdata[vecsize*numpoints-1] = X_(numpoints-1)[vecsize-3]
   //
  
   if(exdata){
      fNED1 = numpoints;
      fNED2 = vecsize;
      fEXDA = exdata;
   }
}


//______________________________________________________________________________
void TFumili::SetFitMethod(const char *name)
{
   // ret fit method (chisquare or loglikelihood)
   
   if (!strcmp(name,"H1FitChisquare"))    SetFCN(H1FitChisquareFumili);
   if (!strcmp(name,"H1FitLikelihood"))   SetFCN(H1FitLikelihoodFumili);
   if (!strcmp(name,"GraphFitChisquare")) SetFCN(GraphFitChisquareFumili);
}


//______________________________________________________________________________
Int_t TFumili::SetParameter(Int_t ipar,const char *parname,Double_t value,Double_t verr,Double_t vlow, Double_t vhigh) {
   // Sets for prameter number ipar initial parameter value, 
   // name parname, initial error verr and limits vlow and vhigh
   // If vlow = vhigh but not equil to zero, parameter will be fixed.
   // If vlow = vhigh = 0, parameter is released and its limits are discarded
   //
   if (ipar<0 || ipar>=fNpar) return -1;
   fANames[ipar] = parname;
   fA[ipar] = value; 
   fParamError[ipar] = verr;
   if(vlow<vhigh) {
      fAMN[ipar] = vlow;
      fAMX[ipar] = vhigh;
   } else {
      if(vhigh<vlow) {
         fAMN[ipar] = vhigh;
         fAMX[ipar] = vlow;
      }
      if(vhigh==vlow) {
         if(vhigh==0.) {
            ReleaseParameter(ipar);
            fAMN[ipar] = gMINDOUBLE;
            fAMX[ipar] = gMAXDOUBLE;
         }
         if(vlow!=0) FixParameter(ipar);
      }
   }
   return 0;
}

//______________________________________________________________________________
Int_t TFumili::SGZ()
{
   //  Evaluates objective function ( chi-square ), gradients and   
   //  Z-matrix using data provided by user via TFumili::SetData
   //
   fS = 0.;
   Int_t i,j,l,k2=1,k1,ki=0;
   Double_t *x  = new Double_t[fNED2];
   Double_t *df = new Double_t[fNpar];
   Int_t nx = fNED2-2;
   for (l=0;l<fNED1;l++) { // cycle on all exp. points
      k1 = k2;
      if (fLogLike) {
         fNumericDerivatives = kTRUE;
         nx  = fNED2;
         k1 -= 2;
      }
  
      for (i=0;i<nx;i++){
         ki  += 1+i;
         x[i] = fEXDA[ki];
      }
      //  Double_t y = ARITHM(df,x);
      Double_t y = EvalTFN(df,x);
      if(fNumericDerivatives) Derivatives(df,x);
      Double_t sig=1.;
      if(fLogLike) { // Likelihood method
         if(y>0.) {
            fS = fS - log(y);
            y  = -y;
            sig= y;
         } else { // 
            delete [] x;
            delete [] df;
            fS = 1e10;
            return -1; // indflg[0] = 1;
         }
      } else { // Chi2 method
         sig = fEXDA[k2]; // sigma of experimental point
         y = y - fEXDA[k1-1]; // f(x_i) - F_i
         fS = fS + (y*y/(sig*sig))*.5; // simple chi2/2
      }
      Int_t n = 0;
      for (i=0;i<fNpar;i++) {
         if (fPL0[i]>0){
            df[n]   = df[i]/sig; // left only non-fixed param derivatives div by Sig
            fGr[i] += df[n]*(y/sig);
            n++;
         }
      }
      l = 0;
      for (i=0;i<n;i++)
         for (j=0;j<=i;j++) 
            fZ[l++] += df[i]*df[j];
      k2 += fNED2;
   }
 
   delete[] df;
   delete[] x;
   return 1;
}


//______________________________________________________________________________
void TFumili::FitChisquare(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag)
{
   //  Minimization function for H1s using a Chisquare method
   //  Default method (function evaluated at center of bin)
   //  for each point the cache contains the following info
   //    -1D : bc,e,xc  (bin content, error, x of center of bin)
   //    -2D : bc,e,xc,yc
   //    -3D : bc,e,xc,yc,zc

   Foption_t fitOption = GetFitOption();
   if (fitOption.Integral) {
      FitChisquareI(npar,gin,f,u,flag);
      return;
   }
   Double_t cu,eu,fu,fsum;
   Double_t x[3];
   Double_t *zik=0;
   Double_t *pl0=0;

   TH1 *hfit = (TH1*)GetObjectFit();
   TF1 *f1   = (TF1*)GetUserFunc();
   Int_t nd  = hfit->GetDimension();
   Int_t j;
   
   npar = f1->GetNpar();
   SetParNumber(npar);
   if(flag == 9) return;
   zik = GetZ();
   pl0 = GetPL0();

   Double_t *df = new Double_t[npar];
   f1->InitArgs(x,u);
   f = 0;
   
   Int_t npfit = 0;
   Double_t *cache = fCache;
   for (Int_t i=0;i<fNpoints;i++) {
      if (nd > 2) x[2]     = cache[4];
      if (nd > 1) x[1]     = cache[3];
      x[0]     = cache[2];
      cu  = cache[0];
      TF1::RejectPoint(kFALSE);
      fu = f1->EvalPar(x,u);
      if (TF1::RejectedPoint()) {cache += fPointSize; continue;}
      eu = cache[1];
      Derivatives(df,x);
      Int_t n = 0;
      fsum = (fu-cu)/eu;
      if (flag!=1) {
         for (j=0;j<npar;j++) {
            if (pl0[j]>0) {
               df[n] = df[j]/eu; 
               // left only non-fixed param derivatives / by Sigma
               gin[j] += df[n]*fsum;
               n++;
            }
         }
         Int_t l = 0;
         for (j=0;j<n;j++)
            for (Int_t k=0;k<=j;k++) 
               zik[l++] += df[j]*df[k];
      }
      f += .5*fsum*fsum;
      npfit++;
      cache += fPointSize;
   }
   f1->SetNumberFitPoints(npfit);
   delete [] df;
}


//______________________________________________________________________________
void TFumili::FitChisquareI(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag)
{
   //  Minimization function for H1s using a Chisquare method
   //  The "I"ntegral method is used
   //  for each point the cache contains the following info
   //    -1D : bc,e,xc,xw  (bin content, error, x of center of bin, x bin width of bin)
   //    -2D : bc,e,xc,xw,yc,yw
   //    -3D : bc,e,xc,xw,yc,yw,zc,zw

   Double_t cu,eu,fu,fsum;
   Double_t x[3];
   Double_t *zik=0;
   Double_t *pl0=0;

   TH1 *hfit = (TH1*)GetObjectFit();
   TF1 *f1   = (TF1*)GetUserFunc();
   Int_t nd  = hfit->GetDimension();
   Int_t j;
   
   f1->InitArgs(x,u);
   npar = f1->GetNpar();
   SetParNumber(npar);
   if(flag == 9) return;
   zik = GetZ();
   pl0 = GetPL0();

   Double_t *df=new Double_t[npar];
   f = 0;
   
   Int_t npfit = 0;
   Double_t *cache = fCache;
   for (Int_t i=0;i<fNpoints;i++) {
      cu  = cache[0];
      TF1::RejectPoint(kFALSE);
      f1->SetParameters(u);
      if (nd < 2) {
         fu = f1->Integral(cache[2] - 0.5*cache[3],cache[2] + 0.5*cache[3],u)/cache[3];
      } else if (nd < 3) {
         fu = f1->Integral(cache[2] - 0.5*cache[3],cache[2] + 0.5*cache[3],cache[4] - 0.5*cache[5],cache[4] + 0.5*cache[5])/(cache[3]*cache[5]);
      } else {
         fu = f1->Integral(cache[2] - 0.5*cache[3],cache[2] + 0.5*cache[3],cache[4] - 0.5*cache[5],cache[4] + 0.5*cache[5],cache[6] - 0.5*cache[7],cache[6] + 0.5*cache[7])/(cache[3]*cache[5]*cache[7]);
      }
      if (TF1::RejectedPoint()) {cache += fPointSize; continue;}
      eu = cache[1];
      Derivatives(df,x);
      Int_t n = 0;
      fsum = (fu-cu)/eu;
      if (flag!=1) {
         for (j=0;j<npar;j++) {
            if (pl0[j]>0){
               df[n] = df[j]/eu; 
               // left only non-fixed param derivatives / by Sigma
               gin[j] += df[n]*fsum;
               n++;
            }
         }
         Int_t l = 0;
         for (j=0;j<n;j++)
            for (Int_t k=0;k<=j;k++) 
               zik[l++] += df[j]*df[k];
      }
      f += .5*fsum*fsum;
      npfit++;
      cache += fPointSize;
   }
   f1->SetNumberFitPoints(npfit);
   delete[] df;
}


//______________________________________________________________________________
void TFumili::FitLikelihood(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag)
{
   //  Minimization function for H1s using a Likelihood method*-*-*-*-*-*
   //     Basically, it forms the likelihood by determining the Poisson
   //     probability that given a number of entries in a particular bin,
   //     the fit would predict it's value.  This is then done for each bin,
   //     and the sum of the logs is taken as the likelihood.
   //  Default method (function evaluated at center of bin)
   //  for each point the cache contains the following info
   //    -1D : bc,e,xc  (bin content, error, x of center of bin)
   //    -2D : bc,e,xc,yc
   //    -3D : bc,e,xc,yc,zc

   Foption_t fitOption = GetFitOption();
   if (fitOption.Integral) {
      FitLikelihoodI(npar,gin,f,u,flag);
      return;
   }
   Double_t cu,fu,fobs,fsub;
   Double_t dersum[100];
   Double_t x[3];
   Int_t icu;

   TH1 *hfit = (TH1*)GetObjectFit();
   TF1 *f1   = (TF1*)GetUserFunc();
   Int_t nd = hfit->GetDimension();
   Int_t j;
   Double_t *zik = GetZ();
   Double_t *pl0 = GetPL0();

   npar = f1->GetNpar();
   SetParNumber(npar);
   if(flag == 9) return;
   Double_t *df=new Double_t[npar];
   if (flag == 2) for (j=0;j<npar;j++) dersum[j] = gin[j] = 0;
   f1->InitArgs(x,u);
   f = 0;
   
   Int_t npfit = 0;
   Double_t *cache = fCache;
   for (Int_t i=0;i<fNpoints;i++) {
      if (nd > 2) x[2] = cache[4];
      if (nd > 1) x[1] = cache[3];
      x[0]     = cache[2];
      cu  = cache[0];
      TF1::RejectPoint(kFALSE);
      fu = f1->EvalPar(x,u);
      if (TF1::RejectedPoint()) {cache += fPointSize; continue;}
      if (flag == 2) {
         for (j=0;j<npar;j++) {
            dersum[j] += 1; //should be the derivative
            //grad[j] += dersum[j]*(fu-cu)/eu; dersum[j] = 0;
         }
      }
      if (fu < 1.e-9) fu = 1.e-9;
      icu  = Int_t(cu);
      fsub = -fu +icu*TMath::Log(fu);
      fobs = GetSumLog(icu);
      fsub -= fobs;
      Derivatives(df,x);
      int n=0;
      // Here we need gradients of Log likelihood function
      // 
      for (j=0;j<npar;j++) {
         if (pl0[j]>0){
            df[n]   = df[j]*(icu/fu-1); 
            gin[j] -= df[n];
            n++;
         }
      }
      Int_t l = 0;
      // Z-matrix here - production of first derivatives  
      //  of log-likelihood function
      for (j=0;j<n;j++)
         for (Int_t k=0;k<=j;k++) 
            zik[l++] += df[j]*df[k];
            
      f -= fsub;
      npfit++;
      cache += fPointSize;
   }
   f *= 2;
   f1->SetNumberFitPoints(npfit);
   delete[] df;
}


//______________________________________________________________________________
void TFumili::FitLikelihoodI(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag)
{
   //  Minimization function for H1s using a Likelihood method*-*-*-*-*-*
   //     Basically, it forms the likelihood by determining the Poisson
   //     probability that given a number of entries in a particular bin,
   //     the fit would predict it's value.  This is then done for each bin,
   //     and the sum of the logs is taken as the likelihood.
   //  The "I"ntegral method is used
   //  for each point the cache contains the following info
   //    -1D : bc,e,xc,xw  (bin content, error, x of center of bin, x bin width of bin)
   //    -2D : bc,e,xc,xw,yc,yw
   //    -3D : bc,e,xc,xw,yc,yw,zc,zw

   Double_t cu,fu,fobs,fsub;
   Double_t dersum[100];
   Double_t x[3];
   Int_t icu;

   TH1 *hfit = (TH1*)GetObjectFit();
   TF1 *f1   = (TF1*)GetUserFunc();
   Int_t nd = hfit->GetDimension();
   Int_t j;
   Double_t *zik = GetZ();
   Double_t *pl0 = GetPL0();

   Double_t *df=new Double_t[npar];

   npar = f1->GetNpar();
   SetParNumber(npar);
   if(flag == 9) return;
   if (flag == 2) for (j=0;j<npar;j++) dersum[j] = gin[j] = 0;
   f1->InitArgs(x,u);
   f = 0;
   
   Int_t npfit = 0;
   Double_t *cache = fCache;
   for (Int_t i=0;i<fNpoints;i++) {
      if (nd > 2) x[2] = cache[4];
      if (nd > 1) x[1] = cache[3];
      x[0]     = cache[2];
      cu  = cache[0];
      TF1::RejectPoint(kFALSE);
      if (nd < 2) {
         fu = f1->Integral(cache[2] - 0.5*cache[3],cache[2] + 0.5*cache[3],u)/cache[3];
      } else if (nd < 3) {
         fu = f1->Integral(cache[2] - 0.5*cache[3],cache[2] + 0.5*cache[3],cache[4] - 0.5*cache[5],cache[4] + 0.5*cache[5])/(cache[3]*cache[5]);
      } else {
         fu = f1->Integral(cache[2] - 0.5*cache[3],cache[2] + 0.5*cache[3],cache[4] - 0.5*cache[5],cache[4] + 0.5*cache[5],cache[6] - 0.5*cache[7],cache[6] + 0.5*cache[7])/(cache[3]*cache[5]*cache[7]);
      }
      if (TF1::RejectedPoint()) {cache += fPointSize; continue;}
      if (flag == 2) {
         for (j=0;j<npar;j++) {
            dersum[j] += 1; //should be the derivative
            //grad[j] += dersum[j]*(fu-cu)/eu; dersum[j] = 0;
         }
      }
      if (fu < 1.e-9) fu = 1.e-9;
      icu  = Int_t(cu);
      fsub = -fu +icu*TMath::Log(fu);
      fobs = GetSumLog(icu);
      fsub -= fobs;
      Derivatives(df,x);
      int n=0;
      // Here we need gradients of Log likelihood function
      // 
      for (j=0;j<npar;j++) {
         if (pl0[j]>0){
            df[n]   = df[j]*(icu/fu-1); 
            gin[j] -= df[n];
            n++;
         }
      }
      Int_t l = 0;
      // Z-matrix here - production of first derivatives  
      //  of log-likelihood function
      for (j=0;j<n;j++)
         for (Int_t k=0;k<=j;k++) 
            zik[l++] += df[j]*df[k];
            
      f -= fsub;
      npfit++;
      cache += fPointSize;
   }
   f *= 2;
   f1->SetNumberFitPoints(npfit);
   delete[] df;
}


//______________________________________________________________________________
//
//  STATIC functions
//______________________________________________________________________________

//______________________________________________________________________________
void H1FitChisquareFumili(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag)
{
//           Minimization function for H1s using a Chisquare method
//           ======================================================

   TFumili *hFitter = (TFumili*)TVirtualFitter::GetFitter();
   hFitter->FitChisquare(npar, gin, f, u, flag);
}

//______________________________________________________________________________
void H1FitLikelihoodFumili(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag)
{
//   -*-*-*-*Minimization function for H1s using a Likelihood method*-*-*-*-*-*
//           =======================================================
//     Basically, it forms the likelihood by determining the Poisson
//     probability that given a number of entries in a particular bin,
//     the fit would predict it's value.  This is then done for each bin,
//     and the sum of the logs is taken as the likelihood.
//     PDF:  P=exp(-f(x_i))/[F_i]!*(f(x_i))^[F_i]
//    where F_i - experimental value, f(x_i) - expected theoretical value
//    [F_i] - integer part of F_i.
//    drawback is that if F_i>Int_t - GetSumLog will fail
//    for big F_i is faster to use Euler's Gamma-function


   TFumili *hFitter = (TFumili*)TVirtualFitter::GetFitter();
   hFitter->FitLikelihood(npar, gin, f, u, flag);
}



//______________________________________________________________________________
void GraphFitChisquareFumili(Int_t &npar, Double_t * gin, Double_t &f,
                       Double_t *u, Int_t flag)
{
//*-*-*-*-*-*Minimization function for Graphs using a Chisquare method*-*-*-*-*
//*-*        =========================================================
//
// In case of a TGraphErrors object, ex, the error along x,  is projected
// along the y-direction by calculating the function at the points x-exlow and
// x+exhigh.
//
// The chisquare is computed as the sum of the quantity below at each point:
//
//                     (y - f(x))**2
//         -----------------------------------
//         ey**2 + (0.5*(exl + exh)*f'(x))**2
//
// where x and y are the point coordinates and f'(x) is the derivative of function f(x).
// This method to approximate the uncertainty in y because of the errors in x, is called
// "effective variance" method.
// The improvement, compared to the previously used  method (f(x+ exhigh) - f(x-exlow))/2
// is of (error of x)**2 order. 
//  NOTE:
//  1) By using the "effective variance" method a simple linear regression
//      becomes a non-linear case , which takes several iterations
//      instead of 0 as in the linear case .
//
//  2) The effective variance technique assumes that there is no correlation 
//      between the x and y coordinate .
//
// In case the function lies below (above) the data point, ey is ey_low (ey_high).

   Double_t cu,eu,exl,exh,ey,eux,fu,fsum;
   Double_t x[1];
   Int_t i, bin, npfits=0;

   TFumili *grFitter = (TFumili*)TVirtualFitter::GetFitter();
   TGraph *gr     = (TGraph*)grFitter->GetObjectFit();
   TF1 *f1   = (TF1*)grFitter->GetUserFunc();
   Foption_t fitOption = grFitter->GetFitOption();
   
   Int_t n        = gr->GetN();
   Double_t *gx   = gr->GetX();
   Double_t *gy   = gr->GetY();
   npar           = f1->GetNpar();

   grFitter->SetParNumber(npar);

   if(flag == 9) return;
   Double_t *zik = grFitter->GetZ();
   Double_t *pl0 = grFitter->GetPL0();
   Double_t *df  = new Double_t[npar];


   f1->InitArgs(x,u);
   f      = 0;
   for (bin=0;bin<n;bin++) {
      x[0] = gx[bin];
      if (!f1->IsInside(x)) continue;
      cu   = gy[bin];
      TF1::RejectPoint(kFALSE);
      fu   = f1->EvalPar(x,u);
      if (TF1::RejectedPoint()) continue;
      npfits++;
      Double_t eusq=1.;
      if (fitOption.W1) {
         eu = 1.;
      } else {
         exh  = gr->GetErrorXhigh(bin);
         exl  = gr->GetErrorXlow(bin);
         ey  = gr->GetErrorY(bin);
         if (exl < 0) exl = 0;
         if (exh < 0) exh = 0;
         if (ey < 0)  ey  = 0;
         if (exh > 0 && exl > 0) {
//          "Effective variance" method for projecting errors
            eux = 0.5*(exl + exh)*f1->Derivative(x[0], u);
         } else
            eux = 0.;
         eu = ey*ey+eux*eux;
         if (eu <= 0) eu = 1;
         eusq = TMath::Sqrt(eu);
      }
      grFitter->Derivatives(df,x);
      n = 0;
      fsum = (fu-cu)/eusq;
      for (i=0;i<npar;i++) {
         if (pl0[i]>0){
            df[n] = df[i]/eusq; 
            // left only non-fixed param derivatives / by Sigma
            gin[i] += df[n]*fsum;
            n++;
         }
      }
      Int_t l = 0;
      for (i=0;i<n;i++)
         for (Int_t j=0;j<=i;j++) 
            zik[l++] += df[i]*df[j];
      f += .5*fsum*fsum;

   }
   delete [] df;
   f1->SetNumberFitPoints(npfits);
}


 TFumili.cxx:1
 TFumili.cxx:2
 TFumili.cxx:3
 TFumili.cxx:4
 TFumili.cxx:5
 TFumili.cxx:6
 TFumili.cxx:7
 TFumili.cxx:8
 TFumili.cxx:9
 TFumili.cxx:10
 TFumili.cxx:11
 TFumili.cxx:12
 TFumili.cxx:13
 TFumili.cxx:14
 TFumili.cxx:15
 TFumili.cxx:16
 TFumili.cxx:17
 TFumili.cxx:18
 TFumili.cxx:19
 TFumili.cxx:20
 TFumili.cxx:21
 TFumili.cxx:22
 TFumili.cxx:23
 TFumili.cxx:24
 TFumili.cxx:25
 TFumili.cxx:26
 TFumili.cxx:27
 TFumili.cxx:28
 TFumili.cxx:29
 TFumili.cxx:30
 TFumili.cxx:31
 TFumili.cxx:32
 TFumili.cxx:33
 TFumili.cxx:34
 TFumili.cxx:35
 TFumili.cxx:36
 TFumili.cxx:37
 TFumili.cxx:38
 TFumili.cxx:39
 TFumili.cxx:40
 TFumili.cxx:41
 TFumili.cxx:42
 TFumili.cxx:43
 TFumili.cxx:44
 TFumili.cxx:45
 TFumili.cxx:46
 TFumili.cxx:47
 TFumili.cxx:48
 TFumili.cxx:49
 TFumili.cxx:50
 TFumili.cxx:51
 TFumili.cxx:52
 TFumili.cxx:53
 TFumili.cxx:54
 TFumili.cxx:55
 TFumili.cxx:56
 TFumili.cxx:57
 TFumili.cxx:58
 TFumili.cxx:59
 TFumili.cxx:60
 TFumili.cxx:61
 TFumili.cxx:62
 TFumili.cxx:63
 TFumili.cxx:64
 TFumili.cxx:65
 TFumili.cxx:66
 TFumili.cxx:67
 TFumili.cxx:68
 TFumili.cxx:69
 TFumili.cxx:70
 TFumili.cxx:71
 TFumili.cxx:72
 TFumili.cxx:73
 TFumili.cxx:74
 TFumili.cxx:75
 TFumili.cxx:76
 TFumili.cxx:77
 TFumili.cxx:78
 TFumili.cxx:79
 TFumili.cxx:80
 TFumili.cxx:81
 TFumili.cxx:82
 TFumili.cxx:83
 TFumili.cxx:84
 TFumili.cxx:85
 TFumili.cxx:86
 TFumili.cxx:87
 TFumili.cxx:88
 TFumili.cxx:89
 TFumili.cxx:90
 TFumili.cxx:91
 TFumili.cxx:92
 TFumili.cxx:93
 TFumili.cxx:94
 TFumili.cxx:95
 TFumili.cxx:96
 TFumili.cxx:97
 TFumili.cxx:98
 TFumili.cxx:99
 TFumili.cxx:100
 TFumili.cxx:101
 TFumili.cxx:102
 TFumili.cxx:103
 TFumili.cxx:104
 TFumili.cxx:105
 TFumili.cxx:106
 TFumili.cxx:107
 TFumili.cxx:108
 TFumili.cxx:109
 TFumili.cxx:110
 TFumili.cxx:111
 TFumili.cxx:112
 TFumili.cxx:113
 TFumili.cxx:114
 TFumili.cxx:115
 TFumili.cxx:116
 TFumili.cxx:117
 TFumili.cxx:118
 TFumili.cxx:119
 TFumili.cxx:120
 TFumili.cxx:121
 TFumili.cxx:122
 TFumili.cxx:123
 TFumili.cxx:124
 TFumili.cxx:125
 TFumili.cxx:126
 TFumili.cxx:127
 TFumili.cxx:128
 TFumili.cxx:129
 TFumili.cxx:130
 TFumili.cxx:131
 TFumili.cxx:132
 TFumili.cxx:133
 TFumili.cxx:134
 TFumili.cxx:135
 TFumili.cxx:136
 TFumili.cxx:137
 TFumili.cxx:138
 TFumili.cxx:139
 TFumili.cxx:140
 TFumili.cxx:141
 TFumili.cxx:142
 TFumili.cxx:143
 TFumili.cxx:144
 TFumili.cxx:145
 TFumili.cxx:146
 TFumili.cxx:147
 TFumili.cxx:148
 TFumili.cxx:149
 TFumili.cxx:150
 TFumili.cxx:151
 TFumili.cxx:152
 TFumili.cxx:153
 TFumili.cxx:154
 TFumili.cxx:155
 TFumili.cxx:156
 TFumili.cxx:157
 TFumili.cxx:158
 TFumili.cxx:159
 TFumili.cxx:160
 TFumili.cxx:161
 TFumili.cxx:162
 TFumili.cxx:163
 TFumili.cxx:164
 TFumili.cxx:165
 TFumili.cxx:166
 TFumili.cxx:167
 TFumili.cxx:168
 TFumili.cxx:169
 TFumili.cxx:170
 TFumili.cxx:171
 TFumili.cxx:172
 TFumili.cxx:173
 TFumili.cxx:174
 TFumili.cxx:175
 TFumili.cxx:176
 TFumili.cxx:177
 TFumili.cxx:178
 TFumili.cxx:179
 TFumili.cxx:180
 TFumili.cxx:181
 TFumili.cxx:182
 TFumili.cxx:183
 TFumili.cxx:184
 TFumili.cxx:185
 TFumili.cxx:186
 TFumili.cxx:187
 TFumili.cxx:188
 TFumili.cxx:189
 TFumili.cxx:190
 TFumili.cxx:191
 TFumili.cxx:192
 TFumili.cxx:193
 TFumili.cxx:194
 TFumili.cxx:195
 TFumili.cxx:196
 TFumili.cxx:197
 TFumili.cxx:198
 TFumili.cxx:199
 TFumili.cxx:200
 TFumili.cxx:201
 TFumili.cxx:202
 TFumili.cxx:203
 TFumili.cxx:204
 TFumili.cxx:205
 TFumili.cxx:206
 TFumili.cxx:207
 TFumili.cxx:208
 TFumili.cxx:209
 TFumili.cxx:210
 TFumili.cxx:211
 TFumili.cxx:212
 TFumili.cxx:213
 TFumili.cxx:214
 TFumili.cxx:215
 TFumili.cxx:216
 TFumili.cxx:217
 TFumili.cxx:218
 TFumili.cxx:219
 TFumili.cxx:220
 TFumili.cxx:221
 TFumili.cxx:222
 TFumili.cxx:223
 TFumili.cxx:224
 TFumili.cxx:225
 TFumili.cxx:226
 TFumili.cxx:227
 TFumili.cxx:228
 TFumili.cxx:229
 TFumili.cxx:230
 TFumili.cxx:231
 TFumili.cxx:232
 TFumili.cxx:233
 TFumili.cxx:234
 TFumili.cxx:235
 TFumili.cxx:236
 TFumili.cxx:237
 TFumili.cxx:238
 TFumili.cxx:239
 TFumili.cxx:240
 TFumili.cxx:241
 TFumili.cxx:242
 TFumili.cxx:243
 TFumili.cxx:244
 TFumili.cxx:245
 TFumili.cxx:246
 TFumili.cxx:247
 TFumili.cxx:248
 TFumili.cxx:249
 TFumili.cxx:250
 TFumili.cxx:251
 TFumili.cxx:252
 TFumili.cxx:253
 TFumili.cxx:254
 TFumili.cxx:255
 TFumili.cxx:256
 TFumili.cxx:257
 TFumili.cxx:258
 TFumili.cxx:259
 TFumili.cxx:260
 TFumili.cxx:261
 TFumili.cxx:262
 TFumili.cxx:263
 TFumili.cxx:264
 TFumili.cxx:265
 TFumili.cxx:266
 TFumili.cxx:267
 TFumili.cxx:268
 TFumili.cxx:269
 TFumili.cxx:270
 TFumili.cxx:271
 TFumili.cxx:272
 TFumili.cxx:273
 TFumili.cxx:274
 TFumili.cxx:275
 TFumili.cxx:276
 TFumili.cxx:277
 TFumili.cxx:278
 TFumili.cxx:279
 TFumili.cxx:280
 TFumili.cxx:281
 TFumili.cxx:282
 TFumili.cxx:283
 TFumili.cxx:284
 TFumili.cxx:285
 TFumili.cxx:286
 TFumili.cxx:287
 TFumili.cxx:288
 TFumili.cxx:289
 TFumili.cxx:290
 TFumili.cxx:291
 TFumili.cxx:292
 TFumili.cxx:293
 TFumili.cxx:294
 TFumili.cxx:295
 TFumili.cxx:296
 TFumili.cxx:297
 TFumili.cxx:298
 TFumili.cxx:299
 TFumili.cxx:300
 TFumili.cxx:301
 TFumili.cxx:302
 TFumili.cxx:303
 TFumili.cxx:304
 TFumili.cxx:305
 TFumili.cxx:306
 TFumili.cxx:307
 TFumili.cxx:308
 TFumili.cxx:309
 TFumili.cxx:310
 TFumili.cxx:311
 TFumili.cxx:312
 TFumili.cxx:313
 TFumili.cxx:314
 TFumili.cxx:315
 TFumili.cxx:316
 TFumili.cxx:317
 TFumili.cxx:318
 TFumili.cxx:319
 TFumili.cxx:320
 TFumili.cxx:321
 TFumili.cxx:322
 TFumili.cxx:323
 TFumili.cxx:324
 TFumili.cxx:325
 TFumili.cxx:326
 TFumili.cxx:327
 TFumili.cxx:328
 TFumili.cxx:329
 TFumili.cxx:330
 TFumili.cxx:331
 TFumili.cxx:332
 TFumili.cxx:333
 TFumili.cxx:334
 TFumili.cxx:335
 TFumili.cxx:336
 TFumili.cxx:337
 TFumili.cxx:338
 TFumili.cxx:339
 TFumili.cxx:340
 TFumili.cxx:341
 TFumili.cxx:342
 TFumili.cxx:343
 TFumili.cxx:344
 TFumili.cxx:345
 TFumili.cxx:346
 TFumili.cxx:347
 TFumili.cxx:348
 TFumili.cxx:349
 TFumili.cxx:350
 TFumili.cxx:351
 TFumili.cxx:352
 TFumili.cxx:353
 TFumili.cxx:354
 TFumili.cxx:355
 TFumili.cxx:356
 TFumili.cxx:357
 TFumili.cxx:358
 TFumili.cxx:359
 TFumili.cxx:360
 TFumili.cxx:361
 TFumili.cxx:362
 TFumili.cxx:363
 TFumili.cxx:364
 TFumili.cxx:365
 TFumili.cxx:366
 TFumili.cxx:367
 TFumili.cxx:368
 TFumili.cxx:369
 TFumili.cxx:370
 TFumili.cxx:371
 TFumili.cxx:372
 TFumili.cxx:373
 TFumili.cxx:374
 TFumili.cxx:375
 TFumili.cxx:376
 TFumili.cxx:377
 TFumili.cxx:378
 TFumili.cxx:379
 TFumili.cxx:380
 TFumili.cxx:381
 TFumili.cxx:382
 TFumili.cxx:383
 TFumili.cxx:384
 TFumili.cxx:385
 TFumili.cxx:386
 TFumili.cxx:387
 TFumili.cxx:388
 TFumili.cxx:389
 TFumili.cxx:390
 TFumili.cxx:391
 TFumili.cxx:392
 TFumili.cxx:393
 TFumili.cxx:394
 TFumili.cxx:395
 TFumili.cxx:396
 TFumili.cxx:397
 TFumili.cxx:398
 TFumili.cxx:399
 TFumili.cxx:400
 TFumili.cxx:401
 TFumili.cxx:402
 TFumili.cxx:403
 TFumili.cxx:404
 TFumili.cxx:405
 TFumili.cxx:406
 TFumili.cxx:407
 TFumili.cxx:408
 TFumili.cxx:409
 TFumili.cxx:410
 TFumili.cxx:411
 TFumili.cxx:412
 TFumili.cxx:413
 TFumili.cxx:414
 TFumili.cxx:415
 TFumili.cxx:416
 TFumili.cxx:417
 TFumili.cxx:418
 TFumili.cxx:419
 TFumili.cxx:420
 TFumili.cxx:421
 TFumili.cxx:422
 TFumili.cxx:423
 TFumili.cxx:424
 TFumili.cxx:425
 TFumili.cxx:426
 TFumili.cxx:427
 TFumili.cxx:428
 TFumili.cxx:429
 TFumili.cxx:430
 TFumili.cxx:431
 TFumili.cxx:432
 TFumili.cxx:433
 TFumili.cxx:434
 TFumili.cxx:435
 TFumili.cxx:436
 TFumili.cxx:437
 TFumili.cxx:438
 TFumili.cxx:439
 TFumili.cxx:440
 TFumili.cxx:441
 TFumili.cxx:442
 TFumili.cxx:443
 TFumili.cxx:444
 TFumili.cxx:445
 TFumili.cxx:446
 TFumili.cxx:447
 TFumili.cxx:448
 TFumili.cxx:449
 TFumili.cxx:450
 TFumili.cxx:451
 TFumili.cxx:452
 TFumili.cxx:453
 TFumili.cxx:454
 TFumili.cxx:455
 TFumili.cxx:456
 TFumili.cxx:457
 TFumili.cxx:458
 TFumili.cxx:459
 TFumili.cxx:460
 TFumili.cxx:461
 TFumili.cxx:462
 TFumili.cxx:463
 TFumili.cxx:464
 TFumili.cxx:465
 TFumili.cxx:466
 TFumili.cxx:467
 TFumili.cxx:468
 TFumili.cxx:469
 TFumili.cxx:470
 TFumili.cxx:471
 TFumili.cxx:472
 TFumili.cxx:473
 TFumili.cxx:474
 TFumili.cxx:475
 TFumili.cxx:476
 TFumili.cxx:477
 TFumili.cxx:478
 TFumili.cxx:479
 TFumili.cxx:480
 TFumili.cxx:481
 TFumili.cxx:482
 TFumili.cxx:483
 TFumili.cxx:484
 TFumili.cxx:485
 TFumili.cxx:486
 TFumili.cxx:487
 TFumili.cxx:488
 TFumili.cxx:489
 TFumili.cxx:490
 TFumili.cxx:491
 TFumili.cxx:492
 TFumili.cxx:493
 TFumili.cxx:494
 TFumili.cxx:495
 TFumili.cxx:496
 TFumili.cxx:497
 TFumili.cxx:498
 TFumili.cxx:499
 TFumili.cxx:500
 TFumili.cxx:501
 TFumili.cxx:502
 TFumili.cxx:503
 TFumili.cxx:504
 TFumili.cxx:505
 TFumili.cxx:506
 TFumili.cxx:507
 TFumili.cxx:508
 TFumili.cxx:509
 TFumili.cxx:510
 TFumili.cxx:511
 TFumili.cxx:512
 TFumili.cxx:513
 TFumili.cxx:514
 TFumili.cxx:515
 TFumili.cxx:516
 TFumili.cxx:517
 TFumili.cxx:518
 TFumili.cxx:519
 TFumili.cxx:520
 TFumili.cxx:521
 TFumili.cxx:522
 TFumili.cxx:523
 TFumili.cxx:524
 TFumili.cxx:525
 TFumili.cxx:526
 TFumili.cxx:527
 TFumili.cxx:528
 TFumili.cxx:529
 TFumili.cxx:530
 TFumili.cxx:531
 TFumili.cxx:532
 TFumili.cxx:533
 TFumili.cxx:534
 TFumili.cxx:535
 TFumili.cxx:536
 TFumili.cxx:537
 TFumili.cxx:538
 TFumili.cxx:539
 TFumili.cxx:540
 TFumili.cxx:541
 TFumili.cxx:542
 TFumili.cxx:543
 TFumili.cxx:544
 TFumili.cxx:545
 TFumili.cxx:546
 TFumili.cxx:547
 TFumili.cxx:548
 TFumili.cxx:549
 TFumili.cxx:550
 TFumili.cxx:551
 TFumili.cxx:552
 TFumili.cxx:553
 TFumili.cxx:554
 TFumili.cxx:555
 TFumili.cxx:556
 TFumili.cxx:557
 TFumili.cxx:558
 TFumili.cxx:559
 TFumili.cxx:560
 TFumili.cxx:561
 TFumili.cxx:562
 TFumili.cxx:563
 TFumili.cxx:564
 TFumili.cxx:565
 TFumili.cxx:566
 TFumili.cxx:567
 TFumili.cxx:568
 TFumili.cxx:569
 TFumili.cxx:570
 TFumili.cxx:571
 TFumili.cxx:572
 TFumili.cxx:573
 TFumili.cxx:574
 TFumili.cxx:575
 TFumili.cxx:576
 TFumili.cxx:577
 TFumili.cxx:578
 TFumili.cxx:579
 TFumili.cxx:580
 TFumili.cxx:581
 TFumili.cxx:582
 TFumili.cxx:583
 TFumili.cxx:584
 TFumili.cxx:585
 TFumili.cxx:586
 TFumili.cxx:587
 TFumili.cxx:588
 TFumili.cxx:589
 TFumili.cxx:590
 TFumili.cxx:591
 TFumili.cxx:592
 TFumili.cxx:593
 TFumili.cxx:594
 TFumili.cxx:595
 TFumili.cxx:596
 TFumili.cxx:597
 TFumili.cxx:598
 TFumili.cxx:599
 TFumili.cxx:600
 TFumili.cxx:601
 TFumili.cxx:602
 TFumili.cxx:603
 TFumili.cxx:604
 TFumili.cxx:605
 TFumili.cxx:606
 TFumili.cxx:607
 TFumili.cxx:608
 TFumili.cxx:609
 TFumili.cxx:610
 TFumili.cxx:611
 TFumili.cxx:612
 TFumili.cxx:613
 TFumili.cxx:614
 TFumili.cxx:615
 TFumili.cxx:616
 TFumili.cxx:617
 TFumili.cxx:618
 TFumili.cxx:619
 TFumili.cxx:620
 TFumili.cxx:621
 TFumili.cxx:622
 TFumili.cxx:623
 TFumili.cxx:624
 TFumili.cxx:625
 TFumili.cxx:626
 TFumili.cxx:627
 TFumili.cxx:628
 TFumili.cxx:629
 TFumili.cxx:630
 TFumili.cxx:631
 TFumili.cxx:632
 TFumili.cxx:633
 TFumili.cxx:634
 TFumili.cxx:635
 TFumili.cxx:636
 TFumili.cxx:637
 TFumili.cxx:638
 TFumili.cxx:639
 TFumili.cxx:640
 TFumili.cxx:641
 TFumili.cxx:642
 TFumili.cxx:643
 TFumili.cxx:644
 TFumili.cxx:645
 TFumili.cxx:646
 TFumili.cxx:647
 TFumili.cxx:648
 TFumili.cxx:649
 TFumili.cxx:650
 TFumili.cxx:651
 TFumili.cxx:652
 TFumili.cxx:653
 TFumili.cxx:654
 TFumili.cxx:655
 TFumili.cxx:656
 TFumili.cxx:657
 TFumili.cxx:658
 TFumili.cxx:659
 TFumili.cxx:660
 TFumili.cxx:661
 TFumili.cxx:662
 TFumili.cxx:663
 TFumili.cxx:664
 TFumili.cxx:665
 TFumili.cxx:666
 TFumili.cxx:667
 TFumili.cxx:668
 TFumili.cxx:669
 TFumili.cxx:670
 TFumili.cxx:671
 TFumili.cxx:672
 TFumili.cxx:673
 TFumili.cxx:674
 TFumili.cxx:675
 TFumili.cxx:676
 TFumili.cxx:677
 TFumili.cxx:678
 TFumili.cxx:679
 TFumili.cxx:680
 TFumili.cxx:681
 TFumili.cxx:682
 TFumili.cxx:683
 TFumili.cxx:684
 TFumili.cxx:685
 TFumili.cxx:686
 TFumili.cxx:687
 TFumili.cxx:688
 TFumili.cxx:689
 TFumili.cxx:690
 TFumili.cxx:691
 TFumili.cxx:692
 TFumili.cxx:693
 TFumili.cxx:694
 TFumili.cxx:695
 TFumili.cxx:696
 TFumili.cxx:697
 TFumili.cxx:698
 TFumili.cxx:699
 TFumili.cxx:700
 TFumili.cxx:701
 TFumili.cxx:702
 TFumili.cxx:703
 TFumili.cxx:704
 TFumili.cxx:705
 TFumili.cxx:706
 TFumili.cxx:707
 TFumili.cxx:708
 TFumili.cxx:709
 TFumili.cxx:710
 TFumili.cxx:711
 TFumili.cxx:712
 TFumili.cxx:713
 TFumili.cxx:714
 TFumili.cxx:715
 TFumili.cxx:716
 TFumili.cxx:717
 TFumili.cxx:718
 TFumili.cxx:719
 TFumili.cxx:720
 TFumili.cxx:721
 TFumili.cxx:722
 TFumili.cxx:723
 TFumili.cxx:724
 TFumili.cxx:725
 TFumili.cxx:726
 TFumili.cxx:727
 TFumili.cxx:728
 TFumili.cxx:729
 TFumili.cxx:730
 TFumili.cxx:731
 TFumili.cxx:732
 TFumili.cxx:733
 TFumili.cxx:734
 TFumili.cxx:735
 TFumili.cxx:736
 TFumili.cxx:737
 TFumili.cxx:738
 TFumili.cxx:739
 TFumili.cxx:740
 TFumili.cxx:741
 TFumili.cxx:742
 TFumili.cxx:743
 TFumili.cxx:744
 TFumili.cxx:745
 TFumili.cxx:746
 TFumili.cxx:747
 TFumili.cxx:748
 TFumili.cxx:749
 TFumili.cxx:750
 TFumili.cxx:751
 TFumili.cxx:752
 TFumili.cxx:753
 TFumili.cxx:754
 TFumili.cxx:755
 TFumili.cxx:756
 TFumili.cxx:757
 TFumili.cxx:758
 TFumili.cxx:759
 TFumili.cxx:760
 TFumili.cxx:761
 TFumili.cxx:762
 TFumili.cxx:763
 TFumili.cxx:764
 TFumili.cxx:765
 TFumili.cxx:766
 TFumili.cxx:767
 TFumili.cxx:768
 TFumili.cxx:769
 TFumili.cxx:770
 TFumili.cxx:771
 TFumili.cxx:772
 TFumili.cxx:773
 TFumili.cxx:774
 TFumili.cxx:775
 TFumili.cxx:776
 TFumili.cxx:777
 TFumili.cxx:778
 TFumili.cxx:779
 TFumili.cxx:780
 TFumili.cxx:781
 TFumili.cxx:782
 TFumili.cxx:783
 TFumili.cxx:784
 TFumili.cxx:785
 TFumili.cxx:786
 TFumili.cxx:787
 TFumili.cxx:788
 TFumili.cxx:789
 TFumili.cxx:790
 TFumili.cxx:791
 TFumili.cxx:792
 TFumili.cxx:793
 TFumili.cxx:794
 TFumili.cxx:795
 TFumili.cxx:796
 TFumili.cxx:797
 TFumili.cxx:798
 TFumili.cxx:799
 TFumili.cxx:800
 TFumili.cxx:801
 TFumili.cxx:802
 TFumili.cxx:803
 TFumili.cxx:804
 TFumili.cxx:805
 TFumili.cxx:806
 TFumili.cxx:807
 TFumili.cxx:808
 TFumili.cxx:809
 TFumili.cxx:810
 TFumili.cxx:811
 TFumili.cxx:812
 TFumili.cxx:813
 TFumili.cxx:814
 TFumili.cxx:815
 TFumili.cxx:816
 TFumili.cxx:817
 TFumili.cxx:818
 TFumili.cxx:819
 TFumili.cxx:820
 TFumili.cxx:821
 TFumili.cxx:822
 TFumili.cxx:823
 TFumili.cxx:824
 TFumili.cxx:825
 TFumili.cxx:826
 TFumili.cxx:827
 TFumili.cxx:828
 TFumili.cxx:829
 TFumili.cxx:830
 TFumili.cxx:831
 TFumili.cxx:832
 TFumili.cxx:833
 TFumili.cxx:834
 TFumili.cxx:835
 TFumili.cxx:836
 TFumili.cxx:837
 TFumili.cxx:838
 TFumili.cxx:839
 TFumili.cxx:840
 TFumili.cxx:841
 TFumili.cxx:842
 TFumili.cxx:843
 TFumili.cxx:844
 TFumili.cxx:845
 TFumili.cxx:846
 TFumili.cxx:847
 TFumili.cxx:848
 TFumili.cxx:849
 TFumili.cxx:850
 TFumili.cxx:851
 TFumili.cxx:852
 TFumili.cxx:853
 TFumili.cxx:854
 TFumili.cxx:855
 TFumili.cxx:856
 TFumili.cxx:857
 TFumili.cxx:858
 TFumili.cxx:859
 TFumili.cxx:860
 TFumili.cxx:861
 TFumili.cxx:862
 TFumili.cxx:863
 TFumili.cxx:864
 TFumili.cxx:865
 TFumili.cxx:866
 TFumili.cxx:867
 TFumili.cxx:868
 TFumili.cxx:869
 TFumili.cxx:870
 TFumili.cxx:871
 TFumili.cxx:872
 TFumili.cxx:873
 TFumili.cxx:874
 TFumili.cxx:875
 TFumili.cxx:876
 TFumili.cxx:877
 TFumili.cxx:878
 TFumili.cxx:879
 TFumili.cxx:880
 TFumili.cxx:881
 TFumili.cxx:882
 TFumili.cxx:883
 TFumili.cxx:884
 TFumili.cxx:885
 TFumili.cxx:886
 TFumili.cxx:887
 TFumili.cxx:888
 TFumili.cxx:889
 TFumili.cxx:890
 TFumili.cxx:891
 TFumili.cxx:892
 TFumili.cxx:893
 TFumili.cxx:894
 TFumili.cxx:895
 TFumili.cxx:896
 TFumili.cxx:897
 TFumili.cxx:898
 TFumili.cxx:899
 TFumili.cxx:900
 TFumili.cxx:901
 TFumili.cxx:902
 TFumili.cxx:903
 TFumili.cxx:904
 TFumili.cxx:905
 TFumili.cxx:906
 TFumili.cxx:907
 TFumili.cxx:908
 TFumili.cxx:909
 TFumili.cxx:910
 TFumili.cxx:911
 TFumili.cxx:912
 TFumili.cxx:913
 TFumili.cxx:914
 TFumili.cxx:915
 TFumili.cxx:916
 TFumili.cxx:917
 TFumili.cxx:918
 TFumili.cxx:919
 TFumili.cxx:920
 TFumili.cxx:921
 TFumili.cxx:922
 TFumili.cxx:923
 TFumili.cxx:924
 TFumili.cxx:925
 TFumili.cxx:926
 TFumili.cxx:927
 TFumili.cxx:928
 TFumili.cxx:929
 TFumili.cxx:930
 TFumili.cxx:931
 TFumili.cxx:932
 TFumili.cxx:933
 TFumili.cxx:934
 TFumili.cxx:935
 TFumili.cxx:936
 TFumili.cxx:937
 TFumili.cxx:938
 TFumili.cxx:939
 TFumili.cxx:940
 TFumili.cxx:941
 TFumili.cxx:942
 TFumili.cxx:943
 TFumili.cxx:944
 TFumili.cxx:945
 TFumili.cxx:946
 TFumili.cxx:947
 TFumili.cxx:948
 TFumili.cxx:949
 TFumili.cxx:950
 TFumili.cxx:951
 TFumili.cxx:952
 TFumili.cxx:953
 TFumili.cxx:954
 TFumili.cxx:955
 TFumili.cxx:956
 TFumili.cxx:957
 TFumili.cxx:958
 TFumili.cxx:959
 TFumili.cxx:960
 TFumili.cxx:961
 TFumili.cxx:962
 TFumili.cxx:963
 TFumili.cxx:964
 TFumili.cxx:965
 TFumili.cxx:966
 TFumili.cxx:967
 TFumili.cxx:968
 TFumili.cxx:969
 TFumili.cxx:970
 TFumili.cxx:971
 TFumili.cxx:972
 TFumili.cxx:973
 TFumili.cxx:974
 TFumili.cxx:975
 TFumili.cxx:976
 TFumili.cxx:977
 TFumili.cxx:978
 TFumili.cxx:979
 TFumili.cxx:980
 TFumili.cxx:981
 TFumili.cxx:982
 TFumili.cxx:983
 TFumili.cxx:984
 TFumili.cxx:985
 TFumili.cxx:986
 TFumili.cxx:987
 TFumili.cxx:988
 TFumili.cxx:989
 TFumili.cxx:990
 TFumili.cxx:991
 TFumili.cxx:992
 TFumili.cxx:993
 TFumili.cxx:994
 TFumili.cxx:995
 TFumili.cxx:996
 TFumili.cxx:997
 TFumili.cxx:998
 TFumili.cxx:999
 TFumili.cxx:1000
 TFumili.cxx:1001
 TFumili.cxx:1002
 TFumili.cxx:1003
 TFumili.cxx:1004
 TFumili.cxx:1005
 TFumili.cxx:1006
 TFumili.cxx:1007
 TFumili.cxx:1008
 TFumili.cxx:1009
 TFumili.cxx:1010
 TFumili.cxx:1011
 TFumili.cxx:1012
 TFumili.cxx:1013
 TFumili.cxx:1014
 TFumili.cxx:1015
 TFumili.cxx:1016
 TFumili.cxx:1017
 TFumili.cxx:1018
 TFumili.cxx:1019
 TFumili.cxx:1020
 TFumili.cxx:1021
 TFumili.cxx:1022
 TFumili.cxx:1023
 TFumili.cxx:1024
 TFumili.cxx:1025
 TFumili.cxx:1026
 TFumili.cxx:1027
 TFumili.cxx:1028
 TFumili.cxx:1029
 TFumili.cxx:1030
 TFumili.cxx:1031
 TFumili.cxx:1032
 TFumili.cxx:1033
 TFumili.cxx:1034
 TFumili.cxx:1035
 TFumili.cxx:1036
 TFumili.cxx:1037
 TFumili.cxx:1038
 TFumili.cxx:1039
 TFumili.cxx:1040
 TFumili.cxx:1041
 TFumili.cxx:1042
 TFumili.cxx:1043
 TFumili.cxx:1044
 TFumili.cxx:1045
 TFumili.cxx:1046
 TFumili.cxx:1047
 TFumili.cxx:1048
 TFumili.cxx:1049
 TFumili.cxx:1050
 TFumili.cxx:1051
 TFumili.cxx:1052
 TFumili.cxx:1053
 TFumili.cxx:1054
 TFumili.cxx:1055
 TFumili.cxx:1056
 TFumili.cxx:1057
 TFumili.cxx:1058
 TFumili.cxx:1059
 TFumili.cxx:1060
 TFumili.cxx:1061
 TFumili.cxx:1062
 TFumili.cxx:1063
 TFumili.cxx:1064
 TFumili.cxx:1065
 TFumili.cxx:1066
 TFumili.cxx:1067
 TFumili.cxx:1068
 TFumili.cxx:1069
 TFumili.cxx:1070
 TFumili.cxx:1071
 TFumili.cxx:1072
 TFumili.cxx:1073
 TFumili.cxx:1074
 TFumili.cxx:1075
 TFumili.cxx:1076
 TFumili.cxx:1077
 TFumili.cxx:1078
 TFumili.cxx:1079
 TFumili.cxx:1080
 TFumili.cxx:1081
 TFumili.cxx:1082
 TFumili.cxx:1083
 TFumili.cxx:1084
 TFumili.cxx:1085
 TFumili.cxx:1086
 TFumili.cxx:1087
 TFumili.cxx:1088
 TFumili.cxx:1089
 TFumili.cxx:1090
 TFumili.cxx:1091
 TFumili.cxx:1092
 TFumili.cxx:1093
 TFumili.cxx:1094
 TFumili.cxx:1095
 TFumili.cxx:1096
 TFumili.cxx:1097
 TFumili.cxx:1098
 TFumili.cxx:1099
 TFumili.cxx:1100
 TFumili.cxx:1101
 TFumili.cxx:1102
 TFumili.cxx:1103
 TFumili.cxx:1104
 TFumili.cxx:1105
 TFumili.cxx:1106
 TFumili.cxx:1107
 TFumili.cxx:1108
 TFumili.cxx:1109
 TFumili.cxx:1110
 TFumili.cxx:1111
 TFumili.cxx:1112
 TFumili.cxx:1113
 TFumili.cxx:1114
 TFumili.cxx:1115
 TFumili.cxx:1116
 TFumili.cxx:1117
 TFumili.cxx:1118
 TFumili.cxx:1119
 TFumili.cxx:1120
 TFumili.cxx:1121
 TFumili.cxx:1122
 TFumili.cxx:1123
 TFumili.cxx:1124
 TFumili.cxx:1125
 TFumili.cxx:1126
 TFumili.cxx:1127
 TFumili.cxx:1128
 TFumili.cxx:1129
 TFumili.cxx:1130
 TFumili.cxx:1131
 TFumili.cxx:1132
 TFumili.cxx:1133
 TFumili.cxx:1134
 TFumili.cxx:1135
 TFumili.cxx:1136
 TFumili.cxx:1137
 TFumili.cxx:1138
 TFumili.cxx:1139
 TFumili.cxx:1140
 TFumili.cxx:1141
 TFumili.cxx:1142
 TFumili.cxx:1143
 TFumili.cxx:1144
 TFumili.cxx:1145
 TFumili.cxx:1146
 TFumili.cxx:1147
 TFumili.cxx:1148
 TFumili.cxx:1149
 TFumili.cxx:1150
 TFumili.cxx:1151
 TFumili.cxx:1152
 TFumili.cxx:1153
 TFumili.cxx:1154
 TFumili.cxx:1155
 TFumili.cxx:1156
 TFumili.cxx:1157
 TFumili.cxx:1158
 TFumili.cxx:1159
 TFumili.cxx:1160
 TFumili.cxx:1161
 TFumili.cxx:1162
 TFumili.cxx:1163
 TFumili.cxx:1164
 TFumili.cxx:1165
 TFumili.cxx:1166
 TFumili.cxx:1167
 TFumili.cxx:1168
 TFumili.cxx:1169
 TFumili.cxx:1170
 TFumili.cxx:1171
 TFumili.cxx:1172
 TFumili.cxx:1173
 TFumili.cxx:1174
 TFumili.cxx:1175
 TFumili.cxx:1176
 TFumili.cxx:1177
 TFumili.cxx:1178
 TFumili.cxx:1179
 TFumili.cxx:1180
 TFumili.cxx:1181
 TFumili.cxx:1182
 TFumili.cxx:1183
 TFumili.cxx:1184
 TFumili.cxx:1185
 TFumili.cxx:1186
 TFumili.cxx:1187
 TFumili.cxx:1188
 TFumili.cxx:1189
 TFumili.cxx:1190
 TFumili.cxx:1191
 TFumili.cxx:1192
 TFumili.cxx:1193
 TFumili.cxx:1194
 TFumili.cxx:1195
 TFumili.cxx:1196
 TFumili.cxx:1197
 TFumili.cxx:1198
 TFumili.cxx:1199
 TFumili.cxx:1200
 TFumili.cxx:1201
 TFumili.cxx:1202
 TFumili.cxx:1203
 TFumili.cxx:1204
 TFumili.cxx:1205
 TFumili.cxx:1206
 TFumili.cxx:1207
 TFumili.cxx:1208
 TFumili.cxx:1209
 TFumili.cxx:1210
 TFumili.cxx:1211
 TFumili.cxx:1212
 TFumili.cxx:1213
 TFumili.cxx:1214
 TFumili.cxx:1215
 TFumili.cxx:1216
 TFumili.cxx:1217
 TFumili.cxx:1218
 TFumili.cxx:1219
 TFumili.cxx:1220
 TFumili.cxx:1221
 TFumili.cxx:1222
 TFumili.cxx:1223
 TFumili.cxx:1224
 TFumili.cxx:1225
 TFumili.cxx:1226
 TFumili.cxx:1227
 TFumili.cxx:1228
 TFumili.cxx:1229
 TFumili.cxx:1230
 TFumili.cxx:1231
 TFumili.cxx:1232
 TFumili.cxx:1233
 TFumili.cxx:1234
 TFumili.cxx:1235
 TFumili.cxx:1236
 TFumili.cxx:1237
 TFumili.cxx:1238
 TFumili.cxx:1239
 TFumili.cxx:1240
 TFumili.cxx:1241
 TFumili.cxx:1242
 TFumili.cxx:1243
 TFumili.cxx:1244
 TFumili.cxx:1245
 TFumili.cxx:1246
 TFumili.cxx:1247
 TFumili.cxx:1248
 TFumili.cxx:1249
 TFumili.cxx:1250
 TFumili.cxx:1251
 TFumili.cxx:1252
 TFumili.cxx:1253
 TFumili.cxx:1254
 TFumili.cxx:1255
 TFumili.cxx:1256
 TFumili.cxx:1257
 TFumili.cxx:1258
 TFumili.cxx:1259
 TFumili.cxx:1260
 TFumili.cxx:1261
 TFumili.cxx:1262
 TFumili.cxx:1263
 TFumili.cxx:1264
 TFumili.cxx:1265
 TFumili.cxx:1266
 TFumili.cxx:1267
 TFumili.cxx:1268
 TFumili.cxx:1269
 TFumili.cxx:1270
 TFumili.cxx:1271
 TFumili.cxx:1272
 TFumili.cxx:1273
 TFumili.cxx:1274
 TFumili.cxx:1275
 TFumili.cxx:1276
 TFumili.cxx:1277
 TFumili.cxx:1278
 TFumili.cxx:1279
 TFumili.cxx:1280
 TFumili.cxx:1281
 TFumili.cxx:1282
 TFumili.cxx:1283
 TFumili.cxx:1284
 TFumili.cxx:1285
 TFumili.cxx:1286
 TFumili.cxx:1287
 TFumili.cxx:1288
 TFumili.cxx:1289
 TFumili.cxx:1290
 TFumili.cxx:1291
 TFumili.cxx:1292
 TFumili.cxx:1293
 TFumili.cxx:1294
 TFumili.cxx:1295
 TFumili.cxx:1296
 TFumili.cxx:1297
 TFumili.cxx:1298
 TFumili.cxx:1299
 TFumili.cxx:1300
 TFumili.cxx:1301
 TFumili.cxx:1302
 TFumili.cxx:1303
 TFumili.cxx:1304
 TFumili.cxx:1305
 TFumili.cxx:1306
 TFumili.cxx:1307
 TFumili.cxx:1308
 TFumili.cxx:1309
 TFumili.cxx:1310
 TFumili.cxx:1311
 TFumili.cxx:1312
 TFumili.cxx:1313
 TFumili.cxx:1314
 TFumili.cxx:1315
 TFumili.cxx:1316
 TFumili.cxx:1317
 TFumili.cxx:1318
 TFumili.cxx:1319
 TFumili.cxx:1320
 TFumili.cxx:1321
 TFumili.cxx:1322
 TFumili.cxx:1323
 TFumili.cxx:1324
 TFumili.cxx:1325
 TFumili.cxx:1326
 TFumili.cxx:1327
 TFumili.cxx:1328
 TFumili.cxx:1329
 TFumili.cxx:1330
 TFumili.cxx:1331
 TFumili.cxx:1332
 TFumili.cxx:1333
 TFumili.cxx:1334
 TFumili.cxx:1335
 TFumili.cxx:1336
 TFumili.cxx:1337
 TFumili.cxx:1338
 TFumili.cxx:1339
 TFumili.cxx:1340
 TFumili.cxx:1341
 TFumili.cxx:1342
 TFumili.cxx:1343
 TFumili.cxx:1344
 TFumili.cxx:1345
 TFumili.cxx:1346
 TFumili.cxx:1347
 TFumili.cxx:1348
 TFumili.cxx:1349
 TFumili.cxx:1350
 TFumili.cxx:1351
 TFumili.cxx:1352
 TFumili.cxx:1353
 TFumili.cxx:1354
 TFumili.cxx:1355
 TFumili.cxx:1356
 TFumili.cxx:1357
 TFumili.cxx:1358
 TFumili.cxx:1359
 TFumili.cxx:1360
 TFumili.cxx:1361
 TFumili.cxx:1362
 TFumili.cxx:1363
 TFumili.cxx:1364
 TFumili.cxx:1365
 TFumili.cxx:1366
 TFumili.cxx:1367
 TFumili.cxx:1368
 TFumili.cxx:1369
 TFumili.cxx:1370
 TFumili.cxx:1371
 TFumili.cxx:1372
 TFumili.cxx:1373
 TFumili.cxx:1374
 TFumili.cxx:1375
 TFumili.cxx:1376
 TFumili.cxx:1377
 TFumili.cxx:1378
 TFumili.cxx:1379
 TFumili.cxx:1380
 TFumili.cxx:1381
 TFumili.cxx:1382
 TFumili.cxx:1383
 TFumili.cxx:1384
 TFumili.cxx:1385
 TFumili.cxx:1386
 TFumili.cxx:1387
 TFumili.cxx:1388
 TFumili.cxx:1389
 TFumili.cxx:1390
 TFumili.cxx:1391
 TFumili.cxx:1392
 TFumili.cxx:1393
 TFumili.cxx:1394
 TFumili.cxx:1395
 TFumili.cxx:1396
 TFumili.cxx:1397
 TFumili.cxx:1398
 TFumili.cxx:1399
 TFumili.cxx:1400
 TFumili.cxx:1401
 TFumili.cxx:1402
 TFumili.cxx:1403
 TFumili.cxx:1404
 TFumili.cxx:1405
 TFumili.cxx:1406
 TFumili.cxx:1407
 TFumili.cxx:1408
 TFumili.cxx:1409
 TFumili.cxx:1410
 TFumili.cxx:1411
 TFumili.cxx:1412
 TFumili.cxx:1413
 TFumili.cxx:1414
 TFumili.cxx:1415
 TFumili.cxx:1416
 TFumili.cxx:1417
 TFumili.cxx:1418
 TFumili.cxx:1419
 TFumili.cxx:1420
 TFumili.cxx:1421
 TFumili.cxx:1422
 TFumili.cxx:1423
 TFumili.cxx:1424
 TFumili.cxx:1425
 TFumili.cxx:1426
 TFumili.cxx:1427
 TFumili.cxx:1428
 TFumili.cxx:1429
 TFumili.cxx:1430
 TFumili.cxx:1431
 TFumili.cxx:1432
 TFumili.cxx:1433
 TFumili.cxx:1434
 TFumili.cxx:1435
 TFumili.cxx:1436
 TFumili.cxx:1437
 TFumili.cxx:1438
 TFumili.cxx:1439
 TFumili.cxx:1440
 TFumili.cxx:1441
 TFumili.cxx:1442
 TFumili.cxx:1443
 TFumili.cxx:1444
 TFumili.cxx:1445
 TFumili.cxx:1446
 TFumili.cxx:1447
 TFumili.cxx:1448
 TFumili.cxx:1449
 TFumili.cxx:1450
 TFumili.cxx:1451
 TFumili.cxx:1452
 TFumili.cxx:1453
 TFumili.cxx:1454
 TFumili.cxx:1455
 TFumili.cxx:1456
 TFumili.cxx:1457
 TFumili.cxx:1458
 TFumili.cxx:1459
 TFumili.cxx:1460
 TFumili.cxx:1461
 TFumili.cxx:1462
 TFumili.cxx:1463
 TFumili.cxx:1464
 TFumili.cxx:1465
 TFumili.cxx:1466
 TFumili.cxx:1467
 TFumili.cxx:1468
 TFumili.cxx:1469
 TFumili.cxx:1470
 TFumili.cxx:1471
 TFumili.cxx:1472
 TFumili.cxx:1473
 TFumili.cxx:1474
 TFumili.cxx:1475
 TFumili.cxx:1476
 TFumili.cxx:1477
 TFumili.cxx:1478
 TFumili.cxx:1479
 TFumili.cxx:1480
 TFumili.cxx:1481
 TFumili.cxx:1482
 TFumili.cxx:1483
 TFumili.cxx:1484
 TFumili.cxx:1485
 TFumili.cxx:1486
 TFumili.cxx:1487
 TFumili.cxx:1488
 TFumili.cxx:1489
 TFumili.cxx:1490
 TFumili.cxx:1491
 TFumili.cxx:1492
 TFumili.cxx:1493
 TFumili.cxx:1494
 TFumili.cxx:1495
 TFumili.cxx:1496
 TFumili.cxx:1497
 TFumili.cxx:1498
 TFumili.cxx:1499
 TFumili.cxx:1500
 TFumili.cxx:1501
 TFumili.cxx:1502
 TFumili.cxx:1503
 TFumili.cxx:1504
 TFumili.cxx:1505
 TFumili.cxx:1506
 TFumili.cxx:1507
 TFumili.cxx:1508
 TFumili.cxx:1509
 TFumili.cxx:1510
 TFumili.cxx:1511
 TFumili.cxx:1512
 TFumili.cxx:1513
 TFumili.cxx:1514
 TFumili.cxx:1515
 TFumili.cxx:1516
 TFumili.cxx:1517
 TFumili.cxx:1518
 TFumili.cxx:1519
 TFumili.cxx:1520
 TFumili.cxx:1521
 TFumili.cxx:1522
 TFumili.cxx:1523
 TFumili.cxx:1524
 TFumili.cxx:1525
 TFumili.cxx:1526
 TFumili.cxx:1527
 TFumili.cxx:1528
 TFumili.cxx:1529
 TFumili.cxx:1530
 TFumili.cxx:1531
 TFumili.cxx:1532
 TFumili.cxx:1533
 TFumili.cxx:1534
 TFumili.cxx:1535
 TFumili.cxx:1536
 TFumili.cxx:1537
 TFumili.cxx:1538
 TFumili.cxx:1539
 TFumili.cxx:1540
 TFumili.cxx:1541
 TFumili.cxx:1542
 TFumili.cxx:1543
 TFumili.cxx:1544
 TFumili.cxx:1545
 TFumili.cxx:1546
 TFumili.cxx:1547
 TFumili.cxx:1548
 TFumili.cxx:1549
 TFumili.cxx:1550
 TFumili.cxx:1551
 TFumili.cxx:1552
 TFumili.cxx:1553
 TFumili.cxx:1554
 TFumili.cxx:1555
 TFumili.cxx:1556
 TFumili.cxx:1557
 TFumili.cxx:1558
 TFumili.cxx:1559
 TFumili.cxx:1560
 TFumili.cxx:1561
 TFumili.cxx:1562
 TFumili.cxx:1563
 TFumili.cxx:1564
 TFumili.cxx:1565
 TFumili.cxx:1566
 TFumili.cxx:1567
 TFumili.cxx:1568
 TFumili.cxx:1569
 TFumili.cxx:1570
 TFumili.cxx:1571
 TFumili.cxx:1572
 TFumili.cxx:1573
 TFumili.cxx:1574
 TFumili.cxx:1575
 TFumili.cxx:1576
 TFumili.cxx:1577
 TFumili.cxx:1578
 TFumili.cxx:1579
 TFumili.cxx:1580
 TFumili.cxx:1581
 TFumili.cxx:1582
 TFumili.cxx:1583
 TFumili.cxx:1584
 TFumili.cxx:1585
 TFumili.cxx:1586
 TFumili.cxx:1587
 TFumili.cxx:1588
 TFumili.cxx:1589
 TFumili.cxx:1590
 TFumili.cxx:1591
 TFumili.cxx:1592
 TFumili.cxx:1593
 TFumili.cxx:1594
 TFumili.cxx:1595
 TFumili.cxx:1596
 TFumili.cxx:1597
 TFumili.cxx:1598
 TFumili.cxx:1599
 TFumili.cxx:1600
 TFumili.cxx:1601
 TFumili.cxx:1602
 TFumili.cxx:1603
 TFumili.cxx:1604
 TFumili.cxx:1605
 TFumili.cxx:1606
 TFumili.cxx:1607
 TFumili.cxx:1608
 TFumili.cxx:1609
 TFumili.cxx:1610
 TFumili.cxx:1611
 TFumili.cxx:1612
 TFumili.cxx:1613
 TFumili.cxx:1614
 TFumili.cxx:1615
 TFumili.cxx:1616
 TFumili.cxx:1617
 TFumili.cxx:1618
 TFumili.cxx:1619
 TFumili.cxx:1620
 TFumili.cxx:1621
 TFumili.cxx:1622
 TFumili.cxx:1623
 TFumili.cxx:1624
 TFumili.cxx:1625
 TFumili.cxx:1626
 TFumili.cxx:1627
 TFumili.cxx:1628
 TFumili.cxx:1629
 TFumili.cxx:1630
 TFumili.cxx:1631
 TFumili.cxx:1632
 TFumili.cxx:1633
 TFumili.cxx:1634
 TFumili.cxx:1635
 TFumili.cxx:1636
 TFumili.cxx:1637
 TFumili.cxx:1638
 TFumili.cxx:1639
 TFumili.cxx:1640
 TFumili.cxx:1641
 TFumili.cxx:1642
 TFumili.cxx:1643
 TFumili.cxx:1644
 TFumili.cxx:1645
 TFumili.cxx:1646
 TFumili.cxx:1647
 TFumili.cxx:1648
 TFumili.cxx:1649
 TFumili.cxx:1650
 TFumili.cxx:1651
 TFumili.cxx:1652
 TFumili.cxx:1653
 TFumili.cxx:1654
 TFumili.cxx:1655
 TFumili.cxx:1656
 TFumili.cxx:1657
 TFumili.cxx:1658
 TFumili.cxx:1659
 TFumili.cxx:1660
 TFumili.cxx:1661
 TFumili.cxx:1662
 TFumili.cxx:1663
 TFumili.cxx:1664
 TFumili.cxx:1665
 TFumili.cxx:1666
 TFumili.cxx:1667
 TFumili.cxx:1668
 TFumili.cxx:1669
 TFumili.cxx:1670
 TFumili.cxx:1671
 TFumili.cxx:1672
 TFumili.cxx:1673
 TFumili.cxx:1674
 TFumili.cxx:1675
 TFumili.cxx:1676
 TFumili.cxx:1677
 TFumili.cxx:1678
 TFumili.cxx:1679
 TFumili.cxx:1680
 TFumili.cxx:1681
 TFumili.cxx:1682
 TFumili.cxx:1683
 TFumili.cxx:1684
 TFumili.cxx:1685
 TFumili.cxx:1686
 TFumili.cxx:1687
 TFumili.cxx:1688
 TFumili.cxx:1689
 TFumili.cxx:1690
 TFumili.cxx:1691
 TFumili.cxx:1692
 TFumili.cxx:1693
 TFumili.cxx:1694
 TFumili.cxx:1695
 TFumili.cxx:1696
 TFumili.cxx:1697
 TFumili.cxx:1698
 TFumili.cxx:1699
 TFumili.cxx:1700
 TFumili.cxx:1701
 TFumili.cxx:1702
 TFumili.cxx:1703
 TFumili.cxx:1704
 TFumili.cxx:1705
 TFumili.cxx:1706
 TFumili.cxx:1707
 TFumili.cxx:1708
 TFumili.cxx:1709
 TFumili.cxx:1710
 TFumili.cxx:1711
 TFumili.cxx:1712
 TFumili.cxx:1713
 TFumili.cxx:1714
 TFumili.cxx:1715
 TFumili.cxx:1716
 TFumili.cxx:1717
 TFumili.cxx:1718
 TFumili.cxx:1719
 TFumili.cxx:1720
 TFumili.cxx:1721
 TFumili.cxx:1722
 TFumili.cxx:1723
 TFumili.cxx:1724
 TFumili.cxx:1725
 TFumili.cxx:1726
 TFumili.cxx:1727
 TFumili.cxx:1728
 TFumili.cxx:1729
 TFumili.cxx:1730
 TFumili.cxx:1731
 TFumili.cxx:1732
 TFumili.cxx:1733
 TFumili.cxx:1734
 TFumili.cxx:1735
 TFumili.cxx:1736
 TFumili.cxx:1737
 TFumili.cxx:1738
 TFumili.cxx:1739
 TFumili.cxx:1740
 TFumili.cxx:1741
 TFumili.cxx:1742
 TFumili.cxx:1743
 TFumili.cxx:1744
 TFumili.cxx:1745
 TFumili.cxx:1746
 TFumili.cxx:1747
 TFumili.cxx:1748
 TFumili.cxx:1749
 TFumili.cxx:1750
 TFumili.cxx:1751
 TFumili.cxx:1752
 TFumili.cxx:1753
 TFumili.cxx:1754
 TFumili.cxx:1755
 TFumili.cxx:1756
 TFumili.cxx:1757
 TFumili.cxx:1758
 TFumili.cxx:1759
 TFumili.cxx:1760
 TFumili.cxx:1761
 TFumili.cxx:1762
 TFumili.cxx:1763
 TFumili.cxx:1764
 TFumili.cxx:1765
 TFumili.cxx:1766
 TFumili.cxx:1767
 TFumili.cxx:1768
 TFumili.cxx:1769
 TFumili.cxx:1770
 TFumili.cxx:1771
 TFumili.cxx:1772
 TFumili.cxx:1773
 TFumili.cxx:1774
 TFumili.cxx:1775
 TFumili.cxx:1776
 TFumili.cxx:1777
 TFumili.cxx:1778
 TFumili.cxx:1779
 TFumili.cxx:1780
 TFumili.cxx:1781
 TFumili.cxx:1782
 TFumili.cxx:1783
 TFumili.cxx:1784
 TFumili.cxx:1785
 TFumili.cxx:1786
 TFumili.cxx:1787
 TFumili.cxx:1788
 TFumili.cxx:1789
 TFumili.cxx:1790
 TFumili.cxx:1791
 TFumili.cxx:1792
 TFumili.cxx:1793
 TFumili.cxx:1794
 TFumili.cxx:1795
 TFumili.cxx:1796
 TFumili.cxx:1797
 TFumili.cxx:1798
 TFumili.cxx:1799
 TFumili.cxx:1800
 TFumili.cxx:1801
 TFumili.cxx:1802
 TFumili.cxx:1803
 TFumili.cxx:1804
 TFumili.cxx:1805
 TFumili.cxx:1806
 TFumili.cxx:1807
 TFumili.cxx:1808
 TFumili.cxx:1809
 TFumili.cxx:1810
 TFumili.cxx:1811
 TFumili.cxx:1812
 TFumili.cxx:1813
 TFumili.cxx:1814
 TFumili.cxx:1815
 TFumili.cxx:1816
 TFumili.cxx:1817
 TFumili.cxx:1818
 TFumili.cxx:1819
 TFumili.cxx:1820
 TFumili.cxx:1821
 TFumili.cxx:1822
 TFumili.cxx:1823
 TFumili.cxx:1824
 TFumili.cxx:1825
 TFumili.cxx:1826
 TFumili.cxx:1827
 TFumili.cxx:1828
 TFumili.cxx:1829
 TFumili.cxx:1830
 TFumili.cxx:1831
 TFumili.cxx:1832
 TFumili.cxx:1833
 TFumili.cxx:1834
 TFumili.cxx:1835
 TFumili.cxx:1836
 TFumili.cxx:1837
 TFumili.cxx:1838
 TFumili.cxx:1839
 TFumili.cxx:1840
 TFumili.cxx:1841
 TFumili.cxx:1842
 TFumili.cxx:1843
 TFumili.cxx:1844
 TFumili.cxx:1845
 TFumili.cxx:1846
 TFumili.cxx:1847
 TFumili.cxx:1848
 TFumili.cxx:1849
 TFumili.cxx:1850
 TFumili.cxx:1851
 TFumili.cxx:1852
 TFumili.cxx:1853
 TFumili.cxx:1854
 TFumili.cxx:1855
 TFumili.cxx:1856
 TFumili.cxx:1857
 TFumili.cxx:1858
 TFumili.cxx:1859
 TFumili.cxx:1860
 TFumili.cxx:1861
 TFumili.cxx:1862
 TFumili.cxx:1863
 TFumili.cxx:1864
 TFumili.cxx:1865
 TFumili.cxx:1866
 TFumili.cxx:1867
 TFumili.cxx:1868
 TFumili.cxx:1869
 TFumili.cxx:1870
 TFumili.cxx:1871
 TFumili.cxx:1872
 TFumili.cxx:1873
 TFumili.cxx:1874
 TFumili.cxx:1875
 TFumili.cxx:1876
 TFumili.cxx:1877
 TFumili.cxx:1878
 TFumili.cxx:1879
 TFumili.cxx:1880
 TFumili.cxx:1881
 TFumili.cxx:1882
 TFumili.cxx:1883
 TFumili.cxx:1884
 TFumili.cxx:1885
 TFumili.cxx:1886
 TFumili.cxx:1887
 TFumili.cxx:1888
 TFumili.cxx:1889
 TFumili.cxx:1890
 TFumili.cxx:1891
 TFumili.cxx:1892
 TFumili.cxx:1893
 TFumili.cxx:1894
 TFumili.cxx:1895
 TFumili.cxx:1896
 TFumili.cxx:1897
 TFumili.cxx:1898
 TFumili.cxx:1899
 TFumili.cxx:1900
 TFumili.cxx:1901
 TFumili.cxx:1902
 TFumili.cxx:1903
 TFumili.cxx:1904
 TFumili.cxx:1905
 TFumili.cxx:1906
 TFumili.cxx:1907
 TFumili.cxx:1908
 TFumili.cxx:1909
 TFumili.cxx:1910
 TFumili.cxx:1911
 TFumili.cxx:1912
 TFumili.cxx:1913
 TFumili.cxx:1914
 TFumili.cxx:1915
 TFumili.cxx:1916
 TFumili.cxx:1917
 TFumili.cxx:1918
 TFumili.cxx:1919
 TFumili.cxx:1920
 TFumili.cxx:1921
 TFumili.cxx:1922
 TFumili.cxx:1923
 TFumili.cxx:1924
 TFumili.cxx:1925
 TFumili.cxx:1926
 TFumili.cxx:1927
 TFumili.cxx:1928
 TFumili.cxx:1929
 TFumili.cxx:1930
 TFumili.cxx:1931
 TFumili.cxx:1932
 TFumili.cxx:1933
 TFumili.cxx:1934
 TFumili.cxx:1935
 TFumili.cxx:1936
 TFumili.cxx:1937
 TFumili.cxx:1938
 TFumili.cxx:1939
 TFumili.cxx:1940
 TFumili.cxx:1941
 TFumili.cxx:1942
 TFumili.cxx:1943
 TFumili.cxx:1944
 TFumili.cxx:1945
 TFumili.cxx:1946
 TFumili.cxx:1947
 TFumili.cxx:1948
 TFumili.cxx:1949
 TFumili.cxx:1950
 TFumili.cxx:1951
 TFumili.cxx:1952
 TFumili.cxx:1953
 TFumili.cxx:1954
 TFumili.cxx:1955
 TFumili.cxx:1956
 TFumili.cxx:1957
 TFumili.cxx:1958
 TFumili.cxx:1959
 TFumili.cxx:1960
 TFumili.cxx:1961
 TFumili.cxx:1962
 TFumili.cxx:1963
 TFumili.cxx:1964
 TFumili.cxx:1965
 TFumili.cxx:1966
 TFumili.cxx:1967
 TFumili.cxx:1968
 TFumili.cxx:1969
 TFumili.cxx:1970
 TFumili.cxx:1971
 TFumili.cxx:1972
 TFumili.cxx:1973
 TFumili.cxx:1974
 TFumili.cxx:1975
 TFumili.cxx:1976
 TFumili.cxx:1977
 TFumili.cxx:1978
 TFumili.cxx:1979
 TFumili.cxx:1980
 TFumili.cxx:1981
 TFumili.cxx:1982
 TFumili.cxx:1983
 TFumili.cxx:1984
 TFumili.cxx:1985
 TFumili.cxx:1986
 TFumili.cxx:1987
 TFumili.cxx:1988
 TFumili.cxx:1989
 TFumili.cxx:1990
 TFumili.cxx:1991
 TFumili.cxx:1992
 TFumili.cxx:1993
 TFumili.cxx:1994
 TFumili.cxx:1995
 TFumili.cxx:1996
 TFumili.cxx:1997
 TFumili.cxx:1998
 TFumili.cxx:1999
 TFumili.cxx:2000
 TFumili.cxx:2001
 TFumili.cxx:2002
 TFumili.cxx:2003
 TFumili.cxx:2004
 TFumili.cxx:2005
 TFumili.cxx:2006
 TFumili.cxx:2007
 TFumili.cxx:2008
 TFumili.cxx:2009
 TFumili.cxx:2010
 TFumili.cxx:2011
 TFumili.cxx:2012
 TFumili.cxx:2013
 TFumili.cxx:2014
 TFumili.cxx:2015
 TFumili.cxx:2016
 TFumili.cxx:2017
 TFumili.cxx:2018
 TFumili.cxx:2019
 TFumili.cxx:2020
 TFumili.cxx:2021
 TFumili.cxx:2022
 TFumili.cxx:2023
 TFumili.cxx:2024
 TFumili.cxx:2025
 TFumili.cxx:2026
 TFumili.cxx:2027
 TFumili.cxx:2028
 TFumili.cxx:2029
 TFumili.cxx:2030
 TFumili.cxx:2031
 TFumili.cxx:2032
 TFumili.cxx:2033
 TFumili.cxx:2034
 TFumili.cxx:2035
 TFumili.cxx:2036
 TFumili.cxx:2037
 TFumili.cxx:2038
 TFumili.cxx:2039
 TFumili.cxx:2040
 TFumili.cxx:2041
 TFumili.cxx:2042
 TFumili.cxx:2043
 TFumili.cxx:2044
 TFumili.cxx:2045
 TFumili.cxx:2046
 TFumili.cxx:2047
 TFumili.cxx:2048
 TFumili.cxx:2049
 TFumili.cxx:2050
 TFumili.cxx:2051
 TFumili.cxx:2052
 TFumili.cxx:2053
 TFumili.cxx:2054
 TFumili.cxx:2055
 TFumili.cxx:2056
 TFumili.cxx:2057
 TFumili.cxx:2058
 TFumili.cxx:2059
 TFumili.cxx:2060
 TFumili.cxx:2061
 TFumili.cxx:2062
 TFumili.cxx:2063
 TFumili.cxx:2064
 TFumili.cxx:2065
 TFumili.cxx:2066
 TFumili.cxx:2067
 TFumili.cxx:2068
 TFumili.cxx:2069
 TFumili.cxx:2070
 TFumili.cxx:2071
 TFumili.cxx:2072
 TFumili.cxx:2073
 TFumili.cxx:2074
 TFumili.cxx:2075
 TFumili.cxx:2076
 TFumili.cxx:2077
 TFumili.cxx:2078
 TFumili.cxx:2079
 TFumili.cxx:2080
 TFumili.cxx:2081
 TFumili.cxx:2082
 TFumili.cxx:2083
 TFumili.cxx:2084
 TFumili.cxx:2085
 TFumili.cxx:2086
 TFumili.cxx:2087
 TFumili.cxx:2088
 TFumili.cxx:2089
 TFumili.cxx:2090
 TFumili.cxx:2091
 TFumili.cxx:2092
 TFumili.cxx:2093
 TFumili.cxx:2094
 TFumili.cxx:2095
 TFumili.cxx:2096
 TFumili.cxx:2097
 TFumili.cxx:2098
 TFumili.cxx:2099
 TFumili.cxx:2100
 TFumili.cxx:2101
 TFumili.cxx:2102
 TFumili.cxx:2103
 TFumili.cxx:2104
 TFumili.cxx:2105
 TFumili.cxx:2106
 TFumili.cxx:2107
 TFumili.cxx:2108
 TFumili.cxx:2109
 TFumili.cxx:2110
 TFumili.cxx:2111
 TFumili.cxx:2112
 TFumili.cxx:2113
 TFumili.cxx:2114
 TFumili.cxx:2115
 TFumili.cxx:2116
 TFumili.cxx:2117
 TFumili.cxx:2118
 TFumili.cxx:2119
 TFumili.cxx:2120
 TFumili.cxx:2121
 TFumili.cxx:2122
 TFumili.cxx:2123
 TFumili.cxx:2124
 TFumili.cxx:2125
 TFumili.cxx:2126
 TFumili.cxx:2127
 TFumili.cxx:2128
 TFumili.cxx:2129
 TFumili.cxx:2130
 TFumili.cxx:2131
 TFumili.cxx:2132
 TFumili.cxx:2133
 TFumili.cxx:2134
 TFumili.cxx:2135
 TFumili.cxx:2136
 TFumili.cxx:2137
 TFumili.cxx:2138
 TFumili.cxx:2139
 TFumili.cxx:2140
 TFumili.cxx:2141
 TFumili.cxx:2142
 TFumili.cxx:2143
 TFumili.cxx:2144
 TFumili.cxx:2145
 TFumili.cxx:2146
 TFumili.cxx:2147
 TFumili.cxx:2148
 TFumili.cxx:2149
 TFumili.cxx:2150
 TFumili.cxx:2151
 TFumili.cxx:2152
 TFumili.cxx:2153
 TFumili.cxx:2154
 TFumili.cxx:2155
 TFumili.cxx:2156
 TFumili.cxx:2157
 TFumili.cxx:2158
 TFumili.cxx:2159
 TFumili.cxx:2160
 TFumili.cxx:2161
 TFumili.cxx:2162
 TFumili.cxx:2163
 TFumili.cxx:2164
 TFumili.cxx:2165
 TFumili.cxx:2166
 TFumili.cxx:2167
 TFumili.cxx:2168
 TFumili.cxx:2169
 TFumili.cxx:2170
 TFumili.cxx:2171
 TFumili.cxx:2172
 TFumili.cxx:2173
 TFumili.cxx:2174
 TFumili.cxx:2175
 TFumili.cxx:2176
 TFumili.cxx:2177
 TFumili.cxx:2178
 TFumili.cxx:2179
 TFumili.cxx:2180
 TFumili.cxx:2181
 TFumili.cxx:2182
 TFumili.cxx:2183
 TFumili.cxx:2184
 TFumili.cxx:2185
 TFumili.cxx:2186
 TFumili.cxx:2187
 TFumili.cxx:2188
 TFumili.cxx:2189
 TFumili.cxx:2190
 TFumili.cxx:2191
 TFumili.cxx:2192
 TFumili.cxx:2193
 TFumili.cxx:2194
 TFumili.cxx:2195
 TFumili.cxx:2196
 TFumili.cxx:2197
 TFumili.cxx:2198
 TFumili.cxx:2199
 TFumili.cxx:2200
 TFumili.cxx:2201
 TFumili.cxx:2202
 TFumili.cxx:2203
 TFumili.cxx:2204