PEARL Procedures rev-distro-1.6.0-0-gcf1399e-dirty Igor procedures for the analysis of PEARL data

- Main Page
- Related Pages
- Packages
- Namespaces
- Files

- File List
- Globals

pearl-fitfuncs.ipf
Go to the documentation of this file.
1 #pragma rtGlobals=3// Use modern global access method and strict wave access.
2 #pragma IgorVersion = 6.2
3 #pragma ModuleName = PearlFitFuncs
4 #pragma version = 1.01
5 #include "mm-physconst", version >= 1.05
6
7 // various fit functions for photoelectron spectroscopy
8
9 // $Id$
10 // author: matthias.muntwiler@psi.ch
11 // Copyright (c) 2013-14 Paul Scherrer Institut
12
13 // Licensed under the Apache License, Version 2.0 (the "License");
14 // you may not use this file except in compliance with the License.
15 // You may obtain a copy of the License at
16 // http://www.apache.org/licenses/LICENSE-2.0
17
18 //------------------------------------------------------------------------------
19 // Doniach-Sunjic fit functions
20 //------------------------------------------------------------------------------
21
22 threadsafe variable DoniachSunjic(variable x, variable amp, variable pos, variable sing, variable fwhm){
23 // Doniach-Sunjic line shape
24 // [S. Doniach, M. Sunjic, J. Phys. C 3 (1970) 285]
25 variable x// independent variable
26 variable amp// amplitude
27 variable pos// position
28 variable sing// singularity index (0 <= sing < 1)
29 variable fwhm// width
30
31 variable nom, denom
32 nom = cos(pi * sing / 2 + (1 - sing) * atan((x - pos) / fwhm * 2))
33 denom = ((x - pos)^2 + fwhm^2 / 4)^((1 - sing) / 2)
34
35 return amp * nom / denom * fwhm / 2
36 };
37
38 threadsafe variable ds1_bg(wave w, variable x){
39 // Doniach-Sunjic fit function
40 // 0 <= sing < 1
41 wave w// coefficients - see below
42 variable x// independent variable
43
44 //CurveFitDialog/ These comments were created by the Curve Fitting dialog. Altering them will
45 //CurveFitDialog/ make the function less convenient to work with in the Curve Fitting dialog.
46 //CurveFitDialog/ Equation:
47 //CurveFitDialog/ f(x) = DoniachSunjic(x, amp, pos, sing, fwhm) + bg
48 //CurveFitDialog/ End of Equation
49 //CurveFitDialog/ Independent Variables 1
50 //CurveFitDialog/ x
51 //CurveFitDialog/ Coefficients 5
52 //CurveFitDialog/ w[0] = bg
53 //CurveFitDialog/ w[1] = amp
54 //CurveFitDialog/ w[2] = pos
55 //CurveFitDialog/ w[3] = sing
56 //CurveFitDialog/ w[4] = FWHM
57
58 return DoniachSunjic(x, w[1], w[2], w[3], w[4]) + w[0]
59 };
60
61 threadsafe variable ds2_bg(wave w, variable x){
62 Wave w
63 Variable x
64
65 //CurveFitDialog/ These comments were created by the Curve Fitting dialog. Altering them will
66 //CurveFitDialog/ make the function less convenient to work with in the Curve Fitting dialog.
67 //CurveFitDialog/ Equation:
68 //CurveFitDialog/ f(x) = w_0+( w_1*cos(pi*w_3/2+(1-w_3)*atan((x-w_2)/w_4)))/(((x-w_2)^2+w_4^2)^((1-w_3)/2)) +(w_5*cos(pi*w_7/2+(1-w_7)*atan((x-(w_6))/w_8)))/(((x-w_6)^2+w_8^2)^((1-w_7)/2))
69 //CurveFitDialog/ End of Equation
70 //CurveFitDialog/ Independent Variables 1
71 //CurveFitDialog/ x
72 //CurveFitDialog/ Coefficients 9
73 //CurveFitDialog/ w[0] = bg
74 //CurveFitDialog/ w[1] = amp1
75 //CurveFitDialog/ w[2] = pos1
76 //CurveFitDialog/ w[3] = sing1
77 //CurveFitDialog/ w[4] = wid1
78 //CurveFitDialog/ w[5] = amp2
79 //CurveFitDialog/ w[6] = pos2
80 //CurveFitDialog/ w[7] = sing2
81 //CurveFitDialog/ w[8] = wid2
82
83 variable ds1 = DoniachSunjic(x, w[1], w[2], w[3], w[4])
84 variable ds2 = DoniachSunjic(x, w[5], w[6], w[7], w[8])
85
86 return w[0] + ds1 + ds2
87 };
88
89 variable ds4_bg(wave w, variable x){
90 Wave w
91 Variable x
92
93 //CurveFitDialog/ These comments were created by the Curve Fitting dialog. Altering them will
94 //CurveFitDialog/ make the function less convenient to work with in the Curve Fitting dialog.
95 //CurveFitDialog/ Equation:
96 //CurveFitDialog/ f(x) = w_0+( w_1*cos(pi*w_3/2+(1-w_3)*atan((x-w_2)/w_4)))/(((x-w_2)^2+w_4^2)^((1-w_3)/2)) +(w_5*cos(pi*w_7/2+(1-w_7)*atan((x-(w_6))/w_8)))/(((x-w_6)^2+w_8^2)^((1-w_7)/2)) +( w_9*cos(pi*w_11/2+(1-w_11)*atan((x-w_10)/w_12)))/(((x-w_10)^2+w_12^2)^((1-w_11)/2)) +( w_13*cos(pi*w_15/2+(1-w_15)*atan((x-w_14)/w_16)))/(((x-w_14)^2+w_16^2)^((1-w_15)/2))
97 //CurveFitDialog/ End of Equation
98 //CurveFitDialog/ Independent Variables 1
99 //CurveFitDialog/ x
100 //CurveFitDialog/ Coefficients 17
101 //CurveFitDialog/ w[0] = w_0
102 //CurveFitDialog/ w[1] = w_11
103 //CurveFitDialog/ w[2] = w_12
104 //CurveFitDialog/ w[3] = w_13
105 //CurveFitDialog/ w[4] = w_14
106 //CurveFitDialog/ w[5] = w_21
107 //CurveFitDialog/ w[6] = w_22
108 //CurveFitDialog/ w[7] = w_23
109 //CurveFitDialog/ w[8] = w_24
110 //CurveFitDialog/ w[9] = w_31
111 //CurveFitDialog/ w[10] = w_32
112 //CurveFitDialog/ w[11] = w_33
113 //CurveFitDialog/ w[12] = w_34
114 //CurveFitDialog/ w[13] = w_41
115 //CurveFitDialog/ w[14] = w_42
116 //CurveFitDialog/ w[15] = w_43
117 //CurveFitDialog/ w[16] = w_44
118 Variable ds1, ds2, ds3, ds4
119 ds1=( w[1]*cos(pi*w[3]/2+(1-w[3])*atan((x-w[2])/w[4])))/(((x-w[2])^2+w[4]^2)^((1-w[3])/2))
120 ds2=( w[5]*cos(pi*w[7]/2+(1-w[7])*atan((x-w[6])/w[8])))/(((x-w[6])^2+w[8]^2)^((1-w[7])/2))
121 ds3=( w[9]*cos(pi*w[11]/2+(1-w[11])*atan((x-w[10])/w[12])))/(((x-w[10])^2+w[12]^2)^((1-w[11])/2))
122 ds4=( w[13]*cos(pi*w[15]/2+(1-w[15])*atan((x-w[14])/w[16])))/(((x-w[14])^2+w[16]^2)^((1-w[15])/2))
123
124
125 return w[0]+ds1+ds2+ds3+ds4
126
127
128 };
129
130 variable ds6_bg(wave w, variable x){
131 Wave w
132 Variable x
133
134 //CurveFitDialog/ These comments were created by the Curve Fitting dialog. Altering them will
135 //CurveFitDialog/ make the function less convenient to work with in the Curve Fitting dialog.
136 //CurveFitDialog/ Equation:
137 //CurveFitDialog/
138 //CurveFitDialog/ Variable g, ds1, ds2, ds3, ds4, ds5, ds6
139 //CurveFitDialog/ ds1=( w_11*cos(pi*w_13/2+(1-w_13)*atan((x-w_12)/w_14)))/(((x-w_12)^2+w_14^2)^((1-w_13)/2))
140 //CurveFitDialog/ ds2=( w_21*cos(pi*w_23/2+(1-w_23)*atan((x-w_22)/w_24)))/(((x-w_22)^2+w_24^2)^((1-w_23)/2))
141 //CurveFitDialog/ ds3=( w_31*cos(pi*w_33/2+(1-w_33)*atan((x-w_32)/w_34)))/(((x-w_32)^2+w_34^2)^((1-w_33)/2))
142 //CurveFitDialog/ ds4=( w_41*cos(pi*w_43/2+(1-w_43)*atan((x-w_42)/w_44)))/(((x-w_42)^2+w_44^2)^((1-w_43)/2))
143 //CurveFitDialog/ ds5=( w_51*cos(pi*w_53/2+(1-w_53)*atan((x-w_52)/w_54)))/(((x-w_52)^2+w_54^2)^((1-w_53)/2))
144 //CurveFitDialog/ ds6=( w_61*cos(pi*w_63/2+(1-w_63)*atan((x-w_62)/w_64)))/(((x-w_62)^2+w_64^2)^((1-w_63)/2))
145 //CurveFitDialog/
146 //CurveFitDialog/ f(x) =w_0+ds1+ds2+ds3+ds4+ds5+ds6
147 //CurveFitDialog/
148 //CurveFitDialog/ End of Equation
149 //CurveFitDialog/ Independent Variables 1
150 //CurveFitDialog/ x
151 //CurveFitDialog/ Coefficients 25
152 //CurveFitDialog/ w[0] = w_0
153 //CurveFitDialog/ w[1] = w_11
154 //CurveFitDialog/ w[2] = w_12
155 //CurveFitDialog/ w[3] = w_13
156 //CurveFitDialog/ w[4] = w_14
157 //CurveFitDialog/ w[5] = w_21
158 //CurveFitDialog/ w[6] = w_22
159 //CurveFitDialog/ w[7] = w_23
160 //CurveFitDialog/ w[8] = w_24
161 //CurveFitDialog/ w[9] = w_31
162 //CurveFitDialog/ w[10] = w_32
163 //CurveFitDialog/ w[11] = w_33
164 //CurveFitDialog/ w[12] = w_34
165 //CurveFitDialog/ w[13] = w_41
166 //CurveFitDialog/ w[14] = w_42
167 //CurveFitDialog/ w[15] = w_43
168 //CurveFitDialog/ w[16] = w_44
169 //CurveFitDialog/ w[17] = w_51
170 //CurveFitDialog/ w[18] = w_52
171 //CurveFitDialog/ w[19] = w_53
172 //CurveFitDialog/ w[20] = w_54
173 //CurveFitDialog/ w[21] = w_61
174 //CurveFitDialog/ w[22] = w_62
175 //CurveFitDialog/ w[23] = w_63
176 //CurveFitDialog/ w[24] = w_64
177
178
179 Variable ds1, ds2, ds3, ds4, ds5, ds6
180 ds1=( w[1]*cos(pi*w[3]/2+(1-w[3])*atan((x-w[2])/w[4])))/(((x-w[2])^2+w[4]^2)^((1-w[3])/2))
181 ds2=( w[5]*cos(pi*w[7]/2+(1-w[7])*atan((x-w[6])/w[8])))/(((x-w[6])^2+w[8]^2)^((1-w[7])/2))
182 ds3=( w[9]*cos(pi*w[11]/2+(1-w[11])*atan((x-w[10])/w[12])))/(((x-w[10])^2+w[12]^2)^((1-w[11])/2))
183 ds4=( w[13]*cos(pi*w[15]/2+(1-w[15])*atan((x-w[14])/w[16])))/(((x-w[14])^2+w[16]^2)^((1-w[15])/2))
184 ds5=( w[17]*cos(pi*w[19]/2+(1-w[19])*atan((x-w[18])/w[20])))/(((x-w[18])^2+w[20]^2)^((1-w[19])/2))
185 ds6=( w[21]*cos(pi*w[23]/2+(1-w[23])*atan((x-w[22])/w[24])))/(((x-w[22])^2+w[24]^2)^((1-w[23])/2))
186
187 return w[0]+ds1+ds2+ds3+ds4+ds5+ds6
188
189 };
190
191 struct DoniachSunjicStruct{;
192 // data structure for DoniachSunjicBroadS structural function fit
193
194 // waves populated by the FuncFit operation
195 wave pw;
196 wave yw;
197 wave xw;
198
199 // convolution parameters to be set upon creation of the structure
200 variable precision;
201 variable oversampling;
202
203 // auxiliary fields used internally by DoniachSunjicBroadS
204 // do not touch these
205 wave xdw;
206 wave model;
207 wave broadening;
208 wave convolution;
209 };
210
211 //------------------------------------------------------------------------------
212 threadsafe variable DoniachSunjicBroadS(DoniachSunjicStruct* s){
213 //------------------------------------------------------------------------------
214 // Two-peak (bulk + surface) Doniach-Sunjic line shape with Gaussian broadening (convolution).
215 // Hold parameter 5 at 0 to fit just one peak.
216
217 // Structural fit function for efficient fitting in procedures.
218 // Calculating the convolution requires auxiliary waves and additional, non-fitting parameters.
219 // To eliminate the time-consuming overhead of creating and killing the auxiliary waves,
220 // these waves are held in the fitting structure.
221
222 // Caution: The function on its own is thread-safe.
223 // However, since FuncFit uses the same structure in all threads, the fitting cannot run in parallel.
224 // Set /NTHR=1.
225
226 // See also Fit_DoniachSunjicBroad (example), DoniachSunjicBroad (conventional fit function)
227 Struct DoniachSunjicStruct &s
228
229 // pw[0] = bulk amplitude
230 // pw[1] = bulk position
231 // pw[2] = Lorentzian FWHM
232 // pw[3] = Donjach-Sunjic singularity index (0..1)
233 // pw[4] = surface shift
234 // pw[5] = surface/bulk ratio
235 // pw[6] = Gaussian FWHM
236 // pw[7] = constant background
237 // pw[8] = linear background
238
239 wave xw = s.xw
240 wave yw = s.yw
241 wave pw = s.pw
242
243 variable precision = s.precision
244 variable oversampling = s.oversampling
245
246 if (WaveExists(s.xdw))
247 wave xdw = s.xdw
248 wave model = s.model
249 wave broadening = s.broadening
250 wave convolution = s.convolution
251 else
252 make /n=0 /free xdw, model, broadening, convolution
253 redimension /d xdw, model, broadening, convolution
254 wave fs.xdw = xdw
255 wave fs.model = model
256 wave fs.broadening = broadening
257 wave fs.convolution = convolution
258 endif
259
260 // calculate wave spacing based on minimum spacing of desired x points
261 differentiate /p xw /d=xdw
262 xdw = abs(xdw)
263 variable xd = wavemin(xdw) / oversampling
264
265 // calculate broadening wave size based on width and precision variable
266 variable x0b = pw[6] * precision
267 variable nb = 2 * floor(x0b / xd) + 1
268
269 // calculate model size based on desired range for yw
270 variable x0m = max(abs(wavemax(xw) - pw[1]), abs(wavemin(xw) - pw[1])) + x0b
271 variable nm = 2 * floor(x0m / xd) + 1
272 nb = min(nb, nm * 10)// limit wave size to avoid runtime errors for unphysically large parameter
273
274 // create and calculate initial waves, normalize exponential
275 redimension /n=(nb) broadening
276 redimension /n=(nm) model
277 setscale/i x -x0b, x0b, "", broadening
278 setscale/i x -x0m, x0m, "", model
279
280 broadening = exp( - (x / pw[6])^2 * 4 * ln(2))
281 variable nrm = area(broadening)
282 broadening /= nrm
283 model = DoniachSunjic(x, 1, 0, pw[3], pw[2])// bulk
284 model += DoniachSunjic(x, pw[5], pw[4], pw[3], pw[2])// surface
285
286 // calculate the convolution
287 Convolve /a broadening, model
288 variable scale = pw[0] / wavemax(model)
289 model *= scale
290
291 // prepare output
292 nm = numpnts(model)
293 x0m = xd * (nm - 1) / 2
294 setscale/i x -x0m, x0m, "", model
295
296 yw = model(xw[p] - pw[1]) + pw[7] + pw[8] * xw[p]
297 yw = numtype(yw) ? 0 : yw
298 };
299
300 //------------------------------------------------------------------------------
301 variable DoniachSunjicBroad(wave pw, wave yw, wave xw){
302 //------------------------------------------------------------------------------
303 // Two-peak (bulk + surface) Doniach-Sunjic line shape with Gaussian broadening (convolution).
304 // Hold parameter 5 at 0 to fit just one peak.
305 // Conventional fit function for use with the curve-fitting dialog.
306 // Compared to DoniachSunjicBroadS this function incurs extra overhead
307 // because auxiliary waves are created and killed between function calls.
308 // See also DoniachSunjicBroadS (optimized structural fit function)
309 Wave pw
310 Wave yw
311 Wave xw
312
313 // pw[0] = bulk amplitude
314 // pw[1] = bulk position
315 // pw[2] = Lorentzian FWHM
316 // pw[3] = Donjach-Sunjic singularity index (0..1)
317 // pw[4] = surface shift
318 // pw[5] = surface/bulk ratio
319 // pw[6] = Gaussian FWHM
320 // pw[7] = constant background
321 // pw[8] = linear background
322
323 // set up data structure
324 struct DoniachSunjicStruct fs
325 fs.precision = 5
326 fs.oversampling = 4
327
328 wave fs.pw = pw
329 wave fs.xw = xw
330 wave fs.yw = yw
331
332 // create temporary calculation waves in a global folder
333 dfref df = root:packages:pearl_fitfuncs:doniach_sunjic
334 if (DataFolderRefStatus(df) == 0)
335 newdatafolder root:packages:pearl_fitfuncs:doniach_sunjic
336 dfref df = root:packages:pearl_fitfuncs:doniach_sunjic
337 endif
338
339 wave /z /sdfr=df fs.xdw = xdw
340 wave /z /sdfr=df fs.model = model
341 wave /z /sdfr=df fs.broadening = broadening
342 wave /z /sdfr=df fs.convolution = convolution
343
344 if (WaveExists(fs.xdw) == 0)
345 dfref savedf = GetDataFolderDFR()
346 setdatafolder df
347 make /n=0 /d xdw, model, broadening, convolution
348 wave fs.xdw = xdw
349 wave fs.model = model
350 wave fs.broadening = broadening
351 wave fs.convolution = convolution
352 setdatafolder savedf
353 endif
354
355 // calculate
356 DoniachSunjicBroadS(fs)
357
358 yw = fs.yw
359 };
360
361 //------------------------------------------------------------------------------
362 variable Calc_DoniachSunjicBroad(wave pw, wave yw){
363 //------------------------------------------------------------------------------
364 // Calculate the DoniachSunjicBroadS line shape
365 Wave pw// coefficient wave
366 Wave yw// output wave, correct x-scaling required on input
367
368 struct DoniachSunjicStruct fs
369 fs.precision = 5
370 fs.oversampling = 4
371
372 duplicate /free pw, fs.pw
373 duplicate /free yw, fs.xw
374 fs.xw = x
375 duplicate /free yw, fs.yw
376
377 DoniachSunjicBroadS(fs)
378
379 yw = fs.yw
380 };
381
382 //------------------------------------------------------------------------------
383 variable Fit_DoniachSunjicBroad(wave pw, wave yw, wave xw, wave ww){
384 //------------------------------------------------------------------------------
385 // Fit the DoniachSunjicBroadS line shape.
386 // The function applies constraints which assume that the energy scale is in eV.
387 // Returns chi^2.
388 wave pw// coefficient wave- pre-load it with initial guess
389 wave yw
390 wave /z xw
391 wave /z ww// weights (standard deviation)
392
393 struct DoniachSunjicStruct fs
394 fs.precision = 5
395 fs.oversampling = 4
396
397 duplicate /free pw, fs.pw
398 if (WaveExists(xw))
399 duplicate /free xw, fs.xw
400 else
401 duplicate /free yw, fs.xw
402 fs.xw = x
403 endif
404 duplicate /free yw, fs.yw
405
406 variable v_chisq = nan
407 variable V_FitMaxIters = 100
408 make /n=1 /t /free constraints = {"K0 >= 0", "K2 > 0", "K2 < 10", "K3 >= 0", "K3 < 1", "K4 >= -10", "K4 <= 10", "K5 >= 0", "K5 <= 1", "K6 >= 0", "K6 < 10"}
409 // note: only single thread allowed
410 FuncFit /NTHR=1 DoniachSunjicBroadS, pw, yw /X=xw /D /STRC=fs /C=constraints /NWOK /I=1 /W=ww
411
412 return v_chisq
413 };
414
415 //------------------------------------------------------------------------------
416 // peak-specific fit functions
417 //------------------------------------------------------------------------------
418
419 variable Au4f(wave w, variable x){
420 // fit function for a nitrogen 1s-pi* absorption spectrum
421 // modelled as multiple Voigt shapes on a constant background
422 // similar to the Igor VoigtFit function
423 // but with a constant gaussian width (instrumental broadening) for all peaks
424 // gaussian and lorentzian widths are specified as FWHM
425 wave w// parameters
426 // w[0] constant background
427 // w[1] linear background
428 // w[2] global gaussian FWHM
429 // w[3 + 0 + (n-1) * 3] peak n area
430 // w[3 + 1 + (n-1) * 3] peak n position
431 // w[3 + 2 + (n-1) * 3] peak n lorentzian FWHM
432 // length of wave defines number of peaks
433
434 // for compatibility with older code the linear background term can be omitted.
435 // if the number of parameters divides by 3, the linear background term is added,
436 // otherwise only the constant background.
437 variable x
438
439 variable np = 15
440 variable ip, ip0
441
442 variable bg = w[0]
443 variable v = bg
444 if (mod(np, 3) == 0)
445 v += w[1] * x
446 ip0 = 3
447 else
448 ip0 = 2
449 endif
450
451 variable vc1, vc2, vc3, vc4
452 vc2 = 2 * sqrt(ln(2)) / w[ip0-1]
453 for (ip = ip0; ip < np; ip += 3)
454 vc1 = w[ip] / sqrt(pi) * vc2
455 vc3 = w[ip+1]
456 vc4 = vc2 * w[ip+2] / 2
457 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
458 endfor
459
460 return v
461
462 };
463
464 variable Au4f_2p2(wave w, variable x){
465 // Au 4f 5/2 and 7/2 2-component Voigt fit with a common gaussian width
466 // gaussian and lorentzian widths are specified as FWHM
467 wave w// parameters
468 // w[0] constant background
469 // w[1] linear background
470 // w[2] global gaussian FWHM
471 // w[3] 5/2 bulk area
472 // w[4] 5/2 bulk position
473 // w[5] 5/2 lorentzian FWHM
474 // w[6] 7/2 bulk area
475 // w[7] 7/2 bulk position
476 // w[8] 7/2 lorentzian FWHM
477 // w[9] surface/bulk area ratio
478 // w[10] surface core level shift
479 variable x
480
481 variable bg = w[0] + w[1] * x
482 variable v = bg
483
484 variable vc1// amplitude
485 variable vc2// width
486 variable vc3// position
487 variable vc4// shape
488 vc2 = 2 * sqrt(ln(2)) / w[2]
489
490 // 5/2 bulk
491 vc1 = w[3] / sqrt(pi) * vc2
492 vc3 = w[4]
493 vc4 = vc2 * w[5] / 2
494 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
495
496 // 5/2 surface
497 vc1 = w[3] / sqrt(pi) * vc2 * w[9]
498 vc3 = w[4] + w[10]
499 vc4 = vc2 * w[5] / 2
500 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
501
502 // 7/2 bulk
503 vc1 = w[6] / sqrt(pi) * vc2
504 vc3 = w[7]
505 vc4 = vc2 * w[8] / 2
506 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
507
508 // 7/2 surface
509 vc1 = w[6] / sqrt(pi) * vc2 * w[9]
510 vc3 = w[7] + w[10]
511 vc4 = vc2 * w[8] / 2
512 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
513
514 return v
515
516 };
517
518 variable ShowComponents_Au4f_2p2(wave coef_wave, wave fit_wave){
519 wave coef_wave
520 wave fit_wave
521
522 duplicate /free coef_wave, coef1, coef2
523 coef1[9] = 0
524 coef2[3] *= coef_wave[9]
525 coef2[4] += coef_wave[10]
526 coef2[6] *= coef_wave[9]
527 coef2[7] += coef_wave[10]
528 coef2[9] = 0
529
530 string s_fit_wave = NameOfWave(fit_wave)
531 string s_fit_p1 = s_fit_wave + "_p1"
532 string s_fit_p2 = s_fit_wave + "_p2"
533 duplicate /o fit_wave, $(s_fit_p1) /wave=fit_p1
534 duplicate /o fit_wave, $(s_fit_p2) /wave=fit_p2
535
536 fit_p1 = Au4f_2p2(coef1, x)
537 fit_p2 = Au4f_2p2(coef2, x)
538
539 string traces = TraceNameList("", ";", 1)
540 if ((WhichListItem(s_fit_wave, traces, ";") >= 0) && (WhichListItem(s_fit_p1, traces, ";") < 0))
541 appendtograph fit_p1, fit_p2
542 ModifyGraph lstyle($s_fit_p1)=2
543 ModifyGraph lstyle($s_fit_p2)=2
544 ModifyGraph rgb($s_fit_p1)=(0,0,65280)
545 ModifyGraph rgb($s_fit_p2)=(0,0,65280)
546 endif
547 };
548
549 variable Au4f_2p3(wave w, variable x){
550 // Au 4f 5/2 and 7/2 3-component Voigt fit with a common gaussian width
551 // gaussian and lorentzian widths are specified as FWHM
552 wave w// parameters
553 // w[0] constant background
554 // w[1] linear background
555 // w[2] global gaussian FWHM
556 // w[3] 5/2 bulk area
557 // w[4] 5/2 bulk position
558 // w[5] 5/2 lorentzian FWHM
559 // w[6] 7/2 bulk area
560 // w[7] 7/2 bulk position
561 // w[8] 7/2 lorentzian FWHM
562 // w[9] surface/bulk area ratio
563 // w[10] surface core level shift
564 // w[11] 2nd layer/bulk area ratio
565 // w[12] 2nd layer core level shift
566 variable x
567
568 variable bg = w[0] + w[1] * x
569 variable v = bg
570
571 variable vc1// amplitude
572 variable vc2// width
573 variable vc3// position
574 variable vc4// shape
575 vc2 = 2 * sqrt(ln(2)) / w[2]
576
577 // 5/2 bulk
578 vc1 = w[3] / sqrt(pi) * vc2
579 vc3 = w[4]
580 vc4 = vc2 * w[5] / 2
581 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
582
583 // 5/2 surface
584 vc1 = w[3] / sqrt(pi) * vc2 * w[9]
585 vc3 = w[4] + w[10]
586 vc4 = vc2 * w[5] / 2
587 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
588
589 // 5/2 2nd layer
590 vc1 = w[3] / sqrt(pi) * vc2 * w[11]
591 vc3 = w[4] + w[12]
592 vc4 = vc2 * w[5] / 2
593 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
594
595 // 7/2 bulk
596 vc1 = w[6] / sqrt(pi) * vc2
597 vc3 = w[7]
598 vc4 = vc2 * w[8] / 2
599 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
600
601 // 7/2 surface
602 vc1 = w[6] / sqrt(pi) * vc2 * w[9]
603 vc3 = w[7] + w[10]
604 vc4 = vc2 * w[8] / 2
605 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
606
607 // 7/2 2nd layer
608 vc1 = w[6] / sqrt(pi) * vc2 * w[11]
609 vc3 = w[7] + w[12]
610 vc4 = vc2 * w[8] / 2
611 v += vc1 * VoigtFunc(vc2 * (x - vc3), vc4)
612
613 return v
614
615 };
616
617 variable ShowComponents_Au4f_2p3(wave coef_wave, wave fit_wave){
618 wave coef_wave
619 wave fit_wave
620
621 duplicate /free coef_wave, coef1, coef2, coef3
622 coef1[9] = 0
623 coef1[11] = 0
624
625 coef2[3] *= coef_wave[9]
626 coef2[4] += coef_wave[10]
627 coef2[6] *= coef_wave[9]
628 coef2[7] += coef_wave[10]
629 coef2[9] = 0
630 coef2[11] = 0
631
632 coef3[3] *= coef_wave[11]
633 coef3[4] += coef_wave[12]
634 coef3[6] *= coef_wave[11]
635 coef3[7] += coef_wave[12]
636 coef3[9] = 0
637 coef3[11] = 0
638
639 string s_fit_wave = NameOfWave(fit_wave)
640 string s_fit_p1 = s_fit_wave + "_p1"
641 string s_fit_p2 = s_fit_wave + "_p2"
642 string s_fit_p3 = s_fit_wave + "_p3"
643 duplicate /o fit_wave, $(s_fit_p1) /wave=fit_p1
644 duplicate /o fit_wave, $(s_fit_p2) /wave=fit_p2
645 duplicate /o fit_wave, $(s_fit_p3) /wave=fit_p3
646
647 fit_p1 = Au4f_2p2(coef1, x)
648 fit_p2 = Au4f_2p2(coef2, x)
649 fit_p3 = Au4f_2p2(coef3, x)
650
651 string traces = TraceNameList("", ";", 1)
652 if ((WhichListItem(s_fit_wave, traces, ";") >= 0) && (WhichListItem(s_fit_p1, traces, ";") < 0))
653 appendtograph fit_p1, fit_p2, fit_p3
654 ModifyGraph lstyle($s_fit_p1)=2
655 ModifyGraph lstyle($s_fit_p2)=2
656 ModifyGraph lstyle($s_fit_p3)=2
657 ModifyGraph rgb($s_fit_p1)=(0,0,65280)
658 ModifyGraph rgb($s_fit_p2)=(0,0,65280)
659 ModifyGraph rgb($s_fit_p3)=(0,0,65280)
660 endif
661 };
662
672 variable FermiGaussConv(wave pw, wave yw, wave xw){
673 WAVE pw, yw, xw
674
675 // half width of temporary gaussian wave is pw[5] multiplied by this factor (may be fractional)
676 variable precision_g = 5
677 variable oversampling = 4
678
679 // calculate wave spacing based on minimum spacing of desired x points
680 duplicate /free xw, xdw
681 differentiate /p xw /d=xdw
682 xdw = abs(xdw)
683 variable xd = wavemin(xdw) / oversampling
684
685 // calculate gausswave size based on pw[5] and precision variable
686 variable x0g = abs(pw[5]) * precision_g
687 variable ng = 2 * floor(x0g / xd) + 1
688
689 // calculate fermiwave size based on desired range for yw
690 variable emax = wavemax(xw)
691 variable emin = wavemin(xw)
692 variable x0f = max(abs(emax - pw[3]), abs(emin - pw[3])) + x0g
693 variable ne = 2 * floor(x0f / xd) + 1
694
695 // create and calculate initial waves, normalize exponential
696 make /d /n=(ng) /free gausswave
697 make /d /n=(ne) /free fermiwave
698 setscale/i x -x0g, x0g, "", gausswave
699 setscale/i x -x0f, x0f, "", fermiwave
700
701 gausswave = exp( - (x / pw[5] )^2 )
702 fermiwave = 1 / (exp( x / (kBoltzmann * pw[4])) + 1.0 )
703
704 // calculate the convolution
705 duplicate /free fermiwave, resultwave
706 Convolve /a gausswave, resultwave
707 variable rmax = wavemax(resultwave)
708 resultwave /= rmax
709
710 // prepare output
711 ng = numpnts(resultwave)
712 x0g = xd * (ng - 1) / 2
713 setscale/i x -x0g, x0g, "", resultwave
714
715 yw = pw[2] * resultwave(xw[p] - pw[3]) + pw[0] + pw[1] * xw[p]
716 };
717
DoniachSunjicBroad
variable DoniachSunjicBroad(wave pw, wave yw, wave xw)
Definition: pearl-fitfuncs.ipf:301
DoniachSunjicStruct::convolution
wave convolution
Definition: pearl-fitfuncs.ipf:208
ds4_bg
variable ds4_bg(wave w, variable x)
Definition: pearl-fitfuncs.ipf:89
DoniachSunjicStruct
Definition: pearl-fitfuncs.ipf:191
FermiGaussConv
variable FermiGaussConv(wave pw, wave yw, wave xw)
convolution of Fermi-Dirac distribution and a Gaussian.
Definition: pearl-fitfuncs.ipf:672
DoniachSunjicBroadS
threadsafe variable DoniachSunjicBroadS(DoniachSunjicStruct *s)
Definition: pearl-fitfuncs.ipf:212
DoniachSunjicStruct::xw
wave xw
Definition: pearl-fitfuncs.ipf:197
Au4f
variable Au4f(wave w, variable x)
Definition: pearl-fitfuncs.ipf:419
Au4f_2p3
variable Au4f_2p3(wave w, variable x)
Definition: pearl-fitfuncs.ipf:549
Au4f_2p2
variable Au4f_2p2(wave w, variable x)
Definition: pearl-fitfuncs.ipf:464
ShowComponents_Au4f_2p2
variable ShowComponents_Au4f_2p2(wave coef_wave, wave fit_wave)
Definition: pearl-fitfuncs.ipf:518
Calc_DoniachSunjicBroad
variable Calc_DoniachSunjicBroad(wave pw, wave yw)
Definition: pearl-fitfuncs.ipf:362
ShowComponents_Au4f_2p3
variable ShowComponents_Au4f_2p3(wave coef_wave, wave fit_wave)
Definition: pearl-fitfuncs.ipf:617
Fit_DoniachSunjicBroad
variable Fit_DoniachSunjicBroad(wave pw, wave yw, wave xw, wave ww)
Definition: pearl-fitfuncs.ipf:383
DoniachSunjicStruct::yw
wave yw
Definition: pearl-fitfuncs.ipf:196
DoniachSunjicStruct::oversampling
variable oversampling
Definition: pearl-fitfuncs.ipf:201
DoniachSunjicStruct::precision
variable precision
Definition: pearl-fitfuncs.ipf:200
ds2_bg
threadsafe variable ds2_bg(wave w, variable x)
Definition: pearl-fitfuncs.ipf:61
DoniachSunjicStruct::pw
wave pw
Definition: pearl-fitfuncs.ipf:191
version
version
Definition: pearl-fitfuncs.ipf:5
DoniachSunjicStruct::xdw
wave xdw
Definition: pearl-fitfuncs.ipf:205
ds6_bg
variable ds6_bg(wave w, variable x)
Definition: pearl-fitfuncs.ipf:130
DoniachSunjicStruct::broadening
wave broadening
Definition: pearl-fitfuncs.ipf:207
ds1_bg
threadsafe variable ds1_bg(wave w, variable x)
Definition: pearl-fitfuncs.ipf:38
DoniachSunjicStruct::model
wave model
Definition: pearl-fitfuncs.ipf:206