Newer
Older
Master_thesis / testfit.py
  1. # Example of a fitting program
  2. # ============================
  3. #
  4. # The fitting function fcn is a simple chisquare function
  5. # The data consists of 5 data points (arrays x,y,z) + the errors in errorsz
  6. # More details on the various functions or parameters for these functions
  7. # can be obtained in an interactive ROOT session with:
  8. # Root > TMinuit *minuit = new TMinuit(10);
  9. # Root > minuit->mnhelp("*",0) to see the list of possible keywords
  10. # Root > minuit->mnhelp("SET",0) explains most parameters
  11.  
  12.  
  13. from ROOT import *
  14. from array import array;
  15.  
  16. Error = 0;
  17. z = array( 'f', ( 1., 0.96, 0.89, 0.85, 0.78 ) )
  18. errorz = array( 'f', 5*[0.01] )
  19.  
  20. x = array( 'f', ( 1.5751, 1.5825, 1.6069, 1.6339, 1.6706 ) )
  21. y = array( 'f', ( 1.0642, 0.97685, 1.13168, 1.128654, 1.44016 ) )
  22.  
  23. ncount = 0
  24.  
  25. ##______________________________________________________________________________
  26. def testfit():
  27.  
  28. gMinuit = TMinuit(5)
  29. gMinuit.SetFCN( fcn )
  30.  
  31. arglist = array( 'd', 10*[0.] )
  32. ierflg = Long(0)
  33.  
  34. arglist[0] = 1
  35. gMinuit.mnexcm( "SET ERR", arglist, 1, ierflg )
  36.  
  37. # Set starting values and step sizes for parameters
  38. vstart = array( 'd', ( 3, 1, 0.1, 0.01 ) )
  39. step = array( 'd', ( 0.1, 0.1, 0.01, 0.001 ) )
  40. gMinuit.mnparm( 0, "a1", vstart[0], step[0], 0, 0, ierflg )
  41. gMinuit.mnparm( 1, "a2", vstart[1], step[1], 0, 0, ierflg )
  42. gMinuit.mnparm( 2, "a3", vstart[2], step[2], 0, 0, ierflg )
  43. gMinuit.mnparm( 3, "a4", vstart[3], step[3], 0, 0, ierflg )
  44.  
  45. # Now ready for minimization step
  46. arglist[0] = 500
  47. arglist[1] = 1.
  48. gMinuit.mnexcm( "MIGRAD", arglist, 2, ierflg )
  49.  
  50. # Print results
  51. amin, edm, errdef = 0.18, 0.19, 0.20
  52. nvpar, nparx, icstat = 1983, 1984, 1985
  53. gMinuit.mnstat( amin, edm, errdef, nvpar, nparx, icstat )
  54. gMinuit.mnprin( 3, amin )
  55.  
  56.  
  57. ##______________________________________________________________________________
  58. def fcn( npar, gin, f, par, iflag ):
  59. global ncount
  60. nbins = 5
  61.  
  62. # calculate chisquare
  63. chisq, delta = 0., 0.
  64. for i in range(nbins):
  65. delta = (z[i]-func(x[i],y[i],par))/errorz[i]
  66. chisq += delta*delta
  67.  
  68. f[0] = chisq
  69. ncount += 1
  70.  
  71. def func( x, y, par ):
  72. value = ( (par[0]*par[0])/(x*x)-1)/ ( par[1]+par[2]*y-par[3]*y*y)
  73. return value
  74.  
  75.  
  76. ##______________________________________________________________________________
  77. if __name__ == '__main__':
  78. testfit()