1 | #pragma rtGlobals=1 // Use modern global access method. |
---|
2 | #pragma IgorVersion=6.0 |
---|
3 | |
---|
4 | //////////////////////////////////////////////// |
---|
5 | // |
---|
6 | // this function is for the form factor of a with some |
---|
7 | // number of shells around a central core (currently 1-2-3-4) |
---|
8 | // |
---|
9 | // monodisperse and polydisperse (and smeared) versions are included |
---|
10 | // - for the polydisperse models, only polydispersity of the core is taken |
---|
11 | // into account, and dPolyOne numerically. for a Schulz distribution, this |
---|
12 | // should be possible to do analytically, whith a great savings in computation |
---|
13 | // time. |
---|
14 | // |
---|
15 | // It may also be useful to think of scenarios where the layers as well are |
---|
16 | // polydisperse - to break up the very regular spacing of the layers, which |
---|
17 | // is not a very natural structure. |
---|
18 | // |
---|
19 | // 03 MAR 04 SRK |
---|
20 | // 07 AUG 08 AJJ - redone for new version of software |
---|
21 | // 08 AUG 08 AJJ - bug fixed three and four shell models. |
---|
22 | //////////////////////////////////////////////// |
---|
23 | |
---|
24 | #include "Core_and_NShells_v40" |
---|
25 | |
---|
26 | //this macro sets up all the necessary parameters and waves that are |
---|
27 | //needed to calculate the model function. |
---|
28 | // |
---|
29 | Proc PlotPolyOneShell(num,qmin,qmax) |
---|
30 | Variable num=200, qmin=0.001, qmax=0.7 |
---|
31 | Prompt num "Enter number of data points for model: " |
---|
32 | Prompt qmin "Enter minimum q-value (^-1) for model: " |
---|
33 | Prompt qmax "Enter maximum q-value (^-1) for model: " |
---|
34 | // |
---|
35 | Make/O/D/n=(num) xwave_PolyOneShell, ywave_PolyOneShell |
---|
36 | xwave_PolyOneShell = alog(log(qmin) + x*((log(qmax)-log(qmin))/num)) |
---|
37 | Make/O/D coef_PolyOneShell = {1.,60,0.1,6.4e-6,10,1e-6,6.4e-6,0.001} |
---|
38 | make/o/t parameters_PolyOneShell = {"scale","core radius (A)","Core Polydispersity(0,1)","Core SLD (A-2)","Shell thickness (A)","Shell SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"} |
---|
39 | Edit parameters_PolyOneShell, coef_PolyOneShell |
---|
40 | |
---|
41 | Variable/G root:g_PolyOneShell |
---|
42 | g_PolyOneShell := PolyOneShell(coef_PolyOneShell, ywave_PolyOneShell, xwave_PolyOneShell) |
---|
43 | Display ywave_PolyOneShell vs xwave_PolyOneShell |
---|
44 | ModifyGraph marker=29, msize=2, mode=4 |
---|
45 | ModifyGraph log=1,grid=1,mirror=2 |
---|
46 | Label bottom "q (\\S-1\\M) " |
---|
47 | Label left "I(q) (cm\\S-1\\M)" |
---|
48 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
49 | |
---|
50 | AddModelToStrings("PolyOneShell","coef_PolyOneShell","PolyOneShell") |
---|
51 | // |
---|
52 | End |
---|
53 | |
---|
54 | Proc PlotPolyTwoShell(num,qmin,qmax) |
---|
55 | Variable num=200, qmin=0.001, qmax=0.7 |
---|
56 | Prompt num "Enter number of data points for model: " |
---|
57 | Prompt qmin "Enter minimum q-value (^-1) for model: " |
---|
58 | Prompt qmax "Enter maximum q-value (^-1) for model: " |
---|
59 | // |
---|
60 | Make/O/D/n=(num) xwave_PolyTwoShell, ywave_PolyTwoShell |
---|
61 | xwave_PolyTwoShell = alog(log(qmin) + x*((log(qmax)-log(qmin))/num)) |
---|
62 | Make/O/D coef_PolyTwoShell = {1.,60,0.1,6.4e-6,10,1e-6,10,2e-6,6.4e-6,0.001} |
---|
63 | make/o/t parameters_PolyTwoShell = {"scale","core radius (A)","Core Polydispersity(0,1)","Core SLD (A-2)","Shell 1 thickness","Shell 1 SLD (A-2)","Shell 2 thickness","Shell 2 SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"} |
---|
64 | Edit parameters_PolyTwoShell, coef_PolyTwoShell |
---|
65 | |
---|
66 | Variable/G root:g_PolyTwoShell |
---|
67 | g_PolyTwoShell := PolyTwoShell(coef_PolyTwoShell, ywave_PolyTwoShell, xwave_PolyTwoShell) |
---|
68 | Display ywave_PolyTwoShell vs xwave_PolyTwoShell |
---|
69 | ModifyGraph marker=29, msize=2, mode=4 |
---|
70 | ModifyGraph log=1,grid=1,mirror=2 |
---|
71 | Label bottom "q (\\S-1\\M) " |
---|
72 | Label left "I(q) (cm\\S-1\\M)" |
---|
73 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
74 | |
---|
75 | AddModelToStrings("PolyTwoShell","coef_PolyTwoShell","PolyTwoShell") |
---|
76 | // |
---|
77 | End |
---|
78 | |
---|
79 | Proc PlotPolyThreeShell(num,qmin,qmax) |
---|
80 | Variable num=200, qmin=0.001, qmax=0.7 |
---|
81 | Prompt num "Enter number of data points for model: " |
---|
82 | Prompt qmin "Enter minimum q-value (^-1) for model: " |
---|
83 | Prompt qmax "Enter maximum q-value (^-1) for model: " |
---|
84 | // |
---|
85 | Make/O/D/n=(num) xwave_PolyThreeShell, ywave_PolyThreeShell |
---|
86 | xwave_PolyThreeShell = alog(log(qmin) + x*((log(qmax)-log(qmin))/num)) |
---|
87 | Make/O/D coef_PolyThreeShell ={1.,60,0.1,6.4e-6,10,1e-6,10,2e-6,10,3e-6,6.4e-6,0.001} |
---|
88 | make/o/t parameters_PolyThreeShell = {"scale","core radius (A)","Core Polydispersity(0,1)","Core SLD (A-2)","Shell 1 thickness","Shell 1 SLD (A-2)","Shell 2 thickness","Shell 2 SLD (A-2)","Shell 3 thickness","Shell 3 SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"} |
---|
89 | Edit parameters_PolyThreeShell, coef_PolyThreeShell |
---|
90 | |
---|
91 | Variable/G root:g_PolyThreeShell |
---|
92 | g_PolyThreeShell := PolyThreeShell(coef_PolyThreeShell, ywave_PolyThreeShell, xwave_PolyThreeShell) |
---|
93 | Display ywave_PolyThreeShell vs xwave_PolyThreeShell |
---|
94 | ModifyGraph marker=29, msize=2, mode=4 |
---|
95 | ModifyGraph log=1,grid=1,mirror=2 |
---|
96 | Label bottom "q (\\S-1\\M) " |
---|
97 | Label left "I(q) (cm\\S-1\\M)" |
---|
98 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
99 | |
---|
100 | AddModelToStrings("PolyThreeShell","coef_PolyThreeShell","PolyThreeShell") |
---|
101 | // |
---|
102 | End |
---|
103 | |
---|
104 | Proc PlotPolyFourShell(num,qmin,qmax) |
---|
105 | Variable num=200, qmin=0.001, qmax=0.7 |
---|
106 | Prompt num "Enter number of data points for model: " |
---|
107 | Prompt qmin "Enter minimum q-value (^-1) for model: " |
---|
108 | Prompt qmax "Enter maximum q-value (^-1) for model: " |
---|
109 | // |
---|
110 | Make/O/D/n=(num) xwave_PolyFourShell, ywave_PolyFourShell |
---|
111 | xwave_PolyFourShell = alog(log(qmin) + x*((log(qmax)-log(qmin))/num)) |
---|
112 | Make/O/D coef_PolyFourShell ={1.,60,0.1,6.4e-6,10,1e-6,10,2e-6,5,3e-6,10,4e-6,6.4e-6,0.001} |
---|
113 | make/o/t parameters_PolyFourShell = {"scale","core radius (A)","Core Polydispersity(0,1)","Core SLD (A-2)","Shell 1 thickness","Shell 1 SLD (A-2)","Shell 2 thickness","Shell 2 SLD (A-2)","Shell 3 thickness","Shell 3 SLD (A-2)","Shell 4 thickness","Shell 4 SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"} |
---|
114 | Edit parameters_PolyFourShell, coef_PolyFourShell |
---|
115 | |
---|
116 | Variable/G root:g_PolyFourShell |
---|
117 | g_PolyFourShell := PolyFourShell(coef_PolyFourShell, ywave_PolyFourShell, xwave_PolyFourShell) |
---|
118 | Display ywave_PolyFourShell vs xwave_PolyFourShell |
---|
119 | ModifyGraph marker=29, msize=2, mode=4 |
---|
120 | ModifyGraph log=1,grid=1,mirror=2 |
---|
121 | Label bottom "q (\\S-1\\M) " |
---|
122 | Label left "I(q) (cm\\S-1\\M)" |
---|
123 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
124 | |
---|
125 | AddModelToStrings("PolyFourShell","coef_PolyFourShell","PolyFourShell") |
---|
126 | // |
---|
127 | End |
---|
128 | |
---|
129 | |
---|
130 | |
---|
131 | |
---|
132 | // |
---|
133 | //this macro sets up all the necessary parameters and waves that are |
---|
134 | //needed to calculate the smeared model function. |
---|
135 | // |
---|
136 | //no input parameters are necessary, it MUST use the experimental q-values |
---|
137 | // from the experimental data read in from an AVE/QSIG data file |
---|
138 | //////////////////////////////////////////////////// |
---|
139 | // - sets up a dependency to a wrapper, not the actual SmearedModelFunction |
---|
140 | Proc PlotSmearedPolyOneShell(str) |
---|
141 | String str |
---|
142 | Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4) |
---|
143 | |
---|
144 | // if any of the resolution waves are missing => abort |
---|
145 | if(ResolutionWavesMissingDF(str)) //updated to NOT use global strings (in GaussUtils) |
---|
146 | Abort |
---|
147 | endif |
---|
148 | |
---|
149 | SetDataFolder $("root:"+str) |
---|
150 | |
---|
151 | // Setup parameter table for model function |
---|
152 | Make/O/D smear_coef_PolyOneShell = {1.,60,0.1,6.4e-6,10,1e-6,6.4e-6,0.001} |
---|
153 | make/o/t smear_parameters_PolyOneShell = {"scale","core radius (A)","Core Polydispersity(0,1)","Core SLD (A-2)","Shell thickness (A)","Shell SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"} |
---|
154 | Edit smear_parameters_PolyOneShell,smear_coef_PolyOneShell //display parameters in a table |
---|
155 | |
---|
156 | // output smeared intensity wave, dimensions are identical to experimental QSIG values |
---|
157 | // make extra copy of experimental q-values for easy plotting |
---|
158 | Duplicate/O $(str+"_q") smeared_PolyOneShell,smeared_qvals |
---|
159 | SetScale d,0,0,"1/cm",smeared_PolyOneShell |
---|
160 | |
---|
161 | Variable/G gs_PolyOneShell=0 |
---|
162 | gs_PolyOneShell := fSmearedPolyOneShell(smear_coef_PolyOneShell,smeared_PolyOneShell,smeared_qvals) //this wrapper fills the STRUCT |
---|
163 | |
---|
164 | Display smeared_PolyOneShell vs smeared_qvals |
---|
165 | ModifyGraph log=1,marker=29,msize=2,mode=4 |
---|
166 | Label bottom "q (\\S-1\\M)" |
---|
167 | Label left "I(q) (cm\\S-1\\M)" |
---|
168 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
169 | |
---|
170 | SetDataFolder root: |
---|
171 | AddModelToStrings("SmearedPolyOneShell","smear_coef_PolyOneShell","PolyOneShell") |
---|
172 | End |
---|
173 | |
---|
174 | Proc PlotSmearedPolyTwoShell(str) |
---|
175 | String str |
---|
176 | Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4) |
---|
177 | |
---|
178 | // if any of the resolution waves are missing => abort |
---|
179 | if(ResolutionWavesMissingDF(str)) //updated to NOT use global strings (in GaussUtils) |
---|
180 | Abort |
---|
181 | endif |
---|
182 | |
---|
183 | SetDataFolder $("root:"+str) |
---|
184 | |
---|
185 | // Setup parameter table for model function |
---|
186 | Make/O/D smear_coef_PolyTwoShell = {1.,60,0.1,6.4e-6,10,1e-6,10,2e-6,6.4e-6,0.001} |
---|
187 | make/o/t smear_parameters_PolyTwoShell = {"scale","core radius (A)","Core Polydispersity(0,1)","Core SLD (A-2)","Shell 1 thickness","Shell 1 SLD (A-2)","Shell 2 thickness","Shell 2 SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"} |
---|
188 | Edit smear_parameters_PolyTwoShell,smear_coef_PolyTwoShell //display parameters in a table |
---|
189 | |
---|
190 | // output smeared intensity wave, dimensions are identical to experimental QSIG values |
---|
191 | // make extra copy of experimental q-values for easy plotting |
---|
192 | Duplicate/O $(str+"_q") smeared_PolyTwoShell,smeared_qvals |
---|
193 | SetScale d,0,0,"1/cm",smeared_PolyTwoShell |
---|
194 | |
---|
195 | Variable/G gs_PolyTwoShell=0 |
---|
196 | gs_PolyTwoShell := fSmearedPolyTwoShell(smear_coef_PolyTwoShell,smeared_PolyTwoShell,smeared_qvals) //this wrapper fills the STRUCT |
---|
197 | |
---|
198 | Display smeared_PolyTwoShell vs smeared_qvals |
---|
199 | ModifyGraph log=1,marker=29,msize=2,mode=4 |
---|
200 | Label bottom "q (\\S-1\\M)" |
---|
201 | Label left "I(q) (cm\\S-1\\M)" |
---|
202 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
203 | |
---|
204 | SetDataFolder root: |
---|
205 | AddModelToStrings("SmearedPolyTwoShell","smear_coef_PolyTwoShell","PolyTwoShell") |
---|
206 | End |
---|
207 | |
---|
208 | |
---|
209 | Proc PlotSmearedPolyThreeShell(str) |
---|
210 | String str |
---|
211 | Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4) |
---|
212 | |
---|
213 | // if any of the resolution waves are missing => abort |
---|
214 | if(ResolutionWavesMissingDF(str)) //updated to NOT use global strings (in GaussUtils) |
---|
215 | Abort |
---|
216 | endif |
---|
217 | |
---|
218 | SetDataFolder $("root:"+str) |
---|
219 | |
---|
220 | // Setup parameter table for model function |
---|
221 | Make/O/D smear_coef_PolyThreeShell = {1.,60,0.1,6.4e-6,10,1e-6,10,2e-6,5,3e-6,6.4e-6,0.001} |
---|
222 | make/o/t smear_parameters_PolyThreeShell = {"scale","core radius (A)","Core Polydispersity(0,1)","Core SLD (A-2)","Shell 1 thickness","Shell 1 SLD (A-2)","Shell 2 thickness","Shell 2 SLD (A-2)","Shell 3 thickness","Shell 3 SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"} |
---|
223 | Edit smear_parameters_PolyThreeShell,smear_coef_PolyThreeShell //display parameters in a table |
---|
224 | |
---|
225 | // output smeared intensity wave, dimensions are identical to experimental QSIG values |
---|
226 | // make extra copy of experimental q-values for easy plotting |
---|
227 | Duplicate/O $(str+"_q") smeared_PolyThreeShell,smeared_qvals |
---|
228 | SetScale d,0,0,"1/cm",smeared_PolyThreeShell |
---|
229 | |
---|
230 | Variable/G gs_PolyThreeShell=0 |
---|
231 | gs_PolyThreeShell := fSmearedPolyThreeShell(smear_coef_PolyThreeShell,smeared_PolyThreeShell,smeared_qvals) //this wrapper fills the STRUCT |
---|
232 | |
---|
233 | Display smeared_PolyThreeShell vs smeared_qvals |
---|
234 | ModifyGraph log=1,marker=29,msize=2,mode=4 |
---|
235 | Label bottom "q (\\S-1\\M)" |
---|
236 | Label left "I(q) (cm\\S-1\\M)" |
---|
237 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
238 | |
---|
239 | SetDataFolder root: |
---|
240 | AddModelToStrings("SmearedPolyThreeShell","smear_coef_PolyThreeShell","PolyThreeShell") |
---|
241 | End |
---|
242 | |
---|
243 | Proc PlotSmearedPolyFourShell(str) |
---|
244 | String str |
---|
245 | Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4) |
---|
246 | |
---|
247 | // if any of the resolution waves are missing => abort |
---|
248 | if(ResolutionWavesMissingDF(str)) //updated to NOT use global strings (in GaussUtils) |
---|
249 | Abort |
---|
250 | endif |
---|
251 | |
---|
252 | SetDataFolder $("root:"+str) |
---|
253 | |
---|
254 | // Setup parameter table for model function |
---|
255 | Make/O/D smear_coef_PolyFourShell = {1.,60,0.1,6.4e-6,10,1e-6,10,2e-6,5,3e-6,10,4e-6,6.4e-6,0.001} |
---|
256 | make/o/t smear_parameters_PolyFourShell = {"scale","core radius (A)","Core Polydispersity(0,1)","Core SLD (A-2)","Shell 1 thickness","Shell 1 SLD (A-2)","Shell 2 thickness","Shell 2 SLD (A-2)","Shell 3 thickness","Shell 3 SLD (A-2)","Shell 4 thickness","Shell 4 SLD (A-2)","Solvent SLD (A-2)","bkg (cm-1)"} |
---|
257 | Edit smear_parameters_PolyFourShell,smear_coef_PolyFourShell //display parameters in a table |
---|
258 | |
---|
259 | // output smeared intensity wave, dimensions are identical to experimental QSIG values |
---|
260 | // make extra copy of experimental q-values for easy plotting |
---|
261 | Duplicate/O $(str+"_q") smeared_PolyFourShell,smeared_qvals |
---|
262 | SetScale d,0,0,"1/cm",smeared_PolyFourShell |
---|
263 | |
---|
264 | Variable/G gs_PolyFourShell=0 |
---|
265 | gs_PolyFourShell := fSmearedPolyFourShell(smear_coef_PolyFourShell,smeared_PolyFourShell,smeared_qvals) //this wrapper fills the STRUCT |
---|
266 | |
---|
267 | Display smeared_PolyFourShell vs smeared_qvals |
---|
268 | ModifyGraph log=1,marker=29,msize=2,mode=4 |
---|
269 | Label bottom "q (\\S-1\\M)" |
---|
270 | Label left "I(q) (cm\\S-1\\M)" |
---|
271 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
272 | |
---|
273 | SetDataFolder root: |
---|
274 | AddModelToStrings("SmearedPolyFourShell","smear_coef_PolyFourShell","PolyFourShell") |
---|
275 | End |
---|
276 | |
---|
277 | |
---|
278 | // nothing to change here |
---|
279 | // |
---|
280 | //AAO version, uses XOP if available |
---|
281 | // simply calls the original single point calculation with |
---|
282 | // a wave assignment (this will behave nicely if given point ranges) |
---|
283 | Function PolyOneShell(cw,yw,xw) : FitFunc |
---|
284 | Wave cw,yw,xw |
---|
285 | |
---|
286 | #if exists("PolyOneShellX") |
---|
287 | yw = PolyOneShellX(cw,xw) |
---|
288 | #else |
---|
289 | yw = fPolyOneShell(cw,xw) |
---|
290 | #endif |
---|
291 | return(0) |
---|
292 | End |
---|
293 | |
---|
294 | Function PolyTwoShell(cw,yw,xw) : FitFunc |
---|
295 | Wave cw,yw,xw |
---|
296 | |
---|
297 | #if exists("PolyTwoShellX") |
---|
298 | yw = PolyTwoShellX(cw,xw) |
---|
299 | #else |
---|
300 | yw = fPolyTwoShell(cw,xw) |
---|
301 | #endif |
---|
302 | return(0) |
---|
303 | End |
---|
304 | |
---|
305 | Function PolyThreeShell(cw,yw,xw) : FitFunc |
---|
306 | Wave cw,yw,xw |
---|
307 | |
---|
308 | #if exists("PolyThreeShellX") |
---|
309 | yw =PolyThreeShellX(cw,xw) |
---|
310 | #else |
---|
311 | yw = fPolyThreeShell(cw,xw) |
---|
312 | #endif |
---|
313 | return(0) |
---|
314 | End |
---|
315 | |
---|
316 | Function PolyFourShell(cw,yw,xw) : FitFunc |
---|
317 | Wave cw,yw,xw |
---|
318 | |
---|
319 | #if exists("PolyFourShellX") |
---|
320 | yw = PolyFourShellX(cw,xw) |
---|
321 | #else |
---|
322 | yw = fPolyFourShell(cw,xw) |
---|
323 | #endif |
---|
324 | return(0) |
---|
325 | End |
---|
326 | |
---|
327 | |
---|
328 | // |
---|
329 | // unsmeared model calculation |
---|
330 | // |
---|
331 | Function fPolyOneShell(w,x) : FitFunc |
---|
332 | Wave w |
---|
333 | Variable x |
---|
334 | |
---|
335 | Variable scale,rcore,thick,rhocore,rhoshel,rhosolv,bkg,pd,zz |
---|
336 | scale = w[0] |
---|
337 | rcore = w[1] |
---|
338 | pd = w[2] |
---|
339 | rhocore = w[3] |
---|
340 | thick = w[4] |
---|
341 | rhoshel = w[5] |
---|
342 | rhosolv = w[6] |
---|
343 | bkg = w[7] |
---|
344 | |
---|
345 | zz = (1/pd)^2-1 //polydispersity of the core only |
---|
346 | // |
---|
347 | // local variables |
---|
348 | Variable nord,ii,va,vb,contr,vcyl,nden,summ,yyy,zi,qq |
---|
349 | Variable answer,zp1,zp2,zp3,vpoly |
---|
350 | String weightStr,zStr |
---|
351 | |
---|
352 | //select number of gauss points by setting nord=20 or76 points |
---|
353 | NVAR/Z gNord=gNord |
---|
354 | if(! NVAR_Exists(gNord) ) |
---|
355 | nord=76 //use 76 pts as default |
---|
356 | else |
---|
357 | if( (gNord == 20) || (gNord ==76) ) |
---|
358 | nord = gNord // should only allow 20 or 76 points |
---|
359 | else |
---|
360 | abort "global value gNord in SchulzSpheres must be either 20 or 76" |
---|
361 | endif |
---|
362 | endif |
---|
363 | |
---|
364 | weightStr = "gauss"+num2str(nord)+"wt" |
---|
365 | zStr = "gauss"+num2str(nord)+"z" |
---|
366 | |
---|
367 | if (WaveExists($weightStr) == 0) // wave reference is not valid, |
---|
368 | Make/D/N=(nord) $weightStr,$zStr |
---|
369 | Wave gauWt = $weightStr |
---|
370 | Wave gauZ = $zStr // wave references to pass |
---|
371 | if(nord==20) |
---|
372 | Make20GaussPoints(gauWt,gauZ) |
---|
373 | else |
---|
374 | Make76GaussPoints(gauWt,gauZ) |
---|
375 | endif |
---|
376 | else |
---|
377 | if(exists(weightStr) > 1) |
---|
378 | Abort "wave name is already in use" //executed only if name is in use elsewhere |
---|
379 | endif |
---|
380 | Wave gauWt = $weightStr |
---|
381 | Wave gauZ = $zStr // create the wave references |
---|
382 | endif |
---|
383 | |
---|
384 | // set up the integration end points and weights |
---|
385 | // limits are technically 0-inf, but wisely choose non-zero region of distribution |
---|
386 | Variable range=8 //multiples of the std. dev. from the mean |
---|
387 | va = rcore*(1-range*pd) |
---|
388 | if (va<0) |
---|
389 | va=0 //otherwise numerical error when pd >= 0.3, making a<0 |
---|
390 | endif |
---|
391 | If(pd>0.3) |
---|
392 | range = range + (pd-0.3)*18 //stretch upper range to account for skewed tail |
---|
393 | Endif |
---|
394 | vb = rcore*(1+range*pd) // is this far enough past avg radius? |
---|
395 | |
---|
396 | //temp set scale=1 and bkg=0 for quadrature calc |
---|
397 | Make/O/D/N=7 temp_1sf |
---|
398 | temp_1sf[0] = 1 |
---|
399 | temp_1sf[1] = w[1] //the core radius will be changed in the loop |
---|
400 | temp_1sf[2] = w[3] |
---|
401 | temp_1sf[3] = w[4] |
---|
402 | temp_1sf[4] = w[5] |
---|
403 | temp_1sf[5] = w[6] |
---|
404 | temp_1sf[6] = 0 |
---|
405 | |
---|
406 | // evaluate at Gauss points |
---|
407 | summ = 0.0 // initialize integral |
---|
408 | for(ii=0;ii<nord;ii+=1) |
---|
409 | zi = ( gauZ[ii]*(vb-va) + vb + va )/2.0 |
---|
410 | temp_1sf[1] = zi |
---|
411 | yyy = gauWt[ii] * Schulz_Point_Nsf(zi,rcore,zz) *fOneShell(temp_1sf,x) |
---|
412 | //un-normalize by volume |
---|
413 | yyy *= 4*pi/3*(zi+thick)^3 |
---|
414 | summ += yyy |
---|
415 | endfor |
---|
416 | // calculate value of integral to return |
---|
417 | answer = (vb-va)/2.0*summ |
---|
418 | |
---|
419 | //re-normalize by the average volume |
---|
420 | zp1 = zz + 1. |
---|
421 | zp2 = zz + 2. |
---|
422 | zp3 = zz + 3. |
---|
423 | vpoly = 4*Pi/3*zp3*zp2/zp1/zp1*(rcore+thick)^3 |
---|
424 | answer /= vpoly |
---|
425 | //scale |
---|
426 | answer *= scale |
---|
427 | // add in the background |
---|
428 | answer += bkg |
---|
429 | |
---|
430 | Return (answer) |
---|
431 | End |
---|
432 | |
---|
433 | |
---|
434 | Function fPolyTwoShell(w,x) : FitFunc |
---|
435 | Wave w |
---|
436 | Variable x |
---|
437 | |
---|
438 | Variable scale,rcore,rhocore,rhosolv,bkg,pd,zz |
---|
439 | Variable thick1,thick2 |
---|
440 | Variable rhoshel1,rhoshel2 |
---|
441 | scale = w[0] |
---|
442 | rcore = w[1] |
---|
443 | pd = w[2] |
---|
444 | rhocore = w[3] |
---|
445 | thick1 = w[4] |
---|
446 | rhoshel1 = w[5] |
---|
447 | thick2 = w[6] |
---|
448 | rhoshel2 = w[7] |
---|
449 | rhosolv = w[8] |
---|
450 | bkg = w[9] |
---|
451 | |
---|
452 | zz = (1/pd)^2-1 //polydispersity of the core only |
---|
453 | // |
---|
454 | // local variables |
---|
455 | Variable nord,ii,va,vb,contr,vcyl,nden,summ,yyy,zi,qq |
---|
456 | Variable answer,zp1,zp2,zp3,vpoly |
---|
457 | String weightStr,zStr |
---|
458 | |
---|
459 | //select number of gauss points by setting nord=20 or76 points |
---|
460 | NVAR/Z gNord=gNord |
---|
461 | if(! NVAR_Exists(gNord) ) |
---|
462 | nord=76 //use 76 pts as default |
---|
463 | else |
---|
464 | if( (gNord == 20) || (gNord ==76) ) |
---|
465 | nord = gNord // should only allow 20 or 76 points |
---|
466 | else |
---|
467 | abort "global value gNord in SchulzSpheres must be either 20 or 76" |
---|
468 | endif |
---|
469 | endif |
---|
470 | |
---|
471 | weightStr = "gauss"+num2str(nord)+"wt" |
---|
472 | zStr = "gauss"+num2str(nord)+"z" |
---|
473 | |
---|
474 | if (WaveExists($weightStr) == 0) // wave reference is not valid, |
---|
475 | Make/D/N=(nord) $weightStr,$zStr |
---|
476 | Wave gauWt = $weightStr |
---|
477 | Wave gauZ = $zStr // wave references to pass |
---|
478 | if(nord==20) |
---|
479 | Make20GaussPoints(gauWt,gauZ) |
---|
480 | else |
---|
481 | Make76GaussPoints(gauWt,gauZ) |
---|
482 | endif |
---|
483 | else |
---|
484 | if(exists(weightStr) > 1) |
---|
485 | Abort "wave name is already in use" //executed only if name is in use elsewhere |
---|
486 | endif |
---|
487 | Wave gauWt = $weightStr |
---|
488 | Wave gauZ = $zStr // create the wave references |
---|
489 | endif |
---|
490 | |
---|
491 | // set up the integration end points and weights |
---|
492 | // limits are technically 0-inf, but wisely choose non-zero region of distribution |
---|
493 | Variable range=8 //multiples of the std. dev. from the mean |
---|
494 | va = rcore*(1-range*pd) |
---|
495 | if (va<0) |
---|
496 | va=0 //otherwise numerical error when pd >= 0.3, making a<0 |
---|
497 | endif |
---|
498 | If(pd>0.3) |
---|
499 | range = range + (pd-0.3)*18 //stretch upper range to account for skewed tail |
---|
500 | Endif |
---|
501 | vb = rcore*(1+range*pd) // is this far enough past avg radius? |
---|
502 | |
---|
503 | //temp set scale=1 and bkg=0 for quadrature calc |
---|
504 | Make/O/D/N=9 temp_2sf |
---|
505 | temp_2sf[0] = 1 |
---|
506 | temp_2sf[1] = w[1] //the core radius will be changed in the loop |
---|
507 | temp_2sf[2] = w[3] |
---|
508 | temp_2sf[3] = w[4] |
---|
509 | temp_2sf[4] = w[5] |
---|
510 | temp_2sf[5] = w[6] |
---|
511 | temp_2sf[6] = w[7] |
---|
512 | temp_2sf[7] = w[8] |
---|
513 | temp_2sf[8] = 0 |
---|
514 | |
---|
515 | // evaluate at Gauss points |
---|
516 | summ = 0.0 // initialize integral |
---|
517 | for(ii=0;ii<nord;ii+=1) |
---|
518 | zi = ( gauZ[ii]*(vb-va) + vb + va )/2.0 |
---|
519 | temp_2sf[1] = zi |
---|
520 | yyy = gauWt[ii] * Schulz_Point_Nsf(zi,rcore,zz) * fTwoShell(temp_2sf,x) |
---|
521 | //un-normalize by volume |
---|
522 | yyy *= 4*pi/3*(zi+thick1+thick2)^3 |
---|
523 | summ += yyy |
---|
524 | endfor |
---|
525 | // calculate value of integral to return |
---|
526 | answer = (vb-va)/2.0*summ |
---|
527 | |
---|
528 | //re-normalize by the average volume |
---|
529 | zp1 = zz + 1. |
---|
530 | zp2 = zz + 2. |
---|
531 | zp3 = zz + 3. |
---|
532 | vpoly = 4*Pi/3*zp3*zp2/zp1/zp1*(rcore+thick1+thick2)^3 |
---|
533 | answer /= vpoly |
---|
534 | //scale |
---|
535 | answer *= scale |
---|
536 | // add in the background |
---|
537 | answer += bkg |
---|
538 | |
---|
539 | Return (answer) |
---|
540 | End |
---|
541 | |
---|
542 | Function fPolyThreeShell(w,x) : FitFunc |
---|
543 | Wave w |
---|
544 | Variable x |
---|
545 | |
---|
546 | Variable scale,rcore,rhocore,rhosolv,bkg,pd,zz |
---|
547 | Variable thick1,thick2,thick3 |
---|
548 | Variable rhoshel1,rhoshel2,rhoshel3 |
---|
549 | scale = w[0] |
---|
550 | rcore = w[1] |
---|
551 | pd = w[2] |
---|
552 | rhocore = w[3] |
---|
553 | thick1 = w[4] |
---|
554 | rhoshel1 = w[5] |
---|
555 | thick2 = w[6] |
---|
556 | rhoshel2 = w[7] |
---|
557 | thick3 = w[8] |
---|
558 | rhoshel3 = w[9] |
---|
559 | rhosolv = w[10] |
---|
560 | bkg = w[11] |
---|
561 | |
---|
562 | zz = (1/pd)^2-1 //polydispersity of the core only |
---|
563 | // |
---|
564 | // local variables |
---|
565 | Variable nord,ii,va,vb,contr,vcyl,nden,summ,yyy,zi,qq |
---|
566 | Variable answer,zp1,zp2,zp3,vpoly |
---|
567 | String weightStr,zStr |
---|
568 | |
---|
569 | //select number of gauss points by setting nord=20 or76 points |
---|
570 | NVAR/Z gNord=gNord |
---|
571 | if(! NVAR_Exists(gNord) ) |
---|
572 | nord=76 //use 76 pts as default |
---|
573 | else |
---|
574 | if( (gNord == 20) || (gNord ==76) ) |
---|
575 | nord = gNord // should only allow 20 or 76 points |
---|
576 | else |
---|
577 | abort "global value gNord in SchulzSpheres must be either 20 or 76" |
---|
578 | endif |
---|
579 | endif |
---|
580 | |
---|
581 | weightStr = "gauss"+num2str(nord)+"wt" |
---|
582 | zStr = "gauss"+num2str(nord)+"z" |
---|
583 | |
---|
584 | if (WaveExists($weightStr) == 0) // wave reference is not valid, |
---|
585 | Make/D/N=(nord) $weightStr,$zStr |
---|
586 | Wave gauWt = $weightStr |
---|
587 | Wave gauZ = $zStr // wave references to pass |
---|
588 | if(nord==20) |
---|
589 | Make20GaussPoints(gauWt,gauZ) |
---|
590 | else |
---|
591 | Make76GaussPoints(gauWt,gauZ) |
---|
592 | endif |
---|
593 | else |
---|
594 | if(exists(weightStr) > 1) |
---|
595 | Abort "wave name is already in use" //executed only if name is in use elsewhere |
---|
596 | endif |
---|
597 | Wave gauWt = $weightStr |
---|
598 | Wave gauZ = $zStr // create the wave references |
---|
599 | endif |
---|
600 | |
---|
601 | // set up the integration end points and weights |
---|
602 | // limits are technically 0-inf, but wisely choose non-zero region of distribution |
---|
603 | Variable range=8 //multiples of the std. dev. from the mean |
---|
604 | va = rcore*(1-range*pd) |
---|
605 | if (va<0) |
---|
606 | va=0 //otherwise numerical error when pd >= 0.3, making a<0 |
---|
607 | endif |
---|
608 | If(pd>0.3) |
---|
609 | range = range + (pd-0.3)*18 //stretch upper range to account for skewed tail |
---|
610 | Endif |
---|
611 | vb = rcore*(1+range*pd) // is this far enough past avg radius? |
---|
612 | |
---|
613 | //temp set scale=1 and bkg=0 for quadrature calc |
---|
614 | Make/O/D/N=11 temp_3sf |
---|
615 | temp_3sf[0] = 1 |
---|
616 | temp_3sf[1] = w[1] //the core radius will be changed in the loop |
---|
617 | temp_3sf[2] = w[3] |
---|
618 | temp_3sf[3] = w[4] |
---|
619 | temp_3sf[4] = w[5] |
---|
620 | temp_3sf[5] = w[6] |
---|
621 | temp_3sf[6] = w[7] |
---|
622 | temp_3sf[7] = w[8] |
---|
623 | temp_3sf[8] = w[9] |
---|
624 | temp_3sf[9] = w[10] |
---|
625 | temp_3sf[10] = 0 |
---|
626 | |
---|
627 | // evaluate at Gauss points |
---|
628 | summ = 0.0 // initialize integral |
---|
629 | for(ii=0;ii<nord;ii+=1) |
---|
630 | zi = ( gauZ[ii]*(vb-va) + vb + va )/2.0 |
---|
631 | temp_3sf[1] = zi |
---|
632 | yyy = gauWt[ii] * Schulz_Point_Nsf(zi,rcore,zz) * fThreeShell(temp_3sf,x) |
---|
633 | //un-normalize by volume |
---|
634 | yyy *= 4*pi/3*(zi+thick1+thick2+thick3)^3 |
---|
635 | summ += yyy |
---|
636 | endfor |
---|
637 | // calculate value of integral to return |
---|
638 | answer = (vb-va)/2.0*summ |
---|
639 | |
---|
640 | //re-normalize by the average volume |
---|
641 | zp1 = zz + 1. |
---|
642 | zp2 = zz + 2. |
---|
643 | zp3 = zz + 3. |
---|
644 | vpoly = 4*Pi/3*zp3*zp2/zp1/zp1*(rcore+thick1+thick2+thick3)^3 |
---|
645 | answer /= vpoly |
---|
646 | //scale |
---|
647 | answer *= scale |
---|
648 | // add in the background |
---|
649 | answer += bkg |
---|
650 | |
---|
651 | Return (answer) |
---|
652 | End |
---|
653 | |
---|
654 | |
---|
655 | Function fPolyFourShell(w,x) : FitFunc |
---|
656 | Wave w |
---|
657 | Variable x |
---|
658 | |
---|
659 | Variable scale,rcore,rhocore,rhosolv,bkg,pd,zz |
---|
660 | Variable thick1,thick2,thick3,thick4 |
---|
661 | Variable rhoshel1,rhoshel2,rhoshel3,rhoshel4 |
---|
662 | scale = w[0] |
---|
663 | rcore = w[1] |
---|
664 | pd = w[2] |
---|
665 | rhocore = w[3] |
---|
666 | thick1 = w[4] |
---|
667 | rhoshel1 = w[5] |
---|
668 | thick2 = w[6] |
---|
669 | rhoshel2 = w[7] |
---|
670 | thick3 = w[8] |
---|
671 | rhoshel3 = w[9] |
---|
672 | thick4 = w[10] |
---|
673 | rhoshel4 = w[11] |
---|
674 | rhosolv = w[12] |
---|
675 | bkg = w[13] |
---|
676 | |
---|
677 | zz = (1/pd)^2-1 //polydispersity of the core only |
---|
678 | // |
---|
679 | // local variables |
---|
680 | Variable nord,ii,va,vb,contr,vcyl,nden,summ,yyy,zi,qq |
---|
681 | Variable answer,zp1,zp2,zp3,vpoly |
---|
682 | String weightStr,zStr |
---|
683 | |
---|
684 | //select number of gauss points by setting nord=20 or76 points |
---|
685 | NVAR/Z gNord=gNord |
---|
686 | if(! NVAR_Exists(gNord) ) |
---|
687 | nord=76 //use 76 pts as default |
---|
688 | else |
---|
689 | if( (gNord == 20) || (gNord ==76) ) |
---|
690 | nord = gNord // should only allow 20 or 76 points |
---|
691 | else |
---|
692 | abort "global value gNord in SchulzSpheres must be either 20 or 76" |
---|
693 | endif |
---|
694 | endif |
---|
695 | |
---|
696 | weightStr = "gauss"+num2str(nord)+"wt" |
---|
697 | zStr = "gauss"+num2str(nord)+"z" |
---|
698 | |
---|
699 | if (WaveExists($weightStr) == 0) // wave reference is not valid, |
---|
700 | Make/D/N=(nord) $weightStr,$zStr |
---|
701 | Wave gauWt = $weightStr |
---|
702 | Wave gauZ = $zStr // wave references to pass |
---|
703 | if(nord==20) |
---|
704 | Make20GaussPoints(gauWt,gauZ) |
---|
705 | else |
---|
706 | Make76GaussPoints(gauWt,gauZ) |
---|
707 | endif |
---|
708 | else |
---|
709 | if(exists(weightStr) > 1) |
---|
710 | Abort "wave name is already in use" //executed only if name is in use elsewhere |
---|
711 | endif |
---|
712 | Wave gauWt = $weightStr |
---|
713 | Wave gauZ = $zStr // create the wave references |
---|
714 | endif |
---|
715 | |
---|
716 | // set up the integration end points and weights |
---|
717 | // limits are technically 0-inf, but wisely choose non-zero region of distribution |
---|
718 | Variable range=8 //multiples of the std. dev. from the mean |
---|
719 | va = rcore*(1-range*pd) |
---|
720 | if (va<0) |
---|
721 | va=0 //otherwise numerical error when pd >= 0.3, making a<0 |
---|
722 | endif |
---|
723 | If(pd>0.3) |
---|
724 | range = range + (pd-0.3)*18 //stretch upper range to account for skewed tail |
---|
725 | Endif |
---|
726 | vb = rcore*(1+range*pd) // is this far enough past avg radius? |
---|
727 | |
---|
728 | //temp set scale=1 and bkg=0 for quadrature calc |
---|
729 | Make/O/D/N=13 temp_4sf |
---|
730 | temp_4sf[0] = 1 |
---|
731 | temp_4sf[1] = w[1] //the core radius will be changed in the loop |
---|
732 | temp_4sf[2] = w[3] |
---|
733 | temp_4sf[3] = w[4] |
---|
734 | temp_4sf[4] = w[5] |
---|
735 | temp_4sf[5] = w[6] |
---|
736 | temp_4sf[6] = w[7] |
---|
737 | temp_4sf[7] = w[8] |
---|
738 | temp_4sf[8] = w[9] |
---|
739 | temp_4sf[9] = w[10] |
---|
740 | temp_4sf[10] = w[11] |
---|
741 | temp_4sf[11] = w[12] |
---|
742 | temp_4sf[12] = 0 |
---|
743 | |
---|
744 | // evaluate at Gauss points |
---|
745 | summ = 0.0 // initialize integral |
---|
746 | for(ii=0;ii<nord;ii+=1) |
---|
747 | zi = ( gauZ[ii]*(vb-va) + vb + va )/2.0 |
---|
748 | temp_4sf[1] = zi |
---|
749 | yyy = gauWt[ii] * Schulz_Point_Nsf(zi,rcore,zz) * fFourShell(temp_4sf,x) |
---|
750 | //un-normalize by volume |
---|
751 | yyy *= 4*pi/3*(zi+thick1+thick2+thick3+thick4)^3 |
---|
752 | summ += yyy |
---|
753 | endfor |
---|
754 | // calculate value of integral to return |
---|
755 | answer = (vb-va)/2.0*summ |
---|
756 | |
---|
757 | //re-normalize by the average volume |
---|
758 | zp1 = zz + 1. |
---|
759 | zp2 = zz + 2. |
---|
760 | zp3 = zz + 3. |
---|
761 | vpoly = 4*Pi/3*zp3*zp2/zp1/zp1*(rcore+thick1+thick2+thick3+thick4)^3 |
---|
762 | answer /= vpoly |
---|
763 | //scale |
---|
764 | answer *= scale |
---|
765 | // add in the background |
---|
766 | answer += bkg |
---|
767 | |
---|
768 | Return (answer) |
---|
769 | End |
---|
770 | |
---|
771 | |
---|
772 | |
---|
773 | |
---|
774 | /////////////////////////////////////////////////////////////// |
---|
775 | // smeared model calculation |
---|
776 | // |
---|
777 | // you don't need to do anything with this function, as long as |
---|
778 | // your PolyOneShell works correctly, you get the resolution-smeared |
---|
779 | // version for free. |
---|
780 | // |
---|
781 | // this is all there is to the smeared model calculation! |
---|
782 | Function SmearedPolyOneShell(s) : FitFunc |
---|
783 | Struct ResSmearAAOStruct &s |
---|
784 | |
---|
785 | // the name of your unsmeared model (AAO) is the first argument |
---|
786 | Smear_Model_20(PolyOneShell,s.coefW,s.xW,s.yW,s.resW) |
---|
787 | |
---|
788 | return(0) |
---|
789 | End |
---|
790 | |
---|
791 | Function SmearedPolyTwoShell(s) : FitFunc |
---|
792 | Struct ResSmearAAOStruct &s |
---|
793 | |
---|
794 | // the name of your unsmeared model (AAO) is the first argument |
---|
795 | Smear_Model_20(PolyTwoShell,s.coefW,s.xW,s.yW,s.resW) |
---|
796 | |
---|
797 | return(0) |
---|
798 | End |
---|
799 | |
---|
800 | Function SmearedPolyThreeShell(s) : FitFunc |
---|
801 | Struct ResSmearAAOStruct &s |
---|
802 | |
---|
803 | // the name of your unsmeared model (AAO) is the first argument |
---|
804 | Smear_Model_20(PolyThreeShell,s.coefW,s.xW,s.yW,s.resW) |
---|
805 | |
---|
806 | return(0) |
---|
807 | End |
---|
808 | |
---|
809 | Function SmearedPolyFourShell(s) : FitFunc |
---|
810 | Struct ResSmearAAOStruct &s |
---|
811 | |
---|
812 | // the name of your unsmeared model (AAO) is the first argument |
---|
813 | Smear_Model_20(PolyFourShell,s.coefW,s.xW,s.yW,s.resW) |
---|
814 | |
---|
815 | return(0) |
---|
816 | End |
---|
817 | |
---|
818 | |
---|
819 | |
---|
820 | /////////////////////////////////////////////////////////////// |
---|
821 | |
---|
822 | |
---|
823 | // nothing to change here |
---|
824 | // |
---|
825 | //wrapper to calculate the smeared model as an AAO-Struct |
---|
826 | // fills the struct and calls the ususal function with the STRUCT parameter |
---|
827 | // |
---|
828 | // used only for the dependency, not for fitting |
---|
829 | // |
---|
830 | Function fSmearedPolyOneShell(coefW,yW,xW) |
---|
831 | Wave coefW,yW,xW |
---|
832 | |
---|
833 | String str = getWavesDataFolder(yW,0) |
---|
834 | String DF="root:"+str+":" |
---|
835 | |
---|
836 | WAVE resW = $(DF+str+"_res") |
---|
837 | |
---|
838 | STRUCT ResSmearAAOStruct fs |
---|
839 | WAVE fs.coefW = coefW |
---|
840 | WAVE fs.yW = yW |
---|
841 | WAVE fs.xW = xW |
---|
842 | WAVE fs.resW = resW |
---|
843 | |
---|
844 | Variable err |
---|
845 | err = SmearedPolyOneShell(fs) |
---|
846 | |
---|
847 | return (0) |
---|
848 | End |
---|
849 | |
---|
850 | Function fSmearedPolyTwoShell(coefW,yW,xW) |
---|
851 | Wave coefW,yW,xW |
---|
852 | |
---|
853 | String str = getWavesDataFolder(yW,0) |
---|
854 | String DF="root:"+str+":" |
---|
855 | |
---|
856 | WAVE resW = $(DF+str+"_res") |
---|
857 | |
---|
858 | STRUCT ResSmearAAOStruct fs |
---|
859 | WAVE fs.coefW = coefW |
---|
860 | WAVE fs.yW = yW |
---|
861 | WAVE fs.xW = xW |
---|
862 | WAVE fs.resW = resW |
---|
863 | |
---|
864 | Variable err |
---|
865 | err = SmearedPolyTwoShell(fs) |
---|
866 | |
---|
867 | return (0) |
---|
868 | End |
---|
869 | |
---|
870 | Function fSmearedPolyThreeShell(coefW,yW,xW) |
---|
871 | Wave coefW,yW,xW |
---|
872 | |
---|
873 | String str = getWavesDataFolder(yW,0) |
---|
874 | String DF="root:"+str+":" |
---|
875 | |
---|
876 | WAVE resW = $(DF+str+"_res") |
---|
877 | |
---|
878 | STRUCT ResSmearAAOStruct fs |
---|
879 | WAVE fs.coefW = coefW |
---|
880 | WAVE fs.yW = yW |
---|
881 | WAVE fs.xW = xW |
---|
882 | WAVE fs.resW = resW |
---|
883 | |
---|
884 | Variable err |
---|
885 | err = SmearedPolyThreeShell(fs) |
---|
886 | |
---|
887 | return (0) |
---|
888 | End |
---|
889 | |
---|
890 | Function fSmearedPolyFourShell(coefW,yW,xW) |
---|
891 | Wave coefW,yW,xW |
---|
892 | |
---|
893 | String str = getWavesDataFolder(yW,0) |
---|
894 | String DF="root:"+str+":" |
---|
895 | |
---|
896 | WAVE resW = $(DF+str+"_res") |
---|
897 | |
---|
898 | STRUCT ResSmearAAOStruct fs |
---|
899 | WAVE fs.coefW = coefW |
---|
900 | WAVE fs.yW = yW |
---|
901 | WAVE fs.xW = xW |
---|
902 | WAVE fs.resW = resW |
---|
903 | |
---|
904 | Variable err |
---|
905 | err = SmearedPolyFourShell(fs) |
---|
906 | |
---|
907 | return (0) |
---|
908 | End |
---|
909 | |
---|
910 | //calculate the normalized distribution |
---|
911 | Static Function Schulz_Point_Nsf(x,avg,zz) |
---|
912 | Variable x,avg,zz |
---|
913 | |
---|
914 | Variable dr |
---|
915 | |
---|
916 | dr = zz*ln(x) - gammln(zz+1)+(zz+1)*ln((zz+1)/avg)-(x/avg*(zz+1)) |
---|
917 | return (exp(dr)) |
---|
918 | End |
---|