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