[BACK]Return to fit.dem CVS log [TXT][DIR] Up to [local] / OpenXM_contrib / gnuplot / demo

Annotation of OpenXM_contrib/gnuplot/demo/fit.dem, Revision 1.1.1.1

1.1       maekawa     1: #
                      2: # $Id: fit.dem,v 1.9 1997/05/27 01:29:34 drd Exp $
                      3: #
                      4:
                      5: print "Some examples how data fitting using nonlinear least squares fit"
                      6: print "can be done."
                      7: print ""
                      8: pause -1 "first plotting the pure data set  (-> return)"
                      9:
                     10: set title 'data for first fit demo'
                     11: plot 'lcdemo.dat'
                     12: set xlabel "Temperature T  [deg Cels.]"
                     13: set ylabel "Density [g/cm3]"
                     14:
                     15: print "now fitting a straight line to the data :-)"
                     16: print "only as a demo without physical meaning"
                     17: load 'line.fnc'
                     18: y0 = 0.0
                     19: m = 0.0
                     20: show variables
                     21: pause -1 "first a plot with all parameters set to zero  (-> return)"
                     22: set title 'all fit params set to 0'
                     23: plot 'lcdemo.dat', l(x)
                     24: pause -1 "now start fitting...  (-> return)"
                     25: fit l(x) 'lcdemo.dat' via y0, m
                     26: pause -1 "now look at the result (-> return)"
                     27: set title 'unweighted fit'
                     28: plot 'lcdemo.dat', l(x)
                     29:
                     30: pause -1 "see the influence of weights for single data points   (-> return)"
                     31: fit l(x) 'lcdemo.dat' using 1:2:3 via y0, m
                     32: pause -1 "now look at the result (-> return)"
                     33: set title 'fit weighted towards low temperatures'
                     34: plot 'lcdemo.dat', l(x)
                     35:
                     36: pause -1 "now prefer the high temperature data   (-> return)"
                     37: fit l(x) 'lcdemo.dat' using 1:2:4 via y0, m
                     38: pause -1 "now look at the result (-> return)"
                     39: set title 'bias to high-temperates'
                     40: plot 'lcdemo.dat', l(x)
                     41:
                     42: pause  0 "now use real single-measurement errors to reach such a result (-> return)"
                     43: pause  0 "(look at the file lcdemo.dat and compare the columns to see the difference)"
                     44: pause -1 "(-> return)"
                     45: set title 'data with experimental errors'
                     46: plot 'lcdemo.dat' using 1:2:5 with errorbars
                     47: fit l(x) 'lcdemo.dat' using 1:2:5 via y0, m
                     48: pause -1 "now look at the result (-> return)"
                     49: set title 'fit weighted by experimental errors'
                     50: plot 'lcdemo.dat' using 1:2:5 with errorbars, l(x)
                     51:
                     52: print "It's time now to try a more realistic model function"
                     53: load 'density.fnc'
                     54: show functions
                     55: print "density(x) is a function which shall fit the whole temperature"
                     56: print "range using a ?: expression. It contains 6 model parameters which
                     57: print "will all be varied. Now take the start parameters out of the"
                     58: pause -1 "file 'start.par' and plot the function    (-> return)"
                     59: load 'start.par'
                     60: set title 'initial parameters for realistic model function'
                     61: plot 'lcdemo.dat', density(x)
                     62: fit density(x) 'lcdemo.dat' via 'start.par'
                     63: pause -1 "now look at the result (-> return)"
                     64: set title 'fitted to realistic model function'
                     65: plot 'lcdemo.dat', density(x)
                     66:
                     67: print  "looks already rather nice? We will do now the following: set"
                     68: print  "the epsilon limit higher so that we need more iteration steps"
                     69: print  "to convergence. During fitting please hit ctrl-C. You will be asked"
                     70: print  "Stop, Continue, Execute: Try everything. You may define a script"
                     71: print  "using the FIT_SCRIPT environment variable. An example would be"
                     72: print  "'FIT_SCRIPT=plot nonsense.dat'. Normally you don't need to set"
                     73: print  "FIT_SCRIPT since it defaults to 'replot'. Please note that FIT_SCRIPT"
                     74: print  "cannot be set from inside gnuplot."
                     75: print  ""
                     76: pause -1  "(-> return)"
                     77: FIT_LIMIT = 1e-10
                     78: fit density(x) 'lcdemo.dat' via 'start.par'
                     79: pause -1 "now look at the result (-> return)"
                     80: set title 'fit with more iterations'
                     81: plot 'lcdemo.dat', density(x)
                     82:
                     83: FIT_LIMIT = 1e-5
                     84: print "\nNow a brief demonstration of 3d fitting."
                     85: print "hemisphr.dat contains random points on a hemisphere of"
                     86: print "radius 1, but we let fit figure this out for us."
                     87: print "It takes many iterations, so we limit FIT_MAXITER to 50."
                     88: #HBB: made this a lot harder: also fit the center of the sphere
                     89: #h(x,y) = sqrt(r*r - (x-x0)**2 - (y-y0)**2) + z0
                     90: #HBB 970522: distort the function, so it won't fit exactly:
                     91: h(x,y) = sqrt(r*r - (abs(x-x0))**2.2 - (abs(y-y0))**1.8) + z0
                     92: x0 = 0.1
                     93: y0 = 0.2
                     94: z0 = 0.3
                     95: r=0.5
                     96: FIT_MAXITER=50
                     97: set title 'the scattered points, and the initial parameter'
                     98: splot 'hemisphr.dat' using 1:2:3, h(x,y)
                     99: pause -1 "(-> return)"
                    100:
                    101: # we *must* provide 4 columns for a 3d fit. We fake errors=1
                    102: fit h(x,y) 'hemisphr.dat' using 1:2:3:(1) via r, x0, y0, z0
                    103: set title 'the scattered points, fitted curve'
                    104: splot 'hemisphr.dat' using 1:2:3, h(x,y)
                    105: print "\n\nNotice, however, that this would converge much faster when"
                    106: print "fitted in a more appropriate co-ordinate system:"
                    107: print "fit r 'hemisphr.dat' using 0:($1*$1+$2*$2+$3*$3) via r"
                    108: print "where we are fitting f(x)=r to the radii calculated as the data"
                    109: print "is read from the file. No x value is required in this case.
                    110: pause -1 "(This is left as an excercise for the user). (-> return)"
                    111: FIT_MAXITER=0   # no limit : we cannot delete the variable once set
                    112:
                    113: print "\n\nNow an example how to fit multi-branch functions\n"
                    114: print  "The model consists of two branches, the first describing longitudinal"
                    115: print  "sound velocity as function of propagation direction (upper data),"
                    116: print  "the second describing transverse sound velocity (lower data).\n"
                    117: print  "The model uses these data in order to fit elastic stiffnesses"
                    118: print  "which occur differently in both branches.\n"
                    119: pause -1 "(-> return)"
                    120: load 'hexa.fnc'
                    121: load 'sound.par'
                    122: set title 'sound data, and model with initial parameters'
                    123: plot 'soundvel.dat', vlong(x), vtrans(x)
                    124: # Must provide an error estimate for a 3d fit. Use constant 1
                    125: fit f(x,y) 'soundvel.dat' using 1:-2:2:(1) via 'sound.par'
                    126: #create soundfit.par, reading from sound.par and updating values
                    127: update 'sound.par' 'soundfit.par'
                    128: print  ""
                    129: pause -1 "(-> return)"
                    130: set title 'pseudo-3d multi-branch fit to velocity data'
                    131: plot 'soundvel.dat', vlong(x), vtrans(x)
                    132: print  "Look at the file 'hexa.fnc' to see how the branches are realized"
                    133: print  "using the data index as a pseudo-3d fit"
                    134: print  ""
                    135: print  "Next we only use every fifth data point for fitting by using the"
                    136: print  "'every' keyword. Look at the fitting-speed increase and at"
                    137: print  "fitting result."
                    138: print  ""
                    139: pause -1 "(-> return)"
                    140: load 'sound.par'
                    141: fit f(x,y) 'soundvel.dat' every 5 using 1:-2:2:(1) via 'sound.par'
                    142: set title 'fitted only every 5th data point'
                    143: plot 'soundvel.dat', vlong(x), vtrans(x)
                    144: print  "When you compare the results (see 'fit.log') you remark that"
                    145: print  "the uncertainties in the fitted constants have become larger,"
                    146: print  "the quality of the plot is only slightly affected."
                    147: print  ""
                    148: print  "By marking some parameters as '# FIXED' in the parameter file"
                    149: print  "you fit only the others (c44 and c13 fixed here)."
                    150: print  ""
                    151: pause -1 "(-> return)"
                    152: load 'sound2.par'
                    153: set title 'initial parameters'
                    154: plot 'soundvel.dat', vlong(x), vtrans(x)
                    155: fit f(x,y) 'soundvel.dat' using 1:-2:2:(1) via 'sound2.par'
                    156: set title 'fit with c44 and c13 fixed'
                    157: plot 'soundvel.dat', vlong(x), vtrans(x)
                    158: print  "This has the same effect as specifying only the real free"
                    159: print  "parameters by the 'via' syntax."
                    160: print  ""
                    161: print  "fit f(x) 'soundvel.dat' via c33, c11, phi0"
                    162: print  ""
                    163: pause -1 "(-> return)"
                    164: load 'sound.par'
                    165: set title 'initial parameters'
                    166: plot 'soundvel.dat', vlong(x), vtrans(x)
                    167: fit f(x,y) 'soundvel.dat' using 1:-2:2:(1) via c33, c11, phi0
                    168: set title 'fit via c33,c11,phi0'
                    169: plot 'soundvel.dat', vlong(x), vtrans(x)
                    170:
                    171: print  "Here comes an example of a very complex function..."
                    172: print  ""
                    173: pause -1 "first plotting the pure data set  (-> return)"
                    174:
                    175: set xlabel "Delta [degrees]"
                    176: set ylabel "Reflectivity"
                    177: set title 'raw data'
                    178: #HBB 970522: here and below, use the error column present in moli3.dat:
                    179: plot 'moli3.dat' w e
                    180:
                    181: print "now fitting the model function to the data"
                    182: load 'reflect.fnc'
                    183:
                    184: #HBB 970522: Changed initial values to something sensible, i.e.
                    185: #  something an experienced user of fit would actually use.
                    186: #  FIT_LIMIT is also raised, to ensure a better fit.
                    187: eta = 1.2e-4
                    188: tc = 1.8e-3
                    189: FIT_LIMIT=1e-10
                    190:
                    191: show variables
                    192: show functions
                    193: pause -1 "first a plot with all parameters set to initial values  (-> return)"
                    194: set title 'initial parameters'
                    195: plot 'moli3.dat' w e, R(x)
                    196: pause -1 "now start fitting...  (-> return)"
                    197: fit R(x) 'moli3.dat' u 1:2:3 via eta, tc
                    198: pause -1 "now look at the result (-> return)"
                    199: set title 'fitted parameters'
                    200: replot
                    201:
                    202: #HBB 970522: added comment on result of last fit.
                    203: print "Looking at the plot of the resulting fit curve, you can see"
                    204: print "that this function doesn't really fit this set of data points."
                    205: print "This would normally be a reason to check for measurement problems"
                    206: print "not yet accounted for, and maybe even re-think the theoretic"
                    207: print "prediction in use."
                    208: print ""
                    209:
                    210: print  "You can have a look at all previous fit results by looking into"
                    211: print  "the file 'fit.log' or whatever you defined the env-variable 'FIT_LOGFILE'."
                    212: print  "Remember that this file will always be appended, so remove it"
                    213: print  "from time to time!"
                    214: print  ""
                    215: pause -1 "Done with fitting demo  (-> return)"

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>