1 | #pragma rtGlobals=1 // Use modern global access method. |
---|
2 | #pragma IgorVersion = 6.0 |
---|
3 | |
---|
4 | |
---|
5 | //////////////////////////////////////////////////// |
---|
6 | // |
---|
7 | // Calculates the scattering from a binary mixture of |
---|
8 | // hard spheres |
---|
9 | // |
---|
10 | // there are some typographical errors in Ashcroft/Langreth's paper |
---|
11 | // Physical Review, v. 156 (1967) 685-692 |
---|
12 | // |
---|
13 | // Errata on Phys. Rev. 166 (1968) 934. |
---|
14 | // |
---|
15 | //(A5) - the entire term should be multiplied by 1/2 |
---|
16 | // |
---|
17 | //final equation for beta12 should be (1+a) rather than (1-a) |
---|
18 | // |
---|
19 | // |
---|
20 | //Definitions are consistent with notation in the paper: |
---|
21 | // |
---|
22 | // phi is total volume fraction |
---|
23 | // nf2 (x) is number density ratio as defined in paper |
---|
24 | // aa = alpha as defined in paper |
---|
25 | // r2 is the radius of the LARGER sphere (angstroms) |
---|
26 | // Sij are the partial structure factor output arrays |
---|
27 | // |
---|
28 | // S. Kline 15 JUL 2004 |
---|
29 | // see: bhs.c and ashcroft.f |
---|
30 | // |
---|
31 | //////////////////////////////////////////////////// |
---|
32 | |
---|
33 | //this macro sets up all the necessary parameters and waves that are |
---|
34 | //needed to calculate the model function. |
---|
35 | // |
---|
36 | // larger sphere radius(angstroms) = guess[0] |
---|
37 | // smaller sphere radius (A) = guess[1] |
---|
38 | // volume fraction of larger spheres = guess[2] |
---|
39 | // volume fraction of small spheres = guess[3] |
---|
40 | // size ratio, alpha(0<a<1) = derived |
---|
41 | // SLD(A-2) of larger particle = guess[4] |
---|
42 | // SLD(A-2) of smaller particle = guess[5] |
---|
43 | // SLD(A-2) of the solvent = guess[6] |
---|
44 | //background = guess[7] |
---|
45 | |
---|
46 | Proc PlotBinaryHS(num,qmin,qmax) |
---|
47 | Variable num=256, qmin=.001, qmax=.7 |
---|
48 | Prompt num "Enter number of data points for model: " |
---|
49 | Prompt qmin "Enter minimum q-value (A-1) for model: " |
---|
50 | Prompt qmax "Enter maximum q-value (A-1) for model: " |
---|
51 | // |
---|
52 | Make/O/D/n=(num) xwave_BinaryHS, ywave_BinaryHS |
---|
53 | xwave_BinaryHS = alog(log(qmin) + x*((log(qmax)-log(qmin))/num)) |
---|
54 | Make/O/D coef_BinaryHS = {100,25,0.2,0.1,3.5e-6,0.5e-6,6.36e-6,0.001} //CH#2 |
---|
55 | make/o/t parameters_BinaryHS = {"large radius","small radius","volume fraction large spheres","volume fraction small spheres","large sphere SLD","small sphere SLD","solvent SLD","Incoherent Bgd (cm-1)"} |
---|
56 | Edit parameters_BinaryHS, coef_BinaryHS |
---|
57 | ModifyTable width(parameters_BinaryHS)=160 |
---|
58 | ModifyTable width(coef_BinaryHS)=90 |
---|
59 | |
---|
60 | Variable/G root:g_BinaryHS |
---|
61 | g_BinaryHS := BinaryHS(coef_BinaryHS, ywave_BinaryHS,xwave_BinaryHS) |
---|
62 | Display ywave_BinaryHS vs xwave_BinaryHS |
---|
63 | ModifyGraph marker=29, msize=2, mode=4 |
---|
64 | ModifyGraph log=1 |
---|
65 | Label bottom "q (A\\S-1\\M) " |
---|
66 | Label left "I(q) (cm\\S-1\\M)" |
---|
67 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
68 | |
---|
69 | AddModelToStrings("BinaryHS","coef_BinaryHS","BinaryHS") |
---|
70 | End |
---|
71 | |
---|
72 | // - sets up a dependency to a wrapper, not the actual SmearedModelFunction |
---|
73 | Proc PlotSmearedBinaryHS(str) |
---|
74 | String str |
---|
75 | Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4) |
---|
76 | |
---|
77 | // if any of the resolution waves are missing => abort |
---|
78 | if(ResolutionWavesMissingDF(str)) //updated to NOT use global strings (in GaussUtils) |
---|
79 | Abort |
---|
80 | endif |
---|
81 | |
---|
82 | SetDataFolder $("root:"+str) |
---|
83 | |
---|
84 | // Setup parameter table for model function |
---|
85 | Make/O/D smear_coef_BinaryHS = {100,25,0.2,0.1,3.5e-6,0.5e-6,6.36e-6,0.001} //CH#4 |
---|
86 | make/o/t smear_parameters_BinaryHS = {"large radius","small radius","volume fraction large spheres","volume fraction small spheres","large sphere SLD","small sphere SLD","solvent SLD","Incoherent Bgd (cm-1)"} |
---|
87 | Edit smear_parameters_BinaryHS,smear_coef_BinaryHS //display parameters in a table |
---|
88 | ModifyTable width(smear_parameters_BinaryHS)=160 |
---|
89 | ModifyTable width(smear_coef_BinaryHS)=90 |
---|
90 | // output smeared intensity wave, dimensions are identical to experimental QSIG values |
---|
91 | // make extra copy of experimental q-values for easy plotting |
---|
92 | Duplicate/O $(str+"_q") smeared_BinaryHS,smeared_qvals // |
---|
93 | SetScale d,0,0,"1/cm",smeared_BinaryHS // |
---|
94 | |
---|
95 | Variable/G gs_BinaryHS=0 |
---|
96 | gs_BinaryHS := fSmearedBinaryHS(smear_coef_BinaryHS,smeared_BinaryHS,smeared_qvals) //this wrapper fills the STRUCT |
---|
97 | |
---|
98 | Display smeared_BinaryHS vs smeared_qvals // |
---|
99 | ModifyGraph log=1,marker=29,msize=2,mode=4 |
---|
100 | Label bottom "q (A\\S-1\\M)" |
---|
101 | Label left "I(q) (cm\\S-1\\M)" |
---|
102 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
103 | |
---|
104 | SetDataFolder root: |
---|
105 | AddModelToStrings("SmearedBinaryHS","smear_coef_BinaryHS","BinaryHS") |
---|
106 | End |
---|
107 | |
---|
108 | |
---|
109 | |
---|
110 | //AAO version, uses XOP if available |
---|
111 | // simply calls the original single point calculation with |
---|
112 | // a wave assignment (this will behave nicely if given point ranges) |
---|
113 | Function BinaryHS(cw,yw,xw) : FitFunc |
---|
114 | Wave cw,yw,xw |
---|
115 | |
---|
116 | #if exists("BinaryHSX") |
---|
117 | yw = BinaryHSX(cw,xw) |
---|
118 | #else |
---|
119 | yw = fBinaryHS(cw,xw) |
---|
120 | #endif |
---|
121 | return(0) |
---|
122 | End |
---|
123 | |
---|
124 | |
---|
125 | //CH#1 |
---|
126 | // you should write your function to calculate the intensity |
---|
127 | // for a single q-value (that's the input parameter x) |
---|
128 | // based on the wave (array) of parameters that you send it (w) |
---|
129 | // |
---|
130 | Function fBinaryHS(w,x) : FitFunc |
---|
131 | Wave w |
---|
132 | Variable x |
---|
133 | // Input (fitting) variables are: |
---|
134 | // larger sphere radius(angstroms) = guess[0] |
---|
135 | // smaller sphere radius (A) = w[1] |
---|
136 | // number fraction of larger spheres = guess[2] |
---|
137 | // total volume fraction of spheres = guess[3] |
---|
138 | // size ratio, alpha(0<a<1) = derived |
---|
139 | // SLD(A-2) of larger particle = guess[4] |
---|
140 | // SLD(A-2) of smaller particle = guess[5] |
---|
141 | // SLD(A-2) of the solvent = guess[6] |
---|
142 | //background = guess[7] |
---|
143 | |
---|
144 | // give them nice names |
---|
145 | Variable r2,r1,nf2,phi,aa,rho2,rho1,rhos,inten,bgd |
---|
146 | Variable err,psf11,psf12,psf22 |
---|
147 | Variable phi1,phi2,phr,a3 |
---|
148 | |
---|
149 | r2 = w[0] |
---|
150 | r1 = w[1] |
---|
151 | phi2 = w[2] |
---|
152 | phi1 = w[3] |
---|
153 | rho2 = w[4] |
---|
154 | rho1 = w[5] |
---|
155 | rhos = w[6] |
---|
156 | bgd = w[7] |
---|
157 | |
---|
158 | phi = w[2] + w[3] //total volume fraction |
---|
159 | aa = r1/r2 //alpha(0<a<1) |
---|
160 | |
---|
161 | //calculate the number fraction of larger spheres (eqn 2 in reference) |
---|
162 | a3=aa^3 |
---|
163 | phr=phi2/phi |
---|
164 | nf2 = phr*a3/(1-phr+phr*a3) |
---|
165 | // calculate the PSF's here |
---|
166 | |
---|
167 | err = Ashcroft(x,r2,nf2,aa,phi,psf11,psf22,psf12) |
---|
168 | |
---|
169 | // /* do form factor calculations */ |
---|
170 | Variable v1,v2,n1,n2,qr1,qr2,b1,b2 |
---|
171 | v1 = 4.0*PI/3.0*r1*r1*r1 |
---|
172 | v2 = 4.0*PI/3.0*r2*r2*r2 |
---|
173 | // a3 = aa*aa*aa |
---|
174 | // phi1 = phi*(1.0-nf2)*a3/(nf2+(1.0-nf2)*a3) |
---|
175 | // phi2 = phi - phi1 |
---|
176 | n1 = phi1/v1 |
---|
177 | n2 = phi2/v2 |
---|
178 | |
---|
179 | qr1 = r1*x |
---|
180 | qr2 = r2*x |
---|
181 | b1 = r1*r1*r1*(rho1-rhos)*BHSbfunc(qr1) |
---|
182 | b2 = r2*r2*r2*(rho2-rhos)*BHSbfunc(qr2) |
---|
183 | inten = n1*b1*b1*psf11 |
---|
184 | inten += sqrt(n1*n2)*2.0*b1*b2*psf12 |
---|
185 | inten += n2*b2*b2*psf22 |
---|
186 | ///* convert I(1/A) to (1/cm) */ |
---|
187 | inten *= 1.0e8 |
---|
188 | |
---|
189 | inten += bgd |
---|
190 | Return (inten) |
---|
191 | End |
---|
192 | |
---|
193 | //AAO version, uses XOP if available |
---|
194 | // simply calls the original single point calculation with |
---|
195 | // a wave assignment (this will behave nicely if given point ranges) |
---|
196 | Function BinaryHS_PSF11(cw,yw,xw) : FitFunc |
---|
197 | Wave cw,yw,xw |
---|
198 | |
---|
199 | #if exists("BinaryHS_PSF11X") |
---|
200 | yw = BinaryHS_PSF11X(cw,xw) |
---|
201 | #else |
---|
202 | yw = fBinaryHS_PSF11(cw,xw) |
---|
203 | #endif |
---|
204 | return(0) |
---|
205 | End |
---|
206 | |
---|
207 | Function fBinaryHS_PSF11(w,x) : FitFunc |
---|
208 | Wave w |
---|
209 | Variable x |
---|
210 | // Input (fitting) variables are: |
---|
211 | // larger sphere radius(angstroms) = guess[0] |
---|
212 | // smaller sphere radius (A) = w[1] |
---|
213 | // number fraction of larger spheres = guess[2] |
---|
214 | // total volume fraction of spheres = guess[3] |
---|
215 | // size ratio, alpha(0<a<1) = derived |
---|
216 | // SLD(A-2) of larger particle = guess[4] |
---|
217 | // SLD(A-2) of smaller particle = guess[5] |
---|
218 | // SLD(A-2) of the solvent = guess[6] |
---|
219 | //background = guess[7] |
---|
220 | |
---|
221 | // give them nice names |
---|
222 | Variable r2,r1,nf2,phi,aa,rho2,rho1,rhos,inten,bgd |
---|
223 | Variable err,psf11,psf12,psf22 |
---|
224 | Variable phi1,phi2,a3,phr |
---|
225 | |
---|
226 | r2 = w[0] |
---|
227 | r1 = w[1] |
---|
228 | phi2 = w[2] |
---|
229 | phi1 = w[3] |
---|
230 | rho2 = w[4] |
---|
231 | rho1 = w[5] |
---|
232 | rhos = w[6] |
---|
233 | bgd = w[7] |
---|
234 | |
---|
235 | phi = w[2] + w[3] //total volume fraction |
---|
236 | aa = r1/r2 //alpha(0<a<1) |
---|
237 | |
---|
238 | //calculate the number fraction of larger spheres (eqn 2 in reference) |
---|
239 | a3=aa^3 |
---|
240 | phr=phi2/phi |
---|
241 | nf2 = phr*a3/(1-phr+phr*a3) |
---|
242 | |
---|
243 | // calculate the PSF's here |
---|
244 | |
---|
245 | err = Ashcroft(x,r2,nf2,aa,phi,psf11,psf22,psf12) |
---|
246 | return(psf11) |
---|
247 | End |
---|
248 | |
---|
249 | //AAO version, uses XOP if available |
---|
250 | // simply calls the original single point calculation with |
---|
251 | // a wave assignment (this will behave nicely if given point ranges) |
---|
252 | Function BinaryHS_PSF12(cw,yw,xw) : FitFunc |
---|
253 | Wave cw,yw,xw |
---|
254 | |
---|
255 | #if exists("BinaryHS_PSF12X") |
---|
256 | yw = BinaryHS_PSF12X(cw,xw) |
---|
257 | #else |
---|
258 | yw = fBinaryHS_PSF12(cw,xw) |
---|
259 | #endif |
---|
260 | return(0) |
---|
261 | End |
---|
262 | |
---|
263 | Function fBinaryHS_PSF12(w,x) : FitFunc |
---|
264 | Wave w |
---|
265 | Variable x |
---|
266 | // Input (fitting) variables are: |
---|
267 | // larger sphere radius(angstroms) = guess[0] |
---|
268 | // smaller sphere radius (A) = w[1] |
---|
269 | // number fraction of larger spheres = guess[2] |
---|
270 | // total volume fraction of spheres = guess[3] |
---|
271 | // size ratio, alpha(0<a<1) = derived |
---|
272 | // SLD(A-2) of larger particle = guess[4] |
---|
273 | // SLD(A-2) of smaller particle = guess[5] |
---|
274 | // SLD(A-2) of the solvent = guess[6] |
---|
275 | //background = guess[7] |
---|
276 | |
---|
277 | // give them nice names |
---|
278 | Variable r2,r1,nf2,phi,aa,rho2,rho1,rhos,inten,bgd |
---|
279 | Variable err,psf11,psf12,psf22 |
---|
280 | Variable phi1,phi2,a3,phr |
---|
281 | |
---|
282 | r2 = w[0] |
---|
283 | r1 = w[1] |
---|
284 | phi2 = w[2] |
---|
285 | phi1 = w[3] |
---|
286 | rho2 = w[4] |
---|
287 | rho1 = w[5] |
---|
288 | rhos = w[6] |
---|
289 | bgd = w[7] |
---|
290 | |
---|
291 | phi = w[2] + w[3] //total volume fraction |
---|
292 | aa = r1/r2 //alpha(0<a<1) |
---|
293 | |
---|
294 | //calculate the number fraction of larger spheres (eqn 2 in reference) |
---|
295 | a3=aa^3 |
---|
296 | phr=phi2/phi |
---|
297 | nf2 = phr*a3/(1-phr+phr*a3) |
---|
298 | |
---|
299 | // calculate the PSF's here |
---|
300 | |
---|
301 | err = Ashcroft(x,r2,nf2,aa,phi,psf11,psf22,psf12) |
---|
302 | return(psf12) |
---|
303 | End |
---|
304 | |
---|
305 | //AAO version, uses XOP if available |
---|
306 | // simply calls the original single point calculation with |
---|
307 | // a wave assignment (this will behave nicely if given point ranges) |
---|
308 | Function BinaryHS_PSF22(cw,yw,xw) : FitFunc |
---|
309 | Wave cw,yw,xw |
---|
310 | |
---|
311 | #if exists("BinaryHS_PSF22X") |
---|
312 | yw = BinaryHS_PSF22X(cw,xw) |
---|
313 | #else |
---|
314 | yw = fBinaryHS_PSF22(cw,xw) |
---|
315 | #endif |
---|
316 | return(0) |
---|
317 | End |
---|
318 | |
---|
319 | Function fBinaryHS_PSF22(w,x) : FitFunc |
---|
320 | Wave w |
---|
321 | Variable x |
---|
322 | // Input (fitting) variables are: |
---|
323 | // larger sphere radius(angstroms) = guess[0] |
---|
324 | // smaller sphere radius (A) = w[1] |
---|
325 | // number fraction of larger spheres = guess[2] |
---|
326 | // total volume fraction of spheres = guess[3] |
---|
327 | // size ratio, alpha(0<a<1) = derived |
---|
328 | // SLD(A-2) of larger particle = guess[4] |
---|
329 | // SLD(A-2) of smaller particle = guess[5] |
---|
330 | // SLD(A-2) of the solvent = guess[6] |
---|
331 | //background = guess[7] |
---|
332 | |
---|
333 | // give them nice names |
---|
334 | Variable r2,r1,nf2,phi,aa,rho2,rho1,rhos,inten,bgd |
---|
335 | Variable err,psf11,psf12,psf22 |
---|
336 | Variable phi1,phi2,phr,a3 |
---|
337 | |
---|
338 | r2 = w[0] |
---|
339 | r1 = w[1] |
---|
340 | phi2 = w[2] |
---|
341 | phi1 = w[3] |
---|
342 | rho2 = w[4] |
---|
343 | rho1 = w[5] |
---|
344 | rhos = w[6] |
---|
345 | bgd = w[7] |
---|
346 | |
---|
347 | phi = w[2] + w[3] //total volume fraction |
---|
348 | aa = r1/r2 //alpha(0<a<1) |
---|
349 | |
---|
350 | //calculate the number fraction of larger spheres (eqn 2 in reference) |
---|
351 | a3=aa^3 |
---|
352 | phr=phi2/phi |
---|
353 | nf2 = phr*a3/(1-phr+phr*a3) |
---|
354 | // calculate the PSF's here |
---|
355 | |
---|
356 | err = Ashcroft(x,r2,nf2,aa,phi,psf11,psf22,psf12) |
---|
357 | return(psf22) |
---|
358 | End |
---|
359 | |
---|
360 | Function BHSbfunc(qr) |
---|
361 | Variable qr |
---|
362 | |
---|
363 | Variable ans |
---|
364 | |
---|
365 | ans = 4.0*pi*(sin(qr)-qr*cos(qr))/qr/qr/qr |
---|
366 | return(ans) |
---|
367 | End |
---|
368 | |
---|
369 | |
---|
370 | Function Ashcroft(qval,r2,nf2,aa,phi,s11,s22,s12) |
---|
371 | Variable qval,r2,nf2,aa,phi,&s11,&s22,&s12 |
---|
372 | |
---|
373 | // CALCULATE CONSTANT TERMS |
---|
374 | Variable s1,s2,v,a3,v1,v2,g11,g12,g22,wmv,wmv3,wmv4 |
---|
375 | Variable a1,a2i,a2,b1,b2,b12,gm1,gm12 |
---|
376 | Variable err,yy,ay,ay2,ay3,t1,t2,t3,f11,y2,y3,tt1,tt2,tt3 |
---|
377 | Variable c11,c22,c12,f12,f22,ttt1,ttt2,ttt3,ttt4,yl,y13 |
---|
378 | Variable t21,t22,t23,t31,t32,t33,t41,t42,yl3,wma3,y1 |
---|
379 | |
---|
380 | s2 = 2.0*r2 |
---|
381 | s1 = aa*s2 |
---|
382 | v = phi |
---|
383 | a3 = aa*aa*aa |
---|
384 | V1=((1.-nf2)*A3/(nf2+(1.-nf2)*A3))*V |
---|
385 | V2=(nf2/(nf2+(1.-nf2)*A3))*V |
---|
386 | G11=((1.+.5*V)+1.5*V2*(aa-1.))/(1.-V)/(1.-V) |
---|
387 | G22=((1.+.5*V)+1.5*V1*(1./aa-1.))/(1.-V)/(1.-v) |
---|
388 | G12=((1.+.5*V)+1.5*(1.-aa)*(V1-V2)/(1.+aa))/(1.-V)/(1.-v) |
---|
389 | wmv = 1/(1.-v) |
---|
390 | wmv3 = wmv*wmv*wmv |
---|
391 | wmv4 = wmv*wmv3 |
---|
392 | A1=3.*wmv4*((V1+A3*V2)*(1.+V+V*v)-3.*V1*V2*(1.-aa)*(1.-aa)*(1.+V1+aa*(1.+V2))) + ((V1+A3*V2)*(1.+2.*V)+(1.+V+V*v)-3.*V1*V2*(1.-aa)*(1.-aa)-3.*V2*(1.-aa)*(1.-aa)*(1.+V1+aa*(1.+V2)))*wmv3 |
---|
393 | A2I=((V1+A3*V2)*(1.+V+V*v)-3.*V1*V2*(1.-aa)*(1.-aa)*(1.+V1+aa*(1.+V2)))*3*wmv4 + ((V1+A3*V2)*(1.+2.*V)+A3*(1.+V+V*v)-3.*V1*V2*(1.-aa)*(1.-aa)*aa-3.*V1*(1.-aa)*(1.-aa)*(1.+V1+aa*(1.+V2)))*wmv3 |
---|
394 | A2=A2I/a3 |
---|
395 | B1=-6.*(V1*G11*g11+.25*V2*(1.+aa)*(1.+aa)*aa*G12*g12) |
---|
396 | B2=-6.*(V2*G22*g22+.25*V1/A3*(1.+aa)*(1.+aa)*G12*g12) |
---|
397 | B12=-3.*aa*(1.+aa)*(V1*G11/aa/aa+V2*G22)*G12 |
---|
398 | GM1=(V1*A1+A3*V2*A2)*.5 |
---|
399 | GM12=2.*GM1*(1.-aa)/aa |
---|
400 | //C |
---|
401 | //C CALCULATE THE DIRECT CORRELATION FUNCTIONS AND PRINT RESULTS |
---|
402 | //C |
---|
403 | // DO 20 J=1,npts |
---|
404 | |
---|
405 | yy=qval*s2 |
---|
406 | //c calculate direct correlation functions |
---|
407 | //c ----c11 |
---|
408 | AY=aa*yy |
---|
409 | ay2 = ay*ay |
---|
410 | ay3 = ay*ay*ay |
---|
411 | T1=A1*(SIN(AY)-AY*COS(AY)) |
---|
412 | T2=B1*(2.*AY*sin(AY)-(AY2-2.)*cos(AY)-2.)/AY |
---|
413 | T3=GM1*((4.*AY*ay2-24.*AY)*sin(AY)-(AY2*ay2-12.*AY2+24.)*cos(AY)+24.)/AY3 |
---|
414 | F11=24.*V1*(T1+T2+T3)/AY3 |
---|
415 | |
---|
416 | //c ------c22 |
---|
417 | y2=yy*yy |
---|
418 | y3=yy*y2 |
---|
419 | TT1=A2*(sin(yy)-yy*cos(yy)) |
---|
420 | TT2=B2*(2.*yy*sin(yy)-(Y2-2.)*cos(yy)-2.)/yy |
---|
421 | TT3=GM1*((4.*Y3-24.*yy)*sin(yy)-(Y2*y2-12.*Y2+24.)*cos(yy)+24.)/ay3 |
---|
422 | F22=24.*V2*(TT1+TT2+TT3)/Y3 |
---|
423 | |
---|
424 | //c -----c12 |
---|
425 | YL=.5*yy*(1.-aa) |
---|
426 | yl3=yl*yl*yl |
---|
427 | wma3 = (1.-aa)*(1.-aa)*(1.-aa) |
---|
428 | Y1=aa*yy |
---|
429 | y13 = y1*y1*y1 |
---|
430 | TTT1=3.*wma3*V*sqrt(nf2)*sqrt(1.-nf2)*A1*(sin(YL)-YL*cos(YL))/((nf2+(1.-nf2)*A3)*YL3) |
---|
431 | T21=B12*(2.*Y1*cos(Y1)+(Y1^2-2.)*sin(Y1)) |
---|
432 | T22=GM12*((3.*Y1*y1-6.)*cos(Y1)+(Y1^3-6.*Y1)*sin(Y1)+6.)/Y1 |
---|
433 | T23=GM1*((4.*Y13-24.*Y1)*cos(Y1)+(Y13*y1-12.*Y1*y1+24.)*sin(Y1))/(Y1*y1) |
---|
434 | T31=B12*(2.*Y1*sin(Y1)-(Y1^2-2.)*cos(Y1)-2.) |
---|
435 | T32=GM12*((3.*Y1^2-6.)*sin(Y1)-(Y1^3-6.*Y1)*cos(Y1))/Y1 |
---|
436 | T33=GM1*((4.*Y13-24.*Y1)*sin(Y1)-(Y13*y1-12.*Y1*y1+24.)*cos(Y1)+24.)/(y1*y1) |
---|
437 | T41=cos(YL)*((sin(Y1)-Y1*cos(Y1))/(Y1*y1) + (1.-aa)/(2.*aa)*(1.-cos(Y1))/Y1) |
---|
438 | T42=sin(YL)*((cos(Y1)+Y1*sin(Y1)-1.)/(Y1*y1) + (1.-aa)/(2.*aa)*sin(Y1)/Y1) |
---|
439 | TTT2=sin(YL)*(T21+T22+T23)/(y13*y1) |
---|
440 | TTT3=cos(YL)*(T31+T32+T33)/(y13*y1) |
---|
441 | TTT4=A1*(T41+T42)/Y1 |
---|
442 | F12=TTT1+24.*V*sqrt(nf2)*sqrt(1.-nf2)*A3*(TTT2+TTT3+TTT4)/(nf2+(1.-nf2)*A3) |
---|
443 | |
---|
444 | C11=F11 |
---|
445 | C22=F22 |
---|
446 | C12=F12 |
---|
447 | S11=1./(1.+C11-(C12)*c12/(1.+C22)) |
---|
448 | S22=1./(1.+C22-(C12)*c12/(1.+C11)) |
---|
449 | S12=-C12/((1.+C11)*(1.+C22)-(C12)*(c12)) |
---|
450 | |
---|
451 | return(err) |
---|
452 | End |
---|
453 | |
---|
454 | |
---|
455 | //wrapper to calculate the smeared model as an AAO-Struct |
---|
456 | // fills the struct and calls the ususal function with the STRUCT parameter |
---|
457 | // |
---|
458 | // used only for the dependency, not for fitting |
---|
459 | // |
---|
460 | Function fSmearedBinaryHS(coefW,yW,xW) |
---|
461 | Wave coefW,yW,xW |
---|
462 | |
---|
463 | String str = getWavesDataFolder(yW,0) |
---|
464 | String DF="root:"+str+":" |
---|
465 | |
---|
466 | WAVE resW = $(DF+str+"_res") |
---|
467 | |
---|
468 | STRUCT ResSmearAAOStruct fs |
---|
469 | WAVE fs.coefW = coefW |
---|
470 | WAVE fs.yW = yW |
---|
471 | WAVE fs.xW = xW |
---|
472 | WAVE fs.resW = resW |
---|
473 | |
---|
474 | Variable err |
---|
475 | err = SmearedBinaryHS(fs) |
---|
476 | |
---|
477 | return (0) |
---|
478 | End |
---|
479 | |
---|
480 | // this is all there is to the smeared calculation! |
---|
481 | Function SmearedBinaryHS(s) :FitFunc |
---|
482 | Struct ResSmearAAOStruct &s |
---|
483 | |
---|
484 | // the name of your unsmeared model (AAO) is the first argument |
---|
485 | Smear_Model_20(BinaryHS,s.coefW,s.xW,s.yW,s.resW) |
---|
486 | |
---|
487 | return(0) |
---|
488 | End |
---|
489 | |
---|
490 | |
---|
491 | |
---|
492 | Macro Plot_BinaryHS_PSF() |
---|
493 | |
---|
494 | if(Exists("coef_BinaryHS") != 1) |
---|
495 | abort "You need to plot the unsmeared model first to create the coefficient table" |
---|
496 | Endif |
---|
497 | |
---|
498 | Make/O/D/n=(numpnts(xwave_BinaryHS)) psf11_BinaryHS,psf12_BinaryHS,psf22_BinaryHS,QD2_BinaryHS |
---|
499 | |
---|
500 | Variable/G root:g_psf11,root:g_psf12,root:g_psf22,root:g_QD2 |
---|
501 | |
---|
502 | g_psf11 := BinaryHS_psf11(coef_BinaryHS, psf11_BinaryHS, xwave_BinaryHS) |
---|
503 | g_psf12 := BinaryHS_psf12(coef_BinaryHS, psf12_BinaryHS, xwave_BinaryHS) |
---|
504 | g_psf22 := BinaryHS_psf22(coef_BinaryHS, psf22_BinaryHS, xwave_BinaryHS) |
---|
505 | QD2_BinaryHS := xwave_BinaryHS*coef_BinaryHS[0]*2 |
---|
506 | // Display psf11_BinaryHS vs xwave_BinaryHS |
---|
507 | // AppendtoGraph psf12_BinaryHS vs xwave_BinaryHS |
---|
508 | // AppendtoGraph psf22_BinaryHS vs xwave_BinaryHS |
---|
509 | |
---|
510 | Display psf11_BinaryHS vs QD2_BinaryHS |
---|
511 | AppendtoGraph psf12_BinaryHS vs QD2_BinaryHS |
---|
512 | AppendtoGraph psf22_BinaryHS vs QD2_BinaryHS |
---|
513 | |
---|
514 | ModifyGraph marker=19, msize=2, mode=4 |
---|
515 | ModifyGraph lsize=2,rgb(psf12_BinaryHS)=(2,39321,1) |
---|
516 | ModifyGraph rgb(psf22_BinaryHS)=(0,0,65535) |
---|
517 | ModifyGraph log=0,grid=1,mirror=2 |
---|
518 | SetAxis bottom 0,30 |
---|
519 | Label bottom "q*LargeDiameter" |
---|
520 | Label left "Sij(q)" |
---|
521 | Legend |
---|
522 | // |
---|
523 | End |
---|
524 | |
---|
525 | |
---|
526 | //useful for finding the parameters that duplicate the |
---|
527 | //figures in the original reference (uses the same notation) |
---|
528 | //automatically changes the coefficient wave |
---|
529 | Macro Duplicate_AL_Parameters(eta,xx,alpha,Rlarge) |
---|
530 | Variable eta=0.45,xx=0.4,alpha=0.7,Rlarge=100 |
---|
531 | |
---|
532 | Variable r1,phi1,phi2,a3 |
---|
533 | r1 = alpha*Rlarge |
---|
534 | a3 = alpha*alpha*alpha |
---|
535 | phi1 = eta*(1.0-xx)*a3/(xx+(1.0-xx)*a3) //eqn [2] |
---|
536 | phi2 = eta - phi1 |
---|
537 | |
---|
538 | Print "phi (larger) = ",phi2 |
---|
539 | Print "phi (smaller) = ",phi1 |
---|
540 | Print "Radius (smaller) = ",r1 |
---|
541 | |
---|
542 | if(Exists("coef_BinaryHS") != 1) |
---|
543 | abort "You need to plot the unsmeared model first to create the coefficient table" |
---|
544 | Endif |
---|
545 | |
---|
546 | coef_BinaryHS[2] = phi2 |
---|
547 | coef_BinaryHS[3] = phi1 |
---|
548 | coef_BinaryHS[1] = r1 |
---|
549 | End |
---|
550 | |
---|
551 | |
---|
552 | //calculates number fractions of each population based on the |
---|
553 | //coef_BinaryHS parameters |
---|
554 | Macro Calculate_BHS_Parameters() |
---|
555 | |
---|
556 | if(exists("coef_BinaryHS") != 1) |
---|
557 | Abort "You must plot the unsmeared BHS model first to create the coefficient wave" |
---|
558 | endif |
---|
559 | Variable r1,r2,phi1,phi2,aa //same notation as paper - r2 is LARGER |
---|
560 | Variable a3,xx,phi |
---|
561 | r1 = coef_BinaryHS[1] |
---|
562 | r2 = coef_BinaryHS[0] |
---|
563 | phi1 = coef_BinaryHS[3] |
---|
564 | phi2 = coef_BinaryHS[2] |
---|
565 | |
---|
566 | phi = phi1+phi2 |
---|
567 | aa = r1/r2 |
---|
568 | a3 = aa^3 |
---|
569 | |
---|
570 | xx = phi2/phi*a3 |
---|
571 | xx /= (1-(phi2/phi)+(phi2/phi)*a3) |
---|
572 | |
---|
573 | Print "Number fraction (larger) = ",xx |
---|
574 | Print "Number fraction (smaller) = ",1-xx |
---|
575 | |
---|
576 | End |
---|